Fix frozen atom for ZimmermanGroup#33

ms860309 · Jun 5, 2021 · 936c730 · 936c730
1 parent 4f59168
commit 936c730
Show file tree

Hide file tree

Showing 4 changed files with 9 additions and 10 deletions.
diff --git a/pygsm/growing_string_methods/main_gsm.py b/pygsm/growing_string_methods/main_gsm.py
@@ -512,7 +512,6 @@ def set_stage(self,totalgrad,sumgradrms, ts_cgradq,ts_gradrms,fp):
                 if self.optimizer[self.TSnode].options['DMAX']>0.1:
                     self.optimizer[self.TSnode].options['DMAX']=0.1
                 self.optimizer[self.TSnode] = eigenvector_follow(self.optimizer[self.TSnode].options.copy())
-                print(type(self.optimizer[self.TSnode]))
                 self.optimizer[self.TSnode].options['SCALEQN'] = 1.
                 self.nhessreset=10  # are these used??? TODO 
                 self.hessrcount=0   # are these used?!  TODO
@@ -620,7 +619,7 @@ def reparameterize(self,ic_reparam_steps=8,n0=0,nconstraints=1):
         Reparameterize the string
         '''
         if self.interp_method == 'DLC':
-            print('reparameterizing')
+            # print('reparameterizing')
             self.ic_reparam(nodes=self.nodes,energies=self.energies,climbing=(self.climb or self.find),ic_reparam_steps=ic_reparam_steps,NUM_CORE=self.mp_cores)
         return
 

diff --git a/pygsm/level_of_theories/base_lot.py b/pygsm/level_of_theories/base_lot.py
@@ -379,7 +379,7 @@ def get_gradient(self,coords,multiplicity,state,frozen_atoms=None):
         if Gradient.value is not None:
             if frozen_atoms is not None:
                 for a in frozen_atoms:
-                    Gradient.value[a,:]=0.
+                    Gradient.value[a, :]= 0.
             if Gradient.unit=="Hartree/Bohr":
                 return Gradient.value *units.ANGSTROM_TO_AU  #Ha/bohr*bohr/ang=Ha/ang
             elif Gradient.unit=="kcal/mol/Angstrom":

diff --git a/pygsm/optimizers/_linesearch.py b/pygsm/optimizers/_linesearch.py
@@ -4,15 +4,15 @@
 
 #TODO remove unecessary arguments: nconstraints, xp, ,...
 
-def NoLineSearch(n, x, fx, g, d, step, xp, constraint_step, parameters,molecule,verbose=False):
+def NoLineSearch(n, x, fx, g, d, step, xp, constraint_step, parameters, molecule, verbose=False):
 
     x = x + d * step  + constraint_step  # 
-    xyz = molecule.coord_obj.newCartesian(molecule.xyz, x-xp,verbose=verbose)
+    xyz = molecule.coord_obj.newCartesian(molecule.xyz, x-xp, frozen_atoms = molecule.frozen_atoms, verbose = verbose)
 
     # use these so molecule xyz doesn't change
     print(" evaluate fx in linesearch")
     fx = molecule.PES.get_energy(xyz)
-    gx = molecule.PES.get_gradient(xyz)
+    gx = molecule.PES.get_gradient(xyz, frozen_atoms = molecule.frozen_atoms)
     g = molecule.coord_obj.calcGrad(xyz,gx)
 
     #print(" [INFO]end line evaluate fx = %5.4f step = %1.2f." %(fx, step))
@@ -24,7 +24,7 @@ def NoLineSearch(n, x, fx, g, d, step, xp, constraint_step, parameters,molecule,
 # TODO might be wise to add to backtrack a condition that says if 
 # the number of iterations was many and the energy increased
 # just return the initial point
-def backtrack(nconstraints, x, fx, g, d, step, xp,constraint_step, parameters,molecule, verbose=False):
+def backtrack(nconstraints, x, fx, g, d, step, xp,constraint_step, parameters, molecule, verbose=False):
     print(" In backtrack")
 
     # n is the non-constrained
@@ -64,13 +64,13 @@ def backtrack(nconstraints, x, fx, g, d, step, xp,constraint_step, parameters,mo
         x = xp
         x = x + d * step  + constraint_step 
         xyzp = molecule.xyz.copy()
-        xyz = molecule.coord_obj.newCartesian(molecule.xyz, x-xp,verbose=verbose)
+        xyz = molecule.coord_obj.newCartesian(molecule.xyz, x-xp, frozen_atoms=molecule.frozen_atoms, verbose=verbose)
         # Evaluate the function and gradient values. 
         # use these so molecule xyz doesn't change
         fx = molecule.PES.get_energy(xyz)
 
         #print('new fx %11.9f' % fx)
-        gx = molecule.PES.get_gradient(xyz,frozen_atoms=molecule.frozen_atoms)
+        gx = molecule.PES.get_gradient(xyz, frozen_atoms=molecule.frozen_atoms)
         g = molecule.coord_obj.calcGrad(xyz,gx)
         #g = molecule.gradient
         width = 1.

diff --git a/pygsm/potential_energy_surfaces/pes.py b/pygsm/potential_energy_surfaces/pes.py
@@ -316,7 +316,7 @@ def get_gradient(self,xyz,frozen_atoms=None):
             for i in self.RESTRAINTS:
                 a= i[0]
                 force=i[1]
-                grada[a] += force*(xyz[a] - self.reference_xyz[a])
+                grad[a] += force*(xyz[a] - self.reference_xyz[a])
 
         grad = np.reshape(grad,(-1,1))
         return grad  #Ha/ang