本文為某計算機本科生的分子動力學學習筆記，在gpt4的輔助下，非體系化地整理相關生物、化學、統計力學知識。所有生成內容經過檢查和調整，均直接代表本人觀點。有科學性錯誤的話歡迎指教。
什麼是分子動力學
定義： 分子動力學是一門結合物理，數學和化學的綜合技術。分子動力學是一套分子模擬方法，該方法主要是依靠牛頓力學來模擬分子體系的運動，以在由分子體系的不同狀態構成的系統中抽取樣本，從而計算體系的構型積分，並以構型積分的結果為基礎進一步計算體系的熱力學量和其他宏觀性質。
自己的理解：使用計算機模擬分子系統，探索其行為、微觀和宏觀性質。
兩類方法
一類基於第一性原理，使用薛定諤方程計算所有原子的波函式隨時間變化，能夠提供非常精確的模擬結果，因為它不依賴於任何經驗引數。然而，由於需要解決電子的量子力學問題，這種方法計算成本非常高，通常限制於較小的系統或者較短的時間尺度。
一類基於統計力學原理，使用經驗勢（如Lennard-Jones勢或者更復雜的多體勢）來描述原子或分子間的相互作用。在這種方法中，系統的宏觀性質是透過對大量微觀狀態的統計來預測的。
本文主要學習第二種
基本框架
1. 系統初始化
- 選擇模型系統：確定要模擬的分子或原子組成的系統。
- 定義初始位置：為系統中的所有粒子（原子或分子）指定初始位置，通常這是透過實驗資料、隨機分佈或晶體結構來確定。
- 分配初始速度：根據溫度給粒子分配初始速度。這些速度通常是從玻爾茲曼分佈中隨機選取的。
2. 勢能函式的選擇
- 選擇適當的勢能函式：這是描述粒子間相互作用的數學表示式。勢能函式可以是簡單的兩體相互作用，如Lennard-Jones勢，也可以是複雜的多體相互作用，如反應場勢或者鍵合勢。
3. 力的計算
- 計算力：在每一步模擬中，根據勢能函式計算作用在每個粒子上的力。這些力是勢能對粒子座標的負梯度。
4. 積分方程的求解
- 時間積分：使用數值積分演算法，如Verlet演算法或其變種（比如速度Verlet演算法），來更新粒子的位置和速度。這些演算法允許從當前粒子的位置和速度計算出一小時間步後的位置和速度。
玻爾茲曼分佈
玻爾茲曼分佈描述了在熱平衡狀態下，系統中的粒子在不同能量水平上的分佈情況。這一分佈是由奧地利物理學家路德維希·玻爾茲曼在19世紀末提出的，它揭示了溫度、粒子能量與粒子數之間的關係。
數學表達
玻爾茲曼分佈可以用以下公式表示：
其中：
- \(P(E)\) 是粒子具有能量\(E\)的機率密度。
- \(E\) 是粒子的能量。
- \(k_B\) 是玻爾茲曼常數，其值約為\(1.38 \times 10^{-23} \, \text{J/K}\)。
- \(T\) 是系統的絕對溫度（單位是開爾文）。
在分子動力學中
在分子動力學模擬中，需要初始化系統，初始化的系統裡的粒子的分佈就需要滿足玻爾茲曼分佈。
自由能
在熱力學和統計物理中，自由能是一個描述系統能量狀態以及系統在給定條件下進行自發過程的能力的重要物理量。自由能有兩種常見的形式：亥姆霍茲自由能（Helmholtz free energy）和吉布斯自由能（Gibbs free energy），它們分別適用於不同的熱力學條件。
亥姆霍茲自由能 (A)
亥姆霍茲自由能（\(A\)）是在恆溫、恆體積條件下系統的自由能，它定義為：
其中：
- \(A\) 是亥姆霍茲自由能，
- \(U\) 是系統的內能，
- \(T\) 是絕對溫度，
- \(S\) 是系統的熵。
亥姆霍茲自由能的減少表示在恆溫恆體積下系統自發進行的過程。
吉布斯自由能 (G)
吉布斯自由能（\(G\)）是在恆溫、恆壓條件下系統的自由能，它定義為：
或者根據\(H = U + PV\)（\(H\)是焓，\(P\)是壓力，\(V\)是體積）也可以寫作：
其中：
- \(G\) 是吉布斯自由能，
- \(H\) 是系統的焓，
- \(P\) 是壓力，
- \(V\) 是體積，
- \(T\) 和 \(S\) 分別是系統的絕對溫度和熵。
吉布斯自由能的減少表示在恆溫恆壓下系統自發進行的過程。
在分子動力學中
由克勞修斯定理知孤立系統的熵只增不減，也就是自由能只降不增。分子動力學認為，自由能降到最低時系統達到平衡態。
熵
這部分建議閱讀張老師的部落格：https://www.cnblogs.com/manuscript-of-nomad/p/17154461.html#_label2_1
系綜理論
系綜理論是統計物理學的一個核心概念，它提供了一種方法來處理和理解大量粒子系統的宏觀性質。在經典物理學中，理想情況下，如果我們知道系統中所有粒子的初始位置和速度，我們可以使用牛頓的運動定律精確地預測系統的未來狀態。然而，對於包含大量粒子的系統（如氣體、液體或固體），這種方法是不切實際的。系綜理論透過使用機率和統計方法來解決這個問題，它不是關注單個粒子的行為，而是研究系統的宏觀性質。
系綜的基本概念
系綜是指一大群虛擬的、宏觀上處於相同條件下的微觀系統的集合。每個微觀系統都是獨立的，並且具有相同的宏觀條件（如體積、能量、粒子數等）。透過研究這個集合的統計性質，我們可以得到系統的熱力學性質。
主要型別的系綜
系綜理論中有三種基本的系綜，它們分別對應於不同的物理條件：
-
微正則系綜（Microcanonical Ensemble）：適用於能量、體積和粒子數都恆定的孤立系統。在這個系綜中，所有微觀狀態具有相同的能量，因此每個狀態出現的機率相等。微正則系綜的研究重點是熵和溫度的關係。
-
正則系綜（Canonical Ensemble）：適用於體積、粒子數和溫度恆定的封閉系統。這種系統與一個大熱庫處於熱平衡。在正則系綜中，不同微觀狀態具有不同的能量，狀態出現的機率遵循玻爾茲曼分佈。正則系綜的研究重點是系統的自由能。
-
巨正則系綜（Grand Canonical Ensemble）：適用於溫度、體積和化學勢恆定的開放系統。這種系統可以與外界交換粒子。在巨正則系綜中，除了能量分佈，粒子數也是變數，系統的性質可以透過巨配分函式來描述。
在系綜理論中，相空間是一個非常重要的概念，它為描述和分析系統的微觀狀態提供了一個強大的框架。相空間是一個抽象的多維空間，其中的每一點代表系統的一個可能的微觀狀態。對於一個由粒子組成的物理系統，每個粒子的位置和動量（或速度）共同定義了系統的微觀狀態。
相空間的定義
對於一個單一粒子的簡單系統，相空間由粒子的位置和動量座標構成。例如，如果我們考慮一個自由移動的單粒子在三維空間中，其相空間將由六個維度構成：三個位置座標\((x, y, z)\)和三個與之對應的動量座標\((p_x, p_y, p_z)\)。因此，這個單粒子的相空間是一個六維空間，每個點\((x, y, z, p_x, p_y, p_z)\)代表粒子的一個可能的微觀狀態。
對於一個包含\(N\)個粒子的系統，相空間的維度將是\(6N\)，因為每個粒子都有三個位置座標和三個動量座標。因此，相空間的每一點都代表了整個系統的一個可能狀態，即所有粒子的位置和動量的一個特定組合。
能量平面
在經典力學中，一個系統的哈密頓量（Hamiltonian）\(H\)是描述系統能量的函式，它通常依賴於系統的廣義座標\(q_i\)和與之共軛的廣義動量\(p_i\)。對於一個不受外力作用、孤立的系統，哈密頓量是守恆的，即不隨時間變化。
在相空間中，系統的狀態可以由一個點表示，該點的位置由廣義座標和廣義動量的值確定。隨著時間的推移，這個點在相空間中的軌跡（稱為相軌跡）將遵循哈密頓方程：
這裡，\(dq_i/dt\)和\(dp_i/dt\)分別是廣義座標\(q_i\)和廣義動量\(p_i\)隨時間的變化率。
對於一個哈密頓量保持不變的孤立系統，其相軌跡將被限制在相空間中能量為常數的曲面上，即滿足\(H(q_i, p_i) = E\)的曲面，其中\(E\)是系統的總能量。這個曲面被稱為能量曲面或能量殼層。
配分函式
配分函式在統計物理學的系綜理論中扮演著核心角色，它是連線微觀物理狀態與宏觀熱力學性質的橋樑。配分函式的具體形式取決於所研究的系綜型別（微正則、正則或巨正則系綜），但在所有情況下，它都提供了一種計算系統宏觀性質的方法。
正則系綜的配分函式
在正則系綜中，系統與一個大熱庫處於熱平衡，能夠交換能量，但粒子數和體積保持不變。這種情況下的配分函式被稱為正則配分函式，通常用\(Z\)表示。
對於一個由\(N\)個粒子組成的系統，處於溫度\(T\)下，其正則配分函式定義為：
其中，求和遍及系統所有可能的微觀狀態\(i\)，\(E_i\)是系統在第\(i\)個微觀狀態下的能量，\(\beta = \frac{1}{k_B T}\)，\(k_B\)是玻爾茲曼常數，\(T\)是絕對溫度。
正則配分函式\(Z\)的重要性在於，它可以用來計算系統的各種熱力學量。例如，系統的自由能\(F\)、內能\(U\)、熵\(S\)、熱容\(C\)等都可以透過\(Z\)及其對溫度的導數來表達：
- 亥姆霍茲自由能：\(F = -k_B T \ln Z\)
- 內能：\(U = -\frac{\partial \ln Z}{\partial \beta}\)
- 熵：\(S = -\left(\frac{\partial F}{\partial T}\right)_V\)
- 熱容：\(C_V = \left(\frac{\partial U}{\partial T}\right)_V\)
巨正則系綜的配分函式
在巨正則系綜中，系統不僅能夠與熱庫交換能量，還能與粒子庫交換粒子。這種情況下的配分函式被稱為巨正則配分函式，通常用\(\Xi\)表示。
巨正則配分函式定義為：
其中，第一個求和遍及所有可能的粒子數\(N\)，第二個求和遍及在給定\(N\)下所有可能的微觀狀態\(i\)，\(E_{i,N}\)是系統在有\(N\)個粒子且處於第\(i\)個微觀狀態下的能量，\(\mu\)是化學勢。
巨正則配分函式\(\Xi\)同樣可以用來計算系統的各種熱力學量，如系統的壓強、化學勢、平均粒子數等。
在分子動力學模擬中
通常使用正則系綜（恆溫恆壓）。給定相空間中的一個點，我們就能計算它的能量，進而找到能量平面的最低點，也就是穩定態。能量是動能+勢能，因為我們都是恆溫模擬，所以相空間內所有狀態總動能是不變的，我們主要關心勢能。計算給定狀態的總勢能的方法被稱為力場。
力場
分子動力學模擬的力場定義了模擬中分子間和分子內部各原子之間的相互作用力。力場的選擇和引數化對模擬結果的準確性和可靠性有著決定性影響。
力場的組成
力場通常由以下幾個部分組成：
-
鍵合相互作用：這些相互作用描述了原子間的化學鍵。它們通常包括鍵長（bond lengths）、鍵角（bond angles）和二面角（dihedral angles）的勢能函式。這些勢能函式通常採用簡單的數學形式，如諧振子勢（對於鍵長和鍵角）和週期勢（對於二面角）。
-
非鍵合相互作用：這部分描述了不透過化學鍵直接相連的原子之間的相互作用，主要包括範德華力（Van der Waals forces）和庫侖力（Coulomb forces）。範德華力通常使用Lennard-Jones勢來描述，而庫侖力則描述了帶電粒子之間的靜電相互作用。
勢能函式
力場透過勢能函式（Potential Energy Function）來描述這些相互作用。總勢能\(U\)是各種相互作用勢能的總和：
其中，\(U_{\text{bond}}\)、\(U_{\text{angle}}\)、\(U_{\text{dihedral}}\)分別代表鍵合相互作用中的鍵長、鍵角和二面角的勢能，而\(U_{\text{nonbonded}}\)則包括所有非鍵合相互作用的勢能。
力場的引數化
力場的引數化是指為勢能函式中的各個項確定引數（如力常數、平衡鍵長、平衡鍵角等）的過程。這些引數通常透過實驗資料和量子化學計算來確定。準確的引數化對於模擬的成功至關重要，因為它直接影響到模擬所預測的分子結構、動力學和熱力學性質的準確性。
粗粒化力場
全原子力場計算系統內所有原子的相互作用，進而求出精確的勢能。而粗粒化力場使用一些策略將幾個臨近的原子打包成一個大粒子，然後使用人工製造的力場來求勢能。這個勢能也是人工的，它並不代表自然界真實的勢能，它為研究目的服務。
Upside-md方法
upside-md是jumper的博士畢業論文。jumper也就是Alphafold2的一作，AF2是目前最具革命性的AI模型之一，被認為可以準確預測蛋白質的三級結構，在結構生物、製藥等領域有革命性的影響。
upside-md是用於研究蛋白質動力學的模型，它的核心思想是抓住蛋白的主要矛盾（主鏈結構），將側鏈打包成有方向的珠子。整個蛋白被數學表示為主鏈座標和側鏈chi角的序列，大大減少了粒子數。在給定主鏈座標的情況下，側鏈chi角的所有狀態構成一個相空間。
使用引數化的側鏈-側鏈相互作用、側鏈-主鏈相互作用力場計算相空間內狀態的勢能、進而求改主鏈座標配置的自由能。這個引數使用機器學習方法最大似然方法求得，使用電鏡測定的真實蛋白質結構作為訓練集，訓練目標是使相空間內，真實狀態存在的機率儘量大（也就是真實狀態能量儘量小，虛假狀態能量儘量高）。
認為自由能是降低的，因此自由能對主鏈原子座標的偏導數，就是該原子受到的力。然後執行朗之萬動力學（保證系統始終滿足玻爾茲曼分佈），迭代多次後收斂到穩定狀態。
該模型的計算量非常小，可以使用單核CPU在幾天內收斂一個蛋白。可以用於蛋白的從頭摺疊（只給氨基酸序列，預測蛋白結構），蛋白的動力學分析（如多個平衡態之間變化的分析），蛋白質結合分析（蛋白-蛋白互作用、蛋白-小分子結合。