In this paper, we describe the first steps towards fully non-perturbative
cosmology. We explain why the conventional methods used by cosmologists based
on the ADM formulation are generally inadequate for this purpose and why it is
advantageous instead to adapt the harmonic formulation pioneered and utilized
in mathematical and numerical relativity. Here we focus on using this approach
to evaluating the linear mode stability in homogeneous and nearly homogeneous
backgrounds and devising a valid scheme and diagnostics for numerical
computation. We also briefly touch on the relevance of these methods for
extracting cosmological observables from non-perturbative simulations.Comment: 43 page
