Canonical structure of the LLG equation for exponential updates in micromagnetism

Technical Summary: Canonical Structure of the LLG Equation for Exponential Updates in Micromagnetism

Problem Statement
The numerical integration of the Landau-Lifshitz-Gilbert (LLG) equation in micromagnetism faces a fundamental challenge: preserving the unit length constraint ($\|m\|=1$) of the normalized magnetization vector $m$ during time evolution. Standard algorithms, such as the implicit Euler scheme, often fail to satisfy this constraint, leading to non-physical results where the magnetization vector deviates from the unit sphere. While post-processing techniques like renormalization or penalty methods can enforce the constraint, they introduce artificial energy contributions or fail to match the converged solution residuals, thereby compromising the physical accuracy of the simulation. Conversely, methods utilizing spherical coordinates reduce degrees of freedom but suffer from singularities at the poles, which can cause non-converging solution procedures.

Methodology
The authors propose an exponential update algorithm grounded in the geometric integration of the LLG equation on spheres. The core of the methodology involves three main steps:

1. Canonical Reformulation: The LLG equation is reformulated into a canonical tensor-algebraic structure suitable for exponential integration. By utilizing the skew-symmetric matrix representation of the cross product and applying the Grassmann identity, the evolution equation is expressed as $\partial_t m = W \cdot m$, where $W$ is a skew-symmetric matrix dependent on the effective field $H_{eff}$ and the magnetization $m$.
2. Exponential Integration: The solution is advanced in time using the matrix exponential: $m_{n+1} = \exp(W \Delta t) \cdot m_n$. This approach inherently preserves the length of the vector because the exponential of a skew-symmetric matrix belongs to the special orthogonal group $SO(3)$, ensuring that the update is a pure rotation.
3. Efficient Computation: To avoid the computational cost of infinite series expansions, the authors apply the Cayley-Hamilton theorem. By analyzing the invariants of the skew-symmetric matrix $W$ (where the first and third invariants vanish), the matrix exponential is reduced to a finite tensor representation involving only the identity matrix, $W$, and $W^2$. The coefficients are derived as trigonometric functions ($\sin(w)/w$ and $(1-\cos(w))/w^2$) of the magnitude of the vector associated with $W$, with specific case differentiations to handle singularities when $w \to 0$.

The proposed method is compared against the implicit Backward Euler method and the implicit Midpoint Rule. The Backward Euler is tested both in its standard form and with post-step renormalization.

Key Contributions

• Canonical Structure: The paper establishes a specific canonical structure for the LLG equation tailored for exponential update algorithms, a formulation the authors note has not been previously detailed in the context of micromagnetism.
• Constraint Preservation: The derived algorithm is a priori length-preserving. Unlike renormalization schemes, it does not require artificial corrections or penalty terms, thereby maintaining the physical integrity of the system's energy.
• Algorithmic Efficiency: By leveraging the specific algebraic properties of skew-symmetric matrices, the authors derive an efficient update scheme that avoids the singularities associated with spherical coordinates and the instability of standard Euler methods.

Results
Numerical experiments were conducted on a single magnetic element representing the soft magnetic material Fe$_{81.3}$Ga$_{18.7}$ under an external magnetic field. The study evaluated precessional switching and step-size stability:

• Length Conservation: The implicit Backward Euler method failed to preserve the unit length, with the vector magnitude decreasing significantly as the time step increased. While renormalization restored the length, it resulted in a solution that did not match the converged residuals of the underlying system. In contrast, the exponential update algorithm maintained $\|m\|=1$ exactly, even for coarse time discretizations.
• Positional Accuracy: In simulations of a full 360° precession, the Backward Euler and renormalized Backward Euler methods failed to return the magnetization vector to its initial position, exhibiting phase errors that grew with larger time steps. The exponential update algorithm demonstrated exact rotational behavior, returning the vector to the initial configuration regardless of the time step size.
• Stability: The exponential method proved robust against large time increments, whereas other methods required significantly smaller steps to achieve comparable accuracy, leading to higher computational costs for the same level of precision.

Significance
The authors claim that the proposed exponential update algorithm offers a superior alternative for micromagnetic simulations by strictly adhering to the geometric properties of the system. The primary significance lies in the ability to utilize larger time steps without sacrificing the physical constraint of unit magnetization length or the accuracy of the vector's trajectory on the sphere. This leads to significantly shorter computing times while maintaining high-quality solutions, addressing the trade-offs inherent in existing constraint-enforcement methods. The work builds upon and extends previous investigations into geometric integrators and quaternion-based schemes, providing a rigorous tensor-algebraic foundation for these approaches.