For all methods the number of strains is

If you know in which direction the atoms relax and if there is only one
atom that moves in each strain type, one can use `frozen_ions=.TRUE.`,
and compute the energy at several atomic positions choosing the one
that minimizes the energy. This approach computes more energy configurations
but does not relax them and should provide results identical to
`frozen_ions=.FALSE.`. It can be used only with `energy`
and `energy_std` algorithms, but should work also with the
free energy in `what=`elastic_constants_geo'`.

At finite pressure all methods give the elastic constants that
relate linearly stress and strain.

The input variables that control this option are:

frozen_ions: if .TRUE. the elastic constants are calculated keeping the ions frozen in the strained positions. Default: logical .FALSE. ngeo_strain: the number of strained configurations used to calculate each derivative. Default: integer 4 ('standard' and 'advanced'), 6 ('energy') elastic_algorithm: 'standard', 'advanced', 'energy_std' or 'energy'. See discussion above. Default: character 'standard' delta_epsilon: the interval of strain values between two geometries. To avoid a zero strain geometry that might have a different symmetry ngeo_strain must be even. Default: real 0.005 epsilon_0: a minimum strain. For small strains the ionic relaxation routine requires a very small threshold to give the correct internal relaxations and sometimes fail to converge. In this case you can increase delta_epsilon, but if delta_epsilon becomes too large you can reach the nonlinear regime. In this case you can use a small delta_epsilon and a minimum strain (To be used only for difficult systems). Default: real 0.0 poly_degree: degree of the polynomial used to interpolate stress or energy. ngeo_strain must be larger than poly_degree+1. Default: 3 ('standard', 'advanced', 2 if ngeo_strain < 6), 4 ('energy', 3 if ngeo_strain < 6). stype(1:strain_types): For each strain type for which this variable is .TRUE. it is possible to compute the energy by moving one atom in one direction and minimizing the energy with respect to the ionic position. This option requires frozen_ions =.TRUE. and should give the same result of frozen_ions=.FALSE. if the atom that moves and its direction is chosen correctly. Default: logical .FALSE.. move_at(1:strain_types): The atom that moves in each strain type. Default: integer none (it must be set if stype(istrain)=.TRUE.). atom_dir(3,1:strain_types): The versor of the atomic movement in Cartesian coordinates (the vector can be given with any module. It is normalized by the code). Default: real none (it must be set if stype(istrain)=.TRUE.). atom_step(1:strain_types): the step for each atomic displacement. Default : real 0.01 a.u.. nmove Number of atomic positions sampled for minimizing the energy at each strain type. Default: integer 5. lcm_ec If .TRUE. the center of mass of the cell is conserved in the displacement introduced with previous variables. Default: logical .TRUE. fl_el_cons: the name of the file that contains the elastic constants. Default: character(len=*) 'output_el_con.dat'

The three algorithms are equivalent only at convergence both with
**k**-point sampling and with the kinetic energy cut-off, but
large differences between the elastic constants obtained with the
`standard` and `advanced` algorithms might point to
insufficient **k**-point sampling. Large differences between the
elastic constants obtained with the `energy_std` or `energy`
algorithms with respect
to the other two might point to insufficient kinetic energy cut-off.

Number of tasks for this option: `ngeo_strain` times the number of
independent strains.

Using the elastic constants tensor the code can calculate and print a few auxiliary quantities: the bulk modulus, the poly-crystalline averages of the Young modulus, of the shear modulus, and of the Poisson ratio. Both the Voigt and the Reuss averages are printed together with the Hill average. The Voigt-Reuss-Hill average of the shear modulus and of the bulk modulus are used to compute average sound velocities. The average of the Poisson ratio and the bulk modulus allow the estimation of the Debye temperature. The Debye temperature is calculated also with the exact formula evaluating the average sound velocity from the angular average of the sound velocities calculated for each propagation direction solving the Christoffel wave equation. The exact Debye temperature is used within the Debye model to calculate the Debye's vibrational energy, free energy, entropy, and constant strain heat capacity. These quantities are plotted in a postscript file as a function of temperature.