ON GENERAL RELATIVISTIC UNIFORMLY ROTATING WHITE DWARFS

The relevance of rotation in enhancing the maximum stable mass of a white dwarf (WD) has been discussed for many years, both for uniform rotation (see, e.g., James 1964; Anand 1965; Roxburgh & Durney 1966; Monaghan 1966; Geroyannis & Hadjopoulos 1989) and differential rotation (see, e.g., Ostriker & Bodenheimer 1968; Ostriker & Tassoul 1969; Tassoul & Ostriker 1970; Durisen 1975). Newtonian gravity and post-Newtonian approximation have mainly been used to compute the structure of the star, with the exception of the work of Arutyunyan et al. (1971), where rotating white dwarfs (RWDs) were computed in full general relativity (GR). From the microscopical point of view, the equation of state (EOS) of cold WD matter has been assumed to be either a microscopically uniform degenerate electron fluid as used by Chandrasekhar (1931) in his classic work, or to have a polytropic form.

However, as first shown by Salpeter (1961) in the Newtonian case and then by Rotondo et al. (2011a, 2011b) in GR, a detailed description of the EOS taking into account the effects of the Coulomb interaction is essential for the determination of the maximum stable mass of non-rotating WDs. The specific microphysics of the ion–electron system forming a Coulomb lattice, together with the detailed computation of the inverse β-decays and the pycnonuclear reaction rates, play a fundamental role.

A new EOS taking into account the finite size of the nucleus, the Coulomb interactions, and the electroweak equilibrium in a self-consistent relativistic fashion has recently been obtained by Rotondo et al. (2011b). This relativistic Feynman–Metropolis–Teller (RFMT) EOS generalizes both the Chandrasekhar (1931) and Salpeter (1961) works using a full treatment of the Coulomb interaction given through the solution of a relativistic Thomas–Fermi model. This leads to a more accurate calculation of the energy and pressure of the Wigner–Seitz cells, and hence a more accurate EOS. It has been shown how the Salpeter EOS overestimates at high densities and underestimates at low densities the electron pressure. The application of this new EOS to the structure of non-rotating ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe was recently done in Rotondo et al. (2011a). The new mass–radius relations generalize the works of Chandrasekhar (1931) and Hamada & Salpeter (1961); smaller maximum masses and a larger minimum radii are obtained. Both GR and inverse β-decay can be relevant for the instability of non-rotating WDs depending on the nuclear composition, as we can see from Table 1, which summarizes some of the results of Rotondo et al. (2011a).

Here, we extend the previous results of Rotondo et al. (2011a) for uniformly RWDs at zero temperatures that obey the RFMT EOS. We use Hartle's approach (Hartle 1967) to solve the Einstein equations accurately up to a second-order approximation of the angular velocity of the star. We calculate the mass M, equatorial R_eq and polar R_p radii, angular momentum J, eccentricity , and quadrupole moment Q, as a function of the central density ρ_c and rotation angular velocity Ω of the WD. We also construct RWD models for the Chandrasekhar and Salpeter EOS and compare and contrast the differences with the RFMT ones.

We analyze in detail the stability of RWDs from both the microscopic and macroscopic point of view in Section 3. Besides the inverse β-decay instability, we also study the limits to the matter density imposed by zero-temperature pycnonuclear fusion reactions using up-to-date theoretical models (Gasques et al. 2005; Yakovlev et al. 2006). We calculate the mass-shedding limit as well as the secular axisymmetric instability boundary.

The general structure and stability boundaries of ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe WDs are discussed in Section 4. From the maximally rotating models (mass-shedding sequence), we calculate in Section 5 the maximum mass of uniformly rotating ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe WDs for the Chandrasekhar, Salpeter, and RFMT EOS, and compare the results with the existing values in the literature. We calculate the minimum(maximum) rotation period(frequency) of an RWD for the above nuclear compositions, taking into account both inverse β-decay and pycnonuclear restrictions to the density; see Section 6.

We discuss in Section 7 the axisymmetric instabilities found in this work. A comparison of Newtonian and general relativistic WDs presented in Appendix C show that this is indeed a general relativistic effect. Furthermore, we estimate in Appendix D the accuracy of the "slow" rotation approximation (power-series solutions up to order Ω²) for the determination of the maximally rotating sequence of WDs. In this manner, we calculate the rotation to gravitational energy ratio and the deviations from spherical symmetry.

In addition, we construct in Section 8 constant rest-mass evolution tracks of RWDs at fixed chemical composition and show that RWDs may experience both spin-up and spin-down epochs while losing angular momentum, depending on their initial mass and rotation period.

Finally, in Section 9, we outline some of the astrophysical implications of the results presented in this work, which we summarize in Section 10.

Hartle (1967) was the first to describe the structure of rotating objects approximately up to second-order terms in the angular velocity of the star Ω, within GR. In this "slow" rotation approximation, the solution of the Einstein equations in the exterior vacuum can be written in analytic closed form in terms of the total mass M, angular momentum J, and quadrupole moment Q of the star (see Appendix A). The interior metric is constructed by solving numerically a system of ordinary differential equations for the perturbation functions (see Hartle 1967; Hartle & Thorne 1968, for details).

The spacetime geometry up to the order of Ω², with an appropriate choice of coordinates, in geometrical units c = G = 1, is described by (Hartle 1967)

$$\begin{eqnarray} &ds^{2}= \left\lbrace e^{\nu (r)}[1+2 h_0(r)+2 h_2(r)P_2(\cos \theta)] -\omega ^2 r^2 \sin ^2\theta \right\rbrace dt^{2} \nonumber \\ &+2 \omega r^2 \sin ^2\theta dt d\phi -e^{\lambda (r)}\left[\!1+2 \frac{m_0(r)+m_2(r)P_2(\cos \theta)}{r-M^{J=0}(r)}\!\right]dr^2 \nonumber \\ &-r^2\left[1+2k_2(r)P_2(\cos \theta)\right](d\theta ^2+\sin ^2\theta d\phi ^2), \end{eqnarray} \tag{ 1 }$$

where P₂(cos θ) is the second order Legendre polynomial, e^ν(r) and e^λ(r) = [1 − 2M^{J = 0}(r)/r]⁻¹, and M^{J = 0}(r) are the metric functions and mass of the corresponding static (non-rotating) solution with the same central density as the rotating one. The angular velocity of the local inertial frames ω(r), proportional to Ω, as well as the functions h₀, h₂, m₀, m₂, k₂, proportional to Ω², must be calculated from the Einstein equations (see Hartle 1967; Hartle & Thorne 1968, for details); their analytic expressions in the vacuum case can be found in Appendix A.

The parameters M, J, and Q, are then obtained for a given EOS from the matching procedure between the internal and external solutions at the surface of the rotating star. The total mass is defined by M = M^{J ≠ 0} = M^{J = 0} + δM, where M^{J = 0} is the mass of a static (non-rotating) WD with the same central density as M^{J ≠ 0}, and δM is the contribution to the mass due to rotation.

3.1. The Mass-shedding Limit

The velocity of particles on the equator of the star cannot exceed the Keplerian velocity of "free" particles, computed at the same location. At this limit, particles on the star's surface remain bound to the star only because of a balance between gravity and centrifugal forces. The evolution of a star rotating at this Keplerian rate is accompanied by a loss of mass, thus becoming unstable (see, e.g., Stergioulas 2003, for details). A procedure to obtain the maximum possible angular velocity of the star before reaching this limit was developed, e.g., by Friedman et al. (1986). However, in practice, it is less complicated to compute the mass-shedding (or Keplerian) angular velocity of a rotating star, Ω^{J ≠ 0}_K, by calculating the orbital angular velocity of a test particle in the external field of the star and corotating with it at its equatorial radius, r = R_eq.

For the Hartle–Thorne (HT) external solution, the Keplerian angular velocity can be written as (see, e.g., Torok et al. 2008, and Appendix A.2 for details)

$$\begin{equation} \Omega ^{J \ne 0}_{K}= \left(G \frac{M}{R_{{\rm eq}}^3}\right)^{1/2}\left[1- j F_{1}(R_{{\rm eq}})+j^2F_{2}(R_{{\rm eq}})+q F_{3}(R_{{\rm eq}})\right], \end{equation} \tag{ 2 }$$

where j = cJ/(GM²) and q = c⁴Q/(G²M³) are the dimensionless angular momentum and quadrupole moment, and the functions F_i(r) are defined in Appendix A.2. Thus, the numerical value of Ω^{J ≠ 0}_K can be computed by gradually increasing the value of the angular velocity of the star, Ω, until it reaches the value Ω^{J ≠ 0}_K expressed by Equation (2).

It is important to analyze the accuracy of the slow rotation approximations, e.g., accurate up to second order in the rotation expansion parameter, for the description of maximally rotating stars as WDs and neutron stars (NSs). We have performed in Appendix D a scrutiny of the actual physical request made by the slow rotation regime. Based on this analysis, we have checked that the accuracy of the slow rotation approximation increases with the density of the WD, and that the mass-shedding (Keplerian) sequence of RWDs can be accurately described by the Ω² approximation within an error smaller than the one found for rapidly rotating NSs, ≲ 6%.

3.2. The Turning-point Criterion and Secular Axisymmetric Instability

During a period when the central density increases, the mass of the non-rotating star is limited by the first maximum of the M–ρ_c curve, i.e., the turning point given by the maximum mass, ∂M/∂ρ_c = 0, marks the secular instability point and also coincides with the dynamical instability point if the perturbation obeys the same EOS as that of the equilibrium configuration (see, e.g., Shapiro & Teukolsky 1983, for details). The situation, however, is much more complicated in the case of rotating stars; the determination of axisymmetric dynamical instability points to finding the perturbed solutions with zero frequency modes, that is, perturbed configurations whose energy (mass) is the same as the unperturbed (equilibrium) one, at second order. However, Friedman et al. (1988) formulated, based on the works of Sorkin (1981, 1982), a turning-point method to locate the points where secular instability sets in for uniformly rotating relativistic stars: along a sequence of rotating stars with fixed angular momentum and increasing central density, the onset of secular axisymmetric instability is given by

$$\begin{equation} \left(\frac{\partial M(\rho _c, J)}{\partial \rho _c}\right)_{J}=0\,. \end{equation} \tag{ 3 }$$

Thus, the configurations on the right side of the maximum mass of a J-constant sequence are secularly unstable. After the secular instability sets in, the configuration evolves quasi-stationarily until it reaches a point of dynamical instability where gravitational collapse should take place (see Stergioulas 2003 and references therein). The secular instability boundary thus separates stable from unstable stars. It is worth stressing here that the turning point of a constant J sequence is a sufficient but not a necessary condition for secular instability, and therefore it establishes an absolute upper bound for the mass (at constant J). We construct the boundary given by the turning points of constant angular momentum sequences as given by Equation (3). The question of whether or not dynamically unstable RWDs can exist on the left side of the turning-point boundary remains an interesting problem and deserves further attention in view of the very recent results obtained by Takami et al. (2011) for some models of rapidly rotating NSs.

3.3. Inverse β-decay Instability

It is known that a WD might become unstable against the inverse β-decay process (Z, A) → (Z − 1, A) through the capture of energetic electrons. In order to trigger such a process, the electron Fermi energy (with the rest mass subtracted off) must be larger than the mass difference between the initial (Z, A) and final (Z − 1, A) nucleus. We denote this threshold energy as ^β_Z. Usually, it is satisfied by ^β_{Z − 1} < ^β_Z, and therefore the initial nucleus undergoes two successive decays, i.e., (Z, A) → (Z − 1, A) → (Z − 2, A) (see, e.g., Salpeter 1961; Shapiro & Teukolsky 1983). Some of the possible decay channels in WDs with the corresponding known experimental threshold energies ^β_Z are listed in Table 2. The electrons in the WD may eventually reach the threshold energy to trigger a given decay at some critical density ρ^β_crit. Since the electrons are responsible for the internal pressure of the WD, configurations with ρ > ρ^β_crit become unstable due to the softening of the EOS as a result of the electron capture process (see Salpeter 1961 for details). In Table 2, the critical density ρ^β_crit given by the RFMT EOS is shown corresponding to each threshold energy ^β_Z; see Rotondo et al. (2011a) for details.

Table 2. Onset for the Inverse β-decay of ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe

Decay	βZ	ρβcrit
	(MeV)	(g cm−3)
4He → 3 H + n → 4n	20.596	1.39 × 1011
12C → 12B → 12Be	13.370	3.97 × 1010
16O → 16N → 16C	10.419	1.94 × 1010
56Fe → 56Mn → 56Cr	3.695	1.18 × 109

Notes. The experimental values of the threshold energies ^β_Z have been taken from Table 1 of Audi et al. (2003); see also Wapstra & Bos (1977) and Shapiro & Teukolsky (1983). The corresponding critical density ρ^β_crit are for the RFMT EOS (see Rotondo et al. 2011a).

Download table as:  ASCIITypeset image

3.4. Pycnonuclear Fusion Reactions

In our WD model, we assume a unique nuclear composition (Z, A) throughout the star. We have just seen that inverse β-decay imposes a limit to the density of the WD over which the current nuclear composition changes from (Z, A) to (Z − 1, A). There is an additional limit to the nuclear composition of a WD. Nuclear reactions proceed when the nuclei in the lattice overcome the Coulomb barrier. In the present case of zero temperatures T = 0, the Coulomb barrier can be overcome because of the zero-point energy of the nuclei (see, e.g., Shapiro & Teukolsky 1983)

$$\begin{equation} E_p=\hbar \omega _p,\quad \omega _p=\left({\frac{4 \pi e^2 Z^2 \rho }{A^2 M_u^2}}\right)^{1/2}, \end{equation} \tag{ 4 }$$

where e is the fundamental charge and M_u = 1.6605 × 10⁻²⁴ g is the atomic mass unit.

Based on the pycnonuclear rates computed by Zeldovich (1958) and Cameron (1959), Salpeter (1961) estimated that in a time of 0.1 Myr, ¹H is converted into ⁴He at ρ ∼ 5 × 10⁴ g cm⁻³, ⁴He into ¹²C at ρ ∼ 8 × 10⁸ g cm⁻³, and ¹²C into ²⁴Mg at ρ ∼ 6 × 10⁹ g cm⁻³. The threshold density for the pycnonuclear fusion of ¹⁶O occurs, for the same reaction time of 0.1 Myr, at ρ ∼ 3 × 10¹¹ g cm⁻³, and for 10 Gyr at ∼10¹¹ g cm⁻³. These densities are much higher than the corresponding density for inverse β-decay of ¹⁶O, ρ ∼ 1.9 × 10¹⁰ g cm⁻³ (see Table 2). The same argument applies to heavier compositions, e.g., ⁵⁶Fe, so that pycnonuclear reactions are not important for those heavier than ¹²C in WDs.

It is important to analyze the case of ⁴He WDs in detail. At densities ρ_pyc ∼ 8 × 10⁸ g cm⁻³ an ⁴He WD should have a mass M ∼ 1.35 M_☉ (see, e.g., Figure 3 in Rotondo et al. 2011a). However, the mass of ⁴He WDs is constrained to lower values from their previous thermonuclear evolution: a cold star with mass >0.5 M_☉ has already burned an appreciable part of its helium content at earlier stages. Thus, WDs of M > 0.5 M_☉ with ⁴He cores are very unlikely (see Hamada & Salpeter 1961, for details). It should be stressed that ⁴He WDs with M ≲ 0.5 M_☉ have central densities ρ ∼ 10⁶ g cm⁻³ (Rotondo et al. 2011a) and at such densities pycnonuclear reaction times are longer than 10 Gyr, and hence are unimportant. However, we construct in this work ⁴He RWDs configurations all the way up to their inverse β-decay limiting density for the sake of completeness, keeping in mind that the theoretical ⁴He WDs configurations with M ≳ 0.5 M_☉ could actually not be present in any astrophysical system.

From the above discussion, we conclude that pycnonuclear reactions can be relevant only for ¹²C WDs. It is important to stress here that the reason the pycnonuclear reaction time, τ^{C + C}_pyc, determines the lifetime of a ¹²C WD is that reaction times τ^{C + C}_pyc < 10 Gyr are achieved at densities ∼10¹⁰ g cm⁻³, lower than the inverse β decay threshold density of ²⁴Mg, ²⁴Mg→ ²⁴Na→ ²⁴Ne, ρ ∼ 3.2 × 10⁹ g cm⁻³ (see, e.g., Salpeter 1961; Shapiro & Teukolsky 1983). Thus, the pycnonuclear ¹²C+¹²C fusion produces unstable ²⁴Mg that almost instantaneously decays due to electron captures, and so the WD becomes unstable as we discussed in Section 3.1.

However, the pycnonuclear reaction rates are not known with precision due to theoretical and experimental uncertainties. Hamada & Salpeter (1961) had already pointed out in their work that the above pycnonuclear density thresholds are reliable only within a factor of three or four. The uncertainties are related to the precise knowledge of the Coulomb tunneling in the high-density, low-temperature regime relevant to astrophysical systems, e.g., WDs and NSs, as well as with the precise structure of the lattice; impurities, crystal imperfections, and the inhomogeneities of the local electron distribution and finite temperature effects, also affect the reaction rates. The energies for which the so-called astrophysical S-factors are known from experiments are larger with respect to the energies found in WD and NS crusts, and therefore the value of the S-factors have to be obtained theoretically from the extrapolation of experimental values using appropriate nuclear models, which at the same time are poorly constrained. A detailed comparison between the different theoretical methods and approximations used for the computation of the pycnonuclear reaction rates can be found in Gasques et al. (2005) and Yakovlev et al. (2006).

The S-factors have been computed in Gasques et al. (2005) and Yakovlev et al. (2006) using up-to-date nuclear models. Following these works, we have computed the pycnonuclear reaction times for C+C fusion as a function of the density as given by Equation (B3), τ^{C + C}_pyc, which we show in Figure 1; we refer to Appendix B for details.

Zoom In Zoom Out Reset image size

We determine that for τ^{C + C}_pyc = 10 Gyr, ρ_pyc ∼ 9.26 × 10⁹ g cm⁻³, while for τ^{C + C}_pyc = 0.1 Myr, ρ_pyc ∼ 1.59 × 10¹⁰ g cm⁻³, to be compared with the value ρ ∼ 6 × 10⁹ g cm⁻³ estimated by Salpeter (1961). In order to compare the threshold densities for inverse β-decay and pycnonuclear fusion rates, we shall indicate in our mass–density and mass–radius relations the above two density values corresponding to these two lifetimes. It is important to stress that the computation of the pycnonuclear reactions rates is subject to theoretical and experimental uncertainties (see Gasques et al. 2005, for details). For instance, Hamada & Salpeter (1961) stated that these pycnonuclear critical densities are reliable within a factor of three or four. If three times larger, then the above value of ρ_pyc for τ^{C + C}_pyc = 0.1 Myr becomes ρ_pyc ∼ 4.8 × 10¹⁰ g cm⁻³, larger than the inverse β-decay threshold density ρ^C_β ∼ 3.97 × 10¹⁰ g cm⁻³ (see Table 2). As we will see in Section 7, the turning-point construction leads to an axisymmetric instability boundary in the density range ρ^{C, J = 0}_crit = 2.12 × 10¹⁰ < ρ < ρ^C_β g cm⁻³ in a specific range of angular velocities. This range of densities is particularly close to the above values of ρ_pyc, which suggests a possible competition between different instabilities at high densities.

The structure of uniformly RWDs has been studied by several authors (see, e.g., James 1964; Anand 1965; Roxburgh & Durney 1966; Monaghan 1966; Geroyannis & Hadjopoulos 1989). The issue of the stability of both uniformly and differentially rotating WDs has been studied as well (see, e.g., Ostriker & Bodenheimer 1968; Ostriker & Tassoul 1969; Tassoul & Ostriker 1970; Durisen 1975). All of the above computations were carried out within Newtonian gravity or at the post-Newtonian approximation. The EOS of cold WD matter has been assumed to be either a microscopically uniform degenerate electron fluid, which we refer to hereafter as Chandrasekhar EOS (Chandrasekhar 1931), or a polytropic EOS. However, microscopic screening caused by Coulomb interactions as well as the process of inverse β-decay of the composing nuclei cannot be properly studied within such an EOS (see Rotondo et al. 2011a, 2011b, for details).

The role of general relativistic effects, shown in Rotondo et al. (2011a), has been neglected in all of the above preceding literature. The only exception to this rule, to our knowledge, is the work of Arutyunyan et al. (1971), who investigated uniformly RWDs for the Chandrasekhar EOS within GR. They use an Ω² approximation following a method developed by Sedrakyan & Chubaryan (1968), independently of the work of Hartle (1967). A detailed comparison of our results with the ones of Arutyunyan et al. (1971) can be found in Appendix C.

In Figures 2 and 3, we show the mass-central density relation and the mass–radius relation of general relativistic rotating ¹²C and ¹⁶O WDs. We explicitly show the boundaries of mass-shedding, secular axisymmetric instability, inverse β-decay, and pycnonculear reactions.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Turning now to the rotation properties, in Figure 4, we show the J–M plane, specifically focusing on RWDs with masses larger than the maximum non-rotating mass, hereafter referred to as super-Chandrasekhar WDs (SCWDs). It becomes clear from this diagram that SCWDs can be stable only by virtue of their non-zero angular momentum: the lower half of the stability line of Figure 4, from J = 0 at M/M^{J = 0}_max all the way up to the value of J at M^{J ≠ 0}_max ∼ 1.06M^{J = 0}_max, determines the critical (minimum) angular momentum under which SCWDs become unstable. The upper half of the stability line determines, instead, the maximum angular momentum that SCWDs can have.

Zoom In Zoom Out Reset image size

The maximum masses of rotating WDs belong to the Keplerian sequence (see Figures 2–4) and can be expressed as

$$\begin{equation} M_{{\rm max}}^{J\ne 0}=k\,M_{{\rm max}}^{J=0}, \end{equation} \tag{ 5 }$$

where M^{J = 0}_max is the maximum stable mass of non-rotating WDs and k is a numerical factor that depends on the chemical composition, see Table 3 for details. For ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe RWDs, we found M^{J ≠ 0}_max ∼ 1.500, 1.474, 1.467, 1.202 M_☉, respectively.

Table 3. Properties of Uniformly Rotating General Relativistic ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe WDs.

Composition	$\rho _{M^{J\ne 0}_{{\rm max}}}$	k	MJ = 0max/M☉	$R_{M^{J=0}_{{\rm max}}}$	Pmin	RPminp	RPmineq	$(T/|W|)^{P_{{\rm min}}}$	$\epsilon ^{P_{{\rm min}}}$	$j^{P_{{\rm min}}}$	$q^{P_{{\rm min}}}$
4He	5.46 × 109	1.0646	1.40906	1163	0.284	564	736	0.0163	0.642	1.004	526
12C	6.95 × 109	1.0632	1.38603	1051	0.501	817	1071	0.0181	0.647	1.287	1330
16O	7.68 × 109	1.0626	1.38024	1076	0.687	1005	1323	0.0194	0.651	1.489	2263
56Fe	1.18 × 109	1.0864	1.10618	2181	2.195	2000	2686	0.0278	0.667	2.879	23702

Notes. $\rho _{M_{{\rm max}}^{J\ne 0}}$ is the central density in g cm−3 corresponding to the rotating maximum mass MJ ≠ 0max; k is the dimensionless factor used to express the rotating maximum mass MJ ≠ 0max as a function of the non-rotating maximum mass MJ = 0max of WDs, in solar masses, obtained in Rotondo et al. (2011a), as defined in Equation (5); the corresponding minimum radius is $R_{M^{J=0}_{{\rm max}}}$, in km; Pmin is the minimum rotation period in seconds. We recall that the configuration with Pmin is obtained for a WD rotating at the mass-shedding limit with a central density equal to the critical density for inverse β-decay (see Table 2 and the right panel of Figure 6). The polar RPminp and equatorial RPmineq radii of the configuration with Pmin are also given in km. The quantity $(T/|W|)^{P_{{\rm min}}}$ is the ratio between the kinetic and binding energies, the parameter $\epsilon ^{P_{{\rm min}}}$ is the eccentricity of the star, rotating at Pmin. Finally, $j^{P_{{\rm min}}}$ and $q^{P_{{\rm min}}}$ are the dimensionless angular momentum and quadrupole moment of WDs, respectively.

Download table as:  ASCIITypeset image

In Table 4, we compare the properties of the configuration with maximum mass using different EOSs, namely, Chandrasekhar μ = 2 (see, e.g., Boshkayev et al. 2011), Salpeter, and RFMT EOS. A comparison with classical results obtained with different treatments and EOSs can be found in Appendix C.

It is worth mentioning that the maximum mass of RWDs is not associated with a critical maximum density for gravitational collapse. This is in contrast with the non-rotating case, where the configuration of maximum mass (turning point) corresponds to a critical maximum density over which the WD is unstable against gravitational collapse.

The angular momentum J along the mass-shedding sequence is not constant and thus the turning-point criterion (Equation (3)) does not apply to this sequence. Therefore, the configuration of maximum rotating mass (Equation (5)) does not separate stable from secular axisymmetrically unstable WDs. We have also verified that none of the RWDs belonging to the mass-shedding sequence are a turning point of some J = constant sequence, and therefore they are indeed secularly stable. We therefore extend the Keplerian sequence all the way up to the critical density for inverse β decay, ρ^β_crit; see Table 2 and Figure 2.

The minimum rotation period P_min of WDs is obtained for a configuration rotating at the Keplerian angular velocity, at the critical inverse β-decay density; i.e., this is the configuration lying at the crossing point between the mass-shedding and inverse β-decay boundaries, see Figures 2 and 4. For ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe RWDs, we found the minimum rotation periods ∼0.28, 0.50, 0.69, and 2.19 s, respectively (see Table 3 for details). In Table 4, we compare the properties of the configuration with the minimum rotation period using different EOS, namely Chandrasekhar μ = 2, Salpeter, and RFMT EOS.

In the case of ¹²C WDs, the minimum period, 0.50 s, has to be compared with the value obtained by assuming as critical density the threshold for pycnonuclear reactions. Assuming lifetimes τ^{C + C}_pyc = 10 Gyr and 0.1 Myr, corresponding to critical densities ρ_pyc ∼ 9.26 × 10⁹ g cm⁻³ and ρ_pyc ∼ 1.59 × 10¹⁰ g cm⁻³, we obtain minimum periods P^pyc_min = 0.95 and 0.75 s, respectively.

It is interesting to compare and contrast some classical results with the ones presented in this work. Using a post-Newtonian approximation, Roxburgh & Durney (1966) analyzed the problem of the dynamical stability of maximally rotating RWDs, i.e., WDs rotating at the mass-shedding limit. The result was a minimum polar radius of 363 km, assuming the Chandrasekhar EOS with μ = 2. The Roxburgh critical radius is rather small with respect to our minimum polar radii, see Table 3. It is clear that such a small radius would lead to a configuration with the central density over the limit established by inverse β-decay: the average density obtained for the Roxburgh's critical configuration is ∼1.47 × 10¹⁰ g cm⁻³, assuming the maximum mass of 1.48 M_☉ obtained in the same work (see Table 6). A configuration with this mean density will certainly have a central density larger than the inverse β-decay density of ¹²C and ¹⁶O, 3.97 × 10¹⁰ g cm⁻³ and 1.94 × 10¹⁰ g cm⁻³, respectively (see Table 2). The rotation period of the WD at the point of dynamical instability of Roxburgh certainly must be shorter than the minimum values presented here.

The above comparison is in line with the fact that we did not find any turning point that crosses the mass-shedding sequence (see Figures 2 and 3). Presumably, ignoring the limits imposed by inverse β-decay and pycnonuclear reactions, the boundary determined by the turning points could cross the Keplerian sequence at some higher density. Such a configuration should have a central density very similar to the one found by Roxburgh & Durney (1966).

Arutyunyan et al. (1971) did not consider the problem of the minimum rotation period of a WD. However, they showed their results for a range of central densities covering the range of interest of our analysis. Thus, we have interpolated their numerical values of the rotation period of WDs in the Keplerian sequence and calculated the precise values at the inverse β-decay threshold for ⁴He, ¹²C, and ¹⁶O that have μ = 2 and therefore, in principle, are comparable to the Chandrasekhar EOS results with the same mean molecular weight. We thus obtained minimum periods ∼0.31, 0.55, 0.77 s, in agreement with our results (see Table 5).

It is important to stress that although it is possible to compare the results using the Chandrasekhar EOS μ = 2 with the ones obtained for the RFMT EOS, both qualitative and quantitative differences exist between the two treatments. In the former, a universal mass–density and mass–radius relation is obtained by assuming μ = 2 while, in reality, the configurations of equilibrium depend on the specific values of Z and A in a non-trivial way. For instance, ⁴He, ¹²C, and ¹⁶O have μ = 2, but the configurations of equilibrium are rather different. This fact was emphasized by Hamada & Salpeter (1961) in the Newtonian case, and further in GR by Rotondo et al. (2011a) for non-rotating configurations. In Figure 5, we present a comparison of the mass–density and mass–radius for the universal Chandrasekhar μ = 2 and the RFMT EOS for specific nuclear compositions.

Zoom In Zoom Out Reset image size

Regarding the stability of rotating WDs, Ostriker & Bodenheimer (1968), Ostriker & Tassoul (1969), and Durisen (1975) showed that uniformly rotating Newtonian polytropes and WDs described by the uniform degenerate electron fluid EOS are axisymmetrically stable at any rotation rate. In clear contrast with these results, we have shown here that uniformly RWDs can indeed be secularly axisymmetrically unstable as can be seen from Figures 2–4 (green boundary). We have constructed in Appendix C Newtonian RWDs for the Chandrasekhar EOS and compare the differences with the general relativistic counterpart. Apart from the quantitative differences for the determination of the mass at high densities, the absence of turning points in the Newtonian mass–density relation can be seen in Figure 7 (left panel). This can be understood from the fact that the maximum stable mass of non-rotating WDs is, in the Newtonian case, formally reached at an infinite central density. We should then expect that turning points will only appear from a post-Newtonian approximation, where the critical mass is shifted to finite densities (see, e.g., Roxburgh & Durney 1966, for the calculation of dynamical instability for post-Newtonian RWDs obeying the Chandrasekhar EOS).

In this respect, Figure 4 is of particular astrophysical relevance. Configurations lying in the filled region are stable against mass-shedding, inverse β-decay, and secular axisymmetric instabilities. RWDs with masses smaller than the maximum non-rotating mass (sub-Chandrasekhar WDs), i.e., M^{J ≠ 0} < M^{J = 0}_max, can have angular momenta ranging from a maximum at the mass-shedding limit all the way down to the non-rotating limit J = 0. SCWDs, however, are stabilized due to rotation, and therefore there exists a minimum angular momentum, J_min > 0, to guarantee their stability. We have shown above that secular axisymmetric instability is relevant for the determination of this minimum angular momentum of SCWDs (see green boundary in Figure 4). In this respect, it is interesting to note that from our results it turns out that SCWDs with light chemical compositions such as ⁴He and ¹²C are unstable against axisymmetric, inverse β-decay, and mass-shedding instabilities. On the contrary, in SCWDs with heavier chemical compositions, such as ¹⁶O and ⁵⁶Fe, the secular axisymmetric instability does not take place; see Figure 4. The existence of the new boundary due to secular axisymmetric instability is a critical issue for the evolution of SCWDs, since their lifetime might be reduced depending on their initial mass and angular momentum.

From the quantitative point of view, we have found that axisymmetric instability sets in for ¹²C SCWDs in the range of masses M^{J = 0}_max < M ≲ 1.397 M_☉, for some specific range of rotation periods ≳ 1.24 s. We can express the minimum rotation period that an SCWD with a mass M within the above mass range can have through the fitting formula

$$\begin{equation} P_{{\rm axi}}=0.062 \left(\frac{M-M^{J=0}_{{\rm max}}}{M_\odot }\right)^{-0.67}\,\,{\rm s}, \end{equation} \tag{ 6 }$$

where M^{J = 0}_max is the maximum mass of general relativistic non-rotating ¹²C WDs, M^{J = 0}_max ≈ 1.386 M_☉ (see Table 1 and Rotondo et al. 2011a). Thus, Equation (6) describes the rotation periods of the configurations along the green-dotted boundary in Figures 2–4. Correspondingly, the central density along this instability boundary varies from the critical density of static ¹²C WDs, ρ^{C, J = 0}_crit = 2.12 × 10¹⁰ g cm⁻³ (see Table 1), up to the inverse β-decay density, ρ^C_β = 3.97 × 10¹⁰ g cm⁻³ (see Table 2).

It is important to note that at the lower edge of the density range for axisymmetric instability, ρ^{C, J = 0}_crit, the timescale of C+C pycnonuclear reactions is τ^{C + C}_pyc ≈ 339 yr (see Figure 1). It then becomes of interest to compare this timescale with that corresponding to the secular axisymmetric instability that sets in at the same density.

The growing time of the secular instability is given by the dissipation time that can be driven either by gravitational radiation or viscosity (Chandrasekhar 1970). However, the gravitational radiation reaction is expected to drive secular instabilities for systems with rotational to gravitational energy ratio T/|W| ∼ 0.14, the bifurcation point between McClaurin spheroids and Jacobi ellipsoids (see Chandrasekhar 1970, for details). Therefore, we expect gravitational radiation to become important only for differentially rotating WDs, which can attain more mass and more angular momentum (Ostriker & Bodenheimer 1968). In the present case of general relativistic uniformly RWDs, only the viscosity timescale τ_v is relevant. A rotating star that becomes secularly unstable first evolves with a characteristic time τ_v and eventually reaches a point of dynamical instability, thus collapsing within a time of $\tau _{{\rm dyn}}\approx \Omega ^{-1}_K \sim \sqrt{R^3/G M} \lesssim 1$ s, where R is the radius of the star (see, e.g., Stergioulas 2003).

The viscosity timescale can be estimated as τ_v = R²ρ/η (see, e.g., Lindblom 1987), where ρ and η are the density and viscosity of the star. The viscosity of a WD assuming degenerate relativistic electrons is given by (Durisen 1973)

$$\begin{equation} \eta _{{\rm fluid}} =4.74\times 10^{-2} \frac{H_\Gamma (Z)}{Z} \rho ^{5/3}\left[ \left(\frac{\rho }{2\times 10^6}\right)^{2/3}+1 \right]^{-1}, \end{equation} \tag{ 7 }$$

where $H_\Gamma (Z)$ is a slowly varying dimensionless constant that depends on the atomic number Z and the Coulomb to thermal energy ratio

$$\begin{equation} \Gamma =\frac{e^2 Z^2}{k_B T}\left(\frac{4 \pi }{3}\frac{\rho }{2 Z Mu}\right)^{1/3}, \end{equation} \tag{ 8 }$$

where k_B is the Boltzmann constant and A ≃ 2Z has been used.

Expression (7) is valid for values of Γ smaller than the critical value for crystallization Γ_cry. The critical Γ_cry is not well constrained but its value should be of the order of Γ_cry ∼ 100 (see, e.g., Durisen 1973; Shapiro & Teukolsky 1983). The critical value Γ_cry defines a crystallization temperature T_cry under which the system behaves as a solid. For Γ_cry ∼ 100, we have T_cry ≈ 8 × 10⁷[ρ/(10¹⁰ g cm⁻³)]^1/3 K, for Z = 6. When Γ approaches Γ_cry, the viscosity can increase drastically to values close to (van Horn 1969; Durisen 1973)

$$\begin{equation} \eta _{{\rm cry}} =4.0\times 10^{-2} \left(\frac{Z}{7}\right)^{2/3} \rho ^{5/6}\exp [0.1(\Gamma -\Gamma _{{\rm cry}})]\,. \end{equation} \tag{ 9 }$$

For instance, we find that at densities ρ^{C, J = 0}_crit and assuming a central temperature of T ≳ 0.5T_cry with T_cry ≈ 10⁸ K, the viscous timescale is in the range 10 ≲ τ_v ≲ 1000 Myr, where the upper limit is obtained using Equation (7) and the lower limit with Equation (9). These timescales are longer than the pycnonuclear reaction timescale τ^{C + C}_pyc = 339 yr at the same density. So if the pycnonuclear reaction rates are accurate, it would imply that pycnonuclear reactions are more important to restricting the stability of RWDs with respect to the secular instability. However, we have to keep in mind that, as discussed in Section 3.4, the pycnonuclear critical densities are subject to theoretical and experimental uncertainties, which could in principle shift them to higher values. For instance, a possible shift of the density for pycnonuclear instability with timescales τ^{C + C}_pyc ∼ 1 Myr to higher values ρ^{C + C}_pyc > ρ^{C, J = 0}_crit, would suggest an interesting competition between secular and pycnonuclear instability in the density range ρ^{C, J = 0}_crit < ρ < ρ^C_β.

It is known that at constant rest mass M₀, entropy S, and chemical composition (Z, A), the spin evolution of an RWD is given by (see Shapiro et al. 1990, for details)

$$\begin{equation} \dot{\Omega } = \frac{\dot{E}}{\Omega } \left(\frac{\partial \Omega }{\partial J} \right)_{M_0,S,Z,A}, \end{equation} \tag{ 10 }$$

where $\dot{\Omega } \equiv d\Omega /dt$ and $\dot{E} \equiv dE/dt$, and E is the energy of the star.

Thus, if an RWD is losing energy by some mechanism during its evolution, that is $\dot{E}<0$, then the change of the angular velocity Ω in time depends on the sign of ∂Ω/∂J; RWDs that evolve along a track with ∂Ω/∂J > 0 will spin down ($\dot{\Omega }<0$), and the ones following tracks with ∂Ω/∂J < 0 will spin up ($\dot{\Omega }>0$).

In Figure 6, we show, in the left panel, the Ω = constant and J = constant sequences in the mass-central density diagram and, in the right panel, contours of constant rest mass in the Ω–J plane.

Zoom In Zoom Out Reset image size

The sign of ∂Ω/∂J can be analyzed from the left panel plot of Figure 6 by joining two consecutive J = constant sequences with a horizontal line and taking into account the fact that J decreases from left to right and from top to bottom. Instead, the angular velocity Ω decreases from right to left and from top to bottom for SCWDs and, for sub-Chandrasekhar WDs, from left to right and from top to bottom. We note that in the SCWDs region, Ω = constant sequences satisfy ∂Ω/∂ρ_c < 0 while, in the sub-Chandrasekhar region, both ∂Ω/∂ρ_c < 0 and ∂Ω/∂ρ_c > 0 appear (see minima). SCWDs can only either spin up by angular momentum loss or spin down by gaining angular momentum. In the latter case, the RWD becomes decompressed with time, increasing its radius and moment of inertia, and then SCWDs following this evolution track will end at the mass-shedding limit (see Figure 6). Some evolutionary tracks of sub-Chandrasekhar WDs and SCWDs are shown in the right panel of Figure 6. It is appropriate to recall here that Shapiro et al. (1990) showed that spin-up behavior by angular momentum loss occurs for rapidly rotating Newtonian polytropes if the polytropic index is very close to n = 3, namely, for an adiabatic index of Γ ≈ 4/3. It was explicitly shown by Geroyannis & Papasotiriou (2000) that these conditions are achieved only by super-Chandrasekhar polytropes.

Besides the confirmation of the above known result for SCWDs in the general relativistic case, here we also report the presence of minima ∂Ω/∂ρ_c = 0 for some sub-Chandrasekhar masses (see, e.g., the evolution track of the RWD with M = 1.38 M_☉ in the right panel of Figure 6), which raises the possibility that sub-Chandrasekhar WDs can experience, by angular momentum loss, not only the intuitively spin-down evolution, but also spin-up epochs.

It is appropriate to analyze the astrophysical consequences of the general relativistic RWDs presented in this work.

Most of the observed magnetic WDs are massive; for instance, REJ 0317-853 with M ∼ 1.35 M_☉ and B ∼ (1.7–6.6) × 10⁸ G (see, e.g., Barstow et al. 1995; Külebi et al. 2010); PG 1658+441 with M ∼ 1.31 M_☉ and B ∼ 2.3 × 10⁶ G (see, e.g., Liebert et al. 1983; Schmidt et al. 1992); and PG 1031+234 with the highest magnetic field ∼10⁹ G (see, e.g., Schmidt et al. 1986; Külebi et al. 2009). However, they are generally found to be slow rotators (see, e.g., Wickramasinghe & Ferrario 2000). It is worth mentioning that it has recently been shown by García-Berro et al. (2012) that such magnetic WDs can be indeed the result of the merger of double degenerate binaries; the misalignment of the final magnetic dipole moment of the newly born RWD with the rotation axis of the star depends on the difference of the masses of the WD components of the binary.

The precise computations of the evolution of the rotation period have to account for the actual value at each time of the moment of inertia, and the equatorial and polar radii of the WD. Whether magnetic and gravitational radiation braking can or can not explain the current relatively long rotation periods of some observed magnetic WDs is an important issue that deserves the appropriate attention and will be addressed elsewhere.

Magnetic braking of SCWDs has recently been invoked as a possible mechanism to explain the delayed time distribution of type Ia supernovae (SNe; see Ilkov & Soker 2012, for details): a type Ia SN explosion is delayed for a time that is typical of the spin-down timescale τ_B due to magnetic braking, provided that the result of the merging process of a WD binary system is a magnetic SCWD rather than a sub-Chandrasekhar one. The characteristic timescale τ_B of a SCWD has been estimated to be 10⁷ ≲ τ_B ≲ 10¹⁰ yr for magnetic fields comprised in the range 10⁶ ≲ B ≲ 10⁸ G. A constant moment of inertia ∼10⁴⁹ g cm² and a fixed critical(maximum) rotation angular velocity,

$$\begin{equation} \Omega _{\rm crit}\sim 0.7 \Omega ^{J=0}_{\rm K}=0.7 \left({\frac{G M^{J=0}}{R^3_{M^{J=0}}}}\right)^{1/2}, \end{equation} \tag{ 11 }$$

have been adopted (Ilkov & Soker 2012).

It is important to recall here that, as discussed in Section 8, SCWDs spin up by angular momentum loss, and therefore the reference to a "spin-down" timescale for them is just historical. SCWDs then evolve toward the mass-shedding limit, which in this case determines the critical angular velocity for rotational instability.

If we express Ω^{J ≠ 0}_K in terms of Ω^{J = 0}_K (see Appendix A.2), taking into account the values of j and q from the numerical integration, then we find for RWDs that the Keplerian angular velocity can be written as

$$\begin{equation} \Omega _{K}^{J\ne 0}=\sigma \Omega _{K}^{J=0}, \end{equation} \tag{ 12 }$$

where the coefficient σ varies in the interval [0.78,0.75] in the range of central densities [10⁵, 10¹¹] g cm⁻³. It is important to mention that the above range of σ remains approximately the same, independently of the chemical composition of the WD. However, the actual numerical value of the critical angular velocity, Ω^{J ≠ 0}_K, is different for different compositions owing to the dependence of the mass–radius relation of non-rotating WDs on (Z, A).

Furthermore, as we have shown, the evolution track followed by an SCWD depends strongly on the initial conditions of mass and angular momentum, as well as on chemical composition and the evolution of the moment of inertia (see Figure 6 and Section 8 for details). It is clear that the assumption of a fixed moment of inertia I ∼ 10⁴⁹ g cm² leads to a spin-down timescale that only depends on the magnetic field strength. A detailed computation will lead to a strong dependence on the mass of the SCWD, resulting in a two-parameter family of delayed times τ_B(M, B). Detailed calculations of the braking-down of the lifetime of SCWDs due to magnetic dipole radiation are then needed to shed light on this important matter. Theoretical work along these lines is currently in progress and the results will be presented in a forthcoming publication.

Massive, fast rotating and highly magnetized WDs have been proposed as an alternative scenario of soft gamma-ray repeaters (SGRs) and anomalous X-ray pulsars (AXPs), see Malheiro et al. (2012) for details. Within such a scenario, the range of minimum rotation periods of massive WDs found in this work, 0.3 ≲ P_min ≲ 2.2 s, depending on the nuclear composition (see Table 5), implies the rotational stability of SGRs and AXPs, which possess observed rotation periods 2 ≲ P ≲ 12 s. The relatively long minimum period of ⁵⁶Fe RWDs, ∼2.2 s, implies that RWDs describing SGRs and AXPs have to be composed of nuclear compositions lighter than ⁵⁶Fe, e.g., ¹²C or ¹⁶O.

We have calculated the properties of uniformly RWDs within the framework of GR using the Hartle formalism and our new EOS for cold WD matter based on the RFMT treatment (Rotondo et al. 2011b), which generalizes previous approaches including the EOS of Salpeter (1961). A detailed comparison with RWDs described by the Chandrasekhar and the Salpeter EOS has been performed.

We constructed the region of stability of RWDs taking into account the mass-shedding limit, secular axisymmetric instability, inverse β-decay, and pycnonuclear reaction lifetimes. The latter have been computed using the updated theoretical models of Gasques et al. (2005) and Yakovlev et al. (2006). We found that the minimum rotation periods for ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe RWDs are ∼0.3, 0.5, 0.7, and 2.2 s, respectively (see Table 5). For ¹²C WDs, the minimum period 0.5 s needs to be compared with the values P^pyc_min = 0.75 and 0.95 s, obtained by assuming the critical density to be the threshold for pycnonuclear reactions for lifetimes τ^{C + C}_pyc = 0.1 Myr and 10 Gyr, respectively. For the same chemical compositions, the maximum masses are ∼1.500, 1.474, 1.467, and 1.202 M_☉ (see Table 4). These results and additional properties of RWDs can be found in Table 3.

We have presented a new instability boundary of general relativistic SCWDs, over which they become axisymmetrically unstable. We have expressed the range of masses and rotation periods where this occurs through a fitting formula given by Equation (6). A comparison with Newtonian RWDs in Appendix C leads to the conclusion that this new boundary of instability for uniformly rotating WDs is a general relativistic effect.

We showed that, by losing angular momentum, sub-Chandrasekhar RWDs can experience both spin-up and spin-down epochs, while SCWDs can only spin up. These results are particularly important for the evolution of WDs whose masses approach, either from above or from below, the maximum non-rotating mass. The knowledge of the actual values of the mass, radii, and moment of inertia of massive RWDs are relevant for the computation of delay collapse times in the models of type Ia SN explosions. A careful analysis of all of the possible instability boundaries such as the one presented here have to be taken into account during the evolution of the WD at pre-SN stages.

We have indicated specific astrophysical systems where the results of this work are relevant: for instance, the long rotation periods of observed massive magnetic WDs, the delayed collapse of SCWDs as progenitors of type Ia SNe, and the alternative scenario for SGRs and AXPs based on massive RWDs.

We thank the anonymous referees for the many comments and suggestions that improved the presentation of our results. J.A.R. is grateful to Enrique García-Berro and Noam Soker for helpful discussions and remarks on the properties of magnetic WDs resulting from WD mergers and on the relevance of this work for the delayed collapse of super-Chandrasekhar WDs.

A.1. The Hartle–Thorne Vacuum Solution

The HT metric given by Equation (1) can be written in an analytic closed form in the exterior vacuum case in terms of the total mass M, angular momentum J, and quadrupole moment Q of the rotating star. The angular velocity of local inertial frames ω(r), proportional to Ω, and the functions h₀, h₂, m₀, m₂, k₂, proportional to Ω², are derived from the Einstein equations (see Hartle 1967; Hartle & Thorne 1968, for details). Thus, the metric can be then written as

$$\begin{eqnarray} ds^2&=&\left(1-\frac{2{ M }}{r}\right)\left[1+2k_1P_2(\cos \theta)+2\left(1-\frac{2{ M}}{r}\right)^{-1}\right. \nonumber\\ &&\left.\times\frac{J^{2}}{r^{4}}(2\cos ^2\theta -1)\right]dt^2+\frac{4J}{r}\sin ^2\theta dt d\phi \nonumber\\ &&-\left(1-\frac{2{ M}}{r}\right)^{-1}\left[1-2\left(k_1-\frac{6 J^{2}}{r^4}\right)P_2(\cos \theta)\right.\nonumber\\ &&-\left.2\left(1-\frac{2{ M}}{r}\right)^{-1}\frac{J^{2}}{r^4}\right]dr^2-r^2\nonumber\\ &&\times [1-2k_2P_2(\cos \theta)](d\theta ^2+\sin ^2\theta d\phi ^2), \end{eqnarray} \tag{ A1 }$$

where

$$\begin{eqnarray*} k_1&=&\frac{J^{2}}{{ M}r^3}\left(1+\frac{{ M}}{r}\right)+\frac{5}{8}\frac{Q-J^{2}/{ M}}{{ M}^3}Q_2^2\left(\frac{r}{{ M}}-1\right),\nonumber\\ k_2&=&k_1+\frac{J^{2}}{r^4}+\frac{5}{4}\frac{Q-J^{2}/{ M}}{{ M}^2r}\left(1-\frac{2{ M}}{r}\right)^{-1/2}Q_2^1\left(\frac{r}{ M}-1\right)\, \end{eqnarray*}$$

and

$$\begin{eqnarray} Q_{2}^{1}(x)&=&(x^{2}-1)^{1/2}\left[\frac{3x}{2}\ln \frac{x+1}{x-1}-\frac{3x^{2}-2}{x^{2}-1}\right],\nonumber\\ Q_{2}^{2}(x)&=&(x^{2}-1)\left[\frac{3}{2}\ln \frac{x+1}{x-1}-\frac{3x^{3}-5x}{(x^{2}-1)^2}\right] \end{eqnarray} \tag{ A2 }$$

are the associated Legendre functions of the second kind, where x = r/M − 1, and P₂(cos θ) = (1/2)(3cos ²θ − 1) is the Legendre polynomial. The constants M, J, and Q are the total mass, angular momentum, and mass quadrupole moment of the rotating object, respectively. This form of the metric corrects some misprints of the original paper by Hartle & Thorne (1968) (see also Berti et al. 2005 and Bini et al. 2009). The precise numerical values of M, J, and Q are calculated from the matching procedure of the exterior and interior metrics at the surface of the star.

The total mass of a rotating configuration is defined as M = M^{J ≠ 0} = M^{J = 0} + δM, where M^{J = 0} is the mass of the non-rotating configuration and δM is the change in mass of the rotating configuration from the non-rotating configuration with the same central density. It should be stressed that in the terms involving J² and Q, the total mass M can be substituted by M^{J = 0} since δM is already a second-order term in the angular velocity.

A.2. Angular Velocity of Equatorial Circular Orbits

The four-velocity u of a test particle on a circular orbit in the equatorial plane of axisymmetric stationary spacetime can be parameterized by the constant angular velocity Ω with respect to an observer at infinity,

$$\begin{equation} u=\Gamma [\partial _t+\Omega \partial _{\phi }], \end{equation} \tag{ A3 }$$

where Γ is a normalization factor which assures that u^αu_α = 1. From normalization and geodesics conditions, we obtain the following expressions for Γ and Ω = u^ϕ/u^t

$$\begin{eqnarray} \Gamma =&\pm (g_{tt}+2\Omega g_{t\phi }+\Omega ^2 g_{\phi \phi })^{-1/2},\nonumber\\ g_{tt,r}&+2\Omega g_{t\phi,r}+\Omega ^2 g_{\phi \phi,r}=0, \end{eqnarray} \tag{ A4 }$$

hence, Ω, the solution of (A4)₂, is given by

$$\begin{equation} \Omega _{\pm {\rm orb}}(r)=\frac{u^{\phi }}{u^{t}}=\frac{-g_{t\phi,r}\pm {((g_{t\phi,r})^2-g_{tt,r}g_{\phi \phi,r}})^{1/2}}{g_{\phi \phi,r}}, \end{equation} \tag{ A5 }$$

where (+ / −) stands for corotating/counterrotating orbits, u^ϕ and u^t are the angular and time components of the four-velocity, and a colon stands for the partial derivative with respect to the corresponding coordinate. In our case, one needs to consider only corotating orbits (omitting the plus sign in Ω_+orb(r) = Ω_orb(r)) to determine the mass shedding (Keplerian) angular velocity on the surface of the WD. For the HT external solution Equation (A1), we have

$$\begin{equation} \Omega _{{\rm orb}}(r)=\left({\frac{M}{r^3}}\right)^{1/2}[1- j F_{1}(r)+j^2F_{2}(r)+q F_{3}(r)], \end{equation} \tag{ A6 }$$

where j = J/M² and q = Q/M³ are the dimensionless angular momentum and quadrupole moment,

$${\begin{eqnarray*} F_{1}&=&\left(\frac{M}{r}\right)^{3/2},\nonumber\\ F_{2}&=&\frac{48{M}^7-80{M}^6r+4{M}^5r^2-18{M}^4r^3+40{M}^3r^4+10{M}^2r^5+15{M}r^6-15r^7}{16{M}^2r^4(r-2{M})}+F, \nonumber\\ F_{3}&=&\frac{6{M}^4-8{M}^3r-2{M}^2r^2-3{M}r^3+3r^4}{16{M}^2r(r-2{M})/5}-F,\nonumber\\ F&=&\frac{15(r^3-2{M^3})}{32{M}^3}\ln \frac{r}{r-2{M}}. \end{eqnarray*}}$$

The mass shedding limiting angular velocity of a rotating star is the Keplerian angular velocity evaluated at the equator (r = R_eq), i.e.,

$$\begin{equation} \Omega _K^{J\ne 0}=\Omega _{{\rm orb}}(r=R_{{\rm eq}}). \end{equation} \tag{ A7 }$$

In the static case, i.e., when j = 0 and hence q = 0 and δM = 0, the well-known Schwarzschild solution and the orbital angular velocity for a test particle Ω^{J = 0}_ms on the surface (r = R) of the WD is given by

$$\begin{equation} \Omega _K^{J=0}=\left({\frac{M^{J=0}}{R_{M^{J=0}}^3}}\right)^{1/2}. \end{equation} \tag{ A8 }$$

A.3. Weak Field Limit

Let us estimate the values of j and q, recovering physical units with c and G. The dimensionless angular momentum is

$$\begin{equation} j=\frac{cJ}{GM^2}=\frac{c}{G}\frac{\alpha MR^2 \Omega }{M^2}=\alpha \left(\frac{\Omega R}{c}\right)\left(\frac{GM}{c^2R}\right)^{-1}, \end{equation} \tag{ A9 }$$

where we have used the fact that J = IΩ, with I = αMR², and α ∼ 0.1 from our numerical integrations. For massive and fast rotating WDs, we have (ΩR)/c ∼ 10⁻² and (GM)/(c²R) ∼ 10⁻³, so j ∼ 1.

The dimensionless quadrupole moment q is

$$\begin{equation} q=\frac{c^4}{G^2}\frac{Q}{M^3}=\frac{c^4}{G^2}\frac{\beta MR^2}{M^3}=\beta \left(\frac{GM}{c^2R}\right)^{-2}, \end{equation} \tag{ A10 }$$

where we have expressed the mass quadrupole moment Q in terms of the mass and radius of the WD, Q = βMR², where β ∼ 10⁻², so we have q ∼ 10⁴.

The large values of j and q might arouse some suspicion concerning the products jF₁, j²F₂, and qF₃ as real correction factors in Equation (A6). It is easy to check this in the weak field limit M/r ≪ 1, where the functions F_i can be expanded as a power series

$$\begin{eqnarray*} F_{1}&=&\left(\frac{M}{r}\right)^{3/2},\nonumber\\ F_{2}&\approx& \frac{1}{2}\left(\frac{M}{r}\right)^3-\frac{117}{28}\left(\frac{M}{r}\right)^4-6\left(\frac{M}{r}\right)^5-\cdots,\nonumber\\ F_{3}&\approx& \frac{3}{4}\left(\frac{M}{r}\right)^2+\frac{5}{4}\left(\frac{M}{r}\right)^3+\frac{75}{28}\left(\frac{M}{r}\right)^4+6\left(\frac{M}{r}\right)^5+\cdots \end{eqnarray*}$$

so evaluating at r = R,

$$\begin{equation} j F_1=\alpha \left(\frac{\Omega R}{c}\right)\left(\frac{GM}{c^2R}\right)^{1/2}, \quad j^2F_2=\frac{\alpha }{2}\left(\frac{\Omega R}{c}\right)\left(\frac{GM}{c^2R}\right)^{2}, \end{equation} \tag{ A11 }$$

so we finally have jF₁ ∼ 10^−9/2, j²F₂ ∼ 10⁻⁹, and qF₃ ∼ 10⁻². We can therefore see that the products are indeed correction factors and, in addition, that the effect due to the quadrupolar deformation is larger than the frame-dragging effect.

The theoretical framework for the determination of the pycnonuclear reaction rates was developed by Salpeter & van Horn (1969). The number of reactions per unit volume per unit time can be written as

$$\begin{eqnarray} R_{{\rm pyc}}&=&Z^4 A\rho S(E_p) 3.90\times 10^{46} \lambda ^{7/4}\nonumber\\ &&\times \exp (-2.638/\sqrt{\lambda })\,\,{\rm cm}^{-3}\,{\rm s}^{-1},\nonumber\\ \lambda &=&\frac{1}{Z^2 A^{4/3}}\left(\frac{\rho }{1.3574\times 10^{11}\,{\rm g}\,{\rm cm}^{-3}}\right)^{1/3}, \end{eqnarray} \tag{ B1 }$$

where S are astrophysical factors in units of MeV barns (1 barn = 10⁻²⁴ cm²) that have to be evaluated at the energy E_p given by Equation (4).

For the S-factors, we adopt the results of Gasques et al. (2005) calculated with the NL2 nuclear model parameterization. For center of mass energies E ⩾ 19.8 MeV, the S-factors can be fitted by

$$\begin{equation} S(E)=5.15\times 10^{16}\exp\! \left[\!-0.428 E\,{-}\,\frac{3 E^0.308}{1+e^{0.613 (8-E)}}\!\right]{\rm MeV\,barn}, \end{equation} \tag{ B2 }$$

which is appropriate for the ranges of the zero-point energies at high densities. For instance, ¹²C nuclei at ρ = 10¹⁰ g cm⁻³ have a zero-point oscillation energy of E_p ∼ 34 keV.

All the nuclei (Z, A) at a given density ρ will fuse in a time τ_pyc given by

$$\begin{equation} \tau _{{\rm pyc}}=\frac{n_N}{R_{{\rm pyc}}}=\frac{\rho }{A M_u R_{{\rm pyc}}}, \end{equation} \tag{ B3 }$$

where n_N = ρ/(AM_u) is the ion density. Gasques et al. (2005) estimated that the S-factors (B2) are uncertain within a factor ∼3.5; it is clear from the above equation that for a given lifetime τ_pyc, such uncertainties are also reflected in the determination of the density threshold.

We have constructed solutions of the Newtonian equilibrium equations for RWDs that are accurate up to the order Ω², following the procedure of Hartle (1967). In Figure 7 (left panel), we compare these Newtonian configurations with general relativistic RWDs for the Chandrasekhar EOS with μ = 2. We can clearly see the differences between the two mass–density relations toward the high-density region, as expected. A most remarkable difference is the existence of an axisymmetric instability boundary in the general relativistic case, which is absent in its Newtonian counterpart.

Zoom In Zoom Out Reset image size

To our knowledge, the only previous work on RWDs within GR is that of Arutyunyan et al. (1971). A method to compute RWDs configurations accurately up to second order in Ω was developed by two authors (see Sedrakyan & Chubaryan 1968, for details), independently of the work of Hartle (1967). In Arutyunyan et al. (1971), RWDs were computed for the Chandrasekhar EOS with μ = 2.

In Figure 7 (right panel), we show the mass-central density relation obtained with their method along with the ones constructed in this work for the same EOS. We note here that the results are different even at the level of static configurations, and since the methods are based on the construction of rotating configurations from seed static ones, those differences extrapolate to the corresponding rotating objects. This fact is to be added to the possible additional difference arising from the different manner used to approach the order Ω² in the approximation scheme. The differences between the two equilibrium configurations are evident.

Turning now to the problem of the maximum mass of an RWD, in Table 6, we present the previous results obtained in Newtonian, post-Newtonian approach, and GR by several authors. Depending on their method, approach, treatment, theory, and numerical code, the authors showed different results. These maximum masses of RWDs are to be compared with the ones found in this work and presented in Table 4 for the Chandrasekhar μ = 2, Salpeter, and RFMT EOS.

Table 4. The Maximum Rotating Mass of General Relativistic Uniformly Rotating ⁴He, ¹²C, ¹⁶O and ⁵⁶Fe WDs for Different EOS

Nuclear Composition	EOS	$\rho _{M^{J\ne 0}_{{\rm max}}}$	$R^{M^{J\ne 0}_{{\rm max}}}_{p}$	$R^{M^{J\ne 0}_{{\rm max}}}_{{\rm eq}}$	MJ ≠ 0max/M☉	$P^{M^{J\ne 0}_{{\rm max}}}$)
		(g cm−3)	(km)	(km)		(s)
μ = 2	Chandrasekhar	1.07 × 1010	1198.91	1583.47	1.5159	0.884
	Salpeter	1.07 × 1010	1193.08	1575.94	1.4996	0.883
4He	RFMT	5.46 × 109	1458.58	1932.59	1.5001	1.199
	Salpeter	1.08 × 1010	1183.99	1564.16	1.4833	0.878
12C	RFMT	6.95 × 109	1349.15	1785.98	1.4736	1.074
	Salpeter	1.09 × 1010	1178.88	1556.68	1.4773	0.875
16O	RFMT	7.68 × 109	1308.09	1730.65	1.4667	1.027
	Salpeter	1.14 × 109	2002.43	2693.17	1.2050	2.202
56Fe	RFMT	1.18 × 109	2000.11	2686.06	1.2017	2.195

Notes. $\rho _{M^{J\ne 0}_{{\rm max}}}$, $R^{M^{J\ne 0}_{{\rm max}}}_{p}$, $R^{M^{J\ne 0}_{{\rm max}}}_{{\rm eq}}$, and $P^{M^{J\ne 0}_{{\rm max}}}$ are central density, polar and equatorial radii, and rotation period of the configuration with the maximum mass, M^{J ≠ 0}_max.

Download table as:  ASCIITypeset image

In his classic work, Hartle (1967) described the slow rotation regime by requesting that fractional changes in pressure, energy density, and gravitational field due to the rotation of the star are all much smaller with respect to a non-rotating star with the same central density. From a dimensional analysis, such a condition implies

$$\begin{equation} \Omega ^2\ll \left(\frac{c}{R}\right)^2 \frac{GM^{J=0}}{c^2R}, \end{equation} \tag{ D1 }$$

where M^{J = 0} is the mass of the unperturbed configuration and R its radius. The expression on the right is the only multiplicative combination of M, R, G, and c, and in the Newtonian limit coincides with the critical Keplerian angular velocity Ω^{J = 0}_K given by Equation (A8). For unperturbed configurations with (GM)/(c²R) < 1, the condition (D1) implies ΩR/c ≪ 1. Namely, every particle must move at non-relativistic velocities if the perturbation to the original geometry have to be small in terms of percentage. Equation (D1) can also be written as

$$\begin{equation} \Omega \ll \Omega ^{J=0}_{K}, \end{equation} \tag{ D2 }$$

which is the reason why it is often believed that the slow rotation approximation is not suitable for describing stars rotating at their mass-shedding value.

Let us discuss this point more carefully. It is clear that the request that the contribution of rotation to pressure, energy density, and gravitational field be small can be summarized in a single expression, Equation (D1), since all of them are quantitatively given by the ratio between the rotational and the gravitational energy of the star. The rotational energy is T ∼ MR²Ω² and the gravitational energy is |W| ∼ GM²/R = (GM/c²R)Mc², hence the condition T/|W| ≪ 1 leads to Equation (D1) or Equation (D2). Now, we will discuss the above condition for realistic values of the rotational and gravitational energy of a rotating star, abandoning the assumption of either fiducial or order-of-magnitude calculations. We show below that the actual limiting angular velocity on the right-hand side of the condition (D2) has to be higher than the Keplerian value.

We can write the gravitational binding energy of the star as |W| = γGM²/R and the rotational kinetic energy as T = (1/2)IΩ² = (1/2)αMR²Ω², where the constants γ and α are structure constants that depends on the density and pressure distribution inside the star. According to the slow rotation approximation, T/|W| ≪ 1, namely,

$$\begin{eqnarray} \frac{T}{|W|}&=&\frac{\alpha M R^2\Omega ^2/2}{\gamma G M^2/R}=\left(\frac{\alpha }{2\gamma }\right)\left(\frac{GM}{R^3}\right)^{-1}\nonumber\\ \Omega ^2&=&\left(\frac{\alpha }{2\gamma }\right)\left(\frac{\Omega }{\Omega _{K}^{J=0}}\right)^2 \ll 1, \end{eqnarray} \tag{ D3 }$$

which can be rewritten analogously to Equation (D2) as

$$\begin{equation} \Omega \ll \left({\frac{2\gamma }{\alpha }}\right)^{1/2} \Omega _{K}^{J=0}. \end{equation} \tag{ D4 }$$

Now, we check that the ratio of the structural constants is larger than unity. Let us first consider the simplest example of a constant density sphere. In this case, α = 2/5 and γ = 3/5, so (2γ/α)^1/2 ≈ 1.73, and the condition (D4) is Ω ≪ 1.73Ω^{J = 0}_K. If we now consider a more realistic density profile, for instance, a polytrope of index n = 3, we have (see, e.g., Shapiro & Teukolsky 1983)

$$\begin{equation} |W|=\frac{3}{5-n}\frac{GM^2}{R}=\frac{3}{2}\frac{GM^2}{R},\quad T=\frac{1}{2}I\Omega ^2=\frac{1}{2}\frac{2}{3}M\langle r^2\rangle \Omega ^2, \end{equation} \tag{ D5 }$$

where 〈r²〉 = 0.11303R². Therefore, in this case we have γ = 3/2 and α = 0.075, and so Equation (D4) becomes Ω ≪ 6.32Ω^{J = 0}_K. This is not surprising since T/|W| → 0.025 when Ω → Ω^{J = 0}_K.

The above analysis has been performed while assuming spherical symmetry. When deviations from the spherical shape are taken into account, the ratio T/|W| turns out to be even smaller than the previous estimates based on spherical polytropes. Since the equatorial radius satisfies R_eq > R, at mass-shedding we will have Ω < Ω^{J = 0}_K. In fact, in the Roche model, the mass-shedding angular velocity is Ω^{J ≠ 0}_K = (2/3)^3/2Ω^{J = 0}_K ≈ 0.544Ω^{J = 0}_K, corresponding to a rotational to gravitational energy ratio of T/|W| ≈ 0.0074 (see, e.g., Shapiro & Teukolsky 1983).

In our RWDs, we have determined that the mass-shedding angular velocity satisfies Ω^{J ≠ 0}_K ≈ 0.75Ω_K^{J = 0} at any density, see Equation (12). Accordingly, we show in the left panel of Figure 8 the ratio T/|W| for RWDs as a function of the central density for the Keplerian sequence. For an increasing central density, T/|W| decreases. On the right panel, we have plotted the eccentricity versus the central density. For increasing central density, the eccentricity decreases, so RWDs become less oblate at higher densities.

Zoom In Zoom Out Reset image size

Now, we turn to evaluate more specifically the deviations from the spherical symmetry. The expansion of the radial coordinate of a rotating configuration r(R, θ) in powers of the angular velocity is written as (Hartle 1967)

$$\begin{equation} r=R+\xi (R,\theta)+O(\Omega ^4), \end{equation} \tag{ D6 }$$

where ξ is the difference in the radial coordinate, r, between a point located at the polar angle θ on the surface of constant density ρ(R) in the rotating configuration, and the point located at the same polar angle on the same constant density surface in the non-rotating configuration. In the slow rotation regime, the fractional displacement of the surfaces of constant density due to the rotation have to be small, namely, ξ(R, θ)/R ≪ 1, where ξ(R, θ) = ξ₀(R) + ξ₂(R)P₂(cos θ) and ξ₀(R) and ξ₂(R) are functions of R proportional to Ω². In the right panel of Figure 9, the difference in the radial coordinate over static radius versus the central density is shown. Here, we see the same tendency as in the case of the eccentricity that these differences are decreasing with an increasing central density. In the left panel, the rotation parameter ΩR/c versus the central density is shown. Here, with an increasing central density, the rotation parameter increases. Thus, for higher densities, the system becomes less oblate, smaller in size with a larger rotation parameter, i.e., a higher angular velocity.

Zoom In Zoom Out Reset image size

In order to estimate the accuracy of the slow rotation approximation for RWDs, based on the above results, it is useful to compare all the above numbers with the known results for NSs. For instance, we notice that in NSs ΩR/c ∼ 10⁻¹, ξ(R, 0)/R ∼ 10⁻², and ξ(R, π/2)/R ∼ 10⁻¹ (see, e.g., Berti et al. 2005), to be compared with the corresponding values of RWDs shown in Figure 9, ΩR/c ≲ 10⁻², ξ(R, 0)/R ∼ 10⁻², and ξ(R, π/2)/R ∼ 10⁻¹. Weber & Glendenning (1992) calculate the accuracy of the Hartle's second-order approximation and found that the mass of maximally rotating NSs is accurate within an error ≲ 4%; Benhar et al. (2005) found that the inclusion of third-order expansion Ω³ improved the mass-shedding limit numerical values by less than 1% for NSs obeying different EOS. On the other hand, it is known that the ratio T/|W| in the case of NSs is as large as ∼0.1 in the Keplerian sequence (see, e.g., Tables 1–5 of Berti & Stergioulas 2004). Since RWDs have T/|W| and ΩR/c smaller than NSs, and δR/R = ξ/R at least of the same order (see left panel of Figure 8), we expect that the description of the structure of RWDs up to the mass-shedding limit within the Hartle's approach to have at least the same accuracy as in the case of NSs.

Table 5. The Minimum Rotation Period of General Relativistic Rotating ⁴He, ¹²C, ¹⁶O, and ⁵⁶Fe WDs

Nuclear Composition	EOS	ρβcrit	$R^{P_{{\rm min}}}_{p}$	$R^{P_{{\rm min}}}_{{\rm eq}}$	$M^{J\ne 0}_{P_{{\rm min}}}/M_\odot$	Pmin
		(g cm−3)	(km)	(km)		(s)
μ = 2	Chandra	1.37 × 1011	562.79	734.54	1.4963	0.281
	Salpeter	1.37 × 1011	560.41	731.51	1.4803	0.281
4He	RFMT	1.39 × 1011	563.71	735.55	1.4623	0.285
	Salpeter	3.88 × 1010	815.98	1070.87	1.4775	0.498
12C	RFMT	3.97 × 1010	816.55	1071.10	1.4618	0.501
	Salpeter	1.89 × 1010	1005.62	1324.43	1.4761	0.686
16O	RFMT	1.94 × 1010	1005.03	1323.04	1.4630	0.687
	Salpeter	1.14 × 109	2002.43	2693.17	1.2050	2.202
56Fe	RFMT	1.18 × 109	2000.11	2686.06	1.2018	2.195

Notes. ρ^β_crit is the critical density for inverse β decay. $M^{J\ne 0}_{P_{{\rm min}}}$, $R^{P_{{\rm min}}}_{p}$, and $R^{P_{{\rm min}}}_{e}$ are the mass, polar, and equatorial radii corresponding to the configuration with minimum rotation period, P_min.

Download table as:  ASCIITypeset image