分子動(dòng)力學(xué)簡(jiǎn)介(簡(jiǎn)明)

1、 什么是分子模擬 分子模擬是在分子模型的基礎(chǔ)上用計(jì)算機(jī)做實(shí)驗(yàn)，“計(jì)算機(jī)實(shí)驗(yàn)” 通過模擬微觀粒子的運(yùn)動(dòng)來計(jì)算宏觀性質(zhì)溫度壓力黏度傳遞性質(zhì)表面張力分子間的作用模型牛頓力學(xué)量子力學(xué)統(tǒng)計(jì)力學(xué)等 分子模擬具有理論和實(shí)驗(yàn)的雙重性質(zhì) 分子模擬不能完全取代實(shí)驗(yàn)理論實(shí)驗(yàn)?zāi)M理論的正確性模擬參數(shù)的正確性模擬方法的選擇理論的更新與自然界相比的準(zhǔn)確程度尺度（米）時(shí)間（秒）10-910-710-510-310-1510-910-610-31電子模擬（量化計(jì)算，DFT）分子模擬（分子動(dòng)力學(xué)，蒙特卡洛）顆粒方法流體力學(xué)原子結(jié)構(gòu)薛定諤方程模擬電子云能量性質(zhì)，化學(xué)鍵等信息 量子化學(xué)計(jì)算一般處理幾個(gè)到幾十個(gè)原子常見軟件：GAUS

2、SIAN，NWCHEM等 密度泛函（DFT）可以算到上百個(gè)原子上百個(gè)原子常見軟件：VASP 分子水平的模擬 以分子的運(yùn)動(dòng)為主要模擬對(duì)象 采用經(jīng)驗(yàn)性的分子間作用函數(shù)模擬微粒之間的作用 一般情況下不考慮電子轉(zhuǎn)移效應(yīng)，因而不能準(zhǔn)確模擬不能準(zhǔn)確模擬化學(xué)成鍵作用 發(fā)展最早 1950s，Alder，勞倫斯利物默實(shí)驗(yàn)室，分子動(dòng)力學(xué)模擬32個(gè)原子 1950s，Metropolis，洛斯阿洛莫斯實(shí)驗(yàn)室，蒙特卡洛模擬32個(gè)原子 分子級(jí)別的模擬應(yīng)用的領(lǐng)域很廣 廣泛應(yīng)用于化學(xué)，物理，生物，化工，材料，機(jī)械，治藥等領(lǐng)域 簡(jiǎn)單易學(xué) 蒙特卡洛是一種優(yōu)化方法通過蒙特卡洛算法來尋求能量最優(yōu)點(diǎn)隨機(jī)方法 通過系綜平均來求取宏觀性質(zhì)

3、 模擬的是平衡狀態(tài)平衡狀態(tài)，不涉及時(shí)間效應(yīng)（KMC除外） 優(yōu)點(diǎn)是可以跨越時(shí)間因素，缺點(diǎn)是得不到有關(guān)時(shí)間信息的性質(zhì) 可以模擬平衡狀態(tài)，也可以模擬中間狀態(tài) 可以獲得有關(guān)時(shí)間的信息 受時(shí)間的限制，無法模擬緩慢過程分子體系（幾百幾億）求解牛頓運(yùn)動(dòng)方程宏觀性質(zhì)( )dU rm rdr= -& & 同時(shí)考慮原子核的運(yùn)動(dòng)（牛頓力學(xué)）和電子的運(yùn)動(dòng)（量子力學(xué)） 能同時(shí)準(zhǔn)確模擬物理作用和化學(xué)鍵作用 目前來說CPMD可以處理的體系還很?。◣资畟€(gè)原子） 將分子基團(tuán)（幾個(gè)或者幾十上百個(gè)原子）當(dāng)成單個(gè)的微粒來處理 微粒之間的作用也是通過類似于分子動(dòng)力學(xué)的未能函數(shù)來描述 可以模擬更長(zhǎng)的時(shí)間跨度電子原子核原子

4、量子級(jí)別模擬分子級(jí)別模擬CG級(jí)別模擬分子動(dòng)力學(xué)對(duì)勢(shì)能函數(shù)的依賴性：所有從分子動(dòng)力學(xué)計(jì)算出來得到的宏觀性質(zhì)最終都取決于勢(shì)能模型( )dU rm rdr= -& &分子動(dòng)力學(xué)的核心：牛頓運(yùn)動(dòng)方程勢(shì)能（位能）模型：( )( )U rf r=rUr例：甲烷，某些惰性氣體質(zhì)點(diǎn)處理Ur方阱模型Ur階梯模型鍵的振動(dòng)鍵角扭矩分子內(nèi)部各原子（基團(tuán)）之間的范德華力、靜電力一般要計(jì)算1-4（相隔超過兩個(gè)鍵的原子或基團(tuán)對(duì)）15432qqq分子之間的范德華力分子之間的靜電力CCCHHHHHHHH 10根鍵長(zhǎng)作用 18個(gè)鍵角作用 8個(gè)扭矩作用 27個(gè)范德華力作用 27個(gè)靜電作用20()1exp() eu

5、lDa ll=- Morse類鍵長(zhǎng)模型能量阱深參數(shù)鍵長(zhǎng)平衡鍵長(zhǎng)參數(shù) 胡克類鍵長(zhǎng)模型20()()2ku lll=-鍵長(zhǎng)平衡鍵長(zhǎng)參數(shù) 胡克類鍵角模型20()()2ku lqq=-鍵角平衡鍵角參數(shù)0( ) 1cos()2NnnVunwwg=+-扭矩障礙參數(shù)扭動(dòng)360度所經(jīng)過的能量最低點(diǎn)的次數(shù)126126( )4U rrrsse驏=-桫Lennard-Jones模型Ur 混合規(guī)則ABAB1()2A BABsx ss=+A BABez e e=通常都取1 一般情況下只考慮點(diǎn)電荷之間的作用力 不考慮極化作用所帶來的長(zhǎng)程項(xiàng)的作用114ABNNi jijijqqUrpe=邋分子A的點(diǎn)電荷總數(shù)分子B的點(diǎn)電荷總數(shù)

6、 用有限的微觀分子體系模擬實(shí)際宏觀體系的必要手段3214342121計(jì)算周期邊界條件下兩個(gè)微粒之間的作用宏觀體系微觀體系主盒影子影子影子 Ewald方法 Particle Mesh Ewald 方法（PME) 計(jì)算量很大BAC( , , )( , )( , )( , )U A B CU A BU A CU B C+ 處理方法(2)(3)( )( , ,)totalijijikjkUUrUr r r=+邋邋 通過量子化學(xué)模擬回歸得到點(diǎn)電荷范德華力鍵長(zhǎng)、鍵角、扭矩力 實(shí)驗(yàn)數(shù)據(jù)回歸鍵長(zhǎng)鍵角范德華力初始化能量?jī)?yōu)化平衡數(shù)據(jù)產(chǎn)出避免局部分子重疊，并不是動(dòng)力學(xué)模擬 根據(jù)所有分子的當(dāng)前坐標(biāo)計(jì)算個(gè)分子的受力（位

7、能函數(shù)） 根據(jù)受力更新分子的坐標(biāo) 在此過程中收集用來計(jì)算宏觀性質(zhì)的有關(guān)信息讀入模型參數(shù)，模擬控制參數(shù)去除某些可能存在的原子重疊去除某些嚴(yán)重扭曲的鍵長(zhǎng)、鍵角、扭矩等方法最速下降法牛頓拉夫森方法其他一般優(yōu)化幾千到幾萬步 Verlet法簡(jiǎn)單易行，但是有精度損失Leap-frog法，Verlet法的變種 速度Verlet法和Verlet法相比，可以同時(shí)獲得位置、速度、加速度，而且沒有精度損失Beeman算法，速度計(jì)算精度更高，但是計(jì)算量大預(yù)測(cè)-校正法算法選擇精度和運(yùn)算量的折衷體系原子體系剛性分子體系非剛性分子，剛性鍵非剛性分子涉及到的分子運(yùn)動(dòng)平動(dòng)平動(dòng)、轉(zhuǎn)動(dòng)平動(dòng)、轉(zhuǎn)動(dòng)、扭矩平動(dòng)、轉(zhuǎn)動(dòng)、扭矩、振動(dòng)建議步長(zhǎng)

8、10 fs5 fs2 fs-1 fs 不同性質(zhì)的作用力采取不同的步長(zhǎng) 最常見的方法：r-RESPA 成功用于多種體系，簡(jiǎn)單流體，有機(jī)分子體系，蛋白質(zhì)等 和普通速度Verlet法相比，對(duì)于復(fù)雜分子體系（高分子，生物分子等），r-RESPA可以大大加快運(yùn)算速度，步長(zhǎng)可以達(dá)到1-2 fs 為保證原子各自運(yùn)動(dòng)時(shí)分子的整體性而添加的分子內(nèi)部的約束條件算法 最常見的是針對(duì)氫原子，因?yàn)闅湓雍芨叩恼駝?dòng)頻率 Verlet算法：SHAKE 速度Verlet算法：RATTLE NVE，最簡(jiǎn)單的分子動(dòng)力學(xué)E=K+U恒定體積，而且和外界沒有能量或者物質(zhì)交換，只有動(dòng)能和勢(shì)能之間的轉(zhuǎn)換通過監(jiān)視各能量項(xiàng)的變化，可以檢查程序

9、是否基本正確21/2(3)NiicBimKNNk T=-p p體系總動(dòng)能和體系溫度的關(guān)系原子總數(shù)體系總的約束條件數(shù)自由度模擬開始時(shí)，體系中各微粒的初始速度也是可以根據(jù)這個(gè)關(guān)系計(jì)算微粒 i 的動(dòng)量微粒 i 的質(zhì)量BPVN k TW=+壓力與分子間的作用勢(shì)能的關(guān)系理想氣體分子間的作用造成的剩余壓力1111( )33iNNiiiiiiWU= -籽=邋r rrrr frrr f維利系數(shù)W的計(jì)算擴(kuò)散系數(shù)二元徑向分布函數(shù)2( )()NNijij iVg rNd=-邋rrrr1( )(0)6iiDtt=-rrrrxtx通過斜率求取D實(shí)際位置，不是PBC后的位置 速度調(diào)節(jié) 簡(jiǎn)單scale v=v*(T0/T(

10、t)1/2 Berendsen法 增加調(diào)節(jié)參數(shù) 重新隨機(jī)設(shè)置速度每隔一定步數(shù)，按照控制溫度隨機(jī)抽取一組高斯速度分布熱量交換Nose-Hoover法將交換熱源也當(dāng)成體系的一部分進(jìn)行積分產(chǎn)生嚴(yán)格的NVT系綜控制體積活塞類似于溫度控制通過調(diào)節(jié)體系的體積來調(diào)節(jié)壓力壓力的控制比溫度更難，因?yàn)閴毫湍芰康奈⒎至浚Γ┏烧龋恳徊降牟▌?dòng)更大 氣體的模擬 小分子體系，不需要復(fù)雜的勢(shì)能模型 幾百到幾千個(gè)分子，分子分布稀疏，大部分是短程作用 一般用一臺(tái)微機(jī)就可以處理，計(jì)算時(shí)間幾分鐘幾小時(shí) 簡(jiǎn)單的液體，不涉及太多的界面性質(zhì) 小分子體系，勢(shì)能模型不是很復(fù)雜 幾百個(gè)分子，可能涉及到靜電作用，可能需要長(zhǎng)程校正 用微機(jī)也

11、可以處理，計(jì)算時(shí)間一般幾小時(shí)幾天 必要性 體系越來越大 模擬時(shí)間越來越長(zhǎng) 解決辦法 制造更快的處理器 并行計(jì)算機(jī) 例子：50000原子的生物體系，1ns模擬 單個(gè)處理器：12天 16個(gè)并行處理器：1天或者 Message Passing Interface 90年代初制定和完善的一套并行語法 支持Fortran, C, C+ 簡(jiǎn)單易學(xué) 并行效率需要1小時(shí)需要1/2小時(shí)完美的并行效率 處理器的速度遠(yuǎn)遠(yuǎn)超過數(shù)據(jù)傳輸?shù)乃俣?，大量的時(shí)間花在處理器之間的信息傳遞上了 CPU的速度幾乎是幾何級(jí)數(shù)增長(zhǎng) 內(nèi)存的速度是代數(shù)級(jí)數(shù)增長(zhǎng) 加快數(shù)據(jù)傳輸，盡量減少花在數(shù)據(jù)傳輸上的時(shí)間 數(shù)據(jù)傳輸硬件上的進(jìn)步 算法上做文章

12、每一個(gè)處理器負(fù)責(zé)處理一部分原子 每一步計(jì)算每一個(gè)處理器都要接受其它處理器負(fù)責(zé)處理的原子的相關(guān)信息 信息傳輸量大，使用的處理器越多并行效率效率越低 一般適合處理5-10萬左右微粒的體系N=20n=15n=610n=1620n=11151201481291961613155271741131018實(shí)際情況 按照體系的實(shí)際物理位置按區(qū)域劃分每個(gè)處理器的處理范圍 每一步計(jì)算每一個(gè)處理器只需要和相鄰的處理器交換信息 數(shù)據(jù)傳輸量小，并行效率高，適合處理大型體系（超過10萬微粒） 算法比較復(fù)雜（邊界的處理）N=201201481291961613155271741131018 主要針對(duì)與生物和化學(xué)軟材料體系

13、 優(yōu)點(diǎn) 程序設(shè)計(jì)水平高，計(jì)算效率高，號(hào)稱可以有效并行到上千個(gè)處理器 兼容多種輸入和輸出文件格式，有很好的分析輔助軟件VMD 有很好的維護(hù)服務(wù) 不需安裝 免費(fèi) 缺點(diǎn) 萬一需要自己安裝的話比較麻煩/Research/namd/ 主要針對(duì)生物體系，也適當(dāng)兼容一般化學(xué)分子 優(yōu)點(diǎn) 有很好的內(nèi)置勢(shì)能模型 自定義新模型和新分子很方便 有很完善的維護(hù)網(wǎng)站 缺點(diǎn) 計(jì)算效率不高（收斂到16個(gè)處理器），運(yùn)算速度慢 $400 主要針對(duì)生物體系，也包含部分化學(xué)體系 優(yōu)點(diǎn) 勢(shì)能模型更新很快 自定義新模型比較方便 維護(hù)服務(wù)很好 缺點(diǎn) 運(yùn)

14、算速度慢，計(jì)算效率低 $600/ 一般性分子動(dòng)力學(xué)軟件，對(duì)生物體系略有偏重 優(yōu)點(diǎn) 支持多種模型 免費(fèi) 缺點(diǎn) 仍在開發(fā)中，某些方面還不完善/tinker/ 一般性分子模擬軟件 優(yōu)點(diǎn) 兼容當(dāng)前大多數(shù)的勢(shì)能模型 編程水平高，計(jì)算效率高（比NAMD差，強(qiáng)于其他所有類似軟件） 可以模擬軟材料和固體物理系統(tǒng) 免費(fèi) 缺點(diǎn) 維護(hù)差/sjplimp/lammps.html 一般性分子模擬軟件 優(yōu)點(diǎn) 界面友好 計(jì)算效率高（有兩個(gè)版本供選擇，適合于不同大小的體系） 維護(hù)服務(wù)很好 缺點(diǎn) 兼

15、容性不好 100英鎊http:/www.cse.clrc.ac.uk/msi/software/DL_POLY/ 主要針對(duì)生物體系，也適當(dāng)照顧一般化學(xué)體系 優(yōu)點(diǎn) 算法好，計(jì)算效率高 界面友好 維護(hù)服務(wù)好 免費(fèi)軟件 缺點(diǎn) 兼容性不好/納米顆粒一般是人工合成物質(zhì)納米管，納米球，納米棒等等納米顆粒在生物醫(yī)學(xué)方面的應(yīng)用(PNAS, v99, p6645, 2002)藥物輸送，基因輸送 納米顆粒的生物安全性 (Science, v300,p243, 2003)C60在魚類腦部的累計(jì) (Environmental Health Perspectives, v112

16、, p1058, 2004)C60導(dǎo)致老鼠的肺部病癥 (Toxicol. Sci. v77, p117, p126, 2004)用分子動(dòng)力學(xué)研究C60和基因分子的作用MD引擎 NAMD靜電力-PME恒溫 (300 K)衡壓 (1 bar)時(shí)間步長(zhǎng) = 2 fs離子-中和體系模擬時(shí)間 4-20 nsDNA模型 AMBER 99 雙鏈, 12-20 堿基對(duì) 單鏈, 12-20 堿基C60, sp2 碳原子用LJ 模型模擬溶劑 4600到6000水分子 TIP3P模型水，離子水，離子水，離子水，離子三個(gè)方向施加周期邊界條件1 ns 的動(dòng)畫起始結(jié)構(gòu) C60和DNA能相互吸引，形成穩(wěn)定的雜化結(jié)構(gòu) 二者的

17、最初作用一般發(fā)生在1-2納秒之內(nèi)C60和B-DNA之間的距離和時(shí)間的函數(shù)C60和A-DNA之間的距離與時(shí)間的函數(shù) C60和DNA分子之間典型的結(jié)合能大約在15-40 kcal/mol, 取決于二者的相對(duì)作用位置和結(jié)合點(diǎn) C60和DNA分子之間的結(jié)合能遠(yuǎn)大于熱波動(dòng)kT (0.6 kcal/mol) 同時(shí)也遠(yuǎn)大于在同等模擬環(huán)境下C60和C60之間的結(jié)合能 (7.5 kcal/mol)C60吸附于DNA側(cè)面C60吸附于DNA尾部 C60和B-DNA形成穩(wěn)定的結(jié)合體，一般的結(jié)合能在20-40kcal/mol 比較穩(wěn)定的結(jié)合部位在DNA的尾部或者側(cè)面（minor grooves）C60進(jìn)入DNA尾部 C

18、60和A型DNA的側(cè)面和尾部結(jié)合，而且，C60可以進(jìn)入A型DNA的尾部如果沒有C60吸附，尾部的堿基對(duì)形成氫鍵的原子對(duì)之間的距離 穿入DNA尾部的C60破壞了DNA尾部的堿基對(duì)的氫鍵C60進(jìn)入DNA尾部后，尾部的堿基對(duì)形成氫鍵的原子對(duì)之間的距離變形的角度 吸附在A型DNA側(cè)面的DNA使得DNA堿基對(duì)之間的相對(duì)結(jié)構(gòu)變形。本應(yīng)平行的堿基對(duì)變形形成了一個(gè)大約20度的角 C60和單鏈DNA之間的結(jié)合導(dǎo)致單鏈DNA的嚴(yán)重變形 C60被單鏈DNA的某些部分包圍，疏水的堿基平面和C60緊密接觸 C60和單鏈DNA的結(jié)合能大約在20到30 kcal/mol C60和DNA之間結(jié)合的驅(qū)動(dòng)力是疏水作用DNA親水部位的水密度分布C60和DNA結(jié)合部位的水密度分布C60進(jìn)入損傷部