Nonlinear massive gravity theory of geometric origin

Thibault Damour<sup>1,*</sup> and Tamanna Jain<sup>2,1,3,†</sup>

<sup>1</sup>Institut des Hautes Etudes Scientifiques, 91440 Bures-sur-Yvette, France

<sup>2</sup>LPENS, Département de physique, Ecole normale supérieure, Université PSL, Sorbonne Université,
Université Paris Cité, CNRS, 75005, Paris

<sup>3</sup>Department of Applied Mathematics and Theoretical Physics, University of Cambridge,
Wilberforce Road, CB3 0WA Cambridge, United Kingdom

<figure><img src="image_1.png" /></figure>

(Received 24 January 2025; accepted 7 March 2025; published 1 April 2025)

We study the number of propagating degrees of freedom, at nonlinear order, in torsion gravity theories, a class of modified theories of gravity that include a propagating torsion in addition to the metric. We focus on a three-parameter subfamily of theories (“torsion bigravity”) that contains, at linear order, only two physical excitations: a massless spin-2 one (with 2 degrees of freedom) and a massive spin-2 one (with 5 degrees of freedom). We study the dynamics of the massive spin-2 field in the limit where the torsion field decouples from the metric. The number of degrees of freedom of the torsion field is found to change, at nonlinear order, from five to nine.

DOI: 10.1103/PhysRevD.111.084001

# I. INTRODUCTION

The predictions of the theory of general relativity (GR) have been so far in outstanding agreement with gravitational observations and experiments over a very large range of scales and field strengths (for reviews, see, e.g., Chap. 21 in [1] or the book [2]). Nevertheless, since its inception, there has been a constant effort to study possible extensions of GR. The search for alternatives to GR is motivated by open questions both at ultraviolet scales (e.g., [3]) and at infrared scales (e.g., [4–6]). The existence of consistent alternative theories of gravity is also useful for suggesting tests of GR [7].

In the geometric approach to gravity, modifying GR can stem from extending the geometrical structure of spacetime. As early as 1923, Cartan [8,9] suggested to extend GR by formally keeping the Einstein-Hilbert action, but in a first-order formalism, considering the metric and a general (metric-preserving) “torsionfull” connection as independent fields (“Einstein-Cartan theory”). Cartan imagined that the torsion could be sourced by some material spin density. This idea was later made concrete by Weyl [10] who studied the coupling of (Dirac) spinor fields to gravity within a first-order formalism and pointed

out that, in this approach, the spinors provide a spin-density source $\tau^{\mu\nu\lambda}$ for the torsion proportional to $\bar{\psi}\gamma^{[\mu}\gamma^\nu\gamma^\lambda]\psi$. He also pointed out that in Einstein-Cartan theory the torsion does not propagate (i.e., is not a dynamical field), but is constrained to be proportional to its source $\tau^{\mu\nu\lambda}$, thereby leading (in a second-order formalism for the metric) to a contact term proportional to $\tau^{\mu\nu\lambda}\tau_{\mu\nu\lambda}$. The ideas of Cartan and Weyl were further extended and developed by Sciama, Kibble, and others (see, e.g., [11] for a review).

In this work, we study a class of geometric theories of gravitation, introduced in Refs. [12–15], which generalize the Einstein-Cartan theory in having a “dynamical” torsion field that propagates away from its source $\tau^{\mu\nu\lambda}$. We focus on the theories defined in Refs. [12–15] that contain (around flat-space backgrounds), in addition to Einstein-like massless spin-2 excitations, (even-parity) massive spin-2 excitations (embodied in the torsion field) as well as (pseudoscalar) massive spin-0 excitations. Such theories will hereafter be referred to as “torsion gravity” theories. Below, we shall focus on the subclass of torsion gravity theories, called “torsion bigravity” theories, which contain only massive and massless spin-2 excitations.

Torsion gravity theories were recently revived with the hope that the geometrically rooted massive spin-2 field they contain might define a theoretically healthy, and phenomenologically interesting, modification of gravity both in cosmological situations [16–20] and in strong-field regimes [21–24]. The hope was that a geometrically well-motivated origin for massive spin-2 excitations (through a dynamical torsion) would alleviate the problematic issues that arise when the extra massive spin-2 field is linked to the addition of a second metric tensor $f_{\mu\nu}$, besides $g_{\mu\nu}$, as in generic

Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.

<footnote>*Contact author: damour@ihes.fr</footnote>

<footnote>†Contact author: tj317@cam.ac.uk</footnote>
“bimetric gravity” theories. Let us recall that there are several problematic issues associated with the usual (Fierz-Pauli-type [25]) descriptions of massive spin-2 fields, notably: the Van Dam–Veltman discontinuity [26], the Vainshtein radius issue [27], and the additional, ghostlike, sixth degree of freedom (d.o.f.) [28].

Generic bimetric gravity theories (see, e.g., Ref. [29]), involving some general potential $V(g_{\mu\nu}, f_{\mu\nu})$, suffer from the arising, at the nonlinear level, of the unwanted, additional, Boulware-Deser ghostlike d.o.f., see, e.g., [30]. The Boulware-Deser ghost was avoided by generalizing the special potential of the “ghost-free” massive gravity theory of Ref. [31] to a bimetric gravity framework [32]. Note, however, that recent works (see [33] and references therein) suggest that the theory of massive spin-2 excitations described by the ghost-free action of Ref. [31] is only valid up to an energy scale which cannot be parametrically larger than the mass of the spin-2 field.

By contrast with bimetric gravity theories, it was shown in Refs. [16,17,19–24] that torsion gravity theories have the good feature of keeping the same number of d.o.f. in many situations, involving, notably, dynamical (but linear) torsion excitations around generic nonlinear Einstein backgrounds or fully nonlinear spherically symmetric (but stationary) torsion gravity solutions. In addition, it was also found [22] that the massive spin-2 torsion excitation does not introduce a Vainshtein radius in static, spherically symmetric nonlinear solutions. On the other hand, the study of time-dependent torsion perturbations suggested that torsion gravity theories are only valid up to curvature scales of order of the squared spin-2 mass [23,24].]

Motivated by the good features of torsion gravity theories found in previous studies, the aim of the present work is to investigate the number of propagating d.o.f. of torsion gravity theories in situations which have not yet been considered. Namely, we consider (vacuum) solutions which are both time dependent and nonlinear in the torsion field. To simplify our investigation, we focus on the nonlinear dynamics of the massive spin-2 (torsion) excitation in a well-defined limit where it decouples from the massless spin-2 (metric) one.

The paper is organized as follows. In Sec. II, we review the formalism and present the field equations for general torsion gravity theories along with reviewing the specialized three-parameter subfamily of torsion bigravity theories. The decoupling limit of the massive spin-2 field from the massless spin-2 field is presented in Sec. III. We first demonstrate the general, all-order structure of the field equations in the decoupling limit and then explicitly derive the field equations and investigate the constraints (degrees of freedom) at the linear and the next-to-linear order in Sec. IV. To explicitly understand the number of degrees of freedom for the field equations, we present several alternative action formulations of the decoupling limit in Sec. V. We then investigate the general structure of the constraints

in a resulting first-order (Hamiltonian) formalism in Sec. VI. The number of d.o.f. is derived at the linear level in Sec. VII and, finally, at the nonlinear level in Sec. VIII. Brief conclusions are presented in Sec. IX.

## II. FORMALISM OF TORSION GRAVITY

In this section, we discuss the formalism of torsion gravity theories following the notations of Ref. [20], as updated in Ref. [21]. The Latin indices $i, j... = 0, 1, 2, 3$ are used to denote vierbein indices linked to a local orthonormal frame ($e^i_\mu$) while Greek indices $\mu, \nu... = 0, 1, 2, 3$ are used to denote the spacetime indices related to the coordinate frame $x^\mu$. The frame indices are moved by the Minkowski metric $\eta_{ij}$ whereas the spacetime indices are moved by the metric $g_{\mu\nu} = \eta_{ij}e^i_\mu e^j_\nu$, where $e^i_\mu$ is the coframe, such that $e^i_\mu e^\mu_j = \delta^i_j$. The metric signature used is $(-, +, +, +)$.

The most general action with propagating torsion [expressed in terms of the coframe $e^i_\mu$ and a general SO (3,1) connection $A_{ij\mu}$] having as spectrum a massless spin-2 field, a massive spin-2 field, and a massive pseudoscalar field, is $S_{\text{TG}} + S_{\text{matter}}$, where the field part of the action reads

$$S_{\text{TG}}[e^i_\mu, A_{ijk}] = \int d^4 x \sqrt{-g} \mathcal{L}_{\text{TG}}[e, \partial e, \partial^2 e, A, \partial A], \quad (2.1)$$

with

$$ \begin{aligned} \mathcal{L}_{\text{TG}} &= c_R R + c_F F + c_2 + c_3 F^{ij} F_{ij} + c_4 F^{ij} F_{ji} \\ &\quad - \frac{1}{3} (c_3 + c_4) F^2 + c_6 (e \cdot F)^2 \\ &= c_R R + c_F F + c_2 + c_{F^2} \left( F_{(ij)} F^{(ij)} - \frac{1}{3} F^2 \right) \\ &\quad + c_{34} F_{[ij]} F^{[ij]} + c_6 (e \cdot F)^2. \end{aligned} \quad (2.2) $$

The second form of the latter action uses the symmetric $(F_{(ij)})$ and antisymmetric $(F_{[ij]})$ parts of $F_{ij} = F_{(ij)} + F_{[ij]}$, and the notation

$$c_{34} \equiv c_3 - c_4, \quad c_{F^2} \equiv c_3 + c_4. \quad (2.3)$$

Here, $c_R, c_F, c_2, c_3, c_4$, and $c_6$ are coupling constants, $R$ is the usual Einsteinian Ricci scalar of the metric $g_{\mu\nu}$, while $F_{ijkl}$, $F_{ij}$, and $F$ are the (frame components of the) curvature tensor, Ricci tensor ($F_{ij} = \eta^{kl}F_{kilj}$), and curvature scalar ($F = \eta^{ij}F_{ij}$) defined by the (torsionfull) connection $A_{ij\mu}$. The metric-preserving character of the SO (3,1) connection $A_{ijk} \equiv A_{ij\mu}e^i_\mu$ is encapsulated in the algebraic fact that $A_{ijk} = -A_{jik}$. The term $e \cdot F \equiv \epsilon_{ijkl}F^{ijkl}$ involves the Levi-Civita antisymmetric tensor $\epsilon_{ijkl}$. The explicit expression of the curvature tensor $F_{ijkl}$
and its corresponding avatars $F_{ij}$ and $F$ are

$$F_{ijkl} = e_k^\mu e_l^\nu (\partial_\mu A_{ij\nu} - \partial_\nu A_{ij\mu} + \eta^{mn} A_{im\mu} A_{nj\nu} - \eta^{mn} A_{im\nu} A_{nj\mu}), \quad (2.4)$$

$$F_{ij} = \eta^{kl} F_{ikjl}, \quad F = \eta^{ij} F_{ij}. \quad (2.5)$$

The curvature tensor $F_{ijkl}$, like the metric Riemann tensor $R_{ijkl}$, is antisymmetric under the exchanges $i \leftrightarrow j$ and $k \leftrightarrow l$, however, contrary to $R_{ijkl}$, it is not symmetric under the exchange $ij \leftrightarrow kl$.

The torsion gravity field equation obtained by varying the vierbein $e_\mu^i$ in $S_{\text{TG}} + S_{\text{matter}}$ reads

$$c_R \left( R_{ij} - \frac{1}{2}\eta_{ij}R \right) + c_F \left( F_{ij} - \frac{1}{2}\eta_{ij}F \right) - \frac{1}{2}c_2\eta_{ij} + c_3(F_i^k F_{kj} + F^{kl} F_{kilj}) + c_4(F_i^k F_{kj} + F^{lk} F_{kilj}) - \frac{2}{3}(c_3 + c_4)FF_{ij} + 2c_6\epsilon^{klmi}F_{klmj}(\epsilon \cdot F) - \frac{1}{2}\eta_{ij}L^{(2)} = \mathcal{T}_{ij}, \quad (2.6)$$

where $\mathcal{T}_{ij}$ is the gravitational matter source, and

$$L^{(2)} = c_3 F^{ij} F_{ij} + c_4 F^{ij} F_{ji} - \frac{1}{3}(c_3 + c_4)F^2 + c_6(\epsilon \cdot F)^2. \quad (2.7)$$

The field equation corresponding to varying the connection $A_{ijk}$ reads

$$\begin{aligned} & \left[ \eta_{ik} \left( D^m P_{jm} - \frac{2}{3} D_j P \right) - D_i P_{jk} \right] \\ & - \left[ \eta_{jk} \left( D^m P_{im} - \frac{2}{3} D_i P \right) - D_j P_{ik} \right] \\ & + 4c_6 \epsilon_{ijkm} D^m (\epsilon \cdot F) + H_{ijk} + S_{ijk} = \mathcal{T}_{ijk}, \end{aligned} \quad (2.8)$$

where $\mathcal{T}_{ijk}$ is the material torsion source, and where

$$H_{ijk} = c_F \left( K_{ikj} - K_{jki} - K_{il}^l \eta_{jk} + K_{jl}^l \eta_{ik} \right), \quad (2.9)$$

$$S_{ijk} = \frac{1}{c_F} H_{mnk} \left( \eta_{im} P_{jn} - \eta_{jm} P_{in} - \frac{2}{3} \eta_{im} \eta_{jn} P + 2c_6 \epsilon_{ijmn} (\epsilon \cdot F) \right), \quad (2.10)$$

$$P_{ij} = c_3 F_{ij} + c_4 F_{ji}, \quad P \equiv \eta_{ij} P_{ij}. \quad (2.11)$$

Here, $D_i$ denotes the covariant derivative with respect to the connection $A_{ijk}$: e.g., $D_i V^j = e_i^\mu (\partial_\mu V^j + A_{k\mu}^j V^k)$; and $K_{ijk}$

is the contorsion tensor, i.e., the difference between the connection $A_{ijk}$ and the metric-induced (Levi-Civita) connection $\omega_{ijk} \equiv \omega_{ij\mu} e_\mu^i$,

$$K_{ijk} = A_{ijk} - \omega_{ijk}. \quad (2.12)$$

The field equations (2.6) and (2.8) of torsion gravity generally depend on the six coupling parameters: $c_F$, $c_R$, $c_2$, $c_3$, $c_4$, and $c_6$.

In the present work, we follow Ref. [21] in focusing on the three-parameter subfamily of torsion bigravity theories where

$$c_2 = 0, \quad c_3 = c_4, \quad \text{and} \quad c_6 = 0, \quad (2.13)$$

such that, in the notation of Eq. (2.3), we have

$$c_{34} = 0, \quad c_3 = c_4 = \frac{1}{2} c_F^2. \quad (2.14)$$

These theories have only a massless spin-2 and a massive spin-2 in their spectrum. The corresponding torsion bigravity (TBG) Lagrangian takes the simple form

$$\mathcal{L}_{\text{TBG}} = c_R R + c_F F + c_F^2 \left( F_{(ij)} F^{(ij)} - \frac{1}{3} F^2 \right). \quad (2.15)$$

The corresponding field equations simplify. In particular, Eq. (2.8) now reads

$$\begin{aligned} & \left[ \eta_{ik} \left( D_m P_{jm} - \frac{2}{3} D_j P \right) - D_i P_{jk} \right] \\ & - \left[ \eta_{jk} \left( D_m P_{im} - \frac{2}{3} D_i P \right) - D_j P_{ik} \right] + H_{ijk} + S_{ijk} = \mathcal{T}_{ijk}, \end{aligned} \quad (2.16)$$

where $P_{ij} = c_F^2 F_{(ij)}$, and $H_{ijk}$ and $S_{ijk}$ are as defined above (with $c_6 = 0$). Note that the latter field equation only involves the symmetric part $F_{(ij)}$ of $F_{ij}$. In this work, we focus on vacuum field equations, in absence of material sources (i.e., with $\mathcal{T}_{ij} = 0$ and $\mathcal{T}_{ijk} = 0$).

Our restriction to a subfamily of torsion gravity theories will be of sufficient generality for our purpose because we shall show below that the subfamily (2.15), which showed healthy properties in previously considered situations, becomes unhealthy when considered at the nonlinear level.

### III. DECOUPLING GRAVITY FROM THE TORSION DYNAMICS

It will be convenient here to further simplify our task by focusing on the nonlinear dynamics of the massive spin-2 excitation contained in the connection $A_{ijk}$ by considering a limit where the dynamics of $A_{ijk}$ decouples from the dynamics of the metric $g_{\mu\nu}$. [Here, we assume that the
vierbein $e_\mu^i$ has been (locally) algebraically expressed in terms of the metric.] To this effect let us write the metric as

$$g_{\mu\nu} = \eta_{\mu\nu} + \frac{1}{\sqrt{c_R}} \hat{h}_{\mu\nu}. \qquad (3.1)$$

Inserting the latter parametrization of the metric in the Einsteinian part of the action density $L_g = \sqrt{-g} c_R R(g)$ yields (after an integration by parts) an action density of the form

$$L_g = (\partial\hat{h})^2 + \frac{1}{\sqrt{c_R}} \hat{h}(\partial\hat{h})^2 + \frac{1}{c_R} \hat{h}^2 (\partial\hat{h})^2 + \dots \qquad (3.2)$$

Therefore, in the limit $c_R \to \infty$, the metric perturbation $\hat{h}_{\mu\nu}$ (which also couples to the matter and to other fields proportionally to $1/\sqrt{c_R}$) decouples as a free massless spin-2 excitation.

As a consequence, in the limit $c_R \to \infty$ (henceforth referred to as the “decoupling limit”), the dynamics of the connection $A_{ijk}$ [equal in this limit to the contorsion tensor $K_{ijk}$ because $g_{\mu\nu} = \eta_{\mu\nu} + \mathcal{O}(\frac{1}{\sqrt{c_R}})$, $e_\mu^i = \delta_\mu^i + \mathcal{O}(\frac{1}{\sqrt{c_R}})$, and $\omega_{ijk} = \mathcal{O}(\frac{1}{\sqrt{c_R}})$] takes place within the flat background spacetime $\eta_{\mu\nu}$ (endowed with a flat vierbein field $e_\mu^i = \delta_\mu^i$ with corresponding Lorentzian coordinates $x^i = x^\mu$).

Below, we study the resulting decoupled (and sourceless) field equations for $A_{ijk}$. The crucial point for our study is that we are retaining in this limit all the nonlinearities in $A_{ijk}$ entailed by the action

$$\mathcal{L}_{\text{TBG}}^{\text{decoupled}}(A) = c_F F(A) + c_{F^2} \left( F_{(ij)}(A) F^{(ij)}(A) - \frac{1}{3} (F(A))^2 \right), \qquad (3.3)$$

where $F(A) \sim \partial A + AA$, namely

$$F_{jl} = \eta^{ik} (\partial_k A_{ijl} - \partial_l A_{ijk} + \eta^{mn} A_{imk} A_{njl} - \eta^{mn} A_{iml} A_{njk}), \\ F = \eta^{jl} F_{jl}. \qquad (3.4)$$

Here, all partial derivatives are with respect to the Lorentzian coordinates $x^i$.

The action $\mathcal{L}_{\text{TBG}}^{\text{decoupled}}(A)$ is “quartic” in $A_{ijk}$. It has a formal similarity with a Yang-Mills action, though it differs from what would be the normal Yang-Mills action for a SO(3,1) connection.

## IV. STRUCTURE OF THE EQUATIONS OF MOTION OF TORSION BIGRAVITY IN THE DECOUPLED LIMIT

A direct derivation of the field equations from the action (3.3) would yield second-order equations for the 24 components of the connection $A_{ijk} = -A_{jik}$. However, following Refs. [17,19], the field equations for $A_{ijk}$ can be conveniently reexpressed as field equations for the symmetric part $F_{(ij)}$ of the Ricci tensor of $A_{ijk}$.

Indeed, taking into account that $D \sim \partial + A$, $P_{ij} \sim c_{F^2} F_{(ij)}$, $H_{ijk} \sim c_F A_{ijk}$, and $S_{ijk} \sim \frac{1}{c_F} A_{ijk} \cdot F_{(lm)}$, the torsion field equations (2.16) have the following structure (indicating only the number of indices):

$$c_{F^2}(\partial_1 + A_3)F_{(2)} + c_F A_3 + c_{F^2}A_3 F_{(2)} = 0. \qquad (4.1)$$

Let us introduce the notation

$$\kappa^2 \equiv \frac{c_F}{c_{F^2}}, \qquad (4.2)$$

to denote the squared mass of the massive spin-2 excitation. [More precisely, Eq. (4.2) defines the squared mass in the limit where $c_F \to 0$. When $c_R$ is finite, we would have $\kappa^2 = \frac{c_F}{c_{F^2}} (1 + \frac{c_F}{c_R})$.] The structure of the above field equation then reads

$$\partial_1 F_{(2)} + (\kappa^2 + F_{(2)}) A_3 = 0. \qquad (4.3)$$

This equation can be formally solved in terms of $A_3$ to yield

$$A_3 = -(\kappa^2 + F_{(2)})^{-1} \partial_1 F_{(2)}, \qquad (4.4)$$

where the formal inverse $(\kappa^2 + F_{(2)})^{-1}$ is perturbatively defined as an all-order expansion in powers of $\frac{F_{(2)}}{\kappa^2}$. We shall not study here the singularities that might appear in this inverse when considering strong fields $\frac{F_{(2)}}{\kappa^2} = \mathcal{O}(1)$. (The occurrence of similar singularities, leading to a lower bound on the mass $\kappa$, were pointed out in Refs. [23,24].)

Formally inserting the expression (4.4) for $A_3$ in the definition $F_{(2)} \sim \partial A_3 + A_3 A_3$ of $F_{(2)}$ then leads to a non-linear (Fierz-Pauli-like [19]) second-order equation for the dimensionless symmetric tensor field

$$\hat{F}_{(ij)} \equiv \frac{F_{(ij)}}{\kappa^2}, \qquad (4.5)$$

having the symbolic structure

$$\partial \left( \frac{\partial \hat{F}}{1 + \hat{F}} \right) + \left( \frac{\partial \hat{F}}{1 + \hat{F}} \right)^2 + \kappa^2 \hat{F} = 0. \qquad (4.6)$$

We have computed the explicit form of this derived field equation for $\hat{F}_{(ij)}$ up to the quadratic order in $\hat{F}_{(ij)}$. Let us
display our results when using the variable $u_{ij} = u_{(ij)}$, defined as

$$u_{ij} = \frac{1}{\kappa^2} \left( F_{(ij)} - \frac{1}{6}\eta_{ij}F \right). \qquad (4.7)$$

For brevity, we do not henceforth explicitly indicate that $u_{ij}$ is a symmetric tensor. We have followed here Refs. [17,19] in defining $u_{ij}$ so that it satisfies, at the linear level, the usual Fierz-Pauli equation. We shall show below how to derive Eq. (4.6) from an action.

### A. Linear order equation for $u_{ij}$

The first step to obtain the linear equation for $u_{ij}$ is to invert the (linear) relation (2.11) between $H_{ijk}$ and $K_{ijk} = A_{ijk}$. An easy computation yields

$$A_{ijk} = \frac{1}{2c_F} \left[ H_{ijk} - H_{jki} - H_{kij} - \frac{1}{2}\eta_{jk} \left( 2H_{il}^l - H_{li}^l \right) + \frac{1}{2}\eta_{ik} \left( 2H_{jl}^l - H_{lj}^l \right) \right]. \qquad (4.8)$$

Inserting Eq. (2.16) in this expression, then yields the explicit form of the equation sketched in Eq. (4.3) above. Proceeding as just indicated and working at the linear level in $u_{ij}$ yields the linear equation

$$U_{ij} = 0, \qquad (4.9)$$

where

$$U_{ij} = \partial_k \partial^k u_{ij} + \partial_i \partial_j u - \partial_i \partial_k u_j^k - \partial_j \partial_k u_i^k + \eta_{ij} \left( \partial_k \partial_l u^{kl} - \partial_k \partial^k u \right) - \kappa^2 (u_{ij} - \eta_{ij} u). \qquad (4.10)$$

with an effective quadratically nonlinear source given by

$$\begin{aligned} U_{ij}^{\text{quad}}(u) = & \partial_i u_{kl} \partial_j u^{kl} - \partial_i u \partial_j u + u^{kl} \partial_i \partial_j u_{kl} - u \partial_i \partial_j u + \partial^k u_{ik} \partial_j u + \partial^k u_{jk} \partial_i u \\ & + \frac{1}{2} (u_{jk} \partial^k \partial_i u + u_{ik} \partial^k \partial_j u) - \partial_k u \partial^k u_{ij} - \partial^k u_{ik} \partial^l u_{jl} + \frac{1}{2} (\partial_i u_{jk} \partial_l u^{kl} + \partial_j u_{ik} \partial_l u^{kl} - u^{kl} \partial_l \partial_i u_{jk} - u^{kl} \partial_l \partial_j u_{ik}) \\ & + u \partial^l \partial_i u_{jl} + u \partial^l \partial_j u_{il} - u_{jk} \partial_l \partial_i u^{kl} - u_{ik} \partial_l \partial_j u^{kl} + u_{ij} \partial_l \partial_k u^{lk} - u \partial_l \partial^l u_{ij} + u_j^k \partial_l \partial^l u_{ik} + u_i^k \partial_l \partial^l u_{jk} - u_{ij} \partial_l \partial^l u \\ & + \partial^k u_{jl} \partial^l u_{ik} + \partial^k u_{jl} \partial^k u_{il} - \frac{1}{2} (u_{jk} \partial^l \partial^k u_{il} + u_{ik} \partial^l \partial^k u_{jl}) - \frac{3}{2} (\partial_j u^{kl} \partial_l u_{ik} + \partial_i u^{kl} \partial_l u_{jk}) \\ & - \eta_{ij} \left[ \frac{1}{2} (\partial^n u^{kl} \partial_n u_{kl} + \partial_k u^{kl} \partial^n u_{nl} - \partial^l u \partial_l u - \partial_l u_{kn} \partial^n u^{kl}) - u \partial_k \partial^k u - u^{kl} \partial^n \partial_l u_{nk} + u \partial_n \partial_l u^{nl} + u^{kl} \partial_n \partial^n u_{kl} \right]. \end{aligned} \qquad (4.15)$$

The equation $U_{ij} = 0$ is the usual form of the (linear) Fierz-Pauli equation [25]. Let us recall that one way to prove that this equation describes the five d.o.f. of a massive spin-2 field is to take the divergence and the trace of $U_{ij} = 0$. These yield the constraints

$$\partial^j U_{ij} \equiv \kappa^2 (\partial_i u - \partial^j u_{ij}) = 0, \qquad (4.11)$$

$$\eta^{ij} U_{ij} + \frac{2}{\kappa^2} \partial^{ij} U_{ij} \equiv -3\kappa^2 u = 0. \qquad (4.12)$$

These $4 + 1$ constraints (last equations on the right) involve only $u$ and $\partial u$ and are equivalent to the following five constraints:

$$\partial^j u_{ij} = 0, \qquad u = 0, \qquad (4.13)$$

which ensure the correct number (five) of d.o.f. at the linear level.

The main aim of our present study is to investigate whether this number of d.o.f. is preserved at the non-linear level.

### B. Next-to-linear order

As a first attempt to investigate the healthiness of the nonlinear massive gravity theory described by the sketchy equation (4.6) we explicitly computed the terms quadratic in $F_{(ij)}$ (i.e., in $u_{ij}$) in Eq. (4.6) by following the procedure explained above. We obtained

$$\partial_k \partial^k u_{ij} + \partial_i \partial_j u - \partial_i \partial_k u_j^k - \partial_j \partial_k u_i^k + \eta_{ij} (\partial_k \partial_l u^{kl} - \partial_k \partial^k u) - \kappa^2 (u_{ij} - \eta_{ij} u) = U_{ij}^{\text{quad}}(u) + \mathcal{O}(u^3), \qquad (4.14)$$
If the nonlinear massive gravity theory described by Eq. (4.14) had the same number of d.o.f. as its linear counterpart, one would expect Eq. (4.14) to give rise again to five quadratically deformed first-order constraints of the type

$$ \partial^j u_{ij} = \mathcal{O}_i^{\text{quad}}(u, \partial u), \quad u = \mathcal{O}^{\text{quad}}(u, \partial u). \quad (4.16) $$

One would expect that such quadratically deformed first-order constraints can be revealed by replacing $\mathcal{U}_{ij}(u)$ on the left-hand sides of Eqs. (4.11) and (4.12) by its quadratically deformed avatar $\mathcal{U}_{ij}(u) - \mathcal{U}_{ij}^{\text{quad}}(u)$. If that were the case, the quadratic quantities

$$ C_i^{\text{quad}} \equiv \partial^j \mathcal{U}_{ij}^{\text{quad}}, \quad (4.17) $$

$$ C^{\text{quad}} \equiv \eta^{ij} \mathcal{U}_{ij}^{\text{quad}} + \frac{2}{\kappa^2} \partial^{ij} \mathcal{U}_{ij}^{\text{quad}} \quad (4.18) $$

would be expected to vanish when working on linear shell $[\mathcal{U}_{ij}(u) = 0]$ and imposing the linear-level constraints given in Eq. (4.13).

However, we found that the so-defined quantities $C_i^{\text{quad}}$, $C^{\text{quad}}$ do not vanish when imposing Eq. (4.13), but, actually, involve up to three derivatives of $u$.

This result suggests that the (second-order) nonlinear massive gravity equations satisfied by $u_{ij}$ propagate more than five d.o.f. or exhibit some worse inconsistency. However, this reasoning does not allow one to understand the real number of d.o.f. of the exact field equations for $u_{ij}$ whose structure was indicated in Eq. (4.6). In addition, the reasoning leading to the derivation of Eq. (4.6) did not allow us to prove that this second-order field equation for $F_{(ij)}$ (or $u_{ij}$) could be derived from an action. In the following sections, we shall therefore take a different tack (directly based on action principles).

## V. ACTION FORMULATIONS OF TORSION BIGRAVITY

In the decoupling limit $c_R \to \infty$ [with the Eq. (3.1) assumption] considered here, the torsion bigravity theory is originally defined by a second-order action depending only on the connection $A_3 = A_{ijk} = A_{[ij]k}$, say [see, Eq. (3.3)],

$$ \mathcal{L}_{\text{TBG}}[A_3] = c_F^2 \left( F_{(ij)}(A, \partial A) F^{(ij)}(A, \partial A) - \frac{1}{3} F^2(A, \partial A) \right) + c_F F(A, \partial A). \quad (5.1) $$

Let us now show that torsion bigravity also admits two other useful action formulations: (1) a first-order Lagrangian $\mathcal{L}_1(b_2, A_3)$ involving both $A_3$ and a symmetric two-tensor $b_2 = b_{(ij)}$ (equivalent to $F_{(ij)}$ or $u_{ij}$); and (2) a second-order Lagrangian $\mathcal{L}_2(b_2)$ involving only $b_2 = b_{(ij)}$

[and leading to the second-order equation of motion given in Eq. (4.6)].

These alternative action formulations are derived from the original action Eq. (5.1) by introducing both a Lagrange multiplier $b_{ij} = b_{(ij)}$ and an auxiliary independent symmetric field $\mathcal{F}_{ij}$.

Indeed, we can construct a Lagrangian depending on the three fields $[b_2, A_3, \mathcal{F}_2]$ such that $b_{ij}$ acts as a Lagrange multiplier imposing the constraint $\mathcal{F}_{ij} = F_{ij}(A)$. A Lagrangian satisfying these requirements is

$$ \mathcal{L}_0[b_2, A_3, \mathcal{F}_2] = c_b b^{ij} (\mathcal{F}_{ij} - F_{ij}(A, \partial A)) + c_F \mathcal{F} + c_{F^2} \left( \mathcal{F}_{ij} \mathcal{F}^{ij} - \frac{1}{3} \mathcal{F}^2 \right), \quad (5.2) $$

where $c_b$ is a constant (to be chosen below), and $\mathcal{F} \equiv \mathcal{F}_{ij}\eta^{ij}$. The variation of $\mathcal{L}_0$ with respect to $b_{ij}$ gives the constraint $\mathcal{F}_{ij} = F_{ij}(A)$. We can use this constraint to eliminate both $b_{ij}$ and $\mathcal{F}_{ij}$ from the Lagrangian. This then reduces $\mathcal{L}_0[b_2, A_3, \mathcal{F}_2]$ to the original Lagrangian $\mathcal{L}_{\text{TBG}[A_3]}$.

On the other hand, as the action $\mathcal{L}_0[b_2, A_3, \mathcal{F}_2]$, Eq. (5.2), is quadratic in the independent auxiliary field $\mathcal{F}_{ij}$, we can instead “integrate out” $\mathcal{F}_{ij}$ by solving its equation of motion $\frac{\delta\mathcal{L}_0}{\delta\mathcal{F}_{ij}} = 0$. The latter equation of motion reads

$$ \frac{\delta \mathcal{L}_0}{\delta \mathcal{F}_{ij}} = c_b b_{ij} + c_F \eta_{ij} + 2 c_{F^2} \left( \mathcal{F}_{ij} - \frac{1}{3} \eta_{ij} \mathcal{F} \right) = 0. \quad (5.3) $$

It is convenient to choose

$$ c_b = -c_F, \quad (5.4) $$

to ensure that $b_{ij} \to \eta_{ij}$ when $\mathcal{F}_{ij} \to 0$. Doing this, the field equation for $\mathcal{F}_{ij}$ becomes (recalling $\mathcal{F} \equiv \eta^{ij}\mathcal{F}_{ij}$)

$$ \mathcal{F}_{ij} - \frac{1}{3} \eta_{ij} \mathcal{F} = \frac{\kappa^2}{2} (b_{ij} - \eta_{ij}). \quad (5.5) $$

Defining for notational convenience

$$ \bar{b}_{ij} \equiv b_{ij} - \eta_{ij}, \quad (5.6) $$

the above equation yields (with $\bar{b} \equiv \eta^{ij}\bar{b}_{ij}$)

$$ \mathcal{F}_{ij} = \frac{\kappa^2}{2} (\bar{b}_{ij} - \eta_{ij} \bar{b}), \quad \mathcal{F} = -\frac{3\kappa^2}{2} \bar{b}. \quad (5.7) $$

Inserting Eq. (5.7) in Eq. (5.2), we then obtain the following first-order action involving $b_2$ and $A_3$:

$$ \mathcal{L}_1[b_2, A_3] = c_F b^{ij} F_{ij}(A) - \frac{c_F \kappa^2}{4} (\bar{b}_{ij} \bar{b}^{ij} - \bar{b}^2). \quad (5.8) $$
This first-order action will be used below to investigate the d.o.f. of torsion bigravity in the decoupling limit.

The Lagrangian given in Eq. (5.8) is linear in $F_2(A) = \partial A_3 + A_3 A_3$, and therefore quadratic in $A_3$. We can therefore integrate out $A_3$ from the Lagrangian $\mathcal{L}_1[b_2, A_3]$, Eq. (5.8), by solving its equation of motion $\frac{\delta \mathcal{L}_1}{\delta A_3} = 0$.

Using the explicit expression of $F_2(A) = \partial A_3 + A_3 A_3$, namely

$$F_{ij}(A) = \frac{1}{2}(\partial_j A_{il} - \partial^l A_{ilj} + A_{imj} A^{ml} - A_i^{ml} A_{mlj} + i \leftrightarrow j), \quad (5.9)$$

we find that the field equation for $A_3$ reads

$$\frac{\delta \mathcal{L}_1}{\delta A_{ijk}} = \frac{1}{2} c_F [\beta_{ijk} - (b_2 * A_3)_{ijk}] = 0, \quad (5.10)$$

where we define the following rank-3 tensors having the same symmetries as $A_{ijk} = A_{[ij]k}$,

$$\beta_{ijk} \equiv \partial_i b_{jk} - \eta_{ik} \partial^l b_{jl} - (i \leftrightarrow j), \quad (5.11)$$

and

$$(b_2 * A_3)_{ijk} \equiv b_{ik} A_{jl}^l - b_{il} A_{jk}^l + b_{kl} A_{i^l}^j + \eta_{ik} b^{lm} A_{jlm} - (i \leftrightarrow j). \quad (5.12)$$

The field equation of $A_3$ is then linear (and algebraic) in $A_3$ and has the structure

$$(b_2 * A_3)_{ijk} = \beta_{ijk}, \quad (5.13)$$

where $\beta_{ijk}$ is the (curl-like) linear combination of the first derivatives of $b_{ij}$ defined in Eq. (5.11).

The structure of the field equation (5.13) is reminiscent of the Lyapunov matrix equation

$$BX + XB = C, \quad (5.14)$$

where $X$ is an unknown matrix (playing the role of $A_3$) while the matrix $B$ plays the role of $b_2$ and the right-hand side matrix $C$ plays the role of $\beta_3 \sim \partial_1 b_2$. The solution of this Lyapunov toy model is easily obtained if the matrix $B$ is diagonalizable. Indeed, when working in a frame where $B_{ij} = b_i \delta_{ij}$, the equation for $X_{ij}$ reads $(b_i + b_j)X_{ij} = C_{ij}$, which is solved as

$$X_{ij} = \frac{C_{ij}}{(b_i + b_j)}, \quad (5.15)$$

if there are no vanishing denominators.

We can morally think of the matrix denominators $(b_i + b_j)^{-1}$ as defining some new “inverse” $B^{-1}$ of the

matrix $B$. In fact, if the matrix $B$ is close to the unit matrix, say $B_{ij} = \delta_{ij} + \bar{b}_{ij}$ (which is similar to our situation if we work in the Euclidean signature), we can expand this special-frame solution in powers of the eigenvalues $\bar{b}_i$ of $\bar{b}_{ij}$ and thereby easily rewrite the general-frame solution $X_{ij}$ in terms of a series of terms made of various contractions of $C_{ij}$ with several factors involving $\bar{b}_{ij}$. Suppressing indices, signs, and factors, we have the structure

$$X_2 \sim C_2 + \bar{b}_2 * C_2 + \bar{b}_2 * \bar{b}_2 * C_2 + \dots \quad (5.16)$$

It is easily seen that the obtention of such a perturbative solution to the Lyapunov equation (5.14) can be extended to the linear $A_3$ equation (5.13) when expanding its solution in powers of $\bar{b}_{ij} = b_{ij} - \eta_{ij}$. Namely, the solution for $A_3$ will have the perturbative structure [valid for $\bar{b}_{ij} < O(1)$]

$$A_3 \sim \beta_3 + \bar{b}_2 * \beta_3 + \bar{b}_2 * \bar{b}_2 * \beta_3 + \dots \quad (5.17)$$

Actually, the successive terms in this expansion can be straightforwardly derived by solving the explicit form of the $A_3$ equation, namely,

$$\beta_{ijk} = [b_{ik}A_{jl}^l - b_{il}A_{jk}^l + b_{kl}A_{i^l}^j + \eta_{ik}b^{lm}A_{jlm} - (i \leftrightarrow j)], \quad (5.18)$$

in an expansion in powers of $\bar{b}_{ij} = b_{ij} - \eta_{ij}$.

The first term, i.e. the linear order (LO), in this expansion is obtained by replacing $b_{ij} \to \eta_{ij}$ on the right-hand side of Eq. (5.18). This yields

$$A_{ijk}^{LO} = \frac{1}{2}(\beta_{ijk} - \beta_{kij} + \beta_{jki}), \quad (5.19)$$

which is actually of order $A^{LO} \sim \partial \bar{b}$. We also explicitly computed the next-to-linear order (NLO) term in the expansion

$$A_3[\bar{b}_2, \partial \bar{b}_2] = A_3^{LO} + A_3^{NLO} + \dots \sim \partial \bar{b}_2 + \bar{b}_2 \partial \bar{b}_2 + \dots \quad (5.20)$$

We can morally think that this solution is the massive gravity analog of the solution expressing the Levi-Civita connection in terms of the (inverse of the) metric $g_2$ and of its first derivatives $\partial g_2$. We can also formally write its structure as

$$A_3[\bar{b}_2, \partial \bar{b}_2] \sim (b_2)^{-1} \beta_3 \sim (b_2)^{-1} \partial_1 b_2. \quad (5.21)$$

Integrating out $A_3$ from the action $\mathcal{L}_1[b_2, A_3]$, Eq. (5.8), simply means replacing $A_3$ in $\mathcal{L}_1[b_2, A_3]$ (after an integration by parts to shuffle the derivatives from $A_3$ to $b_2$) by the solution Eq. (5.21). This yields a second-order Lagrangian for torsion bigravity depending only on $b_2$, and having the structure
$$ \frac{\mathcal{L}_2[b]}{c_F} \sim \partial b * b^{-1} * \partial b + b * b^{-1} * \partial b * b^{-1} * \partial b + \kappa^2 \bar{b}^2. \qquad (5.22) $$

The kinetic terms in this Lagrangian are defined by a power series expansion in $\bar{b}_{ij}$, while the mass terms are quadratic in $\bar{b}_{ij}$ and exactly given by the last term in Eq. (5.8). Namely,

$$ \frac{\mathcal{L}_2^{\text{mass term}}[b]}{c_F} = -\frac{1}{4} \kappa^2 (\bar{b}_{ij} \bar{b}^{ij} - \bar{b}^2). $$

Note in passing that the present reasoning shows that the nonlinear equation (4.6) can indeed be derived from an action, namely Eq. (5.22).

Forgetting that it was obtained by decoupling a dynamical torsionfull connection from an associated dynamical metric, the second-order action, Eq. (5.22), considered as a Poincaré-invariant action in Minkowski spacetime for $b_{ij}$ defines a specific nonlinear massive gravity theory, which is quite different from the usually considered massive gravity theories. Indeed, the latter [29,34] define the kinetic term of $b_{ij}$ by means of the Einstein-Hilbert action, and then look for a nonlinear deformation of the Fierz-Pauli

mass term having good properties. By contrast the Einstein-Cartan origin of the dynamics of

$$ b_{ij} = \eta_{ij} + \frac{2}{\kappa^2} \left( \mathcal{F}_{ij} - \frac{1}{3} \mathcal{F} \eta_{ij} \right) \qquad (5.23) $$

has led to the nonlinear action $\mathcal{L}_2[b]$, Eq. (5.22), which features a mass term of the usual quadratic Fierz-Pauli form, but with specific, non-Einsteinian kinetic terms.

We have explicitly computed the second-order action $\frac{\mathcal{L}_2[b]}{c_F}$ up to cubic order in $\bar{b}_{ij}$, i.e.,

$$ \frac{\mathcal{L}_2[b]}{c_F} = \frac{\mathcal{L}_2^{\text{LO}}[b]}{c_F} + \frac{\mathcal{L}_2^{\text{NLO}}[b]}{c_F} + \dots $$

The leading-order action (quadratic in $\bar{b}_{ij}$) is

$$ \frac{\mathcal{L}_2^{\text{LO}}[b]}{c_F} = \frac{1}{8} \partial_k \bar{b}_i^i \partial^k \bar{b}_j^j + \frac{1}{4} \partial_j \bar{b}^{jk} \partial^i \bar{b}_{ik} + \frac{1}{4} \partial_k \bar{b}_{ij} \partial_i \bar{b}^{jk} \\ - \frac{1}{4} \partial_i \bar{b}^{jk} \partial^i \bar{b}_{jk} - \frac{1}{4} \kappa^2 (\bar{b}_{ij} \bar{b}^{ij} - \bar{b}_i^i \bar{b}_j^j), \qquad (5.24) $$

while the next-to-leading-order one (cubic in $\bar{b}_{ij}$) reads

$$ \frac{\mathcal{L}_2^{\text{NLO}}[b]}{c_F} = \frac{1}{8} \bar{b}_{ij} \left( \beta^{jkl} \beta^{kli} - \frac{1}{2} \beta^{klj} \beta^{kli} - \beta^{jki} \beta_{kl}^l \right) + \frac{1}{8} \bar{b}_{ij} \beta_{kl}^l \partial^i \bar{b}^{jk} + \frac{1}{8} \bar{b}^{ij} \beta_{jk}^k \partial^l \bar{b}_{il} - \frac{1}{8} \bar{b}_i^i \beta_{jk}^k \partial_l \bar{b}^{lj} \\ + \frac{1}{8} \bar{b}_{ij} \beta^{ikj} \partial_k \bar{b}_l^l - \frac{1}{4} \bar{b}^{ij} \beta_{ikj} \partial_l \bar{b}^{lk} - \frac{1}{2} \bar{b}_{ij} \beta_{l}^{ki} \partial_k \bar{b}^{lj} - \frac{1}{8} \bar{b}_{ij} \beta^{kl} \partial_k \bar{b}^{lj}. \qquad (5.25) $$

We have checked that the equations of motion for $\bar{b}_{ij}$ derived from the above Lagrangian agree with the equations of motion (4.6). To check the agreement, one must take into account the relation between $\bar{b}_{ij}$ and $F_{ij} = \mathcal{F}_{ij}$ given by Eq. (5.7). In particular, the latter relation implies that $\bar{b}_{ij}$ is related to the Fierz-Pauli usual variable $u_{ij}$ via

$$ \bar{b}_{ij} = u_{ij} - \frac{1}{2} \eta_{ij} u. \qquad (5.26) $$

Note in passing the curious fact that the Fierz-Pauli mass term has the same expression in terms of $\bar{b}_{ij}$ and in terms of $u_{ij}$, because

$$ \bar{b}_{ij} \bar{b}^{ij} - \bar{b}^2 \equiv u_{ij} u^{ij} - u^2. \qquad (5.27) $$

## VI. STRUCTURE OF THE FIRST-ORDER ACTION

In this section, we finally tackle in a precise manner the issue of the number of d.o.f. of torsion bigravity at the nonlinear level by using the first-order Lagrangian $\mathcal{L}_1[b, A]$, Eq. (5.8), namely,

$$ \frac{\mathcal{L}_1[b_2, A_3]}{c_F} = b^{ij} F_{ij}(A) - \frac{\kappa^2}{4} (\bar{b}_{ij} \bar{b}^{ij} - \bar{b}^2). \qquad (6.1) $$

In order to reach clear conclusions, we shall focus here on the particular case where [after having performed the field transformations leading to Eq. (6.1)] all the variables entering the action $\mathcal{L}_1[b_2, A_3]$ are homogeneous in space and only depend on time, $t = x^0$.

Remembering that $F[A] = \partial A + AA$, we see from Eq. (6.1), that, if we symbolically denote the variables entering the first-order action as $q = A_3$ and $p = \bar{b}_{ij}$, we have (when all variables depend only on time and after discarding the total-derivative term $\eta^{ij} F_{ij}^{\text{lin}}[A] \sim \dot{A}$)

$$ \frac{\mathcal{L}_1[p, q]}{c_F} = p \dot{q} - pqq - \kappa^2 p^2. \qquad (6.2) $$

In other words, we have a Hamiltonian action $p\dot{q} - H(p, q)$ with a cubic Hamiltonian $H(p, q) = pqq + \kappa^2 p^2$.

This would be a trivial problem if this action did not involve many different constraints. Our task will be to
delineate the precise constrained structure of the first-order action (6.1).

In the following, we find it convenient to work with the metriclike Fierz-Pauli variable $u_{ij}$, instead of $\bar{b}_{ij}$. To avoid confusion with the various notations used above, we shall actually henceforth denote this Fierz-Pauli variable as $h_{ij}$. In other words, we henceforth parametrize $b_{ij}$ as

$$ 
\begin{aligned}
b_{ij} &\equiv \eta_{ij} + \bar{b}_{ij} \\
&\equiv \eta_{ij} + h_{ij} - \frac{1}{2}\eta_{ij}h_k^k.
\end{aligned} 
\quad (6.3) $$

Let us first consider the kinetic $p\dot{q}$ component of the first-order Lagrangian $\mathcal{L}_{c_F}^{\mathcal{L}}$, namely

$$ p\dot{q} = b_{00}\tilde{F}_{00} - b_{0a}\tilde{F}_{0a} + b_{ab}\tilde{F}_{ab}, \quad (6.4) $$

where $\tilde{F}$ denotes the $\dot{q} \sim \dot{A}$ part of the Ricci tensor $F[A] = \partial A + AA$. We find [using our $(-, +, +, +)$ signature]

$$ \tilde{F}_{00} = \partial_0 A_{0a}{}^a, \quad (6.5) $$

$$ \tilde{F}_{0a} = \frac{1}{2}\partial_0 A_{ac}{}^c, \quad (6.6) $$

$$ \tilde{F}_{ab} = \frac{1}{2}(\partial_0 A_{a0b} + \partial_0 A_{b0a}). \quad (6.7) $$

Here and henceforth, the indices $a, b, c$ denote spatial indices, $a, b, c = 1, 2, 3$ (moved by the three-dimensional Euclidean metric $\delta_{ab}$), in contrast to the indices $i, j, k, \dots$ used above to denote spacetime frame indices $i, j, k, \dots = 0, 1, 2, 3$ (moved by the four-dimensional Minkowski metric $\eta_{ij}$). (Remember that, in our decoupling limit, the frame indices can be identified with Lorentz-coordinate indices $\mu = i$.)

Substituting the above equations in Eq. (6.4), we obtain

$$ 
\begin{aligned}
p\dot{q} &= \frac{1}{2}h_{00}(\tilde{F}_{00} + \tilde{F}_a{}^a) + \frac{1}{2}\bar{h}(\tilde{F}_{00} - \tilde{F}_a{}^a) - h_{0a}\tilde{F}_0{}^a + \tilde{F}_{ab}h^{ab} \\
&= \bar{h}\partial_0 A_{0a}{}^a - \frac{1}{2}h_{0a}\partial_0 A_{ab}{}^b - \frac{1}{2}h_{ab}(\partial_0 A_{0ab} + \partial_0 A_{0ba}) \\
&= -\frac{1}{2}h_{0a}\partial_0 A_{ab}{}^b - h_{ab}(\partial_0 A_{0(ab)} - \delta^{ab}\partial_0 A_{0c}{}^c),
\end{aligned} 
\quad (6.8) $$

where the bar in $\bar{h}$ denotes a trace only on the space indices: $\bar{h} \equiv \delta_{ab}h^{ab}$. Here we have also used

$$ \partial_0 A_{0(ab)} = \frac{1}{2}(\partial_0 A_{0ab} + \partial_0 A_{0ba}). \quad (6.9) $$

The above equation, Eq. (6.8), shows that the “momentum” $p_{00} = h_{00}$ has disappeared in the $p\dot{q}$ kinetic terms. The momentum $p_{00} = h_{00}$ can then only appear in the $pqq \sim hAA$ terms or the $p^2$ ones. Evidently, it can only appear

linearly in the $pqq \sim hAA$ terms. Concerning the $p^2$ ones, it is well known that the (quadratic) Fierz-Pauli mass term is such that $p_{00} = h_{00}$ only appears linearly. Indeed, from Eq. (6.1),

$$ 
\begin{aligned}
\kappa^2 p^2 &= \frac{\kappa^2}{4}(\bar{b}_{ij}\bar{b}^{ij} - \bar{b}^2) \\
&= \frac{\kappa^2}{4}(h_{ij}h^{ij} - (\eta^{ij}h_{ij})^2) \\
&= \frac{\kappa^2}{4}(h_{ab}h^{ab} - \bar{h}^2 - 2h_{0a}^2 + 2h_{00}\bar{h}).
\end{aligned} 
\quad (6.10) $$

Hence, we conclude that in torsion bigavity $h_{00}$ is an exact Lagrange multiplier whose equation of motion yields a scalar constraint equation. This situation is quite different from the usual massive gravity theories where the lapse appears linearly in the $p\dot{q}$ terms but nonlinearly in the mass term and where the existence of a scalar constraint in ghost-free theories is rather hidden [31,35–37].

Let us now continue our exploration of the constrained nature of the first-order action, Eq. (6.1). Equation (6.8) shows that among a priori possible 24 kinetic terms $\dot{q} = \partial_0 A_{ijk}$, i.e., $\partial_0 A_{0(ab)}$, $\partial_0 A_{0|ab}$, $\partial_0 A_{a00}$, $\partial_0 A_{[ab]0}$, $\partial_0 A_{[ab]c}$, only nine combinations, namely $\partial_0 A_{0(ab)}$ and $\partial_0 A_{ab}{}^b$, enter the first-order Lagrangian $\mathcal{L}_1$. Let us define as follows the corresponding nine $q$ variables entering $p\dot{q}$ terms:

$$ 
\begin{aligned}
q_a &\equiv -A_{ab}{}^b, \\
q_{(ab)} &\equiv -A_{0(ab)} + \delta_{ab}A_{0c}{}^c.
\end{aligned} 
\quad (6.11) $$

In addition to the nine “active” $q$'s, Eq. (6.11), entering the kinetic terms, the first-order action contains 15 remaining $q$-type variables which enter the Hamiltonian $H(p, q)$ without participating in the kinetic terms. These remaining 15$q$-type variables happen to include two different types of variables. Six of them turn out to appear linearly in the action (and therefore to be Lagrange multipliers), while nine of them appear quadratically in the action.

The six Lagrange-multiplier $q$ variables (decorated with a superscript $L$, for Lagrange) are

$$ Q_a^L = A_{a00}, \quad (6.12) $$

$$ Q_{[ab]}^L = A_{[ab]0} = \frac{1}{2}(A_{ab0} - A_{ba0}). \quad (6.13) $$

Because of the antisymmetry over $ab$ in $Q_{[ab]}^L$ there are only three different two-index Lagrange multipliers.

The remaining nine $q$-type variables entering the action quadratically are

$$ \check{Q}_{[ab]} = A_{0[ab]} = \frac{1}{2}(A_{0ab} - A_{0ba}), \quad (6.14) $$

$$ \check{Q}_{abc} = A_{[ab]c} - \frac{1}{2}(A_{ad}{}^d \delta_{bc} - A_{bd}{}^d \delta_{ac}), \quad (6.15) $$
where $\check{Q}_{abc}$ is antisymmetric in $(a \leftrightarrow b)$ and traceless. The ab antisymmetry of $\check{Q}_{[ab]}$ means that they contain only three independent components. The algebraic constraints on $\check{Q}_{abc}$ mean that these contain only six independent components. To derive the equations of motion of $\check{Q}_{abc}$, we will express it in terms of its symmetric dual tensor $\check{Q}^*_{(ab)}$ as

$$ \check{Q}_{abc} = \epsilon_{abd} \check{Q}^*_{(cd)}. \quad (6.16) $$

At this stage, the first-order Lagrangian can be written as

$$ \frac{\mathcal{L}_1}{c_F} = h_{ab} \dot{q}^{ab} + \frac{1}{2} h_{0a} \dot{q}_a - \mathcal{H}(h_{ab}, h_{00}, h_{0a}, q_{ab}, q_a, Q^L, \check{Q}), \quad (6.17) $$

where we recall that $\mathcal{H} = pqq + p^2 = hAA + hh$.

As announced, this first-order action involves seven Lagrange multipliers, namely $h_{00}$, $Q_a^L$, and $Q_{[ab]}^L$. This means that the Hamiltonian part of the action can be written as

$$ -\mathcal{H} = h_{00}C[\text{vars}] + Q_a^L C_a[\text{vars}] + \frac{1}{2} Q_{[ab]}^L C_{[ab]}[\text{vars}] - \mathcal{H}^{\text{remain}}[\text{vars}], \quad (6.18) $$

where $[\text{vars}]$ denote the following (initial) list of variables:

$$ [\text{vars}] = [h_{ab}, q_{ab}, h_{0a}, q_a, \check{Q}_{[ab]}, \check{Q}_{abc}]. \quad (6.19) $$

Note that the seven Lagrange multipliers imply seven corresponding (first-class) constraints, namely

$$ \begin{aligned} C[\text{vars}] &= 0, \\ C_a[\text{vars}] &= 0, \\ C_{[ab]}[\text{vars}] &= 0. \end{aligned} \quad (6.20) $$

However, the impact of these seven constraints depend on their precise form. In addition, the fact that the six variables $\check{Q}_{abc}$ enter only algebraically (and quadratically) in the action means that they will also give rise to six (second-class) constraints, namely

$$ \frac{\partial \mathcal{H}}{\partial \check{Q}^*_{(ab)}} = 0. \quad (6.21) $$

The latter additional constraints are best taken into account after having dealt with the seven Lagrange-multiplier constraints, Eq. (6.20).

To clarify the difference between the linear and the nonlinear dynamics of torsion bigravity, we will consider first the leading-order (linear) dynamics.

## VII. TORSION BIGRAVITY AT THE LINEAR LEVEL

At the leading-order (linear equations of motion), the Hamiltonian has the general form (6.18) with

$$ C^{LO} = -\frac{\kappa^2}{2} \bar{h}, \quad (7.1) $$

$$ C_a^{LO} = 2q_a, \quad (7.2) $$

$$ C_{[ab]}^{LO} = 4\check{Q}_{[ab]}, \quad (7.3) $$

and

$$ -\mathcal{H}^{\text{remain}} = -\frac{1}{2}q_a^2 + \frac{1}{2}(q_a^a)^2 - q_{ab}q^{ab} + \check{Q}_{ab}\check{Q}^{ab} + \check{Q}_{abc}\check{Q}^{cba} - \frac{\kappa^2}{4}(h_{ab}h^{ab} - \bar{h}^2 - 2h_{0a}^2). \quad (7.4) $$

Let us recall that the constraint imposed by any Lagrange-multiplier term, say $q^L C(q')$, has a double effect: (i) it eliminates some variable (possibly to be chosen among several) entering the corresponding constraint equation, $C(q') = 0$, and (ii) solving the latter constraint furthermore eliminates $q^L$ from the action.

At the linear level, the term $Q_a^L C_a^{LO}$ imposes that $q_a = 0$ and also eliminates $Q_a^L$. Then the term $Q_{[ab]}^L C_{[ab]}^{LO}$ imposes that $\check{Q}_{[ab]} = 0$ and also eliminates $Q_{[ab]}^L$. Finally, the scalar term $h_{00}C^{LO}$ sets the spatial trace $\bar{h}$ of $h_{ab}$ to zero and also eliminates $h_{00}$. Those eliminations imply also a simplification of the kinetic terms of the action. In view of the constraint on $\bar{h} = \delta^{ab}h_{ab}$, it will be henceforth convenient to generally decompose $h_{ab}$ in trace-free part $(h_{\langle ab\rangle} \equiv h_{ab} - \frac{1}{3}h_c^c\delta_{ab})$, and trace part, namely

$$ h_{ab} = h_{\langle ab\rangle} + \frac{1}{3}\bar{h}\delta_{ab}, \quad (7.5) $$

and similarly, with $\bar{q} \equiv \delta^{ab}q_{ab}$,

$$ q_{ab} = q_{\langle ab\rangle} + \frac{1}{3}\bar{q}\delta_{ab}. \quad (7.6) $$

With this notation the kinetic terms of the action become

$$ h_{\langle ab\rangle} \dot{q}^{\langle ab\rangle} + \frac{1}{3} \bar{h} \dot{\bar{q}} + \frac{1}{2} h_{0a} \dot{q}_a. \quad (7.7) $$

The linear constraints $\bar{h} = 0$ and $q_a = 0$ then eliminate two of the kinetic terms in the action and lead to a reduced (but not yet fully reduced) action of the form

$$ \frac{\mathcal{L}_1^{\text{red}}}{c_F} = h_{\langle ab\rangle} \dot{q}^{\langle ab\rangle} - \mathcal{H}^{\text{red}}[h_{\langle ab\rangle}, q_{\langle ab\rangle}, \bar{q}, h_{0a}, \check{Q}_{abc}], \quad (7.8) $$
where

$$ -\mathcal{H}^{\text{red}}[h_{\langle ab\rangle}, q_{\langle ab\rangle}, \bar{q}, h_{0a}, \check{Q}_{abc}] = \frac{1}{6}\bar{q}^2 - q_{\langle ab\rangle}^2 + \check{Q}_{abc}\check{Q}_{cba} - \frac{\kappa^2}{4}(-2h_{0a}^2 + h_{\langle ab\rangle}^2). \quad (7.9) $$

At this stage, the variables $\bar{q}$, $h_{0a}$, and $\check{Q}_{abc}$ enter only algebraically in the action, so that their respective equations of motion,

$$ \begin{aligned} \frac{\partial \mathcal{H}^{\text{red}}}{\partial \bar{q}} &= 0, \\ \frac{\partial \mathcal{H}^{\text{red}}}{\partial h_{0a}} &= 0, \\ \frac{\partial \mathcal{H}^{\text{red}}}{\partial \check{Q}_{ab}^*} &= 0, \end{aligned} \quad (7.10) $$

yield the additional constraints

$$ \bar{q} = 0, \quad h_{0a} = 0, \quad \text{and} \quad \check{Q}_{ab}^* = 0. \quad (7.11) $$

Finally, the fully reduced linear action involves only $h_{\langle ab\rangle}$, and $q_{\langle ab\rangle}$ and reads

$$ \frac{\mathcal{L}_1^{\text{full red}}}{c_F} = h_{\langle ab\rangle} \dot{q}^{\langle ab\rangle} - \mathcal{H}^{\text{full red}}[h_{\langle ab\rangle}, q_{\langle ab\rangle}], \quad (7.12) $$

with

$$ -\mathcal{H}^{\text{full red}}[h_{\langle ab\rangle}, q_{\langle ab\rangle}] = -\frac{\kappa^2}{4} h_{\langle ab\rangle}^2 - q_{\langle ab\rangle}^2. \quad (7.13) $$

The above quadratic-in-fields Hamiltonian action explicitly displays the presence of the five d.o.f. (i.e., five $q$'s and five $p$'s) of a linear Fierz-Pauli (spatially homogeneous) massive spin-2 field, namely $h_{\langle ab\rangle}$ and its canonical conjugate $q_{\langle ab\rangle}$.

In the next section, we consider the nonlinear dynamics of torsion bigravity.

## VIII. TORSION BIGRAVITY AT THE NONLINEAR LEVEL

At the nonlinear level, we still have seven first-class constraints associated with the seven Lagrange-multiplier terms $h_{00}C + Q_a^L C_a + \frac{1}{2} Q_{[ab]}^L C_{[ab]}$, but the corresponding constraints are less potent in eliminating dynamical variables than in the linear case, notably because they do not eliminate kinetic terms. One way to see what are the remaining number of d.o.f. is to proceed as follows.

Having decomposed as above $h_{ab}$ and $q_{ab}$ in trace-free and trace parts, and using for later convenience an integration by parts, we can rewrite the action so that $q_a$ appears undifferentiated,

$$ \frac{\mathcal{L}_1}{c_F} = h_{\langle ab\rangle} \dot{q}^{\langle ab\rangle} + \frac{1}{3} \bar{h} \dot{q} - \frac{1}{2} q_a \dot{h}_{0a} - \mathcal{H}. \quad (8.1) $$

Introducing the following new list of variables

$$ [\text{vars}'] = [\bar{h}, q_a, \check{Q}_{[ab]}; \bar{q}, h_{0a}, h_{\langle ab\rangle}, q_{\langle ab\rangle}, \check{Q}_{abc}], \quad (8.2) $$

the Hamiltonian part of the action has the structure

$$ -\mathcal{H} = h_{00}C[\text{vars}'] + Q_a^L C_a[\text{vars}'] + \frac{1}{2} Q_{[ab]}^L C_{[ab]}[\text{vars}'] - \mathcal{H}^{\text{remain}}[\text{vars}'] \quad (8.3) $$

We separated the list of variables $[\text{vars}']$, Eq. (8.2), in two blocks; the first block consists of the seven variables

$$ [\bar{h}, q_a, \check{Q}_{[ab]}]. $$

This is because we have seen above that, in the linear approximation, the seven Lagrange constraints are respectively proportional to $C^{LO} \propto \bar{h}$, $C_a^{LO} \propto q_a$, $C_{[ab]}^{LO} \propto \check{Q}_{[ab]}$, and thereby respectively eliminate these seven variables $[\bar{h}, q_a, \check{Q}_{[ab]}]$. At the nonlinear level, taking into account the $(1+p)qq + p^2 \sim (1+h)AA + h^2$ structure of the Hamiltonian, the constraints $C_a$ and $C_{[ab]}$ are deformed but remain linear in $q_a$ and $\check{Q}_{[ab]}$ (with coefficients depending on $h_{\langle ab\rangle}$ and $\bar{h}$). One can then solve the linear system $C_a = 0$ and $C_{[ab]} = 0$ to get $q_a$ and $\check{Q}_{[ab]}$ as functions of $\bar{h}$, $h_{\langle ab\rangle}$, $q_{\langle ab\rangle}$, $h_{0a}$, $\bar{q}$, $\check{Q}_{abc}$. Inserting the latter expressions in the scalar constraint

$$ 0 = C = -\kappa^2 \bar{h} + \text{nonlinear terms}, \quad (8.4) $$

then yields a nonlinear equation for $\bar{h}$ which can, however, be perturbatively solved to express $\bar{h}$ as a function of the rest of the list of variables of Eq. (8.2), namely

$$ \bar{q}, h_{0a}, h_{\langle ab\rangle}, q_{\langle ab\rangle}, \check{Q}_{abc}. \quad (8.5) $$

The function $\bar{h}[\text{vars}']$ can be obtained at any desired order as a power series in the remaining variables.

At this stage one has eliminated $q_a$, $\check{Q}_{[ab]}$, and $\bar{h}$, as well as their corresponding Lagrange multipliers, $Q_a^L$, $Q_{[ab]}^L$, and $h_{00}$. This leads to a reduced action having the structure
$$ \frac{\mathcal{L}_1^{\text{red}}}{c_F} = h_{\langle ab \rangle} \dot{q}^{\langle ab \rangle} + \frac{1}{3} \dot{q} \bar{h}^{\text{sol}}[\text{vars}''; \check{Q}_{abc}] \\ - \frac{1}{2} \dot{h}_{0a} q_a^{\text{sol}}[\text{vars}''; \check{Q}_{abc}] - \mathcal{H}^{\text{red}}[\text{vars}''; \check{Q}_{abc}]. \quad (8.6) $$

Here, the list of variables [vars''; $\check{Q}_{abc}$], entering both the solutions obtained by solving the Lagrange constraints, $\bar{h}^{\text{sol}}[\text{vars}''; \check{Q}_{abc}]$ and $q_a^{\text{sol}}[\text{vars}''; \check{Q}_{abc}]$, as well as the non-kinetic part of the action, $-\mathcal{H}^{\text{red}}[\text{vars}''; \check{Q}_{abc}]$ comprises, as indicated, two different sets of variables. First, the reduced list

$$ [\text{vars}''] = [\bar{q}, h_{0a}, h_{\langle ab \rangle}, q_{\langle ab \rangle}], \quad (8.7) $$

and second, the six variables $\check{Q}_{abc}$. The reason why we separated $\check{Q}_{abc}$ from the list [vars''], is that the latter variables enter as dynamical variables in the kinetic terms of the action, while $\check{Q}_{abc}$, only appears algebraically in the action.

At this last stage we can further eliminate the six variables $\check{Q}_{abc}$ by solving their (reduced) equation of motion, namely (using the restricted components of $\check{Q}_{abc}$ to compute the derivative)

$$ 0 = c_F^{-1} \frac{\partial \mathcal{L}_1^{\text{red}}}{\partial \check{Q}_{abc}}. \quad (8.8) $$

This last constraint can be perturbatively solved for $\check{Q}_{abc}$ as a power series in the remaining dynamical variables because we have seen above that, at the linear level, $\check{Q}_{abc}$ appears quadratically in the Hamiltonian as a simple term $\propto \check{Q}_{abc} \check{Q}_{bca}$. Therefore, at the nonlinear level, the constraint (8.8) will be of the form

$$ 0 = \check{Q}_{bca} + \text{nonlinear terms.} \quad (8.9) $$

Note, however, the following crucial new feature of the latter constraint at the nonlinear level: In the action (8.6), the six variables $\check{Q}_{abc}$ appeared not only in the reduced Hamiltonian $-\mathcal{H}^{\text{red}}[\text{vars}; \check{Q}_{abc}]$ but also in the kinetic terms through the couplings

$$ + \frac{1}{3} \dot{q} \bar{h}^{\text{sol}}[\text{vars}''; \check{Q}_{abc}] - \frac{1}{2} \dot{h}_{0a} q_a^{\text{sol}}[\text{vars}''; \check{Q}_{abc}]. \quad (8.10) $$

As a consequence, the solution $\check{Q}_{abc}^{\text{sol}}$ of the constraint Eq. (8.8) will depend not only on the list of variables [vars''], Eq. (8.7), but also on $\dot{q}$ and $\dot{h}_{0a}$. The elimination of $\check{Q}_{abc}$ then leads to a final reduced action having the structure

$$ \frac{\mathcal{L}_1^{\text{red}'} }{c_F} = h_{\langle ab \rangle} \dot{q}^{\langle ab \rangle} + \frac{1}{3} \dot{q} \bar{h}^{\text{sol}}[\text{vars}''; \dot{q}, \dot{h}_{0a}] \\ - \frac{1}{2} \dot{h}_{0a} q_a^{\text{sol}}[\text{vars}''; \dot{q}, \dot{h}_{0a}] - \mathcal{H}^{\text{red}}[\text{vars}''; \dot{q}, \dot{h}_{0a}]. \quad (8.11) $$

This structure makes it clear that, at the nonlinear level, the reduced torsion bigravity action involves the time derivatives $\dot{q}$ and $\dot{h}_{0a}$ no longer in a linear way ($\dot{a} la p \dot{q}$) but in a nonlinear way.

By explicitly performing the reduction process indicated above, we have perturbatively computed the final reduced action, Eq. (8.11), as a power series in all the dynamical variables, starting at the quadratic level (corresponding to the linear dynamics), say

$$ \frac{\mathcal{L}_1^{\text{red}'}}{c_F} = \epsilon^2 + \epsilon^3 + \epsilon^4 + \dots, \quad (8.12) $$

where $\epsilon$ is a counting parameter keeping track of the order in all the variables entering the final reduced action

$$ \epsilon \sim h_{\langle ab \rangle}, q_{\langle ab \rangle}, \dot{q}_{\langle ab \rangle}; h_{0a}, \dot{h}_{0a}; \bar{q}, \dot{q}. \quad (8.13) $$

In particular, we found that at order $\mathcal{O}(\epsilon^4)$, i.e., when considering terms quartic in the fields, there arise contributions quadratic in $\dot{h}_{0a}$ with coefficients involving $h_{\langle ab \rangle}$, e.g., terms of the form

$$ \sim h_{\langle ac \rangle} h_{\langle cb \rangle} \dot{q} \dot{h}_{0a} \dot{h}_{0b}. \quad (8.14) $$

The presence of these terms means that the three variables $\dot{h}_{0a}$ which could be eliminated at the leading-order, $\mathcal{O}(\epsilon^2)$, introduce three new d.o.f. (beyond the normal five d.o.f. described by the canonical variables $h_{\langle ab \rangle}, q_{\langle ab \rangle}$) at the quartic-in-field level.

At the quartic level, the variable $\bar{q}$ appears only algebraically (but nonlinearly) in the reduced action. This implies that, at this level, the field equation of $\bar{q}$ would algebraically determine it in terms of the other variables, so that it would not represent a new d.o.f.

However, starting at the quintic level, $\mathcal{O}(\epsilon^5)$, there are terms in the reduced action that couple $\dot{q}$ to the time derivatives of $h_{0a}$, say terms of the structure

$$ \sim h_{\langle ac \rangle} h_{\langle cb \rangle} \dot{q} \dot{h}_{0a} \dot{h}_{0b}. \quad (8.15) $$

When varying these terms with respect to $h_{0a}$ we see that the second time derivative of $\bar{q}$ will appear in the equations of motion. This shows that one now needs to give both $\bar{q}$ and $\dot{q}$ as initial data. In other words, $\bar{q}$ enters the dynamics as a full additional d.o.f. at the quintic level. In addition, at the sextic level, $\mathcal{O}(\epsilon^6)$, there are terms in the reduced action that involve the square of $\dot{q}$, say terms of the form

$$ \sim h_{\langle ac \rangle} h_{\langle cb \rangle} (\dot{q})^2 \dot{h}_{0a} \dot{h}_{0b}. \quad (8.16) $$
The presence of such terms in the nonlinear action confirms that $\bar{q}$ provides a full additional d.o.f. in the nonlinear dynamics of torsion bigravity.

## IX. CONCLUSIONS

We studied the nonlinear regime of torsion gravity, a three-parameter class of modified gravity theories involving a propagating (con)torsion $K_{j\mu}^i = A_{j\mu}^i - \omega_{j\mu}^i$ in addition to the metric $g_{\mu\nu}$ (with associated Levi-Civita connection $\omega_{j\mu}^i$). We focused on the formal weak-gravity limit $16\pi G \equiv \frac{1}{c_R} \to 0$ where $g_{\mu\nu} = \eta_{\mu\nu} + \frac{1}{\sqrt{c_R}}\hat{h}_{\mu\nu}$ decouples both from matter and from the dynamics of $A_{j\mu}^i$ and defines a Minkowski background $\eta_{\mu\nu} = \eta_{ij}$. (Here, both Greek indices $\mu, \nu, ...$ and Latin indices from the second part of the alphabet, $i, j, k, ...$ take four values: 0,1,2,3.)

We studied the Poincaré-invariant nonlinear dynamics of the contorsion field $K_{[ij]k} = A_{[ij]k}$ in the decoupling limit. We showed that it could be described by several different actions, notably: (i) a second-order action $\mathcal{L}_{TBG}[A_3]$, Eq. (5.1), involving only $A_{ijk}$; (ii) a first-order action $\mathcal{L}_1[b_2, A_3]$, Eq. (5.8), involving both $A_{ijk}$ and a (Fierz-Pauli-like) symmetric 2-tensor $b_{ij}$; and (iii) a second-order action $\mathcal{L}_2[b_2]$, Eq. (5.22), involving only the Fierz-Pauli-like tensor $b_{ij}$. The latter field is related to the Ricci tensor $\mathcal{F}_{ij} = F_{ij}[A_3]$ of the torsionfull connection $A_{ijk}$ via

$$b_{ij} = \eta_{ij} + \frac{2}{\kappa^2} \left( \mathcal{F}_{ij} - \frac{1}{3} \mathcal{F} \eta_{ij} \right). \qquad (9.1)$$

We also parametrize $b_{ij}$ by a metric-perturbation-like tensor $h_{ij}$ defined as

$$b_{ij} = \eta_{ij} + h_{ij} - \frac{1}{2} \eta_{ij} h_k^k. \qquad (9.2)$$

The second-order action $\mathcal{L}_2[b_2]$ defines a massive gravity theory which differs from the usually considered ones in having a quadratic Fierz-Pauli mass term, but a non-Einsteinian kinetic term. We found convenient to study the nonlinear dynamics of the massive gravity theory $\mathcal{L}_2[b_2]$ by means of the first-order action $\mathcal{L}_1[b_2, A_3]$, Eq. (5.8).

We focused on “spatially homogeneous” solutions of the latter torsion-inspired massive gravity theory. Our final result is that the dynamics defined by $\mathcal{L}_2[b_2]$ involves (in the spatially homogeneous case) nine d.o.f. (instead of the five d.o.f. of a normal, Fierz-Pauli massive spin-2 field). Our analysis used the first-order action $\mathcal{L}_1[b_2, A_3]$ to display and solve the constraints implied by this action; see Sec. VIII.

We define (with spatial indices $a, b, c, ... = 1, 2, 3$) the torsion-related variable

$$q_{(ab)} = -A_{0(ab)} + \delta_{ab} A_{0c}{}^c, \qquad (9.3)$$

and its spatial trace $\bar{q} = \delta_{ab}q_{(ab)}$, together with its trace-free part $q_{(ab)} \equiv q_{(ab)} - \frac{1}{3}\bar{q}\delta_{ab}$. Similarly, $h_{(ab)} \equiv h_{ab} - \frac{1}{3}\bar{h}\delta_{ab}$, with $\bar{h} = \delta_{ab}h_{ab}$.

The final reduced form of the action, after Legendre transforming the kinetic term $h_{(ab)}\dot{q}_{(ab)}$ (which is the only occurrence of $\dot{q}_{(ab)}$ in the reduced action) can be described by a second-order Lagrangian with the structure

$$\mathcal{L}(h_{(ab)}, \dot{h}_{(ab)}; h_{0a}, \dot{h}_{0a}; \bar{q}, \dot{\bar{q}}), \qquad (9.4)$$

where all the time-derivative terms appear nonlinearly, though at different orders in nonlinearity: $\dot{h}_{(ab)}$ enters at the leading (Fierz-Pauli), quadratic order

$$\frac{\mathcal{L}^{LO}}{c_F} = \frac{1}{4} \left( (\dot{h}_{(ab)})^2 - \kappa^2 (h_{(ab)})^2 \right), \qquad (9.5)$$

while $\dot{h}_{0a}$ enters bilinearly at the quartic order, Eq. (8.14), and $\dot{\bar{q}}$ enters bilinearly at the sextic order, Eq. (8.16).

The general solution of the corresponding second-order equations of motion then requires that one gives initial data for the 18 variables

$$h_{(ab)}, \dot{h}_{(ab)}; h_{0a}, \dot{h}_{0a}; \bar{q}, \dot{\bar{q}}, \qquad (9.6)$$

corresponding to the nine d.o.f. $h_{(ab)}, h_{0a}, \bar{q}$.

It is possible that the appearance of the four new d.o.f. (here parametrized by $h_{0a}$ and $\bar{q}$) would arise in a different way, and at lower perturbative orders, when considering generic, inhomogeneous solutions of torsion bigravity (notably via the presence of spatial gradients in the constraints). However, on the one hand, our simplifying homogeneity ansatz is consistent,<fnref n="1" /> and, on the other hand, it suffices to conclude that torsion gravity theories do not define a healthy theory containing, in all nonlinear situations, only the hoped-for 2 + 5 d.o.f. of massless spin-2 and massive spin-2 excitations.

Our negative conclusion concerning the healthiness of torsion bigravity (and thereby of the more general six-parameter torsion gravity theories) is not in conflict with the positive results obtained in previous works Refs. [16,17,19–24] because the latter results had studied less general dynamical configurations. Our proof above has shown that the four additional d.o.f. were rather hidden, because they only start to arise at the quartic order in torsion variables and in time-dependent situations. Let us finally mention that our negative conclusion extends the negative conclusions of the Hamiltonian analysis of Ref. [38] to a case

<footnote n="1">¹As a side check, we applied the approach of homogeneous solutions to the usual massive gravity theories (containing either a “bad” nonlinear Fierz-Pauli mass term or a “good” the de Rham-Gabadze-Tolley one) and we confirmed that this approach predicts in all cases the correct counting of d.o.f.</footnote>
that was not considered there, namely what would be called
in their nomenclature the simple spin-2⁺ case.

ACKNOWLEDGMENTS

T. J. thanks the hospitality and the stimulating environ-
ment of the Institut des Hautes Etudes Scientifiques. T. J. is
supported by a LabEx Junior Research Chair Fellowship.

The present research was also partly supported by the
“2021 Balzan Prize for Gravitation: Physical and
Astrophysical Aspects” awarded to Thibault Damour.

DATA AVAILABILITY

No data were created or analyzed in this study.

[1] S. Navas et al. (Particle Data Group), Review of particle physics, Phys. Rev. D **110**, 030001 (2024).

[2] C. M. Will, *Theory and Experiment in Gravitational Physics* (Cambridge University Press, Cambridge, England, 2018).

[3] C. Kiefer, Quantum gravity—an unfinished revolution, arXiv:2302.13047.

[4] A. Joyce, B. Jain, J. Khoury, and M. Trodden, Beyond the cosmological standard model, Phys. Rep. **568**, 1 (2015).

[5] N. Frusciante and L. Perenon, Effective field theory of dark energy: A review, Phys. Rep. **857**, 1 (2020).

[6] P. Salucci, N. Turini, and C. Di Paolo, Paradigms and scenarios for the dark matter phenomenon, Universe **6**, 118 (2020).

[7] E. Berti et al., Testing general relativity with present and future astrophysical observations, Classical Quantum Gravity **32**, 243001 (2015).

[8] E. Cartan, Sur les variétés à connexion affine et la théorie de la relativité généralisée. (première partie), Ann. Sci. Ecole Norm. Sup. **40**, 325 (1923).

[9] E. Cartan, Sur les variétés à connexion affine et la théorie de la relativité généralisée. (première partie) (Suite), Ann. Sci. Ecole Norm. Sup. **41**, 1 (1924).

[10] H. Weyl, A remark on the coupling of gravitation and electron, Phys. Rev. **77**, 699 (1950).

[11] F. W. Hehl, P. Von Der Heyde, G. D. Kerlick, and J. M. Nester, General relativity with spin and torsion: Foundations and prospects, Rev. Mod. Phys. **48**, 393 (1976).

[12] E. Sezgin and P. van Nieuwenhuizen, New ghost free gravity Lagrangians with propagating torsion, Phys. Rev. D **21**, 3269 (1980).

[13] K. Hayashi and T. Shirafuji, Gravity from Poincare gauge theory of the fundamental particles. 1. Linear and quadratic Lagrangians, Prog. Theor. Phys. **64**, 866 (1980); **65**, 2079(E) (1981).

[14] K. Hayashi and T. Shirafuji, Gravity from poincare gauge theory of the fundamental particles. 4. Mass and energy of particle spectrum, Prog. Theor. Phys. **64**, 2222 (1980).

[15] E. Sezgin, Class of ghost free gravity Lagrangians with massive or massless propagating torsion, Phys. Rev. D **24**, 1677 (1981).

[16] V. P. Nair, S. Randjbar-Daemi, and V. Rubakov, Massive spin-2 fields of geometric origin in curved spacetimes, Phys. Rev. D **80**, 104031 (2009).

[17] V. Nikiforova, S. Randjbar-Daemi, and V. Rubakov, Infrared modified gravity with dynamical torsion, Phys. Rev. D **80**, 124050 (2009).

[18] P Baekler, F. W. Hehl, J. M. Nester, Poincare gauge theory of gravity: Friedman cosmology with even and odd parity modes. Analytic part, Phys. Rev. D **83**, 024001 (2011).

[19] C. Deffayet and S. Randjbar-Daemi, Non linear Fierz-Pauli theory from torsion and bigravity, Phys. Rev. D **84**, 044053 (2011).

[20] V. Nikiforova and T. Damour, Infrared modified gravity with propagating torsion: Instability of torsionfull de Sitter-like solutions, Phys. Rev. D **97**, 124014 (2018).

[21] T. Damour and V. Nikiforova, Spherically symmetric solutions in torsion bigravity, Phys. Rev. D **100**, 024065 (2019).

[22] V. Nikiforova, Absence of a Vainshtein radius in torsion bigravity, Phys. Rev. D **101**, 064017 (2020).

[23] V. Nikiforova, Spherically symmetric perturbations of a Schwarzschild black hole in torsion bigravity, Phys. Rev. D **104**, 024032 (2021).

[24] V. Nikiforova, Perturbations of a Schwarzschild black hole in torsion bigravity, Phys. Rev. D **105**, 064060 (2022).

[25] M. Fierz and W. Pauli, On relativistic wave equations for particles of arbitrary spin in an electromagnetic field, Proc. R. Soc. A **173**, 211 (1939).

[26] H. van Dam and M. J. G. Veltman, Massive and massless Yang-Mills and gravitational fields, Nucl. Phys. **B22**, 397 (1970).

[27] A. I. Vainshtein, To the problem of nonvanishing gravitation mass, Phys. Lett. **39B**, 393 (1972).

[28] D. G. Boulware and S. Deser, Can gravitation have a finite range?, Phys. Rev. D **6**, 3368 (1972).

[29] T. Damour and I. I. Kogan, Effective Lagrangians and universality classes of nonlinear bigravity, Phys. Rev. D **66**, 104024 (2002).

[30] E. Babichev, C. Deffayet, and R. Ziour, The Vainshtein mechanism in the decoupling limit of massive gravity, J. High Energy Phys. 05 (2009) 098.

[31] C. de Rham, G. Gabadadze, and A. J. Tolley, Resummation of massive gravity, Phys. Rev. Lett. **106**, 231101 (2011).

[32] S. F. Hassan and R. A. Rosen, Bimetric gravity from ghost-free massive gravity, J. High Energy Phys. 02 (2012) 126.

[33] B. Bellazzini, G. Isabella, S. Ricossa, and F. Riva, Massive gravity is not positive, Phys. Rev. D **109**, 024051 (2024).

[34] C. de Rham, Massive gravity, Living Rev. Relativity **17**, 7 (2014).
[35] C. de Rham, G. Gabadadze, and A. J. Tolley, Ghost free massive gravity in the Stückelberg language, Phys. Lett. B **711**, 190 (2012).

[36] S. F. Hassan and R. A. Rosen, Confirmation of the secondary constraint and absence of ghost in massive gravity and bimetric gravity, J. High Energy Phys. 04 (2012) 123.

[37] S. F. Hassan and R. A. Rosen, Resolving the ghost problem in non-linear massive gravity, Phys. Rev. Lett. **108**, 041101 (2012).

[38] H.-J. Yo and J. M. Nester, Hamiltonian analysis of Poincare gauge theory: Higher spin modes, Int. J. Mod. Phys. D **11**, 747 (2002).