Exact diagonalization to resolve ground state degeneracies I am studying a perturbed Toric Code model that is not analytically solvable. 
On a torus the ground state degeneracy of the unperturbed model is 4. 
Once we turn on the perturbation there is a change in the ground state degeneracy. 
I would like to detect this change in ground state degeneracy numerically using exact diagonalization techniques. 
On my computer I have stored the action of the Hamiltonian on a set of basis states. So if you give me some state $\left|\psi\right\rangle$ I can give you 
$\hat{H}\left|\psi\right\rangle$ in terms of the basis states. Now I used this information to compute the spectrum using the Lanczos algorithm and the Jacobi-Davidson algorithm. While the spectrum itself is correctly reproduced neither of these algorithms reproduces the degeneracy of the ground state correctly - not even in the unperturbed case. 
Hence the following question: What are common exact diagonalization algorithms for this type of many-body system that correctly resolve the degeneracy of the ground state? 
I am looking forward to your responses!
 A: What implementations of the two algorithms are you using? In my understanding, the Jacobi-Davidson algorithm is supposed to be good at dealing with degeneracy.
In my personal experience, using the ARPACK package is usually good enough even for highly degenerate eigenvalue problems (~10 fold or even more), as long as a reasonably large workspace is provided.
ARPACK uses what's called the Implicitly Restarted Lanczos Method. All it requires is the application of the matrix on any given vector, which is exactly what you have.
A: You must use re-orthogonalization procedure (Gram-Schmidt) at every Lanczos Step in order to correctly capture this degeneracy. I have found that degeneracy will not be captured correctly if you do not orthogonalize your Lanczos vectors. I have used Lanczos on 2D antiferromagnetic Kitaev model that has 4-fold degenerate ground state, and Lanczos correctly capture this only after the implementation of Gram-Schmidt procedure. 
While this will help to find the correct degeneracy in the eigen-values of the Hamiltonian, Lanczos still has difficulty with the actual state. Lanczos will end up returning the corresponding 4-degenerate states as some arbitrary linear combination. I do not know how to fix this issue. 
