Physical test scenario 2 - single spin precession-nutation motion dynamics with time-dependent driving field

[AdamsMP92]

Motivation

A general physical validation test for LLG dynamics can be constructed by prescribing an analytical magnetisation trajectory and deriving the time-dependent field required to generate it.

This provides an exact test of:

• time-dependent Zeeman fields,
• simultaneous precession and damping,
• the Gilbert-to-Landau–Lifshitz conversion,
• signs and cross-product ordering,
• and the numerical time integration.

Manufactured trajectory

Consider the normalised magnetisation

$$\mathbf{m}(t) = \mathbf{e}_r(\theta(t), \phi(t)),$$

and

$$\dot{\mathbf{m}}(t) = \dot{\theta} \mathbf{e}_\theta(\theta(t), \phi(t)) + \dot{\phi}\sin\theta \mathbf{e}_\phi(\theta(t), \phi(t)),$$

with

$$\begin{align} \theta(t) &= \theta_0+\Delta\theta\sin(\omega_n t) \\\ \phi(t) &= \omega_p t \end{align}$$

Here, $\omega_p$ is the azimuthal precession frequency, $\omega_n$ is the polar oscillation frequency, and $\Delta\theta$ is the polar oscillation amplitude.

The explicit Landau–Lifshitz equation is

$$\dot{\mathbf{m}} = -\gamma\,\mathbf{m}\times\mathbf{H}_{\mathrm{eff}} - \beta\,\mathbf{m}\times \left( \mathbf{m}\times\mathbf{H}_{\mathrm{eff}} \right),$$

where the coefficients corresponding to Ubermag's Gilbert convention are

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

Only the field component perpendicular to $\mathbf{m}$ contributes to the dynamics. For a prescribed trajectory, this field is

$$\mathbf{H}_{\mathrm{eff}}^\perp(t) = \frac{ \beta\,\dot{\mathbf{m}}(t) + \gamma\,\mathbf{m}(t)\times\dot{\mathbf{m}}(t) }{ \gamma^2+\beta^2 }.$$

Therefore, the prescribed $\mathbf{m}(t)$ is an exact solution of the LLG equation driven by $\mathbf{H}_{\mathrm{eff}}^\perp(t)$.

Analytical reference implementation

import numpy as np


def manufactured_m(t, theta0, delta_theta, omega_n, omega_p):
    """Return the prescribed normalised magnetisation trajectory."""
    t = np.asarray(t, dtype=float)

    theta = theta0 + delta_theta * np.sin(omega_n * t)
    phi = omega_p * t

    return np.column_stack(
        (
            np.sin(theta) * np.cos(phi),
            np.sin(theta) * np.sin(phi),
            np.cos(theta),
        )
    )


def manufactured_dm_dt(t, theta0, delta_theta, omega_n, omega_p):
    """Return the exact time derivative of the trajectory."""
    t = np.asarray(t, dtype=float)

    theta = theta0 + delta_theta * np.sin(omega_n * t)
    theta_dot = delta_theta * omega_n * np.cos(omega_n * t)

    phi = omega_p * t
    phi_dot = omega_p

    e_theta = np.column_stack(
        (
            np.cos(theta) * np.cos(phi),
            np.cos(theta) * np.sin(phi),
            -np.sin(theta),
        )
    )

    e_phi = np.column_stack(
        (
            -np.sin(phi),
            np.cos(phi),
            np.zeros_like(t),
        )
    )

    return (
        theta_dot[:, np.newaxis] * e_theta
        + (np.sin(theta) * phi_dot)[:, np.newaxis] * e_phi
    )


def manufactured_field(
    t,
    theta0,
    delta_theta,
    omega_n,
    omega_p,
    alpha,
    gamma0,
):
    """Return the transverse field generating the trajectory."""
    t = np.atleast_1d(np.asarray(t, dtype=float))

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

    m = manufactured_m(
        t,
        theta0,
        delta_theta,
        omega_n,
        omega_p,
    )

    dm_dt = manufactured_dm_dt(
        t,
        theta0,
        delta_theta,
        omega_n,
        omega_p,
    )

    return (
        beta * dm_dt
        + gamma * np.cross(m, dm_dt)
    ) / (gamma**2 + beta**2)

Time-dependent Zeeman field

mm.Zeeman can represent the manufactured field using a time-dependent matrix acting on a fixed reference field.

Choose

H_scale = 1e5
H_reference = (H_scale, 0, 0)

and define a matrix whose first column is the manufactured field divided by H_scale:

def field_matrix(t):
    H = manufactured_field(
        t=[t],
        theta0=theta0,
        delta_theta=delta_theta,
        omega_n=omega_n,
        omega_p=omega_p,
        alpha=alpha,
        gamma0=gamma0,
    )[0]

    hx, hy, hz = H / H_scale

    return [
        hx, 0.0, 0.0,
        hy, 0.0, 0.0,
        hz, 0.0, 0.0,
    ]

The model can then be defined as

system.energy = mm.Zeeman(
    H=H_reference,
    func=field_matrix,
    dt=field_dt,
)

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

Proposed test setup

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

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

theta0 = np.pi / 2
delta_theta = 0.5

omega_n = 2 * np.pi * 0.5e9
omega_p = 2 * np.pi * 2e9

H_scale = 1e5
field_dt = 1e-13

m0 = (
    np.sin(theta0),
    0.0,
    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="manufactured_llg_trajectory")

system.energy = mm.Zeeman(
    H=(H_scale, 0, 0),
    func=field_matrix,
    dt=field_dt,
)

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

system.m = df.Field(
    mesh,
    nvdim=3,
    value=m0,
    norm=Ms,
)

A natural simulation interval is one complete polar oscillation,

$$T=\frac{2\pi}{\omega_n}.$$

Primary assertion

The simulated spatially averaged magnetisation should be normalised and compared directly with the manufactured trajectory:

m_exact = manufactured_m(
    t=times,
    theta0=theta0,
    delta_theta=delta_theta,
    omega_n=omega_n,
    omega_p=omega_p,
)

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

The tolerances should account for both the time-integration error and the temporal discretisation of the field callable.

Expected value

Unlike a constant-field relaxation test, this scenario exercises a field whose direction and magnitude both vary in time.

The test would validate the complete dynamical pathway from the Ubermag model to the calculator output using a nontrivial exact trajectory. The same construction can be reused to generate further manufactured LLG test cases from arbitrary differentiable trajectories on the unit sphere.