Vortex spectroscopy in the vortex glass: A real-space numerical approach

A method is presented to solve the Bogoliubov–de Gennes equations with arbitrary distributions of vortices. The real-space Green's function approach based on Chebyshev polynomials is complemented by a gauge transformation which allows one to treat finite as well as infinite, ordered as well as disordered vortex configurations. This tool gives unprecedented access to vortex lattices at very low magnetic fields and glassy phases. After describing in detail the method and its implementation, we use it to address a series of problems related to $d$-wave superconductivity on the square lattice. We first study the continuity of the vortex-core energy spectrum and its evolution from the quantum regime to the semiclassical limit; we investigate the effect of the band structure on the vortex by following the self-consistent solution through a Lifshitz transition; we then study the evolution from the vortex lattice to the isolated-vortex limit with decreasing field and show that a new emerging length scale controls this transition; finally, we perform a statistical study of the vortex-core local density of states in the presence of positional disorder in the vortex lattice. The calculations reveal a number of qualitative differences between the properties of vortices in the quantum and semiclassical regimes.

Figure 1

Local phase and nonlocal phase corrections for an isolated vortex. The big and small dots indicate the vortex center and lattice sites, respectively. (a) Local phase and branch cut. The phase is $+\pi$ on the cut and zero at the vortex center. (b) Core correction for a vortex sitting on a lattice site in a square lattice. (c) and (d) Long-range correction for bonds touching the branch cut.

Figure 2

(a) Function ${\mathrm{\Phi }}_{\widehat{\mathit{y}}}\left(\mathit{r}\right)$ for an ideal square vortex lattice oriented along the (11) direction. The vortex centers are located on the nodes of the underlying square lattice. An expanded view of a vortex core is displayed in the inset. (b) Function ${\mathrm{\Phi }}_{\widehat{\mathit{x}}}\left(\mathit{r}\right)$ for a disordered vortex configuration embedded in a square vortex lattice oriented along the (10) direction. The vortex centers are in the plaquettes of the underlying square lattice.

Figure 3

(a) Diamond-like system of size $M=3$ with $1+2M(1+M)$ sites on the square lattice. This shape is optimal if the range of the Hamiltonian extends to the first, third, or fifth neighbors. We typically use $M=500$, corresponding to a system of 501 001 sites. (b) Square-like system of size $M=2$ with ${(2M+1)}^2$ sites. This shape is optimal if the range of the Hamiltonian extends to the second or fourth neighbors. We typically use $M=350$ in this case, corresponding to 491 401 sites. A possible sequential numbering of the sites is indicated.

Figure 4

(a) DOS calculated with $\mathfrak{a}=1.1{\mathfrak{a}}_N,\mathfrak{b}={\mathfrak{b}}_N$ (black, blue shade) and $\mathfrak{a}={\mathfrak{a}}_{\mathrm{\Delta }},\mathfrak{b}=0$ (red). The dashed curve is the image of the DOS by particle-hole symmetry. The inset shows the loss of resolution in the gap region for $\mathfrak{a}={\mathfrak{a}}_{\mathrm{\Delta }}$. The other parameters are $M=350$ and $N=4M$. (b) DOS in the gap region for various system sizes and expansion orders. (c) DOS at the lower band edge showing Gibbs oscillations and their removal by various kernels (see Ref. [44]).

Figure 5

Convergence of the self-consistent gap as a function of (a) system size and (b) expansion order and truncation kernel. Inside the gray shaded regions the convergence is better than 0.1%. The calculations are done with $\mathfrak{a}=12|t_1|$ and $\mathfrak{b}=0$ for the same model as in Fig. 4.

Figure 6

LDOS at the vortex center for the half-filled square lattice with nearest-neighbor hopping $t_1$, second-neighbor hopping $t_2=0$, and various values of ${\mathrm{\Delta }}_0$ for (a) $d$-wave pairing and (b) $s$-wave pairing. The thick horizontal bars denote the energy resolution (which is identical in all panels). The dashed vertical lines in (b) show multiples of ${\mathrm{\Delta }}_0/t_1$, indicating the typical interlevel spacing for bound states; the long (short) solid vertical lines mark the actual bound-state energies with (without) weight at the vortex center. All calculations use the system geometry of Fig. 3 for $M=1000,N=4M$ with the Jackson kernel, $\mathfrak{a}=8t_1$, and $\mathfrak{b}=0$.

Figure 7

Fermi surface and low-energy LDOS [58] for a $d$-wave vortex at three values of ${\mathrm{\Delta }}_0/t$ on the square lattice with (a) $t_1=t,t_2=\mu =0$ and (b) $t_2=t,t_1=\mu =0$. The vortex-core size is $t/{\mathrm{\Delta }}_0$ in units of the lattice parameter, corresponding to a single pixel in each image. The color scale is logarithmic going from minimum LDOS (black) to maximum LDOS (white).

Figure 8

(a) Evolution of the Fermi surface with varying chemical potential $\mu$ for a square lattice with $t_1=-250$ meV and $t_2=75$ meV. (b) Evolution with $\mu$ of the average Fermi velocity (left scale) and nearest-neighbor interaction needed to produce a $d$-wave gap of 40 meV on the Fermi surface (right scale). The velocity is calculated for a lattice parameter $a=3.8\AA$. The vertical line indicates the Lifshitz transition. (c) Evolution of the vortex parameters ${\xi }_0,{\xi }_1$, and ${\xi }_c={\xi }_{1}W({\xi }_0/{\xi }_1)$ (left scale). The length $\hslash v_{\mathrm{F}}/\left(\pi \mathrm{\Delta }\right)$ is also indicated for comparison. The squares (right scale) show the maximum amplitude of the $s^{*}$-wave component. The surface plots show, for three values of the chemical potential and for each node of the microscopic square lattice (black lines), the self-consistent order parameter modulus (blue), the difference between the model (8) and the self-consistent data (orange), and the $s^{*}$-wave component in one quarter of the field of view (green). The lattice grid also marks the zero of the vertical scale, such that positive (negative) values of the difference shown in orange appears above (below) the lattice.

Figure 9

(a) Spatial and (b), (c) spectral distributions of the LDOS for open hole-like Fermi surface $\left(\mu =0\right)$, nearly square Fermi surface ($\mu =-275$ meV), and closed electron-like Fermi surface ($\mu =-500$ meV). The first two columns in (a) show the same data on logarithmic (log) and linear (lin) scales. The scale for the other maps is linear. All maps show the same region of size $41a\times 41a$. (b) and (c) display spectral traces running from the core along the (10) and (11) directions, respectively; spectra are shifted vertically. The solid lines show the LDOS for the self-consistent order parameter; the red-dashed lines in (b) show the LDOS calculated with the isotropic model (8) using the fitted parameters; the blue-dashed lines in (c) show the LDOS calculated using (8) with the fixed values ${\xi }_0=1.3a$ and ${\xi }_1=14a$. Only one quarter of the dashed spectra is shown for clarity. All calculation were performed with $M=350,N=4M,\mathfrak{a}=2.6$ eV, and $\mathfrak{b}=0.3\mathrm{eV}-\mu$.

Figure 10

Insets in (a). Parameters ${\xi }_0$ and ${\xi }_1$ as a function of intervortex distance $d$ in the $0^{\circ }$ (black squares) and ${45}^{\circ }$ (red diamonds) vortex lattices. The solid lines show exponential behavior with the characteristic lengths indicated in parentheses. (b) LDOS in the core and (c) in-between vortices for increasing field (top to bottom, from $d/a=100$ to $d/a=10$, shifted vertically) in the $0^{\circ }$ (black solid lines) and ${45}^{\circ }$ (red dashed lines) vortex lattice. The main graph (a) shows how the correlation between the spectra for the $0^{\circ }$ and ${45}^{\circ }$ orientations approaches unity with decreasing field in the core (filled circles) and in-between vortices (empty circles). The magnetic field $B={\mathrm{\Phi }}_0/d^2$ is shown for a lattice parameter $a=3.8\AA$. The self-consistent calculations were performed with $M=100,N=10000,\mathfrak{a}=3.1$ eV, $\mathfrak{b}=0$; the LDOS calculations with $M=500,N=4M,\mathfrak{a}=2.6$ eV, $\mathfrak{b}=0.55$ eV.

Figure 11

Spatial evolution of the LDOS along the antinodal (filled black symbols) and nodal (empty red symbols) directions at various energies for an isolated vortex. The LDOS is normalized to its value in the absence of vortex. The vertical lines indicate the characteristic lengths ${\xi }_0,{\xi }_1$, and $\ell$.

Figure 12

(a) and (b) Average vortex-core spectra in the disordered $0^{\circ }$ (solid black) and ${45}^{\circ }$ (dashed red) vortex lattices for (a) $d=50a$ and (b) $d=30a$ and for three disorder strengths. In each case, the average is calculated from a distribution of 3630 spectra. The standard deviation of the distribution is displayed as gray area and light red lines for the $0^{\circ }$ and ${45}^{\circ }$ orientations, respectively. (c) and (d) Correlation between the average spectra for both orientations as a function of disorder strength. Linear and exponential behaviors are shown by solid lines on top of the calculation (circles).

Figure 13

(a) Zero-energy LDOS over a $200a\times 200a$ region in a disordered ${45}^{\circ }$ vortex lattice with $d=50a$ and $\eta =0.3$. The energy-dependent LDOS calculated at the vortex centers is displayed next to each vortex. (b) Zoom of the central region showing spectral traces across two vortex cores at the positions indicated in the insets. Red crosses mark the vortex centers (phase singularity points). The distance between the two centers is $22a$. (c) Superfluid current density (color scale and black arrows) in the same region as (b) and resultant Lorentz force on each vortex (green arrows). The force is the average of $\mathit{j}\times \mathit{B}$ inside the dashed circles of radius $11a$.