Segmentation fault in applying ASR

[Yr-ren]

Dear developers,
Thanks for your effort for developing and releasing this program. Recently, I have analyzed the anharmonic effect in 2D materials monolayer NbTe2. Things are going well. After calculating 1400 structures, the convergence is successful！However, there are some problems in imposing the acoustic sum rule (ASR). Error displays Segmentation fault. After checking, the possible reason is that the array is out of bounds. Have you ever had a problem like this？The supercell is 4x4x1. What can I do about that?
Then I try to use the NbTe2.Hessian.dyn (I know the frequency coming from the SSCHA free energy Hessian is not a true physical quantity) instead of the harmonic matrix in calculation of phonon spectrum by QE. For the reason, in the tutorial, it writes "The matrices PbTe.Hessian.dyn#q are a generalization of the standard harmonic dynamical matrices that include quantum and thermal effects on a static level." The static level cann't be improved by convenient QE calculation in q2r.x and matdyn.x. I guess? So it's wrong? The spectrum is meaningless.
By the way, in the PbTe spectral functions tutorial, in the part for applying ASR, dyn = CC.Phonons.Phonons("PbTe.SSCHA.dyn",3) seems that 3 should be 8. And I don't get the reason *2.0 in the " d3 = np.load("d3_realspace_sym.npy")*2.0". It will be greatly appreciated if you could spare time to answer my questions.
Regards,
Ren

[mesonepigreco]

Hi Ren,
Can you provide a minimal input file that reproduces the Seg. Fault error you are facing? (Just the with python script and the dynamical matrices)
In this way, we can check whether there is an error on the input (so we can insert a more meaningful error message in the code) or it is a real bug of the code that we need to fix.

Reguarding the other questions:

dyn = CC.Phonons.Phonons("PbTe.SSCHA.dyn",3) seems that 3 should be 8.

The 3 indicates the number of irreducible q points by symmetry (i.e. the number of dynamical matrix files produced by quantum espresso). PbTe in the cubic Fm-3m structure has 48 symmetry operations, and the 8 q points of the 2x2x2 supercell reduce, once symmetries are applied, to 3 points. So, the quantum espresso dynamical matrices are 3 independent (loaded from the files).

I don't get the reason *2.0 in the " d3 = np.load("d3_realspace_sym.npy")*2.0

This is a unit conversion. I understand that it not so clear, however, the FORTRAN subroutines of the SSCHA code works in Hartree atomic units, while the rest of the code uses in Rydberg atomic units (this is a consequence of merging a part of an old FORTRAN implementation of SSCHA with some subroutines from quantum espresso and CellConstructor). The file d3_realspace_sym.npy is expressed in Hartree units and needs to be converted in Rydberg before further processing (1 Ha = 2 Ry).
We will eliminate this ugly behavior in a future release of the code.

I try to use the NbTe2.Hessian.dyn (I know the frequency coming from the SSCHA free energy Hessian is not a true physical quantity) instead of the harmonic matrix in calculation of phonon spectrum by QE. For the reason, in the tutorial, it writes "The matrices PbTe.Hessian.dyn#q are a generalization of the standard harmonic dynamical matrices that include quantum and thermal effects on a static level." The static level cann't be improved by convenient QE calculation in q2r.x and matdyn.x. I guess? So it's wrong? The spectrum is meaningless.

You tried to use espresso q2r.x and matdyn.x to plot the phonon dispersion from the NbTe2.Hessian.dyn file and you get a meaningless spectrum, correct?
In principle, the procedure you used should give the correct phonon spectrum once the SSCHA calculation is well converged. So probably there is a problem in the SSCHA calculation and the original dynamical matrices have some errors (which may also explain the Segmentation Fault you get when applying the acoustic sum rule).
If you can attach your dynamical matrices to the discussion (or send us privately if you prefer to keep them confidential) we may help you further with this issue.

Bests,
Lorenzo

[Yr-ren]

Dear Lorenzo,
Thanks for your reply. Through your explanation, I have a better understanding of the whole process. This is my python script and the dynamical matrices.
for.py.txt
dyn_end_population1_1.txt
dyn_end_population1_2.txt
dyn_end_population1_3.txt
dyn_end_population1_4.txt
It also needs d3_realspace_sym.npy to load, I sent it to your email.
In fact, through the QE calculation, I get a reasonable result. At 300K, the structure is temperature dependent. This is the result.

But I'm not sure if it has a physical meaning, because I don't fully understand the difference between dynamic and static calculation. Is it reasonable to calculate with QE? The result is phonon spectrum? Or Hessian frequency? Or nothing.
In addition, in the tutorial, other loads are dyn = CC.Phonons.Phonons("PbTe.SSCHA.dyn",8) . Since I haven't calculated PbTe, I'm not sure how many irreducible points in BZ, but I think all loads have the same number?
Bests,
Ren