Dear Kwant community, (I don't know if here is a good place to ask this question. Let me know if it is inadequate. ) I calculate the eigen vectors of a close system containing a large s-wave superconductor and a small non-SC region. The spin is conserved, i.e. eigen states can be grouped into pure spin up and down. Can I infer the ground state wavefunction by the eigen vectors I obtained from kwant? The non-SC region is like a quantum dot. As calculated in [1], if there's a single orbital level and it couples to SC lead, the state in the quantum dot will become ground state |->=-v*|up,down>+u|0>, excited states |up> and |down>, and the highest excited state |+>=u|up,down>+v*|0>. However that is an effective model. I'm interested in how the ground state |-> looks like in tight binding model. Here is my procedure/reasoning: kwant can calculate BdG Hamiltonian on the lattice. The entire system has a ground state |g>, which is a product state of: |g> = (u1|psi1,up>|psi1,down> + v1|0>)*(u2|psi2,up>|psi2,down> + v2|0> )* ... . Here |0> is the state with no electron. The eigen vectors I obtained from kwant should be interpreted as excitations on |g>. Now I find a eigen vector localized in non-SC region. It is spin-up, with spin-up electron part [e1, e2, e3, ...] and hole part [h1, h2, h3, ...] at lattice points 1,2,3,... . These coefficients {e1,h1, e2,h2, e3, h3, ...} can be used to construct a fermionics operator gamma_up. Its Kramer's partner gamma_down would have spin-down part [e1, e2, e3, ...] and hole part [-h1, -h2, -h3, ...] at each lattice point. By applying gamma_up and gamma_down on |0> I can get one of the states in |g>. It turns out to be [sum_i {e_i h_i } + sum_i,j {e_i e_j a_i^dagger a_j^dagger}] |0>, where i,j are lattice point index, up to a normalization factor. Just as a reference, the onsite and hopping Hamiltonian are: def onsite(): (4 * p.t - p.mu + p.pot[x] ) * pauli.szs0 + p.delta[x]* pauli.sxs0 def hop(): -p.t * pauli.szs0 p.pot[x] is potential. p.delta[x] is pairing potential. pauli.szs0 is the tensor product of sigma_z (e-h basis) and identity (spin-basis). Thanks, Chien [1] "Self-consistent description of Andreev bound states in Josephson quantum dot devices", Phys. Rev. B 79, 224521 (2009)