Single compact stars

5.2 Single compact stars

The computation of the structure of stationary compact stars dates back to 1939 with the famous solution of Tolman–Oppenheimer–Volkoff. During the last few years, the need for accurate models has become more pressing, especially with the coming online of gravitational wave detectors, which could help to probe the interior of such compact stars. Isolated stars in rotation are essentially axisymmetric, but some physical effects can induce a symmetry breaking that could lead to the emission of gravitational waves. In the following, we will review some computations that aim at including some of these effects, such as spontaneous symmetry breaking, the inclusion of magnetic fields, the effects of exotic dense matter (mainly with strange quarks) and the influence of an interior composed of two different superfluids.

5.2.1 Formalisms

The first computation of models of relativistic rotating stars in general relativity, by means of spectral methods, is presented in [41 ]. The equations are solved in spherical coordinates (see Section 3.2). Doing so, the fields only depend on the azimuthal angle $𝜃$ and on the radius $r$ . The fields are expanded in terms of spherical harmonics with respect to the angle and in Chebyshev polynomials with respect to $r$ . The use of spherical harmonics is a natural way of dealing with the coordinate singularity on the $z$ -axis. In [41 ] the whole space is divided into two spherical domains, the outer one extending up to infinity by making use of the compactification in $1∕r$ seen in Section 3.1.2. With such a setting, Einstein’s equations reduce to a set of four elliptic equations with sources extending up to infinity that are solved using a version of the algorithm based on matching with the homogeneous solutions (presented in Section 2.6.4) for each spherical harmonic. The system is complete once a description of the matter is given. The simplest choice is to consider a polytropic fluid, with or without a magnetic field. The system is solved by iteration.

In [41 ], a particular emphasis is put on various methods to measure the accuracy of the code. For nonrotating stars, the error is found to decrease exponentially, as the number of coefficients increases (see Figures 5 and 6 of [41 ]). However, for fast-rotating configurations, the error only decays as a power law (see Figure 7 of [41 ]). This comes from the fact that quantities like the energy density are no longer $𝒞∞$ across the star’s surface. Nevertheless, the results are in good agreement (to 0.1%) with those obtained by other numerical methods, as can be seen in [161 ].

Spectral convergence can be recovered by using surface-adapted coordinates as first done in [36 ]. A regular mapping of the numerical coordinates to the physical ones is introduced, so that the surface of the star lies at the boundary between two domains (see Section 3.1.1). For polytropes with $γ < 2$ , this is sufficient to recover spectral convergence (see Figures 5 and 6 of [38 ]). However, for $γ > 2$ , some quantities are still diverging at the surface but the convergence can be made closer and closer to the spectral one by analytically regularizing the density (see Section IV of [38 ]). Doing so, the error decreases as a power law, but the decrease can be made arbitrary fast at the cost of increasing the number of operations and so the computational time.

Up until 2006, neutron stars were computed using quasi-isotropic coordinates. However, in order to use these configurations as initial data for evolutionary codes, it may be useful to allow for other choices. Among the possible gauges, Dirac’s is one of the most promising [37 ]. In [134 ] models of rotating neutron stars in the Dirac gauge are computed for both polytropic and realistic equations of state. Contrary to the quasi-isotropic coordinates, the use of this new gauge implies that one must solve one tensor-like Poisson equation. Configurations obtained with the two different formalisms are shown to be in very good agreement.

5.2.2 Rotating neutron star models

Even before adapted mappings were available, interesting results could be obtained. In [184 , 185 ], models of rotating neutron stars with various equations of state have been computed. Among the most surprising findings is the existence of supra-massive stars. These stars do not connect to the nonrotating limit. Indeed, their high mass can only be supported by the presence of a centrifugal force. One of the remarkable features of such stars is the fact that they actually spin up when they lose angular momentum, contrary to what is observed for normal stars. This effect can also be seen for neutron stars containing hyperons and, thus, a softer equation of state [238 ]. Let us mention that, in this case, the stability analysis of the configurations required the great precision that spectral methods with adapted coordinates could provide.

It is known that isolated pulsars spin down due to magnetic braking. As the rotational frequency decreases, it is possible that the star will encounter a transition from one state of matter to another. Stationary rotating models have been used to determine the properties of such transitions [231 , 232 , 233 ]. A puzzling result is that the amount of energy released in a first-order phase transition does not depend on the orbital velocity of the star and is the same as for nonrotating ones. This is shown to be the case for both weak [232 ] and strong [233 ] first-order transitions.

5.2.3 Spontaneous symmetry breaking

It is known that stars can undergo spontaneous symmetry breaking when rotating fast enough. When such a phenomenon occurs, triaxial configurations are formed that are potential emitters of gravitational waves. The departure from axisymmetry is studied in two papers by the Meudon group [34 , 35 ]. The idea of the method is to start from an axisymmetric neutron star configuration and to follow the growth or decay of triaxial instabilities. This work reaffirms well-established results from the Newtonian regime and presents the first results in general relativity for various equations of states. For a few of them, the frequency at which symmetry breaking occurs lies in the frequency band of the LIGO and Virgo detectors.

In 2002, this work was extended [89 ] by making use of surface-fitting coordinates. This enabled the authors to obtain results in the incompressible case by properly dealing with discontinuities lying at the surface of the star. Classical results in the incompressible case are, thus, recovered and it is found that the inclusion of relativity has only a moderate effect. Indeed, the critical ratio between the kinetic energy and the absolute gravitational energy $T ∕|W |$ at which the triaxial instability occurs is only 30% larger for relativistic stars, with respect to their classical counterparts.

If relativistic effects only slightly stabilize stars, the same is not true for differential rotation. Indeed, in [183 ], the authors study various rotation profiles and equations of state using the same technique as in [34 , 35 ] to determine the onset of instability. It appears that the critical value of $T∕ |W |$ can be almost twice as high as for uniformly rotating stars.

5.2.4 Configurations with magnetic fields

Even if magnetic fields are strong in neutron stars, the structure of the object is not affected until it reaches huge values, on the order of at least 10¹⁰ T. In [32 ], models of rapidly-rotating stars with poloidal fields are constructed for various equations of state. The magnetic effects are taken into account consistently by solving the appropriate Maxwell equations, as well as by means of spectral methods. The maximum mass of highly-magnetized neutrons stars is found to be higher by 13 – 29% than for nonmagnetized stars. The magnetic field induces an additional pressure, which can help to support more massive stars, thus explaining this increase.

The presence of a magnetic field can also lead to the deformation of a neutron star. Such deformation could lead to the formation of a triaxial configuration, which would then emit gravitational waves. In [36 ] the emitted signal is computed. Typically the system radiates at two frequencies: $Ω$ and $2Ω$ where $Ω$ is the angular velocity of the star.

In more recent work by the Meudon group [160 ], magnetized configurations have been computed using coordinates matched to the surface of the star, thus making the computation much more accurate. Gyromagnetic ratios of rapidly-rotating neutron stars of various equations of state are obtained. The limit of a ratio of $g = 2$ value for a charged black hole is never reached.

5.2.5 Strange stars

It is possible that the fundamental state of nuclear matter is not the ordinary matter but rather a plasma of deconfined quarks $u$ , $d$ and $s$ , called strange matter. If this is the case, neutron stars should rather be called strange stars. The main difference between these two types of compact stars is that strange ones are somewhat smaller and thus more compact. In particular, they support higher rotation rates. There is a strong density jump at the surface of a strange star and surface-fitting coordinates are required in order to deal with it.

Fast rotating models of strange stars are computed in [103 , 88 ]. Due to higher compactness, it is found that strange stars can rotate significantly faster than their neutron star counterparts. The attained $T ∕|W |$ can be twice as large. As in the neutron star case, supermassive configurations that spin up with angular momentum loss are found. The influence of strange matter on the emission of gravitational waves is studied in [90 ] where viscosity effects and triaxial instabilities are carefully taken into account.

It is believed that millisecond pulsars have been spun up by accreting matter from a companion. However, the details of this mechanism depend on the nature of the compact object. In [237 ], the differences between accretion onto a neutron star and onto a strange star are investigated, using 2D stationary models computed by spectral methods.

5.2.6 Quasiperiodic oscillations

Quasiperiodic oscillations (QPOs) are observed in the kHz regime and are believed to be the signature of matter falling into a compact object. In the standard picture, the frequency of the QPOs is that of the last stable orbit around the compact object. Let us mention that the presence of a last stable orbit around an extended body is not an effect of relativity but can also be seen in the Newtonian regime, as shown in [234 ].

The precise structure of the accreting object has a great influence on the QPO. In a series of papers [235 , 92 , 3 , 27 ], comparisons are made between observations and various compact star models that could account for QPOs.

Using a multidomain approach, strange stars with a crust can also be computed [236 ], one domain describing the interior of the star and another the crust. It is shown that the presence of the crust could change the value of the QPO by up to 20%.

5.2.7 More complex configurations

In this section, objects in more exotic configurations are presented. This is an illustration of both the complexity of neutron stars physics and the ability of spectral methods to deal with complicated systems.

The observation of glitches in isolated pulsars is consistent with the presence of a superfluid interior. The simplest model considers two fluids, one composed of neutrons and the other of protons and electrons, both components being superfluids. However, these two components could have different velocities, in particular different rotation rates. Such configurations are computed in [178 ]. A multidomain setting is crucial to be able to follow the two different fluids because the components do not have the same shape. Among the various results obtained, let us mention the confirmation of the existence of prolate-oblate configurations.

Neutron stars are usually assumed to be at zero temperature. However, this approximation is known to no longer be true for newborn neutron stars just after the supernova. Models of newborn neutron stars in uniform rotations are constructed in [106 ] using an extension of the code developed in [41 ]. Various hypothesis about the interior (different lepton numbers, isothermal versus isentropic) are considered. Sequences of fixed baryon mass and angular momentum are constructed. Such sequences approximate the evolution of the protoneutron star into a cold neutron star. The results have been extended to differentially-rotating protoneutron stars in [107 ].

The effect of finite temperature is also investigated in [226 ]. The authors found that newborn neutron stars are unlikely to undergo bar mode instability, but that the secular ones could take place and result in a significant emission of gravitational waves. Another interesting result of [226 ] is the existence of toroidal-like configurations, which appear for a broad range of parameters and before the mass-shedding limit. In such cases, the surface of the star is highly deformed and surface-fitting coordinates are required.

Axisymmetric rotating neutron stars have also been computed by a code developed by Ansorg et al. [9 , 10 ] . This code is based on Lewis–Papapetrou coordinates $(ρ,ξ)$ , which are closely related to the usual cylindrical coordinates. Typically, space is decomposed into two domains: one for the interior of the star and another for the exterior, which extends up to spatial infinity. Compactification of space and surface-fitting mappings are used. Both variables are expanded in terms of Chebyshev polynomials. Instead of solving the equations harmonic by harmonic and then iterating, as is done by the Meudon group, the equations are written with a collocation method (see Section 2.5.3) and solved as one single system. The price to pay is that the size of the system is somewhat larger (i.e. in $2 m$ , $m$ being the total number of coefficients for each coordinate). The system is solved by means of the Newton–Raphson method. At each step, the linear system is solved using iterative techniques with preconditioning. With this method, impressive accuracy is reached.

The coordinates used in [9 , 10 ] are more general than the ones used by the Meudon group, especially with respect to their surface-fitting capabilities. They can account for more complicated configurations and, in particular, highly-distorted matter distribution can be obtained. This is shown in [12 , 11 ], where relativistic axisymmetric toroidal configurations of matter, known as the Dyson rings, are computed. Such rings are obtained up to the mass-shedding limit. Transition to the limit of an extreme Kerr black hole is also discussed.