2.  力の計算コードを追記する。

$$\begin{array}{rcl}& & \underset{i}{\mathrm{\Lambda }}\underset{j>i}{\mathrm{\Lambda }}\{V+=\epsilon f_2(r_{ij});\\ & & {\mathit{F}}_i+=-\epsilon \frac{\partial f_2(r_{ij})}{\partial r_{ij}}\frac{{\mathit{r}}_i-{\mathit{r}}_j}{r_{ij}};\\ & & {\mathit{F}}_j+=-\epsilon \frac{\partial f_2(r_{ij})}{\partial r_{ij}}\frac{{\mathit{r}}_j-{\mathit{r}}_i}{r_{ij}}\}\\ & & +\underset{i}{\mathrm{\Lambda }}\underset{j\ne i}{\mathrm{\Lambda }}\underset{\stackrel{k>j}{k\ne i}}{\mathrm{\Lambda }}\{V+=\epsilon h(r_{ij},r_{ik},{\theta }_{jik});\\ & & {\mathit{F}}_i+=-\epsilon \frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial {\mathit{r}}_i};\\ & & {\mathit{F}}_j+=-\epsilon \frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial {\mathit{r}}_j};\\ & & {\mathit{F}}_k+=-\epsilon \frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial {\mathit{r}}_k}\}\end{array}$$

３体項の微分を実行すると

$$\begin{array}{rcl}& =& \underset{i}{\mathrm{\Lambda }}\underset{j>i}{\mathrm{\Lambda }}\{V+=\epsilon f_2(r_{ij});\\ & & {\mathit{F}}_i+=-\epsilon \frac{\partial f_2(r_{ij})}{\partial r_{ij}}\frac{{\mathit{r}}_i-{\mathit{r}}_j}{r_{ij}};\\ & & {\mathit{F}}_j+=-\epsilon \frac{\partial f_2(r_{ij})}{\partial r_{ij}}\frac{{\mathit{r}}_j-{\mathit{r}}_i}{r_{ij}}\}\\ & & +\underset{i}{\mathrm{\Lambda }}\underset{j\ne i}{\mathrm{\Lambda }}\underset{\stackrel{k>j}{k\ne i}}{\mathrm{\Lambda }}\{V+=\epsilon h(r_{ij},r_{ik},{\theta }_{jik});\\ & & {\mathit{F}}_i+=-\epsilon A_{jik}\frac{{\mathit{r}}_i-{\mathit{r}}_j}{r_{ij}}-\epsilon B_{jik}\frac{{\mathit{r}}_i-{\mathit{r}}_k}{r_{ik}};\\ & & {\mathit{F}}_j+=-\epsilon A_{jik}\frac{{\mathit{r}}_j-{\mathit{r}}_i}{r_{ij}}-\epsilon D_{jik}\frac{{\mathit{r}}_j-{\mathit{r}}_k}{r_{jk}};\\ & & {\mathit{F}}_k+=-\epsilon B_{jik}\frac{{\mathit{r}}_k-{\mathit{r}}_i}{r_{ij}}-\epsilon D_{jik}\frac{{\mathit{r}}_k-{\mathit{r}}_j}{r_{jk}}\}\end{array}$$

ここで

$$\begin{array}{rcl}{\displaystyle A_{jik}}& =& \frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial r_{ij}}+\frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial \cos {\theta }_{jik}}\frac{r_{ij}–r_{ik}\cos {\theta }_{jik}}{r_{ij}{r}_{ik}}\\ B_{jik}& =& \frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial r_{ik}}+\frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial \cos {\theta }_{jik}}\frac{r_{ik}–r_{ij}\cos {\theta }_{jik}}{r_{ij}{r}_{ik}}\\ D_{jik}& =& -\frac{\partial h(r_{ij},r_{ik},{\theta }_{jik})}{\partial \cos {\theta }_{jik}}\frac{r_{jk}}{r_{ij}{r}_{ik}}\end{array}$$

である。

ポイントは、注目している項に関与しているすべての原子について力を計算することである。２体項では$f_2(r_{ij})$によって生じる力を原子$i$と$j$について計算する。３体項では$h(r_{ij},r_{ik},{\theta }_{jik})$から生じる力を原子$i,j,k$について計算する。$C_{ijk}$がなくなっているのが気になるかもしれないが、$C_{ijk}$は$i$と$j$が入れ替わった$A_{jik}$に他ならず、ここで計算しなくても$h(r_{ji},r_{jk},{\theta }_{ijk})$の項で計算されるから問題ない。

１つの原子に注目し、その原子にかかる力をすべて列挙しようとすると、３体項の場合$h(r_{ij},r_{ik},{\theta }_{jik}),h(r_{ji},r_{jk},{\theta }_{ijk}),h(r_{ki},r_{kj},{\theta }_{kji})$といくつもの項を考えなければならず、見落としや重複のリスクが常につきまとう。ターム指向であれば、$h(r_{ij},r_{ik},{\theta }_{jik})$一つに注目すればよいので、ほぼ機械的に間違いなく記述できる。