Peeling in generalized harmonic gauge

The concept of null infinity plays a very important role in general relativity (GR) and astrophysics [1–4]. From the point of view of both numerical and mathematical relativity, it is necessary in finding solutions to long-standing problems such as the weak cosmic censorship conjecture, global stability analysis of spacetimes and their (non) peeling properties. We will come back to the latter later in this introduction, as it is the main topic that this work aims to address. From an astrophysics point of view, null infinity is also a vital concept, since gravitational radiation is non-localizable and is hence only well-defined at null infinity. Waveforms computed in numerical relativity are however evaluated at large radius and extrapolated to infinity. There are various approaches to avoid this, and instead to include null infinity in the computational domain. For instance, the work of Hübner [5, 6] and Frauendiener [7] makes use of the conformal Einstein field equations (CEFE) introduced by Friedrich [8, 9], in turn building upon Penrose's idea of bringing null infinity to a finite coordinate distance by means of a conformal compactification [1]. Although the CEFE provide a geometric approach to the problem of the inclusion of null infinity, the standard methods of numerical relativity that have proven to work well for the strong field region of spacetimes of physical interest cannot be trivially lifted over. In particular, this approach has not yet been used for compact binary evolutions.

Chief alternative approaches to the problem are Cauchy-characteristic-matching [10], and the use of initial data on a hyperboloidal slice. Such slices are not Cauchy hypersurfaces, as their domain of dependence does not contain the entire spacetime, rather they are spacelike hypersurfaces which intersect null infinity. At the PDE level, this presents technical difficulties, the main issue being the appearance of formally singular terms. Numerical treatments of this in spherical symmetry are [11–14]. Further hyperboloidal formulations of GR (with successful numerics) can be found in [15–18].

It is thus desirable to find formulations that mollify the formally singular terms. A recent proposal is the dual frame approach [19]. This strategy consists of decoupling the coordinates from the tensor basis and carefully choosing each, thus allowing one to write the Einstein field equations (EFE) in generalized harmonic gauge (GHG) and then solve the resulting equations in hyperboloidal coordinates [20–23]. Important ingredients in this method are the coordinate lightspeed condition and the weak null condition [21, 24]. The former is the requirement that derivatives of the radial coordinate lightspeeds have a certain fall-off near null infinity, while the latter is expected to be a sufficient condition on the non-linearities of a quasilinear wave equation for establishing small data global existence. Although it has not been shown in full generality that the weak null condition implies small data global existence, recent work by Keir [25] proved that if a system of quasilinear wave equations satisfies the hierarchical weak null condition, then small data global existence is guaranteed. In [26], the conditions required on the nonlinearities were relaxed even further.

It has been shown in [20] that the dual foliation formulation of GR in GHG, together with hyperboloidal coordinates, can be used to avoid the worst type of formally singular terms. In harmonic gauge, however, even the simplest choices of initial data give rise to terms with decay near null infinity of the type O(R⁻¹ log R), with R a suitably defined radial coordinate. Such terms are problematic in numerical evolutions, which eventually see logarithms as an obstruction to rapid convergence. In [22], the authors used a toy model composed of wave equations with non-linearities of the same kind as those present in the EFE to show that these logarithmically divergent terms can be explicitly regularized by a non-linear change of variables. This toy model is called the good–bad–ugly model as it splits the evolved fields into three categories according to their fall-off near null infinity, and it is known to satisfy the weak null condition.

Exploiting conformal transformations and under the assumption of the smoothness of null infinity, it can be shown, under very general conditions, that the components of the Weyl tensor of spacetimes associated with an isolated source fall off with certain negative integer powers of r, an affine parameter along outgoing null geodesics [2]. This property, known as peeling, implies in particular that the far field of any source of gravitational waves behaves locally as a plane wave. The question of whether physically relevant spacetimes satisfy the peeling property was a subject of extensive debate. Whereas solutions of interest have been shown to allow for a smooth null infinity, namely using hyperboloidal initial data [27–29], and asymptotically flat initial data [30, 31], it is now clear that several physically motivated constructions violate peeling, and hence cannot have smooth null infinity [32, 33]. Rigorous evidence (from the PDE perspective) of logarithmic terms in the asymptotics of solutions to the EFE in harmonic gauge can be found in [34] and in the proof of the stability of the Minkowski spacetime (with polyhomogeneous initial data and developments) of [35].

In [36], a heuristic method to find asymptotic solutions to the good–bad–ugly system was laid out. The non-linearities in this model are known to mimic those present in the EFE. Formal expansions are derived in which terms proportional to the logarithm of the radial coordinate appear at every order in the bad field, from the second order onward in the ugly field but never in the good field. The model was then generalized to wave operators built from an asymptotically flat metric and it was shown that it admits polyhomogeneous asymptotic solutions. Finally the authors define stratified null forms, a generalization of standard null forms, which capture the behavior of different types of fields, and demonstrate that the addition of such terms to the original system bears no qualitative influence on the type of asymptotic solutions found.

In this work we build upon [36] by pursuing the same strategy to find polyhomogeneous expansions of the EFE near null infinity. The question addressed is whether or not, even within a favorable class of initial data, the peeling property is manifest in harmonic or GHG. In section 2 we review the formalism we use, which is very similar to that used in [36]. We furthermore write the derivatives of metric components in terms of components of the transition tensor between two covariant derivatives. In section 3 we state the theorem shown in [36] and generalize its proof to include GR in GHG. In section 4 we derive the EFE as a system of 10 coupled non-linear wave equations for our variables of choice. Section 5 is dedicated to showing that the EFE, in Cartesian harmonic gauge and with a particular constraint addition, constitutes a good–bad–ugly system. Building on previous results, we then integrate the equations for the first few orders to find the functional form of the leading terms in the polyhomogeneous expansion. Informally, we integrate the equations to find the powers of log R that can be present in each term with general initial data. In section 6 we show that there are choices of initial data which prevent the Weyl scalars from peeling by introducing a log R/R³ term in Ψ₂. In section 7 we explore the interplay between gauge, constraint addition and decay. Inspired by the general strategy of [37], we show that there is a particular choice of gauge source function which recovers the peeling property by preventing any powers log R from appearing in the first orders of the expansion. In section 8, we show that by applying a change of coordinates to this good gauge choice, the obstruction to peeling that we saw earlier reappears. We conclude in section 9.

Representation of the metric: abstract tensor indices are represented by Latin letters and coordinate indices are represented by Greek letters from the beginning of the alphabet. We assume a Lorentzian metric g_ab with Levi-Civita connection ∇; indices are raised and lowered with the spacetime metric g_ab and it is inverse. We introduce an asymptotically Cartesian coordinate system ${X}^{\underline{\alpha }}=(T,{X}^{\underline{i}})$, with the respective vector and covector bases ${\partial }_{\underline{\alpha }}$ and $\mathrm{d}{X}^{\underline{\alpha }}$. The sense in which ${X}^{\underline{\alpha }}=(T,{X}^{\underline{i}})$ is asymptotically Cartesian will be clarified later. It is useful to define the vector field ${T}^{a}{:=}{\partial }_{T}^{a}$, with the associated directional derivative denoted ∇_T = T^a ∇_a . Directional derivatives along other vector fields will be defined similarly. The flat covariant derivative associated to ${X}^{\underline{\alpha }}$ is denoted $\stackrel{\circ}{\nabla }$ and has the defining property that ${\stackrel{\circ}{\nabla }}_{a}{\partial }_{\underline{\alpha }}^{b}=0$. We also introduce shell coordinates ${X}^{{\underline{\alpha }}^{\prime }}=({T}^{\prime },{X}^{{\underline{i}}^{\prime }})=(T,R,{\theta }^{A})$, with radial coordinate given by ${R}^{2}={({X}^{\underline{1}})}^{2}+{({X}^{\underline{2}})}^{2}+{({X}^{\underline{3}})}^{2}$ and the respective vector and covector bases ${\partial }_{{\underline{\alpha }}^{\prime }}$ and $\mathrm{d}{X}^{{\underline{\alpha }}^{\prime }}$. The flat covariant derivative associated to the shell coordinates has the defining property that ${\stackrel{\bullet }{\nabla }}_{b}{\partial }_{{\underline{\alpha }}^{\prime }}^{a}=0$. The flat covariant derivatives $\stackrel{\circ}{\nabla }$ and $\stackrel{\bullet }{\nabla }$ reduce to partial derivatives when working in the associated coordinate bases. These definitions are based on the notation used in xAct [38] to represent partial derivatives using abstract index notation. Given the two flat covariant derivatives $\stackrel{\circ}{\nabla }$ and $\stackrel{\bullet }{\nabla }$, we define the Christoffel transition tensor

$$\begin{equation}{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{a}{}^{\!b}{}_{\!\!c}{v}^{c}={\stackrel{\circ}{\nabla }}_{a}{v}^{b}-{\stackrel{\bullet }{\nabla }}_{a}{v}^{b},\end{equation} \tag{ 1 }$$

v^a being an arbitrary vector field. Similarly,

$$\begin{equation}{\Gamma}{[\nabla ,\stackrel{\bullet }{\nabla }]}_{a}{}^{\!b}{}_{\!\!c}{v}^{c}={\nabla }_{a}{v}^{b}-{\stackrel{\bullet }{\nabla }}_{a}{v}^{b},\end{equation} \tag{ 2 }$$

which will be represented with the shorthand ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}{:=}{\Gamma}{[\nabla ,\stackrel{\bullet }{\nabla }]}_{a}{}^{\!b}{}_{\!\!c}$. It follows that

$$\begin{equation}{\Gamma}{[\nabla ,\stackrel{\circ}{\nabla }]}_{b}{}^{\!a}{}_{\!\!c}={\Gamma}{[\nabla ,\stackrel{\bullet }{\nabla }]}_{b}{}^{\!a}{}_{\!\!c}-{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{b}{}^{\!a}{}_{\!\!c}.\end{equation} \tag{ 3 }$$

For clarity, we show these coordinates and their covariant derivatives in the following table:

Coordinates	Cov. der.
${X}^{\underline{\alpha }}=(T,{X}^{\underline{i}})$	$\stackrel{\circ}{\nabla }$
${X}^{{\underline{\alpha }}^{\prime }}=({T}^{\prime },{X}^{{\underline{i}}^{\prime }})$	$\stackrel{\bullet }{\nabla }$

${\underline{\psi }}^{a}$

$$\begin{align}\hfill {\psi }^{a}& ={\partial }_{T}^{a}+{\mathcal{C}}_{+}^{R}{\partial }_{R}^{a},\hfill \\ \hfill {\underline{\psi }}^{a}& ={\partial }_{T}^{a}+{\mathcal{C}}_{-}^{R}{\partial }_{R}^{a},\hfill \end{align} \tag{ 4 }$$

where ${\mathcal{C}}_{+}^{R}$ and ${\mathcal{C}}_{-}^{R}$ are radial coordinate light-speeds, the values for which are determined by the condition that ψ^a and ${\underline{\psi }}^{a}$ are null with respect to g_ab . It is useful to define the following covectors,

$$\begin{equation}{\sigma }_{a}={e}^{-\varphi }{\psi }_{a},\qquad {\underline{\sigma }}_{a}={e}^{-\varphi }{\underline{\psi }}_{a},\end{equation} \tag{ 5 }$$

with φ specified by the condition ${\sigma }_{a}{\partial }_{R}^{a}=-{\underline{\sigma }}_{a}{\partial }_{R}^{a}=1$. The covectors may be decomposed in the following way:

$$\begin{align}\hfill {\sigma }_{a}& =-{\mathcal{C}}_{+}^{R}{\nabla }_{a}T+{\nabla }_{a}R+{\mathcal{C}}_{A}^{+}{\nabla }_{a}{\theta }^{A},\hfill \\ \hfill {\underline{\sigma }}_{a}& ={\mathcal{C}}_{-}^{R}{\nabla }_{a}T-{\nabla }_{a}R+{\mathcal{C}}_{A}^{-}{\nabla }_{a}{\theta }^{A}.\hfill \end{align} \tag{ 6 }$$

It is straightforward to show that the null vectors satisfy the relations

$$\begin{align}\hfill & {\sigma }_{a}{\psi }^{a}={\underline{\sigma }}_{a}{\underline{\psi }}^{a}=0,\hfill \\ \hfill & {\sigma }_{a}{\underline{\psi }}^{a}={\underline{\sigma }}_{a}{\psi }^{a}=-\tau ,\hfill \end{align} \tag{ 7 }$$

with $\tau {:=}{\mathcal{C}}_{+}^{R}-{\mathcal{C}}_{-}^{R}$.

We perform the following decomposition of the inverse spacetime metric,

$$\begin{equation}{g}^{ab}=-2{\tau }^{-1}{e}^{-\varphi }\,{\psi }^{(a}{\underline{\psi }}^{b)}+{\not{\hskip-3pt {g}}}^{\,ab},\end{equation} \tag{ 8 }$$

with the normalization of first term chosen so that ${\not{\hskip-3pt g}}^{\,ab}{\sigma }_{b}={\not{\hskip-3pt g}}^{\,ab}{\underline{\sigma }}_{b}=0$. We emphasize that while ${{\not{\hskip-3pt g}}_{\,a}}^{b}$ serves as a projection operator orthogonal to σ_a and ${\underline{\sigma }}_{b}$, ${\not{\hskip-3pt g}}^{\,ab}$ is not the inverse induced metric on level sets of T and R, as ∇_a T or ∇_a R are not in the kernel of ${\not{\hskip-3pt g}}^{\,ab}$. We obtain the following decomposition for the metric,

$$\begin{equation}{g}_{ab}=-2{\tau }^{-1}{e}^{\varphi }\,{\sigma }_{(a}{\underline{\sigma }}_{b)}+{\not{\hskip-3pt g}}_{\,ab}.\end{equation} \tag{ 9 }$$

Let $\mathcal{S}$ be a Cauchy surface defined by a constant value of the coordinate T. It will be convenient to make a conformal rescaling of this induced metric on $\mathcal{S}$. We define

$$\begin{equation}{q}_{ab}={e}^{-{\epsilon}}{R}^{-2}{\not{\hskip-3pt g}}_{\,ab},\qquad {({q}^{-1})}^{ab}={e}^{{\epsilon}}{R}^{2}{\not{\hskip-3pt g}}^{\,ab}\end{equation} \tag{ 10 }$$

with

$$\begin{equation}{\epsilon}=(\mathrm{ln}\vert \not{\hskip-3pt g}\vert -\mathrm{ln}\vert \!\stackrel{\bullet }{\not{\hskip-3pt { g}}}\vert )/2\end{equation} \tag{ 11 }$$

where $\vert \!\stackrel{\bullet }{\not{\hskip-3pt { g}}}\vert$ denotes the determinant of the metric of ${\mathbb{S}}^{2}$ of radius R embedded in Minkowski spacetime in Shell coordinates. For future reference we also define $\stackrel{\bullet }{{\epsilon}}\equiv \frac{1}{2}\mathrm{ln}\,\vert \!\stackrel{\bullet }{\not{\hskip-3pt { g}}}\vert$. Notice that, as a consequence of the latter, the determinant of q is that of the metric of the unit ${\mathbb{S}}^{2}$ in Shell coordinates, so it is a fixed function of the coordinates θ^A . Although this construction is general for any coordinatization of ${\mathbb{S}}^{2}$, for conciseness, in the rest of the paper, we use standard spherical polar coordinates θ¹ = θ and θ² = ϕ so that the line element on ${\mathbb{S}}^{2}$ reads dΩ² = dθ² + sin² θdϕ². We have used the notation ${({q}^{-1})}^{ab}$ to emphasize that the metric $\not{\hskip-3pt g}$ and not q is used to raise and lower indices on tensors on $\mathcal{S}$. For instance ${q}^{ab}={\not{\hskip-3pt g}}^{\,ac}{\not{\hskip-3pt g}}^{\,bd}{q}_{cd}\ne {({q}^{-1})}^{ab}$. Finally, it is worth introducing a special parameterization of the angular part of ${({q}^{-1})}^{ab}$ in Shell coordinates, in order to find wave equations for each of its two independent metric functions. These independent components encode the two degrees of freedom of gravitational waves usually denoted in the linear theory as h₊ and h_×.

$$\begin{equation}{({q}^{-1})}^{AB}=\left[\begin{matrix}\hfill {e}^{-{h}_{+}}\mathrm{cosh}\,{h}_{\times }\hfill & \hfill \frac{\mathrm{sinh}\,{h}_{\times }}{\mathrm{sin}\,\theta }\hfill \\ \hfill \frac{\mathrm{sinh}\,{h}_{\times }}{\mathrm{sin}\,\theta }\hfill & \hfill \frac{{e}^{{h}_{+}}\mathrm{cosh}\,{h}_{\times }}{\mathrm{sin}\,{\theta }^{2}}\hfill \end{matrix}\right].\end{equation} \tag{ 12 }$$

Thus the ten independent variables we use to represent the metric are,

$$\begin{equation}\varphi ,\qquad {\mathcal{C}}_{\pm }^{R},\qquad {\mathcal{C}}_{A}^{\pm },\qquad {\epsilon},\qquad {h}_{+},\qquad {h}_{\times }.\end{equation} \tag{ 13 }$$

Tensors ${T}_{cd}^{\prime }$ projected with ${{\not{\hskip-3pt g}}_{\,a}}^{b}$ are ${\not{\hskip-3pt T}}_{ab}\equiv {{\not{\hskip-3pt g}}_{\,a}}^{c}{{\not{\hskip-3pt g}}_{\,b}}^{d}{T}_{cd}^{\prime }$. The projected covariant derivative is denoted $\not{\hskip-3pt D}$, and is given explicitly below for vectors satisfying ${v}^{a}={{\not{\hskip-3pt g}}_{\,b}}^{a}{v}^{b}$,

$$\begin{equation}{\not{\hskip-3pt D}}_{b}{v}^{a}{:=}{{\not{\hskip-3pt g}}_{\,c}}^{a}{{\not{\hskip-3pt g}}_{\,b}}^{d}{\nabla }_{d}{v}^{c}.\end{equation} \tag{ 14 }$$

The generalization for higher rank tensors is straightforward. We define the projected covariant derivative $\stackrel{\circ}{\not{\hskip-3pt D}}$ similarly, with ∇_d replaced with ${\stackrel{\circ}{\nabla }}_{d}$.

A brief remark on terminology—since we deal with polyhomogeneous expansions, which include terms of the form R⁻ⁿ (log R)^m , it is appropriate to clarify what is meant by the term 'order' in this article. Throughout, 'order n' will refer to terms proportional to R⁻ⁿ .

Decomposition of the connection: using the Koszul formula,

$$\begin{equation}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}=\frac{1}{2}{\not{\hskip-3pt g}}^{\,bd}({\stackrel{\bullet }{\nabla }}_{a}{\not{\hskip-3pt g}}_{\,db}+{\stackrel{\bullet }{\nabla }}_{b}{\not{\hskip-3pt g}}_{\,ad}-{\stackrel{\bullet }{\nabla }}_{d}{\not{\hskip-3pt g}}_{\,ab}),\end{equation} \tag{ 15 }$$

we can find expressions for each component of ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}$. Here we present a rearranged version of these relations, writing derivatives of metric functions in terms of components of the connection and not the other way around. The reason for us to do this is that writing the EFE in terms of ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}$ makes them significantly more readable. We find that,

$$\begin{align}\hfill & {\nabla }_{a}{\mathcal{C}}_{+}^{R}=-{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\sigma }{}_{\!\!\psi },\quad {\nabla }_{a}{\mathcal{C}}_{-}^{R}={\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\underline{\sigma }}{}_{\!\!\underline{\psi }},\hfill \\ \hfill & {\nabla }_{a}{\mathcal{C}}_{A}^{+}=2{\not{\hskip-3pt g}}_{\,bA}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!(b\sigma )}-\frac{{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\sigma }{}_{\!\!\psi },\hfill \\ \hfill & {\nabla }_{a}{\mathcal{C}}_{A}^{-}=2{\not{\hskip-3pt g}}_{\,bA}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!(b\underline{\sigma })}-\frac{{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\underline{\sigma }}{}_{\!\!\underline{\psi }},\hfill \\ \hfill & {\nabla }_{a}\varphi ={\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}({\delta }_{b}^{c}-{{\not{\hskip-3pt g}}_{\,b}}^{c}),{\nabla }_{a}{\not{\hskip-3pt g}}^{\,AB}=-2{{\not{\hskip-3pt g}}_{\,b}}^{A}{{\not{\hskip-3pt g}}_{\,c}}^{B}{\stackrel{\bullet }{{\Gamma}}}^{(bc)}{}_{\!\!a},\hfill \\ \hfill & {\nabla }_{a}({\epsilon}+\stackrel{\bullet }{{\epsilon}})={{\not{\hskip-3pt g}}_{\,c}}^{b}{\stackrel{\bullet }{{\Gamma}}}_{b}{}^{\!c}{}_{\!\!a}.\hfill \end{align} \tag{ 16 }$$

Note that these covariant derivatives are interchangeable with $\stackrel{\bullet }{\nabla }$ and $\stackrel{\circ}{\nabla }$ even when they act upon ${\mathcal{C}}_{A}^{\pm }$, as A is not a tensorial index, but a label to designate one of the angular coordinates, θ¹ = θ or θ² = ϕ. When writing the asymptotic system near null infinity, we will have to use the version of these relations that gives the components of the transition tensor ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}$ in terms of first derivatives of metric functions. In order not to overload the reader with long expressions, and because they contain basically the same information as (16), we do not present those relations here.

D'Alembert operators: throughout this work we will use three different wave operators, each associated to a different covariant derivative. Thus we define □ as,

$$\begin{equation}\square \phi ={g}^{ab}{\nabla }_{a}{\nabla }_{b}\phi ,\end{equation} \tag{ 17 }$$

and we define $\stackrel{\circ}{\square }$ and $\stackrel{\bullet }{\square }$ analogously, associated with $\stackrel{\circ}{\nabla }$ and $\stackrel{\bullet }{\nabla }$, respectively. Note that when $\stackrel{\circ}{\square }$ acts on a scalar function, one can straightforwardly change to $\stackrel{\bullet }{\square }$ through,

$$\begin{equation}\stackrel{\circ}{\square }\phi =\stackrel{\bullet }{\square }\phi -{g}^{bc}{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{b}{}^{\!a}{}_{\!\!c}{\nabla }_{a}\phi .\end{equation} \tag{ 18 }$$

Moreover, in [36] it was shown that the Cartesian wave operator can be expanded as,

$$\begin{align}\hfill {\mathcal{C}}_{+}^{R}\stackrel{\circ}{\square }\phi & =-2{e}^{-\varphi }{\nabla }_{\psi }{\nabla }_{T}\phi +{\nabla }_{T}\phi ({\not{\hskip-3pt g}}^{\,ab}{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{a}{}^{\!\sigma }{}_{\!\!b}+{X}_{T})\hfill \\ \hfill & \quad +{\nabla }_{\psi }\phi {X}_{\psi }-\frac{2{e}^{-\varphi }{\mathcal{C}}_{-}^{R}}{\tau }{\nabla }_{\psi }^{2}\phi +{\mathcal{C}}_{+}^{R}{/\hskip-3pt {\Delta}}\phi ,\hfill \end{align} \tag{ 19 }$$

where X_T and X_ψ are,

$$\begin{align}\hfill \tau {X}_{T}& {:=}{\mathcal{C}}_{A}{/{\hskip-3pt \!D}}^{A}{\mathcal{C}}_{+}^{R}-\tau {/{\hskip-3pt D}}^{A}{\mathcal{C}}_{A}^{+}-\frac{2{e}^{-\varphi }{\mathcal{C}}_{-}^{R}}{{\mathcal{C}}_{+}^{R}}{\nabla }_{\psi }{\mathcal{C}}_{+}^{R},\hfill \\ \hfill \tau {X}_{\psi }& {:=}\frac{{\mathcal{C}}_{A}}{\tau }\left({\mathcal{C}}_{-}^{R}{{\hskip-3pt D}}^{A}{\mathcal{C}}_{+}^{R}-{\mathcal{C}}_{+}^{R}{/{\hskip-3pt D}}^{A}{\mathcal{C}}_{-}^{R}\right)-{\mathcal{C}}_{-}^{R}{{\hskip-3pt D}}^{A}{\mathcal{C}}_{A}^{+}-{\mathcal{C}}_{+}^{R}{{/\hskip-3pt \!D}}^{A}{\mathcal{C}}_{A}^{-}+{\mathcal{C}}_{-}^{R}{{\hskip-3pt g}}^{\,ab}{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{a}{}^{\!\sigma }{}_{\!\!b}\hfill \\ \hfill & \quad +{\mathcal{C}}_{+}^{R}{{/\hskip-3pt g}}^{\,ab}{\Gamma}{[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]}_{a}{}^{\!\underline{\sigma }}{}_{\!\!b}+\frac{2{e}^{-\varphi }{\mathcal{C}}_{-}^{R}}{{\mathcal{C}}_{+}^{R}}{\nabla }_{\psi }{\mathcal{C}}_{+}^{R},\hfill \end{align} \tag{ 20 }$$

and ${\mathcal{C}}_{A}{:=}{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}$. Most importantly, it was shown that, to leading order, this operator behaves as,

$$\begin{equation}\stackrel{\circ}{\square }\phi \simeq -2{\nabla }_{\psi }{\nabla }_{T}\phi -\frac{2}{R}{\nabla }_{T}\phi .\end{equation} \tag{ 21 }$$

Stratified null forms are defined as terms that involve products of up to one derivative of the evolved fields and fall-off faster than R⁻² close to null infinity. Henceforth ${\mathcal{N}}_{\phi }$, where ϕ is the evolved field that ${\mathcal{N}}_{\phi }$ is associated to, will denominate an arbitrary linear combination of stratified null forms.

The good–bad–ugly system: we introduce the following model,

$$\begin{align}\hfill & \stackrel{\circ}{\square }g={\mathcal{N}}_{g},\hfill \\ \hfill & \stackrel{\circ}{\square }b={({\nabla }_{T}g)}^{2}+{\mathcal{N}}_{b},\hfill \\ \hfill & \stackrel{\circ}{\square }u=\frac{2}{R}{\nabla }_{T}u+{\mathcal{N}}_{u},\hfill \end{align} \tag{ 22 }$$

where g, b and u stand for good, bad and ugly fields, respectively. The leading order of the particular case where the metric g^ab is the Minkowski metric and ${\mathcal{N}}_{\phi }=0$ was studied in detail in [22]. The main reason for the interest in this particular model lies in the fact that its non-linearities are known to mimic those present in the EFE. Previous work [36] has shown that this system admits polyhomogeneous expansions near null infinity. In this work we want to build directly upon that result, by applying it to GR and finding the functional form of terms beyond first order. Then we aim to shed light on a very special interplay between the choice of gauge and constraint addition by showing that peeling is violated for certain choices and satisfied in others. Additionally, this will allow us to know where terms with factors of log R may appear and thus builds towards a full-blown regularization of the EFE in GHG in hyperboloidal slices within the dual foliation framework. Finally, we expect that this method is also applicable to the Maxwell and Yang–Mills equations, as these can also be written as non-linear wave equations in Lorenz gauge.

3.1. Assumptions

To proceed, it is appropriate to outline assumptions on the evolved fields in (22) that allow us to formally equate terms of the same order and identify a hierarchy of equations which are satisfied by the fields g, b and u order-by-order close to null infinity. We also discuss corresponding assumptions for the metric functions.

Evolved fields: consider a null tetrad $\left\{\psi ,\underline{\psi },{X}_{1},{X}_{2}\right\}$, with X₁ and X₂ orthogonal to ψ^a and ${\underline{\psi }}^{a}$ and normalized so that ${g}_{ab}{X}_{A}^{a}{\psi }^{b}={g}_{ab}{X}_{A}^{a}{\underline{\psi }}^{b}=0$ and ${g}_{ab}{X}_{A}^{a}{X}_{B}^{b}={\delta }_{AB}$, with A = 1, 2. Now let ω_g,b represent a good or a bad field or any of its first derivatives. Following the insight of [25] we assume the following behavior at null infinity,

$$\begin{equation}{\omega }_{g,b}={o}^{+}({R}^{-n}){\Rightarrow}\begin{cases}{\nabla }_{\psi }{\omega }_{g,b}={o}^{+}({R}^{-n-1})\quad \hfill \\ {\nabla }_{\underline{\psi }}{\omega }_{g,b}={o}^{+}({R}^{-n})\quad \hfill \\ {\nabla }_{{X}_{A}}{\omega }_{g,b}={o}^{+}({R}^{-n-1})\hfill \end{cases} .\end{equation} \tag{ 23 }$$

The notation f = o⁺(h) refers to the condition

$$\begin{equation}\exists {\epsilon} > 0:\underset{R\to \infty }{\mathrm{lim}}\frac{f}{h{R}^{-{\epsilon}}}=0,\end{equation} \tag{ 24 }$$

which can be informally stated as the condition that f falls-off faster than h¹⁺ as R goes to infinity, which is a faster falloff than f = o(h); more precisely, o⁺(h) = o(hR⁻ ). As discussed in [36] this condition will be needed to ensure that the error terms remain small when integrated. We employ the condition o⁺(R⁻ⁿ ) instead of O(R⁻ⁿ⁻¹) since it was shown previously (for instance [22]) that the system (22) admits asymptotic solutions proportional to R⁻¹ log R, and the O notation naively excludes such solutions. Notice that for fields satisfying a system of the form (22), certain derivatives improve the fall off of the argument, but others derivatives do not; this motivates the terminology employed in the previous literature, namely that directional derivatives corresponding to ψ^a and X_A (those tangent to outgoing null-cones) are termed good derivatives, and directional derivatives corresponding to ${\underline{\psi }}^{a}$ (transverse to outgoing null-cones) are termed bad derivatives.

Now let ω_u represent the field u or its first derivative. As discussed in [22], if $\stackrel{\circ}{\square }$ is constructed from the Minkowski metric, the derivatives of the u fields have different asymptotic behavior; we assume the following,

$$\begin{equation}{\omega }_{u}={o}^{+}({R}^{-n}){\Rightarrow}\begin{cases}{\nabla }_{\psi }{\omega }_{u}={o}^{+}({R}^{-n-1})\quad \hfill \\ {\nabla }_{\underline{\psi }}{\omega }_{u}={o}^{+}({R}^{-n-1})\quad \hfill \\ {\nabla }_{{X}_{A}}{\omega }_{u}={o}^{+}({R}^{-n-1})\quad \hfill \end{cases} .\end{equation} \tag{ 25 }$$

We seek solutions which decay near null infinity, so we restrict ourselves to initial data with that property, i.e.,

$$\begin{equation}g={o}^{+}(1),\qquad b={o}^{+}(1),\qquad u={o}^{+}(1).\end{equation} \tag{ 26 }$$

To allow for a nonzero ADM mass and linear momentum, the initial data is chosen to decay at spacelike infinity (as opposed to null infinity) in the following manner,

$$\begin{align}\hfill & \phi {\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{\phi ,n}}{{R}^{n}},\hfill \\ \hfill & {\nabla }_{T}\phi {\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2}),\hfill \end{align} \tag{ 27 }$$

where m_ϕ,n are scalar functions that are independent of T and R. Though this is not the most general choice which permits nontrivial ADM mass and linear momentum, it is general enough for a very large class of spacetimes of interest.

Metric functions: we now turn to the metric functions. We require that the metric functions are written as,

$$\begin{align}\hfill & \varphi ={\gamma }_{1},\qquad {\mathcal{C}}_{\pm }^{R}=\pm 1+{\gamma }_{2}^{\pm },\qquad {\mathcal{C}}_{A}^{\pm }=R{\gamma }_{3}^{\pm },\hfill \\ \hfill & {\epsilon}={\gamma }_{4},\qquad {h}_{+}={\gamma }_{5},\qquad {h}_{\times }={\gamma }_{6},\hfill \end{align} \tag{ 28 }$$

where the γ = γ(g, b, u) are assumed to be analytic in a neighborhood of null infinity ${\mathcal{I}}^{+}$. Since we are interested in metrics that asymptote to the Minkowski metric near ${\mathcal{I}}^{+}$, we require that

$$\begin{equation}\gamma (g,b,u){\vert }_{{\mathcal{I}}^{+}}=0,\end{equation} \tag{ 29 }$$

or that the γ functions vanish as one approaches asymptotic infinity. In fact, we will see that upon an appropriate choice of gauge and constraint addition, in GR the γ functions turn out to be themselves either good, bad or ugly. The metric variables in this formulation were chosen for this purpose.

3.2. Earlier work

In [36] the following theorem was shown,

Theorem 1. The good–bad–ugly system defined as (22) where g_ab is an asymptotically flat metric, admits a polyhomogeneous expansion near null infinity of the type,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{n,k}({\psi }^{\ast })}{{R}^{n}}\hfill \\ \hfill & b=\frac{{B}_{1,0}({\psi }^{\ast })+\mathrm{log}\,R{B}_{1,1}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{b}}\frac{{(\mathrm{log}\,R)}^{k}{B}_{n,k}({\psi }^{\ast })}{{R}^{n}}\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{n,k}({\psi }^{\ast })}{{R}^{n}},\hfill \end{align} \tag{ 30 }$$

where ψ* means that this scalar function does not vary along integral curves of ψ, Φ_n,k with Φ ∈ {G, B, U} are coefficients where n is their associated power of R⁻¹ and k their associated power of log R, and with initial data on $\mathcal{S}$ of the type,

$$\begin{equation}\begin{cases}g{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{g,n}}{{R}^{n}}\quad \hfill \\ b{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{b,n}}{{R}^{n}}\quad \hfill \\ u{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{u,n}}{{R}^{n}}\quad \hfill \end{cases} ,\quad \begin{cases}{\nabla }_{T}g{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \\ {\nabla }_{T}b{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \\ {\nabla }_{T}u{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \end{cases} ,\end{equation} \tag{ 31 }$$

where m_ϕ,n are scalar functions that are independent of T and R. This is valid outside a compact ball centered at R = 0.

The assumptions on initial data could be relaxed if we wished only to build formal solutions down to a finite order in R⁻¹.

Throughout the article, all functions denoted by m_ϕ,n for any field ϕ and any integer n will be assumed to be independent of T and R. We will show that in GR the metric functions (13) can be separated into good, bad and ugly fields, as defined in (22). Moreover, we will see that there is an interplay between gauge choice and constraint addition that interferes with the asymptotic system and gives rise to different combinations of these fields. However there are three generalizations that we have to do to the above-mentioned theorem in order to include GR with the gauge choices of interest. First, our good–bad–ugly system should be allowed to include arbitrary numbers of good, bad and ugly fields. In particular, in GR, there will be 10 independent metric functions, each with its own nonlinear wave equation. Second, the leading term on the rhs of the bad equation must be allowed to be $-1/2{({\nabla }_{T}{g}_{1})}^{2}-1/2{({\nabla }_{T}{g}_{2})}^{2}$, where g₁ and g₂ are good fields. And third, GR forces us to extend our conception of an ugly equation to include a slightly larger class of equations whose asymptotic expansions behave as uglies to leading order, but differ in the decay rate of terms beyond first order.

3.3. Generalization of theorem 1

Let there be any number of evolved fields, each satisfying one of the following three kinds of nonlinear wave equation,

$$\begin{align}\hfill & \stackrel{\circ}{\square }g={\mathcal{N}}_{g},\hfill \\ \hfill & \stackrel{\circ}{\square }b=\frac{\beta }{{R}^{2}}+{\mathcal{N}}_{b},\hfill \\ \hfill & \stackrel{\circ}{\square }u=\frac{2p}{R}{\nabla }_{T}u+{\mathcal{N}}_{u},\hfill \end{align} \tag{ 32 }$$

where $-2\beta ={R}^{2}{({\nabla }_{T}{g}_{1})}^{2}+{R}^{2}{({\nabla }_{T}{g}_{2})}^{2}$, g₁ and g₂ satisfy an equation of the first type, and p is a natural number. To prove a more general version of theorem 1 we employ a method almost exactly the same that considered in [36]. For this reason we leave out the common details and focus exclusively on the differences.

Motivation for induction hypothesis: asymptotic flatness, together with the fact that the γ functions (see (28)) are analytic functions of the evolved fields at null infinity, allows us to Taylor expand them around g = b = u = 0, because the fields are assumed to have decay near null infinity and find that γ = o⁺(1) and,

$$\begin{equation}{\omega }_{\gamma }={o}^{+}({R}^{-n}){\Rightarrow}\begin{cases}{\nabla }_{\psi }{\omega }_{\gamma }={o}^{+}({R}^{-n-1})\quad \hfill \\ {\nabla }_{{X}_{A}}{\omega }_{\gamma }={o}^{+}({R}^{-n-1})\quad \hfill \end{cases} ,\end{equation} \tag{ 33 }$$

where ω_γ is any γ function or any first derivative of it. Let us rescale g and b as,

$$\begin{equation}{\mathcal{G}}_{1}=Rg,\qquad {\mathcal{B}}_{1}=Rb,\end{equation} \tag{ 34 }$$

and focus on expansions of u which satisfy,

$$\begin{equation}u=\frac{{m}_{u,1}}{R}+\frac{{\mathcal{U}}_{2}}{{R}^{2}},\end{equation} \tag{ 35 }$$

with ${\mathcal{U}}_{2}={o}^{+}(R)$. Plugging (34) and (35) in (32) we get,

$$\begin{align}\hfill & {\nabla }_{\psi }{\nabla }_{T}{\mathcal{G}}_{1}\simeq 0,\hfill \\ \hfill & {\nabla }_{\psi }{\nabla }_{T}{\mathcal{B}}_{1}\simeq -\frac{1}{R}\beta ,\hfill \\ \hfill & {\nabla }_{\psi }({R}^{p-1}{\nabla }_{T}{\mathcal{U}}_{2})\simeq {R}^{p-2}{{\Omega}}_{1}^{u},\hfill \end{align} \tag{ 36 }$$

where ${{\Omega}}_{1}^{u}={o}^{+}(R)$ cannot contain ${\mathcal{U}}_{2}$ or any derivatives thereof. For the first two equations we have,

$$\begin{align}\hfill & {\mathcal{G}}_{1}\simeq {G}_{1,0}({\psi }^{\ast })\hfill \\ \hfill & {\mathcal{B}}_{1}\simeq {B}_{1,0}({\psi }^{\ast })+{B}_{1,1}({\psi }^{\ast })\mathrm{log}\,R.\hfill \end{align} \tag{ 37 }$$

The third equation requires closer attention. To leading order, ${{\Omega}}_{1}^{u}$ is allowed to contain logs, but it is not allowed to contain R⁻¹, so we can write,

$$\begin{equation}{{\Omega}}_{1}^{u}=\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}{{\Omega}}_{1,i}^{u}({\psi }^{\ast }).\end{equation} \tag{ 38 }$$

In order to integrate the last equation in (36) along ψ, we have to integrate (log R)ⁱ R^p−2, which gives one of two results depending on p,

$$\begin{equation}\int \frac{{(\mathrm{log}\,R)}^{i}}{{R}^{-p+2}}\mathrm{d}R=\begin{cases}\sum\limits _{j=0}^{i}-\frac{{(\mathrm{log}\,R)}^{j}}{{(-p+1)}^{i-j+1}{R}^{-p+1}}\frac{i!}{j!},\,p\ne 1\quad \hfill \\ \frac{{(\mathrm{log}\,R)}^{i+1}}{i+1},\,p=1\quad \hfill \end{cases} .\end{equation} \tag{ 39 }$$

The main difference between these two cases is that if p = 1, the maximum power of log R goes up by one, whereas if p ≠ 1 it does not. Since we are allowing our null forms to have any finite power of log R (38), this does not make a big difference for the general case. However, more structure on these null forms will allow us to find the power of log R at each order. More on this point will be said later. For now, for either case we can write,

$$\begin{equation}{\mathcal{U}}_{2}\simeq \sum\limits _{i=0}^{{N}_{2}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}{U}_{2,i}({\psi }^{\ast })+\frac{1}{{R}^{p-1}}\int {\dot{u}}_{2}({\psi }^{\ast })\mathrm{d}T,\end{equation} \tag{ 40 }$$

where the second term can be incorporated in the first for p = 1, while for p > 1 we choose solutions with ${\dot{u}}_{2}({\psi }^{\ast })=0$, which basically amounts to pushing the second term into a higher order term in the expansion. This seems to suggest that g, b and u are polyhomogeneous functions where each term can have up to ${N}_{n}^{\phi }$ powers of log R, where ${N}_{n}^{\phi }$ is a finite number which depends on which field ϕ it refers to and on the power of R⁻¹ in the expansion, n. Formally, we therefore conjecture,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{n,k}({\psi }^{\ast })}{{R}^{n}},\hfill \\ \hfill & b=\frac{{B}_{1,0}({\psi }^{\ast })+\mathrm{log}\,R{B}_{1,1}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{b}}\frac{{(\mathrm{log}\,R)}^{k}{B}_{n,k}({\psi }^{\ast })}{{R}^{n}},\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{n=2}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{n,k}({\psi }^{\ast })}{{R}^{n}}.\hfill \end{align} \tag{ 41 }$$

We proceed by induction. We already know that to first order in g and u, log R terms are not allowed, and the conjecture (41) incorporates this property by construction. Truncating at n = 1, we have seen

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\frac{{\mathcal{G}}_{2}}{{R}^{2}},\hfill \\ \hfill & b=\frac{{B}_{1,0}({\psi }^{\ast })+{B}_{1,1}({\psi }^{\ast })\mathrm{log}\,R}{R}+\frac{{\mathcal{B}}_{2}}{{R}^{2}},\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\frac{{\mathcal{U}}_{2}}{{R}^{2}},\hfill \end{align} \tag{ 42 }$$

with ${\mathcal{G}}_{2}={o}^{+}(R)$, ${\mathcal{B}}_{2}={o}^{+}(R)$ and ${\mathcal{U}}_{2}={o}^{+}(R)$, so in order to show (41), we have to show that if we can write the evolved fields as,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{m=2}^{n-1}\,\sum\limits _{k=0}^{{N}_{n-1}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{G}}_{n}}{{R}^{n}},\hfill \\ \hfill & b=\frac{{B}_{1}}{R}+\sum\limits _{m=2}^{n-1}\,\sum\limits _{k=0}^{{N}_{n-1}^{b}}\frac{{(\mathrm{log}\,R)}^{k}{B}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{B}}_{n}}{{R}^{n}},\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{m=2}^{n-1}\,\sum\limits _{k=0}^{{N}_{n-1}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{U}}_{n}}{{R}^{n}},\hfill \end{align} \tag{ 43 }$$

where B₁ = B_1,0(ψ*) + log RB_1,1(ψ*), ${\mathcal{G}}_{n}={o}^{+}(R)$, ${\mathcal{B}}_{n}={o}^{+}(R)$ and ${\mathcal{U}}_{n}={o}^{+}(R)$, then we can also write them as,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{m=2}^{n}\,\sum\limits _{k=0}^{{N}_{n}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{G}}_{n+1}}{{R}^{n+1}},\hfill \\ \hfill & b=\frac{{B}_{1}}{R}+\sum\limits _{m=2}^{n}\,\sum\limits _{k=0}^{{N}_{n}^{b}}\frac{{(\mathrm{log}\,R)}^{k}{B}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{B}}_{n+1}}{{R}^{n+1}},\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{m=2}^{n}\,\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{U}}_{n+1}}{{R}^{n+1}},\hfill \end{align} \tag{ 44 }$$

where ${\mathcal{G}}_{n+1}={o}^{+}(R)$, ${\mathcal{B}}_{n+1}={o}^{+}(R)$ and ${\mathcal{U}}_{n+1}={o}^{+}(R)$. The cases of g and b go through in exactly the same way as in [36], so we will focus our attention on u.

Induction proof: assuming (43), we take the ugly equation in (32) and formally equate terms proportional to R⁻ⁿ⁻¹. Putting all terms with ${\mathcal{U}}_{n}$ on the lhs and all the rest on the rhs we get,

$$\begin{equation}{\nabla }_{\psi }({R}^{p+1-n}{\nabla }_{T}{\mathcal{U}}_{n})\simeq {R}^{p-n}{{\Omega}}_{n-1}^{u}.\end{equation} \tag{ 45 }$$

where on the rhs ${{\Omega}}_{n-1}^{u}$ depends on the functions {G_m,k , B_m,k , U_m,k , m_u,1}, for m ∈ [1, n − 1] and $k\in [0,{N}_{n}^{\phi }]$, and their derivatives. Also, ${{\Omega}}_{0}^{u}{:=}0$. We can split ${{\Omega}}_{n-1}^{u}$ in the following way,

$$\begin{equation}{{\Omega}}_{n-1}^{u}=\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}{{\Omega}}_{n-1,i}^{u}({\psi }^{\ast }).\end{equation} \tag{ 46 }$$

It is worth noting that the specific form of ${{\Omega}}_{n-1,i}^{u}$ has no influence on the proof of our hypothesis, as long as it is possible to write (46). We can now integrate (45) in order to get the asymptotic behavior of ${\mathcal{U}}_{n}$ in terms of {G_m,k , B_m,k , U_m,k , m_u,1}. Plugging (46) into (45),

$$\begin{equation*}{\nabla }_{\psi }({R}^{p+1-n}{\nabla }_{T}{\mathcal{U}}_{n})\simeq {R}^{p-n}\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}{{\Omega}}_{n-1,i}^{u}({\psi }^{\ast }),\end{equation*}$$

and integrating it along integral curves of ψ^a and ${\partial }_{T}^{a}$ gives, once more, two different results. For n ≠ p + 1,

$$\begin{align}\hfill {\mathcal{U}}_{n}\simeq & \sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}\sum\limits _{j=i}^{{N}_{n}^{{{\Omega}}^{u}}}-\frac{1}{{(n-1-p)}^{j-i+1}}\frac{j!}{i!}{\int }_{{T}_{0}}^{T}{{\Omega}}_{n-1,j}^{U}\mathrm{d}{T}^{\prime }\,+{m}_{u,n}\hfill \\ \hfill & =\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}{(\mathrm{log}\,R)}^{i}{U}_{n,i}({\psi }^{\ast }),\hfill \end{align} \tag{ 47 }$$

whereas for n = p + 1,

$$\begin{align}\hfill {\mathcal{U}}_{n}\simeq \,& \sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}}\frac{{(\mathrm{log}\,R)}^{i+1}}{i+1}{\int }_{{T}_{0}}^{T}{{\Omega}}_{n-1,i}^{U}\mathrm{d}{T}^{\prime }+{m}_{u,n}\hfill \\ \hfill & =\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{u}}+1}{(\mathrm{log}\,R)}^{i}{U}_{n,i}({\psi }^{\ast }).\hfill \end{align} \tag{ 48 }$$

By induction,

$$\begin{equation}u=\frac{{m}_{u,1}}{R}+\sum\limits _{n=2}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{n,k}({\psi }^{\ast })}{{R}^{n}}.\end{equation} \tag{ 49 }$$

This concludes the proof. These results can be packaged in the following theorem.

Theorem 2. Let there be any number of evolved fields, each satisfying one of the following three wave equations,

$$\begin{align}\hfill & \stackrel{\circ}{\square }g={\mathcal{N}}_{g},\hfill \\ \hfill & \stackrel{\circ}{\square }b=-\frac{1}{2}{({\nabla }_{T}{g}_{1})}^{2}-\frac{1}{2}{({\nabla }_{T}{g}_{2})}^{2}+{\mathcal{N}}_{b},\hfill \\ \hfill & \stackrel{\circ}{\square }u=\frac{2p}{R}{\nabla }_{T}u+{\mathcal{N}}_{u},\hfill \end{align} \tag{ 50 }$$

where g_ab is an asymptotically flat metric, g₁ and g₂ satisfy the first equation in (32), and p may be different for different functions. Fields that satisfy the equations for g, b and u admit polyhomogeneous expansions near null infinity of the types,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{n,k}({\psi }^{\ast })}{{R}^{n}}\hfill \\ \hfill & b=\frac{{B}_{1}}{R}+\sum\limits _{n=2}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{b}}\frac{{(\mathrm{log}\,R)}^{k}{B}_{n,k}({\psi }^{\ast })}{{R}^{n}}\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{n=2}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{n,k}({\psi }^{\ast })}{{R}^{n}},\hfill \end{align} \tag{ 51 }$$

respectively, with initial data on $\mathcal{S}$ of the type,

$$\begin{equation}\begin{cases}g{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{g,n}}{{R}^{n}}\quad \hfill \\ b{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{b,n}}{{R}^{n}}\quad \hfill \\ u{\vert }_{\mathcal{S}}=\sum\limits _{n=1}^{\infty }\frac{{m}_{u,n}}{{R}^{n}}\quad \hfill \end{cases} ,\quad \begin{cases}{\nabla }_{T}g{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \\ {\nabla }_{T}b{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \\ {\nabla }_{T}u{\vert }_{\mathcal{S}}={O}_{\mathcal{S}}({R}^{-2})\quad \hfill \end{cases} ,\end{equation} \tag{ 52 }$$

where m_ϕ,n are scalar functions that are independent of T and R. This is valid outside a compact ball centered at R = 0.

Remark 1. Consider a system in which all of the evolved fields satisfy either the good equation or the ugly one. We know that the good equation cannot create logs, it can only inherit them though coupling with the other equations. Since there is no bad field in such a system, the only logs that may appear are the ones created by the ugly equation. Let us also assume that all the ugly fields have the same p. In that case, for all n < p + 1, there can be no logs in (46), or in other words, ${N}_{n}^{{{\Omega}}^{u}}=0$. This implies that up to order p + 1,

$$\begin{equation}{\mathcal{U}}_{n}\simeq {U}_{n,0}({\psi }^{\ast }).\end{equation} \tag{ 53 }$$

However, when we get to order n = p + 1, because the integration generates one more power of log R than what already existed, in our system this is the first order at which a log may appear. For n > p + 1, integration along ψ no longer generates higher powers of logs, so any increase in ${N}_{n}^{u}$ may only come from nonlinearities. So our system admits the following polyhomogeneous expansion,

$$\begin{align}\hfill & g=\frac{{G}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{n=2}^{p}\frac{{G}_{n,0}({\psi }^{\ast })}{{R}^{n}}+\sum\limits _{n=p+1}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{g}}\frac{{(\mathrm{log}\,R)}^{k}{G}_{n,k}({\psi }^{\ast })}{{R}^{n}}\hfill \\ \hfill & u=\frac{{m}_{u,1}}{R}+\sum\limits _{n=2}^{p}\frac{{U}_{n,0}({\psi }^{\ast })}{{R}^{n}}+\sum\limits _{n=p+1}^{\infty }\,\sum\limits _{k=0}^{{N}_{n}^{u}}\frac{{(\mathrm{log}\,R)}^{k}{U}_{n,k}({\psi }^{\ast })}{{R}^{n}}.\hfill \end{align} \tag{ 54 }$$

In other words, in a system with the properties described above, the polyhomogeneous expansions can have no logs up to order p + 1. This result will be critical to our analysis of peeling.

In this section we derive the field equations of GR in GHG using our preferred variables. In this work we are primarily concerned with the asymptotic properties of solutions, which we obtain by brute-force integration, with many terms turning out to be irrelevant. Therefore we do not focus on obtaining the cleanest possible geometric derivation.

A concise way to find some of the wave equations is by computing the wave equation for the null covector σ_a and commuting derivatives to get some component of the Ricci tensor. From (5) we get,

$$\begin{align*}\hfill \square {\sigma }_{a}& =-\square {\mathcal{C}}_{+}^{R}{\nabla }_{a}T+\square {\mathcal{C}}_{A}^{+}{\nabla }_{a}{\theta }^{A}-2{g}^{bc}{\nabla }_{b}{\mathcal{C}}_{+}^{R}{\nabla }_{c}{\nabla }_{a}T+2{g}^{bc}{\nabla }_{b}{\mathcal{C}}_{A}^{+}{\nabla }_{c}{\nabla }_{a}{\theta }^{A}\hfill \\ \hfill & \quad -{\mathcal{C}}_{+}^{R}{g}^{bc}{\nabla }_{b}{\nabla }_{c}{\nabla }_{a}T+{g}^{bc}{\nabla }_{b}{\nabla }_{c}{\nabla }_{a}R+{\mathcal{C}}_{A}^{+}{g}^{bc}{\nabla }_{b}{\nabla }_{c}{\nabla }_{a}{\theta }^{A}.\hfill \end{align*}$$

Commuting the Levi-Civita covariant derivatives ∇ and using the Ricci identity, reveals

$$\begin{align}\hfill \square {\sigma }_{a}& ={R}_{a\sigma }-\square {\mathcal{C}}_{+}^{R}{\nabla }_{a}T+\square {\mathcal{C}}_{A}^{+}{\nabla }_{a}{\theta }^{A}-2{g}^{bc}{\nabla }_{b}{\mathcal{C}}_{+}^{R}{\nabla }_{c}{\nabla }_{a}T\hfill \\ \hfill & \quad +2{g}^{bc}{\nabla }_{b}{\mathcal{C}}_{A}^{+}{\nabla }_{c}{\nabla }_{a}{\theta }^{A}-{\mathcal{C}}_{+}^{R}{\nabla }_{a}\square T+{\nabla }_{a}\square R+{\mathcal{C}}_{A}^{+}{\nabla }_{a}\square {\theta }^{A}.\hfill \end{align} \tag{ 55 }$$

Replacing the first derivatives of metric functions with (16), using the fact that,

$$\begin{equation}\square T=-{\stackrel{\bullet }{{\Gamma}}}^{a}{\nabla }_{a}T,\end{equation} \tag{ 56 }$$

$$\begin{equation}\square R=-{\stackrel{\bullet }{{\Gamma}}}^{a}{\nabla }_{a}R,\end{equation} \tag{ 57 }$$

$$\begin{equation}\square {\theta }^{A}=-{\stackrel{\bullet }{{\Gamma}}}^{a}{\nabla }_{a}{\theta }^{A},\end{equation} \tag{ 58 }$$

and writing ∇_a T in terms of the null covectors with (6),

$$\begin{equation}\tau {\nabla }_{a}T=-{\sigma }_{a}-{\underline{\sigma }}_{a}+{\mathcal{C}}_{A}^{+}{\nabla }_{a}{\theta }^{A}+{\mathcal{C}}_{A}^{-}{\nabla }_{a}{\theta }^{A},\end{equation} \tag{ 59 }$$

we get the expression,

$$\begin{align}\hfill \square {\sigma }_{a}& ={R}_{a\sigma }-\square {\mathcal{C}}_{+}^{R}{\nabla }_{a}T+\square {\mathcal{C}}_{A}^{+}{\nabla }_{a}{\theta }^{A}-{\sigma }_{b}{\nabla }_{a}{\stackrel{\bullet }{{\Gamma}}}^{b}+{\stackrel{\bullet }{{\Gamma}}}^{b}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\sigma }{}_{\!\!b}\hfill \\ \hfill & \quad +\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!b}({\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\sigma b}+{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\underline{\sigma }b})-2{\stackrel{\bullet }{{\Gamma}}}_{\mathit{bAa}}{\stackrel{\bullet }{{\Gamma}}}^{(\sigma A)b}.\hfill \end{align} \tag{ 60 }$$

We can contract this equation with ψ^a , ${\underline{\psi }}^{a}$ and ${{{/\hskip-3pt \!g}}_{\,b}}^{a}$ to obtain wave equations for ${\mathcal{C}}_{+}^{R}$, φ and ${\mathcal{C}}_{A}^{+}$, respectively, and we can contract the analogous equation for $\square {\underline{\sigma }}_{a}$ with ${\underline{\psi }}^{a}$ and ${{\not{\hskip-3pt g}}_{\,b}}^{a}$ to obtain wave equations for ${\mathcal{C}}_{-}^{R}$ and ${\mathcal{C}}_{A}^{-}$. First we introduce the reduced Ricci tensor

$$\begin{align}\hfill {\mathcal{R}}_{ab}{:=}{R}_{ab}-{\nabla }_{(a}{Z}_{b)}+{W}_{ab}\quad \text{with}\quad {Z}^{a}{:=}{\stackrel{\bullet }{{\Gamma}}}^{a}+{F}^{a}& \hfill & \hfill \end{align} \tag{ 61 }$$

where ${\stackrel{\bullet }{{\Gamma}}}^{a}{:=}{g}^{bc}{\stackrel{\bullet }{{\Gamma}}}_{b}{}^{\!a}{}_{\!\!c}$, F^a are the gauge source functions and W_ab denotes a generic constraint addition where by generic we mean any expression for which the propagation of the constraints construction holds. It can be shown, see for instance [21], that as long as W_ab is a homogeneous expression in Z^a , i.e. Z^a = 0 ⟹ W_ab = 0, the constraints propagate. In other words, if Z^a = 0 and ∇_a Z^b = 0 on a spacelike hypersurface ${\Sigma}\subset \mathcal{M}$ then Z^a = 0 in the future domain of dependence ${\mathcal{D}}^{+}({\Sigma})$. For conciseness we also introduce the tensor ${\tilde{\mathcal{R}}}_{ab}$ defined as,

$$\begin{equation}{\tilde{\mathcal{R}}}_{ab}{:=}{\mathcal{R}}_{ab}+{\nabla }_{(a}{F}_{b)}-{W}_{ab}.\end{equation} \tag{ 62 }$$

Now we contract (60) with ψ^a to obtain,

$$\begin{align}\hfill {\psi }^{a}\square {\sigma }_{a}& ={e}^{-\varphi }{\tilde{\mathcal{R}}}_{\psi \psi }-\stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}+\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!a}({\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma a}+{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\underline{\sigma }a})-4{\not{\hskip-3pt g}}_{\,bc}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!b}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}^{(\sigma c)a}.\hfill \end{align} \tag{ 63 }$$

Note that all derivatives of ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!c}$ canceled and also that we have replaced $\square {\mathcal{C}}_{+}^{R}$ with $\stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}$. On the other hand, a direct computation yields,

$$\begin{equation}{\psi }^{a}\square {\sigma }_{a}=-{e}^{-\varphi }{\not{\hskip-3pt g}}_{bc}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!b}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!ca}.\end{equation} \tag{ 64 }$$

We can now equate the lhs of equations (63) and (64), and solve for $\stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}$. This gives,

$$\begin{align}\hfill \stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}& ={e}^{-\varphi }{\tilde{\mathcal{R}}}_{\psi \psi }+\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma a}-{e}^{-\varphi }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!b}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\psi b}{}^{\!a}{}_{}\!\!-2{{{/\hskip-3pt g}}_{a}}^{c}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!a}{}_{\!\!b}{\stackrel{\bullet }{{\Gamma}}}_{c}{}^{\!\sigma b}.\hfill \end{align} \tag{ 65 }$$

If we contract (60) instead with ${\underline{\psi }}^{a}$ we get,

$$\begin{align}\hfill {\underline{\psi }}^{a}\square {\sigma }_{a}& ={e}^{-\varphi }{\tilde{\mathcal{R}}}_{\underline{\psi }\psi }-{e}^{-\varphi }{\psi }^{a}{\nabla }_{\underline{\psi }}{\stackrel{\bullet }{{\Gamma}}}_{a}+{e}^{-\varphi }{\underline{\psi }}^{b}{\psi }^{d}{\nabla }_{(d}{\stackrel{\bullet }{{\Gamma}}}_{b)}-\square {\mathcal{C}}_{+}^{R}\hfill \\ \hfill & \quad +\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!a}\left(\stackrel{\bullet }{{\Gamma}}{\underline{\psi }}^{\,\sigma a}+\stackrel{\bullet }{{\Gamma}}{\underline{\psi }}^{\underline{\sigma }a}\right)-2{{/\hskip-3pt g}}_{de}\stackrel{\bullet }{{\Gamma}}\underline{\psi }{}^{\!d}{}_{\!\!c}\left({\stackrel{\bullet }{{\Gamma}}}^{\sigma ec}+{\stackrel{\bullet }{{\Gamma}}}^{e\underline{\sigma }c}\right)\hfill \end{align} \tag{ 66 }$$

Computing the lhs without commutation gives,

$$\begin{align}\hfill {\underline{\psi }}^{a}\square {\sigma }_{a}& =-{e}^{-\varphi }{\tilde{\mathcal{R}}}_{\underline{\psi }\psi }+{e}^{-\varphi }{\underline{\psi }}^{a}{\nabla }_{\psi }{\stackrel{\bullet }{{\Gamma}}}_{a}-{e}^{-\varphi }{\underline{\psi }}^{b}{\psi }^{d}{\nabla }_{(d}{\stackrel{\bullet }{{\Gamma}}}_{b)}\hfill \\ \hfill & \quad \times \frac{1}{\tau }\left({\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\sigma }{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\sigma a}{}_{\!\!}+2{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\sigma }{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma a}-{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\underline{\sigma }}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\underline{\sigma }a}{}_{\!\!}\right)-\square {\mathcal{C}}_{+}^{R}+\tau \square \varphi -{\stackrel{\bullet }{{\Gamma}}}^{a}{\stackrel{\bullet }{{\Gamma}}}_{\psi a}^{\underline{\sigma }},\hfill \end{align} \tag{ 67 }$$

where, in order to eliminate terms like ${g}^{cd}{\nabla }_{c}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!a}{}_{\!\!bd}$, we have used the following identity,

$$\begin{equation}{g}^{cd}{\nabla }_{c}\stackrel{\bullet }{{\Gamma}}{}^{\!a}{}_{\!\!bd}={\stackrel{\bullet }{{\Gamma}}}^{\mathit{acd}}{\stackrel{\bullet }{{\Gamma}}}_{\mathit{cbd}}-{{R}_{b}}^{a}+{\nabla }_{b}{{\Gamma}}^{a}.\end{equation} \tag{ 68 }$$

Putting the two sides together cancels all derivatives of ${\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!a}$. Solving for $\stackrel{\bullet }{\square }\varphi$ finally gives,

$$\begin{align}\hfill \stackrel{\bullet }{\square }\varphi & =\frac{2}{\tau }{{\tilde{\mathcal{R}}}_{\psi }}^{\underline{\sigma }}-\frac{2}{{\tau }^{2}}\left[{{{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}}^{}{}^{ a}}_{(\psi }{{{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi })}}^{}{}^{ \sigma }}_{a}+{{{\stackrel{\bullet }{{\Gamma}}}_{\psi }}^{}{}^{ a}}_{(\underline{\psi }}{{{\stackrel{\bullet }{{\Gamma}}}_{\psi )}}^{}{}^{\!\underline{\sigma }}}_{a}\right]\hfill \\ \hfill & \quad +{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!db}\stackrel{\bullet }{{\Gamma}}{}^{ ac}{}_{\!e}\left[2{{/\hskip-3pt \!g}}_{cb}({\delta }_{d}^{e}-{{{/\hskip-3pt \!g}}_{d}}^{e})+{{/\hskip-3pt \!g}}_{cd}({\delta }_{b}^{e}-{{{/\hskip-3pt \!g}}_{b}}^{e})\right].\hfill \end{align} \tag{ 69 }$$

The contraction of (60) with ${{\not{\hskip-3pt g}}_{b}}^{a}$ yields,

$$\begin{align}\hfill {{\not{\hskip-3pt g}}_{\,A}}^{a}\square {\sigma }_{a}& ={e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,A}}^{a}{\tilde{\mathcal{R}}}_{a\psi }+\frac{{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}}{\tau }\stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}+\stackrel{\bullet }{\square }{\mathcal{C}}_{A}^{+}+\frac{2}{\tau }{{\not{\hskip-3pt g}}_{A}}^{a}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!b}\left({\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\sigma b}+{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!\underline{\sigma }b}\right)\hfill \\ \hfill & \quad -4{{\not{\hskip-3pt g}}_{\,A}}^{a}{\not{\hskip-3pt g}}_{de}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!d}{}_{\!\!b}\!{\stackrel{\bullet }{{\Gamma}}}^{(\sigma e)b}-{e}^{-\varphi }{{\not{\hskip-3pt g}}_{A}}^{a}{\psi }^{b}{\nabla }_{[a}{\stackrel{\bullet }{{\Gamma}}}_{c]}-{e}^{-\varphi }{\not{\hskip-3pt g}}_{aA}{\stackrel{\bullet }{{\Gamma}}}^{c}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!c}.\hfill \end{align} \tag{ 70 }$$

On the other hand, the lhs gives,

$$\begin{align}\hfill {{\not{\hskip-3pt g}}_{\,A}}^{a}\square {\sigma }_{a}& =-{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,A}}^{a}{\tilde{\mathcal{R}}}_{a\psi }-{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,A}}^{a}{\psi }^{b}{\nabla }_{[a}{\stackrel{\bullet }{{\Gamma}}}_{c]}-{e}^{-\varphi }{\not{\hskip-3pt g}}_{\,aA}{\stackrel{\bullet }{{\Gamma}}}^{c}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!a}{}_{\!\!c}\hfill \\ \hfill & \quad +\frac{2}{\tau }{\not{\hskip-3pt g}}_{\,Ab}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!b}{}_{\!\!a}\left({\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\sigma a}-{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma a}\right)+2{e}^{-\varphi }{\not{\hskip-3pt g}}_{\,Ab}{{\not{\hskip-3pt g}}_{\,c}}^{e}{\stackrel{\bullet }{{\Gamma}}}_{a}{}^{\!b}{}_{\!\!e}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!ca}.\hfill \end{align} \tag{ 71 }$$

Putting (70) and (71) together yields,

$$\begin{align}\hfill \stackrel{\bullet }{\square }{\mathcal{C}}_{A}^{+}& =-{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,A}}^{a}{\tilde{\mathcal{R}}}_{a\psi }+\frac{{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}}{\tau }\stackrel{\bullet }{\square }{\mathcal{C}}_{+}^{R}+\frac{2}{\tau }{e}^{-\varphi }{\stackrel{\bullet }{{\Gamma}}}_{\psi A}{}^{\!a}({\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\sigma }{}_{\!\!a}-{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!a})\hfill \\ \hfill & \quad -\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\sigma }{}_{\!\!b}({\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!\sigma b}+{\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!\underline{\sigma }b})+2{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,f}}^{e}{\stackrel{\bullet }{{\Gamma}}}_{aAe}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!fa}-4{\not{\hskip-3pt g}}_{\,de}{\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!d}{}_{\!\!b}{\stackrel{\bullet }{{\Gamma}}}^{(\sigma e)b}.\hfill \end{align} \tag{ 72 }$$

A completely analogous calculation using the contractions ${\underline{\psi }}^{a}\square {\underline{\sigma }}_{a}$ and ${{\not{\hskip-3pt g}}_{\,A}}^{a}\square {\underline{\sigma }}_{a}$ allows us to find wave equations for ${\mathcal{C}}_{-}^{R}$ and ${\mathcal{C}}_{A}^{-}$, respectively.

$$\begin{align}\hfill \stackrel{\bullet }{\square }{\mathcal{C}}_{-}^{R}& =-{e}^{-\varphi }{\tilde{\mathcal{R}}}_{\underline{\psi }\underline{\psi }}-\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\underline{\sigma }}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\underline{\sigma }a}+{e}^{-\varphi }{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!b}{}_{\!\!a}{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }b}{}^{\!a}+2{{\not{\hskip-3pt g}}_{\,A}}^{c}{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!a}{}_{\!\!b}{\stackrel{\bullet }{{\Gamma}}}_{c}{}^{\!\underline{\sigma }b},\hfill \end{align} \tag{ 73 }$$

and,

$$\begin{align}\hfill \stackrel{\bullet }{\square }{\mathcal{C}}_{A}^{-}& =-{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,A}}^{a}{\tilde{\mathcal{R}}}_{a\underline{\psi }}-\frac{{\mathcal{C}}_{A}^{+}+{\mathcal{C}}_{A}^{-}}{\tau }\stackrel{\bullet }{\square }{\mathcal{C}}_{-}^{R}+\frac{2}{\tau }{e}^{-\varphi }{\stackrel{\bullet }{{\Gamma}}}_{\psi A}{}^{\!a}({\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!\underline{\sigma }}{}_{\!\!a}-{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\underline{\sigma }}{}_{\!\!a})\hfill \\ \hfill & \quad -\frac{2}{\tau }{\stackrel{\bullet }{{\Gamma}}}_{\underline{\psi }}{}^{\!\underline{\sigma }}{}_{\!\!b}({\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!\underline{\sigma }b}+{\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!\sigma b})+2{e}^{-\varphi }{{\not{\hskip-3pt g}}_{\,f}}^{e}{\stackrel{\bullet }{{\Gamma}}}_{aAe}{\stackrel{\bullet }{{\Gamma}}}_{\psi }{}^{\!fa}-4{\not{\hskip-3pt g}}_{\,de}{\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!d}{}_{\!\!b}{\stackrel{\bullet }{{\Gamma}}}^{(\underline{\sigma }e)b}.\hfill \end{align} \tag{ 74 }$$

To obtain the final 3 equations associated to the angular sector of the metric g_ab : , h₊ and h_×, we recall the second Cartan structure equation, which in our notation reads,

$$\begin{equation}{{R}_{\mathit{abc}}}^{d}=2{g}^{ef}\stackrel{\bullet }{{\Gamma}}[a\vert e_{}\!\!{\stackrel{\bullet }{{\Gamma}}}_{f\vert b]c}-2{\nabla }_{[a}\stackrel{\bullet }{{\Gamma}}{}^{\!d}{}_{b}]c.\end{equation} \tag{ 75 }$$

Tracing and projecting equation (75) appropriately gives,

$$\begin{equation}{\not{\hskip-3pt g}}^{\,ab}{R}_{ab}={\stackrel{\bullet }{{\Gamma}}}^{\mathit{abd}}{\stackrel{\bullet }{{\Gamma}}}_{bd}{}^{\!c}{\not{\hskip-3pt g}}_{ac}-{\stackrel{\bullet }{{\Gamma}}}^{a}{\stackrel{\bullet }{{\Gamma}}}^{b}{}_{\!a}{}^{\!\!c}{\not{\hskip-3pt g}}_{bc}+{\not{\hskip-3pt g}}^{\,ab}{\nabla }_{b}{\stackrel{\bullet }{{\Gamma}}}_{a}-{\not{\hskip-3pt g}}^{\,ab}{\nabla }_{c}{\stackrel{\bullet }{{\Gamma}}}_{ab}{}^{\!c}.\end{equation} \tag{ 76 }$$

A calculation using (16) renders the following wave equation for ,

$$\begin{equation}\stackrel{\bullet }{\square }({\epsilon}+\stackrel{\bullet }{{\epsilon}})=-\tilde{{/\hskip-3pt \mathcal{R}}}-{\stackrel{\bullet }{{\Gamma}}}^{\mathit{bac}}{\stackrel{\bullet }{{\Gamma}}}_{e}{}^{\!d}{}_{\!\!c}\left[{{/\hskip-3pt g}}_{ad}({\delta }_{b}^{e}-{{\not{\hskip-3pt g}}_{b}}^{e})-{\not{\hskip-3pt g}}_{bd}{{{/\hskip-3pt g}}_{\,A}}^{e}\right],\end{equation} \tag{ 77 }$$

To get the remaining two wave equations we trace (75) again, raise its indices and use (62) to obtain

$$\begin{equation}{\tilde{\mathcal{R}}}^{ab}=-{\nabla }_{c}{\stackrel{\bullet }{{\Gamma}}}^{(ab)c}-{\stackrel{\bullet }{{\Gamma}}}^{c}{\stackrel{\bullet }{{\Gamma}}}^{(ab)}{}_{\!c}+{{\stackrel{\bullet }{{\Gamma}}}_{d}{}^{\!c}}^{(a}{\stackrel{\bullet }{{\Gamma}}}_{c}{}^{\!b)d}.\end{equation} \tag{ 78 }$$

We need to eliminate the first term on the rhs and to do that we take a Shell derivative of the inverse metric,

$$\begin{equation}{\stackrel{\bullet }{\nabla }}_{c}{g}^{ab}=-{g}^{ad}{g}^{be}{\stackrel{\bullet }{\nabla }}_{c}{g}_{de}=-2{\stackrel{\bullet }{{\Gamma}}}^{(ab)}{}_{\!c},\end{equation} \tag{ 79 }$$

then we take a divergence on the downstairs index with the Levi-Civita covariant derivative and change that derivative to a $\stackrel{\bullet }{\nabla }$ with the appropriate transition tensor,

$$\begin{align}\hfill -2{\nabla }_{c}{\stackrel{\bullet }{{\Gamma}}}^{(ab)c}& =\stackrel{\bullet }{\square }{g}^{ab}+2{\stackrel{\bullet }{{\Gamma}}}^{c}{\stackrel{\bullet }{{\Gamma}}}^{(ab)}{}_{\!c}-{\stackrel{\bullet }{2{\Gamma}}}^{ca}{}_{\!d}{\stackrel{\bullet }{{\Gamma}}}^{(db)}{}_{\!c}-2{\stackrel{\bullet }{{\Gamma}}}^{cb}{}_{\!d}{\stackrel{\bullet }{{\Gamma}}}^{(da)}{}_{\!c}.\hfill \end{align} \tag{ 80 }$$

We can now replace (80) in (78) to get,

$$\begin{equation}{\tilde{\mathcal{R}}}^{ab}=\frac{1}{2}\stackrel{\bullet }{\square }{g}^{ab}+{\stackrel{\bullet }{{\Gamma}}}_{c}{}^{\!a}{}_{\!\!d}{\stackrel{\bullet }{{\Gamma}}}^{\mathit{cbd}}.\end{equation} \tag{ 81 }$$

On the other hand, a rather long calculation using (10) gives,

$$\begin{align}\hfill \stackrel{\bullet }{\square }{\not{\hskip-3pt g}}^{\,AB}& =\frac{{e}^{-{\epsilon}}}{{R}^{2}}\stackrel{\bullet }{\square }{({q}^{-1})}^{AB}-{\not{\hskip-3pt g}}^{\,AB}\left[\stackrel{\bullet }{\square }({\epsilon}+\stackrel{\bullet }{{\epsilon}})-{\not{\hskip-3pt g}}^{\,\theta \theta }\frac{1+{\cos }^{2}\theta }{{\mathrm{sin}}^{2}\theta }+2\,\mathrm{cot}\,\theta {{\not{\hskip-3pt g}}_{\,e}}^{d}{\not{\hskip-3pt g}}^{\,c\theta }{\stackrel{\bullet }{{\Gamma}}}_{d}{}^{\!e}{}_{\!\!c}\right]\hfill \\ \hfill & \quad -4{{\not{\hskip-3pt g}}_{\,f}}^{A}{{\not{\hskip-3pt g}}_{\,e}}^{B}{\stackrel{\bullet }{{\Gamma}}}^{(ef)}{}_{\!c}\left[\mathrm{cot}\,\theta {\not{\hskip-3pt g}}^{\,c\theta }-{{\not{\hskip-3pt g}}_{\,e}}^{d}{\stackrel{\bullet }{{\Gamma}}}_{d}{}^{\!ec}\right].\hfill \end{align} \tag{ 82 }$$

Taking only the angular components of equation (81), we can substitute the first term on the rhs with (82) and project the whole equation with,

$$\begin{equation*}{{{{\not{\hskip-3pt S}}_{a}}^{c}}_{b}}^{d}={{\not{\hskip-3pt g}}_{a}}^{(c}{{\not{\hskip-3pt g}}_{b}}^{d)}-\frac{1}{2}{\not{\hskip-3pt g}}_{ab}{\not{\hskip-3pt g}}^{\,cd}.\end{equation*}$$

Finally, for the two remaining metric functions h₊ and h_×, we get the concise expression,

$$\begin{equation}{{{{\not{\hskip-3pt S}}_{A}}^{c}}_{B}}^{d}\stackrel{\bullet }{\square }{q}^{AB}=2{R}^{2}{e}^{{\epsilon}}{{{{\not{\hskip-3pt S}}_{A}}^{c}}_{B}}^{d}\left[{\tilde{\mathcal{R}}}^{AB}+{\stackrel{\bullet }{{\Gamma}}}_{A}{}^{\!e}{}_{\!\!f}{\stackrel{\bullet }{{\Gamma}}}^{\mathit{eBf}}-{\stackrel{\bullet }{2{\Gamma}}}^{(AB)}{}_{\!c}\left({{\not{\hskip-3pt g}}_{f}}^{e}{\stackrel{\bullet }{{\Gamma}}}_{e}{}^{\!fc}-\mathrm{cot}\,c{\not{\hskip-3pt g}}^{\,c\theta }\right)\right].\end{equation} \tag{ 83 }$$

With this we conclude the derivation of the reduced EFE, a system of non-linear wave equations, one for each of the 10 metric components, (65), (69), (72)–(74), (77) and (83). The next step in this analysis is to see how to identify each field as good, bad or ugly and we explore that in the next section.

The aim of this section is to find the functional form of the first few orders of the polyhomogeneous expansions for our metric variables near null infinity in GR. In other words, we want to use the tools of the last sections to discover the maximum power of log R at each order n for each metric function ϕ. To do that we begin by finding the terms that can contribute to leading order, then we verify that each field can be classified as a good, a bad or an ugly. Finally, we use theorem 2 to show that these fields admit a polyhomogeneous expansion and we integrate the wave equations to find the maximum power of log R, ${N}_{n}^{\phi }$, that is allowed. To illustrate the fact that the information about the functional form of a field ϕ to some order n is contained in the set of integers ${N}_{n}^{\phi }$, we consider an example: let ϕ have the functional form,

$$\begin{equation*}{N}_{1}^{\phi }=0,\qquad {N}_{2}^{\phi }=2,\end{equation*}$$

then we can write its first 2 orders as,

$$\begin{equation*}\phi \simeq \frac{{{\Phi}}_{1,0}}{R}+\frac{{{\Phi}}_{2,0}+\mathrm{log}\,R{{\Phi}}_{2,1}+{(\mathrm{log}\,R)}^{2}{{\Phi}}_{2,2}}{{R}^{2}},\end{equation*}$$

where the functions Φ_n,i = Φ_n,i (ψ*).

Cartesian harmonic gauge: observe that up to this point the gauge source functions in ${\tilde{\mathcal{R}}}_{ab}$ have not been specified. Let,

$$\begin{equation}{F}^{a}={\stackrel{\circ}{F}}^{a},\end{equation} \tag{ 84 }$$

where,

$$\begin{equation}{\stackrel{\circ}{F}}^{a}={g}^{bc}{\Gamma}[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]{}^{a}_{bc}.\end{equation} \tag{ 85 }$$

Notice that this choice implies, using (3) and (61), that ${g}^{bc}{\Gamma}[\nabla ,\stackrel{\circ}{\nabla }]{}_{b}{}^{a}{}_{c}=0$. In other words, ${F}^{a}={\stackrel{\circ}{F}}^{a}$ implies that the Cartesian coordinates ${X}^{\underline{\mu }}$ are harmonic—recall that ${\Gamma}[\stackrel{\circ}{\nabla },\stackrel{\bullet }{\nabla }]{}_{a}{}^{b}{}_{c}$ are given functions of the coordinates. For clarity, we write down the four components of ${\stackrel{\circ}{F}}^{a}$ explicitly,

$$\begin{align}\hfill & {\stackrel{\circ}{F}}^{\sigma }=\frac{2{e}^{-{\epsilon}}}{R}\,\mathrm{cosh}\,{h}_{+}\mathrm{cosh}\,{h}_{\times }+{\mathcal{C}}_{A}^{+}{\stackrel{\circ}{F}}^{A},\hfill \\ \hfill & {\stackrel{\circ}{F}}^{\underline{\sigma }}=-\frac{2{e}^{-{\epsilon}}}{R}\,\mathrm{cosh}\,{h}_{+}\mathrm{cosh}\,{h}_{\times }+{\mathcal{C}}_{A}^{-}{\stackrel{\circ}{F}}^{A},\hfill \\ \hfill & {\stackrel{\circ}{F}}^{\theta }=\frac{{e}^{-{\epsilon}}\mathrm{cot}\,\theta }{{R}^{2}}{e}^{{h}_{+}}\mathrm{cosh}\,{h}_{\times }-\frac{2}{R}{\not{\hskip-3pt g}}^{\,R\theta },\hfill \\ \hfill & {\stackrel{\circ}{F}}^{\phi }=-\frac{2{e}^{-{\epsilon}}\mathrm{cot}\,\theta }{\mathrm{sin}\,\theta {R}^{2}}\mathrm{sinh}\,{h}_{\times }-\frac{2}{R}{\not{\hskip-3pt g}}^{\,R\phi }.\hfill \end{align} \tag{ 86 }$$

Constraint addition: we want to be able to separate all of the 10 metric functions into goods, bads and uglies in the same way as in [36]. As it turns out it is possible to write four of the metric functions as uglies,

$$\begin{equation}\stackrel{\circ}{\square }\phi =\frac{2}{R}{\nabla }_{T}\phi +{\mathcal{N}}_{\phi },\end{equation} \tag{ 87 }$$

through constraint addition, where ${\mathcal{N}}_{\phi }$ is a stratified null form and ϕ is a metric function. To see this, we need to compute the constraints Z^a to leading order and add specific constraint terms encoded in W_ab to some of the wave equations. Asymptotically, the constraints are the following,

$$\begin{align}\hfill & {Z}^{\sigma }=2{\nabla }_{T}{\mathcal{C}}_{+}^{R}+{o}^{+}(1),\hfill \\ \hfill & {Z}^{\underline{\sigma }}=-2{\nabla }_{T}{\epsilon}+{o}^{+}(1),\hfill \\ \hfill & {Z}^{\theta }=-\frac{1}{{R}^{2}}{\nabla }_{T}{\mathcal{C}}_{\theta }^{+}+{o}^{+}({R}^{-2}),\hfill \\ \hfill & {Z}^{\phi }=-\frac{1}{{R}^{2}\,{\mathrm{sin}}^{2}\theta }{\nabla }_{T}{\mathcal{C}}_{\phi }^{+}+{o}^{+}({R}^{-2}).\hfill \end{align} \tag{ 88 }$$

As each of the components of Z^a involves a ∇_T derivative of a specific metric function, we can choose the components of W_ab in order to turn ${\mathcal{C}}_{+}^{R}$, and ${\mathcal{C}}_{A}^{+}$ into ugly fields. One choice of components that does this, and this choice is highly non-unique, is the following,

$$\begin{align}\hfill & {W}_{\psi \psi }=-\frac{1}{2}{({Z}^{\sigma })}^{2}+\frac{1}{R}{Z}^{\sigma },\hfill \\ \hfill & \not{\hskip-3pt W}=-\frac{1}{R}{Z}^{\underline{\sigma }},\hfill \\ \hfill & {W}_{\psi A}=\frac{2}{R}{Z}^{A},\hfill \end{align} \tag{ 89 }$$

with all remaining components set to zero. Note that the indices in the last equation are not incorrect, since A here is not a tensorial index, but a label.

First order asymptotic system: in order to classify all wave equations as good, bad or ugly, we need to rescale ${\mathcal{C}}_{A}^{\pm }$ to account for the R factor in (28). Therefore we define,

$$\begin{equation}R{\hat{\mathcal{C}}}_{A}^{\pm }={\mathcal{C}}_{A}^{\pm }.\end{equation} \tag{ 90 }$$

Using Cartesian harmonic gauge, the constraint additions (89) and the wave equations (65), (69), (72)–(74), (77) and (83), the system takes the following form,

$$\begin{align}\hfill & \stackrel{\circ}{\square }\varphi ={\mathcal{N}}_{\varphi },\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{+}^{R}=\frac{2}{R}{\nabla }_{T}{\mathcal{C}}_{+}^{R}+{\mathcal{N}}_{{\mathcal{C}}_{+}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{-}^{R}=-\frac{1}{2}{({\nabla }_{T}{h}_{+})}^{2}-\frac{1}{2}{({\nabla }_{T}{h}_{\times })}^{2}+{\mathcal{N}}_{{\mathcal{C}}_{-}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{+}=\frac{2}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{+}+{\mathcal{N}}_{{\mathcal{C}}_{A}^{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{-}=\frac{4}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{-}+{\mathcal{N}}_{{\mathcal{C}}_{A}^{-}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\epsilon}=\frac{2}{R}{\nabla }_{T}{\epsilon}+{\mathcal{N}}_{{\epsilon}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{+}={\mathcal{N}}_{{h}_{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{\times }={\mathcal{N}}_{{h}_{\times }}.\hfill \end{align} \tag{ 91 }$$

We can easily see that the fields φ, h₊ and h_× satisfy the good equation, whereas ${\mathcal{C}}_{-}^{R}$ satisfies the bad one, ${\mathcal{C}}_{+}^{R}$, ${\hat{\mathcal{C}}}_{A}^{+}$ and satisfy the ugly equation with p = 1 and ${\hat{\mathcal{C}}}_{A}^{-}$ the ugly equation with p = 2. This already gives us the functional form of the fields to first order in R⁻¹,

$$\begin{align}\hfill & {N}_{1}^{\varphi }={N}_{1}^{{\mathcal{C}}_{+}^{R}}={N}_{1}^{{\mathcal{C}}_{A}^{\pm }}={N}_{1}^{{\epsilon}}={N}_{1}^{{h}_{+}}={N}_{1}^{{h}_{\times }}=0,\hfill \\ \hfill & {N}_{1}^{{\mathcal{C}}_{-}^{R}}=1.\hfill \end{align} \tag{ 92 }$$

Second and third order: we can now apply theorem 2 to assert that GR in Cartesian harmonic gauge and with constraint additions (89) admits a polyhomogeneous expansion of the type (32). In order to find the functional form for second order, we plug (32) into the wave equations and formally equate terms proportional to R⁻³, we put all terms containing Φ₂ on the lhs and all the rest we collect on the rhs and name it ${{\Omega}}_{1}^{\phi }$. Here ${{\Phi}}_{n}{:=}{\sum }_{k=0}^{{N}_{n}^{\phi }}{(\mathrm{log}\,R)}^{k}{{\Phi}}_{n,k}({\psi }^{\ast })$, where Φ_n,k = {G_n,k , B_n,k , U_n,k }. Note that we are not interested in the exact dependence of ${{\Omega}}_{1}^{\phi }$ on the fields, but we are interested in its functional form. Integrating the resulting equations we find that the functional form to second order is given by,

$$\begin{align}\hfill & {N}_{2}^{\varphi }={N}_{2}^{{\mathcal{C}}_{+}^{R}}={N}_{2}^{{\mathcal{C}}_{A}^{+}}={N}_{2}^{{\epsilon}}={N}_{1}^{{h}_{+}}={N}_{2}^{{h}_{\times }}=1,\hfill \\ \hfill & {N}_{2}^{{\mathcal{C}}_{-}^{R}}={N}_{2}^{{\mathcal{C}}_{A}^{-}}=2.\hfill \end{align} \tag{ 93 }$$

Equating terms of order R⁻⁴ and following exactly the same procedure we find the functional form to third order,

$$\begin{align}\hfill & {N}_{3}^{{\mathcal{C}}_{+}^{R}}=1,\hfill \\ \hfill & {N}_{3}^{\varphi }={N}_{3}^{{\epsilon}}={N}_{3}^{{\mathcal{C}}_{A}^{+}}={N}_{3}^{{h}_{+}}={N}_{3}^{{h}_{\times }}=2,\hfill \\ \hfill & {N}_{3}^{{\mathcal{C}}_{-}^{R}}={N}_{3}^{{\mathcal{C}}_{A}^{-}}=3.\hfill \end{align} \tag{ 94 }$$

Given the increasing complexity and computational time these expansions demand, we stop at third order, since that is already enough to show the result of the next section.

The full calculations needed to arrive at these results are given in the Mathematica notebooks that accompany the paper. They rely heavily on the xAct [39] package.

As stated above, smooth null infinity implies the satisfaction of peeling, i.e. that the components of the Weyl curvature tensor fall off with certain negative integer powers of radius. We will see here that our choice of gauge (84), together with our constraint addition (89) gives rise to a violation of peeling by introducing powers of log R in the leading order terms of some of the components. The Weyl tensor has 10 independent components which can be described concisely by the Weyl scalars, five complex numbers that result from the contraction of the Weyl tensor with unit null vectors. Following the conventions of [40], we define the null tetrad (with normalization ${l}_{a}{n}^{a}=-{m}_{a}{\bar{m}}^{a}=-1$):

$$\begin{align}\hfill & {l}^{a}=\frac{{\psi }^{a}}{\sqrt{\tau \,{\text{e}}^{\varphi }}},\hfill \\ \hfill & {n}^{a}=\frac{{\underline{\psi }}^{a}}{\sqrt{\tau \,{\text{e}}^{\varphi }}},\hfill \end{align} \tag{ 95 }$$

$$\begin{align}\hfill & {m}^{a}=\frac{1}{\sqrt{2}}\left({e}_{(2)}^{a}+i{e}_{(3)}^{a}\right),\hfill \\ \hfill & {\bar{m}}^{a}=\frac{1}{\sqrt{2}}\left({e}_{(2)}^{a}-i{e}_{(3)}^{a}\right),\hfill \end{align} \tag{ 96 }$$

where ${e}_{(2)}^{a}$ and ${e}_{(3)}^{a}$ are real unit vectors orthogonal to l^a + n^a and l^a − n^a . From the Weyl tensor C_abcd , one may construct the Weyl scalars Ψ₄, Ψ₃, Ψ₂, Ψ₁ and Ψ₀, which are defined as follows:

$$\begin{align}\hfill & {{\Psi}}_{4}={C}_{\mathit{abcd}}{n}^{a}{\bar{m}}^{b}{n}^{c}{\bar{m}}^{d},\hfill \\ \hfill & {{\Psi}}_{3}={C}_{\mathit{abcd}}{l}^{a}{n}^{b}{\bar{m}}^{c}{n}^{d},\hfill \\ \hfill & {{\Psi}}_{2}={C}_{\mathit{abcd}}{l}^{a}{m}^{b}{\bar{m}}^{c}{n}^{d},\hfill \\ \hfill & {{\Psi}}_{1}={C}_{\mathit{abcd}}{l}^{a}{n}^{b}{l}^{c}{m}^{d},\hfill \\ \hfill & {{\Psi}}_{0}={C}_{\mathit{abcd}}{l}^{a}{m}^{b}{l}^{c}{m}^{d}.\hfill \end{align} \tag{ 97 }$$

We say that the Weyl scalars Ψ_N satisfy the peeling property if, for some appropriately defined radial parameter r, the asymptotic behavior of the Weyl scalars has the form [2]:

$$\begin{equation}{{\Psi}}_{N}\sim \frac{1}{{r}^{5-N}}.\end{equation} \tag{ 98 }$$

It is worth mentioning that in the classical peeling result in [2] the coordinates used are (u, r, θ^A ) where u satisfies the eikonal equation and l^a = g^ab ∇_b u are tangent to a family of null geodesics. Additionally, the freedom in rotating the m and n legs of the tetrad is fixed by demanding that they are parallel propagated along ${l}^{a}={(\frac{\partial }{\partial r})}^{a}$ where r is an affine parameter along the null geodesics. Nonetheless, in other classical studies such as in [4] the areal radial coordinate (also called luminosity parameter) $\tilde{r}$ has been used to perform similar asymptotic expansions of the curvature tensors and the metric. The relation between Newman–Penrose's r and Bondi–Sach's $\tilde{r}$ has been derived in [41]—see also [42] for a discussion of the Bondi–Sachs set up. However, in our case we will not use r nor $\tilde{r}$ but work with the radial-shell coordinate R as defined in section 2. This is a subtle point since, in order to fully compare an asymptotic expansion for the Weyl tensor to that of the classical peeling theorem of [2] one needs not only to find the transformation between the coordinates but also between the frames. The precise relation between our set up (coordinates and frame) and that of [2] will be analyzed elsewhere. We compute the Weyl scalars Ψ_N using the polyhomogeneous expansions for each of the metric functions in the Cartesian harmonic gauge, including terms up to order R⁻³ and logarithmic terms indicated by (92)–(94). We find that Ψ₄ and Ψ₃ satisfy the peeling property (here, the parenthetical subscripts (n, k) denoting the order of coefficients in the polyhomogeneous expansion):

$$\begin{align}\hfill {{\Psi}}_{4}& =-\frac{\left({\ddot{h}}_{+(1,0)}+i{\ddot{h}}_{\times (1,0)}\right)}{R}+o({R}^{-1}),\hfill \\ \hfill {{\Psi}}_{3}& =\frac{1}{4{R}^{2}}\left[-{\partial }_{\theta }{\dot{h}}_{+(1,0)}-i{\partial }_{\theta }{\dot{h}}_{\times (1,0)}+{\partial }_{\theta }{\dot{\varepsilon }}_{(1,0)}-{\partial }_{\theta }{\dot{\varphi }}_{(1,0)}\right.\hfill \\ \hfill & \quad -{\ddot{C}}_{\theta (1,0)}^{+}-i\mathrm{csc}(\theta )\left({\partial }_{\phi }{\dot{h}}_{+(1,0)}+i{\partial }_{\phi }{\dot{h}}_{\times (1,0)}+{\partial }_{\phi }{\dot{\varepsilon }}_{(1,0)}\right.\hfill \\ \hfill & \quad \left.\left.-{\partial }_{\phi }{\dot{\varphi }}_{(1,0)}-{\ddot{C}}_{\phi (1,0)}^{+}\right)-2\,\mathrm{cot}(\theta )\left({\dot{h}}_{+(1,0)}+i{\dot{h}}_{\times (1,0)}\right)\right]+o({R}^{-2}).\hfill \end{align} \tag{ 99 }$$

To compute Ψ₂, it is convenient to introduce the following contractions of the Ricci tensor, which vanish on vacuum solutions:

$$\begin{equation}\mathcal{H}={R}_{ab}({l}^{a}+{n}^{a})({l}^{b}+{n}^{b}),\end{equation} \tag{ 100 }$$

$$\begin{equation}\mathcal{M}={R}_{ab}({l}^{a}+{n}^{a})({l}^{b}-{n}^{b}).\end{equation} \tag{ 101 }$$

Recall that ${R}_{ab}={\mathcal{R}}_{ab}+{\nabla }_{a}{Z}_{b}-{W}_{ab}$, and if the reduced vacuum EFE are imposed—via the wave equations of section 4—one has that ${\mathcal{R}}_{ab}=0$. Thus, requiring $\mathcal{H}=0$ and $\mathcal{M}=0$ is guaranteed by assuming that the Z_a constraints are satisfied. Due to the construction for the propagation of constraints, roughly speaking one can view $\mathcal{H}$ and $\mathcal{M}$ as a subset of the Hamiltonian and momentum constraints. We do not impose that any of the constraints vanish but, inline with the free-evolution philosophy, given a solution to the reduced field equations, we are free to add and subtract any amount of those constraints to the equations when calculating quantities after the fact. Upon adding the following combination of $\mathcal{H}$ and $\mathcal{M}$, we obtain an expression of the form

$$\begin{align}\hfill {{\Psi}}_{2}^{\prime }& ={{\Psi}}_{2}-\frac{\mathcal{H}+\mathcal{M}}{12}\hfill \\ \hfill & =\frac{{{\Psi}}_{2(3,0)}}{12{R}^{3}}+\frac{\mathrm{log}\,R}{4{R}^{3}}{{\Psi}}_{2(3,1)}+o({R}^{-3}),\hfill \end{align} \tag{ 102 }$$

where:

$$\begin{align}\hfill {{\Psi}}_{2(3,0)}& =-3{\mathcal{C}}_{-(1,0)}^{R}+2{\mathcal{C}}_{-(1,1)}^{R}+3{\mathcal{C}}_{+(1,0)}^{R}+2{h}_{+(1,0)}-6{\dot{\varepsilon }}_{(2,0)}+4{\dot{\varepsilon }}_{(2,1)}-{\partial }_{\theta }^{2}{\mathcal{C}}_{+(1,0)}^{R}\hfill \\ \hfill & \quad -{\partial }_{\theta }^{2}{h}_{+(1,0)}+{\partial }_{\theta }^{2}{\varepsilon }_{(1,0)}-2{\partial }_{\theta }^{2}{\varphi }_{(1,0)}+6{\dot{\mathcal{C}}}_{+(2,0)}^{R}+2\,\mathrm{cot}(\theta )\mathrm{csc}(\theta ){\partial }_{\phi }{h}_{\times (1,0)}\hfill \\ \hfill & \quad -3\left({h}_{+(1,0)}-i{h}_{\times (1,0)}\right)\left({\dot{h}}_{+(1,0)}+i{\dot{h}}_{\times (1,0)}\right)-\mathrm{cot}(\theta )\left[{\partial }_{\theta }{\mathcal{C}}_{+(1,0)}^{R}+3{\partial }_{\theta }{h}_{+(1,0)}\right.\hfill \\ \hfill & \quad \left.-{\partial }_{\theta }{\varepsilon }_{(1,0)}+2{\partial }_{\theta }{\varphi }_{(1,0)}\right]+{\mathrm{csc}}^{2}(\theta )\left[{\partial }_{\phi }^{2}{h}_{+(1,0)}-{\partial }_{\phi }^{2}{\mathcal{C}}_{+(1,0)}^{R}+{\partial }_{\phi }^{2}{\varepsilon }_{(1,0)}\right.\hfill \\ \hfill & \quad \left.-2{\partial }_{\phi }^{2}{\varphi }_{(1,0)}+2\,\mathrm{sin}(\theta ){\partial }_{\theta }{\partial }_{\phi }{h}_{\times (1,0)}\right],\hfill \\ \hfill {{\Psi}}_{2(3,1)}& =-{\mathcal{C}}_{-(1,1)}^{R}+2{\dot{\mathcal{C}}}_{+(2,1)}^{R}-2{\dot{\varepsilon }}_{(2,1)}.\hfill \end{align} \tag{ 103 }$$

Since the coefficient ${\mathcal{C}}_{-(1,1)}^{R}$ appears in Ψ_2(3,1), the Weyl scalar will fall off as log R/R³. It follows that for solutions of the asymptotic system in Cartesian harmonic coordinates, and with this particular constraint addition (89), the Weyl scalar Ψ₂ fails to satisfy (98). We can then conclude that the vacuum EFE in GHG do not, in general, satisfy the peeling property (even when working within our restricted class of initial data). However we shall see in the next section that it is possible to pick the gauge choice and constraint addition carefully in order to recover this property.

As previously said, the asymptotic behavior of solutions of GR near null infinity is affected by a special interplay between the gauge and the constraint addition we choose to work with. We have seen in the previous section that our first choice results in a violation of the peeling property of the Weyl scalars, even within our special class of initial data. In this section we will see however that it is possible, with the same type of initial data, to find a gauge and constraint addition in which peeling is manifest. This is done by finding that there is at least one choice which makes 2 of the metric functions behave as goods and 8 as uglies. Then, by remark 1, one can make the first logs appear only at arbitrarily high order, so that the Weyl scalars cannot possibly contain log R to leading order.

Gauge choice: we choose a gauge that is Cartesian harmonic to leading order, and add a higher order correction ${\check{F}}^{a}$,

$$\begin{equation}{F}^{a}={\stackrel{\circ}{F}}^{a}+{\check{F}}^{a}.\end{equation} \tag{ 104 }$$

Before choosing ${\check{F}}^{a}$ we will compute the first order asymptotic system to find the most convenient choice for our purposes.

Constraint addition: for the same reason that led us to keeping ${\check{F}}^{a}$ free, we want to keep some freedom in the constraint addition to choose a posteriori. We parameterize that freedom by the four natural numbers p₁, p₂ and p_A , and we make the following constraint addition,

$$\begin{align}\hfill & {W}_{\psi \psi }=-\frac{1}{2}{({Z}^{\sigma })}^{2}+\frac{{p}_{1}}{R}{Z}^{\sigma },\hfill \\ \hfill & \not{\hskip-3pt W}=-\frac{{p}_{2}}{R}{Z}^{\underline{\sigma }},\hfill \\ \hfill & {W}_{\psi A}=\frac{3-{p}_{A}}{R}{Z}^{A}.\hfill \end{align} \tag{ 105 }$$

GHG in good–bad–ugly form: our choice of gauge, together with constraint addition (105), gives the following system,

$$\begin{align}\hfill & \stackrel{\circ}{\square }\varphi ={\nabla }_{T}{\check{F}}^{\sigma }+{\mathcal{N}}_{\varphi },\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{+}^{R}=\frac{2{p}_{1}}{R}{\nabla }_{T}{\mathcal{C}}_{+}^{R}+{\mathcal{N}}_{{\mathcal{C}}_{+}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{-}^{R}=-\frac{1}{2}{({\nabla }_{T}{h}_{+})}^{2}-\frac{1}{2}{({\nabla }_{T}{h}_{\times })}^{2}-2{\nabla }_{T}{\check{F}}^{\underline{\sigma }}+{\mathcal{N}}_{{\mathcal{C}}_{-}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{+}=\frac{2{p}_{A}}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{+}+{\mathcal{N}}_{{\hat{\mathcal{C}}}_{A}^{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{-}=\frac{4}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{-}-2R{\nabla }_{T}{\check{F}}^{A}+{\mathcal{N}}_{{\hat{\mathcal{C}}}_{A}^{-}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\epsilon}=\frac{2{p}_{2}}{R}{\nabla }_{T}{\epsilon}+{\mathcal{N}}_{{\epsilon}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{+}={\mathcal{N}}_{{h}_{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{\times }={\mathcal{N}}_{{h}_{\times }}.\hfill \end{align} \tag{ 106 }$$

A close look shows that h₊ and h_× satisfy good equations regardless of the freedom we still have. So we are left with eight equations and eight degrees of freedom given by the scalar functions ${\check{F}}^{a}$ and the choice of constraint addition, paramaterized here by the natural numbers p₁, p₂ and p_A , all of which may influence the asymptotic solution space. We want to recover the peeling property of the Weyl scalars, so our criteria in choosing our free functions will be to have as few logs in the first few orders as possible. In fact, there is a particular choice which makes the expansion of the metric functions log-free to arbitrarily high order. As was stated in remark 1, the polyhomogeneous expansion of a field which satisfies an ugly equation of the type (32), provided that no logs are inherited through coupling with other equations, will have no logs up to order p. Therefore our strategy here is to set p₁ = p₂ = p_A = p and choose p large enough to get rid of all logs in the ugly fields that may contribute to the leading decay in the Weyl scalars. However we must guarantee that no logs are inherited through coupling, and we will do that by carefully choosing ${\check{F}}^{a}$ in order to turn φ, ${\mathcal{C}}_{-}^{R}$ and ${\hat{\mathcal{C}}}_{ A}^{-}$ into uglies with no logs to order p. Therefore, we want ${\check{F}}^{a}$ to satisfy the following conditions,

$$\begin{align}\hfill & {\nabla }_{T}{\check{F}}^{\sigma }\simeq \frac{2p}{R}{\nabla }_{T}\varphi ,\hfill \\ \hfill & {\nabla }_{T}{\check{F}}^{\underline{\sigma }}\simeq -{({\nabla }_{T}{h}_{+})}^{2}-{({\nabla }_{T}{h}_{\times })}^{2}-\frac{p}{R}{\nabla }_{T}{\mathcal{C}}_{-}^{R},\hfill \\ \hfill & {\nabla }_{T}{\check{F}}^{A}\simeq -\frac{p-2}{{R}^{2}}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{-},\hfill \end{align} \tag{ 107 }$$

so that the system turns into,

$$\begin{align}\hfill & \stackrel{\circ}{\square }\varphi =\frac{2p}{R}{\nabla }_{T}\varphi +{\mathcal{N}}_{\varphi },\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{+}^{R}=\frac{2p}{R}{\nabla }_{T}{\mathcal{C}}_{+}^{R}+{\mathcal{N}}_{{\mathcal{C}}_{+}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\mathcal{C}}_{-}^{R}=\frac{2p}{R}{\nabla }_{T}{\mathcal{C}}_{-}^{R}+{\mathcal{N}}_{{\mathcal{C}}_{-}^{R}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{+}=\frac{2p}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{+}+{\mathcal{N}}_{{\hat{\mathcal{C}}}_{A}^{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\hat{\mathcal{C}}}_{A}^{-}=\frac{2p}{R}{\nabla }_{T}{\hat{\mathcal{C}}}_{A}^{-}+{\mathcal{N}}_{{\hat{\mathcal{C}}}_{A}^{-}},\hfill \\ \hfill & \stackrel{\circ}{\square }{\epsilon}=\frac{2p}{R}{\nabla }_{T}{\epsilon}+{\mathcal{N}}_{{\epsilon}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{+}={\mathcal{N}}_{{h}_{+}},\hfill \\ \hfill & \stackrel{\circ}{\square }{h}_{\times }={\mathcal{N}}_{{h}_{\times }}.\hfill \end{align} \tag{ 108 }$$

With this we would be able to write the 10 metric functions as 2 goods and 8 uglies. Moreover, p has not yet been chosen so by remark 1, we are allowed to have logs only from an arbitrarily high order onward. Note that (107) is only an asymptotic condition, meaning that we might be interfering with terms that contribute to sub-leading orders. However, our goal here is to guarantee that our equations are only either good or ugly, and that is determined solely by the first order asymptotic system. In the following we explain how (107) can be achieved.

Making sure hyperbolicity is maintained: there is a subtlety in guaranteeing that the gauge source functions satisfy (107). In order to write ${\check{F}}^{a}$ in terms of metric functions, we can simply integrate the first and last equations to get,

$$\begin{equation}{\check{F}}^{\sigma }=\frac{2p}{R}({\text{e}}^{\varphi }-1),\qquad {\check{F}}^{A}=\frac{2-p}{{R}^{2}}{\hat{\mathcal{C}}}_{A}^{-}.\end{equation} \tag{ 109 }$$

This simple-minded approach to the second equation fails, as it does not integrate cleanly, and instead yields an expression for ${\check{F}}^{\underline{\sigma }}$ involving an integral of derivatives of h₊ and h_×. To overcome the issue we take inspiration from the gauge driver approach of [37]. We begin by separating ${\check{F}}^{\underline{\sigma }}$ into two pieces, an explicit one and a gauge driver ${\check{F}}_{1}^{\underline{\sigma }}$, as follows:

$$\begin{equation}{\check{F}}^{\underline{\sigma }}=\frac{1}{R}{\check{F}}_{1}^{\underline{\sigma }}-\frac{p}{R}(1+{\mathcal{C}}_{-}^{R}),\end{equation} \tag{ 110 }$$

so that, according to (107), ${\check{F}}_{1}^{\underline{\sigma }}$ must now satisfy,

$$\begin{equation}\frac{1}{R}{\nabla }_{T}{\check{F}}_{1}^{\underline{\sigma }}\simeq -{({\nabla }_{T}{h}_{+})}^{2}-{({\nabla }_{T}{h}_{\times })}^{2}{:=}\frac{1}{{R}^{2}}H.\end{equation} \tag{ 111 }$$

Since h₊ and h_× are goods, we can write,

$$\begin{equation}H\simeq {H}_{1,0}({\psi }^{\ast }).\end{equation} \tag{ 112 }$$

Our strategy consists of insisting that ${\check{F}}_{1}^{\underline{\sigma }}$ satisfies a wave equation whose solution has the property that (111). Then, second derivatives of ${\check{F}}_{1}^{\underline{\sigma }}$ in wave equations for metric functions can be treated as a variable to be evolved instead of as functions of the metric themselves. We impose the following gauge-driving evolution equation

$$\begin{equation}\stackrel{\circ}{\square }{\check{F}}_{1}^{\underline{\sigma }}=\frac{2q}{R}\left({\nabla }_{T}{\check{F}}_{1}^{\underline{\sigma }}-\frac{1}{R}H\right)+{\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}},\end{equation} \tag{ 113 }$$

where q is a natural number deliberately left free for the moment. The reason for the specific form of this equation will become clear later in this section, but morally one can already think of the first term on the rhs as suppressing the part of the radiation field of ${\check{F}}_{1}^{\underline{\sigma }}$ that would naturally appear from initial data, and the second term forcing the radiation field equal to a desired target. We included a combination of stratified null forms ${\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}$ because, as they do not contribute to leading order, they can be chosen in order to make the wave equation as simple as possible. We shall make that choice later. Rescaling our problematic function as ${f}_{1}{:=}R{\check{F}}_{1}^{\underline{\sigma }}$ then gives,

$$\begin{equation}-\frac{2}{R}{\nabla }_{\psi }{\nabla }_{T}{f}_{1}=\frac{2q}{{R}^{2}}\left({\nabla }_{T}{f}_{1}-H\right)+{\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}^{\prime },\end{equation} \tag{ 114 }$$

where ${\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}^{\prime }$ is just ${\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}$ minus the stratified null forms coming from the lhs. We now add 2R⁻¹∇_ψ H on both sides,

$$\begin{align}\hfill -\frac{2}{R}{\nabla }_{\psi }\left({\nabla }_{T}{f}_{1}-H\right)& =\frac{2q}{{R}^{2}}\left({\nabla }_{T}{f}_{1}-H\right)+\frac{2}{R}{\nabla }_{\psi }H+{\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}^{\prime }.\hfill \end{align} \tag{ 115 }$$

As we left ${\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}$ free, we can still choose it so that $\frac{2}{R}{\nabla }_{\psi }H=-{\mathcal{N}}_{{\check{F}}_{1}^{\underline{\sigma }}}^{\prime }$, since ∇_ψ H is necessarily a stratified null form (112). This is not an essential step, but it makes the resulting equation more tractable. Therefore we get the simple equation,

$$\begin{equation}{\nabla }_{\psi }\left({\nabla }_{T}{f}_{1}-H\right)=-\frac{q}{R}\left({\nabla }_{T}{f}_{1}-H\right),\end{equation} \tag{ 116 }$$

that we can integrate to find,

$$\begin{equation}{\nabla }_{T}{f}_{1}-H=\frac{1}{{R}^{q}}\alpha ({\psi }^{\ast }).\end{equation} \tag{ 117 }$$

From (117) we see that as long as q is a natural number, the condition (111) is satisfied, and we can turn ${\mathcal{C}}_{-}^{R}$ into an ugly while ensuring that hyperbolicity is maintained. However, by introducing an eleventh non-linear wave equation, we may be introducing logs in the system, which will then appear in the EFE through coupling. In fact, the wave equation we forced ${\check{F}}_{1}^{\underline{\sigma }}$ to satisfy, namely (113) does not fit directly into our conception of good, bad or ugly. For that reason, we have to show that ${\check{F}}_{1}^{\underline{\sigma }}$ itself admits a polyhomogeneous expansion in which logs only appear at arbitrarily high order. We rescale ${f}_{1}{:=}R{\check{F}}_{1}^{\underline{\sigma }}$ and assume that f₁ satisfies the conditions (23), so we can write,

$$\begin{equation}{\nabla }_{\psi }\left({R}^{q}{\nabla }_{T}{f}_{1}\right)\simeq q{R}^{q-1}{H}_{1,0},\end{equation} \tag{ 118 }$$

Integrating we get,

$$\begin{equation}{f}_{1}\simeq {F}_{1,0}({\psi }^{\ast }),\end{equation} \tag{ 119 }$$

where F_1,0(ψ*) is constant along ψ^a . Now we assume that,

$$\begin{equation}{\check{F}}_{1}^{\underline{\sigma }}=\frac{{F}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{m=2}^{n-1}\sum\limits _{k=0}^{{N}_{n-1}^{F}}\frac{{(\mathrm{log}\,R)}^{k}{F}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{F}}_{n}}{{R}^{n}},\end{equation} \tag{ 120 }$$

where ${\mathcal{F}}_{n}={o}^{+}(R)$, and must demonstrate,

$$\begin{equation}{\check{F}}_{1}^{\underline{\sigma }}=\frac{{F}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{m=2}^{n}\sum\limits _{k=0}^{{N}_{n}^{F}}\frac{{(\mathrm{log}\,R)}^{k}{F}_{m,k}({\psi }^{\ast })}{{R}^{m}}+\frac{{\mathcal{F}}_{n+1}}{{R}^{n+1}}.\end{equation} \tag{ 121 }$$

Plugging (120) into (113) and collecting terms proportional to R⁻ⁿ⁻¹ we get,

$$\begin{equation}{\nabla }_{\psi }({R}^{q+1-n}{\nabla }_{T}{\mathcal{F}}_{n})\simeq {R}^{q-n}{{\Omega}}_{n-1}^{F},\end{equation} \tag{ 122 }$$

where,

$$\begin{equation}{{\Omega}}_{n-1}^{F}=\sum\limits _{i=0}^{{N}_{n}^{{{\Omega}}^{F}}}{(\mathrm{log}\,R)}^{i}{{\Omega}}_{n-1,i}^{F}({\psi }^{\ast }).\end{equation} \tag{ 123 }$$

Integrating (122) we find that (121) and so, by induction, ${\check{F}}_{1}^{\underline{\sigma }}$ admits a polyhomogeneous expansion near null infinity of the type,

$$\begin{equation}{\check{F}}_{1}^{\underline{\sigma }}=\frac{{F}_{1,0}({\psi }^{\ast })}{R}+\sum\limits _{m=2}^{\infty }\sum\limits _{k=0}^{{N}_{n}^{F}}\frac{{(\mathrm{log}\,R)}^{k}{F}_{m,k}({\psi }^{\ast })}{{R}^{m}}.\end{equation} \tag{ 124 }$$

To this point we have allowed ${\check{F}}_{1}^{\underline{\sigma }}$ to have any power of log R in terms beyond first order. However, an argument completely analogous to remark 1 can be made here in order to show that there can be no logs up to order n = q as long as no logs are inherited through coupling. With this in mind, we choose q = p and notice that we thus avoid the danger of introducing logs in the first orders of the asymptotic expansions. By insisting that ${\check{F}}_{1}^{\underline{\sigma }}$ satisfies the wave equation (113), we have prevented our special choice of gauge from interfering with hyperbolicity. Finally, we have effectively shown that a choice of gauge and constraint addition allows us to have 2 goods and 8 uglies with no logs up to arbitrarily high order. We conclude this section with a paragraph on how this result recovers the peeling property.

Recovering peeling: as we saw in section 6, the peeling property is satisfied iff,

$$\begin{equation}{{\Psi}}_{N}\sim \frac{1}{{R}^{5-N}},\end{equation} \tag{ 125 }$$

where N ∈ {0, 1, 2, 3, 4}. We want to use our freedom in the choice of the natural number p to make sure that logs only appear in the polyhomogeneous expansion at an order that is high enough that they cannot possibly appear in the leading orders of the Weyl scalars. These are components of the Weyl tensor and hence they are non-linear combinations of the metric and its first and second derivatives. Since no derivatives can ever decrease the order of the given term in R⁻¹, we need only know the functional form of the metric components to order R⁻⁵ in order to determine whether the peeling property is satisfied. We need to make sure that the metric components do not contain any logs to order n = 5, or rather,

$$\begin{equation}{g}_{ab}=\sum\limits _{n=1}^{5}\frac{{G}_{ab,n}({\psi }^{\ast })}{{R}^{n}}+\frac{{\mathcal{G}}_{6}}{{R}^{6}},\end{equation} \tag{ 126 }$$

where ${\mathcal{G}}_{6}={o}^{+}(R)$. Now in order to choose p, we must take into account that certain metric components contain factors of R and R² attached to metric functions that are themselves good, bad or ugly, as can be seen from (10) and (28). That said, without actually computing the dependence of the Weyl scalar on the metric functions, we can guarantee that they will not contain any logs to leading order if we require that the metric functions do not contain any logs to order n = 7. This is achieved by choosing p ⩾ 7, in accordance with remark 1. So, by carefully choosing the gauge and constraint addition, we have prevented the existence of any terms proportional to log R in the polyhomogeneous expansions of all the metric functions up to order n = 7, thereby recovering, within our class of initial data, the peeling property of the gravitational field in the EFE-GHG setup.

In the previous sections we saw that the Cartesian harmonic gauge, together with the constraint addition that turns 4 of the metric functions into uglies, generates terms proportional to log R at first order in ${\mathcal{C}}_{-}^{R}$ and at second order in some of the other functions. These logs in turn, give rise to a violation of peeling. Then it was shown that by changing the gauge to a very specific one and modifying the constraint additions slightly to contain the natural number p, we can delay the appearance of those logs to an arbitrarily high order. This ensures that peeling is recovered if p is high enough. Here we will not be concerned with peeling, but with how these logs are generated. So for simplicity we fix p = 1 so that the constraint additions in either case are exactly the same, and focus on the log that appears at leading order, i.e. the one in ${\mathcal{C}}_{-}^{R}$. In that case, the EFE only differ by the choice of gauge, so we should, in principle, be able to generate and get rid of the logs just by a particular coordinate change. To illustrate this, we begin with a coordinate system {T, R, θ^A } which we assume generates no logs to leading order. Then we change one coordinate, say R, to something different, {T, ρ, θ^A }, and check that the leading log in ${\mathcal{C}}_{-}^{R}$ is recovered. The relation between the new ${\mathcal{C}}_{-}^{\rho }$ and the old ${\mathcal{C}}_{-}^{R}$ is given by,

$$\begin{equation}{\mathcal{C}}_{-}^{\rho }={\mathcal{C}}_{-}^{R}{{J}^{\rho }}_{R}+{{J}^{\rho }}_{T}={\nabla }_{\underline{\psi }}\rho ,\end{equation} \tag{ 127 }$$

which may be obtained by comparing expressions for the covector basis σ_a , ${\underline{\sigma }}_{a}$ between the two sets of coordinates. The wave equation that ρ satisfies is,

$$\begin{equation}\square \rho -{F}^{\rho }=0,\end{equation} \tag{ 128 }$$

and the most straightforward way to find the asymptotics of ρ through (128) is by taking advantage of the analysis we already have of the good, bad and ugly equations. However, the variable ρ is expected to grow like R, so it is not exactly a variable like the g, b or u. So we choose to work with ρ/R instead and change (128) accordingly,

$$\begin{equation}R\square \frac{\rho }{R}=-2{g}^{aR}{\nabla }_{a}\frac{\rho }{R}+\frac{\rho }{R}{\stackrel{\bullet }{{\Gamma}}}^{R}+{F}^{\rho }.\end{equation} \tag{ 129 }$$

This can be written as,

$$\begin{equation}\stackrel{\circ}{\square }\frac{\rho }{R}=\frac{2}{R}{\nabla }_{T}\frac{\rho }{R}+{\mathcal{N}}_{\frac{\rho }{R}},\end{equation} \tag{ 130 }$$

which implies that $\frac{\rho }{R}$ is an ugly and admits a polyhomogeneous expansion of the form,

$$\begin{equation}\frac{\rho }{R}=1+\frac{m({\theta }^{A})}{R}+{o}^{+}({R}^{-1}),\end{equation} \tag{ 131 }$$

which, for convenience, we will write as,

$$\begin{equation}\rho =R+m({\theta }^{A})+\frac{{\Phi}}{R},\end{equation} \tag{ 132 }$$

by defining the quantity Φ accordingly. If we plug this in (128) and use constraint addition, we get,

$$\begin{equation}\square \frac{{\Phi}}{R}={F}^{\rho }-{F}^{R}-{/\hskip-3pt \!\!{\Delta}}m={o}^{+}({R}^{-1}),\end{equation} \tag{ 133 }$$

which implies,

$$\begin{equation}-\frac{1}{R}{\nabla }_{\psi }{\nabla }_{\underline{\psi }}{\Phi}={F}^{\rho }-{F}^{R}-{/\hskip-3pt \!\!{\Delta}}m+{o}^{+}({R}^{-2}).\end{equation} \tag{ 134 }$$

Upon integration along integral curves of ψ^a , the only way for ${\nabla }_{\underline{\psi }}{\Phi}$ not to have a log is if,

$$\begin{equation}{F}^{\rho }-{F}^{R}-{/\hskip-3pt \!\!{\Delta}}m={o}^{+}({R}^{-2}),\end{equation} \tag{ 135 }$$

which, in general, is not true. In fact, one of the second order terms in F^R is,

$$\begin{equation}\frac{1}{R}{\check{F}}_{1}^{\underline{\sigma }}\simeq -{\int }_{{T}_{0}}^{T}\left[{({\nabla }_{T}{h}_{+})}^{2}+{({\nabla }_{T}{h}_{\times })}^{2}\right]\mathrm{d}{T}^{\prime },\end{equation} \tag{ 136 }$$

precisely the term that gives rise to the log at first order in R⁻¹ in the polyhomogeneous expansion of ${\mathcal{C}}_{-}^{R}$ in the Cartesian harmonic gauge. So we see that if ${F}^{\rho }-{/\hskip-3pt {\Delta}}m$, both determined by our choice of coordinates, do not contain a term symmetric to (136), a log is unavoidable. Therefore, in general,

$$\begin{equation*}\rho =R+m({\theta }^{A})+\frac{{{\Phi}}_{1,0}({\psi }^{\ast })+{{\Phi}}_{1,1}({\psi }^{\ast })\mathrm{log}\,R}{R}+{o}^{+}({R}^{-1}),\end{equation*}$$

which, together with (127) gives,

$$\begin{equation}{\mathcal{C}}_{-}^{\rho }={\mathcal{C}}_{-}^{R}+\frac{{\nabla }_{\underline{\psi }}{{\Phi}}_{1,0}+{\nabla }_{\underline{\psi }}{{\Phi}}_{1,1}\mathrm{log}\,R}{R}+{o}^{+}({R}^{-1}),\end{equation} \tag{ 137 }$$

thus recovering a log R to first order in R⁻¹ as we wanted.

Consider initial data that is, in some sense, asymptotically flat and that generates a spacetime with a piece of future null infinity. Introduce coordinates and a null-frame appropriately adapted to asymptotic flatness in a neighborhood of null infinity. If the Weyl scalars decay like Ψ_N = O(r^N−5), with r a radial coordinate naturally constructed from the adapted coordinates, then the Weyl tensor is said to peel. Whether or not the Weyl tensor peels is intimately tied both to the specific notion of asymptotic flatness imposed on initial data, and to the coordinates and frame employed. Described as such, the notion of peeling appears rather pedestrian. It could be that even within a strict class of initial data one could conclude that the Weyl tensor both peels and does not peel, depending on the choice of coordinates and frame. A refined formulation, for which we have much sympathy, would be to say that the Weyl tensor peels if any set of coordinates and frame can be found in which the aforementioned decay is manifest. In practice it is often the more pedestrian notion of peeling that is of concern, since that is directly related to smoothness at null-infinity within the gauge employed. Besides which, demonstrating that no coordinates exist in which peeling holds is a much less concrete task than examining specific cases of interest.

A fundamental ingredient in the small-data global existence result in harmonic gauge [43] was to understand and control terms in the metric components of the asymptotic form log R/R near null infinity. Since such terms decay more slowly than solutions to the wave equation and appear in the presence of any gravitational wave content, the questions arise: could the peeling property generically be satisfied when using harmonic gauge, even with a very strict notion of asymptotic flatness? Assuming that this is not the case, could peeling be recovered for the same type of initial data by using a suitable GHG?

The main achievement of this work was to address these questions by exploiting the formal asymptotic expansions of [36]. We found that there are choices of gauge and constraint addition which give rise to a violation of the peeling property. For this it was necessary to generalize the last theorem of [36] to include systems of any number of good, bad and ugly equations with stratified null forms, as well as to generalize the notion of 'ugly'. This was to establish that such a system admits formal polyhomogeneous asymptotic solutions near null infinity for a broad class of initial data. Since in general the metric components ${g}_{\underline{\mu \nu }}$ with respect to the GHG coordinate basis ${\partial }_{\underline{\mu }}$ does not fall into our good–bad–ugly categorization, these methods cannot be directly and cleanly applied to the standard formulation of the EFE in GHG (in which ${g}_{\underline{\mu \nu }}$ are the variables). However, exploiting that this categorization and the first order system discussed in [21, 24] are based on the privileged role played by null directions, we have introduced a reformulation of the EFE in GHG where the variables are written in terms of the coordinate lightspeeds (encoding the non-trivial components of the outgoing and incoming null vectors). One point to notice here is that in [21, 24] the null directions are those of the Minkowski spacetime while in the formulation discussed here, the null directions are the (true) null directions in $(\mathcal{M},g)$. This special set of variables has the property that, the evolution equations in the EFE in GHG, upon a choice of gauge and constraint addition, have the structure of a good–bad–ugly system.

This was carefully seen in section 5 by considering the Cartesian harmonic gauge and then adding constraints in order to maximize the number of uglies. Analyzing the asymptotic system via computer algebra, it was possible to use the theorem shown in section 3 to find the maximum power of log R that each variable may have at each order up to 3rd in R⁻¹. Computing the leading order terms in the Weyl scalars, we were able to assert that, at least within the formal expansions framework of [36] the peeling property is not, in general, satisfied by solutions in this gauge, as these contain leading terms proportional to log R, instead of only negative integer powers of R.

Hence, physically relevant spacetimes evolved (say numerically using the EFE in GHG), will not, in general, possess a smooth null infinity in that gauge. We saw, however that it is possible to choose the gauge and constraint addition carefully to guarantee that terms proportional to log R do not appear in the few first orders of the metric components, hence recovering peeling within the same class of initial data. Interestingly, in contrast to the order-by-order approach used for harmonic gauge in the post-Minkowskian setting [44] the generalized harmonic strategy we apply suppresses log R terms down to a finite order all at once, and may therefore be of use in that context.

It is furthermore interesting that the violation of peeling identified in this paper occurs for Ψ₂ and is of the form log R/R³. This formally coincides with the violation of peeling reported in [45]. But unlike the present work, the result of [45] was given in the Newman–Penrose gauge [2]. Establishing a relation between the gauges used in this paper and other gauges such as the Newman–Penrose gauge or the F-gauge of [46] is of deep interest but requires a delicate analysis that will be done elsewhere. In the framework of the formal expansions of solutions to the CEFEs discussed in [45], see also [46], the violation of peeling can be related to a regularity condition at the level of the initial data in terms of derivatives of the Bach tensor of the initial conformal metric evaluated at spatial infinity. The main focus of this paper has been on the evolution equations and the choice of gauge. However, the relation between the class of initial data, solutions to the Hamiltonian and momentum constraints, and the terms appearing in our asymptotic expansions, and ultimately in those responsible for the violations of peeling, is of particular interest and deserves a full study in itself, so is left for future work. It should be stressed that the asymptotic expansions derived here as well as those reported in [45] are formal and the relation between the expansions and actual solutions is yet to be established. Nonetheless, the proof of the global stability of the Minkowski spacetime of [47] predicts a decay for the Weyl scalars which, despite not being polyhomogeneous, does not follow the fall-off predicted in the classical peeling theorem [2], and hence shows the failure of peeling for generic initial data. The classical peeling result in the setup of [47] can be retrieved by restricting the initial data as shown in [48]. See also [35] for a proof of the nonlinear stability of the Minkowski spacetime with polyhomogeneous initial data leading to polyhomogeneous spacetime metric developments. Perhaps one interpretation of the results here is that the violation of peeling induced by the log R/R terms in the metric are, in a sense, pure gauge. The expectation is that within a larger class of initial data, peeling will fail in the more refined sense, that is, in such a way that it cannot be restored by a simple change of gauge.

In this paper, the expected decay of the Weyl scalars according to the classical peeling theorem was recovered in section 7 by finding a particular choice that turns our system of 10 non-linear wave equations into 2 goods and 8 uglies with no logs up to arbitrarily high order. In the process, we have built a method to find the functional form of polyhomogeneous asymptotic solutions to the EFE in GHG near null infinity. In other words, we now know the powers of log R that may appear in the metric at each order in R⁻¹. This knowledge may allow us to explicitly regularize formally singular terms appearing in a hyperboloidal initial value formulation of the EFE in GHG [19–21, 23] by a non-linear change of variables like the one in [22]. We expect this to be an important step towards a full regularization of the EFE in GHG in the hyperboloidal setup, and that this might result. The full details of this regularization will be presented elsewhere.

The authors wish to thank Alex Vañó-Viñuales for helpful discussions. Many of our derivations were performed in xAct [38] for Mathematica. The notebooks are available at [49]. MD acknowledges support from FCT (Portugal) program PD/BD/135511/2018, DH acknowledges support from the FCT (Portugal) IF Program IF/00577/2015, PTDC/MAT- APL/30043/2017. JF acknowledges support from FCT (Portugal) programs PTDC/MAT-APL/30043/2017, UIDB/00099/2020. EG gratefully acknowledges support from the FCT (Portugal) 2020.03845.CEECIND and from the European Union's H2020 ERC Consolidator Grant 'Matter and Strong-Field Gravity: New Frontiers in Einstein's Theory', Grant Agreement No. MaGRaTh-646597 and the PO FEDER-FSE Bourgogne 2014/2020 program/EIPHI Graduate School (Contract ANR-17-EURE-0002) as part of the ISA 2019 project.

The data that support the findings of this study are openly available at the following URL/DOI: https://centra.tecnico.ulisboa.pt/∼hilditch/Peeling_GHG.