On the linear stability of polytropic fluid spheres in $R^{2}$ gravity

Einstein’s general relativity (GR) is a generally recognized theory of gravitation that is successfully applied on different spatial and time scales. In particular, it is used to model processes which took place in the early Universe. It is, however, obvious that to describe the very early stages of the evolution of the Universe, it is already needed to take into account quantum effects. To do this, it is necessary to have a full theory of quantum gravity, which is absent at present. For this reason, starting from 1960’s, it was being suggested to use some effective description of quantum effects in strong gravitational fields by a change of the classical Einstein gravitational Lagrangian $\sim R$ by various modified Lagrangians containing different curvature invariants. In the simplest case this can be some function $f(R)$ of the scalar curvature $R$. Such modified gravity theories (MGTs) have been widely applied to model various cosmological aspects of the early and present Universe (for a general review on the subject, see, e.g., Refs. Nojiri:2010wj ; Nojiri:2017ncd ).

On the other hand, in considering processes and objects on relatively small scales comparable to sizes of galaxies and even of stars, the effects of modification of gravity can also play a significant role. For example, within $f(R)$ gravity, there have been constructed models of relativistic stars rel_star_f_R_1 ; rel_star_f_R_2 , wormholes wh_f_R_1 ; wh_f_R_2 , and neutron stars Cooney:2009rr ; Arapoglu:2010rz ; Orellana:2013gn ; Alavirad:2013paa ; Astashenok:2013vza ; Ganguly:2013taa . It was shown that the modification of gravity can affect, in particular, a number of important physical characteristics of neutron stars which may be verified observationally. One of the problems of interest is determining the mass-radius (or the mass-central density of matter) relations, which have been obtained, for example, in Refs. Astashenok:2014pua ; Astashenok:2014dja ; Capozziello:2015yza ; Bakirova:2016ffk ; Astashenok:2017dpo ; Folomeev:2018ioy ; Feola:2019zqg ; Astashenok:2020isy . It was shown there that, as in GR, in MGT, for some central density of neutron matter, the mass-central density curves possess maxima whose location depends on the physical properties of the specific matter. In GR, the presence of such a maximum indicates the fact that there is a transition from stable systems (located to the left of the maximum) to unstable configurations (located to the right of the maximum). This fact is established by studying the behavior of matter and metric perturbations using the variational approach Chandrasekhar:1964zz .

In MGTs, different types of perturbations have been repeatedly studied as well (see, e.g., Refs. Blazquez-Salcedo:2018pxo ; Blazquez-Salcedo:2020ibb and references therein). In considering these problems, a transition from $f(R)$ gravity to a scalar-tensor theory is usually performed. Correspondingly, studies of the perturbations are carried out not in the Jordan frame but in the Einstein frame. However, the question of the physical equivalence between these two frames is still under discussion, and it cannot be regarded as completely solved Capozziello:2010sc ; Kamenshchik:2014waa ; Kamenshchik:2016gcy ; Bahamonde:2016wmz ; Ruf:2017xon . Here, the following potential difficulties may be noted: (i)  In performing the transition from $f(R)$ gravity to a scalar-tensor theory, there may, in general, occur some undesirable consequences (singularities, fixed points, etc.); as a result, the equivalence between the frames can be violated. (ii) Objects that are stars in one frame may represent some other configurations in another. (iii) In constructing perturbation theory, the equivalence between the frames can be lost in view of the approximate nature of such a theory.

In this connection, it may be of some interest to study the stability directly in the Jordan frame, and this is the goal of the present paper. For the sake of simplicity, we will work within $R^{2}$ gravity, where linear radial perturbations of a strongly gravitating system supported by a polytropic fluid will be investigated. For this purpose, in Sec. II, we first derive the general equations for $f(R)$ gravity. Using these equations, in Sec. III, we numerically find static solutions describing equilibrium configurations, on the background of which the behavior of matter and spacetime perturbations is studied in Sec. IV.

We consider modified gravity with the action [the metric signature is $(+,-,-,-)$]

where $G$ is the Newtonian gravitational constant, $f(R)$ is an arbitrary nonlinear function of $R$, and $S_{m}$ denotes the action of matter.

For our purposes, we represent the function $f(R)$ in the form

where $h(R)$ is a new arbitrary function of $R$ and $\alpha$ is an arbitrary constant. When $\alpha=0$, one recovers Einstein’s general relativity. The corresponding field equations can be derived by varying action (1) with respect to the metric, yielding

Here $G_{i}^{k}\equiv R_{i}^{k}-\frac{1}{2}\delta_{i}^{k}R$ is the Einstein tensor, $h_{R}\equiv dh/dR$, and the semicolon denotes the covariant derivative.

To obtain the modified Einstein equations and the equation for the fluid, we choose the spherically symmetric metric in the form

where $\nu$ and $\lambda$ are in general functions of $r,x^{0}$, and $x^{0}=c\,t$ is the time coordinate.

As a matter source in the field equations, we take an isotropic fluid with the energy-momentum tensor

where $\varepsilon$ is the fluid energy density and $p$ is the pressure.

The trace of Eq. (3) yields the equation for the scalar curvature

where $T$ is the trace of the energy-momentum tensor (5).

Using Eqs. (4) and (5), the $(^{t}_{t})$, $(^{r}_{r})$, $(^{\theta}_{\theta})$, and $(^{r}_{t})$ components of Eq. (3) can be written as

where the dot and prime denote differentiation with respect to $x^{0}$ and $r$, respectively. In turn, Eq. (6) yields

Finally, the $i=r$ component of the law of conservation of energy and momentum, $T^{k}_{i;k}=0$, gives

Consider that the equilibrium systems described above are perturbed in such a way that spherical symmetry is maintained. In obtaining the equations for the perturbations, we will neglect all quantities which are of the second and higher order. The components of the four-velocity in the metric (4) are given by Chandrasekhar:1964zz

with the three-velocity $\mathpzc{v}=dr/dx^{0}\ll 1$. The index 0 in the metric functions indicates the static, zeroth order solution of the gravitational equations.

Now we consider perturbations of the static solutions of the form

where the index 0 refers to the static solutions, the index $p$ indicates the perturbation, and $y$ denotes one of the functions $\lambda,\nu,\varepsilon,p$ or the scalar curvature $R$. We will use the variational approach of Ref. Chandrasekhar:1964zz when one introduces a “Lagrangian displacement” $\zeta$ with respect to $x^{0}$, $\mathpzc{v}=\partial\zeta/\partial x^{0}$. Then, substituting the expansions (24) in Eqs. (7)-(12) and seeking solutions in a harmonic form

[for convenience, we hereafter drop the tilde sign on $\tilde{y}_{p}(r)$], one can obtain the following set of equations for the perturbations $\theta_{p},\Sigma_{p}$, and $\lambda_{p}$:

where $s_{1}=1+\alpha\left(2\Sigma_{0}+x\Sigma_{0}^{\prime}\right),s_{2}=1+\sigma(n+1)\theta_{0}$, and the dimensionless frequency $\bar{\omega}=L\omega$. Here, Eq. (IV) follows from the conservation law (12), Eq. (IV) follows from the equation for the scalar curvature (11), Eq. (IV) is the $(^{t}_{t})$ component (7). In deriving Eq. (IV), we have used the expression (here $\psi=\zeta/L$ is the dimensionless Lagrangian displacement)

[which follows from the $(^{r}_{t})$ component (10)], expressing from it $\psi$ and eliminating it from (IV).

For this set of equations, we choose the following boundary conditions near the center $x=0$:

where the expansion coefficients $\theta_{pc}$ and $\Sigma_{pc}$ are arbitrary and $\theta_{p2},\lambda_{p2}$, and $\Sigma_{p2}$ can be found from Eqs. (IV)-(IV).

The set of equations (IV)-(IV) together with the boundary conditions (28) defines an eigenvalue problem for $\bar{\omega}^{2}$. The question of stability is therefore reduced to a study of the possible values of $\bar{\omega}^{2}$. If any of the values of $\bar{\omega}^{2}$ are found to be negative, then the perturbations will increase and the configurations in question will be unstable against radial oscillations.

The choice of eigenvalues of $\bar{\omega}^{2}$ is carried out such that we have asymptotically decaying solutions for the perturbations $\Sigma_{p}$ and $\lambda_{p}$. In doing so, it is necessary to ensure the following properties of the solutions: (i) The function $\theta_{p}$ must be finite (though not necessarily zero) at the boundary of the star. This is sufficient to ensure that the perturbation of the fluid pressure $p_{p}\sim\theta_{0}^{n}\theta_{p}$ meets the condition $p_{p}=0$ at the edge of the star where $\theta_{0}=0$ [see, e.g., Eq. (60) in Ref. Chandrasekhar:1964zz ]. (ii) The function $\lambda_{p}$ must be nodeless; this corresponds to a zero mode of the solution (on this point, see, e.g., Ref. Gleiser:1988ih ).

In this connection, it is useful to write out the asymptotic behavior of the solutions.

Examples of solutions for the perturbations $\theta_{p},\Sigma_{p}$, and $\lambda_{p}$ are shown in Fig. 2. The procedure for determining the eigenvalues of $\bar{\omega}^{2}$ is as follows:

1. (1)

   In the case of GR (when $\bar{\alpha}=0$), there are only two equations (IV) and (IV) for the functions $\theta_{p}$ and $\lambda_{p}$. Considering that these equations are linear, the rescaling of the central $\theta_{pc}\to\beta\theta_{pc}$ results in the corresponding rescaling of the metric perturbation $\lambda_{p}\to\beta\lambda_{p}$, but the qualitative behavior of the solutions remains unchanged. That is, it is possible to take any central $\theta_{pc}\ll 1$, and the eigenfrequency $\bar{\omega}^{2}$ will not change.

2. (2)

   In the case of $R^{2}$ gravity (when $\bar{\alpha}\neq 0$), all three equations (IV)-(IV) need to be solved. In doing so, there are two arbitrary central values $\theta_{pc}$ and $\Sigma_{pc}$. Numerical calculations indicate that to ensure regular asymptotically decaying solutions for the perturbations $\Sigma_{p}$ and $\lambda_{p}$ one has to adjust the values both of the eigenfrequency $\bar{\omega}^{2}$ and of one of these two arbitrary parameters. That is, either $\theta_{pc}$ or $\Sigma_{pc}$ is an eigenparameter of the problem. The corresponding numerical values of these parameters are given in Table 1 for several values of $\sigma$. It is seen from the table and Fig. 1 that the square of the lowest eigenfrequency $\bar{\omega}^{2}$ is positive to the left of the maximum and negative to the right of it. That is, as in GR, in $R^{2}$ gravity under consideration the transition from stable to unstable systems occurs strictly at the maximum of the mass.

Summarizing the results obtained, within $R^{2}$ gravity, we have examined the question of stability of compact configurations supported by a polytropic fluid against linear radial perturbations. In contrast to the studies performed earlier in the literature, here the calculations have been carried out in the Jordan frame to avoid the potential difficulties related to the conformal transformation to the Einstein frame (see Introduction). In doing so, we regard the scalar curvature as a dynamical variable for which the behavior of the corresponding perturbation modes is studied as well. As a result, it is shown that, as in GR, within the framework of $R^{2}$ gravity the transition from stable to unstable configurations takes place at the point of the maximum of the curve mass-central density of the fluid. One may expect that similar results will also be obtained for another, more realistic EoSs of matter, including those that are used in constructing models of neutron stars (see, e.g., Refs. Astashenok:2017dpo ; Folomeev:2018ioy ; Astashenok:2020isy ).

The authors are very grateful to S. Odintsov for fruitful discussions and comments. We gratefully acknowledge support provided by Grant No. BR05236322 in Fundamental Research in Natural Sciences by the Ministry of Education and Science of the Republic of Kazakhstan. We are also grateful to the Research Group Linkage Programme of the Alexander von Humboldt Foundation for the support of this research.