Interiors of Terrestrial Planets in Metric-Affine Gravity

Abstract

Using a semiempirical approach, we show that modified gravity affects the internal properties of terrestrial planets, such as their physical characteristics of a core, mantle, and core–mantle boundary. We also apply these findings for modeling a two-layer exoplanet in Palatini $f\left(R\right)$ gravity.

1. Introduction

Discoveries of exoplanets in the Milky Way Galaxy [1,2] and in the Whirlpool Galaxy [3], as well as growing observational datasets of those objects provided by the current and future missions [4,5,6,7,8], have increased the need for theoretical tools which allow us to describe the planets’ interiors and eventual habitable properties on the basis of those data. A common approach is to extrapolate the Preliminary Reference Earth Model (PREM) [9] and its later improvements [10,11,12] (see more at [13]). Therefore, although an Earth-like planet should have at least six differently composed layers, one usually considers two [14]: iron core and silicate mantle, as they have the biggest impact on the observed properties, such as the planet’s mass, radius, and polar moment of inertia. However, a very different composition of the rocky planets may also be possible, as argued in [15], such as quartz-rich mantles, in comparison to the Solar System ones, whose mantles are mainly made of silicates. Clearly, such findings call for more research in planetary physics, not only from an observational point of view, but also in theoretical modeling.

Regarding the planet’s modeling based on PREM, we are still improving our knowledge on the deepest zones of the Earth, as well as the instrumentation and methods used becoming ameliorated, allowing us to obtain more accurate data regarding the planet’s interior. For instance, recent seismic observation [16] has revealed the existence of a liquid/mushy region of the inner core instead of the solid one, as has been believed so far. On the other hand, a new generation of the neutrino telescopes will be settled to provide information on matter density inside the planet, and on characteristics and abundances of light elements in the outer core [17,18,19,20]. In addition, in laboratories, with the use of lasers [21], the high pressures and temperatures, that is, the extreme conditions of the Earth’s core, are recreated in order to understand the properties and behavior of iron, which is the main element of planets’ cores. All those revelations make the research regarding planets’ modeling relevant, especially agreeing with the fact that various models of gravity predict different layers’ structure in comparison to the Newtonian model [22], commonly used in planetary science. Therefore, knowing the planet’s profile with high accuracy, that is, the number of differently composed layers and their thicknesses, might be another tool to test theories of gravity (see the details on the method in [23]).

As already mentioned, some extensions of Einstein’s theory of gravity may impact the internal structure of the rocky planets, as well as their properties [22,23,24]. This is so, as such theories modify the nonrelativistic hydrostatic equilibrium equations [25] and others, crucial for stellar and planetary modeling. For instance, in the Schwarzschild criterion, which is used to constitute a type of the energy transport through an astrophysical object, there appears an additional term making the star more or less stable with respect to convective processes [26]; energy production in a stellar core is also affected [27,28,29,30], as well as stars’ evolution [31], or cooling processes of substellar objects [32]. Therefore, modified gravity theories proposed to provide some explanations of dark matter and dark energy phenomena [33,34,35,36,37,38], spacetime singularities [39], extreme masses of compact objects [40,41,42,43,44,45], or to unify all four interactions into a single theory [46,47], also impact modeling of gravitational systems for which full relativistic description is not necessary.

One of such theories we are interested in is a subclass of the so-called Ricci-based theories [48], that is, Palatini $f\left(R\right)$ gravity. The main geometric property of these theories is that the metric and connection are considered as independent objects in comparison to most extensions of Einstein’s theory. In Ricci-based gravities, the connection is not coupled to the matter fields, assuming that we are dealing with metric theories, that is, the particles are moving along geodesics distinguished by the metric; moreover, in such a formulation, the connection is not dynamical. However, their most important feature is related to their vacuum dynamics, as it provides the same dynamical equations as general relativity ones with a cosmological constant [49,50,51], providing that those proposals pass the Solar System tests [52] and gravitational waves’ observations as the waves are moving with the speed of light in those theories. However, the difference is clear when one deals with matter fields—Ricci-based gravities then introduce terms which depend on energy density, modifying the structural equations [53].

In this work, we focus on a gravitational model which introduces a quadratic Ricci scalar term, and it is considered in the Palatini approach. Since those terms contribute to the structural equations of spherical symmetric low-temperature spheres, such a modification will have an influence on internal properties of the planet. Therefore, using an analytical method allowing us to obtain the core and core–mantle boundary values of pressure from given masses and radii of transiting exoplanets, we demonstrate that those values will differ in modified gravity. Moreover, we also use them to model an exoplanet interior.

2. Simple Model of Small Rocky Planets in Palatini Gravity

In this section, we recall the hydrostatic equilibrium equations for a cold, low-mass spherical symmetric object. Our terrestrial planets, with masses from the range $M_p\in (0.1-10)M_{\oplus }$, where $M_{\oplus }$ is the Earth’s mass, and core mass fraction (CMF), defined as

$$\mathrm{CMF}=\frac{M_{\mathrm{core}}}{M_p},$$

(1)

not exceeding ∼0.4 of the total planet’s mass1, will be modeled as a two-layer planet, that is, consisting of an iron core and a silicate mantle. Then, using the semiempirical expression relating the CMF with the radius and mass of a transiting exoplanet, we derive the planet’s internal characteristics, such as core pressure and density, their boundary values between the core and mantle, and the mantle’s ones.

2.1. Planets’ Structure Equations

Nonrelativistic hydrostatic equilibrium equations for the quadratic Starobinski model2

$$f\left(\mathcal{R}\right)=\mathcal{R}+\beta {\mathcal{R}}^2$$

(2)

considered in the Palatini approach are given by [53,56]

$$\begin{array}{cc}\hfill p^{\prime }\left(r\right)& =-\frac{Gm\rho }{r^2}\left(1-\beta c^2{\kappa }^2(5\rho -2r{\rho }^{\prime })\right),\hfill \end{array}$$

(3)

$$\begin{array}{cc}\hfill m\left(r\right)& ={\int }_0^{r}4\pi {\tilde{r}}^2\rho \left(\tilde{r}\right)d\tilde{r},\hfill \end{array}$$

(4)

where prime denotes the derivative with respect to the radial coordinate. The matter part of the full relativistic field equation is described as a perfect fluid in this approach, $T_{\mu \nu }=(\rho +p)u_{\mu }{u}_{\nu }+p{g}_{\mu \nu }$, where $p=p\left(r\right)$ and $\rho =\rho \left(r\right)$ are pressure and energy density, respectively, while $u^{\mu }$ is a normalized 4-vector, representing an observer comoving with the fluid. Let us notice that the different numerical factors appearing in the modification term in Equation (3) are the results of the considered assumptions; for example, in [56], the equations were obtained by assuming the conformal invariance of the standard polytropic equation of state for the quadratic model demonstrated in [57], while the equations derived in [53] are more general, without adopting any equation of state. In this work, we also use some polytropic equations of state, however it differs slightly with respect to the common one, as explained later in the text.

Our small rocky exoplanet is modeled as a cold sphere consisting of two different layers. As already mentioned, the material they are made of is iron in the core and silicate in the mantle, whose equations of state are given by the Birch equation of state [58,59], working well when temperatures can be considered uniform but less than 2000 K, and when pressure is below 200 GPa. However, in order to be able to consider more massive objects than the terrestrial planets of the Solar System, one has to take into account the electron degeneracy, as the internal pressure can be $p\gtrsim {10}^4\mathrm{GPa}$. The usual procedure is to match this equation of state with the Thomas–Fermi–Dirac one [60,61,62,63], which also qualifies to describe density-dependent correlation energy [64] which appears because of the interactions between electrons when they obey the Pauli exclusion principle and move in the Coulomb field of the nuclei. Such a hybrid equation of state is very well approximated by a modified polytropic equation of state (SKHM) of the form [14]

$$\rho \left(p\right)={\rho }_0+c{p}^n,$$

(5)

whose best-fit parameters ${\rho }_0$, c, and n for iron and silicate (Mg, Fe)SiO${}_3$ are provided in

Table 1. Best-fit parameters for the SKHM equation of state (Equation (5)) obtained in [14].

Material	${\mathit{\rho }}_0$ (kg m${}^{-3}$)	c (kg m${}^{-3}$ Pa${}^{-\mathit{n}})$	n
Fe($\alpha$)	8300	0.00349	0.528
(Mg, Fe)SiO${}_3$	4260	0.00127	0.549

. Because solids and liquids are incompressible at the low-pressure regimes, the additional term ${\rho }_0$ is present to include this effect. Equation of state constructed in such a way is valid up to $p={10}^7\mathrm{GPa}$, giving the maximal value of the central pressure possible in our analysis.

Moreover, to explore the model with the described features, one needs to establish the initial and boundary conditions. In previous works, we have used the shooting method in order to find the initial values of the core’s densities as well as between the layers’ ones [22,23]. This demonstrated that modified gravity can have a significant impact on those values, and this is a result of different physical assumptions to, for example, Newtonian physics. Therefore, even slight modification to the standard hydrostatic equilibrium equation will have an effect on the internal structure. Keeping this in mind, we have restudied a simple but reasonable method [65] used to obtain the internal characteristic of a distant planet, whose mass and radius can be found by the use of the transit observation techniques [66]. Therefore, for the given total mass of the planet and its radius, we derive the central pressure, its value on the core–mantle boundary (CMB), and the mantle one. This will show that modified gravity indeed affects them.

2.2. Internal Structure of Palatini Planets

There is only one planet whose interior structure and materials, that is, equations of state, are known: Earth3. The many-layers structure, their thickness, and equations of state are given by seismic data, that is, PREM [9]. Since some planets of our Solar System and exoplanets are alike dense and possess similar other characteristics, one usually extrapolates the Earth’s model to describe them. Therefore, extrapolating the Earth’s model, one may derive the semiempirical expression for the core mass fraction (CMF) which carries the information on the core–mantle boundary, often used in numerical procedures and simulations of very distant planets, whose mass $M_p$ and radius $R_p$ are given by the transit. Such a relation between CMF and observed radius and mass was given in [68]:

$$\mathrm{CMF}=\frac{1}{0.21}\left[1.07-\left(\frac{R_p}{R_{\oplus }}\right)/{\left(\frac{M_p}{M_{\oplus }}\right)}^{0.27}\right],$$

(6)

where $R_{\oplus }$ and $M_{\oplus }$ are Earth’s radius and mass, respectively. Furthermore, CMF can be also used to obtain the approximated value for the core radius fraction (CRF), defined as

$$\mathrm{CRF}=\frac{R_{\mathrm{core}}}{R_p},$$

(7)

which is also suitable for numerical analysis [65]:

$$\mathrm{CRF}\approx \sqrt{}\mathrm{CMF}.$$

(8)

Using these two values, that is, CMF and CRF, we derive the core’s and mantle’s pressure, as well as its boundary value, for an exoplanet of the mass $M_p$ and radius $R_p$.

Let us firstly use the definition of local gravity, usually defined as

$$g=\frac{Gm\left(r\right)}{r^2},$$

(9)

to rewrite Equation (3) in a more suitable form for further purposes:

$$p^{\prime }\left(r\right)=-g\rho \left(1-\alpha \left[\frac{14g+g^{\prime }r-2g^{″}{r}^2}{4\pi Gr}\right]\right),$$

(10)

where we have defined a new parameter $\alpha ={\kappa }^2{c}^2\beta$ for the further convenience. Using the mass equation (Equation (4)), together with the expression for the local gravity (9), it can be transformed into

$$\frac{dp}{dm}=-\frac{g^2}{4\pi G}\frac{d\ln \left(m\right)}{dm}\sigma ,$$

(11)

where $\sigma =1-\alpha \left[\frac{14g+g^{\prime }r-2g^{″}{r}^2}{4\pi Gr}\right]$ while $\ln \left(m\right)$ is the natural logarithm of m. Assuming that the surface pressure is zero, we integrate the above equation from the surface inward, such that

$${\int }_{\mathrm{surface}}^{\mathrm{interior}}dp=-\frac{1}{4\pi G}{\int }_{M_p}^{\mathrm{mass\; enclosed\; inside}}{g}^{2}d\ln \left(m\right)\sigma .$$

(12)

Before proceeding further, let us define the surface gravity $g_s$ as a local gravity on the planet’s surface with mass $M_p$ and radius $R_p$

$$g_s:=\frac{G{M}_p}{R_p^2},$$

(13)

while the so-called typical pressure $p_{\mathrm{typ}}$ is defined as

$$p_{\mathrm{typ}}:=\frac{g_s^2}{4\pi G}=\frac{G{M}_p^2}{4\pi R_p^4}.$$

(14)

As the local gravity of the mantle can be assumed to be a constant [65], we may integrate Equation (12) to obtain the pressure of the mantle:

$$\begin{array}{c}\hfill p_{\mathrm{mantle}}=2p_{\mathrm{typ}}\ln \left(\frac{R_p}{r}\right)\left[1+\alpha \frac{7g_s}{\pi G}\frac{M_p}{R_p}\left(\frac{1}{\sqrt{M_p}}-\frac{1}{\sqrt{m}}\right)\right],\end{array}$$

(15)

where we use the planet’s characteristics defined before. In particular, the pressure on the core–mass boundary (CMB) can be obtained by inserting $r\to R_{\mathrm{core}}$ and $m\to M_{\mathrm{core}}$ such that

$$\begin{array}{c}\hfill p_{\mathrm{CMB}}=p_{\mathrm{typ}}\ln \left(\frac{1}{\mathrm{CMF}}\right)\left[1+\alpha \frac{7g_s\sqrt{M_p}}{\pi G{R}_p}\left(1-\frac{1}{\sqrt{\mathrm{CMF}}}\right)\right],\end{array}$$

(16)

where we have used Equations (7) and (8).

On the other hand, since in our model the core density ${\rho }_{\mathrm{core}}$ can be assumed to be a constant value, the core mass is given as $M_{\mathrm{core}}=\frac{4}{3}\pi R_{\mathrm{core}}^3{\rho }_{\mathrm{core}}$. Therefore, the hydrostatic equilibrium Equation (3) can be written with the use of Equation (4) as

$$\begin{array}{c}\hfill \frac{d{p}_{\mathrm{core}}}{dr}=-g{\rho }_{\mathrm{core}}\left(1-\alpha \left[\frac{9m^{\prime }}{4\pi r^2}-\frac{m^{″}}{2\pi r}\right]\right)=-\frac{3r{p}_{\mathrm{typ}}}{R_{\mathrm{core}}^2}\left[1-\alpha \frac{15g_s}{4\pi G{R}_{\mathrm{core}}}\right].\end{array}$$

(17)

Integrating the above equation results in

$$p_{\mathrm{core}}\left(r\right)=p_0-\frac{3}{2}{p}_{\mathrm{typ}}{\left(\frac{r}{R_{\mathrm{core}}}\right)}^2\left(1-15\alpha \frac{g_s}{4\pi G{R}_{\mathrm{core}}}\right),$$

(18)

where $p_0$ is the central pressure which can be determined by matching the above $p_{\mathrm{core}}$ at CMB with the pressure on the boundary (Equation (16)):

$$\begin{array}{cc}\hfill p_0& =p_{\mathrm{CMB}}+\frac{3}{2}{p}_{\mathrm{typ}}\left(1-15\alpha \frac{g_s}{4\pi G{R}_{\mathrm{core}}}\right)\hfill \\ \hfill & =p_{\mathrm{typ}}\left(\frac{3}{2}\left[1-15\alpha \frac{g_s}{4\pi G{R}_{\mathrm{core}}}\right]+\ln \left(\frac{1}{\mathrm{CMF}}\right)\left[1+\alpha \frac{7g_s\sqrt{M_p}}{\pi G{R}_p}\left(1-\frac{1}{\sqrt{\mathrm{CMF}}}\right)\right]\right).\hfill \end{array}$$

(19)

The above result allows us to find an approximated value of the central pressure for a given terrestrial exoplanet whose mass and radius are provided by the transit observations. The effect of modified gravity is clearly present; therefore, in the next section, we numerically solve the structural equations with the use of those findings.

3. Numerical Solutions

In order to compare models of different values of the Starobinsky parameter $\beta$, we have introduced earlier a dimension-full parameter $\alpha =c^2{\kappa }^2\beta$, which allows one to write the formulas in a more convenient way. We chose four values of the parameter, $\alpha \in \{0,{10}^{-15},{10}^{-14},{10}^{-13}\}\times {\mathrm{m}}^3/\mathrm{kg}$, that is, $\beta \in \{0,{10}^{11},{10}^{12},{10}^{13}\}\times {\mathrm{m}}^2$, neglecting the possibility of negative values of the parameter4. Let us comment that with the current experiments performed in the Solar System, one is not able to distinguish between Palatini $f\left(R\right)$ gravity and GR [52]. Upper bounds on the absolute value of the parameter $\beta$ in the Palatini approach have been determined to be $\beta \sim {10}^{12}$ cm${}^2$, when one investigates weak-field limit of the theory [69], or $\beta \sim {10}^9$ cm${}^2$, if gravitational forces become as strong as electric forces [70]. It must be noted that the absolute bounds differ from the ones obtained for metric counterpart of the theory [69]. Another issue concerns the critical value of the parameter. It can be shown that, when $\beta$ takes negative values, there also exists an upper bound on the value of the parameter depending on the energy density of matter; any values above the critical one lead to unphysical results [23].

Having established the range of the parameter, we focus on numerically solving Equations (3) and (4), supplemented with the equation of state, Equation (5). In order to determine the exact density profile for a planet of given mass and radius, and chemical composition, one needs to use the shooting method, i.e., find a value of internal density such that, at the surface of the planet, the radius and the mass coincide with the desired values (the surface, if defined by $p\left(R_p\right)=0$). The fact that the masses and radii of the planets we examine are fixed by the transit observations provides the possibility to determine the core density and core size with its mass, as well as to plot the density profiles. As one can see in Figure 1, all curves denoting solutions for different values of $\alpha$ end at the same point; what changes is the size of the core. This allows us to compare CMFs and CRFs obtained from the quasi-empirical Formula (6) (which is constant once the mass and the radius of the planet are given) to the numerical findings.

As far as the pressure is concerned, we simply calculate it for one planet, Kepler-78 b5, using Formulas (15) and (18), as well as the exact value of $R_{\mathrm{core}}$ determined in the previous numerical step. The results are shown in Figure 2, illustrating the effects of modified gravity on pressure within the exoplanet. The analytical solutions are then compared with numerical ones, to determine how good the approximations are. The results are shown in Figure 3 for two values of $\alpha$.

4. Conclusions

Previous studies regarding terrestrial planets in modified gravity [22] revealed that extensions of Einstein’s gravity alter the internal structure of those objects, providing a possibility to test such theories with the use of seismic data [23]. Therefore, the physical quantities, such as core pressure and energy density, as well as their boundary values between layers, should also be affected, which would have an impact on the way we model distant planets, where seismology cannot be applied. This fact forces us to look for methods allowing us to find those values, when only the observed characteristics, such as mass and radius of a transiting exoplanet, are available. In this work, we wanted to check if such methods are model-independent.

As clearly demonstrated, the methods can indeed depend on the applied theory of gravity. For this analysis, we considered quadratic modification to the general relativity’s Lagrangian Equation (2), considered in the Palatini approach; however, our conclusions are valid for other theories of gravity which modify the nonrelativistic limits of their field equations.

• Density profiles, as already noticed in our previous works, can significantly differ in modified gravity with respect to the Newtonian model. We observe not only lower values of central density and on the core–mantle boundary, but also the cores of the given exoplanets are bigger; that is, the cores are less dense in the case of Palatini gravity. Therefore, the observed transiting planets can have different structure for the same masses and radii than the one predicted in the usual way, and can affect the planet’s polar moment of inertia. The fact that internal structure of planets is affected by modifications of gravitational interaction is to be expected, since Equations (3) and (4), allowing one to compute the density profiles, change. This entails the fact that modifications of gravity introduce additional degeneracy when trying to determine planets’ internal composition by looking at the mass–radius relationship [23]. The values of internal pressure and core radius, giving the same total mass and radius, depend on the parameter $\alpha$. Therefore, this fact alone does not allow us to constrain alternative gravity models. What actually could help in distinguishing between different models would be collecting seismic data from Solar System planets, and investigating their density profiles. For example, Earth’s mass and radius are well known, as well as its internal composition, so, after having developed a more realistic model taking into account modifications of gravity, it will be possible to place a stringent constraint on values of $\alpha$.
• A similar situation happens when we plot the pressure curves obtained in this work: its central values decrease in modified gravity; however, when we approach the planet’s surface, the mantles do not differ much. This result derives from the fact that the additional term in Equation (15) for the pressure in the mantle is small, and smaller than the extra term appearing in the analogous equation for the core (Equation (18)).
• We also compared the numerical solutions for the pressure obtained from Equations (3) and (4) to the ones resulting from the analytical approach (which are approximated solutions). As one can see, the pressure drops roughly, similar to $-r^2$ in the core, and then changes in a linear way in the mantle in the case of both numerical and analytical solution (although it is less pronounced for larger values of the parameter $\alpha$). One notices that in the case of Newtonian gravity ($\alpha =0$), analytical (approximated) solution tends to provide smaller values than the numerical one. However, in the case of modified gravity, the effect is the reverse—approximated analytical solution provides larger values than the numerical one. This can be explained in the following way: the analytical approximation does not take into account the effect of modification of gravity in the CMF Formula (6), so it stays constant for various values of the parameter $\alpha$ (as it depends of the mass and radius of the planet only, and these values do not change). On the other hand, the numerical method suggests that the size of the core and its mass grow in modified gravity, and hence the CMF must change. This combined effect of change in $\alpha$ and CMF/CRF results in a bigger drop in internal pressure.
• Moreover, as already mentioned in the previous point, our numerical analysis revealed that the equation for the semiempirical CMF used in that work also depends on modified gravity. This is not a surprise, remembering the fact that for finding that relation, one uses the PREM model, which is based on Newtonian gravity.

Although our studies presented in this paper are based on crude methods and assumptions, such as spherical symmetric, nonrotating planets, their two-layers structure, and constant values for the mantle’s characteristics, it is evident that alternative theories of gravity do impact the planets’ descriptions and modeling. Improving our analytical and numerical methods, that is, taking into account the missing ingredients mainly related to more realistic planet’s geometry, should also manifest similar results. The work along these lines is currently underway.