第6章分子動力學方法

1、- 28 - 第6章 分子動力學方法 第6章 分子動力學方法經(jīng)典分子動力學方法無疑是材料，尤其是大分子體系和大體系模擬有效的方法之一。分子動力學可以用于NPT，NVE，NVT等不同系綜的計算，是一種基于牛頓力學確定論的熱力學計算方法。與蒙特卡羅法相比在宏觀性質(zhì)計算上具有更高的準確度和有效性，可以廣泛應用于物理，化學，生物，材料，醫(yī)學等各個領(lǐng)域。本章在介紹分子動力學的基本概念的基礎(chǔ)上，簡單介紹了分子動力學的基本思想，勢函數(shù)分類和基本方程。然后介紹了分子動力學的常用系綜和典型的NPT，NVE，NVT系綜基本方程。結(jié)合材料建模中的基本簡化方法和技巧，闡述了邊界條件和時間積分的數(shù)值處理技巧。最后，利用

2、統(tǒng)計力學的基本概念給出分子動力學的計算信息的解析方式。并且結(jié)合Materials Explore軟件計算分析了CNT的幾何結(jié)構(gòu)穩(wěn)定性。6.1 引言分子動力學方法(Molecular Dynamics, MD)方法是一種按該體系內(nèi)部的內(nèi)稟動力學規(guī)律來計算并確定位形的變化的確定性模擬方法。首先需要在給定的外界條件下建立一組粒子的運動方程，然后通過直接對系統(tǒng)中的一個個粒子運動方程進行數(shù)值求解，得到每個時刻各個分子的坐標與動量，即在相空間的運動軌跡，再利用統(tǒng)計力學方法得到多體系統(tǒng)的靜態(tài)和動態(tài)特性，從而獲得系統(tǒng)的宏觀性質(zhì)?？梢钥闯?，分子動力學方法中不存在任何隨機因素，這個也是分子動力學方法和后文要提到的

3、蒙特卡洛方法的區(qū)別之一。在分子動力學方法的處理過程中，方程組的建立是通過對物理體系的微觀數(shù)學描述給出的。在這個微觀的物理體系中，每個分子都各自服從經(jīng)典的牛頓力學定律(或者是拉格朗日方程)。每個分子運動的內(nèi)稟動力學是用理論力學上的哈密頓量或者拉格朗日函數(shù)來描述，也可以直接用牛頓運動方程來描述。確定性方法是實現(xiàn)玻爾茲曼的統(tǒng)計力學途徑。這種方法可以處理與時間有關(guān)的過程，因而可以處理非平衡態(tài)問題。但是分子動力學方法的計算機程序相對蒙特卡羅較復雜，其計算成本較高。分子動力學方法發(fā)展歷史改革經(jīng)歷了近60年，分子動力學方法是20世紀50年代后期由Alder B J 和Wainwright T E創(chuàng)造發(fā)展的。

4、Alder和Wainwright在1957年利用分子動力學模擬，發(fā)現(xiàn)了“剛性球組成的集合系統(tǒng)會發(fā)生由其液相到結(jié)晶相的相轉(zhuǎn)變”，后來人們稱這種相變?yōu)锳lder相變。其結(jié)果表明，不具有引力的粒子系統(tǒng)也具有凝聚態(tài)。到20世紀70年代，產(chǎn)生了剛性體系的動力學方法，成功地被應用于水和氮等分子性溶液體系的處理；1972年，Less A W和 Edwards S F 等發(fā)展了該方法并擴展到了存在速度梯度，即處于非平衡狀態(tài)的系統(tǒng)。之后，此方法被 Ginan M J等推廣到了具有溫度梯度的非平衡系統(tǒng)，從而構(gòu)造并形成了所謂的非平衡分子動力學方法體系。進人20世紀80年代之后，出現(xiàn)了在分子內(nèi)部對一部分自由度施加約束

5、條件的新的分子動力學方法，從而使 分子動力學方法可適用于類似蛋白質(zhì)等生物大分子的解析與設(shè)計。分子動力學方法真正作為材料科學領(lǐng)域的一個重要研究方法，開始于由Andersen, Parrinello和Rahman等創(chuàng)立恒壓分子動力學方法和Nse等完成恒溫分子動力學方法的建立及在應用方面的成功。后來，針對勢函數(shù)模型化比較困難的半導體和金屬等，1985 年人們又提出了將基于密度泛函理論的電子論和分子動力學方法有機統(tǒng)一起來的所謂Car-Parrinello 方法，亦即第一性原理分子動力學方法。這樣，分子動力學的方法進一步得到發(fā)展和完善，它不僅可以處理半導體和金屬的間題，同時還可應用于處理有機物和化學反應

6、。關(guān)于Car-Parrinello分子動力學方法將在第7章重點學習討論。本章重點討論經(jīng)典分子動力學方法的基本原理和計算方法。6.2 分子動力學的計算框架6.2.1 基本思想分子動力學方法是沿用牛頓運動方程或拉格朗日方程、哈密頓方程，通過考察粒子的運動來研究多粒子系統(tǒng)的物理性質(zhì)。在處理孤立粒子或原子團簇(不考慮與其他粒子的相互作用)時，單純地對牛頓運動方程進行積分求解就行了；而在處理凝聚系統(tǒng)時，可以認為將所考慮的N個粒子放入具有一定體積的容器中，在這樣組成的封閉系中，建立直角坐標(x, y, z)，每個粒子的位置就由三個坐標分量決定，通過求解3N個聯(lián)立方程組就可以得到保守系的總能量。分子動力學的

7、基本思想在于通過設(shè)定原子之間的相互作用(勢函數(shù))和相關(guān)的系綜(亦即作用對象和條件)，確定其基本的模擬范疇。在給定的牛頓方程、拉格朗日方程或者是哈密頓方程進行時間的迭代，在達到指定的力學收斂條件之后得到一個最終的位置坐標，然后通過相關(guān)的位形信息和運動的速度和加速度等信息通過熱力統(tǒng)計學的方法，給出模擬體系的統(tǒng)計信息，主要包括復雜的光學，電學和一些其他的物理信息。其基本思想，總結(jié)如圖6-1所示。圖 6-1 分子動力學的一般思想6.2.2 計算流程總體來說，分子動力學的基本計算有如下四個步驟：1. 模型的設(shè)定，也就是勢函數(shù)的選取勢函數(shù)的研究和物理系統(tǒng)中對物質(zhì)的描述研究息息相關(guān)。分子動力學的模擬最早使用

8、是硬球勢，即小于臨界值時無窮大，大于等于臨界值時為零。常用的是Lennard-Jones勢函數(shù)，還有嵌入原子EAM勢函數(shù)，不同的物質(zhì)狀態(tài)描述用不同的勢函數(shù)。模型勢函數(shù)一旦確定，就可以根據(jù)物理學規(guī)律求得模擬中的守恒量。關(guān)于勢函數(shù)選取將在6.4節(jié)進一步討論。2. 給定初始條件運動方程的求解需要知道粒子的初始位置和速度，不同的算法要求不同的初始條件。如Verlet算法需要兩組坐標來啟動計算，一組零時刻的坐標，一組是前進一個時間步的坐標或者一組零時刻的速度值。一般意義上講系統(tǒng)的初始條件不可能知道，實際上也不需要精確選擇代求系統(tǒng)的初始條件，因為當模擬時間足夠長時，系統(tǒng)就會忘掉初始條件(對于無記憶的體系而

9、言)。當然，合理的初始條件可以加快系統(tǒng)趨于平衡的時間和步程，獲得好的精度。常用的初始條件可以選擇為，令初始位置在差分劃分網(wǎng)格的格子上，初始速度則從玻爾茲曼分布隨機抽樣得到；令初始位置隨機的偏離差分劃分網(wǎng)格的格子上，初始速度為零；令初始位置隨機的偏離差分劃分網(wǎng)格的格子上，初始速度也是從玻爾茲曼分布隨機抽樣得到。圖6-2 分子動力學計算流程3. 趨于平衡計算在邊界條件和初始條件給定后就可以解運動方程，進行分子動力學模擬。但這樣計算出的系統(tǒng)是不會具有所要求的系統(tǒng)能量，并且這個狀態(tài)本身也還不是一個平衡態(tài)。為使得系統(tǒng)平衡，模擬中設(shè)計一個趨衡過程，即在這個過程中，增加或者從系統(tǒng)中移出能量，直到持續(xù)給出確定

10、的能量值，稱此時系統(tǒng)已經(jīng)達到平衡，達到平衡的時間稱為弛豫時間。 分子動力學中，積分步程的大小選擇十分重要，這決定了模擬所需要的時間。為了減小誤差，步長要小，但過小系統(tǒng)模擬的弛豫時間就長。因此根據(jù)經(jīng)驗選擇適當?shù)牟介L。如，對一個具有幾百個氬氣分子的體系，采用Lennard-Jones勢函數(shù)，發(fā)現(xiàn)取h為0.01量級，可以得到很好的結(jié)果。這里選擇的h是沒有量綱的，實際上這樣選擇的h對應的時間在10-14s的量級。如果模擬1000步，系統(tǒng)達到平衡，弛豫時間只有10-11s。4. 宏觀物理量的計算實際計算宏觀的物理量往往是在模擬的最后階段進行的，它是沿相空間軌跡求平均來計算得到的，根據(jù)各態(tài)歷經(jīng)假說，時間平

11、均代替系綜平均。關(guān)于各態(tài)歷經(jīng)假說將在第8章的蒙特卡洛方法中深入討論。6.2.3 經(jīng)典分子動力學中近似處理實際的計算體系由于其外在條件的復雜性，往往需要對計算的系統(tǒng)進行一定的近似處理，在經(jīng)典分子動力學計算中，引進以下近似處理：v 經(jīng)典粒子相互作用，不考慮電子相互作用量子效應。v 力的作用形式，由參數(shù)可調(diào)的相互作用勢函數(shù)決定，并經(jīng)過實驗測定來驗證。v 模擬體系與實際體系相差較大，一般需要采用周期邊界來擴展計算體系。時間平均是在有限時間內(nèi)完成。【練習與思考】6-1. 查找文獻，根據(jù)上述的分子動力學的處理流程圖，編寫實現(xiàn)分子動力學的簡單程序，可參考Daan F和Berend S編著的分子模擬一書。6.

12、3 分子動力學的系綜系綜(ensemble)是指在一定的宏觀條件下(約束條件)，大量性質(zhì)和結(jié)構(gòu)完全相同的、處于各種運動狀態(tài)的、各自獨立的系統(tǒng)集合，或稱為統(tǒng)計系綜。系綜是用統(tǒng)計方法描述熱力學系統(tǒng)的統(tǒng)計規(guī)律性時引入的一個基本概念；系綜是統(tǒng)計理論的一種表述方式，系綜理論使統(tǒng)計物理成為普遍的微觀統(tǒng)計理論；系綜并不是實際的物體，構(gòu)成系綜的系統(tǒng)才是實際物體。 約束條件是由一組外加宏觀參量來表示。在平衡統(tǒng)計力學范疇下，可以用來處理穩(wěn)定系綜。6.3.1 常用系綜分類 1. 正則系綜(canonical ensemble)全稱應為“宏觀正則系綜”，簡寫為NVT，即表示具有確定的粒子數(shù)(N)、體積(V)、溫度(T

13、)。正則系綜是分子動力學方法模擬處理的典型代表。假定N個粒子處在體積為V的盒子內(nèi)，將其埋入溫度恒為T的熱浴中。此時，總能量(E)和系統(tǒng)壓強(P)可能在某一平均值附近起伏變化。平衡體系代表封閉系統(tǒng)，與大熱源熱接觸平衡的恒溫系統(tǒng)。正則系綜的特征函數(shù)是亥姆霍茲自由能。 2. 微正則系綜(micro-canonical ensemble)簡寫為NVE，即表示具有確定的粒子數(shù)(N)、體積(V)、總能量(E)。微正則系綜廣泛被應用在分子動力學模擬中。假定N個粒子處在體積為V的盒子內(nèi)，并固定總能量(E)。此時，系綜的溫度(T)和系統(tǒng)壓強(P)可能在某一平均值附近起伏變化。平衡體系為孤立系統(tǒng)，與外界即無能量交

14、換，也無粒子交換。微正則系綜的特征函數(shù)是熵。3. 等溫等壓(constant-pressure, constant-temperature)簡寫為NPT，即表示具有確定的粒子數(shù)(N)、壓強(P)、溫度(T)。其總能量(E)和系統(tǒng)體積(V)可能存在起伏。體系是可移動系統(tǒng)壁情況下的恒溫熱浴。特征函數(shù)是吉布斯自由能。 4. 等壓等焓(constant-pressure, constant-enthalpy)簡寫為NPH，即表示具有確定的粒子數(shù)(N)、壓強(P)、焓(H)。由于，故在該系綜下進行模擬時要保持壓力與焓值為固定，其調(diào)節(jié)技術(shù)的實現(xiàn)也有一定的難度，這種系綜在實際的分子動力學模擬中已經(jīng)很少遇到了

15、。 5. 巨正則系綜(grand canonical ensemble)簡寫為VT，即表示具有確定的粒體積(V)、溫度(T)和化學勢()。巨正則系綜通常是蒙特卡羅模擬的對象和手段。此時、系統(tǒng)能量(E)、壓強(P)和粒子數(shù)(N)會在某一平均值附近有一個起伏。體系是一個開放系統(tǒng)，與大熱源大粒子源熱接觸平衡而具有恒定的溫度(T)。特征函數(shù)是馬休(Massieu)函數(shù)J(，V，T)。對于不同的系綜，由于其選擇的運動方程的求解方法亦有不同。針對不同的計算體系具體的細節(jié)比較復雜，下文將以NEV和NTP這兩類簡單的系綜運用兩類不同的條件方法下求解的不同方程做簡要的介紹。6.3.2 NEV系綜基本方程NEV系

16、綜是一類較為常見的系綜。對于該系綜，經(jīng)典分子動力學方法一般采用差分近似法求解牛頓運動方程，并追蹤系統(tǒng)的時間變化?？紤]由N個粒子構(gòu)成的系統(tǒng)，假設(shè)第i個粒子在時刻位置，速度為受到的作用力為，為求出該粒子(i)在時刻的位置，經(jīng)常使用的方法是Verlet算法(對于該算法的具體細節(jié)在下文將有詳細的推導)。根據(jù)系統(tǒng)的動能及統(tǒng)計物理學中的能量均分定理，則系統(tǒng)的平衡溫度可表達為(6.1)根據(jù)此方法，在晶體中粒子的配置與排布、晶形變化等結(jié)構(gòu)相變也可以用計算機模擬進行研究。根據(jù)物理化學的基本氣體的基本方程可以知道，體系的NEV系綜下的計算體系主要變化的量有兩個壓力和溫度。下文就從恒溫和恒壓兩個條件下對NEV系綜的

17、基本方法做簡單的推導。1. 恒溫方法Nóse S和Hoover W G引入與熱源相關(guān)的參數(shù)，嘗試構(gòu)造恒溫狀態(tài)，亦即具有確定N，V，T值的系統(tǒng)，對此可設(shè)想為與大熱源接觸而達到平衡的系統(tǒng)。由于熱源很大，交換能量不會改變熱源的溫度。在兩者建立平衡后，系統(tǒng)將與熱源具有相同的溫度，系統(tǒng)與熱源合起來構(gòu)成一個復合系統(tǒng)，如圖6-3所示。這樣，系統(tǒng)中微觀粒子的動力學方程：(6.2)(6.3)(6.4)式中為系統(tǒng)的勢函數(shù)，為有效質(zhì)量，表示熱力學摩爾系數(shù)。式(6.2)至式(6.4)同往常的分子動力學方法的區(qū)別體現(xiàn)在式(6.3)中，即增加了與熱源的相互作用相關(guān)的并與力的量綱相同的一項()。與熱源相關(guān)的變化參

18、數(shù)的運動方程表明，當系統(tǒng)的總能量大于時，是增加的，從而顯示出使粒子速度減小那樣的作用；反之則顯示出使粒子速度增大；是表示與溫度控制有關(guān)的常數(shù)。圖6-3 采取擴展系方法恒溫分子動力學方法原理示意圖該方法的特征是系統(tǒng)處于微觀狀態(tài)的分布服從正則分布，若在推導式(6.2)、式(6.3)、式(6.4)中，引人由公式定義的變量，則系統(tǒng)的哈密頓量(即能量)可以表示為(6.5)式中第二項和第三項對應于與熱源相聯(lián)系的部分。這說明，由于系統(tǒng)與熱源存在熱接觸，二者可以交換能量，因此粒子的能量發(fā)生變化(H0H)，同時也說明系統(tǒng)可能的微觀狀態(tài)可具有不同的能量值。2. 恒壓方法為了進行壓力的控制，考慮如圖6-4所示的方法

19、，即利用活塞調(diào)控系統(tǒng)的體積變化。對于質(zhì)量一定的理想氣體系統(tǒng)，體積變大，其壓力就變小，反之亦然。利用這樣的原理，可以通過改變體積使壓力達到某一個預定值。就實際的模擬而言，可控制體積大小均勻(三維)地變化，將粒子的坐標、動量表示成如下形式(6.6)從而將粒子的坐標、動量與對分子動力學單元(邊長)歸一化的參變量和聯(lián)系起來。應用歸一化的參變量并根據(jù)正則方程可得到粒子運動方程。圖6-4 恒壓分子動力學方法原理示意圖(實際上，活塞為三維運動)粒子系和活塞組成復合系統(tǒng)，其哈密頓量可由下式給出(6.7)(6.8)式中，是外壓；是體積對應的共軛動量；是決定于體積變化速度的因子。若將由上述哈密頓量導出的運動方程用

20、粒子的實際坐標直接寫出來，則有(6.9)(6.10)(6.11)式(6.9)中右邊第一項是根據(jù)維里(Virial)定理求出的系統(tǒng)的瞬間壓P。在此處理中，特別重要的是式(6.9)右邊第一項與統(tǒng)計力學壓力表達式相同。6.3.3 NTP 系綜質(zhì)點系的基本方程前面已討論了關(guān)于NEV系綜的質(zhì)點系方程。鑒于NTV和NHP系綜的基本方程包含于NTP系綜的基本方程之中，在此只對NTP系綜的基本方程作一些說明。NTP系綜是在宏觀上具有確定的粒子數(shù)N、溫度T和壓力P的系綜。在此系綜中，T和P可作為外部參數(shù)給出。那么，怎樣才能滿足N，T，P恒定的條件呢？顯然，粒子數(shù)N保持一定是最簡單的，因為只要使所考察的基本單元內(nèi)

21、的粒子數(shù)恒定就行了。1. 壓力恒定體系為簡單起見，以由圓柱容器內(nèi)充滿氣體和調(diào)節(jié)活塞構(gòu)成的系統(tǒng)為例，如圖6-5所示。設(shè)活塞的重量為W，則在平衡狀態(tài)下，由活塞的重量產(chǎn)生的外部壓力()與由氣體產(chǎn)生的內(nèi)壓平衡。如果外壓變大，為了保持平衡，則氣體部分的體積將收縮，從而內(nèi)壓變大。圖6-5 活塞與氣體容器構(gòu)成的示意圖因此，把體積作為系統(tǒng)的動力學參數(shù)建立其運動方程，就好像維持壓力恒定。其方程可寫為(6.12)式中，為體積，是由維里(Virial)定理求出的內(nèi)壓，是對應于的假想質(zhì)量。圖6-5為活塞與氣體容器構(gòu)成的示意圖。將上述想法推廣到一般情況下的晶格系統(tǒng)，如圖6-6所示，對一般的晶體，其基本單元可由三個平移矢

22、量 , 決定其形狀和大小。將的三個分量寫成矩陣，該矩陣規(guī)定了晶體基本單元的形狀和大小，這樣也成為原子的實際坐標和晶格矢量之間的變換矩陣，亦即圖6-6 基于單元形狀示意圖在式(6.15)中，體積作為動力學參數(shù)，而在一般晶體系統(tǒng)中，應該作為動力學參數(shù)。若給出此動力學系統(tǒng)的哈密頓量，則有(6.13)式中，是與格矢相應的原子共扼動量，為數(shù)值行列式，中是多粒子體系的勢函數(shù)，是中與體積V相應的共扼動量，是基本單元的體積，且由給出，為相應的質(zhì)量，為外部壓強，代表各向異性應力張量。2. 溫度恒定體系接下來討論滿足溫度恒定的基本方法。恒溫體系比定壓體系要抽象，簡單地給出視覺化圖像是困難的。勢能在原子坐標和動量之

23、外引入了附加自由度，并由下式用把現(xiàn)實系統(tǒng)與具有能量確定的(微正則系統(tǒng))假想系統(tǒng)聯(lián)系起來。(6.14)式中，為現(xiàn)實系統(tǒng)的時間，為假想系統(tǒng)的時間。該假想(力學)系統(tǒng)的哈密頓量可以寫成(6.15)式中，是在假想系統(tǒng)中原子的動量，是對應的共扼動量，是由公式給出的假想系統(tǒng)的自由度數(shù)目，二是熱浴的溫度，為玻爾茲曼常數(shù)。顯然，式(6.18)中右邊第一項代表粒子的動能，第二項為粒子系統(tǒng)的勢能，第三項是熱浴的假想動能，第四項是熱浴的勢能；因為在假想系統(tǒng)中，力學系綜是微正則分布的，所以其配分函數(shù)變成具有總能量 E 的狀態(tài)總數(shù)(6.16)(6.17)若將此配分函數(shù)積分變換到現(xiàn)實系統(tǒng)中，則可得到在恒定溫度T下的正則分

24、布(也稱為玻爾茲曼分布)的配分函數(shù) (6.18)式中，是由下式給出的現(xiàn)實系統(tǒng)的哈密頓量(6.19)是假想系統(tǒng)的守恒量，而對真實系統(tǒng)來講，由于存在熱交換將不是守恒量。為了判斷計算結(jié)果是否合理，只要簡單地檢驗是否保持恒定(守恒)就行了。3. NTP系綜到此，已經(jīng)導出了壓力恒定系綜和溫度恒定系綜的哈密頓量，將兩個系綜的哈密頓量組合起來就可以得到NTP系綜的哈密頓量，其模型如圖6-7所示。NTP系綜的哈密頓量為 (6.20)圖6-7 NTP系統(tǒng)示意圖若根據(jù)此哈密頓量導出正則方程式，并用速度置換動量，則可得到各變數(shù)的運動方程(6.21)(6.22)(6.23)式中，變量上方的圓點表示對時間的微分，而分別

25、由下列各式給出(6.24)(6.25)(6.26)4. NPT系綜分子動力學求解若在三維周期邊界條件下求解上述方程，則可模擬固體的結(jié)構(gòu)相變、溶解現(xiàn)象和晶化過程。各變量的運動方程式(6.21)，式(6.22)和式(6.23)是非線性常微分方程，不能用貝魯勒(Berrele)法進行求解。因此，對這些方程通常用牛頓迭代方法求解。而問題的關(guān)鍵在于確定式(6.22) 和式(6.23)中的假想質(zhì)量。所考察基本單元內(nèi)的粒子體系存在著由粒子作用勢以及粒子數(shù)與粒子質(zhì)量共同決定的固有壓力和溫度等的起伏變化。若想再現(xiàn)這些波動周期那樣給出和，則由體積變化或溫度變化引起的粒子系統(tǒng)對平衡態(tài)的偏離將迅速回復而變成穩(wěn)定狀態(tài)。

26、因此，若設(shè)此周期分別為和，則和可以表示為(6.27)(6.28)式中，為溫度，為玻爾茲曼常數(shù)，為原子數(shù)，為基本單元的線尺度，為壓縮率，和分別根據(jù)對NEV系綜的計算結(jié)果而定?！揪毩暸c思考】6-2. 推演NEV系綜的牛頓方程，哈密頓方程，說明各個物理量的意義，在那些實際的計算體系中得到應用。6-3. 根據(jù)NPT系綜的推導方法，對NPH系綜的哈密頓方程做簡單的推導。6-4. 詳細的分析不同的系綜的使用范疇，給出各個系綜的使用條件，列表進行簡單的比較。6.4 原子勢函數(shù)和分子力場構(gòu)造分子動力學中存在著兩個基本的概念，一是在前面提到的系綜的問題，它所確定的是模擬系統(tǒng)的大小，規(guī)格和相關(guān)的模擬條件等信息；二

27、是本節(jié)需要講述的勢函數(shù)。簡單的講，勢函數(shù)就是用來描述原子和原子之間相互作用的。隨著勢函數(shù)的不斷發(fā)展，現(xiàn)有的勢函數(shù)主要分成了兩個大的部分。即是經(jīng)典的勢函數(shù)和電子理論的勢函數(shù)。如圖6-8所示，經(jīng)典分子動力學的基本的勢函數(shù)的一個簡單分類。圖6-8 經(jīng)典的分子動力學相互作用勢分類作用勢的選擇與動力學計算的關(guān)系極為密切，選擇不同的作用勢，體系的勢能面會有不同的形狀，動力學計算所得的分子運動和分子內(nèi)部運動的軌跡也會不同，進而影響到抽樣的結(jié)果和抽樣結(jié)果的勢能計算，在計算宏觀體積和微觀成分關(guān)系的時候主要采用剛球模型的二體勢，計算系統(tǒng)能量，熵等關(guān)系時早期多采用Lennard-Jones、Morse勢等二體勢模型

28、，對于金屬計算，主要采用Morse勢，但是由于通過實驗擬合的對勢容易導致柯西關(guān)系，與實驗不符，因此在后來的模擬中有人提出采用EAM等多體勢模型，或者采用第一性原理計算結(jié)果通過一定的物理方法來擬合二體勢函數(shù)。但是相對于二體勢模型，多體勢往往缺乏明確的表達式，參量很多，模擬收斂速度很慢，給應用帶來很大的困難，因此在一般應用中，通過第一性原理計算結(jié)果擬合勢函數(shù)的Lennard-Jones，Morse等勢模型的應用仍然非常廣泛。6.4.1 經(jīng)典理論的原子勢函數(shù)在分子動力學的框架之下，早期的勢函數(shù)是主要以描述原子和原子之間的相互作用而產(chǎn)生和發(fā)展起來的。早在1903年的時候，Mie G就研究了兩個粒子之間

29、的相互作用勢，他給出的勢函數(shù)是由兩項構(gòu)成的，一項代表原子之間的相互吸引，另一項則代表著原子之間的相互排斥。這就是非常著名的Lennard-Jones勢函數(shù)的原始思考。對于其他有關(guān)勢函數(shù)的發(fā)展在這里就不做過多的贅述了。從圖6-9可以清楚的知道原子之間的相互作用可以分為四類：對勢、對泛函勢、組合勢、組合泛函勢。下面就對相關(guān)的勢函數(shù)的基本的發(fā)展做一個簡單的總結(jié)。表6-1 原子中相互作用勢函數(shù)勢函數(shù)作用因素典型代表對勢僅僅決定于原子的坐標L-J勢、BMH勢、GK勢對泛函勢多體相互作用EAM勢有效介質(zhì)理論組合勢共價鍵的相互作用SW勢組合泛函勢Abell-Tersoff勢CheliKowsky-Phill

30、ips勢6.4.2 分子力場分子力場根據(jù)量子力學的波恩奧本海默近似，一個分子的能量可以近似看作構(gòu)成分子的各個原子的空間坐標的函數(shù)，簡單地講就是分子的能量隨分子構(gòu)型的變化而變化，而描述這種分子能量和分子結(jié)構(gòu)之間關(guān)系的就是分子力場函數(shù)。分子力場函數(shù)為來自實驗結(jié)果的經(jīng)驗公式，可以講對分子能量的模擬比較粗糙，但是相比于精確的量子力學從頭計算方法，分子力場方法的計算量要小數(shù)十倍，而且在適當?shù)姆秶鷥?nèi)，分子力場方法的計算精度與量子化學計算相差無幾，因此對大分子復雜體系而言，分子力場方法是一套行之有效的方法。以分子力場為基礎(chǔ)的分子力學計算方法在分子動力學、蒙特卡羅方法、分子對接等分子模擬方法中有著廣泛的應用。

31、 1. 力場主要討論的對象及其構(gòu)成：以下參數(shù)構(gòu)成一套力場函數(shù)體系需要有一套聯(lián)系分子能量和構(gòu)型的函數(shù)，還需要給出各種不同原子在不同成鍵狀況下的物理參數(shù)，比如正常的鍵長、鍵角、二面角等，這些力場參數(shù)多來自實驗或者量子化學計算。 v 鍵伸縮能：構(gòu)成分子的各個化學鍵在鍵軸方向上的伸縮運動所引起的能量變化；v 鍵角彎曲能：鍵角變化引起的分子能量變化；v 二面角扭曲能：單鍵旋轉(zhuǎn)引起分子骨架扭曲所產(chǎn)生的能量變化； v 非鍵相互作用：包括范德華力、靜電相互作用等與能量有關(guān)的非鍵相互作用； v 交叉能量項：上述作用之間耦合引起的能量變化。2. 常用力場函數(shù)和分類不同的分子力場會選取不同的函數(shù)形式來描述上述能量與

32、體系構(gòu)型之間的關(guān)系。到目前，不同的研究小組設(shè)計了很多適用于不同體系的力場函數(shù)，根據(jù)他們選擇的函數(shù)和力場參數(shù)，可以分為以下幾類：力場CHARMM力場CVFF力場MMX力場AMBER力場COMPASS力場MMF94力場CFF力場UFF力場Dreiding力場ESFF力場第一代力場第二代力場第三代力場圖6-9 分子力場的發(fā)展和分類v 傳統(tǒng)力場(第一代力場)a) AMBER力場：由Kollman課題組開發(fā)的力場，是目前使用比較廣泛的一種力場，適合處理生物大分子；b) CHARMM力場：由Karplus課題組開發(fā)，對小分子體系到溶劑化的大分子體系都有很好的擬合；c) CVFF力場：CVFF力場是一個可以

33、用于無機體系計算的力場；d) MMX力場：MMX力場包括MM2和MM3，是目前應用最為廣泛的一種力場，主要針對有機小分子。v 第二代力場 第二代的勢能函數(shù)形式比傳統(tǒng)力場要更加復雜，涉及的力場參數(shù)更多，計算量也更大，當然也相應地更加準確；a) CFF力場CFF力場是一個力場家族，包括了CFF91、PCFF、CFF95等很多力場，可以進行從有機小分子、生物大分子到分子篩等諸多體系的計算；b) COMPASS力場由MSI公司開發(fā)的力場，擅長進行高分子體系的計算；c) MMF94力場Hagler開發(fā)的力場，是目前最準確的力場之一。v 通用力場(第一代力場)通用力場也叫基于規(guī)則的力場，它所應用的力場參數(shù)

34、是基于原子性質(zhì)計算所得，用戶可以通過自主設(shè)定一系列分子作為訓練集來生成合用的力場參數(shù)；a) ESFF力場MSI公司開發(fā)的力場，可以進行有機、無機分子的計算；b) UFF力場可以計算周期表上所有元素的參數(shù)；c) Dreiding力場適用于有機小分子、大分子、主族元素的計算。3. 分子力場的計算對于分子之間的相互作用通常需要將計算體系的所有粒子納入考慮范疇之內(nèi)，那么分子之間的相互作用可以作為長程作用來處理。在長程相互作用中又分為周期性長程相互作用力和非周期性相互作用力。周期性的長程作用力通?？梢岳脽o網(wǎng)格算法，包含Ewald求和算法等；而對于非周期性的相互作用通常使用快速多極子方法計算。下面針對這

35、兩個最為常見的算法做簡單的敘述。v Ewald求和算法Ewald求和算法是Ewald于1921年提出，最初是以計算離子晶體的能量。此方法選定一個計算系統(tǒng)的盒子，盒子中的粒子除了和盒子內(nèi)部的粒子相互作用之外，還會與其鏡像系統(tǒng)中的粒子發(fā)生作用，并且鏡像系統(tǒng)和盒子內(nèi)部的系統(tǒng)完全相同。圖顯示的是計算的相互作用的系統(tǒng)，中央的黑色的區(qū)域表示的是計算的區(qū)域，周圍的鏡像隨著系統(tǒng)的距離逐漸變淡，那么這樣的一個極限作用區(qū)域就像一個球，這個球稱為Ewald球。假設(shè)作用的系統(tǒng)的邊長是，含有個相互作用的粒子，其鏡像的作用位置可以表示為()；其中都是正整數(shù)。v 快速多極子法(The Fast Multipole Meth

36、od)在一個開放的系統(tǒng)中，不是對所謂的晶格進行求和，而是對要求計算的整個體系的粒子對進行求和。那么這樣的一個系統(tǒng)不同位置粒子的靜電能可以簡單的表述為(6.29)因此這個體系的靜電作用和前文所講述的Lennard-Jones作用項具有相同的形式?？焖俣鄻O子法將需要計算的體系分為多個體積相同的格胞，由格胞內(nèi)部的所有的原子計算每一個格胞的多極(電荷，偶極，四極)。而格胞和格胞之間的相互作用是利用中央多極展開的方法進行計算。因此中央多極展開的方法的條件是格胞之間的距離必須是大于一定距離的，因為大于一個格胞以上的以多極展開的方式計算相互作用的力，而格胞內(nèi)部的原子對是以成對原子的作用力來計算的?！揪毩暸c思

37、考】6-5. 查找文獻，推導Lennard-Jones勢函數(shù)的公式。6-6. 推導EAM勢函數(shù)的基本公式，并且實現(xiàn)簡單的算法，編寫基本設(shè)計程序。6-7. 以對勢中的Lennard-Jones勢函數(shù)為例，寫出對勢中勢函數(shù)的基本算法，并且進行基本的程序設(shè)計。6.5 數(shù)值方法與編程技巧6.5.1 邊界條件隨著材料的低維化和納米化，材料的表面和邊界對于其性質(zhì)的影響愈來愈大。在處理物質(zhì)表面、界面，以及塊體物質(zhì)時，對其周期邊界條件問題的考慮就變得非常必要了。對于周期性邊界條件可以分為一維、二維及三維的情況。介紹薄膜(可使用二維周期邊界條件)和擴大到半無限的表面情況，然后討論液體與晶體、晶體與晶體之間形成的

38、界面等，最后學習塊狀物質(zhì)的三維周期邊界條件。1. 薄膜的邊界條件(二維周期邊界條件)對于薄膜情況，如圖6-10所示賦予其周期性邊界條件?？梢哉J為薄膜在x-y平面內(nèi)無限擴展，存在周期邊界條件，而在 z方向受到限制(不賦予周期邊界條件)。事實上，二維周期性邊界條件也適用于在z方向不具有真正意義上的自由度的二維物質(zhì)。但就物質(zhì)科學的研究范疇而言，多數(shù)情況下是指在某種襯底上生成的原子或分子的薄膜。此處討論的邊界條件，除針對薄膜以外，也可近似地用于下面將要講到的在 z 的負方向半無限地擴展的表面問題，例如表面重構(gòu)等。但是，當在z方向存在異質(zhì)原子層時，原來塊體物質(zhì)在異質(zhì)原子層附近將有近于表面的性質(zhì)，這種由塊

39、體向表面的過渡特征取決于原子層的厚度，一般要求原子層厚度在力的中斷距離以上(作用力程以上)。在處理實際問題時，只要嘗試進行幾次改變原子層厚度的模擬和計算，并檢驗計算結(jié)果是否變化就行了。圖6-10 薄膜周期性邊界條件示意圖2. 表面邊界條件(半無限的情況)對于表面存在表面重構(gòu)的情況，是與薄膜的情況不同的，其在 z 軸方向擴展至半無限。處理此類系統(tǒng)之邊界條件時，只對x-y平面和 z 的負方向賦予相應的周期性邊界條件，其示意圖如圖6-11所示。但是，在將基本單元原樣地復制中，必須考慮到這樣的一種情況，那就是在距薄膜結(jié)構(gòu)的基本單元最下層的相鄰底層被配置了基本單元表面層的復制品，從而導致了在本來應具有塊

40、體物質(zhì)性質(zhì)的基本單元的最下層的結(jié)構(gòu)中出現(xiàn)了表面層的強烈影響。所以將基本單元對z軸反轉(zhuǎn)(即鏡像映射)進行復制后使用。然而，即使在這種情況下，考慮到基本單元在z方向的厚度對計算結(jié)果有影響，所以應當進行與前面所述同樣的檢驗。圖6-11 采用映射方法半無限表面的周期性邊界條件圖6-12 固定層法半無限表面的周期性邊界條件在利用分子動力學方法對分子束外延的晶體生長模擬時，要充分認識到分子束外延 (MBE)是一個原子數(shù)不斷增加的非平衡系。所謂MBE技術(shù)就由束源產(chǎn)生的原子蒸發(fā)到襯底上，從而形成晶體。因為由束源出來的原子具有動量，所以碰撞襯底，襯底的重心就變得下移，為了克服這種現(xiàn)象的影響，對分子束外延即在z軸

41、方向不賦予周期性邊界條件，而采用通過調(diào)整配置使原子位置坐標不變的所謂固定原子層方法，然后使用邊界條件，如圖6-12所示。固定(原子)層還可理解為具有無限大質(zhì)量的原子層。來自固定層中原子作用力的貢獻可作為外力處理。為了基本單元最下層具有塊體物質(zhì)的性質(zhì)，固定層中的原子層數(shù)目必須足夠多。采用鏡像映射的方法可自然地處理體系伴隨溫度變化而產(chǎn)生的膨脹和收縮。但就使用固定層的方法而言，基本單元發(fā)生變化的同時，必須將固定層的晶格常數(shù)作為外部參量變化來考慮?！揪毩暸c思考】6-8. 通過對比映射方法和固定層法的不同，理解在分子動力學中的二維邊界條件的設(shè)置方法，思考在線形的一維約束條件的邊界的處理方法。6-9. 根

42、據(jù)界面的邊界條件的設(shè)置方法，分析討論在分子動力學中如何構(gòu)造線缺陷模型，做要的圖示和語言表述，說明界面的條件的選取方法。6-10. 根據(jù)上述的邊界的處理方法，試思考塊體材料的邊界條件的實現(xiàn)方法，并且給出簡單的實例。6.5.2 時間積分處理計算機模擬方法的基點是利用現(xiàn)代計算機高速和精確的優(yōu)點，對幾百個以至上千分子的運動方程進行數(shù)值積分。有許多不同的積分方法，它們的效率和方便程度各異?；締栴}就是用有限差分法來對二階常微分方程進行積分。而對于分子動力學而言其基本的工作也就是對運動方程進行時間上的積分。這樣需要將運動方程離散化為有限的差分方程，即是將時間離散化為有限的格點。主要的方法是利用泰勒展開采取

43、不同的截斷半徑。而相鄰的格點之間的距離就是所謂的時間步長，這樣在知道了原始的位置坐標之后，在知道某些物理量對于時間t的微分量之后就可以很容易確定下一個時間體系的坐標。這里很容易出現(xiàn)所謂的兩類誤差，即是截斷誤差和舍入誤差。為此我們需要在實際的模擬中，應該根據(jù)不同的問題來確定所需要的精度，為此發(fā)展了幾種時間積分的方法。1. Verlet方法Verlet法經(jīng)常用于NEV系綜的分子動力學微分方程的數(shù)值積分求解。若將在和的位置坐標分別用時刻t的位置坐標進行泰勒展開，則有(6.30)(6.31)為簡單起見，以后均省略有關(guān)原子標識符號I，V(t)為速度，為了消去速度項，將式(6.30)和式(6.31)相加，

44、并將r()移至右邊，同時利用，則有(6.32)上式(6.32)是差分公式，其中()的奇次項全被抵消，同時忽略及其他高次項。差分精度(誤差)為0()。若由式(6-36)減去式(6-37)則可得到速度表達式(6.33)因為是消去了關(guān)于的偶數(shù)方項，所以速度的誤差為的量級。上述式(6.32)右邊第一、二項均是的量級項，而第三項是量級項。這樣若將大的絕對值項和小的絕對項混合，則隨著重復計算次數(shù)的增多，其誤差發(fā)生積累。為防止誤差積累，Hockney提出了所謂蛙躍(Leap-frog)法。2. 蛙躍法首先，利用速度將二階常微分方程轉(zhuǎn)換為二個一階常微分方程。亦即(6.34)速度的微分由時刻，和的差分給出(6.

45、35)根據(jù)此式，在時刻的速度給出如下(6.36)另一方面，利用時刻和時刻的差分可給出位置矢量的微分(6.37)根據(jù)上式可給出時刻的位置矢量表達式(6.38)上式(6.38)所表示的意義是速度和位置矢量的變化依賴步長。利用時刻的位置坐標可計算出力，利用此力的結(jié)果求出時刻()的速度，進而，利用此速度求出在()時刻的位置坐標，這樣就完成了第一步的計算。正如以上所述，如果使用Verlet法或蛙躍法，則可求解關(guān)于NEV系綜的線性常微分方程。但是，在此之外的NTV，NHP，NTP等系綜給出的是非線性常微分方程。所以用Verlet法或蛙躍法不能進行求解。關(guān)于這些非線性常微分方程的求解要采用下面介紹的所謂Ge

46、ar方法。Hockey提出的Leap-frog算法是Verlet算法的變化，這種方法涉及半時間間隔的速度。這種算法與Verlet算法相比有兩個優(yōu)點：(1)包括顯速度項；(2)收斂速度快，計算量小。這種算法明顯的缺陷是位置和速度不同步。3. Gear方法因為原子軌跡是時間的連續(xù)函數(shù)，所以可將()時刻的位置、速度、加速度等對時刻泰勒展開(6.39)若使用給出的新的位置矢量，則可計算()時刻的力F()和加速度()。根據(jù)加速度可評價預測的加速度的誤差，亦即(6.40)若將此差值加于預測項上，則可得到修正項(6.41)(6.42)4. Velocity-Verlet算法這種算法可以同時給出位置、速度和加

47、速度，并且不犧牲精度，給出了顯速度項，計算量適中，目前應用比較廣泛。(6.43)5. Gear的預測-校正算法這種算法分三步來完成。首先，根據(jù)泰勒展開，預測新的位置、速度和加速度。然后，根據(jù)新的計算的力計算加速度。最后這個加速度再由與泰勒級數(shù)展開式中的加速度進行比較，兩者之差在校正步里用來校正位置和速度項。這種方法的缺點就是占有計算機的內(nèi)存大。【練習與思考】6-11. 查找文獻，總結(jié)分子動力學中現(xiàn)有的時間算法，列表說明各種不同的時間算法優(yōu)勢。6-12. 推導出Verlet算法的基本方程，并編程實現(xiàn)。6-13. Hockey提出的Leap-frog算法是Verlet算法的變化，這種方法涉及半時間

48、間隔的速度，查找原始文獻，推導出Leap-frog算法的基本方程，結(jié)合練習6-12中的問題，修改程序?qū)崿F(xiàn)Leap-frog算法。6-14. 根據(jù)文中的敘述，詳細的推導出Gear算法的基本方程，并編程實現(xiàn)。6-15. 查閱文獻，分析討論時間可逆算法劉維公式在分子動力學中應用。6.6 計算結(jié)果的解析方法總體來說分子動力學(MD)目前能夠討論的問題已經(jīng)可以涉及到比較廣闊的范疇，從基本的模擬體系的位形，到現(xiàn)在的各種材料的光譜學性能的計算。但是總體來說MD的計算結(jié)果僅僅分了兩個層次。第一個是計算模擬直接的輸出，基本上包含了在平衡態(tài)的體系位形信息和體系運動的基本動力學信息；第二個是基于現(xiàn)代統(tǒng)計熱力學的平均

49、統(tǒng)計得到的一些宏觀體系的綜合信息，包含了光學，電學，磁學和光譜學等等。但是，在材料科學中，不僅要很好地知道原子水平上的三維結(jié)構(gòu)，而且還必須計算出表征宏觀物性(熱力學性質(zhì)、光學性質(zhì)等)的物理量值。為了靈活應用分子動力學方法，下面對與這些物性相關(guān)的物理參數(shù)值的計算方法作介紹。6.6.1 宏觀參量的統(tǒng)計平均由分子動力學基本方程出發(fā)進行計算的輸出結(jié)果是關(guān)于原子的位置坐標和速度。物性參數(shù)值可根據(jù)位置坐標和速度通過統(tǒng)計處理求出。在統(tǒng)計物理中，可利用系綜微觀量的統(tǒng)計平均值來計算物性參數(shù)值，亦即： (6.44)式中，ens是系綜(Ensemble)的英文縮寫，是分布函數(shù)(幾率密度)。幾率密度因所考察的系綜不同

50、而不同。例如，對于微正則系綜(NEV系)，可表達為(6.45)而對于正則系綜(NTV系)，等溫、等壓系綜(NPT系)和巨正則系綜，分別由下列各式給出：(6.46)系綜之間的演變可由相應的公式聯(lián)系起來。首先，配分函數(shù)Z和熱力學勢，分別由下列式子給出(6.47)式中，在(,)中，設(shè)定，則可有以下組合：NTV和NEV系綜之間的轉(zhuǎn)變NTP和NTV系綜之間的轉(zhuǎn)變VT和NTV系綜之間的轉(zhuǎn)變在上面最后一個演變例子中，關(guān)于N的積分應理解為對N0，1，的求和。在NEV系綜情況下，表示所有狀態(tài)數(shù)W，熱力學勢為(-TS)；對于NTP和VT系綜，其熱力學勢分別為N和(-PV)。此時，物理量在系綜中的平均值的變換可由下

51、式給出(6.48)在熱力學極限情況下，上述各系綜的平均是相等的(6.49)到此，容易發(fā)現(xiàn)，有了上述這些分布函數(shù)，代入式(6.50)進行積分運算，理論上講，可以求出各種各樣的物理量。然而，在通常情況下，式(6.50)不能給出真正的解析解。在理想的剛體球勢情況下，可以進行解析求解，而對于連續(xù)函數(shù)勢，只能是數(shù)值求解。不難看出，這些解始終都是近似解。從歷史的角度來說，正是為了超越上述積分方程帶來的限制，人們引人了所謂蒙特卡羅方法和分子動力學方法。在分子動力學方法中，人們使用了所謂的時間平均等于系綜平均的各態(tài)歷經(jīng)假設(shè)，亦即(6.49)在熱力學統(tǒng)計物理學中，各態(tài)歷經(jīng)假設(shè)沒有嚴格地證明，它的正確性是由它的推

52、論被實驗結(jié)果所證實而被肯定。式(6.59)中由下式給出 (6.50)在氣體動力論中，原子的速度分布為Maxwell分布。在統(tǒng)計力學中，瞬間的熱力學性質(zhì)可以就其平均值做Gaussian分布描述。所以，對任意一個狀態(tài)函數(shù)X(t)，其瞬間的值是作無規(guī)地起伏變化，而此狀態(tài)值X(t)則依據(jù)它的平均值定義了其分布函數(shù)。但是要使得性質(zhì)分布函數(shù)X(t)的取樣點為可幸賴的，則模擬的向空間長度要多長？在此，先來分辨三種形式的統(tǒng)計平均值：軌跡平均(trajectory average)、取樣平均(sample average)、取樣平均的期望值(the expected value of the sample av

53、erage)。軌跡平均：考慮一個獨立系統(tǒng)的平衡相空間軌跡，并且想象在軌跡上取個不同的取樣點，這些取樣點是均勻分布在軌跡在線，相隔時間的相鄰點。由于軌跡是代表所有的平衡分布，而且很小，所以是一個非常大的值。對一個與時間有關(guān)的性質(zhì)物理量，時間平均值可取為軌跡平均值是我們想要知道的量。對于在任何時刻的相對于平均值的變異數(shù)為(6.51)取樣平均：雖然我們想要知道系統(tǒng)性質(zhì)的平均值(也就是軌跡平均值)，但是卻又不希望對所有軌跡上的做計算。所以找了一個取代的方法：在軌跡上取一小段，并在此段上取M點做計算(M)。也就是在一小段軌跡上取連續(xù)相鄰的M個點，使成一個集合，此集合稱為樣本。這個樣本有它自己的平均值(6

54、.52)而在軌跡上任意時刻的，可針對此平均值做變異數(shù)：(6.53)但是由軌跡中某段所得的平均值是否可代表整個軌跡的平均值？這個問題等介紹完第三個平均值后就可討論。取樣平均的期望值：自整個軌跡上取個片段軌跡，每段軌跡都包含M個點，任意其中一個取樣段(例如第k個)，都有它自己的平均值，所以共有個平均值。取這些個平均值再做平均稱為平均值的期望值，則這些取樣段平均值會對做性質(zhì)函數(shù)分布，稱為取樣分布(sampling distribution)。對任意一個取樣段的平均值相對平均值的期望值的變異數(shù)為(6.54)由以上三種取樣(統(tǒng)計)方法，可得到取樣理論的一個重要結(jié)果：期望值的變異數(shù)與軌跡變異數(shù)S成比例關(guān)系

55、(6.55)所以，若對整個軌跡取樣，則M=，因此，任一小段的平均值(近似期望值)是等于軌跡平均值。另外，若M，則(6.56)取任一段的平均值去估計軌跡平均值得準確度是與任一小段上的取樣點數(shù)M有關(guān)，而與全軌跡長度的無關(guān)。直覺上，我們以為在取軌跡中的一段做平均值時，只有此段長度愈長(也就是)，其平均值才會愈接近全軌跡的平均值。但由上式卻看到，當非常大時，是不重要的。這也就是為什么取1500人做調(diào)查的結(jié)果，可代表數(shù)千萬人的分布。將取樣理論應用在分子動力學：在分子動力學中，常常針對時間及粒子數(shù)目做性質(zhì)函數(shù)的平均值。通常所謂的巨觀尺度，就時間而言，數(shù)量級為秒；就粒子數(shù)目而言，數(shù)量級為。但是在做分子動力模

56、擬時，我們卻可以僅做皮秒及數(shù)百個分子的計算即可?！揪毩暸c思考】6-16. 搜索利用分子動力學方法計算之后的位形顯示軟件，并且列表比較其優(yōu)勢。6-17. 根據(jù)文中提到的配分函數(shù)概念，利用熱力統(tǒng)計學的基礎(chǔ)知識，完成定壓比熱容CP，定容比熱容Cv，熱膨脹系數(shù)，等溫壓縮率參數(shù)的公式推到，并利用程序來實現(xiàn)。6-18. 根據(jù)文獻查找結(jié)果，編程完成一般分子動力學方法中輸運系數(shù)的計算。6-19. 根據(jù)擴散系數(shù)的基本方程，通過查找參數(shù)模擬在Si體系中三個的H原子的擴散的路徑，并編程實現(xiàn)。6.7 分子動力學計算實例6.7.1 CNT的幾何特性自1991年被發(fā)現(xiàn)以來，碳納米管因其所特有的優(yōu)異物理、化學和機械性能及其巨大的潛在應用價值而得到了廣泛的關(guān)注。例如根據(jù)碳納米管獨特的一維管狀結(jié)構(gòu)和極大的長徑比可將其制成高效的傳質(zhì)單元；超薄的針形尖端使碳納米管可用做掃描隧道顯微鏡及原子力顯微鏡的探針；碳納米管具有較強的場發(fā)射性能，是一種優(yōu)良的場效應晶體管材料，Duan等最近在研究中發(fā)現(xiàn)，用金屬銫對單壁碳納米管進行修飾能加強場發(fā)射性能；碳納米管兼具金屬或半導體導電性，是一種理想的電極材料；碳納米管的表面積較大，可應用于儲氫、儲能以及吸附劑等領(lǐng)域。如上所述，碳納米管的發(fā)現(xiàn)進一步拓寬了人們對納米材料的認知領(lǐng)域，隨著制備與