How to verify/falsify the existence of localised edge states numerically?

I have to consider a Hamiltonian given in second quantized form in real space

$$H = \sum c_k^\dagger h_{kl} c_l \, ,$$

describing fermions on a 2d hypercubic lattice. The concrete form of the matrix h is not important, but I should mention that the matrix contains also a random part of Gaussian distributed hopping elements. As far as I can see this means that Fourier transformation is not of much use.

Now I have to find out whether there exist localized zero modes on the edges of the lattice. Ok I thought this is easy, I just have to calculate the eigenvalues of the matrix h using mathematica. Then, I have to check whether there are small eigenvalues (relative to the order of magnitude of the remaining ones). If this is the case, then I consider the eigenvectors of h corresponding to these eigenvalues and check if the entries of the vectors corresponding to the sites at an edge are much larger than in the bulk etc. it's a localized edge state.

As a check I considered the matrix for a 1d Su Schrieffer Heeger chain and expected to verify with this method that there are edge states at the ends (for suitable staggering). But the method described above didn't yield the expected result of localized states at the ends. Now I spend a lot of time on this and am very confused and am not sure anymore if what I did made any sense.

So therefore my (maybe stupid xD) question: how can I check reliably if there are localized edge states? And it's important that I can be very sure about it, so every hint how I can check the consistency additionally would help me a lot.

• Try to plot the band structures and visually check. Feb 11 at 9:50

Your method is correct. I just followed your approach and successfully observed the localized zero modes at the end of the chain. Following is my Mathematica code.

First generate a list of alternating hopping amplitudes with random fluctuation:

L = 20;
ts = 0.5 + Boole@EvenQ@Range[L-1] + RandomVariate[NormalDistribution[0, 0.1], L-1];
ListLinePlot@ts

Then construct the Hamiltonian and diagonalize:

H=Symmetrize[SparseArray[Band[{1,2}]->-ts,{L,L}],Symmetric[{1,2}]];
Eigenvalues@H
ListLinePlot[Eigenvectors[H][[-1]],PlotRange->All]

You should find that the last two eigen values are very close to zero compared to the others. They correspond to the zero modes. As we plot out the last wave vector, it should be localized at both ends of the chain (left and right ends will hybridize because the chain is of finite size).

A more systematic way to diagnose whether a state $\psi$ is localized is to look at its inverse participation ratio (IPR), which is defined as

$$\text{IPR}=\frac{\left(\sum_{i}|\psi_i|^2\right)^2}{\sum_{i}|\psi_i|^4},$$

where $\psi_i$ is the wave amplitude on site $i$ (the $i$th component of the eigen vector in question). The IPR is a simple method to quantify localization:

• If the particle is localized on only one site, then $\text{IPR} = 1$.
• In contrast, if a particle is evenly distributed over $L$ sites, then $\text{IPR} = L$.

To find the localized states, we first calculate the IPR for all eigen vectors, and pick out those with small IPR. We can actually plot each eigen vector as a point at $(\text{IPR},E)$ where $E$ is the corresponding eigen energy.

ListPlot[{1/Norm[#2^2]^2, #1} & @@@ Thread@Eigensystem@H, PlotRange -> All]

Then in the plot, it is easy to identify zero energy states (points on the horizontal axis), and the localized states (left-most points). It should be found that there are indeed two points on the horizontal axis to the left-most of the plot, which correspond to the two end modes.

Strictly speaking, a localized state is well defined only on an infinite system. Therefore, the natural idea is to change the size of the system. Localized states respond differently to the size-changing than the extended states. For example, the center-of-mass of the wave function is located at a constant distance to the edge as long as the system size is large enough, while that of an extended state will be linearly proportional to the system size. Plot the $E$.vs.$X_c$ for two different sizes, and an edge state will emerge as a fixed pair of ($E$, $X_c$). This is what we did in our paper http://journals.aps.org/pra/abstract/10.1103/PhysRevA.87.023613, to filter out the bound states out of the extended states. But i am not sure it will work for your problem.