I know of several ways to calculated the eigen states of 1D potentials (i.e. DVR, Crank–Nicolson, etc). However I wonder what is the most efficient way to do the same for a N-Dimensional potential? How would this change if the potential has some kind of symmetry ?
This is really too big a question to ask in a single post, as there are at least two steps that need to be addressed: discretization and eigensolves.
The first step is that you need to discretize the partial differential equation that underlies your problem. There are of course many ways of doing this, the finite element and finite differences method being two popular ones.
The second one being to find the eigenvalues and eigenvectors of this discrete problem. This is much harder than in 1d because the matrices are typically much larger, and methods are consequently more expensive. There are, again, any number of methods available for this (and they are available in high-quality libraries).
If you want to see some of this in action, here is a tutorial program of the deal.II library that shows you some of this: http://dealii.org/developer/doxygen/deal.II/step_36.html (Disclaimer: I'm one of the authors of this library.)
The main difference is in the dimensionality, that is reflected in the matrices. In 1D, e.g, using (FEM) Finite Element Method or (FDM) Finite Difference Methods the Stiffness (Hamiltonian) matrix is tridiagonal. In 2D, the Stiffness matrix is pentadiagonal for the FDM and octadiagonal for the FEM with linear or bilinear elements, using structured meshes. For Mass (overlap) matrices, we have identity matrix for FDM and octadiagonal for FEM. Depending on the solver the storage of the matrix is also a factor to take into account. Besides these changes in dimensions, the procedure is similar: you look the first $m$ eigenvalues/eigenvectors.
Normally, potentials are defined over non-bounded domains and you will need to do one of two things:
- "truncate" it, make it's value infinite outside some region. Or, equivalently make the function zero (artificial Dirichlet BCs); or
- Use some type of BCs that approximate the unboundedness of the domain, like absorbing boundary conditions or infinite elements.
Regarding the symmetries, you can have smaller meshes with these symmetries using the right BCs. For example, a Quantum Harmonic Oscillator have symmetric and anti-symmetric modes. For the symmetric ones we have $f'(x=0)=0$, for the anti- symmetric cases $f(x=0)=0$.
TL;DR, I would suggest FDM if you are interested in simple cases, and FEM, since it can handle boundary conditions more naturally.