Skip to content

Commit

Permalink
Change the order of if-statements in critical algo
Browse files Browse the repository at this point in the history
First change is to make sure, if there is LinAlg error then
to discontinue.
If the density is zero, then discontinue
if the gradient isn't zero, then discontinue.
  • Loading branch information
Ali-Tehrani committed Aug 4, 2021
1 parent d6dbf0f commit 4ee8e95
Showing 1 changed file with 16 additions and 11 deletions.
27 changes: 16 additions & 11 deletions chemtools/topology/critical.py
Original file line number Diff line number Diff line change
Expand Up @@ -142,19 +142,23 @@ def find_critical_points(self):
(self._coords is not None and index < len(self._coords)):
try:
coord = self._root_vector_func(point.copy())
# add critical point if it is new and it doesn't contain any
if not (
np.any([np.linalg.norm(coord - cp.coordinate) < 1.e-3 for cp in self.cps])
and not np.any(np.isnan(coord))
):
dens = self.func(coord)
# skip critical point if its density value is zero
if np.abs(dens) > 1.e-4:
grad = self.grad(coord)
# Make sure gradient is zero, as it's a critical point.
if np.all(np.abs(grad) < 1e-4):
# compute rank & signature of critical point
eigenvals, eigenvecs = np.linalg.eigh(self.hess(coord))
cp = CriticalPoint(coord, eigenvals, eigenvecs, 1e-4)
self._cps.setdefault((cp.rank[0], cp.signature[0]), []).append(cp)
except np.linalg.LinAlgError as _:
continue
# add critical point if it is new
if not np.any([np.linalg.norm(coord - cp.coordinate) < 1.e-3 for cp in self.cps]):
dens = self.func(coord)
grad = self.grad(coord)
# skip critical point if its dens & grad are zero
if abs(dens) < 1.e-4 and np.all(abs(grad) < 1.e-4):
continue
# compute rank & signature of critical point
eigenvals, eigenvecs = np.linalg.eigh(self.hess(coord))
cp = CriticalPoint(coord, eigenvals, eigenvecs, 1e-4)
self._cps.setdefault((cp.rank[0], cp.signature[0]), []).append(cp)
# check Poincare–Hopf equation
if not self.poincare_hopf_equation:
warnings.warn("Poincare–Hopf equation is not satisfied.", RuntimeWarning)
Expand All @@ -178,6 +182,7 @@ def _root_vector_func(self, guess, maxiter=5000):
guess = guess - np.dot(np.linalg.inv(hess), grad[:, np.newaxis]).flatten()
norm = np.linalg.norm(grad, axis=-1)
niter += 1
# Sometimes guess returns [np.nan, np.nan, np.nan]
return guess

@staticmethod
Expand Down

0 comments on commit 4ee8e95

Please sign in to comment.