Weiss model of spontaneous magnetization

from Statistical Physics using Mathematica
© James J. Kelly, 1996-2002

Introduction

As the temperature is reduced below a critical temperature (Subscript[T, c]), many systems experience a continuous phase transition characterized by the sudden onset of long-range order.  Order-disorder transitions are usually described in terms of an order parameter, m, and an ordering field, h.  A fundamental problem in the study of continuous phase transitions is to formulate a qualitative description of thermodynamic properties that are either discontinuous or divergent at the critical point.  It is customary to represent such properties in terms of critical exponents.  When the ordering field vanishes, the order parameter typically behaves as

\(m[T, 
        0] ∝ \((T\_c - T)\)\^β\ \ \ \ \ \ \ T ≲ T\_c\)
m[T, 0] = 
      0                       T ≥ 
        (Subscript[T, c])

At the critical temperature, the order parameter typically scales with ordering field according to

m[(Subscript[T, c]), 
        h]∝(h^(1/δ))

Similarly, the response function χ = ∂(m/∂)h typically diverges at the critical temperature according to

\(χ[T, 
        h] ∝ \((T\_c - T)\)\^\(-\(γ\_ < \)\)\ \ \ \ \ \ \ \ \ \ \ T ≲ T\_c\)
\(χ[T, 
        h] ∝ \((T - T\_c)\)\^\(-\(γ\_ > \)\)\ \ \ \ \ \ \ \ \ \ \ T ≳ T\_c\)

Finally, the heat capacity for vanishing field is often discontinuous for an order-disorder phase transition.  In this notebook we study a mean-field model for spontaneous magnetization in paramagnetic and ferromagnetic systems, for which the average magnetization serves as the order parameter and the external magnetic field as the ordering field.  Experimentally, one finds that the critical exponents are β ≈ 0.30 - 0.36, δ ≈ 4.2 - 4.8, \(TraditionalForm\`\(γ\_ < \) ≈ 1.0 - 1.2\), and \(TraditionalForm\`\(γ\_ > \) ≈ 1.2 - 1.4\).

The magnetic contribution to the Hamiltonian for a system of spin-(1/2) particles in an external magnetic field takes the form

\(H = \(1\/2\) \(∑\+\(i, j\)\(ε\_\(i, 
                  j\)\) σ\&⇀\_i·σ\&⇀\_j\)\ \  - \ μ \(∑\+i σ\&⇀\_i·B\&⇀\)\)

where μ is the atomic magnetic moment and where the spin-spin interaction (Subscript[ε, i, j]) depends upon separation.  The spin-spin interaction range is typically small enough that one can restrict the summation to nearest neighbors and take (Subscript[ε, i, j]) to be constant, obtaining thereby the Heisenberg model.  Often it is sufficient to replace the scalar product by the single term \(TraditionalForm\`σ\&⇀\_i·σ\&⇀\_j -> \(σ\_i\) σ\_j\), where (Subscript[σ, i]) = ±1 is twice the z-component of the spin for atom i; this approximation is known as the Ising model.  Finally, the Weiss model assumes that the tendency for a particular spin to be found in the up orientation is governed by the average magnetization of the entire sample, with essentially random fluctuations among nearby spins cancelling.  Thus, we write

\(H ≈ \(-\(∑\+i\((\(q\/2\) ε  σ\&_ + μ B)\) σ\_i\ \ ⟹\ \ U\)\) ≈ \(-\((\(q\/2\) ε  σ\&_ + μ B)\)\) σ\&_\)

where the order parameter \(TraditionalForm\`σ\&_\) is twice the average spin per atom, q is the coordination number representing the number of nearest neighbors, and ε is a spin-spin coupling parameter whose sign is positive for ferromagnetism.  The factor of (1/2)in the spin-spin contribution prevents double counting for pair-wise interactions.  Although the spins are correlated through the constraint upon the average magnetization, we can focus upon a single atom and treat the remaining system as a reservoir which establishes the temperature and mean field for this system.  The probability for finding any particular atom within this mean field in the orientation σ is then determined by the Boltzmann factor

\(P[σ]\  ∝ \ 
      Exp[\(-β\) \((μ B + 
              q ε  σ\&_)\) σ]\)

where the double-counting factor of (1/2) is not needed in the energy for a single atom interacting with the mean field.  Thus, we obtain an equation for the average magnetization that is very similar to that for paramagnetism, except that the spin orientation must satisfy a self-consistency relationship of the form

\(σ\&_ = 
      Tanh[β \((μ B + 
              q ε  σ\&_)\)]\)

obtained by the replacement \(TraditionalForm\`μ B -> μ B +  q ε  σ\&_\).

In the absence of an external magnetic field, the self-consistency equation reduces to

\(σ\&_ = 
      Tanh[β q ε  σ\&_]\)

which can be solved graphically by looking for intersections between the line \(TraditionalForm\`y = σ\&_\) and the curve \(TraditionalForm\`y = 
        Tanh[β q ε  σ\&_]\) that represents states with spontaneous magnetization.  A typical graphical solution is shown in the figure below, which demonstrates that there are no nontrivial solutions unless \(TraditionalForm\`β q ε > 1\).  Therefore, the critical temperature is (Subscript[T, c]) = q ε/(Subscript[k, B]).

It is useful to express the self-consistency condition in the form

\(σ\&_ = Tanh[\(b + σ\&_\)\/τ]\)

where τ = T/(Subscript[T, c]) = (Subscript[k, B]) 
            T/q ε is the temperature relative to the critical temperature and b = μ B/q ε describes the relative importance of external and internal magnetic energies.  Similarly, it is convenient to measure energies in units of \(TraditionalForm\`𝒩 q ε\) and entropies in units of 𝒩 (Subscript[k, B]), where 𝒩 is the number of atoms, to obtain dimensionless expressions.  In the absence of an external field the two nontrivial solutions are degenerate and the orientation of spontaneous magnetization is arbitrary.  The presence of an external field breaks the degeneracy

and for sufficiently large b only a single solution dominated by the external field remains

For ferromagnetism with ε > 0 it is sufficient to consider only the τ > 0, b > 0 sector, but antiferromagnetic models with ε < 0 (here obtained using τ < 0) are also of interest.

In the following sections we develop the thermodynamic properties of the Weiss model in more detail.  Computation of the heat capacity and magnetic susceptibility are complicated by the need for numerical solution of the self-consistency condition, but can be accomplished with relative ease using Mathematica.

Initialization

Defaults,  packages, and symbols

ClearAll["Global`*"];
Off[General::spell, General::spell1]

$DefaultFont = {"Times", 12};
$TextStyle = {FontFamily -> "Times", FontSize -> 12, 
          FontSlant -> "Italic"};

Needs["Utilities`Notation`"];
Needs["Graphics`Master`"]

Symbolize[(Subscript[k, B])];Symbolize[(Subscript[C, B])];Symbolize[\(σ\&_\)]

Numerical solution of self-consistency condition

The self-consistency condition can be expressed in the form

eq[τ_, b_] := m == Tanh[((b + m)/τ)]

If we limit ourselves to ferromagnetism, the following function provides the relevant solution for the magnetic alignment and provides protection against numerical errors for extreme values.

\(σ\&_[τ_?NumericQ, b_?NumericQ] := 
      Module[{y}, 
        y = Which[τ == 0, 1, b == 0, 
            If[τ < 1, m /. FindRoot[Evaluate[eq[τ, 0]], {m, 0.5}], 
              0], True, m /. FindRoot[Evaluate[eq[τ, b]], {m, 0.5}]]; 
        Which[y ≤ 0, 0, y ≥ 1, 1, True, 
          y]]\)

Simplification rules

stirling  =  Log[x_!] -> x Log[x] - x;

It will be useful to develop a set of rules for using the self-consistency equation to simplify expressions.

\(trigrules = {
Tanh[\(b + σ\&_[τ, b]\)\/τ] -> σ\&_[τ, 
              b], 
Coth[\(b + σ\&_[τ, b]\)\/τ] -> 
            1\/σ\&_[τ, b], 
Cosh[\(b + σ\&_[τ, b]\)\/τ] -> 
            Cosh[ArcTanh[σ\&_[τ, 
                  b]]], 
Sinh[\(b + σ\&_[τ, b]\)\/τ] -> 
            Sinh[ArcTanh[σ\&_[τ, 
                  b]]], 
Sech[\(b + σ\&_[τ, b]\)\/τ] -> 
            Sech[ArcTanh[σ\&_[τ, 
                  b]]], 
Csch[\(b + σ\&_[τ, b]\)\/τ] -> 
            Csch[ArcTanh[σ\&_[τ, b]]]
} // 
        Simplify\)

Derivatives of \(TraditionalForm\`σ\&_\) can now be expressed in algebraic form as follows.

\(eq2[τ_, b_] := σ\&_[τ, b] == 
        Tanh[\(b + σ\&_[τ, b]\)\/τ]\)

D[eq2[τ, b], τ] /. trigrules

σrules=Join[Solve[D[eq2[τ, b], τ] /. trigrules,Derivative[1, 0][\(σ\&_\)][τ, b]][[1]],Solve[D[eq2[τ, b], b] /. trigrules,Derivative[0, 1][\(σ\&_\)][τ, b]][[1]]]//FullSimplify

Thermodynamic functions

Magnetization

The figure below demonstrates that the magnetization of the sample approaches saturation either for strong fields or low temperatures.  Experimentally one finds that this model provides a good description of the magnetization of most ferromagnetic materials, although it is difficult to calculate ε or (Subscript[T, c]) from first principles.

\(\(Plot3D[σ\&_[τ, b], {τ, 0, 5}, {b, 0, 5}, 
        PlotPoints -> 50, PlotLabel -> "magnetization", 
        AxesLabel -> {\*""T/(Subscript[T, c])"", \*""μB/(Subscript[k, B])(Subscript[T, c])"", \*""\!\(σ\&_\)""}];\)\)

For sufficiently low temperatures, the spin-spin interaction produces spontaneous magnetization in the absence of an external magnetic field.  This internal order decreases with temperature and vanishes at and above the critical temperature.

\(\(plot1 = 
        Plot[σ\&_[τ, 0], {τ, 0, 1}, 
          PlotRange -> {{0, 1.2}, {0, 1.2}}, Frame -> True, 
          FrameLabel -> {\*""T/(Subscript[T, c])"", \*""\!\(σ\&_\)""}, PlotLabel -> "order parameter", 
          Epilog -> Text["B=0", {1, 1}]];\)\)

The magnetic susceptibility

\(χ[τ_, 
        b_] = ∂\_b σ\&_[τ, b] /. σrules // 
        Simplify\)

diverges at τ = 1 when b -> 0 and flattens out as the external field increases.  The peak in the magnetic susceptibility decreases in amplitude and moves to higher temperatures as b increases.

\(\(Show[
        GraphicsArray[{{Plot[χ[τ, 0], {τ, 0, 2}, 
                PlotRange -> {All, {0, 20}}, Frame -> True, 
                FrameLabel -> {\*""T/(Subscript[T, c])"", "χ/Nμ"}, PlotLabel -> "magnetic susceptibility", 
                Epilog -> Text["B = 0", {0.2, 15}], 
                DisplayFunction -> Identity], 
              Plot[Evaluate[
                  Table[χ[τ, b], {b, 0.1, 1, 0.1}]], {τ, 0, 
                  6}, PlotRange -> All, Frame -> True, 
                FrameLabel -> {\*""T/(Subscript[T, c])"", "χ/Nμ"}, PlotLabel -> "magnetic susceptibility", 
                Epilog -> Text["B = {0.1,1,0.1}", {4, 2}], 
                DisplayFunction -> Identity]}}, 
          GraphicsSpacing -> {\(-0.1\), 0.0}]];\)\)

In the immediate vicinity of the critical temperature, spontaneous magnetization has the following behavior.

\(rule1 = \(Solve[σ\&_ == 
              Normal[Series[
                  Tanh[σ\&_\/τ], {σ\&_, 0, 
                    3}]], σ\&_]\)[[3]] // Simplify[#, 0 < τ < 1] &\)

Thus, the mean-field model predicts a critical exponent (Subscript[β, mf]) = (1/2) that should be compared with experimental values typically in the range 0.30 - 0.36.  Similarly, the magnetization at the critical temperature can be obtained by expressing the self-consistency equation in the form

\(τ = 
      1\ \ ⟹\ \ ArcTanh[σ\&_] == 
        b + σ\&_\)

and expanding the lhs around \(TraditionalForm\`σ\&_ = 0\) to obtain

\(Solve[
      b + σ\&_ == 
        Normal[Series[
            ArcTanh[σ\&_], {σ\&_, 0, 
              3}]], σ\&_]\)

Thus, we find  \(TraditionalForm\`σ\&_ ≈ \((3  b)\)\^\(1/3\)\ ⟹\ \ χ ≈ \((3  b)\)\^\(\(-2\)/3\)\) and obtain another critical exponent (Subscript[δ, mf]) = 3 for comparison with experimental values in the range 4.2 - 4.8.  Although neither prediction is perfect, the mean field gives results for the magnetization near the critical point that are qualitatively correct.

For temperatures well above critical, the magnetization

\(τlarge = \(Solve[σ\&_ == 
            Normal[Series[
                Tanh[\(b + σ\&_\)\/τ], {τ, ∞, 
                  1}]], σ\&_]\)[[1]]\)

reduces to the Curie-Weiss law, which is very similar to the corresponding result for ideal paramagnetism except that the denominator is replaced by τ -> τ - 1.  In the opposite extreme of very low temperatures

\(Tanh[\(b + σ\&_\)\/τ] = \(1 - Exp[\(-2\) \(b + σ\&_\)\/τ]\)\/\(1 + Exp[\(-2\) \(b + σ\&_\)\/τ]\) ≈ 
        1 - 2  Exp[\(-2\) \(b + σ\&_\)\/τ]\ \ ⟹\ \ σ\&_ ≈ 
        1 - 2  Exp[\(-2\) \(b + 1\)\/τ]\)

\(\(τsmall = σ\&_ -> 
          1 - 2  Exp[\(-2\) \(b + 1\)\/τ];\)\)

we obtain an exponential behavior that is characteristic of systems with a finite energy gap, but is not entirely realistic.  Just as the Einstein model produces exponential while the more realistic Debye model of lattice vibrations produces power law low-temperature behavior, collective excitations of the spin field lead to power-law behavior in ferromagnetic systems also.  The Weiss model describes the long-range order but neglects short-range correlations between neighboring spins.  Consider a classical model of the magnetic moments.  If one spin is tilted slightly, the spin-spin interaction will tend to tilt the neighbors, while the interaction with the mean field tends to restore those spins to the original direction.  Hence, there is the possibility of propagating waves.  The energy spectrum for these collective excitations is continuous; hence, the magnetization is expected to have power-law behavior at low temperatures.  Quantized magnetic waves are known as magnons.

Energy

The average energy per atom is simply the mean value of the model Hamiltonian

\(\(U[τ_, 
          b_] = \(-\((b + \(1\/2\) σ\&_[τ, 
                    b])\)\) σ\&_[τ, b];\)\)

such that the specific heat (in units of (Subscript[k, B])) becomes

(Subscript[C, B])[τ_, b_] = (∂_τ) U[τ, b] /. σrules // 
        Simplify

\(\(Plot3D[C\_B[τ, b], {τ, 10\^\(-4\), 5}, {b, 0, 5}, 
        PlotRange -> All, PlotPoints -> 50, 
        AxesLabel -> {\*""T/(Subscript[T, c])"", \*""μB/(Subscript[k, B])(Subscript[T, c])"", \*""C/(Subscript[Nk, B])""}, 
        PlotLabel -> "magnetic heat capacity"];\)\)

In the absence of an external magnetic field, the magnetic heat capacity has a finite discontinuity at the critical temperature characteristic of second-order phase transitions.  The heat capacity vanishes at low temperature because the system is completely ordered and above the critical temperature because it is completely disordered.

\(\(plot2 = 
        Plot[C\_B[τ, 0], {τ, 0, 1.2}, PlotRange -> All, 
          Frame -> True, 
          FrameLabel -> {\*""T/(Subscript[T, c])"", \*""C/(Subscript[Nk, B])""}, PlotLabel -> "magnetic heat capacity", 
          Epilog -> Text["B = 0", {0.2, 1.2}]];\)\)

Expanding the zero-field heat capacity around the critical temperature

\(Series[\(C\_B[τ, 
              b] /. {σ\&_[τ, b] -> σ\&_, 
              b -> 0}\) /. rule1, {τ, 1, 1}] // Simplify\)

but the nonmagnetic contributions would have to be included before comparing with experimental measurements.  The heat capacity for very low temperatures has the exponential behavior characteristic of systems with a gap in the excitation spectrum, but the magnon contribution is expected to provide a power-law temperature dependence.

Entropy

Assuming that the spins are arranged randomly, constrained only by the net magnetization, the specific entropy is related to the multiplicity by

\(S[τ_, 
        b_] = \(\(\((\(1\/𝒩\) 
                    Log[\(𝒩!\)\/\(\(\(N\_+\)!\) \(\(N\_-\)!\)\)] // PowerExpand)\) /. 
              stirling\) /. {\(N\_+\) -> \(𝒩\/2\) \((1 + σ\&_[τ, 
                        b])\), \(N\_-\) -> \(𝒩\/2\) \((1 - σ\&_[τ, b])\)} // PowerExpand\) // Simplify\)

In the high temperature limit where spins are found in either orientation with equal probability, the entropy approaches 𝒩 (Subscript[k, B]) Log[2].  In the absence of an external magnetic field, the entropy saturates at the critical temperature and remains constant thereafter.  Hence, the phase transition for spontaneous magnetization is continuous, but the discontinuity in the derivative of entropy with respect to temperature is related to the finite discontinuity in specific heat.  This phase transition is classified as second-order because the second-derivative of free energy with respect to temperature is discontinuous.  However, the presence of an external magnetic field removes this discontinuity and obscures the order-disorder phase transition.  As B increases, tending to order the system, higher temperatures are required to saturate the entropy.

\(\(Plot[Evaluate[Table[S[τ, b], {b, 0, 1, 0.2}]], {τ, 0, 5}, 
        PlotRange -> All, Frame -> True, 
        FrameLabel -> {\*""T/(Subscript[T, c])"", \*""S/(Subscript[Nk, B])""}, PlotLabel -> "magnetic entropy", 
        Epilog -> Text["B = {0,1,0.2}", {3, 0.2}]];\)\)

Ginzburg-Landau description of spontaneous symmetry breaking

The Ginzburg-Landau theory of phase transitions treats the free energy as a thermodynamic potential in which the order parameter serves as the dynamical variable.  The equilibrium state of the system is found by minimizing the free energy with respect to the order parameter.  The internal energy for the Weiss model can be expressed as a function of magnetization as

U = -(b + (m/2)) m

while the combinatoric calculation gives entropy

S = (((1 + m))/2) Log[(2/(1 + m))] + (((1 - m))/2) 
          Log[(2/(1 - m))]

Thus, the free energy becomes

Clear[F]

F[τ_, 
          b_] = -(b + (m/2)) 
            m - τ ((((1 + m))/2) 
                  Log[(2/(1 + m))] + (((1 - m))/2) 
                  Log[(2/(1 - m))]);

It is important to recognize that, unlike our earlier formulation, these energy and entropy functionals are not limited to equilibrium states, but may be evaluated for any state in which the magnetization is uniform.  In equilibrium the system will be found with the magnetization \(TraditionalForm\`m\&_\) which minimizes the free energy.  It is instructive to visualize the dependence of free energy upon spontaneous magnetization (b -> 0), for several temperatures.

plotF[τ_, b_, min_, max_, opts___Rule] := 
      Plot[F[τ, b], {m, min, max}, opts, Frame -> True, 
        FrameLabel -> {"m", "F"}, 
        PlotLabel -> "τ = " <> 
            ToString[τ] <> ",   b = " <> ToString[b]]

Show[
        GraphicsArray[{
{plotF[2, 0, -1, 1, 
                DisplayFunction -> Identity], 
              plotF[0.95, 0, -1, 1, 
                PlotRange -> {Automatic, {-0.662, -0.65}}, 
                DisplayFunction -> Identity]}, 
{plotF[
                1.2, 0, -1, 1, DisplayFunction -> Identity], 
              plotF[0.8, 0, -1, 1, 
                DisplayFunction -> Identity]}, 
{plotF[
                1, 0, -1, 1, DisplayFunction -> Identity], 
              plotF[0.5, 0, -1, 1, DisplayFunction -> Identity]}}], 
        GraphicsSpacing -> {-0.1, 0}];

For temperatures well above critical the free energy displays a qualitatively parabolic dependence about a minimum at \(TraditionalForm\`m\&_ = 0\).  The bottom of the potential well  flattens as the temperature decreases until both first and second derivatives vanish at the critical temperature.  Although there is still a minimum at \(TraditionalForm\`m\&_ = 0\), one expects very strong fluctuations because the potential energy is so flat.  Strong fluctuations in the order parameter are a characteristic signature of a phase transition.  As the temperature is reduced below the critical temperature, two minima symmetrically placed about 0 represent spontaneous magnetization.  Either minimum can be chosen in the absence of an external field, and one expects a random choice determined by fluctuations or by stray fields that may not have been eliminated perfectly.  The degree of spontaneous magnetization increases as the sample is cooled, approaching saturation as T -> 0.  Although there is always a minimum with |m| slightly less than unity for finite temperatures, it can be challenging to find numerically.

In the absence of an external field, the free energy is independent of the sign of the magnetization

F[τ, 0] == (F[τ, 0] /. m -> -m)

More generally, if we had kept the full three-dimensional dependence of energy upon spin, we would have found that the free energy is independent of the orientation of the magnetization vector.  Nevertheless, as the sample is cooled through the phase transition, spontaneous magnetization with a definite direction in space develops even though the energy functional remains rotationally invariant.  This phenomenon is known as spontaneous symmetry breaking — the equilibrium state does not share the symmetry of the Hamiltonian.  One often finds that a continuous phase transition occurs as long-range order minimizing the free energy spontaneously breaks a symmetry of the underlying dynamics.  Rotational invariance is broken by spontaneous magnetization, rotational and translational invariance for crystallization, and gauge invariance for superfluidity.  This phenomenon also plays fundamental roles in quantum field theory and particle physics.

External fields that couple linearly to the order parameter also break symmetry and hence destroy the continuous phase transition.  The figure below illustrates the effect of an external field for a temperature just below critical.  Although it is difficult to discern the symmetric minima in the curve for b = 0 on this scale, we know that such minima occur at the dashed lines.  Increasing the external field skews the free-energy curves and shifts the minimum with the same sign toward larger magnetization.

\(\(Plot[Evaluate[Table[F[0.9, b], {b, 0, 0.5, 0.1}]], {m, \(-1\), 1}, 
        Frame -> True, FrameLabel -> {"m", "F"}, 
        PlotLabel -> "τ=0.9, {b,0,0.5}", 
        Epilog -> {Dashing[{0.02, 0.02}], 
            Line[{{σ\&_[0.9, 0], \(-1\)}, {σ\&_[0.9, 0], 0}}], 
            Line[{{\(-σ\&_[0.9, 0]\), \(-1\)}, {\(-σ\&_[0.9, 
                      0]\), 0}}]}];\)\)

For the sake of completeness, we demonstrate that this approach provides the same magnetization equation as the self-consistency approach above.  Minimizing the free energy gives

eq3[τ_, b_] = 0 == D[F[τ, b], m] // Simplify

The simplest way to compare this equation with our earlier version is to solve both for τ and to compare temperatures for given m and b.

τ1 = τ /. Solve[
            eq[τ, 
              b], τ][[1]]

τ2 = τ /. Solve[
              eq3[τ, 
                b], τ][[1]] // 
        Simplify

FullSimplify[τ1 - τ2, 0 < m < 1]

Therefore, both methods give the same dependence of equilibrium magnetization upon τ and b.