Simulating quantum network of harmonic oscillators Let's say that I have a system of $n$ particles $p_1,\ldots,p_n\in\mathbb{R}^3$ (where $n$ here is on the order of 10,000).  Furthermore, suppose we have a graph $G=(V,E)$ describing some network, where the set of vertices $V$ is the set of particles and the set of edges $E$ satisfies $|E|\ll n^2$.  Each edge in $G$ will represent a spring between the corresponding pair of particles holding them together.  I will also give each particle $p_i$ its own scalar potential $V_i$.
The wavefunction of this system will be of very large dimension:  $\Psi(x_1,\ldots,x_n;t)$.  We can write the Hamiltonian as $$H=\sum_{i=1}^n\left(-\frac{h^2}{2m} \nabla_i^2+V_i(x_i)\right)+\sum_{(p_i,p_j)\in E} c\|x_i-x_j\|^2$$
Obviously I cannot step the Schrodinger equation for this huge system forward in time using any standard discrete simulation method, nor can I even write down a discretized version of $\Psi$ due to its large dimensionality.
Are there techniques for computing approximate low-energy states of this system?  What approximations are reasonable here?  Can one find likely positions of each particle $p_i$?
 A: The way you wrote it is probably not what you intended to write down--- you probably want $||x_i - x_j - C_{ij}||^2$ instead of $||x_i-x_j||^2$, so that the springs have a nonzero stationary length, so that the masses have a shape in the absence of external forces. If you don't do this, the classical least-energy solution is all the particles are exactly on top of each other.
But I will solve it as you wrote it (it's not much harder the other way, but it's more annoying to figure out the stationary configuration). In this case, you write the spring term as
$$ \sum A_{ij} x_i x_j $$
and you diagonalize A by a rotation (this is possible since A is symmetric).
$$ y_i = R_{ij} x_j $$
Where the sum on j is implicit, and R is a rotation matrix. You can compute the R using any diagonalization algorithm, they are all easy in the symmetric case. Then in terms of y, you have decoupled oscillators 
$$ \sum \lambda_i y_i^2 $$
Where the $\lambda_i$ are the eigenvalues of A, and each one is a normal harmonic oscillator that you solve independently. To find the expected values of the x's, you write the x's in terms of the y's.
When you have a symmetry, like an abstract lattice of x's where there are only finitely many classes of points not related by lattice symmetry to each other, you can write the explicit solution for even an infinite lattice, using Fourier theory. This is the starting point for quantum field theory.
