Integrators for the Landau-Lifshitz-Gilbert Equation

One of the key features of magnum.fe is the time integration of the Landau-Lifshitz-Gilbert equation (LLG). The LLG describes the motion of the magnetization configuration in micromagnetics:

\[\partial_t \vec{m} = - \gamma ( \vec{m} \times \vec{H}_\text{eff} ) + \alpha ( \vec{m} \times \partial_t \vec{m} )\]

where the effective field \(\vec{H}_\text{eff}\) is given by the variational derivative of the total free energy

\[\vec{H}_\text{eff} = - \frac{1}{\mu_0 M_\text{s}} \frac{\delta U(\vec{m})}{\delta \vec{m}}.\]

For a detailed discussion on the LLG, see [Brown1963]. magnum.fe implements different algorithms for the integration of the LLG. For most situations the LLGCvode integrator will provide by far the best performance.

All LLG integrators are configured with a list of LLGTerm objects. Each of these objects represents either an effective field contribution, e.g. the demagnetization field, or a more general contribution to the LLG, e.g. the spin-torque term.

# create LLG integrator that incorporates the demagnetization field
# and the exchange field as effective field
llg = LLGCvode([DemagField(), ExchangeField()])

By default all integrators solve the LLG on the ‘magnetic’ region. This default can be changed by passing a region parameter to the initializer.

# solve LLG on region named 'some_region_name'
llg = LLGCvode([DemagField(), ExchangeField()], region = 'some_region_name')

The integration is then performed stepwise by calling the step method on the solver object.

# perform timestep with dt = 1ps
llg.step(state, 1e-12)

Alternatively you can use the step method of the state object. This method does not only call the step method of the integrator, but also increases the time state.t appropriately.

# perform timestep with dt = 1ps and increase state.t
state.step(llg, 1e-12)

Depending on the type of LLG integrator, different additional options might be passed to the initializer. Please refer to the following class reference for details.

LLGCvode

class magnumfe.LLGCvode(terms=[], prec_solver='LU', region='magnetic', normtol=0.005, no_precession=False, **kwargs)

Solver for the Landau-Lifshitz-Gilbert equation using Sundials CVODE for integration. The right-hand-side of the LLG is evaluated nodewise.

If the BDF method is used, the linear solver is preconditioned as described in [Suess2002] with

\[\mat{P} = \mat{I} - \phi \frac{\delta}{\delta \vec{m}} \left[ - \gamma' ( \vec{m} \times \vec{H}_\text{prec}(\vec{m}) ) - \alpha' ( \vec{m} \times (\vec{m} \times \vec{H}_\text{prec}(\vec{m}) ) \right]\]

where \(\phi\) is the scalar parameter in the Newton matrix \(\mat{M} = \mat{I} - \phi \mat{J}\) and \(\vec{H}_\text{prec}\) contains all effective-field terms but the demagnetization field.

Arguments
terms (LLGTerm)
Effective Field terms to be considered.
prec_solver (string)
Solver for preconditioning system, either ‘LU’ or ‘GMRES’, defaults to ‘LU’
region (string)
The region where the LLG is solved, defaults to ‘magnetic’
normtol (float)
Threshold for the renormalization of the magnetization
rtol (float)
relative tolerance of CVODE integrator, defaults to 1e-5
atol (float)
absolute tolerance of CVODE integrator, defaults to 1e-5
**kwargs
Additional options for the CVODE solver (e.g. order, maxl)
reset()

Reset the CVode solver. This method should be called after hard resets of input values, e.g. after calling set on an external field object.

step(state, dt)

Calculate \(\vec{m}(t+\Delta t)\) for a given timestep.

Arguments
state (State)
The magnetization configuration.
dt (float)
The time-step size.

LLGAlougesProject

class magnumfe.LLGAlougesProject(terms=[], region='magnetic')

Solver for the Landau-Lifshitz-Gilbert equation according to Alouges [Alouges2008]. In this two-step algorithm \(\partial_t \vec{m} = \vec{v}\) is computed as

\[\int_\Omega ( \alpha \vec{v} + \vec{m} \times \vec{v}) \cdot \vec{w} \dx + \frac{2 A_\text{ex} \gamma}{\mu_0 M_\text{s}} \int_\Omega \nabla (\vec{m} + \frac{1}{2} k \vec{v}) \cdot \nabla \vec{w} \dx + \gamma \int_\Omega \vec{H}_\text{eff} \cdot \vec{w} \dx = 0 \quad \forall \quad \vec{w} \in T_{\vec{m}}\]

where both the trial and test functions are restricted to \(T_{\vec{m}}\) which is the tangent space of \(\vec{m}\) [Abert2013b]. As suggested in [Goldenits2012] the exchange field is integrated implicitly whereas all other contributions to the effective field \(\vec{H}_\text{eff}\) are integrated explicitly. This class implements the tangent-space restriction through a projection method

\[(\mat{I} - \mat{B}^T \mat{B}) \mat{A} \vec{v} = (\mat{I} - \mat{B}^T \mat{B}) \vec{b}\]

where \(\mat{A}\) and \(\vec{b}\) are the discretized bilinear and linear part of the weak formulation above and \(\mat{B}\) is a matrix that delivers the node-wise scalar product with \(\vec{m}\) if applied to a vector field (e.g. \(\vec{v}\)).

The integration step itself is performed by application of a node-wise renormalization as proposed by Alouges in [Alouges2008]:

\[\vec{m}_i(t+k) = \frac{\vec{m}_i(t) + k \vec{v}_i}{|\vec{m}_i(t) + k \vec{v}_i|}\]
Arguments
terms (LLGTerm)
Effective Field terms to be considered.
region (string)
The region where the LLG is solved, defaults to ‘magnetic’
step(state, dt)

Calculate \(\vec{m}(t+\Delta t)\) for a given timestep.

Arguments
state (State)
The magnetization configuration.
dt (float)
The time-step size.

LLGAlougesLagrange

class magnumfe.LLGAlougesLagrange(terms=[], region='magnetic')

Solver for the Landau-Lifshitz-Gilbert equation according to Alouges [Alouges2008]. In this two-step algorithm \(\partial_t \vec{m} = \vec{v}\) is computed as

\[\int_\Omega ( \alpha \vec{v} + \vec{m} \times \vec{v}) \cdot \vec{w} \dx + \frac{2 A_\text{ex} \gamma}{\mu_0 M_\text{s}} \int_\Omega \nabla (\vec{m} + \frac{1}{2} k \vec{v}) \cdot \nabla \vec{w} \dx + \gamma \int_\Omega \vec{H}_\text{eff} \cdot \vec{w} \dx = 0 \quad \forall \quad \vec{w} \in T_{\vec{m}}\]

where both the trial and test functions are restricted to \(T_{\vec{m}}\) which is the tangent space of \(\vec{m}\) [Abert2013b]. As suggested in [Goldenits2012] the exchange field is integrated implicitly whereas all other contributions to the effective field \(\vec{H}_\text{eff}\) are integrated explicitly. This class implements the tangent-space restriction through a Lagrange multiplier ansatz as proposed in [Goldenits2012].

\[\begin{split}\begin{pmatrix} \mat{A} & \mat{B} \\ \mat{B}^T & 0 \end{pmatrix} \begin{pmatrix} \vec{v} \\ \lambda \end{pmatrix} = \begin{pmatrix} \vec{b} \\ 0 \end{pmatrix}\end{split}\]

where \(\mat{A}\) and \(\vec{b}\) are the discretized bilinear and linear part of the weak formulation above and \(\mat{B}\) is a matrix that delivers the node-wise scalar product with \(\vec{m}\) if applied to a vector field (e.g. \(\vec{v}\)).

The integration step itself is performed by application of a node-wise renormalization as proposed by Alouges in [1]:

\[\vec{m}_i(t+k) = \frac{\vec{m}_i(t) + k \vec{v}_i}{|\vec{m}_i(t) + k \vec{v}_i|}\]
Arguments
terms (LLGTerm)
Effective Field terms to be considered.
region (string)
The region where the LLG is solved, defaults to ‘magnetic’
step(state, dt)

Calculate \(\vec{m}(t+\Delta t)\) for a given timestep.

Arguments
state (State)
The simulation state.
dt (float)
The time-step size.

LLGStochastic

class magnumfe.LLGStochastic(terms=[], region='magnetic', normtol=0.005, seed=None)

Solver for the Landau-Lifshitz-Gilbert equation using semi-implicit integrator suitable for stochastic integration. The right-hand-side of the LLG is evaluated nodewise.

Stochastic LLG equation:

\[\partial_t \vec{m} = - \frac{\gamma}{1+\alpha^2} \vec{m} \times \left( \vec{h}_\text{eff} + \vec{h}_\text{th}\right) + \frac{\alpha \gamma}{1+\alpha^2} \vec{m} \times \vec{m} \times \left( \vec{h}_\text{eff} + \vec{h}_\text{th} \right)\]

where the fluctuating magnetic fields \(\vec{h}_\text{th}\) are uncorrelated Gaussian white noises interpreted in the sense of Stratonovich and

\[\begin{split}\langle \vec{h}_\text{th}^{i,l}(t) \rangle &= 0 \\ \langle \vec{h}_\text{th}^{i,m}(t) \; \vec{h}_\text{th}^{j,n}(t_0) \rangle &= 2D \; \delta_{ij} \; \delta_{mn} \; \delta(t-t_0)\end{split}\]

with \(D=\frac{\alpha \, k_b \, T}{\gamma \, \mu_0 \, m_s \, V}\), \(i,j = 1, \dots, N\) and \(m,n \in {x,y,z}\). Integration of the thermal field over a timeinterval \(\Delta t\) leads to a fluctuating field \(\Delta W\), with the following properties:

\[\begin{split}\langle \Delta W^{i,m}(t) \rangle &= 0 \\ \langle \Delta W^{i,m}(t) \; \Delta W^{j,n}(t_0) \rangle &= 2D \; \delta_{ij} \; \delta_{mn} \; \delta(t-t_0) \; \Delta t\end{split}\]

Using a normally distributed random variable with zero mean and unit variance \(\eta\) allows to write \(\Delta W = \sqrt{2 \, D \, \Delta t} \, \eta\)

The stochastic PDE is often split into a deterministic and a stochastic part in order to distinguish the different temporal scaling of the contributing fields:

\[\begin{split}\Delta \vec{m} &= A(\vec{m}, t) \; \Delta t + B_k(\vec{m}, t) \; \Delta W_k \\ &= A(\vec{m}, t) \; \Delta t + B_k(\vec{m}, t) \; \sqrt{2 \, D \, \Delta t} \; \eta_k\end{split}\]

Due to the stochastic nature of the thermal field the second term contributes only with \(\sqrt{\Delta t}\) to the change in the magnetization, thereas the deterministic term contributes with \(\Delta t\). This discrepancy can be hided by moving the different scaling into the thermal field.

Finally the semiimplicit integration method reads as follows:

\[\begin{split}\partial_t \vec{m} &= \left[ A(\vec{m}, t) + B_k(\vec{m}, t) \; \sqrt{\frac{2 \, D}{\Delta t}} \; \eta_k \right] \\ \vec{m}^0_{t+\Delta t/2} &= \vec{m}_{t} + \frac{\vec{m}_{t}-\vec{m}_{t-\Delta t}}{2} \\ \vec{m}^p_{t+\Delta t/2} &= \vec{m}_{t} + \partial_t \vec{m} \left(\vec{m}^{p-1}_{t+\Delta t/2} \; , \; t + \Delta t / 2\right) \; \frac{\Delta t}{2} \\ \vec{m}_{t+\Delta t} &= \vec{m}_{t} + \partial_t \vec{m} \left(\vec{m}^p_{t+\Delta t/2} \; , \; t + \Delta t / 2\right) \; \Delta t\end{split}\]

where \(p\) is the index of the fixpoint iteration, used to solve the implicit (and nonlinear) equation.

Arguments
terms (LLGTerm)
Effective Field terms to be considered.
region (string)
The region where the LLG is solved, defaults to ‘magnetic’
normtol (float)
Threshold for the renormalization of the magnetization
seed (int)
Seed for the random number generator (defaults to ‘None’ for random seed)
step(state, dt)

Calculate \(\vec{m}(t+\Delta t)\) for a given timestep.

Arguments
state (State)
The magnetization configuration.
dt (float)
The time-step size.
[Brown1963]Brown Jr, W. F. Micromagnetics, 1963. Interscience, New York.
[Suess2002]Suess, D., Tsiantos, V., Schrefl, T., Fidler, J., Scholz, W., Forster, H., ... & Miles, J. J. (2002). Time resolved micromagnetics using a preconditioned time integration method. Journal of Magnetism and Magnetic Materials, 248(2), 298-311.
[Alouges2008](1, 2, 3) Alouges, F. (2008). A new finite element scheme for Landau-Lifchitz equations. Discrete Contin. Dyn. Syst. Ser. S, 1(2), 187-196.
[Abert2013b](1, 2) Abert, C., Exl, L., Bruckner, F., Drews, A., & Suess, D. (2013). magnum. fe: A micromagnetic finite-element simulation code based on FEniCS. Journal of Magnetism and Magnetic Materials, 345, 29-35.
[Goldenits2012](1, 2, 3) Goldenits, P., Hrkac, G., Praetorius, D., & Suess, D. (2012, February). An effective integrator for the Landau-Lifshitz-Gilbert equation. In Proceedings of Mathmod 2012 Conference.