Chemija · Pamoka

Molekulinė dinamika, Monte Karlo metodas ir būsenų atranka

Nagrinėkite laikinį integravimą, pusiausvyrinių būsenų atranką, būsenų svorius ir koreliuotus stebinius pagal nurodytus modelius.

Grįžti į chemijos planą

Baigę pamoką, galėsite

Atlikti vieną greičių Verlet žingsnį ir perskaičiuoti modelio laiką.
Į Metropolis–Hastings priėmimo skaičiavimą įtraukti atvirkštinius siūlymus.
Sujungti mikroskopines energijas ir degeneracijas jų nesuskaičiuojant du kartus.
Apskaičiuoti baigtinės imties vidurkio dispersiją pagal nurodytą kovariacijos modelį.

Molekulinė dinamika seka nurodytą kitimo laike modelį

Klasikiniame molekulinės dinamikos modelyje padėtys ir greičiai kinta pagal nurodytas jėgas ir integravimo taisykles. Pradinės sąlygos, ribos ir bet koks temperatūros ar slėgio valdymas yra modeliavimo dalis. Laiko žingsnis turi tikti atvaizduojamiems judesiams ir skaitiniam metodui.

Stabiliai atrodanti trajektorija savaime nepatvirtina tikslios dinamikos ar reprezentatyvios pusiausvyrinių būsenų atrankos. Tvermę ar kitą tikėtiną elgseną reikia tikrinti pagal tikrąjį ansamblį ir algoritmą.

Laiko žingsnis ir atrankos žingsnis reiškia skirtingus dalykus

Greičių Verlet metodas atnaujina padėtį pagal esamą greitį ir pagreitį, tada greičiui apskaičiuoti naudoja abiejų žingsnio galų pagreičius. Jo laiko prieaugis priklauso dinaminiam modeliui. Pakeitus prieaugį, pasikeičia skaitinė paklaida ir laikas, atitinkantis fiksuotą žingsnių skaičių. Toliau pateiktas vienas osciliatoriaus žingsnis nepatvirtina ilgų trajektorijų patikimumo.

Monte Karlo siūlymas, priešingai, keičia būseną pagal algoritmą. Dviejų būsenų pavyzdyje siūlymai yra asimetriniai, todėl atvirkštinio siūlymo tikimybė įtraukiama į priėmimo santykį. Detalusis balansas patikrina, kad nurodytas tikslinis skirstinys yra stacionarus. Tinkamas svarbių būsenų pasiekiamumas ir pakankama atranka išlieka atskirais reikalavimais.

Atmetus siūlymą, išlaikoma esama būsena. Pašalinus tokius siūlymus iš vidurkio, pakeičiama atrankos apskaita. Šiame pavyzdyje nepateikta, kaip siūlymų skaičių paversti fizikiniu molekulių laiku.

Monte Karlo metodu galima atrinkti būsenas neatvaizduojant fizikinio laiko

Monte Karlo procedūra siūlo pokyčius ir juos priima arba atmeta pagal nurodytą taisyklę. Tinkama pusiausvyrinė schema atitinkamomis sąlygomis siekia apibrėžto skirstinio. Siūlymo tikimybės ir priėmimo taisyklės turi derėti tarpusavyje.

Monte Karlo žingsnių skaičius paprastai skiriasi nuo praėjusio fizikinio laiko. Atmestas siūlymas yra atrankos algoritmo dalis, o ne įrodymas, kad tikras molekulinis įvykis neįvyko. Patikrinkite, ar procedūra pasiekia svarbias būsenas, užuot apsiribojusi viena sritimi.

Užpildai ir įvertintas vidurkis turi turėti teisingus svorius

Energijos eilutė gali reikšti vieną mikroskopinę būseną arba jų grupę. Nurodytame dviejų grupių pavyzdyje keturios lygiavertės B būsenos suteikia keturis vienodus Bolcmano dėmenis. Kai šie dėmenys jau įtraukti į makrobūsenos laisvąją energiją, dar kartą dauginant iš tos pačios degeneracijos jie būtų suskaičiuoti du kartus.

Stebimojo dydžio vidurkiui tuomet reikia įvertinti atrankos neapibrėžtį. Teigiama gretimų stebinių koreliacija padidina vidurkio dispersiją, palyginti su nepriklausomais stebiniais, kurių atskirų verčių dispersija tokia pati. Toliau pateikta baigtinė suma įtraukia tikrąjį stebinių porų skaičių esant kiekvienam atstumui.

Pagal dispersiją lygiavertis efektyvusis skaičius apibendrina vieną nurodytą kovariacijos dėsnį. Net ilga seka gali neapimti svarbios būsenos, o dispersijos formulė, paremta stacionarumo prielaida, neįrodo stacionarumo. Išsaugokite pradines sąlygas, pusiausvyros nusistovėjimo pasirinkimus, ribas, stebinius ir transformacijas, taikytus tikram modeliavimui įvertinti.

Koreliacija padidina vidurkio standartinę paklaidą

Standartinė paklaida / kJ mol⁻¹

Stebinių skaičius N

Nurodyti koreliuoti stebiniai
Nepriklausomas palyginimo atvejis
— Baigtinės imties kovariacijos modelio rezultatas

Pagrindiniai duomenys

Stebinių skaičius N	Standartinė paklaida / kJ mol⁻¹
20	0.54884552846152 · Nurodyti koreliuoti stebiniai
80	0.287174838304302 · Nurodyti koreliuoti stebiniai
160	0.20452782855411442 · Nurodyti koreliuoti stebiniai
20	0.2683281572999747 · Nepriklausomas palyginimo atvejis
80	0.13416407864998736 · Nepriklausomas palyginimo atvejis
160	0.09486832980505137 · Nepriklausomas palyginimo atvejis

Originalus tikslus stacionarios kovariacijos modelis su σ=1.20 kJ mol⁻¹ ir ρ=0.650. Apskritimai ir brūkšninė kreivė atitinka baigtinės koreliuotų stebinių imties vidurkio dispersiją; kvadratai atitinka nepriklausomus stebinius su ta pačia atskirų verčių dispersija. Neatvaizduojama jokia tikra trajektorija ar intervalo aprėptis.

Ansamblio vidurkiui reikia teisingų būsenų svorių

Užpildams apskaičiuoti reikia būsenų energijų arba tinkamų laisvųjų energijų, degeneracijų ir temperatūros susitarimo. Didesnės energijos būsena vis tiek gali prisidėti, o kelios būsenos kartu gali nusverti vieną mažesnės energijos būseną.

Atskirkite mikroskopinę energiją nuo makroskopinės laisvosios energijos, į kurią jau įtrauktas būsenos vidinis kartotinumas. Pakartotinai įskaičiavus degeneraciją, šis indėlis gali būti įtrauktas du kartus. Nurodykite, ką reiškia kiekviena modelio eilutė.

Daug išsaugotų kadrų nebūtinai reiškia daug nepriklausomų imties elementų

Gretimi trajektorijos kadrai gali būti koreliuoti. Ilgesnė atranka ar daugiau išsaugotų kadrų gali apibūdinti seką, tačiau nepriklausomos informacijos apie pusiausvyrinį vidurkį kiekis nebūtinai padidėja proporcingai. Reti perėjimai ir nesusisiekiančios sritys gali riboti aprėptį.

Numatytam įverčiui taikykite tinkamus pusiausvyros nusistovėjimo, konvergavimo, koreliacijos ir kartotinių bandymų patikrinimus. Toliau pateiktuose originaliuose atrankos pavyzdžiuose naudojami išgalvoti skirstiniai ir stebiniai; jie neįrodo tikro molekulinio modeliavimo konvergavimo.

Išnagrinėti pavyzdžiai

Vienas integravimo žingsnis ir jo laiko apskaita

Simboliniam neutraliam osciliatoriui naudokite vienmatį greičių Verlet metodą: m=4.00×10⁻²⁶ kg, k=8.00 N m⁻¹, x₀=1.00×10⁻¹⁰ m, v₀=−120 m s⁻¹ ir Δt=2.00 fs. Nėra nei termostato, nei barostato. Apskaičiuokite vieną žingsnį ir laiką, priskirtą 60,000 tokių žingsnių.

a₀=−kx₀/m=−2.00×10¹⁶ m s⁻². x₁=x₀+Δtv₀+Δt²a₀/2=9.972×10⁻¹¹ m=0.09972 nm.
a₁=−kx₁/m. v₁=v₀+Δt(a₀+a₁)/2=−159.944 m s⁻¹. Pradinė bendroji energija mv₀²/2+kx₀²/2 yra 4.028800×10⁻²⁰ J; energija po vieno žingsnio yra 4.028795526272×10⁻²⁰ J.
Šio žingsnio pokytis yra −4.473728×10⁻²⁶ J. Jis patikrina skaitinį rezultatą šiame nurodytame modelyje; jis neįrodo nei ilgalaikio stabilumo, nei tikslios molekulinės dinamikos.
60,000×2.00 fs=120 ps modelio laiko. Šis skaičiavimas nepriklauso nuo kadrų išsaugojimo dažnio ir nepriskiria fizikinio laiko nesusijusiam Monte Karlo siūlymų skaičiui.

x₁=0.09972 nm, v₁=−159.944 m s⁻¹; šiame kitimo laike modelyje 60,000 žingsnių atitinka 120 ps.

Asimetriniam siūlymui reikia atvirkštinės tikimybės

Dviem simbolinėms tų pačių medžiagų konfigūracijoms duotos tikslinės tikimybės π_A=0.750 ir π_B=0.250. Siūlykite A→B su tikimybe 0.600 ir B→A su 0.200; likusi siūlymo tikimybės dalis palieka kiekvieną būseną tokią pačią. Taikykite Metropolis–Hastings metodą.

Priėmimo tikimybė α_AB=min[1,π_Bq_BA/(π_Aq_AB)]=min[1,0.250×0.200/(0.750×0.600)]=1/9. α_BA=1.
Tikrosios nediagonalinės perėjimo tikimybės yra P_AB=0.600/9=0.0666667 ir P_BA=0.200. Pasilikimo tikimybės yra 0.9333333 ir 0.800, įskaitant atmestus siūlymus.
Pusiausvyriniai srautai yra π_AP_AB=0.0500 ir π_BP_BA=0.0500. Jų lygybė patikrina detalųjį balansą ir nurodytą invariantinį tikslinį skirstinį šioje susietoje baigtinėje grandinėje.
Neatsižvelgiant į nevienodus atvirkštinius siūlymus, būtų gaunamas kitas tikslinis skirstinys. Atmetimai išlaiko esamą būseną ir turi būti įtraukti į atrankos vidurkį. Šis pavyzdys nepriskiria fizikinio laiko ir nepatvirtina tikros molekulinės sistemos konvergavimo.

Priėmimo tikimybės 1/9 ir vienetas; abu tikslinėmis tikimybėmis pasverti perėjimų srautai yra 0.0500 per žingsnį.

Kelios didesnės energijos būsenos gali nusverti vieną mažesnės energijos būseną

Esant 310.0 K, makrobūsenai A priskirkite vieną mikroskopinę būseną, kurios energija 0 kJ mol⁻¹, o makrobūsenai B – keturias būsenas, kurių kiekvienos energija 3.40 kJ mol⁻¹. Šioje supaprastintoje apskaitoje kitų būsenų nėra. Raskite B užpildą.

Naudokite molinę dujų konstantą R=N_Ak_B=8.314462618… J mol⁻¹ K⁻¹, kur N_A yra Avogadro konstanta, o k_B – Bolcmano konstanta. Santykiniai svoriai yra w_A=1 ir w_B=4exp[−3.40/(0.008314462618…×310.0)]=1.069486053.
P_B=w_B/(1+w_B)=0.516788239 ir P_A=0.483211761. Vidutinė nurodyta molinė energija yra 3.40P_B=1.75708001 kJ mol⁻¹.
Šioms būsenų grupėms ΔG_macro=3.40−RTln4=−0.173150719 kJ mol⁻¹. Kiekviena B mikrobūsena turi didesnę energiją, tačiau visos B grupės užpildas yra šiek tiek didesnis.
Jei įvesties laisvoji energija jau apimtų šį vidinį kartotinumą, pakartotinai dauginant iš keturių jis būtų įskaičiuotas dukart. Šios nurodytos būsenos nėra išmatuoti molekuliniai konformerai.

P_B=0.516788239; keturgubas kartotinumas pakeičia makrobūsenų palyginimą.

Koreliuoti kadrai pakeičia vidurkio neapibrėžtį

Duotas stacionarus stebimasis dydis su N=80 vienodais intervalais gautų stebinių, atskirų verčių standartiniu nuokrypiu σ=1.20 kJ mol⁻¹ ir kovariacija Cov(Y_i,Y_j)=σ²ρ^|i−j|, ρ=0.650. Apskaičiuokite imties vidurkio standartinę paklaidą pagal šį tiksliai nurodytą kovariacijos dėsnį.

N stebinių vidurkiui Var(mean)=σ²/N[1+2Σ_(k=1)^(N−1)(1−k/N)ρ^k]. Baigtinės sumos daugiklis yra 4.581632653.
SE=1.20√(4.581632653/80)=0.287174838 kJ mol⁻¹. Jei 80 stebinių būtų nepriklausomi, SE būtų 1.20/√80=0.134164079 kJ mol⁻¹.
Apibrėžkite pagal dispersiją lygiavertį efektyvųjį skaičių N_eff=N/factor=17.4610245. Jis yra trupmeninis ir aprašo šį vidurkio dispersijų palyginimą; tai nėra garantuotas nepriklausomų kadrų skaičius.
Kovariacija yra nurodyta, o ne įvertinta. Tiksli jos dispersija nepatvirtina pusiausvyros nusistovėjimo, retų būsenų aprėpties, intervalo aprėpties ar tikros trajektorijos konvergavimo.

Modeliu pagrįsta SE=0.287174838 kJ mol⁻¹ ir pagal dispersiją lygiavertis N_eff=17.4610245.

Pabandykite patys

Naudokite dt = 4.00 fs su ta pačia pradine būsena ir 60000 žingsnių praėjusiam laikui apskaičiuoti. Raskite būseną po vieno žingsnio ir praėjusį laiką.
Rodyti atsakymą ir paaiškinimą
x₁=0.09936 nm, v₁=−199.744 m s⁻¹; energijos pokytis −4.0828928×10⁻²⁵ J; 60,000 žingsnių duoda 240 ps.
Integravimo būsena ir skaitinė energijos paklaida kinta keičiant dt; praėjęs modelio laikas yra žingsnių skaičius, padaugintas iš dt.
Pakeiskite tik siūlymo B -> A tikimybę į 0.500, kai selfB = 0.500. Raskite abi priėmimo tikimybes, nediagonalinius perėjimus ir pusiausvyrinius srautus.
Rodyti atsakymą ir paaiškinimą
α_AB=5/18, α_BA=1; P_AB=1/6, P_BA=1/2; abu srautai 0.125.
Atvirkštinio siūlymo tikimybė priklauso priėmimo santykiui. Pakoregavus priėmimo tikimybę, detalusis balansas išlieka.
Pakeiskite tik B degeneraciją iš 4 į 2. Raskite P_B ir makrobūsenų laisvųjų energijų skirtumą esant 310.0 K.
Rodyti atsakymą ir paaiškinimą
P_B=0.348425122; ΔG_macro=+1.61342464 kJ mol⁻¹.
Kiekvienos mikroskopinės B būsenos energija išlieka ta pati, o būsenų skaičius pasikeičia.
Padvigubinkite tik stebinių skaičių iki 160, palikdami tą patį nurodytą kovariacijos dėsnį ir atskirų verčių sigma. Raskite N_eff ir standartinę paklaidą.
Rodyti atsakymą ir paaiškinimą
N_eff=34.4237102; SE=0.204527829 kJ mol⁻¹.
Didesnis stebinių skaičius sumažina vidurkio dispersiją, tačiau koreliuoti stebiniai tebėra mažiau informatyvūs nei toks pat nepriklausomai atrinktų verčių skaičius.

Ką verta įsidėmėti

Integravimo laikas nustatomas pagal nurodytą dinaminį laiko žingsnį.
Pusiausvyrinei atrankai reikia suderinamų siūlymo ir priėmimo tikimybių.
Kartotinumas turi būti įtrauktas vieną kartą pagal apibrėžtą būsenų pagrindą.
Koreliuoti stebiniai ir nevisiška būsenų aprėptis riboja skirtingus modeliavimo aspektus.

Šaltiniai

Mokomoji medžiaga ir žinynai, kuriais remiasi pamokos paaiškinimai ir pavyzdžiai.