4.2 Relax the system

Once your system is in the GS, depending on how you have prepared the starting atomic configuration:

1. if you have set the atomic positions "by hand" and/or from a classical code, check the forces on atoms, and if they are large ( $\sim$ 0.1 ÷ 1.0 atomic units), you should perform an ionic minimization, otherwise the system could break up during the dynamics.
2. if you have taken the positions from a previous run or a previous ab-initio simulation, check the forces, and if they are too small ( $\sim$ 10^-4 atomic units), this means that atoms are already in equilibrium positions and, even if left free, they will not move. Then you need to randomize positions a little bit (see below).

Let us consider case 1). There are different strategies to relax the system, but the most used are again steepest-descent or damped-dynamics for ions and electrons. You could also mix electronic and ionic minimization scheme freely, i.e. ions in steepest-descent and electron in with damped-dynamics or vice versa.

(a)	suppose we want to perform steepest-descent for ions. Then we should specify the following section for ions: &ions ion_dynamics = 'sd' / Change also the ionic masses to accelerate the minimization: ATOMIC_SPECIES C 2.0d0 c_blyp_gia.pp H 2.00d0 h.ps while leaving other input parameters unchanged. Note that if the forces are really high (> 1.0 atomic units), you should always use steepest descent for the first ($\sim$ 100 relaxation steps.
(b)	As the system approaches the equilibrium positions, the steepest descent scheme slows down, so is better to switch to damped dynamics: &ions ion_dynamics = 'damp', ion_damping = 0.2, ion_velocities = 'zero' / A value of ion_damping around 0.05 is good for many systems. It is also better to specify to restart with zero ionic and electronic velocities, since we have changed the masses. Change further the ionic masses to accelerate the minimization: ATOMIC_SPECIES C 0.1d0 c_blyp_gia.pp H 0.1d0 h.ps
(c)	when the system is really close to the equilibrium, the damped dynamics slow down too, especially because, since we are moving electron and ions together, the ionic forces are not properly correct, then it is often better to perform a ionic step every N electronic steps, or to move ions only when electron are in their GS (within the chosen threshold). This can be specified by adding, in the ionic section, the ion_nstepe parameter, then the &IONS namelist become as follows: &ions ion_dynamics = 'damp', ion_damping = 0.2, ion_velocities = 'zero', ion_nstepe = 10 / Then we specify in the &CONTROL namelist: etot_conv_thr = 1.d-6, ekin_conv_thr = 1.d-5, forc_conv_thr = 1.d-3 As a result, the code checks every 10 electronic steps whether the electronic system satisfies the two thresholds etot_conv_thr, ekin_conv_thr: if it does, the ions are advanced by one step. The process thus continues until the forces become smaller than forc_conv_thr. Note that to fully relax the system you need many runs, and different strategies, that you should mix and change in order to speed-up the convergence. The process is not automatic, but is strongly based on experience, and trial and error. Remember also that the convergence to the equilibrium positions depends on the energy threshold for the electronic GS, in fact correct forces (required to move ions toward the minimum) are obtained only when electrons are in their GS. Then a small threshold on forces could not be satisfied, if you do not require an even smaller threshold on total energy.

Let us now move to case 2: randomization of positions.

If you have relaxed the system or if the starting system is already in the equilibrium positions, then you need to displace ions from the equilibrium positions, otherwise they will not move in a dynamics simulation. After the randomization you should bring electrons on the GS again, in order to start a dynamic with the correct forces and with electrons in the GS. Then you should switch off the ionic dynamics and activate the randomization for each species, specifying the amplitude of the randomization itself. This could be done with the following &IONS namelist:

          &ions
            ion_dynamics = 'none',
            tranp(1) = .TRUE.,
            tranp(2) = .TRUE.,
            amprp(1) = 0.01
            amprp(2) = 0.01
          /

In this way a random displacement (of max 0.01 a.u.) is added to atoms of species 1 and 2. All other input parameters could remain the same. Note that the difference in the total energy (etot) between relaxed and randomized positions can be used to estimate the temperature that will be reached by the system. In fact, starting with zero ionic velocities, all the difference is potential energy, but in a dynamics simulation, the energy will be equipartitioned between kinetic and potential, then to estimate the temperature take the difference in energy (de), convert it in Kelvin, divide for the number of atoms and multiply by 2/3. Randomization could be useful also while we are relaxing the system, especially when we suspect that the ions are in a local minimum or in an energy plateau.