{'fraction_non_alphanumeric': 0.05787496761937656, 'fraction_numerical': 0.03091270183921941, 'mean_word_length': 3.8850574712643677, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 25, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Dynamical relativistic liquid bodies I: constraint propagation

We introduce a new wave formulation for the relativistic Euler equations with vacuum boundary conditions that consists of a system of non-linear wave equations in divergence form with a combination of acoustic and Dirichlet boundary conditions. We show that solutions of our new wave formulation determine solutions of the relativistic Euler equations that satisfy the vacuum boundary conditions provided the initial data is chosen to make a specific set of constraints vanish on the initial hypersurface. Moreover, we prove that these constraints propagate. This article is the first step of a two step strategy to establish the local-in-time existence and uniqueness of solutions to the relativistic Euler equations representing dynamical liquid bodies in vacuum.


Introduction
Over the past two decades, a number of results that guarantee the local-in-time existence and uniqueness of solutions to the (non-relativistic) Euler equations that represent dynamical fluid bodies in vacuum have been established [7,9,10,18,25,29,30,46,47]. Recently, the first step towards extending these existence results to the relativistic setting have been taken. Specifically, a priori estimates for solutions to the relativistic Euler equations that satisfy the vacuum boundary conditions have been established for both liquids [41] and gases [21,24]. Unlike most initial boundary value problems where well-known approximation schemes can be used to obtain local-in-time existence and uniqueness results from a priori estimates, it is highly non-trivial to obtain existence from a priori estimates for dynamical fluid bodies in vacuum, whether relativistic or not. The main reason for the difficulty is the presence of the free fluid-matter vacuum boundary, which make it necessary to exploit much of the structure of the Euler equations in order to derive a priori estimates. This makes the use of approximation methods problematic since any approximate equation would have to possess all of the essential structure of the Euler equations used to derive estimates, and to find such approximations has proved to be very difficult.
The a priori estimates from [41] were derived using a wave formulation of the Euler equations consisting of a fully non-linear system of wave equations in divergence form together with non-linear acoustic boundary conditions. This system of wave equations and acoustic boundary conditions were obtained by differentiating the Lagrangian representation of the Euler equations and vacuum boundary conditions in time and adding constraints that vanish identically on solutions. A priori estimates, without derivative loss, were then established using an existence and uniqueness theory that was developed for linear systems of wave equations with acoustic boundary conditions together with Sobolev-Moser type inequalities to handle the non-linear estimates. This approach to deriving a priori estimates suggests a two step strategy to obtain the local-in-time existence, without derivative loss, of solutions to the relativistic Euler equations that satisfy the vacuum boundary conditions. The first step is to show that the constraints used to derive the wave formulation propagate; that is, to show that if the constraints, when evaluated on a solution of the wave formulation, vanish on the initial hypersurface, then they must vanish identically everywhere on the world tube defined by support of the solution. The second step is to establish the local-in-time existence and uniqueness of solutions to the wave formulation, which would follow from a standard iteration argument using the linear theory and Sobolev-Moser inequalities developed in [41]. Such solutions would, by step one, then determine solutions of the relativistic Euler equation that satisfy the vacuum boundary conditions thereby establishing the local-in-time existence of solutions representing dynamical relativistic liquid bodies.
The main purpose of this article is to carry out the first step of the above strategy. For technical reasons, we do not use the wave formulation from [41], but instead, we consider a related version that differs by a choice of constraints. This new wave formulation involves an additional scalar field that solves a wave equation with Dirichlet boundary conditions. We will not address the second step in this article. It is carried out in the separate article [38], where we prove a local-in-time existence and uniqueness theorem that establishes the existence of solutions corresponding to dynamical relativistic liquid bodies.
1.1. Related and prior work. In the non-relativistic setting, a number of different approaches have been used to establish the local-in-time existence and uniqueness of solutions to the Euler equations that satisfy vacuum boundary conditions. Important early work was carried out by S. Wu who, in the articles [46,47], solved the water waves problem by establishing the local-in-time existence of solutions for an irrotational incompressible liquid in vacuum. This work improved on the earlier results [11,37,49], where existence for water waves was established under restrictions on the initial data. Wu's results were later generalized, using a Nash-Moser scheme combined with extensions to earlier a priori estimates derived in [6], by H. Lindblad to allow for vorticity in [30]. This work was subsequently extended to compressible liquids in [29].
Due to the reliance on Nash-Moser, Lindblad's existence results involve derivative loss. By using an approximation based on a parabolic regularization that reduces in the limit of vanishing viscosity to the Euler equations, the authors of [9] were able to establish, without derivative loss, a local-in-time existence result for incompressible fluid bodies, which they later generalized to compressible gaseous and liquid bodies in [10] and [7], respectively. Existence for compressible gaseous bodies was also established using a different approach in [25]. For other related results in the non-relativistic setting, which includes other approaches to a priori estimates and existence on small and large time scales, see the works [1,2,4,8,16,20,19,22,23,31,32,33,34,36,43,48,50] and references cited therein.
In the relativistic setting, much less is known. For gaseous relativistic bodies, the only existence result in the most physically interesting case where the square of the sound speed goes to zero like the distance to the boundary that we are aware of is [40], which is applicable to 2 spacetime dimensions. However, based on earlier work by Makino [35] in the non-relativistic setting, Rendall established the existence of solutions to the Einstein-Euler equations representing self-gravitating gaseous bodies that are undergoing collapse [42]. For relativistic liquids, a local-in-time existence result involving derivative loss has been established in [44] using a symmetric hyperbolic formulation in conjunction with a Nash-Moser scheme.
We note that the use of constraints to establish the existence of solutions has a long history in Mathematical Relativity with, perhaps, the most well known and important application being the proof of the existence and geometric uniqueness of solutions to the Einstein equations, which was first established in [12]; see also [3,5,15,27,28,42] for related work when boundaries are present. We would also like to add that the work presented here and in [41] was inspired by the constraint propagation approach to the relativistic fluid body problem from [14].
1.2. Initial boundary value problem for relativistic liquid bodies. In order to define the initial boundary value problem (IBVP) for a relativistic fluid, we first need to introduce some geometric structure starting with a 4-dimensional manifold 1 M equipped with a smooth Lorentzian metric g = g µν dx µ dx ν (1.1) of signature (−, +, +, +). In the following, we let ∇ µ denote the Levi-Civita connection of g µν and Ω 0 ⊂ M be a bounded, connected 2 spacelike hypersurface with smooth boundary ∂Ω 0 . The manifold Ω 0 defines the initial hypersurface where we specify initial data for the fluid. The proper energy density ρ of the fluid is initially non-zero on Ω 0 and vanishes outside 3 . The initial hypersurface Ω 0 forms the "bottom" of the world tube Ω T defined by the motion of the fluid body through spacetime, which is diffeomorphic to the cylinder [0, T ] × Ω 0 . We let Γ T denote the timelike boundary of Ω T , which is diffeomorphic to [0, T ] × ∂Ω 0 . By our conventions, Γ 0 = ∂Ω 0 . The motion of the fluid body is governed by the relativistic Euler equations given by where T µν = (ρ + p)v µ v ν + pg µν is the stress energy tensor, v µ is the fluid 4-velocity normalized by In this article, we, for simplicity, restrict our considerations to the physical spacetime dimension of n = 4. However, the results presented in this article are valid for all spacetime dimensions n ≥ 2. 2 There is no need for Ω 0 to be either connected or bounded; non-connected Ω 0 correspond to multiple fluid bodies, while non-bounded components represent unbounded fluid bodies. 3 Here, by vanishing outside, we mean there exists a spacelike hypersurace Σ ⊂ M that properly contains Ω 0 and ρ vanishes in Σ \ Ω 0 .
p is the pressure, and ρ, as above, is the proper energy density of the fluid. Projecting (1.2) into the subspaces parallel and orthogonal to v µ yields the following well-known form of the relativistic Euler equations In this article, we restrict our attention to fluids with a barotropic equation of state of the form for some positive constants 0 < ρ 0 < ρ 1 , and 0 < s 0 < s 1 < 1. Since the square of the sound speed is given by the assumption (1.7) implies that 0 < s 2 0 ≤ s 2 ≤ s 2 1 < 1, or in other words, that the sound speed is bounded away from zero and strictly less than the speed of light.
The boundary of the world-tube Ω T , which separates the fluid body from the vacuum region, is defined by the vanishing of the pressure, i.e. p| ΓT = 0. By our assumption (1.6), this means that the proper energy density does not vanish at the boundary, and hence, there is a jump in the proper energy density across Γ T . Fluids of this type are referred to as liquids. In addition to the vanishing of the pressure, the condition v| ΓT ∈ T Γ T must be satisfied to ensure that no fluid moves across Γ T . These two conditions form the vacuum boundary conditions satisfied by freely evolving fluid bodies. Collecting these boundary conditions together with the evolution equations (1.3)-(1.4), the complete Initial Boundary Value Problem (IBVP) for a relativistic liquid body is: where (ρ,ṽ µ ) is the initial data, and n ν is the outward pointing unit conormal to Γ T .

1.3.
Overview. We fix our notation and conventions used throughout this article in Section 2 and in Appendix A, where a number of definitions and formulas from differential geometry are collected. In Sections 3.1 and 3.3, we define the primary fields and constraints, respectively, that will be used in our wave formulation of the relativistic liquid body IBVP. The Eulerian representation of our wave formulation, which includes the freedom to add constraints, is introduced in Section 3.4. We then state and prove Theorem 4.1 in Section 4. Informally, this theorem guarantees the constraints, when evaluated on solutions of our wave formulation, vanish in Ω T provided they vanish initially on Ω 0 , i.e. they propagate, and moreover, that solutions to our wave formulation for which constraints vanish correspond to solutions of the relativistic liquid IBVP. We then, in Section 5, make a particular choice of the constraints that appear in the evolution equations and boundary conditions that define our wave formulation in order to bring the total system into form that is favorable for establishing the existence and uniqueness of solutions. In Section 6, we introduce Lagrangian coordinates and express system of equations and boundary conditions from Section 5 in these coordinates. In the final section, Section 7, we briefly discuss how the linear existence and uniqueness theory for systems of wave equations with acoustic boundary conditions from [41] can be applied to the linearization of the Lagrangian representation of our wave formulation from Section 6, and moreover, how this will be used in the article [38] to establish the local-in-time existence and uniqueness of solutions to the relativistic liquid IBVP.

Preliminaries
In this section, we fix our notation and conventions that we will employ throughout this article; see also Appendix A where we collect a number of definitions and formulas from differential geometry.
2.1. Indexing conventions. We will need to index various objects. The conventions that will employ are as follows: to denote partial derivatives with respect to the coordinates (x µ ), and ∂ to denote spacetime gradients so that ∂f = (∂ µ f ) for scalar fields f = f (x µ ). More generally, for k ∈ Z ≥0 , we use ∂ k f = (∂ µ1 ∂ µ2 · · · ∂ µ k f ) to denote the set of all partial derivatives of order k, and similarly, ∂ k to denote the collection of derivatives of order k that are tangent to the boundary Γ T . We also let the collection of partial derivatives of order less than or equal to k.

2.3.
Raising and lowering indices. We lower and raise spacetime coordinate indices without comment using the metric g µν and its inverse g µν , respectively, while frame indices will be lowered and raised, again without comment, using the frame metric γ ij and its inverse γ ij , respectively; see (A.3) for a definition of the frame metric. We will have occasion to raise or lower indices using metrics other than g µν or γ ij . In these situations, we will be explicit about this type of operation. For example, given a metric m µν and a 1-form λ µ , we would define the raised version using m explicitly by setting λ µ = m µν λ ν .

Constraint terms.
To help encode the freedom to add constraints to evolution equations and boundary conditions, we reserve upper case Fraktur letters, e.. R, S, T, possibly endowed with spacetime indices, e.g. R ν , to denote maps that depend linearly on a set of constraints Z. More precisely, if Z is R N -valued, then is added to a boundary condition.

The Eulerian wave formulation
3.1. Primary and auxiliary fields. The primary fields for our wave formulation consist of a scalar field ζ satisfying ζ > 0 and a future pointing, timelike 1-formθ 0 =θ 0 µ dx µ . We use the primary fields to define a timelike 1-form, again future pointing, by which we complete to a coframe by introducing spacelike 1-forms Along with these 1-forms, we introduce a collection of scalar fields σ i k j . The set {θ I , σ i k j } defines the auxiliary fields that will evolve via simple transport equations. For latter use, we introduce a number of additional geometric fields beginning with the frame dual to θ i . Following our notation from Appendix A, we use γ ij and ω i k j to denote the associated frame metric and connection coefficients, respectively; see the formulas (A.3) and (A.5). Finally, we define a future pointing, timelike vector field by

3.2.
Recovering ρ and v µ . The fluid 4-velocity v µ will be shown to be recoverable from the primary fields {ζ,θ 0 µ } by normalizing the vector field ξ µ to get , (3.4) where in obtaining the second equality we used the fact that ζ > 0 and γ 00 = g(θ 0 , θ 0 ) < 0. Recovering the proper energy density is more complicated. The first step is to define the pressure as a solution p = p(λ) of the initial value problem where λ 0 > 0 and p 0 ≥ 0. To be definite, we set p 0 = 0 and λ = 1.
From standard ODE theory, we see that p = p(λ) is smooth for λ > 0, while from the IVP (3.5)-(3.6), it follows that p is strictly increasing, which in turn, implies the invertibility of the map R ≥1 ∋ λ → p(λ) ∈ R ≥0 . We will show that we can then use this map to recover the proper energy density from scalar field ζ by setting λ = ζ to give ρ = ρ p(ζ) .
3.3. Constraints. In this section, we define the constraints that will be essential for the definition of our wave formulation of the relativistic Euler equations with vacuum boundary conditions. We separate the constraints into bulk and boundary constraints, which are to be interpreted as being associated to Ω T and Γ T , respectively.
We collect together the following bulk constraints where we have setγ 00 = g(θ 0 ,θ 0 ) (3.21) and we are using R µ and Q µ to denote constraint terms in line with our notation set out in Section 2.5.
Boundary constraints: and as in above, n is the outward pointing unit conormal to Γ T .
3.4. Eulerian IBVP. The formulation (1.9)-(1.13) of the vacuum IBVP for the relativistic Euler equations is commonly referred to as the Eulerian representation. In this representation, the matter-vacuum boundary is free, or in other words, dynamical. This terminology is useful for distinguishing this form of the IBVP from the Lagrangian representation where the boundary is fixed. We will continue to use the Eulerian terminology for the wave formulation of the relativistic vacuum IBVP that we introduce in this section since the boundary is also free. Later, in Section 6, we will consider the Lagrangian representation where the boundary is fixed.
3.4.1. Eulerian evolution equations. Before stating the evolution equations, we first define a number of tensors that will be used repeatedly throughout this article starting with where m µν is a Riemannian metric and a αβ is a Lorentzian metric 5 that are defined by and , (3.26) respectively. We use W αµ to define and Here, we are using R µν to denote the Ricci tensor of the metric g µν , to denote the coordinate components of the 2-form F defined by (3.15), i.e. F = 1 2 F αβ dx α ∧ dx β , and h αβ is as defined previously (1.5). Using the above definitions, the evolution equations for our wave formulation are given by where Remark 3.1. The term R µ (∂ |1| χ, E) in (3.32) encodes the available freedom to add multiples of E λ and χ and its first derivatives 6 to the evolution equation ∇ v E µ = 0 in Ω T . This freedom provides flexibility in choosing the form of the equations of motion, which we will take advantage of in Section 5 to bring the equations of motion into a form that is useful for establishing the local-in-time existence and uniqueness of solutions. We further note that it is clear from the equations of motion (3.34)-(3.35) that we are free to add terms proportional to L v θ I ν and v(σ i k j ) and their derivatives to evolution equations (3.32)-(3.33).

3.4.2.
Eulerian boundary conditions. Before stating the boundary conditions for our wave formulation, we define There is nothing stopping us from adding higher order derivatives of χ to the evolution equation. The only thing that would change in the analysis below is that the class of solutions that we are dealing with would have to have enough regularity for the derivatives of χ to make sense.
with ν ναβγ = |g|ǫ ναβγ , (3.38) and ǫ, κ ∈ R constants to be fixed later. With B µ so defined, the boundary conditions are given by (3.44) encodes the available freedom to add multiples of B λ , χ h, k and their indicated derivatives 7 to the boundary conditions ∇ v B µ = 0 in Γ T . This freedom provides flexibility in choosing the form of the boundary conditions, which we will take advantage of in Section 5 to bring the boundary conditions into a form that is useful for establishing the local-in-time existence and uniqueness of solutions. We also note that it is clear from the equations of motion (3.34)-(3.35) and the boundary condition (3.45) that we are free to add terms proportional to L v θ I ν and v(σ i k j ) and their derivatives to the boundary conditions (3.44)-(3.45), and ζ − 1 and its derivatives tangential to Γ T to (3.44).

3.4.3.
Projections. The operators h µ ν and p µ ν define projections into the subspaces g-orthogonal to ξ µ and N µ , respectively, and hence satisfy Moreover, these projections commute, that is since it is clear from the definitions (3.3) and (3.39) that the vector fields N ν and ξ ν are g-orthogonal.
It is worthwhile noting that we can raise the lower index of h µ λ using either g νλ or m νλ since h µ λ v λ = 0, so that where h µν = g νλ h λ ν by our raising and lowering conventions. From the definition (3.42), it is then clear by (3.47)-(3.50) that Π µ ν defines a projection operator satisfying We denote the complementary projection operator by and setΠ µν = m µλΠν λ .
(3.54) 7 Again, there is nothing stopping us from adding higher order derivatives of χ and k to the boundary conditions.

Initial conditions.
In general, solutions of the IBVP problem consisting of the evolution equations (3.32)-(3.35) and boundary conditions (3.44)-(3.46) will not correspond to solutions of the relativistic Euler equations with vacuum boundary conditions given by (1.9)-(1.12). As we establish in Theorem 4.1 below, a solution {ζ,θ 0 µ , θ I µ , σ i k j } of (3.32)-(3.35) and (3.44)-(3.46) will only determine a solution, via the formulas (3.8)-(3.9), of (1.9)-(1.12) if the following conditions on the initial data are satisfied: and The constraints on the initial data (3.55)-(3.70) do not involve any constraints on the choice of initial data for the fluid, that is ρ and v µ or equivalentlyθ 0 µ and ζ, beyond the usual compatibility conditions that arise from the IBVP for the relativistic Euler equations with vacuum boundary conditions. The constraints that are unrelated to compatibility conditions for the physical fields involve the auxiliary fields, which are not physical and we are free to choose as we like; see Section 4.2 of [41] for details on how to choose the initial data for the auxiliary fields so that the above constraints are satisfied.

Constraint propagation
The precise relationship between solutions to the evolution equations and boundary conditions defined in the previous section and solutions to the relativistic Euler equations with vacuum boundary conditions is given in the following theorem. The main content of this theorem is that it guarantees that the constraints from Section 3.3 propagate for solutions of the evolution equations and boundary conditions defined in Sections 3.4.1 and 3.4.2 provided that the initial condition from Section 3.4.4 are satisfied. and The assumption (4.1) and the boundary condition (3.46) imply that v is timelike and tangent to the boundary Γ T , which in turn, implies that the set Ω T is invariant under the flow of v. From this fact, the transport equation (4.3) and the initial condition (3.56), it follows via the uniqueness of solutions to transport equations 8 that Propagation of e K in Ω T : From the definition (3.14) of e K , we have that in Ω T , and so by the uniqueness of solutions to transport equations since e K vanishes on the initial hypersurface, see (3.59). We note, with the help of the Cartan structure equations (A.17), that (4.7) is equivalent to We also observe that wherex 0 defines a coordinate on the interval [0, T ] and the (x Λ ) define "spatial" coordinates on Ω 0 . Without loss of generality, we may assumex 3 is a defining coordinate for the boundary [0, T ] × ∂Ω 0 so that dx 3 is an outward pointing conormal on the boundary [0, T ] × ∂Ω 0 . Lettingθ 3 = φ * θ 3 denote the pull-back of θ 3 by φ, we have from the choice of initial data, see (3.64), that for some positive function r > 0. Pulling back the evolution equation (3.34) for I = 3, we see from (4.10) and (A.9) that Lvθ 3 But this shows thatθ 3 defines an outward pointing conormal to [0, T ] × ∂Ω 0 , and hence, that θ 3 is an outward pointing conormal to Γ T and k =θ 3 − n = 0 in Γ T (4.13) by (3.23). Propagation of h in Γ T : By definition (3.2), θ j µ is the inverse of e µ j , and consequently, ǫ µαβγ e α 0 e β 1 e γ 2 , (e = (e µ j )), by cofactor expansion, which is equivalent to (4.14) by (3.38). Since ξ = e 0 by (4.5), we can, using the boundary condition (3.45), write (4.14) as 9 where N ν is defined by (3.39). But, we also have that where in deriving the last equality we used (3.17). Combining (4.14) and (4.16) yields Next, from (3.19) and (3.45), we have that Using this together with (4.5), we find, after a short calculation, that , we see, with the help of (3.1), (A.3) and the boundary condition (3.45), that applying ∇ ξ to −γ 00 ξ α yields where h αβ is as defined previously by (1.5). Then since by the anti-symmetry in the indicies of the volume form ν µαβγ , we can use (4.20) to write (4.19) as where s µν γ is as defined previously by (3.43). From the identity [e 0 , e j ] β = − d θ k (e 0 , e j )e β k , it then follows, with the help of the definitions (3.13)- (3.15), that which, in turn, implies that [e 0 , e j ] β = F (e 0 , e j )e β 0 (4.23) via (4.6)-(4.7). Substituting (4.23) into (4.21) gives while applying ∇ ξ to (4.17), yields, with the help of (4.24), where in deriving the last equality we used (4.4), (4.6) and (4.7). Rearranging (4.27), we find, by (3.1), (3.3) and (4.25), that Contracting the above expression with θ 3ν = θ 3 µ g µν , we find, via the anti-symmetry of F αβ and (4.4), that (4.28) But, by (3.39), (4.17) and (4.18), we observe that which we note, using the definition (3.36) and (4.4), can be used to write (4.28) as for some function a ∈ C 0 (Γ T ). Next, we observe, by (3.45), (4.4), (4.5), (4.6), (4.7), (4.13) and (4.26), that the boundary condition (3.44) can be expressed as for some collection of functions b µ ℓ ∈ C 0 (Γ T ), ℓ = 0, 1, and c µ lν ∈ C 0 (Γ T ). Contracting (4.31) with θ 3 µ , we see, with the help of (4.30), that h satisfies for some collection of function d ℓ ∈ C 0 (Γ T ), ℓ = 0, 1, and e µ ∈ C 0 (Γ T ). From the choice of initial data, see (3.62), (3.67) and (3.70), it then follows immediately from the transport equations (4.31) and (4.32) that and Propagation of j in Ω T : Contracting (4.34) withθ 0 ν shows, with the help of (3.45), (4.4) and (4.13), that the constraint j satisfies the boundary condition where we note that To derive an evolution equation for j , we first observe from (4.4)-(4.7) and (3.32) that E µ satisfies From this equation and the initial conditions (3.69), we then deduce that E µ = 0 in Ω T , and hence, by (3.27), that in Ω T . Contracting this equation withθ 0ν while using (4.4) shows that But solutions to the IBVP defined by (3.63), (3.68), (4.35) and (4.37) are unique by (4.36) and Theorem 2.2 of [26], and therefore, we conclude that ν follows directly from (1.5) and (4.38), we can, with the help of (3.33) and (4.38), write (4.37) as we calculate where a µν is the acoustic metric defined by (3.26). From this, we see that which, we observe, using (3.31) and (4.6), can be written as Using (A.20) to express the last term in (4.42) as where C β ων is defined by (3.30) and∇ is the Levi-Civita connection of the acoustic metric a αβ , we see from (4.39), (4.42) and (4.43) that d g − δ a F = 0 in Ω T . (4.44) Applying the codifferential δâ to to this expression, we see, with the help of (A.15), that g satisfies the wave equation (4.45) To complete the proof that the constraint g propagates, we need to show that g satisfies an appropriate boundary condition. We determine the boundary condition by first noting, see (A.12), that * g (θ 1 ∧θ 2 ∧θ 3 ) is g-orthogonal to the co-frame fields θ I . Since θ 0 is orthogonal to the θ I by (4.4), θ 0 must be proportional to * g (θ 1 ∧ θ 2 ∧ θ 3 ), and so q * g θ 0 = θ 1 ∧ θ 2 ∧ θ 3 in Ω T (4.46) for some function q. To find q, we wedge the above expression with θ 0 to get qθ 0 ∧ * g θ 0 = θ 0 ∧ θ 1 ∧ θ 2 ∧ θ 3 and then use the formulas (A.11) and (A.10) to obtain q = √ − det(γ kl ) γ 00 . Substituting this into (4.46) where in deriving the last two equalities we used the relations (4.6) and (4.7) along with the fact that while we see that In order to verify that the constraint F propagates, we need to show that F satisfies an appropriate boundary condition. To show this, we begin by noting that the equality is an immediate consequence of the relations (4.29), (4.33) and (4.40), and the boundary conditions (3.45)-(3.46).
Note: From this point until we establish the propagation of the constraint F , we will raise and lower spacetime indices with the acoustic metric a µν . After that, we will return to our standard convention of raising and lowering spacetime indices with the metric g µν .
From (4.4), (4.13) and (4.53), it is then clear that Π µ ν , see (3.42), projects onto the a-orthogonal subspace to the span of the vector fields n µ and v µ . Next, we observe from (3.45), (4.26) and (4.33) that (4.25) can be written as From this result and (4.27), it then follows, with the help of (4.13) and (4.34), that Furthermore, we observe, using (4.4)-(4.5) and (A.16)-(A.18), that ∇ ξ θ 0 ν is given by where we have set Since the frame components e k µ satisfy n µ e k µ = 0 in Γ T by (4.13), the vector fields e k must be tangent to Γ T , and consequently, e k (γ 00 ) = 0 in Γ T by (4.4), (4.40) and the boundary condition (3.45). Substituting this into (4.55) yields which, in turn, implies that since ω 0 0 3 = − i e3 ∇ e0 θ 0 > 0 in Γ T by assumption. Adding |∇ ξ θ 0 | −1 m ∇ ξ θ 0 ν to both sides of (4.58), we see that Applying ∇ ξ to this expression, we find after a short calculation that Using this result and (4.57), it is not difficult to see, with the help of (4.53), that we can write (4.54) as for an appropriate choice of functions {b K µ } ⊂ C 0 (Γ T ), where we have defined We then use (4.59) to expressn µ F µα a αβ F νβ ξ ν on the boundary Γ T as follows: We now claim that (4.60) defines the required boundary condition needed to show that the constraint F propagates. To see this, we introduce a time function τ that foliates Ω T and choose a future pointing timelike vector field τ ν such that τ ν ∂ ν τ = 1 in Ω T and τ νn ν = 0 in Γ T . We further define the level sets Since F satisfies the Maxwell's equations in Ω T by (4.51)-(4.52), we obtain the integral identity T µν = F µα F α ν − 1 2 a νµ F αβ F αβ is the stress energy tensor of F , and we have set It is important note that the matrix the 2-by-2 symmetric matrix Π J K = θ J α Π αβ θ K β is positive definite. Using the Divergence Theorem, we can express the boundary integral in (4.61) as (4.62) Noting that −τ ν ξ ν > 0, by virtue ofτ ν and ξ ν both being future pointing timelike vector fields, and |F | 2 m τ ν F νµ ξ ν by Lemma B.2, it is clear that the energy estimate is a direct consequence of (4.61)-(4.62). But this implies by Gronwall's inequality, and so we conclude that since F and Y vanish at t = 0 in Ω 0 and Γ 0 , respectively, by (3.60) and (4.56). For use below, we note, with the help of the Cartan structure equations (A.17) and (4.6), that (4.63) is equivalent to where s 2 = s 2 (−γ 00 ) − 1 2 , while ω 00J + ω 0J0 = 0 in Ω T (4.66) follows from applying e 0 to (4.4). We also note from (4.6), (4.8) and (4.64) that a short calculation shows that the above equations can be written as which, in turn, are easily seen to be equivalent to where

Choice of constraints
With the goal of establishing the local-in-time existence of solutions to guide us, we will, in this section, make particular choices for the constraints that appear in the evolution equations (3.32) and boundary conditions (3.44). The reason for these particular choices will be discussed in more detail in the following sections. Before proceeding, we introduce some notation. Given a vector field A µ , we define the projections A µ := Π µ ω A ω and A µ ⊥ :=Π µ ω A ω . Furthermore, given a tensor field B µν , we define vector fields B iν by B iν := θ i µ B µν , which, using our notation above, we can write the projected components as Lemma 5.1. Suppose ǫ ∈ [0, 1). Then there exists remainder terms R µ and Q µ such that 3) and all other quantities are as previously defined.
Proof. We begin establishing (5.2) by computing where here and below, in line with our conventions, we use Q and R to denote generic remainder terms that may change from line to line. Observing that we can write (5.11) as From the definitions (3.37) and (5.10), we further observe that the last two terms of the above expression may be expressed as (5.14) Noticing that , (by (4.15) and (3.51)) holds for any one form λ µ , we can rewrite (5.14) as Combining this with (5.13) and recalling the definitions (3.21) and (3.17) ofγ 00 and h, respectively, gives Using the commutator formula (A. 19) and the definition (3.24) of W αµ , we further observe that ∇ v W αµ is given by Next, we notice that where in deriving the last equality we used (4.15). But where (γ i j )=(γ i j ) −1 , β i =γ i j γ j 3 and α=γ 33 −γ i j γ i 3 γ j 3 , and so withp µν given by (5.9). The importance of (5.19) is that it allows us to express projected derivatives p µν ∂ ν on the boundary Γ T explicitly in terms of the vector fields v µ and e µ I , which will be shown to be tangent to the boundary Γ T in Section 6. Substituting (5.19) into (5.15) and multiplying the resulting expression by 1 1−ǫ Π τ µ +Π τ µ , we find, with the help of (3.42) and (3.51), that To proceed, we examine the term ∇ v ℓ µ in more detail by first noting that   25) where in deriving the last equality, we used the commutator formula (A. 19). It is then easy to verify using (5.10), (5.16), (5.21), (5.22), (5.24), (5.25) and the identity Turning to the bulk relation, we observe that For the wave formulation to be useful for establishing the existence of solutions, the only derivative of W αµ that can appear outside the principal term (5.27). We begin with L ei W iµ , which we can write as follows: Substituting (5.28) and (5.23) into (5.27) gives which is easily seen to be equivalent to (5.2).
The following corollary is then an immediate consequence of the above lemma and Theorem 4.1. (5.29) and the quadruple {ζ,θ 0 µ , θ J µ , σ i k j } satisfies the initial conditions (3.55)-(3.70) and
Since the Lagrangian coordinates are defined via the flow of the fluid velocity v, it follows that the components of the pullbackv = φ * v are given bȳ Substituting this into the transformation law J µ νv ν = v µ shows that φ µ satisfies where we note, see (3.21), thatγ 00 = g µνθ0 µθ 0 ν . In the Lagrangian representation, the map φ = (φ µ ) is treated as an unknown and (6.3) is viewed as an evolution equation for φ.
Pulling back the evolution equations (5.35) and (5.36) using the map φ, we see, with the help of formulas (6.2) and (A.9), and the naturalness property φ * L v = Lvφ * of Lie derivatives, that By (4.9), it is clear that φ satisfies φ(Ω 0 ) = Ω 0 from which it follows that ∈ Ω 0 , follows from our choice of initial data, see (3.11) and (3.56). The above two results together with (6.3), (6.4) and the transformation lawθ I µ = J ν µθ I ν yield the explicit representations of the unique solution to the evolution equations (5.35) and (5.36). We also recall, see (4.12), that defines an outward pointing conormal to the boundary [0, T ] × ∂Ω 0 by our choice of initial data.
(i) Since (6.5) represents the unique solution of the evolution equations (5.35) and (5.36) given the particular choice of initial data, we consider them solved and remove these equations from our system of equations. Furthermore, from the transformation lawθ I µ = J ν µ θ I ν and (6.5), it is clear that we can express θ I ν as which, in particular, shows that the components θ I ν (x 0 ,x Σ ) are determined completely in terms of the initial data θ I Λ (0,x Σ ) and the derivatives of φ µ , since by definition (J µ ν ) = (∂ µ φ ν ) −1 . (ii) The boundary condition (5.32) is automatically satisfied sinceθ 3 defines an outward pointing conormal to the boundary [0, T ] × ∂Ω 0 andv µθ3 µ = 0 follows immediately from (6.2) and (6.6). We, therefore, consider the boundary condition (5.32) as satisfied and do not consider it further.
The systems (6.18)-(6.24) along with the relations (6.5) defines the Lagarangian representation of the system (5.30)-(5.36). The unknowns to solve for are {φ µ , ζ,θ 0 ν , ψ ν , Ψ ν }. All other quantities in the system can be computed in terms of these fields and the time-independent frame fieldsθ I µ and functions σ i k j defined by (6.5). In particular, the Lagrangian representation of the metric and metric related quantities, i.e. Christoffel symbols and curvature, are to be viewed as functions of φ µ . To see this, we observe that if g αβ (x λ ) are the components of the Lorentzian metric, in the Eulerian representation, with associated Christoffel symbols Γ γ αβ (x λ ) and Curvature tensor R αβγ σ (x λ ), then in the Lagrangian representation, we have that Moreover, the Lagrangian representation of the components of any covariant derivative can be reduced to partial derivatives with respect to the Lagrangian coordinates using computations like the following: (i) The evolution equations (6.23) and (6.24) forθ 0 ν and ψ ν are easily seen to be equivalent to ψ ν = ∇ vθ 0 ν and Ψ ν = ∇v ∇vθ 0 ν |∇vθ 0 |m . So, even though we treat ψ ν and Ψ ν as independent variables, they are, by the choice of evolution equations, just labels for the derivatives ∇ vθ 0 ν and ∇v ∇vθ 0 ν |∇vθ 0 |m , respectively. We take this point of view because it turns out to be convenient for establishing the local-in-time existence of solutions via energy methods; see the following section for further discussion.
The importance of (6.34) and (6.35) is that they allow us to write the evolution equation (6.18) and boundary condition (6.19) as where This form of the equations is crucial for establishing energy estimates, see [41]. However, this form is not enough by itself to derive energy estimates, the following coercive type of elliptic estimate must also hold: there exists positive constants c 0 , c 1 such that , it is not difficult to verify, with the help of the evolution equations (6.22) that the boundary the coefficients and source terms {S µνγ , Q µν , P µν ,∂ |1| 0 Π λ ν , G µ } that appear on the right hand side of the boundary condition (6.19) involve only the fields { ∂ |1| φ µ , ∂|1|θ0 µ }, where in line with our previous notation, we are using ∂ to denote derivatives tangential to the boundary [0, T ] × ∂Ω 0 .

Remarks on local-in-time existence and uniqueness
In this section, we briefly remark on the properties of the system (6.18)-(6.24) that will allow us to use the existence results and energy estimates for systems of wave equations with variable coefficients from the articles [26,41] to prove the local-in-time existence and uniqueness of solutions to (6.18)-(6.24); the complete existence proof will be presented in the separate article [38]. The first step in the existence proof from [38] uses a contraction mapping argument that effectively reduces the problem of establishing localin-time existence of solutions to the non-linear system to that of establishing the existence of solutions to the linear system with variable coefficients that arises from "freezing" the non-principal coefficients of (6.18)-(6.24). Since this step is standard, we will not discuss it further here.
7.1. The modified system. Our strategy to obtain the existences of solutions with energy estimates for systems of the form (7.1)-(7.5) is to use the linear existence theory for systems of wave equations with acoustic boundary conditions that was developed in Section 7 of [41]. However, at the moment, the system (7.1)-(7.5) is not in the required form to apply this theory. So instead, we will consider a modified system and obtain solutions to that system first, and subsequently show these solutions determine solutions to the original system (7.1)-(7.5). As will be clear from Proposition 7.4 below, the modified system is obtained by performing an orthogonal transformation Ψ → EΨ followed by differentiating the evolution equation (7.1) and boundary conditions (7.2) in time, i.e. with respect tox 0 .
In order to define the modified system, we first introduced an orthogonal transformation E via the following orthonormal change of basis. Letting where m αβ is the inverse of the frozen version of the Lagrangian representation of the positive definite metric m αβ defined by (3.25) and we are using the notation m(X, Y ) = m αβ X α Y β . Setting we can use the Gram-Schmidt algorithm to complete {E µ 0 , E µ 3 } to an orthonormal basis We then define the matrices E = (E µ i ) and Θ = (Θ i µ ), which, by (7.9), are inverses of each other, that is ΘE = 1I, (7.11) where 1I is the 4 × 4 identity matrix.
Remark 7.1. It follows from the Gram-Schmidt process that E and Θ can be viewed as matrix valued maps that depend analytically on the vectors v andN . This allows us to write for some maps E(X, Y ) and Θ(X, Y ) that are analytic in their arguments.
Proof. From the definition (6.29) of Q µν and (7.12), we see that where q is defined as above by (7.25). Multiplying this expression of the left and right by Θ tr and Θ, respectively, we find, with the help of (7.11) and (7.13)-(7.14), that Θ tr QΘ = qΘ tr mΘEΠΘ = qP, which establishes the first formula (7.23). To establish the second formula (7.24), we first observe from (5.7), (6.30) and (7.11) that Multiplying this expression of the left and right by Θ tr and Θ, respectively, we get where in deriving this expression we used (7.14) and (7.16)-(7.17), and the proof is complete.
Proof. Writing I as and then differentiating with respect tox 0 gives Using the evolution equation (7.32) to replace∂ 0 ψ with λΘΨ + βψ in (7.37), we see, with the help of the identity E∂ β Θ = −∂ β EΘ (7.38) that follows from differentiating (7.11), that we can write (7.37) as where Y α and H are as defined above by (7.33) and (7.34), respectively. But the right hand side of this expression vanishes by virtue of the evolution equation (7.28) showing that∂ 0 I = 0 in [0, T ] × Ω 0 .
Acknowledgements: Part of the work for this article was completed at the Institut Henri Poincaré during the Mathematical Relativity trimester in 2015. I thank the Institute for its support and hospitality. This work was partially supported by the Australian Research Council grants FT1210045 and DP170100630.
Appendix A. Differential geometry formulas In this appendix, we collect together some useful formulas from differential geometry that will be used throughout this article. In the following, we let g = g µν dx µ dx ν denote a smooth Lorentzian metric on a four dimensional manifold M , and we use ∇ to denote the Levi-Civita connection of this metric. We also use the indexing conventions from Section 2.1. Given a local frame e j = e µ j ∂ µ (A.1) on M , we denote the dual frame by θ j = θ j µ dx µ (θ j µ ) := (e µ j ) −1 , (A.2) and use the notation γ ij = g(e i , e j ) = g µν e µ i e ν j and γ ij = g(θ i , θ j ) = g µν θ i µ θ j ν (A. 3) for the frame metric and its inverse, respectively.
A.2. Volume form. We use ν to denote the volume form of the metric g, which is given locally by where the components are computed using ν µαβγ = |g|ǫ µαβγ .
Here, ǫ ναβγ denotes the completely anti-symmetric symbol and we employ the standard notation for the negative of the determinant of the metric g µν . The volume form is also given locally in terms of the coframe by ν = |γ|θ 0 ∧ θ 1 ∧ θ 2 ∧ θ 3 , (A.10) where |γ| = − det(γ kl ).
A.7. Covariant derivatives and changes of metrics. Letting∇ denote the Levi-Civita connection of another Lorentian metricĝ =ĝ µν dx µ dx ν on M, the covariant derivatives with respect to the metrics g andĝ are related via the formulâ Appendix B. Maxwell's equations As in the introduction, let Ω 0 ⊂ M be a bounded, connected spacelike hypersurface with smooth boundary ∂Ω 0 , and Ω T be a timelike cylinder diffeomorphic to [0, T ] × Ω 0 . We use Γ T to denote the timelike boundary of Ω T , which is diffeomorphic to [0, T ] × ∂Ω 0 , and note that Γ 0 = ∂Ω 0 . We denote the outward unit conormal to Γ T by n = n ν dx ν , which we arbitrarily extend to all of M, and we let Ω T ∼ = {T } × Ω 0 and Γ T = ∂Ω T ∼ = {T } × ∂Ω 0 denote the "top" of the spacetime cylinder and its boundary, respectively. We further assume that τ = τ µ ∂ µ and ξ = ξ µ ∂ µ are timelike, future pointing C 1 vector fields on Ω T that satisfy τ (x), ξ(x) ∈ T x Γ T for all x ∈ Γ T , τ is normal to Ω 0 and Ω T , and τ is unit length.
Maxwell's equations on the world tube Ω T are given by where F = 1 2 F µν dx µ ∧ dx ν is the electromagnetic tensor. We recall that the stress-energy tensor T µν of the electromagnetic field is defined by T µν = 2F µα F ν α − 1 2 g µν F αβ F αβ . Integrating this expression over Ω T leads to the following well-known integral relation, which will be needed in the proof of Theorem 4.1.
Integrating this expression over Ω T , we find using the Divergence Theorem that where in deriving this we have used the fact that n µ ξ µ = 0 in Γ T .
In the proof of Theorem 4.1, we will also need the inequality from the following lemma, which is used in literature to show that the electromagnetic stress-energy tensor satisfies the Dominant Energy Condition. Before stating the lemma, we first denote the unit-normalized version of ξ µ by v µ = (−g(ξ, ξ)) −1/2 ξ µ and we let h µν = g µν + v µ v ν denote the induced positive definite metric on the subspace g-orthogonal to v µ . We also define a positive definite metric by m µν = g µν + 2v µ v ν . Proof. Starting from the standard decomposition, see [17,Ch. 13], of the electromagnetic tensor in terms of the electric and magnetic fields relative to v µ , which are defined by respectively, a straightforward calculation, see [17,Ch. 15], shows that F µν F µν = 2(B µ B µ − E µ E µ ). Using E µ v µ = B µ v µ = 0, (B.7) we can then write F µν F µν as F µν F µν = 2(|B| 2 m − |E| 2 m ). From this and (B.7), it is then not difficult to verify that |F | 2 m = m αβ m µν F αµ F βν = 2(|B| 2 m + |E| 2 m ). (B.8) Next, we recall that the energy density relative to v µ is given by, see [17,Ch. 15], the formula v µ T µν v ν = |E| 2 m + |B| 2 m . (B.9) Using this, we compute Furthermore, since v µ and τ µ are both future pointing and timelike by assumption, we have that v µ τ µ > 0 and τ µ τ µ = |τ | 2 h − (−v µ τ µ ) 2 < 0, from which it follows that (−τ λ v λ ) − |τ | h ≥ c > 0 in Ω T for some positive constant c. The proof now follows from this inequality and (B.11).