17. Structure relaxation calculation
In this section, we show the rudiments of relaxation calculations.
17.1. Example: Si crystal bond length
To perform a structure relaxation, set the keyword
info.calculationType to relaxation. Let’s consider a silicon
primitive cell in which one atom is out of place, a calculation which
boils down to computing the Si crystal bond length. Create a text file
si_real_rlx containing the following keyword definitions.
info.savepath = 'results/si_real_rlx'
info.calculationType = 'relaxation'
atom.element = [1 1]
atom.fracxyz = [0 0 0;0.26 0.27 0.28]
domain.latvec = 10.25*[0.0 0.5 0.5;0.5 0.0 0.5;0.5 0.5 0.0]
domain.lowres = 0.4
kpoint.gridn = [3 3 3]
element(1).species = 'Si'
element(1).path = './Si_TM_LDA.mat'
functional.list = {'XC_LDA_X','XC_LDA_C_PW'}
sr.method = 'qn'
sr.tol = 2e-4
sr.relaxedstructure = 'si2_relaxed.xyz'
eigensolver.emptyBand = 1
Note that keywords setting the relaxation engine and other relaxation
related parameters begin with sr.
We described the rationale behind the above input file.
info.calculationTypeis set to ’relaxation’ so that RESCU optimizes the structure;atom.fracxyzcontains the atoms’ fractional coordinates. Here, the second atom is displaced by [0.01 0.02 0.03] from it’s equilibrium position at [0.25 0.25 0.25];sr.methodis a string determining the algorithm used to predict the equilibrium positions. The allowed values are ’sd’, ’cg’ and ’qn’, which stand for steepest descent, conjugate gradient and quasi-Newton respectively. The default value is ’qn’. One may also use ’fire’ to use the FIRE algorithm, a variation of damped molecular dynamics;sr.tolsets the tolerance for the structure relaxation engine in atomic units. More precisely, RESCU stops the structure optimization when all force coordinates fall below the tolerance. The default value is 1e-4, which corresponds to roughly 0.005 eV/Å;sr.relaxedstructureis a string containing the path to the file which contains the relaxed atomic coordinates. The default value is./relaxedstructure.xyz;eigensolver.emptyBandcontrols the number of unoccupied bands included in a calculation. Although the ground state density includes only occupied bands, the inclusion of unoccupied bands is advisable (e.g. for numerical stability) or crucial (e.g. for metallic materials). The default value depends on system size, but it is rather ‘safe’, i.e. rather large. Since relaxation calculations are computationally demanding, it is best to use the minimal amount of bands possible. For a simple silicon crystal, a single extra band will do.
Various additional parameters are described by the function
inputDescription. One parameter worth mentioning is
sr.moveableAtomList, which controls which atom coordinates are
included in the relaxation’s engine degrees of freedom. When
sr.moveableAtomList is not set, RESCU relaxes all atoms except the
first one, taking advantage of translational invariance to eliminate
three degrees of freedom.
The calculation is executed as follows.
rescu -i si_real_rlx
During the relaxation procedure, RESCU performs several self-consistent
calculations whose output goes according to info.outfile as usual.
The relaxation procedure itself may be monitored via the file
relaxationprocess.txt, which looks as follows. If the positions or
forces vary too wildly, it is thus possible to stop the execution and
prevent waste of computational resources.
-----------------------------------------------------------------------------
Relaxation Iteration: 000
Total energy = -7.905e+00, normi(force) = +2.559e-02
-----------------------------------------------------------------------------
Atomi Spec: xCrd yCrd zCrd xFrc yFrc zFrc
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
1 Si: +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00
2 Si: +2.818e+00 +2.767e+00 +2.716e+00 -2.559e-02 -1.781e-02 -8.477e-03
-----------------------------------------------------------------------------
...
-----------------------------------------------------------------------------
Line search iteration: 001
Total energy = -7.912e+00, normi(force) = +1.712e-04
sk.force = -1.834e-04 (norm(sk) = +5.559e-03, alpha = +1.000e+00)
-----------------------------------------------------------------------------
Atomi Spec: xCrd yCrd zCrd xFrc yFrc zFrc
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
1 Si: +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00 +0.000e+00
2 Si: +2.563e+00 +2.562e+00 +2.563e+00 -1.581e-04 -1.747e-05 -1.712e-04
-----------------------------------------------------------------------------
Structure relaxation calculation finished!
Total time: 9e+01 seconds
-----------------------------------------------------------------------------
At each iteration of the relaxation engine, RESCU prints the iteration
number, the total energy and norm of the forces in atomic units. The
forces associated with frozen degrees of freedom (see section
1.2) are set to zero. The norm of the forces is the
infinity norm, which is essentially the force coordinate with the
largest magnitude. RESCU also prints the atom index, the atom species,
the coordinates and the forces (in Cartesian coordinates in atomic
units) for every atom. The relaxation engines proceed by so-called line
searches (except for the FIRE algorithm). For each iteration of a line
search, RESCU prints similar information, but it will additionally print
the force projected on the line search direction (norm(sk)) and the
current line search step (alpha).
When the relaxation procedure is completed, RESCU prints the resulting atomic coordinates to an XYZ file. For instance,
-50pt-50pt
2
AtomType X Y Z
Si +0.000000000000000e+00 +0.000000000000000e+00 +0.000000000000000e+00
Si +1.356583441021671e+00 +1.356077852198500e+00 +1.356628976586721e+00
The atomic positions are Cartesian coordinates in units of Angstrom.
17.2. Constrained relaxation
It is often advantageous to choose which degrees of freedom should be
relaxed as opposed to which atoms. In RESCU, the three Cartesian
coordinates of every atom are the independent degrees of freedom. Which
is to be relax is again specified via the keyword
sr.moveableAtomList, only that this time it is a
\(n_{\text{atom}}\times 3\) Boolean array. It is usually
inconvenient to write such an array in the main input file. Like the
atomic coordinates, which are specified in a file given by atom.xyz,
the active degrees of freedom may be specified in an XYZ file given by
sr.moveableAtomList. RESCU reads columns 8 to 10 in that file and
fixes the value of the degrees of freedom marked 0. Note that the first
column contains the element label, the columns 2 to 4 the atomic
coordinates and the columns 5 to 7 the initial magnetic moment, and
hence the need for 10 columns. Depending on the spin treatment, the
initial magnetic moment may be ignored, but it must be included
nevertheless. For example, the following xyz-file will minimize
\(f_x^1\), \(f_y^2\) and \(f_z^3\) (where the superscript
denotes the atom index)
8
#element x y z sx sy sz fx fy fz
Si 0.26 0 0 0 0 0 1 0 0
Si 0 4.9052 5.1052 0 0 0 0 1 0
Si 5.1052 0 5.3052 0 0 0 0 0 1
Si 5.1052 5.1052 0 0 0 0 0 0 0
Si 2.5526 2.5526 2.5526 0 0 0 0 0 0
Si 2.5526 7.6578 7.6578 0 0 0 0 0 0
Si 7.6578 2.5526 7.6578 0 0 0 0 0 0
Si 7.6578 7.6578 2.5526 0 0 0 0 0 0