\begin{eqnarray}
E & = & \sum_i E_i = \frac{1}{2} \sum_{i \ne j} V_{ij} \textrm{ ,} \\
V_{ij} & = & f_C(r_{ij}) [ f_R(r_{ij}) + b_{ij} f_A(r_{ij}) ] \textrm{ .}
\begin{eqnarray}
E & = & \sum_i E_i = \frac{1}{2} \sum_{i \ne j} V_{ij} \textrm{ ,} \\
V_{ij} & = & f_C(r_{ij}) [ f_R(r_{ij}) + b_{ij} f_A(r_{ij}) ] \textrm{ .}
\zeta_{ij} & = & \sum_{k \ne i,j} f_C (r_{ik}) \omega_{ik} g(\theta_{ijk}) \textrm{ ,}\\
g(\theta_{ijk}) & = & 1 + c_i^2/d_i^2 - c_i^2/[d_i^2 + (h_i - \cos \theta_{ijk})^2] \textrm{ .}
\end{eqnarray}
\zeta_{ij} & = & \sum_{k \ne i,j} f_C (r_{ik}) \omega_{ik} g(\theta_{ijk}) \textrm{ ,}\\
g(\theta_{ijk}) & = & 1 + c_i^2/d_i^2 - c_i^2/[d_i^2 + (h_i - \cos \theta_{ijk})^2] \textrm{ .}
\end{eqnarray}
with respect to ${\bf r}_i$ are necessary to compute the forces for atom $i$.
According to this, for every triple $(ijk)$ the derivatives of the three
potential contributions, denoted by $V_{ijk}$, $V_{jik}$ and $V_{jki}$
with respect to ${\bf r}_i$ are necessary to compute the forces for atom $i$.
According to this, for every triple $(ijk)$ the derivatives of the three
potential contributions, denoted by $V_{ijk}$, $V_{jik}$ and $V_{jki}$
derivatives for each $(ijk)$ triple.
The $V_{jik}$ and $V_{jki}$ potential and its derivatives will be calculated
in subsequent loops anyways.
derivatives for each $(ijk)$ triple.
The $V_{jik}$ and $V_{jki}$ potential and its derivatives will be calculated
in subsequent loops anyways.
the force contributions for atom $j$ and $k$ due to the $V_{ijk}$ contribution
have to be considered by calculating the derivatives of $V_{ijk}$
with respect to ${\bf r}_j$ and ${\bf r}_k$
the force contributions for atom $j$ and $k$ due to the $V_{ijk}$ contribution
have to be considered by calculating the derivatives of $V_{ijk}$
with respect to ${\bf r}_j$ and ${\bf r}_k$
keeping in mind that all the necessary force contributions for atom $i$
are calculated and added in subsequent loops.
keeping in mind that all the necessary force contributions for atom $i$
are calculated and added in subsequent loops.
\nabla_{{\bf r}_j} f_A(r_{ij}) &=& - \nabla_{{\bf r}_i} f_A(r_{ij}) \\
\nabla_{{\bf r}_j} f_C(r_{ij}) &=& - \nabla_{{\bf r}_i} f_C(r_{ij})
\end{eqnarray}
\nabla_{{\bf r}_j} f_A(r_{ij}) &=& - \nabla_{{\bf r}_i} f_A(r_{ij}) \\
\nabla_{{\bf r}_j} f_C(r_{ij}) &=& - \nabla_{{\bf r}_i} f_C(r_{ij})
\end{eqnarray}
-The pair contributions .... easy.
-Now having a look at $b_{ij}$.
+The pair contributions are, thus, easily obtained.
+The contribution of the bond order term is given by:
\begin{eqnarray}
\nabla_{{\bf r}_j}\cos\theta_{ijk} &=&
\nabla_{{\bf r}_j}\Big(\frac{{\bf r}_{ij}{\bf }r_{ik}}{r_{ij}r_{ik}}\Big)
\begin{eqnarray}
\nabla_{{\bf r}_j}\cos\theta_{ijk} &=&
\nabla_{{\bf r}_j}\Big(\frac{{\bf r}_{ij}{\bf }r_{ik}}{r_{ij}r_{ik}}\Big)
\begin{eqnarray}
\nabla_{{\bf r}_k}\zeta_{ij} &=&
g(\theta_{ijk})\nabla_{{\bf r}_k}f_C(r_{ik}) +
f_C(r_{ik})\nabla_{{\bf r}_k}g(\theta_{ijk}) \\
\nabla_{{\bf r}_k}f_C(r_{ik}) &=& - \nabla_{{\bf r}_i}f_C(r_{ik}) \\
\nabla_{{\bf r}_k}\cos\theta_{ijk} &=&
\begin{eqnarray}
\nabla_{{\bf r}_k}\zeta_{ij} &=&
g(\theta_{ijk})\nabla_{{\bf r}_k}f_C(r_{ik}) +
f_C(r_{ik})\nabla_{{\bf r}_k}g(\theta_{ijk}) \\
\nabla_{{\bf r}_k}f_C(r_{ik}) &=& - \nabla_{{\bf r}_i}f_C(r_{ik}) \\
\nabla_{{\bf r}_k}\cos\theta_{ijk} &=&
-The implementation of the force evaluation shown in the following
-is applied to the potential designed by Erhard and Albe.
-There are slight differences comparted to the original potential by Tersoff:
+The implementation of the force evaluation shown in the following is applied to the potential designed by Erhard and Albe \cite{albe_sic_pot}.
+There are slight differences compared to the original potential by Tersoff:
\begin{itemize}
\item Difference in sign of the attractive part.
\item $c$, $d$ and $h$ values depend on atom $k$ in addition to atom $i$.
\begin{itemize}
\item Difference in sign of the attractive part.
\item $c$, $d$ and $h$ values depend on atom $k$ in addition to atom $i$.
\item The exponent of the $b$ term is constantly $-\frac{1}{2}$.
\end{itemize}
These differences actually slightly ease code realization.
\item The exponent of the $b$ term is constantly $-\frac{1}{2}$.
\end{itemize}
These differences actually slightly ease code realization.
\begin{figure}
\renewcommand\labelitemi{}
\renewcommand\labelitemii{}
\renewcommand\labelitemiii{}
\begin{figure}
\renewcommand\labelitemi{}
\renewcommand\labelitemii{}
\renewcommand\labelitemiii{}
\item calculate: $r_{ij}$, $r_{ij}^2$
\item IF $r_{ij} > S_{ij}$ THEN CONTINUE
\item
\item calculate: $r_{ij}$, $r_{ij}^2$
\item IF $r_{ij} > S_{ij}$ THEN CONTINUE
\item