zadorlab / sella

A Python software package for saddle point optimization and minimization of atomic systems.
https://www.ecc-project.org/
Other
72 stars 20 forks source link

Minimal example of N2 with EMT or LJ leads to a crash #45

Open Andrew-S-Rosen opened 4 months ago

Andrew-S-Rosen commented 4 months ago

The following simple example does not seem to work, and I'm not 100% sure why.

Reproduction

pip install git+https://github.com/zadorlab/sella.git
from sella import Sella
from ase.build import molecule
from ase.calculators.emt import EMT
atoms = molecule("N2")

atoms.calc = EMT()
opt = Sella(atoms)
opt.run()

The same error happens with LennardJones as the calculator instead of EMT.  

Traceback

ValueError                               Traceback (most recent call last)Cell In[1], line 8
      6 atoms.calc = EMT()
      7 opt = Sella(atoms)
----> 8 opt.run()

File ~\miniconda\envs\test\Lib\site-packages\ase\optimize\optimize.py:430, in Optimizer.run(self, fmax, steps)
    415 """Run optimizer.
    416
    417 Parameters
   (...)
    427     True if the forces on atoms are converged.
    428 """
    429 self.fmax = fmax
--> 430 return Dynamics.run(self, steps=steps)

File ~\miniconda\envs\test\Lib\site-packages\ase\optimize\optimize.py:275, in Dynamics.run(self, steps)
    257 def run(self, steps=DEFAULT_MAX_STEPS):
    258     """Run dynamics algorithm.
    259
    260     This method will return when the forces on all individual
   (...)
    272         True if the forces on atoms are converged.
    273     """
--> 275     for converged in Dynamics.irun(self, steps=steps):
    276         pass
    277     return converged

File ~\miniconda\envs\test\Lib\site-packages\ase\optimize\optimize.py:246, in Dynamics.irun(self, steps)
    243 # run the algorithm until converged or max_steps reached
    244 while not is_converged and self.nsteps < self.max_steps:
    245     # compute the next step
--> 246     self.step()
    247     self.nsteps += 1
    249     # log the step

File ~\github\sella\sella\optimize\optimize.py:248, in Sella.step(self)
    247 def step(self):
--> 248     s, smag = self._predict_step()
    250     # Determine if we need to call the eigensolver, then step
    251     if self.nsteps_since_diag >= self.diag_every_n:

File ~\github\sella\sella\optimize\optimize.py:227, in Sella._predict_step(self)
    225         self.pes.calculate_hessian()
    226     else:
--> 227         self.pes.diag(**self.diagkwargs)
    228     self.nsteps_since_diag = -1
    229 self.initialized = True

File ~\github\sella\sella\peswrapper.py:269, in PES.diag(self, gamma, threepoint, maxiter)
    266 Hproj = NumericalHessian(self._calc_eg, self.get_x(), self.get_g(),
    267                          self.eta, threepoint, Ufree)
    268 Hc = self.get_Hc()
--> 269 rayleigh_ritz(Hproj - Ufree.T @ Hc @ Ufree, gamma, P, v0=v0,
    270               method=self.eigensolver,
    271               maxiter=maxiter)
    273 # Extract eigensolver iterates
    274 Vs = Hproj.Vs

File ~\github\sella\sella\eigensolvers.py:58, in rayleigh_ritz(A, gamma, P, B, v0, vref, vreftol, method, maxiter)
     56 while True:
     57     Atilde = V.T @ (symmetrize_Y(V, AV, symm=symm))
---> 58     lams, vecs = eigh(Atilde, V.T @ B @ V)
     59     nneg = max(1, np.sum(lams < 0))
     60     # Rotate our subspace V to be diagonal in A.
     61     # This is not strictly necessary but it makes our lives easier later

File ~\miniconda\envs\test\Lib\site-packages\scipy\linalg\_decomp.py:458, in eigh(a, b, lower, eigvals_only, overwrite_a, overwrite_b, type, check_finite, subset_by_index, subset_by_value, driver)
    454 if driver not in drv_str:
    455     raise ValueError('"{}" is unknown. Possible values are "None", "{}".'
    456                      ''.format(driver, '", "'.join(drv_str[1:]))) --> 458 a1 = _asarray_validated(a, check_finite=check_finite)
    459 if len(a1.shape) != 2 or a1.shape[0] != a1.shape[1]:
    460     raise ValueError('expected square "a" matrix')

File ~\miniconda\envs\test\Lib\site-packages\scipy\_lib\_util.py:321, in _asarray_validated(a, check_finite, sparse_ok, objects_ok, mask_ok, as_inexact)
    319         raise ValueError('masked arrays are not supported')
    320 toarray = np.asarray_chkfinite if check_finite else np.asarray
--> 321 a = toarray(a)
    322 if not objects_ok:
    323     if a.dtype is np.dtype('O'):

File ~\miniconda\envs\test\Lib\site-packages\numpy\lib\_function_base_impl.py:649, in asarray_chkfinite(a, dtype, order)
    647 a = asarray(a, dtype=dtype, order=order)
    648 if a.dtype.char in typecodes['AllFloat'] and not np.isfinite(a).all():
--> 649     raise ValueError(
    650         "array must not contain infs or NaNs")
    651 return a

ValueError: array must not contain infs or NaNs