  Making cosmological inferences from the observed galaxy clustering requires
accurate predictions for the mean clustering statistics and their covariances.
Those are affected by cosmic variance -- the statistical noise due to the
finite number of harmonics. The cosmic variance can be suppressed by fixing the
amplitudes of the harmonics instead of drawing them from a Gaussian
distribution predicted by the inflation models. Initial realizations also can
be generated in pairs with 180 degrees flipped phases to further reduce the
variance. Here, we compare the consequences of using paired-and-fixed vs
Gaussian initial conditions on the average dark matter clustering and
covariance matrices predicted from N-body simulations. As in previous studies,
we find no measurable differences between paired-and-fixed and Gaussian
simulations for the average density distribution function, power spectrum and
bispectrum. Yet, the covariances from paired-and-fixed simulations are
suppressed in a complicated scale- and redshift-dependent way. The situation is
particularly problematic on the scales of Baryon Acoustic Oscillations where
the covariance matrix of the power spectrum is lower by only 20% compared to
the Gaussian realizations, implying that there is not much of a reduction of
the cosmic variance. The non-trivial suppression, combined with the fact that
paired-and-fixed covariances are noisier than from Gaussian simulations,
suggests that there is no path towards obtaining accurate covariance matrices
from paired-and-fixed simulations. Because the covariances are crucial for the
observational estimates of galaxy clustering statistics and cosmological
parameters, paired-and-fixed simulations, though useful for some applications,
cannot be used for the production of mock galaxy catalogs.
