From b5824ff87dc95a6c21ab760bc0ac1b415f9bb5fd Mon Sep 17 00:00:00 2001 From: davidcortesortuno Date: Wed, 31 Jul 2024 16:45:32 -0400 Subject: [PATCH] nebm fs: polishing integrator; removed base class dep --- fidimag/common/chain_method_integrators.py | 58 ++++++++++++---------- fidimag/common/nebm_FS.py | 12 +++-- 2 files changed, 42 insertions(+), 28 deletions(-) diff --git a/fidimag/common/chain_method_integrators.py b/fidimag/common/chain_method_integrators.py index 97bbdf1a..3cfc6137 100644 --- a/fidimag/common/chain_method_integrators.py +++ b/fidimag/common/chain_method_integrators.py @@ -42,6 +42,8 @@ def set_step(self, step): self._step = step_choices[step] +# TODO: these integrators are not general anymore, as they rely on the +# structure of the chain method classes class VerletIntegrator(BaseIntegrator): """A quick Verlet integration in Cartesian coordinates @@ -148,19 +150,24 @@ def _step(self, t, y, h, f): # NEBM integrator for F-S algorithm -class FSIntegrator(BaseIntegrator): +# TODO: these integrators are not general anymore, as they rely on the +# structure of the chain method classes +# TODO: type checks +class FSIntegrator(object): """A step integrator considering the action of the band """ - def __init__(self, band, forces, rhs_fun, action_fun, n_images, n_dofs_image, + def __init__(self, ChainObj, + # band, forces, distances, rhs_fun, action_fun, + # n_images, n_dofs_image, maxSteps=1000, maxCreep=5, actionTol=1e-10, forcesTol=1e-6, etaScale=1.0, dEta=2, minEta=0.001, # perturbSeed=42, perturbFactor=0.1, nTrail=10, resetMax=20 ): - super(FSIntegrator, self).__init__(band, rhs_fun) + # super(FSIntegrator, self).__init__(band, rhs_fun) - self.action_fun = action_fun + self.ChainObj = ChainObj # Integration parameters # TODO: move to run function? @@ -172,24 +179,22 @@ def __init__(self, band, forces, rhs_fun, action_fun, n_images, n_dofs_image, self.forcesTol = forcesTol self.actionTol = actionTol self.maxSteps = maxSteps - + self.step = 0 + self.nTrail = nTrail self.i_step = 0 - self.n_images = n_images - self.n_dofs_image = n_dofs_image + + # Chain objects: + self.n_images = self.ChainObj.n_images + self.n_dofs_image = self.ChainObj.n_dofs_image # self.forces_prev = np.zeros_like(band).reshape(n_images, -1) # self.G : - self.forces = forces - self.forces_old = np.zeros_like(forces) + self.forces = self.ChainObj.forces + self.distances = self.ChainObj.distances + self.forces_old = np.zeros_like(self.ChainObj.forces) - # Rename y to band - self.band = self.y - self.band_old = np.zeros_like(self.band) # y -> band - - self.step = 0 - self.nTrail = nTrail - - def run_until(self, t): - pass + # self.band should be just a reference to the band in the ChainObj + self.band = self.ChainObj.band + self.band_old = np.zeros_like(self.band) def run_for(self, n_steps): @@ -202,6 +207,7 @@ def run_for(self, n_steps): self.trailAction = np.zeros(self.nTrail) trailPool = cycle(range(self.nTrail)) # cycle through 0,1,...,(nTrail-1),0,1,... eta = 1.0 + self.i_step = 0 # In __init__: # self.band_last[:] = self.band @@ -216,7 +222,7 @@ def run_for(self, n_steps): self.band[:] = self.band_old # Compute from self.band. Do not update the step at this stage: - # This step updates the forces in the G array of the nebm module, + # This step updates the forces and distances in the G array of the nebm module, # using the current band state self.y # TODO: remove time from chain method rhs # make a specific function to update G?? @@ -243,12 +249,14 @@ def run_for(self, n_steps): self.trailAction - self.rhs(t, self.y) - self.action = self.action_fun() + self.ChainObj.nebm_step(self.band, ensure_zero_extrema=True) + self.action = self.ChainObj.action_fun() self.trailAction[nStart] = self.action nStart = next(trailPool) + self.i_step += 1 + # Getting averages of forces from the INNER images in the band (no extrema) # (forces are given by vector G in the chain method code) # TODO: we might use all band images, not only inner ones @@ -266,7 +274,7 @@ def run_for(self, n_steps): exitFlag = True break # creep loop - if (n_steps >= self.maxSteps): + if (self.i_step >= self.maxSteps): print('Number of steps reached maximum') exitFlag = True break # creep loop @@ -282,7 +290,7 @@ def run_for(self, n_steps): print('') resetCount += 1 bestAction = self.action_old - RefinePath(self.band) # Resets the path to equidistant structures (smoothing kinks?) + self.refine_path(self.ChainObj.distances, self.band) # Resets the path to equidistant structures (smoothing kinks?) # PathChanged[:] = True if resetCount > self.resetMax: @@ -309,13 +317,13 @@ def run_for(self, n_steps): resetCount = 0 # Taken from the string method class - def refine_path(self, distances): + def refine_path(self, distances, band): """ """ new_dist = np.linspace(distances[0], distances[-1], distances.shape[0]) # Restructure the string by interpolating every spin component # print(self.integrator.y[self.n_dofs_image:self.n_dofs_image + 10]) - bandrs = self.band.reshape(self.n_images, self.n_dofs_image) + bandrs = band.reshape(self.n_images, self.n_dofs_image) for i in range(self.n_dofs_image): cs = si.CubicSpline(distances, bandrs[:, i]) diff --git a/fidimag/common/nebm_FS.py b/fidimag/common/nebm_FS.py index 05453f7d..9dd4e648 100644 --- a/fidimag/common/nebm_FS.py +++ b/fidimag/common/nebm_FS.py @@ -324,7 +324,7 @@ def compute_min_action(self): minAction = np.sum(np.abs(self.energies[1:] - self.energies[:-1])) return 2 * (dE + minAction) - def nebm_step(self, y): + def nebm_step(self, y, ensure_zero_extrema=False): self.compute_effective_field_and_energy(y) nebm_clib.project_images(self.gradientE, y, @@ -342,6 +342,12 @@ def nebm_step(self, y): self.n_dofs_image ) + # The effective force at the extreme images should already be zero, but + # we will manually remove any value + if ensure_zero_extrema: + self.G[:self.n_dofs_image] = 0 + self.G[-self.n_dofs_image:] = 0 + # Should be the same if we project the gradient before, instead # of the total force # nebm_clib.project_images(self.G, y, @@ -388,7 +394,7 @@ def compute_distances(self): # else: # return 0 - def step_RHS(self, t, y): + def step_RHS(self, band): """ This function is called on every iteration of the integrators in @@ -400,7 +406,7 @@ def step_RHS(self, t, y): # Update the effective field, energies, spring forces and tangents # using the *y* array - self.nebm_step(y) + self.nebm_step(band) # Now set the RHS of the equation as the effective force on the energy # band, which is stored on the self.G array