Mott transition and magnetism of the triangular-lattice Hubbard model with next-nearest-neighbor hopping

TL;DR

This study uses the variational cluster approximation to map out the phase diagram of the triangular-lattice Hubbard model with next-nearest-neighbor hopping, revealing magnetic and nonmagnetic insulating phases and their relation to the Mott transition.

Contribution

It provides a detailed phase diagram considering both nearest and next-nearest neighbor hopping, highlighting the expansion of nonmagnetic insulating phases with correlation strength.

Findings

Identification of 120° Nél and stripe ordered phases

Discovery of a broad nonmagnetic insulating phase

Analysis of the Mott transition via van Hove singularity

Abstract

The variational cluster approximation is used to study the isotropic triangular-lattice Hubbard model at half filling, taking into account the nearest-neighbor ($t_1$) and next-nearest-neighbor ($t_2$) hopping parameters for magnetic frustrations. We determine the ground-state phase diagram of the model. In the strong correlation regime, the 120$^\circ$ N\'eel and stripe ordered phases appear, and a nonmagnetic insulating phase emerges in between. In the intermediate correlation regime, the nonmagnetic insulating phase expands to a wider parameter region, which goes into a paramagnetic metallic phase in the weak correlation regime. The critical phase boundary of the Mott metal-insulator transition is discussed in terms of the van Hove singularity evident in the calculated density of states and single-particle spectral function.

Equations

$$H=$$

$$+U\sum _i{n}_{i\uparrow }{n}_{i\downarrow }-\mu \sum _{i,\sigma }{n}_{i\sigma },$$

$$H=J_1\sum _{⟨i,j⟩}{\mathbf{S}}_i\cdot {\mathbf{S}}_j+J_2\sum _{⟨⟨i,j⟩⟩}{\mathbf{S}}_i\cdot {\mathbf{S}}_j$$

$$\Omega [{\Sigma }^{\prime }]={\Omega }^{\prime }+\mathrm{Trln}(G_0^{-1}-{\Sigma }^{\prime }{)}^{-1}-\mathrm{Trln}(G_0^{\prime -1}-{\Sigma }^{\prime }{)}^{-1},$$

$$H_{\mathrm{M}}=H_{120^{\circ }}+H_{\mathrm{str}}$$

$$H_{120^{\circ }}=h_{120^{\circ }}^{\prime }\sum _i{\mathbf{e}}_{a_i}\cdot {\mathbf{S}}_i$$

$$H_{\mathrm{str}}=h_{\mathrm{str}}^{\prime }\sum _i{e}^{i{\mathbf{Q}}_{\mathrm{str}}\cdot {\mathbf{r}}_i}{S}_i^z,$$

$$\rho (\omega )=\frac{1}{L}\sum _{\mathbf{k}}A(\mathbf{k},\omega )$$

$$A(\mathbf{k},\omega )=-\frac{1}{\pi }\underset{\eta \to 0}{\lim }\Im {\mathcal{G}}_{\mathrm{CPT}}(\mathbf{k},\omega +i\eta )$$

$${\mathcal{G}}_{\mathrm{CPT}}(\mathbf{k},\omega )=\frac{1}{L_c}\sum _{i,j=1}^{L_c}{\mathcal{G}}_{ij}(\mathbf{k},\omega )e^{-i\mathbf{k}\cdot ({\mathbf{r}}_i-{\mathbf{r}}_j)},$$

$$G_{ij}^{\prime }(\omega )$$

$$+⟨{\psi }_0|c_{j\sigma }^{†}\frac{1}{\omega +H^{\prime }-E_0}{c}_{i\sigma }|{\psi }_0⟩,$$

$${\chi }_0(\mathbf{q})=\frac{1}{L}\sum _{\mathbf{k}}\frac{f({\epsilon }_{\mathbf{k}})-f({\epsilon }_{\mathbf{k}+\mathbf{q}})}{{\epsilon }_{\mathbf{k}+\mathbf{q}}-{\epsilon }_{\mathbf{k}}},$$

Full text

I Introduction

The physics of geometrical frustration in strongly correlated electron systems has long attracted much attention anderson ; lee ; balents . In particular, possible absence of magnetic long-range orders at zero temperature in the Heisenberg and Hubbard models defined on frustrated lattices, or the realization of a spin liquid phase as an exotic state of matter, has been one of the major issues in this field. The Mott metal-insulator transition is also a fundamental phenomenon in the field of strongly correlated electron systems mott ; imada , which has attracted much experimental and theoretical interest as well. As one of the simplest models with geometrical frustration and Mott transition, we therefore study the Hubbard model at half filling defined on the triangular lattice in this paper, where not only the nearest-neighbor hopping parameters but also the next-nearest-neighbor ones are included.

Much effort has so far been devoted in the study of the triangular-lattice Hubbard model with anisotropic nearest-neighbor hopping parameters morita ; sahebsara ; watanabe ; ohashi ; tocchio1 ; tocchio2 ; yamada1 ; laubach ; misumi , which was motivated by experimental findings of possible spin liquid states in some organic Mott insulators such as $\kappa$-(ET)2Cu2(CN)3 komatsu ; shimizu ; kurosaki ; manna and EtMe3Sb[Pd(dmit)2]2 itou ; yamashita . The triangular-lattice Heisenberg model with the anisotropic exchange interactions has also been studied to find a variety of ordered phases such as Néel and spiral orders, as well as the quantum disordered (or spin liquid) phases in-between weihong ; hauke1 ; hauke2 .

However, to the best of our knowledge, the isotropic triangular-lattice Hubbard model with both the nearest-neighbor ($t_{1}$) and next-nearest-neighbor ($t_{2}$) hopping parameters has not yet been addressed, the study of which will therefore provide useful information on the physics of magnetic frustrations and Mott metal-insulator transition in strongly correlated electron systems.

The isotropic Heisenberg model with the nearest-neighbor ($J_{1}$) and next-nearest-neighbor ($J_{2}$) exchange interactions, which may be derived by the second-order perturbation of the above-mentioned Hubbard model in the strong correlation limit, has on the other hand been studied much in detail, mostly from the theoretical point of view ma . In the classical Heisenberg model where the spins are treated as classical vectors, it is known that a single phase transition occurs at $J_{2}/J_{1}=1/8$ between the three-sublattice 120∘ Néel ordered state and an infinitely degenerate four-sublattice magnetically ordered states jolicoeur . This degeneracy is lifted by quantum fluctuations, thereby selecting a two-sublattice stripe ordered state through the so-called “order-by-disorder” mechanism jolicoeur ; chubukov ; deutscher ; lecheminant . One may then expect in the corresponding quantum Heisenberg model that an intermediate phase can appear near the classical critical point at $J_{2}/J_{1}=1/8$, for which many studies have been carried out to predict that the nonmagnetic disordered phase bordered by the 120∘ Néel ordered phase at $J_{2}/J_{1}\simeq 0.05-0.12$ and the stripe ordered phase at $J_{2}/J_{1}\simeq 0.14-0.19$ actually emerges. In particular, recent studies actually predict the emergence of either a gapless or gapped spin liquid phase in this intermediate region manuel ; mishmash ; kaneko ; li2 ; zhu ; hu2 ; iqbal . These results of the Heisenberg model may be compared with those of our Hubbard model in the strong correlation limit (as we will see below).

In this paper, motivated by the above developments in the field, we will study the triangular-lattice Hubbard model at half filling with the isotropic nearest-neighbor and next-nearest-neighbor hopping parameters in its entire interaction strength. We use the variational cluster approximation (VCA), one of the quantum cluster methods based on the self-energy functional theory (SFT) potthoff1 ; potthoff2 ; dahnken ; potthoff3 ; senechal , which enables us to take into account the quantum fluctuations of the model with geometrically frustrated spin degrees of freedom. We thereby calculate the grand potential of the system as a function of the Weiss fields for spontaneous symmetry breakings; here, we take the 120∘ Néel and stripe magnetic orders and evaluate the order parameters and critical interaction strengths. We also calculate the charge gap as well as the density of states (DOS) and single-particle spectral function using the cluster perturbation theory (CPT) senechal and determine the ground-state phase diagram of the model in its entire parameter region.

We will thereby show that in the strong correlation regime the 120∘ Néel and stripe ordered phases appear, and in-between, the nonmagnetic insulating phase caused by the quantum fluctuations in the frustrated spin degrees of freedom emerges, in agreement with the Heisenberg model studies. We will also show that in the intermediate correlation regime the nonmagnetic insulating phase expands to wider parameter regions located around $0\leq t_{2}/t_{1}\lesssim 0.3$ and $0.4\lesssim t_{2}/t_{1}\leq 1$, which go into a paramagnetic metallic phase in the weak correlation regime via the second-order Mott transition. The characteristic behavior of the critical phase boundary of the Mott transition is discussed in terms of the van Hove singularity appearing in the calculated DOS and single-particle spectral function.

The rest of the paper is organized as follows. In Sec. II, we introduce the model and discuss the method of calculation briefly. In Sec. III A, we present our results obtained in the strong correlation regime and compare them with those of the Heisenberg model. In Sec. III B, we present our results obtained in the intermediate to weak correlation regime and discuss the phase diagram of our model. The critical phase boundary of the Mott transition is also discussed. A summary of the paper is given in Sec. IV.

II Model and Method

We consider the triangular-lattice Hubbard model [see Fig. 1(a)] defined by the Hamiltonian

[TABLE]

where $c^{\dagger}_{i\sigma}$ ($c_{i\sigma}$) creates (annihilates) an electron with spin $\sigma$ at site $i$, and $n_{i\sigma}=c^{\dagger}_{i\sigma}c_{i\sigma}$. $\langle i,j\rangle$ indicates the nearest-neighbor bonds with the hopping parameter $t_{1}$ and $\langle\langle i,j\rangle\rangle$ indicates the next-nearest-neighbor bonds with the hopping parameter $t_{2}$. We consider the parameter region $0\leq t_{2}/t_{1}\leq 1$, including two limiting cases, $t_{2}=0$ (isotropic triangular lattice) and $t_{2}=t_{1}$. $U$ is the on-site Coulomb repulsion between two electrons and $\mu$ is the chemical potential maintaining the system at half filling.

In the large-$U$ limit, this model may be mapped onto the triangular-lattice Heisenberg model of spin-1/2 defined by the Hamiltonian

[TABLE]

with the exchange coupling constants of $J_{1}=4t^{2}_{1}/U$ and $J_{2}=4t^{2}_{2}/U$ for the nearest-neighbor and next-nearest-neighbor bonds, respectively. The spin operator is given by $\bm{S}_{i}=\sum_{\alpha\beta}c^{\dagger}_{i\alpha}\bm{\sigma}_{\alpha\beta}c_{i\beta}$/2 with the vector of Pauli matrices $\bm{\sigma}_{\alpha\beta}$. The results obtained for the Hubbard model [Eq. (1)] in the strong correlation regime are compared with those of the Heisenberg model [Eq. (2)].

Let us describe the VCA briefly, which is a many-body variational method based on the SFT, where the grand potential of the system is formulated as a functional of the self-energy potthoff1 ; potthoff2 ; dahnken . The ground state of the original system in the thermodynamic limit can thus be obtained via the calculation of the grand potential $\Omega$ of the system with the exact self-energy. Then, in the VCA, restricting the trial self-energy to the self-energy of the reference system $\Sigma^{\prime}$, we obtain the approximate grand potential as

[TABLE]

where $\Omega^{\prime}$ is the grand potential of the reference system, and $G_{0}$ and $G^{\prime}_{0}$ are the noninteracting Green’s functions of the original and reference systems, respectively. The Hamiltonian of the reference system $H^{\prime}$ is defined below. Note that the short-range correlations within the clusters of the reference system are taken into account exactly. See Refs. potthoff3, ; senechal, for recent reviews of the method.

The advantage of the VCA is that the spontaneous symmetry breaking can be treated within the framework of the theory, where we introduce the Weiss fields as variational parameters. In the present case, the Hamiltonian of the reference system is taken as $H^{\prime}=H+H_{\rm{M}}$ with the Weiss fields

[TABLE]

where $h^{\prime}_{\mathrm{120^{\circ}}}$ and $h^{\prime}_{\mathrm{str}}$ are the strengths of the Weiss fields for the 120∘ Néel and stripe ordered states, respectively. For the Néel order, the unit vectors $\bm{e}_{a_{i}}$ are rotated by 120∘ to each other, where $a_{i}$ ($=1,2,3$) is the sublattice index of site $i$. For the stripe order, the wave vectors can be taken equivalently as either $\bm{Q}_{\rm{str}}=(\pi,\pi/\sqrt{3})$, $(\pi,-\pi/\sqrt{3})$, or $(0,-2\pi/\sqrt{3})$. The variational parameters are optimized on the basis of the variational principle, i.e., $\partial\Omega/\partial h^{\prime}=0$, for each magnetic order, where the solution with $h^{\prime}\neq 0$ corresponds to the ordered state.

In our VCA calculations, we use the 12-site cluster shown in Fig. 2 as the reference system. This is the best appropriate and feasible choice of the reference cluster because we can treat the two-sublattice order (stripe order) with an equal number of up- and down-spin electrons and the three-sublattice order (120∘ Néel order) with an equal number of the three sublattice sites $a_{i}=1$, 2 and 3. The cluster-size and cluster-shape dependences of our results are discussed in Appendix. Note that longer period phases such as the spiral phase mentioned in a different system tocchio2 cannot be treated in the present approach; in our analysis, we fix the pitch angle of the spiral order to be 120∘ (or the three-sublattice of $a_{i}=1,2,3$) even for $t_{2}\neq 0$. The charge orderings discussed in the extended Hubbard model with intersite Coulomb repulsions tocchio3 are also neglected. To our knowledge, no other orders have been predicted in the present triangular-lattice Hubbard and Heisenberg models.

III Results of calculation

III.1 Strong correlation regime

First, let us discuss the strong correlation regime, $U/t_{1}=60$. We calculate the ground-state energies $E=\Omega+\mu$ (per site) and magnetic order parameters $M$ defined as $M_{\rm{120^{\circ}}}=(2/L)\sum_{i}\bm{e}_{a_{i}}\cdot\langle\bm{S}_{i}\rangle$ for the 120∘ Néel order and $M_{\rm{str}}=(2/L)\sum_{i}e^{i\bm{Q}_{\mathrm{str}}\cdot\bm{r}_{i}}\langle S^{z}_{i}\rangle$ for the stripe order, where $\langle\cdots\rangle$ stands for the ground-state expectation value. The results are shown in Fig. 3, where we find three phases: the 120∘ Néel ordered phase around $t_{2}/t_{1}=0$, the stripe ordered phase around $t_{2}/t_{1}=1$, and the nonmagnetic disordered phase in-between.

At $t_{2}/t_{1}=0$, the 120∘ Néel ordered state has the lowest energy and with increasing $t_{2}/t_{1}$ it approaches the energy of the nonmagnetic disordered state gradually. Then, at $t_{2}/t_{1}=0.20$, the 120∘ Néel ordered state disappears continuously. The calculated order parameter $M_{\rm{120^{\circ}}}$ also indicates the continuous (or second-order) phase transition. On the other hand, at $t_{2}/t_{1}=1.0$, the stripe ordered state has the lowest energy and, with decreasing $t_{2}/t_{1}$, the energy of stripe order crosses to that of the nonmagnetic state at $t_{2}/t_{1}=0.50$, indicating the discontinuous (or first-order) transition between the stripe and disordered phases. The calculated order parameter $M_{\rm{str}}$ also disappears discontinuously at $t_{2}/t_{1}=0.50$.

These results may be compared with the previous studies on the $J_{1}$-$J_{2}$ triangular-lattice Heisenberg model manuel ; mishmash ; kaneko ; li2 ; zhu ; hu2 . The transition point between the 120∘ Néel and nonmagnetic phases has been estimated to be $J_{2}/J_{1}=0.05-0.12$, which corresponds to $t_{2}/t_{1}=0.22-0.35$ of our Hubbard model parameters. A reasonable agreement is thus obtained. The transition point between the stripe and nonmagnetic phases has also been estimated to be $J_{2}/J_{1}=0.14-0.19$, which corresponds to $t_{2}/t_{1}=0.37-0.44$ of our Hubbard model parameters. We again find a reasonable agreement with our estimation. The orders of the phase transitions, i.e., the second-order for the 120∘ Néel phase and the first-order for the stripe phase, are also in agreement with the previous study of the Heisenberg model kaneko . We may point out that the strong quantum fluctuations in the frustrated spin degrees of freedom causes this nonmagnetic phase because the classical spin model predicts either the 120∘ Néel or four-sublattice ordered phase without any intermediate nonmagnetic phases jolicoeur ; chubukov ; deutscher ; lecheminant .

III.2 Intermediate to weak correlation regime

Next, let us discuss the intermediate to weak correlation regime $0\leq U/t_{1}\leq 10$. We here calculate the total energies, order parameters, and charge gaps of the model, as well as the grand potential as a function of the Weiss fields, and summarize them as the ground-state phase diagram in the parameter space $(t_{2}/t_{1},U/t_{1})$, as shown in Fig. 4. We find four phases: the 120∘ Néel and stripe ordered phases at large $U/t_{1}$, which are continuous to the phases at $U/t_{1}=60$ discussed above, and nonmagnetic insulating phase in-between, as well as the paramagnetic metallic phase in the weak correlation regime. In the intermediate correlation regime, the nonmagnetic insulating phase expands to wider parameter regions, which are around $0\leq t_{2}/t_{1}\lesssim 0.3$ and around $0.4\lesssim t_{2}/t_{1}\leq 1$. We note that the presence of the nonmagnetic insulating phase around $0\leq t_{2}/t_{1}\lesssim 0.3$ is in agreement with previous studies of the triangular-lattice Hubbard model at $t_{2}/t_{1}=0$ sahebsara2 ; yoshioka1 ; yoshioka2 ; yamada1 ; laubach .

The calculated order parameters of the 120∘ Néel and stripe phases are shown in Fig. 5 as a function of $U/t_{1}$ for several values of $t_{2}/t_{1}$. We find that the transition to the stripe ordered phase is continuous, irrespective of $t_{2}/t_{1}$, up to a large value of $U/t_{1}\sim 30$, but it changes to the discontinuous transition as seen in Fig. 3 at $U/t_{1}=60$. We also find that the transition to the 120∘ Néel ordered phase is discontinuous at $0<t_{2}/t_{1}\leq 0.35$ for $U/t_{1}\simeq 6$ but it is continuous for larger values of $U/t_{1}$. The transition at $U/t_{1}=60$ is also continuous (see Fig. 3). These behaviors are observed also in the calculated Weiss-field dependence of the grand potentials of our model.

The charge gap is evaluated from the total number of electrons as a function of $\mu$ (see Fig. 5) to examine the Mott metal-insulator transition of the system. We find that the transition is continuous (or second-order) and the phase boundary is located around $U/t_{1}\simeq 4-6$, as shown in Fig. 4. We note that the phase boundary decreases (shifts to a lower $U/t_{1}$ side) with increasing $t_{2}/t_{1}$ up to $t_{2}/t_{1}\simeq 0.5$, but it increases for larger values of $t_{2}/t_{1}$. This behavior is in contrast to that of the square-lattice Hubbard model with the next-nearest-neighbor hopping parameters, where a monotonous increase in the critical interaction strength is observed mizusaki ; nevidomskyy ; yamada2 , which is due to the monotonous increase in the band width of the model.

To find out the origin of this behavior, we calculate the DOS $\rho(\omega)$ and single-particle spectral function $A(\bm{k},\omega)$ in the paramagnetic state of the system using the CPT, which are defined as

[TABLE]

with the CPT Green’s function senechal

[TABLE]

where we define the $L_{c}\times L_{c}$ matrices for the cluster of size $L_{c}$ as $\mathcal{G}(\bm{k},\omega)=[G^{\prime-1}(\omega)-V(\bm{k})]^{-1}$ with $V(\bm{k})=G_{0}^{\prime-1}-G_{0}^{-1}$. The exact Green’s function of the reference system $G^{\prime}(\omega)$ is given by

[TABLE]

where $|\psi_{0}\rangle$ and $E_{0}$ are the ground state and ground state energy of $H^{\prime}$.

The calculated results for the DOS and single-particle spectral function of our model are shown in Figs. 6 and 7, respectively. We find that the sharp peak appeared above the Fermi level at $t_{2}/t_{1}=0$, which is caused by the van Hove singularity in the triangular lattice, shifts to the lower energy side with increasing $t_{2}/t_{1}$, and at $t_{2}/t_{1}=0.5$, the peak position coincides with the Fermi level (see Fig. 6). This situation of the high DOS at the Fermi level is energetically unstable fazekas , so that the band gap opens to gain in the band energy in the presence of the Hubbard interaction $U$. With further increasing $t_{2}/t_{1}$, the peak shifts to the higher energy side again. The Hubbard band gap is then the largest at $t_{2}/t_{1}=0.5$ as seen in Figs. 5 and 6. This singularity is also seen in the single-particle spectral function as the presence of the flat-band region around the K point of the Brillouin zone (see Fig. 7). This behavior thus explains why the critical interaction strength becomes small at around $t_{2}/t_{1}\simeq 0.5$.

To confirm the absence of any magnetic instability in our model in the weak correlation regime, we here calculate the generalized susceptibility (or Lindhard function) in the noninteracting limit, which is defined as

[TABLE]

where $\varepsilon_{\bm{k}}$ is the corresponding noninteracting band dispersion and $f(\varepsilon)$ is the Fermi function. The calculated results at temperature $0.01t_{1}$ are shown in Fig. 8, where we find that, in accordance with the absence of significant Fermi-surface nesting features, no singular behaviors actually appear in $\chi_{0}(\bm{q})$, indicating the absence of magnetic long-range orders in the weak correlation limit. This result supports the validity of our phase diagram shown in Fig. 4 in the weak correlation regime.

IV Summary

We have studied the Mott metal-insulator transition and magnetism of the triangular-lattice Hubbard model at half filling in the entire region of the interaction strength, taking into account the next-nearest-neighbor hopping parameters for the effects of magnetic frustrations. We have employed the method of VCA based on the SFT, which has not been used for the present purposes. We have thereby calculated the grand potential of the system as a function of the Weiss fields for the 120∘ Néel and stripe magnetic orders, and have determined the order parameters. We have also calculated the DOS and single-particle spectral function as well as the charge gap of the system. These results have been summarized as the ground-state phase diagram of the system.

We have found four phases: In the strong correlation regime, there appear (i) the 120∘ Néel ordered phase in a wide parameter region around $t_{2}/t_{1}\simeq 0$ and (ii) the stripe ordered phase in a wide parameter region around $t_{2}/t_{1}\simeq 1$, and in-between, (iii) the nonmagnetic insulating phase caused by the quantum fluctuations in the geometrically frustrated spin degrees of freedom emerges. The obtained phase boundaries in the strong correlation limit have been compared with those of the corresponding Heisenberg model to find a reasonable agreement. The orders of the phase transitions of the two magnetically ordered phases have also been determined. In the intermediate correlation regime, the nonmagnetic insulating phase expands to a wider parameter region of $t_{2}/t_{1}$. Then, decreasing the interaction strength further, the system turns into (iv) the paramagnetic metallic phase in the weak correlation regime via the second-order Mott metal-insulator transition. The characteristic behavior of the critical phase boundary of the Mott transition has also been discussed in terms of the shift in the van Hove singularity due to the presence of $t_{2}$, as seen in the calculated DOS and single-particle spectral function.

We suggest that the phase diagram obtained here may contain different types of nonmagnetic insulator (or spin liquid) states depending on the region in the parameter space. The characterization of the states is however beyond the scope of the VCA approach based on the self-energy (or single-particle Green’s function), for which we hope that our results will encourage future studies.

Acknowledgements.

We thank S. Miyakoshi for enlightening discussions and K. Seki for careful reading of our manuscript. This work was supported in part by a Grant-in-Aid for Scientific Research (No. 26400349) from JSPS of Japan. T. K. acknowledges support from the JSPS Research Fellowship for Young Scientists.

Appendix A Cluster dependence of the phase boundaries

In the VCA, or in any other quantum cluster methods, we can in principle calculate the physical quantities in the thermodynamic limit, but the calculated results necessarily depend on the size and shape of the solver cluster. Thus, the choice of the solver cluster is important in the present approach. In the main text, we have chosen the 12-site cluster shown in Fig. 2, which is the best appropriate one because it is first of all computationally feasible and also because it fits with both the three-sublattice 120∘ order and two-sublattice stripe order without introducing unnecessary frustrations. However, it seems instructive to check the cluster-size and cluster-shape dependences of our results presented in the main text.

Here, we choose several clusters [see Fig. 9(b)] that fit either with the 120∘ order or with the stripe order, and we calculate the phase boundaries to check the solver cluster dependence of the ground-state phase diagram shown in Fig. 4. The calculated results for the phase boundaries are shown in Fig. 9(a). We thus find that the phase boundary between the 120∘ Néel ordered and nonmagnetic insulating phases is located around $0\leq t_{2}/t_{1}\lesssim 0.4$ in an intermediate to large $U/t_{1}$ $(\gtrsim 6)$ region and that the phase boundary between the stripe ordered and nonmagnetic insulating phases is located around $0.5\lesssim t_{2}/t_{1}\leq 1$ in an intermediate to large $U/t_{1}$ $(\gtrsim 7)$ region, irrespective of the appropriate choices of the solver cluster. We also find that the phase boundary of the Mott metal-insulator transition is located around $U/t_{1}\simeq 4-6$ with a minimum at $t_{2}/t_{1}\simeq 0.5$, irrespective of the choices of the solver cluster.