Epicyclic oscillations in the Hartle-Thorne external geometry

TL;DR

This paper analyzes the orbital and epicyclic frequencies in the external spacetime of rotating neutron stars described by the Hartle-Thorne metric, with applications to understanding high-frequency QPOs in X-ray binary systems.

Contribution

It provides a systematic analysis of epicyclic motions in Hartle-Thorne spacetime using realistic neutron star models and universal relations, enhancing QPO modeling accuracy.

Findings

Epicyclic frequencies depend on neutron star parameters and quadrupole moments.

Results support models explaining QPOs with frequency ratios near 3:2 or 5:4.

The approach allows direct application across different neutron star configurations.

Abstract

The external Hartle-Thorne geometry, which describes the space-time outside a slowly-rotating compact star, is characterized by the gravitational mass $M$, angular momentum $J$ and quadrupole moment $Q$ of the star and gives a convenient description which, for the rotation frequencies of more than 95 % of known pulsars, is sufficiently accurate for most purposes. We focus here on the motion of particles in these space-times, presenting a detailed systematic analysis of the frequency properties of radial and vertical epicyclic motion and of orbital motion. Our investigation is motivated by X-ray observations of binary systems containing a rotating neutron star which is accreting matter from its binary companion. In these systems, twin high-frequency quasi-periodic oscillations are sometimes observed with a frequency ratio approaching $3:2$ or $5:4$ and these may be explained by models involving the orbital and epicyclic frequencies of quasi-circular geodesic motion. In our analysis, we use realistic equations of state for the stellar matter and proceed in a self-consistent way, following the Hartle-Thorne approach in calculating both the corresponding values of $Q$, $M$ and $J$ for the stellar model and the properties of the surrounding spacetime. Our results are then applied to a range of geodetical models for QPOs. A key feature of our study is that it implements the recently-discovered universal relations among neutron star parameters so that the results can be directly used for models with different masses $M$, radii $R$ and rotational frequencies $f_\mathrm{rot}$.

Tables

Table 1. Table 1: An overview of the orbital models of QPOs, with specification of the upper  ν U subscript 𝜈 U \nu_{\mathrm{U}} and lower  ν L subscript 𝜈 L \nu_{\mathrm{L}} frequencies.

Frequency	RP	TD	WD	RP1	RP2	ER
${\nu }_U$	${\nu }_K$	${\nu }_K+{\nu }_{\mathrm{rad}}$	$2{\nu }_K-{\nu }_{\mathrm{rad}}$	${\nu }_{\theta }$	$2{\nu }_K-{\nu }_{\theta }$	${\nu }_{\theta }$
${\nu }_L$	${\nu }_K-{\nu }_{\mathrm{rad}}$	${\nu }_K$	$2\left({\nu }_K-{\nu }_{\mathrm{rad}}\right)$	${\nu }_K-{\nu }_{\mathrm{rad}}$	${\nu }_K-{\nu }_{\mathrm{rad}}$	${\nu }_{\mathrm{rad}}$

Equations

$$\mathrm{d}{s}^2$$

$$g_{tt}$$

$$g_{rr}$$

$$g_{\theta \theta }$$

$$g_{\varphi \varphi }$$

$$g_{t\varphi }$$

$$F_1(r)$$

$$F_2(r)$$

$$G_1(r)$$

$$L(r)$$

$$A_1(r)$$

$$H_1(r)$$

$$H_2(r)$$

$$A_2(r)$$

$$u$$

$$r_{\mathrm{BL}}$$

$${\theta }_{\mathrm{BL}}$$

$${\nu }_K(r;M,j,q)$$

$$E_1(r)$$

$$E_2(r)$$

$${\nu }_r^2(r;M,j,q)$$

$${\nu }_{\theta }^2(r;M,j,q)$$

$$C_1(r)$$

$$C_2(r)$$

$$C_3(r)$$

$$B_1(r)$$

$$D_1(r)$$

$$D_2(r)$$

$$D_3(r)$$

$$B_2(r)$$

$$r_{\mathrm{ms}}=6M\left[1\mp j\frac{2}{3}\sqrt{\frac{2}{3}}+j^2\left(\frac{251647}{2592}-240\ln \frac{3}{2}\right)+q\left(-\frac{9325}{96}+240\ln \frac{3}{2}\right)\right],$$

$$f_{\mathrm{ms}}=\frac{c^3}{2\pi GM6\sqrt{6}}\left[1+\frac{11j}{6\sqrt{6}}+\frac{1}{864}{j}^2\left(-160583+397710\ln \frac{3}{2}\right)+\frac{5}{32}q\left(1193-2946\ln \frac{3}{2}\right)\right].$$

$$\frac{d{\nu }_r}{dr}=0.$$

$$r_{\max }=8M\left[1-\frac{47}{64\sqrt{2}}j+\frac{-12596207+43683840\ln \frac{4}{3}}{147456}{j}^2-\frac{5}{576}\left(-9835+34128\ln \frac{4}{3}\right)q\right]$$

$${\nu }_{\mathrm{r}}(r_{\max })=\sqrt{\frac{c^3}{64\sqrt{2}\pi GM}\left[1+\frac{15}{16\sqrt{2}}j+\frac{-3972109+13824000\ln \frac{4}{3}}{8192}{j}^2+\frac{5}{16}q\left(1553-5400\ln \frac{4}{3}\right)\right]},$$

$${\nu }_{\mathrm{r}}(r_{\max })=\sqrt{\frac{714.011}{M}(1+0.663j+0.587j^2-0.15q)}.$$

$${\nu }_{\mathrm{P}}^{\mathrm{r}}={\nu }_K-{\nu }_r,$$

$${\nu }_{\mathrm{P}}^{\theta }={\nu }_K-{\nu }_{\theta },$$

$${\left[\frac{{\nu }_{\theta }(r;M,j,q)}{{\nu }_K(r;M,j,q)}\right]}^2=1-jK(r)+j^{2}L(r)-qM(r),$$

$$K(r)$$

$$L(r)$$

$$M(r)$$

Full text

A key feature of our study is that it implements the recently-discovered universal relations among neutron star parameters so that the results can be directly used for models with different masses $M$, radii $R$ and rotational frequencies $f_{\mathrm{rot}}$.

Hartle–Thorne model — QPOs — stars: epicyclic oscillations — stars: pulsars — X-rays: binaries

††journal: ApJ††software: Mathematica (version 11.3, Wolfram Research (2018)) and LORENE/nrotstar (Bonazzola et al. (1998); Gourgoulhon et al. (2016)

1 Introduction

Exterior space-times of rotating compact stars have been studied extensively for many years. Pioneering work was done by Hartle and Thorne who developed a slow-rotation approximation (Hartle, 1967; Hartle & Thorne, 1968), describing the structure of the star itself and also of the surrounding vacuum space-time, constructed as a perturbation of a corresponding spherically-symmetric non-rotating solution, with the perturbation being taken up to second order in the star’s angular velocity $\Omega$. In this approximation the exterior space-time is fully determined by the gravitational mass $M$, angular momentum $J$ and quadrupole moment $Q$ of the rotating star, with the inner boundary of the exterior region being given by the stellar radius $R$.

Solutions for the equivalent problem without the restriction to slow rotation can be constructed using numerical relativity codes such as RNS (Stergioulas & Friedman (1995)) and LORENE/nrotstar (Bonazzola et al. (1998); Gourgoulhon et al. (2016)), but there are certain advantages to working with the Hartle–Thorne method, including that of having an analytic solution for the vacuum metric outside the star. The main topic of this paper concerns the motion of particles around rotating neutron stars, and we use the rotationally induced change in the frequency with which a test particle orbits the star at the marginally stable circular orbit as a key quantity for comparing different approaches. Comparison of the Hartle–Thorne values for this with those coming from the numerical relativity codes shows good agreement for dimensionless spins $j=J/M^{2}\simeq 0.5$111Throughout this paper, we are using units for which $c=G=1$., suggesting applicability of other results obtained using the Hartle–Thorne approximation in most astrophysically relevant situations, although this is something which requires further checking. Comparisons between numerical-relativity and analytical space-times have been discussed by a number of authors (see, for example, Nozawa et al. (1998); Berti & Stergioulas (2004); Berti et al. (2005)). Something not considered in the present paper is the role of the magnetic field of a neutron star for affecting the motion of particles in its vicinity. Discussion of this can be found in articles by Sanabria-Gómez et al. (2010); Bakala et al. (2012); Gutierrez-Ruiz et al. (2013) and others.

Details and overviews of rotating neutron stars and related physics can be found especially in the books by Weber (1999) and Friedman & Stergioulas (2013) and in the review article by Paschalidis & Stergioulas (2017).

Rotating neutron stars are usually observed as pulsars, either isolated or in binary systems. In binary systems they can reach high rotation frequencies due to accretion of matter from the binary companion which increases both the star’s gravitational mass and its angular momentum (and hence its rotational speed). The rotational frequencies can reach hundreds of Hz (the fastest currently known pulsar rotates with frequency 716Hz (Hessels et al., 2006)). Within this work we are interested mostly in neutron stars which are in Low Mass X-ray Binaries, where the neutron star is accompanied by a low-mass ordinary star. Observations of some of these objects exhibit the phenomenon of twin high-frequency quasi-periodic oscillations (QPOs) (van der Klis, 2006).

These QPOs consist of two adjacent peaks in the X-ray power spectrum and have been observed in a number of sources. They often have a lower frequency of $600-800\,\mathrm{Hz}$ and an upper frequency of $900-1200\,\mathrm{Hz}$, and there is clustering of the observed frequency ratios particularly at around $3:2$ but also at around $4:3$ and $5:4$ (Belloni et al., 2007; Török et al., 2008a, b, c; Boutelier et al., 2010). Similar clustering is observed in microquasars (low mass X-ray binaries containing black holes) where the frequencies are fixed at the ratio of $3:2$ and can be explained as resulting from a non-linear resonance between radial and vertical oscillation modes of an accretion disc around the black hole (Török et al., 2005). Therefore, we can anticipate that the $3:2$ ($4:3$ and $5:4$) resonance should play a role also in systems containing neutron stars instead of black holes (Török et al., 2008c). Nevertheless, it is still unclear whether or not the clustering in the distribution of frequency ratios is a real effect (Gilfanov et al., 2003; Méndez, 2006; Boutelier et al., 2010; Montero & Zanotti, 2012; Ribeiro et al., 2017).

There have been several attempts to explain the physical origin of QPOs using simple orbital models based on motion of matter around a central object (e.g. Stella & Vietri (1998, 1999); Wagoner (1999); Stella & Vietri (2001); Abramowicz & Kluźniak (2001); Kluzniak & Abramowicz (2001); Wagoner et al. (2001); Abramowicz et al. (2003b). A large set of references to models for neutron-star QPOs, along with a direct comparison between them, can be found in Török et al. (2010, 2012, 2016)). A more elaborate model which included corrections to the orbital motion due to finite thickness of the accretion disk has been investigated by Török et al. (2016) and was particularly successful in the case of 4U1636-53 where it was shown that the observed frequency pairs could be explained by a torus with its center (maximum of density) located at positions which were variable but always at large enough radius to allow accretion onto the central object via the inner cusp; this model has been called the cusp torus model (see the papers by Rezzolla et al. (2003); Abramowicz et al. (2006); Šrámková et al. (2007); Ingram & Done (2010); Fragile et al. (2016); Török et al. (2016); Mishra et al. (2017); Török et al. (2017a, b); Parthasarathy et al. (2017); de Avellar et al. (2017) and references therein). Within this model the highest observed frequencies correspond to the torus being located close to the marginally stable orbit and being small enough so that the correction due to the finite thickness of the disk could be neglected, enabling one to put constraints on the radius of the central object.

Motivated by the cusp torus model, we here make a detailed study of quasi-circular epicyclic motion in the external Hartle–Thorne geometry, using the formulas obtained previously by Abramowicz et al. (2003a). We begin our analysis by using the Hartle–Thorne approach to calculate neutron-star properties for a range of realistic equations of state for the neutron-star matter, and demonstrate how the influence of the different prescriptions for the equation of state can be hidden if one plots particular combinations of neutron star properties against one another. These results are manifestations of the universal relations which have been discovered rather recently (Yagi & Yunes, 2013a, b; Urbanec et al., 2013; Maselli et al., 2013; Pappas, 2015; Reina et al., 2017) and have been studied extensively under various circumstances (Chakrabarti et al., 2014; Haskell et al., 2014; Doneva et al., 2014a, b; Sham et al., 2015; Doneva et al., 2015; Pani et al., 2015; Silva et al., 2016; Staykov et al., 2016).

We then move on to discuss properties of particle motion around the star, showing the radial profiles for the orbital (Keplerian) frequency of circular motion and the radial and vertical frequencies of epicyclic quasi-circular motion. These profiles are relevant down to the radius of the innermost stable circular geodesic. We then study the frequency ratio of the oscillatory modes related to various geodetical models which have been proposed for explaining twin HF QPOs: the relativistic precession model and its variants (Stella & Vietri, 1998, 1999), an epicyclic resonance model (Török et al., 2005), a tidal disruption model (Kostić et al., 2009) and a warped disc oscillation model (Kato, 2004)222The epicyclic frequencies in the case of fully numerical models for rapidly rotating neutron stars have been studied in (Pappas, 2012) and (Gondek-Rosińska et al., 2014)..

In Section 2 we give a brief introduction to the Hartle–Thorne methodology and then present results for neutron-star models constructed using our set of equations of state. In Section 3, we present general properties of the radial profiles for the orbital frequency, the epicyclic frequencies and the precession frequencies and also describe the behaviour of the maximum of the profile for the radial epicyclic frequency. In Section 4 we show radial profiles for the frequency ratios corresponding to various models for HF QPOs and in Section 5 we compare models calculated using the Hartle–Thorne approach with fully numerical models not restricted to slow rotation. The main results presented within Sections 1-5 are parametrized in terms of the NS mass, angular momentum and quadrupole moment. This allows their convenient application for interpreting neutron-star X-ray variability data within the context of the universal relations. We end with a short Summary in Section 6.

2 The Hartle–Thorne approach for calculating models of rotating compact stars

The Hartle–Thorne approach (Hartle, 1967; Hartle & Thorne, 1968) is a convenient scheme for calculating models of compact stars in slow and rigid rotation and, since its introduction, has been used and discussed by many other authors. In this scheme, deviations away from spherical symmetry are considered as being small perturbations and treated by means of series expansions with terms up to second order in the star’s angular velocity $\Omega$ being retained.

The metric line element takes the form (in standard Schwarzschild coordinates):

[TABLE]

where $\nu$ and $\lambda$ are functions of the radial coordinate $r$, and reduce to those of the spherically symmetric vacuum Schwarzschild solution in the exterior, $\omega(r)$ is a perturbation of order $\Omega$ which represents the dragging of inertial frames, and the perturbation functions $h_{0}(r)$, $h_{2}(r)$, $m_{0}(r)$, $m_{2}(r)$ and $k_{2}(r)$ are of order $\Omega^{2}$ and are functions only of the radial coordinate with the non-spherical angular dependence on the latitude $\theta$ being given by the Legendre polynomial $P_{2}=P_{2}(\theta)=[3\cos^{2}(\theta)-1]/{2}$.

The interior metric is found by solving the Einstein equations with the stress-energy tensor on the right hand side being that for a perfect fluid rotating rigidly with angular velocity $\Omega$, while the exterior solution is found by solving the vacuum Einstein equations, and constants appearing in the calculations are determined by matching these two solutions at the surface of the star.333As pointed out quite recently by Reina & Vera (2015) and Reina (2016), special care is required when doing this for models having a density discontinuity at the surface in order to avoid getting wrong values for rotational corrections to the mass, but this issue does not arise for the models being discussed here. The exterior solution can be expressed in terms of the star’s gravitational mass $M$, angular momentum $J$ and quadrupole moment $Q$ 444Note that throughout we follow the convention of taking $Q$ to be positive for an oblate object, as used by Hartle & Thorne. and it can be useful to rewrite the external metric in the terms of the dimensionless quantities: $j=J/M^{2}$, a dimensionless angular momentum, and $q=Q/M^{3}$, a dimensionless quadrupole moment (Abramowicz et al., 2003a). Quadrupole moments of rotating stars have been under investigation in the past and various definitions can be applied, for more details see e.g. the discussions by Bonazzola & Gourgoulhon (1996) and Pappas & Apostolatos (2012a); Friedman & Stergioulas (2013).

Within our calculations we also introduce another quadrupole parameter $\tilde{q}=q/j^{2}=QM/J^{2}$ which is also dimensionless. Since $q\sim\Omega^{2}$ and $j\sim\Omega$, $\tilde{q}$ is independent of the star’s rotation speed within the Hartle–Thorne approximation. (A similar property has been shown to apply for rapidly rotating neutron stars as well, extending also to appropriately scaled higher order moments (Pappas & Apostolatos, 2012a, b, 2014; Yagi et al., 2014).)

2.1 Models of rotating neutron stars

Here we present results of calculations made using the Hartle–Thorne approach in order to determine the relations satisfied by the main neutron-star properties for a range of equations of state. For each EoS we show a sequence of models rotating at 400Hz (with masses up to the maximal one) and a second similar sequence with zero rotation, for comparison. Models satisfying the slow-rotation criterion will lie between the two (the motivation for focusing on 400 Hz will be explained later). The models of rotating neutron stars are constructed by solving a set of ordinary differential equations for the perturbation functions; the approach used has been extensively described in the literature and we will not repeat that here, see e.g. Hartle & Thorne (1968); Chandrasekhar & Miller (1974); Miller (1977); Weber & Glendenning (1992); Benhar et al. (2005); Paschalidis & Stergioulas (2017).

For our set of equations of state describing nuclear matter we have selected SLy4 (Stone et al., 2003), UBS (Urbanec et al., 2010), APR (Akmal et al., 1998), FPS (Lorenz et al., 1993), Gandolfi (Gandolfi et al., 2010), NRAPR (Steiner et al., 2005) and KDE0v1 (Agrawal et al., 2005). These follow various approaches to nuclear matter theory and lead to neutron star models with different properties. Fig.1 shows mass-radius relations (top left) and mass as a function of the central energy density (top right) for both rotating models (solid lines) and non-rotating models (dashed lines). The different colours correspond to the different equations of state. For the rotating models, the equatorial radius is the quantity plotted as $R$. By comparing the right and left panels, one can see that the more compact stars (having smaller total radius for a given mass) are also having higher central energy densities, as would be expected. From the right panel one can see the mass becoming less sensitive to central density as the maximum is approached. The existence of the maximum mass is important for testing equations of state using astronomical observations: quite recent observations of two solar mass neutron stars (Demorest et al., 2010; Antoniadis et al., 2013) rule out many equations of state and from our selection of nuclear matter models one can see that FPS is not applicable for these two objects.

The bottom panels of Fig.1 show the extent to which particular combinations of neutron star parameters are universally related to the inverse compactness $R/M$ (we are again using geometrical units here, so that all quantities are measured in units of length or its powers). As we will see, making use of this quasi-universal behaviour can be extremely convenient for investigating orbital motion. For given $R$ and $M$, one can immediately obtain the value of $I/M^{3}$ from the relation in the bottom left panel of Fig.1 and hence calculate the dimensionless angular momentum $j=J/M^{2}=2\pi f_{\mathrm{rot}}I/M^{2}$ for particular values of star’s rotational frequency $f_{\mathrm{rot}}$. The dimensionless quadrupole moment $q=Q/M^{3}$ can then be obtained by using the relation in the bottom right panel together with the value of $j$.

This approach enables one to skip making the calculations of neutron-star models for all possible equations of state and central parameters, and to use just global properties of the star (its mass and radius) for calculating all parameters of the external spacetime for the required rotational frequency.

For illustration, in Fig. 2 we demonstrate results from detailed modelling: the left panel shows values of $j$ as a function of the gravitational mass for stars rotating with $f_{\mathrm{rot}}=400$Hz. Within the Hartle–Thorne approximation $j=J/M^{2}$ depends linearly on rotational frequency, and so one can easily find $j$ for other values of $f_{\mathrm{rot}}$ by simple linear scaling. The maximal value of $j$ is that for a star rotating at its mass-shedding limit and, for neutron stars, reaches values of 0.65–0.7 (Lo & Lin, 2011). The right panel shows the values of $\tilde{q}=QM/J^{2}$ as a function of gravitational mass. It can be seen that for neutron stars with $M\geq 1.2\,M_{\odot}$, the values of $\tilde{q}$ are less than 10, and we will use this later as a maximum in our investigation of the motion of particles in the field of rotating neutron stars, to keep us within the limits of astrophysicaly interesting objects.

In this section we have discussed the universal behaviour of the parameters $j$ and $\tilde{q}$ and how they are related to the mass, radius and rotational frequency of neutron star models calculated using realistic equations of state. In the following sections, we will make use of this for studying the motion of particles in the space-time outside rotating neutron stars within the Hartle–Thorne approximation. A similar investigation, but using a different model for the space-time around the rotating neutron stars, has been performed by Pachón et al. (2012) who presented their results in a more condensed form than being done here. Our approach is appropriate for any model of the rotating neutron stars since we plot the results in a way that is not dependent on any particular choice of equation of state.

3 The external Hartle–Thorne geometry and the frequencies of epicyclic motion

Since our main interest here is to investigate free-particle motion external to the surface of rotating neutron stars, we will next discuss this part of the space-time (for $r>R$ to lowest order) in more detail. We should stress that all of the discussion of this section and the next one, apply only to the region outside the surface of the neutron star. Any reference to free particle motion which would seem to occur inside the neutron star needs to be disregarded.

As we have been discussing, the external Hartle–Thorne geometry can be expressed in terms of three parameters, for which we here take the gravitational mass $M$, the dimensionless angular momentum $j$ and the dimensionless quadrupole moment $q$, in the form (Abramowicz et al., 2003a)555Note misprints in the original paper.:

[TABLE]

where

[TABLE]

For $j=0$ and $q=0$ the external Hartle–Thorne geometry reduces to the standard Schwarzschild one. The Kerr geometry taken up to second order in the dimensionless angular momentum $a=Mj$, and using the standard Boyer-Lindquist coordinates, can be obtained from the external Hartle–Thorne geometry by setting $q=j^{2}$ and making the coordinate transformations

[TABLE]

3.1 Radial profiles of the orbital frequency and the epicyclic frequencies

Formulas for the frequencies of circular and epicyclic motion in the external Hartle–Thorne geometry have been calculated in Abramowicz et al. (2003a) and used by many authors, see e.g. Török et al. (2008a, 2014) and Boshkayev et al. (2014, 2015). Here we give the relations for the orbital (Keplerian) frequency and the radial and vertical epicyclic frequencies, which are needed for the twin HF QPO models based on geodesic quasi-circular motion and were presented in Török et al. (2008a).

The Keplerian frequency is given by

[TABLE]

where

[TABLE]

The radial epicyclic frequency $\nu_{r}$ and the vertical epicyclic frequency $\nu_{\theta}$ are given by

[TABLE]

where

[TABLE]

We next demonstrate the dependence of the radial profiles of the orbital and epicyclic frequencies on the parameters $j$ and $q$. For each of the frequency profiles, we choose the values of the spin parameter as $j=0.1,0.2,0.3,0.5$ so as to fully cover the range which might be relevant for application to astrophysically interesting models of neutron stars or strange stars, while staying within the limit of validity of the slow-rotation approximation $(j/j_{\mathrm{max}})^{2}\ll 1$. It has been shown that $j_{\mathrm{max}}\sim 0.65-0.7$ for neutron stars and that it could be larger for strange (quark) stars (Lo & Lin, 2011). For objects which are less compact than these, the value of $j_{\mathrm{max}}$ could be significantly larger.

Values of the quadrupole parameter $q$ are calculated taking $\tilde{q}=q/j^{2}=1,2,3,4,5,10$ for each of the values of $j$. This selection of values for $j$ and $q$ covers the relevant situations and is obtained from modeling of the neutron stars as described in the previous section or, for example, in Urbanec et al. (2013).

In general relativity, stable circular orbital motion is restricted to radii larger than that of the marginally stable orbit which, in the external Hartle–Thorne geometry, is given by

[TABLE]

with the $-$ and $+$ referring to co-rotating and counter-rotating motion respectively. Only outside this orbit can quasi-circular geodesic motion be stable, giving rise to possibly observable quasi-periodic oscillatory effects. Properties of the marginally stable orbit have been discussed in detail in, for example, Török et al. (2014); Cipolletta et al. (2017) and are shown here in Fig. 3. In the left frame, one can see that the radius of the marginally stable orbit $r_{ms}$ always decreases with increasing $j$ near to $j=0$, but can reach a minimum and then become increasing for higher values of $j$ if $q/j^{2}$ is large enough. This sort of behavior can have a relevant impact on the observed distribution of rotational frequencies of QPO sources (Török et al., 2014). The position of the surface of the neutron star can be calculated from the universal relations (Yagi & Yunes, 2006) (given that $q/j^{2}$ is directly related to the inverse compactness of the star $R/M$); results can be seen on the left side of Fig. 3. In the zoomed area, note the curve for $q/j^{2}=5$ where a dot indicates the position of the surface. For the rest of the specified values of $q/j^{2}$, $r_{\mathrm{ms}}$ is situated above the surface for $q/j^{2}<5$ and below the surface for $q/j^{2}>5$, within the range of $j$ shown.

The formula for the orbital frequency at $r_{\mathrm{ms}}$ can be calculated up to second order in the star’s angular velocity as

[TABLE]

The right panel of Figure 3 shows $f_{\mathrm{ms}}$ plotted against $j$ as calculated from this.

In the discussion above, we have used the term “marginally stable orbit” thinking only of the external vacuum region but, in reality, for a given model of the neutron star, this location may not exist in the external region: the equatorial radius of the star $R_{\mathrm{eq}}$ with given parameters $M$, $j$ and $q$ can be larger than the $r_{\mathrm{ms}}$ as calculated in a vacuum exterior space-time with the same parameters. In the remainder of this paper, we will be using the term “innermost stable circular orbit” (ISCO) in a non standard way, to refer either to the orbit at $r_{\mathrm{ms}}$, if it exists outside the star, or otherwise to a surface-skimming orbit at $R_{\mathrm{eq}}$. One should bear in mind that, in practice, a particle on a surface-skimming orbit would be subject to various physical effects which we are not including here.

In order to obtain a complementary point of view, in Figure 4 we show the radial profiles of the orbital and epicyclic frequencies for specified values of $j$, together with the profiles for the local value of the Keplerian frequency $\nu_{\mathrm{K}}$. One can see that the influence of the parameters $j$ and $q/j^{2}$ increases with decreasing radius, and is relatively small at radii $r\sim 10M$ for all of the three frequency profiles; in all three cases, the influence of $q/j^{2}$ increases with increasing $j$. For the orbital frequency, the role of the parameters $j$ and $q$ is the smallest, and it is significantly stronger for the profile of the vertical epicyclic frequency which has a similar behavior to that of the orbital frequency. The most important influence is in the case of the radial epicyclic frequency – the shift of the marginally stable circular orbit, corresponding to the radius where the radial epicyclic frequency vanishes, is significantly shifted even for spin $j=0.1$. For $j=0.3$ it is shifted from $r_{\mathrm{ms}}\sim 5M$ for $q/j^{2}=1$ to $r_{\mathrm{ms}}\sim 5.8M$ for $q/j^{2}=10$, and in the extreme case of $j=0.5$, the shift is from $r_{\mathrm{ms}}\sim 4.4M$ for $q/j^{2}=1$ to $r_{\mathrm{ms}}\sim 6.6M$ for $q/j^{2}=10$. Clearly, the role of the quadrupole moment is most significant for the radial epicyclic motion where it affects also the values at the maxima of radial epicyclic frequency.

One can see that, contrary to the case of the Kerr black hole spacetime, the vertical epicyclic frequency can be larger than the orbital frequency. This phenomenon has been observed also for Newtonian quadrupole gravitational fields (Gondek-Rosińska et al., 2014) and general relativistic solutions with multipole structure of the Manko type (Pappas, 2012). A similar crossing of the radial profiles of the radial and vertical oscillations has been found also for axisymmetric string loops in the Kerr black hole spacetime (Stuchlík & Kološ, 2014, 2015). This type of behavior was discovered by Morsink & Stella (1999) around rapidly rotating neutron stars with a very stiff equation of state and was also found for analytical approximations to these spacetimes (Pappas, 2012; Pappas & Apostolatos, 2013). Recently it has been discussed by Kluźniak & Rosińska (2013) in the context of Newtonian gravity and then further considered by Gondek-Rosińska et al. (2014) for rapidly rotating relativistic strange stars. Its relation to the higher order multipole moments of spacetimes around compact objects can be found in Pappas (2017)).

Clearly, the quadrupole parameter $q/j^{2}$ is playing an important role and this needs to be investigated further. In the next subsections we will look in more detail at the value of the radial epicyclic frequency at its maximum.

3.2 Local extrema of the radial epicyclic frequency

Of the three frequencies describing orbital motion investigated above, the radial epicyclic frequency is the only one that is non-monotonic in $r/M$ for small $j$. We will now investigate further the radius at which $\nu_{\mathrm{r}}$ has its maximum and the values which it takes there. The location of the maximum can be found by solving the equation

[TABLE]

The solution can be found analytically, and is

[TABLE]

The radial epicyclic frequency at $r_{\mathrm{max}}$ is then given by

[TABLE]

or, in a simplified form,

[TABLE]

The left panel of Figure 5 shows the locations of the maximum $r_{\mathrm{max}}$ plotted as a function of $j$, and the right panel shows the values of $\nu_{\mathrm{r}}$ at $r_{\mathrm{max}}$. One can see that $q/j^{2}$ is having an impact rather similar to that in the case of the marginally stable circular orbit, with a minimum appearing in the curves for larger values of $q/j^{2}$ and moving to smaller values of $j$ as $q/j^{2}$ is increased. From the right panel, one sees that the largest values of $q/j^{2}$ may even cause $\nu_{\mathrm{r}}(r_{\mathrm{max}})$ to become a decreasing function of $j$ at high enough rotation speeds.

3.3 Precession frequencies

For matter orbiting around neutron stars, the quantities that may be directly observable could be related to combinations of the epicyclic frequencies, rather than to the frequencies themselves. In particular, the precession frequencies (periastron or nodal) may play an important role.

The periastron precession frequency is given by

[TABLE]

while the nodal precession frequency is given by

[TABLE]

with $\nu_{K}$, $\nu_{r}$ and $\nu_{\theta}$ being given by equations (5), (6) and (7).

For calculation of precession frequencies, we again chose $j=0.1,0.2,0.3,0.5$ and we varied the quadrupole parameter $q$ in the range $q/j^{2}=1,2,3,4,5,10$ for each value of $j$. Results are plotted in Figures 6 and 7. In both cases, the impact of $q/j^{2}$ becomes progressively more significant as $j$ increases. Note that the nodal precession frequency can become negative, since the radial epicyclic frequency can be higher than the orbital frequency. For low values of $j$, this happens only if $q/j^{2}$ is large and then only very near to $r=3M$, while for larger values of $j$ it happens also for smaller $q/j^{2}$.

3.4 Comparison of the Keplerian frequency and the vertical epicyclic frequency

As mentioned previously, within Newtonian theory it has been found by Gondek-Rosińska et al. (2014) that the vertical epicyclic frequency can be larger than the frequency of orbital motion if the central object is sufficiently oblate and we have seen that this can happen also for neutron stars (in general relativity) even if they are not very non-spherical. Here, we calculate the ratio of vertical epicyclic frequency to Keplerian frequency of orbital motion, and investigate when it is larger than one. The square of the ratio of these two frequencies can be written as

[TABLE]

where

[TABLE]

We plot the results in Figure 8, marking the cases where the Keplerian frequency is equal to the vertical epicyclic frequency. These points correspond to situations where the relativistic nodal frequency vanishes and changes sign.

The position where the nodal frequency vanishes needs to be compared with that of the marginally stable circular geodesic for the values of $j$ and $q/j^{2}$ being used, because the effect of nodal frequency switching is relevant only in the regions where the circular geodesic motion is stable against perturbations, i.e., outward of the marginally stable circular orbit at $r_{\mathrm{ms}}$. We have calculated the locations of $r_{\mathrm{ms}}$ using the relation (8), and these are shown as vertical lines in Figure 8 with the line-style corresponding to that for the $\nu_{\theta}/\nu_{\mathrm{K}}$ profile with same value of $q/j^{2}$. We can see some cases (especially for low values of $j$) where $\nu_{\theta}/\nu_{\mathrm{K}}=1$ occurs at a radius below $r_{\mathrm{ms}}$. On the other hand, for higher values of $j$ and $q/j^{2}$, $\nu_{\theta}>\nu_{\mathrm{K}}$ can occur also for orbits that are stable with respect to radial perturbations.

We now investigate the locations where

[TABLE]

Since $\nu_{\theta}$ is always equal to $\nu_{\mathrm{K}}$ in the non-rotating Schwarzschild case, one cannot find the location for the nodal switching point as an expansion around the Schwarzschild value. The values of $q/j^{2}$ satisfying (17) are given by

[TABLE]

and in Figure 9 we plot these against $r/M$ for $j=0.2,0.3,0.4,0.5$. The dotted lines show the value of $r_{\mathrm{ms}}$ for each $j$ and the purple dashed line indicates the stellar surface; we have marked the point on each of the $j$ curves beyond which $r$ is greater than both of those limits.

4 Models of twin HF QPOs and the frequency ratios in the Hartle–Thorne geometry

In the previous section we have investigated the behavior of the orbital, epicyclic and precession frequencies. We have seen that as $j$ increases, the role of the quadrupole term becomes increasingly important. In this section, we focus on particular combinations of frequencies that play key roles for models of QPOs - the twin peaks observed in power spectra of LMXBs. The origin of these peaks is still under investigation, but several candidate models are using frequencies associated with orbital motion around the central objects.

4.1 Frequency ratios of oscillatory modes for the twin HF QPO models

In this section we study the behavior of the frequencies of oscillatory modes given by some current models for the twin HF QPOs observed in atoll and Z-sources, where a neutron star is accreting matter from a low-mass companion (for the differences between atoll and Z-sources see Hasinger & van der Klis (1989)). We consider the specific models presented in Table 1, where we also give the combinations of frequencies that are attributed to the upper and lower observed peaks in each case. Our selection of models is based on Török et al. (2011). None of the models is currently uniquely preferred, because each of them has (different) theoretical problems. The RP (Relativistic Precession) model (Stella & Vietri, 1999) considers relativistic epicyclic motion of blobs at various radii in the inner parts of the accretion disc; The TD (Tidal Disruption) model (Čadež et al., 2008) proposes that the QPOs are generated by a tidal disruption of large accreting inhomogeneities; The WD (Warped Disk) model (Kato, 2001) considers a somewhat exotic disc geometry that causes a doubling of the observed lower QPO frequency; The ER (Epicyclic Resonance) model (Abramowicz & Kluźniak, 2001) involves different combinations of axisymmetric disk-oscillation modes; RP1 (Bursa, 2005) and RP2 (Török et al., 2007) use different combinations of non-axisymmetric disk-oscillation modes. We are considering here only geodesic oscillation models where the oscillatory modes are combinations of the orbital and epicyclic frequencies of near-circular geodesic motion in the equatorial plane.

In Figure 10 we show the radial profiles of the frequency ratio $\nu_{\mathrm{U}}/\nu_{\mathrm{L}}$ for all of the six QPO models being considered here.

As in the previous section, we use the values $j=0.1,0.2,0.3,0.5$, which includes the maximal range of the spin parameter $j$ for the Hartle–Thorne approximation (Urbanec et al., 2013), and for the quadrupole parameter $q$ we use values obtained from taking $q/j^{2}=1,3,5,10$, covering the most astrophysically relevant models for rotating neutron stars. It can be seen that for the RP models, the ratio $\nu_{\mathrm{U}}/\nu_{\mathrm{L}}=3/2$ occurs just slightly outside the location of the marginally stable orbit, meaning that the QPOs will generally be related to the region very close to the inner edge of the accretion disk, but for the TD and WD models the radius is larger and it is even more so for the ER model. However all of the considered models are giving the radius where $\nu_{\mathrm{U}}/\nu_{\mathrm{L}}=3/2$ as being between $r_{\mathrm{ms}}$ and $r=11M$.

We also note that the RP, RP1 and RP2 models give extremely similar profiles, with only very minor differences. The profiles of the WD and TD models are different, but they always cross each other at the frequency ratio $3:2$. For the ER model, the ratio diverges at the marginally stable orbit, since the radial epicyclic frequency has to vanish there, and this is the only case where the ratio decreases with increasing radius. In all of the models, the radial profiles depend only slightly on the quadrupole parameter for small values of $j<0.2$, but the dependence increases with increasing $j$, and for $0.3<j<0.5$ it can be very strong, especially close to the marginally stable orbit.

In Figure 11, we show these results plotted for each of the QPO models separately, using $j=0.1,0.2,0.3,0.5$, as before, and $q/j^{2}=1,2,3,4,5,10$. Here one can appreciate better the role of $j$ and $q/j^{2}$ for varying the position of marginally stable orbit. Note the qualitative difference between the results for the ER model and the others.

4.2 The frequency ratio $3:2$ in the twin HF QPO models

Observational data for QPOs demonstrates the importance of the frequency ratio $3:2$ (Boutelier et al., 2010) and so we have calculated how the location where this ratio is generated varies as a function of $j$. The left panels of Figure 12, show this for each of our QPO models for $q/j^{2}=1,2,3,4,5,10$. Note that the behavior of $r_{3:2}/M$ demonstrates a qualitatively similar behavior to that for the location of the marginally stable orbit: for small values of $q/j^{2}$, the resonance radius decreases with increasing spin, as for Kerr black holes, but for large enough $q/j^{2}$ it reaches a minimum and then increases again.

The right panels of Figure 12 show equivalent results for the ratio 5:4, which are quite similar.

5 Comparison of results from the Hartle–Thorne approach and the Lorene/nrotstar

numerical code

In this section we test the range of validity of the Hartle–Thorne approach by comparing results obtained with it against ones given by the publicly available Lorene/nrotstar code (Gourgoulhon et al., 2016). We compare models having the same gravitational mass for a static non-rotating star $M_{0}$ and then keep constant the value of the central pressure (and the energy or enthalpy density) when we move to rotating models. We use the APR equation of state (Akmal et al., 1998) for making the comparison.

We focus on the quantity $\Delta f_{\mathrm{ISCO}}=f_{\mathrm{ISCO}}-f_{\mathrm{ISCO}}(f=0Hz)$, where $f_{\mathrm{ISCO}}$ is here the frequency at the innermost stable orbit (either the marginally stable circular orbit at $r=r_{\mathrm{ms}}$ or the equatorial surface-skimming orbit at $r=R_{\mathrm{eq}}$ when $R_{\mathrm{eq}}>r_{\mathrm{ms}}$, as discussed previously (see also Pappas (2015) for discussion of when the surface of the neutron star lies above/below $r_{\mathrm{ms}}$).

This quantity is zero in the non-rotating case and its leading order in the slow-rotation expansion goes as $f_{\mathrm{rot}}$ while the order of its first neglected term goes as $\mathcal{O}(f^{3}_{\mathrm{rot}})$, i.e. the relative error is of order $\mathcal{O}(f^{2}_{\mathrm{rot}})$. This is also the case for $\nu^{\theta}_{\mathrm{P}}$, but the quantities $f_{\mathrm{ISCO}}$ and $\nu_{P}^{r}$ have non-zero values for the non-rotating case, i.e. the leading term is of $\mathcal{O}(f^{0}_{\mathrm{rot}})$, and since the first neglected term is again $\mathcal{O}(f^{3}_{\mathrm{rot}})$, this means that the relative error is of order $\mathcal{O}(f^{3}_{\mathrm{rot}})$, one order higher than that for $\nu^{\theta}_{\mathrm{P}}$. Because of this $\nu^{\theta}_{\mathrm{P}}$ is likely to be the quantity most liable to error and so it is sensible to use as our key test quantity $\Delta f_{\mathrm{ISCO}}=f_{\mathrm{ISCO}}-f_{\mathrm{ISCO}}(f=0Hz)$ which has an error of the same order as the most sensitive quantity of interest. This is not a foolproof argument because of different coefficients in the expansions of the quantities concerned, but it should be a reasonable guide.

For the ISCO frequency in the Hartle–Thorne approach, we calculate a gravitational mass $M$ for the rotating star, a dimensionless angular momentum $J/M_{0}^{2}$ and a dimensionless quadrupole moment $Q/M_{0}^{3}$, use these to calculate the radius $r_{\mathrm{ms}}$ from eq. (8) and then find the Keplerian frequency at this radius using eq. (5). If $r_{\mathrm{ms}}$ is less than $R_{\mathrm{eq}}$, then the latter is used instead in eq. (5). For the gravitational mass of the non-rotating star we take $M_{0}=1.2,~{}1.4,~{}1.6,~{}1.8M_{\odot}$. The rotational frequencies are taken up to the frequency where the difference between the Hartle–Thorne and Lorene values for $\Delta f_{\mathrm{ISCO}}$ reaches $\sim 10\%$. The results of the comparison are shown in Figures 13 and 14, where we show also the 5% and 10% intervals around the Lorene values. The positive gradients correspond to $r_{\mathrm{ms}}$ being used, and the negative gradients to $R_{\mathrm{eq}}$. Since we are here dealing with specific models, the location of the stellar surface is directly known.

The two approaches are in excellent agreement for the lower rotation frequencies, as expected; when differences begin to emerge, Hartle–Thorne is giving higher values for $\Delta f_{\mathrm{ISCO}}$ in all cases. We see that for lower masses the difference between the Lorene and Hartle–Thorne values of $\Delta f_{\mathrm{ISCO}}$ are smaller than 10% up to frequencies of about 500Hz, while for higher masses the rotational frequencies can go up to 700Hz (for $M=1.6M_{\odot}$) or 1000Hz (for $M=1.8M_{\odot}$). Since it is generally believed that neutron stars in QPO sources have accreted matter from the companion star they should have rather high masses.

6 Summary

We have examined the behavior of the orbital frequency, the radial epicyclic frequency and the vertical epicyclic frequency in the external Hartle–Thorne geometry for rotating compact stars. The primary parameters of the geometry are the gravitational mass $M$, the angular momentum $J$ and the quadrupole moment $Q$ and we have calculated the corresponding values for these from models of neutron stars constructed using realistic equations of state for nuclear matter. Our approach can be applied for any future equation of state if the universal relations are used as described in Section 2666During the refereeing process and preparation of the revised version of this manuscript, a paper by Luk & Lin (2018) has been published where it is demonstrated that the frequency at the marginally stable orbit around rotating neutron stars is also a universal function of the neutron star properties.. Apart from the individual frequencies, we have also investigated the behavior of combinations of them which are involved in some models for QPOs observed in X-rays from LMXBs.

We have seen that the dimensionless angular momentum $j$ and dimensionless quadrupole moment $q$ play an important role for the behavior of the orbital and epicyclic frequencies, especially for the ratio between the vertical epicyclic frequency and the Keplerian frequency, becoming particularly important for stars rotating with the higher angular velocities. Our investigation was motivated by the observations of QPOs in the X-ray spectra of LMXBs and we have investigated the behavior of the frequencies relevant for several of the models proposed for explaining these.

It may seem redundant to spend time with more models, but at present we do not know which model may be relevant for the sources of QPOs and are not even sure whether only one model will be applicable for all of them. If the radiation does come from the accretion disc and is modulated by its oscillations, different oscillation modes may well be involved for different sources, depending on the geometry of the binary system. For example, if the equatorial plane of the compact object is the same as the orbital plane of the binary system, one may expect radial perturbations of the disc caused by matter falling in from the donor star. On the other hand, if the equatorial plane of the compact star is almost perpendicular to the orbital plane of the binary, one may expect vertical perturbations of the disk instead. This is something which may well clarify in the near future.

Finally we note that the results presented here can be relevant not only for the particular models mentioned in Section 4 but also for a wider class of orbital models. For instance, the maximum of the radial epicyclic frequency determines propagation of some discoseismic modes (Kato et al., 1998; Wagoner et al., 2001) and so there would be an application for making estimates of stellar spin rates based on association of QPOs with oscillations trapped near to the inner edge of thin accretion discs (Kato, 1989, 2001). This is illustrated in Figure 15 where we show the dependence of the maximum value of the radial epicyclic frequency $\nu_{\mathrm{r}}(r_{\mathrm{max}})$ for compact and highly oblate stars. Similar examples can be found for the vertical epicyclic frequency and models dealing with the Lense-Thirring precession effect (Ingram et al., 2009; Rosińska et al., 2014; Wiśniewicz et al., 2015; Ingram et al., 2015; Tsang & Pappas, 2016). The formulas given in Section 3 for determining the epicyclic frequencies (and the ratio between the vertical epicyclic frequency and the Keplerian frequency) can therefore be used in various applications, especially within the context of the universal relations.

GU and ZS acknowledge the Albert Einstein Center for Gravitation and Astrophysics supported by the Czech Science Foundation grant No. 14-37086G. MU and GT were supported by Czech Science Foundation grant No.17-16287S. GU, MU and ZS are grateful for very kind hospitality from the University of Oxford and for fruitful discussions there; the stay in Oxford was partly supported by grant LTC18058. They acknowledge the internal grants of the Silesian University Opava FPF SGS/14,15/2016. Some of the equations of state used in this work were kindly provided to the authors by J.R. Stone. The authors also acknowledge fruitful discussions with Marek A. Abramowicz concerning particular aspects of models of rotating stars and properties of their space-times. We would like to thank the referee for very helpful comments and suggestions and for pointing out several references that we missed in the previous version of our manuscript. ).