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)
+ \nabla_{{\bf r}_j}\Big(\frac{{\bf r}_{ij}{\bf r}_{ik}}{r_{ij}r_{ik}}\Big)
\nonumber \\
&=& \frac{1}{r_{ij}r_{ik}}{\bf r}_{ik} -
\frac{\cos\theta_{ijk}}{r_{ij}^2}{\bf r}_{ij}
\item
\item LOOP k \{
\begin{itemize}
- \item set $ik$-depending values
+ \item set $ik$-dependent values
\item calculate: $r_{ik}$, $r_{ik}^2$
\item IF $r_{ik} > S_{ik}$ THEN CONTINUE
\item calculate: $\theta_{ijk}$, $\cos(\theta_{ijk})$,