# Density of states in coupled chains with off-diagonal disorder

Since the works of Dyson, it has been known that one dimensional systems with off-diagonal disorder may show anomalous behavior in their density of states and their localization properties . This anomalous behavior has reappeared in many reincarnations, ranging from Dyson’s original application of masses coupled with springs, to the random hopping model , quantum XY spin chains , supersymmetric quantum mechanics , one-dimensional classical diffusion in a random medium , narrow-gap semiconductors , and charge density wave materials . In all these systems, it has been shown that both the density of states and the localization length diverge at zero energy. These two divergencies are closely linked by a theorem by Thouless .

Recently, it was found that the localization length of a wire consisting of $`N`$ of such chains shows a remarkable dependence on the parity of $`N`$ : At zero energy, i.e. at the center of the band, the localization length diverges for odd $`N`$ (odd $`N`$ includes the pure one-dimensional case $`N=1`$), while wavefunctions remain localized for even $`N`$. Moreover, it was found that staggering in the hopping parameter (as the consequence of a Peierls instability), may lead to an additional set of delocalized states for all $`N>1`$, regardless of parity . For comparison, in the pure one-dimensional case $`N=1`$ staggering always enhances localization. In this paper we consider the density of states (d.o.s.) for the $`N`$-chain wire with off-diagonal disorder. Unlike in the one-dimensional case, this is a problem that requires separate attention, since the Thouless theorem which links d.o.s. and localization length does not hold for $`N>1`$ coupled chains. Moreover this problem is also of relevance for the random flux model, which is a special example of off-diagonal randomness, and whose localization properties and density of states near $`\epsilon =0`$ are the subject of an ongoing debate .

To be specific, we consider the d.o.s. for the Schrödinger equation

$$\epsilon \psi _n=\psi _n=t_n^{}\psi _{n+1}t_{n1}^{}\psi _{n1},$$
(1)   
where $`\psi _n`$ is an $`N`$-component wavefunction and $`t_n`$ an $`N\times N`$ hopping matrix. The index $`n=1,\mathrm{},L`$ labels the site index along the wire. The length of the wire, measured in units of the lattice spacing, is $`L`$. The reason why this system can display anomalous behavior near zero energy is the existence of a particle-hole or chiral symmetry, which is absent in the presence of on-site disorder: Under a mapping $`\psi _n(1)^n\psi _n`$, the Hamiltonian $``$ changes sign, $``$. As a result, the eigenvalues of $``$ occur in pairs $`\pm \epsilon `$. There are two mechanisms by which the chiral symmetry is known to affect the d.o.s. near $`\epsilon =0`$. First, level repulsion of the eigenvalue $`\epsilon `$ with its mirror image $`\epsilon `$ causes a universal suppression of the d.o.s. near $`\epsilon =0`$ . This suppression appears on the scale of a level spacing, is independent of the geometry, but becomes unimportant in the thermodynamic limit. The second mechanism, for which several different descriptions exist , is special for a (quasi) one-dimensional geometry, and survives in the thermodynamic limit. It is responsible for the divergence of the d.o.s. in the pure one-dimensional case. Below we present a calculation of the d.o.s. in the multi-chain case, combining ideas from the Fokker-Planck approach to localization in multichannel quantum wires  and the calculation of the d.o.s. in the one-dimensional random hopping model .

We first state our main results. We find that the parity dependence that was previously obtained for the localization properties is also present in the density of states. For odd $`N`$, the d.o.s. diverges at zero energy according to

$$\rho (\epsilon )\frac{w^2}{|\epsilon \mathrm{ln}^3(w^2/\epsilon )|},N\text{ odd},$$
(2)   
where $`w`$ is a dimensionless parameter governing the randomness in the $`t_n`$. The bandwidth is chosen as the unit of energy. In Eq. (2) and in the remainder of the paper we assume $`\epsilon >0`$. The d.o.s. for $`\epsilon <0`$ then follows from $`\rho (\epsilon )=\rho (\epsilon )`$. The form of the divergence (2) is independent of whether time-reversal symmetry is broken or not, i.e. of whether the hopping matrices $`t_n`$ are generically real or complex. For even $`N`$, in contrast, the density of states strongly depends on the presence or absence of time-reversal symmetry (labeled by the parameter $`\beta =1`$ or $`2`$, respectively),

$$\rho (\epsilon )|\epsilon /w^2|^{\beta 1}|\mathrm{ln}(w^2/\epsilon )|,N\text{ even}.$$
(3)   
In the presence of time-reversal symmetry, $`\rho (\epsilon )`$ shows a logarithmic divergence as $`\epsilon 0`$, while in the absence of time-reversal symmetry a pseudo-gap is opened, $`\rho (\epsilon )`$ vanishes at $`\epsilon =0`$. Such a strong dependence on time-reversal symmetry is a remarkable result. An effect of comparable magnitude appears in the suppression of the gap by a weak magnetic field in a normal metal in the proximity of a superconductor. In the remainder of this paper we derive the results (2) and (3) assuming a specific statistical model for the hopping matrices $`t_n`$ in Eq. (1). The details of this model are not relevant, the singular behavior of the d.o.s. near zero energy being governed by the fundamental symmetries of the Hamiltonian $``$ only. As an illustration of the general validity of our result, we close with a comparison to numerical simulations.

As we are interested in the d.o.s. of the random hopping model (1) in the thermodynamic limit, the boundary conditions at the two ends of the chain are not important. For convenience we choose hard wall boundary conditions, $`t_0=t_L=0`$, $`L`$ being the length of the chain. We can solve Eq. (1) recursively in terms of a sequence of hermitian $`N\times N`$ matrices $`a_n`$,

$`a_n\psi _{n+1}`$ $`=`$ $`t_n^{}\psi _n,a_n=t_n^{}\left(\epsilon +a_{n1}\right)^1t_n^{}.`$ (4)   
The boundary condition at $`n=0`$ implies $`a_0=0`$. Evaluating the Schrödinger equation at $`n=L`$ then yields that $`\epsilon `$ is an eigenvalue of $``$ if and only if

$$\epsilon \psi _L=t_{L1}^{}\psi _{L1}=a_{L1}\psi _L,$$
(5)   
i.e. if $`a_{L1}`$ has an eigenvalue $`\epsilon `$, or alternatively, in view of Eq. (4), if $`a_L`$ has a diverging eigenvalue.

As a statistical model for the $`t_n`$ that contains all the relevant symmetries, we parameterize the $`t_n`$ in terms of the generators

$$t_n=e^{W_n},$$
(6)   
where $`W_n`$ is a real (complex) matrix for $`\beta =1`$ ($`2`$), and choose the matrices $`W_n`$ from independent Gaussian distributions with mean and variance given by

$`(W_n)_{\mu \nu }^{}(W_n)_{\rho \sigma }^{}`$ $`=`$ $`{\displaystyle \frac{1}{2}}w^2\beta [\delta _{\mu \rho }\delta _{\nu \sigma }{\displaystyle \frac{1\eta }{N}}\delta _{\mu \nu }\delta _{\rho \sigma }],`$ (7)   
$`(W_n)_{\mu \nu }`$ $`=`$ $`{\displaystyle \frac{1}{2}}(1)^n\mathrm{\Delta }\delta _{\mu \nu }.`$ (8)   
Here $`\eta `$ governs the fluctuations of $`\text{tr}W_n`$  and $`\mathrm{\Delta }`$ measures the staggering of the hopping parameter. ($`\mathrm{\Delta }`$ is the gap size that the staggering would induce in the absence of disorder.) We assume that the disorder and staggering are weak, and that the energy is small compared to the bandwidth ($`w^2,\mathrm{\Delta },\epsilon 1`$).

The matrix $`a_L`$ has eigenvalues $`\alpha _\mu `$, that can be parameterized as $`\alpha _\mu \mathrm{tan}(\varphi _\mu /2)`$ ($`\mu =1,\mathrm{},N`$). As we have discussed below Eq. (5), the energy $`\epsilon `$ is an eigenvalue of the Hamiltonian, if and only if there is an angle $`\varphi _\mu `$ with $`\varphi _\mu =\pi `$. For general $`\epsilon `$, however, none of the $`\varphi _\mu `$ will be equal to $`\pi `$. Nevertheless, we can use the angles $`\varphi _\mu `$ to compute the (disorder averaged) density of states. Hereto we first note that $`a_{L+2}=a_L`$ in the absence of disorder, staggering, and for $`\epsilon =0`$. Then, taking disorder, staggering, and a finite energy into account, and considering the length $`L`$ as a fictitious “time”, the angles $`\varphi _\mu `$ perform a Brownian motion on the unit circle, which is such that upon increasing $`L`$, they move around the circle in positive direction. The rate at which the $`\varphi _\mu `$ pass through $`\pi `$ as we increase $`L`$, i.e. their current, equals the number of states per unit length $`N(\epsilon )`$ with energy between $`0`$ and $`\epsilon `$. (This is a generalization of the node-counting theorem  used to compute the d.o.s. for $`N=1`$ .) For comparison, we remark that, in the absence of disorder, the angles $`\varphi _\mu `$ move around at a constant speed $`\epsilon `$, resulting in a constant d.o.s. With disorder, their motion acquires a random (Brownian) component, which dramatically affects their average speed, and hence the density of states, as we shall see below.

For a quantitative description a different parameterization of the eigenvalues $`\alpha _\mu `$ proves to be more convenient,

$$\alpha _\mu =\mathrm{tan}(\varphi _\mu /2)=e^{u_\mu }.$$
(9)   
The variables $`u_\mu `$ are restricted to the two branches $`\text{Im}u_\mu =0`$ and $`\text{Im}u_\mu =\pi `$ in the complex plane, see Fig. 1a. We refer to these as lower and upper branches, respectively. Noting that the $`u_\mu `$ are related to the angles $`\varphi _\mu `$ on the unit circle, we see that a (fictitious) particle with coordinate $`u_\mu `$ that vanishes on one of the branches at $`\pm \mathrm{}`$ reappears at the opposite branch, as indicated by the arrows in Fig. 1a. Upon increasing $`L`$ by two, the $`u_\mu `$ change according to $`u_\mu u_\mu +\delta u_\mu `$, where, to lowest order in $`w`$, $`\epsilon `$, and $`\mathrm{\Delta }`$, the average and variance of the increments $`\delta u_\mu `$ are

$`\delta u_\mu `$ $`=`$ $`2\epsilon \mathrm{cosh}u_\mu +2\mathrm{\Delta }+w^2\beta {\displaystyle \underset{\nu \mu }{}}\mathrm{coth}{\displaystyle \frac{u_\mu u_\nu }{2}},`$ (10)   
$`\delta u_\mu \delta u_\nu `$ $`=`$ $`4w^2\left[\delta _{\mu \nu }(1\eta )/N\right].`$ (11)   
Taking $`L`$ as a continuous variable, their distribution function $`P(u_1,\mathrm{},u_N;L)`$ obeys the Fokker-Planck equation 

$`{\displaystyle \frac{P}{L}}`$ $`=`$ $`w^2{\displaystyle \underset{\mu ,\nu }{}}{\displaystyle \frac{}{u_\mu }}\left(\delta _{\mu \nu }{\displaystyle \frac{1\eta }{N}}\right)J{\displaystyle \frac{}{u_\nu }}J^1P`$ (13)   
$`{\displaystyle \underset{\mu }{}}{\displaystyle \frac{}{u_\mu }}\left(\epsilon \mathrm{cosh}u_\mu +\mathrm{\Delta }\right)P,`$   
$`J`$ $`=`$ $`{\displaystyle \underset{\mu <\nu }{}}\mathrm{sinh}^\beta [(u_\mu u_\nu )/2].`$ (14)   
For large $`L`$ the solution of Eq. (14) acquires a steady state carrying a current $`j(\epsilon )`$ equal to the integrated d.o.s. $`N(\epsilon )`$ in the thermodynamic limit. Unfortunately, except in the case $`N=1`$, where a solution in closed form is possible , it is notoriously problematic to find the steady-state solution of a Fokker-Planck equation of the type (14), due to the lack of detailed balance . Instead, below we present a qualitative analysis of the Brownian motion process described by Eq. (14), illustrating the mechanisms that lead to the anomalous behavior of the d.o.s. near $`\epsilon =0`$. We give a detailed account of the case $`\eta =1`$, followed by a brief discussion of the general case $`\eta 1`$ and a comparison with numerical simulations.

Let us first identify the relevant parameters in the Brownian motion process (14). There are $`N`$ fictitious Brownian “particles” with coordinates $`u_\mu `$ ($`\mu =1,\mathrm{},N`$) and diffusion coefficient $`D=w^2`$. Three forces $`F_w`$, $`F_\mathrm{\Delta }`$, and $`F_\epsilon `$ act on the particles, arising from the presence of disorder, staggering, and energy, respectively. The force $`F_w=(w^2\beta /2)\mathrm{coth}[(u_\mu u_\nu )/2]`$ is a repulsive interaction, with a hard-core on the same branch and a soft core on different branches. (We adopt a convention where the particles have unit mobility.) Staggering causes a constant force field $`F_\mathrm{\Delta }=\mathrm{\Delta }`$ that favors motion to the right on both branches. Finally, $`F_\epsilon =\pm \epsilon \mathrm{cosh}u`$ pushes the particles to the left (right) on the lower (upper) branch and thus causes the nonzero steady-state current. For small energies $`\epsilon w^2`$, motion is governed by the force $`F_\epsilon `$ for large $`|u|`$, $`|\text{Re}u|u_{\mathrm{max}}`$, where

$$u_{\mathrm{max}}=\mathrm{ln}(w^2/\epsilon ),$$
(15)   
while diffusion and the forces $`F_w`$ and $`F_\mathrm{\Delta }`$ are dominant for small $`u`$, $`|\text{Re}u|u_{\mathrm{max}}`$. Following Ref. , we now approximate our model by truncating the branches at $`|\text{Re}u|=u_{\mathrm{max}}`$, see Fig. 1b, and adding a one-way move towards the upper (lower) branches at the end-points, see Fig. 1b. These simplifications don’t alter the functional dependencies of $`\rho (\epsilon )`$ for $`\epsilon 0`$, though they may affect prefactors.

For $`N=1`$, the “particle” with coordinate $`u`$ needs a time $`(2u_{\mathrm{max}})^2/2D`$ to diffuse to the end point of a branch, as one can verify from a solution of the one-dimensional diffusion equation on the line $`u_{\mathrm{max}}<u<u_{\mathrm{max}}`$, with hard wall boundary conditions at $`u=u_{\mathrm{max}}`$ and absorbing wall boundary conditions at $`u=u_{\mathrm{max}}`$ . Hence the current is $`j(\epsilon )=N(\epsilon )D/4u_{\mathrm{max}}^2`$, and after differentiation with respect to $`\epsilon `$, one finds Eq. (2) for the density of states. Diffusion effectively speeds up the particle, explaining the enhancement of the d.o.s. relative to the clean case . For $`N=2`$, the picture is completely different. As a result of their mutual repulsion, the two particles with coordinates $`u_1`$ and $`u_2`$ get trapped near the end points, say at $`u_1=u_{\mathrm{max}}`$ and $`u_2=u_{\mathrm{max}}+i\pi `$. Now the particles have to diffuse out of their traps against their repulsive interaction, until they eventually meet at $`\text{Re}u_1=\text{Re}u_2`$ and the repulsive force $`F_w`$ starts to favor travel (see Figs. 2a and b). Such a process costs a large time, which can be calculated from the diffusion equation for two particles on a line with hard wall boundary conditions at $`\text{Re}u_1=u_{\mathrm{max}}`$ and $`\text{Re}u_2=u_{\mathrm{max}}`$ and absorbing boundary conditions at $`\text{Re}u_1=\text{Re}u_2`$. We find a current $`j(\epsilon )=(4u_{\mathrm{max}}F_w^3/D^2)\mathrm{exp}(2F_wu_{\mathrm{max}}/D)`$, resulting in the d.o.s. (3). (The prefactor $`u_{\mathrm{max}}`$ arises from the degeneracy of the meeting point on the line.) We conclude that the d.o.s. for two coupled chains has only a logarithmic singularity for real hopping disorder, and a pseudogap for complex disorder. The strong $`\beta `$-dependence of the d.o.s. stems from the $`\beta `$-dependence of the interaction force $`F_w`$.

The qualitative behavior of the d.o.s. for general $`N>2`$ depends crucially on the parity of $`N`$ and closely resembles the scenarios we have outlined above for $`N=1`$ and $`N=2`$. If $`N`$ is even, all particles get “trapped” near the ends of the branches, half of them on the lower branch near $`u_{\mathrm{max}}`$, and half of them on the upper branch near $`u_{\mathrm{max}}`$, like in the case $`N=2`$, cf. Fig. 2d. The repulsive interaction force $`F_w`$ is smallest for the two particles that are closest to the origin. These two particles dominate the current, resulting in a d.o.s. of the form (3). If $`N`$ is odd, on the other hand, the picture is like that of the case $`N=1`$, see Fig. 2c. All particles get trapped at the two ends, except for one “free” particle, for which diffusion is not slowed down by the interaction forces (the repulsive forces from the other particles cancel exactly). Hence, for odd $`N`$, the d.o.s. is of the form (2). This even-odd effect is reminiscent of that found for the conductance at $`\epsilon =0`$, where the existence of a delocalized state could be attributed to the existence of a similar “free” transmission eigenvalue .

The effect of staggering $`\mathrm{\Delta }`$ is to add a constant force pointed to the right on both branches. As a function of the staggering strength, the system alternates between behavior corresponding to even and odd $`N`$. Repeating the above analysis, we find that with staggering, the d.o.s. shows the maximum (2) if $`\mathrm{\Delta }=(N+12j)w^2\beta /2`$, $`j=1,\mathrm{},N`$, whereas it shows the minimum (3) if $`\mathrm{\Delta }=(N2j)w^2\beta /2`$, $`j=1,\mathrm{},N1`$. For all other values of $`\mathrm{\Delta }`$, for $`\epsilon 0`$ we have

$$\rho (\epsilon )\epsilon ^{\beta n1},n=\underset{j=1,\mathrm{},N}{\mathrm{min}}\left|\frac{2\mathrm{\Delta }}{w^2\beta }N1+2j\right|.$$
(16)   
The case of arbitrary $`\eta 1`$ is not much different from the case $`\eta =1`$ we considered above; up to prefactors the $`\epsilon `$-dependence of the d.o.s. is not changed. The only exception is the case $`\eta =0`$, $`N=2`$, when the center of mass $`u_1+u_2`$ is pinned. As a result, the degeneracy giving rise to the logarithm in Eq. (3) is lifted, and the logarithmic prefactor vanishes.

We conclude with a comparison to numerical simulations for the d.o.s. in a quantum wire on a square lattice with a width between $`N=1`$ and $`N=4`$, and a length $`L`$ between $`10^3`$ and $`10^5`$. For $`\beta =1`$ and also for $`N=1`$ the hopping amplitudes are taken from a uniform distribution in the interval \[0.5,1.5\], while for $`\beta =2`$ the random flux model  is used, where the randomness is introduced only via the random phases of the hopping amplitudes. Results of an average over $`4\times 10^4`$$`10^6`$ disorder realizations are shown in Fig. 3. The agreement with our theoretical results, Eqs. (2) and (3), is excellent.

We would like to thank C. W. J. Beenakker, D. S. Fisher, and B. I. Halperin for discussions. PWB acknowledges support by the NSF under grant nos. DMR 94-16910, DMR 96-30064, and DMR 97-14725. CM acknowledges a fellowship from the Swiss Nationalfonds. The numerical computations were performed at the Yukawa Institute Computer Facility.