Computational Physics

9.4 Orthonormalisation; conjugate gradient and RM-DIIS techniques 283

Now we may proceed in two ways:

(i) Calculate as a first guess (^) kl= 2 〈ψk(t)|H|ψl(t)〉−μ〈ψ ̇ ̃k(t)|ψ ̇ ̃l(t)〉and
calculate|ψ(t+h)〉by the algorithm(9.50). Then refine the solution thus
obtained by orthonormalisation.
(ii) Calculate the new solution by the algorithm(9.50)withoutthe constraint
term. Then calculate the Lagrange parametersXijby the above algorithm.
The constraint term is then added to the solution.
It is possible to turn the Verlet algorithm into a velocity Verlet form. In that case,
not only the new orbitals are orthonormalised, but also the constraint〈ψ ̇k(t)|ψl(t)〉
is satisfied rigorously at all times. The algorithm requires a bit more work and
storage. Its advantage is that it can be extended more easily within a Nosé–Hoover
thermostatscheme–fordetailsseeRefs.[ 9 , 13 ].Wegivethealgorithmherefor
completeness:
|ψ ̃ ̇k(t+h)〉=|ψ ̇k(t)〉+
h
μ
H|ψk(t)〉; (9.69a)
|ψ ̃k(t+h)〉=|ψk(t)〉+h|ψ ̇ ̃k(t+h)〉; (9.69b)
|ψk(t+h)〉=|ψ ̃k(t+h)〉+

∑

ij

Xij|ψk(t)〉; (9.69c)

|ψ ̇k′(t+h)〉=|ψ ̇ ̃k(t+h)〉+

h

μ

H|ψk(t+h)〉; (9.69d)

|ψ ̇k(t+h)〉=|ψ ̇k′(t+h)〉+

∑

ij

Yij|ψk(t+h)〉. (9.69e)

The matrixXijis the same as used in the standard Verlet procedure above; the matrix
Yijis simply calculated in terms of the matrixCkl=〈ψk(t+h)|ψ ̇l′(t+h)〉as

Ykl=−

Ckl+C†kl

2

. (9.70)

If we include friction in the equations of motion, we are allowed more freedom,
as the only requirement is that the orbits will become stationary by some damping

mechanism. In this case, one can take the (^) klto be diagonal:
(^) kl=εkδkl (9.71)
and the equation of motion for the stationary state leads directly to
εk=〈ψk|H|ψk〉, (9.72)
and this form of the Lagrange parameters is then used throughout the simulation, that
is, even when if the orbitals move. However, this form of the Lagrange multipliers