{'fraction_non_alphanumeric': 0.044805267269277985, 'fraction_numerical': 0.01863631291150541, 'mean_word_length': 4.094656054576464, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 2, '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}}
Optimal inhomogeneity for pairing in Hubbard systems with next-nearest-neighbor hopping

Previous studies have shown that bipartite Hubbard systems with inhomogeneous hopping amplitudes can exhibit higher pair-binding energies than the uniform model. Here we examine whether this result holds for systems with a more generic band structure. To this end, we use exact diagonalization and the density matrix renormalization group method to study the 4x4 Hubbard cluster and the two-leg Hubbard ladder with checkerboard-modulated nearest-neighbor hopping, t, and next-nearest-neighbor (diagonal) hopping, t_d. We find that the strongest pairing continues to occur at an intermediate level of inhomogeneity. While the maximal pair-binding energy is enhanced by a positive t_d/t, it is suppressed and appears at weaker repulsion strengths and smaller hole concentrations when t_d/t is negative. We point out a possible connection between the pairing maximum and the magnetic properties of the system.


I. INTRODUCTION
Consider two identical N -site clusters, each doped away from half filling by M holes and residing in its ground state of energy E 0 (M ). If holes tend to pair and if M is odd it is energetically preferable to move one hole between the clusters in order to achieve a paired state on both. In this sense a positive pair-binding energy to the lowest S = 1 excitation follows ∆ pb it is an indication that the excitation is a result of a dissociation of a singlet hole pair into two separate holes. In one dimension the relation is even more explicit as the opening of a spin gap entails a non-zero amplitude for the superconducting order parameter [1]. Exact diagonalization studies have shown that a number of small Hubbard clusters exhibit pair binding, which reaches a maximum at an intermediate strength of the on-site repulsion [2,3]. Similar behavior was observed using the density-matrix renormalization group (DMRG) method in the two-leg Hubbard ladder, where the binding energy is of the order of the spin gap and where both diminish with doping [4,5]. These findings have inspired searches for superconductivity in two-dimensional systems composed of coupled lower-dimensional building blocks, in hope of harnessing the pairing tendencies of the latter. Such a strategy naturally gives rise to the question of what is the optimal level of inhomogeneity for superconductivity [6,7].
Much of the research into the relationship between inhomogeneity and superconductivity from repulsive interactions has been carried out using the plaquette Hubbard model [8,9]. The model is constructed from 2 × 2 plaquetes with on-site repulsion U and nearest-neighbor hopping t, where neighboring sites on different plaquettes are coupled by hopping t ′ . Exact diagonalization of the 4 × 4 site system [10] has found a substantial maximum of the pair-binding energy at t ′ /t ≈ 0.5, U/t ≈ 8 and low hole doping. A similar pairing maximum, occurring at intermediate inhomogeneity levels t ′ /t ≈ 0.5−0.7 and interaction strengths U/t ≈ 5 − 8, was subsequently found in larger plaquette systems using the contractor renormalization (CORE) [11] and DMRG [12] methods. Furthermore, by calculating the other necessary ingredient for superconductivity, namely phase stiffness, these studies have provided evidence that optimal inhomogeneity likely exists also for the superconducting transition temperature, T c , and not just for the pairing scale. These findings were contested, however, by calculations using the dynamical cluster approximation (DCA) [13] and cellular dynamical mean-field theory (CDMFT) [14], which have obtained a monotonic increase with t ′ /t of both the d-wave pairing interaction and of T c toward a maximum that is exhibited by the homogeneous model. Nevertheless, a recent quantum Monte-Carlo (QMC) study [15] provides support in favor of the CORE and DMRG findings. While the sign problem prevents reliable calculation of T c , it is manageable to low enough temperatures in order to show that for U/t = 4 the pairing vertex is most attractive at t ′ /t ≈ 0.4 . Finally, despite differences in details DCA [16], CDMFT [17] and QMC [18] have all detected enhanced superconductivity in Hubbard models with inhomogeneous charge density due to external potentials.
To date, optimal inhomogeneity for pairing, or more generally for superconductivity, has been demonstrated only on the bipartite square lattice. It is therefore interesting to explore the robustness of the phenomenon to changes in the band structure, not the least because they are present in the cuprate high-temperature superconductors. Specifically, the cuprates are often modeled using a tight-binding band structure of a square lattice with hopping amplitudes that extend beyond the nearestneighbor amplitude t. In particular, it is necessary to include next-nearest-neighbor (diagonal) hopping t d , with t d /t < 0, to account for the observed Fermi surfaces [19], and there are indications that it plays a role in the physics governing T c of the hole-doped systems [20]. The effects of such a term have been investigated both in the context of the Hubbard model [21][22][23][24][25][26][27] and its strong-coupling descendent, the t − J model [28][29][30][31][32]. It appears that different studies agree that various measures of pairing and superconductivity are suppressed in the presence of t d /t < 0, at least for hole doping in the range x < 0.12. The results vary, however, for higher doping levels where calculations using DCA [24] and DMRG [28] continue to find suppression of superconductivity while density matrix embedding theory [27] and variational QMC [29] indicate enhancement of pairing correlations. A similar dichotomy also exists for positive t d /t where the first group of methods finds enhanced superconductivity while the second yields an opposite trend.
Here we study the existence of optimal inhomogeneity for pairing in the plaquette Hubbard model with diagonal hopping. To this end, we calculate the pair-binding energy and the spin gap using exact diagonalization and DMRG. We show that pairing continues to peak at intermediate levels of inhomogeneity but its strength depends on the sign of t d /t. Our results indicate that pairing is enhanced by the presence of t d /t > 0. On the other hand, when t d /t < 0 pairing is suppressed for the higher hole concentrations examined near x = 0.12. It regains strength, however, at lower doping levels, or when U/t is reduced. We note that these effects can not be understood on the level of a single plaquette and speculate on their possible connection to the magnetic properties of the system.

II. MODEL AND RESULTS
We consider the plaquette Hubbard model where c † iσ creates an electron of spin polarization σ =↑, ↓ at site i, and n iσ = c † iσ c iσ . Here i, j and i, j denote nearest-neighbor and next-nearest-neighbor sites, respectively. The hopping amplitudes are modulated as shown in Fig. 1. Neighboring sites within a plaquette are connected by hopping t while next-nearest-neighbors are connected by t d . The corresponding amplitudes across plaquette boundaries are t ′ and t ′ d . For simplicity we restrict ourselves to the case in which the diagonal amplitudes are modulated with the same ratio as the nearestneighbor amplitudes, i.e., t ′ d /t d = t ′ /t. We have studied the model on the 4 × 4 cluster and on the two-leg ladder, depicted in Fig. 1. To obtain an estimate for the finite- 3. The pair-binding energy of the 4 × 4 cluster with periodic boundary conditions at hole doping 1/16 and various t d /t. The results are for U/t = 1(red), 2, · · · , 10(green). depicted in Fig. 2 for the case U/t = 8. We find that the pair-binding energy of the x = 1/16 system is largely insensitive to changes in the boundary conditions over a range of t ′ /t that shrinks with increasing |t d /t|. For |t d /t| ≤ 0.1 (not shown) the two sets of results follow each other and differ by at most 20% as long as 0 ≤ t ′ /t < ∼ 0.8. In particular, a clear maximum in ∆ pb is observed at t ′ /t = 0.5. This maximum is also present in the t d = −0.3t data, but due to the increased sensitivity to the boundary conditions above t ′ /t ≈ 0.65 we are unable to determine whether it constitutes a global pairing maximum. The even larger sensitivity of the t d = −0.5t results precludes reaching a conclusion about the existence of optimal inhomogeneity for pairing in this case. At the same time, the results for positive t d are more robust and ∆ pb (1/16) exhibits a consistent maximum around t ′ /t = 0.3 − 0.5 with a clear enhancement compared to the t d = 0 case for inhomogeneity levels below the maximum. The evolution of ∆ pb with U/t for t d < 0 is presented in Fig. 3. Evidently, within the range of t ′ /t discussed above the position and the magnitude of the maximal binding energy increases with U/t until it reaches a global maximum around U/t ≈ 8 − 10.
Positive pair-binding energy may also be associated with a tendency of the system to phase separate. In order to distinguish between pairing and phase separation one needs to calculate the surface tension between the holerich and hole-poor phases [33]. A cruder way is to look for negative inverse compressibility, as a sign for instability towards phase separation. We, however, always find its discrete version to be positive. Further support for pairing comes from the fact that the spin gap of the system with two doped holes roughly follows ∆ pb (1/16), as presented in Fig. 2. In contrast, ∆ pb (3/16) exhibits large sensitivity to the boundary conditions and we can not determine whether holes pair on the cluster at this higher doping level.
In an effort to substantiate the exact diagonalization study we have used DMRG to calculate ∆ s and ∆ pb of the two-leg plaquette ladder with t d /t = ±0.3. During the calculation we have truncated the density matrix, keeping up to about 3200 states in order to reach low enough truncation errors. The relatively large number of kept states (larger than needed when t d = 0 [12]) meant we could deal with ladders of up to 2 × 32 sites. Our results for the ladder with hole doping x = 1/16 are summarized in Figures 4 and 5, and for x = 1/8 in Figures 6 and 7. All quantities are given in units of t. Concentrating first on the t d /t = −0.3 case, we find general agreement between the DMRG results and those obtained using exact diagonalization. Specifically, a positive pair-binding energy, accompanied by κ −1 > 0, is observed for x = 1/16 and exhibits a peak at an intermediate inhomogeneity level. The peak is robust in the U/t = 4 ladder, where the results change little upon increasing the length of the system. On the other hand, the results for U/t = 8 still show substantial size dependence in going from the 2 × 16 to the 2 × 32 ladder. Nevertheless, the tendency of ∆ pb to grow with the system size and the fact that it shows similar features to ∆ s make it plausible that optimal inhomogeneity for pairing also exists in the thermodynamic limit. Increasing the hole concentration to x = 1/8 leads to a ∆ pb that is indiscernible when U/t = 8. This stands in contrast to the t d = 0 ladder under similar conditions where ∆ pb attains a maximal value of about 0.15t [12]. Optimal pairing reappears upon lowering the interaction strength to U/t = 4 but it is still somewhat weaker than its value when t d = 0 [12]. We therefore conclude that while optimal inhomogeneity continues to exist in the presence of

III. DISCUSSION
What might be the origin of the pairing maximum and its dependence on t d ? One may try to look for the reasons in the Hubbard plaquette itself. Fig. 8 shows, however, that the pair-binding energy of the square is already negative at the U/t values which correspond to the pairing maximum observed in the extended systems. Therefore, the latter is not a single plaquette effect. Nevertheless, ∆ pb of the large and small systems share some characteristics, such as its tendency to decrease as one moves from positive to negative t d /t. It has been suggested in the context of the t − J model [31], that this may be due to the fact that negative t d /t is less favorable for creation of paired states with d-wave symmetry.
It can be shown that for |t d /t| < 1 the ground state of the undoped plaquette, |N h = 0 , is an S = 0 singlet which is odd under π/2 rotations, and approaches the "RVB" state (1/ √ 12)(P † 12 P † 34 − P † 14 P † 23 )|0 at large U/t. Here, P † ij = c † i↑ c † j↓ + c † j↑ c † i↓ creates a singlet electron pair on sites i and j (numbered sequentially around the square). When the interactions are weak the two holes ground state, |N h = 2 , is generated from |N h = 0 by P 12 − P 23 + P 34 − P 41 . However, as U/t is made stronger the pairing operator includes terms that create holes on next-nearest-neighbor sites as well [34,35]. We find that the effect of these terms increases when t d /t turns from positive to negative, see the inset of Fig. 8, but for all values |t d /t| < 1 the state |N h = 2 remains a spin singlet that lies in the identity representation of the rotation group. Therefore, the pair annihilation operator connecting |N h = 0 and |N h = 2 must transform as d x 2 −y 2 . In this sense the inclusion of diagonal hopping does not affect the pairing symmetry on the square. Nevertheless, it does affect the energetics of the pairing process. Interestingly, both |N h = 0 and its energy are independent of t d (in the range considered here). At the same time, the energies of |N h = 2 and of the single hole ground state |N h = 1 increase with U/t in a rate which depends on t d . We find that the suppression of the plaquette ∆ pb at large U/t and t d /t < 0 is largely driven by the slower increase of the single hole energy. In particular, for −1 < t d /t < 0.25 there exists a critical U/t (that increases with t d ) where |N h = 1 turns from a degenerate quartet with S z = ±1/2 and plaquette momentum (0, π) or (π, 0), to a degenerate quartet with constant energy −2t + t d , made of S z = ±1/2 and S z = ±3/2 doublets with momentum (π, π), thereby leading to the break seen in Fig. 8.
As noted above, while some features of the t d dependence of ∆ pb follow the behavior found on the single plaquette level, the existence of a pairing maximum as function of t ′ /t can not be understood from such considerations. Instead, we would like to point out a correlation between the pairing maximum and the magnetic properties of the system. It is well known that in the large U/t limit the half filled Hubbard model maps onto the S = 1/2 Heisenberg model with J = 4t 2 /U [36]. Correspondingly, the plaquette Hubbard model maps onto the plaquette Heisenberg model with J ′ /J = (t ′ /t) 2 . Since the ground state of the uniform model (J ′ = J) exhibits Néel antiferromagnetic (AFM) long-range order, and that of the disconnected system (J ′ = 0) is a product of RVB states on individual plaquettes, one expects that a quantum critical point (QCP) separates the two at an intermediate J ′ /J. Such expectation has been borne out by numerical calculations [37,38] which find a QCP at J ′ /J ≃ 0.55. A recent QMC study of the plaquette Hubbard model [15] provides evidence that this QCP survives at half filling for lower values of the interaction strength. Interestingly, the observed QCP at t ′ /t ≈ 0.5, 0.6 for U/t = 4, 8 resides in the vicinity of the t ′ /t value for which the product of the pairing vertex and the uncorrelated pairing susceptibility is closest to -1, where a superconducting instability would develop. This coincidence joins a related behavior which we have noticed in our CORE study [11] of the model away from half-filling. While no transition to long-range AFM order is observed (nor expected beyond small hole doping levels), the maximal ∆ pb does occur at t ′ /t around which AFM correlations build up from the RVB background. It is therefore possible that the enhanced magnetic fluctuations generated by the inhomogeneity-induced QCP, or its related crossover at finite doping, are the mediator responsible for the enhanced pairing.
The presence of diagonal hopping turns the large-U limit of the half-filled Hubbard model into the J 1 − J 2 Heisenberg model with AFM couplings satisfying J 2 /J 1 = (t d /t) 2 . The frustration introduced by the nextnearest-neighbor AFM coupling, J 2 , causes a sequence of quantum phase transitions where Néel order is first lost at J 2 /J 1 ≈ 0.4, in favor of a nonmagnetic state that is either columnar or plaquette valence-bond-solid, which then gives way to a collinear (striped) magnetic order at J 2 /J 1 ≈ 0.6 [39][40][41]. Variational cluster approximations of the half-filled Hubbard model with positive t d /t point at a similar picture where a nonmagnetic phase in the range 0.7 < ∼ t d /t < ∼ 0.8 separates the magnetically ordered states down to U/t ≈ 5, where it spreads out [42,43]. The fact that these transitions occur at relatively large values of t d /t suggests that they do not play a role in establishing the results presented by us, or in the cuprate superconductors. Nevertheless, the extent of the nonmagnetic phase grows in the plaquette J 1 − J 2 Heisenberg model [44,45]. Ref. 45 predicts that for J 2 /J 1 = 0.1 a transition to a Néel state takes place once the intra-plaquette couplings are about twice the interplaquette ones. This would correspond to a transition in the large-U limit of the half-filled plaquette Hubbard model with |t d /t| = 0.3 at t ′ /t ≈ 0.7, not too far from where we observe the pairing maximum in the system with t d /t = 0.3 and U/t = 8. However, this observation does not explain why the pairing maximum seems to appear at lower values of t ′ /t when t d /t = −0.3. To answer this question and strengthen the conjectured tie between pairing and a magnetic QCP further study of the doped plaquette Hubbard model is called for.