Skip to content

Commit 83a5f43

Browse files
ManexAlkortaManexAlkorta
committed
Added test_free_energy_hessian_dev and remove dummy variables from ensemble.init()
1 parent 934150c commit 83a5f43

13 files changed

Lines changed: 1780 additions & 24 deletions

File tree

‎Modules/Ensemble.py‎

Lines changed: 22 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -252,19 +252,6 @@ def __init__(self, dyn0, T0, supercell = None, **kwargs):
252252
# A flag that memorize if the ensemble has also the stresses
253253
self.has_stress = True
254254

255-
# Store data for calculations of 3FC elements:
256-
self.a = np.zeros( (Nsc * 3), dtype = np.double)
257-
self.new_pol = np.zeros( (Nsc, Nsc * 3, 3), dtype = np.double)
258-
# If the element is sym:
259-
self.ur = np.zeros( (self.N, Nsc * 3))
260-
self.upsilon = np.zeros( (Nsc*3, Nsc * 3))
261-
# Symmetry data for sym 3FC:
262-
self.nsym = 0
263-
self.s_cart = np.zeros( (3, 3, 48) , dtype = np.float64, order = "F")
264-
self.s_inv_cart = np.zeros( (3, 3, 48) , dtype = np.float64, order = "F")
265-
self.irt = np.zeros( (48, Nsc), dtype = np.intc, order = "F")
266-
self.translations_irt = np.zeros( (Nsc, np.prod(self.supercell)), dtype = np.intc, order = "F")
267-
268255
# A flag for each configuration that check if it possess a force and a stress
269256
self.force_computed = None
270257
self.stress_computed = None
@@ -3911,6 +3898,17 @@ def get_free_energy_hessian_dev(self, include_v4 = False, do_scf = True, eps = 1
39113898
nat_sc = int(np.shape(pols)[0] / 3)
39123899
nat = self.dyn_0.structure.N_atoms
39133900

3901+
a = np.zeros( (nat_sc * 3), dtype = np.double)
3902+
new_pol = np.zeros( (nat_sc, nat_sc * 3, 3), dtype = np.double)
3903+
# If the element is sym:
3904+
ur = np.zeros( (self.N, nat_sc * 3))
3905+
upsilon = np.zeros( (nat_sc*3, nat_sc * 3))
3906+
# Symmetry data for sym 3FC:
3907+
s_cart = np.zeros( (3, 3, 48) , dtype = np.float64, order = "F")
3908+
s_inv_cart = np.zeros( (3, 3, 48) , dtype = np.float64, order = "F")
3909+
irt = np.zeros( (48, nat_sc), dtype = np.intc, order = "F")
3910+
translations_irt = np.zeros( (nat_sc, np.prod(self.supercell)), dtype = np.intc, order = "F")
3911+
39143912
# Get the translational modes
39153913
if not self.ignore_small_w:
39163914
trans = CC.Methods.get_translations(pols, super_structure.get_masses_array())
@@ -3933,25 +3931,25 @@ def get_free_energy_hessian_dev(self, include_v4 = False, do_scf = True, eps = 1
39333931
u = self.u_disps.reshape((self.N, nat_sc, 3), order = "C") #/ Bohr
39343932

39353933
log_err = "err_yesrho"
3936-
self.a = SCHAModules.thermodynamic.w_to_a(w, self.current_T)
3934+
a = SCHAModules.thermodynamic.w_to_a(w, self.current_T)
39373935
# Get the polarization vectors in the correct format
39383936
for i in range(nat_sc):
39393937
for j in range(n_modes):
3940-
self.new_pol[i, j, :] = pols[3*i : 3*(i+1), j]
3938+
new_pol[i, j, :] = pols[3*i : 3*(i+1), j]
39413939

39423940
#Calculating rotated displacements and upsilon matrix"
3943-
self.ur, self.upsilon = SCHAModules.get_ur_upsilon_matrices(self.a, self.new_pol, trans, amass, ityp, u)
3941+
ur, upsilon = SCHAModules.get_ur_upsilon_matrices(a, new_pol, trans, amass, ityp, u)
39443942
#Obtaining symmetry data
39453943
qe_sym = CC.symmetries.QE_Symmetry(super_structure)
39463944
qe_sym.SetupFromSPGLIB()
39473945

39483946
# SpaceGroup symmetries data
3949-
self.nsym, self.s_cart, self.s_inv_cart, self.irt = qe_sym.QE_nsym, qe_sym.QE_s_cart, qe_sym.QE_s_inv_cart, qe_sym.QE_irt
3950-
self.s_cart = np.asfortranarray(self.s_cart)
3951-
self.s_inv_cart = np.asfortranarray(self.s_inv_cart)
3947+
s_cart, s_inv_cart, irt = qe_sym.QE_s_cart, qe_sym.QE_s_inv_cart, qe_sym.QE_irt
3948+
s_cart = np.asfortranarray(s_cart)
3949+
s_inv_cart = np.asfortranarray(s_inv_cart)
39523950

39533951
#Translations
3954-
self.translations_irt = qe_sym.QE_translations_irt
3952+
translations_irt = qe_sym.QE_translations_irt
39553953
#self.prepare_sym_third_fc = False
39563954

39573955
wq = np.empty(nat*3, dtype=np.float64)
@@ -4048,7 +4046,7 @@ def get_free_energy_hessian_dev(self, include_v4 = False, do_scf = True, eps = 1
40484046
nat = self.current_dyn.structure.N_atoms
40494047
n_modes = self.current_dyn.structure.N_atoms*mod[0]*mod[1]*mod[2]*3
40504048

4051-
ref_3fc = SCHAModules.module_hess.get_ref3fc(nat, orbit3a, indep_3fc_elem, n_indep_3fc_elem, kernel_3fc, rot_3fc, self.ur, self.upsilon, f, self.rho, log_err, self.s_inv_cart, self.irt, self.translations_irt, True)
4049+
ref_3fc = SCHAModules.module_hess.get_ref3fc(nat, orbit3a, indep_3fc_elem, n_indep_3fc_elem, kernel_3fc, rot_3fc, ur, upsilon, f, self.rho, log_err, s_inv_cart, irt, translations_irt, True)
40524050

40534051
vs_red = np.empty([nrefq2,nat*3,nat*3,nat*3], dtype=np.complex128)
40544052
vs_red = SCHAModules.module_hess.get_ref_vsq(refq2,trs_l,rot_3fc,ref_3fc,mapping_triplet,True)
@@ -4060,15 +4058,15 @@ def get_free_energy_hessian_dev(self, include_v4 = False, do_scf = True, eps = 1
40604058

40614059
orbit4t, orbit4o, norbit_4, indep_4fc_elem, n_indep_4fc_elem, kernel_4fc, rot_4fc, mapping_quadruplet = Classify.recognize_quadruplet(self.current_dyn, mapping, map_uc, verbose)
40624060

4063-
ref_4fc = SCHAModules.module_hess.get_ref4fc(orbit4t, indep_4fc_elem, n_indep_4fc_elem, kernel_4fc, rot_4fc, self.ur, self.upsilon, f, self.rho, log_err, self.s_inv_cart, self.irt, self.translations_irt, True)
4061+
ref_4fc = SCHAModules.module_hess.get_ref4fc(orbit4t, indep_4fc_elem, n_indep_4fc_elem, kernel_4fc, rot_4fc, ur, upsilon, f, self.rho, log_err, s_inv_cart, irt, translations_irt, True)
40644062

40654063
ws_red = np.zeros([nrefq4,nat*3,nat*3,nat*3,nat*3], dtype=np.complex128)
40664064
ws_red = SCHAModules.module_hess.get_ref_wsq(refq4,trs_l,rot_4fc,ref_4fc,mapping_quadruplet,True)
40674065

40684066
degs = qClassify.find_degeneracies(trs_wq)
40694067
Pmn = qClassify.construct_Pmn(mapping, orbitq1a, orbitq1s, trs_polvecs, rot_cart)
40704068

4071-
v_red, ref_3fc = SCHAModules.module_hess.get_v3_red(nat, norbit3, orbit3a, orbit3s, indep_3fc_elem, n_indep_3fc_elem, kernel_3fc, rot_3fc, self.ur, self.upsilon, f, self.rho, log_err, self.s_inv_cart, self.irt, self.translations_irt)
4069+
v_red, ref_3fc = SCHAModules.module_hess.get_v3_red(nat, norbit3, orbit3a, orbit3s, indep_3fc_elem, n_indep_3fc_elem, kernel_3fc, rot_3fc, ur, upsilon, f, self.rho, log_err, s_inv_cart, irt, translations_irt)
40724070
vs = SCHAModules.module_hess.get_all_vsq(trs_l, v_red, map_uc)
40734071

40744072
if do_scf:
@@ -4089,7 +4087,7 @@ def get_free_energy_hessian_dev(self, include_v4 = False, do_scf = True, eps = 1
40894087
for r in range(nr):
40904088
# Translated Single atomic-cartesian index
40914089
# Fortran to Py: -1
4092-
phi_sc_odd[3*(self.translations_irt[nat1,r]-1)+alpha, 3*(self.translations_irt[nat2,r]-1)+beta] = fc9[index]
4090+
phi_sc_odd[3*(translations_irt[nat1,r]-1)+alpha, 3*(translations_irt[nat2,r]-1)+beta] = fc9[index]
40934091
dynq_odd = CC.Phonons.GetDynQFromFCSupercell(phi_sc_odd, np.array(self.current_dyn.q_tot),
40944092
self.current_dyn.structure, super_structure)
40954093
self.convert_units(UNITS_DEFAULT)

‎tests/test_free_energy_hessian_dev/hessianv4_1‎

Lines changed: 149 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,149 @@
1+
Dynamical matrix file
2+
File generated with the CellConstructor by Lorenzo Monacelli
3+
2 4 0 1.8897259890000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
4+
Basis vectors
5+
3.0799999999999996 0.0000000000000000 0.0000000000000000
6+
-1.5399999999999998 2.6673582400000000 0.0000000000000000
7+
0.0000000000000000 0.0000000000000000 5.0499999999999998
8+
1 'Si ' 25597.9124185021937592
9+
2 'C ' 10947.3571678344233078
10+
1 1 -0.0000000000000002 0.0000000000000001 -0.0040643610314026
11+
2 1 -0.0000000000000001 1.7782388266666669 2.5209356389685977
12+
3 2 -0.0000000000000002 0.0000000000000001 1.9032535866761318
13+
4 2 -0.0000000000000001 1.7782388266666669 4.4282535866761314
14+
15+
Dynamical Matrix in cartesian axes
16+
17+
q = ( 0.000000000000 0.000000000000 0.000000000000 )
18+
19+
1 1
20+
0.2224249498463850 0.0000000000000000 -0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
21+
-0.0000000000000000 0.0000000000000000 0.2224249504252410 0.0000000000000000 0.0000000000059596 0.0000000000000000
22+
-0.0000000000000000 0.0000000000000000 0.0000000000059596 0.0000000000000000 0.4271177583537886 0.0000000000000000
23+
1 2
24+
0.0102391881583775 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
25+
0.0000000000000000 0.0000000000000000 0.0102391881730570 0.0000000000000000 0.0000000000005833 0.0000000000000000
26+
0.0000000000000000 0.0000000000000000 -0.0000000000005833 0.0000000000000000 -0.2191138738805857 0.0000000000000000
27+
1 3
28+
-0.0606756614223764 0.0000000000000000 0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
29+
0.0000000000000000 0.0000000000000000 -0.0606756615756377 0.0000000000000000 0.0000000000013264 0.0000000000000000
30+
0.0000000000000000 0.0000000000000000 0.0000000000012314 0.0000000000000000 -0.0514604669025081 0.0000000000000000
31+
1 4
32+
-0.1719884765823874 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
33+
0.0000000000000000 0.0000000000000000 -0.1719884770226615 0.0000000000000000 -0.0000000000078693 0.0000000000000000
34+
0.0000000000000000 0.0000000000000000 -0.0000000000066077 0.0000000000000000 -0.1565434175706992 0.0000000000000000
35+
2 1
36+
0.0102391881583775 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
37+
-0.0000000000000000 0.0000000000000000 0.0102391881730570 0.0000000000000000 -0.0000000000005833 0.0000000000000000
38+
0.0000000000000000 0.0000000000000000 0.0000000000005833 0.0000000000000000 -0.2191138738805857 0.0000000000000000
39+
2 2
40+
0.2224249498463850 0.0000000000000000 -0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
41+
-0.0000000000000000 0.0000000000000000 0.2224249504252410 0.0000000000000000 -0.0000000000059596 0.0000000000000000
42+
-0.0000000000000000 0.0000000000000000 -0.0000000000059596 0.0000000000000000 0.4271177583537886 0.0000000000000000
43+
2 3
44+
-0.1719884765823874 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
45+
0.0000000000000000 0.0000000000000000 -0.1719884770226615 0.0000000000000000 0.0000000000078693 0.0000000000000000
46+
0.0000000000000000 0.0000000000000000 0.0000000000066077 0.0000000000000000 -0.1565434175706992 0.0000000000000000
47+
2 4
48+
-0.0606756614223764 0.0000000000000000 -0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
49+
-0.0000000000000000 0.0000000000000000 -0.0606756615756377 0.0000000000000000 -0.0000000000013264 0.0000000000000000
50+
-0.0000000000000000 0.0000000000000000 -0.0000000000012314 0.0000000000000000 -0.0514604669025081 0.0000000000000000
51+
3 1
52+
-0.0606756614223763 0.0000000000000000 -0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
53+
-0.0000000000000000 0.0000000000000000 -0.0606756615756377 0.0000000000000000 0.0000000000012314 0.0000000000000000
54+
-0.0000000000000000 0.0000000000000000 0.0000000000013264 0.0000000000000000 -0.0514604669025081 0.0000000000000000
55+
3 2
56+
-0.1719884765823874 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
57+
0.0000000000000000 0.0000000000000000 -0.1719884770226615 0.0000000000000000 0.0000000000066077 0.0000000000000000
58+
0.0000000000000000 0.0000000000000000 0.0000000000078693 0.0000000000000000 -0.1565434175706992 0.0000000000000000
59+
3 3
60+
0.2508160735417163 0.0000000000000000 0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
61+
0.0000000000000000 0.0000000000000000 0.2508160741961777 0.0000000000000000 -0.0000000000085173 0.0000000000000000
62+
-0.0000000000000000 0.0000000000000000 -0.0000000000085173 0.0000000000000000 0.3726347100885012 0.0000000000000000
63+
3 4
64+
-0.0181519355369513 0.0000000000000000 0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
65+
0.0000000000000000 0.0000000000000000 -0.0181519355978773 0.0000000000000000 0.0000000000006783 0.0000000000000000
66+
0.0000000000000000 0.0000000000000000 -0.0000000000006783 0.0000000000000000 -0.1646308256152895 0.0000000000000000
67+
4 1
68+
-0.1719884765823874 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
69+
0.0000000000000000 0.0000000000000000 -0.1719884770226615 0.0000000000000000 -0.0000000000066077 0.0000000000000000
70+
0.0000000000000000 0.0000000000000000 -0.0000000000078693 0.0000000000000000 -0.1565434175706992 0.0000000000000000
71+
4 2
72+
-0.0606756614223764 0.0000000000000000 -0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
73+
-0.0000000000000000 0.0000000000000000 -0.0606756615756377 0.0000000000000000 -0.0000000000012314 0.0000000000000000
74+
0.0000000000000000 0.0000000000000000 -0.0000000000013264 0.0000000000000000 -0.0514604669025081 0.0000000000000000
75+
4 3
76+
-0.0181519355369513 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000
77+
0.0000000000000000 0.0000000000000000 -0.0181519355978773 0.0000000000000000 -0.0000000000006783 0.0000000000000000
78+
-0.0000000000000000 0.0000000000000000 0.0000000000006783 0.0000000000000000 -0.1646308256152895 0.0000000000000000
79+
4 4
80+
0.2508160735417163 0.0000000000000000 -0.0000000000000000 0.0000000000000000 -0.0000000000000000 0.0000000000000000
81+
-0.0000000000000000 0.0000000000000000 0.2508160741961777 0.0000000000000000 0.0000000000085173 0.0000000000000000
82+
0.0000000000000000 0.0000000000000000 0.0000000000085173 0.0000000000000000 0.3726347100885012 0.0000000000000000
83+
84+
Diagonalizing the dynamical matrix
85+
86+
q = ( 0.000000000000 0.000000000000 0.000000000000 )
87+
88+
***************************************************************************
89+
freq ( 1) = 0.00000009 [THz] = 0.00000300 [cm-1]
90+
( 0.145757 0.000000 0.282883 0.000000 0.385658 0.000000 )
91+
( 0.145757 0.000000 0.282883 0.000000 0.385658 0.000000 )
92+
( 0.145757 0.000000 0.282883 0.000000 0.385658 0.000000 )
93+
( 0.145757 0.000000 0.282883 0.000000 0.385658 0.000000 )
94+
freq ( 2) = 0.00000018 [THz] = 0.00000609 [cm-1]
95+
( -0.453056 0.000000 -0.047546 0.000000 0.206106 0.000000 )
96+
( -0.453056 0.000000 -0.047546 0.000000 0.206106 0.000000 )
97+
( -0.453056 0.000000 -0.047546 0.000000 0.206106 0.000000 )
98+
( -0.453056 0.000000 -0.047546 0.000000 0.206106 0.000000 )
99+
freq ( 3) = 0.00000048 [THz] = 0.00001597 [cm-1]
100+
( -0.153281 0.000000 0.409532 0.000000 -0.242463 0.000000 )
101+
( -0.153281 0.000000 0.409532 0.000000 -0.242463 0.000000 )
102+
( -0.153281 0.000000 0.409532 0.000000 -0.242463 0.000000 )
103+
( -0.153281 0.000000 0.409532 0.000000 -0.242463 0.000000 )
104+
freq ( 4) = 8.00349639 [THz] = 266.96786497 [cm-1]
105+
( -0.620838 0.000000 0.000000 0.000000 -0.000000 0.000000 )
106+
( 0.620838 0.000000 -0.000000 0.000000 -0.000000 0.000000 )
107+
( 0.338468 0.000000 -0.000000 0.000000 -0.000000 0.000000 )
108+
( -0.338468 0.000000 0.000000 0.000000 -0.000000 0.000000 )
109+
freq ( 5) = 8.00349640 [THz] = 266.96786533 [cm-1]
110+
( -0.000000 0.000000 -0.620838 0.000000 0.000000 0.000000 )
111+
( 0.000000 0.000000 0.620838 0.000000 0.000000 0.000000 )
112+
( 0.000000 0.000000 0.338468 0.000000 -0.000000 0.000000 )
113+
( -0.000000 0.000000 -0.338468 0.000000 -0.000000 0.000000 )
114+
freq ( 6) = 16.01350909 [THz] = 534.15309058 [cm-1]
115+
( -0.000000 0.000000 -0.000000 0.000000 0.661398 0.000000 )
116+
( 0.000000 0.000000 -0.000000 0.000000 -0.661398 0.000000 )
117+
( -0.000000 0.000000 0.000000 0.000000 -0.250106 0.000000 )
118+
( 0.000000 0.000000 0.000000 0.000000 0.250106 0.000000 )
119+
freq ( 7) = 17.07550534 [THz] = 569.57746747 [cm-1]
120+
( 0.160559 0.000000 0.000000 0.000000 -0.000000 0.000000 )
121+
( -0.160559 0.000000 -0.000000 0.000000 -0.000000 0.000000 )
122+
( 0.688637 0.000000 0.000000 0.000000 0.000000 0.000000 )
123+
( -0.688637 0.000000 -0.000000 0.000000 0.000000 0.000000 )
124+
freq ( 8) = 17.07550537 [THz] = 569.57746823 [cm-1]
125+
( -0.000000 0.000000 0.160559 0.000000 -0.000000 0.000000 )
126+
( 0.000000 0.000000 -0.160559 0.000000 -0.000000 0.000000 )
127+
( -0.000000 0.000000 0.688637 0.000000 0.000000 0.000000 )
128+
( 0.000000 0.000000 -0.688637 0.000000 0.000000 0.000000 )
129+
freq ( 9) = 17.13443009 [THz] = 571.54298523 [cm-1]
130+
( -0.000000 0.000000 -0.000000 0.000000 -0.278046 0.000000 )
131+
( 0.000000 0.000000 0.000000 0.000000 -0.278046 0.000000 )
132+
( -0.000000 0.000000 -0.000000 0.000000 0.650147 0.000000 )
133+
( 0.000000 0.000000 0.000000 0.000000 0.650147 0.000000 )
134+
freq ( 10) = 18.12168868 [THz] = 604.47438219 [cm-1]
135+
( -0.278046 0.000000 -0.000000 0.000000 0.000000 0.000000 )
136+
( -0.278046 0.000000 -0.000000 0.000000 0.000000 0.000000 )
137+
( 0.650147 0.000000 0.000000 0.000000 -0.000000 0.000000 )
138+
( 0.650147 0.000000 0.000000 0.000000 -0.000000 0.000000 )
139+
freq ( 11) = 18.12168870 [THz] = 604.47438296 [cm-1]
140+
( -0.000000 0.000000 0.278046 0.000000 0.000000 0.000000 )
141+
( -0.000000 0.000000 0.278046 0.000000 -0.000000 0.000000 )
142+
( 0.000000 0.000000 -0.650147 0.000000 -0.000000 0.000000 )
143+
( 0.000000 0.000000 -0.650147 0.000000 0.000000 0.000000 )
144+
freq ( 12) = 23.40869849 [THz] = 780.83002148 [cm-1]
145+
( -0.000000 0.000000 0.000000 0.000000 0.112887 0.000000 )
146+
( 0.000000 0.000000 0.000000 0.000000 -0.112887 0.000000 )
147+
( 0.000000 0.000000 -0.000000 0.000000 0.698038 0.000000 )
148+
( 0.000000 0.000000 -0.000000 0.000000 -0.698038 0.000000 )
149+
***************************************************************************

0 commit comments

Comments
 (0)