Implications of Symmetry and Pressure in Friedmann Cosmology. I. Formalism

We show that derivation of Friedmann's equations from the Einstein-Hilbert action, paying attention to the requirements of isotropy and homogeneity during the variation, leads to a different interpretation of pressure than what is typically adopted. Our derivation follows if we assume that the unapproximated metric and Einstein tensor have convergent perturbation series representations on a sufficiently large Robertson-Walker coordinate patch. We find the source necessarily averages all pressures, everywhere, including the interiors of compact objects. We demonstrate that our considerations apply (on appropriately restricted spacetime domains) to the Kerr solution, the Schwarzschild constant-density sphere, and the static de-Sitter sphere. From conservation of stress-energy, it follows that material contributing to the averaged pressure must shift locally in energy. We show that these cosmological energy shifts are entirely negligible for non-relativistic material. In relativistic material, however, the effect can be significant. We comment on the implications of this study for the dark energy problem.

This metric encodes the Copernican Principle: on large scales, the universe has no preferred locations or directions. In other words, the zero-order universe is homogeneous and isotropic.
This equation is inconsistent if T µν encodes O(1) positiondependent sources. Gravitational radiation signatures of ultra-compact object mergers, as first reported by Abbott et al. (2016), provide direct evidence that such sources exist. To address this issue, one can instead proceed with the following definition: Isotropy and homogeneity then constrain the O(1) contribution into the form of a perfect fluid T µν (0) = ρ(η) + P(η) u µ u ν + P(η)a −2 η µν .
The functions P(η) and ρ(η) can then be determined from data, and compared to theoretical expectations. To predict P(η) and ρ(η) in the redshift z 40 universe, overwhelming evidence of a hot big bang justifies summing the second moments of a Boltzmann distribution f i for each species i at some temperature T i T µν (0) (η) ∝ i p µ p ν p 0 f i |p|, T i (η) d 3 p.
In the z 40 universe, however, structures may have begun to form and appropriate distribution functions are unknown.
In the spirit of Kaiser (2017, §3), we can try to estimate P(η) and ρ(η) at late times with linear perturbation theory. We begin from a single star: 1. Hypothesis (vacuum): exterior to the star is effectively vacuum.
2. Hypothesis (spherical symmetry): the star is spherical to a very good approximation.
3. ∴ Birkhoff's theorem =⇒ the star is perceived as a point mass to good approximation.
At sufficient distance from the star, the Newtonian monopole limit is recovered. One then invokes superposition and averages over an ensemble of such point masses to produce ρ(η). Performing these averages over distinct regions in space reveals a lower-bound volume V, beyond which the ρ(η) computed at different points agree. The observational success of linear perturbation theory supports the expectation that the rest-mass densities of such ensembles should dominate any kinetic contributions by ∼ 10 5 . Alternatively, following Weinberg (1972, p. 474) or Peebles (1993, p. 297), one applies Birkhoff to a round ball 1 cut out from some presumed RW universe. The interior region is now vacuum and the conclusion becomes that whatever was inside the region can be averaged and reintroduced as a point mass. In this case, the interior motions contribute zero pressure to the exterior universe. Either way, the predicted contribution to T µν (0) (η) is a P(η) ≃ 0 perfect fluid. Observations, however, are consistent with a late-time universe dominated instead by an apparently fixed-density P(η) = −ρ(η) perfect fluid.
At this point, we have built a cosmological model by assuming an RW metric ansatz and then constructing an appropriate stress tensor. One could instead imagine starting with Einstein's equations for the actual matter distribution and ask what must be done to arrive at an RW metric. This alternative approach proceeds as follows, where we summarize aspects of Wetterich (2003, §1) and Ellis (2011, §2). Observe that the matter in the visible universe is distributed uniformly on sufficiently large spatial volumes V. This suggests taking a spatial average of Einstein's Eqns. (1) and defining The objects appearing in this equation are components, so care must be taken to ensure that the averages are covariant. Suppose we perform averages with an RW metric defined 2 by This is reasonable because a spatially averaged quantity is isotropic and homogeneous, by construction. At early times, we expect g µν ≃ g RW µν because gravitational collapse has not yet introduced spatial dependence into the matter distribution. The nonlinearity of the Einstein tensor as a function of the metric, however, implies that in general. Here we have used square braces to indicate that G µν is a differential operator acting on the metric. The discrepancy is expected to grow larger as gravitational collapse forms inhomogeneous, small-scale structures. If one defines the difference then the spatially averaged Einstein's Eqns. (8) become This additional term X µν is called the cosmological backreaction. It has been claimed by Räsänen (2004), Kolb et al. (2006), and others that the cosmological backreaction will manifest as an apparent P(η) = −ρ(η) perfect fluid. This purported ability to explain the accelerated late-time expansion of the universe, by averaging the left side of Einstein's equations, has generated significant interest in the backreaction community.
There is vigorous disagreement, most recently between Green & Wald (2014) and Buchert et al. (2015), about the relevance of cosmological backreaction. We believe that this disagreement can be distilled into a single question: what is the correct way to construct a large-scale homogeneous model of a small-scale inhomogeneous universe within the framework of GR? It is clear that construction of any isotropic and homogeneous global model from Einstein's local differential equations requires removing position dependence. Indeed, both the Birkhoff approach and the backreaction approach take spatial averages. The lack of any unique prescription for which averages to perform, however, is a fundamental ambiguity in our understanding of how to use Einstein's equations to extract physics on the largest scales.
In this paper, we show how to resolve the ambiguity. We develop an internally consistent perturbative treatment of both the metric and the stress tensor directly from the Einstein-Hilbert (EH) action. The main results of this paper are as follows: • Friedmann's acceleration equation and conservation of stress energy without effective source terms of geometric origin, together with; • an explicit expression for P(η), which may not vanish at late times.
These results are to be contrasted with expectations, respectively, of an additional backreaction term and a pressure-free stress tensor via Birkhoff's theorem. The rest of this paper is organized as follows. In §2, we precisely define our gravitational model within the context of GR. We then derive the appropriate equations of motion, taking great care to enforce the RW model symmetries at each step of the calculation. In §3, we determine how various astrophysical sources, including strong ones, contribute to the equations of motion. In §4, we highlight the observational consequences of P(η) 0. In §5, we summarize our results. In Appendix A, we provide a detailed discussion of boundary terms. In Appendix B, we provide additional calculations, which support §3.3. Throughout this paper, we set c ≡ h ≡ 1 and use the (−, +, +, +) signature.

FRIEDMANN'S EQUATIONS FROM THE ACTION
We adopt GR as encoded in the EH action. Everywhere on the manifold M, there exists a well-defined 3 metric g and a well-defined collection of other fields Ψ. These combine to give a well-defined stress-tensor T. We note that the model builder's arbitrary choice for approximating g does not change the actual Ψ. We make these statements precise, beginning with the following Assumption 1 There exists a series representation (this includes the coordinates) for the actual metric which converges on some compact submanifold (with boundary) U ⊂ M. The components h (n) µν are at least twice continuously differentiable (C 2 ) on U.
This representation need not be unique. For the purposes of all subsequent discussions, however, we only require the existence of such convergent representations. If U is not connected, then we take Assumption 1 to mean that there exists 3 i.e. unique and representation independent such a representation on each connected component. The factors of ǫ n encode the relative magnitude of the terms in the expansion. In other words, for fixed µ and ν, each component −1 h (n) µν 1 for all n. This guarantees that the representation is dominated by a convergent power series. Truncation of Equation (14) at N = 1 produces Equation (2). In this sense, Equation (14) generalizes the covariant linear perturbation theory metric ansatz to arbitrary order.
Formal substitution of Equation (14) into Einstein's equations will give a formal power series in ǫ. We do not know, a priori, the relative magnitudes of these terms or if this series converges. This occurs because the Einstein tensor and the stress tensor are proportional and we have not constrained the matter fields Ψ. Since we wish to solve Einstein's equations perturbatively, order-by-order in ǫ on U, we make the following Assumption 2 There exist length scales L n (η, x) such that the Einstein tensor can be expanded as where Ḡ (n) µν 1 and ǫ is the same as in Equation (14).
Following Lin & Segel (1988, §6.3), these length scales are combined, at each order, with length scales set by derivatives of Ψ within the stress tensor to produce dimensionless field equations. The existence of L n is reasonable because the Einstein tensor contains at most second derivatives of the metric, which we have defined to be at least C 2 . The components of the Einstein tensor are thus continuous on a compact set, and thus bounded. Note that L n 0, as we have assumed that Einstein's equations are well-defined everywhere on M, and U ⊂ M. The particular finite values of L n do not matter. Assumption 2 guarantees their existence so that field equations, at each order in the perturbation, can be well-defined.
Taken together, Assumptions 1 and 2 establish sufficient conditions for Einstein's equations to be solved perturbatively on U. They do not assert that low-order approximations are always useful. On the contrary, it may be that infinitely many terms are required to describe a region of spacetime. Our concern is not with utility, but with existence of well-defined representations for g. This will guarantee that equations governing the lowest-order terms remain valid even on submanifolds of U with strong gravity. As is standard in the literature, we will indicate the order N at which field equations are truncated with O ǫ N , even though we write dimensionful field equations. No confusion will arise, as Assumptions 1 and 2 guarantee that a perturbative solution is always possible.

The Action Principle for Gravity
Given Assumptions 1 and 2, we work exclusively on U. In the context of gravity, the action principle demands that where S G is a gravitational action, S M is a matter action, and Ψ is held fixed for the variation δ. By the definition given in Weinberg (1972, §12.2), the symmetric, rank-(2,0), stress tensor T µν enters Einstein's equations through Here, g is the metric determinant. Note that we have flipped indices, relative to Weinberg's definition, and thus incur a minus sign through the metric variation Weinberg's definition encodes the standard assumption that S M comes from a Lagrange scalar density L M , Comparison of Equation (17) with Equation (20) reveals that, under the integral, one may identify The stress tensor components are, therefore, a mixture of the metric and the non-gravitational contributions Ψ(η, x). The non-gravitational contributions Ψ(η, x) have no explicit dependence on g, and thus are unaltered by the model builder's choice of metric model. Of course, the Ψ(η, x) are implicitly constrained by the eventual equations of motion. Any premature constraint, before equations of motion consistent with the action principle are determined, risks introducing inconsistencies into the field equations. Let R be the Ricci scalar, built from the Levi-Civita connection of any metric. The EH gravitational action is We use the representation of the metric given in Equation (14). For the purposes of computation, we now parameterize U. While there is observational evidence that the universe began in a nearly singular hot and dense point, it is not clear how far back in time one can "safely" apply GR. We thus restrict consideration to Here η i is some arbitrary initial time and η f > η i can either be arbitrary, in the case of matter domination, or any time below the asymptotic value in the case of dark energy domination. Note that because U is compact, an action defined over it will never diverge. For the moment, consider a closed 3-ball V ⊂ R 3 so that We will often use V to represent both this compact submanifold and its volume, which is the integral of the Euclidean volume form over this submanifold. We now work only through O(1). From Equation (14), to leading order, the approximation is This gives, for the metric determinant and inverse metric variation, Following Wald (2010, Equation (D.9), p. 446), we may exploit that the flat RW metric is a conformal rescaling of flat space to write where the derivatives are simple partials. Note that we have used Assumption 2 to regard the Ricci scalar of the bare (not multiplied by the scale factor) metric as O(ǫ). Substitution into the definition of the gravitational action S G gives Performing the variation gives plus terms O(ǫ). Integrating by parts twice gives where we have separated off a total divergence term This term becomes a boundary term, which is usually discarded. We must be more careful in our RW setting, however, because we do not have spatial control of the variations δa(η). Detailed consideration of the boundary term can be found in Appendix A, where we show that the term should still be discarded. Consistent with this result, we now set B = 0.

Symmetry of the Variation δa
We now arrive at the crucial step, which reveals the importance of working directly from the EH action. Only variations consistent with the RW model symmetries are permitted. In other words, the variations δa of the scale factor a(η) must be functions of conformal time alone: δa(η). Since the integrand in Equation (31) depends only on η, we have that Here we have used Fubini's theorem to write the compactly supported integral as an iterated integral. Note that the spatial integration is just the comoving volume V, so To compute the variation of the matter action, first note that where T µ µ is the trace of the stress tensor. Substitution into Equation (17), restricted to U, gives Again, because a and δa depend only on η, we have that The trace of the stress tensor cannot be moved through the integral. The trace continues to be position-dependent, even though the metric at zero order is position-independent, because the stress tensor depends on Ψ(η, x). These fields are not explicitly constrained by our (or any) choice of metric model. The quantity in the square brackets, because of the integration over V, is formally independent of position. Note, however, that it could be sensitive to the choice of the region V. This will be elaborated upon in §2.3, where we discuss the physical significance of V. We now combine the variation of the gravitational action Equation (34) with that of the matter action Equation (39). The action principle as stated in Equation (16) Applying the Fundamental Theorem of Variational Calculus then gives the consistent equation of motion Dividing by the constants on the left naturally reveals which is the spatial-slice average of the stress tensor's trace. This expression is manifestly coordinate-invariant: • the trace of any tensor is a general coordinate scalar; and • by Assumption 1, the notion of a spatial slice is geometrically well-defined (e.g. O'neill 1983, §12, p. 342) The spatial integration in this average, and the volume to normalize it, come directly from the action principle. Expanding the derivatives and trace gives Friedmann's acceleration equation expressed in conformal time. Though the trace is invariant, ρ(η, x) and P i (η, x) are the diagonal components of the generic stress tensor in the RW preferred coordinate system. This is the coordinate system in which the metric takes the form of Equation (14). One may introduce the definitions to arrive at the typical presentation of Friedmann's acceleration equation in conformal time. Note, however, that these defined ρ(η) and P(η) must be spatial-slice averages of the unapproximated quantities.

Significance of the Submanifold V
The spatial average contains the arbitrary closed 3-ball V. This enters the construction because we must guarantee that the action integral remains well-defined. This means that Equation (44) represents a continuum of Friedmann models.
We may understand the parameter V as follows. Fix a radius b and pick some event P ∈ U. Compute the spatial averages, at all η, of the stress trace over a round ball of radius b centered at the spatial location of P. The resulting source will determine how a(η) evolves in conformal time.
With the same b, repeat this procedure at a distinct event Q ∈ U, which lies in the same spatial slice as P. If the Q spatial averages do not agree with the P spatial averages, the dynamics for a(η) need not agree. This can happen because matter is not distributed uniformly on small scales. In other words, well-defined field equations for a(η) require a volume large enough so that the averaged density and pressures are position-independent.
Observations (e.g. Hogg et al. 2005;Scrimgeour et al. 2012;Nadathur 2013) suggest that, at the present epoch, there is a radius beyond which the averaged quantities are position-independent. We will assume that V contains a round ball of at least radius b * . Note that the specific value of b * is unimportant for the purposes of our argument. Strong evidence for hierarchical formation of structure tells us that, throughout time, these averaged quantities remain position-independent. In other words, any such V suffices for all earlier epochs.

Resolution of the Averaging Ambiguity
Friedmann's equation, as given in Equation (44) with V interpreted according to §2.3, resolves all ambiguity in our understanding of how to use Einstein's equations to extract physics on the largest scales. Consistent with our motivation of cosmological backreaction in §1, the action includes an integration of degrees of freedom present in the metric over the spatial volume V. At zero order, however, the only such degree of freedom is the scale factor a(η). It is already constrained by the RW model assumption (not the universe) to be position-independent. The required integration produces the spatial-slice volume V, as seen in Equation (34). In other words, given Assumptions 1 and 2, the action principle produces equations of motion without any cosmological backreaction X µν (η).
It is important to note that, given a different starting ansatz, one can arrive at a different cosmological model. This different cosmological model could possibly exhibit a backreaction: additive and position-independent contributions of geometric origin. The presence or absence of such terms is an artifact of the particular model and is logically distinct from the actual metric tensor on U. Our Assumption 1 reduces to the RW metric when truncated at N = 0, and to the typical covariant linear perturbation metric when truncated at N = 1. What we have shown is that, within these frameworks, there is zero cosmological backreaction X µν (η).
Inspection of Equation (44) and Equation (46) reveals that the position-independent pressure P(η) of any Friedmann model must include all contributions from pressures interior to compact objects on V. In other words, given Assumptions 1 and 2, using Birkhoff's theorem to construct a cosmological source to Friedmann's equations is inconsistent with the action principle. In general, use of Birkhoff's theorem in cosmological contexts is only appropriate under very restricted settings, which we derive and discuss in §3.1. We summarize these statements for clarity.
Our Friedmann's equation is an unambiguous and coordinate-invariant consequence of Assumptions 1, 2, and the EH action.
• There is no backreaction. Inhomogenieties in the spatial distribution of matter do not affect the source to Friedmann's equation.
• Pressures interior to all compact objects contribute to Friedmann's equation.
These consequences arise because the RW metric cannot distinguish spatial regions. Perturbatively consistent equations of motion thus require a zero order source without explicit or implicit notions of interior.

Covariant Conservation of Stress-Energy at Zero Order
It is instructive to show how careful application of symmetry considerations also produces the consistent conservation of stress-energy statement. To derive the appropriate relation, note that a general coordinate scalar, such as the matter action, cannot change under coordinate transformations. It is shown in Weinberg (1972, §12.3) that this restriction implies where ǫ λ is an arbitrary (infinitesimal) vector field. Since ǫ λ is arbitrary, the above equation immediately leads to the familiar statement of covariant conservation of stress-energy: Note that we deviate slightly from Weinberg's notation and use ∆ to distinguish the change in S M induced by a coordinate change from any hypothetical change induced by the variational differential δ. We continue to interpret V consistent with §2.3 For clarity, we will continue to explicitly write the time and position dependence of the generic stress tensor. At zero order, the position-independent vector fields ǫ λ (η) are the only ones permissible by RW model symmetries. For consistency, we must restrict ourselves to these vector fields. Note that Equation (48) is true for arbitrary ǫ λ , so it is true for position-independent ǫ λ (η). Thus, we may commute though the spatial integral as before: This is permissible because the preferred RW frame is a function of η alone. We now expand the derivative term and apply Stokes' theorem on the spatial slice Because we are interested in the conservation condition, set λ = 0 so that the vector T k 0 is an energy flow. By Einstein's equations by Assumption 2, so the boundary term must be dropped at this order A standard calculation gives where H η ≡ a −1 da/dη. All of these connection components are dependent on η alone, so that dividing through by V gives which is the appropriate continuum of conservation of stressenergy statements.

Generalizations
The most general RW metric takes the form where γ i j can be the metric for any 3-space of constant curvature and b(τ) is a gauge degree of freedom. In the previous computation, we have fixed the gauge by defining If one keeps the b(τ) degree of freedom, and varies with respect to it, the typical Friedmann energy equation is obtained (see Suzuki et al. 1996, Equation (9)). In this way, one can reproduce Equation (54) entirely within the variational formalism. Non-flat spatial slices lead to the same conclusions with respect to computation of the source terms. The formalism developed in this section has also been applied to the next order in ǫ. The advantage is that source terms to the field equations are generated unambiguously. This has been verified by Croker (2018, §3) at O(ǫ) for the scalar modes, where a lengthy calculation in longitudinal gauge reproduces the standard equations.

APPLICABILITY OF THE METRIC REPRESENTATION TO PHYSICAL SOURCES
In the previous section, we showed that if Eqns. (14) and (15) converge on some region U, then all densities and pressures in U will affect a(η) according to Equation (44). In this section, we show that Equation (44) continues to hold on small scales and inside particular strong sources.
We will proceed in four steps. First, we will determine a spacetime region u ⊂ U where typical strong, local gravity solutions can be used consistently in a cosmological setting. Second, we consider a Kerr BH, which requires extending the results of §2 to additional boundaries. Third, we consider a Schwarzschild constant-density sphere because this solution dominates many solutions of physical interest. Finally, we consider a static de-Sitter sphere, because generalizations of such spheres have been proposed as BH replacements.
= a 2 (η 0 ) If we define a reciprocal timescale then the series expansion of Equation (14) about η 0 gives Now choose ∆η such that the correction is small For any fixed ǫ < 1, the O ∆ηH 0 η term allows the approximation to be cut off at finite N N < For physical reasons, we may choose a spatial volume v ⊂ V, whose light-crossing time is less than ∆η v 1/3 < ∆η.
Note that mathematically, however, the following arguments will not depend on this spatial restriction. We will restrict our attention to η ∈ [η 0 , η 0 + ∆η] and define a spacetime domain Note that u is a spacetime "hockey puck," which encloses a region of the spatial slice V required for the spatial averages in Equation (44). Under a constant rescaling of the representation Equation (61) becomes This is the same metric, so the dynamics are unchanged. In other words, it may be possible to use asymptotically flat models within a universe described by Equation (14). Conclusions drawn from such models, however, are only valid on the restricted spacetime u. In particular, conclusions drawn on u cannot be trivially extended off of u. We summarize for clarity. If the real world satisfies Assumptions 1 and 2, then for short intervals asymptotically flat models are permissible, but only for short intervals. This is another reason why Birkhoff's theorem cannot be used to infer cosmological behavior. Given a sequence of times {η j } and suitable symmetry, one can use Birkhoff's theorem to construct a sequence of Schwarzschild exterior spacetimes characterized uniquely by a sequence of masses {M j }. Without global knowledge of the actual metric, however, no relation between the M j can be established.
We emphasize, however, that the converse relation holds. In other words, it is always correct to proceed globally from Equation (14) to locally, as represented by Equation (67). This follows because u ⊂ U. This means Friedmann's equation remains valid at all events in u, where the metric takes the form given in Equation (67). This means that any source on u, which produces a metric of the form Equation (67), will influence the dynamics of a(η) at η 0 , in the manner given by Equation (44).
For all subsequent discussion, we fix some time η 0 and we fix some ǫ < 1. Because we have the large-scale solution by virtue of Assumptions 1 and 2, we have a well-defined H 0 η (η 0 ). Now choose ∆η and N > 0 so that Equation (63) is satisfied. We now have a well-defined u given by Equation (65). Strong sources defined in u will the contribute to the dynamics of a(η 0 ) through the spatial average over u at η 0 , as required by Friedmann's equations given in Equation (44). We now study some specific sources of interest.

The Kerr Spacetime
LIGO (e.g. Abbott et al. 2016) and the Event Horizon Telescope Collabo (2019) have established the existence of ultracompact objects consistent with the Kerr geometry. In this section, we will establish that the cosmological formalism developed continues to operate near to ultrarelativistic, spinning, sources. Since the Kerr spacetime is asymptotically flat, the considerations of §3.1 necessarily restrict use of the Kerr geometry to timescales ≪ 1/H 0 η . In particular, we determine an appropriate submanifold u ′ ⊂ u defined in Equation (65), where the Kerr spacetime satisfies Assumptions 1 and 2 in the approximate sense of Equation (67). For a convenient reference on the Kerr spacetime, we refer the reader to Visser (2007).
In Kerr-Schild coordinates, the Kerr metric is where R s is the Schwarzschild radius for some fixed mass M, A encodes dimensionful information about the spin, 4 and r(x) is a function of the position coordinates x, y, and z. We will henceforth drop the explicit indication of x dependence for r. The covector ℓ µ is defined to be null with respect to η µν and has the form Since η µν ℓ µ ℓ ν ≡ 0, it follows that which implicitly defines r. Since each term on the right hand side is positive, we see that the components of ℓ µ satisfy |ℓ µ | 1.
It follows immediately that To avoid a naked singularity, we require that which maintains real-valued solutions for the horizon surfaces. By Equation (72), a necessary and sufficient condition to write the metric in the form of Equation (67) is This restriction is saturated when ǫr 4 − R s r 3 + ǫA 2 z 2 = 0.
By Descartes' rule of signs, this equation implicitly defines two (or zero) non-intersecting surfaces r ± (x). As ǫ → 1, r + and r − become, respectively, the outer and inner ergosurfaces of Kerr. Define a submanifold u ′ ⊂ u such that Then u ′ has an additional boundary, which becomes the outer ergosurface as ǫ → 1. We now establish that this additional boundary will not alter Friedmann's equations. As discussed in Appendix A, the dynamically relevant contribution to Friedmann's acceleration equation from this boundary is the second term in Equation (A4). This term continues to vanish at zero order on u ′ because of the RW model symmetries. As shown in Equation (51), an additional term also appears in the conservation equation Since Kerr is a vacuum solution, however, T k 0 = 0 on ∂u ′ and this term also vanishes. Now that we have a suitable domain, we may write Equation (68) in the form of Equation (67): In summary, we have established a domain u ′ where the Kerr solution can be accommodated under Assumptions 1 and 2. This implies that the spatial averages appearing in Friedmann's equation remain valid in the ultrarelativistic vicinity of a Kerr BH. The domain u ′ can be taken to include all but the region enclosed by the outer ergosurface. Consequently, for the domain we have constructed, a Kerr BH contributes nothing to Friedmann's equations. Local observers in u ′ will still perceive a Kerr BH with mass M and spin A for r + < r < v 1/3 and η ∈ [η 0 , η 0 + ∆η]. This is unsurprising, because the cutting procedure we have performed is mathematical. Visser (2007) emphasizes that the inner horizon and the enclosed inner ergoregion are extremely pathological and should not be regarded as physically relevant. Our result suggests that, additionally, the entire region below the outer ergosurface should be replaced with a distinct interior solution.

Cosmological Contribution of Interior Solutions: Typical Astrophysical Sources
In this section, we show that the interior region of Schwarzschild's constant-density sphere satisfies Assumptions 1 and 2. The exterior region of Schwarzschild's constant-density sphere is an A = 0 Kerr solution (i.e. Schwarzschild's BH), which has already been treated. The constant-density solution can be chosen to dominate the actual energy densities for many spherically symmetric, static, sources. We will show that, for such objects with physical radius R > 3GM, all energy densities and pressures in u will influence a(η) according to Equation (44). This will establish that localized pressures, interior to a large class of compact and relativistic objects, contribute to the global Friedmann average.
The functions Φ and Λ are defined in terms of the Schwarzschild radius R s and physical radius R of the object We will show that there exists a natural ǫ such that the coefficients of a Taylor expansion of Equation (80) take the form of Equation (67).

Denote the Taylor expansion coefficients of
Furthermore, extend the notation to negative indices by Then it can be shown that the expansion of Equation (81) in The terms in these infinite sums can be commuted and rearranged because √ 1 − x converges uniformly for |x| < 1.
To guarantee that the coefficient of the linear term remains bounded by 1 on r ∈ [0, R], we must have Since ǫ < 1, this becomes a restriction on the object's physical radius Note that this is the radius of the innermost stable orbit for a photon. 5 Saturating this bound, we find To bound the remaining coefficients of ǫ n , note that for each coefficient in n, ∞ j=0 q j q n− j 2 3 The inequality follows because q n− j < 1 for all j. Since √ 1 − x 1 for x ∈ [0, 1], each coefficient of ǫ n for n 2 is bounded above by 1 on r ∈ [0, R]. This was to be shown. The result for exp(2Λ) follows immediately from use of the binomial expansion, followed by substitution of Equation (87). In Appendix B, we establish that the Schwarzschild constantdensity sphere also satisfies Assumption 2.
In summary, we have established that the Schwarzschild constant-density sphere is of the form Equation (67). In other words, all pressures interior to any astrophysical object bounded by the Schwarzschild constant-density sphere solution will contribute to the cosmological average in Equation (44).

Cosmological Contribution of Interior Solutions: Static de-Sitter Sphere
In this section, we show that the interior of an isolated de-Sitter sphere satisfies Assumptions 1 and 2. The exterior region of such a sphere is an A = 0 Kerr solution (i.e. Schwarzschild's BH), which has already been treated. The isolated de-Sitter sphere is the simplest model of a GEneric Object of Dark Energy (GEODE). Such objects, like the solution of Dymnikova (1992) or the gravastar of Mazur & Mottola (2015), have been proposed as possible BH replacements. Related, but dynamic, GEODEs called "vacuum bubbles," have also been considered as inflationary relics by Berezin et al. (1987). We will demonstrate shortly that the strong negative pressure inside the de-Sitter sphere will influence a(η) according to Equation (44). This will establish that a physically realistic GEODE could contribute to the cosmological P = −ρ.
Consider a de-Sitter patch in static coordinates We have introduced primes on the η ′ and r ′ coordinates for reasons that will become apparent. Denote these coordinates by Ξ ′ . This patch admits a timelike Killing vector field for r ′ < R s . This means that a sphere of Dark Energy with radius R < R s is static, regardless of the coordinate representation. Again, expand in a series, then multiply and divide by ǫ: We see that h (n) µν 1, and so Equation (92) satisfies Equation (67), provided that In other words, Equation (91) has a convergent representation on B 3 (0, R s ǫ). Contrast this situation with the Schwarzschild BH, which converges on v / B 3 (0, R s ǫ −1 ). Because Einstein's equations hold everywhere on U, we know there exist coordinates ζ, which cover r ′ > R s ǫ. This means there exists z > 0 such that is the overlap between the Ξ ′ and ζ charts. We will continue to use the Ξ ′ chart on this region. The relation between the Ξ coordinates of the exterior Schwarzschild spacetime and the Ξ ′ coordinates of the de-Sitter sphere is determined by Mazur & Mottola (2015, Eqns. (5.1-5.3)). They find that Ξ and ζ are the same, and that Ξ ′ satisfies In other words, the spatial slices are unaltered, but time interior to the sphere runs twice as fast. This means that our static approximation is only valid for ∆η/2. Apart from a reduction in v to v/8, nothing changes with respect to the spatial-slice integration. We may integrate Equation (97) to find where C is an integration constant. Note that we may choose C ≡ −η 0 so that a η ′ (η 0 ) = a(η 0 ).
This shows that the unit redefinition used to write Equation (67) remains consistent inside the sphere. We now cut out the r > R s ǫ − z/2 region, leaving a boundary at the cut radius. It again follows that Equation (A4) continues to vanish at zero order. The required conservation condition in Equation (50) continues to hold because T µ ν ∝ δ µ ν , and so the contribution vanishes identically. Consider the domain u ′′ ⊂ u We have established that Equation (67) is satisfied everywhere on u ′′ and that there are no additional boundary contributions to Equation (44). Finally, static de-Sitter space satisfies Assumption 2 because its stress tensor is constant. By adjusting z and ǫ in Equation (100), we conclude that nearly all of the de-Sitter region contributes P = −ρ to Equation (44).

COSMOLOGICAL ENERGY SHIFTS
Friedmann's equation in §2 clarifies that local pressure contributions, which do not vanish upon spatial averaging, affect the scale factor a(t). In §3, we established that the framework of §2 remains valid very near, and often within, ultrarelativistic sources. In this section, we show how any source that contributes to the cosmologically averaged pressure must itself evolve cosmologically.
In the following, we will use the word "object" to refer to bound systems like stars and clusters. How a particular object responds to the scale factor is entirely dependent on the object. To see this, consider a population of generic objects. Let P(η) denote the spatially averaged pressure over the population. Let ρ(η) denote the spatially averaged energy density over the population. Define the equation of state of the contribution to be Conservation of stress-energy, given in Equation (54), describes the temporal evolution of these averaged quantities From the definition of w(η), we may write Switching to the scale factor a as the independent variable gives the separable differential equation Now suppose that w is constant. Then we may integrate Equation (104) to find the standard result, For simplicity, suppose all of the objects have the same comoving energy E. The energy density can then be written in terms of the physical number density N of the object population The number density dilutes with the expansion because objects either comove or belong to structures that comove. Combining Equation (105) and Equation (107) gives This result is again expected. For example, if the objects are photons, w = 1/3 and so E ∝ 1/a, which is the photon redshift.
Our result, however, applies to all objects. This consequence follows from Assumptions 1 and 2, and the EH action. When w(a) 0, the averaged quantities ρ(η) and P(η) evolve cosmologically. Thus, any material that contributes non-vanishing pressure to Equation (44) and Equation (54) must also evolve cosmologically. The evolution of E persists when w(a) varies in time, as can be seen from Equation (104). For compact objects, the particular value of w is strongly dependent on the non-gravitational fields Ψ that define the material. The cosmological energy shift is completely unaffected by the spatial distribution of material in the universe.

Observations
The cosmological evolution of local energies is a necessary consequence of Assumptions 1 and 2, and the EH action. The essential question is now observational: can this effect be measured? Because this section is observational, it is sometimes clearer to use redshift z as a time variable. We thus convert Equation (108) to redshift using that a = (1 + z) −1 . Where appropriate, we use Planck collaboration cosmological parameters from Aghanim et al. (2018): H 0 = 67.4 ± 0.5 km Mpc −1 s −1 , Ω b = (2.24 ± 0.01) × 10 −2 h −2 ≃ 0.05, Ω Λ ≡ 1 − Ω m = 0.685 ± 0.007. We now consider specific model objects relevant for astrophysical observation, using fixed w approximations. We proceed in order of increasingly relativistic w.

Stars
Consider a population of typical stars, at fixed (comoving) coordinate positions. Each typical star will contribute an extremely small positive pressure to the cosmological average. This follows because the pressure is everywhere positive within a star, and so an integral over the stellar pressure cannot vanish. This is true even in simplified stellar models, where fluid packets are radially static. Each packet's contribution to the pressure is non-zero because gravitational momentum flux is not included in T µν . The equation of state w for a typical star can be approximated from the ideal gas law, where k is Boltzmann's constant, T is the temperature of the star, and m p is the proton mass. To develop an upper bound, set T ∼ 10 6 K, representative of core temperatures. Then Thus, the energy (as perceived by RW observers) of Sun-like stars cosmologically evolves as How large is this effect? Suppose a star is produced at z i = 2 and observed at z f = 0. Define the fractional shift in energy as then This change occurs over 10 Gyr from z i = 2 to z f = 0. In other words, it is dominated by other stellar processes and is thus unobservable. This also establishes that the effect is unobservable for any other material with w < w star .
What is the reciprocal effect on the zero-order expansion? Stars contribute ∼ 2Ω b /5 to Friedmann's equation. The cumulative adjustment to ρ(a) from a conservative first light of z i = 40 to z f = 0 is then ∼ −10 −8 . The effect is unobservable.

Galaxy Clusters
Consider a population of clusters, at fixed (comoving) coordinate positions. The velocity dispersion of clusters is roughly bounded by σ ∼ 10 3 km s −1 . This implies a w cluster ∼ 10 −5 by ideal gas arguments. Thus, the energy (as perceived by RW observers) of clusters cosmologically evolves as How large is this effect? A typical length scale for a cluster is ∼ 1 Mpc, giving a light-crossing time of ∼ 3 Myr. In redshift, this could be between z i = 1.001 and z f = 1. During this time, the fractional shift in energy of the cluster is For comparison, during a photon transit completing at z f = 2, the cluster energy shifts less than 10 −8 . In other words, the correction to any integrated Sachs-Wolfe effect through clusters is unobservable. Similarly, corrections to gravitational lensing are unobservable.
What is the reciprocal effect on the zero-order expansion? The energy density in galaxy clusters is certainly bounded above by Ω b . The cumulative adjustment to ρ(a) from a conservative first cluster of z i = 4 to z f = 0 is then ∼ −10 −6 . The effect is unobservable at present and upcoming sensitivities.

GEneric Objects of Dark Energy (GEODEs)
As discussed in §3.4, GEODEs are explicit GR solutions, which schematically resemble the static de-Sitter sphere. Before gravitational-wave observations of ultrarelativistic object mergers, GEODEs were of theoretical interest because they are often free of physical singularities and horizons. In other words, they are regular solutions for gravitationalcollapse remnants, which resolve the BH Information Paradox.
Consider a population of GEODEs positioned at fixed (comoving) coordinate positions. Recall that the equation of state of a de-Sitter sphere is This equation of state is appropriate if the GEODE edge does not contribute significantly to the population-averaged P. GEODE material is maximally relativistic: it saturates the dominant energy condition. The material is also under tension, instead of under compression: the sign of w is inverted from all previously considered cases. This leads to very useful consequences.
The energy (as perceived by RW observers) of the GEODEs cosmologically evolves as Unlike the cases previously considered, this shift is significant and acts to amplify the energy. In other words, GEODEs cosmologically blueshift. The effect is analogous to the photon redshift. How large is this effect? Consider a 3M ⊙ GEODE formed from stellar gravitational collapse at z = 1.5, near the peak in comoving star formation rate density (e.g. Madau & Dickinson 2014). Its mass, observed in a binary merger at z = 0.1, will be This compares favorably with the masses observed by LIGO. A thorough exploration of the observational signatures of a remnant GEODE population, with respect to gravitationalwave observatories, is given by Croker et al. (2019). They find that the scenario is consistent with the LIGO GWTC-1 observed mass function. The cosmological blueshift of a GEODE can achieve the required masses with the standard common envelope binary formation channel. It does not require low-metallicity regions or prohibitively large progenitor stars.
What is the reciprocal effect on the zero-order expansion? Suppose 1% of all stellar material collapsed to Population III GEODEs, instead of traditional black holes, For simplicity, assume that this happens at z ∼ 14: after the Dark Ages but before reionization. The shift in energy of this population, observed at the present day, gives rise to Note that ρ GEODE is essentially constant within Friedmann's equations. This follows because w dS = −1, so that Equation (105) gives no time evolution of the physical density. In other words, the GEODEs dilute in number ∝ a −3 while gaining in mass ∝ a 3 . The result is an apparent cosmological constant, which compares favorably with the observed Ω Λ . A first exploration of the observational signatures of Population III GEODEs, with respect to cosmological observables, is given by Croker (2018). They find that such a scenario is flexible enough to reproduce the observed late-time accelerated expansion, resolving the coincidence problem.

Discussion
We have shown that non-relativistic material exhibits no observable shift. For relativistic objects, however, the energy shift can lead to measurable consequences. Neutron star (NS) material is highly relativistic (e.g. Abbott et al. 2018), with an object-averaged 0.05 w NS 0.1. Given that binary pulsar dynamics are often measured to exquisite precision, it may be possible to measure the cosmological shift in pulsar systems.
Since the cosmological shift for positive w NS appears as an energy loss, one can use a result of Damour & Taylor (1991, Equation (4.1)) to estimate the shift in orbital period decay Using P b as reported by Weisberg et al. (2010, Table 3) for the Hulse-Taylor binary pulsar, we find at the present day. How does this shift compare to current measurements of the Hulse-Taylor system? From Weisberg et al. (2010, Equation (5) and Table 3), the kinematically corrected orbital period decay iṡ The shift in orbital period decay at the present, due to cosmological evolution of NS energies, would appear to be ∼ 1%.
The central value of Equation (125) agrees with the GR radiative loss prediction to 0.2%. So a cosmological effect of the estimated magnitude has not been observed. Since Damour & Taylor (1991) assumed a flat spacetime, however, it follows from Equation (67) that the metric is known only to O(HP b ). In other words, the estimated shift given by Equation (123) is dominated by the error introduced by working under an asymptotically flat approximation within an RW cosmology. An investigation of binary pulsar systems with cosmologically evolving mass is the topic of future work.
The contribution to the cosmologically averaged quantities, and the necessary local evolution, is highly dependent on the object model. An instructive example is given by Dymnikova (1992). Dymnikova's object contains a very nearly Dark Energy interior, but a spatial average over her object gives exactly P = 0. In her model, the "skin" acts like a vacuum vessel: positive pressures in the skin maintain the static character of the object. A population of Dymnikova's objects would behave as cosmological dust, and exhibit no cosmological shift.
At present, no known GEODE solution strongly rotates. As pointed out by J. D. Bjorken (2018, private discussion), vacuum cannot rotate because it has no privileged directions. Given clear evidence of Kerr exterior spacetimes, all of the spin must reside in material outside of, but extremely near, the de-Sitter core. GEODE solutions are often constructed to have a Schwarzschild exterior spacetime and thus cannot evolve cosmologically. Given the considerations of §3.1, and the consequences in §4, it seems essential to have adiabatically evolving object solutions that spatially asymptote to an arbitrary RW geometry. The construction of a realistic GEODE model is an open question that is beyond the scope of this paper.

CONCLUSION
In summary, we have derived the equations of motion for Friedmann cosmology, paying particular attention to the symmetry enforced by isotropy and homogeneity. Working directly from the EH action maintains this symmetry, with respect to the source terms for Friedmann's equations. Our main results, embodied in Equation (44), follow from only two assumptions. We assume that the metric and Einstein tensor have convergent series representations, in RW coordinates.
Contrary to some earlier literature, we find no influence on the source to Friedmann's equation from inhomogeneities in the spatial distribution of matter. Contrary to other earlier literature, we find that every pressure source contributes to Friedmann's equations. Specifically, the spatial average in Equation (44) includes the pressures interior to compact objects. Evolution of the averaged quantities then implies a local energy shift of a −3w in any material that contributes to the averaged pressure. Here, a is the RW scale factor and w ≡ P/ρ is the equation of state determined from the population-averaged pressure P and energy density ρ of the objects. In non-relativistic material, the shift is too small to observe. The shift can be observable in relativistic material, if the equation of state exceeds w ∼ 0.01. This result takes on particular significance when applied to Population III stellarcollapse remnants with P = −ρ interiors. A population of such stellar collapse remnants can shift in energy ∝ a 3 while diluting in number density ∝ 1/a 3 . The population-averaged energy density is then effectively constant and readily produces the observed Ω Λ .
where we ignore the factor 3/8πG for clarity. The contribution along the boundary of U can be broken down into two parts. An "end-cap" contribution from the future and past Cauchy surfaces (entire spatial 3-volumes at an instant of time) and a spatial 2-surface during the time interval considered The intersection of these sets is {η i , η f } × ∂V, but it is two-dimensional and so it is a set of zero measure. In other words, the integral in Equation (A2) is three-dimensional and the integrand must be non-infinite, so this portion of the boundary contributes nothing. Thus, we may write We have used Fubini's theorem to write the integrals as iterated integrals, permissible because U is bounded and a is assumed to be well-behaved. The second term in this sum vanishes identically. This follows immediately because a(η) and δa(η) are functions of time alone (isotropy and homogeneity), so only temporal derivatives survive. This leaves only η k0 terms, but these all vanish. We now expand the Cauchy surface "end-cap" contribution Physically, the usual procedure is to constrain the variations δa(η) at the two endpoints. This leaves a single term Mathematically, if we take complete control over the variations δa(η), we may always write them as mollified by a C ∞ bump on some compact subset of (η i , η f ). This guarantees control of all derivatives at the endpoints, and can be used to destroy this term. If we do not assert this level of control, this final boundary term is just the York boundary term. To see this, we follow Wald (2010, Equation (E.1.18)) where the varied term of the EH action, which gives rise to the York term, is U ∇ α v β g αβ ǫ. (A7) For our RW ansatz, we find for Wald (2010, Equation (E.1.16 By Stokes' theorem, we have where, for our RW ansatz, we find Note that h is Wald's notation for the induced metric on the hypersurface orthogonal to the timelike unit covector field n α . The resulting boundary term is ∂U v β n α g αβ e = 6V∂ 0 (aδa) where the final equality follows because δa vanishes at the temporal endpoints. This term is proportional to Equation (A6), which was to be shown. A term of this form remains because the gravitational Lagrangian, built from the Ricci scalar, contains second-order derivatives of the field variables. The established procedure, given by Wald (2010, Equation (E.1.42)), is to just extend the EH action with the negative of this term. Thus, this final boundary term vanishes by construction. The York term is often omitted from the action when working classically, as it is non-dynamical.

B. THE EINSTEIN TENSOR OF THE CONSTANT-DENSITY SPHERE
In this section, we explicitly show that the Schwarzschild constant-density sphere of §3.3 satisfies Assumption 2 of §2. First, note that the constant-density sphere is static and spherically symmetric. This implies that G µν = 0 when µ ν. By Einstein's equations, it suffices to consider T µν to develop a bound on G µν . First, we require a bound on the density of the sphere. From Equation (87), the minimum physical radius R must satisfy This implies a density ρ < 1 9πGR 2 s . (B18) we have that To demonstrate that G (k) ss (0) 1, we establish that c k , defined as monotonically decreases in k. First, note that all c n (k) are positive. We write the ratio of successive c n (k), at fixed n c n (k + 1) c n (k) = n + 2k 3(k + 1) .
Note that the ratio is unity when n = k + 3. For terms with n < k + 3, individual c n (k) monotonically decrease in k by the ratio test. We now partition off these terms: We need not assert that the finite sum itself is monotonically decreasing in k. We are only bookkeeping individual terms in the sum over n that monotonically decrease in k. As k increases, the number of such terms increases. It suffices to show that the summed tail contribution is monotonically decreasing in k. Once this is established, then all contributions to c k monotonically decrease in k. This will establish that the c k monotonically decrease in k.
We now proceed to show that the tail monotonically decreases in k. First, re-index the infinite series, and introduce notation for the coefficients of this infinite series where we have grouped each multiplicand as a binomial in m and terms decrementing k. Note that there are 2 k grouped terms in its full expansion. For d n (k) where m 3k, 2 k m k is then an upper bound on the product. The tail of this series contribution to Equation (B35) is then We bound the series with the integral and obtain ∞ m=3k 2 k+1 m k 3 (m+2k+3) k! < 2 k+1 3 (2k+3) (ln 3) k+1 .
This expression is also monotonically decreasing in k.
The head of the series in Equation (B36) contains terms d m (k), where m < 3k. We establish that these d m (k) monotonically decrease in m. The ratio of successive terms in m is Each term is largest when the denominator is smallest. This is achieved when j = m = 0, giving an upper bound for the product We conclude the terms monotonically decrease in m for fixed k. The multiplicands in the head d m (k) are largest when m = 0, giving the following upper bound This expression is also monotonically decreasing in k.
Because c k 0 always, this completes the proof that G (k) ss (0) 1 as required by Assumption 2. As a check on our work, central coefficients can be explicitly computed through parametric differentiation with respect to ǫ of Equation (B20) For example, evaluation at k = 1 gives which agrees with Equation (B32) evaluated at k = 1. As an additional check, the bound has been verified numerically through c 60 .

B.2. Peripheral Pressure Coefficient Bound
Now that we have established that the coefficients of the central pressure are bounded, we show that the coefficients for r > 0 remain bounded. This is done readily with series long division. We will show that each coefficient at r > 0 is decreased relative to its value at r = 0 and that coefficients remain bounded below by zero.
To prepare for the algorithm, note that Equation (B20) can be written using uniformly convergent binomial expansions as where the a n are defined as a n ≡ (−1) n 2 3 n 1/2 n .