A picture of pseudogap phase related to charge fluxes

Since the discovery of high temperature superconductivity, with increasingly more understandings in progress, more complexities also emerge, making the field realize the necessity of new physical concepts and pictures, in order to understand the underlying mechanism. Based on our previous understanding of the materials-dependent properties of cuprate superconductors, charge fluxes along different crystallographic directions were considered to be important. In this paper, a new emergent configuration in cuprates is first introduced, called the charge superplane, which exhibits the property of conserving the dynamic charge flux flows. These charge superplanes further redefine the fundamental collective excitation in cuprates, called pQon, that exhibits the novel property of polaron-Quasiparticle duality. Our angle resolved photoemission spectroscopy simulations based on the new picture provide a good agreement with previous ARPES experiments, evolving from underdoped pseudogap region to overdoped regions. Our percolation simulations further provide an explanation for the recently reported temperature independent critical doping around 19%, and a connection to previous scanning tunneling imaging results. Our work provides a new perspective to further understand the origin of pseudogap phase and other related phases in cuprates.


Main
Understanding high temperature superconductivity has caused a long but exciting debate in the community of condensed matter physics since the discovery of La2CuO4 [1]. Earlier works identified the CuO2 plane as the superconducting plane, where strong electronic correlations and spin fluctuations exist [2][3][4]. Although all cuprate superconductors have the same CuO2 building blocks, previous works have pointed out the importance of structural and electronic features out of the CuO2 plane in modulating the materials dependence of superconducting transition temperature Tc [5][6][7][8][9][10][11]. However, recent reports of light-induced superconducting-like behavior near room temperature in YBa2Cu3O6+x (YBCO) [12][13][14] and femtosecond electron-phonon lock-in in iron selenide [15] further suggest the importance of ultrafast lattice-dynamical processes in unconventional superconductors, which most previous models neglect.
In our recent model [16], such ultrafast lattice effect was explicitly considered to couple with electronic effects. Specifically, the higher apical force generated by perturbing the apical oxygen was for the first time correlated with Tc,max. The apical force was naturally connected with the restoring force for lattice oscillations of apical oxygen mode that was found to also couple with the in-plane breathing mode in our density functional theory (DFT) simulations. It suggests that the apical oxygen and the four in-plane oxygens around each Cu ion prefer to oscillate inward toward the Cu ion in-phase, and hence also oscillate outward away from the Cu ion together. Such coupling is accompanied with a strong local charge flux, i.e., a transient directional flow of electron charge. When oxygen ions and charge fluxes around two neighboring Cu ions oscillate cooperatively, DFT simulations based on hybrid functional shows that a small polaron can hop from one Cu site to the neighboring one with vanishingly small diffusion barrier, which is a significantly reduced value compared with the bare polaron hopping barrier without the help of such cooperative fluxes. Note that, first, the apical fluxes here are very small in magnitude compared with the charge of each hole that moves in-plane in this simulation; second, stronger apical force is associated with stronger charge flux among different cuprate families. Our previous work [14] thus suggests a few important pieces of pictures: the cooperative charge fluxes, their coupling with lattice oscillations, their modulation to in-plane transport properties and their correlation with Tc,max in different cuprate families.
One central puzzle of high-temperature superconductivity is the origin of pseudogap phase [17]. For a given underdoped cuprate, the state above Tc but below a characteristic temperature range of T*, which can be determined from many experimental techniques, features a spectral gap of unknown nature, labelled as the pseudogap state or pseudogap phase. Its classification as a true thermodynamic phase, associated with broken long-range symmetries, however, is still a subject of debate [18][19][20][21][22][23][24][25][26][27]. Recent photodestruction experiment suggests, on an ultrafast timescale, an absence of collective electronic order and instead, an existence of ultrafast carrier localization into polaronic states [28]. On the other hand, in angle resolved photoemission spectroscopy (ARPES) measurements, pseudogap phase is known to exhibit a more quasiparticle-like behavior around the nodal direction [29], as with doping relatively sharp peaks appear near the Fermi level, giving the Fermi arc; while simultaneously around the antinodal direction such sharp peaks are completely suppressed [30], giving the pseudogap. Meanwhile, along both nodal and antinodal directions a broad distribution of spectrum weight was observed to extend far below the Fermi level, leading to the discussion of polaron-like behavior [31,32]. There were thus debates about whether the hole in cuprates should be considered as quasiparticle or polaron.
In this paper, we first introduce a new concept of "charge superplane" that governs the collective flows of charge fluxes in cuprates, and eventually serves as a critical emergent configuration to integrate the abovementioned many pieces into one picture. We then derive a natural but amazing result of such charge superplane from our Hamiltonian, in which polaron and quasiparticle can dynamically transform into each other. This correspondingly defines a new collective excitation with polaron-Quasiparticle duality, which is thus coined as pQon. Our simulation of ARPES based on the electronic structure of pQon agrees well with previous ARPES experiments. A percolation simulation considering the hopping of pQons at various doping levels also gives an interpretation of the temperature independent critical doping around 19% reported recently [33], and a possible connection to previous scanning tunneling imaging results on "inhomogeneity" [34][35][36][37]. Figure 1 shows the (100) and (010) superplanes that are perpendicular to the CuO2 planes. They cut through the apical and the in-plane oxygen ions, the copper ions and the apical atoms along either of the two dx 2 -y 2 orbital directions. The most intriguing property of these superplanes can be observed by the fact that they dynamically conserve the charge flow within certain region in cuprates by driving those phonon couplings that prefer such conservation. In order to quantify such charge conservation property, we extend the superplane from a geometric plane illustrated in Fig.1 to a physical layer with thickness of 1 Å. Electron densities within such layer are defined as being inside the superplane area, while electrons in other regions in the cuprate structure are considered as being outside the superplane area. Using this metric, we further calculate different levels of charge flux leakage outside the superplane induced by all pairwise couplings of the selected phonon modes in a supercell ( Fig. 2a-  Several important trends are observed in Fig. 2. First, phonon modes that are more strongly coupled, represented by regions with darker squares in Fig. 2a, c, tend to show stronger charge conservation ability, or less charge leakage, corresponding to regions with lighter squares in Fig. 2b, d. This reflects the charge conservation property of superplanes to drive the coupling of those phonon modes that can better conserve the charge flow within the superplane area. Second, such phenomena show certain materials dependence. YBCO, with lower Tc,max than Hg-1212, shows a weaker such correlation in general (Fig. 2c,d). Hg-1212 shows almost no charge leakage out of the superplanes for those strongly coupled modes, limited leakage for other modes and less probability to induce a charge leakage for a given mode. While YBCO shows obviously more charge leakage and more chance to leak. Thus, the ability of charge superplanes to conserve charge flux flow seems to be also positively correlated with the experimental Tc,max of different hole doped cuprate families. This is complementary to the apical force correlation [16], where stronger flux was associated with higher Tc,max. Third, in the bottom area of Fig. 2a, there are three frozen modes with obvious couplings to many perturbation modes labelled in the dashed box, corresponding to the apical oxygen mode, the breathing and anti-breathing oxygen modes, respectively, in Hg-1212. They show the strongest phonon coupling and the best charge conservation ability to lock flowing charge fluxes inside the superplane area. The trend is less obvious but still there in YBCO. Note that the coupling of apical and breathing modes was also argued to be crucial to the transport property of cuprates [16]. We further replace the apical atom Hg in Hg-1212 with Tl, Bi, Y, and Pb and calculate the phonon coupling of these hypothetical materials as we did in Fig. 2a. Fig. 2e shows that the total phonon coupling strength of Hg-1212 is higher than the hypothetical structures, suggesting that previous experimentalists might have identified one of the best apical atoms for high Tc. In summary, these observations suggest that stronger charge flux and better conservation of their collective flows by superplanes prefer higher Tc,max for a cuprate family.

Polaron-Quasiparticle Duality: pQon
Our vision is that if the charge superplane we identified in cuprates is crucial to the transport property as we discussed in the previous section, then it may also play a critical role here in "distinguishing" the hole representation by polaron or quasiparticle. Because as long as the hole moves, a charge flux is generated, which should be modulated by the charge superplane. Generally speaking, a hole moving along antinodal direction won't obviously violate such charge conservation property of the superplane, as at least the direction of charge movement associated with the direct hole hopping is inside the superplane (Fig. 1, Fig.  3a); while moving along the nodal direction it generates an immediate tendency for charge flux to leak outside the superplane. Considering previous ARPES results, it is then possible that the hole movement along antinodal direction is well protected by the charge superplane to behave somewhat like a polaron in transport, while such polaron-like hole, once hopping along the nodal direction, may show a transient quasiparticle-like behavior due to a lack of superplane protection. Such a picture is consistent with our DFT simulations ( Fig. 3b), where a polaron hopping along the nodal direction will tend to be converted toward a quasiparticle at the early stage of conversion due to the lack of energy barrier, while the antinodal hopping well protects the polaron entity by the immediate energy barrier. It is worth noting that in later stage of the conversion beyond 30% along the nodal direction, an energy barrier builds up to prevent the further conversion to quasiparticle. That is, along the nodal direction a polaron may only transiently show certain character of a quasiparticle. Note, however, a spectroscopy technique that targets the transient nodal conversion processes between 0% to 30%, such as different ultrafast spectroscopies, or measures their time average, such as ARPES, may convert such transient state to manifest a quasiparticle-like spectrum weight by photoexcitation. We now formularize this picture in cuprates and show that a new fundamental collective excitation, called pQon, can be defined in such systems. We use ",$ % and ",$ % to represent creation operators of spatially more delocalized and spatially more localized carrier states, respectively, and propose the following Hamiltonian: where " and 8 are bare kinetic energies of quasiparticle and polaron, respectively, i.e., kinetic energies when they propagate freely. For the quasiparticle, we employ a tight-binding dispersion relation on a square lattice: where , ′ and represent nearest-and next-nearest-neighbor hoppings, and chemical potential, respectively. Here we use the typical values = 0.4 eV and ′ = −0.08 eV [6], and leave to adjust the doping concentration. Considering that polaron is more localized than quasiparticle, it is reasonable to set 8 = −0.1 " to represent the much lower mobility of polaron. Interaction term KL represents the strength of charge transfer between the polaron and the quasiparticle states, which can be induced by coupling to the lattice and/or electronic correlations. We set KL as a d-wave type function to exhibit the d-wave directional preferences: where KL measures the maximum strength of the interaction. KL is zero along the diagonal direction; while it reaches its maximum along the axial direction. Diagonalization of eq (1) yields two renormalized bands: , For the lower renormalized band "$ S , the corresponding operator is: where the coefficients have the following expressions: The operator ",$ / = ",$ ",$ / + ",$ ",$ / represents a carrier with a dual nature of polaron and quasiparticle. We thus define ",$ / as a new collective excitation in this system, coined as pQon.
Propagating pQon will be in a dynamical transformation of polaron-like and quasiparticle-like states, the ratio of which is strongly momentum-dependent due to the d-wave form of KL . We interpret the transformation process as follows: when pQon propagates along the diagonal or nodal direction, it shows alternating transient states of polaron-like and quasiparticle-like behaviors, with coefficient ",$ much larger than ",$ ; while along the axial or antinodal direction, pQon interacts strongly with the superplane that also protects it well, making it behave more like a polaron throughout the propagation, with coefficient ",$ much larger than ",$ .

ARPES simulation based on pQon
Inspired by previous pioneering ARPES simulations [38,39], we further perform ARPES simulations based on the electronic structure of pQon (see Method). Along lines of ^= 1.00 , 0.81 , 0.62 and 0.50 in k-space, we calculate the ARPES spectra, as shown in Fig. 4, which agree well with the previous ARPES experiment [26]. It's clear to see that a gap opens, with the maximum in the antinodal direction, i.e., ^= 1.00 . When ^ moves away from 1.00 , the degeneracy is gradually lifted with the vanishing effect of KL and the recovery of Fermi surface finally. The antinodal gap formation can be attributed to the strong interaction between pQon and charge superplane in axial direction, so that strong spatial localization characteristics are observed in the antinodal region. This comparison with experiment also justifies our anisotropic KL defined in eq (3). We want to emphasize that our simulation has captured several unusual features in the spectra of the pseudogap phase [26]. First, the maximum of the antinodal gap is not at the Fermi momentum, c , where the gap due to BCS pairing is. Instead, it shows a peak at the momentum d slightly away from the Fermi momentum. Second, the depth of the gap is almost doubled in c (0, ) than that of the gap in the antinodal point. We further add the temperature and doping dependences to KL . First, in previous two-component models [40], the coupling of carriers to the lattice would decrease upon increasing temperature, so that itinerant carriers are released exponentially with temperature. Also, Hall measurements show that the number of thermally activated holes increases exponentially versus temperature [41]. Thus, we employ a similar exponential temperature dependence to KL . In terms of doping dependence, since scanning tunneling microscopy (STM) measurements show that the energy scale of pseudogap rises linearly with decreasing doping [36], we make KL also decrease linearly versus doping, so that it terminates at a doping concentration = f : Note that the expression is valid only when < f . Equation 7 suggests that the pQon character and its interaction with superplane start to significantly weaken at a temperature comparable to the scale of f upon heating, corresponding to the close of the pseudogap. Besides, KL becomes zero when doping concentration reaches f , and the carriers in cuprates will lose the pQon character, making the system behave more like a Fermi liquid [42]. Here, f is fixed at 0.33 in our simulation, because transport measurements on La2-xSrxCuO4 found that the resistivity increases as p (pure Fermi liquid characteristics) only when the doping concentration is larger than 0.33 [43]. Within the proposed picture, the phenomenon that the 'nodal' metal emerges in the pseudogap phase from the Fermi-liquid-like background with reducing doping [44], is driven by the increasing strength of KL . However, the gap-like feature does not open immediately when < f , instead, the spectrum weight first weakens around antinodal, with the pseudogap only opening up at a lower doping often below the optimal doping (Fig. 5a-c). That is also to say when the pseudogap closes, KL is still a non-zero value. We ascribe this to the fact that the value of antinodal gap is small enough to be overwhelmed by the lifetime broadening. Another important property of the pseudogap phase is that the gapless region evolves from a Fermi arc with limited length to the non-interacting Fermi surface with increasing temperature [45]. We calculate the ARPES spectra at = 0 , as shown in Fig. 5d-f, which show that with increasing temperature and decreasing KL , as described by Eq. (7), the Fermi arc gradually expands, and eventually the noninteracting Fermi surface recovers, consistent with experimental observations [46]. The growth of the arc with temperature was explained by static disorders, such as charge density wave or spin density wave [21][22][23]; here we show that, instead, it could also be explained by the destruction of the strength measure of charge superplane, KL , with increasing temperature. Besides, the arc is not a real locus of quasi-particle states, rather it forms because real spectral weight scatters into the d-wave-type gap at finite lifetime broadening . Therefore, the difficulty in observing the back side of the arc is then natural, since the arc should be just arc-like in hole doped cuprates, instead of being the front side of a small pocket.

Temperature-invariant critical doping around 19%
There is a very recent report of the temperature-independent critical doping at around = 0.19 by ARPES measurements [33], where the strange metal abruptly reconstructs into a more conventional metal with Bogoliubov quasiparticles. We speculate such a temperature-independent effect is likely to have a geometric origin due to percolation. Using the pQon picture, we describe such critical doping as the transition point beyond which the hopping process of pQons is forced to generate itinerant quasiparticles. This corresponds to when the pQon hopping pathways no longer form a percolated network due to the high pQon concentration. Figure 6 shows the result of our percolation simulation. Figure 6a to Fig. 6h). The white regions are the possible pathways where the pQons can hop and land, i.e., diffuse. Note that our picture here involves only one type of carrier of pQon, and the black/white regions and diffusion pathways are dynamically self-organized by the hopping induced redistribution of pQons. It is thus different from the two-carrier/component model [40] and/or defect-induced disorder or inhomogeneity [47]. The percolation transition is shown in Fig. 6d, where two schemes of 1) including both the nodal and antinodal hoppings and 2) just including the antinodal hopping are compared (see Method). They actually give similar percolation thresholds around 20% and 19%, respectively. Thus, the examples in Fig. a-c correspond to the percolated, the threshold, and the non-percolated pQon hopping pathways, respectively. This suggests that our method provides a percolation threshold robust against minor variation of the settings in the model, which is likely to be responsible for the measured critical doping around 19%. If the pQon diffusion pathways form a percolated network (p < 19%), pQons are well conserved during hopping processes. While beyond 19%, the landing event in a pQon hopping process will be frequently impeded by the unavailable sites due to the lack of a percolated pathway of available sites, leading to an itinerant quasiparticle observed in ARPES. However, these Bogoliubov quasiparticles are in a dynamical conversion with pQons, that is, after certain lifetime they will be converted back to pQons to maintain a constant ratio between the two types of carriers at a given hole doping level beyond 19%.
Furthermore, the simulation suggests a cluster size of around 5 nm, which corresponds well with the scales of the spatial variation of gap size in previous STM images [34][35][36]. We then simulate the pseudocolor images in Fig. 6e-g based on our results in Fig. 6a-c by assuming that a pQon surrounded by more unavailable neighboring sites locally is easier to be converted to a electric signal in the STM measurement, corresponding to a smaller local gap. These pseudocolor images compare well with previous STM images showing spatial inhomogeneities in local gap sizes. The results thus suggest that the STM tip may serve as a local probe with bias voltage to extract the information about spatial variation in pQon hoppings, which is often averaged out in ARPES measurements and only emerges as a collective behavior beyond 19%. In the case of STM, instead, the local bias voltage forces the pQon in the black and white regions to respond differently in a broad doping range. Note that due to the different types of perturbations provided by the local bias voltage in STM and the global photon excitation in ARPES, and different methods in extracting the gap values in STM analysis, a percolation threshold obtained from analyzing STM images may or may not deviate from 19%.

Conclusion
We have identified the charge superplane as a critical emergent configuration that could play a central role in understanding the complicated phenomena in cuprate superconductors. The superplanes cut through the dx 2 -y 2 and dz 2 bonding directions and possess an amazing material-dependent, dynamic chargeconservation property. We have further proposed pQon as the new fundamental collective excitation in cuprates, with polaron-Quasiparticle duality. Our ARPES simulation based on the electronic structure of pQon achieves a good agreement with experimental ARPES measurements in a broad doping range. Our picture also naturally leads to a geometric percolation interpretation of the temperature independent critical doping w = 0.19 reported recently by ARPES, and also previous STM imaging results with spatial "inhomogeneity", all using an intuitive percolation hopping picture of pQons. Our work therefore provides a new perspective to further understand the origin of pseudogap phase and other related phases in cuprates, which may also serve to unify the results from various experimental techniques with different momentum, spatial and time resolutions.

DFT simulation
All DFT calculations were carried out by the Vienna Ab initio Simulation Package (VASP) code, which implements the pseudopotential plane wave band method [48]. The projector augmented wave Perdew-Burke-Ernzerhof (PAW-PBE) functional was utilized for the exchange-correlation energy. A 520 eV planewave energy cutoff was implemented for all calculations. Non-spin polarized DFT+U calculations with U = 8 eV and J = 1.34 eV [49] integrated with PHONOPY package [50] were used for phonon simulations using the supercell approach with finite displacements based on the Hellmann-Feynman theorem. 21 phonon modes out of the total 36 modes in a 2x2x1 supercell of Hg-1212 were selected for the calculation of phonon coupling. 24 modes out of 39 modes in a supercell of YBCO with the same size are selected. The other modes have similar ion movements to the selected modes and were excluded due to computational limit. We expect that adding those similar modes back will not induce any new trend in Fig. 2. The phonon modes were calculated at = (0.5, 0.5, 0). In Fig. 2e, all the hypothetical structures were relaxed, and the same atomic displacements of each mode as in Fig. 2a were used. For each pair of phonons, one is chosen as the "frozen mode" and one as the "perturbation mode". Five calculations were performed on each frozen mode at a given amplitude, with the added amplitude from each perturbation mode scaled at -0.5, -0.25, 0, 0.25, 0.5 of the frozen mode amplitude. The five energies of each such calculation is used to fit a parabola energy profile to further extract the strength of phonon coupling through the shift of the equilibrium position of the fitted curve, following previous method [16]. The charge leakage is defined by the ratio of out-ofsuperplane extra charge flow over the in-superplane extra charge flow due to the phonon coupling. For the polaron-quasiparticle conversion calculation in Fig. 3b, Heyd-Scuseria-Ernzerhof (HSE) exchangecorrelation functional with = 0.25 and a 2x2x1 supercell was used. The polaron was created by HSE relaxation following previous method [16]. The quasiparticle was created by HSE relaxation with symmetry enforced as well as non-perturbed initial ionic positions in a charged supercell. The intermediate images and energies in Fig. 3b were generated by linear interpolation of the ionic positions in initial and final images followed by electronic self-consistent relaxations.

ARPES simulation
The spectral functional was calculated using the (1,1) component of Green function within our pQon model: where is the energy, and represents lifetime broadening (set as a constant of 20 meV in our simulations). Note that both higher energy "$ % and lower energy "$ S bands are involved as poles of the Green's function. The spectral weight is defined by: For Fig. 4, ^ was fixed when plotting the spectral weight versus m and , in which m changes from −0.3 to 0.3 . For fig. 5, was fixed at zero and the spectral weight was plotted versus both m and ^ in the entire Brillouin zone.

Percolation related simulations
A 300*300 square lattice model was used for percolation simulation to model the copper sites in one CuO2 plane. Each grid of the lattice corresponds to a copper site. The hole distribution corresponding to a doping concentration was generated by random sampling while avoiding nearest-neighboring holes. A sampled grid with any occupied nearest neighbors was discarded and resampled, until a grid with no occupied nearest neighbors was obtained. After obtaining a legal sampling, all hole positions were first marked black. For a hole currently sitting at site ( , ), its ability to move in the model is given by the following criteria (Fig.   6h). When considering four quadrants {‚− ƒ " , ƒ " … , ‚ ƒ " , †ƒ " … , ‚ †ƒ " , ‡ƒ " … , ‚ ‡ƒ " ,ˆƒ " …} as its possible moving directions, for the quadrant ‚− ƒ " , ƒ " …, if at least one site of {( , + 2), ( + 1, + 1), ( − 1, + 1)} was occupied, or both sites {( + 1, + 2), ( − 1, + 2)} were occupied by other holes, the hole at ( , ) was considered as not movable to the ‚− ƒ " , ƒ " … quadrant. The same rule was applied to the analysis of other quadrants to determine whether the hole can hop along other directions. If the hole could not move along at least three quadrants, this hole was classified as not movable and its four nearest-neighbor sites and four second-nearest-neighbor sites would be also marked as black, as shown in Fig. 6a-c. The ratio of the size of the largest connected white cluster over the size of all white grids is the indicator of percolation. The percolation threshold is defined as the doping concentration where such ratio swiftly changes from very close to zero to very close to one. Note that for infinitely large lattice (i.e., the thermodynamic limit) the transition should be a step function. For finite-size simulation the corners are slightly rounded as shown in Fig. 6. In this case, the midpoint of the step (ratio = 0.5) was used as the approximated percolation threshold. For pseudocolor image simulation ( Fig. 6e-g), we assume that a pQon that is surrounded by more unavailable neighboring sites locally is easier to be converted to a electric signal in the STM measurement, thus having a smaller gap measured by STM. We thus let the black regions in the percolation simulation have gap value 0 and white regions have gap value 1, and average over every 8*8 (about 1.5Å*1.5Å) region to produce the pseudocolors in Fig. 6e-g. Contributions X. C. and J. D. contributed equally to this work. X. L. conceived the physical picture, designed and supervised the research. X. C. performed the DFT and percolation simulations. X. L. and J. D. wrote the Hamiltonian. J. D. performed the ARPES simulation. All authors analyzed the results and wrote the manuscript.