Starobinsky-Bel-Robinson Gravity

• Starobinsky-Bel-Robinson gravity is a modified four-dimensional theory that combines the Starobinsky R+R² inflationary model with a quartic-curvature correction built from the Bel-Robinson tensor.
• It employs a two-parameter framework defined by the scalaron mass (m) and a dimensionless coupling (β) to control inflation dynamics, black-hole solutions, and ghost-free bounds.
• Observable effects, including exact de Sitter solutions, black-hole shadow deformations, and modified Hawking radiation, highlight its potential for probing high-curvature regimes.

Starobinsky-Bel-Robinson gravity is a four-dimensional effective modified-gravity theory that combines the Starobinsky $R+R^2$ sector with a quartic-curvature correction built from the square of the Bel-Robinson tensor. In the literature under this name, the model is motivated jointly by the phenomenology of Starobinsky inflation and by higher-curvature terms induced by compactification of eleven-dimensional M-theory, and it is treated as a two-parameter framework involving the scalaron mass scale $m$ and a dimensionless coupling $\beta$ controlling the Bel-Robinson sector (Ketov, 2022, Ketov et al., 2022). Its main applications have been inflationary cosmology, Schwarzschild-type and rotating black holes, Hawking radiation, quasinormal spectra, and optical observables in strong gravity.

1. Action, invariants, and field-content organization

The original proposal formulates the theory as Einstein gravity plus the Starobinsky correction plus a Bel-Robinson-tensor-squared term. One representative form is

$$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$

together with the equivalent expression

$$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$

Subsequent cosmological and black-hole analyses usually parameterize the quartic sector explicitly by a positive dimensionless coupling,

$$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$

or, equivalently,

$$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$

These notations encode the same structural idea: the Einstein-Hilbert term, the Starobinsky $R^2$ term, and a quartic invariant built from the Bel-Robinson sector (Ketov, 2022, Ketov et al., 2022).

The relevant curvature densities are

$${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$

and

$$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$

The Bel-Robinson tensor is written in this literature as

$m$0

with the basic identity

$m$1

where $m$2 is identified with the four-dimensional Euler or Gauss-Bonnet density. In the original presentation, the Bel-Robinson tensor is also emphasized as the gravitational analogue of the Maxwell energy-momentum tensor, with the on-shell properties $m$3, $m$4, and $m$5 under the Einstein equations (Ketov, 2022).

A standard technical simplification is that $m$6 vanishes in the backgrounds most often studied. In spatially flat FLRW and in static spherically symmetric geometries, the dynamics therefore reduces effectively to the $m$7 sector. Both the cosmology and black-hole papers exploit this reduction explicitly, and an auxiliary-field rewriting introduces a scalaron $m$8, a dilaton-like field $m$9, and an axion-like field $\beta$0 through $\beta$1, $\beta$2, and $\beta$3 (Ketov et al., 2022).

2. Cosmological dynamics and inflationary structure

On a spatially flat FLRW background,

$\beta$4

the basic invariants become

$\beta$5

while $\beta$6 drops out. The resulting Hubble equation obtained in the SBR model is

$\beta$7

In the Starobinsky limit $\beta$8, this reduces to the standard quasi-de Sitter dynamics,

$\beta$9

with the slow-roll solution

$$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$0

The same analysis yields theoretical and observational bounds on the Bel-Robinson coupling: a ghost-free bound $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$1, a stronger unitarity/causality bound $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$2, and a preferred observational bound $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$3 from the spectral index $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$4 (Ketov et al., 2022).

A distinctive feature of SBR gravity is the existence of an exact isotropic de Sitter solution absent in pure Starobinsky gravity. With the ansatz

$$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$5

the reduced FLRW equations collapse to

$$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$6

This solution exists for $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$7, and the papers stress that its value is set entirely by the Bel-Robinson sector, whereas the $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$8 term does not determine $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{8M_{\rm Pl}^6}T^2 \right],$$9 itself. The same analyses also show that for the standard Starobinsky sign $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$0, the de Sitter fixed point is unstable, because the linearized characteristic equation has at least one positive root. Stability can occur only for sufficiently negative $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$1, specifically $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$2, which no longer reproduces the standard Starobinsky limit as $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$3 (Do, 8 Jul 2025).

The anisotropic sector sharpens this picture. In Bianchi type I, the original SBR model with $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$4 admits no exact anisotropic exponential solution, whereas a modified model with the sign of the $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$5 term flipped can admit such a solution. Even there, the dynamical-system analysis gives an unstable anisotropic fixed point, with characteristic roots including $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$6. The same work concludes that the only stable exact de Sitter solution occurs in the modified model with negative $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$7 coefficient (Do et al., 2023).

A related comparison with the Grisaru-Zanon quartic-curvature correction from closed superstrings shows that, in a flat Friedmann background, the leading slow-roll pieces of the GZ invariant and the Bel-Robinson square coincide: $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$8 The coincidence of the $$S_{\rm SBR}[g_{ij}] = \frac{M_{\rm Pl}^2}{2} \int d^4x\,\sqrt{-g}\left[ R + \frac{1}{6m^2}R^2 + \frac{1}{32M_{\rm Pl}^6}(P_4^2-E_4^2) \right].$$9 and $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$0 terms indicates that the Bel-Robinson sector captures the dominant slow-roll structure of a wider class of quartic superstring corrections, while differing at subleading order (Toyama et al., 2024).

3. Static and rotating black-hole geometries

The best-studied black-hole background in SBR gravity is a Schwarzschild-type solution obtained perturbatively in $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$1. In static spherical symmetry,

$$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$2

with

$$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$3

The correction is very high order in $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$4, so the geometry remains asymptotically flat and reduces to Schwarzschild as $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$5. Solving $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$6 perturbatively gives

$$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$7

so the horizon shifts inward at fixed mass (Delgado et al., 2022).

The same background has also been analyzed from the perturbation-theory side, where it is noted that the numerical prefactor in the metric correction is very large, so even $$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$8 already gives a noticeable deformation of Schwarzschild. This matters for the effective potential governing wave propagation and for eikonal observables near the photon sphere (Bolokhov, 2023).

A rotating extension has been constructed through a Newman-Janis-type procedure. In Boyer-Lindquist form the metric is governed by

$$S_{\rm SBR}=\frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\frac{1}{6m^2} R^2 -\frac{\beta}{8m^6}T^{\mu\nu\lambda\rho}T_{\mu\nu\lambda\rho} \right]$$9

and

$$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$0

This family reduces to the static SBR black hole when $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$1 and to Kerr when $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$2. In the corresponding shadow analysis, the nonrotating solutions retain circular shadows, whereas rotation generates deformed one-dimensional curves with $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$3- and $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$4-dependent distortion (Belhaj et al., 2023).

4. Thermodynamics and Hawking radiation

The thermodynamics of the Schwarzschild-type SBR black hole has been worked out perturbatively from the metric solution. The Hawking temperature obtained from the surface gravity is

$$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$5

so the correction is positive and scales as $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$6. The entropy from Wald’s Noether-charge method is

$$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$7

which contains both a pure Bel-Robinson or quartic-curvature term and an $$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$8-related contribution. The same analysis finds a nonzero effective pressure,

$$S_{\rm SBR}= \frac{M_p^2}{2}\int d^4 x \sqrt{-g}\left[R+\alpha_1 R^2 +\alpha_2 \left({\cal G}^2-{P_4}^2\right)\right],\qquad \alpha_1=\frac{1}{6m^2},\quad \alpha_2=\frac{\beta}{32m^6}.$$9

and a $R^2$0-dependent correction to the evaporation lifetime, with the classical $R^2$1 scaling preserved at leading order (Delgado et al., 2022).

Hawking radiation has also been computed by the Parikh-Wilczek tunneling method in Painlevé coordinates. The emission rate takes the form

$R^2$2

which, at linear order in $R^2$3, implies

$R^2$4

The crucial point is not only the positivity of the correction but its structure: the tunneling calculation produces the same $R^2$5 dependence found in beyond-semiclassical tunneling formulas. The paper explicitly compares the two viewpoints and argues that modifying the classical background geometry first, then performing a semiclassical tunneling computation, can reproduce the same mass dependence usually attributed to higher-order quantum corrections with fixed background (Annamalai et al., 2024).

A separate thermodynamical analysis, including logarithmic entropy corrections and plasma-modified weak lensing, reports that the corrected entropy is only mildly sensitive to $R^2$6, that the Hawking temperature remains positive for the studied values of $R^2$7, and that the energy emission rate decreases as $R^2$8 increases. This identifies the emission spectrum as more sensitive to the SBR coupling than the entropy curve itself (Mustafa et al., 2024).

5. Perturbations, quasinormal spectra, and optical phenomenology

For a neutral, massless test scalar in the Schwarzschild-type SBR background, the radial perturbation equation has the Schrödinger form

$R^2$9

As $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$0 grows, the effective potential develops a negative gap outside the horizon. In the scalar sector, the time-domain analysis finds stability only for

$${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$1

when $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$2; beyond that value a non-oscillatory exponentially growing mode appears. In the stable region, the ringdown can exhibit two stages governed by different modes, and the sixth-order WKB method with Padé improvement fails to reproduce part of the spectrum, including the formal fundamental mode for small $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$3. As a consequence, the usual eikonal relation

$${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$4

reproduces only one branch of the full eikonal spectrum, so the standard quasinormal-mode/null-geodesic correspondence breaks down in this background (Bolokhov, 2023).

The optical sector displays a similarly non-Schwarzschild structure. In the static case the shadow remains circular and shrinks as $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$5 increases; in the rotating case both $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$6 and $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$7 increase the deformation, and for some parameter choices the shadow develops paired cusps. The weak-field deflection angle for the static SBR black hole is

$${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$8

so the SBR correction reduces the bending angle at order $${\cal G}=R^2-4R_{\mu\nu}R^{\mu\nu}+R_{\mu\nu\rho\sigma}R^{\mu\nu\rho\sigma},$$9. The same study gives a shadow-existence estimate $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$0 and an EHT-based estimate $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$1 from M87$$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$2 data (Belhaj et al., 2023).

A separate shadow analysis of the Schwarzschild-like SBR black hole emphasizes near-horizon tidal effects. For radially infalling particles, the radial and angular tidal forces can switch their initial behavior and become compressive and stretching, respectively, before the event horizon is reached, and the geodesic deviation can display an oscillating trend for the chosen initial conditions. The same work reports that the angular diameter of the SBR shadow is smaller than the Schwarzschild one and quotes EHT-based parameter windows

$$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$3

while also noting that, within those narrow windows, the SBR and Schwarzschild geometries are effectively indistinguishable in the plotted observables (Arora et al., 2023).

In plasma environments, weak lensing remains sensitive to the SBR coupling. One study finds that the deflection angle decreases with $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$4, increases with plasma concentration, and satisfies

$$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$5

so uniform plasma produces the largest bending and the largest image magnification among the three plasma models considered (Mustafa et al., 2024).

6. Effective-theory status and adjacent Bel-Robinson interpretations

SBR gravity is consistently presented as an effective high-curvature theory rather than a UV-complete description. At low curvature it reduces to Einstein gravity, while at high curvature the $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$6 and Bel-Robinson sectors become important for early-universe dynamics, strong gravitational fields, black-hole thermodynamics, and Hawking radiation. The proposal is explicitly described as an approximation not expected to remain valid near the Planck scale, and its parameters are treated as renormalized, scale-dependent quantities (Ketov, 2022).

The Bel-Robinson tensor itself has a broader conceptual role than its use in the SBR action. Within general relativity it has been employed as a super-energy object from which one can construct a gravitational entropy; that construction reproduces the Bekenstein-Hawking entropy when integrated over the Schwarzschild interior and ties cosmological entropy growth to the growth of Weyl curvature during structure formation (Clifton et al., 2013). In three-dimensional topologically massive gravity, a Bel-Robinson-type four-index tensor can be covariantly conserved on shell, but the same construction fails for generic quadratic-curvature models, where a derivative-built analogue $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$7 is instead suggested (Deser et al., 2010). These results indicate that the use of the Bel-Robinson tensor in SBR gravity is best understood as an effective invariant in the action, not as the direct importation of a universally conserved super-stress tensor.

A distinct thermodynamic line of work derives a vacuum equation in which Ricci curvature is sourced by Bel-Robinson super-energy if one abandons local Lorentz invariance and introduces a preferred local time direction. That framework is described as a Bel-Robinson analog of Starobinsky-type higher-curvature corrections, but the authors also stress that it is not a standard covariant metric theory derived from an action in the usual sense (Alonso-Serrano et al., 2024). This distinction is important: it separates SBR gravity proper from other “Bel-Robinson” modifications that share the same curvature-squared or super-energy vocabulary without sharing the same dynamical construction.

Taken together, the literature presents Starobinsky-Bel-Robinson gravity as a minimal, geometry-only, superstring-inspired effective theory in which the successful $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$8 inflationary sector is supplemented by a specific quartic curvature correction. Its phenomenology is correspondingly mixed: the Bel-Robinson sector generates exact de Sitter solutions and measurable black-hole deformations, but viable inflation requires very small $$P_4\equiv \frac{1}{2}\sqrt{-g}\,\epsilon_{\mu\nu\rho\sigma} R^{\rho \sigma}{}_{\alpha \beta} R^{\mu\nu\alpha\beta}.$$9, the standard exact de Sitter branch is unstable, part of the quasinormal spectrum evades standard WKB and geodesic intuition, and present black-hole observations constrain the allowed departure from Schwarzschild and Kerr to a narrow regime.