Bond Order Potentials

There is an alternative approach to evaluating the forces to that described above. Rather than trying to differentiate the energy exactly, we can appeal to the Hellmann-Feynman theorem[33,34], which is given in equation 10. The forces are clearly straightforward to evaluate if we can find the density matrix. The bond order potential method is a scheme for evaluating the density matrix. It should be noted that the forces given by equation 10 will only equal the derivatives of the energy once the density matrix is well converged.

To obtain the forces on the atoms, only the elements of the density matrix within the same range as that of the Hamiltonian matrix are required. These can be evaluated from the off-diagonal elements of the Green's function:

$$\rho_{i\alpha,j\beta} = -\frac{1}{\pi}\lim_{\eta\to 0} {\rm Im} \int dE G_{i\alpha,j\beta}(E+{\rm i}\eta)f(x).$$ (37)

It has been shown that in general,

$$\rho_{i\alpha,j\beta} = -\sum^\infty_{n=0} \chi_{0n,n0}^\Lam... ...hi_{0n,(n-1)0}^\Lambda (\delta b_n^\Lambda)_{i\alpha,j\beta},$$ (38)

where the response functions for a given total number of electrons, N_e, and electron temperature[12], T, are defined by

$$\chi_{0m,n0}^\Lambda(N_e,T) = \frac{1}{\pi}\lim_{\eta\to 0} {\rm Im} \int G_{0m}(E+{\rm i}\eta)G_{n0}(E+{\rm i}\eta)f(x)dE.$$ (39)

The Green's functions $G_{0m}^\Lambda(Z)$ are defined along the Lanczos recursion chain and follow the relation:

$$(Z-a_n^\Lambda)G_{nm}^\Lambda(Z)-b_n^\Lambda G_{n-1,m}^\Lambda(Z) -b_{n+1}^\Lambda G_{n+1,m}^\Lambda(Z) = \delta_{n,m}.$$ (40)

The symbol $\Lambda$ denotes the initial vector for the Lanczos recursion algorithm (c.f. $\vert U_0\rangle$ introduced in Section 4.2) in the augmented vector space spanned by the direct product between the atomic orbitals basis and the set of auxiliary vectors (for more information, see references [7] and [12]). In other words, the reference Green's function for the bond order expansion is uniquely specified by the choice of $\Lambda$ matrix as

$$G_{00}^\Lambda = \sum_{i\alpha,j\beta} G_{i\alpha,j\beta}\Lambda_{i\alpha,j\beta}.$$ (41)

This leads naturally to one-site, two-site or even all-site expansions for the Green's functions:

$$G_{i\alpha,j\beta}= \frac{\partial G_{00}^\Lambda}{\partial \Lambda_{i\alpha,j\beta}}.$$ (42)

This holds for any $\Lambda$.

The derivatives of the $\Lambda$-dependent recursion coefficients $(\delta a_n^\Lambda)_{i\alpha,j\beta}$ and $(\delta b_n^\Lambda)_{i\alpha,j\beta}$ are given by

 

$$(\delta a_n^\Lambda)_{i\alpha,j\beta} {{\partial a_n^\Lambda}\over{\partial\Lambda_{i\alpha,j\beta}}} =... ...ambda}{\partial\mu_s^\Lambda} \langle i\alpha\vert\hat{H}^s\vert j\beta\rangle,$$ =

$$(\delta b_n^\Lambda)_{i\alpha,j\beta} {{\partial b_n^\Lambda}\over{\partial\Lambda_{i\alpha,j\beta}}} =... ...ambda}{\partial\mu_s^\Lambda} \langle i\alpha\vert\hat{H}^s\vert j\beta\rangle,$$ = (43)

where $\mu_s^\Lambda=\sum_{i\alpha,j\beta}\Lambda_{i\alpha,j\beta} \langle i\alpha\vert\hat{H}^s\vert j\beta\rangle$. The expressions (43) provide an intuitive scheme for understanding the physical origin of the bond order through the types of atomic orbitals and the local geometry of the atoms. In practice, however, it is much more efficient and numerically stable to use the technique of O-matrices[12] than using equation 43.

Closed form expressions for the factors $\delta a^{i\alpha,j\beta}_n$ and $\delta b^{i\alpha,j\beta}_n$ are given as follows. Let us introduce the orthogonal polynomials, $P^\Lambda_n(x)$, which satisfy:

$$xP^\Lambda_n(x) = a_nP^\Lambda_n(x) + b_nP^\Lambda_{n-1}(x)+b_{n+1}P^\Lambda_{n+1}(x),$$ (44)

where $P^\Lambda_{-1}(x) = 0$ and $P^\Lambda_0(x) = 1$. If the O-matrix is defined as

$$O^{\Lambda,m,n}_{i\alpha,j\beta} = \langle i\alpha\vert P^\Lambda_m(\hat{H})P^\Lambda_n(\hat{H}) \vert j\beta \rangle,$$ (45)

then it can be shown[12] that the derivatives of $a^\Lambda_n$ and $b^\Lambda_n$ can be expressed as:

 

$${{\partial a^\Lambda_n}\over{\partial\Lambda_{i\alpha,j\beta}}} b^\Lambda_{n+1}O^{\Lambda,n+1,n}_{i\alpha,j\beta}- b^\Lambda_nO^{\Lambda,n,n-1}_{i\alpha,j\beta}$$ =

$${{\partial b^\Lambda_n}\over{\partial\Lambda_{i\alpha,j\beta}}} {1\over 2}b^\Lambda_n\left(O^{\Lambda,n,n}_{i\alpha,j\beta}- O^{\Lambda,n-1,n-1}_{i\alpha,j\beta}\right).$$ = (46)

This is a simple result for the derivatives of the recursion coefficients.

To complete the bond order expansion in a finite number of terms for bulk systems, so-called truncators[7,12] have to be included in order to ensure the equivalence of the band energy calculated from the density matrix and calculated from the density of states. They substitute for the derivatives of the last recursion coefficients, $a_N^\Lambda$ and $b_{N+1}^\Lambda$, that remain unknown in the N-level Lanczos recursion. The set of truncators used in the present work has the form:

 

$$\left(\delta a_N^\Lambda\right)_{i\alpha,j\beta} = \left(\frac{\l... ...ft(b_N^\Lambda\right)^2} {b_N^\Lambda}\right)O_{i\alpha,j\beta}^{\Lambda,N,N-1}$$

$$\left(\delta b_{N+1}^\Lambda\right)_{i\alpha,j\beta} = {1\over 2}... ...a - a_N^\Lambda\right)}{b_N^\Lambda}\right) O_{i\alpha,j\beta}^{\Lambda,N,N-1}.$$ (47)

This novel form of truncator is the simplest possible construction, implied by the recursion relation that O-matrices obey:

 

$$b^{\Lambda}_{m+1}O^{\Lambda,m+1,n}_{i\alpha,j\beta} b^{\Lambda}_{n}O^{\Lambda,m,n-1}_{i\alpha,j\beta} + (a^{\Lambda}_n-a^{\Lambda}_m)O^{\Lambda,m,n}_{i\alpha,j\beta}$$ =

$$b^{\Lambda}_{n+1}O^{\Lambda,m,n+1}_{i\alpha,j\beta} - b^{\Lambda}_mO^{\Lambda,m-1,n}_{i\alpha,j\beta}.$$ + (48)

The key features of this method are, then: it provides a rapidly convergent expansion for the density matrix in terms of moments; it requires only modest amounts of memory to implement; there are errors in the forces that follow from dependence on the Hellmann-Feynman theorem. This method has been implemented on a parallel machine, and almost perfect scaling is found, as little information needs to be passed between processors.