Roles of easy-plane and easy-axis XXZ anisotropy and bond alternation on a frustrated ferromagnetic spin-$1/2$ chain

The spin-$1/2$ Heisenberg chain with a ferromagnetic first-neighbor exchange coupling $J_1$ and an antiferromagnetic second-neighbor $J_2$ has a Haldane dimer ground state and an extremely small spin gap. Thus, the ground state is readily altered by perturbations. Here, we investigate the effects of XXZ exchange magnetic anisotropy of both the easy-axis and easy-plane types and an alternation in $J_1$ on the ground state, the spin gap, and magnetic properties of the frustrated ferromagnetic spin-$1/2$ chain. It is found that there are two distinct dimerized spin-gap phases, in one of which the spin gap and the magnetic susceptibility are extremely small around the SU(2) symmetric case and in the other they are moderately large far away from the SU(2) symmetric case. A small alternation in the amplitude of $J_1$ rapidly shortens the pitch of spin correlations towards the four-spin periodicity, as in the limit of $J_1/J_2\to0$. These effects are not sufficient to quantitatively explain overall experimentally observed magnetic properties in the quasi-one-dimensional spin-gapped magnetoelectric cuprate Rb$_2$Cu$_2$Mo$_3$O$_{12}$ that exhibits ferroelectricity stabilized by a magnetic field. Our results are also relevant to Cs$_2$Cu$_2$Mo$_3$O$_{12}$, where the ferromagnetic intrachain and antiferromagnetic interchain order has recently been found, in a single chain level. We also reveal the nature of symmetry-protected topological phase transitions in the model by mapping onto effective spin-1 chain models.

The spin-1/2 Heisenberg chain with a ferromagnetic first-neighbor exchange coupling J1 and an antiferromagnetic second-neighbor J2 has a Haldane dimer ground state and an extremely small spin gap. Thus, the ground state is readily altered by perturbations. Here, we investigate the effects of XXZ exchange magnetic anisotropy of both the easy-axis and easy-plane types and an alternation in J1 on the ground state, the spin gap, and magnetic properties of the frustrated ferromagnetic spin-1/2 chain. It is found that there are two distinct dimerized spin-gap phases, in one of which the spin gap and the magnetic susceptibility are extremely small around the SU(2) symmetric case and in the other they are moderately large far away from the SU(2) symmetric case. A small alternation in the amplitude of J1 rapidly shortens the pitch of spin correlations towards the four-spin periodicity, as in the limit of J1/J2 → 0. These effects are not sufficient to quantitatively explain overall experimentally observed magnetic properties in the quasi-one-dimensional spin-gapped magnetoelectric cuprate Rb2Cu2Mo3O12 that exhibits ferroelectricity stabilized by a magnetic field. Our results are also relevant to Cs2Cu2Mo3O12, where the ferromagnetic intrachain and antiferromagnetic interchain order has recently been found, in a single chain level. We also reveal the nature of symmetry-protected topological phase transitions in the model by mapping onto effective spin-1 chain models.

I. INTRODUCTION
The frustrated spin-1/2 XXZ chain with a ferromagnetic first-neighbor exchange coupling J 1 < 0 and an antiferromagnetic second-neighbor exchange coupling J 2 > 0 has attracted considerable interest both for realizing nontrivial phases, including the vector-spin-chirality ordered phase 1 and the Haldane dimer phase 2 , and for the relevance to quasi-one-dimensional edge-sharing multiferroic cuprates 1 , such as A 2 Cu 2 Mo 3 O 12 (A=Rb, Cs) 3-8 , LiCu 2 O 2 9-12 , SrCuO 2 13 , LiCuVO 4 14-18 , LiCuSbO 4 19 , Li 2 CuZrO 4 20 , and PbCuSO 4 (OH) 2 21,22 . The Hamiltonian is given by with the spin-1/2 operator S i at the site i. The model has been studied intensively and extensively in the Heisenberg (∆ xy = ∆ z = 1) 2,23-32 and easyplane (∆ xy = 1 and 0 < ∆ z < 1) cases 1,2,24,33-38 . For J 1 /J 2 −4, the ground state belongs to a Tomonaga-Luttinger liquid (TLL) phase with a quasi-long-range ferromagnetic order 24 . For J 1 /J 2 −4, on the other hand, the model exhibits various phases, depending on ∆ xy . The ground state belongs to the Haldane dimer phase 2 , which has been labeled as the D + phase 39 , around the SU(2)-symmetric case. The subscript + or -denotes the relative sign of the xy, and z components of the dimer order parameters and is opposite to the parity eigenvalue of the ground state 39 . The spin gap in this D + phase is orders of magnitude smaller 1,2,26,31,32 than in another D + phase realized for antiferromagnetic J 1 32, [40][41][42][43][44][45][46] . With increasing easy-plane exchange anisotropy, namely, decreasing ∆ xy from unity, a gapless vector-chiral (VC) phase appears robustly 1 . In practice, these states are susceptible to other weak perturbations. For instance, an infinitesimally small alternation δ in the amplitude of J 1 , as described by replaces the gapless VC phase with two topologically distinct vector-chiral dimer (VCD + and VCD − ) phases separated by a gapless VC phase boundary 39,47 . Unfrustrated three-dimensional interchain interactions readily leads to a long-range spiral magnetic order, as is the case in multiferroic cuprates 1, 9,10,[14][15][16]21,22 . With a further decrease in ∆ xy , there appears another dimer (D − ) phase, which can have a larger spin gap than in the D + phase. The coexistent phases, i.e., the VCD + and VCD − phases, survive in narrow regions even for δ = 0 1,2 . For more details, see Ref. 2 in the case of δ = 0 and Refs. 39,47 in the case of δ = 0. In contrast to the easy-plane case, the global phase diagram of the Hamiltonian Eq. (1) in the easy-axis (0 < ∆ xy < 1 and ∆ z = 1) case has not been fully clarified yet. It has been known that a fully polarized ferromagnetic phase (FPF) appears when J 1 /J 2 < −4 23,24 . In the Ising limit ∆ xy = 0, a first-order phase transition occurs at J 1 /J 2 = −2 from the FPF phase to an up-up-downdown (UUDD) (· · · ↑↑↓↓ · · · ) antiferromagnetic phase 48 .
A weak inplane exchange interaction, namely, ∆ xy ≪ 1, induces a partially polarized ferromagnetic (PPF) TLL phase between the FPF and UUDD phases 24,48 . Exactdiagonalization calculations 24 and a bosonization analysis 2 have revealed that the leading correlations in this PPF TLL phase are given by multi-magnons or multispinons, as in the associated Heisenberg model under a magnetic field [49][50][51] . However, the phase diagram with an intermediate region of ∆ xy remains open both with and without the bond alternation δ, and its clarification is one of the two main goals of this paper.
Although the J 1 (< 0)-J 2 (> 0) spin-1/2 chain model should certainly be relevant to quasi-one-dimensional edge-sharing cuprates, quantitative comparisons are not necessarily easy. In particular, since the spin gap is extremely small or even vanishes near the SU(2)-symmetric case, as we have mentioned above, various perturbations may critically alter the ground state and/or excitations of the system. For instance, recent experiments combined with theoretical studies 52 have indicated an emergent spin-1 Haldane gap in the quasi-one-dimensional frustrated ferromagnetic spin-1/2 magnet Rb 2 Cu 2 Mo 3 O 12 3 . It has been shown that quantitative explanations of magnetic properties of the compound demand a two-leg ladder model comprising of antiferromagnetically coupled J 1 -J 2 frustrated spin-1/2 chains with moderately large Dzyaloshinskii-Moriya interactions. Namely, the single frustrated ferromagnetic chain is not sufficient. It is the other main goal of this paper to reveal the effects of various perturbations within a single J 1 -J 2 chain on experimentally observable quantities. In particular, we investigate the effects of easy-plane and easy-axis exchange magnetic anisotropy and the bond alternation δ on the spin gap, the periodicity of dominant spin correlations, and the uniform magnetic susceptibility. These results will be useful for direct comparisons with experiments.
The rest of the paper is organized as follows. In Sec. II, we present the global phase diagram of the Hamiltonian H = H XXZ + H δXXZ and the maps of the spin gap, the periodicity of dominant spin correlations, and the uniform transverse magnetic susceptibility. Our results are consistent with the previous results in already known cases with both easy-plane 1,2,39,47 and easy-axis anisotropy 23,24,48 . In Sec. III, we examine in detail phase transitions in the case of easy-axis anisotropy. In Sec. IV, analytical expressions of the spin gap are derived for the FPF phase that appears with easy-axis anisotropy and for the UUDD phase in the Ising limit. We also introduce a mapping onto effective spin-1 XXZ chain models, starting from the strongly dimerized limit. On the basis of this mapping, we elucidate the nature of the symmetry protected topological (SPT) phase transition of the Gaussian universality class between (VC)D ± phases and the continuous phase transition of the Ising universality class between the UUDD phase and the D + phase. Lastly, in Sec. V, we provide discussion and the conclusions, precluding single-chain scenarios for Rb 2 Cu 2 Mo 3 O 12 . Possible relevance to Cs. 2 Cu 2 Mo 3 O 12 is also discussed.
II. GLOBAL GROUND-STATE PHASE DIAGRAM AND MAGNETIC PROPERTIES Figure 1 presents the main results on the global phase diagram, the spin gap ∆ G , the wave number q max of the maximum spin-spin correlation, and the transverse magnetic susceptibility χ x for the Hamiltonian H = H XXZ + H δXXZ . Properties of each phase are summarized in Table I and will be explained below.

A. Fully polarized ferromagnetic (FPF) phase
The FPF phase 48 has the unique order parameter of the uniform magnetization which is pinned to 1/2, where · · · denotes the groundstate expectation value with L being the total number of sites. As we will show in Sec. IV B, the FPF phase emerges for J 1 /J 2 < 2/(1 − δ) in the Ising limit ∆ xy = 0 (∆ z = 1). The lower bound of |J 1 /J 2 | for the FPF phase monotonically increases with increasing ∆ xy for δ = 0, 0.02, 0.1, and 0.2, and reaches J 1 /J 2 = −2[(1 + δ)/(1 − δ) + 1] at ∆ xy = ∆ z = 1 31 . Early studies 23,24,48 have shown that the instability of the FPF state occurs towards multi-magnon bound states, and that the number of magnons forming the bound state increases with increasing ∆ xy . In the FPF phase, the spin gap ∆ G is given by the energy difference between the lowest-energy states with S z T = L 2 and L 2 − 1. As we will show in Sec. IV A, the analytic form of ∆ G near ∆ xy = 0 is given by ∆ z |J 2 |(1 + J 1 /J 2 ) in the Ising limit and gradually decreases with increasing ∆ xy . The value of the wave number k inc of this first-excited state with S z T = L 2 − 1 depends on J 1 /J 2 as (1 − δ 2 ) − 1+δ 2 1−δ 2 < 1, as will be derived in Sec. IV A. The wave number q max at which the spin-spin correlation shows the maximum among α = x, y, z, trivially vanishes, namely, q max = 0, because of the emergence of the uniform magnetization. We compute the ground-state transverse magnetic susceptibility χ x from with the transverse magnetic field h x = 0.001J 2 . It decreases with decreasing ∆ xy in this phase. Nevertheless, it remains as large as 0.37 for δ = 0, 0.02, 0.1, and 0.2 with (J 1 /J 2 , ∆ xy ) = (−3.28, 0.4), as confirmed numerically [ Fig. 1 (d)]. In the limit of J 2 = ∆ xy = δ = 0, the Hamiltonian is equivalent to one-dimensional Ferromagnetic Ising chain, then the zero-field susceptibility of the transverse Ising chain χ x = lim h x →0 S x T Lh x can be obtained exactly as J 1 χ x = 0.5 54-56 .

B. Partially polarized ferromagnetic (PPF) phase
In the PPF phase 23 , the magnetization continuously changes within 0 < M < 1/2 as a function of J 1 /J 2 and ∆ xy . This phase is described as a single channel TLL and thus has no spin gap. The bosonaization analysis 2 reveals that the transverse spin correlation function S + j S − j ′ decays exponentially with respect to the distance, and the longitudinal spin correlation function S z j S z j ′ and the bond nematic correlation 49,50 vanishes for δ = 0 2 , as also confirmed by our numerical calculations. Once we introduce δ > 0, M and D α may appear simultaneously, because the finite δ introduces relevant cosine terms of bosonic fields into the Hamiltonian 2 and can shift the field-locking positions. As shown in Fig. 1 (a), the area of PPF phase is strongly suppressed by increasing the bond alternation δ because the UUDD state has a large energy gain with respect to δ, as we will show in Sec. IV B and Sec. IV C. The PPF phase has q max = 0, as in the FPF phase. It also has a large transverse magnetic susceptibility, for instance, J 2 χ x = 0. 26 The order parameter O uudd characterizing the UUDD phase is defined by For ∆ xy = 0, the UUDD phase emerges for The D + phase is characterized by (D x + D y )D z > 0, and does not have any local magnetic order. In this phase, effective spin-1 degrees of freedom emerge on the bonds with stronger ferromagnetic correlation, forming a valence bond solid state 57 as in the spin S = 1 Heisenberg chain 58,59 . The pattern of the dimer and thus the stronger ferromagnetic bonds has the twofold degeneracy when δ = 0, but is fixed by finite δ that breaks the one-site translational symmetry and doubles the unit cell. In particular, the weaker ferromagnetic bonds are shortrange entangled. This phase appears around the SU(2) symmetric case, as shown in Fig. 1(a). It is robust against the XXZ anisotropy ∆ xy = ∆ z and the bond alternation δ, because of the Z 2 × Z 2 , time-reversal, and bond-center inversion symmetries protecting the topological property of the D + phase 39 .
The phase diagram and the spin gap ∆ G in the SU(2)symmetric case (∆ xy = ∆ z ) have been studied in detail 2,26,31,32 : ∆ G is extremely small for δ = 0 and monotonically increases with increasing |δ|, and the first-order phase transition occurs between the TLL and D + phases The wave number q max of the maximum spin correlation can be incommensurate in the D + phase and evolves from zero towards π/4 with decreasing |J 1 /J 2 | and increasing δ [ Fig. 1 (c)]. It is notable that the transverse susceptibility is tiny only in the D + phase, as shown in Fig. 1 (d), and it merely amounts to J 2 χ x = 7.9 × 10 −2 at most for the parameter line of (J 1 /J 2 , δ) = (−2.2, 0), because the ground state is adiabatically connected to that in the SU(2) case where χ x = 0, without closing the spin gap.

E. Even-parity dimer (D − ) phase
The D − phase 1,2,24,33,38,60 is characterized by (D x + D y )D z < 0, and also does not have any local magnetic order. This phase shares the same topological properties with the D + phase, except that the stronger ferromagnetic bonds are short-range entangled. As in the case of the D + phase, the doubly-degenerated ground states emerge in the D − phases for δ = 0, and the degeneracy is lifted for δ = 0. This phase appears on both the easy-plane and easy-axis sides.
Let us start with the easy-plane side (∆ xy = 1, ∆ z < 1). It has already been shown that the D − phase appears in a rather large region, typically ∆ z 0.6 1 . While this phase appears in −3 J 1 /J 2 −1 in the case of δ = 0, it expands appreciably and occupies a large portion of the phase diagram in the case of δ = 0. (See Fig. 1 (a).) Deeply inside the D − phase, the spin gap ∆ G and the transverse magnetic susceptibility χ x increase up to 0.4J 2 and 0.72/J 2 for (J 1 /J 2 , ∆ z , δ) = (−1.66, 0, 0) and up to 1.4J 2 and 0.88/J 2 for (J 1 /J 2 , ∆ z , δ) = (−4, 0, 0.2), respectively. Thus, the ground state can show a large magnetization by applying a small transverse magnetic field. Now we turn to the easy-axis side (∆ xy = 1, ∆ z < 1). Actually, the presence of the D − phase has for the first time been uncovered in the easy-axis case. As shown in Fig. 1, it appears in a narrow region surrounded by the PPF phase, the D + phase, and the UUDD phase. With increasing |δ|, the region quickly narrows. The spin gap ∆ G is not larger than 0.03J 2 , since the phase is rather narrow and sandwiched by a gapless PPF phase and a critical phase boundary with the D + phase [ Fig. 1 (a)]. Although the spin gap can be small as in the D + phase, the transverse magnetic susceptibility J 2 χ x is typically an order of magnitude larger (∼ 0.2 for (J 1 /J 2 , ∆ xy , δ) = (−2.2, 0.75, 0) ) than in the D + phase.
Namely, on both the easy-plane and easy-axis sides, in general, the observation on the transverse magnetic susceptibility provides a key to discriminate the D − phase from the D + phase in experimental observations.
On both the easy-plane and easy-axis sides, the behavior of the wave number q max of the maximum spin correlation is similar to that in the D + phase: it can be incommensurate and evolves from zero towards π/4 with decreasing |J 1 /J 2 | and increasing δ [ Fig. 1 (c)].

F. Gapless vector-chiral (VC) phase
The gapless VC phase is characterized by a long-range order of the uniform vector-chirality 61 and a single-channel TLL showing a quasi-long-range inplane spiral spin correlation 45 . This phase occupies a wide region of the phase diagram with easy-plane anisotropy (∆ xy = 1, ∆ z < 1) in the case of δ = 0 [ Fig. 1 (a)] 1,2 . Then, the gapless VC phase is immediately replaced with (vector-chiral) dimer phases by introducing nonzero δ 39,47 . The wave number q max of the maximum spin correlations evolves from ∼ 0.7 to π/2 for ∆ z = 0.85 as |J 1 |/J 2 decreases 1,2 , as shown in Fig. 1 (c).

G. Vector-chiral dimer (VCD±) phases
The VCD ± phases have the relative sign ± of the dimer order parameters D x and D z . They appear as coexistent phases of the vector-chirality order and the dimer order on the easy-plane side 1,2,39,47 . In the case of δ = 0, they are restricted to narrow regions sandwiched by the gapless VC phase and the D ± phases. Turning on finite δ = 0 immediately replaces the gapless VC phase with the VCD ± phases, except at the boundary of the VCD ± phases. Then, with increasing δ, the areas of the VCD ± phases gradually decrease. Behaviors of ∆ G , q max , and χ x in the VCD ± phases are similar to those in the D ± phases.

H. Tomonaga-Luttinger liquid (TLL) phase
The TLL phase in the case of easy-plane anisotropy with δ = 0 [ Fig. 1 (a)] has already been investigated intensively 1,2,24,33,36,60,62 and does not have any longrange order, for instance, of M , D α , O uudd , κ, and O N , but shows a quasi-long-range transverse ferromagnetic order and has gapless excitations [ Fig. 1 (b)]. The maximum peak position q max is always zero, as shown in Fig. 1 (c). Reflecting the quasi-long-range order, the transverse magnetic susceptibility χ x diverges [ Fig. 1 (d)]. This TLL phase is immediately replaced with dimer phases by introducing finite δ.

I. Néel (N) phase
The N phase 37 is characterized by the staggered magnetization The spin gap ∆ G is small and not larger than 0.03J 2 , since the phase is surrounded by the gapless VC phase, the TLL phase, and the D − phase [ Fig. 1 (a)]. This phase has a large transverse magnetic susceptibility up to 2 × 10 2 for (J 1 /J 2 , ∆ z ) = (−3.64, 0.85). Because of the Néel order, the maximum spin correlation occurs at q max = π in the z component of spins, while the shortrange spin correlations in the xy component change from commensurate (q = 0) to incommensurate (q > 0) inside this phase 37,63,64 . This N phase is replaced with the (chiral) dimer (VC)D ± phases by finite δ, as is clear from Fig. 1(a).

J. Vector-chiral Neel (VCN) phase
From the bosonization analysis and numerical calculations in Ref. [2 and 47], the vector-chirality order and the Néel order can emerge simultaneously near the phase boundary among VC, N, D ± , and the TLL phases. This VCN phase has an extremely small spin gap ∆ G < 0.03J 2 and a commensurate pitch of q max = π [ Fig. 1 (b) and (c)]. The behavior of χ x in the VCN phase is similar to that in the Néel phase [ Fig. 1 (d)]. We did not observe the VCN phase for δ ≥ 0.02.

III. PHASE TRANSITIONS ON THE EASY-AXIS SIDE
This section is devoted to a numerically precise determination of phase boundaries in the phase diagrams on the easy-axis side [see upper panels of Fig. 1 (a)] and to a clarification of the nature of the phase transitions. For this purpose, we performed iTEBD 53 and/or iDMRG 65 calculations along several typical vertical and horizontal lines in the phase diagrams. In Fig. 2, the results of relevant order parameters are shown as functions of J 1 /J 2 for (∆ xy , δ) = (0.1, 0) [ Fig. 2 (a)], (0.5, 0) [ Fig. 2 Fig. 3 (b)]. In the following, we will explain the results separately for each of the phase transitions.

A. FPF-UUDD transition
It has already been known that a first-order phase transition occurs at J 1 /J 2 = −2 between the FPF and UUDD phases in the case of ∆ xy = δ = 0 48 . The δ dependence of the FPF-UUDD phase transition point is exactly obtained in the Ising limit as from a crossing of the energy levels of the FPF state (· · · ↑↑↑↑ · · · ) and the UUDD state (· · · ↑↑↓↓ · · · ), as we will explain in detail in Sec. IV B. Equation (10) indicates that with an increase in the bond alternation |δ|, the UUDD phase expands while the FFP phase narrows. On the basis of our numerical results shown in the upper panels of Fig. 1(a), this first-order phase transition is likely to survive a finite but small transverse interaction up to ∆ xy ∼ O(δ), though massive calculations are required for confirming the conjecture. With further increasing ∆ xy within ∆ xy ≪ 1 , the direct FFP-UUDD phase transition disappears and the FPF and UUDD phases are intervened by the PPF phase.

B. FPF-PPF transition
Phase transitions between the FPF and PPF phases are observed from a change in the magnetization M in all the panels of Figs. 2 and in Fig. 3 (a). In the case of δ = 0, this transition has been described as a continuous phase transition associated with the condensation of bound multi-magnons 2,23,24,48 . This mechanism also holds when δ = 0. Note also that a similar mechanism holds in the spin S = 1/2 XXZ model around the saturated magnetic field [49][50][51] . For δ = 0, our results on discontinuity in the relevant order parameter O uudd at the PPF-UUDD transition, as shown in Fig. 2(a). As the PPF phase narrows with an increase in |δ|, the direct PPF-UUDD phase transition is also shortened. The D − phase eventually intervenes the PPF and UUDD phases.

D. PPF-D− transition
The PPF-D − phase transition [Figs. 2(c),(e)-(h) and Fig. 3(a)] is also a first-order phase transition. This is confirmed by the discontinuity of the magnetization M at the PPF-D − transition, as shown in Fig. 3(a).

E. D−-UUDD transition
The D − -UUDD phase transition [Figs. 2(c),(e), and (g)] emerges in the case of |δ| > 0. The UUDD order parameter O uudd rapidly drops towards the D − -UUDD phase transition with either a second-order phase transition or a weakly first-order. Note that the D − phase and the UUDD phase are corresponds to Large-D phase and Néel phase in an effective S = 1 Hamiltonian for |δ| ∼ 1, as we will show in Sec. IV C. The Large-D-Néel phase transition in the effective model is known as a first-order phase transition 66 . As long as the analogy holds, this phase transition could be of the first order.

F. D+-UUDD transition
The mechanism of the D + -UUDD phase transition [ Fig. 2(d), (f) and (h), and Fig. 3(b)] is also understood by the analysis of the effective S = 1 Hamiltonian for |δ| ∼ 1 in Sec. IV C. Actually, the D + -UUDD transition corresponds to the Haldne-Néel transition 66 in the effective S = 1 Hamiltonian. Thus, it could belong to the Ising universality class, as long as the analogy holds. Indeed, the order parameter O uudd shows the Ising critical behavior (∆ xy c − ∆ xy ) 1/8 , as shown in the inset of Fig. 3(b).

G. D−-D+ transition
The D − -D + phase transition [ Fig. 2 Fig. 3(a)] has the same mechanism as that on the easyplane side, and is a second-order phase transition 2,47 . This phase transition is associated with the sign change of D x = D y (D z ) on the easy-axis (-plane) side. A bosonization analysis is also available 2,47 . At the transition, there is a critical divergence of the correlation length of the matrix-product state (MPS), as shown in Fig. 3(a). This points to a continuous phase transition. From the analysis of an effective S = 1 Hamiltonian for |δ| ∼ 1, we show in Sec. IV C that this D − -D + phase transition belongs to the Gaussian universality class, and is a canonical symmetry protected topological (SPT) phase transition protected by the time-reversal symmetry, bondinversion symmetry, and Z 2 × Z 2 symmetry 67,68 .

H. FPF-D+ transition
It has been reported that a first-order phase transition occurs between the FPF phase and D + phase at J 1 /J 2 = −2[(1 + δ)/(1 − δ) + 1] in the SU(2) limit 31 . At a glance of Fig. 1(a), it is expected that this direct FPF-D + phase transition survives tiny easy-axis anisotropy, though more intensive calculations are required for the confirmation. With increasing ∆ xy , this transition eventually bifurcates into the FPF-D − and D − -D + transitions.

A. Spin gap in the FPF phase
Here, we assume that the number L of sites in the system is even, L = 2n, when the bond alternation δ is finite. However, the final results on the spin gap holds even when δ = 0 and L is odd.
The FPF ground state for even L is written as The ground state energy is given by The excited states specified by the z component i S z = n − 1 of the total spin and the wave number k ∈ {0, 2π/n, · · · , 2π(n − 1)/n} can be expressed as with coefficients c ik satisfying i |c ik | 2 = 1, where form a set of orthonormal bases. We have also introduced a translation operatorT througĥ Since tha Hamiltonian H commutes withT 2 , we obtain We can compute H|v 2 in the same manner as Eq. (16) and Eq. (17), except that | ↑↓ and | ↓↑ are swapped and T 2 in the third line of Eq. (17) is replaced withT −2 . Using Eqs. (16) and (17) as well as the relation v i |T 2j |v i ′ = δ i,j e −kj for i, i ′ = 1, 2, the matrix elements of H − E F are obtained as This matrix has two eigenvalues The spin gap ∆ G can then be calculated as the lowest eigenvalue min k ε −,k . In addition to two trivial local minima, there may exist another local minimum at an incommensurate wave number k inc in the thermodynamic limit. This solution can be obtained by differentiating ε −,k with respect to k as leading to (26) in the case of δ 2 ≤ 1 and −1 < (1 − δ 2 ) − 1+δ 2 1−δ 2 < 1. The result ∆ G = min(ε −,0 , ε −,π , ε −,±kinc ) is plotted in the FPF phase of Fig. 1(b). In particular, for δ = 0, ε −,−kinc yeilds the spin gap ∆ G for −4 < J 1 /J 2 < 2, while ε −,π for J 1 /J 2 ≤ −4.

B. Spin gap in the UUDD phase in the Ising limit
Here, we assume that the number L of sites in the system is a multiple of 4, L = 4n, when the bond alternation δ is positive. The UUDD state in the Ising limit is written as The ground state energy under the periodic boundary condition reads In the case of J 1 /J 2 > −1 and δ = 0, the lowest excited states within the S z T = 1 manifold are expressed as {T j |σ a r } 0≤j<4n 0≤r<n with |σ a r = | ↑↑↓ ⊗ r l=1 | ↑↑↓↓ ⊗ | ↑ ⊗ n−r−1 whereT j is the operator that translated the spins by j sites. Then, the spin gap ∆ G is given by On the other hand, in the case of −2 < J 1 /J 2 < −1 and δ = 0, the lowest excited states within the S z T = 1 manifold are given by Then, the spin gap is given by Turning on the positive δ, the degeneracy of the two manifolds {T j |σ a r } and {T j |σ b uvw } is lifted. Then, the lowest-energy states within the S z T = 1 manifold are given by {T 2k |s a ℓ } 0≤k<2n ℓ∈{0,1} for J 1 /J 2 > −1/(1 − δ) [ region A in Fig. 4 ] with |s a 1 = | ↑↑↑↓ ⊗ n−1 l=1 | ↑↑↓↓ (34) and Fig. 4 ] with respectively. Then, the spin gap is given by .
(37) Thus, the effect of finite δ introduces a non-monotonic behavior of the spin gap with respect to J 1 /J 2 , and the result of ∆ G is plotted in Fig. 1(b).
The global phase diagram ofH θ forJ xy > 0 has already been investigated by the numerical exact diagonalization with finite-size scaling analyses 66 . It contains the Néel, Haldane, and Large-D phases forJ xy θ > 0 ∧ J z θ > 0. These phases corresponds to the UUDD, D + (VCD + ) and D − (VCD − ) phases for θ = 0 (θ > 0), respectively, in the original spin-1/2 Hamiltonian H. Then, it is reported that the Haldane-Néel phase transition and the Haldane-Large-D phase transition are continuous phase transitions that belong to the Ising universality class (the central charge c = 1/2) and to the Gaussian universality class (c = 1), respectively 66 . Indeed, our current numerical finding on the Ising universality class for the D + -UUDD phase transition on the easy-axis side of the original spin-1/2 HamiltonianH [Sec. III F] agrees with that of the Haldane-Néel phase transition forH θ with θ = 0. Furthermore, our previous numerical finding on the Gaussian universality class for the VCD + -VCD − phase transition in the easy-plane side of H 39 can also be elucidated from that of the Haldane-large-D phase transition forH θ with finite θ. The D + -D − phase transition that most likely appears for small |J 1 /J 2 | in the case of large |δ| in the spin-1/2 Hamiltonian H with the easyplane anisotropy should also belong to the Gaussian universality class from the analogy to the Haldane-Large-D phase transition forH θ with θ = 0.

V. DISCUSSION AND CONCLUSIONS
We have investigated the ground-state properties and the spin gap of the spin-1/2 frustrated spin XXZ chain with the ferromagnetic first-neighbor coupling J 1 and the antiferromagnetic second-neighbor coupling J 2 , both with and without an alternation in J 1 . Using the iTEBD and iDMRG methods, we have numerically completed the ground-state phase diagram as well as the maps of the spin gap, the wave number of the maximum spin correlation, and the transverse magnetic susceptibility on both the easy-plane and easy-axis sides. We have also derived analytic expressions for the spin gap in the fully polarized ferromagnetic phase (FPF) and in the up-updown-down (UUDD) phase near and in the Ising-limit, respectively, in the case with a finite bond alternation δ in the first-neighbor exchange coupling. These maps, summarized in Fig. 1 and Table I, will be useful for understanding low-temperature experimental results on the basis of a single frustrated J 1 -J 2 spin-1/2 chain.
In particular, our numerical results has for the first time uncovered the anisotropic dimer (D − ) phase on the side of easy-axis magnetic anisotropy. We have shown that this phase is surrounded by the partially polarized ferromagnetic (PPF) phase, the up-up-down-down (UUDD) phase, and the (nearly) isotropic dimer phase (D + ). Clear first-order phase transitions have been found at the D − -PPF and PPF-UUDD phase transitions have been found to be clearly of the first order. The D − -UUDD phase transition either is weakly first-order or belongs to the Ising universality class. The D − -D + phase transitions is of the second order and belongs to the Gaussian universality class with the central charge. The nature of D − -UUDD and D − -D + phase transitions has also been explained as the Haldane-Large-D and Haldane-Néel phase transitions, respectively, by mapping the original spin-1/2 model in the strongly dimerized case onto an effective spin-1 XXZ chain model with single-ion anisotropy D. Note that the D − -D + phase transition is a symmetry-protected topological phase transition protected by the time-reversal, bond-inversion, and Z 2 × Z 2 symmetries. This analysis for obtaining the Gaussian universality class holds for the VCD − -VCD + phase transition on the side of easy-plane exchange mag-netic anisotropy, which solves the issue posed from previous numerical results finding 39 .
Our results have crucial implications to experiments on low-temperature magnetic properties of quasi-onedimensional spin-1/2 compounds. Large exchange magnetic anisotropy on both the easy-plane and easy-axis sides can yield an energy gap in the spin excitation above a nonmagnetic ground state in the D − and VCD − phases. However, the ground state necessarily possesses a moderately large transverse magnetic susceptibility. This feature is robust and also holds even when the U(1) spin symmetry is absent. It may happen that while the magnetic susceptibility is strongly suppressed as in the D + and VCD + phases, the experimentally measured spin gap in a dimerized phase is too large for that of the D + and VCD + phases. Now, we explain the extent to which single J 1 -J 2 spin-1/2 chain models can and cannot explain experimental findings of the quasi-1D spin-1/2 chain compound Rb 2 Cu 2 Mo 3 O 12 that hosts a nonmagnetic ground state 52 . It is reasonable to start from the SU(2) symmetric model with no bond alternation δ = 0 and then to consider effects of perturbations. The wave number q max ∼ π/4, which corresponds to the eight-spin periodicity, of the maximum spin correlation has been observed with inelastic neutron-scattering experiments. From this value, the ratio J 1 /J 2 is estimated as −3.6. Then, the magnetic susceptibility of the powder samples can be fit in the temperature range from 300 K down to 40 K by taking |J 1 | = 0.62 meV and J 2 = 0.17 meV. The first spin excitation energy has also been measured as ∼ 0.2 meV in inelastic neutron-scattering experiments. This energy gap is too large for the spin gap in the D + phase for δ = 0, since ∆ G < 0.03J 2 in the D + phase with and without easy-plane anisotropy. Attributing this large energy gap to the bond alternation demands a rather large value δ = 0.1. However, with this value of δ, q max approximates to π/2, in stark contrast to the neutron result of q max ∼ π/4. There is another way to enhance the energy gap while keeping the relation q max ∼ π/4 intact. Increasing easy-plane anisotropy, namely, decreasing ∆ z , VCD ± and D − phases appear. The VCD + phase is also excluded since the spin gap is too small as in the D + phase. In the VCD − and D − phases, the spin gap may increase up to ∆ G = 0.93J 2 , for instance, for (J 1 /J 2 , ∆ z , δ) = (−2.02, 0.0, 0.2). Thus, if we take J 2 = 0.21 meV, we reproduce ∆ G = 0.2 meV keeping the wave number q max ∼ π/4. However, the transverse magnetic susceptibility is as large as 0.7 µ B /Cu·T in the ground state. This contradicts to experimental observations of a suppression of the magnetic susceptibility to 0.005 µ B /Cu·T below 0.08 K and the associated magnetization curve that starts to rise steeply only at 2 T. In general, Dzyaloshinskii-Moriya interactions as well as XYZ exchange anisotropy ∆ x = ∆ xy < ∆ y < ∆ z do not increase the first spin excitation energy as long as the spin gap is already finite in the case of ∆ x = ∆ y = ∆ xy , and thus should not be responsible for enhancing the size of the spin gap. Thus, we are naturally led to the next candidate, namely, an unfrustrated interchain interaction within a pair of J 1 -J 2 spin-1/2 chains. It has been shown that it can enhance the spin gap while keeping q max intact. Indeed, the two-leg ladder model of a pair of J 1 -J 2 frustrated spin-1/2 chains that are coupled by an antiferromagnetic rung interaction J ′ can quantitatively explain overall experimental findings. Namely, choosing J 1 = −9.82 meV, J 2 = 3.03 meV, and J ′ = 1.77 meV, and taking rather large Dzaloshinski-Moriya interactions of the order 3.8 meV, the model accounts all of the magnetic susceptibility, the magnetic field strength for closing the spin gap, and inelastic neutron-scattering and electron spin resonance spectra 52 .
Recent neutron diffraction experiments on this material have shown that the spins order ferromagnetically within single chains and antiferromagnetically between the two coupled chains 69 . It is naturally expected that the effective single spin-1/2 chain model lies in the ferromagnetic phases on the easy-axis side or the TLL phase, which shows a quasilong-range ferromagnetic order, on the easy-plane side. As in Rb 2 Cu 2 Mo 3 O 12 , the two adjacent chains should be coupled antiferromagnetically. Actually, this rung interaction does not allow the full spin polarization, since the fully spin polarized Néel state cannot be an eigenstate of the ladder Hamiltonian. If we take a naive estimate of the ratio J 1 /J 2 to be −2.8 4 , it is likely that the compound possesses easy-axis exchange anisotropy, since with easy-plane anisotropy, the single-chain model shows the vector-chiral order, which is readily driven to a long-range spiral magnetic order by three-dimensional interactions of the order of the Néel temperature 1.85 K 5 . It is also possible that the compound has a moderately large J ′ as in Rb 2 Cu 2 Mo 3 O 12 and a larger |J 1 |/J 2 than 2.8. Inelastic neutron-scattering experiments are required for further quantitative theoretical analyses. error of iTEBD calculation There are two origins of calculation errors in the iTEBD method 53 . One is the step size δτ of the imaginary time in the Suzuki-Trotter (ST) decomposition 71,72 . The other is the finite bond dimension χ of the matrixproduct state (MPS). In the limit where χ → ∞ with δτ being fixed, the dominant factor of the error is the time step size δτ of the ST decomposition. In this limit, the error of the total energy ε(δτ, χ) = E(δτ, χ)−E exact , where E(δτ, χ) is the variational energy of iTEBD calculation and E exact is the exact ground state energy, behaves as ε(δτ, χ → ∞) ∼ (δτ ) n+1 (A1) with n being the order of the ST decomposition.
On the other hand, in the limit of δτ → 0 with χ being fixed, the leading factor of the error becomes finiteχ effects in MPS. From the finite-entanglement scaling 73 , we expect that the asymptotic form of the error becomes where we assume the real correlation length is extremely large, compared to a controllable length scale introduced by the finite-χ effects.
Here, we propose a scaling hypothesis ε(δτ, χ) = χ c2 F (δτ χ c1 ) (A3) to satisfy both asymptotic form Eq. (A1) and Eq. (A2), where the function F has the asymptotes, and apply this scaling analysis to the energy error for H XXZ in Eq. (1) with J 1 /J 2 = −5/3 and ∆ xy = ∆ z = 1 31,32 . Then, we confirm that the energy error ε with several χ are nicely on the universal function with appropriate c 1 and c 2 and employ a condition δτ χ c1 ∼ 4, where two asymptotic lines intersect. This provides a condition for efficient calculation of iTEBD, and the suitable step size is found to be δτ ∼ 0.008/J 2 for χ = 300. Under the conditions, the error of the ground-state energy per site in the Heisenberg case (∆ xy = ∆ z = 1) for J 1 /J 2 = −5/3 is found to be ∼ 5 × 10 −6 J 2 .