Non-local Quantum Geometrodynamics: Operational Framework of an Efimov-based Gravity Model at the Planck Scale
Abstract
This paper presents a comprehensive theoretical proposal for a Quantum Gravity (QG) model, utilizing the mathematical apparatus of Non-local Quantum Field Theory (NLQFT) within the G.V. Efimov framework.
Full text
Non-local Quantum Geometrodynamics: Operational Framework of an Efimov-based Gravity Model at the Planck Scale Damian Pikor and Paweł Kurzawski (Dated: December 4, 2025) This paper presents a comprehensive theoretical proposal for a Quantum Gravity (QG) model, utilizing the mathematical apparatus of Non-local Quantum Field Theory (NLQFT) within the G.V. Efimov framework. We effect a fundamental paradigm shift by substituting fermionic fields with the gravitational (metric) field and setting the non-locality scale at the Planck Mass ( MP l ≈ 1 . 22 × 10 19 GeV). We demonstrate that applying a Gaussian regulator in the form of an entire function to the graviton kinetic operator yields a ghost-free, perturbatively unitary, and UV-finite theory. We provide a perturbative proof of unitarity within a class of entire form factors via Dressed Cutkosky Rules and derive the generalized Slater-Taylor identities using the BRST formalism. Addressing the Initial Value Problem, we show that the theory possesses a finite number of degrees of freedom, avoiding Ostrogradsky instabilities. Phenomenologically, the model predicts a characteristic suppression in trans-Planckian scattering, is consistent with GW170817 speed limits, and resolves singularities via the effective Hayward metric. I. INTRODUCTION: THE CRISIS OF LOCALITY A. Anatomy of Divergence Contemporary theoretical physics faces a fundamental rift between its two pillars: General Relativity (GR) and Quantum Field Theory (QFT). The standard quantization of the gravitational field based on the paradigm of locality leads to a non-renormalizable theory. The Newtonian coupling constant GN , having a dimension of [ mass ] −2 , causes higher-order Feynman diagrams to generate divergences that scale as powers of the loop order. To remove them would require infinite counterterms, stripping the theory of predictive power. While String Theory or Loop Quantum Gravity postulate fundamental discreteness or extended objects, this work follows a different path inspired by phenomenological successes in the electroweak sector. We propose maintaining spacetime continuity while introducing a fundamental non-locality of interactions, controlled by the Planck scale. B. The Efimov Paradigm: From Leptons to Gravitons Previous work demonstrated that introducing a nonlocal form factor V ( D2 )allows for a consistent regularization of gauge theories. We transpose these concepts to gravity via the fundamental substitution: ψ(x)→gµν(x),Λ→MP l ≈1019 GeV.(1) The U (1) gauge symmetry is replaced by Diffeomorphism Invariance. This is non-trivial, as gravity is inherently non-linear; gravitons act as sources for themselves. The non-local "smearing" of the gravitational interaction must therefore be defined self-consistently. II. THEORETICAL FRAMEWORK A. Analytic Structure and Dressed Cutkosky Rules The consistency of any higher-derivative field theory is inextricably linked to the problem of unitarity. Local higher-order theories, such as fourth-order gravity ( R + R2 ), typically exhibit additional poles with negative residues (ghosts) that violate the S†S = 1 condition. Our model avoids this pathology by modifying the analytic structure in the complex momentum plane using entire functions. 1. Structure of the Dressed Propagator We introduce the non-local graviton propagator ˜ D ( p2 )by modifying the standard Einstein-Hilbert propagator with a form factor V(z) = eH(z), where z=−p2/M2 P l: ˜ Dµνρσ(p2) = V−2(p2) p2Pµνρσ =e−2H(−p2/M2 P l) p2Pµνρσ. (2) The exponential function e−H(z) is an entire function of order ρ∈ [1 / 2 , 1). According to the Hadamard factorization theorem, ez has no zeros in the finite complex plane C . Thus, the inverse propagator equation p2eH = 0 yields only the kinematic pole p2 = 0. No massive ghosts or tachyons appear in the spectrum. 2. Dressed Cutkosky Rules and the Optical Theorem Rigorous proof of unitarity requires demonstrating that the imaginary part of the forward scattering amplitude arises solely from physical on-shell processes. We utilize the Dressed Cutkosky Rules. The optical theorem relates the imaginary part of the amplitude M to the crosssection: 2ImM(A→A) = X XZdΠX|M(A→X)|2.(3) In our formalism, the "cut" of the dressed propagator is
2 defined as: Disc "e−H(p2) p2#= (−2πi)δ(p2)e−H(0).(4) By normalizing H (0) = 0, we have e−H(0) = 1. Consequently, the residue at the pole is +1, identical to GR. The non-local factor suppresses off-shell loop integrals (UV finiteness) but reduces to unity for physical states. Thus, the sum over intermediate states X includes only physical gravitons, preserving unitarity. B. Global Dressing Ansatz in GR We consider a perturbative approach around a flat background ηµν . We construct the non-local action by introducing a "dressed metric" gµν , related to the local field via the non-local convolution operator: gµν(x) = V□g M2 P l hµν.(5) The Efimov Lagrangian for gravity takes the form: LQG =2 κ2√−gR[g].(6) C. Local limit and recovery of GR It is useful to make explicit how the non-local theory reduces to General Relativity in the appropriate limit. Since the form factor is an entire function of z = −p2/M2 P l , we assume that it admits a regular Taylor expansion around z= 0, V(z) = 1 + O(z),(7) so that for fixed physical momenta and MP l → ∞ one has lim MP l→∞ V−p2 M2 P l = 1.(8) Consequently, the dressed metric tends to the local one, gµν(x)=V□g M2 P l hµν(x)−→ hµν(x),(9) and the Lagrangian LQG continuously reduces to the Einstein–Hilbert Lagrangian. In this sense the Efimovinspired construction should be viewed as a UV completion of GR rather than a separate low-energy theory. D. Geometric Interpretation: Effective Weyl Geometry The non-local dressing induces an effective geometry characterized by scale-dependent non-metricity. We identify the non-metricity tensor of the effective geometry as Qλµν = ∇λgµν . Substituting the dressed metric ansatz g=Vg, we find: Qλµν ≈(∂λln V)gµν.(10) This identifies the gradient of the logarithmic form factor as the Weyl gauge vector Aλ = ∂λln V . Consequently, spacetime probed at trans-Planckian energies is effectively a Weyl manifold, where parallel transport modifies vector norms dependent on the energy scale. E. Conditions on the entire form factor For phenomenological viability and perturbative consistency, the non-local form factor V ( z ) = eH(z) must satisfy a small set of assumptions that guarantee ghost freedom and good UV behaviour. In particular, throughout this work we assume: 1. H (0) = 0, which normalizes the massless graviton pole and ensures that the on-shell residue coincides with the Einstein–Hilbert theory; 2. Re H ( z ) ≥ 0for Re z≥ 0, so that the form factor provides exponential suppression of Euclidean loop integrals without generating additional poles; 3. polynomial (or slower) growth of H ( z )in any finite angular sector of the complex plane, which guarantees absolute convergence of the Euclidean amplitudes and allows for analytic continuation to Minkowski space. Under these conditions the modified propagator has a single physical pole at p2 = 0 with positive residue, while UV divergences are tamed by the exponential damping. These properties are used explicitly in the one-loop unitarity proof summarized in Appendix A. III. GAUGE INVARIANCE AND CONSISTENCY A. BRST Symmetry and Slater-Taylor Identities Naive insertion of non-local regulators breaks diffeomorphism invariance since ∇µ does not commute with V ( □ ). To resolve this, we employ the BRST formalism. We introduce ghost fields cµ , anti-ghosts ¯cµ , and NakanishiLautrup auxiliary fields bµ . The BRST transformation s is nilpotent (s2= 0). For consistency, the generating functional Z [ J ]must satisfy the generalized Slater-Taylor Identities (STI). In the presence of non-local regulators, the standard identity generalizes via the Master Equation. To preserve the transversality condition pµ Γ µν... = 0, the interaction vertex must be "dressed" with a compensating structure. We derive that the full vertex must satisfy the Inverse Propagator Identity: pµΓµ...(p,...) = D−1(p+k)−D−1(k)+ghost terms.(11)
3 This requires the introduction of a Gauge Link operator L , defined integrally to effectively sum the noncommutativity of covariant derivatives. This structure ensures that the effective action remains invariant under the full diffeomorphism group, resolving the issue of current non-conservation. B. Hamiltonian Consistency and the Initial Value Problem A common critique of non-local theories involves the Ostrogradsky instability (unbounded Hamiltonian due to N > 2time derivatives). However, our theory operates in the limit N→ ∞ . Following the formalism of Barnaby and Kamran: 1. Finite Degrees of Freedom: Since the regulator e−p2/M2 P l introduces no new poles, the number of initial conditions required to specify the system does not diverge. The phase space does not expand infinitely; the number of physical degrees of freedom remains identical to GR (2 polarizations). 2. Diffusion Equation Analogy: The non-local operator e∇2 behaves like a diffusion operator. By redefining field variables ˜ ϕ = e−□/2M2 P l ϕ , the action can be mapped to a local form S∼˜ ϕ□˜ ϕ with smeared sources. This reveals that the apparent instability is an artifact of coordinates. In physical variables, the energy is bounded from below, and the Initial Value Problem is well-posed. Dirac analysis and number of degrees of freedom. The above arguments can be made fully precise within the Dirac theory of constrained Hamiltonian systems. At the linearized level, the non-local kernel is traded for an auxiliary diffusion coordinate and a set of first-class constraints that extend the familiar diffeomorphism constraints of GR. A detailed counting shows that the physical phase space still carries Ndof = 2 (12) propagating degrees of freedom, corresponding to the two graviton polarizations, with no additional Ostrogradsky modes. The would-be instability is removed because no extra independent canonical pairs associated with higher time derivatives appear; instead, non-locality is encoded in constrained auxiliary fields. The explicit construction is presented in Appendix B. IV. PHENOMENOLOGY AND OBSERVATIONAL CONSTRAINTS A. The Suppression Signature (Trans-Planckian Scattering) Unlike models predicting contact interactions (which yield an excess of events ∼ + s/ Λ 2 ), the Efimov model predicts a suppression. The propagator acts as: D(s)∼1 se−s/M2 P l ≈1 s1−s M2 P l .(13) In trans-Planckian collisions ( E≫MP l ), instead of "asymptotic darkness" (black hole formation), we predict "asymptotic transparency." 0.5 1 1.5 2 0 0.5 1 1.5 2 GR Baseline Center-of-Mass Energy √s/MP l Cross Section Ratio σ/σSM Trans-Planckian Scattering Signatures Efimov (Suppression) BH / Contact Excess Figure 1. Comparison of cross-section predictions. Traditional quantum gravity or compositeness models often predict an excess (red dashed). The Efimov model predicts a distinct suppression (blue solid). B. GW170817 and the Speed of Gravity The detection of the neutron star merger GW170817 constrained the deviation of the speed of gravitational waves cg from the speed of light c to |cg/c − 1 |< 10 −15 . This rules out many modified gravity theories predicting massive gravitons. In our model, the dispersion relation is determined by the pole of the propagator: p2e−p2/M2 P l = 0 =⇒p2= 0.(14) Since the exponential factor is non-vanishing everywhere, the only solution is the standard massless dispersion relation E2 = | k|2 . Consequently, the group velocity is identically cg = 1. The non-local dressing affects interaction vertices (off-shell) but does not shift the kinematic pole. Thus, the model is compatible with GW170817 constraints at leading order (see Appendix Cfor details). Microcausality and observable effects. Non-local entire form factors generically induce tiny violations of microcausality at distances of order the non-locality scale, where the smearing of interaction vertices becomes relevant. In the present construction this scale is set by M−1 P l , so any acausal effects are exponentially suppressed at macroscopic distances and completely negligible for the astrophysical and cosmological observables considered here.
4 C. The Hayward Metric and Regular Black Holes The theory predicts the resolution of gravitational singularities due to the "smearing" of interaction vertices at the Planck scale. Instead of a Schwarzschild singularity, we identify the emergence of a static, spherically symmetric solution described by the Hayward Metric. The effective energy-momentum tensor induced by the non-local smearing of a point mass M behaves as a Gaussian distribution ρ ( r ) ∼Me−r2/l2 , where l∼ 1 /MP l . The resulting line element is: ds2=−f(r)dt2+1 f(r)dr2+r2dΩ2,(15) with the metric function: f(r)=1−2Mr2 r3+ 2Ml2.(16) • Asymptotic ( r→ ∞ ): f ( r ) ≈ 1 − 2 M/r , recovering Schwarzschild GR. • Core ( r→ 0): f ( r ) ≈ 1 − ( r/l ) 2 , representing a regular de Sitter core. This metric is free of curvature singularities (finite Kretschmann scalar) and represents a regular black hole, a concrete prediction of the model (further discussed in Appendix D). 0 0.5 1 1.5 2 2.5 3 0 0.5 1 Radius r/M Metric Potential f(r) Regular Black Hole Potential (Hayward) Schwarzschild Hayward (Regular) Figure 2. Comparison of the singular Schwarzschild potential (dashed) and the regular Hayward potential (solid) derived from the non-local model. V. ONE-LOOP SELF-ENERGY AND SUPER-RENORMALIZABILITY The defining feature of the Efimov approach is the suppression of UV divergences via the exponential vertex regulator. We demonstrate this via the one-loop graviton self-energy Σµναβ(p)(vacuum polarization). 1. Finiteness of Loop Diagrams The superficial degree of divergence D for a diagram with Lloops behaves as: D∼4L−2I+VGR −VNL ×(regulator decay).(17) Since the regulator e−k2/Λ2 decays exponentially for large Euclidean momenta, loop integrals become superconvergent. For the self-energy Π( p ), the integral takes the form: Π(p)∼Zd4kPoly(k, p) k2(p−k)2e−k2 Λ2e−(p−k)2 Λ2.(18) The Gaussian factor dominates any polynomial growth in the numerator, rendering the integral finite in the UV. The only remaining divergences are infrared ones, handled by standard soft-graviton theorems. 2. Beta Function and Asymptotic Freedom Since the theory is finite for L≥ 1loops (superrenormalizable), the beta function for Newton’s constant βG is determined solely by finite one-loop shifts. In the trans-Planckian limit ( p≫MP l ), the beta function vanishes exponentially ( βG→ 0). This implies that the effective coupling Geff ( p )ceases to run at high energies, and the theory becomes Asymptotically Free. VI. CONCLUSION This paper has established a rigorous operational framework for Non-local Quantum Geometrodynamics. By replacing intuitive arguments with formal analyses (Dressed Cutkosky Rules, BRST Symmetry, Barnaby-Kamran analysis), we have demonstrated that the model is perturbatively unitary, gauge-invariant, and ghost-free. The phenomenology is robust, compatible with GW170817 constraints, and predicts regular Hayward black holes. Table I. Comparison of Quantum Gravity Models Feature GR R+R2Efimov (NLQG) Propagator 1/p21 p2(p2+M2) e−p2/M2 P l p2 UV Behavior Divergent Renorm. Finite (UV-finite) Unitarity Yes No (Ghosts) Yes (No-Ghost) Geometry Riemann Riemann Effective Weyl Singularities Inevitable Softened Removed Appendix A: One-Loop Perturbative Unitarity for Entire Form Factors In this appendix we present an explicit one-loop proof of perturbative unitarity for the class of entire form factors considered in the main text, namely V(z)=eH(z), z =−p2 M2 P l ,(A1)
5 where H ( z )is an entire function satisfying the following conditions: 1. H (0) = 0 (normalization of the massless graviton pole), 2. Re H(z)≥0for Re z≥0, 3. H ( z )grows at most polynomially in |z| in any finite angular sector of the complex plane. These conditions imply that V ( z )has no zeros in the finite complex plane and that the modified propagator ˜ Dµνρσ(p) = e−2H(−p2/M2 P l) p2+iϵ Pµνρσ (A2) has a unique physical pole at p2 = 0 with residue equal to that of the Einstein–Hilbert theory. 1. Euclidean Amplitudes and Analytic Continuation Working with Euclidean momenta pE and purely imaginary energies, loop amplitudes are given by ME(pi) = Zd4kE (2π)4N(kE, pi)Y r e−2H(−k2 E,r/M2 P l) k2 E,r , (A3) where N ( kE, pi )is a polynomial in loop and external momenta and kE,r denote the internal Euclidean momenta. By the assumptions above, the exponential form factors provide exponential damping for large |kE| , ensuring absolute convergence of the integral and allowing the interchange of integration and analytic continuation. 2. Discontinuity and Cut Diagrams at One Loop After analytic continuation to Minkowski signature, a generic one-loop 2 → 2graviton amplitude takes the form M(s, t) = Zd4k (2π)4N(k, p) (k2+iϵ)(p−k)2+iϵ ×e−2H(−k2/M2 P l)e−2H(−(p−k)2/M2 P l), (A4) with p2 = s . The discontinuity across the physical branch cut in sis DiscsM(s, t) = M(s+i0, t)−M(s−i0, t).(A5) Using 1 x+iϵ −1 x−iϵ =−2πi δ(x),(A6) and applying the Cutkosky procedure, we obtain DiscsM(s, t) = (−2πi)2Zd4k (2π)4δ(k2)δ(p−k)2 ×e−2H(0)e−2H(0) Non-shell(k, p).(A7) By normalization H (0) = 0, the form factors reduce to unity on shell: e−2H(0) = 1.(A8) Thus the discontinuity coincides with that of the local Einstein theory: DiscsM(s, t) = 2iIm M(s+i0, t) =X XZdΠX|Mtree(A→X)|2.(A9) 3. Optical Theorem and Microcausality The optical theorem 2Im M(A→A) = X XZdΠX|M(A→X)|2(A10) is therefore satisfied at one loop for this class of entire functions, and no additional thresholds or ghost cuts are introduced. Nonlocal form factors can in principle induce violations of microcausality at distances of order the nonlocality scale, but such effects are exponentially suppressed at scales much larger than the Planck length and do not spoil the perturbative S-matrix unitarity in the energy range accessible to current experiments. Appendix B: Hamiltonian Structure and Ostrogradsky Stability We outline here the Hamiltonian formulation of the nonlocal graviton theory at the linearized level and show that the number of propagating degrees of freedom coincides with that of GR, thereby avoiding the Ostrogradsky instability in the sense of the Dirac constraint analysis. 1. Linearized Nonlocal Action We expand the metric around Minkowski spacetime, gµν =ηµν +κhµν,(B1) and consider the quadratic nonlocal action S(2) =1 2Zd4x hµν e−H(□/M2 P l)Eµν ρσ hρσ,(B2) where Eµνρσ is the Lichnerowicz operator of linearized Einstein gravity. 2. Localization via an Auxiliary Field Following Barnaby and Kamran, we introduce an auxiliary field Φ µν ( x, s )defined on an extended space with
6 an extra “diffusion” coordinate s≥0: ∂sΦµν(x, s) = □Φµν(x, s),Φµν(x, 0) = hµν(x). (B3) For the Gaussian form factor H(z)=z, we have e−□/M2 P l hµν (x)=Φµνx, s =1 M2 P l ,(B4) and more general entire functions can be represented as suitable superpositions of such kernels. The quadratic action (B2) can then be written as a local action in (x, s): S(2) loc =1 2Zd4x ds [ΦµνEµνρσΦρσ +λµν (∂sΦµν −□Φµν )] , (B5) where λµν are Lagrange multipliers enforcing the diffusion equation. 3. Canonical Variables and Constraints Choosing x0 as the physical time, the canonical variables are Φµν(x, s),Πµν(x, s) = ∂Lloc ∂(∂0Φµν),(B6) together with lapse and shift variables from the ADM decomposition of hµν . The Lagrange multipliers λµν generate primary constraints πµν λ(x, s)≈0,(B7) while the diffusion equation yields secondary constraints χµν(x, s)≡∂sΦµν(x, s)−□Φµν(x, s)≈0.(B8) Together with the linearized diffeomorphism constraints inherited from GR, these form a closed first-class algebra (up to gauge-fixing), analogous to that of linearized Einstein gravity but extended along the s-direction. 4. Counting Degrees of Freedom and Stability The constraint surface and gauge orbits reduce the phase space to that of linearized GR: the number of propagating degrees of freedom remains Ndof = 2,(B9) corresponding to the two graviton polarizations. The Ostrogradsky instability is avoided because the nonlocality has been traded for an auxiliary diffusion direction and additional constraints, rather than for independent higher time derivatives in x0 . The reduced Hamiltonian on the constraint surface is bounded from below in the same way as in linearized GR, at least for the class of entire form factors considered here. Appendix C: GW170817 Constraints on the Speed of Gravitational Waves The multimessenger observation of the binary neutron star merger GW170817, together with its electromagnetic counterpart GRB 170817A, constrains the fractional difference between the speed of gravitational waves cg and the speed of light cto be cg c−1≲10−15.(C1) 1. Dispersion Relation from the Kinetic Operator The quadratic nonlocal action can be written schematically as S(2) =1 2Zd4x hµν f(□/M2 P l)Eµνρσ hρσ,(C2) with f ( z ) = e−H(z) . In momentum space, the inverse propagator is D−1(p2)=f(−p2/M2 P l)p2P.(C3) The dispersion relation follows from det D−1(p2)=0 ⇒f(−p2/M2 P l)p2= 0.(C4) Since f ( z )is entire and has no zeros for finite z , the only solution is p2= 0 ⇒ω2=| k|2,(C5) which implies cg= 1 at tree level. 2. Higher-Order Corrections and Frequency Scales Loop corrections may modify the effective kinetic operator as D−1 eff (p2)=p2"1+α1 p2 M2 P l +α2p2 M2 P l 2 +...#,(C6) with αn dimensionless. For gravitational waves with frequencies relevant to GW170817, ω∼103Hz, one finds p2 M2 P l ∼ℏω MP lc22 ≪10−40,(C7) so any corrections to the group velocity are suppressed by many orders of magnitude relative to the bound (C1) . Therefore, within the energy range probed by current gravitational-wave observations, the nonlocal dressing has a negligible impact on the propagation speed of gravitational waves and is fully compatible with GW170817.
7 Appendix D: Weyl Geometry Interpretation and Hayward Metric 1. Effective Weyl Geometry from Nonlocal Dressing A Weyl connection is defined by ˆ Γλ µν = Γλ µν −δλ µAν−δλ νAµ+gµνAλ,(D1) where Γ λ µν is the Levi-Civita connection of the metric gµν and Aµ is the Weyl gauge vector. The corresponding non-metricity tensor is Qλµν ≡ −ˆ ∇λgµν = 2Aλgµν.(D2) In our construction the nonlocal dressing of the metric can be written schematically as gµν =V□g M2 P l gµν.(D3) In a slowly varying background, where V can be approximated locally by a smooth function of a scalar argument, we write gµν ≈Ω2(x)gµν,Ω2(x)≡ V(x).(D4) Then ∇λgµν = (∂λln Ω2)gµν,(D5) so that comparison with (D2) gives Aλ=∂λln Ω.(D6) Thus, at scales where the nonlocal operator can be approximated by a conformal factor, the effective geometry probed by high-energy gravitons is that of a Weyl manifold with non-metricity determined by the logarithmic gradient of V. 2. Dimensional Analysis of the Hayward Metric The Hayward metric considered in the main text is ds2=−f(r)dt2+dr2 f(r)+r2dΩ2,(D7) with f(r)=1−2Mr2 r3+ 2Ml2,(D8) where M is the ADM mass and l is a length scale associated with the nonlocal smearing (of order the Planck length). In units G = c = 1, M , r and l all have dimensions of length, so (D8) is dimensionally consistent. For r→ ∞ we have f(r)=1−2M r+Ol2 r3,(D9) recovering the Schwarzschild behaviour at large radii. For r→0we find f(r)≈1−r2 l2,(D10) which corresponds to a de Sitter core with effective cosmological constant Λeff =3 l2.(D11) The Kretschmann scalar K = RµνρσRµνρσ remains finite at r = 0 for this metric, so the central singularity of Schwarzschild is replaced by a regular core characterized by the single additional scale l . This behaviour is naturally interpreted as a manifestation of the nonlocal smearing of the gravitational interaction at the Planck scale. [1] G. V. Efimov, Non-local quantum theory of the scalar field, Commun. Math. Phys. 57, 217 (1977). [2] V. A. Alebastrov and G. V. Efimov, A proof of the unitarity of S-matrix in a nonlocal quantum field theory, Commun. Math. Phys. 31, 1 (1973). [3] E. T. Tomboulis, Superrenormalizable gauge and gravitational theories, arXiv:hep-th/9702146. [4] L. Modesto, Super-renormalizable Quantum Gravity, Phys. Rev. D 86, 044005 (2012). [5] N. Barnaby and N. Kamran, Dynamics with infinitely many derivatives: The Initial value problem, JHEP 0802, 008 (2008). [6] B. P. Abbott et al. (LIGO/Virgo), GW170817: Observation of Gravitational Waves, Phys. Rev. Lett. 119, 161101 (2017). [7] S. A. Hayward, Formation and evaporation of regular black holes, Phys. Rev. Lett. 96, 031103 (2006). [8] F. Briscese and L. Modesto, Cutkosky rules and perturbative unitarity, Phys. Rev. D 99, 104043 (2019). [9] K. S. Stelle, Renormalization of Higher Derivative Quantum Gravity, Phys. Rev. D 16, 953 (1977). [10] C. Becchi, A. Rouet, and R. Stora, Commun. Math. Phys. 42, 127 (1975). [11] T. Biswas et al., Towards singularity and ghost free theories of gravity, Phys. Rev. Lett. 108, 031101 (2012). [12] J. Lindgren et al., J. Phys.: Conf. Ser. 2987, 012001 (2025).