《化學(xué)統(tǒng)計熱力學(xué)》課件Molecular Dynamics

1、分子動力學(xué)方法Molecular Dynamics Method7/10/20221Molecular DynamicsMolecular Dynamics與Monte Carlo模擬的區(qū)別：Monte Carlo模擬中，分子坐標(biāo)的改變是完全隨機(jī)的，跟分子所受的力毫無關(guān)系分子力學(xué)中，根據(jù)力的大小和方向來決定移動的方向7/10/20222Molecular DynamicsMolecular Dynamics動力學(xué)與熱力學(xué)的關(guān)系：熱力學(xué)研究的是宏觀平衡態(tài)的性質(zhì)動力學(xué)研究的是態(tài)與態(tài)之間的轉(zhuǎn)化的過程7/10/20223Molecular DynamicsMolecular Dynamics中心問題：

2、解牛頓運(yùn)動方程，得到分子坐標(biāo)和速度是如何隨時間而改變的軌跡：trajectoryN個粒子，這樣的二階微分方程就有N個7/10/20224Molecular Dynamics1. MD的簡單模型硬球模型(Alder and Wainwright in 1957)：粒子作勻速直線運(yùn)動，直到它碰撞到其它粒子粒子間的相互作用只有在粒子的相互距離小于其半徑之和時才存在：Udistance7/10/20225Molecular DynamicsMD的簡單模型，continue硬球模型粒子間的相互作用的小改進(jìn)：Udistance7/10/20226Molecular DynamicsMD的簡單模型，cont

3、inue硬球模型的模擬步驟：按速率分布隨機(jī)給出初始的位置和速度計算碰撞時的位置根據(jù)動量守恒定律計算碰撞后的新速度Analyze Output是否繼續(xù)計算下一次碰撞時間7/10/20227Molecular Dynamics2. 使用連續(xù)勢能函數(shù)的MD作用在分子或粒子上的力其實是時刻在改變的無法求解N個粒子的運(yùn)動方程，因此采用其它法求解運(yùn)動方程F(t)的分析形式很復(fù)雜物理學(xué)到目前為止還為求出三體問題的精確分析解i=1,.,N，共N個二階微分方程!7/10/20228Molecular Dynamics解決辦法：有限差分方法時間t的狀態(tài)（速度，位置）是通過把0-t的積分分為很小的時間間隔t的運(yùn)動的

4、加和，假定粒子在這段很短的時間間隔之內(nèi)受到的力是常數(shù)ttr(t)r(t)r(t+t)7/10/20229Molecular Dynamics有限差分方法(Finite Difference Methods)假定粒子在某時刻附近的位置和速度可以用Taylor階數(shù)展開：7/10/202210Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法：利用t-t的位置r(t-t)和t時刻的位置r(t)和加速度a(t)計算t+t的位置r(t+t)7/10/202211Molecular Dynamics有限差分方法(Finite Differ

5、ence Methods)Verlet算法速度r(t)r(t+t)r(t-t)t-tt+ttr(t+t)-r(t-t)7/10/202212Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法速度也可以計算half-step (t+1/2t)的速度7/10/202213Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的優(yōu)點存儲量?。簉(t-t)，r(t)，a(t)Verlet算法的缺點精度初始位置r(0)和r(-t)：7/10/202214Molecular

6、Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改進(jìn)1leap-frog (Hockney, 1970)，蛙跳蛙跳的由來7/10/202215Molecular Dynamics有限差分方法(Finite Difference Methods)蛙跳算法的優(yōu)點直接包括速度，不再對較大的量(r(t)做差分7/10/202216Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改進(jìn)2velocity Verlet (Swope, 1982)速度也可以這樣算7/10/202217M

7、olecular Dynamics三種算法的比較a. Verlet; b. Leap flog; c. Velocity Verlet7/10/202218Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改進(jìn)3Beemans 算法(Beeman, 1976)7/10/202219Molecular DynamicsPredictor-Corrector integration methods (Gear, 1971)利用Taylor展開式計算下一時刻的位置，速度，加速度及其它高級項把實際計算ac(t+t)的與預(yù)測的a(t+

8、t)作比較7/10/202220Molecular DynamicsPredictor-Corrector integration methods (Gear, 1971)Gear 給出了一套他認(rèn)為最好的系數(shù)這些系數(shù)Taylor階數(shù)的截斷階數(shù)有關(guān)，這里（三階后截斷）7/10/202221Molecular DynamicsMD模擬的一些實際考慮該選擇什么方法上述算法的計算時間并不是最需要考慮的，最耗時間的步驟其實是計算能量，力和加速度需要保持能量和動量守恒不要期望MD得到的運(yùn)動軌跡和真實的軌跡非常一致甚至很接近7/10/202222Molecular Dynamics實際的一些考慮，conti

9、nue時間的步長(t)太?。褐惶剿髁讼嗫臻g的一小部分太大：能量的漲落很大，能量和動量不守恒時間模擬-真實10fs50fs100fs兩個Ar原子間的碰撞7/10/202223Molecular Dynamics實際的一些考慮，continue時間的步長(t)的設(shè)定SystemTypes of MotionSuggested time step (fs)Atomstranslation10Rigid moleculestranslation, rotation5Flexible molecules, rigid bondstranslation, rotation, torsion2Flexibl

10、e molecules, flexible bondstranslation, rotation, torsion, vibration1 or 0.57/10/202224Molecular Dynamics7.2-3 如何開始和進(jìn)行分子動力學(xué)模擬初始位置隨機(jī)確定初始動量Maxwell Boltzmann分布這是一個Gaussian分布，可以用7/10/202225Molecular Dynamics7.2-3 如何開始和進(jìn)行分子動力學(xué)模擬初始動量x使用生成隨機(jī)數(shù)的方法產(chǎn)生：這樣產(chǎn)生的x，它的平均值為零，且2=1，其它平均值不為零， 21的x這樣得到：注意動量守恒和能量守恒7/10/2022

11、26Molecular Dynamics7.2-3 如何開始和進(jìn)行分子動力學(xué)模擬計算相互作用，力和加速度例：Lennard-Jones potential運(yùn)行一段時間，達(dá)到平衡能量漲落變小溫度基本恒定速度和動量分布符合Maxwell-Boltzmann分布7/10/202227Molecular Dynamics7.2-3 如何開始和進(jìn)行分子動力學(xué)模擬計算熱力學(xué)函數(shù)溫度Nc：對體系的限制數(shù)目對孤立體系，1. 體系總的平動為常數(shù)（守恒）；2. 總的角動量也守恒，所以Nc = 6.對于周期性的體系，角動量不守恒， Nc = 3.3N-Nc：總的平動自由度7/10/202228Molecular D

12、ynamics7.2-3 如何開始和進(jìn)行分子動力學(xué)模擬計算熱力學(xué)函數(shù)其它熱力學(xué)量E(pi,qi)：每一步的體系能量，包括動能和勢能G(pi,qi)：要求平均值的物理量N：MD模擬所進(jìn)行的步數(shù)7/10/202229Molecular Dynamics7.2-4 受約束條件限制的動力學(xué)模擬模擬實際體系的困難：剛性分子相互作用于分子的取向有關(guān)柔性分子，質(zhì)心隨分子的構(gòu)象改變；相互作用同構(gòu)象有關(guān)；要考慮振動，模擬的時間步長取決于數(shù)值最大的振動頻率。SHAKE Procedure動力學(xué)模擬中固定鍵長一個完全限制(holonomic)的方法：7/10/202230Molecular DynamicsSHAK

13、E方法Constraint的實質(zhì)就是對分子中的原子加上了額外的力對鍵長固定，即兩原子間的距離在任意時刻滿足7/10/202231Molecular DynamicsSHAKE方法代入Verlet算法中下一時刻位置矢量的表達(dá)式如何求出k呢？7/10/202232Molecular DynamicsSHAKE方法根據(jù)限制條件，共k個SHAKE Procedure:迭代的方法求解，先猜測k-1個的值，根據(jù)兩個原子間距離(|r12|)的限制求出一個的值(12)，再假定一個未知(13)，代入12和其它k-2個猜測的值，求出13，如此循環(huán)往復(fù)，直至收斂。7/10/202233Molecular Dynam

14、icsSHAKE方法SHAKE Procedure:猜測13, 14,.,1k根據(jù)|r13(t)|=d13，代入12 , 14,.,1k，求13求出所有其它步驟收斂與否？根據(jù)|r12(t)|=d12求127/10/202234Molecular DynamicsSHAKE方法對鍵角的固定：轉(zhuǎn)化為對距離的固定7/10/202235Molecular Dynamics7.2-5 恒溫和恒壓下的動力學(xué)模擬上面的MD沒有提溫度和壓力，所以是在NVE(NVEP)系綜下進(jìn)行的。固定溫度7/10/202236Molecular Dynamics控制溫度因此只要在每一時刻（模擬的每一步）將這點的動量（或速度）

15、乘以一個因子缺點：hot solvent, cold solute7/10/202237Molecular Dynamics控制溫度的改進(jìn)方法：將體系與外部恒定溫度的熱源接觸，因此允許體系的溫度有波動。參考：Aderson, H. C. Molecular Dynamics Simulations at Constant Pressure and/or Temperature: J. Chem. Phys. 1980, 72, 2384-93.Nos, S. Molecular Dynamics Method for Simulations in the Canonical Ensemble:

16、 Mol. Phys. 1984, 53, 255-268.Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions: Phys. Rev. A 1985, 31, 1695-1697.7/10/202238Molecular Dynamics固定壓力保持壓力恒定的方式通過改變體積來實現(xiàn)的。把體系的體積乘以一個恰當(dāng)?shù)囊蜃影洋w系同一個外部的體系相作用，這個外部體系被作為“Pressure Bath”p：耦合因子把體積乘以一個因子（相應(yīng)的坐標(biāo)乘以1/3）：物質(zhì)的壓縮率7/10/202239Molecular Dy

17、namics7.2-5 MD和MC的應(yīng)用構(gòu)象分析系統(tǒng)搜索（手動，工作量大，陷入局部最低點）普通分子動力學(xué)（自動，陷入局部最低點）模擬退伙算法（自動，全局最優(yōu)）Kirkpatrick, S.; Gelatt, C. D.; Vecchi, M. P. Optimization by Simulated Annealing: Science 1983, 220, 671-680.遺傳算法（自動，全局最優(yōu)）7/10/202240Molecular Dynamics模擬退火算法在高溫進(jìn)行Monte-Carlo模擬一定的步數(shù)改變（降低）溫度進(jìn)行Monte-Carlo模擬7/10/202241Molecular Dynamics遺傳算法這已經(jīng)不是動力學(xué)模擬了產(chǎn)生一些代表分子不同構(gòu)象的population010111010011123000110111011123基因7/10/202242Molecular Dynamics遺傳算法將基因解碼為構(gòu)型參數(shù)，并計算其適應(yīng)度(fitnes