NN

NN 模型

操作演示

Neural Network(NN) 模型介绍

Neural Network(NN) 中实现了如下 8 种具有平移、旋转和置换不变性的特征类型

        1. 2-body(2b)
        2. 3-body(3b)
        3. 2-body Gaussian(2b gauss)
        4. 3-body Cosine(3b cos)
        5. Moment Tensor Potential(MTP)
        6. Spectral Neighbor Analysis Potential(SNAP)
        7. DP-Chebyshev(dp1)
        8. DP-Gaussian(dp2)

特征（或描述符）是描述原子局部环境的量。它们需要保持平移、旋转和置换对称性。特征通常用作各种回归器（线性模型、神经网络等）的输入，这些回归器输出原子能量和力。由于特征是空间坐标的可微函数，因此可以计算力：

$$\mathbf{F}_i = - \frac{d E_{\text{tot}}}{d \mathbf{R}_i} = - \sum_{j,\alpha} \frac{\partial E_j}{\partial G_{j,\alpha}} \frac{\partial G_{j,\alpha}}{\partial \mathbf{R}_i}$$

其中，$j$ 是在截断半径内的近邻原子的索引，$\alpha$ 是特征的索引。

2-b and 3-b features with piecewise cosine functions (feature 1 & 2)

给定一个中心原子，利用分段余弦函数来描述其局部环境。通过下面的图表，可以大致了解它们的原理。

我们首先定义分段余弦函数，分别用于两体和三体特征。给定内部和外部截断 $R_{\text{inner}}$ 和 $R_{\text{outer}}$，基函数的阶数 $M$，分段函数的宽度 $h$，以及中心原子 $i$ 和近邻原子 $j$ 之间的原子间距 $R_{ij}$，我们定义基函数为

$$\phi_\alpha (R_{ij}) = \begin{cases} \frac{1}{2} \cos\left( \frac{R_{ij} - R_\alpha}{h} \pi \right) + \frac{1}{2} &, |R_{ij} - R_\alpha| < h \\ 0 &, \text{otherwise} \end{cases}$$

其中

$$R_\alpha = R_{\text{inner}} + (\alpha - 1) h,\ \alpha = 1, 2, ..., M$$

中心原子 $i$ 的 两体特征 表达式为

$$G_{\alpha,i} = \sum_{m} \phi_{\alpha}(R_{ij})$$

而 三体特征 表达式为

$$G_{\alpha\beta\gamma,i} = \sum_{j,k} \phi_{\alpha}(R_{ij}) \phi_{\beta}(R_{ik}) \phi_{\gamma}(R_{jk})$$

其中 $\sum_{m}$ 和 $\sum_{m,n}$ 分别表示在中心原子 $i$ 的截断半径 $R_{\text{outer}}$ 内的所有原子的求和。

这两个特征通常是成对使用的。

参考文献：

Huang, Y., Kang, J., Goddard, W. A. & Wang, L.-W. Density functional theory based neural network force fields from energy decompositions. Phys. Rev. B 99, 064103 (2019)

2-b and 3-b Gaussian feature (feature 3 & 4)

这两个特征是 Behler-Parrinello 神经网络中首次使用的特征。给定截断半径 $R_c$ 及中心原子 $i$ 和近邻原子 $j$ 之间的原子间距 $R_{ij}$，定义截断函数 $f_c$

$$f_c(R_{ij}) = \begin{cases} \frac{1}{2} \cos\left( \frac{\pi R_{ij}}{R_c} \right) + \frac{1}{2} &, R_{ij} < R_c \\ 0 &, \text{otherwise} \end{cases}$$

中心原子 $i$ 的 两体高斯 特征定义为

$$G_i = \sum_{j \neq i} e^{(-\eta(R_{ij} - R_s)^2)} f_c(R_{ij})$$

其中 $\eta$ 和 $R_s$ 是用户定义的参数。

中心原子 $i$ 的 三体高斯 特征定义为

$$G_i = 2^{1-\zeta} \sum_{j,k \neq i} (1 + \lambda \cos \theta_{ijk})^\zeta e^{-\eta(R_{ij}^2 + R_{ik}^2 + R_{jk}^2)} f_c(R_{ij}) f_c(R_{ik}) f_c(R_{jk})$$

其中

$$\cos \theta_{ijk} = \frac{\mathbf{R}_{ij} \cdot \mathbf{R}_{ik}}{|\mathbf{R}_{ij}| |\mathbf{R}_{ik}|}$$

$\eta$、$\zeta$ 和 $\lambda = \pm 1$ 是用户定义的参数。

这两个特征通常是成对使用的。

参考文献：

J. Behler and M. Parrinello, Generalized Neural-Network Representation of High Dimensional Potential-Energy Surfaces. Phys. Rev. Lett. 98, 146401 (2007)

Moment Tensor Potential (feature 5)

在 MTP 中，中心原子 $i$ 的局部环境由

$$\mathbf{n}_i = (z_i, z_j, \mathbf{r}_{ij})$$

定义，其中 $z_i$ 是中心原子的原子类型，$z_j$ 是近邻原子的原子类型，$\mathbf{r}_{ij}$ 是近邻原子的相对坐标。接下来，每个原子的能量贡献被展开为

$$E_i(\mathbf{n}_i) = \sum_\alpha c_\alpha B_\alpha(\mathbf{n}_i)$$

其中 $B_\alpha$ 是用户选择的基函数，$c_\alpha$ 是待拟合的参数。

为了构造基函数，我们引入矩张量 $M_{\mu\nu}$ 来定义基函数

$$M_{\mu\nu}(\mathbf{n}_i) = \sum_j f_\mu(|\mathbf{r}_{ij}|, z_i, z_j) \bigotimes_\nu \mathbf{r}_{ij}$$

这些矩张量包含径向和角度部分。径向部分可以展开为

$$f_\mu(|\mathbf{r}_{ij}|, z_i, z_j) = \sum_\beta c^{(\beta)}_{\mu,z_i,z_j} Q^{(\beta)}(|\mathbf{r}_{ij}|)$$

其中 $Q^{(\beta)}(|\mathbf{r}_{ij}|)$ 是径向基函数。具体地，

$$Q^{(\beta)}(|\mathbf{r}_{ij}|) = \begin{cases} \phi^{(\beta)}(|\mathbf{r}_{ij}|) (R_{\text{cut}} - |\mathbf{r}_{ij}|)^2 &, |\mathbf{r}_{ij}| < R_{\text{cut}} \\ 0 &, \text{otherwise} \end{cases}$$

其中 $\phi^{(\beta)}$ 是定义在区间 $[R_{\text{min}}, R_{\text{cut}}]$ 上的多项式（例如切比雪夫多项式）。

角度部分由 $\bigotimes_\nu \mathbf{r}_{ij}$ 给出，它表示对 $\mathbf{r}_{ij}$ 进行 $\nu$ 次张量积，包含了近邻 $\mathbf{n}_i$ 的角度信息。$\nu$ 决定了矩张量的秩。当 $\nu=0$ 时，得到一个常数标量；当 $\nu=1$ 时，得到一个向量（秩-1 张量）；当 $\nu=2$ 时，得到一个矩阵（秩-2 张量）；以此类推。

最后，我们定义矩张量的级数为

$$\text{lev}(M_{\mu\nu}) = 2 + 4\mu + \nu$$

这是一个经验公式。

参考文献：

I.S. Novikov, etal, The MLIP package: moment tensor potential with MPI and active learning. Mach. Learn.: Sci. Technol, 2, 025002 (2021)

Spectral Neighbor Analysis Potential (feature 6)

在 SNAP 中，不使用高斯基函数。因此不计算两个原子局域环境图之间的距离和核函数。它首先定义一个原子局域环境，然后使用球谐函数（或 4D 球，带有旋转矩阵）来展开原子局域环境。然后使用双谱，使其具有旋转不变性。从某种意义上说，它类似于 MTP，但它使用一种特殊的方法来缩并方向指数，使其具有旋转不变性。它通常与线性回归一起使用。

首先，它定义位于 $\mathbf{r}$ 处的中心原子 $i$ 的邻近原子周围的原子局域环境为三维空间中的 $\delta$ 函数之和：

$$\rho(\mathbf{r}) = \delta(\mathbf{r}) + \sum_{\mathbf{r}_{ki} < R_C} f_C(\mathbf{r}_{ki}) \omega_k \delta(\mathbf{r} - \mathbf{r}_{ki})$$

其中 $\mathbf{r}_{ki}$ 是原子 $i$ 的第 $k$ 个近邻的位置，$\omega_k$ 是第 $k$ 个近邻的权重，径向函数 $f_C(\mathbf{r}_{ki})$ 确保每个近邻的贡献在截断半径 $R_C$ 附近平滑地变为零：

$$f_C(\mathbf{r}) = 0.5 \left[ \cos\left( \frac{\pi r}{R_C} \right) + 1 \right]$$

这个原子局域环境函数的角部分可以用球谐函数展开，球谐函数定义在 $l = 0, 1, 2, ...$ 和 $m = -l, -l+1, ..., l-1, l$ 的基础上。径向分布通常由一组径向基函数表示。然而，在这里，径向信息 $\mathbf{r}$ 被映射到 4D 超球面函数 $U^j_{mm'}(\theta_0,\theta,\phi)$ 中，其中所有点（邻近原子）落入 3D 球面（在 4D 空间中），定向（角度）信息由三个角度给出：

$$\mathbf{r} \equiv \begin{pmatrix} x \\ y \\ z \end{pmatrix} \rightarrow \begin{matrix} \phi = \arctan(y/x) \\ \theta = \arccos(z/r) \\ \theta_0 = \frac{3}{4} \pi r / r_c \end{matrix}$$

因此，上述原子局域环境函数可以用这些 4D 超球面函数 $U^j_{mm'}(\theta_0,\theta,\phi)$ 展开，展开系数为 $u^j_{mm'}$：

$$\rho(\mathbf{r}) = \sum_{j=0,\frac{1}{2},1,...}^\infty \sum_{m=-j,-j+1}^{j} \sum_{m'=-j,-j+1}^{j} u^j_{mm'} U^j_{mm'}(\theta_0,\theta,\phi)$$

使用上述原子局域环境函数，可以计算 $u^j_{mm'}$：

$$u^j_{mm'} = U^j_{mm'}(0,0,0) + \sum_{\mathbf{r}_{ki} < R_C} f_C(\mathbf{r}_{ki}) \omega_k U^j_{mm'}(\theta_0(k),\theta(k),\phi(k))$$

其中，$k$ 是邻近原子的索引，$\theta_0(k),\theta(k),\phi(k)$ 是原子 $i$ 的第 $k$ 个近邻的位置矢量的三个角度。注意，$u^j_{mm'}$ 是由于其指数 $m, m'$ 而具有方向依赖性。它们可以与下面的缩并公式（使用三个 $u^j_{mm'}$）缩并：

$$F(j_1,j_2,j) = \sum_{m_1,m_1'=-j_1}^{j} \sum_{m_2,m_2'=-j_2}^{j} \sum_{m,m'=-j}^{j} (u^{j}_{mm'})^* u^{j_1}_{m_1 m_1'} u^{j_2}_{m_2 m_2'} \times C_{j_1 m_1 j_2 m_2}^{j m} C_{j_1 m_1' j_2 m_2'}^{j m'}$$

这里，$C_{j_1 m_1 j_2 m_2}^{j m} C_{j_1 m_1' j_2 m_2'}^{j m}$ 是 Clebsch-Gordan 系数，最终的标量特征是 $F(j_1,j_2,j)$。通过设置不同的 $j_1,j_2,j$，可以产生不同的特征。注意，在这些特征中，没有径向函数索引，而是有三个角动量索引。这是因为我们已经将径向距离信息转换为 3D 球面中的第三个角度信息。

DP-Chebyshev (feature 7)

这个特征类似于 DP 的嵌入网络。它使用切比雪夫多项式作为基础。

首先，我们将 $S(\mathbf{r}_{ij})$ 定义为加权的距离的倒数函数：

$$S(\mathbf{r}) = \frac{f_C(\mathbf{r})}{r}$$ $$f_C(\mathbf{r}) = \begin{cases} 1 &, r < R_{C_2} \\ \frac{1}{2} \cos\left( \pi \frac{r - R_{C_2}}{R_c - R_{C_2}} \right) + \frac{1}{2} &, R_{C_2} \leq r < R_C \\ 0 &, r > R_C \end{cases}$$

这里，$R_{C_2}$ 是一个平滑的截断参数，它允许在由 $R_C$ 定义的局部区域的边界上平滑地将 $\mathbf{r}_i$ 的分量减小到零。这个平滑函数比之前使用的 $R_{C_2}$ 更复杂一些。$S(\mathbf{r}_{ji})$ 减小了远离中心原子 $i$ 的原子的权重。然后，我们定义径向函数 $g_M(s)$ 为深度势能切比雪夫特征中的切比雪夫多项式 $C_M$：

$$g_M(s) = C_M(2R_{\min} S - 1)$$

这里，$R_{\min}$ 是最小 $r$ 的输入。

为了构造这样的特征，我们首先计算四个分量的向量：

$$T_M(k) = \sum_{\mathbf{r}_{ji} < R_C} \hat{X}_{ji}(k) S(\mathbf{r}_{ji}) g_M(S(\mathbf{r}_{ji}))$$

这里，$k = 0,1,2,3$（四分量向量）。它们是由通常的 $x,y,z$ 分量构成的，再加上 $S$ 分量：

$$\{ x_{ji}, y_{ji}, z_{ji}\} \rightarrow \{ S(\mathbf{r}_{ji}), \hat{x}_{ji}, \hat{y}_{ji}, \hat{z}_{ji} \}$$

其中 $\hat{x}_{ji} = \frac{x_{ji}}{r_{ji}}, \hat{y}_{ji} = \frac{y_{ji}}{r_{ji}}, \hat{z}_{ji} = \frac{z_{ji}}{r_{ji}}$ 是 $\mathbf{r}_{ji}$ 的单位向量。

从这些 4D 向量中，我们可以缩并分量索引以得到标量特征：

$$F(M_1,M_2) = \sum_{k=0}^3 T_{M_1}(k) T_{M_2}(k)$$

这里，$M_1$ 也编码了除切比雪夫外的原子类型的数量。因此，如果最大切比雪夫阶数是 $M$，特征的数量是 $M \cdot n_{\text{type}} \cdot (M \cdot n_{\text{type}} + 1) / 2$。我们可以通过设置不同的 $M$ 来产生不同的特征。

DP-Gaussian (feature 8)

这个特征类似于 DP-Chebyshev，但我们使用高斯函数代替切比雪夫多项式，并且位置和宽度参数由用户指定。

类似于 DP-Chebyshev，4D 向量构造如下：

$$T_M(k) = \sum_{\mathbf{r}_{ji} < R_C} \hat{X}_{ji}(k) g_M(\mathbf{r}_{ji})$$ $$\hat{X}(0) = S(\mathbf{r}'), \quad \hat{X}(1) = \frac{x}{r}, \quad \hat{X}(2) = \frac{y}{r}, \quad \hat{X}(3) = \frac{z}{r}$$ $$g_M(\mathbf{r}) = f_C(\mathbf{r}) \cdot \exp\left( -\frac{(r - r_M)}{\omega_M} \right)$$ $$f_C(\mathbf{r}) = \frac{1}{2} \cos\left( \frac{\pi r}{R_C} \right) + \frac{1}{2}$$

缩并过程如下：

$$F(M_1,M_2) = \sum_{k=0}^3 T_{M_1}(k) T_{M_2}(k)$$

这里，$M_1$ 也编码了除切比雪夫外的原子类型的数量。因此，如果最大切比雪夫阶数是 $M$，特征的数量是 $M \cdot n_{\text{type}} \cdot (M \cdot n_{\text{type}} + 1) / 2$。我们可以通过设置不同的 $M$ 来产生不同的特征。