SSCHAcode / python-sscha

The python implementation of the Stochastic Self-Consistent Harmonic Approximation (SSCHA).
GNU General Public License v3.0
55 stars 21 forks source link

vasp-phonopy-sscha problems #164

Open caolz opened 10 months ago

caolz commented 10 months ago

First I follow the README, after this step “Which is your PHONOPY supercell? 3 ### On a 3x3x3 we could get a good harmonic guess”,The first error occurred : Traceback (most recent call last): File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/interface.py", line 298, in run1.matrix_formatting() File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/interface.py", line 125, in matrix_formatting self.dynamical_matrix = self.dynamical_matrix.drop(0, 1) TypeError: DataFrame.drop() takes from 1 to 2 positional arguments but 3 were given

After that I modified lines 125 to 127 of interface.py self.dynamical_matrix = self.dynamical_matrix.drop([0, 1], axis=1)

I was able to continue using it, but after that “You want to use the Sobol sequence for the configurations? y/n” this issue came up.

Traceback (most recent call last): File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/interface.py", line 314, in n_random,n_population,NQIRR = generate(None) File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/generate.py", line 92, in generate ens = sscha.Ensemble.Ensemble(dyn, T, dyn.GetSupercell()) # (((dyn->new_dyn) File "/home/ustb-clz/software/intel/oneapi/intelpython/latest/envs/sscha/lib/python3.10/site-packages/sscha/Ensemble.py", line 152, in init self.dyn_0 = dyn0.Copy() File "/home/ustb-clz/software/intel/oneapi/intelpython/latest/envs/sscha/lib/python3.10/site-packages/sscha/Ensemble.py", line 230, in setattr self.w_0, self.pols_0 = value.DiagonalizeSupercell() File "/home/ustb-clz/software/intel/oneapi/intelpython/latest/envs/sscha/lib/python3.10/site-packages/cellconstructor/Phonons.py", line 3635, in DiagonalizeSupercell c_e_sc = tilde_e_qnu[itau_modes] np.exp(1jphase) / np.sqrt(supercell_size) ValueError: operands could not be broadcast together with shapes (810,) (150,)

I really don't know how to fix this anymore.

diegomartinez2 commented 10 months ago

I do not know if you put something else in the “Which is your PHONOPY supercell?". It should be only one number that is taken to create the NxNxN system. The later error comes from the incompatible shape. I think all revolves around the supercell size input.

caolz commented 10 months ago

Yes, I found that the problem was that the number of dyn/dynq was filled in incorrectly.But now I have a new problem.

Traceback (most recent call last):
  File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/interface.py", line 312, in <module>
    n_random,n_population,NQIRR = generate(None)
  File "/home/ustb-clz/software/vasp-phonopy-sscha-main/vasp-phonopy-sscha/generate.py", line 95, in generate
    ens.generate(n_random, evenodd=True, sobol = Sobolflag)
  File "/home/ustb-clz/software/intel/oneapi/intelpython/latest/envs/sscha/lib/python3.10/site-packages/sscha/Ensemble.py", line 1125, in generate
    structs = self.dyn_0.ExtractRandomStructures(N // 2, self.T0, project_on_vectors = project_on_modes, lock_low_w = self.ignore_small_w, sobol = sobol, sobol_scramble = sobol_scramble, sobol_scatter = sobol_scatter)  # normal Sobol generator****Diegom_test****
  File "/home/ustb-clz/software/intel/oneapi/intelpython/latest/envs/sscha/lib/python3.10/site-packages/cellconstructor/Phonons.py", line 2014, in ExtractRandomStructures
    raise ValueError(ERR_MSG)
ValueError:
    Error, the current matrix is not positive definite.
           I cannot extract a random ensamble.
           If you want to skip this error,
           consider calling the method ForcePositiveDefinite() before extracting the ensemble.

        It could also be a consequence of a sum rule not well imposed.
        Try to run Symmetrize() to force the sum rule.

I do have a structure that has a imaginary frequency, but it's not possible to continue the operation?

diegomartinez2 commented 10 months ago

That's typical in some harmonic calculations to have imaginary frequencies that's why we take higher order approximations. To solve this problem you have to use ForcePositiveDefinite() in your code; that's done in the "You want to positivize your dynamical matrix? y/n" part.

caolz commented 10 months ago

"You want to positivize your dynamical matrix? y/n". This option I get the same error when I choose y or n. Do I need to go and change the source code?

diegomartinez2 commented 10 months ago

If ForcePositiveDefinite() doesn't work I'm not sure what is happening here. In line 84 of vasp-phonopy-sscha/generate.py there is the Symmetrize() call to force the sum rule. Does this happen only on the Sobol sequence or also happend if you do not use Sobol?

caolz commented 10 months ago

I did 4 tests whether I chose y or n and got the same error.

diegomartinez2 commented 10 months ago

So is not a Sobol error and we have the dynamical matrix positive definite and with the sum rule. Are you sure you are using the right dynamical matrix?. I'm also a bit lost here.

caolz commented 10 months ago
import cellconstructor as CC
import cellconstructor.Phonons
import sscha, sscha.Ensemble
import os
DATA_DIR = "data"
DYN_PREFIX = "pop1/dyn/dynq"
n_population = 3
NQIRR = 14
n_random = 300
Sobolflag = True
T = 0
dyn = CC.Phonons.Phonons(DYN_PREFIX, NQIRR)
dyn.ForcePositiveDefinite()
dyn.Symmetrize()
dyn.save_qe("pop1/dyn/dynq")
ens = sscha.Ensemble.Ensemble(dyn, T, dyn.GetSupercell())
ens.generate(n_random, evenodd=True, sobol = Sobolflag)
ens.save(DATA_DIR, n_population)
os.system("mkdir data/position")
os.system("mv data/scf_pop* data/position")

I wrote my own code from generate.py and found that Sobol = Fales gives the same error as before, when Sobol = True it works. But the original code with the same variables would not run.

diegomartinez2 commented 10 months ago

Ok. I need to check what is happening here. At least you can use your code to make calculations? With that code I think you can calculate the first population if you set n_population=1.

caolz commented 10 months ago

Yes, at least the first populations can be calculated.

caolz commented 10 months ago

I found something interesting, the first time I ran the program normally, after the previous error was reported, after commenting out a few lines of interface.py, and then re-running it again, it worked fine.

line 291 #os.system("rm -r pop1")
line 293 #os.system("mkdir dyn pop1")
line 301 #os.system("mv dynq* dyn")
line 302 #os.system("mv dyn pop1")