確定起始構(gòu)型
進(jìn)行分子動(dòng)力學(xué)模擬的第一步是確定起始構(gòu)型,一個(gè)能量較低的起始構(gòu)型是進(jìn)行分子模擬的基礎(chǔ),一般分子的起始構(gòu)型主要來自實(shí)驗(yàn)數(shù)據(jù)或量子化學(xué)計(jì)算。 在確定起始構(gòu)型之后要賦予構(gòu)成分子的各個(gè)原子速度,這一速度是根據(jù)波爾茲曼分布隨機(jī)生成的,由于速度的分布符合波爾茲曼統(tǒng)計(jì),因此在這個(gè)階段,體系的溫度是恒定的。另外,在隨機(jī)生成各個(gè)原子的運(yùn)動(dòng)速度之后須進(jìn)行調(diào)整,使得體系總體在各個(gè)方向上的動(dòng)量之和為零,即保證體系沒有平動(dòng)位移。
進(jìn)入平衡相
由上一步確定的分子組建平衡相,在構(gòu)建平衡相的時(shí)候會(huì)對(duì)構(gòu)型、溫度等參數(shù)加以監(jiān)控。
進(jìn)入生產(chǎn)相
進(jìn)入生產(chǎn)相之后體系中的分子和分子中的原子開始根據(jù)初始速度運(yùn)動(dòng),可以想象其間會(huì)發(fā)生吸引、排斥乃至碰
分子動(dòng)力學(xué)
撞,這時(shí)就根據(jù)牛頓力學(xué)和預(yù)先給定的粒子間相互作用勢(shì)來對(duì)各個(gè)粒子的運(yùn)動(dòng)軌跡進(jìn)行計(jì)算,在這個(gè)過程中,體系總能量不變,但分子內(nèi)部勢(shì)能和動(dòng)能不斷相互轉(zhuǎn)化,從而體系的溫度也不斷變化,在整個(gè)過程中,體系會(huì)遍歷勢(shì)能面上的各個(gè)點(diǎn),計(jì)算的樣本正是在這個(gè)過程中抽取的。步驟
以下是做模擬的一般性步驟,具體的步驟和過程依賴于確定的系統(tǒng)或者是軟件,但這不影響我們把它當(dāng)成一個(gè)入門指南:
1)首先我們需要對(duì)我們所要模擬的系統(tǒng)做一個(gè)簡單的評(píng)估, 三個(gè)問題是我們必須要明確的:
做什么(what to do)為什么做(why to do)怎么做(how to do)
2)選擇合適的模擬工具,大前提是它能夠?qū)崿F(xiàn)你所感興趣的目標(biāo),這需要你非常謹(jǐn)慎的查閱文獻(xiàn),看看別人用這個(gè)工具都做了些什么,有沒有和你相關(guān)的,千萬不要做到一半才發(fā)現(xiàn)原來這個(gè)工具根本就不能實(shí)現(xiàn)你所感興趣的idea,切記! 考慮1:軟件的選擇,這通常和軟件主流使用的力場(chǎng)有關(guān),而軟件本身就具體一定的偏向性,比如說,做蛋白體系,Gromacs,Amber,Namd均可;做DNA, RNA體系,首選肯定是Amber;做界面體系,Dl_POLY比較強(qiáng)大,另外做材料體系,Lammps會(huì)是一個(gè)不錯(cuò)的選擇
考慮2:力場(chǎng)的選擇。力場(chǎng)是來描述體系中最小單元間的相互作用的,是用量化等方法計(jì)算擬合后生成的經(jīng)驗(yàn)式,有人會(huì)嫌它粗糙,但是它確確實(shí)實(shí)給我們模擬大系統(tǒng)提供了可能,只能說關(guān)注的切入點(diǎn)不同罷了。常見的有三類力場(chǎng):全原子力場(chǎng),聯(lián)合力場(chǎng),粗?;?chǎng);當(dāng)然還有所謂第一代,第二代,第三代力場(chǎng)的說法,這里就不一一列舉了。
再次提醒注意:必須選擇適合于我們所關(guān)注體系和我們所感興趣的性質(zhì)及現(xiàn)象的力場(chǎng)。 3)通過實(shí)驗(yàn)數(shù)據(jù)或者是某些工具得到體系內(nèi)的每一個(gè)分子的初始結(jié)構(gòu)坐標(biāo)文件,之后,我們需要按我們的想法把這些分子按照一定的規(guī)則或是隨機(jī)的排列在一起,從而得到整個(gè)系統(tǒng)的初始結(jié)構(gòu),這也是我們模擬的輸入文件。
4)結(jié)構(gòu)輸入文件得到了,我們還需要力場(chǎng)參數(shù)輸入文件,也就是針對(duì)我們系統(tǒng)的力場(chǎng)文件,這通常由所選用的力場(chǎng)決定,比如鍵參數(shù)和非鍵參數(shù)等勢(shì)能函數(shù)的輸入?yún)?shù)。 分子動(dòng)力學(xué)
5)體系的大小通常由你所選用的box大小決定,我們必須對(duì)可行性與合理性做出評(píng)估,從而確定體系的大小,這依賴于具體的體系,這里不細(xì)說了。6)由于初始構(gòu)象可能會(huì)存在兩個(gè)原子挨的太近的情況(稱之為bad contact),所以需要在正式模擬開始的第一步進(jìn)行體系能量最小化,比較常用的能量最小化有兩種,最速下降法和共軛梯度法,最速下降法是快速移除體系內(nèi)應(yīng)力的好方法,但是接近能量極小點(diǎn)時(shí)收斂比較慢,而共軛梯度法在能量極小點(diǎn)附近收斂相對(duì)效率高一些,所有我們一般做能量最小化都是在最速下降法優(yōu)化完之后再用共軛梯度法優(yōu)化,這樣做能有效的保證后續(xù)模擬的進(jìn)行。 7)以平衡態(tài)模擬為例,你需要設(shè)置適當(dāng)?shù)哪M參數(shù),并且保證這些參數(shù)設(shè)置和力場(chǎng)的產(chǎn)生相一致,舉個(gè)簡單的例子,gromos力場(chǎng)是用的范德華勢(shì)雙截?cái)鄟矶ǚ兜氯A參數(shù)的,若你也用gromos力場(chǎng)的話也應(yīng)該用雙截?cái)鄟硖幚矸兜氯A相互作用。常見的模擬思路是,先在NVT下約束住你的溶質(zhì)(劑)做限制性模擬,這是一個(gè)升溫的過程,當(dāng)溫度達(dá)到你的設(shè)定后, 接著做NPT模擬,此過程將調(diào)整體系的壓強(qiáng)進(jìn)而使體系密度收斂。 經(jīng)過一段時(shí)間的平衡模擬,在確定系統(tǒng)弛豫已經(jīng)完全消除之后,就可以開始取數(shù)據(jù)了。如何判斷體系達(dá)到平衡,這個(gè)問題是比較技術(shù)性的問題,簡單的講可以通過以下幾種方式,一,看能量(勢(shì)能,動(dòng)能和總能)是否收斂;二,看系統(tǒng)的壓強(qiáng),密度等等是否收斂;三看系統(tǒng)的RMSD是否達(dá)到你能接受的范圍,等等。
8)運(yùn)行足夠長時(shí)間的模擬以確定我們所感興趣的現(xiàn)象或是性質(zhì)能夠被觀測(cè)到,并且務(wù)必確保此現(xiàn)象出現(xiàn)的可重復(fù)性。
9)數(shù)據(jù)拿到手后,很容易通過一些可視化軟件得到軌跡動(dòng)畫,但這并不能拿來發(fā)文章。真正的工作才剛剛開始——分析數(shù)據(jù),你所感興趣的現(xiàn)象或性質(zhì)只是表面,隱含在它們之中的機(jī)理才是文章中的主題。