The shape of a rapidly rotating polytrope with index unity

TL;DR

This paper evaluates the accuracy of approximate solutions for rapidly rotating, self-gravitating polytropes with index unity, highlighting their limitations in representing true equilibrium configurations.

Contribution

It critically assesses the solutions from a prior study, clarifying their approximate nature and discussing their validity as models of rotating polytropic stars.

Findings

The solutions are only approximate, not exact, solutions of the Euler-Poisson system.

The paper discusses the degree of approximation and its implications for modeling rotating polytropes.

It provides insights into the limitations of existing analytical solutions for such systems.

Abstract

We show that the solutions obtained in the paper `An exact solution for arbitrarily rotating gaseous polytropes with index unity' by Kong, Zhang, and Schubert represent only approximate solutions of the free-boundary Euler-Poisson system of equations describing uniformly rotating, self-gravitating polytropes with index unity. We discuss the quality of such solutions as approximations to the rigidly rotating equilibrium polytropic configurations.

Figures

Equations

$$\frac{\partial P}{\partial \varpi }-\rho {\Omega }^2\varpi$$

$$\frac{\partial P}{\partial z}$$

$$\Delta \Phi =4\pi G\rho .$$

$$\nabla (h+{\Phi }_c+\Phi )=0,$$

$$h+{\Phi }_c+\Phi =C,$$

$${\Phi }_c=-\frac{1}{2}{\Omega }^2{\varpi }^2,$$

$$\nabla (2K\rho +{\Phi }_c+\Phi )=0,$$

$$\Delta \rho +\frac{2\pi G}{K}\rho =\frac{{\Omega }^2}{K}.$$

Full text

1 Introduction

In a recent paper Kong et al. (2015) attacked a long-standing problem of obtaining analytic solutions representing rigidly rotating stationary polytropes with the polytropic index $n=1$ (polytropic exponent $\gamma=1+1/n=2$). The quest for such solutions started with the works of Chandrasekhar (1933a)111This is actually the first of a series of four Chandrasekhar’s papers on distorted polytropes (Chandrasekhar, 1933a, b, c, d). and Kopal (1937, 1939), but the majority of papers on this subject has been published during three decades between 1960 and 1990 (Roberts, 1963; Papoyan et al., 1967; Blinnikov, 1972; Hubbard, 1974; Kozenko, 1975; Cunningham, 1977; Zharkov & Trubitsyn, 1978; Caimmi, 1980; Williams, 1988). A modern account of these works can be found in Horedt (2004).

Kong et al. (2015) assume that the shape of the boundary surface of such polytropes is an ellipsoid of revolution and construct analytic solutions of the corresponding equation for the distribution of the density. In this note we show that they satisfy only approximately the coupled system of Euler and Poisson equations describing equilibrium configurations of rigidly rotating polytropes with $n=1$. We also clarify the difference between solving the free-boundary Euler–Poisson problem and the corresponding equation for the density that is solved in (Kong et al., 2015).

Solutions obtained by Kong et al. (2015) can still be understood as approximations to the rigidly rotating polytropes with $n=1$, although the true shape of the boundary surface of such polytropes deviates from the ellipsoidal one. In general, the shape of an axisymmetric rotating polytropic body ($0<n<5$) can be described only approximately by a biaxial ellipsoid. The exact ellipsoidal description is possible only for the homogeneous, constant density polytrope $n=0$ (cf. Tassoul, 1978, pp. 82–84).

In the following sections we discuss numerical solutions of rigidly rotating $n=1$ polytropes, testing the ellipsoidal approximation. As expected, the deviation from the ellipsoidal shape is significant for rapidly rotating configurations. We also pay special attention to the model that was computed in Kong et al. (2015) assuming observational parameters of $\alpha$ Eri.

2 Rotating polytropes with index unity

Rigidly rotating axisymmetric equilibrium configurations of a self-gravitating perfect fluid are described by the Euler equations

[TABLE]

and the Poisson equation for the gravitational potential $\Phi$

[TABLE]

Here $P$ denotes the pressure, $\rho$ is the mass-density, $\Omega$ is the angular velocity of the fluid, and $G$ is the gravitational constant. We work in cylindrical coordinates $(\varpi,\phi,z)$, assuming that the $z$ axis coincides with the rotation axis of the fluid.

For barotropic equations of state $P=P(\rho)$, one can introduce the specific enthalpy $h$ so that $dh=dP/\rho$. Then, in the region $U$ where $\rho>0$, Eqs. (1) can be written as

[TABLE]

or, equivalently,

[TABLE]

where the centrifugal potential is

[TABLE]

and $C$ denotes an integration constant.

In this paper we deal exclusively with polytropic equations of state of the form $P=K\rho^{1+1/n}$, where $K$ is a constant, and assume $n=1$. This yields $P=K\rho^{2}$, and a linear relation $h=2K\rho$.

A rigidly rotating $n=1$ polytrope is thus described by the Poisson equation (2) and the equation

[TABLE]

that holds in the region $U$, and where $\Phi_{c}$ is given by Eq. (5). It is important to stress that the boundary of the region $U$, at which the density $\rho$ tends to zero, is not known a priori, but has to be established by the procedure of solving Eqs. (2) and (6) itself. In mathematical literature, such a formulation is known as a free-boundary problem. It is a remarkable property of the set of Eqs. (2) and (6) that it allows for the determination of the boundary of $U$. In other words, Eqs. (2) and (6) do not admit regular solutions with a non-negative density $\rho$ for a wrong guess of the shape $U$.

The approximation presented by Kong et al. (2015) (and advertised as an exact solution) is based on the following, well known, observation. By computing the divergence of Eq. (6) and combining the result with Eq. (2), one obtains an equation for the density

[TABLE]

The above equation also holds in the region $U$, and it is assumed that $\rho$ tends to zero at the boundary of $U$. Equation (7) is equivalent to Eq. (8) of Kong et al. (2015), which is written in oblate spheroidal coordinates, and in fact it is Eq. (7) and not the system of Eqs. (2) and (6) that is solved by Kong et al. (2015). In contrast to the system of Eqs. (2) and (6), it admits solutions with the boundary of $U$ being an ellipsoid of revolution.

It is clear from the above derivation that solutions of Eq. (7) and the distribution of the density $\rho$ obtained by solving Eqs. (2) and (6) may differ by a function $f$ satisfying $\Delta f=0$. It is also easy to see that Eq. (7) yields the proper solution, i.e. a solution equal to that of Eqs. (2) and (6), provided that the true shape of the boundary of $U$ is assumed. The fact that this boundary cannot be an ellipsoid of revolution follows from the classic theorem that can be found in the textbook by Tassoul (1978, pp. 82–84). It states that there are no rigidly rotating, barotropic solutions of Eqs. (1) and (2) with an ellipsoidal boundary, except for the configurations of constant density.222The theorem is actually more general: it excludes all barotropic configurations with an ellipsoidal stratification.

3 Numerical solutions

Because guessing of the right shape of the region $U$ in Eq. (6) is difficult (if at all possible), one has to resort to numerical methods or a combination of numerical and analytic schemes (see, e.g., Odrzywołek, 2003).

In this section we give numerical examples of solutions of Eqs. (2) and (6), obtained without assuming that the surface of the configuration is ellipsoidal. Our algorithm follows strictly the scheme proposed by Eriguchi & Müller (1985). In short, it is based on Eq. (4), which is discretised on a two-dimensional grid. The values of the gravitational potential $\Phi$ at the grid nodes are computed using the expansion of the Green function for the Laplace operator in terms of Legendre polynomials. In this way Eqs. (2) and (6) are converted into a system of algebraic equations, which is then solved using standard numerical methods.

The method of Eriguchi & Müller (1985) allows one to obtain whole sequences of models, say with a fixed value of the total mass and different values of the angular momentum. It is important to stress that due to the scaling symmetry of Eqs. (2) and (6), only one such sequence needs to be computed (see, e.g., Aksenov & Blinnikov, 1994).

In Figs. 1–3 we plot the shapes of three of the obtained configurations. Because Kong et al. (2015) discuss a model of $\alpha$ Eri with an ellipsoidal surface, we believe that a comparison of the numerical solution and their analytic approximation for the observed parameters of $\alpha$ Eri would be in order. Such a comparison is shown in the first of the presented figures.

Our numerical solution from Fig. 1 was computed with the observed mass $M=4.9M_{\odot}$, the equatorial radius $R_{e}=12.0R_{\odot}$, and the angular velocity $\Omega=2.9725\times 10^{-5}\mathrm{s}^{-1}$, as quoted by Kong et al. (2015, Table 1). The values of solar parameters used in the code are: $M_{\odot}=1.98855\times 10^{33}$ g, $R_{\odot}=6.95700\times 10^{10}$ cm. The assumed value of the gravitational constant reads in cgs units $G=6.67384\times 10^{-8}\,\mathrm{cm^{3}s^{-2}g^{-1}}$.

Our computed ratio of polar and equatorial radii for the $\alpha$ Eri model is Rp/Re = 0.624, and the polytropic constant $K=1.54\times 10^{16}\,\mathrm{cm}^{5}\,\mathrm{g}^{-1}\,\mathrm{s}^{-2}$. Kong et al. (2015) report the value $K=1.7500\times 10^{16}\,\mathrm{cm}^{5}\,\mathrm{g}^{-1}\,\mathrm{s}^{-2}$ and the value of the eccentricity $\mathcal{E}=0.745356$, which, up to the used precision, is equivalent to $R_{p}/R_{e}=\sqrt{1-\mathcal{E}^{2}}=2/3$. In our computation the value of $K$ is obtained together with the structure of the configuration as one of the unknowns (Eriguchi & Müller, 1985). The obtained value differs by 12% from that of Kong et al. (2015) because of the differences in the density stratification. The obtained shape of the boundary surface is plotted in Fig. 1. For comparison, we also plot an ellipse with the semi-axes that were obtained by Kong et al. (2015) and an ellipse fitted to the obtained shape, i.e., an ellipse with semi-axes equal to the values of $R_{p}$ and $R_{e}$ of our numerical model. There is a clear discrepancy between all these shapes.

The solution presented in Fig. 1 here is relatively close to the limiting configuration characterized by the maximal allowed (critical) rotation. Finding the ratio of $R_{p}/R_{e}$ corresponding to the critical configuration is a subtle numerical issue. James (1964, Table 4) and Hurley & Roberts (1964, Table 1) give the values 0.558 and 0.493, respectively. Cook et al. (1992, Table 7), who consider general-relativistic configurations report a Newtonian-limit polytropic solution with $n=1$ and $R_{p}/R_{e}=0.552$. The shape of an (almost) critically rotating $n=1$ polytrope found in our numerical computations is shown in Fig. 2. In this example we keep the same values of the total mass $M=4.90M_{\odot}$ and the polytropic constant $K$ as before. The remaining parameters of this solution are: $R_{p}=7.42R_{\odot}$, $R_{e}=13.3R_{\odot}$, $R_{p}/R_{e}=0.557$, $\Omega=3.004\times 10^{-5}\,\text{s}^{-1}$.

As expected, the ellipsoidal approximation is much better for slowly rotating polytropes. A comparison between the actual shape and the matching ellipsoid for $\Omega=0.2146\times 10^{-5}\,\text{s}^{-1}$ is shown in Fig. 3. Here again we keep the same values of $M$ and $K$. The remaining parameters read: $R_{p}=8.19R_{\odot}$, $R_{e}=9.47R_{\odot}$, $R_{p}/R_{e}=0.865$. Note that all solutions presented in Figs. 1–3 form a sequence characterized by the same mass and the same equation of state, but different angular velocities.

The numerical grid used in our computations is also depicted in each of Figs. 1-3. We use 31 grid points in the angular direction and 23 in the radial one. The highest order of Legendre polynomials used in the numerical code is 40.

4 Final remarks

Polytropic solutions with index unity serve as suitable approximations of the structure of Jovian planets (Hubbard, 1974). In a recent work Cao & Stevenson (2015) provide an interesting comparison between parameters of a polytropic model and the data for Jupiter. [Note that it is explicitly stated in Cao & Stevenson (2015) that the outer boundary in the polytropic model differs from an exact ellipsoid of revolution.] The values of $n$ close to unity can be probably also assumed for brown dwarfs and some neutron stars. For the star in the radiative equilibrium (like $\alpha$ Eri) the value $n\approx 3$ seems to be more appropriate. This fact was actually admitted by Kong et al. (2015).

It is possible to search for exact solutions representing rotating polytropes of index $n=1$ assuming an ellipsoidal boundary surface and differential rotation. Such solutions were constructed by Cunningham (1977) (He uses the term ‘spheroids’ referring to ellipsoids of revolution).

It is also worth pointing out that because every solution of Eqs. (2) and (6) satisfies Eq. (7), the latter can be actually used to infer important characteristics of solutions of Eqs. (2) and (6). This is also true in nonlinear cases (e.g., for polytropes with $n\neq 1$); results of this type were obtained by Mach & Malec (2011).

Stability of the rotating $n=1$ polytropes was investigated by Bini et al. (2008).

Acknowledgements

We are grateful to Georg Horedt for his comments on the first version of this article. A.O. was supported by the Kosciuszko Foundation.