  We evaluate the covariance matrix of the matter power spectrum using
perturbation theory up to dominant terms at 1-loop order and compare it to
numerical simulations. We decompose the covariance matrix into the disconnected
(Gaussian) part, trispectrum from the modes outside the survey (beat coupling
or super-sample variance), and trispectrum from the modes inside the survey,
and show how the different components contribute to the overall covariance
matrix. We find the agreement with the simulations is at a 10\% level up to $k
\sim 1 h {\rm Mpc^{-1}}$. We show that all the connected components are
dominated by the large-scale modes ($k<0.1 h {\rm Mpc^{-1}}$), regardless of
the value of the wavevectors $k,\, k'$ of the covariance matrix, suggesting
that one must be careful in applying the jackknife or bootstrap methods to the
covariance matrix. We perform an eigenmode decomposition of the connected part
of the covariance matrix, showing that at higher $k$ it is dominated by a
single eigenmode. The full covariance matrix can be approximated as the
disconnected part only, with the connected part being treated as an external
nuisance parameter with a known scale dependence, and a known prior on its
variance for a given survey volume. Finally, we provide a prescription for how
to evaluate the covariance matrix from small box simulations without the need
to simulate large volumes.
