Magnon Torque Spintronics Spin Wave Torque Magnetization Dynamics Llgs Equation

# Magnon Torque Spintronics, Spin-Wave Momentum Transfer, LLGS Dynamics, and All-Magnonic Logic

## 1. Introduction: Spin-Orbit Coupling and Spintronics

Spin-orbit coupling (SOC) is the interaction between the orbital angular momentum of a charge carrier and its spin magnetic moment. Quantitatively, the energy shift is:

In bulk materials, SOC is weak (typically < 10 meV). However, in nanostructures with broken inversion symmetry—such as heterojunctions at interfaces—SOC can become strong enough to compete with kinetic energy and enable new physics.

This has revolutionized spintronics: using electron spin (in addition to charge) to carry information. The key discovery—spin-orbit torque (SOT)—provides a mechanism to switch magnetization electrically without spin-polarized current injection, offering a path to ultra-low-power magnetic memory.

## 2. Rashba Spin-Orbit Interaction

The Rashba effect occurs at interfaces lacking inversion symmetry (e.g., 2D electron gas at a Schottky barrier). The Rashba Hamiltonian:

where:
- α_R is the Rashba coupling strength
- k = (k_x, k_y) is the in-plane wave vector
- σ_x,y are Pauli matrices

### Key Features

Spin texture: At each k, the spin is locked perpendicular to the momentum:

This creates a helical spin texture in k-space—spins precess as carriers move through the lattice.

Band splitting: The energy bands split into two:

where k = |k|. The spin splitting at the Fermi level:

can reach several meV in strong SOC materials (e.g., BiSb, topological insulators).

## 3. Dresselhaus Spin-Orbit Interaction

In materials with cubic inversion symmetry (bulk semiconductors like GaAs), the Dresselhaus interaction dominates:

where β is the Dresselhaus coupling constant.

Unlike Rashba (which is tunable), Dresselhaus is intrinsic to the material. In quasi-2D structures, both Rashba and Dresselhaus contribute.

## 4. Spin-Orbit Torque: The Fieldlike and Anti-Damping Components

When a charge current flows through a material with SOC, the spin texture of the current-carrying electrons exerts a torque on magnetic moments. The spin-orbit torque (SOT) on a magnetic moment m is:

where θ_SH is the spin Hall angle and σ is the spin accumulation at the magnet-SOC interface.

### Anti-Damping (Fieldlike) Torque

The torque can be decomposed into two components:

Anti-damping (damping-like) torque:

This torque is in-plane and acts to precess the magnetization. It can overcome damping and enable coherent switching.

Fieldlike torque:

This acts out-of-plane and is similar to a Zeeman field (∝ m × H_eff).

## 5. Landau-Lifshitz-Gilbert-Slonczewski Equation

The LLGS equation governs magnetization dynamics under applied fields and spin-orbit torques:

where:
- γ is the gyromagnetic ratio (~2.76 × 10¹¹ rad/s/T for electrons)
- α is the Gilbert damping parameter (~0.01–0.1)
- H_eff includes applied field, anisotropy, and exchange interactions
- τ_SOT is the spin-orbit torque

Rearranging to separate the damping term:

## 6. Critical Switching Current for SOT-MRAM

The threshold current density for SOT switching depends on the torque efficiency and damping:

where:
- e is electron charge
- M_s is saturation magnetization
- t is magnetic layer thickness
- ξ = ξ_fl / ξ_ad is the ratio of fieldlike to anti-damping torques
- θ is the switching angle
- θ_SH ≈ 0.05–0.4 is the spin Hall angle (material dependent)

For a thin magnetic film (t ≈ 1 nm) with α = 0.02 and θ_SH = 0.1:

This is 10–100× lower than spin-transfer torque (STT) requirements, enabling sub-10 nm MRAM cells with μA switching currents.

## 7. Python Implementation: SOT Switching Dynamics

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def llgs_equation(m, t, H_eff, gamma, alpha, tau_sot_coeff, current_pulse):
    """
    Solve the Landau-Lifshitz-Gilbert-Slonczewski equation with SOT.
    
    dm/dt = -γ(m × H_eff) + α m × (dm/dt) + τ_SOT
    
    Parameters:
    -----------
    m : array (3,)
        Magnetization unit vector [mx, my, mz]
    t : float
        Time
    H_eff : array (3,)
        Effective magnetic field
    gamma : float
        Gyromagnetic ratio (rad/s/T)
    alpha : float
        Gilbert damping
    tau_sot_coeff : float
        SOT torque coefficient (proportional to current)
    current_pulse : callable
        Current as function of time I(t)
    
    Returns:
    --------
    dm_dt : array (3,)
        Magnetization time derivative
    """
    
    # Normalize magnetization (maintain |m| = 1)
    m = m / np.linalg.norm(m)
    
    # Current-dependent SOT
    I_t = current_pulse(t)
    tau_ad = tau_sot_coeff * I_t  # Anti-damping torque coefficient
    tau_fl = 0.2 * tau_sot_coeff * I_t  # Fieldlike (20% of anti-damping)
    
    # Anti-damping torque: ξ_ad (ĵ × m̂) × m̂
    j_hat = np.array([1, 0, 0])  # Current in x-direction
    tau_ad_vec = tau_ad * np.cross(np.cross(j_hat, m), m)
    
    # Fieldlike torque: ξ_fl (ĵ × m̂)
    tau_fl_vec = tau_fl * np.cross(j_hat, m)
    
    # Total torque
    tau_total = tau_ad_vec + tau_fl_vec
    
    # LLGS equation: dm/dt = -γ(m × H) + α m × (dm/dt) + τ_SOT
    # Rearranging: dm/dt = (1/(1+α²)) * [(-γ m × H) + τ + α m × (-γ m × H) + α m × τ]
    
    precession = -gamma * np.cross(m, H_eff)
    total_torque = precession + tau_total
    
    # Account for damping
    damping_term = alpha * np.cross(m, total_torque)
    
    dm_dt = (total_torque + damping_term) / (1 + alpha**2)
    
    return dm_dt

# Simulation parameters
gamma_param = 1.76e11  # rad/s/T (gyromagnetic ratio)
alpha_param = 0.02  # Gilbert damping
K_u = 4e5  # Uniaxial anisotropy energy (J/m³)
M_s = 1.2e6  # Saturation magnetization (A/m)

# Derived parameters
H_k = 2 * K_u / (mu_0 * M_s) if 'mu_0' not in locals() else K_u / (0.5e-6 * M_s)
tau_sot_base = 0.05  # Normalized SOT strength

# Define current pulse
def current_pulse(t):
    """Current pulse: 1 mA from t=0 to t=2 ns, then off"""
    I_max = 1e-3  # 1 mA
    t_pulse = 2e-9  # 2 ns
    if t < t_pulse:
        return I_max * np.exp(-(t - t_pulse/2)**2 / (0.5e-9)**2)  # Gaussian
    return 0

# Time array and initial magnetization
t = np.linspace(0, 5e-9, 1000)  # 5 ns
m_init = np.array([0.1, 0, 0.99])  # Start slightly off x-axis (in +z direction)
m_init = m_init / np.linalg.norm(m_init)

# Effective field (anisotropy along z)
H_eff_static = np.array([0, 0, H_k])

# Solve ODE
def llgs_rhs(m, t):
    return llgs_equation(m, t, H_eff_static, gamma_param, alpha_param, 
                        tau_sot_base, current_pulse)

solution = odeint(llgs_rhs, m_init, t)

# Extract trajectories
m_x = solution[:, 0]
m_y = solution[:, 1]
m_z = solution[:, 2]

# Compute energy
K_energy = -0.5 * K_u * m_z**2  # Uniaxial anisotropy energy density

# Plotting
fig, axes = plt.subplots(2, 3, figsize=(15, 10))

# Panel A: Applied current pulse
I_pulse = np.array([current_pulse(ti) for ti in t])
axes[0, 0].plot(t*1e9, I_pulse*1e3, 'b-', linewidth=2)
axes[0, 0].set_xlabel('Time (ns)', fontsize=11)
axes[0, 0].set_ylabel('Current (mA)', fontsize=11)
axes[0, 0].set_title('SOT Current Pulse', fontsize=12, fontweight='bold')
axes[0, 0].grid(True, alpha=0.3)

# Panel B: Magnetization components vs time
axes[0, 1].plot(t*1e9, m_x, 'r-', linewidth=2, label='$m_x$')
axes[0, 1].plot(t*1e9, m_y, 'g-', linewidth=2, label='$m_y$')
axes[0, 1].plot(t*1e9, m_z, 'b-', linewidth=2, label='$m_z$')
axes[0, 1].set_xlabel('Time (ns)', fontsize=11)
axes[0, 1].set_ylabel('Magnetization Component', fontsize=11)
axes[0, 1].set_title('Magnetization Precession', fontsize=12, fontweight='bold')
axes[0, 1].legend(fontsize=10)
axes[0, 1].grid(True, alpha=0.3)
axes[0, 1].axhline(y=0, color='k', linestyle='-', alpha=0.2, linewidth=0.5)

# Panel C: Magnitude verification
m_magnitude = np.linalg.norm(solution, axis=1)
axes[0, 2].plot(t*1e9, m_magnitude, 'k-', linewidth=2)
axes[0, 2].set_xlabel('Time (ns)', fontsize=11)
axes[0, 2].set_ylabel('$|\mathbf{m}|$', fontsize=11)
axes[0, 2].set_title('Magnetization Magnitude Conservation', fontsize=12, fontweight='bold')
axes[0, 2].set_ylim([0.99, 1.01])
axes[0, 2].grid(True, alpha=0.3)

# Panel D: 3D trajectory in magnetization space
ax_3d = fig.add_subplot(2, 3, 4, projection='3d')
ax_3d.plot(m_x, m_y, m_z, 'b-', linewidth=2)
ax_3d.scatter([m_init[0]], [m_init[1]], [m_init[2]], color='green', s=100, label='Initial')
ax_3d.scatter([m_x[-1]], [m_y[-1]], [m_z[-1]], color='red', s=100, label='Final')
ax_3d.set_xlabel('$m_x$', fontsize=10)
ax_3d.set_ylabel('$m_y$', fontsize=10)
ax_3d.set_zlabel('$m_z$', fontsize=10)
ax_3d.set_title('Magnetization Trajectory', fontsize=12, fontweight='bold')
ax_3d.legend()

# Panel E: Phase space (m_x vs m_y)
axes[1, 1].plot(m_x, m_y, 'b-', linewidth=2)
axes[1, 1].scatter([m_init[0]], [m_init[1]], color='green', s=100, label='Start')
axes[1, 1].scatter([m_x[-1]], [m_y[-1]], color='red', s=100, label='End')
# Unit circle (for reference if m_z varies)
theta_circle = np.linspace(0, 2*np.pi, 100)
axes[1, 1].plot(np.cos(theta_circle)*0.1, np.sin(theta_circle)*0.1, 'k--', alpha=0.3)
axes[1, 1].set_xlabel('$m_x$', fontsize=11)
axes[1, 1].set_ylabel('$m_y$', fontsize=11)
axes[1, 1].set_title('In-Plane Magnetization Trajectory', fontsize=12, fontweight='bold')
axes[1, 1].legend(fontsize=10)
axes[1, 1].grid(True, alpha=0.3)
axes[1, 1].axis('equal')

# Panel F: Energy evolution
axes[1, 2].plot(t*1e9, K_energy, 'purple', linewidth=2)
axes[1, 2].set_xlabel('Time (ns)', fontsize=11)
axes[1, 2].set_ylabel('Anisotropy Energy (J/m³)', fontsize=11)
axes[1, 2].set_title('Magnetic Energy During Switching', fontsize=12, fontweight='bold')
axes[1, 2].grid(True, alpha=0.3)

# Remove the automatic 3D plot at position (1, 0)
axes[1, 0].axis('off')

plt.tight_layout()
plt.savefig('sot_mram_switching.png', dpi=150, bbox_inches='tight')
plt.show()

print("=== SOT-MRAM Switching Simulation Complete ===")
print(f"Initial magnetization: ({m_init[0]:.3f}, {m_init[1]:.3f}, {m_init[2]:.3f})")
print(f"Final magnetization: ({m_x[-1]:.3f}, {m_y[-1]:.3f}, {m_z[-1]:.3f})")
print(f"Switching angle: {np.arccos(np.dot(m_init, solution[-1]))*180/np.pi:.1f}°")
print(f"Energy change: {K_energy[-1] - K_energy[0]:.2e} J/m³")

## 8. Comparison: STT vs. SOT vs. Field-Driven Switching

SOT dominates for future sub-5nm MRAM due to ultra-low switching current density (J_c ≈ 10⁷ A/cm²).

## 9. Material Platforms for SOT: Heavy Metals and Topological Insulators

Efficient SOT requires large spin Hall angle θ_SH:

• Pt (θ_SH ≈ 0.07): Industry standard, compatible with magnetic layer deposition
• Ta (θ_SH ≈ 0.15): Higher efficiency, but oxidation concerns
• Bi₂Se₃ (topological insulator, θ_SH ≈ 1 predicted): Ultimate efficiency but integration challenges
• Heterostructures: WTe₂, MoS₂ with magnetic proximity effects

## 10. Multilayer Stacks and Interfacial Effects

Optimal SOT-MRAM stack architecture:

Substrate / Capping Layer (SiO₂)
    ↓
Magnetic Layer (CoFe, 0.5–1 nm)
    ↓
Interface (oxidized for damping control)
    ↓
Heavy Metal / SOC Layer (Pt/Ta, 5–10 nm)
    ↓
Back Contact (Cu, W)

Interfacial Dzyaloshinskii-Moriya Interaction (DMI) at the magnetic-SOC interface:

creates skyrmions at certain field/current values, enabling alternative SOT switching mechanisms with even lower power.

## 11. Thermal Effects and Temperature-Dependent Switching

At elevated temperature, thermal fluctuations assist switching:

where β ≈ 0.5–1.0. This thermal assistance can reduce switching current by ~20% at room temperature, but introduces switching noise.

## 12. Endurance and Reliability in SOT-MRAM

Unlike STT (which suffers from current-induced magnetic degradation), SOT benefits from:
- Lower current density → less Joule heating
- Interface-controlled torque → no current through magnetic layer
- Simpler stack → better manufacturability

Projected endurance: >10¹⁵ cycles (vs. 10¹⁰ for FLASH).

## 13. Read-Out: Tunnel Magnetoresistance (TMR) and Related Effects

Information is read via Tunnel Magnetoresistance (TMR):

where P is the spin polarization of ferromagnet and barrier.

Modern junctions achieve TMR > 300% in CoFeB/MgO/CoFeB stacks, enabling low-power reading.

## 14. Challenges and Open Questions

• Symmetry-breaking: Fieldlike torque limits switching efficiency; tuning ξ_fl/ξ_ad is key
• Thermal stability: Energy barrier Δ E = K_u V must exceed k_B T for data retention
• Latency: Sub-nanosecond switching requires careful pulse engineering
• Fabrication: 10–20 nm feature sizes demand tight material and stack control

## 15. Future Directions: Neuromorphic and Quantum Applications

Emerging SOT platforms:
- Neuromorphic computing: Stochastic switching for probabilistic algorithms
- Quantum simulation: SOT-driven magnetization dynamics encode many-body physics
- Reconfigurable devices: SOT for rapid field reprogramming of magnonic circuits
- Spin-orbit qubits: Using spin precession for quantum information processing

SOT-MRAM represents a watershed moment: for the first time, magnetic information can be written with comparable energy to CMOS, bridging the computer architecture energy wall.

## Magnon angular-momentum transport

A magnon is a quantized spin-wave excitation. In a ferromagnet each magnon reduces the spin projection along the equilibrium magnetization by approximately one unit of ℏ. A population imbalance therefore carries angular momentum even when no net charge crosses the device. This is the basis of non-local magnon transduction: a metal injector creates spin accumulation, interfacial exchange launches magnons into an insulating magnet, and a remote detector converts the arriving angular momentum back to voltage or torque.

For a local unit magnetization m and magnon angular-momentum current J_m, the damping-like torque density requested for the device model is

Here M_s is saturation magnetization and g_m collects interface transparency, normalization, and unit conversion. Since the double cross product removes the component parallel to m, only transverse angular momentum changes the magnetization direction. A field-like term proportional to m × J_m may also appear when the spin-mixing conductance has an appreciable imaginary part.

The extended Landau-Lifshitz-Gilbert-Slonczewski equation is

The effective field contains applied, anisotropy, exchange, demagnetizing, and thermal contributions. In SI simulations every term must have the same time-rate units. If J_m is defined as angular momentum per area per time, g_m must include magnetic thickness. If it is defined as an equivalent spin-current density, the familiar factor γℏ/(2eM_s t_F) enters explicitly.

## Thin-film YIG spin-wave dispersion

For an in-plane magnetized YIG film, a useful exchange-dipolar approximation is

Here ω_0 is the field and anisotropy contribution, A is the exchange-field coefficient under the selected unit convention, and ω_M = γμ_0M_s. The group velocity is v_g = dω/dk. Exchange makes the high-k branch increasingly quadratic, while dipolar fields dominate at small k and introduce geometry-dependent nonreciprocity not captured by the scalar approximation.

Propagation length is approximately ℓ_m = v_g/(2αω) when Gilbert damping dominates. Real YIG devices also show two-magnon scattering, surface roughness loss, nonlinear four-magnon processes, and thermal-magnon backgrounds. Reporting only Gilbert α can therefore overestimate the distance available to an all-magnonic interconnect.

## Reflection, transmission, and torque at an interface

Treat a normally incident spin wave as ψ_1 = exp(ik_1x) + r exp(-ik_1x) and the transmitted wave as ψ_2 = t exp(ik_2x). Continuity of dynamic magnetization and exchange spin current gives an impedance-like result

For lossless media, power coefficients must include flux normalization:

The factors N_i represent the modal angular-momentum normalization. Using T = |t|² without the velocity ratio can falsely violate conservation. A finite interfacial exchange J_int modifies the boundary condition and creates frequency-dependent phase delay. Damping makes k complex, so missing flux represents absorption and torque rather than transmission loss alone.

The angular momentum deposited at the interface is proportional to incident flux minus reflected and transmitted flux. For a perfectly reflecting boundary the longitudinal momentum transfer can be large even though no magnon enters the second layer. For magnetization switching, however, the transverse spin component and its dephasing length determine useful torque.

## Switching threshold and stability

Linearize the LLGS equation around a stable state. A damping-like magnon torque compensates Gilbert loss when its rate exceeds α times the precession frequency. A macrospin estimate is

This estimate is a screening tool, not a universal switching law. Thermal activation lowers the observed threshold for long pulses. Spatially nonuniform magnets nucleate domains before coherent rotation. Magnon heating changes M_s and α, while nonlinear frequency shifts detune the source. A credible switching map reports pulse length and error probability alongside current.

## All-magnonic logic and non-local transduction

Phase is a natural logic variable. A majority gate combines three coherent spin waves; constructive or destructive interference sets the output phase. Directional couplers transfer power between adjacent waveguides, while a magnonic crystal supplies stop bands and phase delay. Cascadability requires gain or regeneration because passive interference cannot restore amplitude indefinitely.

Non-local devices separate charge injection from magnetic action. Heavy-metal strips can create magnons through spin Hall torque and detect them through the inverse spin Hall effect. Thermal injection through the spin Seebeck effect is less phase coherent but useful for long-range population transport. Electrical crosstalk, substrate heating, and direct electromagnetic pickup must be measured with nonmagnetic controls.

## Python LLGS, dispersion, and interface solver

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

gamma = 1.760e11
alpha = 2.0e-4
mu0Ms = 0.176
Hbias = 0.06
Hk = 0.012
Aex = 5.0e-6
gm_over_Ms = 2.0e-3

def dispersion(k):
    w0 = gamma*(Hbias + Hk)
    wm = gamma*mu0Ms
    x = w0 + gamma*Aex*k*k
    return np.sqrt(x*(x+wm))

k = np.linspace(0, 3e6, 1200)
omega = dispersion(k)
vg = np.gradient(omega, k)

def interface_coeff(z1, z2, vg1, vg2):
    r = (z1-z2)/(z1+z2)
    t = 2*z1/(z1+z2)
    R = abs(r)**2
    T = (vg2/vg1)*abs(t)**2
    return r, t, R, T

def llgs(t, m, jm):
    m = m/np.linalg.norm(m)
    heff = np.array([Hbias, 0.0, Hk*m[2]])
    precession = -gamma*np.cross(m, heff)
    tau_m = gm_over_Ms*jm*np.cross(m, np.cross(m, np.array([0.,0.,1.])))
    rhs = precession + tau_m
    return (rhs + alpha*np.cross(m, rhs))/(1+alpha**2)

currents = np.linspace(0, 2e14, 90)
final_mz = []
trajectories = {}
for jm in currents:
    sol = solve_ivp(lambda t,m: llgs(t,m,jm), (0, 30e-9),
                    [0.02,0.0,0.9998], max_step=8e-12, rtol=2e-6, atol=1e-8)
    final_mz.append(sol.y[2,-1]/np.linalg.norm(sol.y[:,-1]))
    if jm in currents[[15,45,75]]:
        trajectories[jm] = sol

fig, ax = plt.subplots(1,3,figsize=(14,4))
ax[0].plot(k/1e6, omega/(2*np.pi*1e9))
ax[0].set(xlabel='wave number (1/um)', ylabel='frequency (GHz)')
ax[1].plot(currents, final_mz)
ax[1].set(xlabel='magnon current (arb. SI)', ylabel='final m_z')
for jm,sol in trajectories.items():
    ax[2].plot(sol.t*1e9, sol.y[2], label=f'{jm:.2e}')
ax[2].set(xlabel='time (ns)', ylabel='m_z'); ax[2].legend()
fig.tight_layout(); plt.show()

The numerical coefficient gm_over_Ms is deliberately explicit because magnon-current conventions differ across experiments. Calibrate it from absorbed angular momentum or a measured damping change before interpreting the threshold. Check norm conservation, halve the time step, and compare macrospin results with micromagnetics when the device is larger than an exchange length.

## Experimental validation checklist

• Measure microwave pickup with the magnetic layer saturated in a geometry that suppresses spin-wave excitation.
• Separate coherent magnons from thermal magnons using phase-sensitive detection and temperature controls.
• Normalize transmission by group velocity and modal angular momentum.
• Extract complex reflection, not power alone, because phase controls interference logic.
• Report magnetic thickness, M_s, α, exchange stiffness, bias field, and current convention.
• Check threshold scaling versus pulse duration to distinguish instability from thermal activation.
• Use injector-detector spacing sweeps to extract propagation length.
• Demonstrate cascading with a loaded downstream gate, not only an isolated interference pattern.

Magnon torque becomes useful when the full angular-momentum ledger closes: injected spin, reflected flux, transmitted flux, absorbed torque, and thermal loss must agree. That accounting links elegant spin-wave physics to switchable, cascadable hardware.