Merge pull request #174 from SSCHAcode/asr_1d

ASR in 1D materials

Author: mesonepigreco

Modules/Ensemble.py

@@ -1733,7 +1733,8 @@ def update_weights_fourier(self, new_dynamical_matrix, newT, timer=None):
 
         # Exclude translations
         if not self.ignore_small_w:
-            trans_original = CC.Methods.get_translations(pols_original, super_struct0.get_masses_array())
+            trans_original = super_struct0.get_asr_modes(pols_original)
+            # trans_original = CC.Methods.get_translations(pols_original, super_struct0.get_masses_array())
         else:
             trans_original = np.abs(w_original) < CC.Phonons.__EPSILON_W__
 
@@ -1750,15 +1751,14 @@ def update_weights_fourier(self, new_dynamical_matrix, newT, timer=None):
 
 
         if not self.ignore_small_w:
-            trans_mask = CC.Methods.get_translations(pols, super_structure.get_masses_array())
+            trans_mask = super_structure.get_asr_modes(pols)
+            # trans_mask = CC.Methods.get_translations(pols, super_structure.get_masses_array())
         else:
             trans_mask = np.abs(w_new) < CC.Phonons.__EPSILON_W__
 
 
         # Check if the new dynamical matrix satisfies the sum rule
-        violating_sum_rule = (np.sum(trans_mask.astype(int)) != 3) or (np.sum(trans_original.astype(int)) != 3)
-        if self.ignore_small_w:
-            violating_sum_rule = np.sum(trans_mask.astype(int)) != np.sum(trans_original.astype(int))
+        violating_sum_rule = np.sum(trans_mask.astype(int)) != np.sum(trans_original.astype(int))
 
 
         if violating_sum_rule:
@@ -1955,7 +1955,8 @@ def update_weights(self, new_dynamical_matrix, newT, update_q = False, timer=Non
 
         # Exclude translations
         if not self.ignore_small_w:
-            trans_original = CC.Methods.get_translations(pols_original, super_struct0.get_masses_array())
+            trans_original = super_struct0.get_asr_modes(pols_original)
+            # trans_original = CC.Methods.get_translations(pols_original, super_struct0.get_masses_array())
         else:
             trans_original = np.abs(w_original) < CC.Phonons.__EPSILON_W__
 
@@ -1972,15 +1973,15 @@ def update_weights(self, new_dynamical_matrix, newT, update_q = False, timer=Non
 
 
         if not self.ignore_small_w:
-            trans_mask = CC.Methods.get_translations(pols, super_structure.get_masses_array())
+            trans_mask = super_structure.get_asr_modes(pols)
+            # trans_mask = CC.Methods.get_translations(pols, super_structure.get_masses_array())
         else:
             trans_mask = np.abs(w_new) < CC.Phonons.__EPSILON_W__
 
 
         # Check if the new dynamical matrix satisfies the sum rule
-        violating_sum_rule = (np.sum(trans_mask.astype(int)) != 3) or (np.sum(trans_original.astype(int)) != 3)
-        if self.ignore_small_w:
-            violating_sum_rule = np.sum(trans_mask.astype(int)) != np.sum(trans_original.astype(int))
+        # violating_sum_rule = (np.sum(trans_mask.astype(int)) != 3) or (np.sum(trans_original.astype(int)) != 3)
+        violating_sum_rule = np.sum(trans_mask.astype(int)) != np.sum(trans_original.astype(int))
 
 
         if violating_sum_rule:
@@ -1995,7 +1996,7 @@ def update_weights(self, new_dynamical_matrix, newT, update_q = False, timer=Non
 DETAILS OF ERROR:
     Number of translatinal modes in the original dyn = {}
     Number of translational modes in the target dyn = {}
-    (They should be both 3)
+    (They should be both be the same - 3 in bulk)
 """.format(np.sum(trans_original.astype(int)), np.sum(trans_mask.astype(int)))
 
             print(ERR_MSG)
@@ -2718,7 +2719,8 @@ def get_preconditioned_gradient(self, subtract_sscha = True, return_error = Fals
             w, pols = self.current_dyn.DiagonalizeSupercell()#supercell_dyn.DyagDinQ(0)
 
         if not self.ignore_small_w:
-            trans = CC.Methods.get_translations(pols, super_struct.get_masses_array())
+            trans = super_struct.get_asr_modes(pols)
+            # trans = CC.Methods.get_translations(pols, super_struct.get_masses_array())
         else:
             trans = np.abs(w) < CC.Phonons.__EPSILON_W__
 
@@ -2958,7 +2960,8 @@ def get_stress_tensor(self, offset_stress = None, use_spglib = False):
         wr, pols = self.current_dyn.DiagonalizeSupercell()
 
         if not self.ignore_small_w:
-            trans = ~ CC.Methods.get_translations(pols, super_structure.get_masses_array())
+            trans = ~ super_structure.get_asr_modes(pols)
+            # trans = ~ CC.Methods.get_translations(pols, super_structure.get_masses_array())
         else:
             trans = np.abs(wr) > CC.Phonons.__EPSILON_W__
 
@@ -3386,7 +3389,8 @@ def get_odd_realspace(self):
         w_sc, pols_sc = super_dyn.DyagDinQ(0)
 
         # Remove translations
-        no_trans_mask = ~CC.Methods.get_translations(pols_sc, super_dyn.structure.get_masses_array())
+        no_trans_mask = ~super_dyn.structure.get_asr_modes(pols_sc)
+        # no_trans_mask = ~CC.Methods.get_translations(pols_sc, super_dyn.structure.get_masses_array())
         w_sc = w_sc[no_trans_mask]
         pols_sc = pols_sc[:, no_trans_mask]
 
@@ -3681,7 +3685,8 @@ def get_free_energy_hessian(self, include_v4 = False, get_full_hessian = True, v
 
         # Get the translational modes
         if not self.ignore_small_w:
-            trans = CC.Methods.get_translations(pols, super_structure.get_masses_array())
+            trans = super_structure.get_asr_modes(pols)
+            # trans = CC.Methods.get_translations(pols, super_structure.get_masses_array())
         else:
             trans = np.abs(w) < CC.Phonons.__EPSILON_W__
 

Modules/SchaMinimizer.py

@@ -1504,10 +1504,11 @@ def finalize(self, verbose = 1, timer=None):
             pols = self.ensemble.current_pols.copy()
 
             #w, pols = super_dyn.DyagDinQ(0)
-            trans = CC.Methods.get_translations(pols, super_struct.get_masses_array())
+            trans = super_struct.get_asr_modes(pols)
+            # trans = CC.Methods.get_translations(pols, super_struct.get_masses_array())
 
             for i in range(len(w)):
-                print ("Mode %5d:   freq %16.8f cm-1  | is translation? " % (i+1, w[i] * __RyToCm__), trans[i])
+                print ("Mode %5d:   freq %16.8f cm-1  | is asr? " % (i+1, w[i] * __RyToCm__), trans[i])
 
             print ()
 
@@ -1535,7 +1536,8 @@ def check_imaginary_frequencies(self, timer=None):
 
         # Get translations
         if not self.ensemble.ignore_small_w:
-            trans_mask = ~CC.Methods.get_translations(pols, ss.get_masses_array())
+            trans_mask = ~ss.get_asr_modes(pols)
+            # trans_mask = ~CC.Methods.get_translations(pols, ss.get_masses_array())
         else:
             trans_mask = np.abs(w) > CC.Phonons.__EPSILON_W__
 
@@ -1556,7 +1558,8 @@ def check_imaginary_frequencies(self, timer=None):
             #ss0 = self.ensemble.dyn_0.structure.generate_supercell(self.dyn.GetSupercell())
 
             if not self.ensemble.ignore_small_w:
-                trans_mask = ~CC.Methods.get_translations(pold, ss.get_masses_array())
+                trans_mask = ~ss.get_asr_modes(pold)
+                # trans_mask = ~CC.Methods.get_translations(pold, ss.get_masses_array())
             else:
                 trans_mask = np.abs(wold) > CC.Phonons.__EPSILON_W__
 
@@ -1945,7 +1948,8 @@ def ApplyLambdaTensor(current_dyn, matrix, T = 0):
     w, pols = current_dyn.DyagDinQ(0)
 
     # Get the translations
-    trans = ~CC.Methods.get_translations(pols, current_dyn.structure.get_masses_array())
+    trans = ~current_dyn.structure.get_asr_modes(pols)
+    # trans = ~CC.Methods.get_translations(pols, current_dyn.structure.get_masses_array())
 
     # Restrict only to non translational modes
     w = np.real(w[trans])
@@ -2009,7 +2013,8 @@ def ApplyFCPrecond(current_dyn, matrix, T = 0):
     w, pols = current_dyn.DyagDinQ(0)
 
     # Get the translations
-    trans = ~CC.Methods.get_translations(pols, current_dyn.structure.get_masses_array())
+    trans = ~current_dyn.structure.get_asr_modes(pols)
+    #trans = ~CC.Methods.get_translations(pols, current_dyn.structure.get_masses_array())
 
     # Restrict only to non translational modes
     w = np.real(w[trans])
@@ -2100,7 +2105,8 @@ def GetStructPrecond(current_dyn, ignore_small_w = False, w_pols = None):
 
     # Select translations
     if not ignore_small_w:
-        not_trans = ~CC.Methods.get_translations(pols, mass)
+        not_trans = ~current_dyn.structure.get_asr_modes(pols)
+        # not_trans = ~CC.Methods.get_translations(pols, mass)
     else:
         not_trans = np.abs(w) > CC.Phonons.__EPSILON_W__