# Mott Transition and Phase Diagram of -(Bedt-Ttf)Cu(NCS) Studied by Two-Dimensional Model Derived from Ab initio Method

Mott Transition and Phase Diagram of -(BEDT-TTF)Cu(NCS) Studied by Two-Dimensional Model Derived from Ab initio Method

Hiroshi Shinaoka^{†}^{†}thanks: E-mail address: h.shinaoka@aist.go.jp, Takahiro Misawa, Kazuma Nakamura, and Masatoshi Imada

Nanosystem Research Institute, AIST, Tsukuba 305-8568

CREST, JST, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-8656

Department of Applied Physics, University of Tokyo, Hongo, Bunkyo-ku, Tokyo 113-8656

(Received July 6, 2019)

We present an ab initio analysis for the ground-state properties of a correlated organic compound -(BEDT-TTF)Cu(NCS). First, we derive an effective two-dimensional low-energy model from first principles, having short-ranged transfers and short-ranged Coulomb and exchange interactions. Then, we perform many-variable variational Monte Carlo calculations for this model and draw a ground-state phase diagram as functions of scaling parameters for the onsite and off-site interactions. The phase diagram consists of three phases; a paramagnetic metallic phase, an antiferromagnetic (Mott) insulating phase, and a charge-ordered phase. In the phase diagram, the parameters for the real compound are close to the first-order Mott transition, being consistent with experiments. We show that the off-site Coulomb and exchange interactions affect the phase boundary; (i) they appreciably stabilize the metallic state against the Mott insulating phase and (ii) enhance charge fluctuations in a wide parameter region in the metallic phase. We observe arc-like structure in Fermi surface around the region where the charge fluctuations are enhanced. Possible relevance of the charge fluctuations to the experimentally observed dielectric anomaly in the -BEDT-TTF family compounds is also pointed out.

KEYWORDS: organic conductors, first principles, Hubbard-type low-energy model, variational Monte Carlo method, Mott transition

## 1 Introduction

In organic conductors, we find diverse properties and a plenty of phases including normal metals, superconductors and various types of insulators with antiferromagnetic, charge or (spin) Peierls orders as well as spin-liquid-type nonmagnetic Mott insulators. Although the unit cells of these conductors contain many atoms constituting the molecules with complicated crystal structures, band structures near the Fermi level are in most cases simple with a small number of bands isolated from other bands located away from the Fermi level. These isolated bands originate from molecular orbitals [lowest unoccupied molecular orbitals (LUMO) and highest occupied molecular orbitals (HOMO)]. Because the number of bands near the Fermi level is small, screening of the electron-electron Coulomb interaction is poor. In addition and more importantly, large lattice constants make the overlap of the neighboring molecular orbitals small, leading to a large ratio of the screened electron interaction to the kinetic energy. This is the reason why electron correlations are in general strong in the organic conductors. Because of the complex unit cell and prominent strong correlation effects, the ab initio calculation and clarification of mechanisms of material properties in the organic conductors remain as big challenges.

Among all, a family of compounds, (ET) with a number of choices of anions alternatingly stacked with BEDT-TTF molecules [where BEDT-TTF is bis(ethylenedithio)-tetrathiafulvalene, abbreviated as ET], offers a variety of prototypical behaviors of strongly correlated electron systems with two-dimensional (2D) anisotropies. In particular, the -type compounds characterized by dimerization of ET molecules have served to discoveries of unconventional quantum phases and unexplored concepts at the forefront of condensed matter physics.

Unconventional superconductivity is found in some of these compounds that show metallic behaviors. Namely, the compounds with the anions =Cu[N(CN)]Br (ref. ?) and =Cu(NCS) (refs. ?, ?) abbreviated as -Br and -NCS hereafter, respectively, show superconducting transitions at 10-13K. Under pressure, the compound with Cu(CN) referred to as -CN also shows superconductivity below 2.8 K. However, the driving mechanism of the superconductivity is not fully understood yet.

An unconventional nonmagnetic Mott-insulating phase found near the Mott transition for -CN is another example of such a discovery. In contrast to a naive expectation for a magnetic order at sufficiently low temperatures, no apparent magnetic orders are observed even at a prominently low temperature =0.03 K that is four orders of magnitude lower than the antiferromagnetic spin-exchange interaction 250 K.

The emergence of the quantum spin liquid near the Mott transition on two types of 2D Hubbard models with geometrical frustration effects has already been predicted in earlier numerical studies with the help of essentially exact algorithm of the path integral renormalization group (PIRG). The ground state does not seem to break symmetries so far proposed in literatures because of the geometrical frustration, while the full understanding of the spin liquid phase remains a challenge. Most of later numerical and theoretical studies have also been performed for a simplified single-band 2D Hubbard model based on an empirical estimate of parameters by following extended Hückel calculations. We certainly need a more realistic and ab initio description of -ET compounds to establish the existence of such truly as-yet-unestablished states and to get insights into possible fundamentally new ideas.

Another seminal finding achieved in this family is the unconventional character of the Mott transition found for =Cu[N(CN)]Cl (abbreviated as -Cl) under pressure. The novel universality class observed by the resistivity at this Mott transition is in good agreement with the marginal quantum criticality at the meeting point of the symmetry breaking and topological change. Its significance to physics and theoretical concept of the quantum criticality calls further experimental test and critical examinations based on the comparison with the realistic and first-principles grounds.

In spite of these innovative ideas and findings, the ab initio studies are so far few and most of the studies were performed using empirical models inferred from the Hückel studies. In addition, the strong electron correlation in the organic conductors hardly justifies its naive applications of the standard density functional theory (DFT).

To overcome the limitation and deficiency of the conventional DFT approach, a hybrid first-principles framework, which we call Multi-scale Ab initio scheme for Correlated Electrons (MACE) has been developed and applied to a number of strongly correlated materials. Our general framework consists of (1) ab initio calculations of the global electronic band structure either by the density functional theory or by other many-body theory such as the GW method and (2) a subsequent downfolding procedure by elimination of degrees of freedom far away from the Fermi level, which generates low-energy effective models. It is followed by (3) the procedure to solve the low-energy models by more reliable low-energy solvers such as PIRG, dynamical mean-field theory(DMFT), and the many-variable variational Monte Carlo (mVMC) methods. The accuracy of this three-stage scheme has critically been tested widely against various cases of exciton excitations in semiconductors, phase diagrams with competing orders in transition metal compounds, and iron superconductors. A specific combination of the local density approximation (LDA) for the first step and DMFT for the final step has also been widely applied. They have been favorably compared with available experimental results.

However, the accuracy is not clear for more complex compounds and for more strongly correlated systems. Indeed applications to the organic conductors remain a grand challenge. The ab initio low-energy models have recently been derived for the -ET compounds along the line of the above first and second parts of the three-stage scheme in the present terminology of MACE. In fact, they have shown that the ratio of the typical correlation strength (local screened Coulomb interaction), , to the typical electron transfer is as large as 10 for the downfolded single band models. Furthermore, we have found that the ratio of the typical off-site Coulomb interaction to is as large as . Elucidating electronic structures under such a strong correlation is a challenge at the forefront of research to develop efficient and accurate numerical algorithms.

The purpose of this study is to make a further step to the third procedure following the spirit of MACE; we solve the ab initio effective Hamiltonian of real -ET compounds by using an accurate solver based on a recently developed and improved mVMC method with many variational parameters to see whether the present general framework combined with the mVMC solver offers an accurate framework for the complex and strongly correlated organic conductors. For this purpose, we study -NCS, as a typical compound close to the Mott transition.

This compound is in fact barely metallic with an enhancement of the antiferromagnetic correlations revealed by the nuclear relaxation rate and located close to the antiferromagnetic Mott insulating phase. In fact, above 90 K, the resistivity shows insulating-like increase with decreasing temperatures and has an inflection point around =55 K (ref. ?) coinciding with the peak of with a crossover to a metallic behavior below it. It eventually becomes superconducting below around 10 K.

In this study, we assume normal states (not the superconductor) for the candidate of the metal, while leave it arbitrary for the insulator to study competitions between metals, antiferromagnetic or charge ordered insulators as well as the Mott insulator without symmetry breakings. We examine the phase boundary between the metal and the Mott insulator in a parameter space by taking the relative electron correlation amplitude as a parameter beyond the ab initio value. In practice, we draw a phase diagram as a function of a parameter that monitors the relative strength of the effective interaction to the electron transfer by uniformly scaling the interaction strength, where the ab initio value is given by . It gives us an idea about the relative location of the real material to this phase boundary and hence the relevance of the Mott physics. As we will show in §3, off-site Coulomb and exchange interactions largely stabilize a paramagnetic metal in the region of strong correlation ().

The present work is, to our knowledge, the first ab initio attempt to estimate the Mott transition and its neighboring phases in organic conductors containing a large number of atoms beyond 100 with four complex ET molecules in a unit cell. Our results indicate that -NCS is indeed near the Mott transition within the accuracy of 20%. It also shows that antiferromagnetic insulating state exists even in the metallic phase near the Mott transition, as a metastable excited phase, whose energy is typically about 1 meV 10 K higher than that of the paramagnetic phase. This is also consistent with the above experimental results of the crossover between high-temperature insulating and low-temperature paramagnetic metallic phases. Although the results show overall agreement with the experiments, a closer look of the ground state obtained with the realistic ab initio parameter (i.e., ) becomes an antiferromagnetic Mott insulator in contrast to the metallic phase in the experimental indications. We then discuss possible origins of this discrepancy.

This paper is organized as follows: We briefly summarize in §2 the method of our calculation with ab initio downfolding of the low-energy model for -NCS and the basic framework of the mVMC method we employed. In §3 calculated results are presented. The summary and discussions are given in §4.

## 2 Method

### 2.1 Derivation of low-energy effective model

Here, we describe a derivation of low-energy effective models for the present system. The scheme is based on first principles calculations and an application of the first two stages of the three-stage scheme. The basis of the Hamiltonian is the Wannier function associated with antibonding states of the highest occupied molecular orbitals (HOMOs) of two ET molecules that form a dimer. In the present paper, we restrict our consideration within the single-band models, where we derive the model containing only the degrees of freedom for the antibonding band while the bonding band is traced out in the downfolding. Possible dynamical effects of the bonding degrees of freedom remain as future issues. The explicit form of this Hamiltonian is given in the form of the two-dimensional (2D) single-band extended Hubbard model as

(\theequation) | |||||

where () is a creation (annihilation) operator of an electron with spin in the Wannier orbital localized at the th BEDT-TTF dimer. The parameters are given by

(\theequation) |

with and being the Kohn-Sham Hamiltonian representing an effective one-body potential. The and parameters are screened Coulomb and exchange integrals in the Wannier-orbital basis, respectively, expressed as

(\theequation) | |||||

and

(\theequation) | |||||

with being a 2D screened Coulomb interaction in the low-frequency limit. Note that the Hamiltonian given in eq. (\theequation) can be rewritten as

(\theequation) | |||||

where the onsite Hubbard parameter is given by , and . The spin operator is defined as , where and .

The derivation of for the purely 2D system follows ref. ?, where a new framework of the constrained random-phase approximation (cRPA) was developed for the purpose to derive effective interactions of models defined in lower spatial dimensions. This new scheme is suitable for quasi-low-dimensional materials such as the present system. The cRPA method is originally formulated in the RPA framework with the constraint for the band degree of freedom to eliminate only the degrees of freedom far from the Fermi level in energy. This is called the band downfolding. In the proposed scheme of the supplementary downfolding, however, the concept of the constraint is additionally relaxed to include the screening by the polarization in the other layers/chains even within the target bands. This is formulated in the real space representation and eliminates the degrees of freedom away from the target layer/chain, which results in the low-dimensional model for the target layer/chain. We call it the dimensional downfolding.

Practically, the band+dimensional downfolding is performed in two steps: We first perform the band downfolding to derive the 3D model for small number of bands near the Fermi level. This is followed by the dimensional downfolding in the second step. With this idea, we can naturally derive the low-energy model in any dimensions. In the present case, we use it for the derivation of a 2D model for -SCN.

### 2.2 Multi-variable variational Monte Carlo method

To investigate ground-state properties of the low-energy 2D model, we employ a multi-variable variational Monte Carlo method (mVMC) combined with quantum-number projection and multi-variable optimization. The variational wave function is defined as

(\theequation) |

where , , and are the Gutzwiller factor, the doublon-holon correlation factor, and the Jastrow factor, respectively. These factors are defined as

(\theequation) | |||||

(\theequation) | |||||

(\theequation) |

where , , and are variational parameters. Here, is a many-body operator which is diagonal in the real-space representations. When a doublon (holon) exists at the -th site and holons (doublons) surround at the -th nearest neighbor, gives . Otherwise, gives . In the present study, we take , where is the relative displacement between the sites and . The spin quantum-number projection operator restores the spin-rotational symmetry with the total spin .

For the one-body part , we employ a generalized pairing wave function defined as

(\theequation) |

where and are variational parameters and the total number of the sites, respectively. The number of electrons is denoted by . Throughout this paper, we consider the half-filled case . The paring wave function given in eq. (\theequation) can flexibly describe paramagnetic metals/insulators, antiferromagnetic insulators/metals, charge-ordered insulators/metals, superconducting phases as well as phases with strong quantum fluctuations showing developed spin and/or charge correlations decaying with any power laws as functions of distance. In fact, the mVMC method based on the variational function in eq. (6) can describe paramagnetic metals with developed spin correlations in hole-doped Hubbard models on square lattices.

In principle, it is better to allow variational wavefunctions as much as flexible without imposing any constraint on . However, if we do not impose the constraint, the possible variational parameters increase in proportion to , which increases the computation time for the optimization of the variational parameters enormously for large system sizes. To reduce the computational cost, we restrict to have a sublattice structure as , where is a sublattice index at site . In the study of the Mott transition in §LABEL:Mott_transition, we take (see Fig. LABEL:fig:sblattice-2x2). This assumption is based on the observation of the staggered magnetization with ordering vector of in the frustrated Hubbard model at , and the spin- quantum Heisenberg antiferromagnets corresponding to for . These appear to be relevant to the present parameter region as we will see in the comparison with the exact diagonalization as well. For the study of the charge order in §LABEL:Charge-ordered_insulator (as we describe in detail there), we take more general and to allow the three-sublattice as well as two-sublattice structures (See Fig. LABEL:fig:sblattice-6x2).