Cosmological surveys are probing extremely large volumes at ever increasing precision. The dramatic improvement in statistical errors allows us to detect even very subtle effects in the large-scale distribution of matter, and to thoroughly test our theory of gravity on cosmological distance scales. I will present a numerical framework for cosmological N-body simulations that offers a consistent treatment of the dynamics of spacetime as described by the full set of Einstein's equations. Such a framework can be used to compute higher-order correlations in a relativistic setting, or to predict the gravitational dynamics of exotic models that do not offer a simple Newtonian description.