High-order study of the quantum critical behavior of a frustrated spin-$\frac{1}{2}$ antiferromagnet on a stacked honeycomb bilayer

We study a frustrated spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ Heisenberg antiferromagnet on an $AA$-stacked bilayer honeycomb lattice. In each layer we consider nearest-neighbor (NN), next-nearest-neighbor, and next-next-nearest-neighbor antiferromagnetic (AFM) exchange couplings $J_1,J_2$, and $J_3$, respectively. The two layers are coupled with an AFM NN exchange coupling $J_1^{\perp }\equiv \delta J_1$. The model is studied for arbitrary values of $\delta$ along the line $J_3=J_2\equiv \alpha J_1$ that includes the most highly frustrated point at $\alpha =\frac{1}{2}$, where the classical ground state is macroscopically degenerate. The coupled cluster method is used at high orders of approximation to calculate the magnetic order parameter and the triplet spin gap. We are thereby able to give an accurate description of the quantum phase diagram of the model in the $\alpha \delta$ plane in the window $0\le \alpha \le 1,0\le \delta \le 1$. This includes two AFM phases with Néel and striped order, and an intermediate gapped paramagnetic phase that exhibits various forms of valence-bond crystalline order. We obtain accurate estimations of the two phase boundaries, $\delta ={\delta }_{c_i}\left(\alpha \right)$, or equivalently, $\alpha ={\alpha }_{c_i}\left(\delta \right)$, with $i=1$ (Néel) and 2 (striped). The two boundaries exhibit an “avoided crossing” behavior with both curves being re-entrant. Thus, in this $\alpha \delta$ window, Néel order exists only for values of $\delta$ in the range ${\delta }_{c_1}^{<}\left(\alpha \right)<\delta <{\delta }_{c_1}^{>}\left(\alpha \right)$, with ${\delta }_{c_1}^{<}\left(\alpha \right)=0$ for $\alpha <{\alpha }_{c_1}\left(0\right)\approx 0.46\left(2\right)$ and ${\delta }_{c_1}^{<}\left(\alpha \right)>0$ for ${\alpha }_{c_1}\left(0\right)<\alpha <{\alpha }_1^{>}\approx 0.49\left(1\right)$, and striped order similarly exists only for values of $\delta$ in the range ${\delta }_{c_2}^{<}\left(\alpha \right)<\delta <{\delta }_{c_2}^{>}\left(\alpha \right)$, with ${\delta }_{c_2}^{<}\left(\alpha \right)=0$ for $\alpha >{\alpha }_{c_2}\left(0\right)\approx 0.600\left(5\right)$ and ${\delta }_{c_2}^{<}\left(\alpha \right)>0$ for ${\alpha }_{c_2}\left(0\right)>\alpha >{\alpha }_2^{<}\approx 0.56\left(1\right)$.

Figure 1

The $J_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the honeycomb bilayer lattice, showing (a) the two layers $A$ (red) and $B$ (blue), the nearest-neighbor bonds $(J_1=$ —–; $J_1^{\perp }=$ - - -), and the four sites $\left(1_A,2_A,1_B,2_B\right)$ of the unit cell; (b) the intralayer bonds $(J_1=$ —–; $J_2=$ - - -; $J_3=-·-·-)$ on each monolayer; (c) the triangular Bravais lattice vectors $\mathbf{a}$ and $\mathbf{b}$, and the monolayer Néel state; and (d) one of the three equivalent monolayer striped states. Sites $\left(1_A,2_B\right)$ and $\left(2_A,1_B\right)$ on the two monolayer triangular sublattices are shown by filled and empty circles respectively.

Figure 2

CCM results for the GS magnetic order parameter $M$ vs the scaled interlayer exchange coupling constant, $\delta \equiv J_1^{\perp }/J_1$, for the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the bilayer honeycomb lattice (with $J_3=J_2$ and $J_1>0)$, for three selected values of the intralayer frustration parameter, $\alpha \equiv J_2/J_1$: (a) $\alpha =0$, (b) $\alpha =0.2$, and (c) $\alpha =0.45$. Results based on the Néel state as CCM model state are shown in $\mathrm{LSUB}n$ approximations with $n=2,4,6,8,10$, together with two corresponding $\mathrm{LSUB}\infty \left(i\right)$ extrapolated results using Eq. (23) and the respective data sets $n=\{2,6,10\}$ for $i=1$ and $n=\{4,6,8,10\}$ for $i=2$. In Fig. 2 we also show, for comparison, an $\mathrm{LSUB}\infty \left(2a\right)$ extrapolation based on Eq. (22) and the data set $n=\{4,6,8,10\}$.

Figure 3

CCM results for the GS magnetic order parameter $M$ vs the scaled interlayer exchange coupling constant, $\delta \equiv J_1^{\perp }/J_1$, for the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the bilayer honeycomb lattice (with $J_3=J_2$ and $J_1>0)$, for three selected values of the intralayer frustration parameter, $\alpha \equiv J_2/J_1$: (a) $\alpha =1.0$, (b) $\alpha =0.8$, and (c) $\alpha =0.56$. Results based on the striped state as CCM model state are shown in $\mathrm{LSUB}n$ approximations with $n=2,4,6,8,10$, together with two corresponding $\mathrm{LSUB}\infty \left(i\right)$ extrapolated results using Eq. (23) and the respective data sets $n=\{2,6,10\}$ for $i=1$ and $n=\{4,6,8,10\}$ for $i=2$.

Figure 4

CCM results for the GS magnetic order parameter $M$ vs the scaled interlayer exchange coupling constant, $\delta \equiv J_1^{\perp }/J_1$, for the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the bilayer honeycomb lattice (with $J_3=J_2$ and $J_1>0)$, for a variety of values of the intralayer frustration parameter, $\alpha \equiv J_2/J_1$, using (a) the Néel state and (b) the striped state as the CCM model state. In each case we show extrapolated results, obtained from using Eq. (23) with the corresponding $\mathrm{LSUB}n$ data sets $n=\{2,6,10\}$.

Figure 5

CCM results for the triplet spin gap $\mathrm{\Delta }$ (in units of $J_1)$ vs the scaled interlayer exchange coupling constant, $\delta \equiv J_1^{\perp }/J_1$, for the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the bilayer honeycomb lattice (with $J_3=J_2$ and $J_1>0)$, for three selected values of the intralayer frustration parameter, $\alpha \equiv J_2/J_1$: (a) $\alpha =0$, (b) $\alpha =0.2$, and (c) $\alpha =0.45$. Results based on the Néel state as CCM model state are shown in $\mathrm{LSUB}n$ approximations with $n=2,4,6,8,10$, together with two corresponding $\mathrm{LSUB}\infty \left(i\right)$ extrapolated results using Eq. (24) and the respective data sets $n=\{2,6,10\}$ for $i=1$ and $n=\{4,6,8,10\}$ for $i=2$.

Figure 6

CCM results for the triplet spin gap $\mathrm{\Delta }$ (in units of $J_1)$ vs the scaled interlayer exchange coupling constant, $\delta \equiv J_1^{\perp }/J_1$, for the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the bilayer honeycomb lattice (with $J_3=J_2$ and $J_1>0)$, for three selected values of the intralayer frustration parameter, $\alpha \equiv J_2/J_1$: (a) $\alpha =1.0$, (b) $\alpha =0.8$, and (c) $\alpha =0.56$. Results based on the striped state as CCM model state are shown in $\mathrm{LSUB}n$ approximations with $n=2,4,6,8,10$, together with two corresponding $\mathrm{LSUB}\infty \left(i\right)$ extrapolated results using Eq. (24) and the respective data sets $n=\{2,6,10\}$ for $i=1$ and $n=\{4,6,8,10\}$ for $i=2$.

Figure 7

$T=0$ phase diagram of the spin-$\frac{1}{2}{J}_1\text{−}{J}_2\text{−}{J}_3\text{−}{J}_1^{\perp }$ model on the $AA$-stacked bilayer honeycomb lattice with $J_1>0,J_3=J_2\equiv \alpha J_1>0,J_1^{\perp }\equiv \delta J_1>0$. The leftmost (sky blue) and the rightmost (sea green) regions are the quasiclassical AFM phases with Néel and striped order respectively, while the central (gray) region is a gapped paramagnetic phase. The red cross $(\times )$ and the green plus $(+)$ symbols are points at which the extrapolated GS magnetic order parameter $M$ vanishes for the two respective quasiclassical phases. They represent the respective values ${\alpha }_{c_1}\left(\delta \right),{\alpha }_{c_2}\left(\delta \right)$ and and ${\delta }_{c_1}^{>}\left(\alpha \right),{\delta }_{c_2}^{>}\left(\alpha \right)$ [and also ${\delta }_{c_1}^{<}\left(\alpha \right),{\delta }_{c_2}^{<}\left(\alpha \right)$ for values of $\alpha$ in the range ${\alpha }_{c_1}\left(0\right)<\alpha <{\alpha }_1^{>}$ and ${\alpha }_2^{<}<\alpha <{\alpha }_{c_2}\left(0\right)$, respectively]. In each case the Néel and the striped states are used as CCM model states, and Eq. (23) is used for the extrapolations with the respective $\mathrm{LSUB}n$ data sets $n=\{2,6,10\}$ used as input. (See also Sec. 5 for a discussion of any possible errors in the phase diagram.)