[e-mail]

Pytania – metody symulacyjne

1. Na czym polega metoda dynamiki molekularnej? Podad zwięzłą, ale celną charakterystykę ogólną.

Czym się ona różni od metody Monte Carlo?

2. Co rozumiesz przez pojęcie „pudło symulacyjne”? Co to są i czemu służą „periodyczne warunki

brzegowe”? Optymalny algorytm ich realizacji. Niekonwencjonalne metody budowy układów
pseudo-nieskooczonych.

3. Jakie cechy powinny mied algorytmy całkowania równao ruchu w metodach dynamiki molekularnej?

Wyprowadzid algorytm Verleta i „żabiego skoku”.

4. Jak sprawdzamy formalną (numeryczną) poprawnośd symulacji dynamiczno-molekularnych?

Oszacowad długośd kroku czasowego odpowiedniego do realizacji symulacji dynamiczno-
molekularnych.

5. Na czym polega problem „identyfikacji najbliższych sąsiadów”? Opisad metodę listy Verleta. Jaki

wpływ na czas obliczeo ma wprowadzenie dwóch promieni Verleta (wewnętrznego i
zewnętrznego)?

6. Opisad metodę cel połączonych (Hockney’a). Czemu ona służy? Podad jej wady i zalety.
7. Co to są poprawki dalekozasięgowe? Kiedy je trzeba liczyd? Podad wyrażenie na poprawkę

dalekozasięgową dla energii potencjalnej dla oddziaływao dwuciałowych. Podad przykłady wielkości
fizycznych, dla których konieczne jest uwzględnianie tych poprawek.

8. Omówid metody charakteryzacji uporządkowania bliskiego zasięgu w materiałach symulowanych

numerycznie (radialne i kątowe funkcje rozkładu, metody geometrii stochastycznej).

9. Co oznacza pojecie „termalizacji” i „próbkowania” układu w symulacjach dynamiczno-

molekularnych? Na czym polega, jak jest realizowane i czemu służy w metodach dynamiczno-
molekularnych tzw. „skalowanie prędkości”?

10. Omówid sposób generacji położeo i prędkości początkowych w symulacjach dynamiczno-

molekularnych.

11. Jak w symulacjach dynamiczno-molekularnych obliczamy współczynnik dyfuzji? Jaką informację

niesie jego zależnośd od temperatury? Jak w symulacjach MD obliczamy średnie przemieszczenie
kwadratowe? Jaki charakter musi mied jego zależnośd od czasu, aby współczynnik dyfuzji nie zależał
od czasu?

12. Omówid metodę obliczania oddziaływao kulombowskich w symulacjach dynamiczno-molekularnych

(metoda Ewalda). Skąd wynikają „kłopoty” z liczeniem sił kulombowskich.

13. Co to jest model powłokowy (Kochrana i Finchama)? Przy symulacjach jakich związków konieczne

jest jego stosowanie? Czym różni się przebieg symulacji z modelem powłokowym oddziaływao od
symulacji układów izotropowych cząstek punktowych?

14. Scharakteryzowad symulację układów molekularnych. Co to jest metoda kwaternionowi. Dlaczego

jej używamy? Podad podobieostwa i różnice symulacji układów molekuł sztywnych i molekuł
deformowalnych.

15. Na czym polegają i co mają na celu „metody układów skrępowanych więzami”? Podad przykład

termostatu realizowanego tą metodą.

16. Na czym polega i co mają na celu metody układów rozszerzonych? Podad przykład barostatu i

termostatu realizowanych tą metodą.Na jakiej zasadzie działają termostaty i barostaty należące do
klasy układów rozszerzalnych (o dodatkowe stopnie swobody, ESM)? Omów wybrany przez siebie
przykład barostatu.

17. Omów algorytm Monte Carlo symulacji strukturalnych materiałów w zespole (NpT).

[e-mail]

Ad. 1.

Dynamika molekularna jest metodą numerycznego całkowania równao ruchu układów
wielocząsteczkowych. Każda cząstka tworząca układ podlega klasycznym prawom ruchu, a
makroskopowe parametry opisujące stan układu oblicza się jako średnie po trajektorii w przestrzeni
fazowe.

Opis działania MD:
- W każdym kroku czasowym obliczamy dla każdej cząstki działającą na nią siłę pochodzącą od
pozostałych cząstek.
- Korzystamy z obliczonych sił i znając położenia cząstek w poprzednim kroku obliczamy nowe położenia
i pędy każdej cząstki numerycznie rozwiązując równania ruchu Newtona.
- Po wyznaczeniu parametrów mikroskopowych (czyli nowych prędkości i położeo) w kilku krokach,
możemy obliczyd wielkości makroskopowe (np. temperaturę, energię kinetyczną i całkowitą, ciepło
właściwe, współczynnik dyfuzji, przewodnośd elektryczną i cieplną i lepkośd ...).

Krótko mówiąc metoda dynamiki molekularnej wylicza bezpośrednio położenie molekuł (atomów) na
podstawie równao ruchu, a metoda Monte Carlo wykorzystuje narzędzie statystyczne.

Cechy:
- nie zależy jawnie od czasu,

- siła nie zależy od prędkości 𝑚𝑉

𝑖

= 𝐹

𝑖

- suma energii kinetycznej i potencjalnej pozostaje stała.

Metoda Monte Carlo (MC) polega na przedstawieniu rozwiązania problemu w postaci parametru pewnej
hipotetycznej populacji i używaniu losowych liczb do tworzenia próbki tej populacji, na podstawie której
można dokonad statystycznego oszacowania wartości badanego parametru.

MD- metoda determistyczna (wynikowa),
MC- metoda stochastyczna (losowa)- losuje stałej przestrzeni różne stany.

[e-mail]

Ad. 2.

Pudło symulacyjne- komórka MD, obszar w przestrzeni o stałej gęstości, w którym przeprowadza się
symulację (najczęściej w kształcie sześcianu).

Periodyczne warunki brzegowe – komórka bazowa powtarza się nieskooczenie wiele razy. Istnienie
niepotrzebnych ścian w komórce MD spowodowałoby odbijanie się atomów od ich powierzchni, co
znacząco wpływa na właściwości układu.

W 1D można je sobie wyobrazid:

Ψ(x + L) = Ψ(x)

Praktycznie rzecz ujmując, jeżeli cząstka przechodzi przez jedną ścianę komórki, to pojawia się ona przy
przeciwległej ścianie.

Algorytm realizacji:

𝑕

𝑖𝑗

= 𝑥

𝑗

− 𝑥

𝑖

2

+ 𝑦

𝑗

− 𝑦

𝑖

2

+ 𝑧

𝑗

− 𝑧

𝑖

2

∆𝑥 = 𝑥

𝑗

− 𝑥

𝑖

𝑖𝑓 ∆𝑥 ≥

1
2

𝑡𝑜 ∆𝑥 = ∆𝑥 − 𝑙

𝑖𝑓 ∆𝑥 ≤

1
2

𝑡𝑜 ∆𝑥 = ∆𝑥 + 𝑙

i tak samo dla y i z.

Aby przyspieszyd należy znormalizowad długośd do 1.

∆𝑥 = 𝑥

𝑗

− 𝑥

𝑖

∆𝑥 = ∆𝑥 − min 𝑡(∆𝑥) <- wybiera najbliższą liczbę całkowitą.

Budowa układu pseudo nieskooczonego:

Dla dużych układów mamy:

N

2

↑ tę częśd możemy wstawid do wzoru

- unikamy obliczania odległości każdej cząstki z każdą, czyli trzeba znaleźd najbliższych sąsiadów w
promieniu około 15 Å,

- musimy liczyd siłę co krok, czyli za każdym razem, gdy cząstka się poruszy,

- liczymy zamiast N

2

to N N sąsiadów.

[e-mail]

Ad. 3

Cechy algorytmu:

- powinien byd szybki, gdyż metody MD wymagają bardzo wielu obliczeo numerycznych oraz dokładny,

- metoda musi przynajmniej raz na krok ewaluowad prawą stronę równania 𝑚𝑉

𝑖

= 𝐹

𝑖

= 𝐹

𝑖𝑗

𝑁

𝑗 <𝑖

, i= 1,.,N

Metoda Verleta:

𝑚𝑥

𝑖

= 𝑓

𝑖

𝑥

1

, … , 𝑥

3𝑁

+

𝑥 𝑡 + ∆𝑡 = 𝑥 𝑡 + 𝑥

′

𝑡 ∆𝑡 +

1
2

𝑥

′′

𝑡 ∆𝑡

2

+

1
6

𝑥

′′′

𝑡 ∆𝑡

3

+ 𝜗 ∆𝑡

4

𝑥 𝑡 − ∆𝑡 = 𝑥 𝑡 − 𝑥

′

𝑡 ∆𝑡 +

1
2

𝑥

′′

𝑡 ∆𝑡

2

−

1
6

𝑥

′′′

𝑡 ∆𝑡

3

+ 𝜗 ∆𝑡

4

𝑥 𝑡 + ∆𝑡 + 𝑥 𝑡 − ∆𝑡 = 2𝑥 𝑡 + 𝑥

′′

𝑡 ∆𝑡

2

+ 𝜗 ∆𝑡

4

𝑥 𝑡 + ∆𝑡 = 2𝑥 𝑡 − 𝑥 𝑡 − ∆𝑡 +

𝑓

𝑚

∆𝑡

2

krok aktualny krok poprzedni

𝑉 = 𝑥 𝑡 =

𝑥 𝑡 + ∆𝑡 − 𝑥 𝑡 − ∆𝑡

2𝑥∆𝑡

Metoda „żabiego skoku”:

𝑥 𝑡 + ∆𝑡 = 𝑥 𝑡 + 𝑥

′

𝑡 ∆𝑡 +

1
2

𝑥

′′

𝑡 ∆𝑡

2

…

𝑥 𝑡 + ∆𝑡 = 𝑥 𝑡 + ∆𝑡 𝑥

′

𝑡 +

∆𝑡

2

𝑥

′′

𝑡 + ⋯

𝑥 𝑡 +

∆𝑡

2

= 𝑥 𝑡 +

∆𝑡

2

𝑥 (𝑡)

𝑥 𝑡 + ∆𝑡 = 𝑥 𝑡 + ∆𝑡 𝑥 𝑡 +

∆𝑡

2

−

𝑥 𝑡 +

∆𝑡

2

= 𝑥 𝑡 +

∆𝑡

2

𝑥 𝑡

𝑥 𝑡 −

∆𝑡

2

= 𝑥 𝑡 −

∆𝑡

2

𝑥 𝑡

𝑥 𝑡 +

∆𝑡

2

− 𝑥 𝑡 −

∆𝑡

2

= ∆𝑡

𝑓

𝑚

𝑥 𝑡 +

∆𝑡

2

= 𝑥 𝑡 −

∆𝑡

2

∆𝑡

𝑓

𝑚

krok wstecz

prędkości liczone są w połówkach kroku

[e-mail]

Ad. 4

Zjawiska dynamiczno-molekularne zachodzą w czasie rzędu pikosekund (1ps = 10

-12

s), dlatego aby

uzyskad wiarygodne wyniki symulujemy układy z krokiem czasowym rzędu femtosekund (1fs = 10

-15

s).

Krok jest mniejszy o 100 razy niż zjawiska w symulacji. Poprawnośd symulacji sprawdza się za pomocą
zasad zachowania (energii i pędu).

Numerycznie sprawdzamy też:

- stabilnośd temperatury w większości obliczeo,

- metoda musi przynajmniej raz na krok ewaluowad prawą stronę równania 𝑚𝑉

𝑖

= 𝐹

𝑖

= 𝐹

𝑖𝑗

𝑁

𝑗 <𝑖

, i= 1,.,N

Ad. 5

Identyfikacja najbliższych sąsiadów

- dla małych układów (liczymy każdy z każdym), rząd obliczeo N

2

- dla dużych układów, rząd obliczeo N

Problem polega na tym, w jaki sposób najlepiej ustalid liczbę sąsiadów i jak rozłożyd obliczenia, aby
trwały one jak najkrócej.

Lista Verleta

Aby oszczędzid moc obliczeniową możemy wyliczad potencjał tylko do pewnego promienia r

cut

,

zakładając w ten sposób, że potencjał od pozostałych cząstek jest równy zeru. W tym celu tworzymy listę
najbliższych sąsiadów dla każdej cząstki w układzie. Należy pamiętad, że owa lista może się zmieniad po
pewnej ilości kroków czasowych. Dlatego, co pewną ilośd kroków n, należy uaktualniad listę. Wadą tej
metody jest fakt, że samo stworzenie listy jest czasochłonne (~N

2

) a przechowywanie zajmuje miejsce na

twardych dyskach (~N).

Stosujemy gdy nie jest dużo cząstek; zasada działania:

- liczy się cząstki każde z każdą i zapamiętujemy ich numery:

1 wiersz – ilośd sąsiadów 1 cząstki

…

…

n-ty wiersz – ilośd sąsiadów n-tej cząstki

- można wykonad tysiąc całkowao, lista jest stała,

- błąd będzie się sam kompresował,

- bezpieczne obliczenia dla: promienia wewnętrznego = 15 Å, promieo zewnętrzny (kuli Verleta) = 17 Å
(różnica tych promieni (2 Å) daje tzw. „skórę”, w obrębie której nie musimy uaktualniad co krok listy
najbliższych sąsiadów -> w efekcie szybsze działanie).

[e-mail]

Ad. 6.

Metoda cel połączonych (Hockney’a)

Podstawową komórkę dzielimy na podkomórki (subcells). W każdym kroku czasowym tworzona jest lista
cząstek znajdujących się w każdej podkomórce. Aby uzyskad listę najbliższych sąsiadów operujemy na
podkomórkach, a nie na samych cząsteczkach. Dzięki temu, że układ podkomórek jest stały i znany,
metoda ta jest szybsza od metody Verleta(~N). Wadą tej metody jest trudnośd efektywnego
zaimplementowania metody na procesorach wektorowych (stosowanych w superkomputerach).

- stosujemy gdy mamy dużo cząstek,

- można liczyd w każdym kroku czasowym,

- wektor „mówi”, gdzie kto „mieszka”, więc eliminujemy porównanie „każdy z każdym”.

Ad. 7.

Poprawki dalekozasięgowe – poprawki wynikające z odcięcia potencjału w punkcie r

cut

.

𝑈
𝑁

=

1
2

𝜌 4𝜋𝑟

2

𝑔 𝑟 𝑈 𝑟 𝑑𝑟 =

1
2

𝜌 4𝜋𝑟

2

𝑔 𝑟 𝑈 𝑟 𝑑𝑟 +

1
2

𝜌 4𝜋𝑟

2

𝑔 𝑟 𝑈 𝑟 𝑑𝑟

∞

𝑟

𝑐𝑢𝑡

𝑟

𝑐𝑢𝑡

0

∞

0

z symulacji(<U>

md

)

poprawka dalekozasięgowa

U = U

poprawki

+ <U>

md

𝑔 𝑟 =

1
𝜌

∙

𝑛(𝑟)∆𝑟

4𝜋𝑟

2

∆𝑟

𝑛(𝑟)∆𝑟 - liczba sąsiadów w promieniu (𝑟, 𝑟 + ∆𝑟)

𝜌 =

𝑁

𝑉

Wielkości fizyczne: energia potencjalna, energia kinetyczna.

[e-mail]

Ad. 8.

Radialna funkcja rozkładu to prawdopodobieostwo znalezienia atomu z przedziału (r, r + dr) od
dowolnego atomu odniesienia (wybranego za początek ukł. wspólrz.).

RDF(r) = J(r)=

4𝜋𝜉

2

𝜌(𝜉

𝑟+𝑑𝑟

𝑟

)𝑑𝜉,

gdzie 𝜌(𝜉) jest liczbą atomów w jednostkowej odległości 𝜉 od atomu odniesieni.

Opisuje średnie otoczenie dowolnego atomu rozważanego ciała. Średnia otoczeo wszystkich atomów
ciała.

- otrzymujemy przez rozkład średniej odległości do sąsiadów
- używa się do pierwszych powłok koordynacyjnych,
- histogramujemy odległośd (podaje się zasięg i rozdzielczośd powłoki 0,01 lub 0,02 Å),
- liczona w trakcie samego próbkowania, podajemy zasięg i podziałkę,
- poszczególne dane możemy zapisywad do rożnych plików bo atomy mamy ponumerowane i wiemy
który jest który.


(wykres RDF od r)


Służy do:
- wyznaczania odległości pomiędzy sąsiadami,
- możemy obliczyd liczbę koordynacyjną

Radialne i kątowe funkcje rozkładu służą do identyfikacji najbliższych sąsiadów w przestrzeni. Gdy
struktura jest krystaliczna to odległości miedzy atomami są stałe i otrzymamy w tych miejscach maksima
funkcji rozkładu. Tej metody używa się także w przypadku ciał amorficznych, częściowo wykazujących
strukturę krystaliczną. Można w ten sposób określid w jakim stopniu badany materiał jest amorficzny czy
krystaliczny.

Metoda pierścieni – szukamy zamkniętych łaocuchów (pierścieni) m-kationów i n-anionów, rozkład
długości mówi o stopniu otwartości, układy np. krótkie -> dużo naprężeo i defektów.

(rysunek)

Metoda wielościanów Wonny’a – wielościany pozwalają identyfikowad obszar o danym typie struktury.
Jeśli tworzą się dyslokacje -> obszar atomów krawędziowych zmniejsza się lub zwiększa -> wówczas
szukamy pary (ten atom krawędziowy i leżący we wnętrzu). Stosuje się do struktur gęsto (ciasno)
upakowanych.

[e-mail]

Ad. 9.

Termalizacja – zmiana temperatury w celu zrównoważenia układu, doprowadzenie do równowagi.

Próbkowanie – sposób na obliczenie całki (wykorzystywany w metodach stochastycznych) 𝑓(𝑥) =

𝐼

𝑏−𝑎

Próbkowanie bezpośrednie: 𝐼 =

𝑓 𝑥 𝑑𝑥

𝑏

𝑎

≈

𝑏−𝑎

𝑛

𝑓(𝑥

𝑖

)

𝑛

𝑖

Próbkowanie ważone: 𝐼 =

𝑓 𝑥
𝜇 𝑥

𝜇(𝑥)

𝑏

𝑎

𝑑𝑥 ≈

𝑏−𝑎

𝑛

𝑓 𝑥

𝑖

𝜇 (𝑥

𝑖

)

𝑛

𝑖

𝜇(𝑥) ≈ 𝑓(𝑥)

Próbkowanie jest niezbędne do reprezentacji otrzymanych wartości (z symulacji).

Skalowanie prędkości – stosujemy skalowanie prędkości w celu uzyskania energii cząstek bliskich
energiom zadanym.

- scałkowad równania ruchu dla pewnej liczby kroków czasowych,
- obliczyd energię kinetyczną i potencjalną,
- jeżeli energia nie jest równa energii żądanej, to przeskalowad prędkości,
- powtarzad od początku dopóki układ nie znajdzie się w równowadze.

W praktyce poza doborem algorytmu jest kilka ważnych warunków, których spełnienie jest wymagane
dla powodzenia symulacji. Po pierwsze początkowy model może zawierad zdeformowane fragmenty,
które podczas symulacji będą podlegad gwałtownej relaksacji, co może zniszczyd lokalną strukturę
układu. Dlatego przed symulacją MD konieczne jest zminimalizowanie energii układu. Dodatkowo należy
poddad układ procesowi termalizacji, w którym temperatura będzie stopniowo podnoszona, a prędkości
początkowe losowane wielokrotnie. W ten sposób wszystkie ewentualne naprężenia, będą mogły ulec
rozprężeniu bez nadawania fragmentom układu dużych prędkości względnych, co mogłoby rozbid
lokalną strukturę.

[e-mail]

Ad. 10.

Lattice start – układ periodyczny zapewnia odpowiednią odległośd między środkami mas by nie uzyskały
gigantycznych energii. Stosuje się dla dużych nieregularnych molekuł. Rozpoczyna się tam symulacje ciał
stałych.

Skew start – start skośny, brak regularności, ale zapewniona odległośd pomiędzy środkami mas. Sposób
bardziej zbliżony do rzeczywistości. Jeśli obecny jest więcej niż jeden typ molekuł to ich pozycje są
losowe. Orientacja molekuł jest losowa. Metoda ta działa dobrze dla małych molekuł o regularnych
kształtach. Używa się jej do rozpoczynania symulacji stanów ciekłych.

[e-mail]

Ad. 11.

MSD – średnie przemieszczenie kwadratowe

𝑀𝑆𝐷 =

1

𝑁

𝑟

𝑖

𝑡 − 𝑟

𝑖

0

2

𝑁

𝑖=1

(wykres)

- służy do analizy dyfuzji (jako pastprocessing)

𝐷 =

𝑀𝑆𝐷

6𝑡

6t, bo 6 kierunków (4 w jednym kierunku, ½- jeśli symulacja płaska)

𝐷 𝑡 = 𝐷

0

𝑒

−

𝐸𝑐

𝑘𝑇

- MSD służy do zaobserwowania zmian termicznych

- służy do zaobserwowania przejśd fazowych

(wykres przejść fazowych)

Im wyższa temperatura tym wyższy współczynnik dyfuzji. Współczynnik dyfuzji niesie informacje czy
układ znajduje się w równowadze, a jego zależnośd od temperatury informuje o energii cząstek układu.

Zależnośd współczynnika dyfuzji od temperatury mówi o przejściach fazowych. Gdy jest duży skok
wartości tego współczynnika to znaczy, że jest jakieś przejście fazowe.

[e-mail]

Ad. 12.

Suma Ewalda – każdy ładunek punktowy jest zamieniany na rozkład Gaussa skupiony na tym ładunku.

Metoda Ewalda

W symulacji nie mamy symetrii translacyjnej. Przyjmujemy całe pudło symulacyjne jako jedną komórkę
elementarną.

- próbujemy sumowad ale zbieżnośd jest bardzo powolna -> nie do wyliczenia

(wykres 1)

Odległośd x powinna byd jak największa.

(wykres 2)

Kłopoty:

- 98-99% czasu trwa liczenie sił,

- jeśli dodamy siły kulombowskie to czas liczenia sił wynosi ok. 99,9 %,

- jeżeli mamy oddziaływania kulombowskie to pojawia się stała Madelunga, wkłady do niej są różnych
znaków, szereg jest rozbieżny, suma jest zbieżna, ale moduł jest rozbieżny.

[e-mail]

Ad. 13.

Model powłokowy Kochrana – cząstka punktowa ma swoją masę i wyobrażamy sobie, że taki obiekt
punktowy jest zamknięty w pewnej sferze- powłoce, gdzie jądro to ładunek q1, sferyczna powłoka – q2.

q = q1 + q2

Powłoka jest przytwierdzona na stałe do środka jądra. Powłoka nie ma masy. Cała masa jest skupiona w
jądrze. Występują oddziaływania: jądro-jądro, jądro-powłoka, powłoka-jądro, powłoka-powłoka.

- w pierwszym kroku całkowania liczymy zmianę położeo jąder, ale położenia jąder nie są
równowagowe, bo zostały tam, gdzie były (nie mają masy),

- liczymy położenie równowagowe powłok przy nieruchomych jądrach,

- mamy siły elektrostatyczne i szukamy położenia, które minimalizuje siły kulombowskie,

- wyliczamy siły działające na jądro i znowu ruszamy jadra, itd.

Model Finchama – powłoka przytwierdzona jest sprężynką o stałej k do jądra. Powłoka musi posiadad
małą masę, aby można było obliczyd jej równania ruchu. Mamy dwa razy więcej ciał.

Model powłokowy stosujemy przy cząstkach o małej liczbie atomowej, gdzie masa elektronów ma
wpływ na ruch molekuły.

Model powłokowy w przeciwieostwie do modelu cząstek punktowych uwzględnia masę i ładunek
powłoki. Modele te służą do symulacji cząsteczek o stałym momencie dipolowym (np. woda). Trzeba
także uwzględnid oddziaływanie między powłokami i jądrami.

Ad. 14.

Symulacja układów molekularnych różni się od symulacji cząstek punktowych tym, że zamiast cząstek
mamy molekuły i jeżeli chcemy uzyskad wysoką dokładnośd powinniśmy symulowad zjawiska zachodzące
wewnątrz tej molekuły.

Metoda kwaternionowa

+ nie występuje dzielenie przez 0 (brak sin i cos funkcji),

- równanie nieliniowe, trzeba rozwiązywad iteracyjnie,

Kwaterniony – skooczony ciąg liczb, 4 elementowy (1 rzeczywista i 3 urojone) długości 1 (obroty w
przestrzeni).

Obrót bryły sztywnej można przedstawid za pomocą kwaternionów o długości 1. Przy pomocy
kwaternionów możemy opisad układ mniejszą ilością równao.
Wyznaczamy położenie i prędkości cząsteczek za pomocą tych samych równao, jednak dodatkowo w
modelu molekuł deformowalnych uwzględniamy wewnętrzne wibracje molekuły.
W symulacjach układów molekularnych dodatkowo trzeba uwzględniad drgania zachodzące w
molekułach.

W bryłach sztywnych moment bezwładności się nie zmienia i nie trzeba poprawiad po każdym kroku
całkowania.

[e-mail]

Ad. 15.
Metoda układów skrępowanych więzami (CSM)
- żądamy by całkowita energia kinetyczna była stała,
- upraszcza układ równao eliminując zbędne zmienne niezależne,
- zmniejszenie ilości stopni swobody przyspiesza obliczenia numeryczne.

równanie więzów: 𝑇~𝐸

𝑘

=

1
2

𝑚𝑣

𝑖

2

𝑖

(suma kwadratów prędkości jest stała)

(NVT) – metoda termostatowania
równanie ruchu = równanie Newtona + a ∙ równanie więzów, a- współcz. Lagrange’a

siła

tarcie

- jeśli mamy zachowane równanie więzów mamy stałą temperaturę
(NVp) – metoda barostatowania

chwilowe ciśnienie

𝑝 =

𝐸

𝑘

𝑘𝑇

Mimo modyfikacji T czy p, wartości te są zachowane ściśle.

Ad. 16.
Metoda układów rozszerzonych (ESM) – polega na rozszerzeniu hamiltonianu o dodatkowe czynniki
związane z warunkami fizycznymi ruchy typu ciśnienie, temperatura, czy współczynnik tarcia. Metoda
intuicyjna.

Metoda Andersona (NpH)
- dodano dodatkowy stopieo swobody -> objętośd pudła, trzeba wprowadzid abstrakcyjną masę (ścianki
pudła),
- normalną zmienna jest objętośd pudła, ale mamy jeszcze ciśnienie wewnętrzne i zewnętrzne ->
badamy, które jest większe:

𝑚

𝑡ł𝑜𝑘𝑎

∙ 𝑉 = 𝑝

𝑤𝑒𝑤𝑛𝑒𝑡𝑟𝑧𝑛𝑒

− 𝑝

𝑧𝑒𝑤𝑛 ę𝑡𝑟𝑧𝑛𝑒

𝑝

𝑤𝑒𝑤

> 𝑝

𝑧𝑒𝑤

to rozprężanie, 𝑝

𝑤𝑒𝑤

< 𝑝

𝑧𝑒𝑤

to sprężanie

bok pudła jest funkcją czasu 𝐿 𝑡 =

𝑉(𝑡)

3

<-ściana

- masa tłoka musi byd taka, by ciśnienie nie fluktuowało, czyli, żeby objętośd się nie zmieniała.
Zachowana jest energia całkowita działająca w tłoku i energia związana z całym układem.

Metoda Nose (NVT) Termostat Nosego (1984)
- należy dodad dodatkowy stopieo swobody,
- temperatura będzie całką ruchu,
- układ ma tendencję do ogrzewania

𝑇~ <

𝑚𝑉

2

2

> ~ < 𝑉

2

>

𝑑𝑥

𝑑𝑡

← za duże, to wydłużamy dx,

𝑑𝑥

𝑑𝑡

← za małe, to zmniejszamy dx

wprowadzono czas wirtualny:
𝑑𝑡

′

= 𝑠 𝑡 𝑑𝑡, gdzie s(t) to dodatkowy element swobody związany z masą termiczną

𝑚

𝑡𝑒𝑟𝑚

∙ 𝑠 = 𝑇

𝑤𝑒𝑤

− 𝑇

𝑧𝑒𝑤

NVT – przy stałej l. cząstek, stałej V otrzymujemy stałą T

[e-mail]

Ad. 17.

Monte Carlo – uśredniamy po konfiguracji (po czasie)
- generujemy konfiguracje początkową losową (nie martwię się o temperaturę, zbyt dużymi siłami),
- losowo wybieram 1 cząstkę,
- losowo przesuwam ją w przestrzeni,
- ponownie obliczam jej E

p

, jeśli jest mniejsza od poprzedniej to ją akceptuję jako to

prawdopodobieostwo. Jeśli jest większa to mogę zaakceptowad z prawdopodobieostwem

Boltzmanowskim 𝑒

−

∆𝑈
𝑘𝑇

,

- powtarzamy to z innymi, na początku jest termalizacja układu, przesuwając cząstkę mamy szanse
obniżyd jej temperaturę obniżenia E

p

,

- zaczynamy próbkowanie, gdzie konfiguracje znajdujemy z czynnikiem Boltzmanowskim, średnia po tym
daje nam szukane wielkości.
NpT
- losuję, przesuwam,
- losuje bok pudła,
- liczę E

c

niższą akceptuję, wyższą jako czynnik Boltzmanowski,

- porównujemy energię swobodną.

Zmienna liczba cząstek
- losuję liczbę 1 (przesuo), 2 (wsadź), 3 (usuo) cząstek,
- liczę przyrost E

p

,

- przesunięcia.