# Preconditioning self-consistent iteration

Publication arXiv:0804.2583 describes a method for doing self-consistent iteration without having to diagonalize the Hamiltonian operator at every step.

### IX. PRECONDITIONING

As already mentioned, direct minimisation of the total energy is used to find the converged wavefunctions. The gradient $$g_i$$ of the total energy with respect to the $$i\text{-th}$$ wavefunction $$\left| \Psi_i \right>$$ is given by $$\left| g_i \right> = H \left| \Psi_i \right> - \sum_{j}{\Lambda_{ij} \left| \Psi_j \right>} \,, \tag{26}$$ where $$\Lambda_{ij} = \left< \psi_j \middle| H \middle| \psi_i \right>$$ are the Lagrange multipliers enforcing the orthogonality constraints. Convergence is achieved when the average norm of the residue $$\sideset{}{^{1/2}}{\langle \overline{g_i | g_i} \rangle}$$ is below an user-defined numerical tolerance.

Given the gradient direction at each step, several algorithms can be used to improve convergence. In our method we use either preconditioned steepest-descent algorithm or preconditioned DIIS method (27; 28). These methods work very well to improve the convergence for non-zero gap systems if a good preconditioner is available.

The preconditioning gradient $$\left| \tilde{g}_{i} \right>$$ which approximately points in the direction of the minimum is obtained by solving the linear system of equations obtained by discretizing the equation $$\left(\frac{1}{2} \nabla^2 - \epsilon_i\right) \tilde{g}_i\left(\mathbf{r}\right) = g_i \left(\mathbf{r}\right) \,. \tag{27}$$ The values $$\epsilon_i$$ are approximate eigenvalues obtained by a subspace diagonalization in a minimal basis of atomic pseudopotential orbitals during the generation of the input guess. For isolated systems, the values of the $$\epsilon_i$$ for the occupied states are always negative, therefore the operator of Eq. (27) is positive definite.

Eq. (27) is solved by a preconditioned conjugate gradient (CG) method. The preconditioning is done by using the diagonal elements of the matrix representing the operator $$\frac{1}{2} \nabla^2 − \epsilon_i$$ in a scaling function-wavelet basis. In the initial step we use $$\ell$$ resolution levels of wavelets where $$\ell$$ is typically 4. To do this we have to enlarge the domain where the scaling function part of the gradient is defined to a grid that is a multiple of $$2^\ell .$$ This means that the preconditioned gradient $$\tilde{g}_i$$ will also exist in a domain that is larger than the domain of the wavefunction $$\Psi_i .$$ Nevertheless this approach is useful since it allows us to obtain rapidly a preconditioned gradient that has the correct overall shape. In the following iterations of the conjugate gradient we use only one wavelet level in addition to the scaling functions for preconditioning. In this way we can do the preconditioning exactly in the domain of basis functions that are used to represent the wavefunctions (Eq. 5). A typical number of CG iterations necessary to obtain a meaningful preconditioned gradient is 5.

Question: How is the preconditioning formula (27) derived?

• Welcome to SE.Physics! Since the content being asked about seemed vaguely referenced, I quoted the excerpt from the paper. It turned out a bit longer than I initially intended; please feel free to reduce it if it's excessive. – Nat Jun 19 '19 at 5:26