# Micromagnetic + Atomistic Modeling — lecture notes These notes are a plain-language companion to the 58-slide lecture. They explain what the equations mean physically, how a calculation becomes a simulated experiment, and where interpretation can go wrong. Slide ranges are included so the notes can be used as a talk track rather than as a separate textbook. The central message is simple: **a magnetic model is a controlled way to ask which magnetic states are physically possible, which are stable, how they respond to an experiment, and how quickly they can change.** ## Learning goals By the end of the lecture, students should be able to: 1. Explain why bulk measurements require a grain-scale model for interpretation. 2. Distinguish atomistic, micromagnetic, and bulk descriptions by what they average over. 3. Describe exchange, anisotropy, Zeeman, and demagnetizing energy in words before using equations. 4. Explain how a mesh, material parameters, boundary conditions, and initial state define a micromagnetic problem. 5. Distinguish a local energy minimum from the global minimum. 6. Describe how hysteresis, susceptibility, energy barriers, and relaxation times are simulated. 7. Identify numerical and physical assumptions that must be tested before trusting a result. ## Slides 1–10: Why model magnetic grains? ### The gap between theory and a bulk measurement Rock-magnetic experiments usually return a small set of bulk observables: total moment, susceptibility, remanence, coercivity, or a curve measured while field or temperature changes. Those observables are produced by a population of grains with different sizes, shapes, compositions, stress states, oxidation histories, and magnetic interactions. A bulk measurement therefore does not answer, by itself, questions such as: - Is the remanence carried by single-domain grains, vortices, domain walls, or interacting clusters? - Did the field rotate moments reversibly, move a wall, nucleate a vortex, or cross an energy barrier? - Is a temperature-dependent change caused by **Mₛ**, anisotropy, exchange, chemical ordering, oxidation, or thermal relaxation? - Would an apparently diagnostic hysteresis shape still occur for a different grain geometry? Modeling connects the observable back to possible physical causes. It does not automatically identify one unique cause. Its value is that it makes assumptions explicit and lets us test consequences. ### One material can require several length scales Different questions need different levels of description. - **Electronic-structure or first-principles models** describe electrons and bonding. They are used for site energies, valence, exchange parameters, and crystal properties. - **Atomistic spin models** assign chemical species and spins to lattice sites. They are useful for cation ordering, magnetic ordering, phase transitions, and nanoscale interfaces. - **Micromagnetic models** replace many atomic spins with a continuous magnetization field. They are useful for domains, vortices, domain walls, hysteresis, and stray fields in grains. - **Kinetic or population models** describe transitions among many stable states or many grains over laboratory and geological time. - **Bulk constitutive models** predict ensemble-scale curves without resolving each domain. No scale is universally “more realistic.” A model is useful when its resolution matches the process being asked about and its parameters can be justified. ### A succinct way to introduce multiscale modeling > We do not choose a model by asking how much detail a computer can afford. We choose it by asking which degrees of freedom must be allowed to change for the phenomenon we want to explain. ## Slide 11: Hamiltonian thinking Newton's equations ask how forces change positions and velocities. Hamiltonian thinking describes the same physics through energy and conjugate variables. In magnetic modeling, we often use “Hamiltonian” more broadly for an energy function that assigns a cost to every possible chemical and magnetic configuration. For a state **s**, write an energy **E(s)**. A low-energy state is favored relative to a high-energy state, but “favored” can mean different things: - At zero temperature, a dissipative system tends to move downhill toward a local minimum. - At finite temperature, higher-energy states can still be occupied with Boltzmann probability proportional to **exp(−E/kBT)**. - A state may be locally stable even when a lower state exists elsewhere, because an energy barrier separates them. The energy landscape picture is the thread that connects the entire lecture. Atomistic Monte Carlo samples it. Micromagnetic minimization descends it. Hysteresis follows one sequence of local minima while field changes. Nudged elastic band calculations search for passes between minima. Kinetic models use the barrier heights to estimate transition rates. ### Local and global minima Imagine a mountain landscape in fog. Dropping a ball from one location finds the nearest valley, not necessarily the lowest valley on the continent. A numerical minimizer behaves the same way. Its answer depends on the starting state and on which paths the algorithm can take. This is why one relaxed micromagnetic state is never proof that the grain has only one possible state. ### From multiplicity and entropy to the Boltzmann distribution The Boltzmann factor is easiest to understand by starting with **counting**, not with an assumed exponential law. Suppose a macroscopic condition—fixed energy, volume, and particle number—can be realized by $\Omega$ microscopic arrangements. If the allowed microstates are equally likely, a macrostate that contains more microstates is more probable. Its probability is proportional to its **multiplicity**: $$ P \propto \Omega. $$ For two independent systems, the numbers of possibilities multiply: $$ \Omega_{A+B}=\Omega_A\Omega_B. $$ Why should entropy **add** rather than multiply? More precisely, additivity is a requirement for **separable systems**: systems that are statistically independent and do not have an appreciable interaction energy. Imagine samples $A$ and $B$ separated by an ideal rigid, impermeable wall. Nothing physical about either sample changes if we call them “two systems” or draw a larger boundary around both and call them “one composite system.” A physical state variable must not depend on that bookkeeping choice. The entropy assigned to the larger boundary must therefore contain the contribution already present in $A$, plus the contribution already present in $B$: $$ \boxed{S_{\mathrm{total}}=S_A+S_B}. $$ This is what it means to call entropy an **extensive property**. If we make $n$ independent copies of the same equilibrium sample, all extensive variables scale with the number of copies: $$ S(nU,nV,nN)=nS(U,V,N). $$ The entropy doubles when the amount of independent material doubles. If it instead multiplied or squared, the measured entropy would depend on whether we described the same material as one sample or as many adjacent pieces. The macroscopic definition gives the same result. For a reversible transfer of heat, $$ dS=\frac{\delta Q_{\mathrm{rev}}}{T}. $$ If separable subsystems receive reversible heats $\delta Q_A$ and $\delta Q_B$, the entropy changes are independently measurable and the total change is $$ \begin{aligned} dS_{\mathrm{total}} &=\frac{\delta Q_A}{T_A}+\frac{\delta Q_B}{T_B}\\ &=dS_A+dS_B. \end{aligned} $$ Integrating this relation, with the same convention for the entropy zeros, gives $S_{\mathrm{total}}=S_A+S_B$. Thus additivity is both an empirical thermodynamic property and a consistency requirement for joining or separating independent systems. Now connect this macroscopic requirement to microscopic counting. Let entropy be an increasing function $f$ of multiplicity: $$ S=f(\Omega). $$ For independent systems, every microstate of $A$ can be paired with every microstate of $B$, so $$ \Omega_{\mathrm{total}}=\Omega_A\Omega_B. $$ Evaluate the total entropy in two equivalent ways. First use the multiplicity of the combined system: $$ S_{\mathrm{total}} =f(\Omega_{\mathrm{total}}) =f(\Omega_A\Omega_B). $$ Then use thermodynamic additivity: $$ S_{\mathrm{total}} =S_A+S_B =f(\Omega_A)+f(\Omega_B). $$ Equating the two expressions gives the functional relationship shown in the accompanying derivation: $$ \boxed{ f(\Omega_A\Omega_B) =f(\Omega_A)+f(\Omega_B) }. $$ We can solve this rather than simply asserting the answer. Write $$ \Omega_A=e^x,\qquad \Omega_B=e^y, $$ and define $g(x)=f(e^x)$. Because $\Omega_A\Omega_B=e^{x+y}$, the functional equation becomes $$ g(x+y)=g(x)+g(y). $$ For a continuous function, the solutions of this equation are linear: $$ g(x)=Cx. $$ Substituting $x=\ln\Omega$ back into the definition of $g$ gives $$ \boxed{S=f(\Omega)=C\ln\Omega}. $$ This result is also easy to check with $N$ independent two-state spins. Their multiplicity is $\Omega_N=2^N$, so $$ S_N=C\ln(2^N)=CN\ln2, $$ which grows linearly with $N$, as an extensive quantity should. Any logarithm base would work because changing the base only changes the constant $C$. Thermodynamics conventionally uses the natural logarithm and chooses $C=k_{\mathrm B}$: $$ \boxed{S=k_{\mathrm B}\ln\Omega}. $$ The central logic is therefore: $$ \boxed{ \begin{array}{c} \text{independent systems: }\Omega_{\mathrm{total}}=\Omega_A\Omega_B\\[3pt] \text{extensive entropy: }S_{\mathrm{total}}=S_A+S_B\\[3pt] \Downarrow\\[3pt] S=C\ln\Omega \end{array} }. $$ There is an important qualification for magnetic systems: simple additivity is exact only when the subsystems are statistically independent or their interaction is negligible. If two groups of spins are exchange-coupled, or two grains interact magnetostatically, their states are correlated and the joint multiplicity need not factor as $\Omega_A\Omega_B$. For a correlated probability distribution, $$ S(A,B)=S(A)+S(B)-k_{\mathrm B}I(A:B), $$ where $I(A:B)$ is the mutual information measuring the correlation. In practice, short-range interactions often contribute mainly near the boundary between large subsystems, so the correction becomes small compared with their bulk entropy. Strongly coupled grains or sublattices, however, should be treated as one combined system with a joint Hamiltonian. The logarithm also appears in the entropy of a probability distribution, $$ S=-k_{\mathrm B}\sum_i p_i\ln p_i. $$ The quantity $-\ln p_i$ measures how surprising outcome $i$ is, and the sum is its probability-weighted average. Maximizing this entropy while fixing $\sum_i p_i=1$ and the mean energy $\sum_i p_iE_i=U$ gives the same Boltzmann distribution derived below. In compact form, introduce Lagrange multipliers for those two constraints: $$ \mathcal L= -k_{\mathrm B}\sum_i p_i\ln p_i -\alpha\left(\sum_i p_i-1\right) -\lambda\left(\sum_i p_iE_i-U\right). $$ Setting $\partial\mathcal L/\partial p_i=0$ gives $$ -k_{\mathrm B}(\ln p_i+1)-\alpha-\lambda E_i=0, $$ so $\ln p_i=\text{constant}-(\lambda/k_{\mathrm B})E_i$, and therefore $p_i$ must be exponential in energy. Thermodynamics identifies $\lambda=1/T$, producing $p_i\propto\exp(-E_i/k_{\mathrm B}T)$. Thus the exponential is not an arbitrary curve fit: it is the highest-entropy probability distribution consistent with a fixed mean energy. Now divide the world into a small magnetic **grain** and a much larger thermal **bath**: $$ E_{\mathrm{grain}}+E_{\mathrm{bath}}=E_{\mathrm{total}}. $$ The combined system is isolated, so $E_{\mathrm{total}}$ is fixed. If the grain is in microstate $i$ with energy $E_i$, the bath must contain $$ E_{\mathrm{bath}}=E_{\mathrm{total}}-E_i. $$ The probability of grain state $i$ is proportional to the number of bath microstates compatible with that remaining energy: $$ p_i\propto \Omega_{\mathrm{bath}}(E_{\mathrm{total}}-E_i). $$ Use $S=k_{\mathrm B}\ln\Omega$, and expand the bath entropy around $E_{\mathrm{total}}$. Because the bath is very large, the grain removes only a tiny fraction of its energy: $$ \begin{aligned} S_{\mathrm{bath}}(E_{\mathrm{total}}-E_i) &\approx S_{\mathrm{bath}}(E_{\mathrm{total}}) -E_i\left(\frac{\partial S_{\mathrm{bath}}}{\partial E}\right),\\ \frac{\partial S}{\partial E}&=\frac{1}{T}. \end{aligned} $$ Therefore, $$ \begin{aligned} \Omega_{\mathrm{bath}}(E_{\mathrm{total}}-E_i) &=\exp\!\left(\frac{S_{\mathrm{bath}}}{k_{\mathrm B}}\right)\\ &\approx \underbrace{\exp\!\left[\frac{S_{\mathrm{bath}}(E_{\mathrm{total}})}{k_{\mathrm B}}\right]}_{\text{same constant for every grain state}} \exp\!\left(-\frac{E_i}{k_{\mathrm B}T}\right). \end{aligned} $$ After normalization, $$ \boxed{p_i=\frac{\exp(-E_i/k_{\mathrm B}T)}{Z}}, \qquad Z=\sum_j\exp(-E_j/k_{\mathrm B}T). $$ $Z$, the partition function, makes the probabilities sum to one. If an energy level has $g_i$ distinct microstates with the same energy, its level probability is $g_i\exp(-E_i/k_{\mathrm B}T)/Z$. The ratio of two single-state probabilities is $$ \frac{p_j}{p_i} =\exp\!\left[-\frac{E_j-E_i}{k_{\mathrm B}T}\right]. $$ The plain-language interpretation is: **when the grain takes more energy, the bath has less energy and fewer microscopic arrangements available.** The bath loses multiplicity so rapidly that the grain's probability decreases exponentially with energy. ## Slides 12–20: Atomistic modeling of Fe–Ti oxides ### What the atomistic model keeps explicit An atomistic lattice model can keep track of: - which cation occupies each crystallographic site; - the valence and magnetic moment associated with that cation; - which sites are neighbors; - whether neighboring spins prefer parallel or antiparallel alignment; - how temperature changes the distribution of chemical and magnetic states. The model does not normally calculate every electron from first principles at every Monte Carlo step. Instead, electronic-structure results, experiments, or fitted thermodynamic data are compressed into effective interaction parameters. ### Titanomagnetite chemistry in plain language Magnetite is an inverse spinel with tetrahedral **A** sites and octahedral **B** sites. The iron moments on A and B sublattices are strongly antiferromagnetically coupled. Because the sublattices contain different moment totals, their opposition leaves a net ferrimagnetic moment. In the ideal titanomagnetite solid solution, Ti⁴⁺ has a strong preference for octahedral B sites. Replacing Fe³⁺ with Ti⁴⁺ adds one positive charge, so charge balance is commonly described as reduction of another Fe³⁺ to Fe²⁺. The useful teaching statement is: > Ti substitution changes both chemistry and magnetism: it changes which ions occupy B sites, changes the Fe²⁺/Fe³⁺ ratio, dilutes magnetic exchange pathways, and changes the balance between A- and B-sublattice moments. It is safer not to say “Ti⁴⁺ always replaces Fe³⁺ on B sites” as an absolute law. That is the standard idealized substitution mechanism, but real cation distributions can depend on temperature, oxygen fugacity, non-stoichiometry, quench history, and the degree of inversion. The model should state whether site occupancies are fixed, equilibrated, or allowed to exchange. ### Chemical and magnetic energy terms A useful schematic Hamiltonian is **Etotal = Echemical + Emagnetic.** The **chemical energy** assigns costs or rewards to cations occupying particular sites and to particular neighboring pairs. A positive pair interaction makes that pair unfavorable; a negative interaction makes it favorable. This is a compact way to represent ordering, clustering, and site preference. The **magnetic energy** assigns a coupling **Jij** between neighboring spins. In a common sign convention, **Emagnetic = −Σ Jij si · sj.** Positive **Jij** favors parallel alignment and negative **Jij** favors antiparallel alignment. Authors use different sign conventions, so always read the equation before interpreting the sign. ### Monte Carlo is controlled sampling, not real-time motion At a fixed temperature, atomistic Metropolis Monte Carlo aims to generate configurations $\mathcal C$ with the canonical probability $$ \pi(\mathcal C)=\frac{\exp[-E(\mathcal C)/k_{\mathrm B}T]}{Z}. $$ In one trial: 1. Propose a change, such as swapping two cations or changing a spin. 2. Calculate $\Delta E=E_{\mathrm{new}}-E_{\mathrm{old}}$. 3. Accept the move with $$ \boxed{P_{\mathrm{accept}}=\min\!\left[1,\exp\!\left(-\frac{\Delta E}{k_{\mathrm B}T}\right)\right]}. $$ Every downhill move is accepted. An uphill move is accepted sometimes. For a symmetric proposal rule, this choice satisfies **detailed balance**: $$ \pi(\mathcal C)P(\mathcal C\rightarrow\mathcal C') = \pi(\mathcal C')P(\mathcal C'\rightarrow\mathcal C). $$ Consequently, repeated sampling settles toward the Boltzmann distribution rather than simply falling into the nearest minimum. Temperature sets the energy scale for tolerated fluctuations: - if $\Delta E\gg k_{\mathrm B}T$, an uphill move is almost never accepted; - if $\Delta E\sim k_{\mathrm B}T$, it is accepted often enough to matter; - if $T$ is high compared with the interaction energies, ordering is disrupted. Units must be consistent. If energies are in joules, use $k_{\mathrm B}=1.380649\times10^{-23}\ \mathrm{J\,K^{-1}}$. If they are in electron-volts, use $k_{\mathrm B}=8.617333\times10^{-5}\ \mathrm{eV\,K^{-1}}$. If fitted interaction parameters are already quoted in kelvin, the exponent is simply $-\Delta E_{\mathrm{(K)}}/T$. The thermal bath is implicit; the code does not need to simulate its atoms. Usually $T$ enters through the factor $\beta=1/k_{\mathrm B}T$, while magnetic and chemical order emerge from the accepted fluctuations. Effective interaction parameters may themselves have been fitted for a limited temperature range, however. If $J$, moment magnitudes, lattice constants, or other parameters are made explicitly temperature dependent, document that choice and avoid double-counting thermal effects. The ordering temperature produced by a finite classical-spin model should be checked against experiment or the intended calibration target. Uphill moves are essential. They allow the system to escape shallow traps and produce the correct thermal distribution. A Monte Carlo “sweep” is an accounting unit—roughly one attempted update per degree of freedom—not automatically a second of physical time. For a temperature series, equilibrate at each temperature, discard the initial burn-in, and average measurements over many later configurations. Check multiple random seeds, autocorrelation, system size, and both heating and cooling schedules. To infer a phase transition, do not rely only on a visually changing structure. Track an order parameter, energy, heat capacity, susceptibility, and equilibration history. The move set defines what is allowed to equilibrate. Cation swaps at fixed bulk composition sample chemical ordering; spin flips or rotations sample magnetic order. If real cations are kinetically frozen but the simulation freely swaps them, the calculation represents **equilibrium cation chemistry**, not necessarily the path followed by a cooling rock. Standard Metropolis Monte Carlo samples equilibrium probability but does not usually provide physical time. Converting sweeps to seconds requires a justified dynamical model. For magnetic kinetics, use calculated barriers with Arrhenius rates, kinetic Monte Carlo, stochastic Landau–Lifshitz–Gilbert dynamics, or another method with a physical clock. ### Hemo-ilmenite lamellar magnetism Exsolution creates interfaces between phases with different compositions, magnetic ordering temperatures, and preferred spin structures. Interface moments can arise when one phase terminates or modifies the exchange network of another. The key physical lesson is that a bulk remanence can be controlled by a small volume fraction at interfaces rather than by the dominant mineral volume. When presenting the lamellar examples, separate three ideas: 1. **Chemical architecture:** where Fe- and Ti-rich layers occur. 2. **Magnetic order:** which sublattices or layers are ordered at a given temperature. 3. **Net moment:** which opposing contributions fail to cancel. An atomistic model is valuable because those three questions are coupled at the lattice scale. ## Slides 21–29: From atoms to a continuum grain ### The micromagnetic approximation Micromagnetics replaces many atomic spins in a small volume with one continuous magnetization vector, **M(r) = Mₛ m(r)**, where **Mₛ** is saturation magnetization and **m** is a unit vector giving local direction. Standard micromagnetics assumes the magnitude is approximately fixed and allows the direction to vary smoothly in space. This approximation is excellent when the feature of interest spans many atomic spacings. It is less reliable at an atomically sharp interface, a defect whose chemistry matters explicitly, very near the Curie temperature, or at a singular structure such as the immediate core of a Bloch point. ### Why natural grains are numerically difficult Ideal materials problems often use cubes, ellipsoids, periodic films, and uniform parameters. Natural grains add irregular surfaces, pores, inclusions, compositional zoning, oxidation rims, neighboring grains, and uncertain three-dimensional geometry. These details matter because: - surface orientation controls magnetostatic charge; - constrictions can pin domain walls; - corners can nucleate reversal; - a shell can change exchange, anisotropy, and saturation magnetization; - nearby grains add an interaction field; - a small geometric feature can be smaller than the chosen mesh. The model geometry is therefore part of the physical hypothesis, not merely a drawing around the calculation. ### Finite difference versus finite element **Finite-difference models** divide space into regular rectangular cells. They are efficient, easy to parallelize, and excellent for films and regular geometries. Curved boundaries are represented by stair steps unless special treatments are used. **Finite-element models** divide the magnetic volume into tetrahedra or other irregular elements. They conform more naturally to grain surfaces, internal interfaces, and spatially varying materials. Their mesh construction and magnetostatic calculation are more complex. Neither method is inherently more physical. They discretize the same continuum equations differently. ### What one finite element means One tetrahedron does not represent one atom. It represents a volume over which the magnetization direction is approximated by interpolation from nodal values. A good mesh must be fine enough to resolve the shortest important magnetic length scale and any critical geometric feature. Two useful scales are: - **Magnetostatic exchange length:** **ℓex = √(2A/μ₀Mₛ²)**, which compares exchange with demagnetizing energy. - **Anisotropy exchange length:** approximately **√(A/K)**, which compares exchange with anisotropy. A simple domain-wall-width estimate is **δ ≈ π√(A/Keff)**. Constants and definitions vary by wall type, but the lesson is stable: stronger exchange broadens a wall; stronger anisotropy narrows it. The mesh edge length should be comfortably smaller than the feature being resolved. Convergence must be tested by refining the mesh and checking that the physical conclusion, not just the picture, remains stable. ## Slides 30–37: The micromagnetic energy terms The total energy is commonly written **Etotal = Eexchange + Eanisotropy + Edemag + EZeeman**, with optional magnetoelastic, interfacial, Dzyaloshinskii–Moriya, or other terms when the material requires them. Each term favors a different behavior. Magnetic structure is the compromise. ### Magnetocrystalline anisotropy: the crystal supplies preferred directions Spin–orbit coupling links electron spin to the crystal lattice. Rotating the magnetization relative to the crystal changes the electronic energy slightly. That directional energy is magnetocrystalline anisotropy. For cubic anisotropy, the energy density is commonly written **ecub = K₁(α₁²α₂² + α₂²α₃² + α₃²α₁²) + K₂α₁²α₂²α₃²**, where **α₁, α₂, α₃** are direction cosines of magnetization relative to the crystal axes. Direction cosines are simply the x, y, and z components of the unit vector **m** in the crystal coordinate system. The signs and relative sizes of **K₁** and **K₂** determine which crystallographic directions are easy or hard. For magnetite, the anisotropy constants change strongly with temperature and change character near the isotropic point and Verwey transition. A simulation that changes temperature but leaves **K₁**, **K₂**, **Mₛ**, and **A** fixed is not a complete temperature-dependent material model. **Plain-language explanation:** anisotropy is the crystal saying, “some directions cost less than others.” ### Exchange energy: neighboring directions resist sharp changes The continuum exchange energy is **Eexchange = ∫ A |∇m|² dV.** The gradient **∇m** measures how quickly the magnetization direction changes with position. Uniform magnetization has zero continuum exchange cost. A sharp turn has a large gradient and high cost. Exchange does not always make the whole grain uniform. It competes with other energies. A large grain may reduce demagnetizing energy by curling into a vortex even though curling costs exchange. The vortex core remains finite because exchange prevents the direction from changing infinitely sharply. **Plain-language explanation:** exchange is the magnetic stiffness that makes neighboring directions prefer a smooth conversation rather than an abrupt argument. ### Zeeman energy: the applied field rewards alignment The external-field energy is **EZeeman = −μ₀ ∫ M · Hext dV.** The dot product is largest when magnetization and field are parallel, so the minus sign makes parallel alignment lower in energy. Antiparallel alignment is higher in energy. This term is local and easy to interpret, but the response need not be simple rotation. A changing field may move a wall, shift a vortex core, nucleate a reversed region, or destroy one local minimum so the system jumps into another. **Plain-language explanation:** the applied field tilts the energy landscape toward states with moment along the field. ### Demagnetizing energy: the grain interacts with its own stray field A nonuniform magnetization and a magnetization component normal to a surface create sources of the demagnetizing field. These are often described as effective magnetic charges: - volume charge proportional to **−∇ · M**; - surface charge proportional to **M · n**. They are mathematical tools, not isolated magnetic monopoles. Solving the magnetostatic problem gives the stray field outside the grain and the demagnetizing field **Hd** inside it. A common energy form is **Edemag = −(μ₀/2) ∫ M · Hd dV.** The factor of one half avoids double-counting the self-interaction. Sign conventions differ depending on whether the field energy is written over all space or as a source–field interaction. Demagnetizing energy explains shape anisotropy. A long particle magnetized along its length has less surface charge at the long sides and usually lower magnetostatic energy than the same particle magnetized across its short dimension. Flux-closure states reduce external stray field by curling magnetization inside the grain. **Plain-language explanation:** demagnetizing energy penalizes magnetic poles and rewards flux closure. ### Effective field and equilibrium The energy terms can be converted into an effective field by a functional derivative: **Heff = −(1/μ₀Mₛ) δE/δm.** At a stationary state, the magnetization has no torque from the effective field: **m × Heff = 0.** This does not require every separate energy to be minimal. It requires the combined first-order torque from all terms to vanish. ## Slides 38–43: From an energy model to a simulated experiment ### Define the problem before pressing “run” A reproducible micromagnetic calculation needs: - three-dimensional geometry and mesh; - material regions and parameters **Mₛ**, **A**, anisotropy constants, and crystallographic orientation; - applied field and temperature protocol; - boundary conditions and external magnetic environment; - initial magnetization state; - numerical solver, tolerances, and stopping criteria; - outputs and the rule used to convert them into an experimental observable. Changing the initial state can change the relaxed state. Changing mesh density can change a nucleation field. Changing the field step can skip an instability. These are parts of the result, not bookkeeping details. ### Energy minimization versus magnetization dynamics Two common numerical tasks answer different questions. **Energy minimization** finds a nearby stationary state. It is appropriate when the experiment is treated as quasi-static and the goal is the sequence of local minima. **Landau–Lifshitz–Gilbert dynamics** evolves magnetization in time: **dm/dt = −γ m × Heff + α m × dm/dt.** The first term describes precession around the effective field. The damping term removes energy and moves the magnetization toward alignment with the effective field. Real-time interpretation requires realistic parameters, time steps, excitations, and thermal treatment. Using very large damping can be an efficient route to equilibrium but should not be presented as a literal physical trajectory. ### Simulating a hysteresis loop A quasi-static loop is usually generated as follows: 1. Saturate the model in a large positive field. 2. Reduce the field by one step. 3. Relax to a local minimum. 4. Record total moment and, ideally, the full magnetization state. 5. Use that state as the starting point for the next field step. 6. Continue through negative saturation and back. Hysteresis occurs because the model follows a history-dependent branch of local minima. At a switching field, the current minimum loses stability or the solver crosses a barrier, and the magnetization jumps to another state. The simulated coercivity depends on field step, mesh, defects, thermal activation, sweep rate, and the reversal mode allowed by the geometry. A zero-temperature quasi-static coercivity is not automatically the same as a laboratory coercivity measured over seconds. ### Simulating low-field susceptibility Susceptibility describes how moment changes in response to a small field. In the simplest static calculation, **χ ≈ ΔM/ΔH.** The perturbation must be small enough to remain in the same local basin if the goal is reversible susceptibility. A field that triggers wall motion or switching measures a nonlinear, history-dependent response instead. Temperature-dependent susceptibility requires temperature-dependent material parameters and, when relevant, thermal occupation or transitions among states. Changing only **Mₛ(T)** may reproduce part of a curve while missing changes caused by anisotropy, exchange, blocking, or phase transitions. ### Connecting the internal state to an external observable The solver returns **m(r)**. The instrument may measure something else: - a VSM or SQUID magnetometer measures total moment; - a FORC protocol measures a family of history-dependent moments; - a scanning field sensor measures a component or flux derived from the stray field; - MFM measures a tip-dependent force gradient; - electron or X-ray microscopy measures a projection or contrast mechanism. Always apply the instrument forward model before comparing a simulated state to measured data. Similar internal states can look different to different instruments, and different internal states can produce similar external fields. ## Slides 44–48: Transitions, barriers, and geological time ### A stable state is a basin in energy space Each magnetization configuration is one point in a very high-dimensional space. A local minimum is a basin. To move from state **i** to state **j**, the system must follow a path through intermediate configurations. The highest energy along the best available path is a saddle point. The barrier is **ΔEij = Esaddle − Ei.** The reverse barrier is generally different because it begins from **Ej**. ### Minimum-energy paths and the nudged elastic band method The nudged elastic band method represents a transition by a chain of intermediate images between two known endpoint states. Artificial springs keep the images distributed along the path, while physical forces push the path downhill perpendicular to itself. A climbing-image refinement drives the highest image toward the saddle. The method requires endpoints that are genuine local minima and an initial path that connects the intended mechanisms. Different initial paths can converge to different saddles. Finding one pathway does not prove it is the lowest or only pathway. ### Thermal activation and the Néel–Arrhenius idea The Boltzmann distribution answers an **equilibrium question**: how often should each state be occupied? Arrhenius theory answers a **kinetic question**: how long does it take to cross from one basin to another? Let state $i$ have energy $E_i$, and let the saddle on the relevant path have energy $E^\ddagger_{ij}$. The forward barrier is $$ \Delta E_{ij}=E^\ddagger_{ij}-E_i. $$ The fraction of thermal attempts energetic enough to reach this high-energy transition region carries the Boltzmann penalty $$ \exp\!\left(-\frac{\Delta E_{ij}}{k_{\mathrm B}T}\right). $$ If the system makes $\nu_0$ effective attempts per second, the transition rate is $$ \boxed{k_{ij}=\nu_0 \exp\!\left(-\frac{\Delta E_{ij}}{k_{\mathrm B}T}\right)}, $$ and the mean waiting time is $$ \boxed{\tau_{ij}=\tau_0 \exp\!\left(\frac{\Delta E_{ij}}{k_{\mathrm B}T}\right)}, \qquad \tau_0=\nu_0^{-1}. $$ The prefactor is not a universal constant. It depends on the curvature of the energy basin and saddle, magnetic damping, gyromagnetic dynamics, and the reversal mechanism. The exponential barrier term usually dominates, but the prefactor still needs physical justification. Do not put the state's absolute energy into the Arrhenius exponent. The relevant quantity is the **barrier above the starting state**. For a common saddle, $$ \frac{k_{ij}}{k_{ji}} =\exp\!\left[-\frac{E_j-E_i}{k_{\mathrm B}T}\right], $$ so the rates reproduce the equilibrium Boltzmann population ratio. This is the connection between equilibrium statistical mechanics and thermally activated kinetics. For an ideal, non-interacting, uniaxial single-domain grain at zero field, the simplest estimate is $$ \Delta E\approx KV, \qquad \tau=\tau_0\exp\!\left(\frac{KV}{k_{\mathrm B}T}\right), $$ where $K$ is the effective anisotropy energy density and $V$ is grain volume. Real rock-magnetic barriers can also depend on shape, stress, magnetostatic energy, field direction, defects, interactions, and nonuniform reversal. A field can tilt the landscape and lower one directional barrier while raising the reverse barrier. A blocking temperature is defined by comparing the relaxation time with the observation time $t_{\mathrm{obs}}$: $$ \tau(T_{\mathrm B})\approx t_{\mathrm{obs}}, \qquad \boxed{T_{\mathrm B}\approx \frac{\Delta E} {k_{\mathrm B}\ln(t_{\mathrm{obs}}/\tau_0)}}. $$ “Blocked” is not an intrinsic yes-or-no property. A state is blocked relative to an observation time. The same grain may be stable during a one-second measurement, relax during a day, and be unstable over a million years. The whole chain can be remembered as $$ \boxed{ \text{multiplicity} \rightarrow S=k_{\mathrm B}\ln\Omega \rightarrow p_i\propto e^{-E_i/k_{\mathrm B}T} \rightarrow k_{ij}\propto e^{-\Delta E_{ij}/k_{\mathrm B}T} }. $$ Boltzmann weights describe **where the system spends its equilibrium time**; Arrhenius rates describe **how quickly it moves between metastable states**; Metropolis Monte Carlo is an algorithm designed to sample the first and is not automatically a clock for the second. ### Disconnectivity graphs A disconnectivity graph summarizes an energy landscape. Leaves represent local minima. Branch points show the energy level at which basins first connect. Deep, well-separated branches indicate families of states separated by large barriers; shallow branching indicates easy interconversion. The graph is only as complete as the states and pathways that were found. Missing minima or saddles can change its topology. ### Transition matrices and continuous-time Markov chains Once rates among states are known, place them in a rate matrix **Q**: - off-diagonal **Qij** is the rate from **i** to **j**; - diagonal **Qii** is minus the total rate leaving **i**. The state-probability vector evolves as **dp/dt = pQ** for a row-vector convention. This coarse-grains fast micromagnetic motion into jumps among long-lived minima. It makes long times accessible without integrating nanosecond precession for geological durations. The Markov assumption requires that the system loses memory within each basin before the next transition and that the selected states and rates capture the important dynamics. ## Slides 49–55: Reading the MERRILL hands-on examples ### What the examples are designed to teach The sequence from random initial magnetization to relaxed vortex and then field-driven behavior illustrates a complete workflow: 1. Build or import a finite-element mesh. 2. Assign material parameters and crystal axes. 3. Choose an initial magnetization. 4. Minimize the total energy. 5. Inspect the state using magnetization vectors, helicity, cores, surfaces, or total moment. 6. Change field or another control parameter. 7. Reuse the previous relaxed state to preserve history. 8. Compare different field orientations and check convergence. ### Why use several visualizations of one state? A vector plot can become unreadable in a dense three-dimensional mesh. Derived quantities highlight particular structure: - **magnetization arrows** show local direction; - **helicity or curl-related measures** highlight rotation and vortex structure; - **isosurfaces** identify a core or region meeting a threshold; - **total moment** compresses the state to the bulk quantity measured in a magnetometer; - **stray-field maps** show what an external sensor would observe. No single visualization is the state itself. It is one projection of a large numerical object. ### Field orientation matters Applying field along a cube edge and along a body diagonal changes the symmetry between the field, crystal anisotropy, and geometry. Switching fields and reversal paths can differ even for the same material parameters. In a natural grain, crystal axes may not align with morphological axes, so both orientations must be specified. ### A practical convergence checklist Before treating a result as physical, repeat a representative calculation while changing: - mesh size and mesh quality; - minimizer tolerance and maximum iterations; - field-step size; - initial condition; - solver or damping parameter; - simulation box or magnetostatic boundary treatment; - small geometric perturbations. Report which conclusions are robust. A switching field that moves slightly may still support the same reversal mechanism. A vortex that disappears under modest refinement is not a robust prediction. ## Slides 56–57: Natural grains and stray-field modeling ### A multidomain micrometre-scale grain A one-micrometre magnetite or titanomagnetite grain can support several domains, walls, vortices, closure structures, and surface-related features. Oxidation can create a shell with different **Mₛ**, anisotropy, exchange, or stress. The shell may change the internal state even if its own moment is small. Modeling a natural grain therefore requires choices about: - the outer three-dimensional morphology; - internal phase boundaries and oxidation depth; - spatially varying material constants; - exchange coupling across interfaces; - crystal orientation and stress; - defects or pinning features; - neighboring magnetic sources. A visually similar domain image is not sufficient validation. Compare multiple observables: total moment, switching fields, field dependence, temperature dependence, and external stray field. ### Stray field is a forward problem Given **M(r)**, magnetostatics predicts the field above the specimen. That is a forward calculation. Comparing with microscopy requires further filtering by sensor height, orientation, active area, and measured component. Recovering **M(r)** from one measured field map is an inverse problem and is generally non-unique. Fine internal structures decay rapidly with height, and silent or nearly flux-closed magnetization can produce little external signal. Geometry, mineral maps, multiple heights, multiple field components, and micromagnetic priors reduce the set of plausible sources. ## Slide 58: Caveats and a standard of evidence Micromagnetic models approximate reality in three distinct ways. ### Physical approximation The continuum model omits atomically resolved chemistry and usually treats material parameters as known, local, and homogeneous within each region. Thermal fluctuations, stress, defects, and interfacial physics may be simplified or absent. ### Numerical approximation The geometry is meshed, derivatives are discretized, long-range fields are approximated numerically, and a solver stops at finite tolerance. Mesh and solver tests quantify this part of uncertainty. ### Interpretive approximation A calculation usually represents one or a few grains, whereas an experiment may contain a broad ensemble. Agreement can be non-unique: different combinations of size, shape, composition, and interactions may reproduce the same bulk curve. ### A result is strongest when it passes four tests 1. **Convergence:** the conclusion survives mesh and solver refinement. 2. **Robustness:** the conclusion survives plausible uncertainty in geometry and material parameters. 3. **Validation:** the model reproduces an observable not used merely to choose its parameters. 4. **Discrimination:** competing physical explanations were simulated and shown to predict measurably different outcomes. ### A succinct closing explanation > Micromagnetics does not produce a photograph of what must be inside a grain. It produces physically constrained possibilities. We trust a conclusion when it is numerically converged, robust to uncertain inputs, and able to explain several independent observations better than plausible alternatives. ## Questions to ask students during the lecture - Which energy term is trying to make the magnetization uniform, and which one often favors flux closure? - Why can a vortex have lower total energy than a uniform state even though it has higher exchange energy? - What determines whether a state is stable for a second or for a billion years? - Why might two minimizations of the same grain give different answers? - Which material parameters must change when simulating a temperature-dependent experiment? - How would increasing sensor height change a simulated stray-field map? - What would you vary to test whether a predicted switching field is numerical or physical? - Why is a computed global minimum not necessarily the state measured in a hysteresis experiment? ## Short explanations worth memorizing **Atomistic versus micromagnetic:** Atomistic models keep individual lattice sites and spins; micromagnetics averages many spins into a continuous direction field. **Exchange:** Exchange penalizes rapid changes in magnetization direction, so it supplies magnetic stiffness. **Anisotropy:** Anisotropy makes some directions cheaper than others because spin is coupled to the crystal or shape. **Demagnetizing energy:** Demagnetizing energy penalizes magnetic poles and often drives the magnetization toward flux closure. **Hysteresis:** Hysteresis occurs because the system follows local minima that depend on its history, then jumps when a minimum loses stability or a barrier is crossed. **Thermal relaxation:** Temperature does not simply weaken magnetization; it also gives the system chances to cross energy barriers. **Finite-element mesh:** The mesh is the spatial language in which the model can express magnetic structure. If the language is too coarse, the model cannot describe the correct state. ## Further reading - Brown, W. F. (1963), *Micromagnetics*, the foundational continuum treatment. - Hubert, A. & Schäfer, R. (1998), *Magnetic Domains*, for energy terms, walls, domains, and imaging. - Ó Conbhuí et al. (2018), MERRILL finite-element micromagnetics: - Fabian & Shcherbakov (2018), energy barriers and thermoviscous magnetization: - Nagy et al. (2017), stability and paleomagnetic significance of magnetite vortex states: - Harrison et al. (2025), three-dimensional magnetic structure in a giant magnetofossil: