[Pw_forum] Davidson algorithm

[Paolo Giannozzi]

On Nov 2, 2007, at 21:08 , Nichols A. Romero wrote:

> Can anyone recommend a basic reference (e.g. thesis, personal notes,
> computional article) that explains the Davidson algorithm used in  
> PWSCF
> in more detail?

The basic reference is E. R. Davidson, "The iterative calculation of  
a few
of the lowest eigenvalues and corresponding eigenvectors of large
real-symmetric matrices," Comput. Phys. 17, 87-94 (1975). The only
PWscf-specific stuff I have is the following
---
\section{Iterative diagonalization}

Iterative diagonalization can be used whenever

i) the number of states to be calculated is much smaller than the
dimension of the basis set, and

ii) a reasonable and economical
estimate of the inverse operator $H^{-1}$ is available.

Both conditions are satisfied in practical calculation in a PW
basis set: the number of PW's is usually much larger than the
number of bands, and the Hamiltonian matrix
is dominated by the kinetic energy at large $\G$ (
the Hamiltonian is {\em diagonally dominant}).

Iterative methods are based
on a repeated refinement of a trial solution,
which is stopped when satisfactory convergence is achieved.
The number of iterative steps cannot be predicted in advance.
It depends heavily on the structure of the matrix,
on the type of refinement used, and on the starting point.
A well-known and widely used algorithm is due to
Davidson. In this method,
a set of correction vectors
$|\delta\psi_i\rangle$ to the $M$ trial eigenvectors
$|\psi_i\rangle$
are generated as follows:
\begin{equation}
|\delta\psi_i\rangle = {1 \over D-\epsilon_i}
                        (H-\epsilon_i)|\psi_i\rangle
\end{equation}
where the $\epsilon_i=\langle\psi_i|H|\psi_i\rangle$ are the
trial eigenvalues. The $|\delta\psi_i\rangle$'s
are orthogonalized
and the Hamiltonian is diagonalized (with conventional
techniques) in the subspace spanned by the trial and correction
vectors. A new set of trial eigenvectors is obtained and the procedure
is iterated until convergence is achieved.
A good set of starting trial vectors is supplied by the
eigenvectors at the preceding iteration of the potential.

An important point is the following. The Hamiltonian matrix is never
explicitly required excepted for its diagonal part. Only
$H \psi_i$ products are required,
which can be calculated in a very convenient way
by applying the {\em dual-space technique}.
In fact the kinetic energy is diagonal in \G-space, whereas the local
potential term is diagonal in real space. Using FFT's (see below)
one can go
quickly back and forth  from real to reciprocal space
and perform the products where it is more convenient.
There is still a nonlocal term which appears to require the
storage of the matrix.
The trick is to write $V_{NL}$ in a {\em separable}
form:
\begin{equation}
V_{NL}(\k+\G,\k+\G')= \sum_{\mu=1}^{N_{at}}
\sum_{j=1}^n f^\mu_j(\k+\G) g^\mu_j(\k+\G'),
\end{equation}
where $n$ is a small number and $N_{at}$ is the number of
atoms in the unit cell.
This allows us to perform the products by storing only
the $f$ and $g$ vectors.
---
Paolo Giannozzi, Dept of Physics, University of Udine
via delle Scienze 208, 33100 Udine, Italy
Phone +39-0432-558216, fax +39-0432-558222