Tensor network study of the $m=\frac{1}{2}$ magnetization plateau in the Shastry-Sutherland model at finite temperature

The two-dimensional infinite projected entangled pair state tensor network is evolved in imaginary time with the full update (FU) algorithm to simulate the Shastry-Sutherland model in a magnetic field at finite temperature directly in the thermodynamic limit. We focus on the phase transition into the $m=\frac{1}{2}$ magnetization plateau, which was observed in experiments on ${\mathrm{SrCu}}_2{\left({\mathrm{BO}}_3\right)}_2$. For the largest simulated bond dimension, the early evolution in the high-temperature regime is simulated with the simple update (SU) scheme and then, as the correlation length increases, continued with the FU scheme towards the critical regime. We apply a small symmetry-breaking bias field and then extrapolate towards zero bias using a simple scaling theory in the bias field. The combined SU $+$ FU scheme provides an accurate estimate of the critical temperature, even though the results could not be fully converged in the bond dimension in the vicinity of the transition. The critical temperature estimate is improved with a generalized scaling theory that combines two divergent length scales: One due to the bias, and the other due to the finite bond dimension. The obtained results are consistent with the transition being in the universality class of the two-dimensional classical Ising model. The estimated critical temperature is 3.5(2) K, which is well above the temperature 2.1 K used in the experiments.

Figure 1

The Shastry-Sunderland model of spins 1/2 (represented by blue dots) arranged on a square lattice. The light blue bands indicate antiferromagnetic Heisenberg couplings between nearest- and some next-nearest-neighbor sites with coupling strengths $J^{\prime }$ and $J$, respectively. It can be represented as a checkerboardlike lattice of dimers (indicated with black ovals). We label dimers belonging to the two sublattices as $A$ and $B$. In our simulations, each dimer is combined into a single effective lattice site with a physical dimension $d=4$. This way, we obtain a nearest-neighbor Hamiltonian on an effective square lattice of dimers. The corresponding tensor network ansatz with one tensor per dimer is shown with the black lines that indicate the virtual iPEPS bonds (the physical legs of the tensors are omitted here).

Figure 2

Temperature dependence of (a) the order-parameter $o\left(T\right)$, (b) its temperature derivative $o^{\prime }\left(T\right)$, and (c) the specific-heat $C_V\left(T\right)$. The results are for the smallest simulated symmetry-breaking bias $h_s/h\approx 2.4\times {10}^{-4}$ and the parameters $J^{\prime }/J=0.63$ and $h/J=1.85$ (fixing the units with $J=\hslash =k_B=1$). Different curves correspond to the iPEPS bond dimension $D=6–9$ where the exact results should be recovered for $D\to \infty$. Both $o^{\prime }\left(T\right)$ and $C_V\left(T\right)$ display sharp peaks at temperatures where the order parameter suddenly rises, indicating the vicinity of the critical point. The inset displays a wider range of temperatures.

Figure 3

Extrapolation of the critical temperature based on a simple scaling theory. In panel (a), we plot the position $T^{*}\left(h_s\right)$ of the maximum of $o^{\prime }(T,h_s)$, see Fig. 2 for different values of the bias and $D=6–9$. Each curve is fitted with the scaling ansatz in Eq. (11). Similarly, in (b), we focus on the maxima of the specific-heat $T_{C_V}^{*}\left(h_s\right)$ and fit it with the scaling ansatz in Eq. (12). The obtained values of the critical temperature $T_c$ and critical exponent $1/\tilde{\beta }\delta$ are collected in Table 1.

Figure 4

Scaling of the correlation length and the magnitude of the order-parameter derivative at its peak. In (a), we show a log-log plot of ${\xi }^{*}\left(h_s\right)$ fitted by the scaling ansatz in Eq. (13), whereas in (b), we focus on $o^{\prime *}\left(h_s\right)$ and the corresponding scaling ansatz in Eq. (14). The fitted critical exponents are gathered in Table 2. Here ${\xi }^{*}$ is measured in units of lattice spacing of the effective dimer lattice, see Fig. 1.

Figure 5

The scaling ansatz combining both the effect of a small bias $h_s$ and the finite iPEPS bond dimension $D$. All three curves $T^{*}(h_s,{\xi }_D)$, indicating the position of the maximum of the order-parameter derivative for $D=7–9$, are fitted by the single scaling ansatz in Eq. (19). The fit yields $T_c=0.043\left(2\right),c=0.83\left(4\right),{\xi }_{D=8}/{\xi }_{D=7}=1.09\left(7\right),{\xi }_{D=9}/{\xi }_{D=7}={25}_{-23.6}^{+\infty }$. The error bars correspond to 68% confidence intervals.

Figure 6

A comparison of FU and FU $+$ SU approaches. In (a), we show $T^{*}(h_s/h)$ where the SU $+$ FU results were obtained with ${\beta }_{SU}=0.24$. In (b), we show the comparison of SU $+$ FU results obtained with ${\beta }_{SU}=0.24$ and ${\beta }_{SU}=0.44$. We observe that for $D=7$ both FU and SU $+$ FU give similar results, but for $D=9$ SU $+$ FU gives more regular results. For both $D=7$ and 9, ${\beta }_{SU}=0.24$ and 0.44 give similar results. The fitted $T_c$ and $1/\tilde{\beta }\delta$ can be found in Table 3.

Figure 7

Comparison of the SU approach with the FU and FU $+$ SU results. In (a) we show $T^{*}(h_s/h)$ where the SU results have been obtained with $D=7–10$. They differ significantly from the FU and SU $+$ FU results, changing slowly with increasing $D$. In (b) we show a log-log plot of the correlation length ${\xi }^{*}(h_s/h)$. The SU results for $D=9\mathrm{to}10$ are clearly inconsistent with the expected power-law behavior of ${\xi }^{*}(h_s/h)$.