# Quadratic Rastall Gravity: From Low-mass HESS J1731$-$347 to High-mass PSR
J0952$-$0607 Pulsars

Waleed El Hanafy  Centre for Theoretical Physics, The British University in
Egypt, P.O. Box 43, El Sherouk City, Cairo 11837, Egypt
<EMAIL_ADDRESS>
###### Abstract

Similar to Rastall gravity we introduce matter-geometry nonminimal coupling
which is proportional to the gradient of quadratic curvature invariants. Those
are mimicking the conformal trace anomaly when backreaction of the quantum
fields to a curved spacetime geometry is considered. We consider a static
spherically symmetric stellar structure with anisotropic fluid and Krori-Barua
metric potentials model to examine the theory. Confronting the model with
NICER+XMM-Newton observational constraints on the pulsar PSR J0740$+$6620
quantifies the amount of the nonminimal coupling via a dimensionless parameter
$\epsilon\simeq-0.01$. We verify that the conformal symmetry is broken
everywhere inside the pulsar as the trace anomaly $\Delta>0$, or equivalently
the trace of the stress-energy tensor $\mathfrak{T}<0$, whereas the adiabatic
sound speed does not violate the conjecture conformal upper limit
$v_{r}^{2}/c^{2}=1/3$. The maximum compactness accordingly is
$C_{\text{max}}=0.752$ which is $4\%$ higher than GR. Notably, if the
conformal sound speed constraint is hold, observational data excludes
$\epsilon\geq 0$ up to $\geq 1.6\sigma$. The stellar model is consistent with
the self-bound structure with soft linear equation of state. Investigating
possible connection with MIT bag model of strange quarks sets physical bounds
from microscopic physics which confirm the negative value of the parameter
$\epsilon$. We estimate a radius $R=13.21\pm 0.96$ km of the most massive
observed compact star PSR J0952$-$0607 with $M=2.35\pm 0.17M_{\odot}$.
Finally, we show that the corresponding mass-radius diagram fits well lowest-
mass pulsar HESS J1731$-$347 and highest-mass pulsar PSR J0952$-$0607 ever
observed as well as the intermediate mass range as obtained by NICER and
LIGO/Virgo observations.

pulsars – strange quark star – nonminimal coupling – conformal anomaly –
modified gravity

##  I Introduction

The principle of minimal coupling has been assumed by the general relativity
(GR) theory which leads to a divergence-free energy-momentum constraint,
$\nabla{{}_{\alpha}}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=0$, in a curved
spacetime. However, a more general framework has been argued by Rastall by
giving up the minimal coupling principle and considering matter-geometry
nonminimal coupling instead [1, 2]. In Rastall’s gravity, the nonminimal
coupling is assumed to be of linear order of Ricci scalar $\mathfrak{R}$, i.e.
$\nabla{{}_{\alpha}}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}\propto\partial_{\beta}\mathfrak{R}$.
Clearly, the theory becomes distinguishable from GR where high curvature
associated with dense matter shows up, while it reduces to GR in vacuum. In
this sense, compact stars such as in the case of neutron stars (NSs) are ideal
regimes to test Rastall gravity, which has been shown to have novel features
when applied to compact stellar objects [3, 4]. On the contrary, it has been
argued that Rastall gravity is completely equivalent to GR [5]. However, this
argument has been criticized and shown to be wrong [6] due to
misinterpretation of the matter stress-energy tensor. In practice, the
inequivalence of the two theories has been ensured in several studies on
thermodynamical, cosmological and astrophysical levels, c.f. [7, 8, 9, 10, 11,
12, 13, 14, 15, 16, 17].

Trace anomaly, on the other hand, is important when the backreaction of the
quantum fields to a curved spacetime geometry is considered. It is well known
that the quantum corrections due to trace anomaly are geometric in nature and
can be expressed in terms of Gauss-Bonnet four dimensional topological
invariant (Type A) and the square of the Weyl tensor (Type B) [18]. Those are,
in general, of quadratic order of curvature invariants. At high densities,
$n\sim 40n_{\text{nuc}}$ where $n_{\text{nuc}}\approx 0.16$ fm-3 is the baryon
density number at nuclear saturation, perturbative quantum chromodynamics
(pQCD) is applied, conformal symmetry is hold, i.e. the trace anomaly
vanishes. At finite densities and/or temperatures compatible with NS
densities, few times nuclear saturation density $\rho_{\text{nuc}}\sim
2.7\times 10^{14}$ g/cm3, it has been argued that the conformal symmetry would
be broken, i.e. trace anomaly does not vanish, due to running of the strong
coupling with energy at quantum level [19]. Thermodynamic stability and
causality indeed set a constraining physical range to the trace anomaly,
$\Delta\in[1/3,-2/3]$, but they cannot determine its sign inside an NS.

It is believed that quantum phenomena should play an essential role in
gravitational collapse of massive objects. On the classical level, the
breaking of minimal coupling provides an approach to this problem [20]. In the
present work, we extend the nonminimal coupling frame to include exactly
quadratic curvature invariants, which we call quadratic Rastall gravity (QRG).
Since quadratic nonminimal coupling is motivated by quantum phenomena, we
argue the possibility of combining these two approaches in a single frame by
considering theories such as QRG. We still have the privilege of recovering GR
vacuum solutions, whereas quadratic contribution is more compatible with trace
anomaly when the backreaction of the quantum fields to a curved spacetime
geometry is considered. For sure, different set up is needed if one considers
exactly trace anomaly. However, the present work may provide a step towards
the exact problem.

For the classical matter sector, it is unlikely for a high dense matter as in
the NS core to be represented by the ideally isotropic fluid. Even in low
density regimes on the Newtonian limit, local anisotropy in matter fields
could arise due to anisotropic velocity distribution. At high density regimes,
local anisotropy naturally appears due to solidification [21], pion
condensation [22], super-fluidity [23], strong magnetic fields [24], hyperons
[25] and strong interactions [26]. Therefore, we assume the fluid to be
anisotropic which is the more realistic case in the compact stars. Local
anisotropy and its consequences on the stability and the gravitational
collapse of self-gravitating systems has been investigated by Herrera and
Santos [27] (see also [28, 29, 30, 31]). Several studies assume
gravitationally bound (hadrons) or self bound (quark) stars according to the
matter assumed. In the present study, we do not assume a specific matter type
in prior, since no particular equation of state (EoS) is assumed. Instead, we
assume the metric potentials inside the compact star to follow Krori-Barua
(KB) ansatz [32], which effectively relates pressure and density.

NS is a core of collapsed massive star with unique characteristics: Mass
$M=1.4M_{\odot}$ and radius $R=10$ km for canonical NS, rapid rotation 1.39 ms
and intense magnetic field $10^{12}$ Gauss. This makes the NS an astrophysical
lighthouse where radio beams are emitted from its magnetic poles and received
as pulses due to its spin. For this reason, it is called as a pulsar and in
the case of rapid spin with ms period it is called Millisecond Pulsars (MSPs).
It is widely believed that pulsars are in fact neutron stars while others
suggest quark star model [33, 34]. As a matter of fact, pulsars do not only
provide a perfect laboratory to examine the physical properties of matter at
high densities and compactness at Black Hole (BH) lower limits, but also to
prob possible modified gravity scenarios in general and nonminimal coupling
ones in particular. Fortunately, the availability of pulsar data using
different multi-messenger astrophysics [35] with unprecedented precision
provides an extremely important tool to constrain the parameter space of a
suggested modified gravity.

Radio pulse times of arrival (ToAs) from pulsars are being used for decades to
measure pulsars’ masses via relativistic Shapiro time delay, while measuring
the pulsar radius is relatively not an easy task. Neutron Star Interior
Composition Explorer (NICER) is devoted to measure pulsars radii by tracing
hot-spot regions on the pulsar surface which emit X-rays [36]. NICER
measurements of the radii of the pulsars PSR J0030+0451 and PSR J0740+6620
represent a challenge to theoretical models, since their masses are
respectively $M=1.44^{+0.15}_{-0.14}M_{\odot}$ and $M=2.08\pm 0.07M_{\odot}$
[37, 38] while their radii are almost $R\sim 13$ km. More precisely NICER
measures the pulsar PSR J0030+0451 radius $R=13.02_{-1.06}^{+1.24}$ km [39]
and another independent measurement $R=12.71^{+1.14}_{-1.19}$ km [40], while
the analysis of NICER + X-ray Multi-Mirror (XMM) Newton data of the pulsar PSR
J0740+6620 measures $R=13.7_{-1.5}^{+2.6}$ km [41] and
$R=12.39_{-0.98}^{+1.30}$ km (68% credible level) [42]. In this sense, the
more squeezable models, which suggest the neutrons to squeeze to quarks, seem
to be disfavored by NICER observations. Notably, the latest analysis of PSR
J0030+0451 using an updated NICER + XMM-Newton data (using ST+PDT model)
measures $R=11.71^{+0.88}_{-0.83}$ [43]. If the latest is confirmed, quark
star model could be no longer disfavored. Nevertheless, massive pulsars
$M\simeq 2M_{\odot}$ such as PSR J0740+6620 represents a challenge to the less
squeezable models, since the sound speed should rapidly increase somewhere
inside the stellar core to large values $c_{s}^{2}\sim 0.7c^{2}$. This
requires a problematic nonmonotonical behaviour of the sound speed violating
the conjectured conformal upper bound $c_{s}=c/\sqrt{3}$ [44, 45]. It is noted
that existence of massive pulsars points out a stiff EoS matter where the
sound speed to speed of light ratio becomes large. On the other hand,
gravitational wave signals, GW170817 [46, 47] and GW190425 [48], indicate no
tidal deformability of the measured signals, as confirmed by Laser
Interferometer Gravitational-Wave Observatory (LIGO) and Virgo collaboration,
which points out a soft EoS matter where the speed of sound is not large. It
is to be noted that soft EoS is compatible with the conjectured conformal
upper limit on the sound speed which evidently favors quark star models. On
the other extreme, observations of the low-mass pulsar in supernova remnant
HESS J1731–347, with mass $M=0.77_{-0.17}^{+0.20}$ and radius
$R=10.4^{+0.86}_{-0.78}$, suggest that the pulsar could be a light strange
quark star [49, 50].

We present the organization of the paper in some detail as follows: In Sec.
II, we derive the field equations of QRG. In Sec. III, we set up the stellar
model assuming that the spacetime is described by Static Spherically Symmetric
(SSS) configuration, the metric potentials are given by KB ansatz and the
matter fluid is anisotropic. We also investigate corresponding matter EoS
which has been shown to have a linear pattern and matching conditions which
determine the model parameters. In Sec. IV, we utilize the astrophysical
constraints on the mass and radius of the pulsar PSR J0740+6620 from
NICER+XMM-Newton data to determine the numerical value of the nonminimal
coupling parameter $\epsilon$ as well as the KB model parameters.
Consequently, we determine the radial dependence of the sound speed and the
conformal anomaly inside the pulsar. In addition, we investigate different
constitutional forces which are in hydrostatic equilibrium including the newly
introduced one due to quadratic nonminimal coupling. Similarly, we obtain the
radial dependence of the adiabatic index which helps to distinguish
gravitational bound star from self bound one. Our calculations show that the
pressures and densities fit well with linear model which confirms the linear
pattern of the KB induced EoS. In Sec. V, we discuss some physical features of
the obtained model showing that the observational data excludes $\epsilon\geq
0$, which may indicate some underlying similarity between QRG and conformal
anomaly when quantum fields backreaction is included. We use the Trace Energy
Condition (TEC) to determine the maximum possible compactness
$C=2GMR^{-1}c^{-2}$. Moreover, we investigate in details possible connection
between QRG and MIT bag model of strange quark star which sets further
constraints on the model parameters but from microscopic structure of the
star. In Sec. VI, by assuming that the conformal constraints on the sound
speed to be hold at the core of the pulsar PSR J0952–0607, we estimate the
pulsar radius accordingly. We also provide mass-radius diagram which confirms
the validity of the model to describe compact stars from lightest-to-heaviest
pulsars ever observed including recent constraints on the intermediate mass
range from NICER and LIGO/Virgo collaboration. In Sec. VII, we summarize and
conclude the work. Furthermore, we followed the last section by two appendices
A and B to discuss weak field limit of the field equations and to confront the
model with more pulsars’ observational data of different observational
techniques.

##  II Quadratic Rastall Gravity

We briefly review the main structure and assumptions which lead to GR and
Rastall field equations. Then, we extend the discussion in details to the QRG.
It is well known that the double contraction of Bianchi identity in Riemannian
geometry leads to the divergence-free second order tensor,

$\nabla{{}_{\alpha}}\mathfrak{G}{{}^{\alpha}}{{}_{\beta}}=\nabla{{}_{\alpha}}(\mathfrak{R}{{}^{\alpha}}{{}_{\beta}}-\frac{1}{2}\delta^{\alpha}_{\beta}\mathfrak{R})\equiv
0,$ (1)

where $\nabla_{\alpha}$ denotes the covariant derivative associated with Levi-
Civita linear connection, the second order tensors
$\mathfrak{G}{{}^{\alpha}}{{}_{\beta}}$ and
$\mathfrak{R}{{}^{\alpha}}{{}_{\beta}}$ denote Einstein and Ricci tensors
respectively, whereas the contraction of the latter obtains Ricci scalar
$\mathfrak{R}=g^{\alpha\beta}\mathfrak{R}_{\alpha\beta}$. On the other hand,
the minimal coupling principle (divergence-free stress-energy tensor) has been
introduced to generalize its validity in special relativity to curved
spacetimes. Then, it yields

$\nabla_{\alpha}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=0.$ (2)

This finally led to write the GR field equations

$\mathfrak{G}{{}^{\alpha}}{{}_{\beta}}=\kappa_{E}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}},$
(3)

where $\kappa_{E}$ denotes Einstein coupling constant, which is given in terms
of $G$ is the Newtonian gravitational constant $G$ and the speed of light $c$
by $\kappa_{E}=8\pi G/c^{4}$. The left hand side represents the geometric
sector configuring the spacetime and the right hand side represents the matter
sector whereas matter-geometry nonminimal coupling vanishes.

In Rastall’s gravity, the minimal coupling principle (2) is no longer assumed,
and replaced by [1, 2]

$\nabla_{\alpha}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=\varepsilon\,\partial_{\beta}\mathfrak{R}.$
(4)

We use $\varepsilon\to-\frac{\varepsilon}{\kappa_{R}}$ where the negative sign
is a convention, $\varepsilon$ is dimensionless and $\kappa_{R}$ is Rastall’s
coupling constant. This leads to a modified version of Einstein field
equations which account for matter-geometry nonminimal coupling corresponds to
Rastall assumption (4) as below

$\mathfrak{G}{{}^{\alpha}}{{}_{\beta}}=\mathfrak{R}{{}^{\alpha}}{{}_{\beta}}-\frac{1}{2}\delta^{\alpha}_{\beta}\mathfrak{R}=\kappa_{R}(\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}+\frac{\varepsilon}{\kappa_{R}}\delta^{\alpha}_{\beta}\mathfrak{R}).$
(5)

The contraction of the above field equation,

$\mathfrak{R}=-\frac{\kappa_{R}}{1+4\epsilon}\mathfrak{T},$

puts a restriction $\epsilon\neq-\frac{1}{4}$. The weak field limit of the
above field equations yields [c.f. 1, 51, 52]

$\kappa_{R}=\eta\,\kappa_{E}\text{ and
}\eta=\frac{1+4\varepsilon}{1+6\varepsilon},\qquad(\varepsilon\neq-\frac{1}{6})$
(6)

where $\kappa_{R}$ has the same dimension of $\kappa_{E}$, whereas
$\kappa_{R}=\kappa_{E}$ if and only if $\eta=1$ (equivalently
$\varepsilon=0$).

In the present paper, we take the matter-geometry nonminimal coupling to be
proportional to the gradient of quadratic curvature invariants as follows

$\nabla_{\alpha}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=-\frac{\epsilon\ell^{2}}{\kappa}\,\partial_{\beta}\left(\mathfrak{R}^{2}-\mathfrak{R}_{\mu\nu}\mathfrak{R}^{\mu\nu}\right),$
(7)

where $\epsilon$ is a dimensionless parameter, $\ell$ is a constant of
dimension [L] characterizing the mass of the matter distribution, i.e. it
cannot be fitted arbitrarily to observational data111In astrophysical
applications we replace $\ell$ by the star radius $R$ which is a reasonable
choice., and $\kappa$ acquires the same dimension as $\kappa_{E}$. Thus, the
nonminimal coupling is characterized by the dimensionless
$\epsilon$-parameter, where the corresponding modified field equations are

$\mathfrak{G}{{}^{\alpha}}{{}_{\beta}}=\kappa\left[\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}+\frac{\epsilon\ell^{2}}{\kappa}\delta{{}^{\alpha}_{\beta}}\left(\mathfrak{R}^{2}-\mathfrak{R}_{\mu\nu}\mathfrak{R}^{\mu\nu}\right)\right].$
(8)

The above field equations will be referred to as QRG. Contracting the above
equation gives rise to

$-\mathfrak{R}=\kappa\left[\mathfrak{T}+\frac{4\epsilon\ell^{2}}{\kappa}\left(\mathfrak{R}^{2}-\mathfrak{R}_{\mu\nu}\mathfrak{R}^{\mu\nu}\right)\right],$
(9)

where $\mathfrak{T}:=g_{\alpha\beta}\mathfrak{T}^{\alpha\beta}$ is the trace
of the matter stress-energy tensor. For empty spacetime, i.e.
$\mathfrak{T}_{\alpha\beta}=0$, we obtain

$\left(1+4\epsilon\ell^{2}\mathfrak{R}\right)g^{\mu\nu}\mathfrak{R}_{\mu\nu}=4\epsilon\ell^{2}\mathfrak{R}^{\mu\nu}\mathfrak{R}_{\mu\nu}.$
(10)

Clearly $\mathfrak{R}_{\mu\nu}=0$ is a solution of the above constraint, then
GR vacuum solutions are recovered in QRG. More interestingly, for
nondegenerate $\mathfrak{R}_{\mu\nu}$, another vacuum state can be obtained as
$\mathfrak{R}=\frac{-1}{3\epsilon\ell^{2}}$, where $\epsilon\neq 0$, which
gives rise to de Sitter/Anti-de Sitter vacuum according to the sign of
$\epsilon$. The latter vacuum case does not allow recovering GR as
$\epsilon=0$ is not valid and we drop this case in the present study.
Therefore, we still have the privilege of recovering GR vacuum solutions, as
in the Rastall gravity, where $\mathfrak{R}_{\mu\nu}=0=\mathfrak{R}$ in empty
spacetime. However, the QRG includes only quadratic terms of curvature
invariants and the Newtonian limit of the present theory is expected to be
identical to the Newtonian limit of GR unlike RT, i.e. $\kappa=\kappa_{E}=8\pi
G/c^{4}$, see appendix A. Therefore, all solar system tests are fulfilled in
QRG. Another motivation to consider quadratic corrections is the trace anomaly
when the backreaction of the quantum fields to a curved spacetime geometry is
considered. It is well known that the quantum corrections due to trace anomaly
are geometric in nature and can be expressed in terms of Gauss-Bonnet four
dimensional topological invariant (Type A) and the square of the Weyl tensor
(Type B). Therefore, it is reasonable to consider quadratic corrections of
curvature invariants for matter-geometry nonminimal coupling as in equations
(8). Notably the sign of $\epsilon$-like parameter of the trace anomaly must
be negative, it however remains undetermined in our approach at this stage to
be discussed in Subsec. V.1.

It is convenient to rewrite the field equations (8) in the following form

$\mathfrak{G}_{\alpha\beta}=\mathfrak{R}_{\alpha\beta}-\frac{1}{2}g_{\alpha\beta}\mathfrak{R}=\kappa\widetilde{\mathfrak{T}}_{\alpha\beta},$
(11)

where the tensor $\widetilde{\mathfrak{T}}_{\alpha\beta}$ represents the total
stress-energy tensor

$\widetilde{\mathfrak{T}}_{\alpha\beta}:=\mathfrak{T}_{\alpha\beta}+\frac{\epsilon\ell^{2}}{\kappa}g_{\alpha\beta}\left(\mathfrak{R}^{2}-\mathfrak{R}_{\mu\nu}\mathfrak{R}^{\mu\nu}\right).$
(12)

By virtue of (7), it is obvious that $\widetilde{\mathfrak{T}}_{\alpha\beta}$
a divergence-free, i.e.

$\nabla_{\alpha}\widetilde{\mathfrak{T}}{{}^{\alpha}}{{}_{\beta}}=0.$ (13)

It can be shown that (11) has an equivalent form

$\mathfrak{R}_{\alpha\beta}=\kappa\left(\widetilde{\mathfrak{T}}_{\alpha\beta}-\frac{1}{2}g_{\alpha\beta}\widetilde{\mathfrak{T}}\right),$
(14)

where $\widetilde{\mathfrak{T}}:=g_{\mu\nu}\widetilde{\mathfrak{T}}^{\mu\nu}$
represents the trace of the effective stress-energy tensor. Consequently, we
write the following important relation

$\mathfrak{R}^{2}-\mathfrak{R}_{\mu\nu}\mathfrak{R}^{\mu\nu}=\kappa^{2}\left(\widetilde{\mathfrak{T}}{{}^{2}}-\widetilde{\mathfrak{T}}_{\mu\nu}\widetilde{\mathfrak{T}}^{\mu\nu}\right).$
(15)

Substituting into (12), we obtain

$\mathfrak{T}_{\alpha\beta}=\widetilde{\mathfrak{T}}_{\alpha\beta}-\epsilon\kappa\ell^{2}g_{\alpha\beta}\left(\widetilde{\mathfrak{T}}{{}^{2}}-\widetilde{\mathfrak{T}}_{\mu\nu}\widetilde{\mathfrak{T}}^{\mu\nu}\right).$
(16)

The contracted form of the above equations yields

$\mathfrak{T}=\widetilde{\mathfrak{T}}-4\epsilon\kappa\ell^{2}\left(\widetilde{\mathfrak{T}}{{}^{2}}-\widetilde{\mathfrak{T}}_{\mu\nu}\widetilde{\mathfrak{T}}^{\mu\nu}\right).$
(17)

To summarize, the QRG field equations are given by (8), or equivalently by
(16) where the effective energy-stress tensor on the right hand side is
determined by the geometric sector, i.e. Einstein tensor, as given by (11).
The latter choice of the QRG field equations, namely (16), obtains the matter
sector in the left hand side and the geometric sector in the right hand side.
This form clearly describes how the matter sector, in nonempty spacetime,
would behave in QRG distinguishable from GR as long as $\epsilon\neq 0$. On
another word,
$\mathfrak{T}_{\alpha\beta}=\widetilde{\mathfrak{T}}_{\alpha\beta}$ if and
only if $\epsilon=0$. One natural laboratory to test the present proposal is
compact stars where dense matter and high curvature are expected to manifest
matter-geometry nonminimal coupling scenario. Fortunately, recent precise
measurements of mass and radius of compact stars via several astrophysical
observations could set strict bounds on the nonminimal coupling parameter
$\epsilon$.

##  III Static Spherically Symmetric Model

In this study, we assume that the stellar structure is static spherically
symmetric in 4-dimensional spacetime, where the coordinates are spherical
polar ($0\leq t<\infty,0\leq r<\infty,0\leq\theta<2\pi,0\leq\phi\leq\pi$).
This gives rise to the line element

$ds^{2}=-e^{\alpha}c^{2}dt^{2}+e^{\beta}dr^{2}+r^{2}(d\theta^{2}+\sin^{2}\theta\,d\phi^{2}),$
(18)

where $\alpha\equiv\alpha(r)$ and $\beta\equiv\beta(r)$ are the metric
function with only radial dependence. Since we aim to test the present theory
in high dense medium, it is reasonable to assume anisotropic fluid case as
discussed in the introduction, where the stress-energy tensor in this case in
its standard form reads222We use signature convention $(-,+,+,+)$. [27]

$\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=(\rho
c^{2}+p_{t})v{{}^{\alpha}}v{{}_{\beta}}+p_{t}\delta^{\alpha}_{\beta}+(p_{r}-p_{t})w{{}^{\alpha}}w{{}_{\beta}},$
(19)

where

$v^{\alpha}v_{\alpha}=-c^{2},\,w^{\alpha}w_{\alpha}=1,\,v^{\alpha}w_{\alpha}=0,$

and the four-velocity (time-like) vector field
$v^{\alpha}=(ce^{-\alpha/2},0,0,0)$, the unit vector field(space-like) in the
radial direction $w{{}^{\alpha}}=(0,e^{-\beta/2},0,0,0)$. In this case, the
component of the energy-stress tensor
$\mathfrak{T}{{}^{0}}{{}_{0}}=-\rho(r)c^{2}$ represents the fluid density,
$\mathfrak{T}{{}^{1}}{{}_{1}}=p_{r}(r)$ represents the radial pressure
($w_{\parallel}$) and
$\mathfrak{T}{{}^{2}}{{}_{2}}=\mathfrak{T}{{}^{3}}{{}_{3}}=p_{t}(r)$ represent
tangential pressures ($w_{\perp}$), i.e. the matter tensor reduces to the
diagonal form $\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=diag(-\rho
c^{2},\,p_{r},\,p_{t},\,p_{t})$. Similarly, the effective stress-energy tensor
can be given by
$\widetilde{\mathfrak{T}}{{}^{\alpha}}{{}_{\beta}}=diag(-\tilde{\rho}c^{2},\,\tilde{p}_{r},\,\tilde{p}_{t},\,\tilde{p}_{t})$.
Applying QRG field equations (16) to the spacetime (18) where the matter
sector is as given by (19), we obtain the following

$\displaystyle\rho c^{2}$ $\displaystyle=$
$\displaystyle\tilde{\rho}c^{2}-2\epsilon\kappa\ell^{2}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right],$
$\displaystyle p_{r}$ $\displaystyle=$
$\displaystyle\tilde{p}_{r}+2\epsilon\kappa\ell^{2}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right],$
$\displaystyle p_{t}$ $\displaystyle=$
$\displaystyle\tilde{p}_{t}+2\epsilon\kappa\ell^{2}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right].\qquad$
(20)

One should note that the right hand sides of the above equations are
completely determined by the geometric sector, i.e.
$\tilde{\rho}c^{2}=\frac{1}{\kappa}\mathfrak{G}{{}_{0}}{{}^{0}}$,
$\tilde{p}_{r}=\frac{1}{\kappa}\mathfrak{G}{{}_{1}}{{}^{1}}$ and
$\tilde{p}_{t}=\frac{1}{\kappa}\mathfrak{G}{{}_{2}}{{}^{2}}=\frac{1}{\kappa}\mathfrak{G}{{}_{3}}{{}^{3}}$.
We express the anisotropy contribution as the pressure difference,
$\delta(r)=p_{t}-p_{r}$, where $\delta=0$ in isotropic case ($p_{t}=p_{r}$).
Using the above equations, it yields

$\delta(r)=p_{t}-p_{r}=\tilde{p}_{t}-\tilde{p}_{r}=\frac{e^{-\beta}}{4\kappa
r^{2}}\left[(2\alpha^{\prime\prime}-\alpha^{\prime}\beta^{\prime}+\alpha^{\prime
2})r^{2}-2(\alpha^{\prime}+\beta^{\prime})r+4(e^{\beta}-1)\right],$ (21)

where prime and double prime denote first and second derivatives with respect
to the radial coordinate. We remark that the matter-geometry coupling does not
contribute to the anisotropy once the spherically symmetric spacetime
configuration is assumed. Therefore, deviations from GR due to matter-geometry
coupling cannot be spoiled up with various anisotropic effects. For the
$\epsilon=0$ case, the differential equations (20) coincide with Einstein
field equations of an interior spherically symmetrical spacetime [c.f., 53].
Alternatively, the field equations (20) can be rewritten as

$\displaystyle\tilde{\rho}c^{2}$ $\displaystyle=$
$\displaystyle\frac{e^{-\beta}}{\kappa
r^{2}}(e^{\beta}+\beta^{\prime}r-1)=\rho
c^{2}-\frac{1+2\epsilon\kappa\ell^{2}(3\rho
c^{2}-p_{r}-2p_{t})}{12\epsilon\kappa\ell^{2}}-\frac{4\epsilon\kappa\ell^{2}(p_{r}-p_{t})^{2}-2(p_{r}+2p_{t})}{3\left[2\epsilon\kappa\ell^{2}(3\rho
c^{2}+p_{r}+2p_{t})+1-\sqrt{X(r)}\right]},$ $\displaystyle\tilde{p}_{r}$
$\displaystyle=$ $\displaystyle\frac{e^{-\beta}}{\kappa
r^{2}}(1-e^{\beta}+\alpha^{\prime}r)=p_{r}+\frac{1}{4}\left(\rho
c^{2}-p_{r}-2p_{t}\right)+\frac{1-\sqrt{X(r)}}{24\epsilon\kappa\ell^{2}},$
$\displaystyle\tilde{p}_{t}$ $\displaystyle=$
$\displaystyle\frac{e^{-\beta}}{4\kappa
r}\left[(2\alpha^{\prime\prime}-\alpha^{\prime}\beta^{\prime}+\alpha^{\prime
2})r+2(\alpha^{\prime}-\beta^{\prime})\right]=p_{t}+\frac{1}{4}\left(\rho
c^{2}-p_{r}-2p_{t}\right)+\frac{1-\sqrt{X(r)}}{24\epsilon\kappa\ell^{2}},\qquad$
(22)

where

$X(r)=1+12\epsilon\kappa\ell^{2}(\rho
c^{2}-p_{r}-2p_{t})+36\epsilon^{2}\kappa^{2}\ell^{4}\left[\frac{4}{3}p_{t}(\rho
c^{2}-p_{r}+p_{t})+\frac{1}{3}\rho c^{2}(3\rho
c^{2}+2p_{r})+p_{r}^{2}\right].$ (23)

The above equations (22) allow to write the geometric sector on the left hand
side in terms of the matter sector including nonminimal coupling terms. It is
straightforward to check that the effective density and pressures reduce to
the matter ones when the nonminimal coupling contribution vanishes as
$\epsilon\to 0$.

###  III.1 Krori-Barua ansatz

We assume that the metric potentials inside the stellar model are as defined
by Krori-Barua ansatz [32], i.e.

$\alpha(r)=a_{0}r^{2}/R^{2}+a_{1},\,\beta(r)=a_{2}r^{2}/R^{2},$ (24)

where $R$ denotes the star radius. The dimensionless KB model parameters
{$a_{0},a_{1},a_{2}$} are determined by matching conditions on the stellar
boundary. Although KB potentials have been used widely in the literature, c.f.
[53, 3, 4], it is worth to mention that the final output depends on the
gravitational theory. The KB ansatz could be understood as an alternative of
an equation of state. It is mainly used to obtain a nonsingular interior
behaviour within the stellar configuration. Since the nature of matter inside
compact stars at high densities and pressures is not known, we take the other
approach by assuming KB potentials instead of assuming equations of state.
However, the main differences would be a consequence of the modified gravity
sector which clearly differentiates the present work from [3]. Similarly, one
can use same equation of state of the matter sector, i.e. radiation, cold dark
matter, in cosmological applications while the consequences will be different
according to the modified gravity theory under investigation. It proves
convenient to write the field equations in dimensionless forms, therefore we
introduce the dimensionless radius $0\leq x=r/R\leq 1$ and the dimensionless
quantities

$\bar{\rho}(r)=\frac{\rho(r)}{\rho_{\star}},\,\bar{p}_{r}(r)=\frac{p_{r}(r)}{\rho_{\star}c^{2}},\,\bar{p}_{t}(r)=\frac{p_{t}(r)}{\rho_{\star}c^{2}}\,\text{and}\,\bar{\delta}(r)=\frac{\delta(r)}{\rho_{\star}c^{2}},$
(25)

where $\rho_{\star}=\frac{1}{\kappa c^{2}R^{2}}$ denotes a characteristic
density. Additionally, it is reasonable to take the characteristic length
$\ell=R$. Substituting into the field equations (20) yields

$\displaystyle\bar{\rho}$ $\displaystyle=$
$\displaystyle\frac{e^{-a_{2}x^{2}}}{x^{2}}(e^{a_{2}x^{2}}-1+2a_{2}x^{2})$
$\displaystyle-$
$\displaystyle\frac{2\epsilon}{x^{4}}\left[\left(-a_{0}^{2}(a_{0}-a_{2})^{2}x^{8}-8a_{0}(a_{0}-a_{2})(a_{0}-3a_{2}/4)x^{6}+(24a_{0}a_{2}-16a_{0}^{2}-5a_{2}^{2})x^{4}\right.\right.$
$\displaystyle\left.\left.-2(5a_{0}-3a_{2})x^{2}-1\right)e^{-2a_{2}x^{2}}+2\left(1+2a_{0}(a_{0}-a_{2})x^{4}+(5a_{0}-3a_{2})x^{2}\right)e^{-a_{2}x^{2}}-1\right],$
$\displaystyle\bar{p}_{r}$ $\displaystyle=$
$\displaystyle\frac{e^{-a_{2}x^{2}}}{x^{2}}(1-e^{a_{2}x^{2}}+2a_{0}x^{2})$
$\displaystyle+$
$\displaystyle\frac{2\epsilon}{x^{4}}\left[\left(-a_{0}^{2}(a_{0}-a_{2})^{2}x^{8}-8a_{0}(a_{0}-a_{2})(a_{0}-3a_{2}/4)x^{6}+(24a_{0}a_{2}-16a_{0}^{2}-5a_{2}^{2})x^{4}\right.\right.$
$\displaystyle\left.\left.-2(5a_{0}-3a_{2})x^{2}-1\right)e^{-2a_{2}x^{2}}+2\left(1+2a_{0}(a_{0}-a_{2})x^{4}+(5a_{0}-3a_{2})x^{2}\right)e^{-a_{2}x^{2}}-1\right],$
$\displaystyle\bar{p}_{t}$ $\displaystyle=$ $\displaystyle
e^{-a_{2}x^{2}}(2a_{0}-a_{2}+a_{0}(a_{0}-a_{2})x^{2})$ $\displaystyle+$
$\displaystyle\frac{2\epsilon}{x^{4}}\left[\left(-a_{0}^{2}(a_{0}-a_{2})^{2}x^{8}-8a_{0}(a_{0}-a_{2})(a_{0}-3a_{2}/4)x^{6}+(24a_{0}a_{2}-16a_{0}^{2}-5a_{2}^{2})x^{4}\right.\right.$
$\displaystyle\left.\left.-2(5a_{0}-3a_{2})x^{2}-1\right)e^{-2a_{2}x^{2}}+2\left(1+2a_{0}(a_{0}-a_{2})x^{4}+(5a_{0}-3a_{2})x^{2}\right)e^{-a_{2}x^{2}}-1\right],$

and the anisotropy indicator (21) reads

$\bar{\delta}=\frac{e^{-a_{2}x^{2}}}{x^{2}}\left[e^{a_{2}x^{2}}-1+a_{0}(a_{0}-a_{2})x^{4}-a_{2}x^{2}\right].$
(27)

The mass function which obtains the matter mass content within a radius $r$ is
given by

$\mathfrak{M}(r)=4\pi\int_{\xi=0}^{r}\rho(\xi)\,\xi^{2}d\xi.$ (28)

Recalling the density profile (22), we find

$\mathfrak{M}(x)=\frac{M}{C}\left[x(1-e^{-a_{2}x^{2}})+\epsilon\,\zeta(x)\right].$
(29)

where the compactness parameter $C$ and the function $\zeta(x)$ are

$\displaystyle C$ $\displaystyle=$ $\displaystyle\frac{2GM}{c^{2}R},$
$\displaystyle\zeta(x)$ $\displaystyle=$
$\displaystyle\tfrac{a_{2}^{-\frac{7}{2}}}{2x}\left\\{\right.\left[\right.-\tfrac{1}{4}(a_{0}x^{2}(78+51a_{0}x^{2}+8a_{0}^{2}x^{4})+16)a_{2}^{\frac{7}{2}}+\left(\right.(a_{0}^{2}x^{4}+\tfrac{11}{2}a_{0}x^{2}+\tfrac{103}{16})a_{0}^{2}a_{2}^{\frac{5}{2}}$
(30)
$\displaystyle+(6a_{0}x^{2}+a_{0}^{2}x^{4}+5)a_{2}^{\frac{9}{2}}+(\tfrac{5}{4}a_{0}^{4}x^{2}+\tfrac{33}{8}a_{0}^{3})a_{2}^{\frac{3}{2}}+\tfrac{15\sqrt{a_{2}}a_{0}^{4}}{16}\left.\right)x^{2}\left.\right]e^{-2a_{2}x^{2}}$
$\displaystyle-8a_{2}^{\frac{5}{2}}\left[a_{2}-(a_{0}-a_{2})a_{0}x^{2}\right]e^{-a_{2}x^{2}}+\tfrac{\sqrt{2\pi}15x}{64}\text{erf}(\sqrt{2a_{2}}x)\left[\tfrac{112a_{2}^{4}}{3}-\tfrac{328a_{2}^{3}a_{0}}{15}-\tfrac{22a_{2}a_{0}^{3}}{5}-\tfrac{103a_{0}^{2}a_{2}^{2}}{15}-a_{0}^{4}\right]$
$\displaystyle+4\sqrt{\pi}xa_{2}^{2}(a_{0}+5a_{2})(a_{0}-a_{2})\text{erf}(\sqrt{a_{2}}x)+4a_{2}^{\frac{7}{2}}\left.\right\\}.$

Obviously, the mass function (29) determines $\mathfrak{M}(x\to 0)=0$, and
reduces to the GR version [53] when $\epsilon\to 0$ .

###  III.2 Equation of state

We note that the field equations (20) contain five unknown functions: Three
matter sector functions, namely the density, radial pressure and tangential
one ($\rho,~{}p_{r},~{}p_{t}$), on the other hand there are two geometric
sector functions, namely the metric potentials ($\alpha,\beta$). Therefore, in
principle, two additional constraints are required in order to close the
system. Those constraints could be EoSs $p_{r}(\rho)$ and $p_{t}(\rho)$, we
instead use the KB ansatz (24). However, this geometric approach enables us to
write the field equations (22) as given by the following relations up to
$O(x^{8})$

$\bar{p}_{r}(\bar{\rho})\approx
c_{1}\bar{\rho}+c_{0},\quad\bar{p}_{t}(\bar{\rho})\approx\tilde{c}_{1}\bar{\rho}+\tilde{c}_{0},$
(31)

where

$\displaystyle c_{0}$ $\displaystyle=$
$\displaystyle{64(a_{0}+a_{2})\left[a_{2}^{2}/32+\left(a_{0}^{3}-15a_{0}^{2}a_{2}/4+35a_{0}a_{2}^{2}/8-a_{2}^{3}\right)\epsilon\right]\over\left(32a_{0}^{3}-168a_{0}^{2}a_{2}+260a_{0}a_{2}^{2}-80a_{2}^{3}\right)\epsilon-5a_{2}^{2}},$
$\displaystyle c_{1}$ $\displaystyle=$
$\displaystyle-1-{4a_{2}(a_{0}+a_{2})\over\left(32a_{0}^{3}-168a_{0}^{2}a_{2}+260a_{0}a_{2}^{2}-80a_{2}^{3}\right)\epsilon-5a_{2}^{2}},$
$\displaystyle\tilde{c}_{0}$ $\displaystyle=$
$\displaystyle{-a_{2}^{3}+8a_{0}a_{2}^{2}-6a_{0}^{2}a_{2}+\left(16a_{0}^{4}-8a_{2}a_{0}^{3}-152a_{2}^{2}a_{0}^{2}+324a_{2}^{3}a_{0}-88a_{2}^{4}\right)\epsilon\over\left(32a_{0}^{3}-168a_{0}^{2}a_{2}+260a_{0}a_{2}^{2}-80a_{2}^{3}\right)\epsilon-5a_{2}^{2}},$
$\displaystyle\tilde{c}_{1}$ $\displaystyle=$
$\displaystyle-1-{3(a_{0}+a_{2})^{2}-5a_{0}^{2}\over\left(32a_{0}^{3}-168a_{0}^{2}a_{2}+260a_{0}a_{2}^{2}-80a_{2}^{3}\right)\epsilon-5a_{2}^{2}}.$
(32)

The above relations are of an extreme importance when confront the model with
astrophysical observations. For example, vanishing of the parameter $c_{0}$
indicates that the star would be gravitationally bound (radial pressure at the
stellar surface vanishes as density approaches zero), while significant value
of this parameter indicates that the star would be self-bound (radial pressure
vanishes at significant surface density). This reflects possible composition
of the internal structure of the star core. On the other hand, the parameter
$c_{1}$ is strongly related to the radial speed of sound inside the star
fluid, $c_{s}^{2}=dp_{r}/d\rho$, see Subsec. IV.2, which indicates how much
soft/stiff its EoS. In this sense, stiff EoS with gravitationally bound star
may point out presence of hadron matter, whereas soft EoS with self-bound star
would indicate strange quark matter content. This will be discussed in details
in Subsec. V.3. Fortunately, all parameters in (III.2) are completely
determined by the compactness of the star and the nonminimal coupling
parameter $\epsilon$, which in turn distinguishes QRG from GR, as will be
discussed in the following subsection.

###  III.3 Matching conditions

It is easy to verify that exterior vacuum region in QRG is exactly the GR
Schwarzschild solution, since both theories are identical in vacuum
($\mathfrak{R}_{\alpha\beta}=0$). Therefore, the exterior vacuum solution is
given as

$ds^{2}=-\left(1-\frac{2GM}{c^{2}r}\right)c^{2}dt^{2}+\left(1-\frac{2GM}{c^{2}r}\right)^{-1}dr^{2}+r^{2}(d\theta^{2}+\sin^{2}\theta
d\phi^{2}).$ (33)

Recalling the interior spacetime (18), utilizing KB ansatz (24), we therefore
have the boundary conditions (at $r=R$)

$a_{0}+a_{1}=\ln(1-C)\,\text{and}\,a_{2}=-\ln(1-C).$ (34)

In addition to the vanishing of the radial pressure (22) on the star boundary,
i.e.

$\displaystyle\bar{p}_{r}(x=1)$ $\displaystyle=$ $\displaystyle
a_{0}+\frac{1}{2}(1-e^{a_{2}})+\epsilon\left[\right.e^{a_{2}}+2(2a_{2}-5)a_{0}+2(3a_{2}-2a_{0}^{2})-2+e^{-a_{2}}\left\\{\right.a_{0}^{4}-2\times$
(35)
$\displaystyle(a_{2}-4)a_{0}^{3}+(16-14a_{2}+a_{2}^{2})a_{0}^{2}+2(3a_{2}^{2}-12a_{2}+5)a_{0}+(5a_{2}-6)a_{2}+1\left.\right\\}\left.\right]=0.\qquad$

Solving for the set of constants {$a_{0}$, $a_{1}$, $a_{2}$}, we obtain

$\displaystyle a_{0}$ $\displaystyle=$
$\displaystyle\frac{C}{2(1-C)}-\chi(\epsilon,C),$ $\displaystyle a_{1}$
$\displaystyle=$ $\displaystyle\ln(1-C)-\frac{C}{2(1-C)}+\chi(\epsilon,C),$
$\displaystyle a_{2}$ $\displaystyle=$ $\displaystyle-\ln(1-C),$ (36)

where $\chi$ is known function, but lengthy, fulfills the following limit

$\displaystyle\lim_{\epsilon\to 0}\chi(\epsilon,C)=0,$ (37)

and the GR limit then is recovered [53]. We note that all the physical
quantities within the KB stellar model, i.e. $\bar{\rho}(\epsilon,C)$,
$\bar{p}_{r}(\epsilon,C)$ and $\bar{p}_{t}(\epsilon,C)$, are given in terms of
two parameters. In this case, observational data provides an extremely
important tool to constrain the parameter space {$\epsilon,C$} of the matter-
geometry nonminimal scenario as suggested by the QRG in comparison with the GR
theory. On the other hand, it puts an upper bound on the stellar compactness
and consequently its maximum mass.

##  IV PSR J0740+6620 Astrophysical Constraints

In the following we use NICER measurements of the pulsar PSR J0740+6620 where
the hotspot regions are not complicated as in PSR J0030+045 case. Since the
pulsar PSR J0740+6620 is in a binary system, its mass can be determined with
high accuracy independent from inclination. Using the relativistic Shapiro
time delay of ToAs, its mass has been determined as $M=2.08\pm 0.07M_{\odot}$
[37, 38]. Since PSR J0740+6620 is a massive pulsar, it is a perfect candidate
to examine possible modified gravity scenarios. Recent NICER observations give
well constraints on its radius too, it determines mass its radius
$R=13.7_{-1.5}^{+2.6}$ km and $R=12.39_{-0.98}^{+1.30}$ km as indicated by two
independent analysis of NICER data [41, 42]. It is to be noted that the
combined XMM-Newton and NICER dataset enhances the low count rate of the
latter, which determines $M=2.07\pm 0.11M_{\odot}$ and
$R=12.34^{+1.89}_{-1.67}$ km [54]. In our investigation, we use this
measurement which is at 68% CL agreement with [55]. In this sense, the pulsar
PSR J0740+6620 is a perfect regime to constrain the parameter space of the
QRG, i.e. the set of parameters {$\epsilon$, $a_{0}$, $a_{1}$, $a_{2}$} or
simply {$\epsilon$, $C$}.

###  IV.1 The model parameters

The mean value of the pulsar compactness $C=0.496$. For $\epsilon=-0.01$ and
by using the mass function (29), we estimate the mass of the pulsar PSR
J0740+6620 to be $M=2.09M_{\odot}$ at radius $R=12.34$ km which is in a
perfect agreement with its measured value $M=2.07\pm 0.11M_{\odot}$ and
$R=12.34^{+1.89}_{-1.67}$ as seen in Fig. 1(a). Using the matching conditions
solution (36), we determine the model parameters $a_{0}=0.496$, $a_{1}=-1.181$
and $a_{2}=0.684$. In Fig. 1(b), we show the regularity of the metric
potentials $g_{tt}$ and $g_{rr}$ on the stellar boundary, where KB spacetime
describes the stellar interior while the outer region is determined by
Schwarzschild vacuum solution. Those numerical values will be used in all
figures related to the pulsar PSR J0740+6620. We note that other possibly
choices for $\epsilon>0$ would estimate masses in agreement with observational
data. Nevertheless, it would spoil up other physical constraints. This will be
explained shortly in Section V.

(a) Mass function

(b) Matching conditions

Figure 1: Astrophysical constraints from PSR J0740+6620 observations: 1(a) The
mass function (29) is in a perfect agreement with NICER+XMM constraints on the
mass and radius of PSR J0740+6620, $M=2.07\pm 0.11M_{\odot}$ and
$R=12.34^{+1.89}_{-1.67}$ km [54]. 1(b) The metric potentials $g_{tt}$ and
$g_{rr}$ consistently match the QRG interior of the pulsar PSR J0740+6620 to
the exterior Schwarzschild vacuum where nonminimal coupling vanishes.

At the center of PSR J0740+6620 we calculate $\rho(0)=7.34\times 10^{14}$
[g/cm3] $\approx 2.72\rho_{\text{nuc}}$ and $p_{r}(0)=p_{t}(0)=8.72\times
10^{34}$ [dyn/cm2]. At the surface $\rho(R)=4.19\times 10^{14}$ [g/cm3]
$\approx 1.55\rho_{\text{nuc}}$, $p_{r}(R)=0$ [dyn/cm2] and
$p_{t}(R)=3.27\times 10^{34}$ [dyn/cm2]. Also we calculate the redshift
$Z(r):=-1+1/\sqrt{-g_{tt}}$ at the pulsar center $Z(0)=0.80$ and at its
surface $Z(R)=0.41$. In comparison to GR we obtain $\rho(0)=7.22\times
10^{14}$ [g/cm3] $\approx 2.68\rho_{\text{nuc}}$,
$p_{r}(0)=p_{t}(0)=9.43\times 10^{34}$ [dyn/cm2] and $Z(0)=0.80$. While at the
surface we obtain $\rho(R)=4.17\times 10^{14}$ [g/cm3] $\approx
1.55\rho_{\text{nuc}}$, $p_{r}(R)=0$ [dyn/cm2] and $p_{t}(R)=3.24\times
10^{34}$ [dyn/cm2] and $Z(R)=0.41$. Although both QRG and GR predict almost
same densities, pressures and redshifts some other physical properties of high
importance are significantly different. The next subsection is devoted to
discuss this issue.

###  IV.2 Speed of sound and conformal anomaly

One of the basic quantities which determines the interior of the compact star
is the adiabatic sound speed. Particularly, in the radial direction

$v_{r}^{2}=\frac{dp_{r}}{d\rho}=\frac{\bar{p}^{\prime}_{r}}{\bar{\rho}^{\prime}}.$
(38)

Obviously, it is constrained by causality and thermodynamic stability as
$0\leq v_{r}^{2}\leq 1$. On the other hand, it is strongly related to the EoS
as discussed in Subsec. III.2. It is well known that the sound speed square
$v_{r}^{2}\ll c^{2}$ of non-relativistic matter at densities much lower than
$\rho_{\text{nuc}}$, while it approaches the conformal upper bound
$v_{r}^{2}=c^{2}/3$ at much high densities (ultra-relativistic matter)
$\rho\sim 40\rho_{\text{nuc}}$ where perturbative quantum chromodynamics
(pQCD) approach is applied. At densities few times $\rho_{\text{nuc}}$ which
is compatible with compact star density, two scenarios are suggested: The
conjectured conformal upper limit is broken where the sound speed increases
non-monotonically to a peak and decreases toward the conformal sound speed
limit from below at high densities. The other scenario is to keep the
conformal constraint at all densities where the sound speed increases
monotonically to approach the upper limit $v_{r}^{2}/c^{2}\to 1/3$. The first
scenario is suggested to explain existence of massive pulsars, $M\gtrsim
2M_{\odot}$, where the sound speed peaks up to $v_{r}^{2}\sim 0.7c^{2}$ at
densities $\rho\sim 3-5\rho_{\text{nuc}}$ somewhere inside the star (for an
isotropic perfect fluid). However, the tidal deformability of the observed
gravitational wave signals of NS-NS merging GW170817 and GW190425 sets an
upper bound on the stiffness of the EoS inconsistent with the first scenario.

In fact, we show that the conjectured conformal sound speed can applied
everywhere, according to the present model, inside such a heavy object as in
the case of the compact star PSR J0740+6620. Using (LABEL:eq:Feqs3) and (38),
we plot the radial dependence of adiabatic sound speed within the pulsar PSR
J0740+6620 as shown in Fig. 22(a). Unlike the case when GR is applied, the QRG
predicts a sound speed consistent with the conjecture conformal limit
everywhere inside the pulsar where $v_{r}^{2}/c^{2}\sim 0.28$ at the star
surface and $v_{r}^{2}/c^{2}\sim 0.32$ at its inner core. Remarkably when the
nonminimal coupling scenario is dropped ($\epsilon=0$) the conformal sound
speed is mostly no longer hold for a radial distance $0\leq r\lesssim 9.59$
km, where the sound speed approaches $v_{r}^{2}/c^{2}\simeq 0.37$ at the
stellar center. We interpret the obtained result as follows: The strong
anisotropic conditions ($p_{t}>p_{r}$) induces a repulsive force to oppose the
gravitational collapse and in turn allows for more mass with no need for
higher stiffness (large sound speed). In fact, this force cannot reduce the
sound speed enough as already obtained for the GR case. The nonminimal
coupling, however, contributes as an additional regulating repulsive force
holding the sound speed below the conformal limit. This argument will be
discussed in more detail in the following subsection.

(a) Sound speed

(b) Conformal anomaly

Figure 2: Distinguished features of QRG: 2(a) The radial dependence of the
adiabatic sound speed (38) inside the pulsar PSR J0740+6620. The conformal
upper limit of the sound speed, $v_{r}^{2}/c^{2}\leq 1/3$, holds everywhere
when nonminimal coupling does not vanish. This distinguishes the QRG from GR
in which the conformal limit is violated for radial distance $0\leq r\lesssim
9.59$. 2(b) The radial dependence of the trace anomaly (39) shows that
$\Delta>0$ everywhere which indicates that the conformal symmetry is broken
inside the pulsar PSR J0740+6620. It decreases monotonically toward the center
where $\Delta$ gets higher values than in GR systematically toward the stellar
center.

Another physical quantity which is related to our approach is the conformal
trace anomaly

$\Delta:=\frac{1}{3}\frac{g_{\alpha\beta}\mathfrak{T}^{\alpha\beta}}{\rho
c^{2}}=\frac{1}{3}\left[1-\frac{p_{r}+2p_{t}}{\rho c^{2}}\right],$ (39)

which provides a dimensionless measure of possible deviation from conformal
symmetry. At high densities, pQCD exhibits conformal symmetry where $\Delta\to
0$. At finite densities and/or temperatures, which is compatible with NS
densities, it has been argued that the conformal symmetry is broken due to
running of the strong coupling with energy at quantum level where $\Delta\neq
0$ [19]. In general, thermodynamic stability and causality put a constraint on
the trace anomaly $-2/3\leq\Delta\leq 1/3$, whereas its value and sign inside
NSs are not known. We remind that the matter-geometry nonminimal coupling in
QRG is motivated by the trace anomaly when the backreaction of the quantum
fields to a curved spacetime geometry is considered as suggested by (8). Using
(LABEL:eq:Feqs3) and (39), we plot the radial dependence of trace anomaly
within the pulsar PSR J0740+6620 as shown in Fig. 22(b). Clearly, the
conformal trace anomaly $\Delta$ is positive everywhere inside the pulsar,
which means that the conformal symmetry is broken inside the pulsar. It
approaches the vacuum value $\Delta\to 1/3$ at the pulsar surface and
monotonically decreases toward $\Delta\to 0.2$ at the center of the pulsar.
Similar behaviour has been obtained for typical NS with mass $M\sim
1.4M_{\odot}$ [56]. We note that the trace anomaly at the surface is not
exactly $\Delta(R)=1/3$, since $\rho(R)\neq 0$ and $p_{t}(R)\neq 0$.
Therefore, it may point out that the pulsar PSR J0740+6620 is self-bound.
Similar behaviour is shown when the nonminimal coupling vanishes
($\epsilon=0$), it also shows that the QRG exhibits more positivity of the
trace anomaly than GR. This will affect the maximum compactness allowed by
both theories.

###  IV.3 Hydrostatic equilibrium and thermodynamic stability

The assumption
$\nabla_{\alpha}\widetilde{\mathfrak{T}}{{}^{\alpha}}{{}_{\beta}}=0$
constrains the spherically symmetric star of anisotropic fluid to be in static
gravitational equilibrium according to the quadratic nonminimal coupling
model. This modifies the Tolman–Oppenheimer–Volkoff (TOV) equation to read

$p^{\prime}_{r}=-\frac{1}{2}\alpha^{\prime}(\rho
c^{2}+p_{r})+2\frac{\delta}{r}+\widetilde{B}^{\prime}.$ (40)

where
$\widetilde{B}=2{\epsilon\kappa\ell^{2}}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right]$.
The left hand side of the above equation represents a repulsive force–given
that the radial pressure gradient is negative everywhere inside the star– due
to hydrostatic pressure $p^{\prime}_{r}=:-F_{h}$. The first term of the right
hand side represents the relativistic gravitational force
$F_{g}:=-\frac{1}{2}\alpha^{\prime}(\rho
c^{2}+p_{r})=-\frac{a_{0}r}{R^{2}}(\rho c^{2}+p_{r})$, which is attractive in
its nature as the density and the pressure are positive everywhere inside the
star, and clearly $a_{0}$ must be positive. The second and the third terms,
$F_{a}:=2\delta/r=2(p_{t}-p_{r})/r$ and $F_{QRG}:=\widetilde{B}^{\prime}$,
represent two additional forces due the anisotropy and the quadratic
nonminimal coupling contributions to the hydrostatic equilibrium. Therefore,
we write the balanced forces $F_{h}+F_{g}+F_{a}+F_{QRG}=0$ which are
represented by Fig 33(a) for the interior of the pulsar PSR J0740+6620.

The figure clearly shows that the collapsing (attractive) gravitational force
is balanced by the repulsive forces where the anisotropic (repulsive) force
$F_{a}>0$ due to strong anisotropic condition $p_{t}>p_{r}$, similarly, for
the matter-geometry nonminimal coupling, it introduces another repulsive force
$F_{QRG}$ which plays a key role to allow for higher values of the maximum
compactness in QRG than GR (see Sec. V.2) and also to determine the size of
the star (see Sec. VI).

(a) Hydrostatic Equilibrium

(b) Thermodynamic Stability

Figure 3: Hydrostatic equilibrium and thermodynamic stability of the pulsar
PSR J0740+6620: 3(a) The contributing forces in modified TOV equation of
hydrostatic equilibrium (40). The newly introduced force due to nonminimal
coupling, $F_{\text{QRG}}$, is clearly repulsive which acts in addition to
anisotropic and hydrostatic forces, $F_{a}$ and $F_{h}$, to balance the
collapsing gravitational force, $F_{g}$. 3(b) The radial dependence of the
adiabatic indices (41) verifies that the pulsar is stable since $\gamma>4/3$
and $\Gamma_{r},~{}\Gamma_{t}>\gamma$. The divergence of the adiabatic index
$\Gamma_{r}$ at the pulsar surface characterizes self-bound objects which
allows for very stable massive compact stars even with low sound speed.

Another physical condition on the neutron star is to examine its stability
against adiabatic perturbation. This can be provided via some physical
constraints on the relativistic adiabatic indices

$\gamma=\frac{4}{3}\left(1+\frac{F_{a}}{2|p^{\prime}_{r}|}\right)_{max},\,\Gamma_{r}=\frac{\rho
c^{2}+p_{r}}{p_{r}}(v_{r}/c)^{2},\,\Gamma_{t}=\frac{\rho
c^{2}+p_{t}}{p_{t}}(v_{t}/c)^{2}.$ (41)

For isotropic spheres, i.e. the anisotropy function $\delta=0$, we obtain
$\gamma=4/3$ and $\Gamma_{t}=\Gamma_{r}$. In the general case of anisotropic
fluid, the star is in neutral equilibrium where $\Gamma=\gamma$, while it is
in a stable equilibrium where $\Gamma>\gamma$ everywhere in the star interior.

Using (LABEL:eq:Feqs3) and (41), we plot the radial dependence of adiabatic
indices within the pulsar PSR J0740+6620 as shown in Fig. 33(b). It can be
figured out that $\gamma>4/3$ inside the star which characterizes strong
anisotropy cases with $p_{t}>p_{r}$, in addition $\Gamma_{r}>\gamma$ and
$\Gamma_{t}>\gamma$ which insure that the pulsar is in stable equilibrium
against radial perturbation [57, 58]. Remarkably, as shown in the figure, the
divergence of the radial adiabatic index at the surface of the pulsar PSR
J0740+6620, in general, characterizes self-bound objects and allows for very
stable heavy stars $\simeq 2M_{\odot}$ even if the sound speed is small [34].
This could provide an evidence for quark matter in massive stars such as the
PSR J0740+6620. Although the present approach is completely geometric, we
showed that the KB ansatz (24) can be used to relate radial/tangential
pressures to density in linear relations as given by (31). In the following
subsection we investigate the EoS which describes the interior of the pulsar
PSR J0740+6620.

###  IV.4 KB induced EoS

The EoS of the matter inside compact stars is the ultimate goal of several
studies. For this reason we rewrite (31) in more physical form

$p_{r}(\rho)\approx v_{r}^{2}(\rho-\rho_{1}),\quad p_{t}(\rho)\approx
v_{t}^{2}(\rho-\rho_{2}),$ (42)

where $c_{1}=v_{r}^{2}/c^{2}$ and $c_{0}=-c_{1}(\rho_{1}/\rho_{\star})$, and
similarly $\tilde{c}_{1}=v_{t}^{2}/c^{2}$ and
$\tilde{c}_{0}=-\tilde{c}_{1}(\rho_{2}/\rho_{\star})$. Since equations (42)
are induced by KB ansatz (24), we refer to them as KB EoSs. For PSR J0740+6620
case, we evaluate $c_{0}=-0.395$, $c_{1}=0.321$, $\tilde{c}_{0}=-0.154$,
$\tilde{c}_{1}=0.206$ and $\rho_{\star}=3.52\times 10^{14}$ g/cm3. This
directly gives $v_{r}^{2}=0.321c^{2}$, surface density
$\rho_{s}=\rho_{1}=4.32\times 10^{14}$ g/cm3, $v_{t}^{2}=0.206c^{2}$ and
$\rho_{2}=2.63\times 10^{14}$ g/cm3, which in turn provides the following EoSs

$\displaystyle p_{r}$ $\displaystyle=$ $\displaystyle
0.32c^{2}(\rho-4.32\times 10^{14})[\text{dyn/cm}^{2}],$ $\displaystyle p_{t}$
$\displaystyle=$ $\displaystyle 0.21c^{2}(\rho-2.63\times
10^{14})[\text{dyn/cm}^{2}].$ (43)

We note that the above EoSs are mostly valid at the center of the pulsar as
required by the series approximation as discussed in Subsec. III.2. Therefore,
the sound speed is identical to the exact value as obtained in Subsec. IV.2,
whereas the surface density according to the above calculations is slightly
higher than the exact value $\rho(R)=4.19\times 10^{14}$ [g/cm3] $\approx
1.55\rho_{\text{nuc}}$ as obtained in Subsec. IV.1. On the other hand, we
generate a sequence of exact values of the pressures and density from the
center to the surface of the pulsar PSR J0740+6620, as labeled by blue
asterisk on the plots of Fig. 4, then we give the best fit EoSs

$\displaystyle p_{r}$ $\displaystyle=$ $\displaystyle
0.31c^{2}(\rho-4.22\times 10^{14})[\text{dyn/cm}^{2}],$ $\displaystyle p_{t}$
$\displaystyle=$ $\displaystyle 0.19c^{2}(\rho-2.36\times
10^{14})[\text{dyn/cm}^{2}].$ (44)

Those are in a perfect agreement with equations (IV.4), which proves the
validity of KB EoSs not only at the center of the pulsar but also at every
point inside it. We plot both KB and the best fit EoSs as shown by Fig. 4 to
verify our numerical calculations. Both coincide at the center with slight
deviations toward the surface of the pulsar as expected.

(a) Radial EoS

(b) Tangential EoS

Figure 4: EoSs of the pulsar PSR J0740+6620: A comparison between the KB
EoSs(IV.4) and the best-fit EoSs (IV.4). Both coincide at the stellar center
while slight systematic deviations appear toward the surface, which should be
understood due to KB series approximation (31). We plot the radial EoSs
$p_{r}(\rho)$ in Fig. 4(a) and the tangential EoSs $p_{t}(\rho)$ in Fig. 4(b).

We summarize our findings as follows: The observational constraints on the
mass and radius of the pulsar PSR J0740+6620 provided a perfect laboratory to
estimate the amount of the matter-geometry nonminimal according to the QRG
approach. It determines the dimensionless nonminimal coupling parameter
$\epsilon=-0.01$. Accordingly, we determined the contribution of the
additional repulsive force due to QRG to the hydrostatic equilibrium. This
allows to reduce the maximum sound speed at the center to fulfil its conformal
upper limit, $v_{r}^{2}/c^{2}=1/3$, which significantly differentiates the QRG
from GR. Moreover, we showed that the conformal symmetry is broken everywhere
where the trace anomaly $\Delta$ is always positive, in particular
$0.20\lesssim\Delta\lesssim 0.28$ for the pulsar PSR J0740+6620. Furthermore,
we showed that the pulsar is stable against radial pulsation since the
adiabatic index $\Gamma_{r}>\gamma>4/3$. In conclusion, the divergence of the
radial adiabatic index at the surface and the validity of the conformal sound
speed suggest that the pulsar is self-bound where the stellar interior
consists of quark matter. More physical features of the QRG and its relevance
to strange quark matter is investigated in the following section.

##  V Physical Features

In the present section, we justify our choice of the negative sign of the
dimensionless parameter $\epsilon$ which determines the amount of matter-
geometry nonminimal coupling due to QRG. In addition, we investigate the
maximum compactness allowed by QRG in comparison to the corresponding value in
GR theory. Also, we investigate possible role of matter-geometry nonminimal
coupling relevant to the microscopic structure of the star. We note that the
results included in the present section are valid for QRG with
$\epsilon=-0.01$, in general, independent from observational data.

###  V.1 negative value of the nonminimal coupling parameter

Our choice of the sign of the parameter $\epsilon$ is not yet justified, since
positive values can produce masses compatible with NICER+XMM-Newton
constraints on the pulsar PSR J0740+6620 too. However, one of the direct
consequences of non negative $\epsilon$ is the violation of conformal
constraints on the sound speed upper limit. For example, for $\epsilon=+0.01$,
we obtain $v_{r}^{2}\simeq 0.43c^{2}$ at the stellar center even higher than
it corresponding GR value $v_{r}^{2}\simeq 0.37c^{2}$ ($\epsilon=0$) as shown
in Fig. 2(a). We interpret this result as the nonminimal coupling force in the
hydrostatic TOV equation becomes negative (collapsing) force for $\epsilon>0$.
In turn, it would squeeze the pulsar mass to smaller sizes with stiffer EoS.
Other consequences on the maximum allowed compactness will be mentioned
shortly in the next subsection.

In order to verify our conclusion, we perform the following simple
calculations. If the conformal sound speed is hold within PSR J0740+6620
according to our claim, and by knowing only its measured mass $M=2.07\pm
0.11M_{\odot}$, the estimated radius then varies according to the sign of
$\epsilon$ as provided by Table 1. In the last column, we calculate the
tension level between the estimated radius, for different $\epsilon$ values,
and the measured one, $R=12.34\pm 1.75$ km, in terms of the standard deviation
$\sigma$. Clearly, the conjectured conformal sound speed excludes
$\epsilon\geq 0$ up to $\gtrsim 1.6\sigma$, whereas positive $\epsilon$ is in
a stronger tension, $\sim 2.85\sigma$, with observational data.

Table 1: Estimated radius, $R_{est}$, of the pulsar PSR J0740+6620 with $M=2.07\pm 0.11M_{\odot}$ for different $\epsilon$-values assuming the conformal constraint on the sound speed is hold, by taking $v_{r}^{2}/c^{2}=0.32$ at the stellar core. The KB model parameters {$a_{0}$, $a_{1}$, $a_{2}$} are listed for each case. The last column shows the tension level, in terms of the standard deviation $\sigma$, with the observed radius $R_{obs}=12.34\pm 1.75$ km [54]. $\epsilon$ | $a_{0}$ | $a_{1}$ | $a_{2}$ |  $R_{est}$ [km] | Tension [$\sigma$]
---|---|---|---|---|---
$-0.01$ | $0.49$ | $-1.17$ | $0.68$ | $12.44\pm 0.66$ | $0.05$
$0$ | $0.33$ | $-0.83$ | $0.50$ | $15.46\pm 0.82$ | $1.61$
$+0.01$ | $0.26$ | $-0.67$ | $0.41$ | $18.03\pm 0.96$ | $2.85$

As we mentioned earlier in Sec. II, our choice of the matter-geometry
nonminimal term in QRG is motivated by the trace anomaly where backreaction of
the quantum fields to a curved spacetime geometry is considered. We note that
the quantum corrections due to trace anomaly of Type A and B are in general
proportional to quadratic curvature invariants. In particular, the coefficient
of Euler density contribution–similar to QRG– to the trace anomaly is
explicitly [18]

$\beta=-\frac{1}{360(4\pi)^{2}}\left(N_{s}+11N_{f}+62N_{\nu}\right),$ (45)

where $N_{s}$, $N_{f}$ and $N_{\nu}$ are the number of scalars, Dirac fermions
and vector fields. As noted in [59], the coefficient $\beta$, in general, must
be negative. The results of Table 1 along with the above comments on the trace
anomaly quantum corrections may indicate an underlying similarity between QRG
approach and trace anomaly by realizing the mutual role of $\beta$ and
$\epsilon$ (precisely $\epsilon\ell^{2}$), if any. It should be understood,
however, the trace anomaly needs a completely different setup which may
manifest solutions with different characteristics. We are currently
investigating this issue.

###  V.2 Trace energy condition and maximum compactness

We have shown that the conformal symmetry is broken inside the pulsar PSR
J0740+6620. This has been featured via nonvanishing trace anomaly (39). Indeed
the positivity of the trace anomaly $\Delta$ as obtained in the present study
indicates the trace of the stress-energy tensor is negative everywhere inside
the pulsar, i.e. $\mathfrak{T}=p_{r}+2p_{t}-\rho c^{2}<0$, which we refer to
as the trace energy condition (TEC). It has been argued that this inequality
could be violated in the most massive NS cores [60]. Although this is not
widely accepted, it could be true when the sound speed violates the conformal
upper bound, i.e. $v_{r}^{2}/c^{2}>1/3$, see the given example by Zeldovich
which violates the TEC while $v_{r}\simeq c$ [61]. This may indicate that the
medium is not at its true ground state similar to the case of interference of
sound waves when high-frequencies are amplified and low-frequencies are
damped. In the present study, we verified that $v_{r}^{2}/c^{2}\lesssim 1/3$
and $\Delta\geq 0$, consequently the TEC must hold everywhere inside the
pulsars.

(a) Trace energy condition

(b) Trace-Compactness diagram

Figure 5: Maximum compactness: 5(a) Using the field equations and the KB
parameters (LABEL:eq:Feqs3) and (36) with $\epsilon=-0.01$, we plot the radial
dependence of the trace $\mathfrak{T}=p_{r}+2p_{t}-\rho c^{2}$ in terms of the
characteristic energy-density $\rho_{\star}c^{2}=\kappa^{-1}R^{-2}$. The TEC,
$\mathfrak{T}\leq 0$, is verified everywhere inside the star, whereas the
maximum compactness accordingly is $C_{\text{max}}\simeq 0.752$ where
$\mathfrak{T}=0$ at the stellar center, which is $\sim 4\%$ higher than the
corresponding value in GR ($\epsilon=0$). 5(b) We combine two phase-spaces,
namely $\Delta-C$ and $\mathfrak{T}-C$, by applying QRG with $\epsilon=-0.01$
at the stellar center, which determine the maximum compactness $\simeq 0.752$
according to the constraints $\Delta=\mathfrak{T}=0$. The $\Delta-C$ diagram
shows that $\Delta\to-\infty$ as $C\to 0.96$. Same value of $C$ can be
obtained for a maximally compact object with EoS,
$p_{r}=v_{r}^{2}(\rho-\rho_{s})$, where $v_{r}=c$.

We utilize the TEC to set a restrictive constraint on the maximum compactness
allowed by the QRG scenario. Using the field equations (LABEL:eq:Feqs3) where
the KB parameters are given by (36), we obtain that the maximum compactness
corresponds to QRG with $\epsilon=-0.01$ is $C_{\text{max}}=0.75$ which is
$\sim 4\%$ higher than the GR case [53]. This results is expected as we
pointed out in Subsec. IV.2, where the trace anomaly $\Delta$ is almost the
same in QRG and GR at the surface, while in QRG it systematically gets higher
positive values toward the center than GR. We plot the radial dependence of
the trace of the matter stress-energy tensor $\mathfrak{T}$ for QRG
($\epsilon=-0.01$) in comparison to the GR ($\epsilon=0$) case as shown by
Fig. 5(a). Notably, for $\epsilon=+0.01>0$, the maximum compactness
$C_{\text{max}}\simeq 0.68$ can be realized.

For the QRG, we plot $\Delta-C$ and $\mathfrak{T}-C$ phase-spaces at the
center as shown by Fig. 5(b). The trace anomaly $\Delta\geq 0$ and
$\mathfrak{T}\leq 0$ as long as $C\leq 0.752$ which determines the maximum
compactness allowed by the QRG. If future observations verify existence of
compact stars with $C>0.752$, it would indicate possibility of positiveness of
the matter trace or equivalently negative value of the trace anomaly at the
stellar center. Remarkably, the trace anomaly diverges, $\Delta\to-\infty$, at
the stellar center as $C\simeq 0.96$, which coincides with the maximum
compactness of a maximally compact EoS, $p_{r}=v_{r}^{2}(\rho-\rho_{s})$ with
$v_{r}=c$, where causality still holds.

###  V.3 From macroscopic to microscopic structure

Up to this moment, we obtained some important features of the QRG on the
macroscopic level of compact stellar structure. The theory modifies the TOV
equation by introducing a new repulsive force to the total force budget. This
force allows for more stellar mass within same size, whereas the adiabatic
sound speed does not violate the conformal upper limit. On the other hand, we
confirmed the validity of the TEC, or equivalently the positiveness of the
trace anomaly, which sets an upper limit on the compactness $C\simeq 0.752$.
The model is characterized by vanishing of the surface radial pressure, while
its surface density is above the nuclear saturation density as obtained by the
KB-EoS $p_{r}=v_{r}^{2}(\rho-\rho_{s})$ with $v_{r}\simeq c/\sqrt{3}$, this
makes the star a self-bound object. This is also confirmed by the divergence
of radial adiabatic index at the stellar surface, which allows for very stable
heavy stars $\simeq 2M_{\odot}$ even if the sound speed is small. Although the
model can describe stable stars macroscopically, one must extend
investigations to examine the model on the microscopic level.

In the case of a self-bound state on the microscopic level, it is natural to
consider the case of quark star model. It has been conjectured that quark
matter could be the true ground state of matter at zero pressure as proposed
by Witten [62], see also [63]. Indeed, deconfinement of quarks process is a
strong interaction which converts baryons into two flavour quarks, up ($u$)
and down ($d$) quarks. Whether $ud$ plasma has been produced inside a NS core
where the density approaches (3$-$5)$\times 10^{14}$ g/cm3, or the NS already
contains small fraction of quark matter of cosmological origin; half of $d$
quarks convert into strange ($s$) quarks via equilibrium strangeness weak
interaction

$u+d\rightarrow u+s.$ (46)

The time suggested, depends on the temperature and the EoSs, for this
conversion $10\,\text{ms}\lesssim\tau\lesssim 10\,\text{min}$, cf. [64, 65,
66]. Even a much shorter supersonic conversion $\tau\sim 0.1\,\text{ms}$ has
been suggested, since hydrodynamical instability of subsonic conversion favors
fast conversion [67, 68]. At any case, strange quark matter hypothesis implies
a new minimum of the energy per baryon $\mu_{b}=\rho c^{2}/n_{b}$, at zero
pressure, lower than the energy per baryon of Iron nuclei 56Fe
($<930~{}\text{MeV}$). In a compact stellar objects, the quark chemical
potential is $\mu\simeq 500~{}\text{MeV}$ at most and the temperature is
negligible (zero in practice). We consider quark matter with three $uds$
flavors, whereas the quark masses are assumed to be negligible with respect to
the chemical potential, which is acceptable approximation in the present
context. In this simplest case, one could adopt an MIT bag model in which
quarks are confined to a bag, $B$, of a QCD vacuum, which have been introduced
as a phenomenological constant to compensate the energy difference between
confined ordinary nuclear matter and deconfined quarks333The bag constant is
measured in MeV4 in the natural units system (i.e. $\hbar=c=1$), whereas it is
measured in MeV/fm3 in nuclear physics units, which will be used in this
section. We remark that $B_{\text{natural}}=B_{\text{nuclear}}(\hbar c)^{3}$..
The energy density, for massless and non-interacting quark, is given by

$\displaystyle\rho c^{2}=\rho_{q}c^{2}+B,$ (47)

where the microscopic quark kinetic energy, $\rho_{q}c^{2}$, of an ultra-
relativistic Fermi gas combined with the perfect flavor-color symmetry (3
color, 3 flavour and 2 spin), is given by [69]

$\rho_{q}c^{2}=\frac{9}{4}\pi^{2/3}\hbar c\,n_{b}^{4/3},$ (48)

where $\hbar$ is the reduced Planck’s constant and $n_{b}$ is the baryon
density number. Consequently, the microscopic pressure

$p_{r}=n_{b}^{2}\frac{d}{dn_{b}}\left(\frac{\rho
c^{2}}{n_{b}}\right)=\frac{3}{4}\pi^{2/3}\hbar c\,n_{b}^{4/3}-B.$ (49)

At the stellar surface ($p_{r}=0$), the baryon density number

$n_{b}=\frac{2\sqrt{2}}{\sqrt{\pi}}\left(\frac{B}{3\hbar c}\right)^{3/4},$
(50)

which consequently implies a self-bound state at zero pressure. This simple
case can be described by the MIT bag model EoS

$p_{r}=\frac{1}{3}(\rho c^{2}-4B).$ (51)

Clearly, the bag constant $B$ is related to the surface density (at which the
radial pressure vanishes) by $\rho_{\text{s}}=4B/c^{2}$. In order to
investigate the similarities between QRG and the MIT bag model, we rewrite the
field equations (20) in the form of

$\tilde{\rho}c^{2}=\rho
c^{2}+\widetilde{B},\,\tilde{p}_{r}=p_{r}-\widetilde{B}\,\text{and}\,\tilde{p}_{t}=p_{t}-\widetilde{B},$
(52)

where
$\widetilde{B}=2{\epsilon\kappa\ell^{2}}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right]$.
This leads to the following relation

$\tilde{p}_{r}=\frac{1}{3}(\tilde{\rho}c^{2}-4B_{eff}),$ (53)

where

$\displaystyle B_{eff}$ $\displaystyle=$ $\displaystyle
B+\widetilde{B}=B+2\epsilon\kappa\ell^{2}\left[\tilde{\rho}c^{2}(\tilde{p}_{r}+2\tilde{p}_{t})-\tilde{p}_{t}(\tilde{p}_{t}+2\tilde{p}_{r})\right].$
(54)

At the GR limit, $\epsilon\to 0$, we obtain $B_{eff}=B$. Equation (53) shows
the role of the nonminimal coupling between matter and geometry due to QRG to
effectively mimic bag constant similar to the cases when QCD corrections are
included.

It remains to question the validity of the QRG provided in this study on the
microscopic level in comparison to the physical constraints of the MIT bag
model. We note that the strange quark star would be valid, where at the
stellar surface, the following constraint is fulfilled

$\mu_{b}=\frac{\rho c^{2}}{n_{b}}=\sqrt{2\pi}(3\hbar
c)^{3/4}B^{1/4}<930\,\text{MeV}.$ (55)

This simply requires an upper bound to the bag constant $B<91.34$ MeV/fm3. On
the other hand, for nuclei to be stable with respect to the $ud$ matter
$B>56.64$ MeV/fm3. A more restrictive lower bound on the bag constant is by
requiring that the neutrons should not coagulate into droplets of the $ud$
matter, which leads to $B>58.92$ MeV/fm3, c.f. [69]. Fortunately, in our case,
the bag constant $B=\rho_{s}c^{2}/4$ in KB EoS, namely (42), is completely
determined by the model parameter, as for example for the pulsar PSR
J0740+6620 we have $\rho_{\text{s}}=-(c_{0}/c_{1})\rho_{\star}\approx
4.32\times 10^{14}$ g/cm3, as obtained in the previous section, which directly
gives

$B=-\frac{1}{4}(c_{0}/c_{1})\rho_{\star}c^{2}\simeq 60.58\,\text{MeV/fm}^{3}.$
(56)

This result confirms the validity of the QRG with the microscopic strange
quark star hypothesis.

We Recall the results of Table 1, which obtain estimated radii of the pulsar
PSR J0740+6620 for different values of $\epsilon$ assuming that the
conjectured conformal upper limit on the sound speed is valid. Therefore, it
is straightforward to see that $B=32.67$ MeV/fm3 for $\epsilon=0$ and
$B=21.08$ MeV/fm3 for $\epsilon=+0.01$ much less than the lower bound, which
means that deconfined quarks would be visible in the atomic nuclei. This
provides another evidence to exclude the GR ($\epsilon=0$) and QRG with
$\epsilon>0$ cases, but from microscopic physics.

##  VI Estimation of PSR J0952–0607 Radius

A recent observation, using multicolor light curves, measures the heaviest
pulsar ever-seen PSR J0952–0607 with mass $M=2.35\pm 0.17\,M_{\odot}$ [70].
The rotational frequency of this ms pulsar is 707 Hz. Since its mass near to
$M_{\text{TOV}}$, the size measure of this pulsar provides an important piece
of information to understand the nature at such dense matters. This is the
main aim of the present section.

###  VI.1 Conformal sound speed constraint on the radius

We assume the conformal constraint on the upper bound of the sound speed is
applied everywhere inside the pulsar PSR J0952–0607. Therefore, we take the
possible maximum value at the stellar center, i.e.

$\lim_{r\to 0}v_{r}^{2}/c^{2}=1/3.$ (57)

Using field equations (LABEL:eq:Feqs3), or simply by realizing the density
coefficient $c_{1}$ in the KB–EoS (III.2), we obtain the following constraint

${a_{2}(a_{0}+a_{2})\over
5a_{2}^{2}-\left(32a_{0}^{3}-168a_{0}^{2}a_{2}+260a_{0}a_{2}^{2}-80a_{2}^{3}\right)\epsilon}=1/3.$
(58)

Inserting (36) into the above inequality and by taking the nonminimal coupling
parameter $\epsilon=-0.01$ as constrained by the astrophysical observations of
the pulsar PSR J0740+6620, we calculate the compactness parameter of the
pulsar PSR J0952–0607 as $C=0.53$. Consequently, the corresponding KB model
parameters read $a_{0}=0.56$, $a_{1}=-1.31$, $a_{2}=0.75$. Given its observed
mass $M=2.35\pm 0.17M_{\odot}$, we obtain the corresponding radius $R=13.21\pm
0.96$ km. This result is in agreement with radii measured values of compact
stars with relevant masses, see Table 2.

We accordingly provide some useful physical quantities. At the center
$\rho(0)=6.99\times 10^{14}$ [g/cm3] $\approx 2.59\rho_{\text{nuc}}$,
$p_{r}(0)=p_{t}(0)=9.18\times 10^{34}$ [dyn/cm2] and $Z(0)=0.92$. At the
surface $\rho(R)=3.80\times 10^{14}$ [g/cm3] $\approx 1.41\rho_{\text{nuc}}$,
$p_{r}(R)=0$ [dyn/cm2], $p_{t}(R)=3.38\times 10^{34}$ [dyn/cm2] and
$Z(R)=0.45$. Similar to PSR J0740+6620, we evaluate the KB EoSs of the pulsar
PSR J0952–0607

$\displaystyle p_{r}$ $\displaystyle=$ $\displaystyle
0.333c^{2}(\rho-3.92\times 10^{14})[\text{dyn/cm}^{2}],$ $\displaystyle p_{t}$
$\displaystyle=$ $\displaystyle 0.214c^{2}(\rho-2.22\times
10^{14})[\text{dyn/cm}^{2}].$ (59)

In addition, we obtain the best fit EoSs

$\displaystyle p_{r}$ $\displaystyle=$ $\displaystyle
0.310c^{2}(\rho-3.69\times 10^{14})[\text{dyn/cm}^{2}],$ $\displaystyle p_{t}$
$\displaystyle=$ $\displaystyle 0.194c^{2}(\rho-2.08\times
10^{14})[\text{dyn/cm}^{2}].$ (60)

Clearly both are in a perfect agreement. In the present model, for the pulsar
PSR J0952–0607, we have $c_{1}\simeq 1/3$, $c_{0}\simeq-0.43$ and
$\rho_{\star}\simeq 151.83$ MeV/fm3. This gives
$B=-\frac{1}{4}(c_{0}/c_{1})\rho_{\star}c^{2}\simeq 58.97$ MeV/fm3 just above
the lower bound. This additionally verifies the validity of the QRG on the
stellar microscopic structure.

If we apply same procedure within GR ($\epsilon=0$), by applying the conformal
constraint on the sound speed upper limit at the center of the pulsar PSR
J0952–0607, then we obtain the compactness $C=0.42$ with a corresponding
radius $R=16.41\pm 1.19$ km, which is in a strong tension with the measured
radii values of relevant stellar masses. On the contrary, if we restrict the
radius to acceptable range $\sim 13$ km, one can not hold the sound speed at
its conformal upper bound. We, additionally, determine the KB model parameters
are $a_{0}=0.37$, $a_{1}=-0.92$, $a_{2}=0.55$. At the stellar center, we
obtain $\rho(0)=3.29\times 10^{14}$ [g/cm3] $\approx 1.22\rho_{\text{nuc}}$,
$p_{r}(0)=p_{t}(0)=3.28\times 10^{34}$ [dyn/cm2] and $Z(0)=0.58$. At the
surface, we obtain $\rho(R)=2.11\times 10^{14}$ [g/cm3] $\approx
0.78\rho_{\text{nuc}}$, $p_{r}(R)=0$ [dyn/cm2] and $p_{t}(R)=1.20\times
10^{34}$ [dyn/cm2] and $Z(R)=0.32$. This in turn determines the bag constant
$B=\frac{1}{4}\rho_{s}c^{2}\simeq 30.72$ MeV/fm3 below the physical lower
bound. Much lower values are obtained for QRG with $\epsilon>0$. We interpret
this result by the presence of an additional collapsing force
$F_{\text{QRG}}<0$ (or it vanishes in the GR case), which requires larger
stellar size (less dense medium) for the conformal upper limit on the sound
speed to be hold. This finally leads to smaller surface density
($\rho_{\text{s}}=4B/c^{2}$) violating the lower bound on the bag constant.
This is unlike the QRG with $\epsilon<0$ which induces an additional repulsive
force allowing for a stable stellar configuration at smaller radii (more dense
medium) compatible with quark matter physical constraints. The present results
show the novel role of the matter-geometry nonminimal coupling within the QRG
to obtain a stellar model consistent with macroscopic and microscopic physical
constraints.

###  VI.2 Mass-radius diagram and crustal structure

One of the characterizing features which distinguishes strange quark star from
normal matter star is the mass-radius relation at small masses, $M\lesssim
0.3M_{\odot}$, where the mass decreases monotonically as $M\propto R^{3}$.
This is radically different from normal matter where mass in fact decreases as
radius increases at relevant low masses. The low-mass strange matter core is
dominated by the confining QCD vacuum, whereas gravitational force is
neglected and quark matter is highly incompressible. In this case, the density
is nearly constant and Newtonian gravity can be applied. Therefore, a mass-
radius diagram would be a useful tool to examine the viability of the model
with different observational datasets with a variety of masses.

In our case, we solve the density profile (LABEL:eq:Feqs3), for radius at
compactness values $0\leq C\leq 1$, assuming that the surface density
$\rho_{s}=3.8\times 10^{14}$ g/cm3 relevant to the obtained value for the
pulsar PSR J0952–0607 via KB EoS (VI.1). Then, we plot the corresponding
mass–radius diagram as shown in Fig. 6 (without crust contribution), which is
a typical strange quark matter EoS pattern. The figure includes observational
constraints of nine pulsars, as provided by Table 2, from the lightest pulsar
HESS J1731–347 (gray solid circle) to the heaviest J0952–0607 (blue solid box)
in addition to other observational constraints on the intermediate mass range
from NICER and LIGO/Virgo collaboration. We note that the radii of the pulsars
J0952–0607, J0348+0432 and J1614–2230 are obtained by assuming that the
conjectured conformal upper bound on the sound speed is applied at the center
of the pulsar. Those are indicated by blue solid boxes in Fig. 6. For NICER
measurements of the pulsar J0030+0451 we include the latest analysis of
updated NICER + XMM-Newton data (using ST+PDT model) [43]. We dropped other
older analysis of NICER data of the pulsar J0030+0451 [39, 40], since they are
in tension with GW170817. All other NICER measurements of the pulsars
J0740+6620 and J0437–4715 are indicated by green solid circles whereas
gravitational wave signals are indicated by red solid diamonds. In general,
the mass-radius curve are in agreement at 68% confidence interval with
observational data except for the pulsar J0437–4715 case which indeed has
relatively large size [71]. If upcoming analysis of this pulsar obtains it
radius of $R\sim 12$ km, the present model would be in agreement. More
pulsars’ observational data by different observational techniques are given in
appendix B as provided by Table 3 with a supplementary numerical estimates of
some important physical features for those pulsars, see Table 4.

Table 2: Astrophysical observations of mass and radius of some pulsars. Pulsar | Mass ($M_{\odot}$) | Radius (km) | Ref.
---|---|---|---
J0952–0607 | $2.35\pm 0.17$ |  $13.21\pm 0.96$444The radius value is predicted by the present model assuming that the conformal constraint on the sound speed $v_{r}^{2}/c^{2}\lesssim 1/3$ is hold at the stellar core. | [70]
J0740+6620 | $2.07\pm 0.11$ | $12.34^{+1.89}_{-1.67}$ | [54]
J0348+0432 | $2.01\pm 0.04$ | $12.79\pm 0.25^{\text{\ref{note1}}}$ | [72]
J1614–2230 | $1.908\pm 0.016$ | $12.54\pm 0.11^{\text{\ref{note1}}}$ | [73, 74, 75]
J0030+0451555We use the result of ST+PDT model which estimates radii in line with current observations [43]. | $1.40^{+0.13}_{-0.12}$ | $11.71^{+0.88}_{-0.83}$ | [43]
J0437–4715 | $1.44\pm 0.07$ | $13.6^{+0.9}_{-0.8}$ | [76, 71]
GW170817-1 | $1.45\pm 0.09$ | $11.9\pm 1.4$ | [47]
GW170817-2 | $1.27\pm 0.09$ | $11.9\pm 1.4$ | [47]
HESS J1731–347 | $0.77_{-0.17}^{+0.20}$ | $10.4^{+0.86}_{-0.78}$ | [49]

We briefly discussed the possibility for a NS to have a small seed of quark
matter of cosmological origin or by forming it at the NS center in Sec. V.3.
Since there is no Coulomb barrier in the stellar interior, the weak
interaction $u+d\to u+s$ consumes free neutrons converting the star entirely
to quark matter. In this case, bare strange star is possible since the surface
density is $10^{3}$ larger than neutron drip density ($\simeq 10^{11}$ g/cm3).
However, thin and light crustal structure of normal matter is also possible
where free neutrons at the surface no longer exist, and there is a
sufficiently large Coulomb barrier between crust normal matter and the
interior quark matter [77]. The development of normal crust matter could be a
NS remnant or via accretion. At any case, it would have a crustal structure
with density at neutron drip density $\simeq 4\times 10^{11}$ g/cm3 at the
bottom of the crust where the light and thin crust approximation is most
likely. Let $\rho_{\text{cr}}$ and $p_{\text{cr}}$ be the density and pressure
at the crust bottom, while the average adiabatic index of the otter crust be
$<\gamma>$. Assuming that the solid crust consists of dense matter below the
neutron drip which can be described by BPS model [78], then the crust mass and
thickness are given by [77],

$\displaystyle\Delta M_{\text{cr}}^{\text{stat}}$ $\displaystyle=$
$\displaystyle 8\pi R^{3}\frac{(1-C)}{C}\frac{p_{\text{cr}}}{c^{2}},$ (61)
$\displaystyle\Delta R_{\text{cr}}^{\text{stat}}$ $\displaystyle=$
$\displaystyle\frac{2<\gamma>}{<\gamma>-1}\frac{R}{C}\frac{p_{\text{cr}}}{\rho_{\text{cr}}c^{2}}.$
(62)

Taking the crust density $\rho_{\text{cr}}=4.3\times 10^{11}$ g/cm3 equal to
the neutron drip density and the crust pressure $p_{\text{cr}}=7.8\times
10^{29}$ dyn/cm2 with an average adiabatic index $<\gamma>=1.28$ of normal
matter. In general, the crust mass is negligible and it is thickness varies
according to the interior mass as shown in Fig. 6. For example, the crust in
the case of the massive pulsar J0740+6620 has a negligible mass contribution
$\Delta M_{\text{cr}}^{\text{stat}}\approx 2.1\times 10^{-5}\,M_{\odot}$ and a
thickness $\Delta R_{\text{cr}}^{\text{stat}}\approx 0.46$ km, while in the
case of the light pulsar HESS J1731–347 the crust mass $\Delta
M_{\text{cr}}^{\text{stat}}\approx 4.41\times 10^{-5}\,M_{\odot}$ is still
tiny with a relatively larger thickness $\Delta
R_{\text{cr}}^{\text{stat}}\approx 0.87$ km. Clearly, when the mass is below
$\sim 0.1M_{\odot}$ the thickness significantly changes toward much larger
values as obtained by Fig. 6.

Figure 6: Mass-Radius diagram: We use $\epsilon=-0.01$ and $\rho_{s}=3.8\times
10^{14}$ g/cm3. The colored regions represent excluded values according to
several physical constraints. The diagonal dotted line represents the
conformal constraint on the sound speed $v_{r}^{2}/c^{2}=1/3$. We present
several observational data as provided by Table 2 including the lightest
pulsar HESS J1731–347 and the heaviest J0740+6620 in addition to our estimated
radius value of the pulsar J0952–0607. Clearly the conformal sound speed limit
is hold for all pulsars. In the low mass range $M\propto R^{3}$ which
characterizes strange quark star models. For a bare quark star ignoring the
crust contribution the radius is underestimated by $<$ 1 km. However, the
crust contribution is provided according to equation (61) for two different
values of adiabatic index $\gamma=1.28$ and $\gamma=1.18$.

##  VII Summary and Conclusion

Neutron stars uniquely provide natural laboratory to test physical conditions
beyond terrestrial-based experiments capability. Since matter is so dense and
the spacetime curvature is efficient in the stellar interior, modification of
general relativity such as matter-geometry nonminimal coupling would be
significant. In the present paper, we provide an extension to GR by including
Rastall-like nonminimal coupling but proportional to the gradient of quadratic
curvature invariants. Those are similar to the conformal trace anomaly when
quantum fields backreaction on the curved spacetime is included.

Recent independent observations which measure mass and radius of pulsars by
LIGO/Virgo and NICER provide unprecedented accurate measurements of those
quantities, which in turn helps not only to understand the matter nature at
high pressure and density conditions and possible EoS but also to constrain
the parameter space of suggested modified gravity. Existence of massive NSs
$\gtrsim 2M_{\odot}$ requires stiff EoS with $v_{r}^{2}/c^{2}\sim 0.7$, on the
contrary to tidal deformability of gravitational wave signals which requires
soft EoS. These evidences challenge theoretical models which in general assume
perfect fluid with hadron EoS when GR is applied.

In the present work, we assume static spherical symmetric spacetime, where the
geometric sector of the stellar interior is governed by KB metric potentials
(24). The fluid is anisotropic (19), which is more realistic in practice,
where the anisotropic force becomes repulsive if $p_{t}>p_{r}$ (strong
anisotropy condition). The general relativistic minimal coupling procedure is
relaxed and replaced by QRG which includes quadratic nonminimal coupling (7).
The QRG field equations are given by equations (8).

Accordingly we obtain the density and pressure profiles, namely equations
(LABEL:eq:Feqs3), by applying the QRG to the above mentioned configuration of
the stellar structure model. The model parameters can be entirely written in
terms of the compactness as in the GR case in addition to another
dimensionless parameter $\epsilon$ which characterizes the nonminimal coupling
deviation from GR (in which $\epsilon=0$). Unlike Rastall gravity the QRG
coupling constant is identical to Einstein’s coupling, since the field
equations in GR and QRG are identical in the weak field limit, see appendix A.
Also, both theories are identical in vacuum, since
$\mathfrak{R}=0=\mathfrak{R}_{\mu\nu}$. In effect Schwarzschild vacuum
solution is valid within QRG framework. However, for nondegenerate
$\mathfrak{R}_{\mu\nu}$, another vacuum state can be obtained in QRG where
$\mathfrak{R}=-1/3\epsilon\ell^{2}$, which has been omitted in the present
study as $\epsilon=0$ is not allowed.

We use the Observational data constraints by NICER+XMM-Newton on the pulsar
PSR J0740+6620 to determine the parameter $\epsilon\sim-0.01$. In addition, if
the conformal upper bound on the sound speed is applied at the center, we show
that $\epsilon\geq 0$ values are excluded by NICER at $\gtrsim 1.6\sigma$
level, see Table 1. This may indicate an underlying similarities between QRG
and conformal anomaly which requires a necessarily negative coefficient. For
the pulsar PSR J0740+6620, NICER constraints confirm that the adiabatic radial
sound speed does not violate the conformal upper bound $v_{r}^{2}/c^{2}=1/3$,
if QRG is assumed unlike GR. The trace anomaly is positive which confirms the
negative value of the matter trace inside the star and indicates that the
conformal symmetry is broken inside the star. This in turn sets an upper bound
on the maximum compactness $C_{\text{max}}=0.752$ which is $4\%$ higher than
GR.

The analysis of the contributing forces of the hydrostatic equilibrium
confirms that, in addition to the repulsive hydrostatic and anisotropic
(strong anisotropy condition) forces, there is an additional repulsive force
due to nonminimal coupling with $\epsilon<0$. This regulating force contribute
to keep the sound speed below the conformal constraint unlike the GR case
which requires sound speed above this limit. The divergence of the radial
adiabatic index at the star surface identifies that the star is self-bound.
This feature together with soft EoS (i.e. $v_{r}^{2}/c^{2}\sim 1/3$) suggests
that the star to be a quark star rather a NS with normal matter.

We show that the KB ansatz effectively relates the density and pressures
linearly up to constants, which we call them KB-EoSs (42). The investigation
on the microscopic structure shows that the nonminimal coupling contributes to
the matter density and pressure effectively as a bag constant of the MIT quark
star model. By requiring the conformal sound speed constraint inside the
pulsar PSR J0740+6620 to be applied, the estimated bag constant is consistent
with the physical range $58.92<B=60.58$ MeV/fm${}^{-3}<91.34$ for
$\epsilon=-0.01$. However, the values $\epsilon\geq 0$ are excluded since they
provide small $B$ below the allowed value.

Since the radius of the pulsar PSR J0952–0607 is not yet measured, we apply
the conformal constraint on the sound speed to estimate its radius according
to the QRG correction. Our estimated value, for $\epsilon=-0.01$, is
$R=13.21\pm 0.96$ km in agreement with the measured radii of other pulsars
with comparable masses. On the microscopic level we confirm the validity of
the numerical values of the model parameters with physical constraints from
MIT bag model. Our calculations show that same procedure can be applied in GR
but the estimated radius then is $R=16.41\pm 1.19$ km which is incompatible
with other pulsars within the same mass range, whereas the case gets worse for
$\epsilon>0$ values. This in turn estimates lower surface density which is
related to small bag constant values not valid within its physical range.

We provide a mass-radius diagram corresponding to the suggested QRG with
$\epsilon=-0.01$ and surface density $\rho_{s}=3.8\times 10^{14}$ g/cm3 above
$\rho_{\text{nuc}}$. The diagram fits well with wide range of observational
data including NICER and LIGO/Vigo collaboration. In conclusion, within QRG
same physics can be applied successfully to low-mass, HESS J1731–347, and
high-mass, PSR J0740+6620, compact stars. Also we confirm that nonminimal
coupling can provide a better framework to keep the sound speed below the
conformal upper limit, while macro and micro physical constraints are
satisfactory fulfilled.

###### Acknowledgements.

The author would like to thank A. Awad and A. Golovnev for helpful comments
during the preparation of this paper.

##  Appendix A The Classical Limit

In this section, we investigate the classical limit of the QRG field
equations, where the gravitational fields are assumed to be weak and the
relativistic effects do not contribute due to low speed approximation as in
Poisson’s equation. In order to simplify the presentation, we include only
$\mathfrak{R}^{2}$ term in equation (7), i.e.

$\nabla_{\alpha}\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}=-\frac{\epsilon\ell^{2}}{\kappa}\,\partial_{\beta}\mathfrak{R}^{2}.$
(63)

However, the conclusion will be the same for the full form. Then, equation
(12) reads

$\widetilde{\mathfrak{T}}{{}^{\alpha}}{{}_{\beta}}=\mathfrak{T}{{}^{\alpha}}{{}_{\beta}}+\frac{\epsilon\ell^{2}}{\kappa}\delta^{\alpha}_{\beta}\mathfrak{R}^{2}.$
(64)

Using $\mathfrak{R}=-\kappa\widetilde{\mathfrak{T}}$, the contraction of the
above gives

$\widetilde{\mathfrak{T}}={1\pm\sqrt{1-16\epsilon\kappa\ell^{2}\mathfrak{T}}\over
8\epsilon\kappa\ell^{2}}.$ (65)

At the linear order, we have
$\sqrt{1-16\epsilon\kappa\ell^{2}\mathfrak{T}}\approx
1-8\epsilon\kappa\ell^{2}\mathfrak{T}$, where
$|16\epsilon\kappa\ell^{2}\mathfrak{T}|\ll 1$. This proves that
$\widetilde{\mathfrak{T}}\approx\mathfrak{T}$ or
$\widetilde{\mathfrak{T}}\approx\frac{1}{4\epsilon\ell^{2}\kappa}-\mathfrak{T}$
at the weak field limit. Similar to the full form in Sec. II we drop the
latter solution which allows vacuum solution
$\mathfrak{R}=-\frac{1}{4\epsilon\ell^{2}}$ in the present case, since
$\epsilon=0$ is not valid. Consequently, for the first solution, Ricci tensor
in (14) follows the relation

$\mathfrak{R}_{\alpha\beta}=\kappa\left(\mathfrak{T}_{\alpha\beta}-\frac{1}{2}g_{\alpha\beta}\mathfrak{T}\right),$
(66)

which is identical to the GR theory. For the $tt$-component, we write

$\mathfrak{R}_{tt}=\nabla^{2}\psi=\frac{\kappa}{2}\rho,$ (67)

where $\psi$ denotes the conventional Newtonian potential. The comparison with
Poisson equation, $\nabla^{2}\psi=\frac{4\pi G}{c^{4}}\rho$, allows us to
obtain the QRG coupling constant

$\kappa=\frac{8\pi G}{c^{4}}=\kappa_{E},$ (68)

identical to GR. Similarly, one would obtain the same result if the full QRG,
namely equation (7), is applied. In conclusion, QRG does not modify the
coupling constant unlike Rastall gravity, where Rastall’s correction
($\propto\mathfrak{R}$) spoils up Einstein tensor and finally leads to a
modified coupling constant (6).

##  Appendix B More pulsars’ observational data

In this section, we extend our study to include more pulsars’ observational
data. This is to test the viability of the QRG scenario with larger dataset of
different observational approaches. In Table 3, we give a list of four
categories of pulsars according to their types and observational technique. As
noted before that the radii measurements are not an easy task, so some of the
presented values of radii are in fact model dependent in way or another. We
generated this table, for $\epsilon=-0.01$, by imposing the conformal
constraint on the sound speed to derive the corresponding pulsar radius (if
not observed) given its mass at 68% CL. Those predicted values will be labeled
by “p” letter in Table 3. Consequently, we list the corresponding compactness
and the KB parameters (36). In addition, we provide the estimated mean mass as
given by equation (29) to ensure its agreement with the observed value within
68% CL. We note that since $C<C_{\text{max}}=0.752$, the pulsar is stable and
verifies all physical conditions.

Table 3: The KB model parameters (with $\epsilon=-0.01$) of different types of observed pulsars which generate mass and compactness in agreement with the observed mass and radius values within 1$\sigma$ CL. Pulsar | Ref | observed mass | radius | estimated mass | $C$ | $a_{0}$ | $a_{1}$ | $a_{2}$
---|---|---|---|---|---|---|---|---
|  |  ($M_{\odot}$) | [km] |  ($M_{\odot}$) |  |  |  |
High-mass X-ray binaries
Her X–1 | [79] | $0.85\pm 0.15$ | $7.56\pm 1.33^{\text{p}}$ | $0.82$ | $0.29$ | $0.21$ | $-0.56$ | $0.35$
4U 1538–52 | [80] | $0.87\pm 0.07$ | $7.74\pm 0.62^{\text{p}}$ | $0.84$ | $0.30$ | $0.21$ | $-0.57$ | $0.36$
SMC X–1 | [80] | $1.04\pm 0.09$ | $8.34\pm 0.72^{\text{p}}$ | $1.00$ | $0.34$ | $0.26$ | $-0.67$ | $0.41$
LMC X–4 | [80] | $1.29\pm 0.05$ | $9.10\pm 0.35^{\text{p}}$ | $1.27$ | $0.40$ | $0.34$ | $-0.85$ | $0.51$
Cen X–3 | [80] | $1.49\pm 0.08$ | $9.48\pm 0.51^{\text{p}}$ | $1.47$ | $0.44$ | $0.40$ | $-0.99$ | $0.59$
Vela X–1 | [80] | $1.77\pm 0.08$ | $9.97\pm 0.45^{\text{p}}$ | $1.75$ | $0.50$ | $0.51$ | $-1.20$ | $0.70$
Low-mass X-ray binaries (thermonuclear bursts/quiescence)
4U 1608-52 | [81] | $1.57^{+0.30}_{-0.29}$ | $9.99\pm 1.91^{\text{p}}$ | $1.58$ | $0.46$ | $0.44$ | $-1.06$ | $0.62$
KS 1731-260 | [81] | $1.61^{+0.35}_{-0.37}$ | $10.06\pm 2.31^{\text{p}}$ | $1.62$ | $0.46$ | $0.45$ | $-1.09$ | $0.64$
EXO 1745-248 | [81] | $1.65^{+0.21}_{-0.31}$ | $10.50\pm 1.59^{\text{p}}$ | $1.66$ | $0.46$ | $0.44$ | $-1.06$ | $0.62$
4U 1820-30 | [81] | $1.77^{+0.25}_{-0.28}$ | $11.26\pm 1.59^{\text{p}}$ | $1.78$ | $0.49$ | $0.44$ | $-1.06$ | $0.62$
4U 1724-207 | [81] | $1.81^{+0.25}_{-0.37}$ | $11.52\pm 1.72^{\text{p}}$ | $1.82$ | $0.46$ | $0.44$ | $-1.06$ | $0.62$
SAX J1748.9-2021 | [81] | $1.81^{+0.25}_{-0.37}$ | $11.52\pm 1.91^{\text{p}}$ | $1.82$ | $0.46$ | $0.44$ | $-1.06$ | $0.62$
M13 | [82] | $1.38^{+0.08}_{-0.23}$ | $9.95^{+0.24}_{-0.27}$ | $1.39$ | $0.41$ | $0.35$ | $-0.88$ | $0.53$
X7 | [83] | $1.4$ | $11.1^{+0.8}_{-0.7}$ | $1.41$ | $0.37$ | $0.30$ | $-0.76$ | $0.47$
Millisecond Pulsars
PSR J0030+0451 | [40] | $1.34\pm 0.16$ | $12.71\pm 1.19$ | $1.490$ | $0.381$ | $0.309$ | $-0.788$ | $0.479$
| [39] | $1.44\pm 0.16$ | $13.02\pm 1.24$ | $1.591$ | $0.397$ | $0.331$ | $-0.838$ | $0.506$
| [43] | $1.40^{+0.13}_{-0.12}$ | $11.71^{+0.88}_{-0.83}$ | $1.494$ | $0.395$ | $0.329$ | $-0.832$ | $0.503$
PSR J0437-4715 | [76, 71] | $1.44\pm 0.07$ | $13.6\pm 0.9$ | $1.508$ | $0.350$ | $0.270$ | $-0.700$ | $0.430$
PSR J1614-2230 | [84] | $1.908\pm 0.016$ | $12\pm 1$ | $1.922$ | $0.470$ | $0.477$ | $-1.081$ | $0.634$
PSR J0348+0432 | [72] | $2.01\pm 0.04$ | $12.424\pm 0.517$ | $2.043$ | $0.490$ | $0.485$ | $-1.158$ | $0.673$
PSR J0740+6620 | [41] | $2.07\pm 0.11$ | $14.3\pm 3.5$ | $2.084$ | $0.515$ | $0.536$ | $-1.259$ | $0.723$
| [54] | $2.07\pm 0.11$ | $12.34\pm 1.75$ | $2.088$ | $0.501$ | $0.508$ | $-1.203$ | $0.696$
Gravitational-wave Signals
GW170817-1 | [47] | $1.45\pm 0.09$ | $11.9\pm 1.4$ | $1.502$ | $0.394$ | $0.327$ | $-0.829$ | $0.501$
GW170817-2 | [47] | $1.27\pm 0.09$ | $11.9\pm 1.4$ | $1.320$ | $0.360$ | $0.283$ | $-0.730$ | $0.447$

The four categories in Table 3 are: (i) High-mass X-ray binaries category, it
includes six eclipsing X-ray pulsars in binary systems with high mass
companions, where the X-rays from the pulsar in orbital motion are partly
blocked by the companions. We note that there are two different mass
measurements of the pulsar Her X–1, $M=0.85\pm 0.15M_{\odot}$ using radial
velocity curve [79] and $M=1.07\pm 0.36M_{\odot}$ using X-rays [80]. In the
present study, we take the former mass value which has been obtained by a well
justified technique with more precise measurement, and the availability of the
corresponding radius measurement $R=8.1\pm 0.41$ km using Roche lobe geometry
[85]. In general, the predicted radii are at 68% agreement with [85]. (ii)
Low-mass X-ray binaries category, it includes eight pulsars, mass and radius
can be simultaneously measured by spectroscopic analyzing the spectral lines
obtained during thermonuclear X-ray bursts or in quiescence [81]. Imposing
conformal sound speed bound leads to predicted radii are at 68% agreement with
[86]. For M13 and X7, which are quiescence low-mass X-ray binaries, we used
the observed radii as given by [82, 83]. In practice, (i) and (ii) categories
give less precise measurements of masses in comparison to ToAs analysis which
have been used in the third category. (iii) MSPs category includes five
pulsars, in this category some radii are already measured by NICER and some
are on going. We include the two independent radius measurements of the pulsar
PSR J0030+0451 as obtained by [40] and [39] in addition to its latest updated
mass-radius analysis of NICER data [43]. (iv) The gravitational-wave signals
category includes two signals labeled as GW170817-1 and GW170817-2 as obtained
by LIGO and Virgo collaboration.

Table 4: Estimated values of density, sound speed, trace anomaly, redshift and the bag constant. We note that $B_{60}$ is the value of the bag constant normalized to $B=60$ MeV/fm3, where the allowed range is $0.982<B_{60}<1.53$. Pulsar |  $\rho$ [$\rho_{\text{nuc}}$] | $v_{r}^{2}/c^{2}$ | $\Delta(r)$ | Z(r) | $B_{60}$
---|---|---|---|---|---
| center | surface | center | surface | center | surface | center | surface |
High-mass X-ray binaries
Her X–1 | $3.07$ | $2.31$ | $0.26$ | $0.24$ | $0.27$ | $0.31$ | $0.32$ | $0.19$ | $1.47$
4U 1538–52 | $3.13$ | $2.33$ | $0.26$ | $0.24$ | $0.27$ | $0.31$ | $0.33$ | $0.19$ | $1.49$
SMC X–1 | $3.26$ | $2.32$ | $0.27$ | $0.25$ | $0.26$ | $0.30$ | $0.40$ | $0.23$ | $1.49$
LMC X–4 | $3.54$ | $2.32$ | $0.29$ | $0.26$ | $0.24$ | $0.30$ | $0.53$ | $0.29$ | $1.50$
Cen X–3 | $3.77$ | $2.33$ | $0.30$ | $0.27$ | $0.22$ | $0.29$ | $0.64$ | $0.34$ | $1.51$
Vela X–1 | $4.05$ | $2.29$ | $0.32$ | $0.29$ | $0.20$ | $0.27$ | $0.83$ | $0.42$ | $1.49$
Low-mass X-ray binaries (quiescence/thermonuclear bursts)
4U 1608-52 | $3.77$ | $2.26$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.70$ | $0.37$ | $1.47$
KS 1731-260 | $3.81$ | $2.26$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.73$ | $0.38$ | $1.47$
EXO 1785-248 | $3.41$ | $2.05$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.70$ | $0.37$ | $1.33$
4U 1820-30 | $2.97$ | $1.78$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.70$ | $0.37$ | $1.15$
4U 1724-207 | $2.84$ | $1.70$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.70$ | $0.37$ | $1.10$
SAX J1748.9-2021 | $2.84$ | $1.70$ | $0.31$ | $0.28$ | $0.21$ | $0.28$ | $0.70$ | $0.37$ | $1.10$
M13 | $3.20$ | $2.07$ | $0.29$ | $0.26$ | $0.24$ | $0.29$ | $0.55$ | $0.30$ | $1.34$
X7 | $2.27$ | $1.55$ | $0.28$ | $0.25$ | $0.25$ | $0.30$ | $0.47$ | $0.26$ | $0.99$
Millisecond Pulsars
PSR J0030+0451 | $2.165$ | $1.459$ | $0.283$ | $0.254$ | $0.245$ | $0.298$ | $0.483$ | $0.271$ | $0.921$
| $2.190$ | $1.444$ | $0.288$ | $0.258$ | $0.239$ | $0.296$ | $0.520$ | $0.288$ | $0.911$
| $2.447$ | $1.618$ | $0.288$ | $0.257$ | $0.240$ | $0.296$ | $0.516$ | $0.286$ | $1.021$
PSR J0437-4715 | $1.597$ | $1.121$ | $0.275$ | $0.248$ | $0.255$ | $0.303$ | $0.419$ | $0.240$ | $0.707$
PSR J1614-2230 | $2.658$ | $1.580$ | $0.312$ | $0.277$ | $0.212$ | $0.282$ | $0.717$ | $0.373$ | $0.997$
PSR J0348+0432 | $2.726$ | $1.571$ | $0.319$ | $0.283$ | $0.204$ | $0.277$ | $0.784$ | $0.400$ | $0.991$
PSR J0740+6620 | $3.119$ | $1.727$ | $0.329$ | $0.291$ | $0.192$ | $0.270$ | $0.877$ | $0.435$ | $1.090$
| $2.828$ | $1.600$ | $0.324$ | $0.286$ | $0.199$ | $0.274$ | $0.825$ | $0.416$ | $1.010$
Gravitational-wave Signals
GW170817-1 | $2.399$ | $1.588$ | $0.287$ | $0.257$ | $0.240$ | $0.296$ | $0.514$ | $0.285$ | $1.002$
GW170817-2 | $2.304$ | $1.595$ | $0.278$ | $0.250$ | $0.251$ | $0.301$ | $0.440$ | $0.250$ | $1.006$

Next, in Table 4, we provide corresponding numerical values of some important
physical quantities which identify the physical properties of the pulsar as
estimated by QRG model. Obviously, the surface densities are above
$\rho_{s}>\rho_{\text{nuc}}=2.7\times 10^{14}$ g/cm3 in all cases which is
required for a stable quark star, whereas the core densities are few times
nuclear saturation density, (2–4)$\rho_{\text{nuc}}$. Remarkably, the
conjectured conformal upper limit on the sound speed is fulfilled for all
pulsars even for those which are derived from direct observations. This
confirms the role of QRG to obtain this constraint satisfied. This is unlike
GR which breaks the conformal constraint in some cases using the same KB
matric ansatz [53]. Another interesting feature is the positiveness of the
conformal trace anomaly $\Delta$ everywhere inside the pulsar which confirms
that the conformal symmetry is broken. On the other hand, it indicates that
the stellar matter fulfills the TEC $\mathfrak{T}<0$ everywhere inside the
pulsar. The redshift $Z$ is finite although the fluid is assumed to be
anisotropic. Finally, we provide the corresponding bag constant
$B_{60}=\frac{B}{60~{}\text{MeV/fm}^{3}}$ where $B=\frac{1}{4}\rho_{s}c^{2}$
and $0.982<B_{60}<1.53$ for physical range. If the latter constraint is
fulfilled the strange quark model would be valid in this case. Otherwise, it
is not unless some systematic errors are obtained in the data.

According to our calculations as summarized in Table 4, the strange quark star
model is in tension with NICER measurements of the PSR J0030+0451 radius
$R=12.71\pm 1.19$ [40] and $R=13.02\pm 1.24$ [39], since the corresponding bag
constants are respectively $B_{60}=0.921$ and $B_{60}=0.911$ below the allowed
minimal value. However, the model is in a perfect agreement with the latest
analysis of NICER data of the pulsar PSR J0030+0451 [43] which gives
$R=11.71^{+0.88}_{-0.83}$ relatively smaller than the previous measurements,
and therefore allows for a relatively higher bag constant $B_{60}=1.021$ in
the allowed physical range. Also, in the case of PSR J0437-4715, we obtain
much lower bag constant value $B_{60}=0.707$, if we consider the relatively
much large radius measurement $R=13.6\pm 0.9$ km [71]. Therefore, confirmation
of the updated NICER analysis of the pulsar PSR J0030+0451 [43], would
validate the QRG scenario of strange quark star model. However, confirmation
of the measured radius of PSR J0437-4715 [71] would indicate a significant
tension with the QRG scenario of strange quark star model.

## References

  * Rastall [1972] P. Rastall, Generalization of the einstein theory, Phys. Rev. D 6, 3357 (1972).
  * Rastall [1976] P. Rastall, A Theory of Gravity, Can. J. Phys. 54, 66 (1976).
  * El Hanafy [2022] W. El Hanafy, Impact of Rastall Gravity on Mass, Radius, and Sound Speed of the Pulsar PSR J0740+6620, Astrophys. J. 940, 51 (2022), arXiv:2209.10287 [astro-ph.HE] .
  * El Hanafy and Awad [2023] W. El Hanafy and A. Awad, Implications of the Conformal Constraint on Sound Speed on the Radius of PSR J0952–0607 within Rastall Gravity, Astrophys. J. 951, 144 (2023), arXiv:2305.14953 [astro-ph.HE] .
  * Visser [2018] M. Visser, Rastall gravity is equivalent to Einstein gravity, Phys. Lett. B 782, 83 (2018), arXiv:1711.11500 [gr-qc] .
  * Darabi  _et al._ [2018] F. Darabi, H. Moradpour, I. Licata, Y. Heydarzade, and C. Corda, Einstein and Rastall Theories of Gravitation in Comparison, Eur. Phys. J. C 78, 25 (2018), arXiv:1712.09307 [gr-qc] .
  * Lin and Qian [2020] K. Lin and W.-L. Qian, Cosmic evolution of dark energy in a generalized Rastall gravity, Eur. Phys. J. C 80, 561 (2020), arXiv:2006.03229 [gr-qc] .
  * Li  _et al._ [2019] R. Li, J. Wang, Z. Xu, and X. Guo, Constraining the Rastall parameters in static space–times with galaxy-scale strong gravitational lensing, Mon. Not. Roy. Astron. Soc. 486, 2407 (2019), arXiv:1903.08790 .
  * Ziaie  _et al._ [2020] A. H. Ziaie, H. Moradpour, and H. Shabani, Structure formation in generalized Rastall gravity, European Physical Journal Plus 135, 916 (2020), arXiv:2002.12757 [physics.gen-ph] .
  * Oliveira  _et al._ [2015] A. M. Oliveira, H. E. S. Velten, J. C. Fabris, and L. Casarini, Neutron Stars in Rastall Gravity, Phys. Rev. D 92, 044020 (2015), arXiv:1506.00567 [gr-qc] .
  * Moradpour and Salako [2016] H. Moradpour and I. G. Salako, Thermodynamic analysis of the static spherically symmetric field equations in Rastall theory, Adv. High Energy Phys. 2016, 3492796 (2016), arXiv:1606.06589 [gr-qc] .
  * Bamba  _et al._ [2018] K. Bamba, A. Jawad, S. Rafique, and H. Moradpour, Thermodynamics in Rastall Gravity with Entropy Corrections, Eur. Phys. J. C 78, 986 (2018), arXiv:1812.00753 [physics.gen-ph] .
  * Cruz  _et al._ [2019] M. Cruz, S. Lepe, and G. Morales-Navarrete, A thermodynamics revision of Rastall gravity, Class. Quant. Grav. 36, 225007 (2019), arXiv:1904.11945 [gr-qc] .
  * Fabris  _et al._ [2023] J. C. Fabris, O. F. Piattella, and D. C. Rodrigues, On Rastall gravity formulation as a $f(R,\mathcal{L}_{m})$ and a f(R, T) theory, Eur. Phys. J. Plus 138, 232 (2023), arXiv:2011.10503 [gr-qc] .
  * da Silva  _et al._ [2021] F. M. da Silva, L. C. N. Santos, and C. C. Barros, Rapidly rotating compact stars in Rastall’s gravity, Class. Quant. Grav. 38, 165011 (2021), arXiv:2010.00086 [astro-ph.HE] .
  * Li  _et al._ [2022] Q. Li, C. Ma, Y. Zhang, Z.-W. Lin, and P.-F. Duan, Shadow, absorption and Hawking radiation of a Schwarzschild black hole surrounded by a cloud of strings in Rastall gravity, Eur. Phys. J. C 82, 658 (2022), arXiv:2307.04144 [gr-qc] .
  * Nashed and El Hanafy [2022] G. G. L. Nashed and W. El Hanafy, Non-trivial class of anisotropic compact stellar model in Rastall gravity, Eur. Phys. J. C 82, 679 (2022), arXiv:2208.13814 [gr-qc] .
  * Duff [1994] M. J. Duff, Twenty years of the Weyl anomaly, Class. Quant. Grav. 11, 1387 (1994), arXiv:hep-th/9308075 .
  * Fujimoto  _et al._ [2022] Y. Fujimoto, K. Fukushima, L. D. McLerran, and M. Praszalowicz, Trace Anomaly as Signature of Conformality in Neutron Stars, Phys. Rev. Lett. 129, 252702 (2022), arXiv:2207.06753 [nucl-th] .
  * Goenner [1976] H. Goenner, A metric theory of gravitation without minimal coupling of matter and gravitational field, Zeitschrift für Naturforschung A 31, 1451 (1976).
  * Palmer and Anderson [1974] R. G. Palmer and P. W. Anderson, Corresponding-states approach to nuclear and neutron-star matter, Phys. Rev. D 9, 3281 (1974).
  * Sawyer [1972] R. F. Sawyer, Condensed pi- phase in neutron star matter, Phys. Rev. Lett. 29, 382 (1972).
  * Ramanan [2019] S. Ramanan, Superfluidity in neutron stars, DAE Symp. Nucl. Phys. 64, 45 (2019).
  * Weber  _et al._ [2007] F. Weber, R. Negreiros, P. Rosenfield, and M. Stejner, Pulsars as Astrophysical Laboratories for Nuclear and Particle Physics, Prog. Part. Nucl. Phys. 59, 94 (2007), arXiv:astro-ph/0612054 .
  * Rahmansyah  _et al._ [2020] A. Rahmansyah, A. Sulaksono, A. B. Wahidin, and A. M. Setiawan, Anisotropic neutron stars with hyperons: implication of the recent nuclear matter data and observations of neutron stars, Eur. Phys. J. C 80, 769 (2020).
  * Bowers and Liang [1974] R. L. Bowers and E. Liang, Anisotropic spheres in general relativity, The Astrophysical Journal 188, 657 (1974).
  * Herrera and Santos [1997] L. Herrera and N. O. Santos, Local anisotropy in self-gravitating systems, Physics Reports 286, 53 (1997).
  * Herrera [1992] L. Herrera, Cracking of self-gravitating compact objects, Physics Letters A 165, 206 (1992).
  * Herrera  _et al._ [2008a] L. Herrera, J. Ospino, and A. Di Prisco, All static spherically symmetric anisotropic solutions of Einstein’s equations, Phys. Rev. D 77, 027502 (2008a), arXiv:0712.0713 [gr-qc] .
  * Herrera  _et al._ [2008b] L. Herrera, N. O. Santos, and A. Wang, Shearing Expansion-free Spherical Anisotropic Fluid Evolution, Phys. Rev. D 78, 084026 (2008b), arXiv:0810.1083 [gr-qc] .
  * Herrera  _et al._ [2011] L. Herrera, A. Di Prisco, and J. Ibanez, On the Role of Electric Charge and Cosmological Constant in Structure Scalars, Phys. Rev. D 84, 107501 (2011), arXiv:1111.1071 [gr-qc] .
  * Krori and Barua [1975] K. D. Krori and J. Barua, A singularity-free solution for a charged fluid sphere in general relativity, Journal of Physics A 8, 508 (1975).
  * Bhattacharyya  _et al._ [2016] S. Bhattacharyya, I. Bombaci, D. Logoteta, and A. V. Thampan, Fast spinning strange stars: possible ways to constrain interacting quark matter parameters, Mon. Not. Roy. Astron. Soc. 457, 3101 (2016), arXiv:1601.06120 [astro-ph.HE] .
  * Annala  _et al._ [2020] E. Annala, T. Gorda, A. Kurkela, J. Nättilä, and A. Vuorinen, Evidence for quark-matter cores in massive neutron stars, Nature Phys. 16, 907 (2020), arXiv:1903.09121 [astro-ph.HE] .
  * Mészáros  _et al._ [2019] P. Mészáros, D. B. Fox, C. Hanna, and K. Murase, Multi-Messenger Astrophysics, Nature Rev. Phys. 1, 585 (2019), arXiv:1906.10212 [astro-ph.HE] .
  * Gendreau and Arzoumanian [2017] K. Gendreau and Z. Arzoumanian, Searching for a pulse, Nature Astronomy 1, 895 (2017).
  * Cromartie  _et al._ [2019] H. T. Cromartie _et al._ (NANOGrav), Relativistic Shapiro delay measurements of an extremely massive millisecond pulsar, Nature Astron. 4, 72 (2019), arXiv:1904.06759 [astro-ph.HE] .
  * Fonseca  _et al._ [2021] E. Fonseca _et al._ , Refined Mass and Geometric Measurements of the High-mass PSR J0740+6620, Astrophys. J. Lett. 915 (2021), arXiv:2104.00880 [astro-ph.HE] .
  * Miller  _et al._ [2019] M. C. Miller _et al._ , PSR J0030+0451 Mass and Radius from $NICER$ Data and Implications for the Properties of Neutron Star Matter, Astrophys. J. Lett. 887 (2019), arXiv:1912.05705 [astro-ph.HE] .
  * Raaijmakers  _et al._ [2019] G. Raaijmakers _et al._ , A $NICER$ view of PSR J0030+0451: Implications for the dense matter equation of state, Astrophys. J. Lett. 887 (2019), arXiv:1912.05703 [astro-ph.HE] .
  * Miller  _et al._ [2021] M. C. Miller _et al._ , The Radius of PSR J0740+6620 from NICER and XMM-Newton Data, Astrophys. J. Lett. 918 (2021), arXiv:2105.06979 [astro-ph.HE] .
  * Riley  _et al._ [2021] T. E. Riley _et al._ , A NICER View of the Massive Pulsar PSR J0740+6620 Informed by Radio Timing and XMM-Newton Spectroscopy, Astrophys. J. Lett. 918 (2021), arXiv:2105.06980 [astro-ph.HE] .
  * Vinciguerra  _et al._ [2024] S. Vinciguerra _et al._ , An Updated Mass–Radius Analysis of the 2017–2018 NICER Data Set of PSR J0030+0451, Astrophys. J. 961, 62 (2024), arXiv:2308.09469 [astro-ph.HE] .
  * McLerran and Reddy [2019] L. McLerran and S. Reddy, Quarkyonic Matter and Neutron Stars, Phys. Rev. Lett. 122, 122701 (2019), arXiv:1811.12503 [nucl-th] .
  * Altiparmak  _et al._ [2022] S. Altiparmak, C. Ecker, and L. Rezzolla, On the Sound Speed in Neutron Stars, Astrophys. J. Lett. 939 (2022), arXiv:2203.14974 [astro-ph.HE] .
  * Abbott  _et al._ [2017] B. P. Abbott _et al._ (LIGO Scientific, Virgo), GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral, Phys. Rev. Lett. 119, 161101 (2017), arXiv:1710.05832 [gr-qc] .
  * Abbott  _et al._ [2018] B. P. Abbott _et al._ (LIGO Scientific, Virgo), GW170817: Measurements of neutron star radii and equation of state, Phys. Rev. Lett. 121, 161101 (2018), arXiv:1805.11581 [gr-qc] .
  * Abbott  _et al._ [2020] B. P. Abbott _et al._ (LIGO Scientific, Virgo), GW190425: Observation of a Compact Binary Coalescence with Total Mass $\sim 3.4M_{\odot}$, Astrophys. J. Lett. 892 (2020), arXiv:2001.01761 [astro-ph.HE] .
  * Doroshenko  _et al._ [2022] V. Doroshenko, V. Suleimanov, G. Pühlhofer, and A. Santangelo, A strangely light neutron star within a supernova remnant, Nature Astronomy 6, 1444 (2022).
  * Horvath  _et al._ [2023] J. E. Horvath, L. S. Rocha, L. M. de Sá, P. H. R. S. Moraes, L. G. Barão, M. G. B. de Avellar, A. Bernardo, and R. R. A. Bachega, A light strange star in the remnant HESS J1731$-$347: Minimal consistency checks, Astron. Astrophys. 672 (2023), arXiv:2303.10264 [astro-ph.HE] .
  * Moradpour  _et al._ [2017a] H. Moradpour, A. Bonilla, E. M. C. Abreu, and J. A. Neto, Accelerated cosmos in a nonextensive setup, Phys. Rev. D 96, 123504 (2017a), arXiv:1711.08338 [physics.gen-ph] .
  * Moradpour  _et al._ [2017b] H. Moradpour, N. Sadeghnezhad, and S. H. Hendi, Traversable asymptotically flat wormholes in Rastall gravity, Can. J. Phys. 95, 1257 (2017b), arXiv:1606.00846 [gr-qc] .
  * Roupas and Nashed [2020] Z. Roupas and G. G. L. Nashed, Anisotropic neutron stars modelling: constraints in Krori–Barua spacetime, Eur. Phys. J. C 80, 905 (2020), arXiv:2007.09797 [gr-qc] .
  * Legred  _et al._ [2021] I. Legred, K. Chatziioannou, R. Essick, S. Han, and P. Landry, Impact of the PSR J0740+6620 radius constraint on the properties of high-density matter, Phys. Rev. D 104, 063003 (2021), arXiv:2106.05313 [astro-ph.HE] .
  * Landry  _et al._ [2020] P. Landry, R. Essick, and K. Chatziioannou, Nonparametric constraints on neutron star matter with existing and upcoming gravitational wave and pulsar observations, Phys. Rev. D 101, 123007 (2020), arXiv:2003.04880 [astro-ph.HE] .
  * Ecker and Rezzolla [2022] C. Ecker and L. Rezzolla, Impact of large-mass constraints on the properties of neutron stars, Mon. Not. Roy. Astron. Soc. 519, 2615 (2022), arXiv:2209.08101 [astro-ph.HE] .
  * Heintzmann and Hillebrandt [1975] H. Heintzmann and W. Hillebrandt, Neutron stars with an anisotropic equation of state: mass, redshift and stability, Astronomy and Astrophysics 38, 51 (1975).
  * Chan  _et al._ [1993] R. Chan, L. Herrera, and N. O. Santos, Dynamical instability for radiating anisotropic collapse, Monthly Notices of the Royal Astronomical Society 265, 533 (1993), http://oup.prod.sis.lan/mnras/article-pdf/265/3/533/3807712/mnras265-0533.pdf .
  * Awad [2016] A. Awad, Weyl Anomaly and Initial Singularity Crossing, Phys. Rev. D 93, 084006 (2016), arXiv:1512.06405 [hep-th] .
  * Podkowka  _et al._ [2018] D. M. Podkowka, R. F. P. Mendes, and E. Poisson, Trace of the energy-momentum tensor and macroscopic properties of neutron stars, Phys. Rev. D 98, 064057 (2018), arXiv:1807.01565 [gr-qc] .
  * Zel’dovich [1961] Y. B. Zel’dovich, Equation of state at ultra-high densities and its relativistic limitations, Zhur. Eksptl’. i Teoret. Fiz. 41 (1961), [Engl. transl.: 1962, Sov. Phys. JETP 14, 1143-1147 (1962)].
  * Witten [1984] E. Witten, Cosmic Separation of Phases, Phys. Rev. D 30, 272 (1984).
  * Farhi and Jaffe [1984] E. Farhi and R. L. Jaffe, Strange Matter, Phys. Rev. D 30, 2379 (1984).
  * Olinto [1987] A. V. Olinto, On the conversion of neutron stars into strange stars, Physics Letters B 192, 71 (1987).
  * Doroba [1989] W. Doroba, Strange stars. Are they bare?, Acta Physica Polonica B 20, 967 (1989).
  * Heiselberg and Pethick [1993] H. Heiselberg and C. J. Pethick, Transport and relaxation in degenerate quark plasmas, Phys. Rev. D 48, 2916 (1993).
  * Horvath and Benvenuto [1988] J. E. Horvath and O. G. Benvenuto, On the stability of slow neutron combustion in astrophysical objects, Physics Letters B 213, 516 (1988).
  * Drago  _et al._ [2007] A. Drago, A. Lavagno, and I. Parenti, Burning of an hadronic star into a quark or a hybrid star, Astrophys. J. 659, 1519 (2007), arXiv:astro-ph/0512652 .
  * Haensel  _et al._ [2007] P. Haensel, A. Y. Potekhin, and D. G. Yakovlev,  _Neutron Stars 1: Equation of State and Structure_, 1st ed., Astrophysics and Space Science Library, Vol. 326 (Springer, 2007).
  * Romani  _et al._ [2022] R. W. Romani, D. Kandel, A. V. Filippenko, T. G. Brink, and W. Zheng, PSR J0952$-$0607: The Fastest and Heaviest Known Galactic Neutron Star, Astrophys. J. Lett. 934 (2022), arXiv:2207.05124 [astro-ph.HE] .
  * Gonzalez-Caniulef  _et al._ [2019] D. Gonzalez-Caniulef, S. Guillot, and A. Reisenegger, Neutron star radius measurement from the ultraviolet and soft X-ray thermal emission of PSR J0437$-$4715, Mon. Not. Roy. Astron. Soc. 490, 5848 (2019), arXiv:1904.12114 [astro-ph.HE] .
  * Antoniadis  _et al._ [2013] J. Antoniadis _et al._ , A Massive Pulsar in a Compact Relativistic Binary, Science 340, 6131 (2013), arXiv:1304.6875 [astro-ph.HE] .
  * Demorest  _et al._ [2010] P. Demorest, T. Pennucci, S. Ransom, M. Roberts, and J. Hessels, Shapiro Delay Measurement of A Two Solar Mass Neutron Star, Nature 467, 1081 (2010), arXiv:1010.5788 [astro-ph.HE] .
  * Fonseca  _et al._ [2016] E. Fonseca _et al._ , The NANOGrav Nine-year Data Set: Mass and Geometric Measurements of Binary Millisecond Pulsars, Astrophys. J. 832, 167 (2016), arXiv:1603.00545 [astro-ph.HE] .
  * Arzoumanian  _et al._ [2018a] Z. Arzoumanian _et al._ (NANOGRAV), The NANOGrav 11-year Data Set: Pulsar-timing Constraints On The Stochastic Gravitational-wave Background, Astrophys. J. 859, 47 (2018a), arXiv:1801.02617 [astro-ph.HE] .
  * Reardon  _et al._ [2016] D. J. Reardon _et al._ , Timing analysis for 20 millisecond pulsars in the Parkes Pulsar Timing Array, Mon. Not. Roy. Astron. Soc. 455, 1751 (2016), arXiv:1510.04434 [astro-ph.HE] .
  * Zdunik  _et al._ [2001] J. L. Zdunik, P. Haensel, and E. Gourgoulhon, The Crust of rotating strange quark stars, Astron. Astrophys. 372, 535 (2001), arXiv:astro-ph/0104116 .
  * Baym  _et al._ [1971] G. Baym, C. Pethick, and P. Sutherland, The Ground State of Matter at High Densities: Equation of State and Stellar Models, Astrophys. J. 170, 299 (1971).
  * Abubekerov  _et al._ [2008] M. K. Abubekerov, E. A. Antokhina, A. M. Cherepashchuk, and V. V. Shimanskii, The mass of the compact object in the x-ray binary her x-1/HZ her, Astronomy Reports 52, 379 (2008).
  * Rawls  _et al._ [2011] M. L. Rawls, J. A. Orosz, J. E. McClintock, M. A. P. Torres, C. D. Bailyn, and M. M. Buxton, Refined Neutron-Star Mass Determinations for Six Eclipsing X-Ray Pulsar Binaries, Astrophys. J. 730, 25 (2011), arXiv:1101.2465 [astro-ph.SR] .
  * Özel and Freire [2016] F. Özel and P. Freire, Masses, Radii, and the Equation of State of Neutron Stars, Ann. Rev. Astron. Astrophys. 54, 401 (2016), arXiv:1603.02698 [astro-ph.HE] .
  * Webb and Barret [2007] N. A. Webb and D. Barret, Constraining the equation of state of supra-nuclear dense matter from XMM-Newton observations of neutron stars in globular clusters, Astrophys. J. 671, 727 (2007), arXiv:0708.3816 [astro-ph] .
  * Bogdanov  _et al._ [2016] S. Bogdanov, C. O. Heinke, F. Özel, and T. Güver, Neutron Star Mass-Radius Constraints of the Quiescent Low-mass X-ray Binaries X7 and X5 in the Globular Cluster 47 Tuc, Astrophys. J. 831, 184 (2016), arXiv:1603.01630 [astro-ph.HE] .
  * Arzoumanian  _et al._ [2018b] Z. Arzoumanian _et al._ (NANOGrav), The NANOGrav 11-year Data Set: High-precision timing of 45 Millisecond Pulsars, Astrophys. J. Suppl. 235, 37 (2018b), arXiv:1801.01837 [astro-ph.HE] .
  * Gangopadhyay  _et al._ [2013] 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. Roy. Astron. Soc. 431, 3216 (2013), arXiv:1303.1956 [astro-ph.HE] .
  * Özel  _et al._ [2016] F. Özel, D. Psaltis, T. Guver, G. Baym, C. Heinke, and S. Guillot, The Dense Matter Equation of State from Neutron Star Radius and Mass Measurements, Astrophys. J. 820, 28 (2016), arXiv:1505.05155 [astro-ph.HE] .