Non-Relativistic Gravity and its Coupling to Matter

We study the non-relativistic expansion of general relativity coupled to matter. This is done by expanding the metric and matter fields analytically in powers of $1/c^2$ where $c$ is the speed of light. In order to perform this expansion it is shown to be very convenient to rewrite general relativity in terms of a timelike vielbein and a spatial metric. This expansion can be performed covariantly and off shell. We study the expansion of the Einstein--Hilbert action up to next-to-next-to-leading order. We couple this to different forms of matter: point particles, perfect fluids, scalar fields (including an off-shell derivation of the Schr\"odinger--Newton equation) and electrodynamics (both its electric and magnetic limits). We find that the role of matter is crucial in order to understand the properties of the Newton--Cartan geometry that emerges from the expansion of the metric. It turns out to be the matter that decides what type of clock form is allowed, i.e. whether we have absolute time or a global foliation of constant time hypersurfaces. We end by studying a variety of solutions of non-relativistic gravity coupled to perfect fluids. This includes the Schwarzschild geometry, the Tolman--Oppenheimer--Volkoff solution for a fluid star, the FLRW cosmological solutions and anti-de Sitter spacetimes.


Introduction
Nature is relativistic at a fundamental level but nonetheless it often effectively appears to us as non-relativistic (NR). This typically happens in many-body or condensed matter type systems but it can also be true for gravitational phenomena. General relativity (GR) can often be well-approximated by a theory of non-relativistic gravity such as the post-Newtonian (PN) approximation. Hence, uncovering the mathematical description of non-relativistic geometries, their dynamics and interaction with matter (classical and quantum) is a very relevant subject. In this paper, building on earlier work by [1,2,3,4], we try to systematically build up a geometrical language that allows us to formulate certain gravitational problems in a manifestly non-relativistic manner. This work is mainly foundational and therefore not so much concerned with applications. Nevertheless we will show that the Friedmann equations, the Tolman-Oppenheimer-Volkoff (TOV) fluid star and the usual effects due to the Schwarzschild geometry can all be captured by the theory of non-relativistic gravity described here. More generally we expect this approach to be relevant whenever the gravitational interaction can effectively be treated as instantaneous.

Background and motivation
Recent years have seen a revival in the study of non-relativistic gravity and its formulation in terms of Newton-Cartan (NC) type geometries. NC geometry was originally introduced by Cartan in 1923 [5,6] (see also e.g. [7,8]) to geometrise Newton's laws of gravitation, following the successful use of pseudo-Riemannian geometry in the formulation of Einstein's theory of general relativity. The recent developments in non-relativistic gravity have been spurred in part by modern advances leading towards a more general understanding of nonrelativistic geometry. This includes in particular the discovery of a torsionful generalisation of NC geometry, which allows for a non-closed clock form and was first observed as the boundary geometry in the context of Lifshitz holography [9,10,11]. Besides being relevant in various non-relativistic gravity theories (see e.g. [12,13,14,15,16,17,18]) this geometry, called (type I) torsional Newton-Cartan (TNC) geometry in [4], plays an important role as the background geometry to which non-relativistic field theories naturally couple [19,20,21,22]. It furthermore appears in the context of non-relativistic string theory [23,24,25,26,27] (see also references [28,29,30,31,32] for related theories involving string Newton-Cartan (SNC) geometry).
Importantly, a novel version of TNC geometry (denoted as type II) was uncovered in [4] and shown to arise directly from a careful analysis of the large speed of light expansion of pseudo-Riemannian geometry, as considered also in [3] following earlier work [1,2]. Correspondingly, it was found that type II TNC geometry is the correct framework to describe General Relativity (GR) in the non-relativistic limit. In particular, this geometry allows us to formulate a non-relativistic gravity action [4] (see also [33,34]) in any spacetime dimension 1 . This action has the property that the equations of motion (EOMs) contain the Poisson equation of Newtonian gravity, thus providing for the first time an action principle for Newton's laws of gravity. Moreover, it generalises the latter by allowing for the effects of gravitational time dilation due to strong gravitational fields. The connection with Newtonian gravity follows since type II geometry (just like type I) reduces to standard NC geometry when time is absolute. It was furthermore shown in [33] that the three classical tests of GR, namely perihelion precession, deflection of light and gravitational red-shift, are passed perfectly by this extension of Newtonian gravity since it includes gravitational time dilation effects even though retaining a non-relativistic causal structure.
There are several motivations to study non-relativistic gravity as the dynamical theory of non-relativistic geometry. First and foremost, as already mentioned above, it appears in a 1/c 2 expansion of GR and is thus a relevant limit of a celebrated and well-tested theory. Moreover, it paves the way for a covariant formulation of the post-Newtonian expansion, to any order in principle. Even more so, it generalises this expansion since to a given order in 1/c 2 the theory retains all-order effects in Newton's constant G N . Central to all this, is the appearance of a symmetry principle [4,33], which naturally arises from the 1/c 2 expansion of the Poincaré algebra. The mathematical framework is that of Lie algebra expansions 2 which precisely determines the local symmetry algebra of the "effective geometry" at any given order in 1/c 2 . In this way, type II TNC geometry arises from gauging a novel non-relativistic algebra that differs from the Bargmann algebra, while gauging the Bargmann algebra yields standard NC geometry [43]. For clarity, we emphasise here that there are other non-relativistic gravity theories than the one considered in this paper, arising from gauging type I or related avatars of NC geometry which, though not directly connected to General Relativity, are interesting in their own right from a more general perspective.
Going beyond the classical level, there are even more fundamental reasons to pursue a deeper insight into non-relativistic gravity theories, including the specific one considered in this paper, which arises from GR. One question is whether these theories have their own UV completion in terms of a non-relativistic quantum gravity theory. This in turn is interesting since the construction of such non-relativistic quantum gravity theories could provide an alternate route towards (relativistic) quantum gravity, as opposed to approaching the latter from either the classical GR or quantum field theory perspective. In fact, the question of a UV completion of non-relativistic gravity theories in terms of non-relativistic string theory has recently received a lot of attention [28,23,24,29,44,25,45,31,26,27,32,46,47] 3 (see also [50,51,46,52,53] for a relation to double field theory). Related to this, another relevant application is that classical non-relativistic gravity may have an important role to play in novel types of holographic dualities (see e.g. [16,18]). Finally, it is relevant to mention that for fixed backgrounds, non-relativistic geometry has proven to be useful for understanding aspects such as energy-momentum tensors, Ward identities, hydrodynamics and anomalies in the context of non-relativistic field theories, which are ubiquitous in condensed matter and biological systems (see e.g. [54,55,56,57,58,59]).

Outline and summary of the main results
In this paper we give a comprehensive treatment of non-relativistic gravity and its coupling to matter as it appears from the large speed of light expansion of GR. In particular, we first provide in Section 2 a detailed treatment of various aspects of the 1/c 2 expansion of GR. We will do this both from the geometric point of view, involving expanding the Lorentzian metric, as well as the algebraic perspective, which makes use of a Lie algebra expansion of the Poincaré algebra. From the geometric side, it will be convenient to first write the metric g µν in a certain "pre-non-relativistic" parameterisation as g µν = −c 2 T µ T ν + Π µν where T µ is the timelike vielbein and Π µν a spatial tensor, i.e. a symmetric tensor with signature (0, 1, . . . , 1). One finds that type II TNC geometry arises from the leading order (LO) and the next-to-leading order (NLO) fields in the 1/c 2 expansion of T µ and Π µν . In particular, the LO fields are the timelike vielbein τ µ and symmetric spatial tensor h µν with signature (0, 1, . . . , 1) familiar from standard NC geometry. These are then accompanied by two further gauge fields, m µ and Φ µν , respectively, which appear at NLO in the 1/c 2 expansion of T µ and Π µν .
This set of four spacetime tensors, together with their gauge transformation properties defines type II TNC. The gauge transformations consist of the 1/c 2 expansions of the diffeomorphisms of GR as well as the local Lorentz transformation that acts on T µ and Π µν . The latter lead to local Galilean boosts and their subleading counterparts. The former lead to diffeomorphisms plus gauge transformations (originating from the NLO terms in the 1/c 2 expansion of the diffeomorphisms of GR) acting on the NLO fields m µ and Φ µν . We note that in NC geometry the torsion is determined by the properties of τ µ : dτ = 0 corresponds to zero torsion (absolute time), τ ∧ dτ = 0 twistless torsion (i.e. twistless torsional Newton-Cartan (TTNC) geometry with a foliation in terms of equal time hypersurfaces) and no condition on τ having arbitrary torsion (full TNC geometry).
Just as in [4] we distinguish between type I and type II TNC geometry. Type I is reserved to refer to the more familiar torsional Newton-Cartan geometry that can be viewed as originating from the gauging of the Bargmann algebra [43] while type II is reserved for the version of Newton-Cartan geometry that originates from the gauging of the expansion of the Poincaré algebra as discussed in [4] as well as the present paper. The difference between the two consists of a) the gauge transformation properties of m µ and b) the fact that in type II there is an extra field, namely Φ µν which is not present in type I TNC geometry. We note that when the clock 1-form τ µ is closed the type I and type II gauge transformations of m µ agree and that furthermore the Φ µν field decouples on shell from the equations of motion of Newton-Cartan gravity. This explains why the type II structure was not manifestly present in older approaches to TNC gravity. Nevertheless it is important to understand the structure of the Bianchi identities of TNC gravity as we discuss further in appendix C.2. Type II TNC geometry arises from gauging a non-relativistic algebra that originates from a Lie algebra expansion, which can be viewed as tensoring the Poincaré algebra with the polynomial ring in σ n , with σ ≡ c −2 , truncated at NLO order (quotienting the algebra by removing all levels strictly higher than level 1). The resulting algebra has twice as many generators as the Poincaré algebra and has the Galilei algebra as a subalgebra, corresponding to the LO algebra (level zero). In this algebra the mass generator N (which is the level 1 Hamiltonian) is not central anymore, see equation (2.80) and so the Bargmann algebra is not a subalgebra.
We also display the structure that follows from expanding relativistic Lagrangians in σ, where the Lagrangians can be the Einstein-Hilbert (EH) Lagrangian itself or that of matter coupled to GR. The resulting equations of motion exhibit a cascading structure, such that at any given order the equations of motion include those of the previous order plus a set of new equations for the extra fields that appear at that given order. This is an important feature as it implies that there is a unique Lagrangian describing the geometric fields that arise at a given order in the 1/c 2 expansion and that this Lagrangian can be obtained by computing the corresponding order in the expansion of the relativistic Lagrangian one is interested in.
In order to perform the 1/c 2 expansion of the Einstein-Hilbert Lagrangian it is very useful to express it first in terms of the fields T µ and Π µν . This rewriting involves the choice of a connection that is different from the Levi-Civita connection. This new connection has the property that at leading order in the 1/c 2 expansion it provides us with a useful Newton-Cartan connection. In this way we can rewrite the EH Lagrangian in a form that makes it substantially easier to perform a large speed of light expansion. The general structure and properties at any given order in the 1/c 2 expansion is then studied in Section 3. We work out in detail the expansion of the EH Lagrangian up to next-to-next-to-leading order (NNLO). The leading order action has the property that its equations of motion restrict the clock 1-form to be hypersurface orthogonal. This means that the geometry admits a foliation in terms of equal time hypersurfaces, with Riemannian geometry on these spatial slices, a case that is known in the literature as TTNC geometry [9]. This result is interesting by itself as it shows that the dynamics restricts on shell to non-relativistic geometries that are causally well-defined (at least locally) [60]. One of the main results of the paper is the derivation of the NNLO Lagrangian, which is the Lagrangian that involves the type II TNC fields described above as well as NNLO fields. However we show that, in case we truncate the expansion after the NLO, we are allowed, without loss of generality, to impose the TTNC condition off shell via a Lagrange multiplier. We refer to the resulting Lagrangian as non-relativistic gravity (NRG). It is presented in equation (3.29). We also present the equations of motion that result from this action.
In our earlier paper a Lagrangian on these type II TNC fields was obtained via another method, which we review here including many details that were not given in [4]. This method employs the type II TNC gauge symmetries and constructs the unique two-derivative action respecting this symmetry, starting with the correct kinetic term required for Newton's law of gravitation and then completing the full action. The result is presented in equation (3.67). We also give the resulting form of the equations of motion that follow from that action. We furthermore show that, as expected, the two Lagrangians (3.29) and (3.67) are identical. The difference between the two non-relativistic Lagrangians stems from the fact that slightly different geometric variables are used. Depending on taste and type of application, one can work with either one of them.
In Section 4 we discuss the general properties of matter coupled to type II TNC geometry. We work out the 1/c 2 expansion of the energy-momentum tensor and relate this to the responses, with respect to variations of the fields of the type II TNC geometry, of the Lagrangians one obtains by expanding some matter Lagrangian order by order. We use this to derive the general matter coupled equations of motion for non-relativistic gravity. We also derive the form of the Ward identities resulting from various gauge invariances. These provide the analogue of the conservation of the relativistic energy-momentum tensor in the non-relativistic regime. We end Section 4 by specifying what happens in the simpler case when the clock 1-form τ µ is closed which leads to NC gravity. We reproduce the well-known equations of motion of Newtonian gravity in equation (4.60).
We then proceed with a detailed analysis of the non-relativistic expansion and coupling of various types of well-known relativistic matter systems in Section 5. In this section it will become clear that it is the matter sector that decides whether the geometry must have a closed clock 1-form, i.e. no torsion or whether TTNC geometry is allowed. We start with the simplest case, namely that of a point particle. Already here the expansion exhibits a rich structure, revealing that there are two distinct cases depending on whether one has absolute time or one considers the more general case of TTNC geometry. For both cases the geodesic equation is obtained. We briefly consider the case of adding an electric charge and the role of the Lorentz force. As an illustration of the difference between Lorentzian and Newton-Cartan geometries and the role of geodesics we also include a brief analysis of two-dimensional Rindler spacetime. Here, we will also see an example of the fact that because of the analytic structure of the 1/c 2 expansion, given a relativistic spacetime and two different charts that are related by a diffeomorphism that is not analytic in c, the 1/c 2 expansion of the two charts will give rise to distinct, i.e. non-gauge equivalent charts of two non-relativistic spacetimes.
After this we treat the expansion of perfect fluids and show that there are different regimes depending on how we expand the energy and pressure as a function of 1/c 2 . We then turn to various important field theory examples, namely the case of a complex and a real scalar field as well as electrodynamics. In the case of a complex scalar field we show how we can expand this in 1/c 2 such that we end up with the Lagrangian for the Schrödinger-Newton equation. This novel off shell description of the Schrödinger-Newton system includes fields whose equations of motion tell us that the clock 1-form must be closed. In the case of electrodynamics it is well-known that there are two limits, a magnetic and an electric limit depending on how we expand the gauge connection. We discuss the Lagrangian descriptions of both the magnetic and the electric expansion of Maxwell's theory.
A detailed and more extensive study of solutions is left for the future, but we conclude the paper by presenting some of the simplest solutions of the non-relativistic gravity action in Section 6. In our list of solutions we first discuss two different expansions of the Schwarzschild solution depending on whether we treat the mass parameter as constant or as being of order c 2 when we expand in 1/c 2 , following [3]. In the former case we obtain the well-known NC solution for a massive point particle while in the latter case we find a non-trivial TTNC geometry with spherical symmetry. We then proceed by studying the geodesics in this spherically symmetric TTNC background and we observe that the geodesic equations of motion for orbital motion around a centre are the same as in GR. This can be used to show that we can describe the three classical tests of GR using this non-relativistic perspective [33]. Next we obtain the non-relativistic analogue of the TOV equation. The main result here is that the resulting equations are again the same as in GR. Thus, the physical structure of fluid stars can be correctly described by non-relativistic gravity. We conclude by analysing cosmological FriedmannLematreRobertsonWalker (FLRW) metric solutions for which it is shown that, in parallel with the results above for massive objects, non-relativistic gravity yields the same Friedmann equations as one would obtain from GR. We end the section with some comments about the NR expansion of anti-de Sitter (AdS) spacetime in various coordinates.
The paper is concluded with a discussion and outlook in Section 7. Many of the more technical details are collected in various appendices. In Appendix A we discuss our notation and conventions. Appendix B provides a review of the main elements of type I TNC geometry. Newtonian gravity is presented in C. This appendix also includes a discussion about the null reduction of GR (which is invariant under type I NC gauge transformations) to contrast it with Newtonian gravity. We present an argument showing that Newtonian gravity cannot originate from a type I NC gauge invariant theory. Appendix D contains a collection of many useful identities in TNC geometry without restriction on the type of clock 1-form we use. A large number of variational identities, necessary for obtaining the equations of motion of the non-relativistic Lagrangians discussed in the main text are also collected in Appendix D. The final appendix E contains many useful results that apply to TTNC geometries.

Expansion generalities
In this first section we set up the framework for working with 1/c 2 expansions of field theories and geometry in a systematic way. A very useful so-called "pre-non-relativistic" parameterisation of Lorentzian geometry is defined. This gives a convenient starting point for studying non-relativistic (NR) expansions of the geometry in both vielbein and metric formalisms. We then continue to show that a 1/c 2 expansion of a Cartan connection taking values in the Poincar algebra reveals the underlying NR local symmetry algebra of non-relativistic gravity (NRG). With the tools needed to study 1/c 2 expansions of general Lagrangians then fully developed, we are prepared for tackling the Einstein-Hilbert (EH) Lagrangian and general relativity (GR). Finally, we also study Ward identities (WIs) and equations of motion (EOMs) in "pre-non-relativistic" parameterisation.

Non-relativistic expansions
The distinguishing feature of non-relativistic (NR) physics is that in tangent space the lightcone is flattened out completely because of the causal structure of spacetime. In Lorentzian geometry the slope of the light-cone is 1/c, with c denoting the speed of light. This means that in order to relate this to non-relativistic physics we need to perform an expansion around c = ∞. With c being dimensionful, what is meant more precisely by this statement is that we set c =ĉ/ √ σ, with σ a small dimensionless parameter which in the non-relativistic expansion is expanded around 0. For convenience we choose units in whichĉ = 1. We will consider in this paper the most conventional case of expanding in even powers of 1/c, i.e. in σ = 1/c 2 . Thus, up to an overall factor of c to some power, actions and equations of motion of relativistic gravity and matter theories are studied in a 1/c 2 expansion 4 .
We begin by discussing some general considerations related to the expansions of the fields [34]. Our starting assumption is that, up to an overall power of c which will be factored out, any field φ I (x; σ) (with I a shorthand for all spacetime and/or internal indices) is analytic in σ such that it enjoys a Taylor expansion around σ = 0, where φ I (n) (x) indicates that this is the coefficient of c −n . We will apply this expansion to both the spacetime fields of relativistic gravity as well as other types of relativistic (bosonic) fields that couple to relativistic geometry.
The first main interest of this paper is to consider the expansion of general relativity (GR) itself, and hence we first turn to applying it to the fields characterising a (d + 1)-dimensional Lorentzian manifold, which are taken here to be the relativistic vielbeine E A µ , with spacetime indices µ = 0 . . . d and tangent space indices A = 0 . . . d. Importantly, we need to explicitly choose the overall factors of c in these, such that the remaining fields have the expansion (2.1) starting at order σ 0 . The light-cone structure of the spacetime implies that the timelike vielbeine should scale different with c as compared to the spacelike ones. Thus we write the vielbeine and their inverses as where the flat metric is η AB = diag (−1, 1, . . . , 1). Hence the spatial tangent space indices a, b are raised and lowered with the Kronecker delta. We will denote this way of parameterising the vielbeine as the pre-non-relativistic parameterisation of Lorentzian geometry. The fields T µ , T µ , E a µ , E µ a are assumed to be analytic in σ = 1/c 2 and exhibit the Taylor expansion (2.1). Since the relativistic vielbeine satisfy the completeness relations E µ we have the relations The vielbeine transform under the gauge transformations of general relativity as where Ξ µ is a vector field generating the diffeomorphisms and The factors of c −1 follow from demanding that the local Lorentz transformations respect the appearance of c in (2.2) and (2.3). The inverse vielbeine transform as δE µ We thus find In terms of the geometric fields defined in (2.2), (2.3), the Lorentzian metric and its inverse take the form (2.10) We will define the spatial part of the metric and its inverse as These transform under the gauge transformations as Note that we have the relations (2.14) We will see in Section 2.6 that the pre-non-relativistic form of the metric in (2.9) enables to recast GR in a form that significantly simplifies its non-relativistic expansion.
The main goal at this stage is thus to rewrite GR in terms of fields whose 1/c expansion starts at order c 0 and whose leading order fields are unconstrained. Note that the metric (2.9) does not satisfy these criteria because g µν starts at order c 2 with a leading order term that is constrained to be a product of two 1-forms. We thus need to write the Einstein-Hilbert (EH) Lagrangian in terms of the fields T µ and Π µν . To this end it will prove convenient to work with a different connection than the usual Levi-Civita connection and consider what happens to the curvature of the spacetime with explicit factors of c appearing due to (2.9). In the following we will denote the power in 1/c with an overscript and leave the fields T µ and E a µ unexpanded (see also appendix A for notation conventions). We can write the Christoffel connection Γ ρ µν = 1 2 g ρλ (∂ µ g νλ + ∂ ν g µλ − ∂ λ g µν ) as Here we define with C ρ µν the "pre-non-relativistic" connection 17) and the remaining terms given by

19)
(2) 20) where all the Lie derivatives are with respect to T µ . We emphasise that there is still the implicit c dependence of the background fields, which will be dealt with later. It will be useful in the following to use the torsionful connection C ρ µν in (2.17) for calculating covariant derivatives. We denote this by (C) ∇ µ and the associated Riemann tensor The Ricci tensor R µν associated with the Levi-Civita connection is 1/c 2 -expanded as It is useful to note from (2.19) that S ρ µν T ν = 0 and that the C-connection satisfies One can furthermore derive that the components of the Ricci tensor obey Π µν It follows that the Ricci scalar associated with the Levi-Civita connection takes the following form when expressed in terms of pre-non-relativistic fields and the C-connection This result will be useful when we expand the Einstein-Hilbert Lagrangian in Section 2.6.

Vielbeine
By assumption the fields are taken to be analytic in σ = 1/c 2 and they thus admit a Taylor expansion. For the vielbeine this means that they can be expanded to subleading orders as For the inverse vielbeine we have where m a ≡ e ρ a m ρ as well as the 1/c 2 expansion of the completeness relations (2.4). The leading order fields in (2.36), (2.37) satisfy the completeness relations which are simply the leading order counterparts of (2.14).
Let us now turn our attention towards the gauge transformation of the various fields. The components of the infinitesimal local Lorentz transformation Λ A B and diffeomorphism generating vector Ξ µ are expanded as where Λ a 0 = c −1 Λ a . This ensures that the gauge transformations respect the 1/c 2 expansion of the fields. The interpretation is that ξ µ is a diffeomorphism generating vector field, while ζ µ (being subleading) generates gauge transformations that act on the subleading fields m µ and π a µ . They will be studied in more detail below. In particular, we will see that λ a is a local Galilean (Milne) boost and η a its subleading version, while the parameter λ a b corresponds to a local spatial rotation and σ a b its subleading version.
It is important to realise that because of the analytic structure of the 1/c 2 expansion, given a relativistic spacetime and two different charts that are related by a diffeomorphism that is not analytic in c, the 1/c 2 expansion of the two charts will give rise to distinct, i.e. non-gauge equivalent charts of two non-relativistic spacetimes. We shall give two examples of this in Sections 5.2.1 and 6.6 for the 1/c 2 expansion of 2-dimensional Rindler spacetime and of anti-de Sitter (AdS) spacetime in any dimension.
Expanding the parameters and fields in (2.5)-(2.6) leads to the following transformations deduced from collecting terms order by order in 1/c 2 : For the inverse vielbeine the leading order terms transform as (2.49) In appendix B we give a review of (torsional) Newton-Cartan (NC) geometry 5 . By contrasting this with the 1/c 2 expansion of a Lorentzian geometry we see from (2.44)-(2.47) that there are two noticeable differences. The first one is the transformation rule for m µ . This agrees with equation (B.3) only when τ µ is exact. Secondly the 1/c 2 expansion gives us a new field, namely π a µ . For these reasons we will refer to the geometry originating from the 1/c 2 expansion as type II NC geometry and the version of NC geometry reviewed in appendix B as type I NC geometry.

Metric
In Section 2.2 we found that T µ , T µ are expanded as in equations (2.36) and (2.38), respectively. The fields Π µν , Π µν , defined in (2.11), admit the following 1/c 2 expansions where we defined, in terms of the vielbeine of Section 2.2, the spacetime tensors With these conventions the expansion of the metric (2.9) can be written concisely as where we defined the tensorsh and we have τ µ τ ν Y µν = 0 . (2.61) We will not need any other components of Y µν . Notice that the contribution from the O(c −4 ) term in T µ , the next-to-next-to-leading order (NNLO) field B µ , enters inΦ µν . The leading order (LO) terms in the expansion of the metric and its inverse define the metric structures of the non-relativistic geometry. This is thus given by the clock form via τ µ τ ν and the inverse spatial metric h µν which has rank d and thus one zero eigenvalue. Its kernel is spanned by τ µ . The objects v µ and h µν are not metric tensors because they transform under the Milne boosts with parameter λ a as can be seen from equations (2.46) and (2.48). We have the Milne boost invariant relations The next-to-leading order (NLO) fields m µ , Φ µν are to be thought of as gauge fields in this context. The 1-form m µ is related to the Newtonian potential The NLO field Φ µν is less well-known but also appeared in our previous work [4]. Under a diffeomorphism generated by Ξ µ the metric transforms as δg µν = L Ξ g µν . Expanding the transformation of the metric using Ξ µ as given by (2.41) we find that the fields in the expansions (2.55)-(2.56) transform as Some of the fields appearing in the expansion of the vielbeine transform under Milne boosts with parameter λ a . The combinationsh µν ,v µ andΦ µν that appear in the metric are all Milne boost invariant. In fact they are also invariant under the spatial rotations λ a b and the subleading transformations with parameters η a and σ a b that appear in the transformation of π a µ . Depending on the setting we will either work with the fields that appear in the expansions of T µ and Π µν , i.e. τ µ , m µ , h µν and Φ µν or we will work with the fields that appear in the expansion of the metric, i.e. τ µ ,h µν andΦ µν . The transformations of the first set of fields can be readily found from the results in this and the previous section and we have δΦ µν = L ξ Φ µν + 2λ a τ (µ π a ν) + m (µ e a ν) + 2η a e a (µ τ ν) + 2ΛK µν +∇ µ ζ ν +∇ ν ζ µ , (2.72) where λ µ ≡ e a µ λ a is the Galilean boost parameter which obeys v µ λ µ = 0, η a is the parameter for subleading boosts and where we wrote the subleading diffeomorphisms ζ µ as (2.73) We furthermore defined the torsion vector 74) and the extrinsic curvature It is convenient to define a covariant derivative with respect to the torsionful connectionΓ ρ µν defined as the leading order of C ρ µν (2.17): This is the minimal collection of terms that transforms as an affine connection under diffeomorphisms, the remaining terms in the expansion of the Levi-Civita connection are tensorial. Γ ρ µν is a Newton-Cartan metric compatible connection satisfying the properties (D.2)-(D.5), but it does transform under local Galilean boosts.

Poincaré algebra
It is well-known that the properties of Lorentzian geometry can be understood by starting from a Cartan connection that takes values in the Poincaré algebra. It is therefore natural to study the 1/c 2 expansion from this algebraic point of view using the method of Lie algebra expansions. The latter has been considered e.g. in [37,38,39] and recently been applied to the 1/c 2 expansion of the Poincaré algebra in [33] and subsequently in [40,41,62].
Writing the Poincaré generators as T I = {H, P a , B a , J ab }, where H is the Hamiltonian, P a the spatial momenta, B a the Lorentz boost and J ab the spatial rotations and re-instating all factors of σ in the structure constants, the Poincar algebra becomes [J ab , J cd ] = δ ac J bd − δ bc J ad − δ ad J bc + δ bd J ac . (2.77)

The Cartan connection is
where the boost connection Ω µ a and rotation connection Ω µ ab together form the usual Lorentz connection. Let us schematically write this as A µ = T I A I µ . If we now expand the gauge A I µ then we obtain the new generators T (n) I ≡ T I ⊗ σ n , where n ≥ 0 will be referred to as the level. One then obtains an algebra by expanding in the basis of generators T (n) , with nonzero commutation relations of the form [33] ac . (2.79) We can quotient out all generators with level n > L for some L which amounts to truncating the 1/c 2 expansion. At the lowest level level L = 0 the algebra is isomorphic to the Galilean algebra when identifying H ≡ H (0) , P a ≡ P where G a is the Galilean boost generator. At the next level L = 1 we have furthermore the generators ab . Written out in detail the non-zero commutation relations of the algebra obtained by modding out all levels n > 1 are [J ab , J cd ] = δ ac J bd − δ bc J ad − δ ad J bc + δ bd J ac , where X a denotes P a , T a , G a and B a .
In particular one finds that N = H (1) is not central and the Bargmann algebra is not a subalgebra. This algebra was determined in [4] to be the relevant local symmetry algebra of type II Newton-Cartan geometry. As pointed out in [33], the fact that one does not get the Bargmann algebra, gives a group theoretical perspective on the difference between type I and type II Newton-Cartan geometry.
Let us briefly review how we can obtain type II NC geometry by gauging (2.80). This procedure has previously been studied in [43,12] for other local symmetry algebras and is a powerful way to construct relevant non-relativistic geometries 6 . Consider first a Cartan connection A µ that takes values in the algebra (2.80): A µ = Hτ µ + P a e a µ + N m µ + T a π a µ + G a ω µ a + B a Ω µ a + 1 2 J ab ω µ ab + 1 2 S ab Ω µ ab , (2.81) whose adjoint transformation is given by , Λ] = Hδ Ad τ µ + P a δ Ad e a µ + N δ Ad m µ + T a δ Ad π a µ + . . . , (2.82) where δ Ad denotes the adjoint transformations. Define a new set of transformations via In here R µν (X) denotes a curvature corresponding to generator X, defined by from which the R µν (X)'s can be determined. Without loss of generality we can take for Λ (appearing in the adjoint transformation) the following gauge transformation If we now compute δτ µ etc as defined in (2.83)-(2.86) we reproduce the transformations of τ µ , m µ , e a µ and π a µ given in equations (2.44)- (2.47). This shows that the 1/c 2 expansion of Lorentzian geometry to subleading order can be viewed as the gauging of the level 1 expansion of the Poincar algebra. More generally, this procedure can be used to any order in the 1/c 2 expansion to obtain the relevant geometric fields describing gravity to that particular order.

Lagrangians
We now present the systematics of the 1/c 2 expansion of a given theory at the Lagrangian level. Consider a Lagrangian that is a function of some field φ(x; σ) and its derivatives, i.e. L = L(σ, φ, ∂ µ φ) where we also allow for an explicit dependence on the speed of light. We now want to expand the Lagrangian and φ according to (2.1). Further below we generalise this to a Lagrangian depending on multiple fields. The explicit σ dependence can for example come from the expansion of the background metric or matter fields as well as from parameters appearing in the kinetic or potential terms. Assuming the overall power of the Lagrangian is σ −N/2 = c N , we defineL(σ) = σ N/2 L(σ, φ, ∂ µ φ) such thatL starts at order zero. Now we can Taylor expandL(σ) around σ = 0, i.e.
where the prime denotes differentiation with respect to σ. We have (2.92) Hence we can write an expansion in σ = 1/c 2 as where all the c-dependence is in the prefactors. The task is then to determine the coefficients. These can be found to be given by Hence we see that the equation of motion (EOM) of the subleading field of the subleading action is the EOM of the leading field of the leading action. A very similar calculation gives for the NNLO Lagrangian The second line forms the second variation of the LO Lagrangian and is a quadratic form involving the Hessian of the LO Lagrangian. It can be shown that Combining this with the fact that the EOM of φ (4) gives the EOM of the LO Lagrangian we see that the NNLO Lagrangian reproduces all of the EOM of the NLO Lagrangian. When there is more than one field φ then we also get mixed derivatives of the LO Lagrangian in the second line of (2.96). We can generalise the result by adding an index I to φ as follows at NNLO , (2.98) and similar at pre-leading orders. The expansion can be straightforwardly extended to include higher orders in σ. Doing so we will find relations analogous to (2.97), i.e. lower-order EOM are reproduced when going to higher orders.

Einstein-Hilbert Lagrangian
As a first step towards a 1/c 2 expansion of the Einstein-Hilbert (EH) Lagrangian, we must discuss its dimensionful normalisation. When we write g µν in terms of T µ and Π µν the line element maintains its property that its dimension is L 2 with L = length. Thus T µ dx µ has the dimension of length (if we setĉ = 1) or time (if we keepĉ) and Π µν dx µ dx ν has the dimension of L 2 . For the integration measure we find The measure Ed d xdt has dimensions T L d . Since we have √ −g = cE we take the EH action to be Using the results of Section 2.1, in particular equation (2.35), we find that the EH Lagrangian can now be written as whereL(σ, T, Π, ∂) only depends on fields analytic in σ. The prefactor of c 6 follows from the fact that the Ricci scalar is order c 2 and √ −g is order c. In hereL can be written as This is the form of the EH Lagrangian to which we can apply the results (2.93)-(2.98). At this stage, without the expansion of the fields T µ and Π µν , equation (2.102) is completely equivalent to the Einstein-Hilbert Lagrangian. This form of the Lagrangian is crucial as it establishes a starting point for a non-relativistic expansion of general relativity. In principle one could expand it to any desired order, keeping a manifest non-relativistic symmetry structure at each order. This is expected to be closely related to the usual post-Newtonian expansion of general relativity, except that we here work in a manifestly covariant framework and do not make the assumption that the fields are weak.
We define the variation of the Einstein-Hilbert action as 103) and the coupling to matter through a given matter Lagrangian L mat = L mat (σ, φ, ∂ µ φ) starting at order O(c N ) with variations defined as where we left out the variations of the matter fields. This gives the Einstein field equations of motion in the form The equations of motion satisfy two Ward identities as a consequence of diffeomorphism invariance of the action as well as invariance under local Lorentz transformations acting on T µ and Π µν as expressed in equation (2.12).

Energy-momentum conservation
Diffeomorphism invariance implies the equivalent of the divergencelessness of the usual Einstein tensor, namely Under a local Lorentz transformation the equations transform into each other as Likewise the matter Lagrangian must be invariant under diffeomorphisms and local Lorentz transformations that act simultaneously on the matter fields and on the geometric fields they couple to. On shell this leads to the relativistic conservation law for energy-momentum conservation as derived from diffeomorphism invariance of (2.104) (where we leave out the diffeomorphisms acting on the matter fields which is justified on shell as these are proportional to the matter equations of motion), The 1/c 2 expansion of these results will be studied in Section 4.2.
Assuming that the matter fields are inert under the local Lorentz transformation we find that for the matter currents we have (off shell) The currents E µ mat , E µν mat are related to the usual Hilbert energy-momentum tensor T µν via where the latter equation follows from writing Π µν in terms of the spatial vielbeine and varying those. This is clearly consistent with the Lorentz boost Ward identity (2.111) which is of course nothing other than the symmetry of the Hilbert energy-momentum tensor. The power of c in the first equality of (2.112) is fixed by demanding that the Einstein equation The conservation equation (2.108) is equivalent to that of the Hilbert energy-momentum tensor where the covariant derivative is taken with respect to the Levi-Civita connection.

The Pre-Poisson equation
An interesting combination of the equations of motion comes from a particular rescaling of the metric components: where α and β are numbers. We keep the scalar function ω arbitrary. To find the equation of motion of ω the following results are useful If we now set Expressing the left hand side using the chain rule in terms of variations of the Lagrangian with respect to the T µ and Π µν and using the Einstein equations (2.105) this this can be seen to be equivalent to It will turn out that this equation, when expanded in σ = c −2 , contains important equations like the Poisson equation and the sourcing of Newton-Cartan torsion. This will be shown in Section 3.1.3.

Non-relativistic gravity
In this section we will use the results of the previous section to obtain the action and equations of motion (EOMs) that result from the 1/c 2 expansion of general relativity (GR). In particular we will focus on the theory that governs the dynamics of the leading order (LO) and nextto-leading order (NLO) fields in the expansion of the metric. For definiteness we refer to this as "non-relativistic gravity" (NRG). As was shown already in [4] this includes Newtonian gravity but goes beyond it, as it also includes geometries with gravitational time dilaton 7 . We discuss in this section two distinct methods to obtain non-relativistic gravity. We start with the direct approach which uses the 1/c 2 expansion. Alternatively one can follow a symmetrybased route which uses gauge invariance, from which a unique two-derivative action can be obtained given a kinetic term that is required to include Newtonian gravity. Satisfyingly, we will show that the two methods lead to the same action, though in a slightly different form.

General structure
In this section we want to determine the Lagrangian that arises when we expand the fields in (2.102) using the methods of Section 2.5. This means that we will end up with a theory that is expressed in terms of the fields φ I In the next section we will rederive similar results in terms of the fields τ µ ,h µν ,Φ µν that appear in the expansion of the metric.
The 1/c 2 expansion of the Einstein-Hilbert (EH) Lagrangian will take the form We now define for n ∈ N, including zero, the equations of motion of h αβ and τ α of the N n LO Lagrangian as where the Galilean boost invariant integration measure of both type I and type II Newton-Cartan (NC) geometry is given by These equations of motion will also appear in the 1/c 2 expansion of E µ g , E µν g , defined as the response to the variations of T µ and Π µν , respectively (see (2.103)). We will give the explicit relationship between the two in Section 3.1.3.
When we go beyond leading order (for which n = 0) we encounter subleading fields in the Lagrangian. For example (−4) L NLO depends on both LO fields τ µ , h µν as well as on the NLO fields m µ and Φ µν etc. At order N n LO there are 2(n + 1) fields, each of which has its own equation of motion. However, the equations of motion that appear at order N n−1 LO are all reproduced at the nth order, see for example (2.97). The additional equations appearing at order n that are not already present at order n − 1 involve the nth order subleading fields. Figure 1: Structure of the equations of motion in the 1/c 2 expansion, of which many will enter in the Lagrangian at subleading orders. Because of the way the EH Lagrangian is expanded and the property (2.97) there will only be two new EOMs at each order to solve, the remaining ones being recursively equal to those of the previous order. Notice that when we impose TTNC off shell, all the outermost equations are zero since the LO EOMs are ∝ τ ∧ dτ as explained in Section 3.1.2.
For the NLO fields we define in analogy with (3.2) G α m = 0 because Φ µν , m µ do not appear at the leading order. In particular because of (2.97) we have The structure of the expansion of the equations of motion is summarised in Figure 1.

NNLO Lagrangian: Non-relativistic gravity
Using (2.94) we find the leading order part of the EH Lagrangian L EH given by (2.101) to be With the above conventions we can then write the variation of the leading order Lagrangian: where the leading order equations of motion are As this is a sum of squares it implies that h µν h ρσ τ µρ = 0 and thus that τ ∧ dτ = 0 on shell. This is the Frobenius integrability condition for the existenc of a foliation with normal 1-form τ = N dT where N and T are scalars. This implies that there is a foliation of the Newton-Cartan spacetime in terms of hypersurfaces of absolute simultaneity, foliated by leaves of constant T . The on shell geometry arising from the expansion is thus a twistless torsional Newton-Cartan (TTNC) geometry [9] (albeit of type II but that distinction only affects the gauge fields defined on the geometry described by the LO fields τ µ and h µν ).
This completes the discussion of the LO Lagrangian. We will continue with the NLO Lagrangian which can be obtained from (2.95) generalised to include multiple fields. From (2.95), (2.101) and (2.102) we can see that one first of all needs to compute the derivative of (2.102), which gives where we recall thatŘ µν is defined in (A.8). Using (2.95) this combines with the leading order EOMs contracted with the subleading fields so that we obtain We will refer to this Lagrangian as the Lagrangian of Galilean gravity. This theory was studied in [64] using a first-order formalism. Equation (3.16) can be related to the Lagrangian appearing in that work by a specific choice of the undetermined Lagrange multipliers that appears. That theory also has the scaling properties described in Section 2.6.2. All the leading order equations of motion are included in the NLO Lagrangian as they are obtained by varying with respect to the subleading fields. The τ µ and h µν equations of motion are where the dots denote terms that vanish on shell upon using the m µ and Φ µν equations of motion. These extra terms can easily be calculated from m µ and Φ µν variations of the NNLO Lagrangian (3.20) below using G αβ Φ . The NNLO Lagrangian is found from (2.96). For that we need the second order derivative ofL which reads Consider the general form of the Lagrangian (2.96). Adapted to the case of the EH Lagrangian we have the general result (3.20) The term in square brackets is the second variation of the LO Lagrangian and we used that L LO does not depend on derivatives of h µν . The field ψ µν is the NNLO field in the expansion of Π µν . The term in square brackets is given by which can be thought of as the field strength of m µ (for the torsional U (1) transformation (2.71)). The terms in (3.1.2) can be written as where X ρσ is some tensor. With the help of the results of appendix D.4.3 it can be shown that (for general τ µ ) where we used that (D.35) yields We note that the term in square brackets in (3.25) is symmetric in αβ due to the fact that the antisymmetric part of the Ricci tensor iš The reason we do not need to use X ρσ is as follows. The terms involving the NNLO fields B µ and ψ µν are both also of the form eh µρ h νσ τ µν X ρσ . Furthermore the EOM of the NNLO fields lead to the familiar condition h µρ h νσ τ µν = 0. This means that any variation of eh µρ h νσ τ µν X ρσ that is proportional to h µρ h νσ τ µν = 0 does not contribute on shell. It turns out that the only variation of eh µρ h νσ τ µν X ρσ that contributes on shell is the variation of τ µ except for the special case δτ µ = Ωτ µ for arbitrary Ω. Put another way, terms of the kind h µρ h νσ τ µν X ρσ can be ignored except for variations of the type h µν δτ µ . These variations give us an equation for B µ and so one arrives at the important conclusion that if we only care about EOM for the NLO fields we can ignore the term eh µρ h νσ τ µν X ρσ in the NNLO Lagrangian as well as the terms involving the NNLO fields. Effectively we can set h µρ h νσ τ µν = 0 i.e. impose the TTNC condition τ ∧ dτ = 0 off shell. If we do this we are only allowed to vary τ µ as δτ µ = Ωτ µ .
This procedure gives us what we call the non-relativistic gravity (NRG) Lagrangian: where Φ ≡ −v µ m µ is the Newtonian potential. We added a Lagrange multiplier term to enforce the TTNC condition and also used the identity In order to obtain (3.29), which is one of the central results of this paper, we worked out the variations of the first line of (3.20) using the TTNC condition. In the next subsection it will be shown that this Lagrangian is equivalent to the one given in our previous work [4] which was obtained from gauge symmetry principles. It is clear that the expansion of the Einstein-Hilbert Lagrangian in σ = 1/c 2 can be done systematically using the above framework. There is nothing, except computational complexity, that prevents an expansion to any given order in σ, yielding more and more relativistic corrections to non-relativistic gravity. We will discuss applications of this approach further in Section 7.

Equations of motion
We derive here the equations of motion based on the variations of L NRG with respect to τ µ , h µν , m µ and Φ µν where for the τ µ variation we only consider δτ µ = Ωτ µ with Ω an arbitrary function.
Because TTNC is imposed off shell a number of identities that applies to TTNC geometry has to be used in the process, which we have collected in appendix E. With these identities at hand, the variations that need to be done are straightforward albeit slightly tedious.
The variation of NRG in terms of the LO and NLO fields can be written as We denoted the responses with caligraphic symbols to distinguish them from the variations of the full NNLO Lagrangian in which the TTNC condition has not been imposed, i.e. G φ ≡ Ignoring diffeomorphisms and using TTNC, it follows from (2.64)-(2.69) that h µν and the subleading fields transform as where v µ λ µ = 0. We define the spatial projector P ν µ as The Ward identity for Galilean boost invariance with parameter λ µ is given by We can use this to simplify the process of varying h µν which is by far the most laborious variation. We can write We see that the part in front of the v µ δh µν variation is fixed by the Ward identity (3.36). We can thus without any loss of generality restrict ourselves to the P -projected variation One finds after a bit of work that the equations of motion are given by These equations are somewhat lengthy and perhaps they would be easier to handle in a firstorder formalism at the price of working with more fields. This will be studied further in upcoming work [65]. We will see in Section 3.2.2 that the equations of motion acquire a more compact form if one uses boost-invariant fields. Finally, it can be shown that when dτ = 0 it follows that Φ µν decouples and that there are drastic simplifications as we will be studied further in Section 4.4.
There are three Ward identities (WIs) that result from invariances under the LO and NLO diffeomorhisms. The latter are the gauge transformations with parameters Λ and ζ µ . For the Λ-transformation we find while for ζ µ we obtain An important role will be played by two equations that can be obtained from various contractions of the above equations of motion. The first of these is the combination We note that TTNC implies τ = N dT where N is like a non-relativistic lapse function and that this in turn implies h µν a ν = h µν N −1 ∂ ν N so that the right hand side of the above equation is the Laplacian of the non-relativistic lapse function. When we study matter couplings in the next section this will tell us something about what type of matter sources TTNC torsion. The second equation follows from the following combination Equation (3.46) contains the Laplacian of the Newtonian potential Φ. When dτ = 0 and hence a µ = 0 the field Φ µν decouples from this equation. When coupling this equation to matter we will be able to identify the sources of Newtonian gravity. Given the simplicity of (3.46) in comparison to the h µν equation of motion one might wonder if there exists a simpler way to obtain this result. The answer is affirmative and this conclusion can be obtained by going back to (2.124). To show this we need to expand the Einstein equation in powers of 1/c 2 . We can relate the 1/c 2 expansion of the Einstein equations E µν g , E µ g defined through (2.105) to the EOMs defined through (3.2). When TTNC is imposed off shell we find the non-vanishing orders to give where the Φ + 1 2 h αβ Φ αβ terms originate from expanding √ −g. In deriving this result we expanded E µ g and E µν g as (3.51) The notation follows from (2.103) which has an overall factor of c 6 and this is why we start with (−6) E µ g . Since we work here with off shell TTNC it follows that Similar statements apply to E µν g mutatis mutandis. When expanding the anisotropic scaling equation (2.124) in GR in Section 2.6 we see that at LO it reproduces the TTNC condition while at NLO it gives rise to (3.45). Finally at the NNLO it reproduces (3.46).

Lagrangian
In the previous section we derived the Lagrangian of non-relativistic gravity from the 1/c 2 expansion of the Einstein-Hilbert Lagrangian. In this section we will derive the same result using a different method. Starting with the field content that originates from the 1/c 2 expansion of the metric we will derive a (two-derivative) Lagrangian that has all the gauge invariances associated with the type II NC gauge transformations of these fields.
We will work with manifestly Galilean boost invariant quantities, i.e. τ µ ,h µν ,Φ (and their inversesv µ and h µν ) as well asΦ µν (see equations (2.57)-(2.60) for their definitions). These fields transform as in (2.64)-(2.69). In this section we will mostly be concerned with the Λ = τ µ ζ µ part of the these transformations which can be rewritten as where we defined the boost invariant torsion vector and the extrinsic curvature 8 Furthermore we will work with the Galilean boost invariant connectionΓ ρ µν defined in equation (B.9). We will refer to this as the torsional U (1) gauge transformation to contrast it with the type I Bargmann U (1) gauge transformation with parameter σ.
We will assume that the Lagrangian is at most second order in derivatives and that it has at least the kinetic termv µvνR µν . We will simply add terms until we obtain invariance under all the desired non-relativistic symmetries (2.71)-(2.72) with off shell TTNC condition τ ∧ dτ = 0 imposed. This was the original approach used to find the action presented in our previous work [4].
If we ignore the fieldΦ µν the only difference between type I and type II NC geometry is the type II transformation of m µ under the Λ-transformation which should be contrasted with the type I σ transformation given in (B.3). Note that for zero torsion the type II transformation of m µ reduces to the Bargmann U (1) gauge transformation.
Type I NC geometry can be obtained by null reduction as detailed in appendix C.1. Since from the point of view of the fields τ µ , m µ and h µν the difference between type I and type II is just one term in the transformation rule of m µ we will work with the reasonable assumption that to build a type II invariant action the non-relativistic gravity Lagrangian should contain the termĜ uu . This is the uu component of the 'null uplifted' Einstein tensorĜ M N (see appendix C.1 for more details). This term indeed contains the kinetic termv µvνR µν we would like there to be. The other terms in the Lagrangian can be found by demanding that they cancel the non-invariance ofĜ uu under the type II Λ-transformation.
To this end let us consider the transformations ofĜ uu under the variation of m µ . We start with the transformation of the connectionΓ ρ µν defined in (B.9). Using TTNC throughout we have Using this result as well as equation (C. 16) it follows that for the µν component of the a result that will be useful later. It also implies (up to a total derivative) that when taking where in the second equality we used (C.17), we find In other words since h µν δ mâν = 0 we can write this as We can straightforwardly replaceR uu in (3.62) byĜ uu because for TTNCR is independent of m µ (for TNC it would vary intoR µ u δm µ but for TTNCR µ u = 0). Finally, since for TTNC we also have that δ mĜ µν = 0, which follows essentially from (3.58), we can subtract δ ΛΦµνĜ µν from both sides of (3.62). Putting it all together we obtain the transformation rule with respect to the type II Λ-transformation, It then follows that the Lagrangian given by is invariant under the torsional U (1) transformation (3.52)-(3.53) (after partial integration) by virtue of the Bianchi identity (C.15). We can then writê whereR µν andR are given in (C. 16) and (C.19), respectively. In terms of more intrinsically defined objects this can be rewritten using (up to total derivatives) as follows from (C.20) and (C.21). The Lagrangian (3.64) (with a restored prefactor 1/16πG N ) can finally be written as what we will call the (primed) Non-Relativistic Gravity Lagrangian.
The Lagrangian (3.67) has an additional gauge symmetry which reads In order to show this we need to use two identities. The first one is the N = ν component of the Bianchi identity∇ MĜ M N = 0, which can be written as The second one is the identity The ζ transformations together with the Λ-transformation make the Lagrangian invariant under (2.64)-(2.69). It is interesting to note that, as shown above, the Lagrangian is fixed already by the Λ-invariance, with the ζ transformations appearing as an extra gauge symmetry. This completes the construction of a Lagrangian that is invariant under the type II gauge transformations.

Equations of motion
The details of the variational calculus of the action (3.67) can be found in appendix E. A number of identities that applies to TTNC geometry has to be applied in the process. If we define one finds after a bit of work that the EOMs are given by where we found it convenient to split thev µ variation in two terms h ρµ Gv µ andv µ Gv µ . We only need to consider the variation P α µ P β ν δh µν because this projection is the one where the NNLO fields decouple. This is taken care of by contracting the G h µν δh µν ∈ δL ′ NRG variation with the inverse spatial metrics G αβ h ≡ h µα h νβ G h µν . By taking the trace of G αβ Φ and usingv µ Gv µ we find the equation analogous to (3.46), 3.3 Equality of L NRG and L ′

NRG
The two actions (3.29) and (3.67) are equivalent. To see this one can express all Galilean boost invariant fields τ µ , h µν ,v µ ,h µν ,Φ,Φ µν in terms of the fields τ µ , h µν , v µ , h µν , m µ , Φ µν and B µ , although the latter will drop out when we use TTNC off shell. To make the comparison we also have to change the connection and the Ricci tensor associated to it. The difference between the two connections is (for TTNC geometries) Assuming τ ∧ dτ = 0 we can use some of the results derived in appendix E. In particular the following relations are useful A straightforward calculation then shows that the Lagrangian (3.67) is equal to (3.29). In deriving the actions we see that they arise from two different (but of course closely connected) approaches: The 1/c 2 expansion makes it obvious how the NRG theory is related to Einstein's theory of general relativity. On the other hand, from the perspective of gauging the algebra, the second approach makes it clear what the role of the local symmetry algebra is. We can also relate the equations of motion of the two Lagrangians to each other by changing the basis of the variational calculus as: When TTNC is imposed, B µ and its variation decouple in the projections we are interested in. Inserting these variations in (3.71) and equating with (3.31), we can read off the following relation between the EOMs It can be checked that these relation obey the Galilean boost Ward identity (3.36).
3.4 Comments on imposing τ ∧dτ = 0 and dτ = 0 with a Lagrange multiplier We have seen that the TTNC condition is imposed via the NNLO fields. Alternatively we can enforce this condition with a Lagrange multiplier by adding the term 86) to the Lagrangian where ζ µν = −ζ νµ . When varying the measure or τ µ (in the direction along τ µ ) in this expression we do not find any new contributions to the on shell equations of motion because of the TTNC condition enforced by ζ µν . On the other hand the variation h µν δτ µ leads to an equation of motion for ζ µν which is decoupled from all the other equations.
If we take instead with ζ µν = −ζ νµ unconstrained then the equation of motion of ζ µν enforces a NC geometry with dτ = 0. However the field ζ µν does not decouple from the equations of motion. This is because it now also appears in the equation of motion for Ω defined as δτ µ = Ωτ µ . This is what happens in the 3D CS actions for extended Bargmann algebras where ζ µν = ǫ µνρ ζ ρ in which ζ ρ is associated with the central extension of the 3D Bargmann algebra [66,15,16,18]. We conclude that setting dτ = 0 with a Lagrange multiplier does not lead to the equations of motion that are obtained from the on shell 1/c 2 expansion of GR with dτ put to zero by hand. However, it does provide us with an alternative theory obtained by adding L ′ LM to (3.64) for an unconstrained ζ µν . Since in this theory on shell dτ = 0 we can remove all terms withâ µ since they can be absorbed into the ζ µν term. Doing so leads to the following Lagrangian The interesting feature of this theory is that it is Bargmann U (1) invariant since all the non-invariance is proportional to dτ which can be compensated for by an appropriate transformation of ζ µν . TheΦ µν term can be rewritten using (C.16) andâ µ = 0 tō The term in parenthesis is the Einstein tensor of the Riemannian geometry of the constant time slices. If we are in 2+1 dimensions then this vanishes identically and by writing ζ µν = ǫ µνρ ζ ρ we recover the Chern-Simons (CS) model for the extended Bargmann algebra. Hence (3.88) can be thought of as a novel higher-dimensional generalisation of the CS model. We have thus found two different classes of theories. As it turns out only the one based on the torsional U (1) symmetry has a pure GR origin, while the Bargmann invariant case requires already in 3D to consider GR coupled to a pair of U (1) gauge fields [15]. In 4D the ζ µν can be dualised to another 2-form, B µν say, which has a 1-form gauge symmetry. In other words in 4D we can write 90) which means that there is a gauge symmetry δB µν = ∂ µ Σ ν − ∂ ν Σ µ .

Coupling to matter
In this section we discuss the coupling of matter to the non-relativistic gravity (NRG) action obtained in the previous section. We consider this by expanding a generic matter Lagrangian using the same methods as used for the Einstein-Hilbert (EH) Lagrangian. This will enable us to find the sourced equations of motion in the 1/c 2 expansion, and in particular those of NRG. We also discuss Ward identities (WIs) of the sources that follow from leading and subleading order diffeormphisms as well as the expansion of the relativistic conservation laws of the energy-momentum tensor. In order to make contact with the alternate (or primed) formulation of NRG, we will also discuss the boost-invariant currents that source the equations of motion in that case. Finally, we show how the Poission equation of Newtonian gravity can be obtained from NRG coupled to matter.

Expansion of the matter Lagrangian
T µ m = = Figure 2: Structure of the currents in the 1/c 2 expansion of the matter Lagrangian similar to Figure 1, but with τ ∧ dτ = 0 imposed off shell. For this to be consistent, the leading order currents and those related to them by variational calculus identities must be zero.
Consider any matter Lagrangian L mat = L mat (c 2 , φ, ∂ µ φ) where φ is a generic matter field with the spacetime indices suppressed. Let us suppose that the most leading term in the 1/c 2 expansion is of order c N . Since the Einstein-Hilbert Lagrangian is at most of order c 6 we assume N ≤ 6.
The expansion of the matter Lagrangian is performed using the general methods of Section 2.5 as: At each order n ∈ N including zero, we define matter currents as responses to varying the geometric fields as and similarly for the subleading fields mutatis mutandis as summarised in Figure 2. For n = 0 we have T α m = 0 since the next-to-leading order (NLO) fields do not appear at leading order (LO). Furthermore because of (2.97) we have the relations and similarly for the more subleading fields. These currents are natural to work with as they allow us to write the expansion of matter coupled general relativity to any desired order. The equations of motion (EOMs) with matter couplings at a given order 2m ≥ −6 then becomes and so on for even more subleading fields that we will not consider in this paper.
Notice that G α τ ∝ τ ∧ dτ so it is necessary to have matter with N ≤ 4 such that T α τ = 0 in order to avoid acausal non-relativistic Newton-Cartan (NC) spacetimes. We will see in the next section that N ≤ 4 in all the examples we have encountered.
Because of the above analysis we will for the rest of this section restrict to matter with N ≤ 4 and study how they can source the NRG sector of the 1/c 2 expansion of full general relativity (GR). Recall that the NRG Lagrangian is defined to be the next-to-next-to-leading order (NNLO) Lagrangian with the off shell twistless torsional Newton-Cartan (TTNC) condition τ ∧ dτ = 0. This term appears at order O(c 2 ) and is therefore sourced by (−2) L mat . We will likewise impose TTNC and define the currents as follows We again use the notation that when we are dealing with TTNC off shell, the variations of the matter Lagrangian at order O(c 2 ) are denoted by (suppressing spacetime indices) T φ | τ ∧dτ =0 = T φ for the field φ. If N = 2 the LO matter Lagrangian is of order O(c 2 ). Therefore in that case T µν Φ = T µ m = 0. When N = 4 the LO matter Lagrangian is of order O(c 4 ) and the NLO matter Lagrangian is the one that couples to NRG. In the latter case the responses T µν Φ and T µ m are generically nonzero. The equations of motion of matter coupled NRG are where the left hand sides are given by (3.39)-(3.42). We can also relate the 1/c 2 expansion of the relativistic energy-momentum tensors E µ mat and E µν mat defined in (2.104) to the expansion (4.4)-(4.5) and related equations. Analogously to Section 3.1.3 the relations are In deriving this result we expanded E µ m as and similarly for E µν mat . To explain the notation we refer to (2.104) (where we factored out c N with N = 4). Figure 3: Structure of the Ward identities (WIs): At each order there are LO WIs generated by the LO vector field ξ µ through the Lie derivative L ξ . The subleading vector field ζ µ generates a WI through L ζ , which is equivalent to the LO WI at LO. Similarly L ζ at NNLO generates a WI which is equivalent to the LO WI at NLO. This works similarly for subsubleading vector field Ξ µ (4) and is systematically extended to higher orders in the expansion. Hence, when working at a particular order, energy-momentum conservation of the previous orders in the expansion is always included at that given order.

Ward identities
The matter sector must be invariant under the gauge transformations that act simulta-neously on the geometric objects as well as on the matter fields. Since the latter variations are proportional to the matter equations of motion we can ignore these terms at the expense of being able to only derive on shell Ward identities. We can derive these on shell Ward identities from applying the transformation laws (2.64), (2.69)-(2.72) acting on the geometric fields in (4.10) and requiring invariance of the matter action up to the matter equations of motion. The structure of the tower of WIs that follows from this is summarised in Figure 3.
Diffeomorphism invariance of (4.10) implies that we have the following conservation law 9 : The v ρ and h ρσ projections are: In the first equation we used and in the second equation we used τ µ T µν Φ = 0, a property that will be derived further below. In addition the subleading diffeomorphism gives another conservation equation between the currents. This is the same as the conservation equation at NLO because of (4.4)-(4.5). Explicitly we have The two projections along v ρ and h ρσ , i.e.
give the matter counterparts of the equations (3.43) and (3.44). The v ρ projection agrees with the leading order term in the expansion of (2.109). To show this we used (4.19) and (4.15)-(4.18). Likewise the h ρσ projection agrees with the leading order term in the expansion of (2.110).

Finally the Galilean boost Ward identity is
T µ m e a µ + T µν h e a µ τ ν + T µν Φ (τ µ π a ν + m µ e a ν ) = 0 . This can be read as saying that the spatial components of T µ m are completely determined in terms of the other currents. This is the matter counterpart of (3.36) provided that T µν Φ τ µ = 0. We can show that we always have that T µν Φ τ µ e a ν = 0: This follows from as a Ward identity for the subleading Galilean boosts with parameter η a (see equation (2.72)). Later, in equation (4.53), we will see that when we assume that the order O(c 2 ) matter Lagrangian does not depend on the NNLO field B µ for the case with off shell TTNC that this in fact implies that T µν Φ τ µ = 0. This will be assumed to hold throughout.
Note that the h ρσ projection of the diffeomorphism Ward identity, equation (4.22), only contains h µρ h νσ Φ µν (to see this we need to use that T µν Φ τ µ = 0) which is the part of Φ µν that appears in the NRG Lagrangian. On the other hand the v ρ projection (4.21) contains v µ h νρ Φ µν which does not appear in the NRG Lagrangian. The terms involving v µ Φ µν will therefore drop out from the Ward identity. They must cancel against v µ Φ µν contributions to the currents. More specifically the combination T µ τ − v ρ Φ ρν T µν Φ does not depend on v µ Φ µν . This can be seen by using that the matter Lagrangian is at most of order c 4 so that at NLO, which is c 2 , the field Φ µν appears linearly in a term of the form h µρ h νσ Φ µν X ρσ where X ρσ depends on the matter fields. It is straightforward to see that then

Expansion of the Hilbert energy-momentum tensor
We collect here a few general remarks about the 1/c 2 expansion of the Hilbert energymomentum tensor T µν . From equations (2.113) and (2.114) we can see that for the case N = 4 (which is the highest value N can have for a causal spacetime, referring back to our earlier discussion in the beginning of the section) Hence we expand T µν as where τ µ T µν = 0. More specifically equations (2.113) and (2.114) tell us that T µν m ν , (4.31) T µν , (4.32) T µν m ν , (4.33) where we used (4.19) as well as (4.15)-(4.18). At the leading order the Ward identity derived from the conservation law (2.115) is given by This equation is identical to (4.26). At NLO we find Contracting this with τ ν we find which is the same as (4.25). We used the expansion of the Christoffel connection Γ ρ µν discussed in appendix E.3. To show equality with (4.25) we used equation (3.75).

Boost invariant currents
The matter currents that naturally couple to the boost invariant NRG formulation of the NNLO Lagrangian (3.67) are defined as so that the equations of motion in the presence of matter become

39)
Gv µ = 8πG N Tv µ , (4.40) Like for the geometry part, we can use the variational relations (3.78)-(3.81) to express the boost invariant currents in terms of the currents defined from varying the set of fields τ µ , m µ , h µν , Φ µν : The presence of the NNLO field B σ in (4.43) does not need to concern us as this is in agreement with the fact that only the τ µ variation in the direction of τ µ does not couple to NNLO fields (see the discussion of Section 3.1.3). The τ ρ projection of (4.43) vanishes provided τ ρ T ρσ Φ = 0. A set of currents that have all indices up and that is boost invariant is defined by using as independent variables τ µ ,h µν ,Φ µν , i.e.
Using (3.81) and δh µν = δh µν − 2τ (µ δm ν) − 2m (µ δτ ν) , (4.48) we can relate them to the other currents via T ρσ h e a ρ = T ρσ h e a ρ , (4.51) The variation with respect to the NNLO field B ν is τ µ T µν Φ . Since we will assume that there will be no B µ dependence in the matter Lagrangian with off shell TTNC it follows that One could also formulate the equations of motions by varying the gravity Lagrangian with respect to τ µ ,h µν ,Φ µν , but we shall refrain from doing that.

Newtonian gravity
In Section 3.4 it was shown that we cannot add a Lagrange multiplier to enforce dτ = 0. We will see in the next section that the situation τ ∧ dτ = 0 versus dτ = 0 is decided on by the nature of the 1/c 2 expansion of the matter Lagrangian. This happens via equations of motion imposed by the matter fields. The matter equations of motion may sometimes force dτ = 0. We will see this happening for a massive point particle, certain approximations of perfect fluids and for the Schrödinger field approximation to a massive complex scalar field. In this section we will study the necessary conditions on the matter sector in order that dτ = 0 is compatible with the equations of motion of gravity coupled to matter.
The special properties that the matter currents must satisfy can be understood completely from the scaling equation (3.45) which after using the Einstein equations becomes, This gives an equation for the non-relativistic lapse function N defined as τ = N dT since h µν a ν = h µν N −1 ∂ ν N . It then follows that if there is to be no torsion the matter must satisfy the necessary condition Existence of non-trivial solutions to e −1 ∂ µ (eh µν a ν ) = 0 depends on boundary conditions and topology of the manifold. This is crucial in establishing a twistless torsionful Schwarzschildtype vacuum solution to the EOMs in Section 6.3.
When a µ = 0 the trace equation (3.46) simplifies to By changing the connection to (B.9) all terms on the right hand side combine to form the scalarv µvνR µν . For matter coupled Newtonian gravity what the examples in the next section will show is that in those cases the most leading order in the expansion of the matter Lagrangian is of order c 2 which guarantees that the currents T µ m and T µν Φ are both zero. In that case the equations (4.56) together with the m µ and Φ µν equations of motion are all independent of Φ µν . Hence for matter coupled Newtonian gravity the m µ and Φ µν equations of motion together with (4.56) can be used to solve for the fields τ µ , h µν and m µ . When that happens Φ µν decouples from the other fields τ µ , h µν and m µ . In the next section we will see examples of this. The left hand side of (4.56) contains τ µ T µ τ which, as we will later see, is minus the mass density −ρ that enters in the geometric Poisson equation (C.29). In particular it is not a Bargmann mass as we elaborate on in appendix C.2. When dτ = 0 and the currents T µ m and T µν Φ are both zero the equations of motion of matter coupled Newtonian gravity are simply (after some rewriting) This can equivalently be written as We will see explicit realisations of this equation in the next section.

Examples of matter couplings
In this section we study some canonical examples of matter that can be put on type II Newton-Cartan (NC) geometries and their coupling to non-relativistic gravity (NRG): Point particles, fluids, scalar fields and electrodynamics. We always start with a 1/c 2 expansion of their relativistic parent and derive the equations of motion (EOMs). Their currents and Ward identities (WIs) are studied and we comment on what kind of NRG the theories can source. For the point particles we see that there is basically a branching into the usual non-relativistic (NR) particle and a novel type of particle motion that lives on a twistless torsional Newton-Cartan (TTNC) geometry and couples to torsion. To get a better conceptual understanding of the latter we study the motion on Rindler spacetime. Perfect fluids, which will play a distinguished role later on in Section 6, are then studied with and without an extra U (1) current. We then turn to both real and complex scalar fields and see, among other results, that we can derive Schrdinger-Newton theory as a special case. The 1/c 2 expansion of Maxwell electrodynamics yields novel magnetic and electric theories on (type II) TTNC geometry that we study. Finally we see how we can obtain Galilean electrodynamics (GED) on torsionless NC geometry.

Lagrangian
The proper time particle Lagrangian is In here X µ (λ) are the embedding scalars and λ is the geodesic parameter. Expanding the metric according to (2.55) we obtain the 1/c 2 expansion of the Lagrangian We still need to expand the embedding scalars according to (2.1): This is necessary for otherwise we would overconstrain the equations of motion for X µ . Note that τ µ is a function of X and so we need to expand where all functions τ µ andh µν depend on x µ (λ). The leading order (LO) Lagrangian is After a partial integration the next-to-leading order (NLO) Lagrangian becomes This is the Lagrangian of a particle on type II TNC geometry.
The equations of motion (EOMs) of the LO Lagrangian arė which is correctly reproduced by the EOMs of y µ in the subleading Lagrangian. On a TTNC background this becomesẋ µ a µ = 0 ,ẋ µ τ µ a ρ = 0 . (5.9) Since we assumed that τ µẋ µ = 0 the equations of motion force a µ = 0. On a fixed NC background the action is the same as the standard point particle action on type I TNC geometry [67,43,13]. Further notice that the LO action is order c 2 and that this couples to the next-to-nextto-leading order (NNLO) gravity action. The NLO particle action therefore only backreacts to the N 3 LO gravity action where it will source NNLO fields. This means that we can solve the geodesic equation before we backreact the solution.
If we restrict ourselves to TTNC backgrounds we can rewrite (5.7) to This observation is useful for computing the x µ equation of motion. This turns out to be identical to the type I geodesic equation on a background with dτ = 0 (which is forced upon us by the y µ equation of motion). In a gauge in which τ µẋ µ = 1 this equation is given bÿ where the connection depends on m µ (x, y) which on shell is identical to m µ since a µ = 0. In other words the y µ field has decoupled from the leading and subleading equations of motion, (5.8) and (5.12) respectively. One could thus say that y µ is a Lagrange multiplier for the condition a µ = 0.

Newtonian gravity coupled to point particles
Consider the NRG Lagrangian coupled to the LO point particle Lagrangian. They couple to each other because they both appear at order c 2 . The combined system is what we call Newton-Cartan gravity (NCG) coupled to a point particle, The equations of motion consist of (4.60) where and T µν h = 0 as is easily obtained from the variations of (5.6). In a worldline gauge for which τ µẋ µ = 1 we thus have Besides the NC equations of motion there are additional decoupled equations of motion for the field Φ µν as well as the Lagrange multiplier. The latter can be replaced by the NNLO fields by replacing the gravity part of the above Lagrangian by the NNLO Lagrangian of the expansion of the EH Lagrangian. Finally there is the x µ equation of motion which enforces dτ = 0. This thus gives a complete off shell description of NC gravity with a point particle source.

On shell expansion
Equations (5.8) and (5.12) can also be obtained from an on shell expansion by starting with the relativistic geodesic equationẌ where g µνẊ µẊ ν = −c 2 and by expanding X µ = x µ + . . .. The gauge choice g µνẊ µẊ ν = −c 2 tells us that τ µẋ µ = 1. The leading order term in the expansion of (5.16) tells us that a µ = 0.
Using this result the new leading order expansion of (5.16) gives us (5.12), so in this manner we never needed to work with y µ at the level of the EOM. The τ µ projection of (5.12) is trivially satisfied using that τ µẋ µ = 1. This suggests that we should expand (5.16) to one further subleading order and project that equation with τ µ to find something nontrivial. Indeed doing so leads to d 2 dλ 2 (τ µ y µ ) +K µνẋ µẋν + 2ẋ µ ∂ µΦ = 0 . which is the expression for energy conservation: Consider for example the geometry τ = dt, h µν dx µ dx ν = d x 2 and m = Φ(x)dt with Φ time independent. In this case (5.19) becomes d dλ 1 2˙ x ·˙ x + Φ = 0 , (5.20) which is the classical expression for energy conservation for a particle moving in a timeindependent Newtonian potential.

Fluid description
Instead of using embedding coordinates we can also say that the geodesics correspond to the integral curves of a parallel transported unit normalised vector field U µ , i.e.

Coupling to electrodynamics
We can easily generalise the action to couple the particle to the 1/c 2 expansion of a background electromagnetic potential A µ , whose dynamics we study further in Section 5.5. The expansion of the electromagnetic potential is assumed to be a 1/c 2 expansion to match the orders of the expansion of the point particle Lagrangian (5.2), so that Expanding the usual electric coupling to electrodynamics with electric charge q then yields The resulting total Lagrangian for a massive point particle with mass m and charge q is hence expanded as In order to cancel the leading order term, we see that we need to take A (−2) µ = m q τ µ . In that case also the term that gives coupling to torsion at the next-to-leading order vanishes and we get This is exactly the Lagrangian of a (type I) Newton-Cartan particle coupled to electrodynamics [43,28]. Notice that there is now no y ρ dependence enforcing torsionlessness so that these particles can propagate in torsionful geometry 10 . We return to the 1/c 2 expansion of Maxwell electrodynamics in section 5.5.

TTNC geodesics
In the previous subsection we studied the 1/c 2 expansion of the massive point particle Lagrangian and concluded that this is only consistent on a background with dτ = 0. This begs the question what about point particles moving on a torsionful NC geometry. Consider again the action This action is worldline reparameterisation invariant with respect to δX µ = ξẊ µ . The X µ equation of motion is given by If we fix the worldline reparameterisations by setting where C is any constant, then the geodesic equation becomes X µ + Γ µ νρ (X)Ẋ νẊ ρ = 0 . so only the sign in (5.33) is not automatic. Since we are dealing with a massive point particle we will take it to be timelike. The norm of the timelike tangent vector is In the previous subsection we took T µ (X)Ẋ µ = O(c 0 ) and Π µν (X)Ẋ µ = O(c 0 ). Here we will instead take the following starting point where both are expanded in a series of 1/c 2 . We can achieve this by expanding the embedding scalar as The leading order equation obtained from (5.34) and the expansion of the Christoffel symbols For NC geometry this is automatic but for TTNC geometry this gives τ νẋ ν = 0. In the latter case the only way to keep the tangent vectorẊ µ timelike is for there to be a term at order 1/c. Using This gives where τ µ and h µν are functions of x µ . We conclude that in the large c limit where we took C 2 in (5.33) to be independent of c. Using (5.34), (E.21)-(E.22), (5.38) and (5.41) we find that the leading order geodesic equation is where we used thatΓ µ νρẋ νẋρ =Γ µ νρẋ νẋρ . By contracting (5.34) with τ µ (X) and expanding up to order O(c −1 ) we furthermore obtain d dλ log (N (τ µẏ µ +ẋ µ y ν ∂ ν τ µ )) = 0 , (5.45) where we used τ = N dT . Here T is a time function. This equation exists in any coordinate system, i.e. both N and T are scalar functions of the coordinates. Equations (5.44) and (5.45) can also be obtained from an action with (5.43) appearing as a gauge fixing condition. To see this we will set τ νẋ ν = 0 off shell. The leading term in the expansion of (5.31) is given by We will define the Lagrangian as The Lagrangian is a function of x µ , y µ and their derivatives. The EOMs can be written as The y µ EOM can be written as where we used that for any TTNC geometry h µν a ν = h µν N −1 ∂ ν N . This follows from τ µ = N ∂ µ T . This agrees with equation (5.45) when F = constant which is what was assumed in deriving (5.45). Let the integration constant be a, then we find which we note is manifestly positive as required in (5.43). A useful identity is The x µ equation of motion comes out to bë We will now choose a gauge in which F is constant with F = e a = C 2 . This implies that The x µ equation of motion simplifies tö The last two equations together with determine the geodesics in a TTNC background. Note that m µ and hence the Newtonian potential does not appear. Instead we now have a force that is dictated by minus the gradient of − 1 2 N −2 which plays the role of potential energy. The fact that C 2 > 0 (for massive relativistic point particles moving below the speed of light) means that we only have bound states in this potential field. The last equation is automatically satisfied if one contracts it with τ µ or h µκẋ κ .
The fact that τ µẋ µ = 0 means that we cannot replace the λ geodesic parameter with coordinate time. This makes a particle interpretation challenging. The objects probe only the LO fields τ µ and h µν , which is dictated by local Galilean symmetries and perhaps one should think of these particles as (massless) Galilean particles. In Section 6.3.1 we shall see how the above is realised explicitly for the case of a spherical symmetric Schwarzschild-type twistless torsional Newton-Cartan (TTNC) background. Regardless of conceptual difficulties, one nevertheless obtains the same orbits as from the relativistic geodesic equation in Schwarzschild spacetime. Finally we note that the LO Lagrangian (5.46) is O(c) and so backreactions of this object would require that we include odd powers of 1/c in the metric expansion [61].

Rindler spacetime
To illustrate the difference between Lorentzian and Newton-Cartan geometries and the role of geodesics we make a slight digression and study here the simple case of 2-dimensional Rindler spacetime. In Section 6 we will consider many more examples of solutions of non-relativistic gravity.
Consider the 2D Lorentzian line element Perform the following coordinate transformation where T has dimensions of inverse velocity and R has dimensions of length. We then find This is 2D Rindler spacetime. It corresponds to the left and right wedges of a lightcone with centre at (t, x) = (0, 0). Lines of constant T are straight lines through the origin since ct x = tanh(cT ) , (5.61) and lines of constant R are hyperbolae since In the sense of type II NC geometry the metrics (5.58) and (5.60) give rise to Since the first of these has dτ = 0 and the second has τ ∧ dτ = 0 but dτ = 0 they are clearly not diffeomorphic spacetimes. We thus learn that diffeomorphic spacetimes in Lorentzian geometry can correspond to non-diffeomorphic spacetimes in NC geometry. The reason in this case is because the diffeomorphism (5.59) is not analytic in 1/c. We would like to understand the type II NC limit (5.64) of Rindler spacetime. Since the clock 1-form components vanish at R = 0 we need to check if this is in fact a coordinate singularity. To this end we will perform the same coordinate transformation (5.59) as before. We will take c =ĉ/ǫ whereĉ is numerically equal to the speed of light and ǫ is some dimensionless small quantity. The 1/c expansion then becomes an expansion around ǫ = 0. We will setĉ = 1 (so that τ has dimensions of length), and write We then have the time-like and space-like vielbeins where we write h = ee for the metric on spatial slices. The lines t = ±x correspond to T = ±∞. Since the total lapse of time (i.e. γ τ along some curve γ) is the same in all coordinate systems we can see that future and past infinity correspond to T = ±∞, i.e. t = ±x, with the exception of the origin R = 0 or what is the same (t, x) = (0, 0). The lapse of time along a curve with constant R = R 0 is R 0 T f T i dT from some initial to some final time. Clearly this goes to infinity for T i → −∞ and T f → ∞. That means that in the t, x coordinates we can consider t = ±x to represent the boundaries of spacetime except at the origin.
To understand what happens at the origin consider a straight line t = αx where |α| < 1.
It follows from (5.66) that along such a curve the components of τ become while for those of the 1-form e we get Let us consider first the time-like vielbein τ . We see that the values of the components of τ depend on the direction with which one approaches the origin and secondly when passing through the origin the sign changes. Thus we observe a discontinuous change in τ t and τ x when passing through the origin. For example if we consider a curve along the t = 0 axis then α = 0 so that we jump from τ t = 1 and τ x = 0 to τ t = −1 and τ x = 0. Since we think of τ as the normal to constant time (here T ) hypersurfaces we see that the direction of time is reversed. Going up in the right Rindler wedge is going to the future and in the left Rindler wedge going to the future means going down. This is simple to visualise by drawing a straight line T = constant through the origin. When moving it forward in time on the right means that it is going down on the left. We do not observe such a feature with the h 'metric' because it is quadratic in e and so the sign functions disappear. So as we go through the origin the metrics τ and h depend on the direction through the origin and the sign of τ is flipped. To summarise the type II TNC version of Rindler spacetime can be visualised as the two Rindler wedges of Minkowski spacetime joined at the origin and with the Milne patches of Minkowski spacetime removed entirely. We continue by considering the geodesics in this spacetime. From the expansion of the timelike geodesics we find, using equations (5.55)-(5.57), In terms of the parameter λ we have For finite values of the geodesic parameter λ we reach y T = ±∞. We can also write the solutions as We cannot view these solutions as curves that are entirely described in terms of the spacetime coordinates T and R. The reason is that we can only replace the geodesic parameter by y T but the latter is a worldline scalar y T that is the subleading embedding function for the coordinate T .
If we were to define (note that due to the appearance of y T this is not a coordinate transformation)x = R cosh y T ,t = R sinh y T , (5.72) then for y T 0 = 0 we getx = 1/C and for y T 0 = 0 we get (5.73) These are straight lines with slopes larger than +1 or less than −1, i.e. the timelike geodesics of Minkowksi spacetime. However in the sense of NC geometry that interpretation is lost. From the point of view of type II NC geometry the field y T is just a Lagrange multiplier in the action for a massless Galilean particle. Hence the only real type II geodesics are those for which T = constant. In the sense of an approximation of relativistic geodesics the field y T is of course important and simply captures the NLO effect, correctly reproducing the straight line geodesics of Minkowski spacetime.

Perfect fluids
We continue our study of the 1/c 2 expansion of matter systems by presenting the case of a perfect relativistic fluid.
We expand the normalised fluid velocity according to (2.1) as Then the normalisation condition Next we turn to the relativistic energy-momentum tensor where E and P are the relativistic internal energy and pressure. We will assume that these quantities have an expansion given by The energy-momentum tensor expands according to (4.30), which using (5.74) and the ex-pansion of the inverse metric in (2.56) gives We can use the results of Section 4.2 for the expansion of the relativistic energy-momentum conservation equation. The Ward identity (4.35) becomes while the subleading equation contracted with τ ν (see (4.37)) gives We observe that the quantity P (−4) is a 'pressure' needed to balance the force due to torsion. If P (−4) = E (−4) = P (−2) = 0 we find that we must have a µ = 0 and hence dτ = 0. In that case (4.37) becomes a mass conservation equation with E (−2) the mass density. We set E (−2) = n for mass density and P (0) = P for pressure. Then the conservation equation (4.35) becomes∇ This is the Cauchy stress-mass tensor. Equation (5.87) describes mass-momentum conservation. Contracting (5.87) with τ ν leads to the mass conservation equation∇ µ (nu µ ) = 0. The coupling to Newton-Cartan gravity (4.60) can be found by considering (4.33) and (4.34) which tells us that The subleading conservation equation (4.36) with a µ = 0 can be written as If we contract this with τ ν we obtain (5.91) The relativistic conservation equations for a perfect fluid ∇ µ T µν = 0 when projected with U ν give ∇ µ [(E + P ) U µ ] = U µ ∂ µ P , (5.92) which at subleading order leads tō Combining this with (5.91) we obtain This equation is independent of u µ (2) and for zero pressure it reduces to the fluid description of the non-relativistic particle given in (5.26). It can be shown that upon using (5.87) it is identically satisfied. In other words the energy conservation couples the LO fluid variables to the NLO fluid variable u µ (2) . In order to obtain the standard equations for a massive Galilean (Bargmann) fluid we need to include a relativistic conserved U (1) current J µ = QU µ . We expand the charge density as then for dτ = 0 the leading term in ∇ µ J µ = 0 is∇ µ (nu µ ) = 0 while the subleading term reads∇ Combining this with (5.93) we obtain If we use this to eliminate u µ ∂ µ P from (5.94) we obtain where we defined E = E (0) − Q (0) . Using (5.94) this can be rewritten as Defining the covariant energy-momentum tensor for a non-relativistic Bargmann fluid and extracting the energy current we find the energy conservation equation What we thus see is that in the case of the 1/c 2 expansion of the relativistic perfect fluid without the U (1) current we find mass and momentum conservation at leading order but the energy conservation equation couples to the subleading field u µ (2) . On the other hand, when there is a U (1) current present one can find a limit in which all the usual non-relativistic fluid equations (on a type I NC geometry), i.e. mass, momentum and energy conservation, are obtained, forming a closed set of equations. In the case of the point particle this distinction did not arise because there is no internal energy. We note that various non-relativistic fluids have been studied in the literature, see for example [68,69,54,55,56,57,58,59].

Complex scalar field
Consider the action of a complex scalar field We split the field in terms of its modulus and phase according to φ = 1 √ 2 ϕe ic 2 θ so that where we used (2.9) and (2.11). Next we expand the modulus and phase of φ according to The expansion of the Lagrangian is The leading order Lagrangian is given by so that the ϕ (0) equation of motion tells us that This condition is a sum of squares, so it implies As explained in Section 2.5, this condition will be repeated at any order of the Lagrangian through the equations of motion of the most subleading field in the expansion of φ .
With these comments, we can determine the NLO Lagrangian to be Using that h µν ∂ µ θ (0) = 0 (which now follows from the ϕ (2) EOM), the EOM of Equations (5.110) and (5.112) then imply that so that τ is exact. The LO matter Lagrangian is O(c 4 ) but it does not source gravity at that order due to (5.110). The first sourcing of gravity appears at O(c 2 ). We thus couple (3.29) to (5.111). This means that the coupling of the matter Lagrangian (5.104) to the EH Lagrangian will give rise to Newtonian gravity coupled to a scalar field in the large speed of light expansion. The NNLO Lagrangian is where Y µν is defined in (2.56) and where we added a Lagrange multiplier χ to enforce h µν ∂ µ θ (0) ∂ ν θ (0) and where we ignored all terms quadratic in h µν ∂ µ θ (0) = 0. The field χ is given by ϕ 2 (2) + 2ϕ (0) ϕ (4) . The term ϕ 2 (0) Y µν ∂ µ θ (0) ∂ ν θ (0) only contributes to the θ (0) equation of motion because Y µν is of the form h µρ Y ρ ν + h νρ Y ρ µ as it obeys τ µ τ ν Y µν = 0. Since we will not concern ourselves with NNLO fields we did not write out the Y µν term. The χ and ϕ (2) equations of motion setv µ ∂ µ θ (0) = ±m. We will take the plus sign. The θ (4) equation of motion is automatically satisfied. The θ (2) equation of motion becomes Finally the ϕ (0) equation of motion is whose equations of motion are (5.115) and (5.116). The diffeomorphisms generated by Hence if we define the wavefunction ψ = √ mϕ (0) e iθ (2) then this wavefunction satisfies the Schrödinger equation on a type I NC background with dτ = 0. In other words equations (5.115) and (5.116) can be combined into the complex Schrödinger equation For previous approaches to formulating the Schrdinger equation on Newton-Cartan spacetimes, see also [70,71,72,20,11].

Schrdinger-Newton theory
Let us next consider the coupling of the complex scalar field to gravity. The variation of (5.111) with respect to m µ , Φ µν and h µν all give zero upon using h µν ∂ µ θ (0) = 0. On the other hand the τ variation in the direction of τ gives Consider the NRG Lagrangian coupled to the NLO scalar Lagrangian. They couple to each other because they both appear at order c 2 . The combined system is If we include the next order in the expansion and restrict to EOM containing at most NLO fields we obtain the Schrödinger-Newton theory. This theory is essentially the scalar field analogue of the massive point particle theory of Section 5.1.2. The Lagrangian for the actual Schrödinger equation appears at O(c 0 ) and so the backreaction problem is considerably simplified. One first solves the equations of motion of the Newtonian gravity (4.60) with source given by (5.120), i.e.
This leads to a Newtonian potential that then appears at the next order in the Schrödinger equation giving rise to the well-known Schrödinger-Newton equation.
To see this more explicitly, choose a background that solves (5.122), i.e.
If we take d = 3 we can solve (using a Green's function) for Φ to give The Schrödinger equation (5.119) becomes the Schrödinger-Newton equation: There is an extensive literature on the subject of the Schrödinger-Newton equation [73,74,75,76,77,78]. It appears in many different experimental, numerical and theoretical studies, including applications to quantum interference [79,80] and laboratory tests of quantum gravity aspects [81,82].

Real scalar field
We briefly consider the case of a real Klein-Gordon scalar field with Lagrangian where we assume that the potential V does not depend on c explicitly. Let us expand the scalar field ϕ as where we find In a Kaluza-Klein reduction it would be more natural to give the scalar field Lagrangian the same prefactor as the EH Lagrangian. In that case we have to multiply (5.127) by c 4 16πG N . When we do this the LO scalar field Lagrangian couples to Galilean gravity, i.e. the NLO Lagrangian in the expansion of the EH Lagrangian.

Electrodynamics
Non-relativistic versions of electrodynamics have been investigated in various formulations in the literature [83,84,85,86,87,88]. In this section we will consider the 1/c 2 expansion of Maxwellian electrodynamics.
Consider the Maxwell Lagrangian   The NLO action (5.138) is similar but not the same as the covariantised version of Galilean electrodynamics (GED) presented in [89]. We here have additional couplings toΦ,Φ µν , that were not present in this previous work. The reason is that GED is naturally formulated on a (type I) Newton-Cartan background and obtained through a different non-relativistic limit: To obtain GED one takes a strict c → ∞ limit where a different scaling of the temporal and spatial components of the gauge field is allowed in addition to an extra coupling to a real scalar field. We shall see later that it is also possible to obtain GED via a 1/c 2 expansion.

Magnetic theory
It is useful to decompose the leading order Maxwell field separating the temporal and spatial components since v µǍ µ = 0 . It follows that its field strengthF ρσ is basically the magnetic field strength tensor. On the other hand the temporal componentφ is closely related to the electric potential. The transformation ofφ andǍ µ under diffeomorphisms and gauge transformations are: They both transform under the U (1) gauge transformation, but with temporal and spatial derivatives, respectively. We find that the LO and NLO Lagrangians can be written in these variables as As the O(c 4 ) term is non-zero, spacetime torsion will be sourced in the gravity EOMs (4.11)-(4.14). The LO and NLO equations of motion (5.139)-(5.140) when TTNC is imposed become Since the spatial component of the gauge field dominates the expansion we refer to this as the magnetic limit. The equation of motion forφ is given by the τ ν projection of the latter, which gives (∂ µ − a µ ) eh µρ v σF ρσ + (∂ ρ + a ρ )φ = 0 . (5.148) When τ ∧ dτ = 0 equations (5.146) and (5.148) are the same as in the magnetic limit of Maxwell's equations coupled to TTNC geometry studied in [89]. In addition we still have the equation from the spatial projection of the equation of motion (5.147), which is not present in that work. However, this is the only equation to involve the subleading gauge field A (0) µ , so we do not lose compatibility with the magnetic limit. The reason for the extra equation is that the magnetic limit studied in [89] is not the same kind of non-relativistic expansion as done here. In the previous work we scaled the temporal and spatial components of the gauge field differently and took c → ∞ as a strict limit (on shell), which projects out the extra equation that appears here.

Electric theory
Let us now see how the electric limit of Maxwell's equations studied in [89] fits into the non-relativistic expansions studied here.
First we impose off shell that The LO equation of motion is Contracting (5.154) with τ σ reproduces the LO equation of motion which is thus contained as a NLO equation of motion as expected. For dτ = 0 equation (5.154) agrees with the electric limit of Maxwell's equations coupled to NC geometry as studied in [89]. In addition we have the variation of the NLO action wrt. ϕ which gives This latter equation does not appear in the on shell strict c → ∞ limit of the electric limit of Maxwell's equations. Similar to the previous section the origin of this extra equation can be traced back to the fact that here we are performing an expansion as opposed to taking a limit as in [89].

Galilean electrodynamics
We would like to obtain Galilean electrodynamics (GED) from a 1/c 2 expansion of a relativistic theory. To do this we couple relativistic electrodynamics with gauge field A µ to a massless real free scalar Ψ and study a particular expansion. Consider then the decomposition δǍ µ = L ξǍµ + L ζǍµ − ϕe a µ e λ a L ζ τ λ + e a µ e λ a ∂ λ σ .
The sum of the relativistic Lagrangians (5.127) (with ϕ replaced by Ψ and no potential term) and (5.132) starts at order O(c 2 ) and yields: where we definedF and with TTNC off shell we have h µρ h νσF µν = h µρ h νσF µν . We see that in the leading order Lagrangian the term with two derivatives cancels out, leaving just a single spatial derivative and torsion vector terms.
When dτ = 0 we obtain GED exactly as it appears in [89] with and the subleading scalar χ and the fieldΦ αβ decouple.

Solutions of non-relativistic gravity
In this section we will consider solutions to the non-relativistic gravity (NRG) theory. Here we will explicitly demonstrate that it is much richer than just Newtonian gravity. We see that many of the canonical general relativity (GR) solutions are also exact solutions to NRG. The discussion is initiated by looking at isometries of twistless torsional Newton-Cartan (TTNC) spacetimes and how to do gauge fixing. The 1/c 2 expansion of the Schwarzschild solution is studied first and done in two different ways: a) A weak field expansion related to the PN expansion and b) a strong field expansion that will be an exact torsionful solution of NRG.
For the latter we also study the geodesics, which turn out to be the same as in GR, albeit conceptually different. Next the Tolman-Oppenheimer-Volkoff (TOV) fluid star is studied and we show that the TOV equation can be derived entirely in our NR framework. We then look at cosmological solutions and show that the FriedmannLematreRobertsonWalker (FLRW) spacetime is also an exact solution of NRG. We conclude this section by discussing inequivalent spacetimes that arise from different 1/c 2 expansions of the anti-de Sitter (AdS) spacetime.

Isometries and gauge fixing
The geometric fields τ µ , h µν , m µ and Φ ρσ transform according to studied in Section 2.3 by (2.44), (2.69), (2.70) and (2.71). An isometry is a transformation which leaves the fields unchanged, that is (with τ ∧ dτ = 0): We can also say that these are diffeomorphisms generated by K µ for which there exist λ µ , Λ and ζ µ such that L K τ µ = 0 , (6.5) By fixing diffeomorphisms and Milne boosts we can always write Demanding that δτ i = 0, δh tt = δh ti = 0 leads to the residual gauge transformations (using (2.44) and (2.70)) 14) The nonzero components transform as The residual gauge transformations act on m µ as We have defined the general notion of a Killing vector and discussed a convenient gauge for the leading order (LO) fields. In principle one could next study ansätze that preserve certain symmetries, but instead we will simply discuss a number of solutions to the equations of motion (EOMs).
For the special case of solutions of the LO theory that are also exact solutions of GR, i.e. for which m µ = 0 and Φ µν = 0 the equations of motion (4.11)-(4.14), where the left hand side is given by (3.39)-(3.42), reduce to where K µν = 1 2N ∂ t h µν . We will study solutions of these matter coupled equations in the last two subsections. We refer to [36] for more details and comments about the structure of the equations of motion in the gauge (6.9) and (6.10).

Weak gravity expansion of the Schwarzschild metric
One way to generate solutions to non-relativistic gravity is by considering the 1/c 2 expansion of GR. Therefore, let us consider the Schwarzschild metric with factors of c restored: We can perform two different physically relevant expansions depending on how we treat the mass parameter as a function of c 2 . The first option is to take m constant as we expand. In that case, by comparing to (2.56), (2.56), we can read off the fields in the expansion of the Lorentzian metric as h µν dx µ dx ν = dr 2 + r 2 dΩ S 2 , (6.26) The result is a flat torsionless Newton-Cartan spacetime with non-zero subleading fields m µ and Φ µν . One can verify that this is a vacuum solution of the EOMs (3.39)- (3.42). The solution is expressed in terms of the Newtonian potential Φ ≡ −v µ m µ = −G N m/r. In this case the 1/c 2 expansion does not terminate. The expansion of the temporal vielbein is completely captured by τ µ and m µ (with B µ and further subleading fields equal to zero), while the 1/c 2 expansion of the spatial part of the metric does not terminate. In fact the higher order spatial fields Φ where n ∈ N and Φ (2) µν = Φ µν . When all the fields in the expansion are resummed, we obtain the Lorentzian metric again. Since the torsion is zero, the expansion is really describing weak gravitational fields. From the study of geodesics in Section 5.1 we see that the geodesic equation simply becomes (5.19). In particular at this order in the expansion we do not see any terms that would give rise to the deflection of light etc., i.e. everything agrees with the prediction of Newtonian gravity.

Strong gravity expansion of the Schwarzschild metric
The situation is quite different if we perform an expansion of the Schwarzschild metric (6.23) where we take the mass to be of order c 2 so that M = m/c 2 = constant as was done in [3]. This is a strong gravity expansion of the Schwarzschild metric, i.e. one not captured by Newtonian gravity, but still described as a Newton-Cartan geometry. This provides us with a different approximation of GR as compared to the post-Newtonian expansion. In this case the expansion terminates at NLO and the geometry is described the by the LO fields This is a torsionful Newton-Cartan spacetime which is actually a solution of the equations of motion of the NLO Lagrangian (3.16) in the expansion of the EH Lagrangian (Galilean gravity) as it does not involve the subleading fields. This is a vacuum solution with torsion. In Section 6.4 we will show that this can be viewed as the exterior solution of a fluid star which can be interpreted as a source for the torsion.

Geodesics in static and spherically symmetric backgrounds
Let us consider the results of Section 5.2 and apply them to the case of geodesics in the torsionful geometry (6.29)-(6.32). This will lead to the results reported in [33].
We will start with a slightly more general case than the one in (6.29)-(6.32) and consider a geometry with spherical symmetry which can be written as h µν dx µ dx ν = γ ij dx i dx j = R 2 (r)dr 2 + r 2 dθ 2 + sin 2 θdφ 2 . (6.33) The relevant equations for geodesic motion are given in (5.55)-(5.57). The time component of (5.56) is automatically satisfied because of τ µẋ µ = 0 which impliesṫ = 0. The spatial components obeyẍ We will consider motion in the equatorial plane only, i.e.θ = 0 and θ = π/2. In this case (6.34) reduces to The latter equation can be integrated to r 2φ = l. This has been used in the first equation.
Equation (5.55) becomes which can be rewritten asṙ The λ derivative of this equation gives (6.35). Conversely, integrating (6.35) gives (6.38) with integration constant C 2 . The geodesic equations have an overall scale symmetry which involves rescaling the geodesic parameter λ and thus the angular momentum l as well as N −2 (which in τ can be compensated by rescaling t). This means that the value of C 2 is not important. The only thing that matters is whether it is zero, positive or negative. For timelike geodesics it should be positive. We now specialise to the geometry described by (6.29)-(6.32) where we will call r s = 2G N M the Schwarzschild radius (treated as independent of c). Let us restrict attention to geodesics for whichφ = 0 then after rearranging we find dr dφ This is a well-known equation describing planetary motion (for C 2 > 0) in the Schwarzschild geometry including the effects of the perihelion precession. It also captures the phenomenon of light deflection (for C 2 = 0). For more discussion we refer to [33].

Tolman-Oppenheimer-Volkoff equation
In this section we will show that the Tolman-Oppenheimer-Volkoff (TOV) equation for the hydrostatic equilibrium of a spherically symmetric isotropic body (fluid star) can be derived entirely within the non-relativistic gravity framework. The solution we are after is known to be static, and hence we need K µν = 0. From the equations of motion (6.20)-(6.23) we infer that the sources must obey From the boost Ward identity (4.27) we also learn T µν h P α µ τ ν = 0. In order not to source any subleading orders we can fulfil these conditions if we take a perfect fluid as defined in section 5.3 with only E (−4) and P (−4) nonzero. This can be seen to follow from equations (4.32)-(4.34) We furthermore take for the fluid velocity so that he fluid energy-momentum tensor reads The equations of motion (6.20)-(6.23) then reduce to It can be shown that the 1/c 2 expansion found in [3] agrees with equations (6.43) and (6.44).
The fluid equations of motion are (5.85) and (5.86).
Let us now turn to the most general d = 3 static spherically symmetric ansatz for the spacetime geometry that follows from using the results of Section 6.1 and requiring the relevant isometries, summarized by h µν = diag 0, e −2β(r) , 1/r 2 , 1/(r 2 sin 2 θ) , (6.47) h µν = diag 0, e +2β(r) , r 2 , r 2 sin 2 θ . (6.48) where α(r) and β(r) are arbitrary functions. The same ansatz can be obtained from the 1/c 2 expansion of the corresponding analysis for a Lorentzian metric g µν using Birkhoff's theorem.
Inserting the static spherically symmetric ansatz into (6.43) and (6.44) we find that the equations of motion take the form From these we can solve for P (−4) , E (−4) to find This conservation equation is exactly the same as the one that appears in the relativistic case.
With this one finds that the remaining equations can be rewritten as after reinstating factors of c 2 . This is exactly the relativistic TOV equation for a stellar body of mass of the order of c 2 , i.e. for which M (r) = m(r) c 2 is order c 0 , with m the dimensions of mass. This is the same point of view as taken in [3] and in section 6.3 of this paper.
We thus conclude that the physical structure of stellar bodies can be described completely by non-relativistic (strong) gravity. Its description does not require the principle of relativity.

Cosmological solutions and Friedmann equations
We will next show that FriedmannLematreRobertsonWalker (FLRW) solution solves the LO equations of motion (6.20)-(6.23). Using again the form (6.9), (6.10) for the LO fields, we can show that in this case we must have where σ ij is constant in time and describes a maximally symmetric space in d dimensions.
We assume that the scale factor a is independent of c. It follows that the acceleration and extrinsic curvature satisfy where the dot denotes differentiation with respect to time. The equations of motion (6.20)-(6.23) then reduce to Just like in the previous subsection we will translate this set of equations into conditions on the 1/c 2 expansion of a perfect fluid. Using equations (4.32)-(4.34) it follows that we need to take a perfect fluid as defined in section 5.3 with only E (−4) , P (−4) , E (−2) and P (−2) nonzero. For the fluid velocity we will take again This means that the fluid energy-momentum tensor takes the form The above equations of motion then simplify further to the set of equations where k = −1, 0, 1 depending whether the spatial metric σ ij in γ ij = a 2 (t)σ ij is a maximally symmetric space of constant negative, zero or positive curvature. The third equation follows from the absence of a source for torsion. We see that the sources for the spatial derivatives are P (−4) and E (−4) while the sources for the time derivatives are P (−2) and E (−2) . After resumming the energy E and pressure P according to (5.80) and (5.81), one obtains the Friedmann equations for a d + 1-dimensional cosmological spacetime. It is straightforward to see that for d = 3 one can put these equations in the conventional form The cosmology one obtains from non-relativistic gravity thus agrees with the (relativistic) Friedmann equations obtained from GR. If the spatial curvature k vanishes, then the spacetime can be sourced by a perfect fluid with E (−4) = P (−4) = 0. These equations were derived from an Einstein equation written as G µν = 8πG N c 4 T µν . In the presence of a cosmological constant Λ this is also written as G µν + Λg µν = 8πG N c 4 T ′ µν . If we define T ′ µν = E ′ +P ′ c 2 U µ U ν + P ′ g µν then we have the relations P ′ = P + c 4 8πG N Λ and E ′ = E − c 4 8πG N Λ. In the case of de Sitter spacetime we have a = exp (Ht) where H is the constant Hubble parameter and k = 0. This leads to E (−4) = P (−4) = 0 and The de Sitter radius is c/H. This implies that E ′ = 0 = P ′ as it should.
The leading order fluid conservation equations are given by (5.85) and (5.86), reflecting that the quantities E (−4) , P (−4) (which are nonzero for k = 0) are homogeneous and that energy is conserved. The subleading conservation equations similarly become: and when we rewrite these in terms of E and P , they are equivalent to the conservation equations appearing in GR. For different and more canonical approaches to Newtonian cosmology, see references [90,91,92,93,94,95].
6.6 1/c 2 expansion of AdS spacetimes define l = c H with H independent of c and we treat both signs of l 2 at the same time we find where the upper sign is for AdS and the lower sign is for dS spaces. Expanding this to NLO the resulting NC geometry can be read off as where we left out Φ µν and where we transformed to Cartesian coordinates. Such a NC geometry is known as the Newton-Hooke spacetime. In [96] we showed that such a spacetime can be written in the form of a non-relativistic FLRW geometry with flat spatial slices by a sequence of NC gauge transformations. For the AdS case, i.e. the upper sign in (6.84), one can show using the techniques of [96] that this can be written as where we transformed h µν and m µ using a Galilean boost and abelian gauge transformation and where we furthermore transformed the coordinates. For the dS case we similarly find These should however not be confused with the FLRW spacetimes discussed above as the latter result from a different 1/c 2 expansion. We conclude that, starting with the Lorentzian AdS spacetime one encounters a situation similar to the two different 1/c 2 expansions of the Schwarzschild geometry discussed at the beginning of this section. In that case, the difference depends on how the mass as a function of c 2 is treated. In analogy, we see here that the expansion depends on how we treat the cosmological constant as a function of c 2 .

Discussion and outlook
The main purpose of this paper has been the development of non-relativistic gravity (NRG) as it appears from a large speed of light expansion of general relativity (GR). We have given a detailed introduction to the underlying geometry, which we dubbed type II Newton-Cartan (NC) geometry. We have presented the gauge transformations of the fields and how they can be thought of as arising from the gauging of an algebra that in turn can be obtained from an algebra expansion of the Poincaré algebra. We defined the Lagrangian of NRG to be given by the next-to-next-to-leading order (NNLO) Lagrangian in the 1/c 2 expansion of the Einstein-Hilbert Lagrangian in which we impose the twistless torsional Newton-Cartan (TTNC) condition for a global foliation in terms of constant time slices with the help of a Lagrange multiplier. We derived this Lagrangian using two different methods: by direct 1/c 2 expansion and by using gauge invariance under type II gauge transformations. We have subsequently discussed the resulting equations of motion and the coupling to matter. We have furthermore described some of the main examples of matter couplings, i.e. point particles, perfect fluids, real and complex scalar fields and electrodynamics. Finally, we have presented some of the simplest solutions of non-relativistic gravity (coupled to a perfect fluid) and commented on their physical relevance.

Open problems and future directions
As a first avenue of further analysis, understanding NRG from a Hamiltonian perspective would tell us more about the number of degrees of freedom. This can be achieved by the usual counting of the phase space dimension and constraints per spacetime point. The Hamiltonian perspective would furthermore provide us with natural candidates for the definition of asymptotic charges such as mass, energy, momentum, angular momentum etc. In this light it would for example be interesting to see what would happen with asymptotic symmetry groups in the non-relativistic regime. This might help us in understanding if NRG has the potential to admit a holographic dual description. For example in the case of the AdS/CFT correspondence one could wonder about what happens with the Brown-Henneaux analysis in 3 dimensions [97,98] when we expand in 1/c 2 or what happens to the fluid/gravity correspondence [99,100] in the non-relativistic (NR) expansion. One could also examine how to implement the 1/c 2 expansion in the bulk at the level of the boundary theory of known dualities. More generally, it is interesting to speculate whether there is a relation with the entropic and emergent gravity ideas of [101,102] which are also connected to Newtonian gravity and modifications thereof.
Another perspective on the theory would be provided by performing a detailed analysis of the linearised spectrum (for example around flat NC space). We do not expect that the theory has propagating degrees of freedom, and hence we expect that the gravitational interactions are instantaneous as in Newtonian gravity. Nevertheless, it would be interesting to understand the structure of the propagators and how the theory would behave from a perturbative QFT point of view. Obviously, it would be important to study further 1/c 2 expansions of relativistic solutions in detail. This will teach us more about the conceptual nature of non-relativistic gravity. In particular it would be interesting to see how the 1/c 2 expansion of the Kerr geometry fits into this framework. This is a sufficiently general spacetime to study in order to understand if there is a notion of a non-relativistic black hole. We can then also hopefully shine some more light on the correct interpretation of the geodesics studied in Section 5.2.
It is clear from the analysis presented here that if one were to continue the expansion of the Einstein-Hilbert (EH) Lagrangian beyond NNLO it would quickly become very challenging. We expect that performing the same analysis in first order formalism should be more suited to a higher order expansion. As we have stressed in this paper, we know the underlying symmetry principle at any given order along with the systematics of the expansion of the EH action, but one needs to develop an efficient way to extract results. To this end we plan to pursue the analysis of the 1/c expansion in first order formalism in [65]. In this connection see also the references [35,40,41]. A related point that needs to be addressed is the question about the status of the odd powers of 1/c. These have been discarded in this work as a simplifying assumption, but ultimately we need to understand their physical significance. In this light we refer to [61].
One of the possible applications of a higher order analysis would be to make contact with the post-Newtonian (PN) approximation. In particular it would be of interest to construct a map relating non-relativistic gravity in our formalism to more conventional PN parameterisations since it goes beyond the approximation where the torsion is zero and encodes strong field effects more naturally. Importantly, there may be relevant domains of validity in physical processes, such as the early phase of inspirals of compact objects, where non-relativistic gravity can either give new results or alternative methods to check known results. Examining the two-body problem in non-relativistic gravity will thus be important as well. In another direction, it would be worthwhile to obtain the action and equations of motion at higher order in the 1/c 2 expansion.
It would of course also be important to examine how the non-relativistic action is related to string theory. The current state-of-the-art includes non-relativistic strings that probe type I NC geometry, as well as the closely related string Newton-Cartan (SNC) geometry. It would thus be very interesting to uncover how strings couple to the type II NC geometry and in particular whether the beta-functions of this putative theory reproduce our NRG theory. More generally, it would be interesting to see how branes couple to type II TNC geometry differently from type I TNC geometry as studied in [103,104].
Finally there are of course various other open issues one could consider, for example the coupling of a non-relativistic spinning particle to type II TNC, fermionic matter actions and adding spacetime supersymmetry 11 . With the richness of non-relativistic physics demonstrated so far there are certainly still numerous other interesting studies to be done. X indicates the order of some coefficient of a Laurent/Taylor expansion in 1/c for some object X (c). There is one exception to this rule. When expanding a field φ whose 1/c 2 expansion starts at order c 0 we write instead (A.1)

A.2 Curvature
For any torsionful connection Γ ρ µν with covariant derivative ∇ µ the Riemann tensor R µνσ ρ and torsion tensor T ρ µν are universally defined through The Bianchi identities are The Ricci tensor is also universally defined as We will always work with a connection such that where M is the measure. This implies that R µνρ ρ = 0 , (A. 10) and hence that the antisymmetric part of the Ricci tensor is In this paper we use three different choices of affine connections. The formulae of this appendix apply to all of these three choices.  Table 1: Comparison of notation used in three different papers including the present one. A '-' denotes that the corresponding object has not been defined in the corresponding paper.

A.3 Comparison of notations
The notation in this paper have been streamlined and differs from some of the choices in previous works. To make comparison easier we present in this appendix Table 1 with notations used in two other papers [4,3].

B Review of Torsional Newton-Cartan geometry
Torsional Newton-Cartan (TNC) geometry has been reviewed extensively in the literature. We repeat here the most fundamental aspects, see also references [72,43,9,10,60,107,12]. TNC geometry is characterised by three tensors τ µ , h µν , m µ with h µν symmetric and of signature (0, 1, . . . , 1), subject to the following gauge redundancies where λ µ obeys v µ λ µ = 0 with v µ defined as follows. The inverse of −τ µ τ ν + h µν is given by −v µ v ν + h µν with v µ τ µ = −1, v µ h µν = 0, τ µ h µν = 0 and h µρ h ρν = δ ν µ + v ν τ µ . The inverse objects transform as The parameter λ µ corresponds to local Galilean (or Milne) boosts and the parameter σ to Abelian gauge transformations associated with particle number conservation. TNC geometries only admit degenerate metric structuresgiven by τ µ τ ν and h µν respectively. For example, lower indices can no longer raise be raised at will because contravariant and covariant tensors of the same rank can not be mapped to each other in a one-to-one way. The non-uniqueness in v µ , h µν can be interpreted as the ambiguity due to frames related by local Galliean boost transformations (also sometimes called Milne boosts in the literature).
In addition to τ µ and h µν , one can define the following Galilean boost-invariant spacetime tensorŝ These form a convenient set of variables to describe TNC geometry and they satisfy the completeness relationsv It should be noted thatv µ ,h µν andΦ are not invariant under the Abelian gauge transformation with parameter σ. One can also define an affine connection Γ λ µν so that we may take covariant derivatives. It is natural to require the TNC equivalent of metric compatibility ∇ ρ τ µ = 0, ∇ ρ h µν = 0. The first of these conditions implies that any metric compatible connection must have the same temporal projection of the torsion tensor 2τ ρ Γ ρ [µν] = 2∂ [µ τ ν] . Thus constraints on torsion imply restrictions on the geometric data in contradistinction the case of Riemannian geometry.
The full TNC case is acausal as has been argued in [60], but is still interesting in applications to field theory and holography because the energy current is the response to varying an unconstrained τ µ . On the other hand in the torsionless Newton-Cartan geometry there is a notion of absolute time as τ = 0 implies that all observers agree on the time interval between events. For most purposes we will restrict ourselves to the twistless torsional Newton-Cartan geometry which defines a spacetime foliation whose normal 1-form is τ .
In [108,89] it was shown that any boost invariant TNC connection may be written as where we define a distinguished TNC connection as and C λ µν is a spacetime tensor; a TNC analogue of the "contortion" tensor. ForΓ λ µν the torsion tensor is given byT λ µν = −2v λ ∂ [µ τ ν] . The connectionΓ ρ µν is manifestly boost invariant while, unless dτ = 0, it is in general not invariant under the Abelian gauge transformation with parameter σ.
We will also make use of the non-boost invariant connectioň which has the nice property that it does not contain m µ and is therefore invariant under the σ gauge transformation. It has non-zero torsion given byŤ λ C Null reduction of Einstein gravity

C.1 General properties of null reductions
It is well-known that what we refer in this paper to as type I NC geometry can be obtained from null reduction of a Lorentzian metric with a null isometry, see for example [109,110]. We will denote the null Killing vector by ∂ ∂u . If the d + 1 dimensional NC spacetime has coordinates x µ then the null uplifted Lorentzian geometry has coordinates x M = (u, x µ ).
Any Lorentzian metric with a null isometry and its inverse can be written as where g uu = 0 due to the existence of the null Killing vector ∂ u . The null reduction of the components of the (d + 2)-dimensional Levi-Civita connectionΓ L M N (without any constraints on τ µ ) isΓ We denote the higher-dimensional Levi-Civita connection and associated curvatures with a hat. Note that the null-reduced Levi-Civita connection is equal to the symmetric part of the boost invariant NC connection (B.9). The definitions of extrinsic curvaturesK µν ,â µ can be found in Table 2.
A very useful object is the null reduction of the Ricci tensor, which has the following componentsR whereR µν is the Ricci tensor corresponding to the connection (B.9). The extremely useful property of these expressions is that they transform nicely under the Bargmann U (1) gauge transformation of type I NC geometry. These transformations are easy to derive using the fact that the U (1) corresponds to the following higher dimensional diffeomorphism Infinitesimally this reads δu = −ξ u where ξ u is the u-component of ξ M , the generator of (d + 2)-dimensional diffeomorphisms. Using the tensorial transformation rule ofR M N under a diffeomorphism generated by ξ M = −σδ M u one shows that These transformation rules are fully general and thus true for any TNC geometry. From the higher-dimensional Bianchi identities for the Einstein tensor, i.e.∇ MĜ M N = 0, we can derive two very important results e −1 ∂ µ eĜ µ u = 0 , (C.14) They are true for any TNC geometry. In fact, since these are identities, they are true regardless of which U (1) transformation we assign the fields to have! To derive these results one needs to use the null reduction formulae for the higher dimensional Christoffel connection (C.3)-(C.7). The first of these Bianchi identities, (C.14), is the geometrical counterpart of Bargmann mass conservation. This follows by using the (d + 2)-dimensional Einstein equation and recognisingT µ u as the mass current of the lower-dimensional theory, see e.g. appendix A of [56]. The second identity (C.15) is the geometrical counterpart of energy conservation. The difference between the two is just a raising or lowering of the u index. Furthermore, using an argument similar to the one leading up to (C.13), it can be shown thatĜ µ u is Bargmann U (1) invariant for any TNC geometry whereasĜ µu is not, not even for TTNC geometry.
If we specialise to TTNC geometries the null reduction of the Ricci tensor simplifies tô The higher-dimensional Ricci scalarR for TTNC geometries is given bŷ Furthermore, for TTNC geometries we have thatĜ µ u = 0 andĜ uu is given bŷ We will thus considerĜ uu as a nicely transforming completion of the kinetic term.

C.2 Newtonian gravity vs. null reduced general relativity
In this section we review how to obtain standard Newtonian gravity in the torsionless Newton-Cartan framework. We will also compare this to the null reduction of general relativity (GR). When dτ = 0 we can write the connection (B.9) as Γ ρ µν =Γ ρ µν + K ρ µν , (C. 22) whereΓ ρ µν = −v ρ ∂ µ τ ν + 1 2 h ρσ (∂ µ h νσ + ∂ ν h µσ − ∂ σ h µν ) , (C.23) where F µν = ∂ µ m ν −∂ ν m µ is the Bargmann U (1) curvature. The connectionΓ ρ µν only depends on the U (1) invariant fields τ µ and h µν . This connection is not Galilean boost invariant. Note that the dependence ofΓ ρ µν on m µ is via the addition of a tensor K ρ µν . The sourceless NC equations of motionR µν = 0 can then be written as In NC gravity this should be supplemented with the condition dτ = 0. The left hand side is pure geometric data and the right hand side depends entirely on the "electric" and "magnetic" field strength components of F µν . The divergence of the electric field strength in the third equation, i.e. e −1 ∂ ρ (ev ν h ρσ F νσ ), is what gives rise to Newton's law of gravity when appropriately sourced by a mass density. The null reduction of the Einstein equations of motion, in the absence of sources, also leads toR µν = 0. It does not lead to dτ = 0 on the nose. However one can make the argument that the sourceless null reduced Einstein equations of motion force dτ = 0. This happens in two steps. First of all the equation of motion (C.10) leads to the TTNC condition and secondly the EOMv µR µu = 0 implies that ∂ µ (eh µν a ν ) = 0 . (C.28) Since furthermore (E.2) holds, which means that h µν a ν = N −1 h µν ∂ ν N for some function N defined via τ = N dT , this equation is a Laplacian acting on N . This follows from the fact that e = N √ γ where γ is the metric on the T = constant slices. These spatial slices are described by d-dimensional Riemannian geometry. In particular they are Ricci flat as follows from the equationR µν = 0. In the absence of sources the clock 1-form τ , which defines the foliation of the spacetime, must be everywhere regular. That means that N is a bounded harmonic function on a Ricci flat d-dimensional geometry. Such functions must be constant and so a µ = 0 implying that dτ = 0. Hence for dτ = 0 we observe a Bargmann invariance of the sourceless equations of NC geometry. It is therefore tempting to suggest that the sourceful generalisation should also obey Bargmann invariance. The sourceful NC equations that correspond to Newtonian gravity are given by 12R We will now argue that this equation is not compatible with a Bargmann invariant coupling of NC geometry to matter. The mass current J µ in a Bargmann invariant theory is U (1) invariant and conserved. The only candidate geometrical quantity that obeys the same properties is theĜ µ u component of the Einstein tensor on a background with a null isometry (see (C.14)). Hence the coupling must beĜ µ u ∝ J µ . (C.30) From a null reduction point of view we have of course that J µ =T µ u , whereT M N is the null uplifted energy-momentum tensor. This equation implies that upon contraction with τ µ we obtainR uu ∝T uu = ρ . (C.31) From the form ofR uu in (C.10) we thus see that mass sources τ ∧ dτ = 0. This is in direct conflict with Newtonian gravity described in the previous section because in that theory the notion of mass is compatible with dτ = 0. Hence ρ in (C.29) is not a Bargmann mass density. We thus conclude that Newtonian gravity cannot originate from a Bargmann invariant theory 13 . We can make this a bit more explicit by performing a null reduction of the Einstein-Hilbert Lagrangian, which up to total derivatives yields L = e 16πG N h µνR µν + 1 2 h µνâ µâν −

2Φ
h µρ h νσ τ µν τ ρσ . (C.32) This is not a consistent reduction but the inconsistency is extremely mild in that all of its equations of motion agree with the null reduction of Einstein's equation. It only fails to reproduce theĜ uu ∝T uu equation of motion. The reason for this is simply that the null reduction setsĝ uu = 0 off shell and so we cannot vary this component. Furthermore, thê G uu ∝T uu equation of motion does not impose any constraints on any of the other equations of motion that follow from the null reduced action. The reason is that for any givenΦ, the equations of motionĜ µν ∝T µν together withĜ µu ∝T µu form a closed set. The remaining equationĜ uu ∝T uu , instead of imposing a constraint merely completes the other equations of motion by suppling an equation of motion from whichΦ can be determined. Finally we remark that theT uu component of the null reduced energy-momentum is not a new independent source but a composite object that is formed from various other sources. All of the above implies that we can use the null reduced Einstein-Hilbert (EH) action to study the coupling between geometry and matter for Bargmann invariant theories. On shell one then simply complements the equations of motion withĜ uu ∝T uu to provide an equation for Φ. The point is the when adding matter in a Bargmann invariant manner the fieldΦ couples to the mass density. It then follows that theΦ equation of motion leads to the same conclusion as before, namely that Bargmann mass density sources TNC torsion for which τ ∧ dτ = 0. We conclude that in order to couple matter to NC gravity we cannot use type I NC geometry. In particular, in Section 4 we show that type II NC geometry will lead to a consistent coupling between gravity and matter.

D Torsional Newton-Cartan identities
In this appendix we collect useful identities that hold for type I and type II NC geometries. Throughout this appendix we will not impose any restrictions on the clock 1-form τ .

D.1 Summary of definitions
We present in in Table 2 a summary of the definitions of fields used throughout the paper.

D.2 Identities for covariant derivatives, Riemann and Ricci tensors
The torsion tensor ofΓ λ µν is given by

D.3 Identities for extrinsic curvatures
We have the following useful contractionŝ

E Twistless Torsional NewtonCartan identities
This appendix is very similar to the previous one except that we now assume that τ obeys the hypersurface orthogonality (or TTNC) condition τ ∧ dτ = 0. All of the identities of appendix D of course all apply here as well but there are many simplifications when the TTNC condition is imposed.

(E.4)
If instead of contracting the second Bianchi identity (for κ = λ) with v ν h µσ we contract it with two inverse spatial metrics we obtain (for TTNC),

E.2.1 Basic identities
For TTNC geometry we only need to calculate some projective variations. The τ µ variation along τ µ is best computed by setting δ Ω τ µ ≡ Ωτ µ . (E. 6) for arbitrary Ω. In order to compute the h variation it is sufficient to consider δ P h µν ≡ P ρ µ P σ ν δh ρσ (E.7) where P µ ν ≡ h µλ h λν . The P variations of the other geometric objects are given by: We have some relevant variations that are needed in Section 3.1.3.
(E. 16) With this result it follows that the variation ofv µvνR µν with respect toΦ results in a total derivative.