Charge orders in organic charge-transfer salts

TL;DR

This paper investigates charge and magnetic phases in an extended two-orbital Hubbard model on a triangular lattice, revealing polar charge-ordered states that could explain ferroelectricity in organic salts, independent of magnetism.

Contribution

It introduces a detailed microscopic model showing how ferroelectricity can arise from charge order without magnetic order in organic charge-transfer salts.

Findings

Polar charge-ordered insulating phases are stabilized in the strongly correlated limit.

Ferroelectricity can occur without antiferromagnetic order.

Metallic states with unique charge patterns emerge in weak interaction regimes.

Abstract

Motivated by recent experimental suggestions of charge-order-driven ferroelectricity in organic charge-transfer salts, such as $\kappa$-(BEDT-TTF)$_2$Cu[N(CN)$_2$]Cl, we investigate magnetic and charge-ordered phases that emerge in an extended two-orbital Hubbard model on the anisotropic triangular lattice at $3/4$ filling. This model takes into account the presence of two organic BEDT-TTF molecules, which form a dimer on each site of the lattice, and includes short-range intramolecular and intermolecular interactions and hoppings. By using variational wave functions and quantum Monte Carlo techniques, we find two polar states with charge disproportionation inside the dimer, hinting to ferroelectricity. These charge-ordered insulating phases are stabilized in the strongly correlated limit and their actual charge pattern is determined by the relative strength of intradimer to interdimer couplings. Our results suggest that ferroelectricity is not driven by magnetism, since these polar phases can be stabilized also without antiferromagnetic order and provide a possible microscopic explanation of the experimental observations. In addition, a conventional dimer-Mott state (with uniform density and antiferromagnetic order) and a nonpolar charge-ordered state (with charge-rich and charge-poor dimers forming a checkerboard pattern) can be stabilized in the strong-coupling regime. Finally, when electron-electron interactions are weak, metallic states appear, with either uniform charge distribution or a peculiar $12$-site periodicity that generates honeycomb-like charge order.

Equations

$$\mathcal{H}={\mathcal{H}}_t+{\mathcal{H}}_V+{\mathcal{H}}_U,$$

$${\mathcal{H}}_t$$

$${\mathcal{H}}_V$$

$${\mathcal{H}}_U$$

$$t_{b1}=1,t_{b2}=0.359,t_p=0.539,t_q=0.221,$$

$$E_{\mathrm{polar}}$$

$$E_{\mathrm{pola}{\mathrm{r}}^{\prime }}$$

$$E_{\mathrm{nonpolar}}$$

$$E=U+2V_{b1}+4V_p+4V_q+2V_{b2}.$$

$$|\Psi ⟩=\mathcal{J}|\Phi ⟩.$$

$$\mathcal{J}=\exp \left(-\frac{1}{2}\sum _{i,j,\alpha ,\beta }{v}_{ij}^{\alpha \beta }{n}_i^{\alpha }{n}_j^{\beta }\right),$$

$${\mathcal{H}}_{\mathrm{ins}}={\mathcal{H}}_t+{\mathcal{H}}_{\mathrm{COI}}+{\mathcal{H}}_{\mathrm{AF}},$$

$${\mathcal{H}}_{\mathrm{COI}}$$

$${\mathcal{H}}_{\mathrm{AF}}$$

$$m_{i}^{\alpha}=\left\{\begin{array}[]{ll}m_{1}^{\alpha}&{\rm~{}if~{}}e^{i\bm{Q}\cdot\bm{R}_{i}}\mu^{\alpha}<0,\\ m_{2}^{\alpha}&{\rm~{}if~{}}e^{i\bm{Q}\cdot\bm{R}_{i}}\mu^{\alpha}>0,\end{array}\right.$$

$${\mathcal{H}}_{\mathrm{met}}={\mathcal{H}}_t+{\mathcal{H}}_{\mathrm{COM}},$$

$${\mathcal{H}}_{\mathrm{COM}}=\sum _i[{\mu }_{R(i)}^c{n}_i^c+{\mu }_{R(i)}^f{n}_i^f];$$

$$R(i)=\mathrm{mod}(x_i-y_i,6),$$

$${\mathcal{H}}_{\mathrm{full}}={\mathcal{H}}_t+\sum _i({\mu }_i^c{n}_i^c+{\mu }_i^f{n}_i^f),$$

$$N(q)$$

$$N_{\mathrm{CD}}(q)$$

Full text

1 Introduction

Orbital, charge, and spin degrees of freedom are intertwined in correlated electron systems and the search for unconventional quantum phases emerging from the interplay of these degrees of freedom is a very active field of research in condensed-matter physics. In particular, multiferroicity [1], where magnetism and ferroelectricity coexist, has received a lot of attention in recent years. Conventionally, one can divide multiferroics into two groups: In type-I multiferroics, ferroelectricity and magnetism have different origins [2], while in type-II multiferroics, ferroelectricity occurs only in the magnetically ordered state, where, for example, it is induced by helical magnetic order in geometrically frustrated antiferromagnets [3, 4, 5]. Recently, charge-order-driven ferroelectricity was proposed in organic charge transfer salts [6], such as $\kappa$-(ET)2Cu[N(CN)2]Cl [7] and $\alpha$-(ET)2I3 [8, 9, 10], where ET stands for BEDT-TTF [bis(ethylenedithio)-tetrathiafulvalene]. In the former one, ferroelectric and antiferromagnetic order appear simultaneously and the emergence of charge order is still under debate [11, 12, 13]; instead, the latter one is nonmagnetic and ferroelectricity is observed in the presence of pronounced charge order. These observations have opened a debate about the nature and interplay of charge order, ferroelectricity and magnetism in these materials, which will be at the focus of this study.

The building blocks of the $\kappa$-(ET)2X family (where X indicates a monovalent anion) are strongly-coupled dimers of ET molecules forming a triangular lattice. These materials have been widely studied within the half-filled single-band Hubbard model on the anisotropic triangular lattice, where only a single orbital per dimer is retained [14]. Indeed, because of the strong hybridization between ET molecules belonging to the same dimer, the gap between the bonding and anti-bonding orbitals is large; the former one is fully occupied, while the latter one is half filled, thus justifying a single-band picture. However, this coarse-grained approach cannot explain the possible emergence of ferroelectricity (or multiferroicity) in these materials, which has been suggested to arise from a charge disproportionation within each dimer. In this sense, the minimal model that could capture these features must include two molecular orbitals on each dimer and $3/4$ filling.

In the last decades there have been several attempts to obtain accurate values of the parameters defining microscopic models that would capture the low-energy properties of charge-transfer salts. The hopping integrals between the different molecular orbitals are found to significantly affect the nature of the ground states, as already reported in the first Hartree-Fock studies of correlated models for charge-transfer salts [15, 16, 17]. These considerations motivated a revision of the first estimates of the hopping parameters, that were based on the extended Hückel method [18], by means of density-functional calculations. Here, consistent results for the hopping parameters of the $\kappa$-(ET)2X family have been reported by three independent calculations [19, 20, 21], while slightly different values have been recently proposed [22]. Besides the role of hopping parameters, the importance of Coulomb interactions between different molecules in organic systems has been intensively discussed within ab-initio calculations [19, 23, 24, 25]. More recently, the analysis of various (low-energy) multiorbital models also points to the key role of intermolecular Coulomb interactions in order to describe complex phases relevant for charge-transfer salts. In particular, possible stripe and non-stripe charge orderings [26, 27] and the mutual exclusion of ferroelectricity and magnetism [28] have been discussed for various models with intermolecular interactions. In addition, the existence of a dipolar spin-liquid phase has been suggested [29, 30] (possibly also explaining the dielectric anomaly in $\kappa$-(ET)2Cu2(CN)3 [31]). Furthermore, the two-orbital Hubbard model has been claimed to be relevant for the description of superconductivity in charge-transfer salts [32, 33, 34, 35], including its proximity to charge-ordered phases [36, 37]. In addition, spin and charge fluctuations near the metal-insulator transition in multiorbital models have been analyzed [38].

In this paper, we concentrate on the question of what kind of charge orderings are driven by competing Coulomb interactions and which is their relation to ferroelectricity and magnetism. By using variational Monte Carlo methods, we investigate the phase diagram of an extended two-orbital Hubbard model on the triangular lattice at $3/4$ filling. Our results show that there exist two polar charge-ordered insulating phases, where charge disproportionation occurs within the dimer, and one nonpolar charge-ordered phase, with charge disproportionation between different dimers. All these phases are present also in the absence of magnetic order, indicating that they are not driven by magnetism. When magnetism is also included in the variational wave functions, we find that it coexists with charge order. These results could explain the observed behavior in $\kappa$-(ET)2Cu[N(CN)2]Cl. On the contrary, magnetism is crucial to stabilize the uniform dimer-Mott insulator, which appears in a narrow region between the two polar phases. In this respect, it has been experimentally suggested that a transition between the dimer-Mott insulator and charge-ordered states is a common feature among organic systems [39]. Finally, when intramolecular and intermolecular Coulomb interactions are small and similar in magnitude, a metallic phase emerges, featuring charge order in the form of an effective honeycomb-lattice superstructure.

The paper is organized as follows: In Sec. 2, we present the extended two-band Hubbard model for the organic charge transfer salts and the variational Monte Carlo method to study the phase diagram at zero temperature. In Sec. 3, we show the numerical results and discuss the nature of charge-ordered phases. Finally, in Sec. 4, we draw our conclusions.

2 Model and methods

2.1 The extended two-orbital Hubbard model

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 Fig. 1(a). An equivalent description is given by considering a two-orbital model on the square lattice, see Fig. 1(b). Here, we can define a partition in two sub-lattices ${\cal A}$ and ${\cal B}$, where the ET molecules form horizontal and vertical dimers, respectively. The full Hamiltonian, in this latter description, is given by:

[TABLE]

where

[TABLE]

Here, $c^{\dagger}_{i,\sigma}$ ($f^{\dagger}_{i,\sigma}$) creates an electron with spin $\sigma$ 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 Fig. 1(b). Belonging to the former class, there are terms connecting orbitals within the same dimer ($b_{1}$-type), along $x$ and $y$ nearest-neighbor sites ($q$-type), and along the $x=y$ diagonal ($b_{2}$-type); instead, $p$-type terms connect orbitals at nearest-neighbor sites and belong to the latter class. Accordingly, the noninteracting Hamiltonian ${\cal H}_{t}$ contains four hopping parameters, i.e., $t_{b1}$, $t_{b2}$, $t_{q}$, and $t_{p}$. Similarly, the interacting Hamiltonian ${\cal 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 ${\cal A}$ and ${\cal B}$ sub-lattices is only due to the presence of $p$-type terms. Finally, ${\cal 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$.

As discussed in Ref. [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}$):

[TABLE]

which are based on the results obtained by density-functional-theory calculations on $\kappa$-(ET)2Cu[N(CN)2]Cl [35, 41]. The noninteracting band structure is reported in Fig. 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}$.

2.2 The atomic limit

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 Fig. 3): Two of them are polar charge-order insulators (hereinafter denoted as PCOI and PCOI′), 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:

[TABLE]

where we defined:

[TABLE]

The phase diagram in the $V_{p}{-}V_{q}$ plane is shown in Fig. 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$.

2.3 Variational wave functions

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, 43, 44, 45], which are described as:

[TABLE]

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

[TABLE]

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 $\bm{R}_{i}=(x_{i},y_{i})$. $|\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 [46, 47, 40]. In particular, for the insulating states, we consider:

[TABLE]

where ${\cal H}_{t}$ is the kinetic part of Eq. (2) and

[TABLE]

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

[TABLE]

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 Eq. (14).

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

[TABLE]

where

[TABLE]

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

[TABLE]

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 $|\Phi\rangle$ as the ground state of the tight-binding Hamiltonian with site-dependent chemical potentials $\mu_{i}^{c}$ and $\mu_{i}^{f}$:

[TABLE]

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 Refs. [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 Eqs. (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.

3 Results

3.1 Phase diagram for large $U$

We now investigate how the hopping terms in Eq. (1) modify the phase diagram obtained for the atomic limit, in the region $U\gg t_{b1}$. The schematic phase diagram is shown in Fig. 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′ (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.

In addition to these three states, a uniform dimer-Mott insulator (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, 16, 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′) 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_{\rm CD}(q)$, defined as:

[TABLE]

where $\langle\dots\rangle$ indicates the expectation value over the variational wave function of Eq. (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 $|\bm{q}|\to 0$, $N(q)\propto|\bm{q}|$ for a metal, while $N(q)\propto|\bm{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_{\rm 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}=3t_{b1}$. As shown in Fig. 5, all the phases presented in the phase diagrams are insulating, since $N(q)\propto|\bm{q}|^{2}$ in the limit $|\bm{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_{\rm CD}(q)$, see Figs. 6 and 7. The DMI does not show any Bragg peak either in $N(q)$ or in $N_{\rm CD}(q)$ [Figs. 6(d) and  7(d)], suggesting that no long-range charge order occurs. The nonpolar charge-ordered state shows the Bragg peak at $\bm{Q}=(\pi,\pi)$ in $N(q)$ [Fig. 6(c)], but no sharp peaks in $N_{\rm CD}(q)$ [Fig. 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 $\bm{Q}=(0,0)$ (PCOI′) or $\bm{Q}=(\pi,\pi)$ (PCOI) in $N_{\rm CD}(q)$ [Fig. 7 (a-b)], but no sharp peaks in $N(q)$ [Fig. 6 (a-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 $\bm{Q}=(0,0)$, while each orbital alternates between charge-rich and charge-poor configurations for $\bm{Q}=(\pi,\pi)$, see Fig. 4.

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 ${\cal H}_{\rm AF}$ of Eq. (14), see Fig. 8. The charge patterns are similar to the ones that have been obtained previously with the inclusion of magnetic order (not shown).

3.2 Competition between charge and magnetic orders

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}=3t_{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 Eq. (13)] can be used as a diagnostic to detect polar and nonpolar states. As shown in Fig. 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.

To further investigate the connection among these three phases, we calculate the charge-disproportionation structure factor as function of $V_{p}$, as shown in Fig. 10. When $V_{q}$ ($V_{p}$) is sufficiently large, $N_{\rm CD}(q)$ shows a sharp peak at $\bm{Q}=(0,0)$ ($\bm{Q}=(\pi,\pi)$) [Fig. 10(a) and (b), respectively]. By contrast, when $V_{p}\approx V_{q}$, $N_{\rm CD}(q)$ shows a broad crest along the $q_{x}+q_{y}=2\pi$ direction [Fig. 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_{\rm 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 Fig. 1. It is then natural to expect that correlations do not show any dependence on the transverse direction.

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_{\rm CD}(Q)/L^{2}$ for $L\to\infty$. The results are reported in Fig. 11. For $V_{p}/t_{b1}\lesssim 1.3$ and $V_{p}/t_{b1}\gtrsim 1.7$, we have that $N_{\rm CD}(Q)/L^{2}$ is finite for $\bm{Q}=(0,0)$ and $\bm{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.

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 Fig. 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}|>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′ 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′ one they are connected by $t_{p}$; since $t_{p}>t_{q}$ [see Eq. (5)] the latter case has more charge fluctuations (i.e., it is closer to a metal-insulator transition), thus implying smaller magnetic moments.

3.3 $12$-site ordered metallic phase

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 Fig. 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.

The emergence of the honeycomb-like COM can be easily understood when all the bonds (of $b_{1}$-, $b_{2}$-, $p$-, and $q$-type) are equivalent. In this “isotropic” limit, when considering each molecule as an independent site, the underlying lattice becomes triangular, see Fig. 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{-}B{-}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 Fig. 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 Fig. 15. Both phases are indeed metallic, since $N(q)\propto|\bm{q}|$ for small momenta [Figs. 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_{\rm CD}(q)$, corresponding to the real-space configuration illustrated in Fig. 13(b).

4 Summary and conclusions

By using variational wave functions and quantum Monte Carlo techniques, we have investigated the ground-state phase diagram of an extended two-orbital Hubbard model at $3/4$ filling on the anisotropic triangular lattice, which is relevant for the $\kappa$-(ET)2X family of organic charge-transfer salts. As a representative example, we have chosen the hopping parameters that correspond to $\kappa$-(ET)2Cu[N(CN)2]Cl and varied the interaction terms. For large values of the intramolecular repulsion $U$ and by varying the strength of the competing intermolecular Coulomb interactions, we stabilize two polar and one nonpolar charge-ordered insulating phases, as well as a uniform dimer-Mott insulator. All these phases posses magnetic order, which is mainly determined by the behavior of the spins at the charge-poor molecules (that are effectively at half filling). We have also found that the dimer-Mott insulator is continuously connected to the two polar charge-ordered states: When the anisotropy between the intersite Coulomb interactions $V_{p}$ and $V_{q}$ goes to zero, the Bragg peaks of the two polar phases melt and form a one-dimensional-like structure. For smaller values of the intramolecular interaction $U$, we find a charge-ordered metal, that is similar to the three-sublattice (rich-rich-poor) charge order on the triangular lattice; however, the anisotropy in the intermolecular parameters modify the period of the charge ordering to a $12$-sublattice structure. Although charge-ordered metals in ET organic compounds often show a stripe-like charge pattern [52, 53], the observation of the COM phase would be an intriguing proof for the possibility of stabilizing nontrivial charge orders in metallic phases.

In organic charge transfer salts, the size of the intermolecular Coulomb interactions is expected to be larger when the molecules are closer. In this respect, it is plausible to assume that $V_{b1}$ is the stronger intermolecular Coulomb interaction and that $V_{p}\gtrsim V_{q}$, see Fig. 1. In addition to the fact that the strongest Coulomb interaction is the intramolecular one $U$, most of the compounds should be located at the border between the PCOI and the DMI phases. Since the two phases are continuously connected, a small amount of anisotropy $V_{p}\gtrsim V_{q}$ will lead to a weak charge order, as shown for example in Fig. 9; this fact may explain the difficulty in finding stable charge ordering in $\kappa$-(ET)2Cu[N(CN)2]Cl. Nevertheless, our results indicate that the PCOI phase is polarized, suggesting that charge order is the correct mechanism to induce a finite polarization. Moreover, we observe that magnetism coexists with electronic polarization, as observed in experiments, even if it is not the driving mechanism for it, since polarization occurs also in the absence of magnetic order. In this respect, our study shows that ferroelectricity in organic charge-transfer salts is not driven by magnetism.

Finally, we would like to conclude by mentioning that superconducting pairing correlations (with unconventional pairing symmetries) may be enhanced close to charge-ordered phases in multiorbital Hubbard models [37]. Investigating possible superconductivity (also coexisting with charge ordering) is left for future studies.

The authors would like to thank M. Altmeyer, C. Gros, D. Guterding, A.J. Kim, and S.M. Winter for fruitful discussions. R.K. and R.V. acknowledge the support of the German Science Foundation (Deutsche Forschungsgemeinschaft) through Grant No. SFB/TR49.