The exterior spacetime of relativistic stars in scalar-Gauss-Bonnet gravity

TL;DR

This paper derives an analytical exterior spacetime solution for neutron stars in scalar-Gauss-Bonnet gravity, enabling future observational tests of deviations from general relativity in strong gravitational fields.

Contribution

It provides a novel, equation-of-state independent metric for neutron stars in scalar-Gauss-Bonnet gravity, facilitating strong-field tests of gravity theories.

Findings

Deviations from general relativity are small and higher order.

The derived metric is independent of the neutron star's equation of state.

Potential for future observational constraints on scalar-Gauss-Bonnet gravity.

Abstract

The spacetime around compact objects is an excellent place to study gravity in the strong, nonlinear, dynamical regime where solar system tests cannot account for the effects of large curvature. Understanding the dynamics of this spacetime is important for testing theories of gravity and probing a regime which has not yet been studied with observations. In this paper, we construct an analytical solution for the exterior spacetime of a neutron star in scalar-Gauss-Bonnet gravity that is independent of the equation of state chosen. The aim is to provide a metric that can be used to probe the strong-field regime near a neutron star and create predictions that can be compared with future observations to place possible constraints on the theory. In addition to constructing the metric, we examine a number of physical systems in order to see what deviations exist between our spacetime and that of general relativity. We find these deviations to be small and of higher post-Newtonian order than previous results using black hole solutions. The metric derived here can be used to further the study of scalar-Gauss-Bonnet gravity in the strong field, and allow for constraints on corrections to general relativity with future observations.

Figures

Equations

$$S=S_{\mathrm{EH}}+S_{\phi }+S_{\mathrm{GB}}+S_{\mathrm{m}}$$

$$S_{\mathrm{EH}}\equiv \kappa \int {\mathrm{d}}^{4}x\sqrt{-g}R,$$

$$S_{\phi }\equiv -\frac{1}{2}\int {\mathrm{d}}^{4}x\sqrt{-g}\left[{\nabla }_a\phi {\nabla }^a\phi +2U\left(\phi \right)\right].$$

$$\mathcal{G}\equiv R^2-4R_{ab}{R}^{ab}+R_{abcd}{R}^{abcd}$$

$$S_{\mathrm{GB}}\equiv \int {\mathrm{d}}^{4}x\sqrt{-g}\alpha f(\phi )\mathcal{G},$$

$$S_{\mathrm{GB}}\equiv \int {\mathrm{d}}^{4}x\sqrt{-g}\alpha \phi \mathcal{G}.$$

$$\bar{\alpha }=\frac{\alpha }{\kappa M_{\odot }^2},$$

$$G_{ab}$$

$$\square \phi -U_{,\phi }$$

$${\mathcal{K}}_{ab}$$

$$-4g_{ab}{R}^{cd}{\nabla }_c{\nabla }_d\phi +4R_{acbd}{\nabla }^c{\nabla }^d\phi ,$$

$$T_{ab}^{\phi }={\nabla }_a\phi {\nabla }_b\phi -\frac{1}{2}{g}_{ab}\left[{\nabla }_c\phi {\nabla }^c\phi -2U(\phi )\right].$$

$$T_{\mathrm{m}}^{ab}=\left(\epsilon +p\right)u^a{u}^b+p{g}^{ab},$$

$${\nabla }_a{T}_{\mathrm{m}}^{ab}=0,$$

$$\mathrm{d}{s}^2=-e^{2\tau }\mathrm{d}{t}^2+e^{2\sigma }\mathrm{d}{r}^2+r^2\mathrm{d}{\Omega }^2,$$

$$\vec{z}(r)=\sum _{n=0}^N{\vec{z}}_n(r),$$

$$G_{ab}^0=\frac{1}{2\kappa }\left[\left({\epsilon }_0+p_0\right)u_a{u}_b+p_0{g}_{ab}\right].$$

$$m_0^{\prime }$$

$${\tau }_0^{\prime }$$

$$p_0^{\prime }=\frac{\left({\epsilon }_0+p_0\right)\left(4\pi p_0{r}^3+m_0\right)}{r\left(2m_0-r\right)},$$

$${\square }_0{\phi }_1=-\alpha {\mathcal{G}}_0,$$

$$\frac{e^{-2{\sigma }_0}}{r}\left[{\phi }_1^{\prime \prime }r+{\phi }_1^{\prime }\left({\tau }_0^{\prime }r-{\sigma }_0^{\prime }r+2\right)\right]=-\alpha \mathcal{G}$$

$$\left\{-64\pi \alpha {\phi }_1^{\prime \prime }+\left[64\pi \alpha {\phi }_1^{\prime }+2r\left(2{\tau }_2-2{\sigma }_2+1\right)\right]{\sigma }_0^{\prime }+2r{\sigma }_2^{\prime }-2{\tau }_2+2{\sigma }_2-1\right\}{e}^{-2{\sigma }_0}-192\pi \alpha {\phi }_1^{\prime }{\sigma }_0^{\prime }{e}^{-4{\sigma }_0}$$

$$+64\pi \alpha {\phi }_1^{\prime \prime }{e}^{-4{\sigma }_0}+2{\tau }_2+1=4\pi r^2\left[4{\epsilon }_0{\tau }_2+{\phi }_1^{\prime 2}{e}^{-2{\sigma }_0}+2\left({\epsilon }_0+{\epsilon }_2\right)\right]$$

$$-\left(1+2{\sigma }_2\right)e^{2{\sigma }_0}-192\pi \alpha {\phi }_1^{\prime }{\tau }_0^{\prime }{e}^{-2{\sigma }_0}+\left(64\pi \alpha {\phi }_1^{\prime }+2r\right){\tau }_0^{\prime }+2r{\tau }_2^{\prime }+1=4\pi r^2\left[\left(4p_0{\sigma }_2+2p_0+2p_2\right)e^{2{\sigma }_0}+{\phi }_1^{\prime 2}\right],$$

$$p_2^{\prime }$$

$$-8\alpha {\tau }_0^{\prime \prime }-8\alpha \left({\tau }_0^{\prime 2}-{\tau }_0^{\prime }{\sigma }_0^{\prime }\right)){\phi }_1^{\prime }{e}^{-2{\sigma }_0}$$

$$+8\alpha \left({\tau }_0^{\prime \prime }+{\tau }_0^{\prime 2}-3{\tau }_0^{\prime }{\sigma }_0^{\prime }\right){\phi }_1^{\prime }{e}^{-4{\sigma }_0}]$$

$$-\left(p_0+p_2+{\epsilon }_0+{\epsilon }_2\right){\tau }_0^{\prime }+\left(p_0+{\epsilon }_0\right){\tau }_2^{\prime }+p_0^{\prime }$$

$$e^{2{\tau }_0}=e^{-2{\sigma }_0}=1-\frac{a}{r},$$

$${\mathcal{G}}_0=\frac{12a^2}{r^6},$$

$$r\left(a-r\right){\phi }_1^{\prime \prime }+\left(a-2r\right){\phi }_1^{\prime }=\alpha \frac{12a^2}{r^4},$$

$${\phi }_1$$

$$+\frac{2\alpha }{r}\left(\frac{2}{a}+\frac{1}{r}+\frac{2a}{3r^2}\right)+c_2,$$

$${\phi }_1=-\frac{c_1}{r}-\frac{a{c}_1}{2r^2}-\frac{a^2{c}_1}{3r^3}+\mathcal{O}(r^{-4}),$$

Full text

I Introduction

Einstein’s general relativity (GR) has proved to be an exceptional theory to describe gravitational phenomena in nature. From its early success in explaining the hitherto mysterious advance of the perihelion of Mercury’s orbit around the Sun [1] to its consistency with the gravitational-wave observations of merging black hole (BH) and neutron star (NS) binaries by the LIGO/Virgo Collaboration [2], GR has passed – with flying colors – all experiments it has been confronted with.

Given the continual success of the theory, it is natural to ask: should we consider GR as the final theory of the gravitational interaction? Is it worth the effort to keep developing further tests, seeking glimpses of a more complete theory? Regarding the first question, field-theoretic considerations have shown that GR is non-renormalizable, placing a major obstacle to its quantization, and indicating that the theory must be modified in the ultraviolet regime. Indeed, a generic prediction of the low-energy limits of quantum gravity theories, such as string theory and loop quantum gravity, is that GR ought to be augmented by both additional fields and higher-order curvature scalars. Regarding the second question, GR’s firm place in our vault of fundamental physical theories implies that experimental evidence for a deviation would shake the foundations of this vault. Where should we search for signatures

Where should we search for signatures of beyond-GR phenomenology? Compact objects, NSs and BHs, provide a strong-field arena on which to put GR to the test in a regime beyond the weak fields and low velocities of our Solar System. The prototypical example is radio observations of binary pulsars, which through the detailed and careful monitoring of received pulses can reconstruct the orbital motion of relativistic binaries to stupendous precision [3, 4, 5]. Another example of tests of GR with compact objects is through the observation of electromagnetic radiation emitted by the accretion disks that surround black holes, although these tests are more challenging because of the complex astrophysics in play during such observations  [6, 7]. A final and more recent example is through the observation of the x-ray pulse profile emitted by hot spots on the surface of rapidly rotating stars [8, 9, 10, 11, 12].

All of the tests mentioned probe the exterior spacetime of compact objects in one way or another. Therefore, the construction of spacetimes close to these compact objects are required to place constraints on theories that go beyond GR. While there are many modified theories that attempt to explain anomalies between observations and the theoretical predictions of GR [13, 14], a particularly interesting one is Einstein-dilaton-Gauss-Bonnet (EdGB) gravity. This theory is interesting because it emerges in the low energy limit of heterotic string theory [15], and it agrees well with GR in the weak field region [16]. EdGB gravity modifies the Einstein-Hilbert action through the coupling of the Gauss-Bonnet invariant and a dynamical (dilaton) scalar field [17]. BH solutions in this theory have already been developed [18, 19, 20, 21], but until now, NS solutions had only been obtained numerically [22, 23, 24].

In this paper, we present the first analytical solution of the field equations in the small-dilaton expansion of EdGB gravity [i.e. of scalar-Gauss-Bonnet (sGB) gravity] that represents the exterior spacetime of non-rotating NSs, working in the small coupling approximation. These solutions depend only on the mass of the NS and the strength of the sGB coupling parameter, without any dependance on the dilaton scalar charge or any additional constants of integration. This is contrary to what has been found in other theories containing a scalar field (see. [25, 26]), where a scalar charge depends on integrals over the interior of the source. The absence of this term allows our final analytic solution to be implemented directly, without the need to integrate the interior solution numerically to find the charge term.

With the known analytical exterior metric, we study the properties of this spacetime by considering (timelike) geodesics, and derive sGB corrections to the innermost stable circular orbit (ISCO), to the (circular) orbital frequency and to the epicyclical radial frequencies of perturbed circular orbits. We also consider null geodesics and derive sGB corrections to the visible fraction of a NS hot spot as observed from spatial infinity.

A representative result is shown in Fig. 1, where we present the sGB-corrected mass-radius relation of NS for different values of the sGB coupling parameter and different equations of state (EoSs). Each point in this plane represents a NS solution of a given total gravitational mass and a given total radius, fixed through a numerical integration of a given central density that requires the metric be asymptotically flat at spatial infinity and $C^{1}$ everywhere. Observe that the largest deviations arise in the high compactness regime of the mass-radius relation, where the central densities are highest. This makes sense given that sGB gravity introduces higher curvature corrections to GR, which are bound to be largest when the compactness is as large as possible.

The structure of this paper is as follows. Section II presents the basics of sGB gravity, including its action and its field equations. Section III explains our approach to finding a NS solution in sGB gravity for the exterior spacetime, while Sec. IV focuses on the interior regime and it presents the numerical solution to the interior fields. Section V discusses some astrophysical applications that differ from the known results of GR. Finally, Sec. VI concludes and discusses how this metric may be used in the future. In the remainder of this paper, we use the $(-,+,+,+)$ metric signature and the conventions of [27], as well as units in which $c=1=G$.

II Scalar Gauss-Bonnet gravity

In this section, we present the action of sGB gravity and its field equations. We then introduce the perturbative scheme we will employ to analytically solve the field equations, and we conclude by presenting the perturbatively expanded field equations.

II.1 Action

We start by considering the action of a class of theories that contain modifications proportional to the Gauss-Bonnet invariant, whose taxonomy was described in [28]:

[TABLE]

where $S_{\rm EH}$ is the Einstein-Hilbert action given by

[TABLE]

with $\kappa\equiv(16\pi)^{-1}$, $g$ is the determinant of the metric $g_{ab}$, $R\equiv g^{ab}R_{ab}=g^{ab}\tensor{R}{{}_{a}{}_{c}{}_{b}^{c}}$ is the Ricci scalar (with $R_{ab}$ and $R_{abcd}$ being the Ricci and Riemann tensors respectively) and

[TABLE]

is the action of a canonical scalar field $\varphi$ with potential $U$. The coupling between the scalar field and the Gauss-Bonnet density

[TABLE]

is given by

[TABLE]

where $f(\varphi)$ specifies the functional form of the coupling and $\alpha$ (with dimensions of [length]2) its strength. Finally, $S_{\textrm{m}}$ is the action of matter fields minimally coupled to the metric.

The choice of $f(\varphi)$ defines the particular member in the class of Gauss-Bonnet theories. For example, EdGB gravity is defined via $f(\phi)=e^{\phi}$, with typically a massless dilaton so $U(\varphi)=0$. Other coupling function $f(\phi)$ were also introduced in [29, 30] and in the context of spontaneous black hole scalarization in Refs. [31, 32, 24, 33]. Hereafter, we expand $f(\varphi)$ in a Taylor series $f(\varphi)=f(0)+f_{,\varphi}(0)\varphi+\mathcal{O}(\varphi^{2})$ and work in the so-called decoupling limit of the theory [28]. Since $\mathcal{G}$ is a topological density, the first term in the series yields a boundary term to the action which does not contribute to the equations of motion. In the second term, $f_{,\varphi}(0)$ can be absorbed into the definition of $\alpha$ and we obtain:

[TABLE]

The action in Eq. (1) with $S_{\rm GB}$ given by Eq. (6), is sometimes called decoupled dynamical Gauss-Bonnet gravity [28] or more simply sGB gravity in this paper. This theory is invariant under constant shifts $\varphi\to\varphi+c$ when $U=0$ [34, 35, 36], and thus, it belongs to shift-symmetric Horndeski gravity [37]. In this paper, we will restrict attention to sGB gravity with $U(\varphi)=0$.

Here we here work in the small-coupling approximation in which sGB modifications are small relative to GR predictions. This approximation can be enforced by requiring that $\alpha/\ell^{2}\ll 1$, where $\ell$ is the characteristic length of our system. For isolated NSs, the characteristic length scale is $\ell=\sqrt{R^{3}/M}=R/\sqrt{\mathscr{C}}$, where $R$ is the radius of the star, $M$ its mass and $\mathscr{C}$ its compactness. This length scale suggests the introduction of the dimensionless coupling parameter111Note that this dimensionless coupling parameter is different from that chosen in other work [38], since here we normalize $\alpha$ by $M_{\odot}$ instead of $M$.

[TABLE]

in terms of which the small coupling approximation reduces to $\bar{\alpha}\ll\kappa^{-1}(R/M_{\odot})^{2}\mathscr{C}^{-1}$. For NSs with $R\sim 11$ km and $\mathscr{C}\sim 0.2$, the small coupling approximation then requires that $\alpha\ll 600\;{\rm{km}}^{2}$ or equivalently $\bar{\alpha}\ll 1.5\times 10^{4}$. The small-coupling approximation is well-justified because of current constraints on $\alpha$. Observations of the orbital decay of low-mass x-ray binaries [38] require that $\alpha<9\;{\rm{km}}^{2}$, or $\bar{\alpha}<220$. When making plots, we will here work with $\bar{\alpha}\in(0,30)$, which satisfies both the small coupling approximation and current constraints from low-mass x-ray binary observations.

II.2 Field equations

We can obtain the field equations of the theory by varying the action $S$ with respect to the metric and the scalar field, with the result

[TABLE]

where $G_{ab}$ is the Einstein tensor,

[TABLE]

while the stress-energy for the scalar field is

[TABLE]

As we stated before, we will here choose the scalar field to be massless and not self-interacting, meaning that we can set $U=0=U_{,\varphi}$ in the field equations.

Since we are interested in obtaining NS solutions in this theory, we assume that matter is described by a perfect fluid, whose stress-energy tensor is

[TABLE]

where $u^{a}$ is the four-velocity of the fluid (with pressure $p$ and total energy density $\varepsilon$) subject to the constraint $u^{a}u_{a}=-1$. Due to the diffeomorphism invariance of the theory, $T^{ab}_{\rm m}$ satisfies the conservation law

[TABLE]

as can be verified directly by taking the divergence of the field equations and using the equations of motion for the scalar field.

The EoS of cold nuclear matter characteristic of old NSs can be well approximated by a barotropic EoS, that is $p=p(\varepsilon)$. The large uncertainties on the properties of matter in NS interiors result in a wide variety of competing EoS models [39]. Here, to remain agnostic on which EoS correctly describes NS interiors we consider eight different EoSs, which cover a wide range of underlying nuclear physics models. In increasing order of stiffness we use: FPS [40], SLy [41], WFF1 [42], WFF2 [42], AP4 [40], ENG [43], AP3 [40], and MPA1 [44].

II.3 Perturbative expansion for the metric and fluid variables

Having obtained the field equations, we now present the perturbative scheme that we will use throughout this work. This approach was first introduced in [17] and was used in a number of studies involving BHs [20, 45, 35, 46]. Here, we apply this scheme for the first time to relativistic stars.

Let us consider a static, spherically symmetric star with spacetime described by the line element

[TABLE]

where the metric functions $\tau$ and $\sigma$ contain only radial dependence, and ${\rm d}\Omega^{2}={\rm d}\theta^{2}+\sin^{2}{\rm d}\phi^{2}$ is the line element of the unit two-sphere. The first step in the small-coupling approximation is to expand all variables $\vec{z}\in\{\tau,\sigma,\varphi,\varepsilon,p\}$ in a power series in $\bar{\alpha}$ as follows

[TABLE]

where the subscript $n$ determines the power of $\bar{\alpha}$ associated with $\vec{z}_{n}$, i.e. $\vec{z}_{n}=\mathcal{O}(\bar{\alpha}^{n})$.

With these expansions, we can immediately make a few observations. First, at ${\cal{O}}(\bar{\alpha}^{0})$, the scalar field is everywhere constant, because its source is zero [cf., Eq. (8b)]. We can then exploit shift-symmetry to impose $\varphi_{0}=0$. Second, at ${\cal{O}}(\bar{\alpha}^{1})$, we have $\tau_{1}=\sigma_{1}=0$. This follows from the fact that $\varphi_{0}=0$ and by Eq. (8a), the Einstein equations are identical to those of GR at this order. Furthermore, since the metric is unaffected to this order and there is no direct coupling between $\varphi$ and matter, we also have that $\varepsilon_{1}=p_{1}=0$.

In this paper, we will obtain solutions for all variables $\vec{z}$ up to $\mathcal{O}(\bar{\alpha}^{2})$. From the proceeding discussion, we can outline the steps of the calculation ahead as follows:

at ${\cal{O}}(\bar{\alpha}^{0})$, the problem is identical to GR and we have to calculate $\{p_{0},\varepsilon_{0},\tau_{0},\sigma_{0}\}$; 2. 2.

at ${\cal{O}}(\bar{\alpha}^{1})$, we have to determine $\varphi_{1}$ on the background of a GR star obtained in the previous step; 3. 3.

at ${\cal{O}}(\bar{\alpha}^{2})$, we must take into account the backreaction of the scalar field $\varphi_{1}$ onto the star to calculate $\vec{z}_{2}\in\{p_{2},\varepsilon_{2},\tau_{2},\sigma_{2}\}$. The first two quantities tell us how the fluid is redistributed, while the latter how the spacetime is modified relative to the background GR metric.

The perturbative scheme outlined above could be carried out to higher orders. For instance, at $\mathcal{O}(\alpha^{3})$ we would need to calculate $\varphi_{3}$ using the solutions for $\vec{z}_{2}$. Then, at $\mathcal{O}(\alpha^{4})$, $\varphi_{3}$ would be used to obtain $\vec{z}_{4}$. We here stop our calculations at $\mathcal{O}(\bar{\alpha}^{2})$ because this is the lowest order at which the metric is modified, and therefore, the lowest-order we must have at hand if we want to calculate sGB corrections to astrophysical observables.

II.4 Perturbative expansion of the field equations

At $\mathcal{O}(\bar{\alpha}^{0})$, the GR limit of the field equations give

[TABLE]

As usual [47], it is convenient to introduce a mass function $m_{0}=(r/2)[1-\exp(-2\sigma_{0})]$, and then from the $(t,t)$ and $(r,r)$-components of Eq. (15), we find

[TABLE]

Additionally, we can use the conservation law of Eq. (12) to obtain

[TABLE]

The system of equations (16) and (17) are known as the Tolman-Oppenheimer-Volkoff (TOV) equations [48, 49] and they are valid inside the star. The field equations outside the star can be obtained from the set above through the limits $(\epsilon_{0},p_{0})\to 0$.

At $\mathcal{O}(\bar{\alpha}^{1})$, we have to solve the following equation

[TABLE]

both inside and outside the star, where the d’Alembertian operator and the Gauss-Bonnet curvature invariant are constructed from the metric functions found at $\mathcal{O}(\bar{\alpha}^{0})$, i.e. $\tau_{0}$ and $\sigma_{0}$. Thus, equation (18) can be rewritten explicitly as

[TABLE]

At $\mathcal{O}(\bar{\alpha}^{2})$, the $(t,t)$ and $(r,r)$ components of the field equations yield

[TABLE]

while the conservation law of Eq. (12) gives

[TABLE]

in the stellar interior. The equations in the exterior can be found through the limits $(\varepsilon_{0},\varepsilon_{2},p_{0},p_{2})\to 0$.

III Solutions of the field equations outside the star

In this section, we first solve analytically, in vacuum, the equations presented in Sec. II order by order in $\bar{\alpha}$. The general solutions to these equations will depend on integrations constants. These constants can be fixed by examining the solutions’ asymptotic behavior at spatial infinity and imposing that (i) the spacetime is asymptotically flat and that (ii) the scalar field approaches zero at spatial infinity.

III.1 $\mathcal{O}(\bar{\alpha}^{0})$ equations

At this order, the solutions of Eqs. (16a)–(16b) have the usual Schwarzschild form

[TABLE]

where $a$ is an integration constant which (as we will see shortly) is related with the gravitational mass $M$ of the star. In obtaining this solution, we required that the metric be asymptotically flat near spatial infinity.

III.2 $\mathcal{O}(\bar{\alpha}^{1})$ equations

At this order, we need to consider Eq. (18). To solve it, we first calculate $\mathcal{G}_{0}$ which can easily be found using Eqs. (22) to be

[TABLE]

and in turn Eq. (19) becomes

[TABLE]

where Eq. (22) was used once again.

Equation (24) can be solved analytically to find

[TABLE]

where $c_{1}$ and $c_{2}$ are two integration constants. Requiring that the field vanishes at spatial infinity (i.e. that the cosmological background value of the scalar field is zero), we set $c_{2}=0$. Expanding $\varphi_{1}$ about spatial infinity, we find that

[TABLE]

which shows that $c_{1}$ is the scalar monopole charge. Reference [28] showed that this charge vanishes for all stars, and therefore, we can set $c_{1}=0$. The final expression for the scalar field outside the star is then

[TABLE]

III.3 $\mathcal{O}(\bar{\alpha}^{2})$ equations

At this order, we can substitute $\varphi_{1}$ [cf. Eq. (27)] into Eqs. (20). The resulting system of differential equations can be solved to find

[TABLE]

where $d_{1}$ and $d_{2}$ are integration constants and we defined the dimensionless parameter $\zeta$222Note with this definition of $\zeta$, the condition $\zeta\ll 1$ is not necessarily true. The small coupling approximation requires that $\alpha/\ell^{2}\ll 1$, and $\alpha/\ell^{2}\neq\zeta^{1/2}$. via

[TABLE]

The constants of integration can be determined by studying the asymptotic behavior of the metric functions about spatial infinity. For the $g_{tt}$ metric component we find

[TABLE]

and thus, we set $d_{2}=0$ without loss of generality, as any other choice corresponds to a simple rescaling of the time coordinate $t\to t\exp(d_{2})$. From the $1/r$ term we identify

[TABLE]

as a renormalized mass: the gravitational mass of the star that would be measured by an observer at spatial infinity when performing a Keplerian observation. Decomposing the mass via $M=M_{0}+M_{2}$, we can identify $M_{0}=a/2$ as the gravitational mass of a GR NS, and $M_{2}=\zeta d_{1}M_{0}$ as the sGB correction to it. A similar mass renormalization occurs for black holes [20].

We can now reexpress our exterior solution in terms of the renormalized mass. First, we eliminate $a$ in favor of $M$ in Eq. (31) and substitute the resulting equation into Eq. (28)–(28). The resulting expressions for $\tau_{2}$ and $\sigma_{2}$ can now be inserted in $g_{tt}=-\exp[2(\tau_{0}+\tau_{2})]$ and $g_{rr}=\exp[2(\sigma_{0}+\sigma_{2})]$ and then, after a reexpansion in powers of $\zeta$, we obtain our final expressions for the metric up to $\mathcal{O}(\bar{\alpha}^{2})$:

[TABLE]

The equations (32a)–(32b) are independent of $d_{1}$, and are instead fully determined by the mass $M$ of the star and the strength of the coupling constant (through $\zeta$) only.

For consistency, let us now reexpress the scalar field also in terms of the renormalized mass. The astute reader will notice that to $\mathcal{O}(\bar{\alpha}^{1})$ we can simply replace $a\to 2M$ in Eq. (27) to obtain:

[TABLE]

which is our final expression for the scalar field at $\mathcal{O}(\bar{\alpha}^{1})$. The sGB corrections to $a$ can be ignored in the scalar field, as they would enter at $\mathcal{O}(\bar{\alpha}^{3})$.

III.4 Comparison with black holes spacetimes

Before proceeding with the interior solution, let us compare the solutions obtained above to their counterparts for BHs [20], focusing first on the scalar field solution. The only difference between the calculation performed here and the one carried out for BHs is that in the latter case $\varphi_{1}$ must be regular at the event horizon. This results in a nonzero value of $c_{1}$ that yields

[TABLE]

which is identical to the second term in Eq. (33) with $M$ replaced by the hole’s mass $M_{\bullet}$. In a sense then, $\varphi_{1}$ is equal to $\varphi^{\rm BH}_{1}$ plus a correction that arises from its continuity across the stellar surface. We also observe that the BH limit of the NS solution for $\varphi_{1}$ is discontinuous. One can see this easily by evaluating Eq. (33) at the surface of the star $R_{0}$ and taking the BH limit, $M/R\sim M_{0}/R_{0}\to 1/2$, which is possible for certain anisotropic fluids in GR [50, 51].

Let us now compare the NS and BH exterior solutions for the exterior metric. As in the case of the scalar field, requiring that the metric tensor be regular at the horizon yields

[TABLE]

As in the case of the scalar field, the BH solution contains no logarithmic terms, implying that the BH limit of the NS solution is singular. This is because of the different choice of constants of integration in the NS and BH cases. Notice also that in the BH case the metric differs from GR through terms of $\delta g_{tt}^{\rm BH}={\cal{O}}(M_{\bullet}^{3}/r^{3})$ and $\delta g_{rr}^{\rm BH}={\cal{O}}(M_{\bullet}^{2}/r^{2})$ in the far field, while in the NS case it differs through terms of $\delta g_{tt}={\cal{O}}(M^{7}/r^{7})=\delta g_{rr}$.

III.5 Note on the absence of a scalar charge in the NS solution

In the case of other scalar-tensor theories, the work of [25, 26] shows that the exterior metric depends explicitly on a scalar charge, which in turn depends on the metric and matter terms within the star. If one wished to find the values of this charge, one would have to solve the TOV equations numerically for a specific set of initial conditions and a given equation of state. In our case, the exterior metric does not depend on any scalar charge, but rather it depends only on the mass, radius and coupling constant of the theory.

IV Solutions of the field equations inside the star

For completeness, let us now tackle the problem of solving for the fluid variables, scalar field, and metric components inside the star. This step will inevitably require numerical integrations, for a relationship between pressure $p$ and energy density $\varepsilon$ (i.e. the EoS) must be given and the resulting equations cannot be solved analytically. In this section, we present the numerical scheme and the numerical solutions for the interior fields. We stress however that the exterior solutions found in the previous section do not require these interior numerical solutions.

IV.1 $\mathcal{O}(\bar{\alpha}^{0})$ equations

As we saw in Sec. II.3, to this order we need to solve the TOV equations of GR, i.e. Eqs. (16a)–(16b) and (17). We start by choosing an EoS from our catalog for which, given a central total energy density $\varepsilon_{\rm c}$, gives the corresponding central pressure $p_{\rm c}=p(\varepsilon_{\rm c})$. We can then integrate Eqs. (16a)–(16b) and (17) from $r=0$ up to a point where $p_{0}(R_{0})=0$, which determines the star’s radius $R_{0}$.

In practice, we do this integration starting from a small, finite value of $r_{\rm c}$ and using a series solution valid in this region

[TABLE]

We terminate all integrations at the location where $p_{0}/p_{\rm 0c}=10^{-11}$. The constant $\tau_{0{\rm c}}$ in the series solution of the metric is arbitrary and is fixed a posteriori.

At the star’s surface $R_{0}$ we impose that the metric functions $\tau_{0}$ and $\sigma_{0}$ are continuous, that is:

[TABLE]

We can analytically match Eqs. (16a) and (22) at $R_{0}$ to find

[TABLE]

where $m_{0}(R_{0})$ is the mass of the star enclosed inside the radius $R_{0}$. Furthermore, Eq. (37a) fixes the value of the constant $\tau_{\rm 0c}$. Our final numerical solution for $\tau_{0}$ corresponds to a simple shift $\tau_{0}\to\tau_{0}+\tau_{\rm 0c}$.

The outcome of these integrations can be summarized in a mass-radius relation, shown in Fig. 1. In this figure, the solid lines correspond to various mass-radius curves for the EoSs in our catalog.

IV.2 $\mathcal{O}(\bar{\alpha}^{1})$ equations

At $\mathcal{O}(\bar{\alpha}^{1})$ we only need to solve Eq. (18). From the $\mathcal{O}(\bar{\alpha}^{0})$ solution, we know $a$ and $R_{0}$, which fully determines $\varphi_{1}^{\rm ext}$ and its derivative at $R_{0}$ [cf. Eq. (33)]. This information can be used as initial conditions to integrate Eq. (18) inside the star: we start our integration at $r=R_{0}$ and move in toward $r=0$. In this calculation, it is useful to note that $\mathcal{G}_{0}$ is given by

[TABLE]

inside the star [32], where the functions $m_{0}$, $p_{0}$ and $\varepsilon_{0}$ are all known from the $\mathcal{O}(\bar{\alpha}^{0})$ calculation.

The radial profiles of $\mathcal{G}_{0}$ and $\varphi_{1}$ are shown in Fig. 2 using the SLy EoS with the scalar-Gauss-Bonnet coupling fixed to $\bar{\alpha}=15$. In the top-panel, we see that $\mathcal{G}_{0}$ is mostly negative within the star, except near the surface (indicated by the dashed vertical line) where it changes sign and then matches smoothly to its exterior form, given in Eq. (23). We also observe that $\mathcal{G}_{0}$ has a larger magnitude for stars with larger values of $\varepsilon_{\rm 0c}$. This is can be seen by substituting the expansions of Eqs. (36a)–(36c) into Eq. (39). We find that $\mathcal{G}_{0}$ is negative and nearly constant close to the center of the star at $r\approx 0$, with its magnitude proportional to $\varepsilon_{\rm 0c}$. In the bottom-panel, we see that NSs with larger central energy densities $\varepsilon_{\rm 0c}$ have larger amplitudes of $\varphi_{1}$ at their cores. This is unsurprising given the fact that the source of the scalar, i.e. $\mathcal{G}_{0}$, has a larger magnitude near the stellar center. At the surface, $\varphi_{1}$ connects smoothly with its exterior solution, given by Eq. (III.2) (at this order in $\alpha$). The results for other EoSs are qualitatively the same as the ones shown here.

Let us now investigate how the central values of the scalar field $\varphi_{1}$ vary as a function of both $\varepsilon_{\rm 0c}$ and of $\bar{\alpha}$. This dependence is shown in Fig. 3 for four representative values of $\bar{\alpha}=\{5,10,15,20\}$ covering a range of central energy densities $\varepsilon_{\rm 0c}$ that span stars with masses $0.552\,M_{\odot}$ to $2.03\,M_{\odot}$ using the SLy EoS. We see that for small $\varepsilon_{\rm 0c}$ (i.e. low-mass stars) all values of $\varphi_{\rm 1c}$ converge towards zero regardless of the strength of the coupling. This is can understood by noticing that in this limit $\mathcal{G}_{0}$ is very small and nearly flat (cf. Fig. 2), thus sourcing $\varphi_{1}$ weakly. For larger $\varepsilon_{\rm 0c}$, the situation is different and we see a stronger dependence of the central value of $\varphi_{1}$ on $\bar{\alpha}$. Unsurprisingly, the magnitude of $\varphi_{1}$ is larger at the stellar core the larger the strength of the coupling $\bar{\alpha}$.

IV.3 $\mathcal{O}(\bar{\alpha}^{2})$ equations

At $\mathcal{O}(\bar{\alpha}^{2})$ we need to solve Eqs. (20)–(20b) and (II.4). The boundary conditions are similar to those at $\mathcal{O}(\bar{\alpha}^{0})$. Imposing continuity at the surface gives us

[TABLE]

where $R_{2}$ is the radius of the NS at $\mathcal{O}(\bar{\alpha}^{2})$ given by the condition:

[TABLE]

As in the $\mathcal{O}(\bar{\alpha}^{0})$ integrations, we start from $r_{\rm c}$ and integrate outwards until the point $R_{2}$ where the condition $(p_{2}+p_{0})/p_{\rm 0c}=10^{-11}$ is met. Equation (40) allow us to determine the numerical value of $d_{1}$, which in turn allows us to calculate the renormalized mass $M$ [Eq. (31)] and thereby determine the exterior metric in terms of interior quantities.

When integrating Eqs. (20)–(20b) and (II.4) we need to be careful on how we calculate the perturbed density $\varepsilon_{2}$. To do this, we take our total density $\varepsilon(p)=\varepsilon_{0}+\varepsilon_{2}$, and solve for $\varepsilon_{2}$ as

[TABLE]

where $\varepsilon(p_{0}+p_{2})$ is a spline interpolation of our EoS table. This allows us to eliminate the variable $\varepsilon_{2}$ in favor of the perturbed pressure $p_{2}$. With Eq. (42) and the solutions to all fields up to $\mathcal{O}(\bar{\alpha})$, we can solve our system of equations given by Eqs. (20)–(20b) and (II.4).

The dashed curves in Fig. 1 show the mass-radius relations calculated to $\mathcal{O}(\bar{\alpha}^{2})$ for the various EoSs of our catalog. We see that for a fixed value of $\bar{\alpha}$ the deviations from the GR mass-radius relation occur at larger masses. This is consistent with our previous observations on $\varphi_{1}$, which had a larger magnitude for larger masses. Consequently, these large scalar fields backreact more strongly onto the GR solution, causing larger changes to the mass and the radius. The sGB corrections typically lead to less massive NSs regardless of the EoS considered, a result consistent with those of [22]. For clarity, in Fig. 1 we only showed curves with a fixed $\bar{\alpha}=15$, but how do the mass-radius curves change (for a fixed EoS) as we vary $\bar{\alpha}$? This is shown in Fig. 4 for the SLy EoS. As expected, from our previous discussion of the $\mathcal{O}(\bar{\alpha}^{1})$ results, an increase in $\bar{\alpha}$ causes larger deviations in the mass-radius curve. Indeed, the larger the value of the sGB coupling, the smaller the maximum NS mass that is allowed for a given EoS.

Figure 1 also shows vividly the difficulties of testing modified theories of gravity with masses and radii measurements of NSs. In the absence of a complete understanding of matter in the NS interior, the various competing EoS models predict NSs that cover a wide portion of the $(M,R)$ plane. But this problem could be averted if, in the future, the EoS is constrained through NICER [8, 9, 10] and/or LIGO/VIRGO [52]. Let us imagine, for example, that the SLy EoS is favored by observations. If so, the observation of a $\approx 2M_{\odot}$ NS (see e.g. [53]) would place the stringent constraint $\bar{\alpha}\lesssim 10$ (roughly one order of magnitude more stringent than current bounds), since for larger values a SLy EoS could not predict such a massive NS (see Fig. 4). This constraint would be weaker if the true EoS is stiffer (e.g. MPA1 and AP3), as larger values of $\bar{\alpha}$ would be required to pull the mass-radius curve below $\approx 2M_{\odot}$, but stiffer EoSs are disfavored by recent tidal deformability constraints from the GW170817 gravitational-wave event [54].

Our results shown above are in agreement with those obtained previously in the literature. For example, NS solutions in the full EdGB theory were presented first in Ref. [22]. We compared our final results with those from Ref. [22] (see e.g. Figs. 1-3 in that paper) and found good agreement when the coupling was small (i.e. in the decoupling limit).

V Astrophysical applications

Now that we have a analytic solution for the exterior spacetime of a NS in sGB gravity [see Eqs. (32a)–(32b)], let us explore some astrophysical applications to investigate the physical effects of the corrections on observables.

Probing astrophysical phenomena in the vicinity of NSs naturally requires that one first analyze the geodesic motion of massive test particles and of light in the stellar exterior. Since our metric is static and axisymmetric, we know it possesses a timelike and azimuthal Killing vector, which imply the existence of two conserved quantities: the specific energy $E$ and the specific angular momentum $L$

[TABLE]

where the dots indicate differentiation with respect to proper time. From normalization condition of the four-velocity, $u^{a}u_{a}=\epsilon$, where $\epsilon=(-1\text{ or }0)$ for time-like or null trajectories respectively, we obtain

[TABLE]

which describes the radial motion of the particle in terms of the effective potential

[TABLE]

Because of spherical symmetry, we can set $\theta=\pi/2$ (and therefore $g_{\phi\phi}=r^{2}$) without loss of generality.

V.1 Circular orbits around the star

Let us study the circular motion of massive test particles ($\epsilon=-1$) around a NS with exterior metric given by Eqs. (32a)–(32b). For a circular orbit at $r=r_{\ast}$, the conditions $V_{\rm eff}(r_{\ast})=0$ and $V_{\rm eff}^{\prime}(r_{\ast})=0$ must be satisfied. Using Eq. (45), we can solve for $E$ and $L$, and expand in powers of $\zeta$ to obtain

[TABLE]

where $E_{0}$ and $L_{0}$ are the GR specific energy and angular momentum for circular orbits [47]

[TABLE]

and $E_{2}$ and $L_{2}$ are modifications of $\mathcal{O}(\bar{\alpha}^{2})$. The latter are given by

[TABLE]

and

[TABLE]

We may now make use of Eqs. (44) and (45) along with our circular orbit conditions to find the sGB modifications to the location of the ISCO. Doing so, we find that the ISCO radius is

[TABLE]

which reproduces the well-known GR result when $\zeta=0$. Notice that the sGB correction pushes the ISCO location farther away from the stellar surface (assuming the star is sufficiently compact so that the ISCO is outside the surface in GR in the first place) by a small amount $R_{\rm ISCO}-6M\approx+0.089M\zeta$.

V.2 Modified Kepler’s third law

Now let us derive an expression for the orbital frequency $\Omega_{\phi}={\rm d}\phi/{\rm d}t$ of a massive particle in circular orbit at radius $r_{\ast}$ as measured by an observer at infinity. Using Eqs. (43) and (46a)–(46b) we find

[TABLE]

where $\Omega_{0}^{2}=M/r_{\ast}^{3}$ is the usual GR result. Expanding Eq. (V.2) in the far field limit we find to leading order in $\zeta$

[TABLE]

which is consistent with our expansions of the sGB metric deformation in Sec. III.4. Unlike the case for BHs, where the correction to the frequency occurs as $\mathcal{O}(M^{2}/r^{2})$ [20], the presence of the logarithmic term in Eq. (V.2) gives a small correction. This suggests that weak-field observables will be very poor probes of sGB gravity.

V.3 Quasiperiodic oscillations

Let us now focus on the frequencies of quasiperiodic oscillations (QPOs). There are a number of models which have been proposed as possible causes of QPOs including the relativistic motion of matter [55] and resonance between orbital and epicyclic motion [56]. Regardless of the model in question, it may be interesting to calculate the sGB corrections the QPO frequencies to study what effect, if any, this modification to GR has.

The orbital frequency was already calculated in Eq. (52), so let us now calculate the epicyclic frequency for timelike geodesics. This frequency is determined by a radially-perturbation to the circular orbit equation [Eq. (44)], which yields

[TABLE]

Solving Eq. (53) with the $V_{\rm eff}(r)$ defined in Eq. (45) gives us

[TABLE]

If one were to asymptotically expand this frequency about spatial infinity, one would again find that the sGB corrections are highly suppressed. As we will see below, however, QPOs are sensitive to physics near the ISCO, and in this regime, the sGB corrections are not nearly as suppressed.

In addition to these two frequencies, there is often a third one that is important in QPOs and measures the rate of periastron precession of the orbit. This precession frequency can be found via

[TABLE]

and it is usually important in lower frequency QPOs333Some models treat this frequency as stemming from inhomogeneities near the inner accretion disk boundary, causing a beat frequency [57]. However, this was found to be inconsistent with observations [58]. [59]. With these three frequencies in hand, one could imagine using the observation of QPOs to place constraints on sGB. Figure 5 depicts two of our frequencies against one another (in dimensionless units) and illustrates how there are noticeable deviations from the GR predictions as $\alpha$ increases. Observe that the frequencies approach each other when either of them is small, since here one approaches the weak-field regime described in Eq. (52).

One may present these frequencies in terms of an observable quantity, namely the dimensionless linear orbital velocity $v$, as done in [60, 59]. By introducing the orbital velocity as $v=\left(M\Omega_{\rm sGB}\right)^{1/3}$, we may reexpress the ratio of the precession frequency to the orbital frequency as a series in velocity to obtain

[TABLE]

where the modification to the GR solution again is suppressed by a high power of velocity that is consistent with the expansion of Eq. (52). As before, the largest deviations will then occur for observables that are sensitive to physics near the surface of the NS, i.e. where the orbital velocity is not extremely small.

V.4 Light bending

Let us now consider photon motion in the sGB exterior spacetime, as depicted in Fig. 6. Imagine then a photon leaving the surface of the NS along the unit vector $\vec{k}$, which makes an angle $\gamma$ with the unit vector $\vec{n}$ normal to the star’s surface. The angle $\psi$, between $\vec{n}$ and the line of sight, is an important quantity in astrophysical applications. For instance, when $\gamma=\pi/2$, $\psi=\psi_{\rm crit}$ is the critical angle between the line of sight and the normal to the surface beyond which the photon cannot reach the observer. This allows one to define a visible fraction of the star as

[TABLE]

Moreover, in the context of pulse profile modeling, photons emitted by the hot spot can only reach the observer when the are emitted if emitted with $\cos\psi>\cos\psi_{\rm crit}$ [61].

Let us now derive an expression for $\psi$. We again restrict attention to equatorial orbits, such that $\theta=0$ and $\dot{\theta}=0$, and change notation $\phi\to\psi$ in Eqs. (43) and (44) with $\epsilon=0$. Solving for the fraction ${\rm d}\psi/{\rm d}r$ yields

[TABLE]

Since $E$ and $L$ are constant, we can simplify the above expression through the emission angle $\gamma$, defined via [61]

[TABLE]

The above expression allows us to find a relation between $E$, $L$, and $\gamma$, namely444The ratio of $L/E$ is also called the impact parameter [62], which is denoted as $b$ in Fig. 6.

[TABLE]

where we evaluate the metric functions at the stellar surface. Substituting Eq. (60) into Eq. (58) gives a direct relation between $\psi$ and $\gamma$ for a given $R$, which can be solved to obtain

[TABLE]

The integral in Eq. (61) may not be straightforward to solve, even numerically, but following [63, 64], we can rewrite it in terms of the compactness $M/R$ and a new variable $x=\sqrt{1-R/r}$ to ease the numerical integration.

The results of evaluating Eq. (61) as a function of the compactness are shown in Fig. 7.

Observe that there is a greater deflection of light for NSs of greater compactness. This is apparent even in the GR limit, and it is due to the effects of curvature near compact objects. However, this effect is enhanced in sGB gravity, increasing with larger $\zeta$555The relation between $\zeta$ and $\alpha$ depends on the mass of the NS, which is not specified here. As a reference, for a $1.4\,M_{\odot}$ NS, $\bar{\alpha}=\left(10,20,30\right)$ corresponds to $\zeta\approx\left(0.5,2.1,4.6\right)$., which dictates how strongly the $\mathcal{G}$ correction contributes to the system. For stars with smaller masses and larger radii, there is a negligible change in the deflection of light, regardless of the strength of $\zeta$. The curvature of spacetime near the surface of these NS is simply not large enough even with the quadratic curvature nature of our theory to cause any deviations that may be detectable in future observations.

We may also look at how light bending in sGB gravity compares to light bending in GR, as shown in Fig. 8 for various choices of $\zeta$ values and two fixed compactnesses. As with Fig. 7, there are only tiny deviations when the compactness is small. However, NSs with larger compactnesses do present sGB corrections to light bending that make it stronger relative to GR.

As a final calculation, we can also find the visible fraction of the NS surface, given in Eq. (57). This is shown in Fig. 9. In agreement with our previous results, there is little to be learned about sGB gravity from observations of low compactness stars. However, as the compactness increases, so does the effects of the coupling with the Gauss-Bonnet invariant. Likewise, larger values of the coupling constant lead to larger changes in the visible fraction. In GR, it is known that for NSs with $M/R\approx 0.28$, strong gravitational light bending can make the whole surface of the star visible [65]. The effect of the scalar-Gauss-Bonnet coupling is to reduce the necessary compactness the whole surface of the star to become visible. For instance, when $\zeta=5$, this compactness is 0.264.

VI Conclusions and outlook

In this paper, we obtained an analytical metric that represents the exterior spacetime of a NS in sGB gravity as well as an analytical expression for the scalar field. The metric was derived through a small-coupling perturbative scheme and depends only on the mass of the NS in question and the desired strength of the coupling constant. Our metric is valid to $\mathcal{O}(\bar{\alpha}^{2})$ and we have outlined how higher-order corrections can be obtained. We applied the new spacetime to a sample of astrophysical applications, including the motion of test particle (which is important for instance to model QPOs) and light bending (which is important to model x-ray pulse profiles generate by hot spots at the surface of rotating NSs).

Our work opens the door for number a future studies with NSs in sGB gravity, with the convenience of now being able to treat the metric analytically. One application could be the development of an effective-one-body (EOB) formalism, to model NS binaries in sGB, along the lines of the recent work in scalar-tensor theories [66]. Having the theory expressed in an EOB framework allows one to understand the dynamics of the two-body system, the radiation-reaction components of the system, and knowledge of the gravitational-waveform emitted from a coalescing binary [67].

Another possible use for the exterior metric is in the modeling of x-ray pulse profiles as a possible test bed for sGB gravity. These pulse profiles are generated by the x-ray emission from hot spots on the surface of rotating NS [65] (see [68, 69, 70] for reviews). As the photons propagate from the surface towards the observer, they probe the spacetime around the NS which, in principle, can leave detectable deviations in the observed pulse profile relative to what is predicted in GR. This possibility was recently explored in the context of scalar-tensor theories [71, 11, 12]. In particular, Ref. [12] showed that in principle observations made by NICER can constrain these theories. It would be interesting to see if the same is true in sGB gravity.

A final observation of interest is the absence of any sGB gravity integration constants in the final expression for the exterior metric presented in Eqs. (32a)–(32b). This is rather unexpected because in other theories (such as in scalar tensor theories) the exterior metric does depend on charges that must be computed numerically. In sGB gravity, however, the exterior metric is fully determined in terms of the mass of the star $M$ and the coupling constant of the theory $\alpha$. A deeper physical or mathematical understanding of why this is the case in sGB gravity would be most interesting and will be studied elsewhere.

Acknowledgements.

This work was supported by NSF Grant No. PHY-1250636 and PHY-1759615, as well as NASA grants NNX16AB98G and 80NSSC17M0041. We thank Alejandro Cárdenas-Avendaño, Paolo Pani, Thomas Sotiriou, and Kent Yagi for helpful discussions. Computational efforts were performed on the Hyalite High Performance Computing System, operated and supported by University Information Technology Research Cyberinfrastructure at Montana State University.