Physical test scenario 1 - single spin precession-relaxation dynamics

[AdamsMP92]

Physical model

In Ubermag, the dynamics are defined as

system.dynamics = (
    mm.Precession(gamma0=gamma0)
    + mm.Damping(alpha=alpha)
)

The corresponding Gilbert equation is

$$\frac{d\mathbf{m}}{dt} = -\gamma_0\, \mathbf{m}\times\mathbf{H}_{\mathrm{eff}} + \alpha\, \mathbf{m}\times\frac{d\mathbf{m}}{dt}.$$

Solving explicitly for the time derivative gives

$$\frac{d\mathbf{m}}{dt} = -\frac{\gamma_0}{1+\alpha^2} \mathbf{m}\times\mathbf{H}_{\mathrm{eff}} - \frac{\alpha\gamma_0}{1+\alpha^2} \mathbf{m}\times \left( \mathbf{m}\times\mathbf{H}_{\mathrm{eff}} \right).$$

Thus, in Landau–Lifshitz form, the corresponding coefficients are

$$\gamma = \frac{\gamma_0}{1+\alpha^2}, \qquad \beta = \frac{\alpha\gamma_0}{1+\alpha^2}.$$

Analytical solution

Considering a constant effective field $\mathbf{H}_{\mathrm{eff}}=H\mathbf{e}_z$ and a normalised magnetisation $\lVert\mathbf{m}\rVert=1$, the Landau-Lifshitz equation reduces using a suitable parametrization to a set of two independent and linear ODEs:

$$\mathbf{m}= \begin{bmatrix} m_x\\m_y\\m_z \end{bmatrix} = \begin{bmatrix} \mathrm{sech}(\psi)\cos(\phi) \\\ \mathrm{sech}(\psi)\sin(\phi) \\\ \tanh(\psi) \end{bmatrix} \quad \leadsto \quad \frac{d}{dt} \begin{bmatrix} \psi \\\ \phi \end{bmatrix} = \begin{bmatrix} \beta H \\\ \gamma H \end{bmatrix}$$

with the solution:

$$\begin{bmatrix} \psi(t)\\\phi(t) \end{bmatrix} = \begin{bmatrix} \beta H t \\\ \gamma H t \end{bmatrix} + \begin{bmatrix} \mathrm{artanh}\left(\cos\theta_0\right) \\\ \phi_0 \end{bmatrix}, \quad \mathrm{i.e.} \quad \mathbf{m}_0= \begin{bmatrix} \sin\theta_0 \cos\phi_0 \\\ \sin\theta_0 \sin\phi_0 \\\ \cos \theta_0 \end{bmatrix} .$$

A relaxation time $T_r$ may be defined through the condition $\theta_r = 1 \mathrm{degree}$, where $\cos\theta_r \simeq 0.9998$ and $\psi_r = \mathrm{artanh}(\cos\theta_r)\approx 4.741$, such that

$$T_r = \frac{\mathrm{artanh}(\cos\theta_r) - \mathrm{artanh}(\cos\theta_0)}{\beta H}.$$

Further, the precession period $T_p$ is defined as

$$T_p = \frac{2\pi}{\gamma H}$$

For controlling the number of precessional revolutions $K_p=1, 2, 3, ...$ in the interval $t\in [0, T_r]$, the fraction of the damping and gyromagnetic constant must to be selected to satisfy

$$\frac{\beta}{\gamma} = \frac{1}{2\pi} \frac{\mathrm{artanh}(\cos\theta_r) - \mathrm{artanh}(\cos\theta_0)}{K_p}.$$

Analytical reference implementation

import numpy as np


def analytical_constant_field_llg(
    t,
    theta0,
    phi0,
    H,
    alpha,
    gamma0,
):
    """Return the analytical LLG trajectory for H_eff = H e_z.

    The convention corresponds to the Gilbert parameters used by Ubermag.

    Parameters
    ----------
    t : array-like
        Times in seconds.
    theta0 : float
        Initial polar angle in radians.
    phi0 : float
        Initial azimuthal angle in radians.
    H : float
        Constant effective-field magnitude in A/m.
    alpha : float
        Gilbert damping parameter.
    gamma0 : float
        Gyromagnetic ratio in m/(A s).

    Returns
    -------
    numpy.ndarray
        Array of shape ``(len(t), 3)`` containing the normalised
        Cartesian magnetisation components.
    """
    t = np.asarray(t, dtype=float)

    gamma = gamma0 / (1.0 + alpha**2)
    beta = alpha * gamma0 / (1.0 + alpha**2)

    psi0 = np.arctanh(np.cos(theta0))

    psi = psi0 + beta * H * t
    phi = phi0 + gamma * H * t

    transverse = 1.0 / np.cosh(psi)

    mx = transverse * np.cos(phi)
    my = transverse * np.sin(phi)
    mz = np.tanh(psi)

    return np.column_stack((mx, my, mz))

Proposed test setup

A minimal single-cell system can be used so that the magnetisation remains spatially uniform and no exchange or demagnetisation terms are present.

import numpy as np
import discretisedfield as df
import micromagneticmodel as mm
import oommfc as oc

Ms = 8e5
H = 1e5
alpha = 0.1
gamma0 = mm.consts.gamma0

theta0 = np.deg2rad(60)
phi0 = np.deg2rad(20)

m0 = (
    np.sin(theta0) * np.cos(phi0),
    np.sin(theta0) * np.sin(phi0),
    np.cos(theta0),
)

mesh = df.Mesh(
    p1=(0, 0, 0),
    p2=(1e-9, 1e-9, 1e-9),
    n=(1, 1, 1),
)

system = mm.System(name="constant_field_llg")
system.energy = mm.Zeeman(H=(0, 0, H))
system.dynamics = (
    mm.Precession(gamma0=gamma0)
    + mm.Damping(alpha=alpha)
)
system.m = df.Field(
    mesh,
    nvdim=3,
    value=m0,
    norm=Ms,
)

The system can then be evolved using oc.TimeDriver.

Primary assertion

The simulated spatially averaged magnetisation should be compared directly with the analytical Cartesian solution.

m_exact = analytical_constant_field_llg(
    t=times,
    theta0=theta0,
    phi0=phi0,
    H=H,
    alpha=alpha,
    gamma0=gamma0,
)

np.testing.assert_allclose(
    m_simulated,
    m_exact,
    rtol=...,
    atol=...,
)

If the calculator output contains the physical magnetisation

$$\mathbf{M}=M_s\mathbf{m},$$

rather than the normalised direction $\mathbf{m}$, the simulated values should be divided by $M_s$ before comparison.

Suggested parameterisation

The test could be parametrised over several damping values:

alpha_values = [0.0, 0.1, 0.5]

These cases test different aspects of the implementation:

• alpha = 0.0 tests pure precession.
• alpha = 0.1 tests simultaneous precession and relaxation.
• alpha = 0.5 makes an incorrect treatment of the factor $1+\alpha^2$ clearly visible.

For $\alpha=0.5$,

$$1+\alpha^2=1.25.$$

Omitting this factor would produce a 25% error in the explicit coefficients. Equivalently, the correct precession frequency and damping rate are 20% smaller than the corresponding values without the conversion.

The initial angle should satisfy

$$0<\theta_0<\pi,$$

and should not be too close to either pole.

A nonzero initial azimuthal angle is useful for detecting component swapping and sign errors.

Expected value

This physical test would detect several possible implementation errors:

1. interpreting gamma0 using the wrong LL/Gilbert convention;
2. omitting the factor $1+\alpha^2$;
3. applying the factor $1+\alpha^2$ twice;
4. using the wrong precession sign;
5. using the wrong damping sign;
6. incorrectly coupling mm.Precession and mm.Damping;
7. translating the Ubermag dynamics incorrectly into an external solver configuration;
8. inconsistencies between calculator backends.

The analytical reference is backend-independent and could later be reused to validate additional Ubermag calculators.