  We explore a signature of phase correlations in Fourier modes of dark matter
density fields induced by nonlinear gravitational clustering. We compute the
distribution function of the phase sum of the Fourier modes,
theta_k1+theta_k2+theta_k3, for triangle wavevectors satisfying k1+k2+k3=0, and
compare with the analytic prediction in perturbation theory recently derived by
one of us. Using a series of cosmological N-body simulations, we extensively
examine the time evolution and the dependence on the configuration of triangles
and the sampling volume. Overall we find that the numerical results are
remarkably consistent with the analytic formula from the perturbation theory.
Interestingly the validity of the perturbation theory at a scale k is
determined by P(k)/V_samp, the ratio of the power spectrum P(k) and the
sampling volume V_samp, not by k^3P(k) as in the case of the conventional
cosmological perturbation theory. Consequently this statistics of phase
correlations is sensitive to the size of the sampling volume itself. This
feature does not show up in more conventional cosmological statistics including
the one-point density distribution function and the two-point correlation
functions except as a sample-to-sample variation. Similarly if the sampling
volume size V_samp is fixed, the stronger phase correlation emerges first at
the wavevector where P(k) becomes largest, i.e., in linear regimes according to
the standard cosmological perturbation theory, while the distribution of the
phase sum stays fairly uniform in nonlinear regimes. The above feature can be
naturally understood from the corresponding density structures in real space as
we discuss in detail.
