From this:
Normal mode oscillations. If the Hamiltonian turns out to be a quadratic function of coordinates and momenta for a system of $N$ objects, e.g. $$H=\sum_{ij} M_{ij} q_i q_j + \sum_{ij} M_{ij} p_i p_j$$ then you can simply do a canonical transformation along the eigenvectors of $M_{ij}$ to diagonalize $M_{ij}$, and your system separates into independent harmonic oscillators.
Can someone please elaborate a little bit on this? I don't know how to "simply do a canonical transformation along the eigenvectors of $M_{ij}$ to diagonalize $M_{ij}$".