Detecting quantum critical points in the $t-t'$ Fermi-Hubbard model via complex network theory

TL;DR

This paper introduces a novel approach using complex network theory to detect quantum critical points in the $t-t'$ Hubbard model, linking entanglement patterns to phase transitions in high-temperature superconductors.

Contribution

It demonstrates that analyzing quantum mutual information graphs as complex networks can identify quantum critical points even on small lattices, overcoming finite size limitations.

Findings

Quantum mutual information graphs reveal critical behavior.

The method detects QCPs on small lattices.

Complex network characteristics change near the transition.

Abstract

A considerable success in phenomenological description of high-T$_{\rm c}$ superconductors has been achieved within the paradigm of Quantum Critical Point (QCP) - a parental state of a variety of exotic phases that is characterized by dense entanglement and absence of well-defined quasiparticles. However, the microscopic origin of the critical regime in real materials remains an open question. On the other hand, there is a popular view that a single-band $t-t'$ Hubbard model is the minimal model to catch the main relevant physics of superconducting compounds. Here, we suggest that emergence of the QCP is tightly connected with entanglement in real space and identify its location on the phase diagram of the hole-doped $t-t'$ Hubbard model. To detect the QCP we study a weighted graph of inter-site quantum mutual information within a four-by-four plaquette that is solved by exact diagonalization. We demonstrate that some quantitative characteristics of such a graph, viewed as a complex network, exhibit peculiar behavior around a certain submanifold in the parametric space of the model. This method allows us to overcome difficulties caused by finite size effects and to identify the transition point even on a small lattice, where long-range asymptotics of correlation functions cannot be accessed.

Figures

Equations

$$H=-t\sum _{⟨i,j⟩,\sigma }{c}_{i,\sigma }^{†}{c}_{j,\sigma }-t^{\prime }\sum _{⟨⟨l,k⟩⟩,\sigma }{c}_{l,\sigma }^{†}{c}_{k,\sigma }+h.c.+U\sum _i{n}_{i,\uparrow }{n}_{i,\downarrow },$$

$${\mathcal{I}}_{AB}=S_A+S_B-S_{A\cup B},$$

$$C=\frac{\mathrm{Tr}{\mathcal{I}}^3}{{\sum }_{j\ne i}^N{\sum }_{i=1}^N{\left[{\mathcal{I}}^2\right]}_{ij}},$$

$$Y_i=\frac{{\sum }_{j=1}^N{\left({\mathcal{I}}_{ij}\right)}^2}{{\left({\sum }_{j=1}^N{\mathcal{I}}_{ij}\right)}^2}$$

$$D=\frac{1}{N\left(N-1\right)}\sum _{i,j=1}^N{\mathcal{I}}_{ij},$$

$$D\le \frac{2}{N\left(N-1\right)}\sum _{i=1}^N{S}_i\le \frac{2\ln 4}{N-1}\underset{N\to \infty }{\to }0$$

$$r_{ij}=\frac{{\sum }_{k=1}^N\left({\mathcal{I}}_{ik}-⟨{\mathcal{I}}_i⟩\right)\left({\mathcal{I}}_{jk}-⟨{\mathcal{I}}_j⟩\right)}{\sqrt[]{{\sum }_{k=1}^N{\left({\mathcal{I}}_{ik}-⟨{\mathcal{I}}_i⟩\right)}^2}\sqrt[]{{\sum }_{k=1}^N{\left({\mathcal{I}}_{ik}-⟨{\mathcal{I}}_i⟩\right)}^2}},$$

$$⟨{\mathcal{I}}_i⟩=\frac{1}{N}\sum _{j=1}^N{\mathcal{I}}_{ij}$$

$${\rho }_A(a,a^{\prime })=\frac{1}{Z}\sum _n{e}^{-\beta E_n}{\mathrm{Tr}}_{\bar{A}}|{\psi }_{n,(a,\bar{a})}⟩⟨{\psi }_{n,(a^{\prime },\bar{a})}|.$$

$$G_{i,\sigma }(\omega )=\frac{1}{Z}\sum _{m,n}\frac{|⟨m|c_{i,\sigma }^{†}|n⟩{|}^2}{\omega +E_n-E_m}\left(e^{-\beta E_n}+e^{-\beta E_m}\right).$$

$$A(\omega )=-\frac{1}{\pi }\text{mbox}Im{\mathrm{Tr}}_{i,\sigma }{G}_{i,\sigma }(\omega ).$$

Full text

Introduction

The phenomenon of high-temperature superconductivity (HTSC) still remains very puzzling after more than thirty years since the discovery of superconducting copper-oxide compounds Bednorz . Serious hopes for the understanding of this phenomenon are related to the concept of a quantum critical point (QCP) Sachdev_book ; Sachdev_new , - an exotic state of matter that exhibits scale invariance and lacks long-lived quasiparticles, and thus cannot be described by means of conventional Fermi-liquid theory. Contemporary discussions of observed properties of HTSC are frequently organized around this concept Hussey ; Michon . Precise nature of this critical point is still unclear, – different studies relate it to charge density waves Arpaia , nematic Auvray , or antiferromagnetic fluctuations Wang . However, its theoretical treatment can be conducted universally, though it requires a change of basic mathematical tools: the diagrammatic approach, the main apparatus of quantum many-body theory during the last sixty years AGD ; Mahan , is not really fitted to the description of systems lacking quasiparticles. A paradigmatic shift in studying strongly coupled systems near the QCP has occurred when it was realized that the anti de Sitter/Conformal field theory (holographic) correspondence Maldacena can be used to analyze certain universal phenomenological properties of correlated electronic matter in the regime where the traditional Fermi-liquid picture breaks down Zaanen ; Hartnoll_book . While the number of direct experimental evidences of quantum critical points in high-Tc superconducting materials is limited, assuming its existence and employing the methods of holography allowed to resolve within a relatively short time frame a number of puzzles that remained perplexing for decades Anderson . The correspondence provided an explanation for the linear-$T$ scaling of DC resistivity in the normal state of cuprates Legros (known as strange metals), relating it to general hydrodynamic properties of systems with minimal viscosity proportional to the thermodynamic entropy Davison . It was shown Blake that the Hall angle, - the temperature dependent ratio of the Hall and DC conductivities, $\tanh{\theta_{H}}=\sigma_{xy}/\sigma_{xx}\sim 1/T^{2}$, can be naturally interpreted in terms of a two-constituent quantum liquid, where the regular quasiparticles and the critical sectors give independent contributions to the conductivity, leading to an anti-Matthiessen rule for transport. A new mechanism of the interaction-driven metal-insulator transition that causes anisotropic localization has been suggested Hartnoll , and it appears to be fully in line with the localization of conducting electron gas in two-dimensional CuO planes, while the conductivity in the orthogonal direction is suppressed. Other phenomena, such as the formation of Fermi arcs seen in the angle-resolved photoemission spectra of high-Tc compounds Iliasov , or charge density waves Krikun also fit pretty naturally into the context of quantum criticality.

The main problem of this approach is its purely phenomenological character. It cannot explain by itself why the high-Tc compounds, contrary to the most of interesting condensed matter systems, do not behave as the Fermi liquid but instead are characterized by minimal quantum viscosity and other fancy properties. Such an explanation requires an analysis of electronic structure of specific materials.

Current understanding of high-Tc superconductivity in cuprates assumes a crucial role of strong electron correlations Anderson ; Dagotto ; Millis ; Scalapino , which are taken into account within a particular minimal model that was formulated Andersen1 on the basis of the density functional band structure of cuprates, - the single-band $t-t^{\prime}$ Hubbard model on a square lattice given by the Hamiltonian

[TABLE]

where, the first sum is taken over the pairs $\langle i,j\rangle$ of nearest neighbors, the second one - over the pairs $\langle\langle l,k\rangle\rangle$ of next-to-nearest (diagonal) neighbors, $c_{i,\sigma}$ is the electron annihilation operator, and the on-site occupation operator is $n_{i,\sigma}=c^{\dagger}_{i,\sigma}c_{i,\sigma}$. For convenience, in what follows we express all energies in units of $t$ .

In an attempt to connect the phenomenological and the microscopic levels of description of HTSC, we shall focus on the $t-t^{\prime}$ Fermi-Hubbard model and try to detect QCP on its phase diagram.

Correlation effects beyond the band structure approximation in this model have been thoroughly analyzed with different methods Maier ; Haule ; Khatami ; Gull ; Civelli , and there are a number of good indications that it captures all the relevant features of cuprate superconductors. In a series of papers IKK_parquet ; Metzner1 ; Metzner2 ; Metzner3 ; Metzner4 ; Honerkamp , perturbative renormalization group studies of the model have been conducted, and the emergence of the superconducting order parameter and the competition between superconductivity and antiferromagnetism were demonstrated. In particular Metzner4 , it was argued that the next-to-nearest neighbor hopping $t^{\prime}$ plays a crucial role in the stabilization of superconductivity. A complementary approach is based on the cluster dynamical mean-field studies which consider a 2-by-2 plaquette as an elementary unit LK2000 . Recently Harland1 , it was noticed that this plaquette has a very special electronic structure for the parameters and the electron occupation number typical for the the optimal doping regime in YBa2Cu3O7 ($t^{\prime}=-0.3$, $U\simeq 6$), with an “accidental” degeneracy of many-electron energy levels and formation of the soft fermion mode due to this degeneracy. The pseudogap forms via this mode by a mechanism of the Fano antiresonance, and the superconducting d-wave susceptibility dominates over other instability channels. This behavior was interpreted in terms of formation of a local plaquette valence bond state. On a larger scale, the ground state of the model has been analyzed by means of density matrix renormalization group Devereaux1 (DMRG) (see also Kivelson for the related studies of its cousin, $t-J$-model), and additional arguments in favor of stabilization of superconductivity by the next-to-nearest neighbor hopping were provided. In turn, at temperatures above the superconducting phase transition, determinantal Monte Carlo computations Devereaux2 demonstrated that the DC resistivity exceeds the Mott-Ioffe-Regel limit and scales linearly with temperature.

The search for the QCP in the $t-t^{\prime}$ Hubbard model has been performed within the dynamical cluster approximation Khatami , and its existence has been proven by studying thermodynamics properties of the model at finite temperature and their further extrapolation to $T=0$. However, it is tempting to get a deeper insight into the microscopics of the QCP and demonstrate its emergence due to interactions of electrons at low temperatures.

Since large scale simulations of the fermionic Hubbard model away from half-filling are challenging because of the sign problem, it is natural to ask whether we can extract any information about the tendency to form critical states out of small cluster solutions obtained by means of exact diagonalization. At first, this goal does not seem realistic since studying systems in the critical regime unavoidably requires dealing with long-range correlations, while all the microscopic precursors of the transition on small lattices would be washed out by the finite-size effects. However, it is useful to bear in mind that, in the context of many-body quantum dynamics, the concept of entanglement and the phenomenon of collective emergence go hand in hand. An archetypical example of such relation is the Cooper pairs in the BCS theory of superconductivity: while the ground state wavefunction has a form of a product state of the Cooper pairs, each pair itself is a two-body entangled system. Therefore it is natural to expect that major transitions in phenomenological properties of many-body systems would be reflected in the patterns of entanglement, and quantum criticality should leave its fingerprint on all scales, not only in the deep infrared limit. A nice example of how fruitful this way of thinking can be was given in Refs. Tremblay1 ; Tremblay2 , where entanglement measures were used to determine universality class of the Mott transition in the 2d Hubbard model.

Recently, a novel approach to phase transitions in quantum lattice models based on complex network theory has been suggested Carr1 ; Carr2 . It was noticed that a particular structure that can be computed with relative ease and appears to be very sensitive to reconfigurations of the quantum state is the network of quantum mutual information.

The mutual information between two subsystems $A$ and $B$ of a larger systems is defined as

[TABLE]

where $S_{A}=-\operatorname{Tr}{\rho_{A}\log\rho_{A}}$ is the von Neumann entropy, and $\rho_{A}=\operatorname{Tr}_{\bar{A}}\rho$ is the density matrix of subsystem $A$. Then we can associate a weighted graph with a state of a quantum lattice system, e.g. the Hubbard model, by considering the lattice sites $i=1\dots N$, where $N$ is the number of sites, as nodes of the graph, and the values of pairwise inter-site mutual information $I_{ij}$ play the role of weights on the graph links (see Fig. 1). This representation is appealing for the following reason. Once a wave function on the lattice is known, it is easy to compute the entanglement entropy of a pair of sites and thus the mutual information. At the same time, such a network contains information of quantum correlations which could be very important to understand the dynamics of strongly correlated systems. In the cases of the transverse field Ising and the Bose-Hubbard models in 1d, it was demonstrated that certain characteristics of the mutual information network can be used to detect quantum phase transitions Carr1 ; Carr2 . Namely, behavior of the following functions upon changing parameters of the models has been studied:

• •

Clustering of a weighted graph is defined as

[TABLE]

where $N$ is the total number of sites in the lattice, and $\cal I$ is the $N\times N$ matrix of inter-site mutual information. One can see that this quantity maximizes on graphs with a lot of three-link loops with high weights. For the cases studied in Ref. Carr1 , it was shown that it serves as sensitive detector that exhibits a clear dip at the phase transition point. A natural explanation of this fact is that, at the criticality, one can expect the corresponding network to be scale-free, and for generic scale-free networks clustering is usually quite low scalefreeclustering .

• •

Disparity of a single node in a network is defined as a measure to capture how non-uniformly weights on the links attached to this node are distributed:

[TABLE]

For example, if the node has the same value of mutual information with all the other nodes of the network, its disparity would be $Y_{i}=1/(N-1)$, while if it correlates only with one neighbor, the disparity maximizes as $Y_{i}=1$. Physically speaking, high disparity of a lattice site means that it tends to correlate only with a few other sites, and “factorize out” of the rest of the system. In the context of quantum many-body physics such a behavior would be typical for states that can be nearly decomposed into product states. On the other hand, low disparity means that the site correlates with a large number of degrees of freedom.

• •

Density is an overall characteristic of a network given by

[TABLE]

i.e. it is the averaged fraction of all the weights (mutual information values) of the network. To gain more intuition on what properties of the many-body quantum state it reflects, we shall estimate an upper bound on this measure. If site $i$ of the network is maximally entangled with the rest of the system, its entanglement entropy equals $S_{i}=\ln d=\ln 4$, where $d=4$ is dimension of the local on-site Hilbert space in Hubbard model. On the other hand, mutual information monogamy theorem implies that $2S_{i}\geq\sum\limits_{i,j}{\cal I}_{ij}$ monogamy , leading to

[TABLE]

i.e. the mutual information network is generally sparse even if the system is highly entangled. Note that bound (6) can be saturated in physically very distinct cases. $D$ is maximal if either each single site is maximally entangled with just one partner site, and the state as a whole decomposes into a product of Bell pairs, or if the entanglement between the site and the rest of the system is homogeneously scrambled over all the sites. To distinguish between such configurations one has to refer to the disparity which we defined above.

• •

Pearson correlations measure how much two nodes $i$ and $j$ of a network differ from each other:

[TABLE]

In Ref. Carr1 Pearson correlations of neighboring nodes were shown to develop a cusp around the phase transition point.

For one-dimensional Ising and Bose-Hubbard models Carr1 , this approach to detecting quantum phase transitions points was successfully applied for systems of $\sim 10^{2}$ sites, and was demonstrated to be very robust upon finite-size effects. In the two-dimensional case, we are limited by much smaller system sizes (we perform exact diagonalization for a 4-by-4 plaquette), and should not expect our results to be free from finite-size artifacts. Still, as we shall see in the next section, the network measures exhibit clearly distinguishable features on a submanifold of the $t-t^{\prime}$ Hubbard model parametric space. In particular, this submanifold includes the level-crossing point observed in a 2-by-2 plaquette for the choice of parameters corresponding to YBa2Cu3O7 superconductor Harland1 .

Results

We have computed the complex network measures discussed above across the space of parameters of a 4-by-4 $t-t^{\prime}$ Hubbard plaquette. Within each fixed particle number sector, from $(5,5)$ ($37,5\%$ hole doping) to $(7,7)$ ($12,5\%$ hole doping), we scan over $t^{\prime}$ and $U$, Fig. 2. As an indicative value, we take $t^{\prime}=-0.3$, which is estimated to be the next-neighbor hopping in the Hubbard model of YBaCuO compounds, and search for transition points around it. The temperature is fixed to $1/T=\beta=100$ (all energies are expressed in the units of $|t|$), and the system is studied in the canonical ensemble.

We assume that a transition point is evident if all the measures exhibit some clear features around the same point. Accepting this criterion, we can claim with a high confidence that, for non-periodic boundary-conditions, there is a family of transition points in each sector (except $(7,7)$) forming a nearly perfect straight line in the $t^{\prime}-U$ plane that extends in a certain range of $t^{\prime}$ 111For a more detailed picture of how complexity measures behave at different values of $t^{\prime}$, see Ref. Supplemental , Fig. 3 (for too small $|t^{\prime}|$ the signs of criticality are faded away from the complexity measures). Moreover, for different values of hole doping, all these lines have very similar slope. This can be interpreted as that phase transition occurs on a 2d manifold in the 3d parametric ($U/t$, $t^{\prime}/t$, particle number) space of the model.

We stick to non-periodic boundary conditions for the reason that the mutual information network has a richer structure in that case. If periodic boundary conditions are imposed, all lattice sites are identical, and every site has only five inequivalent connections to others, making the ${\cal I}_{ij}$ matrix highly degenerate. Hence, the corresponding network structures are constrained by symmetries and much less sensitive to variations of the model parameters. Still, we would like to note that, when boundary conditions are changed for periodic ones, all the phase transition lines are smeared out with the only exception of the $(6,6)$ sector which corresponds to the hole doping of $\delta=25\%$. For the latter, only the concrete values of Coulomb repulsion $U$ gets shifted Supplemental . While we do not expect the information network constructed for periodic boundary to be sensitive enough to properly detect phase transitions, it is interesting to note that this single sector where the transition is evident for both choices of b.c. is the same as the one where level-crossing associated with formation of the pseudogap via Fano antiresonance occurs in a 2-by-2 plaquette Harland1 .

At the same time, in the density of states (d.o.s.) the transition point is (almost) invisible. Some minor peculiarity at the quantum critical point is visible in the density of states at $t^{\prime}=-0.3$. Around the transition point identified by means of the complex network theory ($U=7.5$, sector (6,6)) the $\omega=0$ peak in the d.o.s. starts splitting and the pseudogap emerges, see Fig. 4. Further decrease of the hole doping leads to enhancement of the gap. The particular role of $U$ in this transition is less clear, as the d.o.s. profile varies very mildly upon changing $U$, and it would be safer to claim that the spectral properties are not sensitive to the discussed quantum phase transition. Ideologically, this situation is somewhat similar to the Anderson localization in disordered systems which is a clear example of a phenomenon that cannot be detected on the level of the average Green’s functions Disordered .

Discussion.

By associating the quantum state of the $t-t^{\prime}$ Hubbard model with a weighted network of inter-site mutual information, for different values of the next-neighbor hopping $t^{\prime}$, we have found a set of transition lines in the $U-t^{\prime}$ plane of the model parametric space, where characteristics of the network have a clearly distinguishable cusp. Such a behavior was previously shown to be an indication of quantum phase transitions in different one-dimensional models Carr1 ; Carr2 . The modern experimental understanding of the putative QCP in cuprates tells that it must be associated with the emergence of the pseudogap phase Hussey . For example, for YBaCuO compounds the onset of pseudogap was experimentally demonstrated to happen at hole doping $\delta\simeq 22\%$ Sato . The hole doping $\delta=25\%$ is the closest value one can get for a 4-by-4 cluster (the $(6,6)$ sector), and, interestingly, it is exactly the sector where the complex network measures demonstrate the most robust transition features. The particular values of the on-site Coulomb repulsion is affected by the finite size effects, and estimated to be about $U\simeq 7\,-\,8$. At the same time, no peculiarity is seen in the density of states at the transition point, apart from a slight splitting of d.o.s. around the Fermi level, which might be a good indication that the low-order correlation functions that define the spectral and the response properties of the system could be blind to restructuring of many-body quantum states, and does not contain enough information on the role of quantum correlations behind phase transitions in electron systems.

Both the strength and the weakness of the employed approach is that it helps to identify any critical point while being ignorant about its nature. Therefore we can claim the existence of a manifold of QCP in the $t-t^{\prime}$ Hubbard model with a high confidence, but we cannot deduce what order parameter of the corresponding transition is. Still, we tend to relate the observed transition to the critical point discussed in Harland1 , where it was associated with emergence of soft fermion modes.

Within the exact diagonalization approach, we were able to consider only a small cluster, where one unavoidably has to deal with strong finite size effects. It would be interesting to conduct a similar analysis for a larger system. The state-of-the-art DMRG allows to study quasi-2d stripes as large 4-by-64 atomic sites Devereaux1 ; Kivelson ; Devereaux2 , and we hope to apply the complex network approach to such systems in the future.

Methods.

In this section we give the relevant technical details of the calculation of the entanglement measures defined above. The first step is to diagonalize the Hubbard model (1) for a 4-by-4 cluster. This can be done either for a periodic or a non-periodic model. The diagonalization is performed using the Lanczos algorithm with 200 Krylov basis vectors Lanczos . The particle number and spin conservation laws are used so that the diagonalization can be restricted to a sector with a fixed number of up- and down-spins. Those eigenstates with the corresponding eigenvectors are then used to calculate the reduced density matrices for each possible pair of sites as well as for each single site.

The reduced density matrix is computed using its definition that can be symbolically written as:

[TABLE]

Here $a,a^{\prime}$ denote the many-particle (Fock) basis states describing the subsystem $A$ we calculate the density matrix for, $\bar{a}$ stands for the many-particle basis state of the complementary subsystem $\bar{A}$, thus a couple of those $(a,\bar{a})$ denotes a basis Fock state for the whole cluster explicitly split into two parts. As before, $n$ stands for a particular eigenvector, the density matrices for given eigenstates are weighted with the Boltzmann factors corresponding to their energies. In a given sector for a given set of parameters we use the Boltzmann factor cut-off of 1% meaning $e^{(E_{0}-E_{k})\beta}>10^{-2}$, where $E_{0}$ is the ground state energy and $E_{k}$ is the energy of the highest ($k$th) level taken into account. Note that while performing the partial trace over $\bar{A}$ one has to correctly account for the fermionic commutation relations. To this aim one has to effectively change the numeration of sites so that the sites for which we calculate the density matrix stand first. Explicitly it means that each component of an eigenvector, corresponding to a given basis state of the cluster, gets a factor determined as the parity acquired while ”dragging” the occupied sites of $A$ to the beginning past the occupied states of $\bar{A}$. In other words for each basis vector one takes each occupied site from $A$ and for each occupied spin component counts the number of same spin occupied sites from $\bar{A}$ standing before the considered site in the original numeration. Summing up the parities of those numbers for all occupied sites and spins from $A$ one gets the parity that is assigned to a given basis vector with respect to the subsystem $A$. Having multiplied the eigenvector components with the acquired parities one finally performs the partial trace over the complementary subset $\bar{A}$.

Given the reduced density matrix, we first calculate the von Neumann entropy of a given subsystem and then, with Eq. (2), the mutual information for each pair of sites, that serves as the basis for our network.

The $\omega$-dependent Green function is given by:

[TABLE]

Here $m,n$ denote the eigenstates of the system, $i$ and $\sigma$ denote a given site and spin (in the paramagnetic case the answer is spin-independent), $E_{n}$ is the energy of the $n$-th state, and $Z=\sum_{m}e^{-\beta E_{m}}$ is the partition function. Note that $m$ and $n$ necessarily belong to different sectors.

The Green function is then used to calculate the spectral density of states, which is defined as

[TABLE]

To perform numerical computations, delta-peaks in the Green function are broadened with $\delta=\pi/\beta$. Acknowledgements – the authors thank Lincoln Carr for inspiring discussions.

Competing Interests – the Authors declare no Competing Financial or Non-Financial Interests.

Author contributions – A.A.B., M.I.K. and A.I.L. designed the project and directed it with the help of S.B., M.D. and M.H. performed the calculations. A.A.B., M.I.K. and S.B. wrote the manuscript. All authors contributed to discussions.

Funding – the work of A.A.B. was supported by Russian Science Foundation, Grant No. 18-12-00185. M.I.K. and A.I.L. acknowledge a support by European Research Council via Synergy Grant 854843 - FASTCORR. M.H., S.B. and A.I.L. acknowledge support by the Cluster of Excellence ’Advanced Imaging of Matter’ of the Deutsche Forschungsgemeinschaft (DFG) - EXC 2056 - project ID 390715994.

Data Availability – the data that support the findings of this study are available from the corresponding author upon reasonable request.