Magnetization process in a frustrated plaquette dimerized ladder

TL;DR

This paper investigates the magnetization process in a frustrated plaquette dimerized ladder, revealing complex plateau structures through combined numerical and analytical methods, including the emergence of intermediate magnetization plateaus.

Contribution

It introduces a comprehensive analysis of the magnetization process in a frustrated plaquette dimerized ladder using novel combined techniques and effective models.

Findings

Identification of intermediate magnetization plateaus at M=1/4 and 3/4.

Rich plateau structure due to frustration and dimerization.

Validation of effective models with DMRG simulations.

Abstract

The magnetic phase diagram of a plaquette dimerized antiferromagnetic system is studied by using a combination of numerical and analytical techniques. For the strongly frustrated regime, series expansions and bond operators techniques are employed to analyze zero magnetization plateau, whereas low energy effective models are used to study the complete magnetization process. The interplay between frustration and dimerization gives rise to a rich plateaus structure that is captured by effective models and corroborated by numerical density matrix renormalization group simulations, in particular the emergence of intermediate plateaus at M = 1/4 and 3/4 of saturation in the magnetization curve.

Tables

Table 1. Table 1: Eigenvalues corresponding to the 16 eigenstates of the plaquette. The four low energy states involved in the perturbative calculation are labeled as | s 0 ⟩ ket subscript 𝑠 0 |s_{0}\rangle , | t A ⟩ ket subscript 𝑡 𝐴 |t_{A}\rangle , | t B ⟩ ket subscript 𝑡 𝐵 |t_{B}\rangle and | q 0 ⟩ ket subscript 𝑞 0 |q_{0}\rangle .

$S_T$	$S_T^z$		Energy
$0$	$0$	$|s_0⟩$	$-\frac{J}{2}-\frac{J_d}{4}-\frac{1}{2}\sqrt{4J^2-2J{J}_d-4J+J_d^2-2J_d+4}-\frac{1}{2}$
	$0$		$-\frac{J}{2}-\frac{J_d}{4}+\frac{1}{2}\sqrt{4J^2-2J{J}_d-4J+J_d^2-2J_d+4}-\frac{1}{2}$
	$-1$	$|t_A⟩$	$-h-\frac{J}{2}+\frac{J_d}{4}-\frac{1}{2}$
	$-1$	$|t_B⟩$	$-h-\frac{J_d}{4}-\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
	$-1$		$-h-\frac{J_d}{4}+\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
	$0$		$-\frac{J_d}{4}-\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
$1$	$0$		$-\frac{J_d}{4}+\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
	$0$		$-\frac{J}{2}+\frac{J_d}{4}-\frac{1}{2}$
	$1$		$h-\frac{J}{2}+\frac{J_d}{4}-\frac{1}{2}$
	$1$		$h-\frac{J_d}{4}-\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
	$1$		$h-\frac{J_d}{4}+\frac{1}{2}\sqrt{J^2-2J+J_d^2+1}$
	$-2$	$|q_0⟩$	$-2h+\frac{J}{2}+\frac{J_d}{4}+\frac{1}{2}$
	$-1$		$-h+\frac{J}{2}+\frac{J_d}{4}+\frac{1}{2}$
$2$	$0$		$\frac{J}{2}+\frac{J_d}{4}+\frac{1}{2}$
	$1$		$h+\frac{J}{2}+\frac{J_d}{4}+\frac{1}{2}$
	$2$		$2h+\frac{J}{2}+\frac{J_d}{4}+\frac{1}{2}$

Equations

$$N(S-m)\in \mathbb{Z},$$

$$\mathbf{H}$$

$$|s{⟩}_i=\frac{1}{\sqrt{2S+1}}\sum _{m=-S}^S(-1{)}^{m+S}|m,-m{⟩}_i,$$

$$\mathbf{H}$$

$${\mathbf{h}}_J(n)$$

$${\mathbf{h}}_K(n)$$

$$\frac{E_{gs}}{N_s}$$

$$\frac{E_{gs}}{N_s}\ge -\frac{3}{2}{J}_0.$$

$${\vec{\mathbf{L}}}_{1,2}(n)$$

$${\vec{\mathbf{L}}}_{3,4}(n)$$

$${\vec{\mathbf{K}}}_{1,2}(n)$$

$${\vec{\mathbf{K}}}_{3,4}(n)$$

$$S_{\stackrel{1,n}{2,n}}^{\alpha }$$

$$S_{\stackrel{3,n}{4,n}}^{\alpha }$$

$$s_{A,n}^{†}{s}_{A,n}+\sum _{\alpha }{a}_{\alpha ,n}^{†}{a}_{\alpha ,n}^{}$$

$$s_{B,n}^{†}{s}_{B,n}+\sum _{\alpha }{b}_{\alpha ,n}^{†}{b}_{\alpha ,n}^{}$$

$$E=-\frac{3}{4}{J}_0-\frac{3}{2}{s}^2{J}_0+(5-2s^2)\lambda +\frac{3}{N}\sum _k({\omega }_k^{+}+{\omega }_k^{-}),$$

$${\omega }_k^{\pm }=|a|\sqrt{1\pm |g|\frac{s^2}{a}{\beta }_{\delta }(k)},$$

$${\beta }_{\delta }(k)$$

$$d$$

$$a$$

$$|s_0⟩\stackrel{h{c}_{1A}}{\to }|t_A⟩\stackrel{h{c}_{2A}}{\to }|q_0⟩,J_d<\frac{2J}{1+J},$$

$$|s_0⟩\stackrel{h{c}_{1B}}{\to }|t_B⟩\stackrel{h{c}_{2B}}{\to }|q_0⟩,J_d>\frac{2J}{1+J}.$$

$${\mathbf{H}}_{\mu }={\mathbf{H}}_{0,\mu }+{\mathbf{H}}_{int,\mu },$$

$${\mathbf{H}}_{0,\mu }$$

$${\mathbf{H}}_{int,\mu }$$

$${\mathbf{H}}_{\text{eff},\mu }={\mathbf{H}}_{\text{eff},\mu }^{(1)}+{\mathbf{H}}_{\text{eff},\mu }^{(2)}+\cdots ,$$

$${\mathbf{H}}_{\text{eff},\mu }^{(1)}=\sum _{i,j=1,2}|v_i⟩⟨v_i|{\mathbf{H}}_{\text{int},\mu }|v_j⟩⟨v_j|,$$

$${\mathbf{H}}_{\text{eff},\mu }^{(1)}$$

$${\mathbf{S}}_n^z$$

$${\mathbf{S}}_n^{†}$$

$${\mathbf{S}}_n^z$$

$${\mathbf{S}}_n^{†}$$

$$h_{L-F}=\pm (1+\Delta ),$$

$$h_{L-N}=\pm \sinh (g)\sum _{n=-\infty }^{\infty }\frac{(-1{)}^n}{\cosh (ng)},$$

Full text

The interplay between frustration and dimerization gives rise to a rich plateaus structure that is captured by effective models and corroborated by numerical density matrix renormalization group simulations, in particular the emergence of intermediate plateaus at M = 1/4 and 3/4 of saturation in the magnetization curve.

I Introduction

Antiferromagnetic spin ladders are paradigmatic examples of quantum magnetism in low dimensions Dagotto-Rice-Science ; Mila-2016 ; mila-2016-bis ; mila-berry-phase ; laflorencie-2011 ; Totsuka1 ; NH4CuCl3-magnetization ; NH4CuCl3-model ; cabra-2001 ; Mila-magnons-spinons . These systems exhibit a plenty of interesting properties, such as absence of long range order, energy gaps, magnetic disordered phases, among others Mikeska-ch1 . One particular issue of these systems that captured both experimental and theoretical efforts is the existence of magnetization plateaus Mila-capitulo-plateaus . For a quantum spin-S chain, the necessary condition for the occurrence of a magnetization plateau has been established by Oshikawa, Yamanaka and Affleck (OYA)OYA as

[TABLE]

where $N$ is the period of the spin state, $S$ the magnitude of spin and $m$ the magnetization per site in units of $g\mu_{B}$. The OYA criterium has been successfully applied to different one dimensional magnets Lamas-Pujol_2010 ; cabra-2001 and revisited by different techniquesTTH ; totsuka-2015 ; topological-condition ; plat-2015 . In this context several quasi-one-dimensional compounds has been studied in the past years. As an example we mention the compound Cu2Cl${}_{4}\cdot$D8C4SO2, which presents a singlet ground state. In order to describe the properties of this material different models have been proposed Fujisawa-2003 ; experimental_tube_prb ; tube_inelastic_scatering , however the magnetic properties of this compound has not been completely understood. Another system with a remarkable magnetic behavior is the $S=\frac{1}{2}$ zig-zag ladder compound NH4CuCl3. In this material the magnetization curve presents two plateaus at $m=\frac{1}{4}$ and $m=\frac{3}{4}$ of the saturation magnetization, irrespective of the external field direction, but not at $m=0$ and $m=\frac{1}{2}$ NH4CuCl3-magnetization .

In this paper we study an antiferromagnet model in the zig-zag geometry showing that the interplay between geometrical frustration and plaquette dimerization promotes the existence of plateaus at $m=\frac{1}{4}$ and $m=\frac{3}{4}$. Although we do not intend to provide a theoretical description of these materials, our results may be relevant for the discussion of the emergence of rational magnetization plateaus and its magnon and spinon excitations. The model analyzed in this work belongs to a family of frustrated antiferromagnetic systems with a fully dimerized exact ground state Brenig-tetrahedral ; Matera-Lamas_2014 ; Lamas-Matera_2015 . Recently it was shown that by imposing certain exact condition to the couplings, the exact ground state of the system is robust against the inclusion of weak disorder Lamas-Matera_2015 . This property provides to the model interesting magnetic properties when the unit cell is increased in order to consider plaquette dimerization. The exact condition can be imposed both in one or the two different square plaquettes giving a rich plateaus structures, an aspect that will be exploited in this work.

The outline of the paper is as follows. In Section II we study the zero magnetic field phase diagram. Starting form the line where the ground state is exactly determined we use a bond operator formalism to study the excitation gap. In Section III we analyze the effect of an external magnetic field on the model. We study the magnetization process by means of density matrix renormalization group (DMRG) calculations and exploring the parameters space we determine regions where different plateaus are present. In the highly dimerized plaquette regime we derive low energy effective models and studying the excitations we estimate the critical magnetic fields and the plateau weights. The interplay between dimerization and frustration results in a rich plateaus structure that is captured by the low energy model and is consistent with the numerical results obtained through DMRG. Finally in Section IV We present the conclusions. Some relevant formulas and expressions used throughout the work are shown in the Appendix.

II Zero magnetic field phase diagram

II.1 Exact ground state in the absence of magnetic field

We consider the following Heisenberg model on a two legs spin-$S$ ladder

[TABLE]

where $\vec{\mathbf{S}}_{j,n}$ represents the $j$-th spin on the unit-cell $n$. This ladder is schematized in Fig. 1. Following the lines scketched in Ref. Matera-Lamas_2014, we can show that, provided the conditions $J_{d}=2J$ and $K_{d}=2K$, the ground state of the system corresponds to the fully dimerized state $|{\psi}\rangle=\bigotimes_{i=1}^{N_{r}}|s\rangle_{i}$, with

[TABLE]

where the index $i$ label the rungs in the ladder and $N_{r}$ is the number of rungs. The state $|m,-m\rangle_{i}$ is a product state such that ${\bf S}^{z}_{2i-1,n}|m,-m\rangle_{i}=-{\bf S}^{z}_{2i,n}|m,-m\rangle_{i}=m|m,-m\rangle_{i}$ where the index $i=1,2$ corresponds to the two rungs in the unit-cell $n$. This state remains the ground state of the system in the range of couplings $J<J_{0}/2$ and $K<J_{0}/2$ . For this purpose, it is convenient to rewrite the Hamiltonian in the case $J_{d}=2J$, $K_{d}=2K$ as follows

[TABLE]

where $\mathbf{h}_{J}(n)$ and $\mathbf{h}_{K}(n)$ are Hamiltonians of the square plaquettes on unit-cell $n$, given by

[TABLE]

[TABLE]

In this way, is easy to see that $\mathbf{H}$ is a semi-definite positive operator if $J<J_{0}/2$ and $K<J_{0}/2$. In the case $S=1/2$, the energy per square plaquette is bounded from bellow by

[TABLE]

where $l_{1,2}$ ($l_{3,4}$) is the total spin in the rung between sites $1$ and $2$ ($3$ and $4$). It is straightforward to see that the minimum value correspond to $l_{1,2}=l_{3,4}=0$, given the bounding condition

[TABLE]

To conclude, writing the Hamiltonian in terms of local operators in each rung

[TABLE]

is easy to see that if we restrict the couplings to satisfy the condition $J_{d}=2J$ and $K_{d}=2K$, the state $|{\psi}\rangle$ is an eigenstate with eigenvalue $E=-\frac{3}{2}J_{0}N_{s}$. Therefore we have proved that the state $|{\psi}\rangle$ is the ground state of the system along the lines $J_{d}=2J$ and $K_{d}=2K$ as long as the couplings satisfy $J<J_{0}/2$ and $K<J_{0}/2$. Notice that these two last constraints represent a sufficient condition for the existence of this singlet ground state, but the true range of couplings can be larger.

This exact result is a very convenient starting point for dimer expansions. Notice that the existence of a line where the singlet product state is the true ground state of the system allows us to use dimer expansions around a point in the parameter space where the Hamiltonian is maximally frustrated. In the following Section we present some results obtained by different techniques based on dimer expansions.

II.2 Bond Operators appproach

Close to the line in the parameter space where the ground state corresponding to Hamiltonian II.1 is exactly known it is convenient to exploit the description in terms of bosonic operators, the so called bond operators (BO) SachdevBO which label the dimer’s singlet-triplet spectrum. On the exact line the ground state is a direct product of states on each square plaquette. On each plaquette the state corresponds to a product of a singlet between spins located at sites $1,n$ and $2,n$ and another singlet between spins at $3,n$ and $4,n$. Within BO theory the four spins $\vec{S}_{i}$ on each square plaquette are expressed as

{fleqn}

[TABLE]

where $s_{A,n}^{(\dagger)}$and $a_{\alpha,n}^{(\dagger)}$ destroy(create) the singlet and triplet states of the dimer between sites $1,n$ and $2,n$ and greek labels, $\alpha=1,2,3$, refer to the threefold triplet multiplet. Equivalently operators $s_{B,n}^{(\dagger)}$and $b_{\alpha,n}^{(\dagger)}$ act on the dimer between sites $3,n$ and $4,n$. A hard-core constraint

[TABLE]

is implied on each sublattice, which renders the algebra of the r.h.s of Equations. (10) and (11)identical to that of spins. Inserting the BO representation into a spin model leads to an interacting Bose gas. In the BO-MFT, singlets are condensed by $s^{(\dagger)}\rightarrow s\in\,$Re and the constraint in the number of bosons per site is satisfied on average with a global Lagrange multiplier $\lambda$. In this approach terms only up to second order in the BOs are retained and the quadratic Hamiltonian can be diagonalized by standard Bogoliubov transformation leading to an energy $E$ per unit cell of

[TABLE]

with the triplon dispersion

[TABLE]

where

[TABLE]

$a=\frac{J_{0}}{4}-\lambda$, $g=J_{d}-2J$, and $s=\langle s_{A,n}\rangle=\langle s^{\dagger}_{A,n}\rangle=\langle s_{B,n}\rangle=\langle s^{\dagger}_{B,n}\rangle$ is a real parameter and we have restricted the study to the case $K=\delta J$ and $K_{d}=\delta_{d}J_{d}$ with $\delta=\delta_{d}$ ie homogeneous dimerization. The general case will be analyzed with other techniques in the following Sections. In order to obtain the mean field parameters $s$ and $a$ the energy $E$ is extremized, obtaining two self-consistency equations $\partial E/\partial a=0$ and $\partial E/\partial s=0$. These two equations can be combined obtaining

[TABLE]

where $\gamma_{\pm}(k)=\sqrt{1\pm d|g|\beta_{\delta}(k)}$ and $d=s^{2}/a$. These equations can be solved sel-consistently.

We mention in passing, that the trivial limit, i.e. $g=0$, leads to $d=1/J_{0}$, $s=1$, and $\lambda=-\frac{3}{4}J_{0}$, and therefore to a singlet-triplet gap of $\Delta=1$ and a ground state energy of $E/N={-}\frac{3}{2}J_{0}$, which is consistent with two saturated singlets per unit cell.

II.3 Triplet dispersion and energy gap

In this Section we will analyze the triplet dispersion in absence of an external magnetic field in the model given by Eq.(II.1). For this, we will compare the results obtained by applying techniques of different nature, which complement each other and bring a perspective that none of the techniques can provide separately. The methods used here are on the one hand the numerical technique DMRG, and on the other hand strong coupling dimer and plaquette expansions, as the BO method from the previous Section, and the series expansions technique (SE) considered in the Appendix.

The Fig. 2 summarizes the main result of this Section, obtained by the different techniques mentioned previously, for the triplet gap, ie the minimum of the dispersion, as a function of the parameter $g=J_{d}-2J$. The case considered in this figure is that of a homogeneous structure without dimerization, ie $\delta=\delta_{d}=1$. The DMRG results, shown with blue dots in Fig. 2, for $J=0,0.1,0.2,0.3$ and 0.4 from right to left, and were calculated on finite structures of $L=40$ plaquettes (number of sites $N_{s}=160$) with periodic boundary conditions (PBC), maintaining up to $M=520$ states in the computation, which has shown to be enough to achieve the required precision. Unless specified otherwise, in the remainder of the work these will be the DMRG parameters used. All the DMRG computations presented in this work were performed with the open-source code ALPS ALPS .

On the other hand, the SE results are indicated by red line in this figure with the same $J-$parametrization as DMRG, and were obtained from an expansion in dimers to $O(10)$ by using the method of continuous unitary transformations (CUT) SE-CUT . Finally, the BO results are of two types. The first one is the Holstein-Primakoff approximation (BO-HP), which is indicated by green line in the Fig. 2. In the BO-HP the self-consistent parameters in Eqs.(17,18) are solved for $g=0$ (obtaining $s=a=d=1$) and are replaced in the dispersion Eq.(15) which, after the re-scaling $k\rightarrow 2k$, imposed by halving the unit cell (dimer basis), takes the simple form $\omega_{HP}=\sqrt{1-|g|\cos k}$. Note that this expression coincides with the treatment of the model under the random phase approximation (RPA) method link-BO-SE . The second type of BO result is the self-consistent solution of mean-field BO, Eqs.(17 ,18), and is indicated by the orange line in the Fig. 2.

Several comments are in order. First of all note that, apart from the $g-$dependence, exists a residual interaction, represented here through the $J-$parametrization, that is captured by SE and DMRG, but not by BO. That is, BO encapsulates the entire dependence on $J$ and $J_{d}$ in $g$, which is artifact of the mean field approximation.

In addition, as can be observed, there is a very good quantitative agreement between DMRG and SE results, which is maintained until intermediate values of $J$, according to the perturbative character of the SE approximation.

On the other hand, note that the absolute maximum of the gap $(=1)$ at $g=0$, as well as the other maxima in Fig. 2, are non-derivable (peaks), representing a crossing of the dispersion branches at $k=0$ (right) and $k=\pi$ (left) of the peak, which are detected by the SE. This in turn reflects a change from a dimerized chain-type to ladder-type of spectrum, respectively.

Regarding the role of the residual interaction, it shifts the gap peak to the left $(g<0)$. However its effect on the two dispersion branches is very different. While the branch $k=0$ (to the right of the peak) moves along with the peak to the left by increasing $J$, the branch $k=\pi$ on the left tends to collapse on a single curve, ie its dependence on $J$ and $J_{d}$ is captured mainly by $g$. These features are detected by means of SE and numerically corroborated by DMRG as can be seen in Fig. 2.

Finally, with respect to BO, it is worth noting that although the mean field level does not take into account the residual interaction, but only g, it is able to capture some essential aspects of the gap behaviour, not only qualitatively, but also up to a certain quantitative level. Here the two approaches BO describe the gap with different precision on both sides of the gap peak. As can be seen in Fig. 2, the BO-HP (green line) adequately describes the right branch of the gap (in particular the $J=0$ line), while the BO-MFT (orange line) agreement for the left branch is very good, as compared with SE and DMRG results. This tendency of the dispersions to collapse with a single $g-$dependence for $g<0$ could be related to the fact that in that region the gap remains large, ie more adiabatically connected to the exact line $g=0$, as shown in the next Section through a sweep in the $J-J_{d}$ plane of the gap. Finally let us mention that in quantitative terms, the role of the residual interaction and the dependence on $g$ is clarified in the Appendix, where it is explicitly shown that at leading order the dispersions obtained by SE and BO-HP coincide and depend only on $g$.

III Magnetic field phases

III.1 Plateau phases without plaquette dimerization

In this Section we will analyze the effect of an external magnetic field on the model given by Eq.(II.1). To simplify the analysis we will start by studying the case of homogeneous frustrated plaquettes, ie $\delta=\delta_{d}=1$. Although this case has been previously studied in different regions of the space of parameters and employing different techniques Shastry-Sutherland81 ; Tonegawa87 ; Chitra95 ; White96 ; Mikeska96 ; Yokoyama99 , our analysis presents it in a unified way and contextualizes our subsequent study of the general case with frustration and dimerization between-plaquettes of the following Sections.

As it is known the presence of a magnetic field opens the possibility of the emergence of plateau phases in the magnetization curves, where the magnetization is fixed at a certain value, for a finite interval of the magnetic field.

In Fig. 3 we show the critical fields $h_{c}$ delimiting the different possible values of magnetization, determined by DMRG, in function of $J/J_{d}$, along the line $J+J_{d}=1$. The choice of the parameters in this figure is in order to show in a single representation the emergence of intermediate plateaus. As can be observed, apart from the zero field $(M=0)$ and saturation $(M=1)$ plateaus, two additional plateaus at $M=1/3$ and 1/2 are present. Note that the critical line crossing each of this two plateaus represents a small jump that is an artifact due to the open boundary conditions used for the computation.

Additionally, in Fig. 3, the critical line delineating the plateau at $M=0$, obtained by the SE triplet gap closure condition, is presented with red line. The peak of this curve represents the crossing of the gap closing branches at $k=0$ and $k=\pi$. As can be observed, the quantitative agreement between DMRG and SE for the line $M=0$ is excellent throughout the range of parameters shown in the figure.

On the other hand, Fig 4 shows the extension of the plateaus at $M=0$ (a), $1/3$ (b) and $1/2$ (c), in the plane $J-J_{d}$ calculated on finite systems with the DMRG technique. The determination is not intended to be quantitatively accurate, for which an analysis of a finite-size scaling that goes beyond the objectives of the work should be carried out, but rather to determine the regions in the parameters space where the plateaus are most prominent.

Let us first consider Fig 4(a) where the width of $M=0$ plateau is considered. The largest plateau width ($\Delta=1$) is at the origin, corresponding to the limiting case of isolated dimers (see Fig. 1). The maximum width of the plateau follows the line of maximum frustration, $g=0$, indicated in this figure with dotted blue line, up to $J\simeq 0.2$, beyond which deviations towards $g<0$ become increasingly marked. In the lower right inset of Fig. 4 (a) the SE triplet gap is depicted showing the same trend as the DMRG determination. This deviation is consistent with the observed shift of the gap peak in Fig. 2, predicted by SE and DMRG, as well as the tendency to collapse into a single $g-$dependent triplet gap curve for $g<0$.

Another particular case of our model is represented by the point $(J=0,J_{d}=1)$, corresponding to the homogeneous Heisenberg chain, which is gapless. Note that both DMRG and SE predict a tendency toward the gap closure (dark zone) approaching that point. The model remains gapless along the line $J_{d}=1$ up to the point $(J\simeq 0.241,J_{d}=1)$, indicated by a green circle in the Fig. 4 (a)). From this point the system is gapped and for $(J=0.5,J_{d}=1)$ (Majumdar-Gosh point) the model displays the exact dimer-product zero field ground state, which extends all over the line $g=0$ analyzed in Section II.1.

Regarding the $M=1/3$ plateau, the Fig. 4 (b), shows that this plateau emerges in a certain reduced and highly frustrated chain limit of the model, where $J,J_{d}$, are of the order of the unity, reaching values of plateau width of $\approx 0.6$. The nature of this plateau is completely different to the others analyzed in this work, having a classical origin NN-Ising72 , characterized by an ”up-up-down” ordering which stabilizes this magnetization value. This state adiabatically evolves from the Ising chain with first and second neighbors interactions, and survives up to the Heisenberg isotropic limit. This plateau has been extensively analyzed in the works of Refs. Okunishi-Tonegawa03PRB, ; Okunishi-Tonegawa03JPSJ, ; AndreasM1/3, ; Tonegawa04, . In addition, the techniques of strong plaquette expansions that we will use (see next Section) are not specially suitable to treat this plateau, since it emerges in a region where all the couplings are of the same order. For these reasons we will not elaborate more on the same in this work.

Finally the panel (c) of Fig. 4 depicts the extension of the $M=1/2$ plateau determined by DMRG. The emergence of this plateau has been analyzed in terms of a dimerized and frustrated ladder from a numerical point of view Tonegawa-ladder-98 ; Tonegawa04 . The dashed yellow lines indicate estimations of high field $(h>1)$ limit of $M=1/2$ phase plateau, from an effective low energy dimer model description of the system Mila-ladder-98 ; Sen-ladder-99 , indicating a similar trend to our findings.

In the following Sections we will develop a low energy plaquette effective model, which will account for present plateaus (except $M=1/3$, wich would require higher order terms) and predict the emergence of others at $M=1/4$ and $3/4$ due to the interplay between frustration and inter-plaquette dimerization.

III.2 Effective model description of dimerized plaquette Hamiltonian

In the previous Sections we have considered the frustrated homogeneous plaquette case, where the frustration degree is controlled by $g=J_{d}-2J$, being maximal along the line $g=0$, where the model exhibits an exact rung-dimer product state. In this Section we will explore the effect of plaquette dimerization via the parameters $\delta=K/J$ and $\delta_{d}=K_{d}/K$. The interplay between on-plaquette frustration and inter-plaquette dimerization gives rise to a richer plateaus structure that we will analyze in a weak interacting plaquettes regime, starting from the limit of decoupled plaquettes. It is worth mentioning that an equivalent model, although with a different parametrization has been studied analytically using a field theory bosonization approach Cabra-Grynberg-00 . This FLuttstudy, however, is complementary to ours, since it starts from the strong chain limit.

We start by analyzing the Hilbert space corresponding to isolated plaquettes, which is spanned by a basis of sixteen states containing two singlets, nine triplets and five quintuplets, which we identify generically by $|s_{i}\rangle,|t_{j}\rangle$ and $|q_{k}\rangle$, respectively. This states are listed in the Appendix A (Table 1). On each plaquette in the absence of magnetic field $h$ the ground state corresponds to the singlet $|s_{0}\rangle$. When the magnetic field is turned on the energy corresponding to the singlet and a triplet becomes closer until a critical value where a level crossing occurs and a triplet becomes the ground state. By further increasing the magnetic field we obtain a second level crossing between triplet and the quintet state $|q_{0}\rangle$. Here the local structure of the plaquettes couplings becomes important. According to the values of the internal couplings two different types of triplets are involved in the crossings, depending on whether $J_{d}<\frac{2J}{1+J}$ or $J_{d}>\frac{2J}{1+J}$. In the following we refer these two different scenarios as cases A and B, respectively. In the case A (B) the singlet state $|s_{0}\rangle$ is degenerate with a triplet state $|t_{A(B)}\rangle$ in the first transition, at the critical field $hc_{1A(B)}$, whereas $|t_{A(B)}\rangle$ is degenerate with the quintet state $|q_{0}\rangle$ in the second transition at the critical field $hc_{2A(B)}$. These possibilities are represented schematically in eq.(III.2)

[TABLE]

Around these level crossing points we can derive an effective model. Let us write the Hamiltonian as follows

[TABLE]

where the extra index $\mu$ indicates the possible cases mentioned, ie $\mu=\{1A,2A,1B,2B\}$, being

[TABLE]

We construct the effective Hamiltonian via degenerate perturbation theory Totsuka1 ; Mila-ladder-98 ; Sen-ladder-99 ; mila-10 ; 4-tube-3 ,

[TABLE]

where the superscript indicates the order of the perturbation, although for our purpose in this paper we will only consider the first order of the expansion.

Since at the crossing point the ground state of each plaquette is doubly degenerate (except for $J_{d}=\frac{2J}{1+J}$ which is triply degenerate and will not be considered here), the original dimerized plaquette model with sixteen states per site is reduced to an effective spin-1/2 chain model involving the two mentioned degenerate low energy states at the crossover.

The first order effective Hamiltonian around $hc_{\mu}$ is obtained by applying the standard degenerate perturbation theory

[TABLE]

where states $|v_{i}\rangle$ and $|v_{j}\rangle$ span the set $\{|s_{0}\rangle,|t_{A(B)\rangle}\}$ and $\{|q_{0}\rangle,|t_{A(B)\rangle}\}$ for $\mu=1A(B)$ and $\mu=2A(B)$, respectively. In this way is obtained, up to a constant term,

[TABLE]

where effective couplings $J_{xy,\mu},J_{zz,\mu}$ and $\tilde{h}_{\mu}$ are functions of the original plaquette model couplings $J,J_{d},K,K_{d}$ and $h$. The explicit expressions of these functions are available in the Appendix A (Eqs. 33-48). In addition, the pseudo spin$-1/2$ operators in Eq. 24 are projectors of the degenerate plaquette basis

[TABLE]

[TABLE]

where $\mathbf{S}^{\dagger}_{n}=\mathbf{S}^{x}_{n}+i\mathbf{S}^{y}_{n}$.

III.3 Phase diagram of the spin$-1/2$ chain effective Hamiltonian

The effective model given by Eq.(24) corresponds to an XXZ spin$-\frac{1}{2}$ chain and can be solved exactly via the Bethe ansatz bethe . Here we briefly review the main characteristics of the Bethe ansatz solution, considering the different phases present, and their implications in the magnetization process. For simplicity in this discussion we will simplify the notation of Eqs.(24)by removing supra- and sub- indices and indicating the three principal parameters of the model as $J_{xy}$, $J_{zz}$ and the magnetic field $h$. The phase diagram in the plane $(\Delta=J_{zz}/J_{xy},h/J_{xy})$ is shown in the Fig.(5,a). As can be observed there are three phases: Ferromagnetic (F, yellow), Néel (N, blue) and Luttinger Liquid (LL, gray). The F and N phases are gapped, with spin$-1$ magnon and spin$-1/2$ spinon domain wall type of excitations, respectively. On the other hand, the intermediate phase LL is gapless and exhibits quasi-long range order. Regarding critical lines in Fig. 5(a), the straight lines separating LL and F phases are explicitly given in terms of $\Delta$ by

[TABLE]

whereas the border between LL and N phase is obtained by solving chain-in-field

[TABLE]

where $g=\operatorname{arcosh}\Delta$. Note that both, $h_{L-F}$ and $h_{L-N}$, are given in units of $J_{xy}$.

The different phases present in the model translate into distinctive characteristics in the magnetization curves. The gapped phases F and N exhibit a plateau in the magnetization curve, while in the LL phase the magnetization continuously increases with the magnetic field. This dependence is illustrated in Fig.5. In panel (a) of this figure three paths with red, blue and green dotted lines are indicated in the phase diagram, while in Fig. 5(b) the magnetization curves corresponding to these three paths are showed schematically. Note that the vertical axis in panel (a) is represented as horizontal axis in panel (b). The red curve represents the simplest case, where the system is always in the F phase and has a stepwise structure with two plateaus, jumping directly from one to another with the sign change of the magnetic field. The blue magnetization curve represents an intermediate case, in which the two plateaus corresponding to the F phase are connected through a region of continuous growth of the magnetization with the field, whose shape is characteristic of the LL phase. Finally, the green magnetization curve in Fig. 5(b) indicates the more general case, which comprises the three phases. In this case, in addition to the plateaus at the extremes due to the F phase, an additional intermediate plateau indicative of the N phase is present. These three plateaus are connected by LL phases, consistently with the path indicated in green in Fig. 5(a).

Before concluding with this Section we would like to comment on another case where the plateau structure has a simple stepped form. This case is the Ising limit of the model, in which $J_{x,y}\rightarrow 0\,(\Delta\rightarrow\infty)$, ie the right infinite end of Fig. 5(a). Here the magnetization curve exhibits the three plateaus with direct jumps between them, without intermediate LL phase.

III.4 Effective model description of magnetization process

In this Section we will analyze how the effective low-energy spin-1/2 chain model with anisotropy $\Delta$, obtained previously, captures the main aspects of the magnetization process of the original plaquette model in a certain, strong plaquette coupling, regime of the parameters space. The effective model is not only useful to account for the numerical results obtained through DMRG, but also brings a direct physical interpretation of the structure of plateaus that can emerge in the original model. In addition, it provides a tool that allows us to detect some characteristics, such as small intermediate plateaus, which would be difficult to detect by direct numerical scanning.

To apply the results of the effective model to the original plaquette model, note that according to Eq.(III.2), the effective pseudo-spin operators eigenvalues $S_{z}=-1/2(1/2)$ correspond to plaquette $S_{z}=0(1)$, of the singlet $|s_{0}\rangle$ and triplet $|t_{A(B)}\rangle$ eigenvalues of the first level crossing $1A(B)$. Therefore the three possible magnetization plateaus in the effective chain model at $\{-1,0,1\}$ are mapped onto $\{0,1/4,1/2\}$ magnetization plateaus in the plaquette model, normalized to maximum spin per site and plaquette, respectively. Similarly for the second level crossing $2A(B)$ between the triplet $|t_{A(B)}\rangle$ and quintet $|q_{0}\rangle$, the Eq.(III.2) maps $\{-1,0,1\}$ onto $\{1/2,3/4,1\}$, corresponding to effective chain and plaquette models, respectively.

Now we will describe qualitatively the structure of the magnetization curves in terms of the effective models. To this end, the Fig. 6 illustrates different representative plateaus structures in the original model obtained by means of DMRG, showed with blue solid line. In particular, the Fig. 6(a) shows magnetization curves with a single transition between type-F phases at zero magnetic field and saturation, intermediated by a LL phase. Note that the structure of small steps in the LL portion of the curve is due to the finite size of the numerical calculation. In Fig.6(b) an additional plateau emerges at mid-magnetization. However, according to the previous analysis, this plateau is type-F, like the other two present here. The cases of Figs. 6(a,b) are not dimerized, ie $\delta=\delta_{d}=1$, and illustrate a general feature: in absence of dimerization between plaquettes there can only be plateaus at $M=0,1/2$ and $1$ in the magnetization curve, sharing spin$-1$ magnon excitations. The origin of these plateaus is intrinsic to each isolated plaquette, as a consequence of the two (1A(B), 2A(B)) possible level crossings, which are only renormalized by the interaction between plaquettes.

On the other hand, the Fig. 6(c) show a plaquette-dimerized case, where additional plateaus at $M=1/4,3/4$, are present. These plateaus are indicative of the N-type phase which, as we have mentioned, present spinon-like excitations and are generated by the interaction between plaquettes. A necessary (but not sufficient) condition for their emergence is the presence of plaquette dimerization.

Regarding the size of type-N plateaus, we have observed that in general are smaller than those of type-F. For this reason the effective model is very useful for their detection, since numerically they could easily be overlooked, given the quite large number of parameters involved in the original model.

In addition, the effective model not only facilitates the detection of the intermediate N-type plateaus, but also highlights the role of frustration in the emergence of them. In fact it is observed that apart from dimerization, type-N plateaus grow with frustration, being largest along the line of maximum frustration, ie $g=J_{d}-2J=0$ and $g_{k}=K_{d}-2K=0$, whenever this condition can be satisfied, which depends on the crossing levels (A or B) involved. It is interesting to note that in a certain range of validity frustration-induced type-N plateaus are maximal along the maximum frustration line, where in addition the original model (Eq.(1)) exhibits an exact dimer singlet-product zero field ground state. In some sense it suggests that this condition is not only relevant at zero magnetic field, but also play a role in the magnetization process.

Finally, regarding the dotted red and yellow vertical lines showed in the three panels of Fig. 6, they are quantitative estimations of critical fields delimiting the different plateaus based on the dispersion of magnon and spinon excitations of the effective model, and will be discussed in the next Section.

III.5 Magnetic excitations and critical fields

Here we will analyze the excitations on the plateau structures in the effective chain model, which in turn will allow us to obtain the critical fields, ie the borders of the plateaus and compare with the numerical DMRG results from the original model.

First we consider the saturation plateaus at $-(+1)$ in the effective chain model. These correspond to ferromagnetic states, with all the effective spins on each site polarized along the negative and positive direction of the effective field, respectively. The elemental magnon excitation consists of a spin-wave composed by a linear combination of states with a local one-spin inversion, carrying $\Delta S_{z}=+(-1)$. The application of the effective Hamiltonian (Eq.(24)) on this state (subtracting the ground state) gives rise, by a standard calculation, to the magnon (M) dispersion

[TABLE]

Eq.(29) describes all magnon-type of excitations present in the strong-plaquette expansion of the original model, via the mapping given by Eqs.(III.2,III.2), and the effective couplings ((Eqs. 33-48)), namely (i): right (high field) borders of $M=0\quad(+,\mu=1A(B))$ and $M=1/2\quad(+,\mu=2A(B))$; (ii): left (low field) boundaries of $M=1/2\quad(-,\mu=1A(B))$ and $M=1\quad(-,\mu=2A(B))$.

Let us now consider the $M=0$ plateau in the effective model, which is associated to a double degenerate Néel state. Here the elemental spinon excitation consists of a spin-wave composed by a linear combination of domain wall states with two consecutive spin inversions, carrying $\Delta S_{z}=+(-1/2)$. Following a similar procedure as the magnon case we obtain the spinon (S) dispersion

[TABLE]

Equivalently to the previous case, Eq.(30) describes all spinon-type of excitations, namely (i): right (high field) borders of $M=1/4\quad(+,\mu=1A(B))$ and $M=3/4\quad(+,\mu=2A(B))$; (ii): left (low field) boundaries of $M=1/4\quad(-,\mu=1A(B))$ and $M=3/4\quad(-,\mu=2A(B))$.

From the above magnon $\omega^{\pm}_{M,\mu}(k)$ and spinon $\omega^{\pm}_{S,\mu}(k)$ dispersions we can calculate the critical fields, ie the edges and the width of the plateaus present in the model. For this, we impose the condition of gap closure (minimum of the dispersion), which allows to determine the effective fields in terms of the effective couplings and thus the critical fields according to the original couplings of the model. Due to the number of parameters of the model, in general we obtain the critical fields numerically which is much simpler compared to the DMRG or the numerical solution of the Bethe-ansatz equations approaches.

The techniques used in this work to determine the critical fields: DMRG, and effective models involve different types of approximation. DMRG suffers finite size effects and the effective models have a perturbative origin. However, these methods are somewhat complementary and provide consistent results. Although, as we have mentioned, it is not the aim of this work to quantitatively determine the edges of each plateau in the whole parameter space of the model, it is interesting to compare the results of the different techniques in some representative cases. For this purpose, in Fig. 6 we show the critical fields determined by DMRG (full magnetization curve in blue line) and the edges of the different plateaus present by applying the gap closure condition of the corresponding magnon dispersion (red-dashed lines) or spinon dispersions (orange-dashed lines). As can be observed there is a very good agreement between the different techniques, which is maintained up to intermediate couplings (compared with $J_{0}$).

As previously mentioned and as can be seen in the Fig. 6(c) the plateaus at $M=1/4$ and $3/4$ are small compared to the others, so that they could be difficult to detect numerically by sweeping a fairly large parameter space as in this model using DMRG. It is therefore of particular interest to be able to estimate not only the presence but the maximum possible size of these plateaus as well as their interdependence with frustration and dimerization.

From this point of view the effective model offers a direct way of evaluating positions and maximum widths of the plateaus. This will occur when there is no intermediate Luttinger liquid phase between plateaus, ie the limit of the effective Ising model where $J_{xy,\mu}=0$. In this condition the magnon and spinon disperson bands do not propagate (they are flat) and the edges of fields predicted by both merge. From the Eqs (29,30) we see that the effective critical fields at the Ising limit satisfy

[TABLE]

Note in passing that this condition is obtained from the effective model exact line $h_{L-F}$ (Eq.(27)) in the limit $\Delta\rightarrow\infty$.

By solving the Eqs.(31) in terms of the original variables of the model (Eqs. 33-48)) we get the following cases for the widest plateaus $\Delta h$. For case A $(J_{d}<\frac{2J}{1+J})$, the Ising limit is obtained for $g_{k}=K_{d}-2K=0$ with the same $\Delta h=K/2$ for both $M=1/4$ and $M=3/4$ plateaus.

On the other hand, for case B $(J_{d}>\frac{2J}{1+J})$, the plateau widths are not the same, although we can also obtain analytical expressions for both. In the case $M=1/4$ the Ising condition is obtained for both $g=J_{d}-2J=0$ and $g_{k}=K_{d}-2K=0$ although these conditions are not satisfied for $M=3/4$.

IV Conclusions

In the present paper we have studied the magnetization process in a zig-zag quantum antiferromagnet in the presence of both frustration and plaquette dimerization. Due to the complexity of the unit cell, several cases are present giving a rich structure of the magnetization curve. At zero magnetic field the ground state is exactly determined under the conditions $g=J_{d}-2J=0$ and $g_{k}=K_{d}-2K=0$ and corresponds to a product state of spin singlets. Around this highly frustrated line we investigate the $M=0$ plateau by using dimer series expansions, bond operators mean field theory and numerically by means of density matrix renormalization group.

In order to analyze the magnetization process we complement the numerical results with first order low energy effective models starting from the limit of decoupled plaquettes. From a qualitative point of view, the effective models capture the essential features obtained numerically and bring a simple physical interpretation of the emergent structure of plateaus in the original zig-zag model. In addition, the analysis of magnon and spinon excitations in the efective models allows us to obtain estimations of the critical fields bordering the plateaus, which are in a very good quantitative agrement with the numerical computations in a strong plaquette dimerization regime.

Our study also suggest that the combined effect of on-plaquette frustration, controlled by $g$ and $g_{k}$ and inter-plaquette dimerization, via $\delta=K/J$ and $\delta_{d}=K_{d}/J_{d}$, plays a central role in the structure of the magnetization curve. In fact whenever is possible, according to the values of the couplings involved, $M=1/4$ and $3/4$ plateaus are widest along the line of maximum frustration. In this sense we say that these plateaus are induced by frustration, which together with the dimerization between plaquettes, provide a rich interplay between both aspects that are intrinsic to the model. For this reason, the maximally frustrated $g$ and $g_{k}=0$ line is not only fundamental to zero field, but plays a central role beyond zero magnetic field, which is reflected in its influence on the resulting magnetization process.

Acknowledgments

We acknowledge useful discussions with M. Matera. C. A. Lamas is supported by CONICET (PIP 1691) and ANPCyT (PICT 2013-0009)

Appendix A Effective Hamiltonians

As we discuss in the main body of the text, we construct an effective Hamiltonian via degenerate perturbation theory. The energy corresponding to the states on a square plaquette are presented in Table 1. By applying a standard degenerate perturbation theory around the points $hc_{\mu}$ where the lowes eigenvalues become degenerate we obtain, up to a constant term, a first order effective Hamiltonian in the form

[TABLE]

where the couplings are given by

[TABLE]

where

[TABLE]

Appendix B Dimer series expansion

In this Section we display the explicit expressions for ground state energy and triplet dispersion employed in the work, for the case of homogeneous plaquettes ie $\delta=\delta_{d}=1$. The method employed is the continuous unitary transformation (CUT) SE-CUT . Here for brevity we present explicit results to $O(6)$ although the calculations shown have been made using $O(10)$. The ground state energy (per dimer) reads

[TABLE]

where $g=J_{d}-2J$. Note that for $g=0$ only the first term $-\frac{3}{4}$ survives, corresponding to the dimer-product exact ground state analyzed in Section I. Similarly the triplet dispersion is given by $\omega_{SE}=\sum_{n=0}^{\infty}c_{n}\cos(nk),$ where $c_{n}$ coefficients up to O(6) are

[TABLE]

It is interesting to compare previous SE coefficients with the BO-HP expansion $\omega_{HP}=\sqrt{1-g\cos k}=\sum_{n=0}^{\infty}\tilde{c_{n}}\cos(nk)$, which to O(6) reads

[TABLE]

From Eqs(51) and Eqs(52), we note that bond operators technique is only exact at leading order (blue highlighted terms). This connection between both methods illustrates the effect of neglecting interactions beyond the quadratic order in the bosonic model (ref Mikeska).