小知識：分子動力學基本原理及應用

為了微觀模擬體系能夠反映宏觀實驗現象, 需要通過周期性邊界條件對模擬對象體系進行周期性復制, 以避免在實際中并不存在的邊緣效應(edge effects)。原則上，對于任何分子體系的理論研究都需要求解含時薛定諤方程。但在實際工作中，更關注的是原子核的運動軌跡，這樣的軌跡可以利用波恩-奧本海默近似(Born-Oppenheimer approximation)，通過求解經典力學運動方程獲得。Alder和Wainwright曾表示計算機模擬實驗會成為聯系宏觀實驗現象和微觀本質的重要橋梁，在他們首次進行分子動力學模擬實驗之后10 年，法國物理學家Verlet提出了一套牛頓運動方程的積分算法，與此同時提出的還有另一套產生和記錄成對近鄰原子的算法，從而大大簡化了原子間相互作用的計算。這兩種算法至今仍以一些變形的形式在實踐中被廣泛應用^[1,2]。

過去幾十年開發了多種原子級模擬方法，包括晶格靜力學、晶格動力學、蒙特卡羅和分子動力學等。其中，分子動力學特別適用于塑性變形的研究，它通過一些規定的原子間相互作用勢函數的原子相互作用系統的牛頓方程的解，研究變形過程中的實時行為，并包括晶格的非簡諧性、內應力的高度不均勻，以及系統的瞬態響應等方面的影響。

分子動力學主要依靠牛頓力學來模擬分子體系的運動，以在由分子體系的不同狀態構成的系統中抽取樣本，從而計算體系的構型積分，并以構型積分的結果為基礎進一步計算體系的熱力學量和其他宏觀性質。它對原子核和電子構成的多體系統，求解運動方程，是一種能夠解決大量原子組成的系統動力學問題的計算方法，不僅可以直接模擬物質的宏觀演變特性，得出與試驗結果相符合或相近的計算結果，還可以提供微觀結構、粒子運動以及它們和宏觀性質關系的明確圖像，從而為新的理論和概念的發展提供有力的技術支撐。

分子動力學的對象是一個粒子系統，系統中的原子間的相互作用用勢函數來描述，因此，正確選擇勢函數的類型及其參數，對于模擬的結果優劣具有重要作用。勢能函數在大多數情況將描述分子的幾何形變最大程度地簡化為僅僅使用簡諧項和三角函數來實現；而非鍵合原子之間的相互作用，則只采用庫侖相互作用和Lennard-Jones 勢相結合來描述。其中，對于原子間相互作用力的描述通常是經驗或半經驗的，這樣雖然能夠提高計算效率，但無法完全揭示電子鍵合的多體性質，尤其對于缺陷附近與自身結構和化學性有關的復雜自洽變分函數。Daw和Baskws的EAM（Embedded-atom model）勢函數在某種程度上融合了電子鍵合的多體性質，將系統的總勢能表示為：

式中：Fi是原子i的嵌入能函數；ρi是除第i個原子以外所有原子在i處產生的電子云密度之和；Φij是第i個原子與第j個原子之間的對勢作用函數；rij是第i個原子與第j個原子之間的距離^[1]。

勢函數的可靠性主要取決于力場參數的準確性，而力場參數可以通過擬合實驗觀測數據和量子力學從頭算數據得到。目前在生物大分子體系模擬中使用最為廣泛的分子力場是CHARMM力場和AMBER力場，是早期研究生物大分子的分子力場。現有的力場參數仍在不斷優化之中，并且涵蓋的分子類型也在不斷擴大。粗粒化(coarse-grained)模型在計算生物物理研究中越來越引起人們的關注，由于在該模型中定義了粗粒化粒子，對應于全原子模型中的若干原子或原子基團甚至分子，減少了體系中的粒子數，使得模擬的時間和空間尺度得以大幅度提高，但同時也將丟失原子細節信息。基于這種模型的分子動力學模擬適合于研究緩慢的生物現象或者依賴于大組裝體的生物現象。

設計一個基本力場的根本原則就是要單位時間步長(time step)內計算能量的開銷最小化，從而實現模擬尺度的最大化。這一點對于全原子力場尤為重要, 即便對于所謂的粗粒化模型也是一樣。特別地，如果想要進行微秒甚至毫秒級時間尺度的模擬，這一原則極其重要。

圖1顯示了分子動力學的時間和空間尺寸的反比關系，圖中從左至右依次為：(1)水，細胞的基本成分；(2)牛胰蛋白酶抑制體，一種酶，其“呼吸”行為可以在毫秒級時間尺度上進行考察；(3)核糖體，一個復雜的生物器件，可以對基因信息進行解碼并生產蛋白質；(4)紫細菌光合膜片段，擁有2500 萬個原子，圖中顯示的是嵌在磷脂雙分子層上的捕光復合物及光化學反應中心^[2]。

圖1 經典分子動力學的時間和空間的尺度關系

隨著計算機處理器速度的快速增長以及大規模并行計算架構的發展，大規模并行化或專用的架構技術與可擴展分子動力學程序的結合，計算機模擬包含從位錯到基于晶界的變形機制的整個晶粒尺寸范圍，為探索材料體系的研究前沿領域開辟了新的途徑。

例如，William Gon?alves等人通過選用Wolf BKS（van Beest, Kramer and van Santen）勢函數來描述原子之間的相互作用力，使用大規模原子/分子并行模擬器LAMMPS（Large-scale Atomic/Molecular?Massively Parallel Simulator）對二氧化硅氣凝膠的彈性和強度進行了分子動力學方面的研究，他們使用velocity-Verlet算法和1.0 fs時間步長，并且在三個方向上均使用周期性邊界條件^[3]。

圖2為模擬的超過7000000個原子的大體積樣品3D示意圖，以及20nm厚的樣品切片和局部放大圖（藍色為氧原子，紅色為硅原子），圖3（a）是對803nm³氣凝膠樣品進行單軸拉伸試驗，以獲得300K的應力應變曲線，（b-d）是典型韌性斷裂圖像，（e）抗拉強度與樣品體積的對數關系。他們分析認為，為了確保正確評估彈性等機械性能，模擬樣品的尺寸至少為孔徑的8倍，同時，表面積極高的二氧化硅氣凝膠需要相對低的應變率以確保準靜態條件。

圖2 模擬的二氧化硅氣凝膠樣品（超過七百萬個原子）

圖3 單軸拉伸試驗的應力應變曲線（a）、強度與體積的關系（e）和斷裂圖像（b-d）

一般地，臨界晶粒尺寸dc約為20-30nm，對于晶粒尺寸（50-100nm）較大的變形主要由位錯決定；晶粒尺寸小于30nm時主要由GB變形過程主導，降低晶粒尺寸則會導致強度和流動應力減小，即出現“反Hall-Petch效應”。然而，在fcc和bcc金屬中用于模擬GBs的多體和對電位之間的廣泛比較揭示了這些不同力描述預測的行為中幾乎沒有質量差異，這表明多體效應可能不會主導GB行為。

B_ejaud, J. Durinck等人利用分子動力學模擬研究了變形孿晶與納米結構Cu / Ag界面的相互作用，分析了界面結構對孿晶形核、擴展和增厚的影響，并闡述了失配界面位錯網格的作用^[4]。圖4?顯示了肖克利部分位錯網格（由黑線突出顯示），三角形圖案（白色區段）以及界面處的堆垛層錯分布。其中，原子根據中心對稱參數著色，藍色的原子處于完美的FCC環境中，而紅色的原子處于層錯或孿晶斷層中。

圖4 （a）沿界面的Cu和Ag原子的頂視圖：（a.i）COC界面、（a.ii）TO界面，（b）沿X = <011>方向的側視圖顯示：（b.i）在COC界面和情況下，相干區域與固有堆垛層錯（ISF）區域交替、（b.ii）TO界面，在Cu層和Ag層中連續存在雙缺陷區域

圖5為應力應變曲線和孿晶的原子比例與應變的函數關系。通過分析，他們發現界面可以直接或間接地通過Lomer位錯誘導孿晶位錯的成核，以及異相界面結構如何影響機械孿晶過程的不同步驟，從而影響納米結構Cu / Ag中形成的孿晶的大小。通過這種原子尺度方法，為納米級復合材料中的機械孿晶過程提供了一些有用的理論依據。

圖5 （a）應力 - 應變曲線，（b）孿晶的原子比例與應變的函數關系

設計多層材料來調整機械性能是一個熱門話題，同時控制變形機制，因為孿晶允許結合納米層和納米孿晶材料的機械性質。在這方面，這項研究提供了理解雙界面相互作用機制的關鍵，并支持了異相界面促進孿晶的觀點。

對于包含低對稱性hcp 結構金屬的超細尺度層狀復合材料，由于大量的異質界面能夠有效地吸收核輻照引起的空位、間隙原子等缺陷，且hcp金屬自身具有密度低、比強度和比剛度高、導電導熱性好等特點，近年來由Ti、Zr、Mg與其它金屬復合而成的六方系多層材料開始受到人們的關注。但是與晶體結構對稱性高的fcc 和bcc 金屬相比，hcp金屬的室溫塑性變形能力差，制約了相關復合材料的使用^[1]。

除了原子尺度的空間和時間分辨率，分子動力學模擬一方面可以闡述完全表征的理想化納米晶體模型的界面結構、驅動力和原子機制等行為；另一方面可以在非常高的晶界和位錯密度下觀測大塑性變形行為。例如，位錯形核機制、晶界湮滅、納米晶Al中的機械孿晶、減小晶粒尺寸時從位錯到基于晶界的變形機制、剪切帶的觀測及其與凹陷斷裂面的關系。

此外，在實際應用和研究過程中，對于不同的問題描述和選用的理論模型，動力學已發展了很多的理論分支，比如，賓夕法尼亞大學的Jian Han, Spencer L. Thomas等人依賴于晶界動力學的斷開描述總結了多晶材料的晶界動力學概念^[5]，Zheng Ma等人對FeCO₃進行了析出動力學研究^[6]，還有表面/界面動力學^[7]等。

參考文獻

1.Wolf, V. Yamakov, S.R. Phillpot, A. Mukherjee, H. Gleiter，Deformation of nanocrystalline materials by molecular-dynamicssimulation: relationship to experiments??[J]，Acta Materialia 53 (2005) 1–40

2.Cai, Wensheng, Christophe Chipot，Frontiers in High-Performance, Large-Scale Molecular Dynamics.35 Years of Molecular-Dynamics Simulations of Biological Systems?[J]，Acta?Chimica?Sinica，DOI: 10.6023/A12110930

3.William Gon?alves, Julien Morthomas, Patrice Chantrenne, Michel Perez,Genevi_eve Foray ,?Christophe L. Martin，Elasticity and strength of silica aerogels: A molecular dynamics study?on large volumes?[J]，Acta Materialia 145 (2018) 165-174

4.B_ejaud, J. Durinck, S. Brochard，Twin-interface interactions in nanostructured Cu/Ag: Moleculardynamics study?[J]，Acta Materialia 144 (2018) 314-324

5.Jian Han, Spencer L. Thomas, David J. Srolovitz，Grain-boundary kinetics: A unified approach[J]，Progress in Materials Science 98 (2018) 386–476

6.Zheng Ma, Yang Yang, Bruce Brown, Srdjan Nesic, Marc Singer，Investigation of precipitation kinetics of FeCO₃by EQCM?[J]，Corrosion Science 141 (2018) 195–202

7.Ndzondelelo Bingwa, Semakaleng Bewana, Marco Haumann, Reinout Meijboom，Revisiting kinetics of morin oxidation: Surface kinetics analysis[J]，Applied Surface Science 426 (2017) 497–503

本文由材料人科技顧問張博士供稿。

投稿以及內容合作可加編輯微信：cailiaokefu。

歡迎大家到材料人宣傳科技成果并對文獻進行深入解讀，投稿郵箱tougao@cailiaoren.com。