On an issue regarding the convergence of SSCHA calculations

[hellolori]

Dear developers,

I recently conducted a study on the anharmonic effects of a charge-density-wave material, which exhibits the lowest phonon frequency at the M point of a hexagonal lattice, utilizing the SSCHA package. The harmonic phonon dispersion of the material is illustrated below:

I aimed to minimize the free energy of this material at 0 K, employing a 2x2 supercell that contains 28 atoms. However, I encountered a convergence issue, where the structure gradient abruptly increased to an exceedingly large value, as depicted below.

I found a relevant topic in the FAQs section titled "The gradients on my simulations are increasing a lot, why is this happening?", suggesting that the ensemble might have fewer configurations than necessary. Following this advice, I increased the number of configurations per population from 1000 to 5000. Unfortunately, this adjustment did not resolve the issue. The figures below display the results using 2000 and 5000 configurations per population.

Number of configuration: 2000

Number of configuration: 5000

Moreover, I plotted the auxiliary phonon dispersion using the last population of the minimization process (without interpolation with the harmonic dynamical matrix). It is evident that there is a sharp dip in the phonon at the M point, though its frequency remains positive.

When the calculated SSCHA dynamical matrix (using a 2x2 supercell) is interpolated with the harmonic one (on a 6x6 q-grid), imaginary phonons appear, as shown in the subsequent figures.

Could you please provide some advice on how to achieve convergence in the SSCHA calculation?

[mesonepigreco]

Hi, it probably means that your system is unstable at 0K: quantum fluctuations are not enough to melt the CDW (you can see from the frequencies plot that the M frequency goes really low during the minimization).

I suggest using a higher temperature (e.g., 50 K) to stabilize the calculation. This should not change the results much, as except for that frequency, all the others have frequencies with an energy above 50 K, so this temperature should not change much the physics with respect to the result at 0K.
If you want the result at 0 K, you can decrease the temperature progressively after getting the result converged with 50 K

[hellolori]

Dear mesonepigreco,

Thank you very much for your suggestion. Acturally, I have tried to perform the calculation under higher temperatures before, such as 100 and 200 kelvin. However, it still cannot reach convergence.

Here are the results at T=200K, with the number of configurations per population set to 1000. The sudden increase in the structure gradient still occurs.

Here are the auxiliary phonon frequencies at 0 and 200 K, obtained using the dynamical matrix from the last population, using a 2x2 supercell, without interpolation.

And the SSCHA phonon frequencies at 0 and 200 K, respectively, interpolated with the harmonic ones (on a 6x6 q-grid) are shown below.

Though at T=200 kelvin, all the auxiliary phonon frequencies are positive, the SSCHA minimization still cannot converge.

I am very confused. Should I perform the calculation at a higher temperature?

[mesonepigreco]

Very wired: it seems even 200 K is not enough to stabilize the structure. Are you using a force field or DFT for the total energy and forces? What happens if you try to get the CDW structure? Move the atoms in the 2x2 supercell along the CDW instability and relax the structure. You can try to do the SSCHA of the CDW structure to see the phase-transition

[hellolori]

I am using a machine learning potential trained on 1000 configurations generated by SSCHA at corresponding temperatures. Actually, I can get the CDW structure by adding the random displacements on the atoms, followed by the subsequent structure optimization. There is indeed atomic reconstruction in a 2x2 supercell. I have tried to increase the electronic smearing to stabilize the crystal. I found that a critical sigma value above 0.1 Ry makes the phonon frequency at the M points just positive (harmonic phonon calculation), indicating a potentially high CDW transition temperature. This indicates the CDW transition temperature may be very high. Given this, could the lack of convergence in the SSCHA calculations below 200 k I previously described be considered reasonable?

Thank you for your suggestion. However, I still do not understand how to "perform SSCHA on the CDW structure." Could you provide more detailed guidance?