# lmp-plugins

**Repository Path**: z9527567/lmp-plugins

## Basic Information

- **Project Name**: lmp-plugins
- **Description**: No description available
- **Primary Language**: Unknown
- **License**: Not specified
- **Default Branch**: master
- **Homepage**: None
- **GVP Project**: No

## Statistics

- **Stars**: 0
- **Forks**: 0
- **Created**: 2026-04-17
- **Last Updated**: 2026-04-27

## Categories & Tags

**Categories**: Uncategorized

**Tags**: None

## README

# GPU Plugin for LAMMPS - 势函数文档

本项目包含两个自定义势函数，用于模拟 patchy 粒子体系的各向异性相互作用。

---

## 1. PairLJAngleSEF（各向异性 LJ + 偶极-偶极相互作用）

### 1.1 概述

`PairLJAngleSEF` 势函数结合了标准的 Lennard-Jones 相互作用和各向异性的偶极-偶极相互作用。粒子的取向由偶极矩向量 `μ` 表示。该势函数的特点是同时对能量和力在截断距离处进行移位（Shifted Force, SEF）。

### 1.2 总势能表达式

总相互作用势能为：

$$
U_{ij} = U_{\text{LJ}}(r_{ij}) + U_{\text{dd}}(r_{ij}, \theta_i, \theta_j)
$$

其中：
- $r_{ij}$ 是两粒子间距离
- $\theta_i$ 是粒子 $i$ 的偶极矩与两粒子连线的夹角
- $\theta_j$ 是粒子 $j$ 的偶极矩与两粒子连线的夹角

### 1.3 Lennard-Jones 势

$$
U_{\text{LJ}}(r) = 4\epsilon \left[ \left(\frac{\sigma}{r}\right)^{12} - \left(\frac{\sigma}{r}\right)^{6} \right]
$$

力：
$$
\mathbf{F}_{\text{LJ}} = -\frac{1}{r} \frac{\partial U_{\text{LJ}}}{\partial r} \mathbf{r}_{ij} = r^{-2} \left[ 48\epsilon \left(\frac{\sigma}{r}\right)^{12} - 24\epsilon \left(\frac{\sigma}{r}\right)^{6} \right] \mathbf{r}_{ij}
$$

### 1.4 偶极-偶极相互作用

标准偶极-偶极相互作用为：

$$
U_{\text{dd}}^{\text{std}} = \frac{1}{r^{3}} \left[ \boldsymbol{\mu}_i \cdot \boldsymbol{\mu}_j - 3(\boldsymbol{\mu}_i \cdot \hat{\mathbf{r}}_{ij})(\boldsymbol{\mu}_j \cdot \hat{\mathbf{r}}_{ij}) \right]
$$

其中 $\hat{\mathbf{r}}_{ij} = \mathbf{r}_{ij} / r$。

### 1.5 角相关调制函数 $h(\theta)$

本势函数引入角相关调制函数 $h(\theta)$，使得相互作用强度依赖于相对取向：

$$
h(\theta) = (1 - x^2)^3, \quad x = \frac{\theta - \theta_0}{\theta_c}
$$

其中：
- $\theta_0$ 是特征取向角（reference angle）
- $\theta_c$ 是角度截断（angular cutoff），定义 $x \in [-1, 1]$ 的有效范围

当 $|x| \geq 1$ 时，$h(\theta) = 0$。

### 1.6 修正的偶极-偶极势

引入参数 $\Delta$（MyDelta）调制各向异性强度，修正后的势能为：

$$
U_{\text{dd}} = -\Delta \cdot 24\epsilon \left(\frac{\sigma}{r}\right)^{6} \left[ h(\theta_i) + h(\theta_j) - \frac{64}{35\pi}\theta_c \right]
$$

其中常数项 $\frac{64}{35\pi}\theta_c$ 用于在截断处进行移位。

### 1.7 移位力修正（Shifted Force）

在截断距离 $r_c$ 处，能量和力都被移位以确保连续性：

**能量移位：**
$$
U_{\text{shift}}(r) = U(r) - U(r_c) + F(r_c)(r_c - r)
$$

**力移位：**
$$
\mathbf{F}_{\text{shift}}(r) = \mathbf{F}(r) - \mathbf{F}(r_c)
$$

移位后的力在 $r_c$ 处为零。

### 1.8 参数

| 参数 | 符号 | 说明 |
|------|------|------|
| epsilon | $\epsilon$ | LJ 能量参数 |
| sigma | $\sigma$ | LJ 尺寸参数 |
| MyDelta | $\Delta$ | 各向异性调制参数 |
| thetac | $\theta_c$ | 角度截断 |
| theta0 | $\theta_0$ | 特征取向角 |

### 1.9 扭矩

除了平移力，该势函数还计算粒子上的扭矩：

$$
\boldsymbol{\tau}_i = \boldsymbol{\mu}_i \times \mathbf{F}_i
$$

---

## 2. PairBodyNparticlePatchy（Patchy 粒子势）

### 2.1 概述

`PairBodyNparticlePatchy` 势函数用于模拟表面具有 patch（贴点）的球形粒子。Patch 位于粒子表面，通过 body 风格 nparticle 定义。粒子间的相互作用取决于 patch 之间的对齐程度。

### 2.2 中心-中心相互作用

基础的 LJ 势（12-6）：

$$
U_{\text{LJ}}(r) = 4\epsilon \left[ \left(\frac{\sigma}{r}\right)^{12} - \left(\frac{\sigma}{r}\right)^{6} \right]
$$

### 2.3 Patch 选择

对于给定的两粒子 $a$ 和 $b$，选择与连线 $\mathbf{r}_{ab}$ 夹角最小的 patch：

- 在粒子 $a$ 上：选择使 $\cos\psi_{\alpha}$ 最大的 patch $\alpha$
- 在粒子 $b$ 上：选择使 $\cos\psi_{\beta}$ 最大的 patch $\beta$

其中：
$$
\cos\psi_{\alpha} = -\frac{\mathbf{r}_{ab} \cdot \mathbf{d}_\alpha}{r_{ab} \cdot l_p}
$$
$$
\cos\psi_{\beta} = \frac{\mathbf{r}_{ab} \cdot \mathbf{d}_\beta}{r_{ab} \cdot l_p}
$$

这里 $\mathbf{d}_\alpha$ 是 patch $\alpha$ 相对于粒子中心的位移向量，$l_p$ 是 patch 长度。

### 2.4 方向性调制函数 $\Phi$

引入高斯型调制函数：

$$
\Phi = \lambda \exp(-\mu_a^2) \exp(-\mu_b^2)
$$

其中：
$$
\mu_a = \frac{1 - \cos\psi_\alpha}{\delta}, \quad \mu_b = \frac{1 - \cos\psi_\beta}{\delta}
$$

- $\lambda$：调制强度参数
- $\delta$：角度展宽参数

### 2.5 排斥与吸引分支

该势函数在 WCA 截断半径 $r_{\text{WCA}} = 2^{1/6}\sigma$ 处分裂为两支：

**排斥支**（$r < r_{\text{WCA}}$）：
$$
U_{\text{rep}} = U_{\text{LJ}}(r) + \epsilon + U_{ab}^m \cdot \Phi
$$

其中 $U_{ab}^m = -\epsilon - U_{abc}$，用于移位。

**吸引支**（$r \geq r_{\text{WCA}}$）：
$$
U_{\text{att}} = (U_{\text{LJ}}(r) - U_{abc}) \cdot \Phi
$$

其中 $U_{abc} = 4\epsilon\left[\left(\frac{\sigma}{r_c}\right)^{12} - \left(\frac{\sigma}{r_c}\right)^{6}\right]$ 是在截断距离 $r_c$ 处的 LJ 势能值。

### 2.6 力和扭矩

力：
$$
\mathbf{F}_{ab} = F_{\text{LJ}} \cdot \mathbf{r}_{ab} + (U_{\text{LJ}} - U_{abc}) \cdot \nabla \Phi
$$

扭矩：
$$
\boldsymbol{\tau}_a = \mathbf{d}_\alpha \times \left[ -U_{ab}^m \cdot \frac{\partial \Phi}{\partial \mathbf{d}_\alpha} \right]
$$
$$
\boldsymbol{\tau}_b = \mathbf{d}_\beta \times \left[ -U_{ab}^m \cdot \frac{\partial \Phi}{\partial \mathbf{d}_\beta} \right]
$$

### 2.7 AB 因子（异种粒子相互作用）

对于不同类型的粒子，可通过 ABfactor 参数调节 patch 选择：

- ABfactor = 1.0：完全按角度选择
- ABfactor < 1.0：降低方向性，ABfactor → 0 时退化为纯各向同性 LJ

### 2.8 参数

| 参数 | 符号 | 说明 |
|------|------|------|
| epsilon | $\epsilon$ | LJ 能量参数 |
| sigma | $\sigma$ | LJ 尺寸参数 |
| delta | $\delta$ | 角度展宽 |
| lambda | $\lambda$ | 调制强度 |
| ABfactor | $f_{AB}$ | 异种粒子方向性因子 |

---

## 3. 使用方法

### PairLJAngleSEF

```bash
pair_style lj/angle/SEF cutoff
pair_coeff type1 type2 epsilon sigma MyDelta thetac theta0
```

### PairBodyNparticlePatchy

```bash
pair_style body/nparticle/patchy cutoff
pair_coeff type1 type2 epsilon sigma delta lambda ABfactor
```

---

## 4. 参考

本势函数基于以下工作：

- Yuan-Chao Hu, The University of Tokyo
- 修改自 `pair_lj_cut_dipole_cut.cpp`
- 应用于 patchy 粒子体系的建模

LAMMPS 官方文档：https://lammps.sandia.gov/doc/Pair_style.html