《計(jì)算物理學(xué)》課件 第四章 分子動力學(xué)模擬方法

分子動力學(xué)引言分子動力學(xué)基礎(chǔ)知識分子動力學(xué)的基本步驟分子間相互作用勢平衡態(tài)分子動力學(xué)模擬應(yīng)用舉例

分子動力學(xué)方法是一種計(jì)算機(jī)模擬實(shí)驗(yàn)方法。通過模擬實(shí)驗(yàn)或微觀粒子運(yùn)動過程，研究系統(tǒng)宏觀性質(zhì)的微觀機(jī)理等。

什么是分子動力學(xué)方法？

MolecularDynamics1957年至1959年間，Alder和Wainwright首次采用分子動力學(xué)模擬研究了硬球模型系統(tǒng)中無序-有序轉(zhuǎn)變問題，證實(shí)了硬球模型系統(tǒng)的一級相變行為。1964年，Rahman采用分子動力學(xué)方法研究了Lennard-Jones相互作用下液態(tài)Ar的結(jié)構(gòu)和動力學(xué)性質(zhì)。1967~1971年，Verlet模擬了更復(fù)雜的液體——水、溶鹽、聚合物和蛋白質(zhì)等。1974年，Rahman和Stillinger使用分子動力學(xué)模擬研究真實(shí)系統(tǒng)，研究對象為液體水。1977年，McCammon第一個蛋白質(zhì)的分子動力學(xué)研究，研究對象是牛胰島素抑制劑。分子動力學(xué)的簡要發(fā)展歷程分子動力學(xué)實(shí)例蛋白質(zhì)的折疊通過分子動力學(xué)模擬可以研究蛋白質(zhì)的折疊路徑。分子動力學(xué)實(shí)例水化（溶劑化）描述生物大分子周圍的溶劑（水）分子分布。抹香鯨肌紅蛋白的水化水化位置

（藍(lán)，實(shí)驗(yàn)）水化數(shù)密度最大值

（黃，MD模擬）分子動力學(xué)實(shí)例脆性斷裂分子動力學(xué)實(shí)例：晶化的分子動力學(xué)模擬FCC（綠）HCP（紫）BCC（紅）分子動力學(xué)的引入對于大量粒子體系，例如一滴水包含1021個水分子，幾乎不可能計(jì)算所有涉及到的狀態(tài)的物理量，再求其平均。而觀測量確是平均值，如何比較計(jì)算值和觀測值？近似：假設(shè)感興趣的物理量在體系到達(dá)一定尺度后，就基本不隨尺度變化而變化了。因此：對一些列有限的態(tài)做統(tǒng)計(jì)平均。隨機(jī)模擬方法（MC）和確定性模擬方法（MD）MD：由運(yùn)動方程描述，利用統(tǒng)計(jì)方法得到體系的靜態(tài)和動態(tài)性質(zhì)，從而得到體系的宏觀性質(zhì)?？梢钥醋魇求w系在一段時間內(nèi)發(fā)展過程的模擬。不存在任何的隨機(jī)因素

確定性模擬方法分子動力學(xué)簡介分子動力學(xué)是在原子、分子水平上求解多體問題的重要的計(jì)算機(jī)模擬方法。通過求解所有粒子的運(yùn)動方程，分子動力學(xué)方法可以用于模擬與粒子運(yùn)動路徑相關(guān)的基本過程。在分子動力學(xué)中，粒子的運(yùn)動行為是通過經(jīng)典的Newton運(yùn)動方程所描述。系統(tǒng)的所有粒子服從經(jīng)典力學(xué)的運(yùn)動規(guī)律，它的動力學(xué)方程就是從經(jīng)典力學(xué)的運(yùn)動方程——拉格朗日方程和哈密頓方程導(dǎo)出。分子動力學(xué)方法的優(yōu)勢分子動力學(xué)得天獨(dú)厚的優(yōu)勢是能夠模擬分子的運(yùn)動軌道，即能夠提供分子運(yùn)動和變化的最為微觀的、細(xì)致的信息。這對傳統(tǒng)的力學(xué)、熱力學(xué)、光譜學(xué)等實(shí)驗(yàn)方法都是難以辦到的?，F(xiàn)代分析理論和實(shí)驗(yàn)測量手段要準(zhǔn)確地描述物質(zhì)的非平衡態(tài)性質(zhì)是非常困難的。而分子動力學(xué)方法能夠計(jì)算得到運(yùn)行軌道，通過非常直接的計(jì)算公式給出感興趣的非平衡態(tài)物理量。這種特點(diǎn)對現(xiàn)今所有的理論、實(shí)驗(yàn)、計(jì)算手段來講都是獨(dú)一無二的。正是這種不可替代的作用，使得分子動力學(xué)在當(dāng)今自然科學(xué)研究中獲得了廣泛和深入的應(yīng)用。實(shí)際使用的限制在實(shí)際應(yīng)用中，分子動力學(xué)模擬方法和隨機(jī)模擬方法都面臨著兩個基本限制：有限觀測時間的限制有限系統(tǒng)大小的限制通常人們感興趣的是體系在熱力學(xué)極限下（即粒子數(shù)目趨于無窮時）的性質(zhì)。但是計(jì)算機(jī)模擬允許的體系大小要比熱力學(xué)極限小得多，因此可能會出現(xiàn)有限尺寸效應(yīng)。為了減小有限尺寸效應(yīng)，人們往往引入周期性、全反射、漫反射等邊界條件。當(dāng)然邊界條件的引入顯然會影響體系的某些性質(zhì)。數(shù)值求解中的問題對體系的分子運(yùn)動方程組采用計(jì)算機(jī)進(jìn)行數(shù)值求解時，需要將運(yùn)動方程離散化為有限差分方程。常用的求解方法有歐拉法、龍格-庫塔法等。數(shù)值計(jì)算的誤差階數(shù)顯然取決于所采用的數(shù)值求解方法的近似階數(shù)。原則上，只要計(jì)算機(jī)計(jì)算速度足夠大，內(nèi)存足夠多，我們可以使計(jì)算誤差足夠小。分子動力學(xué)的基礎(chǔ)知識運(yùn)動方程的導(dǎo)出求解運(yùn)動的數(shù)值算法：Verlet算法等相互作用勢周期性邊界條件截?cái)嗑嚯x最小像力約定分子動力學(xué)運(yùn)動方程數(shù)值求解系統(tǒng)的動力學(xué)機(jī)制決定運(yùn)動方程的形式。在分子動力學(xué)方法處理過程中，方程組的建立是通過對物理體系的微觀數(shù)學(xué)描述給出的。在這個微觀的物理體系中，每個分子都各自服從經(jīng)典的牛頓力學(xué)。每個分子運(yùn)動的內(nèi)稟動力學(xué)是用理論力學(xué)上的哈密頓量或者拉格朗日量來描述，也可以直接用牛頓運(yùn)動方程來描述。這種方法可以處理與時間有關(guān)的過程，因而可以處理非平衡態(tài)問題。運(yùn)動方程采用分子動力學(xué)方法時，必須對一組分子運(yùn)動微分方程做數(shù)值求解。從計(jì)算數(shù)學(xué)的角度來看，這是個求一個初值問題的微分方程的解。實(shí)際上計(jì)算數(shù)學(xué)為了求解這種問題已經(jīng)發(fā)展了許多的算法。但是并不是所有的這些算法都可以用來解決分子動力學(xué)問題。空間描述在空間描述如何物體的運(yùn)動，如果其本身的大小可以忽略時，就可以將其看作是粒子（或質(zhì)點(diǎn)）。粒子描述：空間位置：r速度：v=dr/dt加速度：a=dv/dt=d2r/dt2若一個系統(tǒng)由N個粒子組成，則粒子描述：空間位置：r1，r2，…

，rN笛卡爾坐標(biāo)系，粒子有3N個自由度設(shè)系統(tǒng)有s個自由度廣義坐標(biāo)：q1，q2，…

，qs廣義速度：

拉格朗日(Lagrange)方程

哈密頓方程粒子運(yùn)動方程的數(shù)值解法

采用分子動力學(xué)方法時，必須對一組分子運(yùn)動微分方程做數(shù)值求解，這是一個微分方程的初值問題。

速度Verlet算法

速度Verlet算法與Verlet算法是等價的。將改寫為代入可得方程（1）。

利用方程

Verlet的速度形式Verlet速度形式的算法比前一種算法好些。它不僅可以在計(jì)算中得到同一時間步上的空間位置和速度，而且數(shù)值計(jì)算的穩(wěn)定性也提高了。一般情況下，對于能量確定的系統(tǒng)不可能給出精確的初始條件。這時需要先給出一個合理的初始條件，然后在模擬過程中逐漸調(diào)節(jié)系統(tǒng)能量達(dá)到給定值。

如果消去速度項(xiàng)，該算法本質(zhì)上等價于Verlet算法。但比Verlet算法優(yōu)越，其好處在于不需要處理兩個大數(shù)的相減，從而盡可能地減小精度丟失。但相對于速度Verlet算法而言，該算法不能得到同一時刻的坐標(biāo)和速度。三種Verlet算法對比Verlet算法速度Verlet算法Leap-Frog算法《Computersimulationofliquids》（1986）Gear預(yù)測校正法預(yù)測校正法是分子動力學(xué)模擬中的常用算法之一，其基本思想是Taylor展開：

根據(jù)預(yù)測值和誤差可以得到各個量的校正值

方程組

可以表示成矩陣形式方程組

可以表示成矩陣形式Gear預(yù)測校正因子

1/65/611/3

19/1203/411/21/12

3/20251/360111/181/61/60保留不同階數(shù)的泰勒展開所得到的的修正系數(shù)如下Gear預(yù)測校正因子

1/65/611/3

19/903/411/21/12

3/16251/360111/181/61/60得到的的修正系數(shù)如下原子間相互作用勢

要進(jìn)行分子動力學(xué)模擬就必須知道原子間的相互作用勢。在分子動力學(xué)模擬中，我們一般采用經(jīng)驗(yàn)勢來代替原子間的相互作用勢，如Lennard-Jones勢、Morse勢、嵌入原子勢（Embedded-atommethodpotential）等。然而采用經(jīng)驗(yàn)勢必然丟失了局域電子結(jié)構(gòu)之間存在的強(qiáng)相關(guān)作用信息，即不能得到原子動力學(xué)過程中的電子性質(zhì)。

由于在模擬中，粒子之間的相互作用是通過勢函數(shù)來計(jì)算的，因此勢函數(shù)的選取非常重要，它決定了計(jì)算量和所計(jì)算模型與真實(shí)系統(tǒng)的近似程度。因此，人民在研究過程中提出了若干勢函數(shù)形式，如今，勢函數(shù)已從最初的兩體勢函數(shù)發(fā)展到多體勢函數(shù)。堿金屬鹽、堿土金屬鹽、硅氧化合物、蛋白質(zhì)、半導(dǎo)體、過渡金屬等多種物質(zhì)相互作用勢函數(shù)的確立，可以極大地提高運(yùn)算精度。根據(jù)所研究對象的不同，我們可以選用不同的勢函數(shù)。

勢函數(shù)發(fā)展簡介1924年，Lennard-Jones發(fā)表了著名的負(fù)冪函數(shù)式的Lennard-Jones勢函數(shù)的解析式。1929年，Morse發(fā)表了指數(shù)式的Morse勢。1931年，Born和Mayer發(fā)表了描述離子晶體中兩離子間相互作用的Born-Mayer勢函數(shù)。1983年至1984年間，Daw和Baskes提出了金屬材料中的嵌入原子勢(Embedded-AtomMethod)的概念和算法。幾乎同時，F(xiàn)innis和Sinclair根據(jù)密度函數(shù)二次矩理論提出形式上與EAM基本相同的經(jīng)驗(yàn)F-S模型。1985年Stillinger和Weber提出了描述半導(dǎo)體材料中原子間相互作用的多體半經(jīng)驗(yàn)勢。37N個相互作用的粒子系統(tǒng)的能量可以寫成單體、二體、三體以及更多體勢之和，即勢函數(shù)簡介因此Lennard-Jones勢、Morse勢都是二體勢函數(shù)EAM勢為多體勢函數(shù)。

39Lennard-Jones勢（1）早期的分子動力學(xué)主要用于模擬惰性氣體和液體，發(fā)展了一些對勢模型。其中應(yīng)用非常普遍的是Lennard-Jones勢。Lennard-Jones勢由JohnEdwardLennard-Jones在1924年提出，其表達(dá)式為：

40第一項(xiàng)為短程排斥力項(xiàng)，這一項(xiàng)對rij的負(fù)導(dǎo)數(shù)為正，即為排斥，冪次高，在rij小的時候起主導(dǎo)作用（左上圖中的藍(lán)線）第二項(xiàng)為遠(yuǎn)程吸引力項(xiàng)，這一項(xiàng)對rij的負(fù)導(dǎo)數(shù)為負(fù)，即為吸引，冪次低，在rij大的時候起主導(dǎo)作用（左上圖中的紅線）-ε為位勢的最小值，其出現(xiàn)在距離為r=21/6σ。α和β分別代表粒子A和B。r=σ時位勢為0。在分子動力學(xué)中，一般使用約化單位。當(dāng)選取原子的質(zhì)量=1時，原子的動量和速度以及力和加速度在數(shù)值上相等。

Lennard-Jones勢（2）41在Lennard-Jones勢作用下，第i個粒子與元胞內(nèi)其它N-1個粒子或其最鄰近的影像粒子的相互作用力在x方向的分量為：

Lennard-Jones勢函數(shù)與力場函數(shù)曲線Lennard-Jones勢（3）Morse勢1929年，Morse提出指數(shù)形式的勢函數(shù)來解決雙原子分子分解的問題。該勢函數(shù)形式如下

LJ和Morse模型勢只考慮了兩體相互作用而沒有考慮多體相互作用，因此無法用來描述原子間的成鍵情況，尤其是電荷密度分布呈現(xiàn)非對稱的體系。

1983年Daw和Baskes在電子的密度泛函理論和有效介質(zhì)理論基礎(chǔ)上發(fā)展了保證一定精度、用以描述金屬體系相互作用的半經(jīng)驗(yàn)的嵌入原子勢（EmbeddedAtomMethod，EAM）。

嵌入原子勢(embeddedatommethod,EAM)（1）EAM模型是一種半經(jīng)驗(yàn)多體勢模型。其中心思想是將原子周圍復(fù)雜的環(huán)境用膠凍模型簡化描述（所謂膠凍指的是在均勻正電荷背景上的均勻電子氣）。在嵌入原子勢方法中，主要的參數(shù)都放在原子的電子密度表示及相關(guān)形式中。這樣就把“原子對間”的性質(zhì)主要?dú)w結(jié)到“原子”的性質(zhì)，大大地簡化了計(jì)算嵌入原子勢（2）得到EAM能量表達(dá)式的過程中，最主要用到了兩個近似。體系的總能量可以表示為嵌入能和二體斥能之和。體系的電荷分布近似為各個原子電荷的簡單線性疊加。

上式第二項(xiàng)是對勢，與原子間距離相關(guān)，系數(shù)1/2避免兩體相互作用的重復(fù)計(jì)算。

Stillinger-Weber勢

1985年Stillinger和Weber提出了描述半導(dǎo)體材料中原子間相互作用的多體半經(jīng)驗(yàn)勢，能夠很好第描述半導(dǎo)體硅的四面體結(jié)構(gòu)。從基本的勢能展開式出發(fā)，假定二體項(xiàng)為

三體項(xiàng)表示為

Stillinger-Weber勢的參數(shù)選擇

分子動力學(xué)模擬中的條件設(shè)置

原胞與邊界條件截?cái)嗑嚯x最小像力約定Verlet近鄰表分子動力學(xué)模擬的基本步驟原胞與邊界條件（1）分子動力學(xué)模擬方法通常用于研究宏觀系統(tǒng)在給定條件下的性質(zhì)，但實(shí)際上不可能對一個近乎無窮大的系統(tǒng)進(jìn)行計(jì)算模擬。因此，分子動力學(xué)模擬需要引入一個體積元，通常稱為原胞，該原胞具有在給定條件下的宏觀系統(tǒng)的性質(zhì)，比如密度。對氣體和液體，如果所占體積足夠大，并且系統(tǒng)處于熱平衡狀態(tài)，這個體積的形狀是無關(guān)緊要的。對于晶體，原胞的形狀是有影響的。為了模擬方便，通常取一個立方體做為分子動力學(xué)模擬的原胞。為了計(jì)算簡便，通常取一個立方體為分子動力學(xué)原胞。設(shè)原胞的尺度為L，則體積為L3。由于引進(jìn)這樣的立方體盒子，將產(chǎn)生六個表面。如果采用這樣的原胞進(jìn)行分子動力學(xué)模擬，模擬中與表面碰撞的粒子將被反射回原胞內(nèi)部。如果原胞內(nèi)的粒子數(shù)目很少，處在原胞邊緣的粒子較多，其受力狀態(tài)與內(nèi)部粒子大為不同。這些表面的存在對系統(tǒng)的任何一種性質(zhì)都會有重大的影響。為了解決以上問題，分子動力學(xué)方法需要對模擬的原胞設(shè)置適當(dāng)?shù)倪吔鐥l件。原胞與邊界條件（2）邊界條件采用分子動力學(xué)方法，必需對被計(jì)算的粒子系統(tǒng)給定適當(dāng)?shù)倪吔鐥l件。這些邊界條件大致可分成四種。自由表面邊界

這種邊界條件通常用于模擬大型自由分子、團(tuán)簇等。因?yàn)?，這樣的系統(tǒng)本身就具有表面，并且表面的存在對系統(tǒng)的性質(zhì)起著重要的作用。

固定邊界(rigidboundary)

在所有要計(jì)算的粒子晶胞之外還要包上幾層結(jié)構(gòu)相同的位置不變的粒子，包層的厚度必須大于粒子間相互作用的力程范圍。包層部分代表了與運(yùn)動粒子起作用的宏觀晶體的那一部分。

比如在薄膜生長的模擬中，通常假設(shè)襯底的原子固定不動。薄膜生長

柔性邊界(flexibleboundary)

這種邊界比固定邊界更接近實(shí)際，它允許邊界上的粒子有微小的移動。一方面反映內(nèi)層粒子的作用力施加到它們身上時的情況，另一方面，由于熱運(yùn)動，邊界粒子也會有熱振動。因此，柔性邊界條件更反映模擬體系的實(shí)際情況。周期性邊界條件（1）在模擬較大的系統(tǒng)時，為了消除表面效應(yīng)或邊界效應(yīng)，通常引入了周期性邊界條件，即令元胞在三個方向周期性地重復(fù)出現(xiàn)，如圖所示。這樣可以消除引入原胞后的表面效應(yīng)，構(gòu)造出一個準(zhǔn)無窮大的體積來更精確地代表宏觀系統(tǒng)，從而將有限的原胞模擬，擴(kuò)展到真實(shí)系統(tǒng)的模擬。二維的周期性邊界條件示意圖。圖中陰影區(qū)為基本原胞，其周圍的8個原胞代表基本原胞在平面上的周期性重復(fù)。周期性邊界條件的數(shù)學(xué)表示形式為：周期性邊界條件（2）

周期性邊界條件就是命令基本原胞完全等同地重復(fù)無窮多次。實(shí)際模擬過程中通過以下操作實(shí)現(xiàn)周期性邊界條件，也就是當(dāng)一個粒子穿過基本原胞的某一個表面時，就讓該粒子也以相同的速度穿過該表面相對的表面，重新進(jìn)入該原胞內(nèi)。截?cái)嗑嚯x（1）一般情況下一個分子動力學(xué)步，求和運(yùn)算遍及所有相鄰粒子。若只進(jìn)行一次求和運(yùn)算，則不存在任何困難。在分子動力學(xué)模擬中，這種求和運(yùn)算將反復(fù)不斷地進(jìn)行，計(jì)算量非常大。如果按照常規(guī)的求和運(yùn)算，分子動力學(xué)模擬中可能絕大多數(shù)時間是用來尋找相鄰粒子、計(jì)算位勢及作用在粒子上的力，極大地消耗了計(jì)算資源?？紤]到兩個粒子之間的距離超過一定范圍后就可以忽略不計(jì)，因此求和不需要遍歷所有相鄰粒子，而只需要考慮有效的力場范圍。通常選取一個適當(dāng)?shù)慕財(cái)嗑嚯x將位勢截?cái)?/p>

截?cái)嗑嚯x（2）

截?cái)嗑嚯x（3）

60

這樣可以保證原子間作用力的連續(xù)性。截?cái)嗑嚯x（4）最小像力約定theminimumimageconvention（1）61由于采用了周期性邊界條件，在計(jì)算某個粒子的受力時，處在位勢截?cái)喟霃椒秶鷥?nèi)的粒子可能都在中心原胞內(nèi)，也有可能部分粒子位于近鄰原胞內(nèi)。因此，在分子動力學(xué)模擬中，每個粒子只與最近鄰粒子或者最近鄰像粒子發(fā)生相互作用，這就是所謂的最小像力約定，如圖所示。在考慮粒子間的相互作用時，通常采用最小像力約定。最小像力約定示意圖。紅色粒子與其他粒子之間的距離需要考察紅色方框中的粒子，其中既包含基本原胞內(nèi)的粒子，也包含鏡像原胞中的粒子。

62

最小像力約定（2）63

對每個原子，都可以構(gòu)造一個近鄰原子列表，每次計(jì)算只需要計(jì)算列表內(nèi)的原子間相互作用，之后再更新列表。提高運(yùn)算效率。Verlet近鄰表（1）

超過一定時間后，中心粒子的近鄰可能發(fā)生了變化，因?yàn)榻徢蛲獾牧Ｗ涌赡苓M(jìn)入球內(nèi)，需要對近鄰列表進(jìn)行更新。通常情況下，每隔10~20步更新一次近鄰列表。

Verlet近鄰表（2）66a

Updateintervalistheaveragenumberofstepsbetweenupdates.b

TimeisCPUtimeperstep,inseconds.《Computersimulationofliquids》（1986）分子動力學(xué)模擬的基本步驟（1）分子動力學(xué)模擬的實(shí)際步驟可以劃分為四步：1、設(shè)定模擬所采用的模型模型的設(shè)定，也就是勢函數(shù)的選取。勢函數(shù)的研究和物理系統(tǒng)上對物質(zhì)的描述研究息息相關(guān)。2、給定初始條件運(yùn)動方程的求解需要知道粒子的初始位置和速度，不同的算法要求不同的初始條件。如：Verlet算法需要兩組坐標(biāo)來啟動計(jì)算，一組零時刻的坐標(biāo)，一組是前進(jìn)一個時間步的坐標(biāo)或者一組零時刻的速度值。673、趨于平衡的計(jì)算過程模擬從一個初始態(tài)出發(fā)，該初始態(tài)通常都不是系統(tǒng)在該宏觀條件下的平衡態(tài)，因此需要模擬一個趨衡過程，即在這個過程中，增加或者從系統(tǒng)中移出能量，直到模擬得到的系統(tǒng)能量值趨于某個定值。稱這時的系統(tǒng)已經(jīng)達(dá)到平衡。4、宏觀物理量的計(jì)算當(dāng)?shù)玫较到y(tǒng)的平衡態(tài)后，可以根據(jù)統(tǒng)計(jì)物理的知識來計(jì)算系統(tǒng)的宏觀物理量，研究系統(tǒng)的物理性質(zhì)。68分子動力學(xué)模擬的基本步驟（2）1、設(shè)定模擬模型設(shè)定模型是分子動力學(xué)模擬的第一步工作。例如在一個分子系統(tǒng)中，假定兩個分子間的相互作用勢為硬球勢，其勢函數(shù)表示為常用的是Lennard-Jones勢函數(shù)、EAM勢函數(shù)等，不同的物質(zhì)狀態(tài)描述用不同的勢函數(shù)。69設(shè)定模擬模型為了計(jì)算方便，取分子動力學(xué)元胞為立方形，元胞的線度大小為L，則其體積為L3。根據(jù)給定密度ρ和指定的單個元胞中的粒子數(shù)N，確定出元胞的L值。采用周期性邊界條件和最小像力約定。根據(jù)實(shí)際需要解決的問題，構(gòu)造初始原胞中粒子的物理結(jié)構(gòu)模型。模型勢函數(shù)一旦確定，就可以根據(jù)經(jīng)典物理學(xué)規(guī)律求得模擬中的守恒量。例如對在微正則系綜的模擬中能量、動量和角動量均為守恒量。702、給定初始條件分子動力學(xué)模擬的過程進(jìn)入對系統(tǒng)微分方程組做數(shù)值求解時，需要知道粒子的初始位置和速度的數(shù)值。不同的算法要求不同的初始條件。一般講，系統(tǒng)的初始條件不可能知道，實(shí)際上也不需要精確選擇代求系統(tǒng)的初始條件，因?yàn)槟M時間足夠長時，系統(tǒng)就會忘掉初始條件。當(dāng)然，合理的初始條件可以加快系統(tǒng)趨于平衡的時間和步伐，獲得好的精度。71三種初始條件常用的初始條件可以選擇為：令初始位置在差分劃分網(wǎng)格的格子上，初始速度則從玻爾茲曼分布隨機(jī)抽樣得到（速度隨機(jī)）令初始位置隨機(jī)的偏離差分劃分網(wǎng)格的格子上，初始速度為零（位置隨機(jī)）令初始位置隨機(jī)的偏離差分劃分網(wǎng)格的格子上，初始速度也是從玻爾茲曼分布隨機(jī)抽樣得。（位置、速度均隨機(jī)）723、趨于平衡計(jì)算按照給出的運(yùn)動方程、邊界條件和初始條件，就可以進(jìn)行分子動力學(xué)模擬計(jì)算。但是，這樣計(jì)算出的系統(tǒng)不會具有所要求的系統(tǒng)能量，并且這個狀態(tài)本身也還不是一個平衡態(tài)。為了使系統(tǒng)達(dá)到平衡，模擬中需要一個趨衡過程。在這個過程中，增加或從系統(tǒng)中移出能量，直到系統(tǒng)具有所要求的能量。然后，再對運(yùn)動方程中的時間向前積分若干步，使系統(tǒng)持續(xù)給出確定能量值。稱這時系統(tǒng)已經(jīng)達(dá)到平衡態(tài)。這段達(dá)到平衡所需的時間稱為弛豫時間。73步長選取分子動力學(xué)計(jì)算的基本思想是：賦予分子體系初始運(yùn)動狀態(tài)之后，利用分子的自然運(yùn)動在相空間中抽取樣本進(jìn)行統(tǒng)計(jì)分子動力學(xué)計(jì)算，時間步長就是抽樣的間隔，因而時間步長的選取對動力學(xué)模擬非常重要。它決定了模擬所需要的時間。時間步長δt太長會造成分子間激烈碰撞，體系數(shù)據(jù)溢出；

為了減小誤差，步長δt必須取得小一些；但是取得太小，系統(tǒng)模擬的弛豫時間就很長，太短的時間步長還會降低模擬過程搜索相空間的能力。因此一般選取的時間步長為體系各個自由度中最短運(yùn)動周期的十分之一。這里需要積累一定的模擬經(jīng)驗(yàn)，選擇適當(dāng)?shù)臅r間步長δt。74分子動力學(xué)模擬中的系綜在經(jīng)典分子動力學(xué)模擬中，需要根據(jù)所研究的問題選擇合適的系綜。76系綜是指在一定的宏觀條件下，將所研究系統(tǒng)的各種微觀狀態(tài)假想成一個大量系統(tǒng)的集合，這些系統(tǒng)與研究的系統(tǒng)的性質(zhì)完全相同，且彼此相互獨(dú)立，各自處于該真實(shí)系統(tǒng)的某一微觀運(yùn)動狀態(tài)。這樣的假象的系統(tǒng)集合全稱為系綜。系綜是用統(tǒng)計(jì)方法描述熱力學(xué)系統(tǒng)的統(tǒng)計(jì)規(guī)律性時引入的一個基本概念；系綜是統(tǒng)計(jì)理論的一種表述方式，系綜理論使統(tǒng)計(jì)物理成為普遍的微觀統(tǒng)計(jì)理論；系綜并不是真實(shí)的系統(tǒng)，構(gòu)成系綜的系統(tǒng)才是真實(shí)的。約束條件是由一組宏觀參量來表示。常用系綜分類（1）微正則系綜

系統(tǒng)具有確定的粒子數(shù)（N）、體積（V）和總能量（E），簡寫為NVE。系統(tǒng)與外界即無能量交換，無粒子交換，即為孤立系統(tǒng)。假定N個粒子處在體積為V的盒子內(nèi)，并固定總能量(E)。此時，系綜的溫度(T)和系統(tǒng)壓強(qiáng)(P)可能在某一平均值附近起伏變化。平衡體系為孤立系統(tǒng)，與外界即無能量交換，也無粒子交換。微正則系綜的特征函數(shù)是熵S(N,V,E)。77常用系綜分類（2）正則系綜(canonicalensemble)全稱應(yīng)為“宏觀正則系綜”，簡寫為NVT，即表示具有確定的粒子數(shù)(N)、體積(V)、溫度(T)。正則系綜是蒙特卡羅和分子動力學(xué)模擬處理的典型代表。

假定N個粒子處在體積為V的盒子內(nèi)，將其埋入溫度恒為T的熱浴中。此時，總能量(E)和系統(tǒng)壓強(qiáng)(P)可能在某一平均值附近起伏變化。平衡體系代表封閉系統(tǒng)，與大熱源熱接觸平衡的恒溫系統(tǒng)。正則系綜的特征函數(shù)是亥姆霍茲自由能F(N,V,T)。78常用系綜分類（3）等溫等壓系綜

系統(tǒng)具有確定的粒子數(shù)（N）、壓強(qiáng)（P）、溫度（T），簡寫為NPT。系統(tǒng)為閉合系統(tǒng)，與大熱源接觸維持恒溫，并且系統(tǒng)壁可移動以維持恒壓狀態(tài)。

其總能量(E)和系統(tǒng)體積(V)可能存在起伏。體系是可移動系統(tǒng)壁情況下的恒溫?zé)嵩　?/p>

特征函數(shù)是吉布斯自由能G(N,P,T)。79常用系綜分類（4）

80系綜調(diào)節(jié)調(diào)溫技術(shù)在NVT系綜或NPT系綜中，即使在NVE系綜模擬的平衡態(tài)中，也經(jīng)常調(diào)整溫度到期望值。

如果希望知道系統(tǒng)的平衡態(tài)性質(zhì)怎樣依賴于溫度，那么就必須在不同的溫度下進(jìn)行模擬。

目前常用的實(shí)現(xiàn)對溫度的調(diào)節(jié)有4種方式：

速度標(biāo)度法、Andersen熱浴、Langevin熱浴、Nose—Hoover熱浴。81系綜調(diào)節(jié)主要是指在分子動力學(xué)模擬中，對溫度和壓強(qiáng)等參量進(jìn)行調(diào)節(jié)，從而實(shí)現(xiàn)正則系綜、巨正則系綜等的模擬。系綜調(diào)節(jié)–調(diào)壓技術(shù)調(diào)壓技術(shù)

在等壓模擬中，可以通過改變模擬原胞的三個方向或一個方向的尺寸來實(shí)現(xiàn)體積的變化。

類似于溫度控制的方法，也有許多方法用于壓力控制。

常用的有以下2種技術(shù)：

Anderson方法、Parrinello-Rahman方法。82溫度調(diào)節(jié)（1）

83溫度調(diào)節(jié)這樣就可以采用速度標(biāo)度法來控制系統(tǒng)的溫度。速度標(biāo)度因子:84

該控溫方法簡單直接，但是可能導(dǎo)致所模擬的系統(tǒng)并不符合等溫系綜。Andersen熱浴法

85Langevin熱浴法該方法也是引入熱池，但是通過虛構(gòu)的熱池粒子引入阻尼和碰撞，使得模擬系統(tǒng)內(nèi)的粒子保持平均動能，從而達(dá)到調(diào)節(jié)溫度的目的。86

以上溫度調(diào)節(jié)算法比較粗糙，沒有嚴(yán)格的理論基礎(chǔ)。研究人員在系統(tǒng)的廣義坐標(biāo)和廣義動量之外，引入一個額外的自由度與熱浴耦合的方法，實(shí)現(xiàn)對溫度的調(diào)節(jié)。該方法也稱為恒溫?cái)U(kuò)展法，其本質(zhì)是通過改變或擴(kuò)展模擬系統(tǒng)的哈密頓函數(shù)或拉格朗日函數(shù)實(shí)現(xiàn)對正則系綜的分子動力學(xué)模擬，例如Nose-Hoover控溫法。Nose-Hoover控溫法Nose-Hoover方法引入了一個額外的自由度s，新的自由度的引入相當(dāng)于引入了一個恒溫?zé)釒?，使得系統(tǒng)與之接觸，允許熱流在系統(tǒng)和熱池之間交換，趨于和熱庫熱平衡。系統(tǒng)的溫度與恒溫?zé)釒煜嗤３趾銣亍?/p>

壓強(qiáng)調(diào)節(jié)對于壓強(qiáng)的調(diào)節(jié)，最直接的方法是通過調(diào)整模擬系統(tǒng)的原胞三個方向或一個方向的尺寸來改變系統(tǒng)的體積，從而實(shí)現(xiàn)對壓強(qiáng)的調(diào)控。但是，該方法存在一個嚴(yán)重問題，如果直接將每個