     Next: How to set up Up: Step-by-step description of calculational Previous: Choice of the k-point

# Total Energy Minimization Schemes

In an electronic structure calculation using a plane-wave basis, the Hilbert space is typically spanned by a huge number of basis functions (up to 105 plane waves). Therefore it would be unwise to attempt to diagonalize the Hamiltonian operator in this high-dimensional space directly. Instead, one uses algorithms which only imply vector operations on the wave function vector (in Hilbert space), rather than matrix operations. The wave functions are gradually improved in an iterative process, until they eventually converge towards the eigenvectors.

• Electronic minimization
The goal is to minimize the total energy with respect to the wave function starting with a trial wave function . The energy minimization scheme is formulated in terms of an equation of motion for the wave function in a fictitious time variable t.

• Steepest Descent
The simplest scheme to iterate the wave functions is the steepest descent approach . It can be derived from a first-order equation of motion imposing the ortho-normality constraint , where is the Kohn-Sham Hamiltonian and are the Lagrange multipliers introduced to account for the ortho-normality constraint.
In the simplest possible discretization of this differential equation, only information from the last step is used, with and . However, it turns out that this discretization scheme is not very efficient.

• Damped Joannopoulos
A more efficient scheme based on a second order equation of motion might also be used (3.1)

where is a damping parameter. The equation of motion is integrated for a step length by the Joannopoulos approach , which iteratively improves the initial wave functions. In this algorithm the new wave function is constructed from the wave functions of the last two iteration steps t and (t-1) where the coefficients are with . The function is defined by with .

• Williams-Soler
Although the damped Joannopoulos algorithm is more efficient than the first order scheme, additional storage for the wave function is needed. Therefore the Williams-Soler algorithm  is recommended whenever storage requirements do not permit to employ the damped Joannopoulos algorithm. The coefficients of this scheme are with . Thus, the damped Joannopoulos scheme contains the Williams-Soler scheme as a limiting case, when . On the other hand, the Williams-Soler scheme itself approaches the steepest descent scheme, if is sufficiently small.

• The choice of and depends on the atomic species and the configuration. The corresponding parameters in the inp.mod  input file are delt and gamma. Typically delt lies between 1 and 40 and gamma is within the range .
• Set up the parameters delt, gamma, delt2 and gamma2 in the file inp.mod. If the improvement in the total energy per   iteration   is   less than the parameters delt2 and gamma2 are used instead, in order to ensure the stability of convergence. Note that the values for delta2 and gamma2 should be smaller than that ones for delta and gamma.

 Parameter Type/Value delt real > 1 Step length of electronic iteration gamma If 2 (see below), damping parameter delt2 see delt c.f. gamma2 see gamma c.f.  real >0 If the total energy varies less than , delt2 and gamma2 replace delt and gamma

• Choose the desired scheme to iterate the wave functions in the input file inp.mod, setting up the parameter according to the following numbers shown below.

 Parameter Value Schemes to iterate the wave functions 0 steepest descent 1 Williams-Soler 2 damped Joannopoulos

• Example

• As an example we consider a total energy calculation for GaAs bulk in the zinc-blende structure. The Williams-Soler minimization scheme is used to iterate the wave functions and different values for the electronic time step delt were used. The results are shown below. Note how the convergence can be accelerated by changing the electronic time step delt. The input files start.inp and inp.mod used for these calculations are shown below.

start.inp

1                      : number of processors
1                      : number of minimal processors
1                      : number of processors per group
2                      : number of species
0                      : excess electrons
4                      : number of empty states
2  0                   : ibrav pgind
10.68 0.0 0.0  1 1 1   : celldm
1                      : number of k-points
0.5 0.5 0.5 1.0        : k-point coordinates, weight
4 4 4                  : fold parameter
.true.                 : k-point coordinates relative or absolute?
8  4.0                 : Ecut [Ry], Ecuti [Ry]
0.1   .true. .f.       : ekt tmetal tdegen
.true. .false. 1       : tmold tband nrho
5  2  1234             : npos  nthm  nseed
873.0 1400.0 10000000 1: T_ion T_init Q nfi_rescale
.true. .false.         : tpsmesh coordwave
1 5 'Arsenic'          : number of atoms, zv, name
1.0 74.92  0.7  3  3   : gauss radius, mass, damping,l_max,l_loc
.t. .t. .f.            : t_init_basis
0.0 0.0 0.0  .f.       : tau0 tford
1 3 'Gallium'          : number of atoms, zv, name
1.0 69.72  0.7  3  3   : gauss radius, mass, damping,l_max,l_loc
.t. .t. .f.            : t_init_basis
0.25 0.25 0.25  .f.    : tau0 tford

inp.mod

-1  100 1000000              : nbeg   iprint timequeue
100 1                        : nomore nomore_init
10.0  0.2                    : delt  gamma
4.0   0.3  0.0001            : delt2 gamma2 eps_chg_dlt
400 2                        : delt_ion nOrder
0.0 1.0                      : pfft_store mesh_accuracy
2 2                          : idyn i_edyn
0 .false.                    : i_xc t_postc
.F. 0.001 .F. 0.002          : trane ampre tranp amprp
.false. .false. .false. 1800 : tfor tdyn tsdp nstepe
.false.                      : tdipol
0.0001 0.0005 0.3            : epsel epsfor epsekinc
0.001 0.001 3                : force_eps max_no_force
2                            : init_basis • We attach a table containing some values to the electronic time step delt, the damping parameter gamma and the minimization method adopted which have been successfully used. Of course this table provides only a few numbers and should be used only as rough guide as everyone has to find her/his "optimum" parameters.
• For an example of a geometry optimization calculation see Section: How to set up a structural relaxation run.      Next: How to set up Up: Step-by-step description of calculational Previous: Choice of the k-point
Matthias Wittenberg
1999-08-06