1、分子動(dòng)力學(xué)模擬入門第1頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一引 言物質(zhì)基本構(gòu)成分子、原子在分子、原子這個(gè)微觀水平上來考察物質(zhì)：多體世界查清楚微觀世界，宏觀就清楚了從微觀考慮問題的現(xiàn)實(shí)可行性從微觀考慮問題的必要性物性的觀測性參數(shù)：熱傳導(dǎo)、溫度、壓力、粘性、. 第2頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一微觀處理的前提已知微觀粒子間的相互作用假設(shè)分子為球，惰性，分子間的作用只取決于分子間的距離分子動(dòng)力學(xué) (Molecular Dynamics，MD)第3頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一MD的應(yīng)用領(lǐng)域：物理、化學(xué)、生物、材料等MD方法能實(shí)時(shí)將

2、分子的動(dòng)態(tài)行為顯示到計(jì)算機(jī)屏幕上, 便于直觀了解體系在一定條件下的演變過程MD含溫度與時(shí)間, 因此還可得到如材料的玻璃化轉(zhuǎn)變溫度、熱容、晶體結(jié)晶過程、輸送過程、膨脹過程、動(dòng)態(tài)弛豫(relax)以及體系在外場作用下的變化過程等水和離子在微小硅孔中的運(yùn)動(dòng)聚乙烯的結(jié)晶第4頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一MD的基本原理用牛頓經(jīng)典力學(xué)計(jì)算許多分子在相空間中的軌跡求解系統(tǒng)中的分子或原子間作用勢能和系統(tǒng)外加約束共同作用的分子或原子的牛頓方程。模擬系統(tǒng)隨時(shí)間推進(jìn)的微觀過程。通過統(tǒng)計(jì)方法得到系統(tǒng)的平衡參數(shù)或輸運(yùn)性質(zhì)計(jì)算程序較為復(fù)雜，占用較多內(nèi)存第5頁，共39頁，2022年，5月20日，9

3、點(diǎn)51分，星期一MD的主要步驟選取要研究的系統(tǒng)及其邊界，選取系統(tǒng)內(nèi)粒子間的作用勢能模型設(shè)定系統(tǒng)中粒子的初始位置和初始動(dòng)量建立模擬算法，計(jì)算粒子間作用力及各粒子的速度和位置當(dāng)體系達(dá)到平衡后，依據(jù)相關(guān)的統(tǒng)計(jì)公式，獲得各宏觀參數(shù)和輸運(yùn)性質(zhì)第6頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一分子間勢能及相互作用N個(gè)粒子系統(tǒng)的總勢能第7頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一剛球模型斥力力心點(diǎn)模型Southerland模型第8頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一分子間勢能及相互作用Lennard-Jones勢能 能量尺度； 長度尺度 為方便，時(shí)常歸一化：記 ；

4、第9頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一分子間勢能及相互作用一些氣體的參數(shù)kB=1.38x1023(J/K): Boltzmann常數(shù)NeonArgonKryponXenonNitrogen(nm)0.2750.34050.3600.4100.370/kB(K)36119.817122195第10頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一分子間勢能及相互作用相互作用標(biāo)量形式：直角坐標(biāo)：至此，各粒子間相互作用已知，可進(jìn)行模擬了第11頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法Euler法和Euler-Cromer方法？不能用：不能保持總

5、能量守恒Verlet算法：速度形式第12頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法Leap-frog算法：x的截?cái)嗾`差為 ，v的截?cái)嗾`差第13頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法邊界條件模擬能力限制，不能模擬大量分子，只能模擬有限空間中的有限個(gè)分子：有限空間邊界固體（剛性）邊界條件不僅僅有分子間的相互作用，還引入了壁面的作用分子量大時(shí)，壁面作用可忽略不計(jì)第14頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法取 ， 前比值為0.20.01。取前值，模擬粗糙；取后值，模擬計(jì)算量太大處理方法：使用周期性邊界條件周期性

6、邊界條件第15頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法兩個(gè)不同粒子在x或y方向上的最大分離距離為a/2最小像約定：兩粒子分離距離最大分離距離，相互作用力可以忽略，而加入其中像粒子之一相互作用力來考慮第16頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法考慮的粒子總數(shù)不變初始條件隨機(jī)初始條件給法之一 要求 大?。簵l件一：規(guī)則給法條件二：隨機(jī)給法第17頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬的數(shù)學(xué)方法 random：隨機(jī)數(shù)產(chǎn)生函數(shù)，產(chǎn)生(0，1)之間的隨機(jī)數(shù)。 方向(按球坐標(biāo)給法)： 分量 ： 第18頁，共39頁，2022年，

7、5月20日，9點(diǎn)51分，星期一模擬微觀量溫度 根據(jù)統(tǒng)計(jì)熱力學(xué)，平衡態(tài)下經(jīng)典系統(tǒng)的能量中的每一個(gè)二次項(xiàng)具有平均值kBT/2,即 注意：上式在系統(tǒng)質(zhì)心速度為0時(shí)適用空間維數(shù)粒子個(gè)數(shù)：取時(shí)間平均第19頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬問題：如何給定系統(tǒng)的初始條件，得到所需要的平衡態(tài)溫度Teq？解決方法之一：速度標(biāo)定法 任給初始條件，模擬到平衡，得到系統(tǒng)平衡態(tài)溫度T。一般TTeq。令 用速度 再模擬直到平衡，若所得溫度仍不等于Teq，再進(jìn)行上述過程第20頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一給定初始條件：xi，vi計(jì)算到平衡態(tài)|T - Teq| f = Te

8、q/ Tvi = vi f 1/d計(jì)算結(jié)束YN第21頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一其它方法：Gaussian熱浴法（約束溫度調(diào)節(jié)方法） 其基本原理在運(yùn)動(dòng)方程中加入“摩擦力”項(xiàng)，并將其與粒子速度聯(lián)系起來。平衡態(tài)時(shí)，系統(tǒng)溫度不變，因此dEk/dt=0第22頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一宏觀性質(zhì)的統(tǒng)計(jì)系統(tǒng)的勢能系統(tǒng)的內(nèi)能系統(tǒng)的總能 E = Ep+Ek系統(tǒng)的溫度第23頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬熱容 定義熱容 計(jì)算系統(tǒng)在溫度T和T+T時(shí)的總能ET、ET +T，E：系統(tǒng)總能第24頁，共39頁，2022年，5月20日，9點(diǎn)

9、51分，星期一模擬壓強(qiáng)對(duì)壁面的壓強(qiáng)t時(shí)間里作用在單位面積壁上的壓力時(shí)刻，速度為時(shí)刻，速度為剛性壁第25頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬粒子速度分布 選速度間隔v，模擬nt個(gè)時(shí)間步，記錄在每個(gè)速度間隔中的粒子數(shù)，最后歸一化。第26頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬氣、液狀態(tài)方程維里定理(Virial Theorem)壓強(qiáng)體積粒子i的位置矢量粒子i所受到的其它粒子的合相互作用力溫度的模擬可得此項(xiàng)在溫度的模擬基礎(chǔ)上再模擬此項(xiàng)第27頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一模擬例：用此可確定高密度氣體和液體狀態(tài)方程(van der W

10、aals方程)中的系數(shù)理想氣體狀態(tài)方程在高密度情況下不可用確定系數(shù)a和b氣體密度第28頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE(1) - Couette Flow If is less than 0, then the two species are immiscible.i, j represent different species第29頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE(1) - Couette Flow ztime第30頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE(1) - Couette Flow z第31頁，

11、共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE(2) - Contact Angle SimulationMass: m1=1, m2=8, m3=0.8L=25.05,W=6.56H=10.29T=1.2i, j=1,2,3, 1 red fluid, 2 - wall, 3 green fluid 第32頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一V第33頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE(3)- Rayleigh-Taylor instability generation(重力場)Dzwinel, W., Alda, W., Pog

12、oda, M., and Yuen, D.A., 2000, Turbulent mixing in the microscale: a 2D molecular dynamics simulation, Physica D, Vol. 137, pp. 157-171.Fig. 1. The snapshots of MD simulations of the RT instability for two particle systems: (A) closed; (B) open (the heavy fluid is coloured inlight grey while the lig

13、hter one is dark grey, the black part of the figure is empty, i.e., it does not contain particles).Fig. 2. The effect of converging geometry obtained by MD simulation of one million particles in the microscale.0.5 m第34頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一CASE (4)Typical translocation eventA 1.4V bias applied

14、 to membrane.20 base-pair fragment of double stranded DNA placed in front of a nanopore.End of DNA nearest to the pore is pulled into the pore by its charged backbone (a,b)System reaches a meta-stable state (c) and translocation halts.Base-pairs start to split. Some freed nucleotides adhere to pore

15、surface.Voltage increased momentarily to drive system out of metastable state.DNA exits pore. One of the bases holds firmly to the pore surface.After 50ns, most of DNA has left pore. Nine of twenty base pairs are split.第35頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一 Bubble nucleation on solid surfacesMaruyama, S, a

16、nd Kimura, T., 2000, A Molecular Dynamics Simulation of a Bubble Nucleation on Solid Surface, Heat and Technology, Vol. 18, pp. 69-73.Molecular transport in dropletsMaruyama, S., Matsumoto, S., and Ogita, A., 1994, Surface Phenomena of Molecular Clusters by Molecular Dynamics Method, Thermal Science and Engineering, Vol. 2, No. 1. Full View Sliced View Void View第36頁，共39頁，2022年，5月20日，9點(diǎn)51分，星期一參考書目D. C. Rapaport.