Bogoliubov–de Gennes Theory in Inhomogeneous Systems

Bogoliubov–de Gennes (BdG) theory turns a spatially varying pairing saddle into a particle–hole-redundant eigenproblem. Given a regulated one-body Hamiltonian, boundary conditions, and either a prescribed or self-consistent pair field, it predicts normalized quasiparticle wavefunctions and local observables such as the tunnelling density of states. Its conclusions remain mean-field conclusions: cutoff dependence, finite size, charging physics, fluctuations, and disorder averaging must be controlled before a computed subgap state is assigned a unique physical identity.

Required background. Nambu–Gor’kov propagators fixes the basis, signs, gauge covariance, and double-counting rule.

Helpful background. Basis construction and symmetry sectors supplies finite-dimensional regulator checks.

The real-space eigenproblem

For a spin-singlet continuum model in the basis $(\psi_\uparrow,\psi_\downarrow^\dagger)^T$, take

$$H_{\mathrm{BdG}}(\mathbf r)= \begin{pmatrix} h_0(\mathbf r)-\mu & \Delta(\mathbf r)\\ \Delta^\ast(\mathbf r)&-[h_0(\mathbf r)-\mu]^\ast \end{pmatrix},$$

with $h_0=[-i\boldsymbol\nabla-e\mathbf A(\mathbf r)]^2/(2m)+V(\mathbf r)$ in this local convention. The eigenvectors satisfy

$$H_{\mathrm{BdG}}\begin{pmatrix}u_n\\v_n\end{pmatrix} =E_n\begin{pmatrix}u_n\\v_n\end{pmatrix}, \qquad \int\mathrm d^d r\,(\lvert u_n\rvert^2+\lvert v_n\rvert^2)=1.$$

For this reduced even-parity spin-singlet block, particle–hole symmetry is represented by $\mathcal C=\tau_yK$. It gives $(v_n^\ast,-u_n^\ast)^T$, up to an overall phase, at $-E_n$. This reduced representation has $\mathcal C^2=-1$; a full spin-orbital Nambu basis has its corresponding particle–hole representation. A numerical spectrum that lacks the required pairing beyond its solver tolerance has a basis, boundary, or assembly error.

The local density of states per spin for this spin-degenerate reduced block, summed only over $E_n>0$, is

$$\rho(\mathbf r,E)=\sum_{E_n>0} \left[ \lvert u_n(\mathbf r)\rvert^2\delta(E-E_n) +\lvert v_n(\mathbf r)\rvert^2\delta(E+E_n) \right].$$

Broadening $\delta(E)\to\eta/[\pi(E^2+\eta^2)]$ is an analysis choice. Report $\eta$, level spacing, and experimental resolution separately; otherwise numerical smearing can be mistaken for a lifetime.

The spin-summed density of states is twice this expression when spin degeneracy is intact. With spin–orbit coupling, magnetism, or a larger orbital basis, compute the electron-sector trace in the full BdG basis instead of inserting that factor by hand.

Self-consistency, gauge covariance, and boundaries

For a local attraction, the regulated self-consistency condition is

$$\Delta(\mathbf r)=g\sum_{0<E_n<E_c} u_n(\mathbf r)v_n^\ast(\mathbf r) \bigl[1-2f(E_n)\bigr].$$

The pair $(g,E_c)$ must be matched to a physical interaction or retained as a declared model regulator. Increasing the basis at fixed bare $g$ generally does not define a continuum limit. Convergence requires varying spatial resolution, box size, $E_c$, and iteration tolerance while holding matched low-energy data fixed.

de Gennes 1999, ch. 5 develops the inhomogeneous quasiparticle equations and their self-consistency conditions.

Under a gauge transformation with phase $\alpha(\mathbf r)$,

$$u_n\mapsto e^{i\alpha}u_n, \quad v_n\mapsto e^{-i\alpha}v_n, \quad \Delta\mapsto e^{2i\alpha}\Delta, \quad \mathbf A\mapsto\mathbf A+\boldsymbol\nabla\alpha/e.$$

The spectrum and $\rho$ are invariant only when hopping phases, vector potential, and pair phase transform together. Setting $\mathbf A=0$ while imposing a spatially varying phase is a gauge choice that can change boundary twists; it is not automatically equivalent to every magnetic configuration.

Hard-wall, periodic, lead self-energy, and open-system boundary conditions answer different questions. A bound state in a closed finite box may broaden or disappear when coupled to a continuum. Conversely, a finite wire always hybridizes nominal end modes with a splitting controlled by length, coherence scale, and oscillatory wavefunction overlap.

Interface and vortex example

Across a slowly varying normal–superconducting interface, $\Delta(x)$ converts electron-like and hole-like components through Andreev reflection. In the Andreev regime $E,\lvert\Delta\rvert\ll E_F$, factoring out the rapid Fermi wavelength yields

$$\begin{pmatrix} -iv_F\partial_s & \Delta(s)\\ \Delta^\ast(s)&iv_F\partial_s \end{pmatrix} \begin{pmatrix}u\\v\end{pmatrix}=E\begin{pmatrix}u\\v\end{pmatrix}.$$

This controlled linearization fails near a band edge, in a low-density region, or when the texture changes on the scale $k_F^{-1}$. For a vortex $\Delta(\mathbf r)=\Delta(r)e^{i\varphi}$, ordinary Caroli–de Gennes–Matricon levels have spacing of order $\Delta^2/E_F$; their near-zero energy is not topological protection. Caroli, de Gennes, and Matricon 1964, pp. 307–309 derives the core spectrum, while Gygi and Schlüter 1991, §§II–III gives a self-consistent vortex calculation.

The structure diagram identifies exactly where a BdG eigenvalue sits in the chapter’s inference chain.

BdG theory maps an inhomogeneous saddle to spectra and local amplitudes after regulator, gauge, and boundary conditions are fixed. A subgap eigenvalue alone does not identify its origin. Original schematic, not to scale.

The paired-matter claim test matrix lists the independent tests needed for vortex, topological, and platform claims.

Numerical stop rules

A BdG result is not ready for interpretation until:

• positive and negative eigenvalues pair within the eigensolver residual;
• each eigenvector and the integrated local spectral weight satisfy normalization sum rules;
• the claimed low-energy feature is stable under basis, mesh, box, cutoff, and broadening changes;
• self-consistent solutions are compared from multiple initial conditions and by free energy, not iteration convergence alone;
• gauge-related inputs give the same gauge-invariant outputs; and
• disorder claims include an ensemble and sample-to-sample distribution rather than one favorable realization.

Exercise

Show particle–hole pairing. For the reduced singlet block $H=\xi\tau_z+\Delta_x\tau_x-\Delta_y\tau_y$, verify that $\tau_yK$ maps an eigenstate at $(\mathbf k,E)$ to one at $(-\mathbf k,-E)$ when $\xi_{-\mathbf k}=\xi_{\mathbf k}$ and $\Delta(-\mathbf k)=\Delta(\mathbf k)$.

Solution

Complex conjugation reverses the sign of the matrix $\tau_y$. Conjugation by $\tau_y$ changes the signs of $\tau_z$ and $\tau_x$ while preserving $\tau_y$, so the stated even-parity conditions give $\tau_yH^\ast(-\mathbf k)\tau_y=-H(\mathbf k)$. Applying $\tau_yK$ to $H(\mathbf k)\lvert u\rangle=E\lvert u\rangle$ therefore yields the partner at $-E$. In components, $\tau_y(u^\ast,v^\ast)^T$ is phase-equivalent to $(v^\ast,-u^\ast)^T$. Odd-parity or larger spin-orbital bases require their corresponding particle–hole operator; reusing a matrix from a different Nambu convention is an error.

References

• Caroli, C., de Gennes, P. G., and Matricon, J. (1964). “Bound fermion states on a vortex line in a type II superconductor.” Physics Letters 9, 307–309. doi:10.1016/0031-9163(64)90375-0.
• de Gennes, P. G. (1999). Superconductivity of Metals and Alloys. Westview Press. Publisher record.
• Gygi, F., and Schlüter, M. (1991). “Self-consistent electronic structure of a vortex line in a type-II superconductor.” Physical Review B 43, 7609–7621. doi:10.1103/PhysRevB.43.7609.