{'fraction_non_alphanumeric': 0.036326851068969684, 'fraction_numerical': 0.010085759149655756, 'mean_word_length': 4.529637669059944, 'pattern_counts': {'":': 1, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Phase diagram of the Hubbard model on the anisotropic triangular lattice

We investigate the Hubbard model on the anisotropic triangular lattice as a suggested effective description of the Mott phase in various triangular organic compounds. Employing the variational cluster approximation and the ladder dual-fermion approach as complementary methods to adequately treat the zero-temperature and the finite-temperature domains, we obtain a consistent picture of the phase diagram as a function of anisotropy and interaction strength. The metal-insulator transition substantially depends on the anisotropy, and so does the nature of magnetism and the emergence of a nonmagnetic insulating phase. We further find that geometric anisotropy significantly influences the thermodynamics of the system. For increased frustration induced by anisotropy, the entropy of the system increases with interaction strength, opening the possibility of adiabatically cooling a frustrated system by an enhancement of electronic correlations.


I. INTRODUCTION
The single crystals of the organic charge transfer salts recently received substantial interest, where fascinating phenomena such as several emergent many-body phases are observed. κ-(BEDT-TTF) 2 Cu[N(CN) 2 ]Cl and κ-(BEDT-TTF) 2 Cu 2 (CN) 3 are two prototypical examples of such organics that show a wide range of similarities, but are also clearly different from each other. Fermi liquid, Mott insulator behavior, and the crossover between the two phases are similarly detected in both materials 1,2 . The metal-insulator transition is found to be associated with unusual critical exponents which, face value, do not fall into any established universality class 3 . What, however, sparks even more substantial interest is the striking difference of these two materials in terms of their magnetic behavior. κ-(BEDT-TTF) 2 Cu[N(CN) 2 ]Cl displays longrange magnetic order at low temperatures, which is in sharp contrast with the apparent spin liquid, or at least strongly magnetically frustrated, behavior of κ-(BEDT-TTF) 2 Cu 2 (CN) 3 . The bulk spin susceptibility 4 of κ-(BEDT-TTF) 2 Cu 2 (CN) 3 reflects this strong frustration, as it shows no indication of long-range antiferromagnetic order at temperatures significantly lower than the scales inferred from the Heisenberg exchange estimated theoretically from high-temperature series expansions 5 . These results highlight the similarity of the two compounds in the charge sector, but also the clear difference regarding the spin degrees of freedom, which also manifests itself in effective field theory description of the problem [6][7][8][9] . Separating charge from spin degrees of freedom, this paves the way for the investigation of effective spinon theories in the Mott insulating phase (with or without gapless spinon modes yielding a potentially unstable spinon Fermi surface 10 ).
These effective field theories, however, can take on dif-ferent perspectives. For example, the thermal conductivity is found to be contributed by the spin-1/2 spinons in a theory of U (1) gauge field coupled with a spinon Fermi surface 6 , but associated with visons in a Z 2 theory 7 . A similarly puzzling situation as in theory likewise exists on the experimental side, where the thermodynamic understanding of the compounds is far from settled. Specific heat 11 and thermal conductivity 12 experiments can be interpreted to yield different conclusions on the nature of the fermionic excitations in κ-(BEDT-TTF) 2 Cu 2 (CN) 3 . On top of all these complications, even if we assume a spin liquid state in the latter compound, is it still debated whether this state would be fully gapped or not 13,14 .
What are the microscopic parameters whose variations impose such a diversity of many-body phenomena in these organic compounds? While the interaction strength might actually be rather comparable in all compounds, a clear difference that catches the eye lies in the anisotropy strength inherently determined by the underlying crystal structure and chemical components. The anisotropy strength can be obtained from ab initio calculations [15][16][17] . Assuming the Hubbard model on a triangular lattice be the correct model for describing the interplay of geometric frustration with electronic correlations in these compounds, it is a natural further step of complexity to consider the lattice anisotropy.
In this paper, we study the triangular Hubbard model with varying anisotropy strengths, devising methods to treat both the zero-and finite-temperature regime. The Hubbard model on an anisotropic triangular lattice is given as: from the limiting values 0 to 1, we effectively change the lattice geometry from square to triangular type. Unless stated otherwise in the paper, we choose the phrasing of small and large anisotropy according to the value of t /t. The most interesting regime for the organic compounds is located around the isotropic triangular limit t /t = 1, where apparently only small variations can yield crucially different scenarios. It is still instructive to see how the phase diagram evolves from the square limit towards the triangular limit, which is why we analyze the complete domain 0 ≤ t /t ≤ 1.1. For t /t > 1.1, the system quickly evolves towards an effectively one-dimensional scenario of weakly coupled chains, quickly rejecting magnetic order 18 . This is neither a relevant regime for the organic compounds nor particularly suited for the methods employed in this work, which is why the regime t /t > 1.1 will not be further addressed.
The paper is organized as follows. Sec. II briefly introduces the methodology of the VCA and the LDFA. The zero and finite temperature phase diagrams are obtained by VCA and LDFA, respectively, providing a complementary and consistent perspective on the anisotropic triangular Hubbard model. Our results are discussed in Sec. III. In particular, as a function of the increasing anisotropy parameter, we find a magnetic phase transition from Néel-AFM to 120 • -AFM order along with a growing metallic regime. Furthermore, already rather nearby the triangular limit, the onset of a non-magnetic insulating regime is found, which is the candidate domain for possible spin liquid states, where the charge degree is frozen without spin ordering. The LDFA additionally offers the possibility of addressing questions of thermodynamics in the Hubbard model. In particular, we find indication for adiabatic cooling caused by the change of frustration as a function of interaction strength, which might be observed in highly tunable scenarios such as triangular optical lattices loaded with ultra-cold fermionic isotopes. In Sec. IV, we conclude that our analysis sets the initial stage for further refined investigations of the many-body phases in the anisotropic triangular Hubbard model, for which we can identify the promising spin liquid candidate regime. Whatever these unconventional phases may be and which effective theories they are best described by, the lattice anisotropy is likely to be a crucial microscopic parameter.

II. METHODOLOGY
In this section, we briefly review the methods we employ in this paper, namely the VCA and the LDFA. These two methods will be subsequently applied for zero and finite temperature, respectively, focusing on the quantum phase diagram and certain thermodynamic properties.

A. T=0: Variational Cluster Approach
The VCA 19 is based on the Self-Energy-Functional Theory (SFT) 20,21 , which provides an efficient numerical technique for studying strongly correlated systems, especially in the presence of different competing orders. VCA simplifies the lattice problem, as defined in Eq. (1), to an exactly-solvable problem defined in a reference system consisting of decoupled finite-size clusters. The thermodynamic limit is recovered by reintroducing the inter-cluster hopping to the decoupled cluster via a non-perturbative variational scheme based on self-energy functional theory. The VCA has been successfully applied to many interesting problems, including the high-T c cuprates 22-26 and topological insulators [27][28][29] . It is a method that is particular suitable for our current study, as the anisotropy induces several phase transitions in our geometrically frustrated system. The competing nature of these phases can be conveniently described by the VCA grand potential.
In the SFT, the grand potential of a system defined by H = H 0 (t) + H 1 (U) is written as a functional of the self-energy Σ: where F [Σ] is the Legendre transform of the Luttinger-Ward functional and G 0 = (ω + µ − t) −1 is the noninteracting Green's function. It can be shown that the functional Ω[Σ] gets stationary at the physical selfenergy, i.e. δΩ [Σ phys ] = 0 20 . As the Luttinger-Ward functional is universal, it has the same interaction dependence for systems with any set of t as long as the interaction U remains unchanged. Note that the functional Ω [Σ] itself is not approximated by any means; we restrict, however, the "parameter" space of possible self-energies to the self-energies of the reference system. Thus, the stationary points are obtained from the self-energy Σ = Σ [t ] of a system defined by H = H 0 (t ) + H 1 (U), which we name the reference system. In terms of the reference system, the VCA grand potential is calculated conveniently as with Ω , Σ , and G denoting the grand potential, the self-energy and the single-particle Green's function of the reference system, respectively. The reference system is chosen so that it can be exactly solved, i.e. Ω , Σ , and G should be readily obtained from, e.g. numerically exact methods. Here, we choose the reference system as a set of decoupled clusters and solve them with open boundary condition using exact diagonalization. In this sense, the short-range correlations within the reference system are fully taken into account in the VCA. The correlation beyond the reference system size will be treated in a mean-field manner via the variational scheme. The choice of the reference system, i.e. the cluster shape and size, has to respect the fact that tuning the anisotropy from 0 to 1 effectively modifies the system geometry from square to triangle. Thus, a reference system which incorporates both the Néel and spiral order in the VCA becomes essential. In this paper, two reference systems are explored with the cluster size of L c = 6 and L c = 12. They are mirrored in our calculation to recover the lattice translation symmetry with a supercluster. 30 The Green's function of this supercluster G consists of the cluster + mirrored cluster is given as, where G 1 is the reference cluster Green's function and G 2 is the Green's function of the mirrored cluster, which is a simple transformation of G 1 (in the simplest case, it is just a copy of G 1 ). The reference and the mirrored clusters are connected through the single-particle hopping t 12 , as indicated by the dotted links in Fig. 1.

B. Finite-T: Ladder Dual-Fermion Approach
For the finite-temperature study, we employ the LDFA 31,32 . LDFA decouples a correlated lattice defined in Eq. (1) into a group of impurities which couple to each other through an effective interaction mediated by auxiliary dual-fermions. The local problem can be well described within the dynamical mean-field theory (DMFT) 33 . The perturbation expansion over the effective interaction of the dual-fermion variables can systematically generate non-local corrections to the DMFT. The basic idea of the dual-fermion approach is schematically shown in Fig. 2, where the lattice problem defined in Fig. 2(a) is decoupled into an impurity problem as in the DMFT ( Fig. 2(b)). While the difference of the lattice and decoupled impurity problems are treated perturbatively in the dual-fermion approach, as shown in Fig. 2

(c).
To be more precise, let's formulate the basic idea in Fig. 2 by starting with the lattice action corresponding to the Hamiltonian in Eq. (1).
dτ is the action of an impurity coupled with a continuum bath which is described dynamically by a hybridization function ∆(iω n ). The second sum in Eq. (5) is the term that we are going to treat perturbatively.
FIG. 2. As a non-local extension of the DMFT, the dual fermion method perturbatively expands the difference of the one particle hopping and the DMFT hybridization function, i.e. k − ∆(iωn), which generates systematic non-local corrections to the DMFT solution. Ideally, with all the expansion terms taken into account, the lattice problem defined in Fig. 2(a) can be exactly solved by the dual fermion method.
If the local action S loc is already a good description of the original system, the second term in the r.h.s. of Eq. (5) will effectively be a small parameter. A expansion over this small term, i.e. k − ∆(iω n ), generates further corrections to S i loc and can be calculated order by order. A convenient way for such an expansion is to rewrite the second term in the r.h.s. of Eq. (5) with a dual variable f, f * through the Gaussian integration. This transformation is a pure mathematical reformulation. After integrating out the c-variables, the original lattice problem described in Eq. (5), now can be equally written as an action that only depends on f -variables: The interaction between the f -variables, i.e. V i in Eq. (6), is the reducible two-particle vertex of the cvariables, which can also be obtained by solving S i loc with an appropriate impurity solver, such as the continuoustime Quantum Monte Carlo method 34 . Eq. (5) and Eq. (6) are two equivalent ways to describe the same problem as no approximation is introduced in the transformation. Thus, the lattice Green's function G k (iω n ) can be equally constructed from this new action, which reads where g d (iω n ) is given as: Now it becomes clear that the problem of solving an interacting many-body problem defined in Eq. (1) is equivalent to solving an Anderson impurity problem, i.e. selfconsistently determine g(iω n ) and ∆(iω n ), and additionally calculating Σ d k (iω n ) from the perturbation expansion of V i [f * , f ] in Eq. (6). In the following calculation, we will consider ladder-type diagrams of Σ d k (iω n ) up to infinite order 32 in the particle-hole channel. The expression of Σ d k (iω n ) from the first two orders expansions can be found in reference 35. Eq. (7) sets up an exact relation for the single-particle Green's function of the original c-variables and the dual f -variables. Similar exact relations can also be found for higher-order correlators. For example, for the spin susceptibility, which is employed to identify different magnetic phases in this work, we have stands for the reducible vertex of the dual variables. The high momentum-resolution of the spin susceptibility calculated from the LDFA is very hard to achieve in other approaches. It is essential in the study of the spin structure at different anisotropy strengths, as the magnetic order changes from Néel-type (with magnetic wave vector (π, π)) to a 120 • -type (with magnetic wave vector (2π/3, 2π/3)). In the spin susceptibility, the corresponding momentum Q (where χ sp Q diverges) directly provides the information of the magnetic ordering, as this momentum is the magnetic wave vector of the corresponding long-range order. Fig. 3 summarizes the zero temperature phase diagram as a function of anisotropy and interaction strength as obtained by VCA calculations. As stated in the introduction, the geometrical frustration is parametrized by the size ratio of the horizontal and diagonal hoppings in Fig. 1. Changing t /t connects two geometrically different lattices, i.e. the unfrustrated square lattice and the frustrated triangular lattice. We expect the phase diagram to rather continously interpolate between the established square and triangular limit, which is indeed what is observed in Fig. 3.

III. RESULTS
In the strongly correlated limit, Eq. (1) is equivalent to the Heisenberg model, and the system develops the Néel-AFM at t /t = 0 36 and the 120 • -AFM at t /t = 1 [37][38][39] . Varying the anisotropy strength, i.e. changing the value of t /t, we observe a transition between these two magnetically ordered states. The Néel-AFM state is surprisingly stable against the geometric frustration, as, among other methods, has also been found in more sophisticated functional renormalization group studies 18 : t /t has to be larger than 0.89 in order to destroy the collinear antiferromagnetic order to establish the 120 • -AFM. In contrast, in the weak-coupling limit, the geometric frustration plays a much more significant role already for small t /t, i.e. it stabilizes the metallic phase in the entire range of t /t. In the square lattice, due to the Fermi surface nesting, the ground state of the system is antiferromagnetically ordered for any finite Coulomb interaction, with the nesting vector being the magnetic vector. (The small offset seen in Fig. 3 is a minor finite size artefact in the VCA.) Due to the reduced symmetry, the system opens up a gap in the charge sector. Thus, along the t /t = 0 line, our anisotropic triangular Hubbard model becomes an antiferromagnetic insulator for arbitrary finite U/t. With the increase of t /t, the metallic state is stabilized and extends up to larger values of U/t. For small t /t, the metal to insulator transition (MIT) is mainly due to the development of the antiferromagnetic order. Accordingly, the boundary for the metal to Néel-AFM transition conincides with the boundary for the MIT transition as obtained in VCA.
With larger anisotropy, the MIT as a function of coupling strength acquires a different character. In Fig. 3, we observe a transition from a metallic "phase" to a nonmagnetic insulating (NMI) phase, which is induced by the charge degrees of freedom in the intermediate coupling region. Approximately for 5.2 < U/t < 9, the system opens up a charge gap without developing longrange magnetic order when t /t > 0.7 (see Fig. 3 for the precise boundaries for the metal-NMI and NMI-magnetic transition.) The NMI phase is the natural regime where one or several kinds of quantum spin-liquid (QSL) phases might be located. The strong geometric frustration completely suppresses the magnetic ordering in this coupling region, thus, only the charge degrees of freedom dominate this phase. Note that since our numerical methods are adjusted to the computation of single-particle quantities, it is impossible for further analyze the specific properties of the NMI phase, which would be indispensable to make concrete predictions for the spin liquid states to be expected.
The appearance of the NMI phase in the intermediate coupling region qualitatively agrees with other theoretical investigations 30,[40][41][42][43][44][45][46][47] . In these works, including ours, the emergence of the NMI phase are consistently shown to be due to the competition of electronic correlation and geometric frustration. However, the size of the NMI phase differs. For VCA, this is partially due to the fact that the lower bound of the NMI phase slightly depends on the cluster size. The metal to NMI phase transition occurs at U/t = 5.4 for L c = 3, U/t = 5.25 for L c = 6, U/t = 6.3 for L c = 12. As electrons gain more mobility in a larger cluster, the kinetic energy of the ground state will become lower. Recently, we discovered this process can be in general more efficiently simulated by introducing another variation parameter, i.e. t, to the VCA procedure 28 : Hopping t describes the itinerancy of a single electron. The variation of t, thus, allows to minimize the kinetic energy, which largely recovers the same physics in a small cluster that woud likewise emerge in a larger cluster. In contrast, the upper bound of the NMI phase, i.e. the NMI to 120 • -AFM phase transition boundary, is determined by the collective behaviors of all electrons in the system, and thus is less affected by the single electron hopping variation. As a result, we find that, by varying t, the lower bound of the NMI phase becomes U/t = 7.2 for the isotropic triangular lattice with L c = 12, while the upper bound U/t = 9.4 is nearly unchanged.
In our VCA calculations, the metal to insulator transition boundary is determined from the opening of the single particle gap ∆ sp . As shown in Fig. 4(a), the MIT transition takes place at U = 5.25 for t /t = 1. The boundary between different magnetic phases as well as the nature of the phase transitions are identified from the comparison of the VCA grand potential in different phases. The stable phase possesses the lower grand potential and a phase transition is characterized by the crossing of the grand potential of different phases. In Fig. 4(b), we show an example of the 120 • -AFM to NMI transition for t /t = 1.0. The two VCA grand potentials cross at U/t = 9.4, which indicates that this transition is of first order. (In a phase transition, if the two grand potentials smoothly change from one to another without any crossing, the transition will be of a second order or even higher.) The stabilization of the metallic state due to geometric frustration can also be seen from the finite-temperature LDFA calculations displayed in Fig. 5(a). In this figure, the Néel temperatures are displayed as a function of interaction for different anisotropy strengths. Strictly speaking, there is no finite-temperature magnetic transition in two dimension, according to the Mermin-Wagner theorem. 49 . Such a transition appears in a method that includes certain implicit IR cutoffs, such as given by the partial mean field character in DMFT and LDFA. Still this information instructive as the magnetic correlations are correctly described in this type of calculations. Furthermore, the finite-temperature magnetically ordered phase can be realized in a slab of multilayer triangular systems. A small but finite interlayer coupling introduces an additional dimension to the system and, therefore, will favor the magnetic symmetry breaking at finite-temperature. In our study, the Néel temperatures are obtained from the extrapolation of the inverse spin susceptibility. Following Eq. (9), we calculate χ s Q (iΩ m = 0) in the entire first Brillouin Zone (BZ) and extrapolate the inverse of the leading value of the temperature-dependent χ Q (iΩ m = 0) to T /t = 0. The Néel temperature is simply given as the temperature where the extrapolated χ −1 Q (iΩ m = 0) becomes zero. Examples of the spin susceptibilities and the extrapolation can be found, for example, in Fig. 5(b) for t /t = 0.6.
As shown in Fig. 5(a), increasing the anisotropy strength greatly suppresses the Néel temperature, especially in the weak-coupling region. At t /t = 0, the square lattice is recovered and the Néel temperature is non-zero for any finite interaction, indicating that long-range magnetic correlations are well established for U/t = 0. In this case, the critical value for the paramagnet to antiferromagnet transition is U c = 0. Increasing t /t results in the suppression of the Néel temperature, as a result U c increases. This coincides with what is obtained from the T = 0 VCA calculations shown in Fig. 3. However, we would not expect an exact quantitative agreement on the critical values U c from these two methods, as both of them are subject to certain approximations. Still, note that both methods are approximation-controlled, i.e. a VCA calculation with an infinite-size reference system will exactly coincide with the calculation of the dualfermion method with all non-local self-energy diagrams included. On the current level of approximation (see Sec. II for more details about the reference system size and the self-energy diagrams considered), Fig. 5(a) and Fig. 3 converge to the same conclusion, i.e. that the enhanced spatial anisotropy stabilizes an extended metallic phase and suppresses the propensity to magnetic ordering.
When the system approaches the antiferromagnetic state, the magnetic correlation, on the one hand, drives the spin susceptibility divergent, on the other hand, it also leaves fingerprints in the single-particle spectra. In the case of the square lattice, with the development of the (π, π)-antiferromagnetism, the effective magnetic unit cell becomes twice of the size of the original unit cell. The single-particle spectra then pick up the new symmetry associated with the magnetic unit cell, which results in a "shadow band" around the Γ point 50 . In the fully isotropic triangular case, the magnetic correlation is of 120 • type. The resulting magnetic unit cell is then three times the size of the original unit cell. The original band will further be folded with respect to the magnetic zone boundary, which would also generate a shadow band around Γ. Thus, detecting the appearance of the shadow band can actually help us to track the magnetic correlations of the system. Comparing dynamical mean field theory (DMFT) and LDFA data is instructive to highlight the additional nonlocal corrections kept in the latter: A detrimental problem of the DMFT lies in the local approximation, making it incapable of describing long-range correlations. Thus, the shadow band induced by the magnetic correlation should be less obvious or even absent from a DMFT calculation. In Fig. 6, we show the single-particle spectra for t /t = 0.8, U/t = 9 and T /t = 0.154. The chosen temperature is slightly above the Néel-temperature, where the magnetic correlations is about to fully develop. Thus, we would expect new folded bands to appear in the single particle spectra. In Fig. 6(b), as indicated in the circle, there is a band with less intensity to develop around Γ. This band is absent in the tight-binding model of Eq. (1) and is a direct result of the band backfolding with respect to the magnetic zone boundary. The new symmetry point induced by the enlarged magnetic unit cell is indicated by the arrow in Fig. 6(b). This "shadow band" locates at a finite binding energy at the Γ point, resulting a slope change of the band above the Fermi level, as indicated by the blue dashed line in Fig. 6(b). This magnetic ordering induced shadow band is not correctly resolved in the DMFT calculation. As displayed in Fig. 6(a), the band close to Γ is not a band folded from the emergent symmetry, i.e. the 120 • symmetry, it is only the remanent of the band between M and K. It gradually approaches to zero binding energy at Γ point, as can be seen by following the blue dashed line in this figure. Another clear difference between the DMFT and the LDFA results lies in the suppression of the quasi-particle peak at the Fermi-level. This is again due to the non-locality missed in the local approximation of the DMFT. The transition between the two types of magnetic order in Fig. 3 is also observed at finite-temperature. In Fig. 7 the momentum dependent susceptibilities χ Q (iΩ m = 0) are shown in the entire first BZ at U/t = 9 and T /t ∼ 0.286, for four representative values of t /t. As discussed in Sec. II, the instability of the paramagnetic solution reveals the formation of a magnetic ordered phase. As shown in Fig. 7(a) and (b), with smaller anisotropic strengths, χ Q (iΩ m = 0) shows a single-peak located at Q = (π, π), the Néel-AFM. The increase of t /t results in peak broadening as for the case of t /t = 0.6 depicted in Fig. 7(b). A further increase of the anisotropy strength splits the peak and generates a double-peak structure of the spin susceptibility, which highlights the evolution from Néel-AFM to the 120 • -AFM. At t /t = 1.0, the valley between the two peaks at (2π/3, 2π/3) and (4π/3, 4π/3) becomes even deeper. It coincides with what we discussed in Fig. 3.
At finite temperatures, the LDFA allows us to analyze the thermodynamics of the anisotropic triangular Hubbard model. In one of our recent works 51 , by using the LDFA, we showed that the thermodynamics, especially the entropy, substantially enhances our understanding of the competing roles of geometric frustration and electronic correlations. We found that geometric frustration favors the effect of "adiabatic cooling": i.e. following a constant entropy curve, increasing the interaction results in the effective decrease of temperature. This is in contrast to the situation in square 52 and cubic 53 lattices, where the entropy is nearly a constant function of interaction in the weak-coupling region. However, "adiabatic cooling" has been theoretically realized in the honeycomb type lattice 54 , which can be viewed as two triangular sublattices. Therefore, we speculate the geometric frustration enhances the decrease of entropy with the increase of interactions, which will be discussed in the following. Here U/t = 6 and for all temperature studied, the system remains metallic.
In Fig. 8, the double occupancy D = n ↑ n ↓ is displayed, whose low temperature behavior reveals more information than just the degree of electron localization. The double occupancy relates to the entropy S through a Maxwell relation, i.e.
As also discussed in reference 51, the negative-entropy slope for temperatures lower than a characteristic temperature T * indicates that the entropy will increase with the increase of interactions. As shown in Fig. 8, for temperatures smaller than T * /t ∼ 0.8, in the isotropic triangular case (i.e. t /t = 1), D decreases upon increasing T 51 . It becomes more intuitive, if we rewrite Eq. (10) as where C v is the specific heat. This immediately implies that keeping the entropy constant, an increase of U results in a decrease of T for T ≤ T * . The influence of the frustration effect becomes clear in Fig. 8: with the reduction of the anisotropy, the negative slope of the double occupancy below T * becomes less obvious. In the unfrustrated limit, i.e. t /t = 0, the double occupancy would resemble that in a square lattice which indicates no "adiabatic cooling" is possible. At a fixed temperature below T * , frustration increases the value of the double occupancy resulting in the enhancement of its negative slope. We conclude the frustration to be the reason for the "adiabatic cooling" in the triangular system. In addition to the change of the slope, we find the characteristic temperature T * becomes slightly smaller with the decrease of t /t (see the solid arrows around T /t ∼ 0.8 in Fig. 8).
In the curves for the t /t = 0.4, 0.6, we also observe a second characteristic temperature T due to the evolution of system geometry from square to triangle, see the dashed arrows around T /t ∼ 0.25.

IV. CONCLUSIONS
We have conducted a detailed single-particle spectra analysis of the anisotropic triangular Hubbard model for zero and finite temperature. Focusing on the role of anisotropy and interaction strength, we have identified the significant features of the phase diagram which exhibits interesting regions such as the magnetic transition regime between Néel-AFM and 120 • -AFM order as well as, in particular, a non-magnetic insulating (NMI) regime setting in at an anistropy value of t /t which quickly grows in range of U/t as a function of anisotropy along with a more and more extended metallic regime for weaker coupling. It is exactly the NMI regime which might prove most relevant for the unconventional organic compounds. Identified as the NMI domain, one or sev-eral spin liquid phases can potentially appear in this window of parameter space. While this is beyond the framework of the current investigation, it will be fascinating to follow up on this NMI regime and to adapt methods capable of calculating multi-particle vertices in order to identify the specific state at a given point in the NMI phase. Since methods such as variational Monte Carlo can serve this goal but are rather time-consuming, our VCA/LDFA study is helpful in that it constrains the interesting parameter window to be scanned.
At finite temperature, in line with our T = 0 study, we find that the anisotropy substantially suppresses the magnetic ordering of the system, even down to zero temperature. The formation of the shadow band at the Γpoint in the LDFA calculations shows that, going beyond DMFT, the LDFA is capable of describing magnetic correlations due to the inclusion of non-local correlations. A characteristic temperature T * is identified in the double occupancy, below which the double occupancy decreases upon increasing temperature. This opens up the possibility of 'adiabatically cooling' the system by increasing the interactions, while keeping the entropy constant. We find that the geometrical anisotropy favors 'adiabatic cooling', i.e. increasing the anisotropy results in a larger negative slope of the double occupancy below T * . This is a fascinating effect which might be detectable in highly tunable condensed matter systems, such as appropriate optical lattices loaded with a cold gas of fermionic isotopes.

Note added.
While finalizing the draft of this manuscript, we became aware of a related VCA study of the anisotropic triangular Hubbard model by A. Yamada 55 . The subset of VCA findings contained in our paper agrees with this study.