# Questions on mean field theory and the enforcement of local spin conserving constraints

Hello Physics StackExchange community,

I've recently been working on a problem that seems like it should be straightforward, but I can't seem to overcome what seems to be a basic obstacle. I'll try to be clear and lay things out from the beginning.

We start with an Ising model of $${N}$$ interacting spins sitting in an external field. At each lattice site $${i}$$, the associated spin can take on two discrete values ($${s_{i}=\pm1}$$), corresponding to up (+) or down (-) spin (we only care about spin/fields in the $${z}$$ direction).

We have the Hamiltonian

$${H=-\frac{1}{2} \sum_{i,j\ne i}^{N,N} J_{ij}s_i s_j - \sum_{i}^{N} h_i s_i, } \qquad (1)$$ where $${J_{ij}}$$ is the interaction strength between spins $${s_i}$$ and $${s_j}$$. For ($${J_{ij}>0}$$) we have an Ising ferromagnet (FM) which favors parallel spins, where $${s_i}$$ and $${s_j}$$ have the same sign. For ($${J_{ij}<0}$$) we've got an antiferromagnet (AFM) where antialignment is preferred. And our second term is just the interactions between the local field $${h_i}$$ and $${s_i}$$.

If we take a mean-field approximation, we obtain the self-consistency equation for the local average spin/magnetization $${m_i}$$ in this system. And we find that for a given site

$${ m_i = \bar s_i = \tanh \left[ \frac{1}{k_BT}\left( \sum_{j \ne i}^{N}J_{ij}m_j + h_i \right) \right], \qquad (2)}$$

in the most general case.

Now this is all fine and good. Right now I'm just playing around with some basic 2D models, so we can apply this to the classical example of an Ising AFM on a square lattice -- to be divided into two square sublattices $${A}$$ and $${B}$$. The $${A}$$ sites have 4 $${B}$$ sites as NNs, and vice versa. This is shown in the figure below In this AFM case, a simplifying assumption would be to take $${J_{ij} = -J}$$, where $${J>0}$$, and consider NN interactions only. Further, we'll assume the external field is uniform so that $${h_i = h}$$ for all sites. So the sublattices have their own basic square crystal geometry/translational invariance, and using $${(2)}$$, we generate two equations -- one for sublattice $${A}$$ and another for $${B}$$, which are coupled and summarized in the equation below.

$${m_{A,B}= \tanh \left[ -\frac{4J}{k_BT} m_{B,A}+\frac{h}{k_BT} \right] \qquad (3)}$$

And this is still all fine and dandy for me. I see that, if I just assign values to the parameters $${\{J,k_BT,h\}}$$ in $${(3)}$$, I can solve the transcendental equation for $${m_A}$$

$${m_{A}= \tanh \left[ -\frac{4J}{k_BT} m_{B}+\frac{h}{k_BT} \right] \\\quad = \tanh \left[ -\frac{4J}{k_BT} \left( \tanh \left[ -\frac{4J}{k_BT} m_{A}+\frac{h}{k_BT} \right] \right)+\frac{h}{k_BT} \right] \qquad (4)}$$

numerically by finding the intersections of the left- and right- hand sides of $${(4)}$$... i.e. we find where we have crossings between the two functions in the following plot. When we have three intersections/solutions, as in the above, this is the AFM phase, where we can look at the free energy, corresponding to the left-most and right-most solutions, and see that these are local minima, while the middle solution corresponds to a local maximum. When we only have one solution, this corresponds to the paramagnetic phase, where the previously described contour collapses into a single free energy minimum. And we can monitor the phase transition by fixing $${\{J,h\}}$$ and looking at the resulting $${m_A}$$ as we tune through $${k_BT}$$.

For example, at zero field $${h=0}$$, we can treat $${\frac{k_BT}{J}}$$ as the only variable, and tune this while monitoring $${m_A}$$. This will look like (we plot the 2/3 available solutions for $${m_A}$$, corresponding to the 2 local free energy minima, on the AFM side)

We can make a similar plot, but of the 2 $${m_B}$$ solutions/free-energy minima as a function of $${\frac{k_BT}{J}}$$. And if we do, we'll see that it will actually be identical, but the corresponding curves will be flipped.

In other words, if I take above figure's upper $${m_A}$$ curve and plot it alongside its corresponding $${m_B}$$ (which is just calculated through $${(3)}$$'s $${m_B}$$ equation) as a function of $${\frac{k_BT}{J}}$$, we'll see that Doing the same for the other $${m_A}$$ solution (which is just the other free-energy minimum; we're ignoring the third AFM solution, which is an unstable free energy maximum) and its corresponding $${m_B}$$ curve, we'll see the same thing but flipped. Both of these cases are therefore equivalent solutions, because the only difference between them is that A and B lattices switch (but these sublattices are just equal with opposite spin magnitudes...swapping them would be like turning the system upside down, which shouldn't change any of the physics)

With this in mind, we don't have to worry about the other solution, which is redundant, and can pick out and focus on the one we have -- the set of $${\{m_A,m_B\}}$$ in the most recent figure for the zero-field case. This just corresponds to $${m_A}$$ having all the spin up at 0 temperature (while the other solution means $${m_A}$$ starts with the spin down)

We clearly see that $${m_A = -m_B}$$, which is a result of the AFM interaction $${-J<0}$$ and the symmetry/equivalence between the sublattices (they're geometrically equivalent and neither sublattice is special, the way this problem is constructed at zero field). And as we increase the temperature, we're going to allow spins to excite and start swapping places, right up to the transition, where the solutions converge to a single (paramagnetic) solution. In this case, it's just 0 -- which checks out because at $${k_BT=0}$$, we start with say $${N/2}$$ spin ups on one sublattice (A in this example, according to the figure) and an equal number of spin downs on an equivalent sublattice (B). So when the thermal energy jumbles the spins and ruins the AFM ordering, each site on either sublattices A or B will have an equal number of either spins up or down available to them, which nets to 0.

Finishing up our discussion of the zero-field case, we can now start turning on and tuning through a nonzero field. Now we should treat $${\{k_BT,J,h\}}$$ as independent parameters. But fixing $${J=1}$$, and tuning through $${h}$$ we see that So in these cases, the paramagnetic solution doesn't drop to 0 right away, but is persistently finite (commensurately with the field strength) until the temperature becomes sufficiently high. This makes sense because when we have a (positive-valued) field, then our previous spin symmetry is broken and one of the spins (up) become favored. Since, for example, all the A sites have up and B sites are down at 0 temperature, but finite field, it's like the two sublattices experience different field-induced shifts in chemical potentials, so their occupations don't evolve as symmetrically with temperature as they did before, which explains why it takes a while for $${m_{A,B} = 0}$$ at $${T=\infty}$$.

OKAY! So I've covered quite a bit of ground introducing all of this basic mean-field stuff in the context of this Ising AFM story (feel free to correct anything I may have misstated). But I think I've set up enough background to start really asking the question that's been bothering me.

Now what I want to do is I want to extend this approach to another kind of system, where we similarly have two sublattices A and B who are no longer equivalent. In fact, there are twice as many A sites as B sites, which is depicted in the figure below. This means that we have twice as many spins up as we do spin down, and at zero temperature we have an AFM configuration of spins up on lattice A and spins down on lattice B. At infinite temperature, any spin is possible, so if we are drawing from a grabbag of 2/3 spin ups (+1) and 1/3 spin down (-1), then our average value should be 1/3. This implies that we should be satisfying

$${2m_A+m_B=1}$$ for $${\{m_A,m_B\} = \{1,-1\} \qquad T=0}$$ and $${\{m_A,m_B\} = \{1/3,1/3\} \qquad T=\infty}$$.

My question is this: how do I impose such a constraint in my mean-field approach? We saw that in the zero field example, $${m_A = - m_B}$$ was automatically satisfied for any temperature, but we lost this constraint as soon as we turned the field on.

Some more points about the problem I'm dealing with: (1) instead of just NN interactions, I'm assuming interactions from all other sites. This can be done with infinite summation techniques (2) the net magnetization of the system is clearly nonzero from having twice as many spins up as we have spins down, which causes problems for such infinite summation techniques. To remedy this, we impose a background neutralizing magnetization field to cancel this net spin per unit volume. With this, my transcendental equations in $${(2)}$$ will be of the form

$${m_{A,B} = \tanh \left[\frac{1}{k_BT} \left( C m_A +D m_B +h \right) \right]}$$

for some coefficients $${\{C,D\}}$$. Now it's clear that if I take $${T=\infty}$$, the solutions for $${m_{A,B}=0}$$ and not 1/3 as we desire. Further, my attempts to shift the solution by doing transformations of $${m_A}$$ to $${m_A-1/3}$$ or other such line equations have proved fruitless in trying to match these spin-conservation-demanding temperature boundary conditions...If I get $${\{m_A,m_B\}=1/3}$$ at infinite temperature, this screws up my values at 0 temperature.

Is there any hope for me in this mean-field problem? Or is local conservation of spin in a unit cell (i.e. +1 overall spin in a 3-site cell of 2 A sites and 1 B site) not an applicable feature in these mean-field self-consistency calculations?

I know this was a lot but if anyone has any insight into this matter of conserving the spin in these mean-field calculations, advice would be greatly appreciated. If anyone would like more details on what I've tried and/or looked at, please let me know. Thanks.