From 53d0deda3405980fa09c8827870874cf5fd66a35 Mon Sep 17 00:00:00 2001 From: davidcorteso Date: Tue, 12 Jul 2016 23:32:33 +0100 Subject: [PATCH 1/2] Fixed VTK saving. The save_vtk option needs to be simplified further --- fidimag/common/neb_spherical.py | 16 +++++++++++++++- 1 file changed, 15 insertions(+), 1 deletion(-) diff --git a/fidimag/common/neb_spherical.py b/fidimag/common/neb_spherical.py index 5cd5457e..adafa421 100755 --- a/fidimag/common/neb_spherical.py +++ b/fidimag/common/neb_spherical.py @@ -548,8 +548,22 @@ def save_vtks(self): self.coords.shape = (self.total_image_num, -1) + # We use Ms from the simulation assuming that all the + # images are the same for i in range(self.total_image_num): - self.vtk.save_vtk(spherical2cartesian(self.coords[i]), step=i, vtkname='m') + # We will try to save for the micromagnetic simulation (Ms) + # or an atomistic simulation (mu_s) + # TODO: maybe this can be done with an: isinstance + try: + self.vtk.save_vtk(spherical2cartesian(self.coords[i]).reshape(-1, 3), + self.sim.Ms, + step=i, + vtkname='m') + except: + self.vtk.save_vtk(spherical2cartesian(self.coords[i]).reshape(-1, 3), + self.sim._mu_s, + step=i, + vtkname='m') self.coords.shape = (-1, ) From b72d42361eb9fff6eabd1503a599723d37a602e0 Mon Sep 17 00:00:00 2001 From: davidcorteso Date: Thu, 28 Jul 2016 17:33:39 +0100 Subject: [PATCH 2/2] Double checked nebm in spherical coordinates. It appears to work much slower than the Cartesian version. The test using two particles is working now for both coordinate systems giving similar results, thus Spherical coordinates are working. Other tests are too slow to implement for now, like the elongated particle, this may indicate that we need to optimise the nebm code somehow --- fidimag/common/neb_spherical.py | 15 ++-- fidimag/common/vtk.py | 1 + ...articles.py => test_two_particles_nebm.py} | 70 +++++++++++++------ 3 files changed, 57 insertions(+), 29 deletions(-) rename tests/{test_two_particles.py => test_two_particles_nebm.py} (66%) diff --git a/fidimag/common/neb_spherical.py b/fidimag/common/neb_spherical.py index adafa421..2157105e 100755 --- a/fidimag/common/neb_spherical.py +++ b/fidimag/common/neb_spherical.py @@ -365,12 +365,11 @@ def init_m(pos): self.coords = np.zeros(2 * self.n * self.total_image_num) self.last_m = np.zeros(self.coords.shape) - # For the effective field we use the extremes (to fit the energies - # array length). We could save some memory if we don't consider them - # (CHECK this in the future) - self.Heff = np.zeros(self.coords.shape) + # For the effective field we do not consider the extremes + # since they are fixed + self.Heff = np.zeros(2 * self.n * self.image_num) # self.Heff = np.zeros(2 * self.n * self.total_image_num) - self.Heff.shape = (self.total_image_num, -1) + self.Heff.shape = (self.image_num, -1) # Tangent components in spherical coordinates self.tangents = np.zeros(2 * self.n * self.image_num) @@ -620,7 +619,7 @@ def compute_effective_field(self, y): # Set the simulation magnetisation to the (i+1)-th image # spin components - self.sim.spin = spherical2cartesian(y[i + 1]) + self.sim.spin[:] = spherical2cartesian(y[i + 1]) # Compute the effective field using Fidimag's methods. # (we use the time=0 since we are only using the simulation @@ -632,7 +631,7 @@ def compute_effective_field(self, y): h = self.sim.field # Save the field components of the (i + 1)-th image # to the effective field array which is in spherical coords - self.Heff[i + 1, :] = cartesian2spherical_field(h, y[i + 1]) + self.Heff[i, :] = cartesian2spherical_field(h, y[i + 1]) # Compute and save the total energy for this image self.energy[i + 1] = self.sim.compute_energy() @@ -688,7 +687,7 @@ def sundials_rhs(self, time, y, ydot): # D_climb = -nabla E + 2 * [nabla E * t] t for i in range(self.image_num): # The eff field has (i + 1) since it considers the extreme images - h = self.Heff[i + 1] + h = self.Heff[i] # Tangents and springs has length: total images -2 (no extremes) t = self.tangents[i] sf = self.springs[i] diff --git a/fidimag/common/vtk.py b/fidimag/common/vtk.py index e0af17eb..c6f69678 100644 --- a/fidimag/common/vtk.py +++ b/fidimag/common/vtk.py @@ -51,6 +51,7 @@ def init_VtkData(self, structure, header): def reset_data(self): self.vtk_data = pyvtk.VtkData(self.structure, self.header) + @ignore_stderr def save_scalar(self, s, name="my_field", step=0): self.vtk_data.cell_data.append(pyvtk.Scalars(s, name)) diff --git a/tests/test_two_particles.py b/tests/test_two_particles_nebm.py similarity index 66% rename from tests/test_two_particles.py rename to tests/test_two_particles_nebm.py index e69126a7..1239d6e3 100644 --- a/tests/test_two_particles.py +++ b/tests/test_two_particles_nebm.py @@ -4,8 +4,9 @@ # FIDIMAG: from fidimag.micro import Sim from fidimag.common import CuboidMesh -from fidimag.micro import UniformExchange, UniaxialAnisotropy, Demag +from fidimag.micro import UniformExchange, UniaxialAnisotropy from fidimag.common.neb_cartesian import NEB_Sundials +from fidimag.common.neb_spherical import NEB_Sundials as NEB_Sundials_spherical import numpy as np # Material Parameters @@ -25,7 +26,7 @@ def two_part(pos): - x, y = pos[0], pos[1] + x = pos[0] if x > 6 or x < 3: return Ms @@ -34,15 +35,16 @@ def two_part(pos): # Finite differences mesh mesh = CuboidMesh(nx=3, - ny=1, - nz=1, - dx=3, dy=3, dz=3, - unit_length=1e-9 - ) + ny=1, + nz=1, + dx=3, dy=3, dz=3, + unit_length=1e-9 + ) # Simulation Function -def relax_neb(k, maxst, simname, init_im, interp, save_every=10000): +def relax_neb(k, maxst, simname, init_im, interp, save_every=10000, + coordinates='Cartesian'): """ Execute a simulation with the NEB function of the FIDIMAG code, for an elongated particle (long cylinder) @@ -78,11 +80,18 @@ def relax_neb(k, maxst, simname, init_im, interp, save_every=10000): # equal to 'the number of initial states specified', minus one. interpolations = interp - neb = NEB_Sundials(sim, - init_images, - interpolations=interpolations, - spring=k, - name=simname) + if coordinates == 'Cartesian': + neb = NEB_Sundials(sim, + init_images, + interpolations=interpolations, + spring=k, + name=simname) + elif coordinates == 'Spherical': + neb = NEB_Sundials_spherical(sim, + init_images, + interpolations=interpolations, + spring=k, + name=simname) neb.relax(max_steps=maxst, save_vtk_steps=save_every, @@ -90,16 +99,17 @@ def relax_neb(k, maxst, simname, init_im, interp, save_every=10000): stopping_dmdt=1e-2) +def mid_m(pos): + if pos[0] > 4: + return (0.5, 0, 0.2) + else: + return (-0.5, 0, 0.2) + + # this test runs for over a minute @pytest.mark.slow def test_energy_barrier_2particles(): # Initial images: we set here a rotation interpolating - def mid_m(pos): - if pos[0] > 4: - return (0.5, 0, 0.2) - else: - return (-0.5, 0, 0.2) - init_im = [(-1, 0, 0), mid_m, (1, 0, 0)] interp = [6, 6] @@ -114,12 +124,27 @@ def mid_m(pos): interp, save_every=5000) + # Relax the same system using spherical coordinates + relax_neb(float(k), 2000, + 'neb_2particles_k{}_10-10int_spherical'.format(k), + init_im, + interp, + save_every=5000, + coordinates='Spherical' + ) + # Get the energies from the last state data = np.loadtxt('neb_2particles_k1e8_10-10int_energy.ndt')[-1][1:] + data_spherical = np.loadtxt( + 'neb_2particles_k1e8_10-10int_spherical_energy.ndt')[-1][1:] ebarrier = np.abs(np.max(data) - np.min(data)) / (1.602e-19) print(ebarrier) + ebarrier_spherical = np.abs(np.max(data_spherical) - + np.min(data_spherical)) / (1.602e-19) + print(ebarrier_spherical) + # Analitically, the energy when a single particle rotates is: # K V cos^2 theta # where theta is the angle of the direction of one particle with respect @@ -133,6 +158,9 @@ def mid_m(pos): assert ebarrier < 0.017 assert ebarrier > 0.005 + assert ebarrier_spherical < 0.017 + assert ebarrier_spherical > 0.005 + -if __name__=='__main__': - test_energy_barrier_2particles() +if __name__ == '__main__': + test_energy_barrier_2particles()