Charge orders in organic charge-transfer salts

In the following, we will consider a model in which every site (i.e. dimer) accommodates two orbitals (hereinafter referred to as c and f), one for each ET molecule. The original lattice is triangular, with hopping and interaction terms depicted in figure 1(a). An equivalent description is given by considering a two-orbital model on the square lattice, see figure 1(b). Here, we can define a partition in two sub-lattices ${ \mathcal A }$ and ${ \mathcal B }$, where the ET molecules form horizontal and vertical dimers, respectively. The full Hamiltonian, in this latter description, is given by:

$$\begin{eqnarray}&&{ \mathcal H }={{ \mathcal H }}_{t}+{{ \mathcal H }}_{V}+{{ \mathcal H }}_{U},\end{eqnarray} \tag{ 1 }$$

where

$$\begin{eqnarray}\,{{ \mathcal H }}_{t} & = & {t}_{b1}\displaystyle \sum _{i,\sigma }{c}_{i,\sigma }^{\dagger }{f}_{i,\sigma }+{t}_{b2}\displaystyle \sum _{i,\sigma }{c}_{i,\sigma }^{\dagger }{f}_{i+x+y,\sigma }\\ & & +{t}_{q}\displaystyle \sum _{i,\sigma }({c}_{i,\sigma }^{\dagger }{f}_{i+x,\sigma }+{c}_{i,\sigma }^{\dagger }{f}_{i+y,\sigma })\\ & & +{t}_{p}\displaystyle \sum _{i\in { \mathcal A },\sigma }({c}_{i,\sigma }^{\dagger }{c}_{i+x,\sigma }+{c}_{i,\sigma }^{\dagger }{c}_{i-y,\sigma }+{f}_{i,\sigma }^{\dagger }{f}_{i-x,\sigma }+{f}_{i,\sigma }^{\dagger }{f}_{i+y,\sigma })+{\rm{h}}.{\rm{c}}.,\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}{{ \mathcal H }}_{V} & = & {V}_{b1}\displaystyle \sum _{i}{n}_{i}^{c}{n}_{i}^{f}+{V}_{b2}\displaystyle \sum _{i}{n}_{i}^{c}{n}_{i+x+y}^{f}\\ & & +{V}_{q}\displaystyle \sum _{i}({n}_{i}^{c}{n}_{i+x}^{f}+{n}_{i}^{c}{n}_{i+y}^{f})\\ & & +{V}_{p}\displaystyle \sum _{i\in { \mathcal A }}({n}_{i}^{c}{n}_{i+x}^{c}+{n}_{i}^{c}{n}_{i-y}^{c}+{n}_{i}^{f}{n}_{i-x}^{f}+{n}_{i}^{f}{n}_{i+y}^{f}),\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&{{ \mathcal H }}_{U}=U\displaystyle \sum _{i}({n}_{i,\uparrow }^{c}{n}_{i,\downarrow }^{c}+{n}_{i,\uparrow }^{f}{n}_{i,\downarrow }^{f}).\,\end{eqnarray} \tag{ 4 }$$

Here, ${c}_{i,\sigma }^{\dagger }$ (${f}_{i,\sigma }^{\dagger }$) creates an electron with spin σ on the c (f) molecular orbital at site i; ${n}_{i}^{c}={\sum }_{\sigma }{n}_{i,\sigma }^{c}$ (${n}_{i}^{f}={\sum }_{\sigma }{n}_{i,\sigma }^{f}$) are the density operator for c (f) electrons at site i. Hopping and interaction terms can be divided into those that connect c and f orbitals and those that connect orbitals of the same kind, see figure 1(b). Belonging to the former class, there are terms connecting orbitals within the same dimer (b₁-type), along x and y nearest-neighbor sites (q-type), and along the x = y diagonal (b₂-type); instead, p-type terms connect orbitals at nearest-neighbor sites and belong to the latter class. Accordingly, the noninteracting Hamiltonian ${{ \mathcal H }}_{t}$ contains four hopping parameters, i.e. t_b1, t_b2, t_q, and t_p. Similarly, the interacting Hamiltonian ${{ \mathcal H }}_{V}$ contains four intermolecular Coulomb interactions, i.e. V_b1, V_b2, V_q, and V_p. Note that the translational symmetry and the consequent partition between ${ \mathcal A }$ and ${ \mathcal B }$ sub-lattices is only due to the presence of p-type terms. Finally, ${{ \mathcal H }}_{U}$ describes the Hubbard-U interaction on each molecule. Our calculations are performed on finite clusters of size ${N}_{s}={L}^{2}$ (where on each site there are two molecules, i.e. orbitals), with periodic–antiperiodic boundary conditions on both directions. The filling factor is fixed to be 3/4.

Zoom In Zoom Out Reset image size

As discussed in [35], this two-orbital model reduces to the single-band Hubbard model (at half filling), when the intradimer hopping is very large (i.e. ${t}_{b1}\gg {t}_{b2},{t}_{p},{t}_{q}$). Furthermore, at ${t}_{b2}=0$ and ${t}_{p}=0$ (or ${t}_{q}=0$) the Hamiltonian reduces to the recently investigated Hubbard model on the honeycomb lattice with anisotropic terms [40].

In this work, we consider the following hopping parameters (in units of t_b1):

$$\begin{eqnarray}&&{t}_{b1}=1,\,{t}_{b2}=0.359,\,{t}_{p}=0.539,\,{t}_{q}=0.221,\end{eqnarray} \tag{ 5 }$$

which are based on the results obtained by density functional theory calculations on κ-(ET)₂Cu[N(CN)₂]Cl [35, 41]. The noninteracting band structure is reported in figure 2. As far as the interaction terms are concerned, for realistic systems, one expects U to be the largest Coulomb repulsion term and V_b1 to be the second largest one, while V_b2 should be comparable to V_p and V_q.

Zoom In Zoom Out Reset image size

We first discuss the possible ground states in the atomic limit, i.e. for ${t}_{b1}={t}_{b2}={t}_{p}={t}_{q}=0$ at 3/4 filling. If the only finite interaction is the intramolecular Hubbard-U term, the ground state is highly degenerate, with all possible charge patterns having N_s doubly-occupied molecules. A further degeneracy arises from the remaining N_s molecules being singly-occupied, where any spin configuration gives the same energy. The charge degeneracy can be lifted by including the intermolecular interactions V_b1, V_b2, V_p, and V_q. We concentrate here on three particular relevant cases that show regular patterns of charge order (see figure 3): two of them are polar charge-order insulators (hereinafter denoted as PCOI and PCOI${}^{{\prime} }$), since there is a charge disproportionation inside each dimer, and one is a nonpolar charge-order insulator (denoted as NPCOI), since the two molecules of the same dimer have the same amount of charge. Their energies per site (i.e. per dimer) can be easily evaluated in the atomic limit:

$$\begin{eqnarray}&&{E}_{\mathrm{polar}}=E+{V}_{q},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{E}_{\mathrm{polar}^{\prime} }=E+{V}_{p},\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{E}_{\mathrm{nonpolar}}=E+\displaystyle \frac{1}{2}({V}_{b1}+{V}_{b2}),\end{eqnarray} \tag{ 8 }$$

where we defined:

$$\begin{eqnarray}&&E=U+2{V}_{b1}+4{V}_{p}+4{V}_{q}+2{V}_{b2}.\end{eqnarray} \tag{ 9 }$$

The phase diagram in the ${V}_{p}\mbox{--}{V}_{q}$ plane is shown in figure 3. Here, V_b1 and V_b2 only modify the phase boundaries between the polar and nonpolar charge-ordered phases. The NPCOI appears when both V_p and V_q dominate over V_b1 and V_b2. Otherwise, the two polar states are stable and the competition between V_p and V_q determines the detailed charge pattern. All three phases are degenerate for ${V}_{p}={V}_{q}=({V}_{b1}+{V}_{b2})/2$.

Zoom In Zoom Out Reset image size

Our numerical results are obtained by means of the variational Monte Carlo method, that is based on the definition of suitable wave functions that approximate the ground-state properties beyond perturbative approaches. We consider Jastrow–Slater wave functions [42–45], which are described as:

$$\begin{eqnarray}&&| {\rm{\Psi }}\rangle ={ \mathcal J }| {\rm{\Phi }}\rangle .\end{eqnarray} \tag{ 10 }$$

Here, ${ \mathcal J }$ is a long-range density–density Jastrow factor given by:

$$\begin{eqnarray}&&{ \mathcal J }=\exp \left(-\displaystyle \frac{1}{2}\displaystyle \sum _{i,j,\alpha ,\beta }{v}_{{ij}}^{\alpha \beta }{n}_{i}^{\alpha }{n}_{j}^{\beta }\right),\end{eqnarray} \tag{ 11 }$$

where ${n}_{i}^{\alpha }$ is the electronic density at site i and orbital $\alpha =c,f$, while ${v}_{{ij}}^{\alpha \beta }$ are translational invariant variational parameters, that are optimized by only imposing translational and inversion symmetry in the square lattice defined by ${{\boldsymbol{R}}}_{i}=({x}_{i},{y}_{i})$. $| {\rm{\Phi }}\rangle$ is a noninteracting fermionic state that is defined as the ground state of an auxiliary Hamiltonian with site-dependent chemical potentials and magnetic order parameters. Such a choice of an auxiliary Hamiltonian allows us to describe both charge and spin orders induced by the intermolecular Coulomb interactions [40, 46, 47]. In particular, for the insulating states, we consider:

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{ins}}={{ \mathcal H }}_{t}+{{ \mathcal H }}_{\mathrm{COI}}+{{ \mathcal H }}_{\mathrm{AF}},\end{eqnarray} \tag{ 12 }$$

where ${{ \mathcal H }}_{t}$ is the kinetic part of equation (2) and

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{COI}}=\displaystyle \sum _{i}{{\rm{e}}}^{{\rm{i}}{\boldsymbol{Q}}\cdot {{\boldsymbol{R}}}_{i}}({\mu }^{c}{n}_{i}^{c}+{\mu }^{f}{n}_{i}^{f}),\end{eqnarray} \tag{ 13 }$$

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{AF}}=\displaystyle \sum _{i}[{m}_{i}^{c}({c}_{i,\uparrow }^{\dagger }{c}_{i,\downarrow }+{c}_{i,\downarrow }^{\dagger }{c}_{i,\uparrow })+{m}_{i}^{f}({f}_{i,\uparrow }^{\dagger }{f}_{i,\downarrow }+{f}_{i,\downarrow }^{\dagger }{f}_{i,\uparrow })].\end{eqnarray} \tag{ 14 }$$

Here, ${\boldsymbol{Q}}=(\pi ,\pi )$ describes the NPCOI (with ${\mu }^{c}={\mu }^{f}$) and the PCOI (with ${\mu }^{c}=-{\mu }^{f}$) phases of figure 3, while ${\boldsymbol{Q}}=(0,0)$ (with ${\mu }^{c}=-{\mu }^{f}$) gives the PCOI${}^{{\prime} }$ phase of figure 3. We optimize the variational magnetic parameters at charge-rich and charge-poor molecular orbitals independently, according to the condition:

$$\begin{eqnarray}{m}_{i}^{\alpha }=\left\{\begin{array}{ll}{m}_{1}^{\alpha } & \,\mathrm{if}\,{{\rm{e}}}^{{\rm{i}}{\boldsymbol{Q}}\cdot {{\boldsymbol{R}}}_{i}}{\mu }^{\alpha }\lt 0,\\ {m}_{2}^{\alpha } & \,\mathrm{if}\,{{\rm{e}}}^{{\rm{i}}{\boldsymbol{Q}}\cdot {{\boldsymbol{R}}}_{i}}{\mu }^{\alpha }\gt 0,\end{array}\right.\end{eqnarray} \tag{ 15 }$$

which implies that ${m}_{1}^{\alpha }$ and ${m}_{2}^{\alpha }$ are associated to the magnetization of the charge-rich and charge-poor molecular orbitals on site i, respectively. In general, incommensurate magnetic order may coexist with commensurate charge order; however, this is beyond the scope of the present paper, and we restrict ourselves to commensurate (and collinear) magnetic order. Notice that, within our variational description based upon an auxiliary Hamiltonian, it is particularly easy to consider nonmagnetic states, which can be described by taking ${m}_{i}^{c}={m}_{i}^{f}=0$ in equation (14).

In order to describe metallic states, we consider the following auxiliary Hamiltonian:

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{met}}={{ \mathcal H }}_{t}+{{ \mathcal H }}_{\mathrm{COM}},\end{eqnarray} \tag{ 16 }$$

where

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{COM}}=\displaystyle \sum _{i}[{\mu }_{R(i)}^{c}{n}_{i}^{c}+{\mu }_{R(i)}^{f}{n}_{i}^{f}];\end{eqnarray} \tag{ 17 }$$

here, the sublattice index R(i) at position ${{\boldsymbol{R}}}_{i}$ is defined as:

$$\begin{eqnarray}&&R(i)=\mathrm{mod}({x}_{i}-{y}_{i},6),\end{eqnarray} \tag{ 18 }$$

which allows a 12-sublattice charge ordering, since there are two orbitals on each site. We neglect in the calculation the possible presence of magnetic order in the metallic states.

In order to exclude the presence of further ordered phases in the explored regions of the phase diagram, we have also employed unbiased wave functions, where different charge orderings can spontaneously emerge. In particular, we constructed a noninteracting wave function $| {\rm{\Phi }}\rangle$ as the ground state of the tight-binding Hamiltonian with site-dependent chemical potentials ${\mu }_{i}^{c}$ and ${\mu }_{i}^{f}$:

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{full}}={{ \mathcal H }}_{t}+\displaystyle \sum _{i}({\mu }_{i}^{c}{n}_{i}^{c}+{\mu }_{i}^{f}{n}_{i}^{f}),\end{eqnarray} \tag{ 19 }$$

where the chemical potentials are variational parameters independently optimized for each site i (several initial configurations of ${\mu }_{i}^{c}$ and ${\mu }_{i}^{f}$ have been chosen, in order to assess the possibility to remain stuck in local minima). Since the number of parameters to be optimized grows as 2N_s, we considered this approach only for N_s = 36. In this case, we have observed that the selected charge orderings are consistent with the states described by the simpler approach above. In addition, we notice that charge order can be also generated by a translationally invariant Jastrow factor, without explicitly breaking the symmetry in the Slater determinant, as shown in [46, 47]. The advantage is that we do not need to assume a priori any type of charge ordering and, if long-range order exists, charge-ordered states should be selected by the optimization of the Jastrow factor. In general, for the chosen set of hoppings and interaction terms, we never found charge orders that cannot be captured by the previous parametrization of equations (13) and (16).

Given the presence of the correlation term (i.e. the Jastrow factor), an analytical evaluation of the variational energy or of any correlation function is impossible on large sizes; nevertheless, a standard Monte Carlo sampling can be employed to obtain all the physical quantities with high accuracy.

We now investigate how the hopping terms in equation (1) modify the phase diagram obtained for the atomic limit, in the region $U\gg {t}_{b1}$. The schematic phase diagram is shown in figure 4. The three charge-ordered phases obtained in the atomic limit are stable also in the presence of finite values of the hopping terms; however, magnetic order is generated from virtual hopping processes involving charge-poor molecules, which form one-dimensional patterns in the lattice and are effectively half filled in all these phases. Antiferromagnetic correlations are then expected along these one-dimensional chains, which are formed by the bonds with hopping t_b1 and t_b2 for the nonpolar state and by the bonds with hopping t_p (t_q) for the PCOI${}^{{\prime} }$ (PCOI). Therefore, in the nonpolar charge-ordered state, the two spins on the molecules of the same dimer have opposite orientations, thus implying that the dimer has no net magnetization. By contrast, the two polar states show ferromagnetic spin correlations within the dimers; here, each dimer contains one charge-rich and one charge-poor molecule, the magnetization being large in the latter one. Moreover, we observe long-range antiferromagnetic order of the magnetic moments of dimers.

Zoom In Zoom Out Reset image size

In addition to these three states, a uniform DMI intrudes between the polar phases. This correlated phase should appear when U is much larger than all the V terms, in the region where V_p and V_q are competing [15–17]. Here, spin correlations are ferromagnetic within each dimer, since there is an average of three electrons per site: two electrons have opposite spins and do not contribute to the magnetic moment, which is fully due to the third electron that is delocalized between the two molecules. Similarly to the two polar charge-ordered states, also here the spins of the dimers possess long-range antiferromagnetic order. We find that the transitions between the DMI and the polar charge-ordered phases (PCOI and PCOI${}^{{\prime} }$) are continuous (see below). Close to the boundaries between nonpolar and polar charge-ordered phases, the DMI state can also be stabilized; however its energy remains higher than the energies of the other phases, indicating that it is a metastable state.

In order to understand the nature of the charge properties of all the insulating phases, we calculate the total charge structure factor N(q) and the charge-disproportionation structure factor ${N}_{\mathrm{CD}}(q)$, defined as:

$$\begin{eqnarray}&&N(q)=\displaystyle \frac{1}{{N}_{s}}\displaystyle \sum _{i,j}\langle ({n}_{i}^{c}+{n}_{i}^{f})({n}_{j}^{c}+{n}_{j}^{f})\rangle {{\rm{e}}}^{{\rm{i}}{\boldsymbol{q}}\cdot ({{\boldsymbol{R}}}_{i}-{{\boldsymbol{R}}}_{j})},\end{eqnarray} \tag{ 20 }$$

$$\begin{eqnarray}&&{N}_{\mathrm{CD}}(q)=\displaystyle \frac{1}{{N}_{s}}\displaystyle \sum _{i,j}\langle ({n}_{i}^{c}-{n}_{i}^{f})({n}_{j}^{c}-{n}_{j}^{f})\rangle {{\rm{e}}}^{{\rm{i}}{\boldsymbol{q}}\cdot ({{\boldsymbol{R}}}_{i}-{{\boldsymbol{R}}}_{j})},\end{eqnarray} \tag{ 21 }$$

where $\langle ...\rangle$ indicates the expectation value over the variational wave function of equation (10). Here, $N(q=0)$ is set to zero. The metallic or insulating character can be assessed by inspecting the small-q limit of the total charge structure factor. Indeed, in the limit $| {\boldsymbol{q}}| \to 0$, $N(q)\propto | {\boldsymbol{q}}|$ for a metal, while $N(q)\propto | {\boldsymbol{q}}{| }^{2}$ for an insulator [48, 49]. In addition, charge order is indicated by the presence of a Bragg peak in N(q) or ${N}_{\mathrm{CD}}(q)$. In the former case, charge order is characterized by charge-rich dimers on one sublattice and charge-poor dimers on the other one, while, in the latter case, charge disproportionation occurs within the dimers.

In the following, we fix the Coulomb interactions to $U/{t}_{b1}=10$, ${V}_{b1}/{t}_{b1}=4$, and ${V}_{b2}/{t}_{b1}=2$, and vary V_p and V_q. Within this choice, in the atomic limit, the polar and nonpolar phases are degenerate for ${V}_{p}={V}_{q}=3{t}_{b1}$. As shown in figure 5, all the phases presented in the phase diagrams are insulating, since $N(q)\propto | {\boldsymbol{q}}{| }^{2}$ in the limit $| {\boldsymbol{q}}| \to 0$, both along the q_y = 0 and the q_x = q_y lines in reciprocal space. Then, each insulating phase can be fully characterized by N(q) and ${N}_{\mathrm{CD}}(q)$, see figures 6 and 7. The DMI does not show any Bragg peak either in N(q) or in ${N}_{\mathrm{CD}}(q)$ (figures 6(d) and 7(d)), suggesting that no long-range charge order occurs. The nonpolar charge-ordered state shows the Bragg peak at ${\boldsymbol{Q}}=(\pi ,\pi )$ in N(q) (figure 6(c)), but no sharp peaks in ${N}_{\mathrm{CD}}(q)$ (figure 7(c)). This implies that staggered charge disproportionation appears between different dimers, but not within the dimers. Finally, the polar charge-ordered states show the Bragg peak at ${\boldsymbol{Q}}=(0,0)$ (PCOI${}^{{\prime} }$) or ${\boldsymbol{Q}}=(\pi ,\pi )$ (PCOI) in ${N}_{\mathrm{CD}}(q)$ (figures 7(a) and b)), but no sharp peaks in N(q) (figures 6(a) and (b)). This fact indicates charge disproportionation within the dimers, while the number of electrons in each dimer is constant. Each orbital has the same number of electrons at each site for ${\boldsymbol{Q}}=(0,0)$, while each orbital alternates between charge-rich and charge-poor configurations for ${\boldsymbol{Q}}=(\pi ,\pi )$, see figure 4.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Remarkably, all polar and nonpolar phases can be stabilized within the variational approach also without considering magnetic order in the Slater determinant. By contrast, the DMI cannot be stabilized without including the ${{ \mathcal H }}_{\mathrm{AF}}$ of equation (14), see figure 8. The charge patterns are similar to the ones that have been obtained previously with the inclusion of magnetic order (not shown).

Zoom In Zoom Out Reset image size

We focus now on the interplay between charge and spin degrees of freedom near the boundary of the polar charge-ordered phases. In particular, we show the numerical results along the ${V}_{p}+{V}_{q}=3{t}_{b1}$ line (still fixing $U/{t}_{b1}=10$, ${V}_{b1}/{t}_{b1}=4$, and ${V}_{b2}/{t}_{b1}=2$), which crosses the two polar and the DMI phases. The absolute value of the difference between the optimized chemical potentials $| {\mu }^{c}-{\mu }^{f}|$ for orbitals c and f (see equation (13)) can be used as a diagnostic to detect polar and nonpolar states. As shown in figure 9, we find that $| {\mu }^{c}-{\mu }^{f}|$ is finite for ${V}_{p}/{t}_{b1}\lesssim 1.3$ and ${V}_{p}/{t}_{b1}\gtrsim 1.7$, while it vanishes in a narrow but finite region for $1.4\lesssim {V}_{p}/{t}_{b1}\lesssim 1.6$, indicating the existence of the DMI. In addition, $| {\mu }^{c}-{\mu }^{f}|$ does not show any evidence of discontinuities at the transition points, strongly suggesting that the three phases are continuously connected. Indeed, near the phase boundary obtained in the atomic limit (${V}_{p}/{t}_{b1}=1.5$) it is not possible to stabilize metastable wave functions with higher energies.

Zoom In Zoom Out Reset image size

To further investigate the connection among these three phases, we calculate the charge-disproportionation structure factor as function of V_p, as shown in figure 10. When V_q (V_p) is sufficiently large, ${N}_{\mathrm{CD}}(q)$ shows a sharp peak at ${\boldsymbol{Q}}=(0,0)$ (${\boldsymbol{Q}}=(\pi ,\pi )$) (figures 10(a) and (b), respectively). By contrast, when ${V}_{p}\approx {V}_{q}$, ${N}_{\mathrm{CD}}(q)$ shows a broad crest along the ${q}_{x}+{q}_{y}=2\pi$ direction (figures 10(c)–(e)). Importantly, there are no divergences in the thermodynamic limit, since the crest remains finite when increasing the size of the cluster. The peculiar one-dimensional-like shape of ${N}_{\mathrm{CD}}(q)$ might be understood in the following simple way: for ${V}_{p}\approx {V}_{q}$, the emergence of charge order is controlled only by V_b1 and V_b2, which define diagonal chains in the lattice, see figure 1. It is then natural to expect that correlations do not show any dependence on the transverse direction.

Zoom In Zoom Out Reset image size

The absence of charge disproportionation for $1.4\lesssim {V}_{p}/{t}_{b1}\lesssim 1.6$ is clearly demonstrated by performing the size scaling of ${N}_{\mathrm{CD}}(Q)/{L}^{2}$ for $L\to \infty$. The results are reported in figure 11. For ${V}_{p}/{t}_{b1}\lesssim 1.3$ and ${V}_{p}/{t}_{b1}\gtrsim 1.7$, we have that ${N}_{\mathrm{CD}}(Q)/{L}^{2}$ is finite for ${\boldsymbol{Q}}=(0,0)$ and ${\boldsymbol{Q}}=(\pi ,\pi )$, respectively. Instead, for $1.4\lesssim {V}_{p}/t\lesssim 1.6$ this quantity goes to zero for $L\to \infty$, suggesting that the DMI is stable in this region.

Zoom In Zoom Out Reset image size

Even if all the charge-ordered phases are present also in the absence of magnetic order, all of them are found to possess stable magnetic order when this possibility is included in the variational state, as shown in figure 12. The DMI shows the same absolute value of the magnetic moment for the orbitals c and f, as expected for a charge uniform state. When the intersite Coulomb interactions become anisotropic (i.e. $| {V}_{p}-{V}_{q}| \gt 0$), magnetic orders for the orbitals c and f start to deviate. This is due to the fact that the charge-rich (charge-poor) molecular orbitals possess a smaller (larger) magnetic moment in the polar charge-ordered phases. Notice that the PCOI state has a larger magnetic order than the PCOI${}^{{\prime} }$ one. This may be due to the anisotropy in the hopping terms. Indeed, in the PCOI state the singly-occupied molecules are connected by t_q, while in the PCOI${}^{{\prime} }$ one they are connected by t_p; since ${t}_{p}\gt {t}_{q}$ (see equation (5)) the latter case has more charge fluctuations (i.e. it is closer to a metal–insulator transition), thus implying smaller magnetic moments.

Zoom In Zoom Out Reset image size

Finally, we focus our attention on the charge-ordered metal (COM) that appears for small values of $U/{t}_{b1}$. Therefore, we now fix ${V}_{b1}/{t}_{b1}=4$, ${V}_{p}/{t}_{b1}=3.5$, and ${V}_{q}/{t}_{b1}=3$ and vary $U/{t}_{b1}$ and ${V}_{b2}/{t}_{b1}$. The phase diagram is shown in figure 13(a). Here, a large metallic phase, with honeycomb-like charge ordering, appears, for relatively small values of the intramolecular interaction. This charge-ordered pattern is similar to the three-sublattice one [46, 50, 51], which has been stabilized on the triangular lattice for intermediate values of the nearest-neighbor interaction. In our case, the periodicity is extended to 12 sites due to the anisotropy of the parameters.

Zoom In Zoom Out Reset image size

The emergence of the honeycomb-like COM can be easily understood when all the bonds (of b₁-, b₂-, p-, and q-type) are equivalent. In this 'isotropic' limit, when considering each molecule as an independent site, the underlying lattice becomes triangular, see figure 13(b). In this case, by decreasing the intramolecular Coulomb interaction U, there is an insulator to metal transition, with a metallic phase below the critical point. Moreover, since intermolecular Coulomb interactions screen the actual value of U, the metallic phase is even more stable when the V's are present in the model. However, the presence of intermolecular Coulomb interactions leads to a spontaneous symmetry breaking in the translational symmetry and to charge disproportionation, that on the triangular lattice it is natural to assume with a three-sublattice ordering $A\mbox{--}B\mbox{--}C$. For an average occupation per site (i.e. molecule) $n=3/2$, the only possible choice to minimize the energy loss due to the intermolecular interactions is then to reduce the electron occupation on one sublattice and increase it on the other two (the limiting case being n_A = 0.5 and ${n}_{B}={n}_{C}=2$).

We investigate now the stability of the COM against a normal metal when decreasing the intermolecular interactions. In this case, we vary all the Coulomb terms together, taking ${V}_{b1}/{t}_{b1}=\alpha$, ${V}_{b2}/{t}_{b1}=0.5\alpha$, ${V}_{p}/{t}_{b1}=0.875\alpha$, and ${V}_{q}/{t}_{b1}=0.75\alpha$, while the intramolecular interaction is fixed to $U/{t}_{b1}=6$. As shown in figure 14, for $\alpha \lesssim 2$ the ground state is found to be a uniform metal, with no charge disproportionation. Charge order appears for $\alpha \approx 2.5$ and is characterized by the rich–rich–poor pattern. A direct comparison between the COM at $\alpha =4$ and the uniform metallic phase is presented in figure 15. Both phases are indeed metallic, since $N(q)\propto | {\boldsymbol{q}}|$ for small momenta (figures 15(a), (b)). On the contrary, the formation of charge order in the COM phase, is signaled by the appearance of strong peaks in both the total charge structure factor N(q) and the structure factor for charge disproportionation ${N}_{\mathrm{CD}}(q)$, corresponding to the real-space configuration illustrated in figure 13(b).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size