Journal of the Physical Society of Japan

FULL PAPERS

Band-Renormalization Eﬀects on Superconductivity and Antiferromagnetismin Two-Dimensional t-J Model
Ryo Sato ∗ and Hisatoshi Yokoyama
Department of Physics, Tohoku University, Sendai 980-8578, Japan
Coexistence or exclusivity of d x2 −y2 -wave superconducting (d-SC) and antiferromagnetic (AF) orders and instability
toward phase separation (PS) are reconsidered for the square-lattice t-J model with the diagonal transfer (t′ ) term nearhalf filling. To reliably treat the strong correlation including the local constraint of no double occupancy, we use a variational Monte Carlo method, in which band-renormalization eﬀect (BRE) is also incorporated. The results are comparedwith those of the corresponding Hubbard model recently obtained [J. Phys. Soc. Jpn. 85, 074701 (2016)]. It is foundthat BRE is very eﬀective for the AF order and greatly modifies the ground-state phase diagram in the t′ -δ space (δ:doping rate), especially for large |t′ /t|, proposed by previous studies; in a wide range of t′ /t (|t′ /t| ≤ 0.5) and for any
underdoped δ (≲ 0.15), the AF order arises for J/t = 0.3 similarly to the Hubbard model. In contrast with the Hubbard
model, however, the d-SC order also arises down to small δ in a wide range of t′ /t (≳ −0.2), and coexists with the AForder even for t′ /t < 0 (type-II AF regime). Furthermore, there is no instability toward PS for any t′ /t and δ for J/t = 0.3.

been known that the d-SC order is predominant and overcomes antiferromagnetism (AF) even for very small doping
rates. Furthermore, the domain of d-SC extends up to sizable
doping rates δ ∼ 0.5 for negative large values of t′ /t with amoderate J/t.21–25) The t-J model seems to have a somewhat
stronger tendency toward d-SC than the Hubbard model.
Because the t-J model works in the strong-correlation limit,
methods applicable to the present subject beyond mean-fieldlevels are restricted. Among others, the VMC method26) has a
crucial advantage, namely, it can treat correlation factors exactly, including the local constraint of no double occupancy,with suﬃciently large systems. Thus, this method has been
often applied to the problems of (a) coexistence and exclusivity of the two orders and of (b) instability toward PS.27–33)
However, the conclusions of (a) and (b) were not necessarily
unified among the studies but seemed rather scattered. Furthermore, the physics underlying them has not been clarified.
The previous paper,16) addressing this subject in the Hubbard
model, showed that vital to these problems is the value of t′ /t,whose importance had been often disregarded. Here, we will
carry out a similar analysis for the t-J model.
In Table I, we summarize the conclusions as to the above
two problems of relevant VMC studies for the t-J model andof recent ones for the Hubbard model. We add a result of theextended Gutzwiller approximation34) as a representative ofanalytic calculations that respect the local constraint. The conclusions are classified according to t′ /t = 0 or ∼ −0.3 (typicalvalues of hole-doped cuprates). For the simple square lattice(t′ /t = 0), all studies unanimously showed the coexistence ofthe two orders, and the values of the critical doping rate δAFare broadly similar. This feature is also common to that of the
Hubbard model. However, the t-J model with J/t ∼ 0.3 hasstability against PS, in contrast with the Hubbard model with
U/t = 12. On the other hand for t′ /t ∼ −0.3, the conclusions
as to whether the two orders coexist or not and the values ofδAF largely depend on the studies. This dependence seems toresult from whether the band renormalization (BR) is properly
introduced or not.
In this paper, we study the interplay of AF and d-SC or-

1. Introduction
In cuprate superconductors (SCs),1) the antiferromagnetic(AF) orders in the insulating parent compounds rapidly vanish on carrier doping and the d x2 −y2 -wave superconducting (dSC) orders arise in wide ranges of doping rate (δ).2) Formerly,
many theoretical studies supported this fundamental feature
as well as other properties of d x2 −y2 -wave superconductivity(d-SC) on the basis of Hubbard-type or t-J-type models.3–12)
However, recent numerical studies on square-lattice Hubbard(t-t′ -U) models using refined techniques which simultaneously consider the two orders13–15) argued that AF long-rangeorders prevail in wide ranges of δ and t′ /t (diagonal hopping)or phase separation (PS) occurs for t′ ∼ 0, and as a result, theranges where pure d-SC orders appear are limited to narrowwindows in the overdoped regime.
In a previous study,16) applying a variational Monte Carlo(VMC) method to a strongly correlated Hubbard (t-t′ -U)model, we obtained basically the same results: A bandrenormalization eﬀect (BRE) plays an essential role for thestabilization of AF orders for large values of |t′ /t| through retrieving the nesting. Consequently, states with AF orders or
PS become remarkably more stable than a pure d-SC state in
the whole underdoped regime and the whole realistic range oft′ /t (Fig. 27 in Ref. 16). This result is serious in that it is obtained using a variation theory, which gives an upper boundof the exact energy. Namely, for realizing a pure d-SC order, breakthrough improvement is required in the trial d-SC
states, notwithstanding such d-SC states already properly describe various aspects of cuprate SCs as mentioned above. It
is important to check whether this perplexing feature is characteristic of the Hubbard model.
In this context, it is significant to study the same subjects
in the t-J model,17) another key model of cuprates. The t-J
model possibly has aspects diﬀerent from the Hubbard model
for realistic values of J/t (∼ 0.2 – 0.5), although the two models are connected for J/t (t/U) → 0 and δ ∼ 0.18–20) In thestudies of pure (single-order AF and d-SC) states,5–8) it has∗ satoryo@cmpt.phys.tohoku.ac.jp

1

J. Phys. Soc. Jpn.

FULL PAPERS

Table I. Results of interplay between the AF and d-SC orders in mixed
states are compared among relevant VMC (and extended Gutzwiller approximation) studies for the square-lattice t-J model and recent VMC studies
for the Hubbard model (lowest three). The upper row of the second column
shows whether (or how) BR is introduced. The total number of BR parameters, if known, is entered in brackets. The lower row indicates a typical valueof J/t (or U/t) studied. The upper row of the third and fourth columns indicates whether the state is coexistent or exclusive (with realized order). The
lower row indicates whether the realized state is homogeneous or phase separated, and the value of δAF , above which the AF orders vanish, is shown inbrackets. For δ > δAF , d-SC appears for all cases shown here.
Literature(Method, Year)

Ref. 34(Ext.GA, ’03)

Ref. 29(VMC, ’99)

Ref. 30(VMC, ’04)

Ref. 31(VMC, ’04)

Refs. 22, 32(VMC,’06,’08)

Ref. 24(VMC, ’09)

Ref. 25(VMC, ’16)

This work(VMC, ’18)

Ref. 35(VMC, ’16)

Ref. 14(m-VMC, ’14)

Ref. 16(VMC, ’16)

BR
(J/t, U/t)
No
(0.3)
No
(0.3)
No
(1/3)
Partially
(0.3)
No
(0.2)
No
(0.3)
Partially (4)(0.3)
Yes (8)
(0.3)
No
(10)
Implicitly(10)
Yes (8)(12)t′ /t = 0
(δAF )
Coexistent— (∼ 0.11)
CoexistentHomo. (∼ 0.1)
CoexistentPS (∼ 0.1)
Coexistent— (∼ 0.1)
CoexistentHomo. (∼ 0.12)——
CoexistentHomo. (∼ 0.14)
CoexistentPS, Homo. (∼ 0.12)
CoexistentPS (∼ 0.18)
CoexistentPS (∼ 0.15)

(a) H

(b) Ȍ

t’

t1

t

t

t3
t2

t4

Fig. 1. (Color online) (a) Hopping processes in Hamiltonian [Eq. (2)] and(b) those corresponding to band-adjusting parameters tη (η = 1 – 4) in trial
wave functions [Eqs. (12) and (18) – (21)]. In both panels, t is the unit.

t′ /t ∼ −0.3
(δAF )——

2.

—

Formulation

After introducing the model Hamiltonian in Sect. 2.1, we
explain the VMC method, especially, the trial wave functionsin Sect. 2.2. Relevant quantities are defined in Sect. 2.3.

Exclusive, AF— (∼ 0.06)
Coexistent— (∼ 0.1)
Coexistent— (∼ 0.1)
CoexistentHomo. (∼ 0.08)
Exclusive, AF
Homo. (∼ 0.20)—

2.1 t-t′ -J model
We consider the t-J model on the square lattice17) with adiagonal transfer term [Fig. 1(a)]H = Hhop + H J = Ht + Ht′ + H J
∑ (
∑(
)
)
= −t
c̃†iσ c̃ jσ + H.c. − t′c̃†iσ c̃ jσ + H.c.
⟨i, j⟩,σ

Exclusive, AF
Homo. (∼ 0.24)
Exclusive, AF
Homo. (∼ 0.25)

+J

∑(

)
1
Si · S j − ñi ñ j ,4⟨i, j⟩(1)(2)(i, j),σ(3)where c̃ jσ = c jσ (1−n j−σ ) is the annihilation operator that actsin the space without double occupation, c jσ is the ordinaryannihilation operator for the Wannier state, ñ jσ = c̃†jσ c̃ jσ , and∑S j = 21 α,β c†jα σαβ c jβ (σ : Pauli matrix). ⟨i, j⟩ and (i, j) in thesummations in Eqs. (2) and (3) indicate the pairs of nearestneighbor and diagonal-neighbor sites, respectively. Thus, thebare band dispersion is()ε̃k = −2t cos k x + cos ky − 4t′ cos k x cos ky .
(4)ders in the t-t′ -J model by comparing pure (AF and d-SC)
states and a mixed state of the two orders, in the frameworkof VMC. Into the trial states, we introduce BR parameters independently for the AF and d-SC parts. We also take accountof nearest-neighbor correlation factors which distinguish all
patterns of bond configurations. In this way, we clarify the following points: (i) How the pure AF state is stabilized by BRE
and exhibits a kind of Lifshitz transition in which the loci of
pocket Fermi surfaces (FSs) switch, as seen in the Hubbard
model.16) (ii) How BRE aﬀects the pure d-SC state. (iii) How
the property of coexistence or exclusivity of the two orders
evolves as t′ /t and δ are varied in the mixed state. And, we
will trace the cause of this property by comparing the behavior
between the pure states and mixed state. (iv) Whether or not
the states are stable against PS. As a summary, we construct aphase diagram in t′ -δ plane using the mixed state (Fig. 21), in
order to compare with that for the Hubbard model.
This paper is organized as follows: In Sect. 2, we explain
the model and method used. In Sect. 3, we discuss BRE on
the pure AF and d-SC states. In Sect. 4, we study the interplay of the two orders using the mixed state. In Sect. 5, werecapitulate this work. In Appendix A, we describe the details of optimizing the trial wave functions. In Appendix B,
we argue that the energy reduction in the AF state is primarily caused by FS renormalization but not by the band form.
Preliminary results of the present study were published in a
conference proceedings.36)

The value of t′ /t depends on the kind of cuprate SCs andis considered −0.1 (LSCO) – −0.3 (YBCO,BSCCO) for holedoped systems and ∼ 0.3 for electron-doped systems (throughelectron-hole transformation).37) According to band calculations,38, 39) the existence of apical oxygen tends to increaset′ /t. We repeat that the value of t′ /t plays a leading role in
the present subject. As for J/t, we mainly discuss a typicalcase for cuprates J/t = 0.3, although we sometimes refer to
J/t dependence. The Hamiltonian Eq. (1) is connected to the
Hubbard Hamiltonian through the strong-coupling expansion
for J/t, t/U → 0 with J = 4t2 /U near half filling.18) Thereby,
J/t = 0.3 corresponds to U/t = 13.3. We use t and the latticespacing as the units of energy and length, respectively.
2.2 Variational Monte Carlo method
A merit of the VMC method26) is the possibility of specifying important factors or physics (e.g. BRE) in the systemof our interest. To this end, it is vital to construct simple
trial wave functions that capture the essence. In this study,
as many-body trial states, we use a Jastrow type, Ψ = PΦ,40)
where P is a two-body correlation factor (projector) and Φ is2

J. Phys. Soc. Jpn.

FULL PAPERS
Table II. Comparison of form of εΛk and choice of {k}occ for band renormalization in mixed states among related studies. The second column (d-SC)
shows the energy dispersion to be optimized in the one-body d-SC part andthe third column (AF) the same but in the AF part. In the fourth column, thetotal number of band parameters are entered. The rightmost column indicates
how to choose the occupied k-points in the AF part.

a one-body (mean-field-type) wave function.
We take the form of P common to all trial states, P =∏PG PJ . Here, PG is the onsite projector PG = j [1−n j↑ n j↓ ]41)
that ensures the above local constraint. In the intersite factor PJ , we consider both charge-density and spin correlations
between nearest-neighbor sites: PJ = Ph PS . For the chargedensity part, we use an interhole correlation,
∏[
]
Ph =
1 − (1 − α) hi h j(5)⟨i, j⟩(6)

In Eqs. (5) and (6), α, β1 and β2 are variational parameters.
Ph is eﬀective, especially, in a highly doped regime. Ps plays
a corrective role for the (often overestimated) AF order. Using Ph and P s , we can assign distinct weights to all diﬀerent
nearest-neighbor bond configurations.
We turn to the one-body part Φ, into which BRE33, 42) is introduced. The forms used in the previous study16) are applied
to the present cases. We start with the normal (paramagnetic)state (Fermi Sea),∏ΦN =c†k,σ |0⟩,
(7)

4

Order of γk

Ref. 16εAFk

8

Order of εAF
k

εAF
k

8

Optimizing Ψmix

εSC
k −µ+

√

∆k
2
2
(εSC
k − µ) + ∆k

,

(15)

with ∆k = ∆d (cos k x − cos ky ). ∆d and µ are variational parameters corresponding to the d-SC gap and chemical potential.
In εSCk , BRE is introduced. Φd is reduced to ΦN for ∆d = 0
and µ = µ0 (µ0 : the non-interacting chemical potential).
The mixed state, which can simultaneously have d-SC and
AF orders,27–29) are constructed by replacing the bare electron
operator c†k,σ in Eq. (14) by the AF quasiparticle operator a†k,σin Eqs. (9) and (10) as28)
 N2


∑
†
†
Φmix =  ϕ(k) a a  |0⟩.
(16)k↑ −k↓kΦmix is reduced to the pure states in certain limits, namely, Φdfor ∆AF = 0 and ΦAF for ∆d = 0 and µ = µ0 .
For introducing BRE into ΦAF , Φd and Φmix , we extend
SC
the band dispersions εAF
k in Eq. (13) and εk in Eq. (15) byincluding tight-binding hopping terms up to three-step processes shown in Fig. 1(b) as,
Λ
Λ
Λ
Λ
εΛ
k = γk + ε1 (k) + ε2 (k) + ε3 (k) + ε4 (k),(17)

with Λ = SC or AF and

In Eqs. (9) and (10), k is taken in the folded AF Brillouin zone
with the lower-band energy dispersion,√
U
EkAF =− γk2 + ∆2AF ,(11)
2γk = −2t(cos k x + cos ky ),(12)
Q = (π, π), sgn(σ) = 1 (−1) for σ =↑ (↓), and
v
u
t
εAF
1
1 − (+) ( )2k.
αk (βk ) = √2
εAF+ ∆2AFkSCεAFk ≡ εk

εSC
k
εSC
k

ϕ(k) ={k}occ , σ

(10)

εSC
k

where

where, a†k,σ is an AF quasiparticle operator given asa†k+Q,σ = −sgn(σ) βk c†k,σ + αk c†k+Q,σ .

Ref. 25k↑ −k↓

where {k}occ indicates the set of occupied wave-number pointsof k ∈ kF . Although there is no explicit variational parameterin ΦN , we have to determine {k}occ variationally as in the AFcase when BRE is considered (see below and Appendix A).
For the (commensurate) AF state, we apply the form ofHartree-Fock (HF) solution for t′ = 0 at half filling, in whichcase the nesting condition is completely satisfied, to the othercases of t′ /t and δ in a sense of variation theory:∏ΦAF =a†k,σ |0⟩,
(8)

(9)

Param.
4k{k}occ , σ

a†k,σ = αk c†k,σ + sgn(σ) βk c†k+Q,σ ,AF
γkdescribe, in Appendix A, how we deal with the optimizationof {k}occ in this study. ΦAF is reduced to ΦN for ∆AF = 0.
The pure d x2 −y2 -wave singlet pairing (BCS) state of a fixed
electron number N is given by44)

 N2
∑

†
†
(14)
Φd =  ϕ(k) c c  |0⟩,⟨i, j⟩[]
× 1 − (1 − β2 ) (ñi↑ ñ j↓ + ñi↓ ñ j↑ ) .

d-SC
εSC
k

Present workwith h j = (1 − ñ j↑ )(1 − ñ j↓ ). For the spin-dependent part, we
use the form,
∏[
]
Ps =
1 − (1 − β1 ) (ñi↑ ñ j↑ + ñi↓ ñ j↓ ){k}occ in AF part
Order of γk

Literature
Ref. 33ΛεΛ
1 (k) = −4t1 cos k x cos ky ,

(18)

Λ
εΛ
2 (k) = −2t2 (cos 2k x + cos 2ky ),

(19)

Λ
εΛ
3 (k) = −4t3 (cos 2k x cos ky + cos k x cos 2ky ),

(20)

Λ
εΛ
4 (k) = −2t4 (cos 3k x + cos 3ky ).

(21)

These band-adjusting (variational) parameters tηΛ /t (η′= 1–4) are independent of the model parameter t /t in H and of
tp /t used for determining {k}occ (see Appendix A). Thus, wehave four band parameters in ΦAF and Φd , and eight in Φmix .
In Φmix , it is crucial to optimize tηSC and tηAF , independently.
Actually, for the Hubbard model,16) only t1SC /t is eﬀective inthe d-SC part, but tηAF /t up to η = 4 become important for thenesting in the AF part. Furthermore, there is room for how tochoose {k}occ in the AF part. In Table II, we compare the ways

(13)

∆AF corresponds to the AF gap in the sense of HF theory.43) In
Eq. (13), allowing for BRE (explained shortly), we write εAFkfor the energy dispersion, instead of γk . The choice of {k}occ ,namely FS, is included in the BR processes of ΨN , ΨAF andΨmix , and is not easy operation for finite-size systems.11) We
3

J. Phys. Soc. Jpn.

FULL PAPERS

of introducing BRE into mixed states among related studies.
In this work (the last row), we take tηAF /t and tηSC /t as eight
independent band parameters, and determine {k}occ so as tominimize the total energy with respect to Ψmix , instead of using some guiding energy dispersion. Therefore, the present
way of introducing BRE comprehends the other ways in Table II as special cases. The independence of tηAF /t and tηSC /t
is, in particular, important; otherwise, qualitatively diﬀerent
results are sometimes derived,25, 33) as shown in Table I above
and Table IV in Ref. 16. We describe the details of optimizing{k}occ in this work in Appendix A.
Variational expectation values of Ô with respect to ΨΛ (=PΦΛ ; Λ = N, AF, d, and mix),⟨Ô⟩Λ =

⟨ΨΛ | Ô |ΨΛ ⟩⟨ΨΛ |ΨΛ ⟩d-SC J/t=0.3 t′/t=-0.3

t1/t t2/t t3/t t4/t

0.2
0.1
0
-0.1
-0.2

L=10

-0.3
0

0.05

0.1

0.15

δ

0.2

t1
t2
t3
t4

L=12

0.25

0.3

Fig. 2. (Color online) The optimized band parameters tηSC ≡ tη (η = 1–4)
in Ψd are depicted as a function of doping rate for t′ /t = −0.3. The values oftη without BR are indicated by dotted lines in coordinate colors.

(22)

are numerically estimated using a VMC method26) for finite
systems of Ns = L × L sites (L = 10 – 24) with the periodicantiperiodic boundary conditions. The number of samples forcomputing ⟨Ô⟩ are typically 2.5 × 105 . Optimization of thevariational parameters except for tp /t (see Appendix) is performed using the stochastic reconfiguration method.45)

BR /t
Table III. Energy improvement owing to band renormalization ∆EΛ(Λ =AF, SC) [Eq. (28)] in the pure states for three typical values of t′ /tand doping rate. J/t = 0.3 and L = 10. For comparison, we show the corresponding data of the Hubbard model with U/t = 12 borrowed from Ref. 16.
The last digit of each has some error.

Statesd-SC

2.3 Physical quantities
Here, we define quantities often referred to. The total energy per site of the state Λ is written as EΛ = ⟨H⟩Λ /Ns . As
the AF order parameter, staggered magnetization,2 ∑ iQ·r j zm=e⟨S j ⟩ ,
(23)
Ns j

AF

d-SC (Hubbrad)is used, where Q is the AF nesting vector (π, π). m becomes
unity for the full moment. As the indicator of d-SC order, we
use the real-space d-SC correlation function for the nearestneighbor-site pairing, Pd ≡ Pd (r∞ ) with r∞ = (L/2, L/2) and
⟨
⟩
1 ∑ ∑
′
(−1)1−δ(τ,τ ) ∆†τ (Ri )∆τ′ (Ri + r) , (24)
Pd (r) =Ns i τ,τ′ =x̂,ŷ

AF (Hubbrad)t′ /t \ δ
−0.3
0.0
+0.3
−0.3
0.0
+0.3

0.0
0.00090.00010.00090.06010.0001
0.05860.08
0.00170.00020.00010.04050.00090.03070.16
0.00380.00050.00010.02510.00100.0087

−0.3
0.0
+0.3
−0.3
0.0
+0.3

0.00190.00000.00010.16530.00580.16530.00040.00020.00030.06220.01020.06030.00010.00050.00050.02200.00200.01083.1 d-wave superconducting state
In a previous study16) for the Hubbard model, it was foundthat BRE in Ψd aﬀects only very slightly physical quantities as well as energy for large U/t and |t′ /t|, compared toΨAF , although the renormalization of t1SC /t in Ψd is neversmall.21, 23, 46) In Fig. 2, we plot δ dependence of the optimized band parameters of the d-SC states for the t-J modelwith t′ /t = −0.3. For δ ∼ 0 and in the overdoped regime(δ ≳ 0.16), t1 /t is greatly renormalized, namely, going awayfrom −0.3; the large BR in the latter regime is absent in the
Hubbard model (See Fig. 3 in Ref. 16). Since t2 /t linearly increases with δ, BR becomes large (t2 /t goes away from 0) in
the overdoped regime. This behavior contrasts with that for
the Hubbard model, where t2 /t substantially vanishes for any
δ. For t′ /t = 0.3 (not shown), more marked BR exists in t1 /t(−0.1 ≲ t1 /t ≲ 0.1) for δ ≲ 0.2, and t2 /t ∼ −0.05 for most δ.
Thus, the renormalization of t1 /t and t2 /t becomes appreciable in Ψd for the t-J model.
In Table III, we show the energy improvement owing to
BRE for the t-J model:

where x̂ (ŷ) denotes the lattice vector in x (y) direction, δ(τ, τ′ )indicates the Kronecker delta, and ∆†τ (Ri ) is the creation operator of a nearest-neighbor singlet pair at site Ri ,√∆†τ (Ri ) = (c†i↑ c†i+τ↓ + c†i+τ↑ c†i↓ )/ 2.
(25)
To evaluate the contribution of AF correlation, we display theq = Q element of the spin structure factor,
1 ∑ iq·(Ri −R j ) ⟨ z z ⟩S (q) =e
SiS j .
(26)
Ns i jFor studying the electronic state, we use the momentum distribution function,1∑ †n(k) =⟨c ckσ ⟩,(27)
2 σ kσin which available k points are shifted by π/L in y direction
owing to the boundary conditions used.

∆EΛBR = EΛ (no BR) − EΛ ,3. Results of Pure States(28)

where EΛ (no BR) is the energy per site of the state Λ [Λ =AF or d] calculated without BR, namely, using the band parameters tηΛ (and tp for AF) fixed at the values in H. Notablefeatures are found in ∆EΛBR : (i) Both for the d-SC and AF

In Sects. 3.1 and 3.2, we discuss noteworthy aspects of thepure d-SC state Ψd and pure AF state ΨAF , respectively.

4

J. Phys. Soc. Jpn.

FULL PAPERS

0.025 t′/t L

8d-SC J/t=0.3-0.4 12-0.3 120.02 -0.3 16-0.2 12-0.1 12120.015 0.0
0.0 16
0.1 120.2 120.01 0.3 12
0.3 16
0.4 120.012

d-SC δ=0.0833 L=12

7

0.01

6

0.005
00

0.05

0.1

0.15

0.2

δ

0.25

t′/t

0.004

2

-0.3
-0.1
0.1
0.3

0.002

10

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.80t′/tL=12 δ=0.1389

J/t

-0.4
-0.3
-0.2
-0.1
0.0
0.1
0.2
0.3
0.4

0.8

n(k)

0.006
4
3

0.3

Fig. 3. (Color online) Doping rate dependence of d-wave SC correlation
function Pd calculated using Ψd for various values of t′ /t (and L).
1

5

Pd

S(π,π)

Pd

0.008

0.6

Fig. 5. (Color online) J/t dependence of the spin correlation function at
q = (π, π) (green) and the d-SC correlation function (blue) calculated using
Ψd are compared for four t′ /t. The figures of the symbols are common to
S (q) and Pd .

0.4

0.4
Pure states

Nor
d-SC

0
(0,0)

L=10
AF d-SC

L=12
AF d-SC

0.3
(π,0)

k

(π,π)

(0,0)

δ

0.2

Fig. 4. (Color online) The momentum distribution function of the pure dSC state is drawn along the path of k: (0, 0)-(π, 0)-(π, π)-(0, 0). Results forvarious values of t′ /t are drawn in blue. For comparison, corresponding results of the normal (paramagnetic) state are added in green. J/t = 0.3.

0.2

d-SC

0.1

0

states, BRE is almost ineﬀective for t′ /t = 0. (ii) Energy improvement in the d-SC state is small for any t′ /t and δ, despite
the above large BR in t1 /t and t2 /t. (iii) Energy improvementin the AF state is one or two orders of magnitude larger thanin the d-SC state for t′ /t = ±0.3. These features of the t-J
model is common to those of the Hubbard model.16) Owing
to the feature (i), the results of previous studies for t′ /t = 0
in Table I are almost identical despite whether BRE is introduced or not. The feature (ii) indicates that various propertiesof Ψd are almost unchanging by BRE. The feature (iii) will be
discussed later in Sect. 3.2. These features are identical withthose of the Hubbard model.
We notice in Table III that δ dependence of ∆EdBR is opposite between the t-J and Hubbard models for t′ /t = −0.3, although its magnitude is not large. This feature stems from theabove contrastive behavior of t1 and t2 , and is probably related
to the prevalence of d-SC order for large negative t′ /t in thet-J model.21–25) This prevalence is confirmed by Pd [Eq. (24)]shown in Fig. 3 and is hardly changed by BRE; the range ofd-SC seems very wide (probably up to δ = 0.5), as comparedto the Hubbard model.11) The dome-like δ dependence of dSC correlation has been well-known since early studies of thet-J model.6)
Shown in Fig. 4 is n(k) for various values of t′ /t. As a feature of d x2 −y2 -wave symmetry, a Fermi point exists in the nodal
direction near (π/2, π/2),9) whereas a gap opens near the antinodal (π, 0) in contrast with the case of ΨN . Note that, as t′ /tis varied, the behavior near the antinodal markedly changes,
but is almost unchanging for the other k. Because the property of d-SC appreciably changes with t′ /t as shown in Fig. 3,

AF

-0.4

-0.2

AF

0

0.2

0.4

t' / tFig. 6. (Color online) t′ -δ phase diagram constructed by comparing energies of pure AF and d-SC states. For large |t′ /t|, the AF domain is likely to
expand to the large-δ direction for technical reasons.48) The bold gray dashedline indicate a similar AF-d-SC boundary determined using the data without
BR for the Hubbard model of U/t = 12.11, 16) The diﬀerence of models does
not aﬀect the aspect for large |t′ /t|.

the electrons with k near antinodal must play a leading rolefor d-SC.16) This feature is common to that of the Hubbard
model.11)
We turn to the relation between the d-SC order and AF
correlation. In Fig. 5, we plot Pd and S (q) [Eq. (26)] with
q = Q = (π, π) as a function of J/t. As the spin exchangeinteraction J/t increases, Pd increases faithfully according to
S (Q) for any t′ /t, meaning that the driving force of d-SC is
the AF exchange interaction. Therefore, the scattering process(or a vector connecting FSs) of Q comes to play a crucial rolefor raising d-SC.17, 47)
These properties are consistent with the results of previous
studies. We will discuss the stability against PS of Ψd together
with the case of ΨAF in Sect. 3.2.
3.2 Antiferromagnetic state
As mentioned in Sect. 3.1 [feature (iii)], BRE is highly effective in reducing the energy of the AF state for large |t′ /t|.
This energy reduction is primarily caused by the FS renormalization (choice of {k}occ ) but not by the band form εAF
k5

J. Phys. Soc. Jpn.

FULL PAPERS1δ

AF L=16

-0.5

-0.55
-0.6

E/t

δ

L

0.0556 12
0.1111 12
0.1200 10m-0.65-0.4 -0.3 -0.2 -0.1

0

(a)
0
1

n(k)

0.7

0.6

0.00000.02780.05560.08330.11110.13890.16670.19440.22220 0.05 0.1 0.15 0.2 0.25 0.3 0.2500
0.2778
δ

0.5

0.2

(b)

(0,0)

(π,0)

k

(π,π)

(0,0)

Fig. 9. (Color online) Doping-rate dependence of the momentum distribution function of ΨAF for J/t = 0.3. In (a), we plot the results of typical t′ /tvalues of type-I (red) and type-II (blue) AF states for three values of doping
rate. The red (blue) arrow indicates the locus of pocket Fermi surface of typeI (type-II) AF. In (b), we plot the results for t′ /t = 0.3 (type-I) in a wide rangeof doping rate up to the PM regime δ > δAF (∼ 0.20, pale color). The arrow
indicates a Fermi surface appearing for the PM state. Shown in the inset is δ
dependence of the staggered magnetization of the same systems.

0.2

(π,π)

0.6

1
0.8
0.6
0.4
0.2
0

(0,0)

δ

0.3

0.4

k

(π,π)

0.4

-0.4
-0.3
-0.2
-0.1
0.1
0.2
0.3
0.4
0.5

0.8

n(k)

(π,0)

AF L=12 t′/t=0.3

0.8

t′/t
AF L=12 δ=0.0833 -0.5

(π,0)

(0,0)

0.9

Fig. 7. (Color online) t′ /t dependence of total energy per site (green, leftaxis) and staggered magnetization [Eq. (23)] (red, right axis) of ΨAF . Threecases of δ (and L) are shown.

0
(0,0)

0.06250.10940.1563t′/t=−0.3

0.6

0.2

t′/t

1

t′/t=0.3

0.8

0.4

m

E/t

-0.45

0.75
0.7
0.65
0.6
0.55
0.5
0.45
0.4
0.35
0.1 0.2 0.3 0.40.3m

AF J/t=0.3

n(k)

-0.4

0.06250.10940.1563(0,0)

Fig. 8. (Color online) The momentum distribution function calculated using metallic (δ > 0) ΨAF is shown along the same path as in Fig. 4. The
cases of various t′ /t are compared. The data of type-I (II) AF state are shownwith red (blue) symbols. The Lifshitz point is tL′ /t ∼ 0. The red (blue) arrow
indicates the loci of pocket Fermi surfaces of the type-I (type-II) AF state.

this insulating AF state changes over to metallic with FSs.
In Fig. 7, we show t′ /t dependence of total energy and staggered magnetization for three doping rates (and L). We findcusps (or small discontinuities in m) in both quantities att′ = tL′ ∼ 0, suggesting some transition. If this anomaly indicates the aforesaid Lifshitz transition, the loci of the FSs
should be switched when t′ passes through tL′ . To confirm
it, we show, in Fig. 8, the momentum distribution function[Eq. (27)] at δ = 0.083 for various values of t′ /t. We find apocket FS (discontinuity) around k = (π, 0) [(π/2, π/2)] fort′ > tL′ [t′ < tL′ ]; the loci of the pocket FS suddenly changes att′ = tL′ . This behavior is basically the same as what was foundfor a strongly correlated Hubbard model.16) We named theformer (latter) state the type-I (type-II) AF state. Related phenomena or behavior in the t-J model as well as the Hubbard
model have been studied since the early days of cuprates51, 52)
and later in various studies53) in connection with the diﬀerence in ARPES spectra between electron-doped54) and holedoped55) cuprates. Furthermore, this diﬀerence greatly aﬀects
whether AF and d-SC orders are coexistent or exclusive.16)
We confirmed that, even if J/t is varied, the above t′ /t dependence of n(k) is preserved as far as J < JPS , where JPS(∼ 2.5t for δ ∼ 0.08 – 0.20) is the transition point to PS (notshown).56) In Fig. 9, we show the δ dependence of n(k). In(a), we find that the loci of pocket FS are unchanging when
δ is varied, regardless of the type of AF, I (t′ /t = 0.3) or II(t′ /t = −0.3). Shown in (b) are the data (t′ /t = 0.3, type-I)
for a wider range of δ up to ∼ 0.28 in the PM area. The lo-itself, as explained in Appendix B. The energy improvementBRowing to BRE ∆EAF
[Eq. (28)] is one or two orders of magnitude larger than ∆EdBR , as shown in Table III. Thereby, thet′ -δ phase diagram constructed within the the pure states isnoticeably changed, as shown in Fig. 6. The bold gray dashedline indicates the boundary between AF and d-SC in the casewithout BRE for the corresponding Hubbard case. The AF
range shrinks and vanishes for large negative t′ /t. In contrast,
by introducing BRE, the AF range expands as |t′ /t| increases
especially in the negative-t′ /t side.48) Incidentally, it is known
for the t-J model with t′ /t = 0 that the stable state is rapidly
switched from AF to d-SC on carrier doping even without
BRE.8) This feature is unchanging by BRE and contrasts with
the Hubbard case (gray line). Anyway, this result with BRE iscontrary to the common knowledge that the AF state is rapidly
destabilized on doping and the d-SC order appears for holedoped cases (t′ /t < 0).
Next, we study the (expected) Lifshitz transition in themetallic (doped) AF states, which corresponds to what was
found in the Hubbard model.16) At half filling, an insulating AF state is realized, because the t-J model is reduced to
the Heisenberg antiferromagnet.49, 50) This AF is preserved atleast for |t′ /t| ≤ 0.5, because the nesting condition is completely retrieved owing to BRE. When carriers are doped,
6

J. Phys. Soc. Jpn.

FULL PAPERS

L
Type-I AF 10
12Type-II AF 10
12

0.35
0.3

Para I
Para II

L
10
12
10
12

3

J/t=0.3

2.5

AF L=10
AF L=12d-SC L=10
d-SC L=12
Hubbard AF L=10
Hubbard AF L=12

2
1.5

0.2

1/χc

δ

0.25

0.15

1
0.5
0

0.1

-0.5

0.05

-1
0

-0.4

-0.2

0

t′/t

0.2

0.4

-1.5
-0.4 -0.3 -0.2 -0.1

Fig. 10. (Color online) t′ -δ diagram according to Fermi-surface topology(occupied k-points) in the optimized metallic AF (PM) states. Red circles(triangles) indicate the type-I AF for L = 10 (L = 12), and similarly blue
symbols the type-II AF. The symbols in pale colors indicate the PM cases(m = 0), in which the occupied k-points are the same as those in the AF
states of corresponding types (and L). The boundary between blue and red(dark and pale colors) indicates tL /t (δAF ).

cus of the pocket FS is preserved until the AF order vanishesat δ = δAF ∼ 0.2 (see m in the inset). The pocket FS of the
AF state changes to a global FS centered at k = (0, 0) of the
PM state at δ = δAF , but {k}occ is identical between ΦAF andΦN . The evolution of n(k) according to δ is similar in type-II
cases.
As a summary of FS topology, we construct a diagram of
the types of AF order in ΨAF and of {k}occ in ΨN in the t′ -δspace (Fig. 10). The Lifshitz point tL′ /t (the boundary between
red and blue) slightly shifts to a smaller value as δ increases.
Incidentally, in the Hubbard model, the Lifshitz point is situated at −0.1 < tL′ /t < 0 for any model-parameter set.16) Notethat the AF order is stable in the whole underdoped regime(δ ≲ 0.16) for any t′ /t.
Finally, we discuss the intrinsic stability against PS of ΨAF .
Following the previous study,16) we judge this property by thesign of charge susceptibility χc ,AF J/t=0.3

(a)

-0.1

EJ/t

-0.2
-0.3
-0.4
t′/t L Ehop/t-0.3 1012
0.0 1012
0.3 1012

Ehop/t

-0.5
-0.6

EJ/t

-0.7

-0.3 1012
0.0 1012
0.3 10
12

-0.8
-0.9

0

0

0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

δd-SC J/t=0.3

(b)

-0.1
-0.2

EJ/t

(29)

For χc > 0 (χc < 0), the state is stable against (unstable toward) PS. We found that δ dependence of the energy of ΨAF
is fitted well by the parabolic form

-0.3
-0.4
t′/t L Ehop/t-0.3 1012
0.0 1012
0.3 1012-0.5Ehop/tE(δ) ≃ c0 + c1 δ + c2 δ2 ,0.1 0.2 0.3 0.4

Fig. 11. (Color online) Inverse charge susceptibilities of pure AF (red) andd-SC (blue) states as functions of t′ /t. For comparison, the corresponding
data of the AF state for the Hubbard model with U/t = 1216) (gray) areadded. The horizontal dash-dotted line (1/χc = 0) is the boundary of stability
against PS.

01∂2 E(δ) E(δ + ∆δ) + E(δ − ∆δ) − 2E(δ)
=
=
.
χc
∂δ2
(∆δ)2

0

t′/t

-0.6

(30)

for δ < δAF (except in the δ → 0 limit). Therefore, we have aunique value χc = c−1
2 in the AF phase. In Fig. 11, the valuesof 1/χc thus estimated are plotted in red as a function of t′ /tfor L = 10 and 12. χc of ΨAF for J/t = 0.3 is positive for any
t′ /t, in contrast with that of the Hubbard model with U/t = 12
(gray symbols), which becomes negative around t′ /t = 0. The
AF state is always stable against PS for J/t = 0.3.
We mention J/t dependence of this property. In Fig. 12(a),we show δ dependence of the two energy elements Ehop =⟨Hhop ⟩/Ns and E J = ⟨H J ⟩/Ns of ΨAF for J/t = 0.3. For any
t′ /t, Ehop /t is downward-convex and E J /t is concave; the con-

EJ/t

-0.7

-0.3 1012
0.0 1012
0.3 10
12

-0.8
-0.9

0

0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

δFig. 12. (Color online) Doping-rate dependence of two energy elements
per site Ehop and E J with respect to pure (a) AF and (b) d-SC states for three
typical values of t′ /t.

7

J. Phys. Soc. Jpn.

FULL PAPERS

0.008

t′/t L

0.006
0.004
0.002

∆E/t

J/t=0.3 t′/t=0.2 δ=0.0833 L=24
AF
d-SC
(a)
(b)

-0.3 1012-0.1 1012
0.1 1012
0.2 10
12

Mixed
(c)

ky

0

-0.002-3-0.004-0.008-0.012 0

-1

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

3

Fig. 15. (Color online) Contour maps of momentum distribution functionn(k) of (a) ΨAF , (b) Ψd , and (c) Ψmix for (t′ /t, δ) = (0.2, ∼ 0.08), a typicalcase of category (1) [EAF (type-I) < Ed ], in the first Brillouin zone. The
arrow in (a) indicate the Q vector connecting two pocket FSs.

-0.006-0.01-2Ed-SC − EAF
0.05

0.1

0.15

δ

0.2

0.25

clusivity of the two orders in Ψmix . Let us focus on the casesof δ = 0.08, because all the above categories (1)–(4) appear.
First, we check whether the two orders coexist or not ineach category. In Fig. 14(a), we depict δ dependence of Pd andm calculated using the pure state Ψd and ΨAF , respectively.
Shown in Fig. 14(b) are the same quantities simultaneously
calculated using the mixed state Ψmix . In (a), both Pd and m at
δ = 0.08 (marked with a ring) are finite for any t′ /t, indicatingthat each order can arise as a single order for any t′ /t if thecounter order does not arise. And the values of both Pd andm are almost constant with respect to t′ /t. In (b), however, Pd
vanishes for t′ /t = −0.3 with m almost unchanging, whereas
Pd is almost unchanging for other values of t′ /t but m tends tobe suppressed, especially, for t′ /t = −0.1. Thus, we find thatin the regime of type-II AF [t′ /t ≲ 0, categories (2) and (4)],the AF and d-SC orders tend to exclude each other, whereas
in the regime of type-I AF [t′ /t ≳ 0, categories (1) and (3)],the two orders are likely to coexist. This is the same tendency
as that of the Hubbard model.16)
Before studying the features of each category, we introduce
common terminologies and an empirical law. If EA < EB (A,B=AF or d-SC) for the pure states, we call the order A (B)the leading (subordinate) order. In all categories (1)–(4), theleading order necessarily arises in the mixed state Ψmix , if theleading order arises in the pure state with the same model
parameter set (t′ /t, δ).
In category (1) [e.g., (t′ /t, δ) = (0.2, 0.08)], the leading andsubordinate orders are the type-I AF and d-SC, respectively.
According to the above law, the AF order of type-I arises inΨmix ; the problem is whether the d-SC order can simultaneously arise or not. Actually, the d-SC order coexists with the
AF order in Ψmix as in Fig. 14(b). Let us consider the mechanism of the coexistence. Figure 15 shows the contour mapsof n(k) for the three kinds of Ψ. As described in Sect. 3.2,
the pure AF state (left panel) has pocket FSs near antinodal[(±π, 0), (0, ±π)] and exhibits a gap in the other region of k.
In the pure d-SC state (middle panel), a SC gap opens except
in the nodal directions (ky = ±k x ), where Fermi points appear
near (±π/2, ±π/2). As mentioned in Sect. 3.1 and in the previous paper,16) the occurence of d-SC as the subordinate order
requires the FSs which are connected with the scattering vector Q = (π, π) in the leading-order (AF) state, especially, nearantinodal owing to the largest d x2 −y2 -wave gap and density ofstate. In category (1), the pocket FSs in Fig. 15(a) meet thisrequirement. Consequently, the d-SC order arises in Ψmix andthe d-SC gap opens at the loci of pocket FS near antinodal.

0.3

Fig. 13. (Color online) Energy diﬀerence between the pure AF and SC
states [Eq. (31)] is depicted as a function of doping rate for severa values oft′ /t. The zero (black line) is the boundary whether the lower energy is given
by ΨAF or Ψd . For positive (negative) ∆E/t, the AF (d-SC) state is morestable. J/t = 0.3. The black ring indicates δ = 0.08 used for explanation.

vexity of total energy, namely, stability against PS stems fromthe hopping energy. Becuase E J is the factor of instability toward PS, ΨAF will become unstable as J/t increases and E J
becomes predominant. We would like to address this subject
again elsewhere.
At this opportunity, we discuss the stability against PS ofother states. In Fig. 11, we also plot 1/χc of Ψd in blue. Because χc is always positive, Ψd is stable against PS in therange shown. This property is the same as that in the Hubbard
model with U/t = 12.16) In Fig. 12(b), we show δ dependenceof the two energy elements. Since both Ehop and E J exhibita tendency similar to that of ΨAF in Fig. 12(a), Ψd will also
phase separates at a large J/t. For the mixed state, we estimate
1/χc for each δ for J/t = 0.3 using the expression of finite
diﬀerences in Eq. (29), because the quadratic fit of E(δ) is deteriorated by the existence of subdivided domains [(1)–(4)],as will be discussed in Sect. 4. Anyway, we find that 1/χc > 0
holds for any δ (not shown), so that Ψmix is also stable againstPS.
4. Mixed State of AF and d-SC Orders
In this section, we study the interplay of AF and d-SC orders in the mixed state Ψmix . We start with the energy diﬀerence between the two pure states ΨAF and Ψd :∆E = Ed − EAF .

(31)

In Fig. 13, δ dependence of ∆E/t for J/t = 0.3 is shown forvarious values of t′ /t (L = 10 and 12). Note that there exist
cases in which d-SC is more stable than AF (∆E/t > 0) in the< 0.16), in contrast with the Hubbard
underdoped regime (δ ∼model with U/t = 12.16) We can divide the δ-t′ space into fourdomains (categories) according to(1) EAF (type-I) < Ed
(2) EAF (type-II) < Ed
(3) Ed < EAF (type-I)(4) Ed < EAF (type-II)when the AF and/or d-SC orders arise. This classification isconvenient to understand the behavior of coexistence or ex8

J. Phys. Soc. Jpn.

FULL PAPERS

1
m

0.025

(b) Mixed J/t=0.3
Pd t′/t L
0.02-0.3 12-0.1 120.1 120.2 12
0.015

0.8
0.6

m
Pd

Pd

(a) AF and d-SC J/t=0.3
Pd t′/t L
0.02-0.3 12-0.1 120.1 120.2 12
0.015

1
m

0.8
0.6

m

0.025

0.01

0.4

0.01

0.4

0.005

0.2

0.005

0.2

00

0.05 0.1 0.15 0.2 0.25 0.30

00

δ

0.05 0.1 0.15 0.2 0.25 0.30

δFig. 14. (Color online) Doping-rate dependence of order parameters of d-SC (Pd , left axis) and AF (m, right axis). In (a), Pd (m) is calculated using the purestate Ψd (ΨAF ). In (b), both Pd and m are calculated simultaneously using the mixed state Ψmix . The rings indicate δ = 0.08 for explanation.

J/t=0.3 t′/t=-0.3 δ=0.0833 L=24
AF
d-SC
(a)
(b)

J/t=0.3 t′/t=-0.1 δ=0.0833 L=24
AF
d-SC
(a)
(b)

Mixed
(c)

Mixed
(c)

ky

-3

-2

-1

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

ky

3

-3

Fig. 16. (Color online) Contour maps similar to Fig. 15 but for (t′ /t, δ) =(−0.3, ∼ 0.08), a typical case of category (2) [EAF (type-II) < Ed ].

J/t=0.3 t′/t=0.1 δ=0.0833 L=24
AF
d-SC
(a)
(b)

-2

-1

0

kx

1

2

3 -3

-2

-1

0

1

2

3 -3

-2

-1

0

kx

1

2

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

3 -3

-2

-1

0

kx

1

2

3

are ready to coexist in Ψmix by using diﬀerent parts of n(k),resulting in a full gap as shown in Fig. 17(c). Because theleading order is d-SC in category (3), the subordinate order
AF is weakened to some extent. This tendency contrasts withthat of category (1).
In category (4), the leading order is d-SC, and the FS inΨAF is type-II. The situation of n(k) in ΨAF and Ψd [Fig. 18]
is similar to that of category (2) [Fig. 16]. Therefore, the twoorders tends to exclude each other, and the subordinate AForder is weakened (δ = 0.08) or removed for high doping rates
in Ψmix as in Fig. 14(b). Consequently, n(k) of Ψmix becomesakin to n(k) of Ψd as seen in Figs. 18(b) and 18(c). A point
diﬀerent from category (2) is that the subordinate order (AF)tends to be excluded not completely in Ψmix . This is probably
because the nesting condition to be satisfied for the AF orderis less strict than the condition for d-SC.
The above aspect of interplay between the AF and d-SC
orders is more clearly revealed in J/t dependence of Pd andm. In the following, we focus on the cases of δ = 0.08. We
start with the type-I regime [t′ /t ≳ 0, categories (1) and (3)].
Shown in Figs. 19(a) is J/t dependence of the energy diﬀerence between the two pure states ∆E [Eq. (31)]. As J/t increases, ∆E/t becomes negative at J0 /t (∼ 0.2 for t′ /t = 0.1,
∼ 0.55 for t′ /t = 0.2), where the stable state switches fromAF to d-SC. Shown in Fig. 19(b) is J/t dependence of Pd andm simultaneously measured using Ψmix . As J/t increases fromzero, both Pd and m rapidly increase and exhibit no anomaly
near J = J0 ; some inflection at J/t ∼ 0.15 not only fort′ /t = 0.1 but for 0.2 probably reflects the rapid variation of
the states. Note that the subordinate order never vanishes bothfor J < J0 (Pd ) and J > J0 (m). Thus, in the type-I regime,

Mixed
(c)

kx

-1

Fig. 18. (Color online) Contour maps similar to Fig. 15 but for (t′ /t, δ) =(−0.1, ∼ 0.08), a typical case of category (4) [Ed < EAF (type-II)].

ky

-3

-2

3

Fig. 17. (Color online) Contour maps similar to Fig. 15 but for (t′ /t, δ) =(0.1, ∼ 0.08), a typical case of category (3) [Ed < EAF (type-I)].

The resultant state Ψmix comes to have both orders and becomes a fully gapped state, as shown in Fig. 15(c).
In category (2) [e.g., (t′ /t, δ) = (−0.3, 0.08)], because theleading order is AF, the type-II AF order arises in Ψmix . In
this case, the loci of FSs in ΨAF are in the nodal directions
as shown in Fig. 16(a); a scattering vector Q to create thed-SC order cannot be placed between two antinodals. As aresult, the d-SC order does not arise, and the resultant stateΨmix remains the pure AF state [Fig. 16(c)]. Namely, the dSC order is excluded. The mechanisms in categories (1) and(2) are basically the same as those for the strongly correlated
Hubbard model.16) Comparing with Ref. 25 (See Tables I andII), we notice that the independence of tηAF and tηSC is crucial
for this exclusivity.
In category (3) [e.g., (t′ /t, δ) = (0.1, 0.08)], the d-SC order necessarily appears, because it is the leading order. In thiscase, the situation of FSs of ΨAF and Ψd [Figs. 17(a) and17(b)] is similar to category (1). Therefore, the two orders
9J. Phys. Soc. Jpn.

Ed-SC − EAF

0.005

0.01 (a)

t′/t
0.2
0.1

0

∆E/t

∆E/t

0

-0.005-0.005(1) (3)(2) (4)
-0.01

-0.01

(2) (4)(1) (3)

δ=0.0833

Mixed-0.015

L=12

TypeI

0.012

0.5

Pd

0.014

0.6

0.01

Ed-SC − EAF
δ=0.0833-0.02
0.016 (b)

0.4

Mixed

L=12

TypeII

0.014

0.6

0.012

0.5

Pd

-0.015-0.02
0.016 (b)

t′/t
-0.3
-0.1

0.005

0.01

0.4

0.008

m

0.008

0.3

0.3
0.006

m

0.01 (a)FULL PAPERS

0.006

t′/t
Pd 0.2-0.3
m-0.3
Pd 0.1
-0.1
0.002
m
-0.1
0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.90

t′/t
Pd 0.20.2
m0.2
Pd 0.1
0.1
0.002
m
0.1
0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.90

0.004

0.004

J/t

J/t

Fig. 19. (Color online) (a) The energy diﬀerence between the two purestates [Eq. (31)] and (b) the order parameters calculated using Ψmix are plotted as functions of J/t for typical values of t′ /t in the type-I AF regime. As
guides, the points where ∆E/t crosses zero (J0 /t) are indicated by the graypentagons and open triangles for t′ /t = 0.1 and 0.2, respectively. The regimesof J < J0 and J > J0 belong to categories (1) and (3), respectively, as shownby arrows for each t′ /t.

Fig. 20. (Color online) The same quantities as in Fig. 19 are plotted but for
the values of t′ /t in the type-II AF regime: t′ /t = −0.1 and −0.3. The regimesof J < J0 and J > J0 belong to categories (2) and (4), respectively.

Para

d-SC
the two orders tend to coexist.
In Fig. 20, ∆E [in (a)], Pd and m [in (b)] are similarly shownfor the type-II regime (t′ /t ≲ 0). In contrast with the type-I
cases, the behavior of both Pd and m anomalously changes
near J = J0 . For J ≲ J0 [category (2)], the subordinate dSC order is completely suppressed by the leading AF order.
For J ≳ J0 , [category (4)] the leading d-SC order suddenly
increases and the subordinate AF order drops near J = J0 . In
this regime, substantial magnitude of m is preserved in ΨAF .
Anyway, in the type-II regime, the two orders tend to exclude
each other.
As a summary, we construct a t′ -δ phase diagram of AF
and d-SC orders for J/t = 0.3 using the mixed state Ψmix
(Fig. 21). Compared to the corresponding phase diagram of
the Hubbard (t-t′ -U) model for U/t = 12 (Fig. 27 in Ref. 16),
the aspect of AF is quite similar, but the aspect of d-SC is different in that its area greatly expands from a part of type-II
area to a large range of type-I area. And that the d-SC order coexists with the AF order even in the type-II area. Another characteristic point is that there is no area where a stateis unstable toward PS. The behavior of χc > 0 in Ψmix (cf.
Sect. 3.2) coincides with previous results.25, 32)
When we consider the correspondence with the Hubbard
model with U/t = 12, the value J/t = 0.3 seems too large inthat the d-SC order excessively appears, but too small in that

AF
Type-II

AF
Type-I

Fig. 21. Phase diagram of AF and d-SC orders in t′ -δ space constructedusing Ψmix for J/t = 0.3. The AF order exists under the bold red curve(guide line), and the d-SC order appears in the area between the two boldblue curves. Solid (open) symbols indicate the boundary points determinedusing the systems of L = 12 (10). In the whole area, all states are stableagainst PS.

PS does not occur. The correspondence is not simple. This ispossibly because the three-site (pair-hopping) and J ′ terms are
disregarded in the model. Anyway, in the context of cuprates,
the main problem remains—why the robust AF order survives
up to high doping rates in theory.

10

J. Phys. Soc. Jpn.

FULL PAPERS-0.3 AF

5. Summary

J/t=0.3 L=12t′/t

δ=0.0278

As an extension of previous work for the Hubbard
model,16) we studied band-renormalization eﬀects (BRE) on
the interplay between AF and d-SC orders in the t-t′ -J model,
which is expected to more favor the d-SC order than the Hubbard model. To reliably treat strong correlation of the t-J
model, a VMC method was used. We mainly discussed thecase of J/t = 0.3 realistic values for cuprate SC. We summarize the main results below.
(i) BR modifies the properties of d-SC only slightly but
greatly stabilizes the AF state, especially, for large valuesof |t′ /t| (Table III). Consequently, the t′ -δ phase diagram islargely modified; the AF order prevails in almost whole underdoped range (Fig. 21), similarly to the Hubbard case.
(ii) The metallic AF order for δ > 0 in ΨAF is classified
into two types according to t′ > tL′ (type-I) or t′ < tL′ (typeII), where tL′ /t (∼ 0) is a Lifshitz transition point (Fig. 10).
As t′ /t decreases, the loci of the pocket Fermi surface switch
at t′ = tL′ from antinodal [∼ (π, 0)] to the nodal directions[∼ (π/2, π/2)]. This distinction of AF type plays a crucial rolefor coexistence/exclusivity of the AF and d-SC orders in Ψmix .
The two orders tend to coexist (exclude each other) in thearea of type-I (II) AF for the same reason as in the Hubbard
model,16) in short, compatibility of the electronic states.
(iii) In contrast with the Hubbard model with U/t = 12,
the pure d-SC state becomes more stable than the pure AF
state in the underdoped regime (δ ≲ 0.16) for small valuesof |t′ /t|, namely, d-SC becomes the leading order. As a result,
it is convenient for discussing coexistence/exclusivity that thet′ -δ space is divided into four domains (categories) according
to whether or not EAF < Ed in addition to whether or nott′ < tL .The area of coexistence comes to appear in the regimeof type-II AF (−0.2 ≲ t′ /t ≲ 0), besides widely in the type-I
regime (t′ /t > 0) (Fig. 21).
(iv) For J/t = 0.3, any state considered here is stableagainst PS (Fig. 11), in contrast with for the Hubbard modelwith U/t = 12, where states with the AF order is unstabletoward PS near t′ /t = 0. However, from analysis of energy
components, both AF and d-SC states probably become unstable toward PS for larger values of J/t.
On the basis of (i) and (iii), we repeat that the problem ofcoexistence/exclusivity largely depends on the value of t′ /tand BRE should be properly introduced for t′ /t ∼ ±0.3 (Table I). Regarding (iv), PS in the t-J model is a long-standing
problem. We will reconsider it in another publication. We concentrated on the interplay between the AF and d-SC orders; it
is important to study interplay among other low-energy statessuch as staggered-flux3, 57, 58) and striped states.59)
The tendency toward predominant AF long-range orders
has been reported not only for the single-band models discussed above but also for the d-p model.60, 61) In addition tothis predominant AF, some results are inconsistent with the
behavior of cuprates: For instance, except for the multilayered systems,62) the AF and d-SC long-range orders do notcoexist and the AF states are always insulating. These points
suggest that the uniform models are possibly insuﬃcient todescribe cuprate SCs. One possibility for reconciling the inconsistency with experiments is that some disorders such as
the impurity potential of carrier dopants, which is inherent incuprates, destabilize the long-range AF order but aﬀect the

-0.35

-0.4
-0.3
-0.2
-0.1
0.0
0.1
0.2
0.3

-0.35

-0.4

E/t

E/t

-0.4

δ=0.0833

-0.45
-0.45
-0.5
-0.5 -0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4 0.5

tp/t

Fig. A·1. (Color online) tp /t (parameter for determining {k}occ ) dependence of total energy per site for the AF state is shown for various valuesof a model parameter t′ /t and two doping rates δ. The system size is L = 12.
The variational parameters except for tp /t are optimized for each model parameter set. The ranges of tp /t that give the minima of E/t are indicated byfilled (for Type-I AF) and empty (for Type-II AF) arrows for each δ.

d-SC order only slightly. Inhomogeniety allows static shortrange AF orders. We will address such subjects elsewhere.63)
Acknowledgments
We thank Kenji Kobayashi, Yuta Toga, Tsutomu Watanabe,and Hidekazu Mukuda for useful discussions and information. This work is supported in part by Grants-in-Aid from the
Ministry of Education, Culture, Sports, Science and Technology, Japan, and JSPS KAKENHI Grant Number JP17J00670.
Appendix A:

Optimization of Occupied k-Points {k}occ
in Normal, AF, and Mixed States

In this Appendix, we describe the details of how we actually optimize ΨN , ΨAF and Ψmix for finite systems in thisstudy. The operation in band renormalization for finite-size
systems is composed of two elements:16) (i) Optimization of
the energy dispersion εk itself, and (ii) optimization of {k}occ(the set of discrete k-points occupied by electrons), which
is usually obtained by filling the k-points with electrons in
the order of small εk optimized in (i).64) These two elements
merge in the thermodynamic limit, but it is convenient to treat
the two elements independently for finite systems for technical reasons.
In optimizing the band parameters tηAF /t for finite systems,
E/t becomes discontinuous as a function of tηAF /t at specificvalues tη(i) /t (i = 1, 2, · · · ), owing to the discrete k-points.
Namely, {k}occ switches from one to another at tη(i) /t. Furthermore, εk sometimes becomes almost constant in the ranges
between two discontinuities, say tη(1) < tη < tη(2) . As a result, it becomes inconvenient to apply ordinary optimization
tools based on the derivative of εk , such as the quasi-Newton
method and the stochastic reconfiguration method.
To overcome this diﬃculty, we follow the next prescription in this study. In the process of (ii), we generate {k}occ
AFindependently of εAFk or tη /t and fix it at a certain configu11

J. Phys. Soc. Jpn.

FULL PAPERS

ration. Then, E/t becomes continuous as a function of every
variational parameter, and the ordinary optimization tools become applicable. Comparing E/t thus optimized for a certain{k}occ with those optimized for other {k}occ , we can determine
the optimized energy and wave function. In this process, it is
not realistic that one checks all possible {k}occ one by one, because the number of {k}occ grows exponentially as the systemsize increases. For the present case, we found that almost all
optimized {k}occ for small L are included in those generatedby εk up to the diagonal hopping. Therefore, we consider only{k}occ generated by()(p)εk = −2t cos k x + cos ky − 4tp cos k x cos ky ,(A·1)

tp/t t1/t t2/t t3/t t4/t

AF J/t=0.3 L=12

0.2
0.1
0

-0.1
tp
t1
t2
t3
t4

-0.2

in a proper range of tp /t, instead of εAF
in Eq. (17) or
k
√
AF
AF 2
2
Ek = U/2 − (εk ) + ∆AF . Here, tp /t is a kind of variational parameter that optimizes {k}occ and is independent ofεAF
k . Thereby, the number of {k}occ to be checked is greatly
reduced. In Fig. A·1, we show the tp /t dependence of variational energy per site for the AF state (L = 12), where all theresidual variational parameters including tηAF /t are optimized.
In the range of tp /t corresponding to a certain {k}occ , E/t isconstant. The optimized {k}occ is determined by finding therange with the lowest E/t. Shown in this figure are the datafor δ ∼ 0.03 and 0.08, and the model parameter t′ /t = −0.4
– +0.3 for each δ. We find for δ = 0.083 that the optimized{k}occ is generated using tp /t = −0.15 [0] (or nearby tp /t in thesame range) for t′ /t < 0 (type-II AF) [t′ /t > 0 (type-I AF)] as
indicated with an open [a solid] arrow.
We apply this scheme to all cases of L = 10 and 12 forΨN , ΨAF , and Ψmix . In some cases, however, it was found thatthe true optimized {k}occ is not generated within the above
scheme, especially, for considerably large |t′ /t| and δ. In suchcases, we search within several plausible {k}occ on the basisof the optimized data for small L. The above scheme, if thecorrect {k}occ is obtained, is being an optimization in a wider
parameter range in the sense that εk and {k}occ are independently optimized; actually, a parameter tp /t is added. However, it remains within the finite-size correction, because thetwo elements merge for L → ∞.
Appendix B:

0.3 δ=0.0833-0.3-0.4 -0.3 -0.2 -0.1 0 0.1 0.2 0.3 0.4

t′/t(Color online) t′ /t dependence of optimized band parameters in

Fig. B·1.
ΨAF . The energy for a fixed t′ /t becomes constant in certain ranges of tp as
explained in Appendix A; the optimized range of tp is shown in green.

(≡ {k}(ε)
occ ) determined using the optimized tη .
On the other hand for t′ /t < 0, in Fig. B·1, the optimized
range of tp /t is switched to a diﬀerent narrow one (∼ −0.14),
whereas the optimized tη are smoothly extended from the(p)′regime of t′ /t > 0. It follows that {k}occ , {k}(ε)
occ for t /t < 0.
As discussed in Sect. 3.2, some properties of ΨAF are criticallydiﬀerent between for t′ /t ≳ 0 (type I) and t′ /t ≲ 0 (type II).
Therefore, such properties are considered to be led by {k}occ
or FS and not by εAF
k itself. This predominance of {k}occ over
εAF
can
be
seen
in
energy reduction. In contrast to the case inkFig. B·1, the optimized tη in Fig. 13(b) in Ref 16 are switched
to diﬀerent values for t′ /t < 0. This is because tη are adjusted(p)
so as to satisfy {k}(ε)
occ = {k}occ . Namely, the optimization of{k}occ takes priority over that of εAF
k in energy minimization.
Thus, we can conclude that the reduction in EAF and some relevant properties are primarily caused by FS renormalization,and not by the band form εAF
k itself.

Choice of {k}occ as Leading Role of BRE1) J. G. Bednorz and K. A. Müller, Z. Phys. B 64, 189 (1986).
2) For the recent status of experimental research, see, for instance, articles
in the special issue of J. Phys. Soc. Jpn. 81 vol. 1 (2012) on “Recent
developments in superconductivity”.
3) P. A. Lee, N. Nagaosa, and X.-G. Wen, Rev. Mod. Phys. 78, 17 (2006).
4) M. Ogata and H. Fukuyama, Rep. Prog. Phys. 71, 036501 (2008).
5) H. Yokoyama and H. Shiba, J. Phys. Soc. Jpn. 56, 3570 (1987).
6) H. Yokoyama and H. Shiba, J. Phys. Soc. Jpn. 57, 2482 (1988).
7) C. Gros, Ann. Phys. (New York) 189, 53 (1989).
8) H. Yokoyama and M. Ogata, J. Phys. Soc. Jpn. 65, 3615 (1996).
9) A. Paramekanti, M. Randeria, and N. Trivedi, Phys. Rev. B 70, 054504(2004).
10) P. W. Anderson, P. A. Lee, M. Randeria, T. M. Rice, N. Trivedi, andF. C. Zhang, J. Phys. Condens. Matter 16, R755 (2004).
11) H. Yokoyama, M. Ogata, Y. Tanaka, K. Kobayashi, and H. Tsuchiura,
J. Phys. Soc. Jpn. 82, 014707 (2013).
12) Many other studies. See, for instance, literature in Refs. 3, 4, 10, and
11.
13) J. Otsuki, H. Hafermann, and A. I. Lichtenstein, Phys. Rev. B 90,
235132 (2014).
14) T. Misawa and M. Imada, Phys. Rev. B 90, 115137 (2014).
15) B.-X. Zheng and G. K.-L. Chan, Phys. Rev. B 93, 035126 (2016).
16) R. Sato and H. Yokoyama, J. Phys. Soc. Jpn. 85, 074701 (2016).

As discussed in Sects. 3.1 and 3.2, BRE plays a crucial role
in reducing EAF for large |t′ /t|. In this Appendix, we argue thatthe choice of {k}occ (not optimizing εAFk ) primarily contributes
to this energy reduction and other properties by analyzing theAFbehavior owing to {k}occ (or tp ) and to εAFk [or tη (= tη ), η =1–4], and by comparing with results of the previous study16)
in which {k}occ and εAF
k depend on each other.
Shown in Fig. B·1 is the t′ /t dependence of the optimized
band parameters tη and the parameter tp , which determines{k}occ , for δ = 0.0833. As described in Appendix A, tp is determined independently of tη in ΨAF . Corresponding results
for the Hubbard model are presented in Fig. 13(b) in Ref. 16,
in which, contrastively, {k}occ is determined according to εAFk(or tη ). Inversely speaking, available ranges of tη are regulatedby {k}occ . The diﬀerence of models is irrelevant here. Comparing these two results for t′ /t > 0, we find that the opitizedvalues of tη are virtually identical. In this regime of t′ /t, the(p)optimized {k}occ (≡ {k}occ ) based on tp is identical with {k}occ
12

J. Phys. Soc. Jpn.

FULL PAPERS

45) S. Sorella, Phys. Rev. B 64, 024512 (2001).
46) K. Kobayashi and H. Yokoyama, Physica C 463-465, 141 (2007);
L. F. Tocchio, F. Becca, and C. Gros, Phys. Rev. B 86, 035102 (2012).
47) D. J. Scalapino, Phys. Rev. 250, 329 (1995).
48) For large values of |t′ /t| and δ, the optimization of {k}occ is possibly
incomplete within the scheme of a single parameter tp described in Appendix A. The true EAF is possibly smaller. As a result, the AF domain
may somewhat expand to the large-δ direction for large |t′ /t|. In contrast, in the mixed state, we determine {k}occ using εAFk (four-parameter
scheme). Therefore, the choice of {k}occ is more correct.
49) N. Trivedi and D. M. Ceperley, Phys. Rev. B 41, 4552 (1990);
K. J. Runge, Phys. Rev. B 45, 12292 (1992).
50) S. Liang, B. Doucot, and P. W. Anderson, Phys. Rev. Lett. 61, 365(1988).
51) B. I. Shraiman and E. D. Siggia, Phys. Rev. Lett. 61, 467 (1988).
52) R. J. Gooding, K. J. E. Vos, and P. W. Leung, Phys. Rev. B 50, 12866(1994).
53) For instance, T. K. Lee, C.-M. Ho, and N. Nagaosa, Phys. Rev. Lett.
90, 067001 (2003); H. Kusunose and T. M. Rice, Phys. Rev. Lett.
91, 186407 (2003); D. Senechal, P.-L. Lavertu, M.-A. Marois, and A.M. S. Tremblay, Phys. Rev. Lett. 94, 156404 (2005).
54) N. P. Armitage, et al., Phys. Rev. Lett. 88, 257001 (2002); N. P. Armitage, P. Fournier, and R. L. Greene, Rev. Mod. Phys. 82, 2421 (2010).
55) For instance, F. Ronning, et al., Science 282, 2067 (1998).

17) P. W. Anderson, Science 235, 1196 (1987); F. C. Zhang and T. M. Rice,Phys. Rev. B 37, 3759(R) (1988).
18) A. B. Harris and R. V. Lange, Phys. Rev. 157, 295 (1967).
19) J. E. Hirsch, Phys. Rev. Lett. 54, 1317 (1985).
20) C. Gros, R. Joynt, and T. M. Rice, Phys. Rev. B 36, 381 (1987).
21) C. T. Shih, T. K. Lee, R. Eder, C.-Y. Mou, and Y. C. Chen, Phys. Rev.
Lett. 92, 227002 (2004).
22) L. Spanu, M. Lugas, F. Becca, and S. Sorella, Phys. Rev. B 77, 024510(2008).
23) T. Watanabe, H. Yokoyama, K. Shigeta, and M. Ogata, New J. Phys. 11,075011 (2009).
24) S. Pathak, V. B. Shenoy, M. Randeria, and N. Trivedi, Phys. Rev. Lett.
102, 017002 (2009).
25) K. Kobayashi, R. Sato, and H. Yokoyama, Phys. Proc. 81, 1 (2016).
26) W. L. McMillan, Phys. Rev. 138, A442 (1965); D. Ceperley,
G. V. Chester, K. H. Kalos, Phys. Rev. B 16, 3081 (1977); H. Yokoyama
and H. Shiba, J. Phys. Soc. Jpn. 56, 1490 (1987); C. J. Umrigar,
K. G. Wilson, and J. W. Wilkins, Phys. Rev. Lett. 60, 1719 (1988).
27) T. K. Lee and S. Feng, Phys. Rev. B 38, 11809 (1988).
28) T. Giamarchi and C. Lhuillier, Phys. Rev. B 43, 12943 (1991).
29) A. Himeda and M. Ogata, Phys. Rev. B 60, R9935 (1999).
30) D. A. Ivanov, Phys. Rev. B 70, 104503 (2004).
31) C. T. Shih, Y. C. Chen, C. P. Chou, and T. K. Lee, Phys. Rev. B 70,
220502 (2004).
32) M. Lugas, L. Spanu, F. Becca, and S. Sorella, Phys. Rev. B 74, 165122(2006).
33) K. Kobayashi and H. Yokoyama, JPS Conf. Proc. 3, 015012 (2014);
Phys. Proc. 58, 22 (2014).
34) M. Ogata and A. Himeda, J. Phys. Soc. Jpn. 72, 374 (2003).
35) L. F. Tocchio, F. Becca, and S. Sorella, Phys. Rev. B 94, 195126 (2016).
36) R. Sato and H. Yokoyama, JPSCP 871, 012014 (2017).
37) For instance, T. Tohyama and S. Maekawa, Supercond. Sci. Technol.
13, R17 (2000).
38) R. Raimondi, J. H. Jeﬀerson, and L. F. Feiner, Phys. Rev. B 53, 8774
(1996); L. F. Feiner, J. H. Jeﬀerson, and R. Raimondi, Phys. Rev. Lett.
76, 4939 (1996).
39) E. Pavarini, I. Dasgupta, T. Saha-Dasupta, O. Jepsen, and O. K. Anderson, Phys. Rev. Lett. 87, 047003 (2001).
40) R. Jastrow, Phys. Rev. 98, 1479 (1955).
41) M. C. Gutzwiller: Phys. Rev. Lett. 10, 159 (1963).
42) A. Himeda and M. Ogata: Phys. Rev. Lett. 85, 4345 (2000).
43) For t′ , 0 or δ > 0, the HF solution deviates from the present formof ΦAF . However, the present form is not inappropriate, because theoptimized band is renormalized so as to satisfy the nesting condition.
Actually, it is confirmed that the optimized energies of the present form
and of the HF solution obtained by numerically diagonalizing HMF coincide with each other in another context by Y. Toga (unpublished).
44) J. P. Bouchaud, A. Georges, and C. Lhuillier, J. Phys. (Paris) 49, 553(1988).

56) R. Sato and H. Yokoyama, to be published in JPSCP , (2018).
57) H. Yokoyama, S. Tamura, and M. Ogata, J. Phys. Soc. Jpn. 85, 124707(2016).
58) K. Kobayashi and H. Yokoyama, JPS Conf. Proc. 871, 012032 (2017);
K. Kobayashi, H. Yokoyama, and Y. Toga, to be published in JPS Conf.
Proc. (2018).
59) K. Kobayashi and H. Yokoyama, Low Temp. Phys. 117, 199 (1999);
A. Himeda, T. Kato and M. Ogata, Phys. Rev. Lett. 88, 117001 (2002);
P. Corboz, T. M. Rice, and M. Troyer, Phys. Rev. Lett. 113, 046402(2014); K. Ido, T. Ohgoe, and M. Imada, Phys. Rev. B 97, 045138(2018).
60) S. Tamura, Ph. D thesis, ‘Variational Monte Carlo studies of d-p
model’, (Tohoku University, 2016). [In Japanese]61) T. Yanagisawa, I. Hase, M. Miyazaki, and K. Yamaji, to be publishedin JPS Conf. Proc. (2018).
62) H. Mukuda, S. Shimizu, A. Iyo, and Y. Kitaoka, J. Phys. Soc. Jpn. 81,
011008 (2012); S. Shimizu, S.-i. Tabata, S. Iwai, H. Mukuda, Y. Kitaoka, P. M. Shirage, H. Kito, and A. Iyo, Phys. Rev. B 85, 024528(2012).
63) H. Yokoyama, R. Sato, and K. Kobayashi, to be published in JPS Conf.
Proc. (2018), and in preparation.
64) The element (ii) is nonexistent for superconducting states, because all
k-points in the Brillouin zone are used.
13