Ferromagnetic Coulomb phase in classical spin ice

Spin ice is a frustrated magnetic system that at low temperatures exhibits a Coulomb phase, a classical spin liquid with topological order and deconfined excitations. This work establishes the presence of a Coulomb phase with coexisting ferromagnetic order in a microscopic model of classical spin ice subject to uniaxial lattice distortion. General theoretical arguments are presented for the presence of such a phase, and its existence is confirmed using Monte Carlo results. This example is used to illustrate generic properties of spin liquids with magnetic order, including deconfinement of monopoles, signatures in the neutron-scattering structure factor, and critical behavior at phase transitions. An analogous phase, a superfluid with spontaneously broken particle-hole symmetry, is demonstrated in a model of hard-core lattice bosons, related to spin ice through the quantum-classical correspondence.

Figure 1

Three configurations of a single tetrahedron, and their energy $E_t$ in the nearest-neighbor Hamiltonian ${\mathcal{H}}_{\text{nn}}$, Eq. (1), with $J_{ij}$ given in Sec. 2b. Pairs of spins situated in the same horizontal plane, indicated with dashed (red) lines, have reduced coupling $J_{ij}=J-3p$, while others have $J_{ij}=J$. All three configurations obey the ice rule, having two spins pointing in and two pointing out. The first two are lowest-energy configurations for a single tetrahedron (since the antiferromagnetically aligned pairs are those with reduced coupling), while the one on the right is one of the remaining four ice-rule configurations whose energy is higher by $4p$. Excitations above the ground state are described by strings of spins flipped relative to a fully polarized configuration, and increase the energy by $4p$ per tetrahedron. When a pair of strings pass through the same tetrahedron, all four spins are inverted and the energy is again minimized; the strings therefore feel an attractive interaction.

Figure 2

Illustration of the four-spin interaction ${\mathcal{H}}_{\text{4s}}$ added to the model to stabilize the ferromagnetic Coulomb phase. The arrows represent spins on the sites of a pyrochlore lattice, a network of corner-sharing tetrahedra. (This configuration obeys the ice rule, with two spins pointing into and two pointing out of each tetrahedron.) The additional interaction couples pairs of tetrahedra at opposite sides of hexagons; one such pair and its hexagon are highlighted.

Figure 3

Magnetization vs temperature, for fixed $V_4/T=0$ (left) and $V_4/T=-0.01$ (right), and $L=24$ $(N_{\text{s}}=16L^3\simeq 2\times {10}^5$ spins). In both cases, there is a jump from saturation to zero magnetization, at a transition temperature indicated with a vertical line. For each temperature, the spontaneous magnetization is found by applying a weak field along the $z$ direction and extrapolating to zero field.

Figure 4

Magnetization vs temperature, for fixed $V_4/T=0.05$ and system size $L=24$. The vertical lines at $T/p\simeq 3.15$ and $3.35$ indicate positions of phase transitions, determined as described in Secs. 5 and 4 respectively. Below the lower-temperature transition, the magnetization remains at its saturation value, apart from small finite-size corrections, while above the higher-temperature transition, it vanishes. The intermediate phase is a ferromagnet with nonzero and continuously varying magnetization.

Figure 5

Variance of energy, $(\langle E^2\rangle -{\langle E\rangle }^2)/N_{\text{s}}$, vs temperature for fixed $V_4/T=0.05$ and various system sizes. The vertical lines indicate the positions of phase transitions in the thermodynamic limit (determined by other means). The double-peak structure, with peak heights at most weakly diverging with $L$, is consistent with a pair of continuous transitions.

Figure 6

Energy and magnetization histograms for $L=16,V_4/T=0.05$, and $T/p=3.31$ (near the higher-temperature transition). The unimodal energy distribution is indicative of a continuous transition, while the two peaks in the magnetization distribution show that this transition is associated with magnetic ordering and breaking of spin-reversal symmetry.

Figure 7

Energy and magnetization histograms for $L=8,V_4/T=-0.01$, and $T/p=6.76$. In this case, the transition is strongly first order, as indicated by the bimodal energy distribution, and occurs directly between the saturated ferromagnet $\left(M_z=\pm M_{\text{sat}}\right)$ and the paramagnet $\left(M_z=0\right)$. (A small system size is necessary to observe distributions with multiple peaks; for larger $L$ the competing states are metastable.)

Figure 8

Ratio of monopole distribution function $G_{\text{m}}$ evaluated at the maximum displacement in the lattice, ${\mathbf{R}}_{\text{max}}$, and at the maximum displacement along the $z$ direction, ${\mathbf{R}}_{\text{max,z}}$. The curves (from bottom to top) have $T/p$ corresponding to the dashed lines (from left to right) in the inset, which reproduces the magnetization curve of Fig. 4. The ratio approaches unity with increasing system size for all temperatures except the lowest, which is within the low-temperature confining phase. No qualitative change is seen across the higher-temperature transition, indicating that the two higher-temperature phases both exhibit deconfinement of monopoles. Explicitly, the curves have, from bottom to top, $T/p=3.125$ (black), $3.175$ (red), $3.226$ (orange), $3.279$ (yellow), $3.333$ (green), $3.390$ (light blue), and $3.448$ (dark blue). In all cases, $V_4/T=0.05$ is fixed.

Figure 9

Effective (dimensionless) interaction $U_{\text{m}}\left(\mathbf{r}\right)/T=-ln{G}_{\text{m}}\left(\mathbf{r}\right)$ between monopoles, for fixed system size $L=32$. For each temperature $T$, the interaction is plotted for separations $\mathbf{r}$ parallel $(\parallel )$ and perpendicular $(\perp )$ to the axis of the applied pressure, and the zero of interaction is chosen as $U_{\text{m}}\left({\mathbf{R}}_{\text{max}}\right)=0$. The lines show least-square fits to the Coulomb form $a-b/\left|\mathbf{r}\right|$ in the region $0.1<\left|\mathbf{r}\right|/L<0.4$, with different parameters $a$ and $b$ for each case. There are deviations from the fit at large separation, because of finite-size effects, and at small separation, because of lattice-scale effects and because of the finite range of the additional interactions. (The parameter $b$ is given by $b_{\parallel }=0.0057,b_{\perp }=0.0036$ for $T/p=3.448$; and $b_{\parallel }=0.0247,b_{\perp }=0.0116$ for $T/p=3.226$.)

Figure 10

Structure factor for (spin-flip) polarized neutron scattering [17] ${\mathcal{S}}_{\text{SF}}\left(\mathbf{Q}\right)$, defined in Eq. (5), for scattering wave vector $\mathbf{Q}$ in the $\left(hh\ell \right)$ plane and incident polarization $\mathbf{P}$ along $\left[1\overline{1}0\right]$. The first three plots are for temperatures above $T_{\text{c>}}$, the fourth is in the intermediate phase, $T_{\text{c<}}<T<T_{\text{c>}}$, and the last is below $T_{\text{c<}}$. Pinch points, characteristic of the dipolar correlations of the Coulomb phase, are visible at all temperatures but the lowest. In the intermediate phase, there are also Bragg peaks at certain reciprocal-lattice vectors, indicating spontaneous magnetization. The system size is $L=16$ and all plots have $V_4/T=0.05$. (Wave vectors are measured in units corresponding to the conventional cubic unit cell.)

Figure 11

Plot of $L\langle M_z^2\rangle$ vs temperature near the higher-temperature transition, for various system sizes $L$. (In each case $V_4/T=0.05$ is fixed.) This quantity has vanishing scaling dimension for the mean-field universality class, consistent with the crossing point for large $L$, at ${(T/p)}_{\text{c>}}=3.3509\left(3\right)$.

Figure 12

Log-log plot of the temperature derivative (in arbitrary units) of $L\langle M_z^2\rangle$, evaluated at $T=T_{\text{c>}}$, vs system size $L$. The slope gives the reciprocal of the correlation-length exponent $\nu$. The (blue) solid line has $\nu =\frac{1}{2}$, as expected for the mean-field universality class, while the (purple) dashed line has the Ising value $\nu =0.63$.

Figure 13

Magnetization near the lower-temperature transition at $T_{\text{c<}}$, indicated with a vertical line, for large system sizes. As $L$ increases, the magnetization approaches saturation below the transition, and a kink develops at $T_{\text{c<}}$ (compare also Fig. 4, for $L=24)$. In this case, the MC simulation is run starting from a state with saturated magnetization, $\langle M_z\rangle =M_{\text{sat}}$; the low temperatures and large system sizes ensure that ergodicity is broken and the order parameter takes a nonzero value.

Figure 14

Order parameter for the lower-temperature transition, $1-\langle M_z\rangle /M_{\text{sat}}$, multiplied by $L$ and plotted vs temperature for various $L$ (using the same symbols as Fig. 13). This quantity has vanishing scaling dimension at the Kasteleyn transition; a crossing is observed at ${(T/p)}_{\text{c<}}=3.15302\left(9\right)$.

Figure 15

Ground-state energy $E_{\text{gs}}$, in the Hilbert-space sector with $N$ particles, for hard-core bosons with particle-hole symmetry, Eq. (A1). Results were obtained using exact diagonalization on a $4\times 4$ square lattice with periodic boundary conditions. The main figure shows $E_{\text{gs}}\left(N\right)$, for the labeled values of $N$, vs nearest-neighbor attraction $V$, both in units of the tunneling strength $t$. The model has an RK point at $V=t$, at which $E_{\text{gs}}$ is independent of $N$. For $V<t$ the ground state of the system occurs for half filling, $N=8$. For $V>t$ there are two generate ground states, with boson density 0 and 1 $(N=0$ and $N=16$ respectively), which spontaneously break particle-hole symmetry. The insets show $E_{\text{gs}}$ vs $N$ for fixed values of $V/t$, indicated by the arrows.

Figure 16

Ground-state energy $E_{\text{gs}}\left(N\right)$ vs nearest-neighbor attraction $V$, as in Fig. 15, with an additional four-body repulsion $V_4=0.3t$, Eq. (A3). As in the case with $V_4=0$, the ground state is at half filling for $V\ll t$ and has density either 0 or 1 for $V\gg t$. In between, there is a regime where the minimum of $E_{\text{gs}}$ crosses over between the extremes, stepping through each intermediate density in turn. In the thermodynamic limit, this is expected to become a phase with continuously varying density, separated from small- and large-$V/t$ phases by continuous quantum phase transitions.