分子动力学模拟基础知识

概要

从分子动力学（Molecular Dynamics, MD）模拟中学习原子间势（Interatomic Potential）是一个复杂的过程，但可以通过以下几个步骤来实现：

1. 选择合适的力场模型：在计算原子间势之前，需要选择一个合适的力场模型，该模型用于描述原子间相互作用。常见的力场模型包括：lj potential（Lennard-Jones 势）、eam potential（Embedded Atom Method 势）、 Morse potential（Morse 势）等。

2. 构建系统：构建一个由原子组成的系统，该系统应该包含足够多的原子以确保系统的统计性质。

3. 设置模拟参数：包括选择合适的时间步长、模拟温度和压力等。

4. 运行分子动力学模拟：使用选择的力场模型和模拟参数运行分子动力学模拟。

5. 计算原子间势：从模拟结果中计算出原子间势，这可以通过计算相邻原子对之间的平均势能来实现。可以使用以下公式计算原子间势：

   U® = <E®> - <E(r=infinity)>

   其中U®表示原子间势，E®表示相邻原子对之间的平均势能，r表示原子间距离，<E(r=infinity)>表示在无穷远处原子对之间的平均势能。

6. 验证原子间势：验证计算出的原子间势是否与实际情况一致，可以通过与实验数据或其他计算方法（如密度泛函法）进行比较来实现。

7. 优化原子间势：如果计算出的原子间势与实际情况存在 discrepancy，则需要进一步优化力场模型。

通过上述步骤，您可以从分子动力学模拟中学习原子间势。需要注意的是，这是一个复杂的过程，需要对物理学和计算机科学有深入的了解。

一、计算原理

分子动力学（Molecular Dynamics, MD）模拟是一种通过计算机模拟分子系统随时间演化的方法，用于研究分子运动和相互作用。在从头计算（ab initio）中学习原子间势（interatomic potentials）的过程中，通常采用量子力学计算（如密度泛函理论，Density Functional Theory, DFT）获得原子间的相互作用信息，然后将这些信息用来构建经典力场（classical force field），并用于MD模拟。以下是详细的计算过程、计算公式和计算结果。

1. 从头计算获取原子间相互作用

首先，使用量子化学方法（例如密度泛函理论-DFT）计算体系中原子间的相互作用能量。计算的步骤包括：

a. 几何优化（Geometry Optimization）

确定分子或固体的平衡几何结构。在DFT计算中，目标是找到体系的基态能量最小的结构。

$$E(\mathbf{R})=\min \left[⟨\Psi (\mathbf{R})|\hat{H}|\Psi (\mathbf{R})⟩\right]$$

其中，R 是原子坐标， $\hat{H}$ 是哈密顿算符，$\Psi (R)$是体系的波函数。

b. 力的计算（Force Calculation）

在平衡几何结构基础上，计算每个原子受到的力：

$${\mathbf{F}}_i=-{\nabla }_{{\mathbf{R}}_i}E(\mathbf{R})$$

2. 构建经典力场

使用DFT计算得到的能量和力信息来参数化经典力场。经典力场通常采用经验势能函数的形式，如Lennard-Jones势、嵌入原子法（Embedded Atom Method, EAM）等。常见的势能函数形式包括：

a. Lennard-Jones势（Lennard-Jones Potential）

$$V(r_{ij})=4ϵ\left[{\left(\frac{\sigma }{r_{ij}}\right)}^{12}-{\left(\frac{\sigma }{r_{ij}}\right)}^6\right]$$

其中，rij​ 是原子i和j之间的距离， ϵ 和 σ 是势能参数。

b. EAM势（Embedded Atom Method）

EAM势能函数结合了双体势（pair potential）和嵌入函数（embedding function）：

$$E=\sum _i{F}_i({\rho }_i)+\frac{1}{2}\sum _{i\ne j}\varphi (r_{ij})$$

其中，Fi​(ρi​) 是嵌入函数， ρi​ 是原子i的电子密度， ϕ(rij​) 是双体势。

3. 分子动力学模拟

使用参数化的力场进行MD模拟。MD模拟的基本步骤如下：

a. 初始化（Initialization）

根据系统的初始条件（如温度、压力和初始位置）初始化原子的速度和位置。

b. 时间积分（Time Integration）

使用数值积分方法（如Verlet算法）计算随时间演化的原子位置和速度：

R i ( t + Δ t ) = R i ( t ) + V i ( t ) Δ t + 1 2 A i ( t ) ( Δ t ) 2 \mathbf{R}_i(t+\Delta t)=\mathbf{R}_i(t)+\mathbf{V}_i(t)\Delta t+\frac12\mathbf{A}_i(t)(\Delta t)^2