Generalised Boundary Terms for Higher Derivative Theories of Gravity

In this paper we wish to find the corresponding Gibbons-Hawking-York term for the most general quadratic in curvature gravity by using Coframe slicing within the Arnowitt-Deser-Misner (ADM) decomposition of spacetime in four dimensions. In order to make sure that the higher derivative gravity is ghost and tachyon free at a perturbative level, one requires infinite covariant derivatives, which yields a generalised covariant infinite derivative theory of gravity. We will be exploring the boundary term for such a covariant infinite derivative theory of gravity.


Introduction
Einstein's General theory of Relativity (GR) has seen tremendous success in matching its predictions with observations in the weak field regime in the infrared (IR) [1], including the recent confirmation of the detection of Gravitational Waves [2]. However the same cannot be said for the ultraviolet (UV), because at small time scales and short distances, the theory breaks down. At a classical level, GR permits a solution with a blackhole singularity and a cosmological singularity [3]. Besides these classical singularities, the GR also begs for a quantum UV completion. For instance, a pure quadratic gravity is renormalizable [4], but being a higher derivative theory, the theory suffers from ghosts, which can not be cured order by order. In a classical context, higher derivatives also bring in an instability known as the Ostrogradsky instability [5]. In any case, gravity being a diffeomorphism invariant theory, one would naturally expect higher order quantum corrections which would allow all possible diffeomorphism invariant terms, such as covariant higher and infinite derivative contributions in the Ricci scalar, Ricci tensor and Weyl, an explicit computation has been performed in [6]. The covariant construction of infinite derivative theories of gravity (IDG), which is parity invariant and torsion-free, was first studied around the Minkowski background in [6], and then generalised to constant curvature backgrounds, such as deSitter (dS) and Anti-deSitter (AdS) backgrounds in [7].
Since these theories contain infinite covariant derivatives, there is no highest momentum operator and as a result their perturbative stability cannot be analysed via Ostrogradsky analysis. One has to understand their tree level propagator and study what are the true dynamical degrees of freedom which propagate in the space time. If one maintains that the dynamics of the original theory of gravity is governed by massless gravitons, then IDG must not contain any other degrees of freedom, other than spin-2 and spin-0 components of the massless graviton [8,9]. In order to prohibit no new states in the propagator, one has to make sure that there are no zeroes in the complex plane, such that both spin-2 and spin-0 components remain massless. IDG will inevitably introduce new poles in the propagator. However, infinite poles corresponding to infinite covariant derivatives can be summed up as the exponent of an entire function [6,7,[10][11][12][13]. By definition, such a modification will not incur any new poles in the propagator and also ensures that the modified propagator has a correct IR limit, where one recovers the predictions of pure Einstein-Hilbert action at large distances from the source and large time scales. This is indeed encouraging as there exists a non-singular blackhole solution, at least at the level of linearized equations of motion for such covariant infinite derivative theory of gravity in a static [6,12,14], and in a time dependent background [15]. In both cases the universal feature is that at short distances there exists no singularity and at large distances the theory behaves like in the case of GR.
These classical results transcend through the quantum domain. In the UV, such modifications soften the UV divergences by weakening the graviton propagator [10,[16][17][18][19]. Since gravitational interactions contain only derivatives, therefore the vertex operators in such theories are exponentially enhanced, which modifies the counting of superficial degree of divergence, and also gives rise to nonlocal interactions, determined by the new scale M ≤ M p ∼ 2.4 × 10 18 GeV in four dimensions. Indeed, with a similar procedure one can now go to any arbitrary dimensions. As such, there is no restriction on dimensionality from perturbative unitarity. For an exponential suppressed graviton propagator, explicit 2-loop computations have been performed for a scalar analogue of a graviton, and it was found that the 2-loop computation yields a UV finite result [20], along with quantum scattering, where the scattering amplitude does not grow with external momenta [21].
Note that such nonlocality is a common thread for many approaches to quantum gravity, such as loop quantum gravity (LQG) [22,23], causal dynamical approach [24], and string theory (ST) [25]. In LQG and causal dynamical approaches the basic formulation is based on nonlocal objects, such as Wilson loops and fluxes coming from the gravitational field. ST introduces nonlocal interactions, where classically strings and branes interact over a region in space. In string field theory (SFT), the appearance of noncommutative geometry [26], p-adic strings [27], zeta strings [28], and strings quantized on a random lattice [29,30] introduce nonlocality. For a review on SFT, see [31]. An important common feature in all these approaches is the presence of an infinite series of higher-derivative terms incorporating the nonlocality in the form of an exponential kinetic correction.
It is also worth mentioning that in recent years there has been a growing interest in infinite-derivative gravitational theories in not only addressing the Big Bang singularity problem [12,[32][33][34][35], but also in finding cosmological inflation and UV completeness of Starobinsky model of inflation [12,36,37], and [38]. In such classes of theory, the gravitational entropy [39,40], of a static and axisymmetric metric receives zero contribution from the infinite derivative sector of the action, when no additional scalar propagating modes are introduced. In this case, the gravitational entropy is strictly given by the area-law arising solely from the contribution of the Einstein-Hilbert action [41].
Irrespective of classical or quantum computations, one of the key features of a covariant action is to have a well-posed boundary condition. In particular, in the Euclidean path integral approach -requiring such an action to be stationary, one also requires all the boundary terms to disappear on any permitted variation. For instance, calculating the black hole entropy using the Euclidean semiclassical approach shows that the entire contribution comes from the boundary term [43][44][45]. It is well known that the variation of the Einstein-Hilbert (EH) action leads to a boundary term that depends not just on the metric, but also on the derivatives of the metric. This is due to the fact that the action itself depends on the metric, along with terms that depend linearly on the second derivatives. Normally, in Lagrangian field theory, such linear second derivative terms can be introduced or eliminated, by adding an appropriate boundary term to the action. In gravity, the fact that the second derivatives arise linearly and also the existence of total derivative indicates that the second derivatives are redundant in the sense that they can be eliminated by integrating by parts, or by adding an appropriate boundary term. Indeed, writing a boundary term for a gravitational action schematically confines the non-covariant terms to the boundary [46]. For the EH action this geometrically transparent, boundary term is given by the Gibbons-Hawking-York (GHY) boundary term [47]. Adding this boundary to the bulk action results in an elimination of the total derivative, as seen for f (R) gravity [48,49].
In the Hamiltonian formalism, obtaining the boundary terms for a gravitational action is vital. This is due to the fact that the boundary term ensures that the path integral for quantum gravity admits correct answers. As a result, in the late 1950s the 3+1 decomposition received a great deal of attention; Richard Arnowitt, Stanley Deser and Charles W. Misner (ADM) showed [50] that upon decomposing spacetime such that for the four dimensional Einstein equation we have three-dimensional surfaces (later to be defined as hypersurfaces) and one fixed time coordinate for each slices. We can therefore formulate and recast the Einstein equations in terms of the Hamiltonian and hence achieve a better insight into GR.
In the ADM decomposition, one foliates the arbitrary region M of the space-time manifold with a family of spacelike hypersurfaces Σ t , one for each instant in time. It has been shown by the authors of [51] that one can decompose a gravitational action, using the ADM formalism and without necessarily moving into the Hamiltonian regime, such that we obtain the total derivative of the gravitational action. Using this powerful technique, one can eliminate this total derivative term by modifying the GHY term appropriately.
The aim of this paper is to find the corresponding GHY boundary term for a covariant IDG. We start by providing a warm up example of how to obtain a boundary term for an infinite derivative, massless scalar field theory. We then in section 3 briefly review the boundary term for EH term and introduce infinite derivative gravity. We then set our preliminaries by discussing the time slicing in section 4, and reviewing how one may obtain the boundary terms in the 3 + 1 formalism in 5. We finally turn our attention to our gravitational action and find the appropriate boundary terms for such a theory in section 6.
2 Warm up exercise: Infinite derivative massless scalar field theory Let us consider the following action of a generic scalar field φ of mass dimension 2: where η µν is the Minkowski metric 1 and n ∈ N >0 . Generalising, we have that n = n i=1 η µ i ν i ∇ µ i ∇ ν i . The aim is to find the total derivative term for the above action. We may vary the scalar field φ as: φ → φ + δφ. Then the variation of the action is given by The term comes with a scale /M 2 , where M is a new scale below M ≤ Mp = (16πG) −1/2 ∼ 2.4 × 10 18 GeV in 4 dimensions. The physical significance of M could be any scale beyond 10 −2 eV, which arises from constraints on studying the 1/r-fall of the Newtonian potential [14]. In our notation, we suppress the scale M in order not to clutter our formulae for the rest of this paper. However, for any physical comparison one has to bring in the scale M along with Mp.
where now X are the 2n total derivatives: where "· · · " in the above equation indicates the intermediate terms.
Let us now consider a more general case where F( ) = ∞ n=0 c n n , where the 'c n 's are dimensionless coefficients. In this case the total derivatives are given by where the superscript ∇ (j) indicate the number of covariant derivatives acting to the right. Therefore, one can always determine the total derivative for any given action, and one can then preserve or eliminate these terms depending on the purpose of the study. In the following sections we wish to address how one can obtain the total derivative for a given gravitational action.

Introducing Infinite Derivative Gravity
The gravitational action is built up of two main components, the bulk part and the boundary part. In the simplest and the most well known case [47], for the Einstein-Hilbert (EH) action, the boundary term are the ones known as Gibbons-Hawking-York (GHY) term. We can write the total EH action in terms of the bulk part and the boundary part simply as where R is the Ricci-scalar, and K is the trace of the extrinsic curvature with K ij ≡ −∇ i n j , M indicates the 4-dimensional region and ∂M denotes the 3-dimensional boundary region. h is the determinant of the induced metric on the hypersurface ∂M and ε = n µ n µ = ±1, where ε is equal to −1 for a spacelike hypersurface, and is equal to +1 for a timelike hypersurface when we take the metric signature is "mostly plus"; i.e. (−, +, +, +). A unit normal n µ can be introduced only if the hypersurface is not null, and n µ is the normal vector to the hypersurface. Indeed, one can derive the boundary term simply by using the variational principle. In this case the action is varied with respect to the metric, and it produces a total-divergent term, which can be eliminated by the variation of S B , [47]. Finding the boundary terms for any action is an indication that the variation principle for the given theory is well posed.
As mentioned earlier on, despite the many successes that the EH action brought in understanding the universe in IR regime, the UV sector of gravity requires corrections to be well behaved. The most general covariant action of gravity, which is quadratic in curvature, can be written as [13], where α is a constant with mass dimension −2 and the 'f in 's are dimensionless coefficients.
For the full equations of motion of such an action, see [52]. Around Minkowski spacetime, the ghost free condition gives a constraint on the form factors F i ( )'s [6,32,52], Around Minkowski spacetime the Weyl contribution vanishes, since the last term in the action can be recast in terms of Weyl, one can take F 3 ( ) = 0, in which case the condition for a ghost-free graviton propagator leads to a particular choice of the form factor [6,12]: where a( ) is an exponential of an entire function, which does not contain any zeroes. A particularly simple class which mimics the stringy gaussian nonlocalities is given by [6,12,31,32,52] a( ) = e − M 2 , where M is the scale of nonlocality. The aim of this paper is to seek the the boundary terms corresponding to S U V , while retaining the Riemann term. A technically challenging question, but the answer will have tremendous impact on various aspects of gauge theory and gravity which should be explored in future.

Time Slicing
Any geometric spacetime can be recast in terms of time like spatial slices, known as hypersurfaces. How these slices are embedded in spacetime, determines the extrinsic curvature of the slices. One of the motivations of time slicing is to evolve the equations of motion from a well-defined set of initial conditions set at a well-defined spacelike hypersurface, see [53,54].

ADM Decomposition
In order to define the decomposition, we first look at the foliation. Suppose that the time orientable spacetime M is foliated by a family of spacelike hypersurfaces Σ t , on which time is a fixed constant t = x 0 . We then define the induced metric on the hypersurface as where the Latin indices run from 1 to 3 2 . The line element is then given by the ADM (Arnowitt-Deser-Misner) decomposition [55]: is the "lapse" function, and g 00 is the "shift" vector. We may then define n µ , the vector normal to the hypersurface, as: where for the ADM metric, given in Eq. (4.1), n µ takes the following form: In the above line element Eq. (4.1), we also have The induced metric of the hypersurface can be related to the 4 dimensional full metric via the completeness relation, where, for a spacelike hypersurface, where ε = −1 for a spacelike hypersurface, and +1 for a timelike hypersurface, and are basis vectors on the hypersurface which allow us to define tangential tensors on the hypersurface 3 . We note 'x's are coordinates on region M, while 'y's are coordinates associated with the hypersurface and we may also keep in mind that, where h ij is the inverse of the induced metric h ij on the hypersurface, see for instance [55]. The change of direction of the normal n as one moves on the hypersurface corresponds to the bending of the hypersurface Σ t which is described by the extrinsic curvature. The extrinsic curvature of spatial slices where time is constant is given by: where D i = e µ i ∇ µ is the intrinsic covariant derivative associated with the induced metric defined on the hypersurface, and e µ i is the appropriate basis vector which is used to transform bulk indices to boundary ones.
Armed with this information, one can write down the Gauss, Codazzi and Ricci equations, see [51]: where in the left hand side of Eq. (4.8) we have the bulk Riemann tensor, but where all indices are now spatial rather than both spatial and temporal, and R ijkl is the Riemann tensor constructed purely out of h ij , i.e. the metric associated with the hypersurface; and £ β is the Lie derivative with respect to shift 4 .

Coframe Slicing
A key feature of the 3+1 decomposition is the free choice of lapse function and shift vector which define the choice of foliation at the end. In this paper we stick to the coframe slicing.
The main advantage for this choice of slicing is the fact that the line element and therefore components of the infinite derivative function in our gravitational action will be simplified greatly. In addition, [56] has shown that such a slicing has a more transparent form of the canonical action principle and Hamiltonian dynamics for gravity. This also leads to a well-posed initial-condition for the evolution of the gravitational constraints in a vacuum by satisfying the Bianchi identities. In order to map the ADM line element into the coframe slicing, we use the convention of [56]. We define where x i and i = 1, 2, 3 is the spatial and t is the time coordinates. 5 The metric in the coframe takes the following form where upon substituting Eq. (4.11) into Eq. (4.12) we recover the original ADM metric given by Eq. (4.1). In this convention, if we take g as the full spacetime metric, we have the following simplifications: The convective derivatives ∂ α with respect to θ α are (4.14) 4 We have £ β Kij ≡ β k D k Kij + K ik Dj β k + K jk Diβ k . 5 We note that in Eq. (4.11), the "i" for θ i is just a superscript not a spatial index.
For time-dependent space tensors T , we can define the following derivative: where £ β is the Lie derivative with respect to the shift vector β i . This is because the off-diagonal components of the coframe metric are zero, i.e., g 0i = g 0i = 0.
We shall see later on how this time slicing helps us to simplify the calculations when the gravitational action contains infinite derivatives.

Extrinsic Curvature
A change in the choice of time slicing results in a change of the evolution of the system. The choice of foliation also has a direct impact on the form of the extrinsic curvature. In this section we wish to give the form of extrinsic curvature K ij in the coframe slicing. This is due to the fact that the definition of the extrinsic curvature is an initial parameter that describes the evolution of the system, therefore is it logical for us to derive the extrinsic curvature in the coframe slicing as we use it throughout the paper. We use [56] to find the general definition for K ij in the coframe metric. In the coframe, where Γ is the ordinary Christoffel symbol and "∧" denotes the exterior or wedge product of vectors θ. By finding the coefficients Cs and subsequently calculating the connection coefficients γ α βγ , one can extract the extrinsic curvature K ij in the coframe setup. We note that the expression for dθ α is the Maurer-Cartan structure equation [60]. It is derived from the canonical 1-form θ on a Lie group G which is the left-invariant g-valued 1-form uniquely determined by θ(ξ) = ξ for all ξ ∈ g.
We can use differential forms (See Appendix A) to calculate the Cs, the coefficients of dθ where now we can write, where k = 1, 2, 3. Now when we insert the Cs from Appendix A, and From the definition of dθ α in Eq. (4.17) and using the antisymmetric properties of the ∧ product, we get Using these properties, we find that C m 0i = ∂β m ∂x i , C m ij = 0 and C 0 ij = 0. Using Eq. (4.16), we obtain that Since from Eq. (4.7) the expression for the extrinsic curvature in coframe slicing is given by: where ∂ 0 is the time derivative and β l is the "shift" in the coframe metric Eq. (4.12).

Riemann Tensor in the Coframe
The fact that we move from the ADM metric into the coframe slicing has the following implication on the form of the components of the Riemann tensor. Essentially, since in the coframe slicing in Eq. (4.12) we have g 0i = g 0i = 0, therefore we also have, from Eq. (4.2), n i = n i = 0 (n 0 and n 0 stay the same as in Eq. (4.2)). Hence the non-vanishing components of the Riemann tensor in the coframe, namely Gauss, Codazzi and Ricci tensor, become: with∂ 0 defined in Eq. (4.15) and D j = e µ j ∇ µ . It can be seen that the Ricci equation, given in Eq. (4.8) is simplified in above due to the definition of Eq. (4.15). We note that Eq. (4.26) is in the coframe slicing, while Eqs. (4.8-4.10) are in the ADM frame only.

D'Alembertian Operator in Coframe
Since we shall be dealing with a higher-derivative theory of gravity, it is therefore helpful to first obtain an expression for the operator in this subsection. To do so, we start off by writing the definition of a single box operator in the coframe, [62], where we note that the Greek indices run from 1 to 4 and the Latin indices run from 1 to 3 (ε = −1 for a spacelike hypersurface). We call the spatial box operator hyp = h ij D i D j , which stands for "hypersurface" as the spatial coordinates are defined on the hypersurface meaning hyp is the projection of the covariant d'Alembertian operator down to the hypersurface, i.e. only the tangential components of the covariant d'Alembertian operator are encapsulated by hyp . Also note that in the coframe slicing g ij = h ij . Generalising this result to the nth power, for our purpose, we get where the f in s are the coefficients of the series.

Generalised Boundary Term
In this section, first we are going to briefly summarise the method of [51] for finding the boundary term. It has been shown that, given a general gravitational action one can introduce two auxiliary fields ̺ µνρσ and ϕ µνρσ , which are independent of each other and of the metric g µν , while they have all the symmetry properties of the Riemann tensor R µνρσ . We can then write down the following equivalent action: The reason we introduce these auxiliary fields is that the second derivatives of the metric appear only linearly in Eq. (5.2). Note that in Eq. (5.2), the terms involving the second derivatives of the metric are not multiplied by terms of the same type, i.e. involving the second derivative of the metric, so when we integrate by parts once, we are left just with the first derivatives of the metric; we cannot eliminate the first derivatives of the metric as well -since in this paper we are keeping the boundary terms. Note that the first derivatives of the metric are actually contained in these boundary terms if we integrate by parts twice, see our toy model scalar field theory example in Eqs. (2.2,2.3) 6 . Therefore, terms which are linear in the metric can be eliminated if we integrate by parts; moreover, the use of the auxiliary fields can prove useful in a future Hamiltonian analysis of the action. From [51], we then decompose the above expression as are equivalent to the components of the Gauss, Codazzi and Ricci equations given in Eq. (4.8), also, where φ ijkl , φ ijk and Ψ ij are spatial tensors evaluated on the hypersurface. The equations of motion for the auxiliary fields ϕ µνρσ and ̺ µνρσ are, respectively given by [51], where R µνρσ is the four-dimensional Riemann tensor. One can start from the action given by Eq. (5.2), insert the equation of motion for ϕ µνρσ and recover the action given by Eq. (5.1). It has been shown by [51] that one can find the total derivative term of the auxiliary action as where K = h ij K ij , with K ij given by Eq. (4.25), and Ψ = h ij Ψ ij , where Ψ ij is given in Eq. (5.10), are spatial tensors evaluated on the hypersurface Σ t and L is the Lagrangian density. In Eq. (5.7), the second term is the total derivative. It has been shown that one may add the following action to the above action to eliminate the total derivative appropriately. Indeed Ψ can be seen as a modification to the GHY term, which depends on the form of the Lagrangian density [51].
where n µ is the normal vector to the hypersurface and the infinitesimal vector field is normal to the boundary ∂M and is proportional to the volume element of ∂M; in above ε µαβγ = √ −g[µ α β γ] is the Levi-Civita tensor and y are coordinates intrinsic to the boundary 7 , and we used Eq. (4.5). Moreover in Eq. (5.8), we have: where f indicates the terms in the Lagrangian density and is built up of tensors ̺ µνρσ , ̺ µν and ̺ as in Eq. (5.2); G is the universal gravitational constant and Ω ij is given in Eq. (5.4). Indeed, the above constraint is extracted from the equation of motion for Ω ij in the Hamiltonian regime [51]. In the next section we are going to use the same approach to find the boundary terms for the most general, covariant quadratic order action of gravity. 7 We shall also mention that Eq.(5.8) is derived from Eq.(5.7) by performing Stokes theorem, that is

Boundary Terms for Finite Derivative Theory of Gravity
In this section we are going to use the 3+1 decomposition and calculate the boundary term of the EH term R, and as prescribed in previous section, as a warm-up exercise. We then move on to our generalised action given in Eq. (3.2). To decompose any given term, we shall write them in terms of their auxiliary field, therefore we have R = ̺, R µν ≡ ̺ µν , and R µνρσ ≡ ̺ µνρσ , where the auxiliary fields ̺, ̺ µν and ̺ µνρσ have all the symmetry properties of the Riemann tensor. We shall also note that the decomposition of the operator in 3+1 formalism in the coframe setup is given by Eq. (4.27).

R
For the Einstein-Hilbert term R, in terms of the auxiliary field ̺ we find in Appendix B.1 where Ω = h ij Ω ij and we used h ij h kl ρ ijkl = ρ, and h ij ρ iνjσ n ν n σ = h ij Ω ij and ̺ ≡ R in the EH action and the right hand side of Eq. (6.1) is the 3 + 1 decomposed form of the Lagrangian and hence ρ and Ω are spatial. We may note that the last term of the expansion on the second line of Eq. (6.1) vanishes due to the symmetry properties of the Riemann tensor. Using Eq. (5.10), and calculating the functional derivative, we find This verifies the result found in [51], and it is clear that upon substituting this result into Eq. (5.8), we recover the well known boundary for the EH action, as K = h ij K ij and Ψ · K ≡ Ψ ij K ij where K ij is given by Eq. (4.25). Hence, where dΣ µ is the normal to the boundary ∂M and is proportional to the volume element of ∂M while n µ is the normal vector to the hypersurface.

R µνρσ R µνρσ
Next, we start off by writing R µνρσ R µνρσ as its auxiliary equivalent ̺ µνρσ ̺ µνρσ to obtain where ̺ µνρσ = δ α µ δ β ν δ γ ρ δ λ σ ̺ αβγλ (where δ α µ is the Kronecker delta). This allowed us to use the completeness relation as given in Eq. (4.4). In Eq. (6.4), we used the antisymmetry properties of the Riemann tensor to eliminate irrelevant terms in the expansion. From Eq. (6.4), we have three types of terms: hhhh, hhhnn, hhnnnn.
The aim is to contract the tensors appearing in Eq. (6.4) and extract those terms which are Ω ij dependent. This is because we only need Ω ij dependent terms to obtain Ψ ij as in Eq. (5.10) and then the boundary as prescribed in Eq. (5.8).
A closer look at the expansion given in Eq. (6.4) leads us to know which term would admit Ω ij type terms. Essentially, as defined in Eq. (5.4), Ω ij = n µ n ν ̺ iµjν , therefore by having two auxiliary field tensors as ̺ αβγλ and ̺ µνρσ in Eq. (6.4) (with symmetries of the Riemann tensor) we may construct Ω ij dependent terms. Henceforth, we can see that in this case the Ω ij dependence comes from the hhnnnn term.
To see this explicitly, note that in order to perform the appropriate contractions in presence of the d'Alembertian operator, we first need to complete the contractions on the left hand side of the operator. We then need to commute the rest of the tensors by using the Leibniz rule to the right hand side of the components of the operator, i.e. the∂ 0 's and the hyp , and only then do we obtain the Ω ij type terms.
We first note that the terms that do not produce Ω ij dependence are not involved in the boundary calculation, however they might form ρ ijkl , ρ ijk , or their contractions. These terms are equivalent to the Gauss and Codazzi equations as shown in Eq. (5.4), and we will address their formation in Appendix B.2. In addition, as we shall see, by performing the Leibniz rule one produces some associated terms, the X ij 's, which appear for example in Eq. (6.5). Again we will keep them only if they are Ω ij dependent, if not we will drop them.
hhnnnn terms: To this end we shall compute the hhnnnn terms, hence we commute the h's and n's onto the right hand side of the in the hhnnnn term of Eq. (6.4): where Ω ij ≡ h ik e k κ h jm e m λ n γ n δ ̺ γκδλ = h ik h jm n γ n δ ̺ γkδm ; we note that X ij 1 only appears because of the presence of the operator.
The term Ω rs X rs 1 will yield X ij 1 when functionally differentiated with respect to Ω ij as in Eq. (5.10). Also note X ij 1 does not have any Ω ij dependence. Similarly for the other X terms which appear later in the paper. We shall note that when we take = 1 in Eq. (6.4), we obtain, where we just contract the indices and we do not need to use the Leibniz rule as we can commute any of the tensors, therefore we do not produce any X ij terms at all 8 . Finally, one can decompose Eq. (6.4) as where "· · · " are terms such as ρ ijkl ρ ijkl , ρ ijk ρ ijk and terms that are not Ω ij dependent and are the results of performing the Leibniz rule (see Appendix B.2). When we take M 2 → ∞, i.e., when we set → 0 (recall that has an associated mass scale /M 2 ), which is also equivalent to considering α → 0 in Eq. (3.2), we recover the EH result.
When → 1, we recover the result for R µνρσ R µνρσ found in [51]. At both limits, → 0 and → 1, the X ij 1 term is not present. To find the boundary term, we use Eq. (5.10) and then Eq. (5.8). We are going to use the Euler-Lagrange equation and drop the total derivatives as a result. We have, Hence the boundary term for R µνρσ R µνρσ is, where K ij is given by Eq. (4.25).

R µν R µν
We start by first performing the 3+1 decomposition of R µν R µν in its auxiliary form ̺ µν ̺ µν , where we have used appropriate contractions to write the Ricci tensor in terms of the Riemann tensor. As before, we then used the completeness relation Eq. (4.4) and used the antisymmetric properties of the Riemann tensor to drop the vanishing terms. We are now set to calculate each term, which we do in more detail in Appendix B.3. Again our aim is to find the Ω ij dependent terms, by looking at the expansion given in Eq. (6.11) and the distribution of the indices, the reader can see that the terms which are Ω ij dependent are those terms which have at least two nns contracted with one of the ̺s such that we form n µ n ν ̺ iµjν .
• hhhnn terms: We start with the hhhnn terms in Eq. (6.11). We calculate the first of these in terms of Ω ik and ρ ik also by moving the 'h's and 'n's onto the right hand side of the , where the contraction is h ij e κ j h kl e λ l h mn e γ m e δ n ̺ γκδλ = h ij h kl ρ jl = ρ ik , and • hhhnn trems: The next hhhnn term in Eq. (6.11) is h mn e λ n n γ n δ ]D a ̺ γκδλ = ρ km Ω km + · · · , (6.14) where we used h kl e κ l h mn e λ n n γ n δ ̺ γκδλ = h kl h mn n γ n δ ̺ γlδn = Ω km and we note that "· · · " are extra terms which do not depend on Ω km .
• hhnnnn terms: The the next term in Eq. (6.11) is of the form hhnnnn: j e λ l n γ n δ ] − D a [e κ j e λ l n γ n δ ]D a ̺ γκδλ = Ω jl Ω jl + Ω jl X 2(b)jl , (6.15) where e κ j e λ l n γ n δ ̺ γκδλ = n γ n δ ̺ γjδl = Ω jl , and • hhnnnn terms: Finally, the last hhnnnn terms in Eq. (6.11) is where we used n κ n λ h mn e γ m e δ n ̺ γκδλ = h mn Ω mn = Ω, and Summarising this result, we can write Eq. (6.11), as where "· · · " are the contractions of ρ ijkl and ρ ijk (see Appendix B.3) and the terms that are the results of performing Leibniz rule, which have no Ω ij dependence. When → 1, we recover the result for R µν R µν found in [51]. At both limits, → 0 and → 1, the X 2 terms are not present. Obtaining the boundary term requires us to extract Ψ ij as it is given in Eq. (5.10). Hence the boundary for R µν R µν is given by, where K ≡ h ij K ij and K ij is given by Eq. (4.25).

R R
We do not need to commute any h's, or n's across the here, we can simply apply Eq. (6.1) to ̺ ̺, the auxiliary equivalent of the R R term: whereupon extracting Ψ ij using Eq. (5.10), and using Eq. (5.8) as in the previous cases, we obtain the boundary term for R R to be where K ≡ h ij K ij and K ij is given by Eq. (4.25). Again when → 1, we recover the result for R 2 found in [51].

Full result
Summarising the results of Eq. (6.10), Eq. (6.20) and Eq. (6.22), altogether we have This result matches with the EH action [51], when we take the limit → 0; that is, we are left with the same expression for boundary as in Eq. (6.3): since the X-type terms are not present when → 0. When → 1, we recover the result for R + α(R 2 + R µν R µν + R µνρσ R µνρσ ) found in [51]; that is, we are left with again the X-type terms are not present when → 1. We should note that the X 1 and X 2 terms are the results of having the covariant d'Alembertian operator so, in the absence of the d'Alembertian operator, one does not produce them at all and hence the result found in [51] is guaranteed.
We may now turn our attention to the R 2 R, R µν 2 R µν and R µνρσ 2 R µνρσ . Here the methodology will remain the same. One first decomposes each term into its 3+1 equivalent. Then one extracts Ψ ij using Eq. (5.10), and then the boundary terms can be obtained using Eq. (5.8). In this case we will have two operators, namely This means that upon expanding to 3 + 1, one performs the Leibniz rule twice and hence obtains eight total derivatives that do not produce any Ω ij s or its contractions that are relevant to the boundary calculations and hence must be dropped.

Generalisation to Infinite Derivative Theory of Gravity
We may now turn our attention to the infinite derivative terms; namely, RF 1 ( )R, R µν F 2 ( )R µν and R µνρσ F 3 ( )R µνρσ . For such cases, we can write down the following relation (see Appendix B.4): where X and Y are tensorial structures such as ̺ µνρσ , ̺ µν , ̺ and their contractions, while D denotes any operators. These operators do not have to be differential operators and indeed this result can be generalised to cover the case where there are different types of operator and a similar (albeit more complicated) structure is recovered. From (6.27), one produces 2n total derivatives, analogous to the scalar toy model case, see Eqs. (2.2,2.3). We can then write the 3 + 1 decompositions for each curvature by generalising Eq. (6.10), Eq. (6.20) and Eq. (6.22) and writing R µνρσ F 3 ( )R µνρσ , R µν F 2 ( )R µν and RF 1 ( )R in terms of their auxiliary equivalents ̺ µνρσ F 3 ( )̺ µνρσ , ̺ µν F 2 ( )̺ µν and ̺F 1 ( )̺. Then where we have dropped the irrelevant terms, as we did before, while X ij 1 and X ij 2 are the analogues of Eqs. (6.6), (6.13), (6.16), (6.18). We now need to use the generalised form of the Euler-Lagrange equations to obtain the Ψ ij in each case: where we have imposed that δΩ ij = 0 on the boundary ∂M. Hence, by using Ω = h ij Ω ij and ρ = h ij ρ ij , we find in Appendix C that: and so using Eq. (5.10), the Ψ ij s are: where we have used Eq. (6.1) in the last line. Finally, we can use Eq. (5.8) and write the boundary terms corresponding to our infinite-derivative action as, where Ω ij = n γ n δ ̺ γiδj , Ω = h ij Ω ij , ρ ij = h km ρ ijkm , ρ = h ij ρ ij , K = h ij K ij and K ij is the extrinsic curvature given by Eq. (4.25). We note that when we decompose the , after we perform the Leibniz rule enough times, we can reconstruct the in its original form, i.e. it is not affected by the use of the coframe. In this way, we can always reconstruct F i ( ). However, the form of the X-type terms will depend on the decomposition and therefore the use of the coframe. In this regard, the X-type terms depend on the coframe but the F i ( ) terms do not.

Conclusion
This paper generalises earlier contributions for finding the boundary term for a higher derivative theory of gravity. Our work has focused on seeking the boundary term or GHY contribution for a covariant infinite derivative theory of gravity, which is quadratic in curvature. Such modifications are inevitable for any covariant construction of gravitational theory. In stringy parlance, such infinite derivative quadratic order curvature correction entails to one-loop in string coupling g s , and all order α ′ corrections. Indeed, here one can possibly make the connection even stronger, if one could find an exact form of the form factors, i.e. F( )'s from closed-string field theory. However, note that our considerations are much more general than any such possible stringy correction, as the analytical construction does not rely either on supersymmetry or in any specific spacetime dimensions. Indeed, we only concentrated on pure massless gravitational degrees of freedom, which is also in stark contrast with the string theory set-up.
Indeed, in this case some novel features distinctively filter through our analysis. Since the bulk action contains nonlocal form factors, F i ( ), the boundary action also contains the nonlocality, as can be seen from our final expression Eq. (6.34). The above expression also has a smooth limit when M → ∞, or → 0, which is the local limit of Eq. (3.5), and our results then reproduce the GHY term, and when F i ( ) → 1, our results coincide with that of [51].
Our study will have implications for finding the Hamiltonian for an infinite derivative theory of gravity, as well as seeking the entanglement entropy, and to understand gravity/field theory correspondence in Anti (de)Sitter spacetimes. Some of these issues will be studied in future publication.

A K ij in the Coframe Metric
In this section we wish to use the approach of [56] and find the general definition for K ij in the coframe metric. Given, where Γ is the ordinary Christoffel symbol, ∧ is the ordinary wedge product and the Cs are coefficients to be found. By comparing the values given in [56] with the ordinary Christoffel symbols, we can see that Now in the coframe metric in Eq. (4.12), .4) and C j0i = g αj C α 0i = g jk C k 0i . In general, for a p-form α and a q-form β, Hence, if p is odd, From Eq. (A.2) and Eq. (A.3) we can see that We get a similar result for dθ 2 and dθ 3 , so we can say that where k = 1, 2, 3. Now from the definition of θ in Eq. (4.17), and Let us point out that Now using Eq. (A.9) and Eq. (A.11), where k = 1, 2, 3. From the definition of dθ α in Eq. (A.2) and using the antisymmetric properties of the ∧ product from Eq. (A.7), and therefore we can then write In order for this to be satisfied, each term must vanish separately as the dx α ∧ dx j are linearly independent and so the coefficient of each must be zero and thus Again, in order for this relation to be satisfied, C 1 12 = C 1 13 = C 1 23 = 0 and C 1 01 = ∂β 1 ∂x 1 , We deduce that C m 0i = ∂β m ∂x i , C m ij = 0 and C 0 ij = 0. Using Eq. (A.2) and that in the coframe Γ 0 Since from Eq. (4.7) Eq. (A.1) becomes noting that the term with four n α s vanishes due to the antisymmetry of the Riemann tensor in the first and last pair of indices (recall that ̺ µνρσ has the same symmetry properties as the Riemann tensor)

B.2 Riemann Tensor
In this section we wish to show the contraction of the rest of the terms in Eq. (6.4) for the sake of completeness. We have, from hhhh, which produced ρ ijkl ρ ijkl and the terms which are the results of Leibniz rule. Next in Eq. (6.4) is, with ρ ijk ≡ n µ ρ ijkµ . Here we produced ρ ijk ρ ijk and the extra terms which are the results of the Leibniz rule. Similarly we can find the contractions for different terms in Eq. (6.4).

B.3 Ricci Tensor
In similar way as we did in the Riemann case we can find all the other contractions in the expansion of Eq. (6.11) which we omitted. They are: with (h jn e κ n )(h kx e λ x )(h ly e γ l e δ y )̺ γκδλ = ρ jk . Above we produced ρ jk ρ jk plus other terms that are results of the Leibniz rule. And, where we used n µ ̺ µk = ρ k and n κ h kl e λ l h mn e γ m e δ n ̺ γκδλ = n κ h kl ̺ κl = ρ k . We produced ρ k ρ k plus other terms that are results of the Leibniz rule. We may also note that one can write, ρ ij ≡ h kl ρ ikjl , ρ ≡ h ik h il ρ ijkl and ρ i ≡ h jk ρ jik .

B.4 Generalisation from to F( )
In Eq. (6.5) for 2 , we have, As a general rule we can write, where X and Y are some tensors and D is some operator. Applying this we can write, where we dropped the irrelevant terms. We moreover can generalise the result of (B.7) and write,

C Functional Differentiation
Given the constraint equation suppose that f = ΩF( )Ω and F( ) = ∞ n=0 f n n , where the coefficients f n are massless 9 . Then, using the generalised Euler-Lagrange equations, we have in the coframe (and 9 Recall that the term comes with an associated scale /M 2 . imposing the condition that δΩ ij = 0 on the boundary ∂M) where we have used that g ij = h ij = 0. Note also that: So we can summarise the results and write, = Ω ij + Ω ij = 2 Ω ij . (C.5) and generalise this to: δ ΩF( )Ω δΩ ij = 2h ij F( )Ω, δ Ω ij F( )Ω ij δΩ ij = 2F( )Ω ij , (C.10)

D Riemann tensor components in ADM gravity
Using the method of [61], we can find the Riemann tensor components. The Christoffel symbols for the ADM metric in Eq. (4.1) are where K ij is the extrinsic curvature given by (4.7) and in the ADM metric, N is the lapse, β i is the shift and h ij is the induced metric on the hypersurface. Now we can find the Riemann tensor components where R ijkl is the Riemann tensor of the induced metric on the hypersurface. Then Relabelling the indices, we obtain that Finally, we have that Hence n µ n ν R µiνj = n 0 n µ R µi0j + n k n µ R µikj where £ β K ij ≡ β k D k K ij + K ik D j β k + K jk D i β k . Therefore overall, we have
Then using the same method as in Eq. (D.2), Finally we have that in the coframe, Hence the non-vanishing components of the Riemann tensor in the coframe, namely the Gauss, Codazzi and Ricci tensor, become: where K ij is the extrinsic curvature of the hypersurface, given in the coframe by Eq. (4.25) and R ijkl is the Riemann tensor of the induced metric on the hypersurface.