Karmarkar−Tolman embedded charged anisotropic stars in f(\mathcal{R}) gravity

Abstract

We investigate various anisotropic spherical distributions of charged celestial bodies within the context of f(\mathcal{R}) gravity, where \mathcal{R} represents the Ricci scalar. The properties of specific charged compact objects are analyzed by using the Karmarkar−Tolman spacetime and three distinct gravitational models. The behavior of the structural parameters is examined via graphical methods. Energy constraints are applied to assess how well the results align with the Karmarkar−Tolman spacetime model. The physical acceptability of the stellar models is evaluated by checking the energy conditions and the equation of state parameter. Additionally we explore the influence of anisotropy on the stability and internal structure of the models. Our findings are compared with predictions from general relativity to highlight the effects of f(\mathcal{R}) gravity on charged compact stars. The obtained results are useful to enhance our understanding of how modified gravity theories affect the properties of compact astrophysical objects.

Graphical abstract

Keywords

${\color{khaki}{ f({\cal{R}}) }}$ gravity / gravitational models / charged compact geometries / stability

1 Introduction

Astrophysical compact star objects, including white dwarfs, massive neutron stars and black holes which originate from gravitational collapse driven by relativistic degeneracy pressure, play a crucial role in astrophysics. These celestial bodies possess exceptionally high densities because of their substantial masses, despite their relatively small volumes. Although their exact properties remain uncertain, compact stars are generally characterized as massive stars with extremely small radii except for black holes, all such objects are classified as degenerate stars [1]. Studying the equilibrium configurations of compact stars within the framework of general relativity GR is of great significance, given their high densities and massive nature. Consequently, research on compact stars incorporating both GR and modified gravitational theories MGTs, remains a vital and evolving field in astrophysics. So modeling these objects requires obtaining exact solutions to the Einstein field equations. In this regard various gravitational theories have been proposed to derive complex and precise solutions to the Einstein field equations [2, 3].

In particular Tolman and Oppenheimer within the context of observational data, investigated realistic models demonstrating that the physical characteristics of these objects are the result of the interaction between internal pressure and gravitational forces ultimately leading to equilibrium structures [4, 5]. This phenomenon, which frequently produces accurate results, is very important for studying the internal configurations of stars. Furthermore after observing strongly magnetized spinning neutrons Baade and Zwicky carried out a comprehensive analysis of compact stars and found that supernovae may develop into smaller compact objects [6, 7]. In astrophysical research, the concept of a perfect fluid has been proposed as a precursor to the formation of compact stars and the emergence of stellar bodies. These objects are typically described as ultra-dense and isotropic exhibiting an apparent spherical symmetry. Although isotropy is often regarded as a favorable property, it is not an intrinsic characteristic of all compact stars. The concept of non zero anisotropy in stellar structures was first introduced by Bowers and Liang [8]. Additionally Ruderman [9], observed that the nuclear density at the core of stellar objects demonstrates anisotropic behavior. He further investigated the possibility that stars could attain extremely high densities, i.e., {10}^{15}\mathrm{g}/{\mathrm{c}\mathrm{m}}^3, when nuclear matter exhibits anisotropic properties.

In the modified theory of gravity an additional function f(\mathcal{R}), which can include higher order derivative terms such as {\mathcal{R}}^2 and others is incorporated into the Einstein−Hilbert EH action. Since these terms solely influence the geometric sector of the EH action, it follows that dark components are, at least theoretically purely geometric in nature [10, 11]. MGTs have emerged as a compelling and active field of study in relativistic astrophysics in recent years. Qadir et al. [12] emphasized the importance of incorporating advanced dynamical frameworks in research suggesting that such approaches could offer insights into resolving dark matter DM related challenges, Various MGTs, such as f(G), f(G,T), f(\mathcal{R}) and f(G,T^2), have been explored, where G denotes the Gauss−Bonnet term and T represents the trace of the energy momentum tensor. Recently several mathematical methodologies integrating these concepts have been introduced to enhance our understanding of cosmic phenomena [13-19]. In these works, it is proposed that MGTs could serve as a tool for uncovering various mysteries of our universe. Theoretical studies suggest that during periods of high density, the pressure inside a star is likely to become anisotropic [9, 20]. In recent years numerous researchers have conducted in depth investigations into the role of pressure anisotropy in shaping the structure of the universe [21-27]. The structural properties of compact matter models can be significantly influenced by pressure anisotropy. To enhance our understanding of cosmic phenomena mathematical techniques that integrate these concepts have been introduced.

It is widely recognized that in locally anisotropic systems, the ratio 2M/R (where M represents mass and R denotes radius) can approach unity [28]. However in isotropic models the Buchdahl limit imposes a constraint, ensuring that [29]. Yazadjiev explored the equations of motion incorporating electric charge to analyze the hydrodynamics of relativistic systems. His findings suggested that the resulting solutions exhibited characteristics of a non perturbative magnetar. The presence of a strong magnetic field ultimately induces anisotropy in the matter distribution. By accounting for the effects of the magnetic field in relativistic interiors, ideal fluid configurations can be utilized to construct precise analytical models. Astashenok et al. [30] examined the impact of strong magnetic fields on neutron star populations within the framework of modified Einstein gravity. Their study concluded that the inclusion of quadratic Gauss−Bonnet terms could support more massive neutron stars with a significant strangeness fraction. Recently, gravitational wave observations of a compact binary coalescence have identified a compact object with a mass ranging between 2.50M_{\odot } and 2.67M_{\odot } [31]. They investigated the upper limits of neutron stars using the causal equation of state EoS, framework by numerically solving the corresponding Tolman−Oppenheimer−Volkoff (TOV) equations. It is well known that a neutron star will eventually collapse into a black hole if its baryon mass exceeds the stability threshold for static configurations. Our understanding of the challenging mass gap region in such systems has advanced following the work of Astashenok et al. [32], which explored the formation of these unique structures within {\mathcal{R}}^2 gravity. Furthermore several other studies have examined the impact of MGTs on the potential modeling of stellar structures [33-37]. In a newly proposed scalar covariant theory as well as in vector and tensor gravity models, variations in \omega, \mu and G are permitted. This theory suggests that since the modified gravitational acceleration law can be derived from the equations of motion, galaxy rotation curves and cluster data can be explained without invoking DM. The concept has been successfully tested within our solar system prompting an investigation into the linear evolution of perturbations in relation to cosmic microwave background measurements. However scalar, vector and tensor gravity theories face several challenges, such as the need to align with observational data issues related to fine tuning and potential violations of the equivalence principle. Additionally these theories exhibit high sensitivity to parameter variations and are inherently complex. As a result various gravitational models distinct from GR have emerged over time offering alternative perspectives on gravitational interactions [15, 16]. Einstein’s GR struggles to account for the late time accelerated expansion of the universe without introducing an unknown dark energy component into the field equations [38-41].

In the TOV limit the maximum mass of neutron stars is a crucial subject in astrophysics, as it determines the threshold separating neutron stars from black holes. A theoretical model predicting neutron star masses considerably exceeding observed values would contradict both empirical data and the fundamental physics of dense matter [4, 42]. The maximum possible mass of a neutron star is dictated by the TOV equation derived from GR for a spherically symmetric object in hydrostatic equilibrium,

{\displaystyle \frac{\mathrm{d}P}{\mathrm{d}r}}=-{\displaystyle \frac{G(ϵ+P)(m+4\pi r^{3}P)}{r(r-2Gm)}},

where P represents pressure, ϵ is the energy density, m(r) is the enclosed mass within a radius r, G is the gravitational constant. To solve the TOV equation an EoS describing neutron abundant matter is required. The maximum neutron star mass is influenced by the stiffness of the EoS. There are two limits, a stiffer EoS supports higher maximum masses and a softer EoS imposes lower masses. Also astronomical observations suggest an upper mass limit for neutron stars. The most massive confirmed neutron stars indicate that realistic EoS models should support a mass range of approximately 2.0-2.3M_{\odot }, but not significantly beyond this [43, 44]. For a neutron star to exceed 3M_{\odot }, an exceptionally stiff EoS would be necessary. However most nuclear physics models, whether derived from microscopic approaches like chiral effective field theory or astrophysical observations such as gravitational wave detections suggest that such an extreme EoS would be unrealistic due to fundamental constraints [45-47]. The study of compact stellar objects within modified theories of gravity has gained significant attention due to limitations of GR in explaining high density astrophysical phenomena. Among various analytical approaches, the Karmarkar embedding condition has emerged as a powerful geometric tool for constructing exact solutions of the gravitational field equations [48, 49]. This condition restricts a static spherically symmetric spacetime to embedding class-I, thereby establishing a functional relationship between the metric potentials and reducing the complexity of the field equations [50-52]. Early applications of the Karmarkar condition were developed within GR to model isotropic and anisotropic compact stars. However, recent works have extended this framework to modified gravity theories, including f(\mathcal{R}),f(\mathcal{R},T),f(\mathcal{R},G) and more generalized matter geometry coupling models. These extensions allow the investigation of deviations from GR at high curvature scales relevant to compact stars [53-55].

In this work, the methodology in Refs. [56-60] is adopted to incorporate f(\mathcal{R}) modifications into a static configuration that represents an anisotropic fluid, using the Karmarkar−Tolman space time. The Karmarkar−Tolman embedding technique offers a systematic way to construct physically realistic stellar models, especially when combined with electric charge and anisotropy, which are essential features of real astrophysical objects. Our work provides a comprehensive modeling of charged anisotropic compact stars within f(\mathcal{R}) gravity using three distinct gravitational models. By incorporating the Karmarkar−Tolman condition, we achieve an exact embedding of the spacetime geometry, ensuring more realistic and stable solutions. A notable advantage of our study is the detailed analysis of structural properties through graphical methods, which enhances the understanding of the physical behavior of these stars under different curvature scenarios. Additionally, the validation of our models through energy conditions and the TOV equilibrium equation further strengthens the physical reliability of the proposed configurations.

Previous studies often dealt with either anisotropic stars or charged configurations separately, or did not consider embedding techniques under f(\mathcal{R}) gravity [61-63]. Our study fills this gap by integrating the Karmarkar−Tolman condition within the f(\mathcal{R}) gravity framework, leading to new and physically meaningful exact solutions. We explore and compare three different functional forms of f(\mathcal{R}), offering a broader perspective on how various curvature corrections affect the structure and stability of compact stars, unlike earlier works that usually focus on a single model. Through extensive graphical analysis of energy conditions, anisotropy, mass to radius relations and charge distributions, our study offers a clear and detailed visualization of the results, reinforcing the stability and physical acceptability of the models. This approach provides more transparency and deeper insight compared to many previous purely analytical studies [64, 65]. Unlike many prior models that ignore the role of charge, we show that introducing charge not only affects the equilibrium and stability of compact stars but also helps in constructing more realistic models that align better with observed astrophysical objects. Thus our work presents an important step forward in the modeling of realistic stellar structures in modified gravity, offering several theoretical and methodological improvements over previous literature.

This manuscript is organized as follows. In Section 2.1, general equations for irrotational spheres with f(\mathcal{R}) gravity corrections are presented. Section 2.2 examines specific solutions and discusses their characteristics to illustrate the approach. Section 2.3 discusses the process of matching with Reissner−Nordström exterior metric. Section 2.4 introduces three different models for f(\mathcal{R}) gravity. The physical aspects of f(\mathcal{R}) gravity models are discussed in details in Section 3. The conclusion will be presented in Section 4.

2 Theoretical framework

2.1 Key field equation in f(\mathcal{R}) gravity

The f(\mathcal{R}) theory of gravity is a generalization of Einstein’s general theory of relativity, in which the gravitational action is modified by an arbitrary function f(\mathcal{R}) of the Ricci scalar \mathcal{R}, instead of being simply proportional to \mathcal{R}. This modification allows for more flexibility in modeling the dynamics of spacetime, especially in contexts like cosmology and astrophysics, where modifications to GR might be necessary to explain phenomena such as the accelerated expansion of the universe or dark energy. The modified function in the f(\mathcal{R}) gravity theory can be expressed as a generalization of Einstein’s GR typically, it appears in the form of [66],

(1)I={\displaystyle \int \mathrm{d}{x}^4\sqrt{-g}\left[f(\mathcal{R})+L_e+L_m\right],}

where L_e and L_m represent the Lagrangians for electromagnetic and matter fields. The field equations derived from the action are more general than the standard Einstein field equations. The modified form reads

(2)f_{\mathcal{R}}{\mathcal{R}}_{\mu \nu }-{\displaystyle \frac{1}{2}}{g}_{\mu \nu }f+(g_{\mu \nu }◻-{\mathrm{\nabla }}_{\mu }{\mathrm{\nabla }}_{\nu })f_{\mathcal{R}}=T_{\mu \nu }^{\mathrm{E}\mathrm{M}}+T_{\alpha \beta }^{\mathrm{m}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}},

where {\mathrm{\nabla }}_{\mu } is a covariant derivative, f_{\mathcal{R}}(\mathcal{R}) is the shorthand notation of \frac{\mathrm{\partial }f(\mathcal{R})}{\mathrm{\partial }\mathcal{R}} operators, T_{\alpha \beta }^{\mathrm{E}\mathrm{M}} is the electromagnetic stress energy tensor and T_{\mu \nu }^{\mathrm{m}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}} is the matter tensor. The explicit expressions of the electromagnetic stress energy and matter tensor read

(3)T_{\mu \nu }^{\mathrm{m}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}}=(\rho +P)U_{\mu }{U}_{\nu }-P_t{g}_{\mu \nu }+(P_r-P_t)V_{\mu }{V}_{\nu },

(4)T_{\alpha \beta }^{\mathrm{E}\mathrm{M}}=F_{\alpha }^{\mu }{F}_{\beta \mu }-{\displaystyle \frac{1}{4}}{g}_{\alpha \beta }{F}^{\mu \nu }{F}_{\mu \nu },

where U_{\mu } is the four velocity of the fluid and U_{\mu }{U}^{\mu }=1, V_{\mu } is the radial unit vector and V_{\mu }{V}^{\mu }=-1, \rho is the energy density, the metric tensor g_{\mu \nu } encodes the geometry of spacetime, F^{\mu \nu } is the electromagnetic field strength tensor, P_r and P_t are radial pressure components and transverse pressure components, respectively. The components of T_{\mu \nu }^{\mathrm{m}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}} and T_{\alpha \beta }^{\mathrm{E}\mathrm{M}} are found out respectively as

(5)\mathrm{M}\mathrm{a}\mathrm{t}\mathrm{t}\mathrm{e}\mathrm{r}:T_0^0=\rho ,T_1^1=P_r,T_2^2=r^2{P}_t,

(6)\mathrm{E}\mathrm{M}:T_0^0={\displaystyle \frac{Q^2}{8\pi r^4}},T_1^1={\displaystyle \frac{Q^2}{8\pi r^4}},T_2^2={\displaystyle \frac{r^2{Q}^2}{8\pi r^4}}.

Since we assume spherical symmetry, the time component of the current depends only on the radial coordinate r. In the electromagnetic field tensor, the only non-vanishing components are F^{01} and F^{10}, with the relation F^{01}=F^{10}, representing the radial electric field. The viability of the adopted f(\mathcal{R}) gravity model has been carefully examined to ensure its physical consistency and compatibility with both strong field and weak field regimes. In modified gravity theories, it is well known that physically acceptable models must satisfy certain fundamental conditions, namely the positivity of the first and second derivatives of the gravitational function with respect to the Ricci scalar, f_{\mathcal{R}}\equiv \mathrm{d}f/\mathrm{d}\mathcal{R}>0 and f_{\mathcal{R}\mathcal{R}}\equiv {\mathrm{d}}^{2}f/\mathrm{d}{\mathcal{R}}^2>0. The condition f_{\mathcal{R}}>0 guarantees a positive effective gravitational coupling and the absence of ghost degrees of freedom, while f_{\mathcal{R}\mathcal{R}}>0 ensures the stability of the scalar degree of freedom and avoids the Dolgov−Kawasaki instability. For the specific form of f(\mathcal{R}) considered in this work these derivatives are evaluated explicitly and it is found that both f_{\mathcal{R}} and f_{\mathcal{R}\mathcal{R}} remain strictly positive throughout the entire stellar interior for the chosen parameter range. This confirms that the model is free from pathological instabilities in the strong gravity regime relevant to compact stars.

2.2 Anisotropic matter distribution in f(\mathcal{R}) gravity

In the framework of anisotropic matter distributions, the dynamics of self gravitating systems, including stars, compact objects and galaxies can be examined within f(\mathcal{R}) gravity to explore deviations from the conventional predictions of GR. Anisotropy is essential for understanding the internal structure of compact stars, where the radial pressure differs from the tangential pressure, a phenomenon that arises in extremely dense astrophysical systems due to intense gravitational fields, phase transitions or electromagnetic influences. The general form for the line element for dynamical spherically symmetric geometry is

(7)\mathrm{d}{s}^2={\mathrm{e}}^{a(r)}\mathrm{d}{t}^2-{\mathrm{e}}^{b(r)}\mathrm{d}{r}^2-r^2(\mathrm{d}{\theta }^2+r^2{\sin }^2\theta \mathrm{d}{\varphi }^2).

For the metric above the Karmarkar condition simplifies to a differential relation between the metric functions a(r) and b(r) [58] in terms of the rank-4 Riemann tensor

(8){\displaystyle \frac{R_{1414}{R}_{2323}}{R_{1212}{R}_{3434}+R_{1224}{R}_{1334}}}=1,

leading to the differential condition on the metric components

(9){\displaystyle \frac{b^{\prime }(r)}{1-{\mathrm{e}}^{b(r)}}}={\displaystyle \frac{a^{\prime }(r)\left[a^{\prime }(r)+2/r\right]}{2}}.

One of the solutions for the a(r) and b(r) is given by

(10){\mathrm{e}}^a=A{\left(1+B{r}^2\right)}^4,{\mathrm{e}}^b=1+64B^{2}AC{r}^2{\left(1+B{r}^2\right)}^2,

which will be used to explore the properties of compact stars under the influence of the Karmarkar and Tolman condition [1, 67]. Here, A, B and C are the arbitrary constants to be determined by using some physical conditions.

In Karmarkar−Tolman spacetime, the metric functions {\mathrm{e}}^{a(r)} and {\mathrm{e}}^{b(r)} serve as relativistic analogs to gravitational potentials. This is a static, spherically symmetric line element commonly used in general relativity to describe stellar interiors. It forms the basis of many solutions to the Einstein field equations, including those incorporating embedding conditions such as the Karmarkar condition. This constraint ensures the spacetime can be embedded in 5D flat space, making the solution physically consistent for modeling relativistic stars. Both functions together replace the single Newtonian potential, providing a complete description of curved spacetime. To embed a 4D Riemannian space in a 5D flat space, the metric must satisfy the Karmarkar condition, which imposes a specific relationship between the metric functions a(r) and b(r). The functions a(r) and b(r) are metric potentials that depend only on the radial coordinate r, reflecting the symmetry and static nature of the configuration. Tolman developed several exact solutions to Einstein’s field equations under the assumption of this form of the metric, using different choices for a(r), b(r), pressure, and energy density. The Karmarkar condition is a geometric constraint that must be satisfied for a four-dimensional Riemannian manifold to be class one embedded in a five-dimensional flat spacetime.

Now, solving the field equations in Eq. (2), we get

(11)\begin{array}{rl}\rho =& \frac{1}{8}\{\frac{4{\mathrm{e}}^{-b(r)}}{r^2}[[2r{b}^{\prime }(r)+{\mathrm{e}}^{b(r)}(r^2\mathcal{R}(r)+2)-2]f^{\prime }(\mathcal{R}(r))+r[f^{″}(\mathcal{R}(r))((r{b}^{\prime }(r)-4){\mathcal{R}}^{\prime }(r)\\ & -2r{\mathcal{R}}^{″}(r))-2r{f}^{(3)}(\mathcal{R}(r)){\mathcal{R}}^{\prime }(r{)}^2]]-4f(\mathcal{R})-\frac{Q(r{)}^2}{\pi r^4}\},\end{array}

(12)\begin{array}{rl}{P}_r=& \frac{1}{8}\{4{\mathrm{e}}^{-b(r)}{r}^{-2}[f^{\prime }(\mathcal{R}(r))(2r{a}^{\prime }(r)-{\mathrm{e}}^{b(r)}(r^2\mathcal{R}(r)+2)+2)\\ & +r(r{a}^{\prime }(r)+4){\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))]+4f(\mathcal{R})+\frac{Q^2(r)}{\pi r^4}\},\end{array}

(13)\begin{array}{rl}{P}_t=& \frac{1}{8}\{-\frac{2{\mathrm{e}}^{-b(r)}}{r}[-2(2(r{f}^{(3)}(\mathcal{R}(r)){\mathcal{R}}^{\prime }(r{)}^2+r{\mathcal{R}}^{″}(r)f^{″}(\mathcal{R}(r))+{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r)))\\ & +r{a}^{″}(r)f^{\prime }(\mathcal{R}(r)))+a^{\prime }(r)((r{b}^{\prime }(r)-2)f^{\prime }(\mathcal{R}(r))-2r{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r)))-r{a}^{\prime }(r{)}^2{f}^{\prime }(\mathcal{R}(r))\\ & +2b^{\prime }(r)(r{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))+f^{\prime }(\mathcal{R}(r)))+2r{\mathrm{e}}^{b(r)}\mathcal{R}(r)f^{\prime }(\mathcal{R}(r))]+4f(\mathcal{R})-\frac{Q(r{)}^2}{\pi r^4}\},\end{array}

where the Ricci scalar \mathcal{R} reads

(14)\begin{array}{rl}& \\ & \mathcal{R}={\displaystyle \frac{{\mathrm{e}}^{-b}}{2r^2}}(4-4{\mathrm{e}}^b+r^2{a}^{\prime 2}-4r{b}^{\prime }+r{a}^{\prime }\left(4-r{b}^{\prime }\right)+2r^2{a}^{\mathrm{\prime }\mathrm{\prime }}).\end{array}

With the definitions of a and b specified in Eq. (10), Eqs. (11)−(13) can be recast to the following forms:

(15)\begin{array}{rl}& \\ \rho =& \frac{1}{8}[\frac{4}{r^2}{\left[64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1\right]}^{-1}\{f^{\prime }(\mathcal{R}(r))[\left(r^2\mathcal{R}(r)+2\right)(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1)\\ & +256A{B}^{2}C{r}^2\left(3B^2{r}^4+4B{r}^2+1\right)(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1{)}^{-1}-2]\\ & +r(f^{″}(\mathcal{R}(r))[4{\mathcal{R}}^{\prime }(r)\left(32A{B}^{2}C{r}^2\left(B^2{r}^4-1\right)-1\right)(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1{)}^{-1})\\ & -2r{\mathcal{R}}^{″}(r)]-2r{f}^{(3)}(\mathcal{R}(r)){\mathcal{R}}^{\prime }(r{)}^2)\}-4f(\mathcal{R})-\frac{Q(r{)}^2}{\pi r^4}],\end{array}

(16)\begin{array}{rl}{P}_r=& \frac{1}{8}[\frac{4}{r^2}\{\left[64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1\right]{\}}^{-1}\{f^{\prime }(\mathcal{R}(r))[-\left(r^2\mathcal{R}(r)+2\right)(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2\\ & +1)+16B{r}^2(B{r}^2+1{)}^{-1}+2]+r[8B{r}^2(B{r}^2+1{)}^{-1}+4]{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))\}\\ & +4f(\mathcal{R})+\frac{Q(r{)}^2}{\pi r^4}],\end{array}

(17)\begin{array}{rl}{P}_t=& \frac{1}{8}\{-\frac{4}{r}(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1{)}^{-1}[r\mathcal{R}(r)\left(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1\right)f^{\prime }(\mathcal{R}(r))\\ & +128A{B}^{2}Cr\left(3B^2{r}^4+4B{r}^2+1\right)\left(r{\mathcal{R}}^{\prime }(r)f^{″}(R(r))+f^{\prime }(R(r))\right)(64A{B}^{2}C{r}^2{\left(B{r}^2+1\right)}^2+1{)}^{-1}\\ & +8Br(B{r}^2+1{)}^{-1}(\left(128A{B}^{3}C{r}^4)\left(B{r}^2+1\right)-1\right)f^{\prime }(\mathcal{R}(r))(64A{B}^{2}C{r}^2){\left(B{r}^2+1\right)}^2+1{)}^{-1}\\ & -r{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r)))-32B^2{r}^3{f}^{\prime }(\mathcal{R}(r)){\left(B{r}^2+1\right)}^{-2}+8Br\left(B{r}^2-1\right)f^{\prime }(\mathcal{R}(r))(B{r}^2+1{)}^{-2}\\ & -2(r{f}^{(3)}(\mathcal{R}(r)){\mathcal{R}}^{\prime }(r{)}^2)-2(r{\mathcal{R}}^{″}(r)f^{″}(\mathcal{R}(r))+{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r)))]+4f(\mathcal{R})-\frac{Q(r{)}^2}{\pi r^4}\}.\end{array}

In order to obtain an explicit form of charge, the following EoS is assumed,

(18)P_r={\displaystyle \frac{1}{3}}(\rho -4B_g),

where B_g represents the bag constant, which is the fixed energy density per unit volume of space. It signifies the inward pressure required to confine quarks within a finite spatial region. Utilizing this EoS, we obtain the explicit function of charge as

(19)\begin{array}{rl}{Q}^2=& \pi r^2\{-{\mathrm{e}}^{-b(r)}[3r^2{a}^{\prime }(r){\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))+6r{a}^{\prime }(r)f^{\prime }(\mathcal{R}(r))-r^2{b}^{\prime }(r){\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))\\ & -2r{b}^{\prime }(r)f^{\prime }(\mathcal{R}(r))+8B_g{r}^2{\mathrm{e}}^{b(r)}-4r^2{\mathrm{e}}^{b(r)}\mathcal{R}(r)f^{\prime }(\mathcal{R}(r))-8{\mathrm{e}}^{b(r)}{f}^{\prime }(\mathcal{R}(r))+4r^2{\mathrm{e}}^{b(r)}f(\mathcal{R}(r))\\ & +2r^2{f}^{(3)}(\mathcal{R}(r)){\mathcal{R}}^{\prime }(r{)}^2+2r^2{\mathcal{R}}^{″}(r)f^{″}(\mathcal{R}(r))+16r{\mathcal{R}}^{\prime }(r)f^{″}(\mathcal{R}(r))+8f^{\prime }(\mathcal{R}(r))]\}.\end{array}

2.3 Matching with Reissner−Nordström exterior metric

The process of matching with the Reissner−Nordström (RN) exterior metric ensures a seamless transition between an interior solution of Einstein’s field equations, representing a charged gravitating body and the well known RN metric, which characterizes the exterior spacetime of a spherically symmetric, charged, non rotating mass. This matching is essential for developing physically realistic models of charged compact objects, such as charged stars or black holes. By enforcing this continuity, the process establishes constraints on the interior matter distribution, pressure profiles and charge distribution, thereby affecting the stability and physical feasibility of models for charged compact stars or other exotic objects within GR [68]. Here, we match the inner metric to the vacuum outer spherically symmetric metric defined by

(20)\begin{array}{rl}\mathrm{d}{s}^2=& \left(1-\frac{2M}{r}+\frac{Q^2}{r^2}\right)\mathrm{d}{t}^2-{\left(1-\frac{2M}{r}+\frac{Q^2}{r^2}\right)}^{-1}\mathrm{d}{r}^2\\ & -r^2\left(\mathrm{d}{\theta }^2+\sin {\theta }^2\mathrm{d}{\phi }^2\right),\end{array}

which leads to

(21)\begin{array}{rl}A& =\frac{{\left(-9MR+5Q^2+4R^2\right)}^4}{256R^2{\left[R(R-2M)+Q^2\right]}^3},\\ B& =\frac{MR-Q^2}{R^2\left[R(4R-9M)+5Q^2\right]},\\ C& =\frac{R^4\left(2MR-Q^2\right)}{4{\left(Q^2-MR\right)}^2}.\end{array}

In what follows, we are interested in three typical compact stars, Her X-1, SAXJ 1808.4-3658, and 4U 1820-30. With the formulae specified in Eq. (21), the values of the A, B and C constants can be readily calculated by using their estimated masses M and radii R [69–70], which are compiled in Table 1.

In the present manuscript, our primary goal is to develop analytically tractable and physically consistent charged anisotropic stellar configurations within f(\mathcal{R}) gravity, emphasizing the role of curvature corrections, anisotropy and electromagnetic fields on equilibrium and stability. As such the study is intended as a theoretical and exploratory framework, rather than a parameter fitting or statistical inference analysis. Nevertheless, we clarify that our models are not purely abstract and are already anchored to observationally inferred compact objects. Specifically, we employ measured mass radius estimates of well studied compact stars, Her X-1, SAX J1808.4-3658, and 4U 1820-30, to fix the integration constants via matching to the exterior Reissner−Nordström spacetime. This ensures that the interior solutions satisfy

{\mathrm{e}}^{\nu (R)}=1-{\displaystyle \frac{2M}{R}}+{\displaystyle \frac{Q^2}{R^2}},{\mathrm{e}}^{-\lambda (R)}=1-{\displaystyle \frac{2M}{R}}+{\displaystyle \frac{Q^2}{R^2}},

where M and R are taken directly from observationally inferred values reported in the literature. Furthermore, we compute key observable quantities that can be directly confronted with astrophysical measurements. A full numerical evaluation of tidal Love numbers and a Bayesian comparison with LIGO−Virgo−KAGRA data is beyond the scope of the present work and will be addressed in future studies. We acknowledge that no explicit statistical inference or error propagation is performed in this study. This is because our aim is not parameter estimation but rather to isolate the impact of modified gravity, anisotropy and charge on stellar structure using an analytically controlled framework. A full statistical analysis would require, a microphysical nuclear equation of state, a large parameter scan, Bayesian likelihood modeling using NICER and gravitational wave datasets. Such an analysis is planned as a natural extension of the present work. While the manuscript does not perform a full statistical fit, it quantitatively incorporates observed masses and radii, predicts observable surface redshifts and yields compactness values consistent with X-ray and gravitational wave constraints.

2.4 Three different f(\mathcal{R}) models

This section will explore three different models that can be used to study the properties of charged compact stars, including aspects like energy conditions, mass and other key physical characteristics. The three models are described as follows.

　● Model 1. First we will take the quadratic correction firstly proposed by Starobinsky [71],

(22)f(\mathcal{R})=\mathcal{R}+\alpha {\mathcal{R}}^2.

　 Here, \alpha is an arbitrary constant and \mathcal{R} represents the Ricci scalar. By setting \alpha =0, we can recover GR. The function f consists of a linear term and a quadratic term, with the parameter \alpha determining the relative contribution of the quadratic component compared to the linear one. The behavior of f varies with different values of \alpha. When \alpha is positive, the quadratic term dominates and grows more rapidly than the linear term for large \mathcal{R}. Conversely, if \alpha is negative, the quadratic term suppresses the overall growth of f.

　● Model 2. We take the following viable model, in which the correction term to GR has a form of exponential like [72]

(23)f(\mathcal{R})=\mathcal{R}+\alpha \mathcal{R}\left({\mathrm{e}}^{-\mathcal{R}/\gamma }-1\right).

　 Here, \alpha and \gamma are arbitrary constants. The function f combines a linear term \mathcal{R} with a modified exponential term. Namely, it has two components, the first is the linear term \mathcal{R}, and the second involves \alpha \mathcal{R} multiplied by the difference between the exponential decay {\mathrm{e}}^{-\mathcal{R}/\gamma } and 1. This form indicates that f depends on \mathcal{R} through both linear and decaying exponential behavior, with \alpha controlling the amplitude and \gamma determining the rate of exponential decay.

　● Model 3. Here, we apply some higher Ricci scalar terms to the quadratic gravity model, one of the cubic correction to GR [73], as like

(24)f(\mathcal{R})=\mathcal{R}+\alpha {\mathcal{R}}^2\left(1+\gamma \mathcal{R}\right).

　 Here, \alpha and \gamma are arbitrary constants, where f depends on the variable \mathcal{R}. The first term, \mathcal{R}, indicates a linear dependence only on \mathcal{R}, while the second term, \alpha {\mathcal{R}}^2(1+\gamma \mathcal{R}), represents a quadratic component that includes an additional factor involving \gamma, making the term effectively cubic in \mathcal{R}. The constant \alpha controls the magnitude of the quadratic contribution, while \gamma adds further complexity by enhancing the dependence on higher powers of \mathcal{R}. This type of function can describe systems where both linear and nonlinear effects influence the behavior of f as \mathcal{R} increases.

The comparison between Model 1, Model 2 and Model 3 is essential for understanding how different curvature corrections modify the internal physics of compact stars. The Model 1 introduces a simple quadratic correction, leading to moderate deviations from general relativity GR. The additional {\mathcal{R}}^2 term enhances the effective pressure and energy density at high curvature central regions, typically resulting in smoother and well behaved profiles. It serves as a baseline to assess how higher order terms influence stellar structure. In Model 2 the exponential correction introduces a scale dependent modification controlled by \gamma. At low curvature, the model closely mimics GR, while at high curvature the exponential suppression alters the growth of effective quantities. This often leads to softer profiles of \rho ,P_r and P_t improve the satisfaction of energy conditions over a wider radial range. It is particularly useful for testing viability against observational constraints. The Model 3 inclusion of higher order curvature terms {\mathcal{R}}^3 significantly enhances nonlinear effects. This model typically produces stronger deviations in the core region, affecting anisotropy and potentially allowing more compact or massive configurations. These models are constructed to pass solar system and weak field tests, remain stable under high curvature conditions and smoothly interpolate between cosmological and astrophysical scales. Our results demonstrate that curvature corrections from these models can significantly influence, compactness bounds, surface redshift and anisotropic force balance. This establishes a meaningful bridge between cosmological modified gravity and compact star physics.

3 Physical aspects of f(\mathcal{R}) gravity models

An interesting way to calculate f(\mathcal{R}) modifications characterized by two independent arbitrary constant parameters \alpha and \gamma. In our work, we consider the values of these parameters are \alpha = 0.03, and \gamma = 0.5 for the f(\mathcal{R}) gravity models [74]. The importance of these values was chosen to maintain the physical consistency of our charge compact star models and to comply with the experimental and observational limits. As for the small value of \alpha, which is equal to 0.03, it aligns with our solar system observations. These include phenomena such as light deflection, perihelion precession and radar echo delays, all of which impose constraints on any deviations from GR in weak gravitational fields. While the \gamma parameter such as \gamma = 0.5, gives the exponential correction in models [74]. The parameters are adjusted carefully to ensure stability, satisfy energy conditions and avoid ghost instabilities across all the models. By selecting values consistent with both cosmological observations and the solar system, our models remain consistent with the well tested predictions of general relativity, while successfully describing the behavior of charge compact stars within the framework of modified f(\mathcal{R}) gravity.

Phenomenological anisotropic fluid models have therefore been extensively employed in the literature to encode these complex microphysical effects in a tractable form (see, e.g., [8–9, 75]). The anisotropy may arise due to, relativistic nuclear interactions and many body correlations at extreme densities, superfluidity and superconductivity of neutrons and protons in the stellar core [21], phase transitions to hyperonic, meson condensed or quark matter phases [76], strong magnetic fields and rotation induced stresses [77], and electromagnetic effects in charged or partially charged stellar configurations. These mechanisms naturally lead to unequal radial and tangential stresses, making pressure anisotropy a physically expected feature of realistic neutron stars rather than a purely mathematical assumption. The higher order curvature terms in f(\mathcal{R}) gravity can be interpreted as effective geometric corrections that mimic unknown microphysical contributions at high density [77–78]. Here, we discuss some physical conditions which are necessary for the interior solution. In the following, we present the anisotropic behavior and stability conditions.

3.1 Energy density and pressure evolutions

The density plot for the strange star candidates, Her X-1, SAX J 1808.4-3658 (SS1) and 4U 1820-30, as shown in Fig. 1, illustrates that as r\to 0, the density \rho reaches its maximum value. In simple terms, \rho behaves as a decreasing function of r, meaning that as r increases, \rho decreases. This pattern signifies the high compactness of the stellar core, confirming that our f(\mathcal{R}) models under investigation remain valid in the outer regions of the core.

As it is well known, the radial pressure vanishes at the surface, i.e., P_r(R)=0. Using the radial pressure condition, we determine the value of B_g for the three relativistic spheres, as presented in Table 2. With Eqs. (15)−(17) and the involved constants A, B and C given in Table 1, we obtain \rho, P_r and P_t. Figures 2 and 3 illustrate the variation in anisotropic radial and transverse pressures, denoted as P_r and P_t, respectively.

The variation of radial derivative of density, radial and transverse pressure are shown in Figs. 4, 5 and 6, respectively. We see that , and for all three models and strange stars. For r=0, we get

(25){\displaystyle \frac{\mathrm{d}\rho }{\mathrm{d}r}}=0,{\displaystyle \frac{\mathrm{d}{P}_r}{\mathrm{d}r}}=0.

As expected, we have the maximum density for the minimum value of r (star core density \rho (0)={\rho }_c), because these are the decreasing function.

3.2 Energy conditions

The general expression for the energy conditions can be derived from the Raychaudhuri equation for expansion [79]. From these conditions, it can be inferred that gravity is inherently attractive, with a positive energy density that cannot exceed the speed of light. The energy conditions are defined as follows.

● Null Energy Condition (NEC):

(26)P_r+\rho \ge 0,P_t+\rho \ge 0.

● Weak Energy Condition (WEC):

(27)\rho \ge 0,\rho +P_r\ge 0,\rho +P_t\ge 0.

● Strong Energy Condition (SEC):

(28)\rho +P_r\ge 0,\rho +P_t\ge 0,\rho +P_r+2P_t\ge 0.

● Dominant Energy Condition (DEC):

(29)\rho \ge |P_r|,\rho \ge |P_t|.

In our study, the NEC (26) and SEC (28) have been thoroughly examined across all three considered f(\mathcal{R}) models, following the framework discussed in Ref. [80]. In fact, the evolution of all the above energy conditions for compact stars is thoroughly satisfied in our viable models, as illustrated graphically for all three strange star candidates in Figs. 7, 8, and 9. This confirms that the matter distribution remains physically acceptable and avoids exotic behavior such as negative energy densities.

3.3 TOV equation

The modeling of stellar structure in scenarios where general relativistic effects play a significant role, such as in neutron stars, relies on solving Einstein’s equation for a spherically symmetric star in static equilibrium. This framework describes the star in terms of energy density, pressure and a quantity analogous to the gravitational potential. The resulting equation can be seen as a natural extension of the principles governing stars in hydrostatic equilibrium [81].

In f(\mathcal{R}) gravity, the matching of an interior stellar solution to an exterior spacetime requires careful consideration of the junction conditions, which can be more restrictive than in General Relativity due to the higher order nature of the field equations. In addition to the continuity of the first fundamental form, namely the metric coefficients g_{tt} and g_{rr} across the stellar boundary r=R, the absence of surface layers generally demands the continuity of the effective radial pressure. Moreover, because the field equations involve derivatives of the Ricci scalar, regular matching typically requires the Ricci scalar and its normal derivative to be continuous at the boundary, i.e., {\mathcal{R}}_{r=R}=0 and {\mathrm{\nabla }}_{\mu }R{n}^{\mu }{]}_{r=R}=0. In the present model, the interior solution is constructed such that the curvature scalar smoothly approaches a constant or vanishes at the stellar surface, which implies that f_R remains continuous across the boundary.

The TOV equation can be expressed in a generalized form as

(30){\displaystyle \frac{\mathrm{d}{P}_r}{\mathrm{d}r}}+{\displaystyle \frac{a^{\prime }(\rho +P_r)}{2}}+{\displaystyle \frac{2(P_r-P_t)}{r}}+{\displaystyle \frac{Q}{4\pi r^4}}\frac{\mathrm{d}Q}{\mathrm{d}r}=0.

Furthermore, it can be written in form of gravitational, hydrostatic, anisotropic and electromagnetic forces (denoted by F_g, F_h, F_a and F_e in order)

(31)F_g+F_h+F_a+F_e=0,

which yields

(32)\begin{array}{rl}& F_g=-\frac{4Br(P_r+\rho )}{1+B{r}^2},F_h=-\frac{\mathrm{d}{P}_r}{\mathrm{d}r},\\ & F_a=\frac{2(P_r-P_t)}{r},F_e=\frac{Q}{4\pi r^4}\frac{\mathrm{d}Q}{\mathrm{d}r}.\end{array}

Using these definitions, we generate plots for three strange compact stars as illustrated in Fig. 10.

Figure 10 shows the variation of the four forces with respect to the radial coordinate r (km). More specifically, these plots illustrate the variation of different force components contributing to the hydrostatic equilibrium of a charge compact stellar object under three different gravitational models. As can be seen from Fig. 10, the gravitational force acts inward throughout the star. It is the dominant attractive force and decreases in magnitude with increasing radius. The anisotropic force is directed inward, indicating negative anisotropy (radial pressure greater than tangential). It slightly modifies the equilibrium but remains smaller in magnitude than gravity. The hydrostatic force is positive and outward directed. It represents the pressure gradient that resists gravitational collapse. It gradually decreases toward the surface. Finally, the electric force is strong near the center due to high charge density, then quickly diminishes. It provides an outward repulsive effect and contributes to the stability of the star. Our results confirm that the combined effects of hydrostatic, electric and anisotropic forces balance the gravitational pull, maintaining equilibrium throughout the stellar interior. This ensure the configuration is stable and physically acceptable under the specified three different models.

3.4 Stability analysis

In this section, we examine the stability of compact star models within the framework of f(\mathcal{R}) theory. To assess the stability of our model, we compute the radial velocity (v_{sr}) and transverse velocity (v_{st}), which are defined via

(33){\displaystyle \frac{\mathrm{d}{P}_r}{\mathrm{d}\rho }}=v_{sr}^2,{\displaystyle \frac{\mathrm{d}{P}_t}{\mathrm{d}\rho }}=v_{st}^2.

For stability, both the radial and transverse speeds must satisfy the condition, 0\le v_{sr}^2\le 1, i.e., the value of v_{sr}^2 is nearly equal to \frac{1}{3} and 0\le v_{st}^2\le 1.

As shown in Fig. 11, the evolution of the radial and transverse sound speeds for all three types of strange star candidates remains within the stability limits, as previously discussed. The speed of sound v_s in dense matter must remain below the speed of light c. Some extreme EoS models require v_s\approx c, which is physically implausible. Observational data indicate a distinct mass gap between the heaviest neutron stars (2.3M_{\odot }) and the lightest black holes (5M_{\odot }). When a neutron star exceeds its maximum mass, such as through mass accretion in a binary system, it undergoes rapid gravitational collapse within milliseconds, forming a black hole. This phenomenon explains the observed absence of compact objects with masses between 2.5-5M_{\odot }. Core collapse supernovae typically yield either neutron stars (\lesssim 2.5M_{\odot }) or black holes (\gtrsim 5M_{\odot }), with little space for intermediate objects. In a similar way, it can be observed from Fig. 12 that , so the stability is attained for compact stars in our considered f(\mathcal{R}) models.

3.5 EoS parameter

The EoS parameter w is essential for understanding the thermodynamic properties of cosmic fluids and their impact on the universe large scale structure and its evolution. In extended gravitational models beyond GR, such as f(\mathcal{R}) gravity, scalar-tensor theories and other modified gravity frameworks, the EoS parameter can vary over time, providing a dynamic view of cosmic acceleration. These models often introduce additional scalar fields that influence w, resulting in its evolution and potential transitions between different values throughout cosmic history. The EoS parameter is an essential tool for describing both the universe expansion and the properties of cosmic fluids. Whether constant or time dependent, its behavior offers crucial understanding of dark energy, the early universe and the long term evolution of the cosmos. By linking theoretical models with observational data, it plays a central role in modern cosmology and astrophysics.

For the anisotropic case, the EoS can be expressed as

(34)P_r=w_r\rho ,P_t=w_t\rho ,

where w_{r,t} are relevant radial and transverse EoS parameters, respectively. The limits are such that and . The graphical representation of the behavior of w_r and w_t is shown Figs. 13 and 14, where and are visible, indicating that the matter content is typical real matter. The comparative analysis across these models is essential for several reasons. It isolates how different functional forms of f(\mathcal{R}) influence equilibrium and stability. It helps identify parameter ranges that yield physically viable stars. It also clarifies whether higher order curvature terms improve or worsen agreement with observational constraints such as mass radius relations. Finally, it provides insight into the interplay between modified gravity, anisotropy and charge which cannot be captured within general relativity alone.

3.6 The measurement of anisotropy

Anisotropy in compact stars, like neutron stars and quark stars, refers to the differences in pressure that occur along various directions within the star. This deviation from isotropic pressure distribution is significant in understanding the internal structure, stability and EoS of these dense objects. Various physical mechanisms, such as strong magnetic fields, phase transitions and anisotropic stress from exotic matter, can induce anisotropy in compact stars. The measurement of anisotropy in compact stars is typically inferred through astrophysical observations and theoretical modeling. Additionally, theoretical models employing GR, the TOV equation and modified gravity frameworks help estimate the impact of anisotropy on mass radius relations and maximum mass constraints. In realistic neutron star interiors, matter is expected to depart from perfect isotropy, particularly at supranuclear densities where extreme physical conditions prevail. In such regimes, p_r=p_t becomes overly restrictive. Naturally arises and plays a crucial role in determining the equilibrium, stability and observable properties of compact stars. Several well established physical mechanisms contribute to anisotropy in neutron star matter. At high densities, strong nuclear interactions and many body effects lead to direction dependent stresses, especially when the mean free path of particles becomes comparable to microscopic interaction scales. Additionally, relativistic motion of particles and superfluidity in the stellar core can induce anisotropic pressure components. The presence of intense magnetic fields, particularly in magnetars, further enhances anisotropy by introducing Lorentz forces that act differently along radial and tangential directions. Coulomb interactions and charge separation effects, especially in charged or partially charged configurations, also contribute to pressure anisotropy.

The anisotropy is defined as

(35)\mathrm{\Delta }={\displaystyle \frac{2}{r}}(P_t-P_r).

Plots of the anisotropy are displayed in Fig. 15. We observe that \mathrm{\Delta }>0, such as P_t>P_r, indicating that the anisotropy is oriented outward.

3.7 Electric field and charge

The equilibrium configuration of a charged compact star is shaped by the interplay between gravitational attraction and electrostatic repulsion. The charge distribution q(r), determined using Gauss’s law in curved spacetime, is depicted in Fig. 16.

The inclusion of electric charge in compact stellar models is motivated by the extreme physical conditions prevailing at supranuclear densities, where additional interactions beyond pure gravity may become significant. Moreover, exotic phases of matter, such as quark matter or superconducting nuclear matter, may naturally support non uniform charge distributions within the stellar interior. The relatively large charge values considered in the model should therefore be interpreted as effective or phenomenological quantities that primarily reside inside the star and are expected to be largely screened by surrounding plasma, rather than as macroscopic charges directly observable at large distances. From a physical standpoint, the magnitude of charge is carefully constrained so that electromagnetic repulsion remains subdominant to gravitational attraction, which is ensured by maintaining a sufficiently small charge to mass ratio Q/M\ll 1. This condition guarantees that the star remains gravitationally bound and dynamically stable. Consequently, the inclusion of electric charge in the present work provides a realistic and physically motivated mechanism that enhances stability and equilibrium without violating known observational or theoretical bounds.

Additionally, understanding the behavior of \sigma (r) is essential for determining whether charge is concentrated near the core or extends towards the outer layers of the star. The charge density \sigma (r), which describes the spatial distribution of charge within the stellar interior, is illustrated in Fig. 17. The charge density \sigma (r), shown in Fig. 17, is obtained from the derivative of q(r) and reveals the manner in which charge is localized within the star. These profiles provide valuable insights into how charge accumulates and varies across different regions of the star.

The presence of charge introduces an electric field that significantly alters the star structure. A strong electric field generates an outward force that counteracts gravitational collapse, thereby contributing to the equilibrium of the system. The electric field E^2(r), derived from the charge function, plays a crucial role in modifying the internal pressure gradients and influencing the overall stability of the configuration. Figure 18 presents the variation of the squared electric field E^2(r) as a function of radial distance. The behavior of E^2(r) provides crucial information on how charge contributes to stabilizing the stellar configuration against gravitational collapse. In Fig. 18, the radial dependence and peak values of E^2(r) highlight its impact on the star structure.

3.8 The comparison between previous work and the present study

In the present work, our primary objective is to construct self consistent relativistic models of charged anisotropic compact stars within the framework of f(\mathcal{R}) gravity, focusing on how modifications to gravity, anisotropy and electromagnetic fields collectively influence the global structure, stability and equilibrium of ultra compact objects. As the analysis is carried out at a phenomenological and effective level, which is standard in relativistic stellar modeling when the exact microphysical equation of state at supranuclear densities remains uncertain. Regarding the origin of anisotropy, it is well established in the literature that pressure anisotropy naturally arises in neutron stars due to several physical mechanisms operating at extreme densities. Importantly, our analysis demonstrates how modified curvature terms in f(\mathcal{R}) gravity, combined with anisotropy and electromagnetic contributions, can alter pressure gradients, enhance stability and modify equilibrium conditions relative to general relativity. These effects provide valuable insight into how deviations from GR may influence the macroscopic properties of neutron stars, independently of the specific microscopic equation of state.

4 Summary

In this work, we explore stable stellar models with an anisotropic matter distribution within the framework of f(\mathcal{R}) gravity, considering a static and spherically symmetric configuration. The stellar metric system is formulated using the Krori−Barua solution, where junction conditions on a hypersurface determine the unknown constants. The study utilized observational data from three different star models, incorporating parameters such as mass and radius. The primary goal is to examine the impact of additional curvature terms in the f(\mathcal{R}) theory on the stability of charged compact stars. To perform this analysis, three specific f(\mathcal{R}) models were considered, which explored the relationship between mass and radius in neutron stars with quadratic and exponential Ricci scalar corrections. Our findings also indicate that the incorporation of an EoS leads to a significant increase in the maximum mass of neutron stars, supporting the plausibility of f(\mathcal{R}) gravity. The radial coordinate r influences physical properties such as energy density and anisotropic stresses. We show that the energy density decreases as the star radius increases, suggesting a maximally compact configuration at the core. A similar trend is observed in the behavior of both tangential and radial pressures. This observation suggests the formation of dense celestial structures. The first derivative of the pressure components transitions from negative to positive, indicating the presence of a local minimum. This highlights the compact nature of the object, a similar pattern observed for \frac{\mathrm{d}{P}_t}{\mathrm{d}r}. The locally anisotropic sphere is structured to satisfy all energy conditions, at a specific radius r. In the TOV equation, the forces like F_g, F_h, F_a and F_e suggest that all charged compact stars will eventually reach a state of hydrodynamic equilibrium. The relation between {\omega }_{r,t} and r within a specified range, the compact internal core of stars in f(\mathcal{R}) gravity is demonstrated by the fact that both {\omega }_r and {\omega }_t remain between 0 and 1. This study presents a systematic and physically self consistent investigation of charged anisotropic compact stellar configurations within the framework of f(\mathcal{R}) gravity, emphasizing the nontrivial interplay among spacetime geometry, matter anisotropy and electromagnetic interactions. Future extensions of this work may include the incorporation of realistic nuclear equations of state, the effects of rotation and strong magnetic fields, and the application of statistical or Bayesian methods to confront the models with astrophysical observations from X-ray missions and gravitational wave detections. Such studies will further clarify the viability of modified gravity theories in the context of compact stars. We hope that the present results will stimulate further research into the role of modified gravity in extreme astrophysical environments.

References

[1] R. C. Tolman , Static solutions of Einstein’s field equations for spheres of fluid, Phys. Rev. 55(4), 364 (1939)

[2] M. Z. Bhatti , B. Almutairi , Z. Yousaf , S. Ijaz , and A. S. Khan , Analytical solutions of planar geometry with zero complexity in extended gravity, Int. J. Geom. Methods Mod. Phys. 22(6), 2550001 (2025)

[3] Z. Yousaf , H. Asad , M. Alshehri , T. Suzuki , and M. Z. Bhatti , Electromagnetic field and spatially hyperbolical spacetime models, Int. J. Geom. Methods Mod. Phys 23(02), 2550145 (2026)

[4] J. R. Oppenheimer and G. M. Volkoff , On massive neutron cores, Phys. Rev. 55(4), 374 (1939)

[5] Z. Yousaf , A. Adeel , M. Rizwan , G. Mustafa , and A. Ali , Construction of complexity-free fuzzy dark matter wormholes, Int. J. Geom. Methods Mod. Phys. 22(12), 2550093 (2025)

[6] W. Baade and F. Zwicky , Remarks on super-novae and cosmic rays, Phys. Rev. 46(1), 76 (1934)

[7] Z. Yousaf , K. Bamba , M. Z. Bhatti , and U. Farwa , Quasi-static evolution of axially and reflection symmetric large-scale configuration, Int. J. Geom. Methods Mod. Phys. 21(8), 2430005 (2024)

[8] R. L. Bowers and E. P. T. Liang , Anisotropic Spheres in General Relativity, Astrophys. J. 188, 657 (1974)

[9] M. Ruderman , Pulsars: Structure and dynamics, Annu. Rev. Astron. Astrophys. 10(1), 427 (1972)

[10] A. R. Athar , M. Ilyas , and B. Masud , Anisotropic strange compact stars in Krori−Barua spacetime under f(R;G) gravity, Int. J. Geom. Methods Mod. Phys. 20(1), 2350003 (2023)

[11] M. Z. Bhatti , M. Y. Khlopov , M. Shrahili , A. Rehman , and Z. Yousaf , Complexity of charged self-gravitating system in f(R;G) theory, Int. J. Geom. Methods Mod. Phys. 22(10), 2540014 (2025)

[12] A. Qadir , H. W. Lee , and K. Y. Kim , Modified relativistic dynamics, Int. J. Mod. Phys. D 26(5), 1741001 (2017)

[13] S. Nojiri and S. D. Odintsov , Modified Gauss−Bonnet theory as gravitational alternative for dark energy, Phys. Lett. B 631(1−2), 1 (2005)

[14] S. Capozziello , S. Nojiri , and S. D. Odintsov , Unified phantom cosmology: Inflation, dark energy and dark matter under the same standard, Phys. Lett. B 632(5−6), 597 (2006)

[15] S. Capozziello and M. De Laurentis , Extended theories of gravity, Phys. Rep. 509, 167 (2011)

[16] S. D. Odintsov,D. Saez-Gomez, f(R; T; R_μνT^μν) gravity phenomenology and ΛCDM universe, Phys. Lett. B 725(4−5), 437 (2013)

[17] Z. Yousaf , K. Bamba , and M. Z. H. Bhatti , Causes of irregular energy density in f(R, T) gravity, Phys. Rev. D 93(12), 124048 (2016)

[18] Z. Yousaf,M. Z. Bhatti,S. Khan,P. K. Sahoo, f(G; T_αβT^αβ) theory and complex cosmological structures, Phys. Dark Univ. 36, 101015 (2022), arXiv:

[19] M. X. Xu , T. Harko , and S. D. Liang , Quantum cosmology of f(R; T) gravity, Eur. Phys. J. C 76(8), 449 (2016)

[20] Z. Yousaf , Viscous fluid cosmology: A path to cosmic acceleration, Phys. Dark Univ. 48, 101884 (2025)

[21] L. Herrera and N. O. Santos , Local anisotropy in self-gravitating systems, Phys. Rep. 286(2), 53 (1997)

[22] M. Sharif and Z. Yousaf , Instability of a dissipative restricted non-static axial collapse with shear viscosity in f(R) gravity, J. Cosmol. Astropart. Phys. 06, 019 (2014)

[23] Z. Yousaf , M. Ilyas , and M. Z. Bhatti , Influence of modification of gravity on spherical wormhole models, Mod. Phys. Lett. A 32(30), 1750163 (2017)

[24] G. G. L. Nashed , S. D. Odintsov , and V. K. Oikonomou , Anisotropic compact stars in higher-order curvature theory, Eur. Phys. J. C 81(6), 528 (2021)

[25] G. G. L. Nashed and S. Capozziello , Anisotropic compact stars in f(R) gravity, Eur. Phys. J. C 81(5), 481 (2021)

[26] G. G. L. Nashed,S. D. Odintsov,V. K. Oikonomou, Anisotropic compact stars in D→4 limit of Gauss−Bonnet gravity, Symmetry (Basel) 14(3), 545 (2022)

[27] D. Müller , V. C. de Andrade , C. Maia , M. J. Reboucas , and A. F. F. Teixeira , Future dynamics in f(R) theories, Eur. Phys. J. C 75(1), 13 (2015)

[28] S. Karmakar , S. Mukherjee , R. Sharma , and S. D. Maharaj , The Role of pressure anisotropy on the maximum mass of cold compact stars, Pramana 68(6), 881 (2007)

[29] S. Yazadjiev , Relativistic models of magnetars: Nonperturbative analytical approach, Phys. Rev. D 85(4), 044030 (2012)

[30] A. V. Astashenok , S. Capozziello , and S. D. Odintsov , Magnetic Neutron Stars in f(R) gravity, Astrophys. Space Sci. 355(2), 333 (2015)

[31] M. Ilyas , Compact stars in f(R; G; T) gravity, Int. J. Mod. Phys. A 36(24), 2150165 (2021)

[32] A. V. Astashenok , S. Capozziello , S. D. Odintsov , and V. K. Oikonomou , Maximum baryon masses for static neutron stars in f(R) gravity, Europhys. Lett. 136(5), 59001 (2021)

[33] M. Z. Bhatti , Z. Yousaf , and M. Ilyas , Existence of wormhole solutions and energy conditions in f(R; T) gravity, J. Astrophys. Astron. 39(6), 69 (2018)

[34] M. Sharif and Z. Yousaf , Instability of meridional axial system in f(R) gravity, Eur. Phys. J. C 75(5), 194 (2015)

[35] M. Sharif and Z. Yousaf , Radiating cylindrical gravitational collapse with structure scalars in f(R) gravity, Astrophys. Space Sci. 357(1), 49 (2015)

[36] M. Sharif and Z. Yousaf , Energy density inhomogeneities with polynomial f(R) cosmology, Astrophys. Space Sci. 352(1), 321 (2014)

[37] M. Sharif and Z. Yousaf , Dynamics of spherical stars with structure scalars and R + ϵRⁿ cosmology, Can. J. Phys. 93(8), 905 (2015)

[38] S. Perlmutter and B. P. Schmidt , Measuring cosmology with supernovae, Lect. Notes Phys. 598, 195 (2003)

[39] C. L. Bennett , M. Bay , M. Halpern , G. Hinshaw , C. Jackson , N. Jarosik , A. Kogut , M. Limon , S. S. Meyer , L. Page , and D. N. Spergel , The microwave anisotropy probe mission, Astrophys J. 583(1), 1 (2003)

[40] E. J. Copeland , M. Sami , and S. Tsujikawa , Dynamics of dark energy, Int. J. Mod. Phys. D 15(11), 1753 (2006)

[41] K. Bamba , S. Capozziello , S. Nojiri , and S. D. Odintsov , Dark energy cosmology: the equivalent description via different theoretical models and cosmography tests, Astrophys. Space Sci. 342(1), 155 (2012)

[42] H. T. Cromartie , E. Fonseca , S. M. Ransom , P. B. Demorest , Z. Arzoumanian , H. Blumer , P. R. Brook , M. E. DeCesar , T. Dolch , J. A. Ellis , and R. D. Ferdman , Relativistic Shapiro delay measurements of an extremely massive millisecond pulsar, Nature Astronomy 4(1), 72 (2020)

[43] T. Naseer , Role of non-conserved gravity theory and electric charge in constructing complexity-free stellar models: A novel approach under non-minimal coupling, Chin. J. Phys. 96, 1212 (2025)

[44] E. Demir , T. Naseer , A. Ashraf , and E. Güdekli , Investigating physical existence of charged stellar models, Chin. J. Phys. 91, 299 (2024)

[45] R. W. Romani,D. Kandel,A. V. Filippenko,T. G. Brink,W. Zheng, PSR J0952-0607: The fastest and heaviest known galactic Neutron Star, Astrophys. J. Lett. 934(2), L17 (2022)

[46] B. P. Abbott,R. Abbott,T. D. Abbott,F. Acernese,K. Ackley,C. Adams,T. Adams,P. Addesso,R. X. Adhikari,V. B. Adya,C. Affeldt,M. Afrough,B. Agarwal,M. Agathos,K. Agatsuma,N. Aggarwal,O. D. Aguiar,L. Aiello,A. Ain,P. Ajith,B. Allen,G. Allen,A. Allocca,P. A. Altin,A. Amato, GW170817: Observation of gravitational waves from a binary Neutron Star inspiral, Phys. Rev. Lett. 119(16), 161101 (2017)

[47] Y. Z. Fan , M. Z. Han , J. L. Jiang , D. S. Shao , and S. P. Tang , Maximum gravitational mass MTOV=2.25−0.07+0.08M_ʘ inferred at about 3% precision with multimessenger data of neutron stars, Phys. Rev. D 109(4), 043052 (2024)

[48] Y. Feng,T. Naseer,A. Ashraf,D. Sofuoglu,I. Abdullayeva, A Brief analysis of isotropic Karmarkar models in modified gravity theory, Chin. J. Phys. 90, 372 (2024)

[49] M. Sharif and T. Naseer , Influence of charge on anisotropic class-one solution in non-minimally coupled gravity, Phys. Scr. 97(12), 125016 (2022)

[50] T. Naseer and J. L. Said , Existence of non-singular stellar solutions within the context of electromagnetic field: A comparison between minimal and non-minimal gravity models, Eur. Phys. J. C 84(8), 808 (2024)

[51] M. Sharif and T. Naseer , Effects of non-minimal matter-geometry coupling on embedding class-one anisotropic solutions, Phys. Scr. 97(5), 055004 (2022)

[52] M. Sharif and T. Naseer , Influence of charge on extended decoupled anisotropic solutions in f(R; T; R_λξT^λξ) gravity, Indian J. Phys. 96(14), 4373 (2022)

[53] Y. Feng,T. Naseer,A. Ashraf,D. Sofuoğlu,I. Abdullayeva, A brief analysis of isotropic Karmarkar models in modified gravity theory, Chin. J. Phys. 90, 372 (2024)

[54] T. Naz,A. Usman,M. F. Shamir, Embedded class-i solution of compact stars in f(R) gravity with Karmarkar condition, arXiv: 2021)

[55] Z. Asghar,M. F. Shamir,A. Usman,A. Malik, Study of embedded class-i fluid spheres in f(R, T) gravity with Karmarkar condition, arXiv: 2023)

[56] Z. Yousaf , M. Sharif , M. Ilyas , and M. Z. Bhatti , Influence of f(R) models on the existence of anisotropic self-gravitating systems, Eur. Phys. J. C 77(10), 691 (2017)

[57] W. U. Rahman and M. Ilyas , Exploring compact anisotropic stars in f(R) gravity within Karmarkar−Tolman space-time framework, Int. J. Geom. Methods Mod. Phys. 23(02), 2550149 (2026)

[58] G. Mustafa , M. F. Shamir , and T. C. Xia , Physically viable solutions of anisotropic spheres in f(R, G) gravity satisfying the Karmarkar condition, Phys. Rev. D 101(10), 104013 (2020)

[59] M. Zubair , H. Javaid , H. Azmat , and E. Gudekli , Relativistic stellar model in f(R, T) gravity with Karmarkar’s condition, New Astron. 88, 101610 (2021)

[60] T. Naseer , M. Sharif , H. Afzal , M. R. Eid , and R. Altuijri , Imprints of electric charge and Rastall’s non-conservation phenomenon on dynamical properties of anisotropic compact stellar models, Phys. Dark Univ. 49, 102024 (2025)

[61] R. Kumar , S. K. Maurya , A. Errehymy , K. Myrzakulov , Z. Umbetova , and V. N. Pathak , Anisotropic model of stellar objects in modified f(R) gravity, Chin. J. Phys 92, 166 (2024)

[62] A. Malik , T. Naz , U. Marwa , P. Bhar , A. Ali , and Z. Yousaf , Discussion of singularity-free embedding stellar structures in f(R) gravity utilizing scalar potential, Int. J. Geom. Methods Mod. Phys. 21(11), 2450186 (2024)

[63] T. Naz and M. F. Shamir , Study of charged stellar models in f(G) gravity with Tolman−Kuchowicz space-time, Int. J. Mod. Phys. A 35(9), 2050040 (2020)

[64] A. Errehymy , Y. Khedif , G. Mustafa , and M. Daoud , Anisotropic stars of class one space-time in f(R, T) gravity under the simplest linear functional of the matter−geometry coupling, Chin. J. Phys. 77, 1502 (2022)

[65] M. Sharif and K. Hassan , Study of charged celestial objects in modified gravity, Chin. J. Phys. 86, 227 (2023)

[66] U. Kasper , Wheeler−De Witt equations for fourth order quantum cosmology, Class. Quantum Gravity 10(5), 869 (1993)

[67] K. Karmarkar , Gravitational metrics of spherical symmetry and class one, Proc. Indian Acad. Sci. Sect. A 27(1), 56 (1948)

[68] G. F. T. Castillo and R. Cartas-Fuentevilla , Scattering by a Reissner−Nordstrom black hole, Phys. Rev. D 54(8), 4886 (1996)

[69] T. Gangopadhyay , S. Ray , X. D. Li , J. Dey , and M. Dey , Strange star equation of state fits the refined mass measurement of 12 pulsars and predicts their radii, Mon. Not. R. Astron. Soc. 431(4), 3216 (2013)

[70] T. Güver , F. Ozel , A. Cabrera-Lavers , and P. Wroblewski , The distance, mass, and radius of the Neutron Star in 4U 1608-52, Astrophys. J. 712(2), 964 (2010)

[71] A. A. Starobinsky, A new type of isotropic cosmological models without singularity, Phys. Lett. B 91(1), 99 (1980)

[72] K. Bamba , S. Nojiri , and S. D. Odintsov , The Universe future in modified gravity theories: Approaching the finite-time future singularity, J. Cosmol. Astropart. Phys. 10, 045 (2008)

[73] A. V. Astashenok , S. Capozziello , and S. D. Odintsov , Nonperturbative models of quark stars in f(R) gravity, Phys. Lett. B 742, 160 (2015)

[74] M. Ilyas , N. A. Shah , F. Khan , and R. Habib , Dynamics of compact objects in higher-order curvature gravity using Finch−Skea spacetime, Eur. Phys. J. C 85(3), 253 (2025)

[75] R. F. Sawyer , Condensed π-phase in neutron star matter, Phys. Rev. Lett. 29(6), 382 (1972)

[76] U. Lombardo,H. J. Schulze, Superfluidity in neutron star matter, Phys. Neutron Star Interiors. 48, 30 (2002), arXiv:

[77] V. Canuto and H. Y. Chiu , Quantum theory of an electron gas in intense magnetic fields, Phys. Rev. 173(5), 1210 (1968)

[78] K. Dev and M. Gleiser , Anisotropic stars: Exact solutions, Gen. Rel. Grav. 34(11), 1793 (2002)

[79] S. W. Hawking,G. F. Ellis, The Large Scale Structure of Space-Time, Cambridge University Press, 2023

[80] W. U. Rahman , M. Ilyas , Z. Yousaf , S. Ullah , F. Khan , and R. Khan , Compact relativistic geometries in f(R, G) gravity, Int. J. Geom. Methods Mod. Phys. 19(8), 2250126 (2022)

[81] A. Smith, Tolman−Oppenheimer-Volkoff (Tov) Stars, University of Texas, 2012

RIGHTS & PERMISSIONS

Higher Education Press