Dynamic magnetization process in the frustrated Shastry-Sutherland system TmB₄

Geometrically frustrated spin systems which exhibit very rich magnetic properties have drawn considerable attentions during the last several decades [1]. For example, various experimental and theoretical explorations have been devoted to the emergence of multi-step magnetization (M) curves in frustrated systems such as the triangular spin-chain system Ca₃Co₂O₆[2–5] and the Shastry-Sutherland (S-S) magnets [6–11]. So far, the step-like magnetization curves observed in Ca₃Co₂O₆ are generally believed to be caused by the non-equilibrium magnetization dynamics [12,13], while those in S-S systems are far from being fully understood and remains to be checked.

The S-S lattice [14] has attracted special attentions since its experimental realization in the compound SrCu₂(BO₃₂) in which a fascinating sequence of magnetization plateaus at fractional values of the saturated magnetization ($M_{\textrm{S}}$) have been reported [15,16]. Most recently, similar magnetic behaviors have been identified in rare-earth tetraborides RB₄ (R = Tb, Dy, Ho, Tm, etc.) with the rare-earth moments located on a lattice which is topologically equivalent to the S-S lattice [8,9,17–20]. For instance, the fractional magnetization plateaus at $M/M_{\textrm{S}}=1/2, 1/7, 1/9$ and $1/11,\ldots$ have been experimentally observed at temperature (T) below 4 K in TmB₄ [9]. In contrast to SrCu₂(BO₃₂) with Cu²⁺ ions carrying a quantum spin 1/2, TmB₄ presents a large total magnetic moment (∼6.0 μB) and can be considered as a classical spin system, triggering an extensive theoretical investigation of classical spin models on the S-S lattice [21–28].

It is experimentally indicated that TmB₄ is of strong easy-axis anisotropy caused by crystal field effects. Based on this point, the magnetization process of the classical Ising model on the S-S lattice was investigated using the tensor renormalization group approach, and a single magnetization plateau at $M/M_{\textrm{S}}=1/3$ was predicted at low T for certain coupling constants [21]. In fact, the ground states of the Ising model on the S-S lattice were investigated most recently and the existence of a single 1/3 plateau was rigorously proved [22,23]. The effect of further-neighbor interactions was suggested to eventually explain the magnetization plateaus in TmB₄, and three different ground states with $M/M_{\textrm{S}}=1/2$ were recognized when the additional third-neighbor interaction was considered. On the other hand, the quantum spin-1/2 Ising-like XXZ model with additional interactions on the S-S lattice was studied using the quantum Monte Carlo method, and the plateau at $M/M_{\textrm{S}}=1/2$ was identified [29–31]. It was believed that the emergence of the $M/M_{\textrm{S}}=1/2$ plateau may be due to the quantum fluctuations and long-range interactions. In our earlier work, the presence of the $M/M_{\textrm{S}}=1/2$ plateau was also confirmed when the additional long-range interactions were taken into account in the classical spin model [32].

So far, the main magnetization plateau at $M/M_{\textrm{S}}=1/2$ and some of other small fractional plateaus in TmB₄ can be reproduced by the consideration of the long-range interactions for certain coupling constants. However, the origins of other small fractional plateaus are still under debate. Generally speaking, a frustrated spin system can be easily trapped into metastable states at low T and is hard to relax to the equilibrium state. The time available experimentally may not be sufficient for the spin rearrangement, even though the energy difference between plateaus may be very small. Thus, it is reasonable to assume that the spins in TmB₄ are easily trapped into a metastable state rather than into the equilibrium one at low T, which, to some extent, can be also favored by the obvious hysteresis loop observed in earlier experiments. Furthermore, the formation and the motion of domain walls due to the non-equilibrium magnetization process may play an important role in the emergence of the magnetization plateaus, which has been verified in the study of the triangular spin-chain compound Ca₃Co₂O₆ [13,33,34]. As a matter of fact, the $M/M_{\textrm{S}}=1/8$ plateau has been predicted in the XXZ model on the S-S lattice for a short relaxation time, indicating that this plateau may arise from non-equilibrium state [31]. Thus, one may question if the non-equilibrium magnetization dynamics is also essential for the small fractional magnetization plateaus in TmB₄. A detailed discussion of this question is definitely crucial for the understanding of the fascinating magnetic properties in such a frustrated system. However, as far as we know, few works on this subject have been reported.

In order to make clear this question, we will study the magnetization dynamics of the classical Ising model on the S-S lattice with additional long-range interactions with the Glauber dynamics which has been successfully used in the study of the magnetic properties of Ca₃Co₂O₆ [13]. This algorithm allows us to investigate the dependence of the magnetization curves on temperature and magnetic-field (h) sweep rate. The experimental magnetization plateaus at small fractional values of $M/M_{\textrm{S}}=1/7$, 1/9, and 1/11 followed by the main magnetization plateau at $M/M_{\textrm{S}}=1/2$ can be reproduced at low T for certain magnetic-field sweep rates. In addition, the hysteresis loop observed in experiments can also be qualitatively explained in our simulation [9,19]. These results clearly demonstrate that the fascinating plateaus in TmB₄ magnetization curve may be closely related to the non-equilibrium magnetization dynamics. Our work provides a new insight into the study of the magnetization process for S-S magnets and other similar frustrated spin systems.

The remainder of this paper is organized as follows: In the second section, the model and the simulation method will be presented and described. The third section is attributed to the simulation results and discussion. At last, the conclusion is presented in the fourth section.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In the presence of the long-range interactions and h, the Hamiltonian can be described as follows:

$$\begin{eqnarray} H &= J_1 \sum\limits_{\left\langle {i,j} \right\rangle _1 } {S_i \cdot S_j } + J_2 \sum\limits_{\left\langle {i,j} \right\rangle _2 } {S_i \cdot S_j } + J_3 \sum\limits_{\left\langle {i,j} \right\rangle _3 } {S_i \cdot S_j } \nonumber\\ &+ J_4 \sum\limits_{\left\langle {i,j} \right\rangle _4 } {S_i \cdot S_j } - h\sum\limits_i {S_i^z } , \end{eqnarray} \tag{ 1 }$$

where $\langle i,j\rangle_{1}, \langle i, j\rangle_{2}, \langle i, j\rangle_{3}$, and $\langle i, j\rangle_{4}$ denote the summations over all pairs on the bonds with J₁, J₂, J₃ and J₄ couplings, respectively, as shown in fig. 1, S_i represents the Ising spin with unit length on site i,h is applied along the +z axis and S_i^z denotes the z component of S_i. J₁ = 1 is the antiferromagnetic (AFM) coupling, the coupling ratio J₂/J₁ = 1 is expected from the crystal structure of TmB₄, similar with earlier estimation [30,31], the AFM $J_{3}=0.15J_{1}$ and the ferromagnetic (FM) $J_{4}=-0.15J_{1}$ are estimated to qualitatively reproduce the experimental results. To investigate the magnetization dynamics of the spin system, the simulation is carried out by a single spin-flip rate in the Glauber form [13,34]. The spins are assumed to interact not only with the neighbors and external magnetic field but also with a heat reservoir, based on the Glauber theory [35]. The probability of a spin flip of the i-th spin per Monte Carlo step (MCs) can be described as

$$\begin{equation} W_i = \frac{\alpha}{2}\left[ {1 - S_i \tanh \left( { - \frac{D}{k_{B} T} + \frac{\mu h}{k_B T}} \right)} \right], \end{equation} \tag{ 2 }$$

with

$$\begin{equation} D = J_1 \sum\limits_{\left\langle {i,j} \right\rangle _1 } {S_j } + J_2 \sum\limits_{\left\langle {i,j} \right\rangle _2 } {S_j } + J_3 \sum\limits_{\left\langle {i,j} \right\rangle _3 } {S_j } + J_4 \sum\limits_{\left\langle {i,j} \right\rangle _4 } {S_j } , \end{equation} \tag{ 3 }$$

where α = 0.5 MCs ⁻¹ is the constant of the interaction of a spin with the heat reservoir, k_B = 1 is the Boltzmann constant, μ = 1 is the magnetic moment of the Tm ion. In earlier work, the Glauber-type form of the spin-flip probability has been discussed in detail and successfully used in the frustrated spin-chain system Ca₃Co₂O₆ [13,34]. Similarly, a reasonable value α = 0.5 MCs⁻¹ is used in our simulation to qualitatively coincide with the experimental results.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Our simulation is performed on an L ×L (L = 96 is chosen unless stated otherwise) lattice with periodic boundary conditions. The simulation is started from the saturated magnetization state under high h, which is in the best accordance with the real case. Then the magnetization curves in the decreasing h at different sweep rates for various temperatures are studied to investigate the magnetization dynamics in detail. In addition, the sweep rate is defined by 1/m MCs⁻¹ which means that the 1/m of h unit is updated per MCs.

Figure 2(a) shows the calculated $M/M_{\textrm{S}}$ as a function of h for different magnetic-field sweep rates at a low temperature T = 0.01. Three steps can be observed for the lowest sweep rate of 1/6000 MCs⁻¹. When h decreases down to ∼4.75, M reaches the plateau at $M/M_{\textrm{S}}=1/2$, and then falls down to the step at M = 0 below h ∼ 2.2. More interestingly, the step at M = 0 decomposes into two substeps (zero and nonzero) separated at h ∼ 1 when the magnetic-field sweep rate is increased. The height of the nonzero substep increases with increasing sweep rate. When the magnetic-field sweep rate increases up to 1/60 MCs⁻¹, the $M/M_{\textrm{S}}=1/9$ plateau reported in experiments can be well reproduced in addition to the major magnetization plateaus at M = 0 and $M/M_{\textrm{S}}=1/2$. Furthermore, the dynamic magnetization curves in response to T at the constant sweep rate of 1/60 MCs⁻¹ are also investigated, and the simulated results are presented in fig. 2(b). It is clearly demonstrated that the nonzero substep is heightened as T decreases. When T falls down to 0.008, a magnetization step at $M/M_{\textrm{S}}=1/7$ is observed at intermediate h range, which is consistent with experimental observation [9]. In addition, some additional narrow plateaus can also be noticeable in our simulation, which deserves to be checked in further experiments.

As stated earlier [32], the ferrimagnetic (FI) state spin arrangement consisting of alternative AFM and FM stripes is more favored than the FM state when h is decreased down to the first critical field, while the Neel state is likely stabilized below the following critical field. For the extremely slow sweep rate (1/6000 MCs⁻¹), the single-domain FI state with the plateau at $M/M_{\textrm{S}}=1/2$ and the Neel state with the M = 0 plateau are, respectively, stabilized below these two critical fields, leading to the three-step magnetization curve which is similar to that obtained by the Monte Carlo simulation. To uncover the origin of the nonzero substep, the specimens of configurations for various plateaus under different h at the magnetic-field sweep rate of 1/60 MCs⁻¹ are presented in fig. 3. It is clearly shown that the Neel state grows at a lot of nucleation centers when h falls down to the second critical field, resulting in the domain formation. A mixed state with the Neel order and domain walls constructed of polarized spin chains is responsible for the emergence of the $M/M_{\textrm{S}}=1/9$ plateau in the magnetization curve, as shown in fig. 3(c). The domains become smaller with increasing magnetic-field sweep rate, leading to the heightening of the nonzero substep. When h is further decreased, the domain walls almost disappear (fig. 3(d)), and the magnetization plateau at M = 0 can be observed. On the other hand, the domain boundary mobility may be greatly decreased as T decreases. Thus, the domains become smaller with decreasing T, making the additional magnetization steps more apparent, as confirmed in our simulation (fig. 2(b)). It is noted that the temperature at which the experimental $M/M_{\textrm{S}}=1/7$ plateau is observed is higher than that of the $M/M_{\textrm{S}}=1/9$ one. The inconsistency between the present theory and the experiment may be due to the fact that the disorder effect caused by the inhomogeneity in realistic materials is completely ignored in our simulation. However, our work clearly indicates that the non-equilibrium magnetization dynamics may play an important role in the appearance of the fractional magnetization steps in TmB₄.

In addition, the hysteresis loop is also studied in our work, and the results can qualitatively reproduce the experimental observations [9,19]. Figure 4(a) shows the hysteresis loop at the magnetic-field sweep rate of 1/60 MCs⁻¹ at T = 0.01. The nonzero substep emerges in both field-decreasing and field-increasing branches of the magnetization curve. The nonzero substep ($M/M_{\textrm{S}}=1/9$) in the field-decreasing branch is much higher than that ($M/M_{\textrm{S}}=1/11$) in the field-increasing branch, in a good agreement with the experimental report [9]. This phenomenon demonstrates that the domain structures can be strongly affected by the initial state at a high magnetic-field sweep rate. In fact, the spin configurations of the simulation reveal that the Neel domains for the nonzero substep in the field-increasing branch are generally larger than those in the field-decreasing branch, resulting in the emergence of the hysteresis loop. More interestingly, it is confirmed in our simulation that the value of plateau magnetization varies between different runs, similar with earlier experimental report [9]. The corresponding results are not shown here for brevity. Furthermore, it has been noticed in earlier experiment that the critical fields in the field-decreasing branch of the magnetization curve are respectively smaller than those in the field-increasing branch [19]. This magnetic behavior can also be well reproduced by the simulation.

As stated earlier, the time required for the spin rearrangement likely exceeds the time available experimentally. Thus, the non-equilibrium magnetization dynamics may be essential for the emergence of the fractional magnetization plateaus in TmB₄. This point has been confirmed in this work in which the non-equilibrium evolution is performed by means of the Glauber dynamics. The fractional magnetization plateaus and the hysteresis loop at low T reported in experiments can be reproduced in our simulation. Thus, our work may provide a new insight into the study of the magnetization process of TmB₄, although not all the experimental results can be excellently explained based on the present theory.

At last, the dependence of the step-like magnetization feature on the lattice size L has been investigated in order to exclude the artificial facts caused by the finite lattice size. Figure 4(b) shows the simulated magnetization curves for different L (L = 60, 96, 120 and 150) at T = 0.01 for the extremely fast magnetic-field sweep rate of 1/60 MCs⁻¹. All the simulated curves for various L are almost merged into one, indicating that the finite-size effect on the magnetization of the system is almost negligible and never affects our conclusion.

In summary, we have investigated the classical Shastry-Sutherland Ising model with long-range interactions employing a Glauber-type form of the spin-flip probability in order to understand the dynamic magnetization process in TmB₄. Besides the main $M/M_{\textrm{S}}=1/2$ plateau, other fractional magnetization plateaus at $M/M_{\textrm{S}}=1/n$ (n = 7, 9, 11) observed in experiments can be reproduced in our simulation of the model at low temperatures for certain magnetic-field sweep rates. In addition, the hysteresis loop can be also well explained in the present theory. It is indicated that the magnetization dynamics may be essential for the emergence of those fractional magnetization plateaus. Thus, the present work may provide a new insight into the understanding of the magnetization process for frustrated S-S magnets and other similar frustrated spin systems.

This work was supported by the Natural Science Foundation of China (11204091, 11274094, 11234005), the National Key Projects for Basic Research of China (2011CB922101), China Postdoctoral Science Foundation funded project (2012T50684, 20100480768, 2011M500088), and the Priority Academic Program Development of Jiangsu Higher Education Institutions, China.