Molekulární dynamika - co to je? technika počítačové simulace, kde je časový vývoj množiny interagujících atomů popsán integrací jejich pohybových rov

PARALELNÍ VÝPOČTY V MOLEKULÁRNÍ DYNAMICE

Vladimír Pelikán, Petr Hora Ústav termomechaniky AV ČR, v.v.i.

Paralelní výpočty v molekulární dynamice

Molekulární dynamika - co to je?

• technika počítačové simulace, kde je časový vývoj množiny interagujících atomů popsán integrací jejich pohybových rovnic. • vychází ze zákonů klasické mechaniky, především Newtonových pohybových rovnic. • jedná se o deterministickou techniku (rozdíl oproti metodě Monte Carlo). Počáteční množina pozic a rychlostí zcela určuje následný časový vývoj. • jedna z metod statistické mechaniky.

1


Historie molekulární dynamiky 1957 Alder a Wainwright První článek referující o MD-simulaci. Fázový diagram systému tuhých koulí. Počítač UNIVAC a IBM 704. 1960 Gibson, Goland, Milgram a Vineyard První použití spojitého potenciálu a časové integrační metody konečných diferencí. Počítač IBM 704, 500 atomů, jeden časový krok trval asi minutu. 1964 Rahman Vlastnosti tekutého argonu s použitím Lennardova-Jonesova potenciálu. Počítač CDC 3600, 864 atomů. 1967 Verlet Fázový diagram argonu s pomocí Lennardovo-Jonesova potenciálu. Poprvé uveden Verletův seznam sousedů. Mimoto byl použit Verletův časový integrační algoritmus. 2


Aplikace molekulární dynamiky Tekutiny. Dostupnost nových realistických interakčních modelů umožňuje studovat nové systémy, elementární i vícesložkové. Pomocí nerovnovážných technik jsou vyšetřovány transportní jevy (např. viskózita) a tepelné toky. Defekty. Pozornost se přesouvá z bodových defektů (vakance, intersticiální porucha) k lineárním (dislokace) a plošným defektům (hranice zrn, vrstevné chyby). Reálnějších modelů je dosaženo díky lepším potenciálům. Lomy. Lomový proces může v závislosti na příslušných parametrech nastat různými způsoby a probíhat různou rychlostí. Technická závažnost je zřejmá a simulace poskytuje proniknutí k podstatě. Povrchy. Simulace stále hraje velkou roli při pochopení jevů jako např. rekonstrukce povrchu, povrchové tavení, povrchová difúze, zdrsňování, atd. Tyto simulace často vyžadují rozsáhlé vzorky a dlouhé simulační časy.

3


Tření. Vývoj atomového silového mikroskopu (AFM) podnítil výzkumy adheze a tření mezi dvěma pevnými látkami. Souhrn makroskopických znalostí je korigován a rozšiřován o mikroskopické základy. Klastry. Klastry, tj. konglomeráty mnoha atomů, představují most mezi molekulárními systémy a pevnou fází a vykazují pozoruhodné vlastnosti. Klastry kovů jsou díky své úloze katalyzátoru v důležitých chemických reakcích (např. katalyzátory v automobilech) zvláště důležité z technologického hlediska. Biomolekuly. Molekulární dynamika umožňuje studovat dynamiku rozsáhlých makromolekul, včetně biologických systémů jako např. proteiny, nukleové kyseliny (DNA, RNA) a membrány. Dynamické jevy mohou hrát klíčovou úlohu při regulačních procesech, které ovlivňují funkční vlastnosti biomolekul. Simulace léků je běžně používaná ve farmaceutickém průmyslu. Elektronové vlastnosti a dynamika. Vývoj Carovy-Parrinellovy metody, při které jsou síly působící na atomy získány řešením problému elektronové struktury místo meziatomového potenciálu, umožňuje plně studovat vlastnosti materiálů včetně jejich dynamiky (a tedy fázových přechodů a jiných teplotně závislých jevů). 4


Základní aparát molekulární dynamiky • Hlavní součástí simulace je model fyzikálního systému. • Pro simulace MD to znamená volbu potenciálu. Potenciál: Funkce pozic atomů V (r1, . . . , rN ), která představuje potenciální energii systému atomů uspořádaných v dané konfiguraci. Nejjednodušší tvar pro V je součet párových interakcí: V (r 1, . . . , r N ) =

XX i

φ(|r i − r j |).

j>i

• Síly jsou odvozeny jako gradienty potenciálu vzhledem k výchylkám atomů: F i = −∇ri V (r 1, . . . , r N ).

5

φLJ (r) = 4ε

σ 12 r

−

σ 6 r

φLJ (r )


Lennardův-Jonesův potenciál 12-6

r

0 √ 6

2σ

−ε

σ

r

6


Postup simulace 1. Zahájení simulace 2. Řízení simulace 3. Získání výsledků ze simulace

7


1. Zahájení simulace 1. Definovat pomyslný prostor, ve kterém bude simulace probíhat. 2. Přiřadit částicím počáteční pozice a rychlosti. • Vytvořit zcela novou množinu počátečních pozic a rychlostí. • Převzít konečné pozice a rychlosti z předchozího běhu MD simulace.

8


2. Řízení simulace Jak lze přivést systém do požadovaného stavu? • Hustota se řídí volbou objemu V simulačního boxu. • Teplotních změn je obvykle dosaženo přeškálováním rychlostí. • V pokročilejších simulacích se volí tlak a měří se objem boxu. Indikátory stavu systému: • stacionární, tj. kolísající okolo fixní hodnoty, • relaxující, tj. kolísající okolo hodnoty, která se pomalu mění. Stavová změna může být: • vyvolána uměle, např. dojde-li k přenastavení parametrů simulace, • samovolná, např. systém prochází fázovou transformací.

9


3. Získání výsledků ze simulace • prostá inspekce pozic atomů, • jednoduché statistické veličiny, • měření teploty tavení, • prostorové korelace, • dynamická analýza atd.

10


Vícečásticový potenciál G.J.Ackland, D.J.Bacon, A.F.Calder, T.Harry: Computer simulation of point defect properties in dilute Fe-Cu alloy using a many-body interatomic potential. Philosophical Magazine A, 1997, Vol. 75, No. 3, 713–732

Energie souboru N atomů je dána výrazem 

1/2

N N N X X 1 X  V (rij ) − φ (rij ) E= 2 i=1 i6=j =1

j6=i=1

V (rij ) − repulzivní potenciál φ (rij ) − kohezivní potenciál

11


Pohybové rovnice Pro každý atom jsou definovány 3 pohybové rovnice: dsn d2sn mn 2 + C + Rsn = Fsn, dt dt

pro s = x, y, z

kde t je čas, mn je hmotnost n-tého atomu, C je koeficient útlumu, Rsn jsou složky gradientu potenciálu a Fsn jsou složky vnější síly působící na n-tý atom. Koeficient útlumu je různý od nuly pouze v případě relaxace.

12


Metoda centrálních diferencí

dsn 1 ∼ (sn(i + 1) − sn(i − 1)) = dt t=ti 2∆t 2

d sn dt2

t=ti

1 ∼ = 2 (sn (i + 1) − 2sn (i) + sn (i − 1)) ∆t

Dosazením této rovnice do pohybových rovnic dostaneme: 2 ∆t2 sn (i + 1) = sn (i − 1) + (sn (i) − sn (i − 1)) + (Fsn − Rsn) , An Anmn kde

C∆t An = 1 + . 2mn

13


Popis paralelní úlohy Úlohy molekulární dynamiky jsou postaveny na neustálém, ale především mnohonásobném vyčíslování relativně jednoduché rekurentní formule, do které vstupují pouze údaje z nejbližšího okolí právě zpracovávaného atomu. • U bcc železa se jedná např. pouze o 14 nejbližších sousedních atomů. → Takovéto úlohy jsou přímo předurčeny k paralelnímu zpracování.

PC

PC

PC

PC

PC

PC

14


Paralelní program • Zdrojová verze programu je napsána v jazyce Fortran 90. • Programový kód využívá systém MPI + OpenMP. Používané funkce MPI: • MPI INIT, MPI FINALIZE, • MPI COMM RANK, MPI COMM SIZE, • MPI SEND, MPI RECV, MPI BCAST, • MPI ISSEND, MPI IRECV, • MPI WAIT.

15


• Hash: metoda indexace buněk (cell index method) M.P.Allen, D.J.Tildesley: Computer Simulation of Liquids. Oxford University Press, New York, 1987 D. Frenkel, B. Smit: Understanding Molecular Simulations. Academic Press, New York, 1996 • Požadavky na paměť: 6×8 + 2×8 = 64 bytů/atom → 6 GB (4 uzlů 2 GB RAM) • Požadavky na diskový prostor: 6×8 = 48 bytů/atom → 4,5 GB

16


Kde se počítalo?

Praha

Plzenˇ Brno

17


Kde se počítalo?

Praha Plzenˇ Brno

NYMPHA 20 uzlu˚ 2x Quad Core Xeon E5472 3 GHz, 12 MB cache 16 GB 300 GB (15 k rpm, SAS) 1 Gbps Ethernet

17


Kde se počítalo?

Praha Plzenˇ Brno

NYMPHA

TARKIL

20 uzlu˚

28 uzlu˚

2x Quad Core Xeon E5472

2x Quad Core Intel Xeon X5570

3 GHz, 12 MB cache

2,93 GHz

16 GB

24 GB

300 GB (15 k rpm, SAS)

2x 300 GB (15 k rpm, SAS)

1 Gbps Ethernet

1 Gbps Ethernet, Infiniband 4x QDR

17


Kde se počítalo?

Praha Plzenˇ Brno

NYMPHA

TARKIL

SKIRIT

20 uzlu˚

28 uzlu˚

35 uzlu˚



2x Dual Core Xeon 5160

3 GHz, 12 MB cache

2,93 GHz

3 GHz, 4 MB cache

16 GB

24 GB

4 GB



73 GB

1 Gbps Ethernet


1 Gbps Ethernet

17


Kde se počítalo?

Praha Plzenˇ Brno

NYMPHA

TARKIL

SKIRIT

20 uzlu˚

28 uzlu˚

35 uzlu˚



2x Dual Core Xeon 5160

3 GHz, 12 MB cache

2,93 GHz

3 GHz, 4 MB cache

16 GB

24 GB

4 GB



73 GB

1 Gbps Ethernet


1 Gbps Ethernet

17

140

400 350 300 240 160 80

120 100

Ts [s]


Závislost doby výpočtu na velikosti vzorku a na počtu procesorů

80 60 40 20 0 0

4

8

12

16 Np

20

24

28

32

18


Vizualizace výsledků Vlastní nástroj pro vizualizaci v MATLABu. Proč? • Existující programy jsou zbytečně složité. • Výstupy našich programů mají specifický formát. • Cizí program je cizí program. Nejdůležitější možnosti vizualizačního programu: • stanovení výřezu, • stanovení pohledu, • volba uvažované mřížky (základní,centrální nebo obě), • volba konkrétního kroku simulace, • filtrace atomů podle energie.

19


Testovací úlohy Klasický experiment → test přesnosti přístrojů Počítačový experiment → test správnosti a přesnosti algoritmů a programů

Testování paralelního programu pro MD simulaci ve 3D bylo prováděno na krystalu α-železa.

Šlo o následující typy ověřovacích experimentů: • povrchová relaxace, • teplotní roztažnost, • Hookeův zákon.

20


MD simulace trhlin ve 3D

1000 × 50 × 1000 + 999 × 49 × 999 ↓ 99 miliónů atomů Použitý potenciál: vícečásticový potenciál Finnisova–Sinclairova typu Integrace pohybových rovnic: metoda centrálních diferencí Integrační krok: 1 × 10−14 s 21


Popis trhliny

Velikost: 199 × 49 × 1 Metoda: vyjmutí centrálních atomů Realizace: interakce přes počáteční rovinu trhliny jsou zakázány

22

a0

0K 2.331

50

16514

ˇ asovy´ krok C

σ A [GPa]

a0

σ A [GPa]


Zatěžování

2.331

300 K

50

16514

ˇ asovy´ krok C

23


Příklad iniciace křehké trhliny při 0 K • iniciace v časovém kroku 23 340 • ve střední rovině krystalu

24


Časový krok: 23300

25



25



25



25



25



25



25


Šíření vln v krystalu bcc železa

kz /ω qSV

SH qP

kx /ω

26


Rychlý celoplošný tah na vnější stěně krychlového vzorku t = 2 ps

t = 4 ps

t = 6 ps

t = 8 ps

t = 8 ps (detail)

27

40 2 ps 4 ps 6 ps 8 ps

35 30 25 |v | [m/s]


Moduly rychlostí jednotlivých atomů na ose y v simulačních časech 2, 4, 6 a 8 ps

20 4s

4s

4s

15 10 5 0 0

50

100 y

150

200

28


Rychlý lokální tah uprostřed vnější stěny vzorku t = 2 ps

t = 4 ps

t = 6 ps

t = 8 ps

t = 8 ps (detail)

29


Lokální tah uprostřed vnější stěny v čase 10 ps

30

0.06

0.04

0.02 vy [m/s]


Průběh vx vers. vy atomu na povrchu vzorku

0

-0.02

-0.04

-0.06

-0.06

-0.04

-0.02

0 0.02 vx [m/s]

0.04

0.06

31


Rychlá celoplošná exploze v centrální rovině vzorku t = 2 ps

t = 4 ps

t = 6 ps

t = 8 ps

t = 8 ps (detail)

32


Rychlá exploze ve středu krychlového vzorku t = 2 ps

t = 4 ps

t = 6 ps

t = 8 ps

t = 8 ps (detail)

33


Děkuji za pozornost.

34


OBSAH Co to je? Historie Aplikace Místo Aparát Lenardův–Jonesův potenciál Simulace – kroky Vícečásticový potenciál Pohybové rovnice Metoda centrálních diferencí Paralelní úloha Paralelní program Kde se počítalo Vizualizace Testovací úlohy Trhliny Vlny

35

Molekulární dynamika - co to je? technika počítačové simulace, kde je časový vývoj množiny interagujících atomů popsán integrací jejich pohybových rov

Recommend Documents