Second-order structural phase transitions, free energy curvature, and temperature-dependent anharmonic phonons in the self-consistent harmonic approximation: theory and stochastic implementation

The self-consistent harmonic approximation is an effective harmonic theory to calculate the free energy of systems with strongly anharmonic atomic vibrations, and its stochastic implementation has proved to be an efficient method to study, from first-principles, the anharmonic properties of solids. The free energy as a function of average atomic positions (centroids) can be used to study quantum or thermal lattice instability. In particular the centroids are order parameters in second-order structural phase transitions such as, e.g., charge-density-waves or ferroelectric instabilities. According to Landau's theory, the knowledge of the second derivative of the free energy (i.e., the curvature) with respect to the centroids in a high-symmetry configuration allows the identification of the phase-transition and of the instability modes. In this work we derive the exact analytic formula for the second derivative of the free energy in the self-consistent harmonic approximation for a generic atomic configuration. The analytic derivative is expressed in terms of the atomic displacements and forces in a form that can be evaluated by a stochastic technique using importance sampling. Our approach is particularly suitable for applications based on first-principles density-functional-theory calculations, where the forces on atoms can be obtained with a negligible computational effort compared to total energy determination. Finally, we propose a dynamical extension of the theory to calculate spectral properties of strongly anharmonic phonons, as probed by inelastic scattering processes. We illustrate our method with a numerical application on a toy model that mimics the ferroelectric transition in rock-salt crystals such as SnTe or GeTe.