Extrasolar planets as a probe of modified gravity

We propose a new method to test modified gravity theories, taking advantage of the available data on extrasolar planets. We computed the deviations from the Kepler third law and use that to constrain gravity theories beyond General Relativity. We investigate gravity models which incorporate three screening mechanisms: the Chameleon, the Symmetron and the Vainshtein. We find that data from exoplanets orbits are very sensitive to the screening mechanisms putting strong constraints in the parameter space for the Chameleon models and the Symmetron, complementary and competitive to other methods, like interferometers and solar system. With the constraints on Vainshtein we are able to work beyond the hypothesis that the crossover scale is of the same order of magnitude than the Hubble radius $r_c \sim H_0^{-1}$, which makes the screening work automatically, testing how strong this hypothesis is and the viability of other scales.


Introduction
Gravity theories beyond General Relativity (GR) are a possible theoretical framework to explain several cosmological problems [1]. In particular, the intriguing present day's cosmic accelerated expansion [2,3]. Theoretical models which predict an extension to General Relativity must, however, comply with strong requirements: first of all the model must have similar predictions to those of the benchmark model ΛCDM at cosmological scales: observational data from both the background evolution and the linear large scale structure formation regime is fully consistent with ΛCDM [4]. Another condition is that the modifications to General Relativity must be suppressed, by physical mechanisms, in the regimes which are well tested, e.g. solar system scales. This requirement is assured via the so called screening mechanisms [5].
Modified gravity models with screening mechanism have been extensively studied in the literature: either focusing on the background cosmology [6,7], large and linear cosmological scales [8][9][10][11] or on astrophysical scales in the nonlinear regime [12], and finally at the small solar system scales using local gravity tests [13,14].
Any weakness in the screening mechanisms should result in appreciable deviations in what we predict from the General Relativity or, in weak field regime, from the Newtonian gravity. This kind of deviations have been used to test modified gravity inside the solar system with, for instance, spectral deviation data from the Cassini space mission, which ensures that the gravitational potential at the Sun surface must deviates less than 10 −5 from the value predicted by the Newtonian gravity Email addresses: mvsantos@protonmail.com (Marcelo Vargas dos Santos), d.f.mota@astro.uio.no (David F. Mota) [13]. Another interesting work investigates how much the gravity may deviates from the Newtonian gravity using the well measured solar system bodies orbits [14].
An important feature in the screening mechanisms is the dependence on the physical properties of the environment as the density field, for instance. Thus we expect that the information from different planetary systems should give more statistical significance once the screening works in different ways for each one. On the other side any significantly deviation should already be measured in Solar System, therefore this deviation must be small even in other planetary system.
In this work we investigate the possibility of using exoplanet data to test and constrain modified gravity models. This is not the first attempt to constrain modified gravity with exoplanets [15], but in there the authors just compared the theoretical prediction with only one measurement, the transiting exoplanet HD209458b "Osiris". This lacks statistical rigor, and in this work, we use more than one hundred data points and propose a statistical method to make a thorough treatment and study of these systems. The data are obtained from exoplanets.org portal [16].

The method
For any gravitational theory the planetary motion is described by the dynamics of a particle under the influence of a central force, i.e. the spatial dependence of the force is only on the distance of the planet to the force center, inside the host star. The relation between the revolution period, T , of a planet in a circular orbit, of radius r, and the absolute value of the gravitational force, F(r), is given by [17] where M p and M are the masses of the planet and the star respectively. For a modified gravity the total force is the sum between the Newtonian force and a fifth force, F(r) = F N (r) + F 5th (r). If the fifth force is null this relation reduces to the third Kepler law where G = 6.67384 × 10 −8 cm 3 g −1 s −2 is the Newton gravitational constant. So the deviation of the square period from the third Kepler law is Therefore, we can use the measured values of ε to constrain modified gravity using a χ 2 given by the sum between the weighted residuals of all the N measurements where ε th is the theoretical prediction for ε (the ratio between the fifth and Newtonian forces predicted by theory), ε obs is the observed value of ε (the deviation of the square period from the third Kepler law computed from the data), and σ i is the standard deviation (computed from the data by error propagation, the theoretical error can be neglected because it is proportional to ε 2 .). θ is a vector of model parameters and x is a vector of the physical properties of the star-planet system, which are: r -the planet orbit radius, R -the star radius, ρ -the star density, Φ S -the surface Newtonian potential. The field is also a function of the galaxy density, ρ g = 10 −24 g/cm 3 .
To compute the credible regions with 95% of confidence level (C.L.), we find the values of χ 2 which delimit the bounds (χ 2 b ), i.e. the value which gives equals to 0.95. Where ∆χ 2 = χ 2 − χ 2 min , and χ 2 min is the minimum value of χ 2 . Fig. 1 (top) shows a comparison between the residuals distribution with a normal distribution of the same mean and standard deviation, which suggests that this is a good approximation for the data distribution. Therefore assuming this distribution to solve 5 we find ∆χ 2 b 5.99 and ∆χ 2 b 8.08 for models with 2 and 3 free parameters, respectively [18].
The values of χ 2 min can me found minimizing the function 4, but this is not a single point for the tested models, there is a degeneracy between the parameters. To avoid this problem we assume χ 2 min equal to the value of χ 2 for which the fifth force is null, i.e. ε th ≡ 0, this assumption does not change the results. For a constant deviation, for example, we find ε = (0.0 +6.0 −6.0 ) × 10 −3 with this assumption and ε = (−0.1 +6.1 −5.9 ) × 10 −3 without it. A shift less than 1% of the confidence interval, and it is reduced in the cases with screening.

The data
Our observational data comes from the website exoplanets.org [16], which has a compilation of all observed exoplanetary systems: there are 2926 planets with well defined orbits. From those we pick for our analysis 177 planets with circular orbits and with measurements of all properties listed above. These systems are typically composed by a star, similar to the Sun, with a mass that varies from 0.5M to 1.5M , and planets like the jovian planets, with masses between 3.2M ⊕ and 600M ⊕ and typically close to the host star, with orbit radius smaller than 0.5AU. This corresponds to short periods, less than a few months.
All the properties are measured by gravity-independent methods, except the orbit radius which use the third Kepler law [19]. However this looks like a circularity it does not affect considerably our analysis, any appreciable deviation in orbit radius such would already have been measured in the Solar System. The Cassini mission [13] measured, by light time delay, the possible deviation from Newtonian gravitational constant in Solar System and the obtained value is very short However the photons may be coupled differently to the field, we expect a correction of the same order in the orbit radius, which is much smaller than σ i . Fig. 1 (bottom) shows that the deviations from the third Kepler law (ε obs ∼ 10 −5 ) are very smaller than the errors (σ ∼ 10 −1 ), suppressing any possible bias. In summary the possibility to measure deviations from Newtonian gravity is not in the measurements per se but with their errors. The relativistic corrections, which are less than 10 −8 , are not appreciable, in solar system, for example, the most affected body is Mercury, which the orbit radius is well determined by the Newtonian gravity, the effect is seen only in the precession of the perihelion. The orbits of all the exoplanets used in this work are circular ones, therefore there is no perihelion to be precessed.

Modified gravity models
The models described above were developed in order to explain the cosmic accelerated expansion, so for that reason most of that was designed in a such way that the screening may be automatically satisfied if the background evolution reproduces the same as the benchmark model ΛCDM. In other words the effective evolution of the Hubble parameter is where H 0 and Ω m0 are the Hubble constant and the present value of the matter density parameter, respectively. Ones of the most confident measures are that made by the Plack mission [4], which are H 0 = (67.6 ± 1.0)km s −1 Mpc −1 and Ω m0 = 0.310 ± 0.008. We use that values in our analysis to reduce the number of free parameters in the chameleon models.
In planetary scales we are interested in the correction to the two-body potential energy to find the correction to the third Kepler law. For a complete description of these spherical solutions see [20,21]. Here we only present a brief description of the theories and show the theoretical prediction for the ratio between the fifth force and the Newtonian force, which is needed for our method as mentioned above.

Yukawa potential
In gravitation, this potential appears when there is a massive propagator for gravitation or a scalar-field coupled to the matter fields, for example. The action of this kind of theory is where M pl = (8πG) −1/2 is the reduced Planck mass, φ is the scalar field, V(φ) is the potential, S m is the matter action, Figure 2: Constraints of the unscreened Yukawa potential in (λ, α = 2β 2 ) parameter space, the blue region is ruled out by 95% of confidence level (95% C.L.). The GR is recovered when β → 0 or λ → 0. g µν = A 2 (φ)g µν is the Jordan frame metric, A(φ) is the conformal factor and Ψ represents the matter fields. Minimizing this with respect to the field we find the Klein-Gordon equation where the effective potential is given by Perturbing the field around the value which minimizes V e f f (φ g ) the equation can be approximated as where (11) are the coupling and the mass of the field, which is related to the force range by λ = /mc, respectively. In general, even β and λ are functions that depend on the parameters of the theory, and can suppress the fifth force in specific conditions. When that happens we say that the model has a screening mechanism. We begin by testing the case which β and λ are constants, i.e. it has no physical mechanism to suppress the fifth force. In this case we have, for a static point star of constant mass (M ) in a background of constant density (ρ g ), that In the non-relativistic case the test particle (planets) motion moves under the influence of a effective gravitational potential given by Φ = Φ N + Φ 5th , where Φ 5th = log A(φ) βφ is the fifth force potential. From the eq. 12 solution of that is a Yukawa correction Taking the gradient we find that the ratio between the forces is This model has two free parameters: β, λ.

Chameleon
This screening mechanism occurs when the field is coupled to matter in such way that the effective range of the field, λ, depends on the local matter density [7,27]: when the density is low, the range is large and the fifth force acts at large scales; when density is high, the range is short, so the fifth force is screened.
As an example, we choose a field coupled to the matter via a linear conformal factor A(φ) = e φ/M 1 + φ/M. The screening happens when we consider the star as a spherically symmetric distribution of matter of constant density ρ embedded in medium, the galaxy, of constant density ρ g . The solution outside the sphere depends both of the field inside and outside the star whereφ =φ(ρ ) andφ g =φ(ρ g ) are the values of the field which minimizes the effective potential V e f f (ρ, φ) = V(φ) + ρφ/M, inside the star and in the galactic medium, and the mass of the field is given by m 2 φ = V ,φφ (φ g ). Therefore Choosing a Ratra-Peebles potential, V(φ) = Λ n+4 /φ n we have This model has the following free parameters: Λ, M, n.

f(R) gravity
This is the simplest modification that can be made to the Einstein-Hilbert Lagrangian, adding a general non linear function of R, L = R + f (R). One positive aspect of f (R) is the equivalence with chameleon theories, where the scalar field is the first derivative of the f (R), φ = f R = d f /dR. In these models the coupling constant is M = 2, and the potential is Therefore where f R = f R (R ) and f Rg = f R (R g ). R and R g are the values of the R in the minimum of the effective potential inside the star and in the galactic medium respectively. In those environments, the mass of the field is high, so gravity is close to General Relativity, and we can make the approximation R g ≈ 8πGρ g and R ≈ 8πGρ . From the effective potential we find that the range is In this work we investigate two different f (R) models. The first one is the model proposed by Hu and Sawicki [28]. For this model we have in the high curvature regime, where f R0 is the present value of f R , n is a positive constant and R 0 = 3H 2 0 (4 − 3Ω m0 ) is the effective present value of the background Ricci scalar. H 0 is the present value of the Hubble parameter and Ω m0 is the matter density parameter.
The second model is the γ gravity model proposed in [29]. In this case we have where R * is a positive constant, which is related to cosmological background by R * = 6nH 2 0 (1 − Ω m0 )/αΓ(1/n). These models have the following free parameters: f R0 , n for Hu-Sawicki; α, n for γ gravity.

Symmetron
The Symmetron is a matter-coupled scalar field which effectively decouples from matter in high-density regions through a symmetry restoration. This is done by introducing a potential on the symmetry breaking form [20,31] Figure 5: Left: Constraints of Symmetron; the blue region is the ruled out region for λ 0 1AU, the green and orange dasehd lines are the bounds for smaller values of λ 0 , are ruled out by 95% C.L. Red line indicates the approximate bounds from the Cassini mission [30]. The two curves in the right side are the values ruled out by the torsion pendulum experiment in [31]. The GR is recovered when λ → ∞ or µ 2 M 2 → ρ g . Right: Constraints of Vainshtein screening, the blue region is ruled out by 95% C.L., the black and the red dashed lines are the bounds for the cubic only and the maximally quartic cases, respectively, note that the deviations are very small. The GR is recovered when β → 0 or r c → ∞. and a quadratic conformal factor, A(ϕ) 1 + ϕ 2 /2M 2 . Which gives the follow effective potential From the spherical solution we get where This model has the following free parameters: µ, M, λ.

Vainshtein Screening
The Vainshtein screening works for higher dimensionally models, like Galileon or the brane world DGP model [32,33]. The action for these models can be written as a brane (B), where the RG and the matter fields are defined, embedded in a fivedimensional bulk spacetime (G) containing a five-dimensional generalization of GR where the (5) index denotes the quantities in the bulk spacetime. This reduces, in the weak field regime, to [34] where L i is a function of i-th order of φ and S p is the action of the movement of point particles. Minimizing and integrating it around the star we find stability conditions imposes that taking the lower limit (α 4 = 0) we find the cubic only case, which the exact solution, and using that F 5th = φ (r)/2, we find where β = α 2 /3M pl is the coupling, r V = 16Φ S R r 2 c /9β 2 1/3 is the Vainshtein radius and r 2 c = α 3 /M pl Λ 3 is the crossover scale. The upper limit (α 4 = 2α 2 3 /3α 2 ) gives the maximally quartic case The free parameters in these models are: r c and β.

Results
The fig. 2 shows the ruled out region for Yukawa potential in (λ, α = 2β 2 ) parameter space. Our constraints are weaker than the ones obtained when using data from the solar system [14]. In there the authors found an upper-bound for the coupling constant of α ≤ 3.1×10 −11 and α ≤ 5.2×10 −11 for λ = 0.15AU and λ = 0.21AU respectively. While here we find α ≤ 1.5 × 10 −2 and α ≤ 1.4 × 10 −2 in 95% (C.L.) for the same values of λ respectively. The reason for our weaker bounds is due to the fact that there are much more, and better, data of orbits around the Sun than the exoplanets data. The consequence is that their constraints are about nine orders of magnitude better than ours. This difficulty will be softened only by increasing the number of exoplanets, since the error decreases with the number of observations, σ ∼ 1/ √ N. However other models, with screening, can be better constrained because the fifth force strength depends strongly on the environment physical properties.
The results for the Chameleon Ratra-Peebles model are displayed in the top panel of the figure 3, the filled regions are ruled out in 95% C.L.. The bounds from other experiments are also included [22][23][24][25], here n was fixed equal to 1 to compare with those experiments. As one can see from the figure, the exoplanets does not give us the best upperbounds for Λ, in particular the value Λ = 2.25meV, which corresponds to the value of the cosmological constant predicted by the Planck satelite [4], is allowed for any value of M. Therefore impose a background evolution close to ΛCDM is enough to ensure that the screening works in planetary scales. Figure 4 (left) shows the results for f (R) models. In bottom panel the plot the excluded region, in a 95% of confidence level, for Hu-Sawicki model in blue. The exoplanet data provides a less restrictive constraint compared to the constraints from clusters of galaxies [26] (also in figure 4). The authors find log 10 | f R0 | ≤ −4.79 for n = 1 in a 95.4% confidence level, while exoplanets gives log 10 | f R0 | ≤ 2.16. Other methods should give better constraints, such as stellar and gaseous rotation curves in dwarf galaxies [35], which gives log | f R0 | ≤ −6 or distance indicators in the nearby universe [36], which gives log | f R0 | ≤ −6.4. Figure 4 (Right) shows the excluded region in parameter space for γ gravity. These constraints are weaker than the ones from cosmological bounds, the upper-bounds for n should be less than 0.25 for any value of α. According the [29] these values should not generate a final de Sitter attractor. This condition can also exclude the allowed region for small values of n (n 0.04). The constraints from exoplanets data allow a large region in the parameter space of these models because the information from background cosmology used already ensures that the screening works very well. Fig. 5 (left) shows the results for Symmetron. For this models the screening is automatically satisfied, if one chooses µ and M such as φ 2 g = 0, however we explore cases which deviate from this choices. The blue region is the ruled out region for λ 0 1AU, the green and orange dasehd lines are the bounds for smaller values of λ 0 . The red line represents the approximated bounds for the Sun surface, data from the Cassini Satelite [13], given by the relation A(φ g )−1 5×10 −4 Φ [30] for λ 0 = 1AU, where Φ 2 × 10 −6 is the surface gravitational potential of the Sun. Notice that our results are more restrictive, probing the power of the method especially because the field is now proportional to the square root of the deviation, i.e. φ g ∝ √ ε. We also include the bounds from a torsion pendulum experiment [31] for comparison. The torsion pendulum gives an upper-bound for M, while the exoplanets gives us a lower-bound. For the Vainshtein screening mechanism, Fig. 5 (right), we allow r c and β to be free parameters. The lower-bound for β is larger as r c is bigger as shown in fig. 5. For r c ∼ 1AU (planetary orbits), β must be higher than 10, in atom scale (r c ∼ 1fm) it reduces to 1.7 × 10 −20 and to 2.5 × 10 −39 for Planck scale (r c ∼ pl ). The choice of a very large values for β or r c (β > 10 or r c > 4 × 10 4 AU, which includes cosmological scales r c ∼ H −1 0 ) is enough to make the screening works in planetary scales, suppressing any modifications to General Relativity. The difference between the cubic only and the maximally quartic cases are very small.

Conclusions and perspectives
In this Letter we show a promising, novel and complementary method to constrain modified gravity using extrasolar planets data. Using the orbits of exoplanets we are able to put tighter constraints on the Chameleon models and the Symmetron. These bounds are comparable and competitive with the ones from other experiments, because exoplanets systems have a wide range of masses for the host star, resulting in different screening scales and effectiveness. For Vainshtein screening we find a relation between the free parameters (β, r c ) beyond the trivial choice r c ∼ H −1 0 .