+If not only pair potentials are considered, three body potentials $U_3$ or multi body potentials $U_n$ can be included.
+Usually these higher order terms are avoided since they are not easy to model and it is rather time consuming to evaluate potentials and forces originating from these many body terms.
+Ordinary pair potentials have a close-packed structure like face-centered cubic (FCC) or hexagonal close-packed (HCP) as a ground state.
+A pair potential is, thus, unable to describe properly elements with other structures than FCC or HCP.
+Silicon and carbon for instance, have a diamand and zincblende structure with four covalently bonded neighbors, which is far from a close-packed structure.
+A three body potential has to be included for these types of elements.
+
+\subsubsection{The Tersoff bond order potential}
+
+Tersoff proposed an empirical interatomic potential for covalent systems \cite{tersoff_si1,tersoff_si2}.
+The Tersoff potential explicitly incorporates the dependence of bond order on local environments, permitting an improved description of covalent materials.
+Due to the covalent character Tersoff restricted the interaction to nearest neighbor atoms accompanied by an increases in computational efficiency for the evaluation of forces and energy based on the short-range potential.
+Tersoff applied the potential to silicon \cite{tersoff_si1,tersoff_si2,tersoff_si3}, carbon \cite{tersoff_c} and also to multicomponent systems like silicon carbide \cite{tersoff_m}.
+The basic idea is that, in real systems, the bond order, i.e. the strength of the bond, depends upon the local environment \cite{abell85}.
+Atoms with many neighbors form weaker bonds than atoms with only a few neighbors.
+Although the bond strength intricately depends on geometry the focus on coordination, i.e. the number of neighbors forming bonds, is well motivated qualitatively from basic chemistry since for every additional formed bond the amount of electron pairs per bond and, thus, the strength of the bonds is decreased.
+If the energy per bond decreases rapidly enough with increasing coordination the most stable structure will be the dimer.
+In the other extreme, if the dependence is weak, the material system will end up in a close-packed structure in order to maximize the number of bonds and likewise minimize the cohesive energy.
+This suggests the bond order to be a monotonously decreasing function with respect to coordination and the equilibrium coordination being determined by the balance of bond strength and number of bonds.
+Based on pseudopotential theory the bond order term $b_{ijk}$ limitting the attractive pair interaction is of the form $b_{ijk}\propto Z^{-\delta}$ where $Z$ is the coordination number and $\delta$ a constant \cite{abell85}, which is $\frac{1}{2}$ in the seond-moment approximation within the tight binding scheme \cite{horsfield96}.
+
+Tersoff incorporated the concept of bond order in a three-body potential formalism.
+The interatomic potential is taken to have the form
+\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{ .}
+\end{eqnarray}
+$E$ is the total energy of the system, constituted either by the sum over the site energies $E_i$ or by the bond energies $V_{ij}$.
+The indices $i$ and $j$ correspond to the atoms of the system with $r_{ij}$ being the distance from atom $i$ to atom $j$.
+The functions $f_R$ and $f_A$ represent a repulsive and an attractive pair potential.
+The repulsive part is due to the orthogonalization energy of overlapped atomic wave functions.
+The attractive part is associated with the bonding.
+\begin{eqnarray}
+f_R(r_{ij}) & = & A_{ij} \exp (- \lambda_{ij} r_{ij} ) \\
+f_A(r_{ij}) & = & -B_{ij} \exp (- \mu_{ij} r_{ij} )
+\end{eqnarray}
+The function $f_C$ is the a cutoff function to limit the range of interaction to nearest neighbors.
+It is designed to have a smooth transition of the potential at distances $R_{ij}$ and $S_{ij}$.
+\begin{equation}
+f_C(r_{ij}) = \left\{
+ \begin{array}{ll}
+ 1, & r_{ij} < R_{ij} \\
+ \frac{1}{2} + \frac{1}{2} \cos \Big[ \pi (r_{ij} - R_{ij})/(S_{ij} - R_{ij}) \Big], & R_{ij} < r_{ij} < S_{ij} \\
+ 0, & r_{ij} > S_{ij}
+ \end{array} \right.
+\label{eq:basics:fc}
+\end{equation}
+As discussed above, $b_{ij}$ represents a measure of the bond order, monotonously decreasing with the coordination of atoms $i$ and $j$.
+It is of the form:
+\begin{eqnarray}
+b_{ij} & = & \chi_{ij} (1 + \beta_i^{n_i} \zeta^{n_i}_{ij})^{-1/2n_i} \\
+\zeta_{ij} & = & \sum_{k \ne i,j} f_C (r_{ik}) \omega_{ik} g(\theta_{ijk}) \\
+g(\theta_{ijk}) & = & 1 + c_i^2/d_i^2 - c_i^2/[d_i^2 + (h_i - \cos \theta_{ijk})^2]
+\end{eqnarray}
+where $\theta_{ijk}$ is the bond angle between bonds $ij$ and $ik$.
+This is illustrated in Figure \ref{img:tersoff_angle}.
+\begin{figure}[t]
+\begin{center}
+\includegraphics[width=8cm]{tersoff_angle.eps}
+\end{center}
+\caption{Angle between bonds of atoms $i,j$ and $i,k$.}
+\label{img:tersoff_angle}
+\end{figure}
+The angular dependence does not give a fixed minimum angle between bonds since the expression is embedded inside the bond order term.
+The relation to the above discussed bond order potential becomes obvious if $\chi=1, \beta=1, n=1, \omega=1$ and $c=0$.
+Parameters with a single subscript correspond to the parameters of the elemental system \cite{tersoff_si3,tersoff_c} while the mixed parameters are obtained by interpolation from the elemental parameters by the arithmetic or geometric mean.
+The elemental parameters were obtained by fit with respect to the cohesive energies of real and hypothetical bulk structures and the bulk modulus and bond length of the diamond structure.
+New parameters for the mixed system are $\chi$, which is used to finetune the strength of heteropolar bonds, and $\omega$, which is set to one for the C-Si interaction but is available as a feature to permit the application of the potential to more drastically different types of atoms in the future.
+
+The force acting on atom $i$ is given by the derivative of the potential energy.
+For a three body potential ($V_{ij} \neq V{ji}$) the derivation is of the form
+\begin{equation}
+\nabla_{{\bf r}_i} E = \frac{1}{2} \big[ \sum_j ( \nabla_{{\bf r}_i} V_{ij} + \nabla_{{\bf r}_i} V_{ji} ) + \sum_k \sum_j \nabla_{{\bf r}_i} V_{jk} \big] \textrm{ .}
+\end{equation}
+The force is then given by
+\begin{equation}
+F^i = - \nabla_{{\bf r}_i} E \textrm{ .}
+\end{equation}
+Details of the Tersoff potential derivative are presented in appendix \ref{app:d_tersoff}.