Electron pairing: from metastable electron pair to bipolaron

Since the discovery of high T_c superconductivity by Bednorz and Müller [1] in 1986, great progresses have been made in the experimental and theoretical investigation of unconventional superconductivity. However, the mechanism of electron pairing in unconventional superconductors remains one of the most challenging and unresolved problems in condensed matter physics [2–4]. Vast experimental evidences have shown that electron pairing and unconventional superconductivity occur in many different materials, such as cuprates [1–4], iron-based superconductors [5–8], and carbon-based superconductors [9, 10], etc. Although there are many different theories for unconventional superconductivity, almost all theories follow the basic idea of the BCS theory [11]. They presume that there is some effective attraction between electrons leading to Cooper pairing which spontaneously condense into a collective non-Fermi liquid state.

We would like to mention some very recent experimental results related to the electron-pairing mechanism in unconventional superconductors. Božović et al reported very impressive and accurate results on superconductivity in high-T_c cuprates [12]. They synthesized atomically perfect thin films and multilayers of cuprates La_2xSr_xCuO₄ (LSCO) and measured the absolute value of the magnetic penetration depth and the phase stiffness with high accuracy in thousands of samples. The large statistics revealed clear trends in the intrinsic properties of the cuprate superconductors. They found that the obtained results disagree with the BCS theory in any variant, i.e. clean or dirty, including the Migdal-Eliashberg theory. Rather, the experimental data indicated small (local) and very light electron pairs with mass on the order of an electron mass. These pairs are preformed well above T_c and at T_c undergo Bose–Einstein condensation [12, 13]. Investigations performed by Zhong et al [14] and by Ren et al [15] challenged the d-wave pairing mechanism in cuprates. With scanning tunneling spectroscopy, they revealed anisotropic and nodeless superconducting gaps in the cuprate superconductors Bi₂Sr₂CaCu₂O_8+δ (Bi2212) and YBa₂Cu₃O_7−x (YBCO). In their paper, Zhong et al [14] affirmed that 'this is contradictory to the nodal d-wave pairing scenario that is often thought to be the most important result in the 30-year study of the HTS mechanism of cuprates'.

Important progresses in recent investigations on the electron pairing mechanism in iron-based superconductors indicate small and preformed Cooper pairs [6–8]. For instance, using Bogoliubov quasiparticle interference imaging, Sprau et al [6] found that the superconducting energy gap in FeSe is extremely anisotropic and nodeless. Their investigation discovered the existence of orbital-selective Cooper pairing in FeSe. Gerber et al [8] combined two time-domain experiments into a coherent lock-in measurement in the terahertz regime and was able to quantify the electron-phonon coupling strength in FeSe. Their study revealed a strong enhancement of the electron-phonon coupling strength in FeSe owing to electron correlations and highlighted the importance of the cooperative interplay between electron-electron and electron-phonon interactions.

In this paper, we present a theory for electron pairing in crystals where we consider the electron-electron correlation, the periodic potential of the crystal lattice, and the electron-phonon interaction. Our theory is different from all previous ones. We obtain preformed small electron pairs in the periodic potential of the crystal. This is essentially an orbital dependent electron-pairing theory in which the exchange-correlation between two electrons occupying the same orbital is decisive for pair formation. However such pairs are metastable unless the electron-phonon interaction is included. The calculations show that the electron-phonon coupling and the polaron effects are responsible for the stabilization of the electron pairs.

Our basic idea for the electron pairing process is that the electron-electron correlation is orbital dependent. Within the jellium model for the electron gas in solids [16], the atomic nuclei that form the periodic lattice are smeared out into a uniform positive charge distribution. Each electron is totally delocalized. Therefore, many electrons 'see' each other with their fluctuation potential at the same time and thus correlate all at once, giving rise to collective screening and oscillation effects. However, in atoms and molecules, significant correlations occur within electron pairs [17–19]. Strong exchange-correlation interaction between two electrons in the same orbital manifests in the shell structure of atoms and also in the covalent and ionic bonding in molecules. Our starting point in this study is to distinguish the exchange-correlation effects between two electrons of opposite spins occupying the same atomic orbital from the average correlation among many electrons in a crystal. This may happen in a crystal but the electrons have to 'feel' the nuclei potential well. This leads to a preliminary condition that the effective Bohr radius of the valence electrons in the crystal has to be comparable or smaller than the lattice constant. For instance, for a cuprate crystal with effective electron mass ${m}^{* }\simeq 5{m}_{0}$ and static dielectric constant ${\epsilon }_{0}\simeq 30$, the effective Bohr radius ${{\rm{a}}}_{B}\simeq 3.2$ Å is smaller than the lattice constant of about 3.8 Å. Because we want to show that the electron-pair correlation in atoms can manifest themselves in electron transport in crystals, our calculations have to start first with the formation of energy bands [20].

In order to find out the electron-pair states in the crystal, we will first establish a simple crystal model to discuss the physical process. We consider a 'hydrogen solid' model with single-electron state of H atom and electron-pair state of the H⁻ ion. We will show that, besides the energy bands from the single-electron energy levels of individual atoms in the crystal, there can exist a metastable electron-pair energy band from the correlated electron pairs of the H⁻ state for the lattice constant λ being larger than the effective Bohr radius a_B [21]. The electron pairs are metastable because the Coulomb repulsion is strong overwhelming the exchange-correlation. In order to stabilize them we have to include the electron-phonon interaction to counterbalance the Coulomb repulsion. Therefore, the electron-phonon interaction is necessitated in a natural way in the electron pairing process.

This paper is organized as follows. In section 2 we calculate and discuss the metastable electron-pair energy band in two-dimensional square-lattice crystals. In section 3 we present the many-particle Hamiltonian consisting of electrons in both the single-electron and electron-pair bands being coupled to the LO-phonons. In sections 4 and 5, we calculate the many-particle exchange-correlation (xc) energies due to electron-electron interactions and polaron energies due to electron-phonon interactions. Then, we show in section 6 the conditions under which the metastable electron pairs can be stabilized in the ground state including many-body effects. Finally, we summarize our work in section 7.

In the calculations, we will use effective Bohr radius ${{\rm{a}}}_{B}={\epsilon }_{0}{{\hslash }}^{2}/{m}^{* }{e}^{2}$ and effective Rydberg ${{\rm{R}}}_{y}={{\hslash }}^{2}/2{m}^{* }{{\rm{a}}}_{B}^{2}$ as the units for length and energy, respectively.

Electronic band structure is of fundamental importance for our understanding of many physical properties of solids. Within the independent electron approximation, the electron states are of Bloch form in the periodic potential of the crystal lattice. The effects of electron-electron interactions are accounted for by an effective potential which repeats this periodicity [20]. In this section, we will show that, besides the energy bands from the single-electron states there can exist a metastable electron-pair energy band depending upon the crystal structure and potential. This metastable electron-pair band originates from two correlated electrons of opposite spins occupying the same atomic orbital.

2.1. Two-electron atoms

Our study will start with the simplest electron-pair system, i.e., two-electron atoms. These helium-like atoms with two electrons of opposite spins occupying the same orbital, e.g. helium atom He and negatively charged hydrogen anion H⁻ have played an important role in the development of theoretical physics in the last century [22]. It is a challenge to determine accurately the correlation energy, even in simple systems such as He atom and H⁻ ion [22–25]. Hylleraas' result for the ground-state energy of He atom obtained in 1929 was −5.80648 R_y [26]. After generations of calculation, very accurate (non-relativistic) ground-state energies [27–32] of two electron atoms have been obtained: −5.807448754068 ⋯ R_y for He and −1.055502033088 ⋯ R_y for H⁻. Recently, using high-precision variational calculations Estienne et al [33] determined the critical nuclear charge Z_c = 0.911 028 224 077 255 73(4) which is the minimum charge required to bind two electrons in a helium-like atom.

On the other hand, the famous experiment on two-electron atoms by Madden and Codling [34] revealed that the simple model based on independent particle picture is inappropriate to characterize a series of doubly excited states because of strong electron-electron correlation [22, 34]. In comparison with the single-electron states of the H atom, the H⁻ ion is a closed-shell system with two strongly correlated electrons. Such an electron-pair state is different in its nature from the single-electron states because of the strong correlation. It should be recognized as a new strongly correlated electronic state.

Negative hydrogen ion H⁻ in two-dimensional (2D) system has also been investigated in the last decades mostly because of the discovery of its counterpart D⁻ center in 2D semiconductor quantum wells [35]. The D⁻ center is a negatively charged shallow donor impurity center in semiconductors, such as a negatively charged Si impurity in a GaAs quantum well. It is an H⁻-like state in solid-state environment but with very different energy and length scales (e.g., in GaAs, the effective Rydberg R_y = 5.9 meV and effective Bohr radius a_B = 98 Å). Therefore, the D⁻ center in semiconductors is considered as an ideal 'laboratory' to study the H⁻ properties, for instance, in high magnetic fields [36]. Earlier variational calculation by Phelps and Bajaj found the energy of the H⁻ in 2D is −4.48 R_y [37]. Further numerical calculations obtained −4.48054 R_y by Ivanov and Schmelcher [38] and −4.4804798 R_y by Ruan et al [39].

In table 1 we compare the ground-state energies of the 2D and 3D H⁻ states. The binding energy E_b is defined as the difference between the energies E_H of the neutral H atom and ${E}_{{{\rm{H}}}^{-}}$ of the H⁻ ion. This is the energy required to remove one of the two electrons from the H⁻ ion to infinity. It is also called electron affinity of the hydrogen atom. One sees that the binding energy E_b of the H⁻ in 2D is almost 9 times larger than that in 3D because electron correlation in 2D is much stronger. In the last column we also give the energy per electron ε_p in the H⁻ state.

Table 1.  The ground-state energies of negative hydrogen ion H⁻ and hydrogen atom H in 2D and 3D. E_b is the binding energy of H⁻. The energy per electron in H⁻ state is given by ${\varepsilon }_{p}={E}_{{{\rm{H}}}^{-}}/2$. Energies are in units of R_y.

	${E}_{{{\rm{H}}}^{-}}$	EH	Eb	εp
3D	−1.055	−1.0	0.055	−0.528
2D	−4.48	−4.0	0.48	−2.24

2.2. 'Hydrogen solid' model with both the single-electron and electron-pair states

In order to explain the so-called Mott insulator and metal-insulator transition [40, 41], Mott considered a hydrogen solid model with the single-electron energy band only, i.e., a simple cubic lattice crystal of one-electron atoms and made the following discussion [40]. For small values of the lattice constant λ, there is a half-filled band in such a crystal and thus it is metallic. If one varies the lattice constant to large values (but not so large as to prevent tunnelling), the Coulomb interaction U for two electrons occupying the same atomic site overcomes the kinetic energy (characterized by the band width W). In this case, each electron should be assigned to its parent atom. The crystal must be nearly the same as a collection of isolated neutral hydrogen atoms and thus it is an insulator. This reveals a competition between potential and kinetic effects. At large λ (small W) the Coulomb repulsion U dominates, the electrons are localized, and the system is insulating [40, 41]. The idea of Mott led to the theoretical model introduced by Hubbard [42]. The Hubbard model traces the insulating behavior to strong Coulomb repulsion between electrons occupying the same orbital. The competition between the kinetic and Coulomb energies gives rise to strong electron-electron correlations. The Hubbard model was proposed originally to describe the transition between conducting and insulating systems. It has also been widely used to study materials with strongly correlated electrons and high-temperature superconductivity [43].

In this paper we will consider both the single-electron and electron-pair states in the 'hydrogen solid' model. Our calculation will be performed for 2D systems because we can obtain more accurate numerical results in 2D. Another reason why 2D systems are more interesting is that electron-electron correlations are stronger. Many unconventional superconductor materials are found to be essentially two dimensional. Figure 1 shows diagrams representing (a) a 2D H atom and (b) a 2D H⁻ ion with their respective energy levels. The nuclear potential V_a(r − R_m) of the atom is represented by the black curves, where R_m is the position of the nucleus. The single-electron levels of a 2D H atom are given by ${\varepsilon }_{i}=-{{\rm{R}}}_{y}/{(i+1/2)}^{2}$ (for i = 0, 1, 2, ...). The energy level of a correlated electron pair in H⁻ ion, i.e., the energy per electron in the ground state, is given by ${\varepsilon }_{p}={E}_{{{\rm{H}}}^{-}}/2=-2.24$ R_y. We remind that a single H⁻ ion is stable.

Zoom In Zoom Out Reset image size

We now consider the following 'hydrogen solid': N atoms are arranged into a simple-lattice crystal at positions R_m (for m = 0, 1, 2, ⋯, N − 1) with lattice constant λ of the order of the Bohr radius a_B as indicated in figure 1(c). The crystal potential for an electron at r_j is given by

$$\begin{eqnarray}&&{V}_{c}({{\bf{r}}}_{j})=\displaystyle \sum _{m=0}^{N-1}{V}_{a}({{\bf{r}}}_{j}-{{\bf{R}}}_{m}).\end{eqnarray} \tag{ 1 }$$

It is known that the single-electron levels ε_i of individual atoms form energy bands E_i(k) in such a crystal [20] as indicated by the horizontal thick-blue lines in figure 1(c). In principle, there is also the possibility that two electrons of opposite spins occupy the same atomic orbital forming an H⁻-like state in the crystal, but it becomes unstable due to the presence of the neighbor atoms. Therefore, the counterpart of the H⁻ state in a crystal has never been investigated. In this section we will show that, though such an electron pair is unstable in a crystal due to the Coulomb repulsion, they may form a metastable energy band (indicated by the horizontal thick-green line in figure 1(c)) depending on the crystal structure and potential. In such a crystal the lattice constant λ should not be so small as to prevent individual atoms to bind two electrons, but not so large as to prevent co-tunnelling of an electron pair between the neighbor unit cells.

In the center of mass and relative coordinates (R, r) of the two electrons at r₁ and r₂, defined by

$$\begin{eqnarray}&&{\bf{R}}=\displaystyle \frac{1}{2}({{\bf{r}}}_{1}+{{\bf{r}}}_{2})\,\,\,\mathrm{and}\,\,\,\,{\bf{r}}={{\bf{r}}}_{1}-{{\bf{r}}}_{2},\end{eqnarray} \tag{ 2 }$$

the wavefunction of an individual electron pair with energy 2ε_p bound to the atom at R_m is given by $\phi ({\bf{R}}-{{\bf{R}}}_{m},{\bf{r}})$. The Schrödinger equation for two electrons in the crystal potential given by equation (1) can be written as,

$$\begin{eqnarray}&&\left[-\displaystyle \frac{1}{2}{{\rm{\nabla }}}_{{\bf{R}}}^{2}-2{{\rm{\nabla }}}_{{\bf{r}}}^{2}+{V}_{c}({{\bf{r}}}_{1})+{V}_{c}({{\bf{r}}}_{2})+\displaystyle \frac{2}{| {\bf{r}}| }\right]{\rm{\Psi }}({\bf{R}},{\bf{r}})\\ &&\quad =\,2{E}_{p}{\rm{\Psi }}({\bf{R}},{\bf{r}}),\end{eqnarray} \tag{ 3 }$$

where ${{\bf{r}}}_{1}={\bf{R}}+\tfrac{{}_{1}}{{}^{2}}{\bf{r}}$ and ${{\bf{r}}}_{2}={\bf{R}}-\tfrac{{}_{1}}{{}^{2}}{\bf{r}}$. The term $2/| {\bf{r}}|$ is the Coulomb repulsion potential between the two electrons. We should bear in mind that, due to electron-electron repulsion, the ground state of this system can be found for $| {{\bf{r}}}_{1}-{{\bf{r}}}_{2}| =| {\bf{r}}| \to \infty$. In other words, the ground state of this two-electron system corresponds to two non-interacting single electrons separated by an infinitely long distance. But we are looking for the quantum states of two-correlated electrons occupying the same orbital in the same unit cell in the crystal with the average separation $\langle r\rangle$ being less than the lattice constant λ. Therefore, the electron-pair states in the crystal are metastable. The calculations in section 2.3 will confirm that such a metastable state does exist in 2D periodic potentials.

If the electron-pair states of two correlated electrons in the crystal can be approximated by a linear combination of the electron-pair wavefunctions ϕ(R − R_m, r ) of single atoms, written as,

$$\begin{eqnarray}&&{\rm{\Psi }}({\bf{R}},{\bf{r}})=\displaystyle \sum _{m}{c}_{m}\phi ({\bf{R}}-{{\bf{R}}}_{m},{\bf{r}})\end{eqnarray} \tag{ 4 }$$

for $\langle r\rangle \lt \lambda$, we can obtain the following homogeneous linear equations,

$$\begin{eqnarray}&&\displaystyle \sum _{m}[{J}_{p}({{\bf{R}}}_{m}-{{\bf{R}}}_{n})+2({E}_{p}-{\varepsilon }_{p}){\alpha }_{p}({{\bf{R}}}_{m}-{{\bf{R}}}_{n})]{c}_{m}=0,\end{eqnarray} \tag{ 5 }$$

for n, m = 0, 1, 2, ..., N − 1, where α_p(R_l ) is the overlap integral

$$\begin{eqnarray}&&{\alpha }_{p}({{\bf{R}}}_{l})=\int d{\bf{R}}\int d{\bf{r}}{\phi }^{* }({\bf{R}}-{{\bf{R}}}_{l},{\bf{r}})\phi ({\bf{R}},{\bf{r}}),\end{eqnarray} \tag{ 6 }$$

with α_p(R_l = 0 ) = 1, and

$$\begin{eqnarray}&&{J}_{p}({{\bf{R}}}_{l})=-\int d{\bf{R}}\int d{\bf{r}}{\phi }^{* }({\bf{R}}-{{\bf{R}}}_{l},{\bf{r}}){\rm{\Delta }}{V}_{l}({\bf{R}},{\bf{r}})\phi ({\bf{R}},{\bf{r}}),\end{eqnarray} \tag{ 7 }$$

with

$$\begin{eqnarray}&&{\rm{\Delta }}{V}_{l}({\bf{R}},{\bf{r}})=\displaystyle \sum _{n\ne l}\left[{V}_{a}({\bf{R}}-{{\bf{R}}}_{n}+\displaystyle \frac{{}_{1}}{{}^{2}}{\bf{r}})+{V}_{a}({\bf{R}}-{{\bf{R}}}_{n}-\displaystyle \frac{{}_{1}}{{}^{2}}{\bf{r}})\right].\end{eqnarray} \tag{ 8 }$$

We observe that equation (5) for c_m depends only on R_l = R_m − R_n . This is an eigenvalue problem of a block circulant matrix [44]. The solution has the following form

$$\begin{eqnarray}&&{c}_{m}=C\cdot {e}^{i{\bf{k}}\cdot {{\bf{R}}}_{m}},\end{eqnarray} \tag{ 9 }$$

where the vector ${\bf{k}}$ should be a reduced wavevector in the first Brillouin zone and C is the normalization constant. We finally obtain the electron-pair wavefunction in the crystal given by

$$\begin{eqnarray}&&{{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}})=\displaystyle \frac{{\displaystyle \sum }_{l}{e}^{i{\bf{k}}\cdot {{\bf{R}}}_{l}}\phi ({\bf{R}}-{{\bf{R}}}_{l},{\bf{r}})}{\sqrt{N\left(1+{\displaystyle \sum }_{l\ne 0}{e}^{i{\bf{k}}\cdot {{\bf{R}}}_{l}}{\alpha }_{p}({{\bf{R}}}_{l})\right)}}.\end{eqnarray} \tag{ 10 }$$

The above wavefunction is a Bloch wavefunction in the center of mass coordinates R because it can be written as

$$\begin{eqnarray}&&{{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}})={{Ce}}^{i{\bf{k}}\cdot {\bf{R}}}\left[\displaystyle \sum _{l}{e}^{-i{\bf{k}}\cdot ({\bf{R}}-{{\bf{R}}}_{l})}\phi ({\bf{R}}-{{\bf{R}}}_{l},{\bf{r}})\right],\end{eqnarray} \tag{ 11 }$$

where the part in the square brackets is a periodic function in the coordinates R with period of the crystal lattice. The dispersion relation of the electron-pair band is given by

$$\begin{eqnarray}&&{E}_{p}({\bf{k}})={\varepsilon }_{p}-\displaystyle \frac{{\displaystyle \sum }_{l}\,{J}_{p}({{\bf{R}}}_{l}){e}^{i{\bf{k}}\cdot {{\bf{R}}}_{l}}}{2{\displaystyle \sum }_{l}\,{\alpha }_{p}({{\bf{R}}}_{l}){e}^{i{\bf{k}}\cdot {{\bf{R}}}_{l}}}.\end{eqnarray} \tag{ 12 }$$

Because the considered electron-pair wavefunction $\phi ({\bf{R}}-{{\bf{R}}}_{m},{\bf{r}})$ of a two-electron atom has essentially the s-symmetry [24, 37], the dispersion relation of the electron-pair band in the 2D square-lattice crystal with lattice constant λ can be approximated as

$$\begin{eqnarray}&&{E}_{p}({\bf{k}})={\varepsilon }_{p}-\displaystyle \frac{1}{2}{J}_{p}(0)-{J}_{p}({{\bf{R}}}_{1})[\cos ({k}_{x}\lambda )+\cos ({k}_{y}\lambda )],\end{eqnarray} \tag{ 13 }$$

where only the nearest-neighbor tunneling term J_p(R₁) is considered. The value of J_p(R₁) determines the bandwidth of the electron-pair band. Because this is a co-tunneling process of two paired electrons between the neighbor sites and the effective mass of the electron pair is twice of a single electron, the electron-pair bandwidth should be much smaller than that of the single-electron band. This dispersion relation will be confirmed in the next section by making numerical calculations of the metastable electron-pair band in 2D periodic potential of a square lattice.

2.3. Metastable electron-pair band in 2D square lattices

For a quantitative demonstration of the metastable electron-pair band in a crystal and its renormalization due to many-body effects, we will use the following 2D periodic potential. For an electron at r_j = (x_j, y_j) in a 2D square lattice with the lattice constant λ, the considered potential is given by

$$\begin{eqnarray}&&{V}_{c}({{\bf{r}}}_{j})={V}_{0}[\cos ({{qx}}_{j})+\cos ({{qy}}_{j})],\end{eqnarray} \tag{ 14 }$$

where q = 2π / λ and V₀ is the amplitude of the crystal potential. Notice that, V₀ is not a measurable quantity, e.g., the amplitude of the crystal potential in figure 1(c) should be infinity. The potential defined in equation (14) with two parameters λ and V₀ will simplify our numerical calculations without losing any essential features of the theory. In this 2D periodic potential, the energy has a continuous spectrum for $E\geqslant 0$. Therefore, two electrons can possibly bind into a pair for E < 0 only. We have calculated the single-electron and metastable electron-pair states in this periodic potential. For the calculation details we refer to [21].

The single-electron states are well known in this potential. The Schrödinger equation for a single electron is given by

$$\begin{eqnarray}&&{H}_{0}{\psi }_{{\bf{k}}+{{\bf{G}}}_{l}}({{\bf{r}}}_{j})={E}_{l}({\bf{k}}){\psi }_{{\bf{k}}+{{\bf{G}}}_{l}}({{\bf{r}}}_{j}),\end{eqnarray} \tag{ 15 }$$

with

$$\begin{eqnarray}&&{H}_{0}({{\bf{r}}}_{j})=-{{\rm{\nabla }}}_{j}^{2}+{V}_{c}({{\bf{r}}}_{j}),\end{eqnarray} \tag{ 16 }$$

where k is the wavevector in the first Brillouin zone, l is the band index, and ${{\bf{G}}}_{l}={l}_{x}q{\bf{i}}+{l}_{y}q{\bf{j}}$ (with ${l}_{x},{l}_{y}=0,\pm 1,\pm 2,\ldots$) is the reciprocal-lattice vector; E_l(k) and ${\psi }_{{\bf{k}}+{{\bf{G}}}_{l}}({\bf{r}})$ are the eigenvalue and eigenfunction, respectively.

When we consider two electrons in this periodic potential, their Hamiltonian is given by

$$\begin{eqnarray}&&H={H}_{0}({{\bf{r}}}_{1})+{H}_{0}({{\bf{r}}}_{2})+\displaystyle \frac{2}{| {{\bf{r}}}_{1}-{{\bf{r}}}_{2}| },\end{eqnarray} \tag{ 17 }$$

where the last term is the Coulomb repulsion potential between the two electrons. In the center of mass R = (X, Y) and relative r = (x, y) coordinates defined in equation (2), the two-electron Hamiltonian becomes

$$\begin{eqnarray}H & = & -\displaystyle \frac{1}{2}{{\rm{\nabla }}}_{{\bf{R}}}^{2}-2{{\rm{\nabla }}}_{{\bf{r}}}^{2}+\displaystyle \frac{2}{r}\\ & & +2{V}_{0}\left[\cos ({qX})\cos \left(\displaystyle \frac{{qx}}{2}\right)+\cos ({qY})\cos \left(\displaystyle \frac{{qy}}{2}\right)\right].\end{eqnarray} \tag{ 18 }$$

This Hamiltonian is periodic in X and Y with period λ. We can choose a Bloch wavefunction in the center-of-mass coordinates for our basis. As to the function in the relative coordinates r = (r, θ), we have to consider the symmetry of the electron-electron Coulomb potential and the periodic potential representing a 2D square lattice. We use the following basis for our wavefunction,

$$\begin{eqnarray}&&{\psi }_{{l}_{x},{l}_{y};n,m}({\bf{R}},{\bf{r}})=\displaystyle \frac{1}{\sqrt{A}}{e}^{i({\bf{k}}+{{\bf{G}}}_{l})\cdot {\bf{R}}}{R}_{n,m}(r){\phi }_{m}(\theta ),\end{eqnarray} \tag{ 19 }$$

with

$$\begin{eqnarray}&&{R}_{n,m}(r)=\beta {c}_{n,m}{(2\beta {\xi }_{n}r)}^{m}{e}^{-\beta {\xi }_{n}r}{L}_{n-m}^{2m}(2\beta {\xi }_{n}r),\end{eqnarray} \tag{ 20 }$$

and

$$\begin{eqnarray}&&{\phi }_{m}(\theta )=\displaystyle \frac{1}{\sqrt{{b}_{m}\pi }}\cos (m\theta ),\end{eqnarray} \tag{ 21 }$$

where n = 0, 1, 2, ⋯, m = 0, 1, 2, ⋯, n, ξ_n = 2/(2n+1), ${c}_{n,m}={[2{\xi }_{n}^{3}(n-m)!/(n+m)!]}^{1/2}$, b₀ = 2, b_m = 1 for m ≥ 1, and ${L}_{n-m}^{2m}(x)$ is the generalized Laguerre polynomial. The function R_{n, m}(r) is taken from the wavefunction of a 2D hydrogen atom [45, 46] with a modification introduced by a dimensionless scaling parameter β. The two-electron wavefunction can be written as

$$\begin{eqnarray}&&{{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}})=\displaystyle \sum _{{l}_{x},{l}_{y}}\displaystyle \sum _{n,m}{a}_{{l}_{x},{l}_{y};n,m}({\bf{k}}){\psi }_{{l}_{x},{l}_{y};n,m}({\bf{R}},{\bf{r}}).\end{eqnarray} \tag{ 22 }$$

Considering the antisymmetry of the electron wavefunctions with spin states, the two-electron wavefunction of the spin singlet state is given by the above expression with the sum over even m only.

Solving the corresponding eigenvalue equation of the two-electron Hamiltonian given by equation (18) with the above basis, we find a metastable electron-pair state of spin-singlet in the 2D square lattice potential. As shown in [21], for fixed period λ, a metastable electron-pair state can be found when V₀ is larger than a certain value. A minimum potential V₀ required for a metastable electron-pair state in the periodic potential corresponds to the critical nuclear charge Z_c for a helium-like two-electron atom. The metastable electron-pair state exist for E < 0 only. This indicates that co-tunneling of the paired electrons occurs in the formation of the electron-pair band in this 2D periodic potential. In the calculations we found that the average separation $\langle r\rangle$ between two electrons in a metastable pair state is always smaller than half the lattice constant, $\langle r\rangle \lt \lambda /2$. The global minimum of the eigenenergy of the two-electron system occurs at $\langle r\rangle \to \infty$ corresponding to two non-interacting single electrons. We also want to emphasize that the parameter β in the wavefunction in equation (20) plays the role of variational parameter to improve the correlation energy of the electron pair. Because this parameter is directly related to the average distance between the two electrons in a pair, it helps us to understand better the metastable electron-pair state. However, the parameter β does not determine the existence of the electron-pair state in the periodic potential.

In the following, we will present numerical results for the band structure and show that the dispersion relation of the electron-pair band in this potential fits very well the expression given by equation (13). Figure 2(a) shows the dispersion relations of the electron states in the 2D crystal potential with λ = 1.3 a_B and V₀ = 15 R_y. The dispersion relation E_p(k) of the spin-singlet metastable electron-pair band is given together with that of the lowest single-electron band ${E}_{0}({\bf{k}})$. The electron-pair band remains above the corresponding single-electron band because the Coulomb repulsion between the two electrons is stronger than their correlation. However, the shape of the dispersion relations of the two bands are very similar. This confirms the dispersion relation of the electron-pair band discussed in previous section within the framework of the tight-binding approach. The similarity is due to the fact that two paired electrons are closely bound in the relative coordinates in real space with an average separation $\langle r\rangle$ less than half of the lattice period. Furthermore, the single-electron and electron-pair states in the relative coordinates are of the same symmetry. The average distance between the two electrons in the case of figure 2(a) is $\langle r\rangle =0.44\lambda$. The energy gap between the electron-pair band and the single-electron band is ${E}_{g}^{(0)}$. The energy difference between the bottoms of the two bands at Γ point is defined as ν₀. The electron pair behaves as a larger particle with both the mass and charge twice of a single electron. Consequently, the tunneling probability of an electron pair to its neighbor site is much smaller than that of a single electron leading to a much narrower electron-pair band.

Zoom In Zoom Out Reset image size

We also find that the dispersion relation of the electron-pair band can indeed be described by equation (13). For instance, the dispersion of the electron-pair band in figure 2(a) is fitted very well by ${E}_{p}({\bf{k}})={E}_{p,0}-{J}_{p,1}[\cos ({k}_{x}\lambda )+\cos ({k}_{y}\lambda )]$, where E_p,0 = −1.35401 ± 0.00002 R_y and J_p,1 = 0.03368 ± 0.00002 R_y. The fitting gives an extremely small fitting parameter χ² = 3 × 10⁻⁸.

In figure 2(b) we plot the electron-pair band as a function of the lattice constant λ together with the two lowest single-electron bands for fixed potential amplitude V₀ = 15 R_y. The two curves for each band indicate the minimum and maximum energies of the band. We see that the electron-pair band appears for λ ≳ 1.22 a_B with a bandwidth $\lesssim 0.2$ R_y. It stays above the lowest single-electron band with a gap ${E}_{g}^{(0)}$ of about 3 to 5 R_y. The bandwidth of the electron-pair band is at least one order of magnitude smaller than that of the single-electron band. Figure 2(c) shows the electron-pair bands as a function of λ for different V₀. We see that the electron-pair band appears for E < 0 because the energy spectrum is continuous for positive energy in the crystal potential given by equation (14). It means that co-tunneling of the paired electrons is required to form the energy band. Therefore, their bandwidth is much smaller than that of the single-electron band. The existence of the metastable electron-pair band is a result of the local confinement in each unit cell, the electron-electron correlation, and the co-tunneling of the electron pair in the crystal. In figure 2(d) we plot the energy gap ${E}_{g}^{(0)}$ together with ν₀. They are important quantities for the renormalization of the electron-pair states.

In the rest of the paper, we will demonstrate that the metastable electron pairs can be stabilized at certain electron densities by including electron-electron and electron-phonon interactions in the crystal. Since the electron pairs are spin singlet, they will be considered as bosonic quasiparticles and mostly distributed at the bottom of the electron-pair band at low temperature. The many-body effects in the crystal renormalize the band structure. If the band renormalization can bring down the bottom of the electron-pair band to the Fermi surface of the single-electron band, the electron pairs at the bottom of the band cannot decay into two single electrons because of the Pauli exclusion principle. In such a case, the electron pairs can be stabilized.

We consider the 2D electron system with two energy bands: a lower single-electron band E_s(k) and a higher metastable electron-pair band ν₀ + E_p(k). The bottom of the single-electron band is taken as reference for energy E = 0 and the bottom of the electron-pair band is at $E={\nu }_{0}$ (see figure 2). Assuming that there are N_t electrons in the system consisting of N_s single electrons and N_p electron pairs with N_t = N_s + 2N_p, the many-particle exchange-correlation (xc) interactions will renormalize the energy bands reducing the energies of both the single electrons and electron pairs. In real materials, a 2D electron system can be found in the interface and surface of bulk materials or in a 2D layer of layered crystal structure such as superconducting cuprate, therefore interaction between the electron system and crystal lattice vibration and polarization affects the electron states. Considering the ionic and polar-covalent characteristics of many superconducting materials, the electron-LO-phonon interactions can be significant and their contribution to the energy band renormalization is important [16, 47]. In this section, we will present the many-particle Hamiltonian including electron-electron (e-e) and electron-phonon (e-ph) interactions assuming that the 2D electron layer is immersed in a 3D phonon field. The Hamiltonian of the considered system is giving by

$$\begin{eqnarray}&&H={H}_{\mathrm{el}}+{H}_{\mathrm{ph}}+{H}_{\mathrm{el} \mbox{-} \mathrm{ph}},\end{eqnarray} \tag{ 23 }$$

where ${H}_{\mathrm{el}}$ is for the electronic part with single electrons and electron pairs, H_ph for 3D LO-phonons, and H_el−ph for e-ph interaction. The Hamiltonian of the electronic part H_el was derived in [48], and is given by

$$\begin{eqnarray}&&{H}_{\mathrm{el}}={H}_{\mathrm{single}}+{H}_{\mathrm{pair}}+{H}_{{\rm{s}} \mbox{-} {\rm{p}}},\end{eqnarray} \tag{ 24 }$$

where H_single is for the single-electron (se) band, H_pair for the electron-pair (ep) band, H_s−p for single-electron-electron-pair (se-ep) interaction. For the N_s electrons in the single-electron band, their Hamiltonian is given by

$$\begin{eqnarray}&&{H}_{\mathrm{single}}=\displaystyle \sum _{{\bf{k}},\sigma }{E}_{s}({\bf{k}}){c}_{{\bf{k}},\sigma }^{\dagger }{c}_{{\bf{k}},\sigma }\\ &&\,+\,\displaystyle \frac{1}{2A}\displaystyle \sum _{{{\bf{k}}}_{1},{{\bf{k}}}_{2},{\bf{q}}}\displaystyle \sum _{\sigma ,{\sigma }^{{\prime} }}{v}_{q}{c}_{{{\bf{k}}}_{1}-{\bf{q}},{\boldsymbol{\sigma }}}^{\dagger }{c}_{{{\bf{k}}}_{2}+{\bf{q}},{{\boldsymbol{\sigma }}}^{{\prime} }}^{\dagger }{c}_{{{\bf{k}}}_{2},{{\boldsymbol{\sigma }}}^{{\prime} }}{c}_{{{\bf{k}}}_{1},{\boldsymbol{\sigma }}},\end{eqnarray} \tag{ 25 }$$

where ${v}_{q}={v}_{{ss}}(q)=2\tfrac{2\pi }{q}$ is the single-electron-single-electron (se-se) Coulomb potential, the operators c_k,σ^† and c_k,σ are creation and annihilation operators, respectively, for a single electron of momentum ℏk and spin σ. They obey the fermion anti-commutation relations $\{{c}_{{\bf{k}},\sigma },{c}_{{{\bf{k}}}^{{\prime} },{\sigma }^{{\prime} }}^{\dagger }\}={\delta }_{{\bf{k}},{{\bf{k}}}^{{\prime} }}{\delta }_{\sigma ,{\sigma }^{{\prime} }}$, $\{{c}_{{\bf{k}},\sigma },{c}_{{{\bf{k}}}^{{\prime} },{\sigma }^{{\prime} }}\}=0$, and $\{{c}_{{\bf{k}},\sigma }^{\dagger },{c}_{{{\bf{k}}}^{{\prime} },{\sigma }^{{\prime} }}^{\dagger }\}=0$. For N_p electron pairs in the electron-pair band, the Hamiltonian is given by

$$\begin{eqnarray}{H}_{\mathrm{pair}} & = & \displaystyle \sum _{{\bf{k}}}2({\nu }_{0}+{E}_{p}({\bf{k}})){b}_{{\bf{k}}}^{\dagger }{b}_{{\bf{k}}}\\ & & +\displaystyle \frac{1}{2A}\displaystyle \sum _{{{\bf{k}}}_{1},{{\bf{k}}}_{2},{\bf{q}}}{v}_{{pp}}(q){b}_{{{\bf{k}}}_{1}-{\bf{q}}}^{\dagger }{b}_{{{\bf{k}}}_{2}+{\bf{q}}}^{\dagger }{b}_{{{\bf{k}}}_{2}}{b}_{{{\bf{k}}}_{1}},\end{eqnarray} \tag{ 26 }$$

where v_pp(q) = 4v_qf_pp(q) is the electron-pair-electron-pair (ep-ep) interaction potential with the form factor f_pp(q) given in [48], the operators b^†_k and b_k are creation and annihilation operators, respectively, for a spin-singlet electron pair of momentum ℏk. They obey the boson commutation relations $[{b}_{{\bf{k}}},{b}_{{{\bf{k}}}^{{\prime} }}^{\dagger }]={\delta }_{{\bf{k}},{{\bf{k}}}^{{\prime} }}$, $[{b}_{{\bf{k}}},{b}_{{{\bf{k}}}^{{\prime} }}]=0$, and $[{b}_{{\bf{k}}}^{\dagger },{b}_{{{\bf{k}}}^{{\prime} }}^{\dagger }]=0$.

The se-ep interband interaction is give by

$$\begin{eqnarray}&&{H}_{{\rm{s}} \mbox{-} {\rm{p}}}={H}_{\mathrm{int}}^{{\rm{s}}}+{H}_{\mathrm{int}}^{{\rm{t}}},\end{eqnarray} \tag{ 27 }$$

with

$$\begin{eqnarray}&&{H}_{\mathrm{int}}^{{\rm{s}}}=\displaystyle \frac{1}{A}\displaystyle \sum _{{\bf{k}},{{\bf{k}}}_{1},{\bf{q}},\sigma }{v}_{{sp}}(q){b}_{{{\bf{k}}}_{1}-{\bf{q}}}^{\dagger }{c}_{{\bf{k}}+{\bf{q}},\sigma }^{\dagger }{b}_{{{\bf{k}}}_{1}}{c}_{{\bf{k}},{\boldsymbol{\sigma }}}\end{eqnarray} \tag{ 28 }$$

and

$$\begin{eqnarray}{H}_{\mathrm{int}}^{{\rm{t}}} & = & \displaystyle \frac{1}{\sqrt{A}}\displaystyle \sum _{{\bf{k}},{\bf{q}}}{v}^{{\rm{t}}}(q)\\ & & \times ({b}_{{\bf{k}}}^{\dagger }{c}_{\tfrac{{\bf{k}}}{2}+{\bf{q}},\uparrow }{c}_{\tfrac{{\bf{k}}}{2}-{\bf{q}},\downarrow }+{b}_{{\bf{k}}}{c}_{\tfrac{{\bf{k}}}{2}-{\bf{q}},\downarrow }^{\dagger }{c}_{\tfrac{{\bf{k}}}{2}+{\bf{q}},\uparrow }^{\dagger }),\end{eqnarray} \tag{ 29 }$$

where v_sp(q) = 2v_q f_sp(q) is the se-ep interband scattering potential (without breaking the electron pair) with form factor f_sp(q) and v^t(q) is the se-ep interaction potential for interband transition (breaking or forming electron pairs). They are given in [48] and f²_sp(q) = f_pp(q). Notice that in the above Hamiltonian, the summation over q does not include q = 0 because it is cancelled out with the background ion-ion interaction and the system is neutral.

As we discussed in the previous section, the single electrons and the paired electrons share the same space in the crystal. The maximum total electron density is two electrons per unit cell (per atom), paired or not. However, when the potential amplitude V₀ of the crystal is larger than a certain value, two electrons in the same unit cell will occupy the same atomic orbital in the relative coordinates forming an electron pair. It means that in this case, each unit cell can be occupied by a single electron or by an electron pair, but not both at the same time. Therefore, for a certain single-electron density n_s = N_s/A (where A is the area of the sample), the maximum electron-pair density n_p = N_p/A in the square-lattice crystal is given by

$$\begin{eqnarray}&&{n}_{p}^{\max }={\lambda }^{-2}-{n}_{s},\end{eqnarray} \tag{ 30 }$$

where λ⁻² is the density of the unit cell of the crystal. Notice that n_s = λ⁻² corresponds to half filling of the single-electron band. The above condition indicates that the single-electron band should be less than half filled if there are any electron pairs in the crystal.

The Hamiltonian of the optical-phonon modes in bulk materials with energy ${\hslash }{\omega }_{{}_{{LO}}}$ and 3D wavevector Q = (q, q_z) is given by

$$\begin{eqnarray}&&{H}_{\mathrm{ph}}=\displaystyle \sum _{{\bf{Q}}}{\hslash }{\omega }_{{}_{{LO}}}{a}_{{\bf{Q}}}^{\dagger }{a}_{{\bf{Q}}},\end{eqnarray} \tag{ 31 }$$

where ${a}_{{\bf{Q}}}^{\dagger }$ (a_Q ) is the creation (annihilation) operator of the LO-phonons. The interaction Hamiltonian of a many-electron system with the LO-phonons is given by [16],

$$\begin{eqnarray}&&{H}_{\mathrm{el} \mbox{-} \mathrm{ph}}=\displaystyle \sum _{j=1}^{{N}_{t}}\displaystyle \sum _{{\bf{Q}}}({V}_{{\bf{Q}}}{a}_{{\bf{Q}}}{e}^{i{\bf{Q}}\cdot {{\bf{r}}}_{j}}+{V}_{{\bf{Q}}}^{* }{a}_{{\bf{Q}}}^{\dagger }{e}^{-i{\bf{Q}}\cdot {{\bf{r}}}_{j}}),\end{eqnarray} \tag{ 32 }$$

with the Fourier coefficient of the e-ph interaction potential

$$\begin{eqnarray}&&{V}_{{\bf{Q}}}=-i{\hslash }{\omega }_{{}_{{LO}}}{\left(\displaystyle \frac{{\hslash }}{2{m}^{* }{\omega }_{{}_{{LO}}}}\right)}^{1/4}\sqrt{\displaystyle \frac{4\pi \alpha }{{{VQ}}^{2}}},\end{eqnarray} \tag{ 33 }$$

where r_j is the position of the electron j with band mass m^*. The Fröhlich electron-phonon coupling constant α is defined by

$$\begin{eqnarray}&&\alpha =\displaystyle \frac{{e}^{2}}{{\hslash }}\sqrt{\displaystyle \frac{{\hslash }}{2{m}^{* }{\omega }_{{}_{{LO}}}}}\left(\displaystyle \frac{1}{{\epsilon }_{\infty }}-\displaystyle \frac{1}{{\epsilon }_{0}}\right),\end{eqnarray} \tag{ 34 }$$

where ₀ and ${\epsilon }_{\infty }$ are the static and high-frequency dielectric constant, respectively.

In the considered two-band system, the gap between the two bands is ${E}_{g}^{(0)}$ in the zero density limit. Although the e-e and e-ph interactions in the crystal reduce the energies of both the single-electron and electron-pair bands, their effects on the electrons in the paired states are larger than those in the unpaired single-electron band. Therefore, we expect that many-body effects will reduce the energy gap between the two bands. We will calculate the exchange-correlation corrections as well as the polaron energies for both single electrons and electron pairs including the screening effects in order to find out whether the electron pairs can be stabilized or not.

In order to obtain the ground-state energy of the present many-particle system with the single electrons, electron pairs and LO-phonons, we will employ the Lee-Low-Pines transformation [49] in dealing with the many-polaron system [49–52]. For weak e-ph coupling, we can assume the ground state of the electron-phonon system as $| \mathrm{GS}\rangle =| {\mathrm{GS}}_{\mathrm{el}}\rangle | {\mathrm{VAC}}_{\mathrm{ph}}\rangle$, where $| {\mathrm{GS}}_{\mathrm{el}}\rangle$ is the ground state of the electronic part and $| {\mathrm{VAC}}_{\mathrm{ph}}\rangle$ is the phonon vacuum state with zero real phonons [49, 50]. The above approximation is valid for weak and intermediate e-ph coupling strength [47, 49, 50, 53, 54] and it allows us to calculate separately the electron exchange-correlation and polaron contributions to the ground-state energy of the system, given by

$$\begin{eqnarray}&&{E}_{\mathrm{GS}}={E}_{\mathrm{GS}}^{(\mathrm{el})}+{E}_{\mathrm{pol}}^{(\mathrm{tot})},\end{eqnarray} \tag{ 35 }$$

where ${E}_{\mathrm{GS}}^{(\mathrm{el})}$ is the ground-state energy of the electronic part without interaction with phonons and ${E}_{\mathrm{pol}}^{(\mathrm{tot})}$ is the total polaron correction due to electron-phonon interaction.

In the next two sections, we will calculate the contributions of the electronic exchange-correlation interaction and electron-phonon interaction to the renormalization of the single-electron end electron-pair energies.

It is known that in both 3D and 2D systems the exchange-correlation energy for band-gap renormalization (BGR) is almost independent of the band characteristics. The many-particle exchange-correlation energy depends only on the inter-particle distance r_s (determined by the particle density) in appropriate rescaled natural units in a universal manner [55–57]. The contribution of the electron-electron interaction to the BGR can be obtained by calculating the average exchange-correlation energies per particle or by calculating the self-energies of the particles involved. The kinetic energy is usually assumed to be unchanged in a renormalization process.

In this section, we will calculate the ground-state energy ${E}_{\mathrm{GS}}^{(\mathrm{el})}$ of the electronic part consisting of the single electrons and electron pairs. The corresponding Hamiltonian is given by equations (24)–(29). Although the electron pairs are metastable, we will treat them in the calculations as if they were stable particles. Only the final results including full many-body corrections will tell us whether they can really be stabilized or not. In the calculations, we will first take the single-electron density n_s and electron-pair density n_p as inputs. The single electrons are considered as fermions and electron pairs as bosons. Therefore, we are dealing with a many-particle system consisting of a boson-fermion mixture [58]. The exchange-correlation energies are obtained as a function of n_s and n_p. We then determine their contributions to the ground-state energy and the band renormalization. Within such a scheme, the se-ep interaction Hamiltonian ${H}_{\mathrm{int}}^{t}$ given in equation (29) will not be invoked explicitly in the calculations. However, transitions between the single-electron and electron-pair bands are permitted because of this term. The many-particle interaction energy in such a two-component system of boson-fermion mixture can be obtained by [59]

$$\begin{eqnarray}&&{E}_{{ij}}={\int }_{0}^{{e}^{2}}\displaystyle \frac{{E}_{{ij}}^{\mathrm{int}}(\xi )}{\xi }d\xi ,\end{eqnarray} \tag{ 36 }$$

where the inter-particle interaction potential ${E}_{{ij}}^{\mathrm{int}}(\xi )$ is given by

$$\begin{eqnarray}&&{E}_{{ij}}^{\mathrm{int}}({e}^{2})=\displaystyle \frac{1}{2A}\displaystyle \sum _{{\bf{q}}}{v}_{{ij}}(q)[{S}_{{ij}}(q)-{\delta }_{{ij}}].\end{eqnarray} \tag{ 37 }$$

The above potential depends on the 'bare' inter-particle potential v_ij(q) and static structure factor S_ij(q), for i, j = s (single electron) and p (electron pair). The static structure factor can be calculated by,

$$\begin{eqnarray}&&{S}_{{ij}}(q)=-\displaystyle \frac{1}{\pi \sqrt{{n}_{i}{n}_{j}}}{\int }_{0}^{\infty }d\omega {\chi }_{{ij}}(q,i\omega ).\end{eqnarray} \tag{ 38 }$$

Within the linear response theory, the density-density response function χ_ij(q, ω) of this two-component plasma is given by [60, 61],

$$\begin{eqnarray}&&\{{[\chi (q,\omega )]}^{{}^{-1}}\}{}_{{ij}}={[{\chi }_{{ii}}^{(0)}(q,\omega )]}^{{}^{-1}}{\delta }_{{ij}}-{\varphi }_{{ij}}(q),\end{eqnarray} \tag{ 39 }$$

where ${\chi }_{{ii}}^{(0)}(q,\omega )$ is the non-interacting polarizability [60–63] of the ith component and φ_ij(q) is the static effective interaction potential. The function ${\chi }_{{ss}}^{(0)}(q,\omega )$ is for a non-interacting 2D electron gas given in [62]. The function ${\chi }_{{pp}}^{(0)}(q,\omega )$ is the polarizability of a non-interacting 2D boson gas of electron pairs with density n_p. The effective potential φ_ij(q) defines a local field correction [60, 64] in terms of the 'bare' potential v_ij(q). Within the random-phase approximation (RPA), the local field correction on the static effective interaction is neglected and therefore ${\varphi }_{{ij}}(q)={v}_{{ij}}(q)$. The potential v_ij(q) has been determined in the previous section given by v_ss(q) = v_q, ${v}_{{pp}}(q)=4{v}_{q}{f}_{{pp}}(q)$, and ${v}_{{sp}}(q)={v}_{{ps}}(q)=2{v}_{q}{f}_{{sp}}(q)$.

For a non-interacting 2D boson (i.e., electron-pair) gas with density n_p, we can assume that all the electron pairs are in the same state at the bottom of the electron-pair band at zero temperature, i.e., in the condensate phase [61, 63]. The polarizability of the non-interacting boson gas is given by,

$$\begin{eqnarray}&&{\chi }_{{pp}}^{(0)}(q,i\omega )=-\displaystyle \frac{2{n}_{p}{\varepsilon }_{q,p}}{{\omega }^{2}+{\varepsilon }_{q,p}^{2}},\end{eqnarray} \tag{ 40 }$$

where ε_q,p = q² /2 .

The contribution of the e-e interaction to the single-electron band renormalization is given by ΔE_s = E_ss + E_sp and to the electron-pair band given by ΔE_p = (E_pp + E_sp)/2. The ground-state energy of the two-band system is given by

$$\begin{eqnarray}&&{E}_{\mathrm{GS}}^{(\mathrm{el})}={E}_{\mathrm{kin}}+{N}_{s}({E}_{{ss}}+{E}_{{sp}})+2{N}_{p}\left({\nu }_{0}+\displaystyle \frac{{E}_{{pp}}+{E}_{{sp}}}{2}\right),\end{eqnarray} \tag{ 41 }$$

where E_kin is the kinetic energy of the many-electron system. We calculate these energies within the RPA. Although the RPA overestimates the exchange-correlation energies ΔE_s and ΔE_p, only the difference between them contributes to the band-gap renormalization. The errors resulting from the RPA should be partially cancelled in the process when determining the condition of stability of the electron pairs. Therefore, we consider the RPA a reasonable approximation for our purpose. As mentioned above, the kinetic energy E_kin will be assumed unchanged in the renormalization. It is given by ${E}_{\mathrm{kin}}={N}_{s}{\varepsilon }_{F}/2$, where ε_F is the Fermi energy of the single-electron band in relation to its band bottom. In a two-dimensional system, the average kinetic energy of a single electron is ε_F /2. The electron pairs have no kinetic energy because they are assumed to be at the bottom of the electron-pair band in the condensate phase.

If we further assume that v_sp(q) = 0, the single electrons and electron pairs become independent in different energy bands. The single-electron system is an usual one-component 2D electron gas which has been widely investigated within different methods and approximations [62, 65, 66]. The system of electron pairs is new. Although the electron pairs are metastable states, we will treat them as stable particles when searching for their many-particle ground state. The static structure factor of a one-component boson system of electron pairs within the RPA is given by,

$$\begin{eqnarray}&&{S}_{{pp}}^{(1)}(q)={[1+2{n}_{p}{v}_{{pp}}(q)/{\varepsilon }_{q,p}]}^{-\tfrac{1}{2}}.\end{eqnarray} \tag{ 42 }$$

Consequently, we obtained the many-electron-pair exchange-correlation energy,

$$\begin{eqnarray}&&{E}_{{pp}}^{(1)}=8{({r}_{s}^{p})}^{{}^{-\tfrac{2}{3}}}I,\end{eqnarray} \tag{ 43 }$$

where ${r}_{s}^{p}={(\pi {n}_{p})}^{-\tfrac{1}{2}}$, and

$$\begin{eqnarray}&&I={\int }_{0}^{\infty }{dx}\left\{\left(\displaystyle \frac{{x}^{3}}{4}\right)\left[\sqrt{1+\displaystyle \frac{8}{{x}^{3}}{f}_{{pp}}(q)}-1\right]-{f}_{{pp}}(q)\right\}.\end{eqnarray} \tag{ 44 }$$

The calculations in [48] showed that f_pp(q) is a monotone decreasing function with f_pp(0) = 1 and ${f}_{{pp}}(\infty )=0$. If we take f_pp(q) ≡ 1 (consequently $\langle r\rangle =0$), i.e., assuming the electron pair as an ideal boson with charge −2e and mass $2{m}^{* }$, the above integral becomes I = I₀ = −1.29 . This is the well known RPA result for an ideal charged 2D boson gas [67]. Notice that the factor 8 in equation (43) is due to the units used here.

In figure 3 we show the many-particle interaction energies E_ss, E_sp, and E_pp within the RPA in the system keeping the same number of electrons in the single-electron and electron-pair bands, i.e., n_s = 2n_p. The energy E_sp = E_ps is positive because of the se-ep repulsion. Figure 3(a) is for the 2D crystal with $\lambda =1.5$ a_B and V₀ = 15 R_y. It shows the density dependence of the interaction energies and the effects of se-ep interaction. The blue-dash and red-dash curves are the exchange-correlation energies ${E}_{{pp}}^{(1)}$ and ${E}_{{ss}}^{(1)}$, respectively, without the se-ep interaction. In this case, the energy ${E}_{{pp}}^{(1)}$ is obtained from equations (43) and (44). We see that the nonzero distance between the two electrons in the pair (i.e., $\langle r\rangle \gt 0$) affects the energy E_pp. In the calculations we found that $\lambda /2\gt \langle r\rangle \gt \lambda /3$. If we assume $\langle r\rangle =0$ for the electron pairs, we obtain ${E}_{{pp}}^{(1)}={E}_{\mathrm{boson}}=-1.29\,\times \,8{({r}_{s}^{p})}^{{}^{-\tfrac{2}{3}}}$ for an ideal charged boson system in 2D indicated by the black-dotted curve.

Zoom In Zoom Out Reset image size

We observe that the se-ep interband interaction not only introduces the energy E_sp but also reduces the energies E_pp and E_ss being evident in the difference between the solid and dashed curves in figure 3(a). The density dependence of the ep-ep interaction energy E_pp is different from that of the se-se interaction E_ss. For example, E_pp is about 7 times larger than E_ss at lower density n_s = 2n_p = 0.001 ${{\rm{a}}}_{B}^{-2}$ and this ratio is reduced to 3 times at higher density n_s = 2n_p = 1.0 ${{\rm{a}}}_{B}^{-2}$. In figure 3(b) we show the many-particle interaction energies for different potential V₀. With increasing V₀ for the same lattice constant λ, the average distance $\langle r\rangle$ between the two electrons in the same pair decreases and, consequently, the form factor f_pp(q) increases (as shown in the inset) and the energy E_pp becomes larger.

In this section, we will study the polaron effects on the N_t electrons in the two-band system interacting with LO-phonons assuming N_s electrons in the single-electron band and N_p electron pairs in the electron-pair band, N_s + 2N_p = N_t. The Hamiltonian of the system is given in equation (23) and the electron-phonon interaction given by equation (32). Considering a boson-fermion mixture (namely, the electron pairs and single electrons) interacting with the phonons, the electron-phonon interaction Hamiltonian in equation (32) can be separated into two parts,

$$\begin{eqnarray}&&{H}_{\mathrm{el} \mbox{-} \mathrm{ph}}={H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{single}}+{H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{pair}},\end{eqnarray} \tag{ 45 }$$

with

$$\begin{eqnarray}&&{H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{single}}=\displaystyle \sum _{j=1}^{{N}_{s}}\displaystyle \sum _{{\bf{Q}}}({V}_{{\bf{Q}}}{a}_{{\bf{Q}}}{e}^{i{\bf{Q}}\cdot {{\bf{r}}}_{j}}+{V}_{{\bf{Q}}}^{* }{a}_{{\bf{Q}}}^{\dagger }{e}^{-i{\bf{Q}}\cdot {{\bf{r}}}_{j}}),\end{eqnarray} \tag{ 46 }$$

for N_s single electrons, and

$$\begin{eqnarray}{H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{pair}} & = & \displaystyle \sum _{j=1}^{{N}_{p}}\displaystyle \sum _{{\bf{Q}}}\left[\displaystyle \sum _{\nu =1}^{2}({V}_{{\bf{Q}}}{a}_{{\bf{Q}}}{e}^{i{\bf{Q}}\cdot {{\bf{r}}}_{j,\nu }}+{V}_{{\bf{Q}},}^{* }{a}_{{\bf{Q}}}^{\dagger }{e}^{-i{\bf{Q}}\cdot {{\bf{r}}}_{j,\nu }})\right]\\ & = & \displaystyle \sum _{j=1}^{{N}_{p}}\displaystyle \sum _{{\bf{Q}}}({V}_{{\bf{Q}}}^{p}({\bf{r}}){a}_{{\bf{Q}}}{e}^{i{\bf{Q}}\cdot {{\bf{R}}}_{j}}+{V}_{{\bf{Q}}}^{{p}^{* }}({\bf{r}}){a}_{{\bf{Q}}}^{\dagger }{e}^{-i{\bf{Q}}\cdot {{\bf{R}}}_{j}})\end{eqnarray} \tag{ 47 }$$

for N_p electron pairs, where ${V}_{{\bf{Q}}}^{p}({\bf{r}})=2{V}_{{\bf{Q}}}\cos \left(\tfrac{{\bf{Q}}\cdot {\bf{r}}}{2}\right)$ and R_j and r are the center-of-mass and relative coordinates of the electron pair, respectively.

In the calculations of the polaron energies, we will ignore the direct participation of the single-electron-electron-pair interaction potential v_sp(q). The potential v_sp(q) has a minor effect on the e-ph interaction. Its influence on the e-ph interaction is mostly indirect through the electronic screening and it is taken into account in the static structure-factor. Within such an approximation—ignoring the potential v_sp(q) in dealing with the e-ph interactions, the Hamiltonian of the whole system defined in equation (23) can be separated into two subsystems. They are the one consisting of single electrons interacting with phonons and the other of electron pairs interacting with phonons. In this way, the contribution of the electron-phonon interaction to the ground-state energy ${E}_{\mathrm{pol}}^{\mathrm{tot}}$ in equation (35) can be calculated as

$$\begin{eqnarray}&&{E}_{\mathrm{pol}}^{(\mathrm{tot})}={N}_{s}{E}_{\mathrm{pol}}^{\mathrm{single}}+{N}_{p}{E}_{\mathrm{pol}}^{\mathrm{pair}},\end{eqnarray} \tag{ 48 }$$

where ${E}_{\mathrm{pol}}^{\mathrm{single}}$ and ${E}_{\mathrm{pol}}^{\mathrm{pair}}$ are the polaron energies of the single electron and electron pair, respectively.

The subsystem composed of single electrons interacting with LO-phonons is

$$\begin{eqnarray}&&{H}_{\mathrm{pol}}^{\mathrm{single}}={H}_{\mathrm{single}}+{H}_{\mathrm{ph}}+{H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{single}},\end{eqnarray} \tag{ 49 }$$

where H_single is given by equation (25), H_ph by equation (31), and ${H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{single}}$ by equation (46). The e-ph coupling leads to a polaron consisting of an electron and a surrounding phonon cloud. When the e-ph interaction is not too strong, the polaron correction to the ground-state energy of an electron gas can be calculated within the Lee-Low-Pines (LLP) unitary transformation method [49, 50]. This method has been used to study polaron gases in bulk materials and also in low-dimensional systems [51, 52]. It is known that the polaron energy obtained from the LLP transformation is exact for $\alpha \to 0$. In the low electron density limit, the polaron energy obtained from the LLP method for e-ph coupling constant α = 6 is 90% of the exact value [47, 54]. Here we are dealing with a polaron gas in which the screening reduces the electron-phonon interaction strength. Therefore, the LLP method should yield a reasonable polaron energy for α < 6.

We calculate the polaron energy within the LLP method for the above Hamiltonian in equation (49), given by,

$$\begin{eqnarray}&&{E}_{\mathrm{pol}}^{\mathrm{single}}=-\displaystyle \sum _{{\bf{q}},{q}_{z}}\displaystyle \frac{| {V}_{{\bf{Q}}}{| }^{2}{S}_{{ss}}^{2}(q)}{{\hslash }{\omega }_{{}_{{LO}}}{S}_{{ss}}(q)+{{\hslash }}^{2}{q}^{2}/2{m}^{* }},\end{eqnarray} \tag{ 50 }$$

where S_ss(q) is the static structure factor of the single-electron gas. In the low electron-density limit, ${S}_{{ss}}(q)=1$. This leads to the well known perturbation result ${E}_{\mathrm{pol}}^{\mathrm{single}}=-(\pi /2)\alpha {\hslash }{\omega }_{{}_{{LO}}}$ for the polaron energy of an electron in 2D coupled with 3D-phonons [51, 68].

Figure 4(a) shows the polaron energy as a function of the single-electron density n_s in the 2D square-lattice crystal with λ = 1.5 a_B coupled with the 3D LO-phonons. The polaron energy in the low-density limit without screening is indicated by the thin-dotted line. It is seen that the screening of the electron gas considered in the structure factor S_ss(q) in equation (50) reduces the polaron effect. At higher electron density n_s = 1.0 ${{\rm{a}}}_{B}^{-2}$, the polaron energy is only about 20% of its low-density value. When there are also electron pairs in the crystal and the se-ep interaction is included in the structure factor S_ss(q) in the screening, the calculation results show that the se-ep interaction reduces the screening and consequently, enhances the polaron energy. The dashed curves in figure 4(a) are obtained with electron-pair density n_p = n_s/2 interacting with single-electrons. The red, green and blue dashed curves are for V₀ = 12, 15, and 18 R_y, respectively. The potential V₀ affects the potentials v_sp and v_pp. Its influences on the polaron energy is indirectly through the structure factor S_ss(q). Therefore, we obtained almost the same value of ${E}_{\mathrm{pol}}^{\mathrm{single}}$ for different V₀.

Zoom In Zoom Out Reset image size

The subsystem consisting of electron pairs and LO-phonons is given by the following Hamiltonian,

$$\begin{eqnarray}&&{H}_{\mathrm{pol}}^{\mathrm{pair}}={H}_{\mathrm{pair}}+{H}_{\mathrm{ph}}+{H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{pair}},\end{eqnarray} \tag{ 51 }$$

where H_pair is given by equation (26), H_ph by equation (31), and ${H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{pair}}$ by equation (47). Following a similar procedure as before and applying the LLP transformation to the electron-pair-phonon interactions, we obtain the polaron contribution to the energy of the electron pairs, given by

$$\begin{eqnarray}&&{E}_{\mathrm{pol}}^{\mathrm{pair}}=-\displaystyle \sum _{{\bf{q}},{q}_{z}}\displaystyle \frac{{| {M}_{\mathrm{pair}}({\bf{Q}})| }^{2}{S}_{{pp}}^{2}(q)}{{\hslash }{\omega }_{{}_{{LO}}}{S}_{{pp}}(q)+{{\hslash }}^{2}{q}^{2}/4{m}^{* }},\end{eqnarray} \tag{ 52 }$$

where S_pp(q) is the static structure factor of the electron pairs. The matrix element is given by

$$\begin{eqnarray}&&{| {M}_{\mathrm{pair}}({\bf{Q}})| }^{2}=| \langle {{\rm{\Psi }}}_{{\bf{k}}-{\bf{q}}}({\bf{R}},{\bf{r}});{\bf{Q}}| {H}_{{\rm{e}} \mbox{-} \mathrm{ph}}^{\mathrm{pair}}| {{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}});0\rangle {| }^{2},\end{eqnarray} \tag{ 53 }$$

where $| {{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}});{\bf{Q}}\rangle =| {{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}})\rangle | {\bf{Q}}\rangle$, with $| {{\rm{\Psi }}}_{{\bf{k}}}({\bf{R}},{\bf{r}})\rangle$ for the electron-pair state and $| {\bf{Q}}\rangle$ for the phonon state. Using the metastable electron-pair wavefunction in equation (22), we obtain

$$\begin{eqnarray}&&{| {M}_{\mathrm{pair}}({\bf{Q}})| }^{2}=4| {V}_{{\bf{Q}}}{| }^{2}{f}_{{pp}}(q),\end{eqnarray} \tag{ 54 }$$

where V_Q is the e-ph interaction potential given by equation (33), f_pp(q) is the form-factor that appears in the pair-pair interaction potential in equation (26).

We now obtain a Fröhlich bipolaron formed by an electron pair coupled with LO-phonons. Figure 4(b) shows the bipolaron energy in the crystals with λ = 1.5 a_B and V₀ = 12, 15, and 18 R_y. The solid curves give the bipolaron energies for different V₀ with only screening of the electron pairs. For larger V₀, the average distance $\langle r\rangle$ between the two electrons in a pair is smaller, i.e., the size of the electron pair is smaller. Therefore, the bipolaron energy is larger. But the screening reduces significantly the bipolaron energy. We also see that, when se-ep interaction is included in the structure factor S_pp(q), the screening of the electron-pair-phonon coupling is reduced enhancing the bipolaron energy as shown by the dashed curves for n_s = 2n_p in the figure.

If we assume ${S}_{{pp}}(q)=1$, we obtain the low-density limit of the bipolaron energy without screening given by the horizontal dotted-dash lines. This is the case of a single bipolaron without screening. If we further assume f_pp(q) = 1, we obtain the upper bound of the bipolaron energy, ${E}_{\mathrm{pol}}^{\mathrm{pair}}=-4\sqrt{2}(\pi /2)\alpha {\hslash }{\omega }_{{}_{{LO}}}$ [69]. This is equivalent to assuming the electron pair as a 'larger single electron' with charge −2e and mass 2m* coupled to the LO-phonons.

As a matter of fact, the effects of electron-phonon interaction and possible bipolaron formation have been extensively studied as the electron pairing mechanism for unconventional superconductivity [47, 70–72]. In order to form a bipolaron in the crystal, a crucial point is that the electron-phonon coupling induced attraction has to overcome the electron-electron Coulomb repulsion [47]. This requires not only a small ratio ${\epsilon }_{\infty }/{\epsilon }_{0}$ of dielectric constants but also a large enough electron-phonon coupling constant α leading to a critical e-ph coupling constant α_c = 2.9 in 2D and α_c = 6.8 in 3D for bipolaron formation.

However, the bipolaron formation mechanism in the present theory is distinct from the traditional bipolaron theory in the literature. In this paper, we show a preformed metastable electron pair due to strong correlation of two electrons occupying the same orbital. The electron-phonon interaction dresses the electron pairs up with a phonon cloud forming bipolarons. Therefore, in the present context we introduce an internal interaction due to orbital dependent electron correlation of the electron pair. The electron-phonon coupling involves in the renormalization of the preformed electron pairs. Notice that polaron effects can be much larger on the electrons in the correlated pairs than on the single electrons, and therefore the bipolaron contribution to the stabilization of the electron pairs overcoming the Coulomb repulsion becomes essential. We will show in the next section that the certain electron pairs can indeed be stabilized as bipolarons.

The condition for stabilization of the electron pairs when including renormalization is that the bottom of the electron-pair band occurs at the Fermi energy of the single-electron band. Using the same energy reference defined in section 3, i.e., taking the bottom of the single-electron band at Γ point before the renormalization as E = 0, the many-particle interactions lower the single-electron band bottom to

$$\begin{eqnarray}&&{E}_{{\rm{s}}}^{{\rm{\Gamma }}}={E}_{{ss}}+{E}_{{sp}}+{E}_{\mathrm{pol}}^{\mathrm{single}},\end{eqnarray} \tag{ 55 }$$

where E_ss and E_sp are determined by equation (36) and ${E}_{\mathrm{pol}}^{\mathrm{single}}$ is given by equation (50). The Fermi energy of the single-electron band can be obtained by ${E}_{F}={E}_{{\rm{s}}}^{{\rm{\Gamma }}}+{\varepsilon }_{F}$ where ε_F is the energy difference between the Fermi energy and the bottom of the band. It is determined by the single-electron density and the density of states of the band. We will not consider the many-body effects on ε_F. The average kinetic energy per single electron in this band is given by ε_F /2 in a two-dimensional system.

Including the band renormalization, the bottom of the electron-pair band is given by

$$\begin{eqnarray}&&{E}_{{\rm{p}}}^{{\rm{\Gamma }}}={\nu }_{0}+\displaystyle \frac{{}_{1}}{{}^{2}}({E}_{{pp}}+{E}_{{sp}}+{E}_{\mathrm{pol}}^{\mathrm{pair}}),\end{eqnarray} \tag{ 56 }$$

where E_pp and E_sp are given by equation (36) and ${E}_{\mathrm{pol}}^{\mathrm{pair}}$ by equation (52). Notice that the energy of an electron pair is given by $2{E}_{{\rm{p}}}^{{\rm{\Gamma }}}$.

Figure 5 shows the energies for the renormalized single-electron and electron-pair bands in the 2D crystal with λ = 1.3 a_B and V₀ = 15 R_y coupled with the 3D LO-phonons. The band structure without many-body corrections was given in figure 2(a). The many-particle interaction energies are obtained for a constant total electron density n_t = n_s + 2n_p = 0.5 ${{\rm{a}}}_{B}^{-2}$ with different e-ph coupling constant α. Because the polaron energy is given by the coupling constant α and the LO-phonon energy ${\hslash }{\omega }_{{LO}}$, the ratio ${{\rm{R}}}_{y}/{\hslash }{\omega }_{{LO}}$ is required in the calculations [73]. This ratio is a material dependent parameter and has been used in the bound-polaron and bound-bipolaron problems [73, 74]. Its value is in the range from 0.5 to 2 for different materials [74, 75]. For the present discussion, we take ${{\rm{R}}}_{y}/{\hslash }{\omega }_{{LO}}$ = 1.

Zoom In Zoom Out Reset image size

Figure 5(a) shows for fixed total electron density n_t = 0.5 ${{\rm{a}}}_{B}^{-2}$ and e-ph coupling constant α = 4, the dependence of the single-electron band bottom ${E}_{{\rm{s}}}^{{\rm{\Gamma }}}$, its Fermi energy E_F, and the electron-pair band bottom ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ on the electron distribution in the two bands n_s and n_p. The total energy per electron E_tot is also given in the figure. It is the average ground-state energy per electron in the system calculated by equations (35), (41), and (48), and is given by

$$\begin{eqnarray}{E}_{\mathrm{tot}} & = & \displaystyle \frac{{E}_{\mathrm{GS}}}{{N}_{t}}=\displaystyle \frac{{E}_{\mathrm{GS}}^{(\mathrm{el})}+{E}_{\mathrm{pol}}^{(\mathrm{tot})}}{{N}_{t}}\\ & = & \displaystyle \frac{{n}_{s}}{{n}_{t}}\left({E}_{{\rm{s}}}^{{\rm{\Gamma }}}+\displaystyle \frac{{\varepsilon }_{F}}{2}\right)+2\displaystyle \frac{{n}_{p}}{{n}_{t}}{E}_{{\rm{p}}}^{{\rm{\Gamma }}},\end{eqnarray} \tag{ 57 }$$

where ${E}_{{\rm{s}}}^{{\rm{\Gamma }}}$ and ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ are given by equations (55) and (56) respectively. It is seen that, though the total electron density n_t is a constant, the above obtained energies depend on the distribution of the electrons between the two bands. Especially the energy of the electron-pair band ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ depends strongly on n_p. In figure 5(a) we see that with increasing the single-electron density n_s, the single-electron band ${E}_{{\rm{s}}}^{{\rm{\Gamma }}}$ and the Fermi energy E_F vary slowly. Because n_t is a constant, increasing n_s means decreasing n_p. At low density n_p, the weak screening in the electron-pair band enhances greatly the bipolaron energy resulting in a lowering of the bottom of the electron-pair band ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ which reaches the Fermi surface. We find that under the considered condition, ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ touches the Fermi surface at n_s = 0.4774 ${a}_{B}^{-2}$ and n_p = 0.0113 ${a}_{B}^{-2}$ as indicated by the black dot. For these densities, the electron pairs (or the bipolarons) are stabilized on the Fermi surface. It means that only 4.5% of the electrons in the system form stable electron pairs in this case. If we look at the total energy E_tot, it tends to reduce the energy of the system at the higher single-electron density side. But the electron pairs on the Fermi surface cannot decay into two single electrons due to the Pauli exclusion principle. Moreover, breaking an electron pair needs a cost to overcome their correlation energy. Therefore, the electron-pair density is stabilized at n_p = 0.0113 ${a}_{B}^{-2}$ in the ground state of the system with α = 4 and n_t = 0.5 ${a}_{B}^{-2}$. Figure 5(b) shows the electron-pair band energy ${E}_{{\rm{p}}}^{{\rm{\Gamma }}}$ and the single-electron band Fermi energy E_F in the same system with n_t = 0.5 ${a}_{B}^{-2}$ but for different e-ph coupling constant α. We see that for α = 2 no stable electron pairs are found. The calculations indicate that for $\alpha \gtrsim 3$, part of the electrons can form electron pairs on the Fermi surface at low electron-pair density. The stabilized electron-pair densities are n_p = 0.001, 0.011, and 0.026 ${{\rm{a}}}_{B}^{-2}$ for α = 3, 4 and 5, respectively.

In order to determine the electron-pair density in the ground-state of the system, we solve the following equation

$$\begin{eqnarray}&&{E}_{{\rm{p}}}^{{\rm{\Gamma }}}({n}_{s},{n}_{p})={E}_{F}({n}_{s},{n}_{p}),\end{eqnarray} \tag{ 58 }$$

as a function of n_s and n_p for fixed λ, V₀ and α. Figures 6(a) and (b) show the density of stabilized electron-pairs or bipolarons in linear and logarithmic scale, respectively, as a function of the single-electron density in the system with λ = 1.3 a_B and V₀ = 15 R_y for α = 2.8, 3.0, 3.5, 4.0, 4.5, and 5.0. We see that for α = 2.8 a very low density of electron pairs become stabilized. With increasing α, more electron pairs appear in the system. For α = 5, their density may reach n_p = 0.032 ${{\rm{a}}}_{B}^{-2}$. This corresponds to about 11% of the electrons in the system are in the electron-pair band. The thick-red curve in the figure is the possible maximum electron-pair density ${n}_{p}^{\max }$ given by equation (30). It is important to notice that the relation between ${n}_{p}^{\max }$ and n_s not only restrict the electron-pair density being less than ${n}_{p}^{\max }$. It also means that if n_s is larger than λ⁻², i.e., the single-electron band filling is more than half, there are no stable electron pairs in the system. For the crystal with λ = 1.3 a_B, the half filling of the single-electron band occurs at n_s = 0.5917 ${{\rm{a}}}_{B}^{-2}$.

Zoom In Zoom Out Reset image size

In figure 6(c) we show the electron-pair density in a different 2D potential with λ = 2.8 a_B and V₀ = 5 R_y. In this case, the band-gap ${E}_{g}^{(0)}$ and also ν₀ are small, as shown in figure 2. The half filling of the single-electron band is at n_s = 0.1276 ${{\rm{a}}}_{B}^{-2}$. We find that for $\alpha \simeq 2$, some electron pairs can already be stabilized with small density. For α = 5, the electron-pair density may reach 0.02 ${{\rm{a}}}_{B}^{-2}$ corresponding to 31% of the total electrons being in the paired state.

The starting point of our work is that the electron-electron correlation is orbital dependent. The exchange-correlation energy of an electron pair occupying the same orbital is larger than the average of the exchange-correlation energy of many electrons. Depending upon the crystal structure and potential, such electron pairs can form a metastable electron-pair band. The metastable electron-pair band is obtained in two different ways. One follows the tight-binding method and the other is similar to the nearly-free electron model. They give consistent results for the dispersion relation of the electron-pair states in the period potential.

Combining the metastable electron-pair band together with the corresponding single-electron band, we constructed a many-particle system consisting of single electrons and metastable electron pairs. When further considering many-body interaction renormalization including single electrons, electron pairs, and optical phonons, we found that the metastable electron pairs can be stabilized. The calculations show that the polaron effects play an essential role in counterbalancing the Coulomb repulsion in the stabilization of the electron pairs. The electron-phonon coupling is enhanced in the strongly correlated electron pairs and leads to the formation of bipolarons. On the other hand, screening affects significantly both the polaron and bipolaron energies manifesting a cooperative interplay of electron-electron and electron-phonon interactions. The obtained ground state of the system consists of polarons in the single-electron band and bipolarons in the electron-pair band sitting on the Fermi surface of the single-electron band.

The numerical calculations presented in this paper are performed for simple potentials and with the s-orbital in 2D square-lattice crystals. However, the physical processes and numerical calculations for electron pairing can be extended to a quasi-2D crystal of a single atomic layer with finite thickness or to 3D systems. In principle, it can also be extended to study the electron pairing in p- and d-orbitals. The obtained results within the present simple model predict light and small electron pairs (bipolarons) in the crystal. In the center-of-mass coordinates, the electron pair is a Bloch wavefunction with an effective mass of twice a single electron. The pair is small and local because the average separation between the two electrons is less than half of the lattice constant and they are localized in the same unit cell when expressed in the relative coordinates. Furthermore, only a fraction of the electrons in the system form pairs.

We have obtained a many-particle ground state with spin-singlet electron pairs in the condensate phase. We expect that these preformed electron pairs at certain densities in coherent state will contribute to superconductivity. Finally, we want to comment on the possible 'superconducting energy gap'. From our calculations we naturally infer such a gap to the transition energy required to break the stabilized electron pairs sitting on the Fermi surface. This transition is determined by the Hamiltonian ${H}_{\mathrm{int}}^{{\rm{t}}}$ in equation (29) where the potential ${v}^{{\rm{t}}}(q)$ was given by equation (11) in [48]. Because the electron pairs are stabilized on the Fermi surface in the center of the Brillouin zone at k = 0, an external excitation has to overcome the internal correlation energy of an electron pair to bring two electrons to the single-electron states above the Fermi surface at $\pm {\bf{k}}\ne 0$. The minimum energy required (or the gap) is primarily determined by the potential ${v}^{{\rm{t}}}(q)$ and the final momenta $\pm {\hslash }{\bf{k}}$ of the single-electron states. It is band structure dependent, anisotropic, and nodeless.

This work was supported by the Brazilian agencies FAPESP and CNPq. GQH would like to thank Prof. Bangfen Zhu for his invaluable support and expert advice.