# Quasi-radial modes of rotating stars in general relativity

## 1 Introduction

Radial oscillations of non-rotating relativistic stars have been studied for over thirty years. Methods for obtaining their spectra have been well established (Bardeen, Thorne & Meltzer 1966 ; see also Chapter 26 of Misner et al. 1973 and the references therein), and have been applied to several equations of state (see for example Meltzer & Thorne 1966). These works were mainly motivated by consideration of stellar stability because general relativistic effects tend to destabilize stellar models .

On the other hand, the effect of rotation on stellar oscillations is less well understood in a general relativistic context. As in the non-axisymmetric mode case, the slow rotation approximation has been the only accessible way for investigating the eigenmode behaviour of rotating stars . Recently, numerical relativistic hydrodynamic codes have been developed by several authors and some numerical simulations of rapidly rotating stars have been carried out. Stergioulas et al.  and Font et al.  have shown that initial small perturbations around an equilibrium star evolved to a superposition of normal mode oscillations (Note that their hydrodynamic simulation is done in the fixed background spacetime. On the other hand, Shibata et al.  have solved the full system of Einstein equations to investigate the dynamical stability of rapidly rotating stars.).

Although the excitation and evolution of these modes in realistic situations should be investigated by time dependent hydrodynamic simulations, it is also important to study the mode behaviour along rotational equilibrium sequences by using linear perturbation theory.

So far we have studied a few sets of non-axisymmetric eigenmodes of rotating stars in general relativity (f-modes by Yoshida & Eriguchi 1997,1999 ; p-modes by Yoshida 1999 (unpublished; presented at the 9th Yukawa International Seminar ’BLACK HOLES AND GRAVITATIONAL WAVES \- New Eyes in the 21st Century-’) ). These results have been obtained within the Cowling approximation in which Euler perturbations of the metric coefficients have been neglected (see McDermott et al. 1983 and Finn 1988 for a definition ; see Lindblom and Splinter 1990 for the accuracy of the method when applied to non-radial modes of spherical stars). Apart from some low order modes, these results are in good agreement with those of the full perturbation theory including metric perturbations (For comparison of the eigenfrequencies for slowly rotating stars, see Yoshida & Kojima 1997 ; for comparison of the neutral points of the CFS instability, see Yoshida & Eriguchi 1999 which compare the results with the one obtained by full computation of Stergioulas & Friedman 1998 and Morsink et al. 1999)

It is therefore natural to expect that the Cowling approximation could also be successfully applied to the quasi-radial modes which are the smooth extensions of the radial modes of spherical stars to rotating stars. In the present paper, we study quasi-radial modes by using the Cowling approximation. Contrary to the expectations, our results indicate that computations with this approximation can not reproduce the relativistic instability of spherical stars. This is plausible because the instability is essentially caused by the loss of balance between gravity and the pressure gradient, and in calculations of it even the small corrections of gravity cannot be neglected. Moreover, the phase cancellation of the perturbed gravitational potentials, which may be effective in the case of non-axisymmetric modes, cannot be expected to happen for radial modes.

See the Appendix of the present paper for the comparison of two methods in the case of radial modes of non-rotating stars.

Although the validity of Cowling approximation for rotating stars is not fully assessed, we here expect that at least a qualitative picture of the eigenmode dependence on stellar rotation could be studied by this approximation.

## 2 Results

The equation of state used here is the polytropic one,

$$p=\kappa \rho ^{1+\frac{1}{N}},ϵ=\rho +Np$$
$`(1)`$   
where $`\rho `$, $`ϵ`$ and $`p`$ are the rest mass density, energy density and pressure of the stellar matter, respectively. Geometrized units, $`c=G=1`$, are adopted in this paper as well as $`M_{}=1`$, following Font et al. . The constant $`N`$ is the polytropic index. The adiabatic exponent of the perturbed matter is assumed to coincide with $`1+1/N`$. The factor $`\kappa `$ is another constant.

Each equilibrium sequence is computed with $`\kappa `$ and $`N`$ fixed.

In the present study polar-like coordinates are used and the metric components are written as:

$$ds^2=e^{2\nu }dt^2+e^{2\alpha }(dr^2+r^2d\theta ^2)+e^{2\beta }r^2\mathrm{sin}^2\theta (d\varphi \omega dt)^2.$$
$`(2)`$   
The rotational axis is located at $`\mathrm{sin}\theta =0`$.

The coordinates used in the actual numerical computation are surface-fitted ones $`(r^{},\theta ^{})`$ which are defined as:

$$r^{}=r/R_s(\theta ),\theta ^{}=\theta ,$$
$`(3)`$   
where $`r=R_s(\theta )`$ is the form of the stellar surface in equilibrium.

The numerical method used here is basically the same as that in Yoshida & Eriguchi  where non-axisymmetric modes were investigated. A minor modification is needed to obtain the axisymmetric modes. In the case of non-axisymmetric modes, the Eulerian variable $`\delta p/(ϵ+p)`$ is explicitly set to zero at the centre of the star ($`\delta p`$ is the Eulerian variation of the pressure). In the case of axisymmetric modes, however, this is not the case since the regularity of the solution requires $`(\delta p)/r`$ to be zero at the stellar centre. Therefore we simply modify the finite-difference scheme at the innermost grid points in our numerical code. Moreover to avoid the coordinate singularity on the rotation axis, points on the axis are excluded from the computational region.

Most of the results shown in this paper are computed with a resolution of 40 uniformly distributed gridpoints in the radial $`r^{}`$ direction and 10 in the angular $`\theta ^{}`$ one. The computational region is a quarter of the meridional section of stars, thus the range of the radial and the angular coordinates are $`0r^{}1`$, $`0\theta ^{}\pi /2`$.

In Figs. 1 and 2 the eigenfrequencies of the axisymmetric modes are plotted against the rotational frequency of the equilibrium model.

The model parameters are tabulated in Table 1.

The eigenfrequency and the rotational frequency are normalized by $`\sqrt{\rho _c/4\pi }`$, where $`\rho _c`$ is the central rest mass density of the models. The sequences ’$`F`$’, ’$`H_1`$’ and ’$`H_2`$’ are the fundamental, first and second overtones of the quasi-radial modes, respectively.

Some of axisymmetric f-modes and overtones of p-modes are also plotted. <sup>1</sup><sup>1</sup>1Our numerical code assumes reflection symmetry of the eigenmodes with respect to the equatorial plane of the equilibrium star. As a result, modes with odd integer $`l`$ cannot be computed. Considering the order of the even $`l`$ modes in the figure, sequences of $`{}_{}{}^{5}\text{f}`$ and $`{}_{}{}^{3}\text{p}_{1}^{}`$ may be located somewhere between the the corresponding p-modes with $`l=2,4`$. These are the continuation of the corresponding f- and p-modes with $`l=L`$ and $`m=0`$ where $`l`$ and $`m`$ are the indices of the ordinary spherical harmonics $`Y_{lm}(\theta ,\phi )`$. The label $`{}_{}{}^{L}\text{p}_{n}^{}`$ refers to a mode corresponding to the p<sub>n</sub>-mode with degree $`l=L`$ and order $`m=0`$ in the non-rotating limit. Similarly, the mode with the label $`{}_{}{}^{L}f`$ is the f-mode with $`l=L,m=0`$ in the non-rotating limit.

Generally the quasi-radial mode sequence encounters other sequences of f- or p-modes. It is seen that the so-called avoided crossing occurs on these sequences (see Fig.3). This seems to be the general relativistic extension of what has been found for the oscillations of rotating Newtonian stars. See for example Clement  and Unno et al.  for the Newtonian cases.

Fig.3 clearly shows the presence of the avoided crossing whereby two eigenfrequency curves approach smoothly, and then depart from each other without crossing. At the point of closest approach, the characteristics of the modes on each sequence exchange. Thus the sequence of a mode with given characteristics has a discontinuity there.

There are several discontinuities along a sequence whose number depends on the mode order as well as on the polytropic parameters of the equilibrium star. For the selected polytropic parameters, discontinuities appear in the rapidly rotating models whose rotational frequencies are over $`80\%`$ of the mass-shedding limit.

## 3 Some remarks on the analysis

### 3.1 Convergence and accuracy

Since we use a finite number of grid points, our results necessarily contain some errors. To see how much the results differ when the number of grid points is changed, we compute the fundamental quasi-radial mode sequence by using different grid numbers $`(M,N)=(40,10),(40,20),(80,10)`$ where $`M`$ is the number of grid points in $`r^{}`$-direction and $`N`$ is that in $`\theta ^{}`$-direction. By extrapolating three eigenfrequencies obtained from these grid numbers, $`\nu _{[40,10]},\nu _{[40,20]},\nu _{[80,10]}`$, we estimate the converged value of the eigenfrequency $`\nu _0`$ in the limit of infinitesimal mesh size to be:

$$\nu _0=2(\nu _{[40,20]}+\nu _{[80,10]})3\nu _{[40,10]}.$$
$`(4)`$   
The ratio of $`\nu _{[M,N]}`$ to $`\nu _0`$ can be treated as a measure of convergence of the results when the grid number is changed (Fig.4). As can be seen from this figure, eigenfrequencies obtained by using $`(M,N)=(40,10)`$ mesh numbers (our standard resolution) agree with the converged values to within 3 percent.

It is noted that the relative error of $`\nu _{[80,10]}`$ is smaller than that of $`\nu _{[40,20]}`$ in the slowly rotating cases, however for rapidly rotating cases, the error of $`\nu _{[80,10]}`$ becomes larger. This is because for the rapidly rotating cases, the deformation of the stars from the spherical configuration is so large that the lack of angular resolution becomes the main source of inaccuracy.

To check the accuracy of our two dimensional (2D) code, we have compared the results obtained here with those obtained with the one dimensional (1D) code described in the Appendix (Table 4). In the 1D code we employ the standard scheme to solve the eigenvalue problem of the linear ordinary differential equation within the Cowling approximation.

As seen from this table, two results agree well to within several percent.

### 3.2 Eigenfunctions and classification of modes

Comparing with the eigenfunctions of non-axisymmetric modes, we notice that the eigenfunctions of quasi-radial modes change their shapes significantly along rotational equilibrium sequences. For example, the radial distributions of the Eulerian pressure perturbation and the radial component of the velocity perturbation change considerably near the equatorial plane of the star. The number of radial nodes of these functions increases as the stellar rotation rate increases. On the other hand, the overall radial dependence of the $`\theta `$-component of the velocity perturbation changes little as the star spins up. Therefore we can use the shape of this function as a tracer of the selected mode along the equilibrium sequence. The nomenclature of mode sequences is based upon the behaviour of the mode of a non-rotating star to which the sequence is continued smoothly.

### 3.3 Realistic neutron star models

In addition to the polytropic case presented here, we have tried to compute quasi-radial modes of realistic neutron stars by using some of the candidate zero temperature equations of state. However it was rather difficult to obtain full sequences of these modes with our method. As the star begins to rotate, the convergence to the quasi-radial modes suddenly becomes much more difficult. This may partly be due to the fact that these modes are sensitive (as compared with the non-radial modes) to the surface condition of the equilibrium star. Unfortunately, in the case of more realistic EOS, the adiabatic exponent is subject to large oscillations, becoming negative in some parts. These large variations decrease considerably the accuracy of our method which then becomes inadequate.

## Acknowledgment

We thank Luciano Rezzolla for detailed comments on the manuscript, and John Miller for his help to improve it. We also thank Kōji Uryū for useful discussions.

## Appendix A Comparison of radial modes in the full theory and in the Cowling approximation

To test the validity and accuracy of the Cowling approximation for a spherical configuration, we here present a comparison between the results obtained by the full perturbation theory and by the Cowling approximation for low order radial modes.

### A.1 Equilibrium model

The space-time of the equilibrium star is characterized by the following metric:

$$ds^2=e^{2\nu (r)}dt^2+e^{2\lambda (r)}dr^2+r^2(d\theta ^2+\mathrm{sin}^2\theta d\phi ^2).$$
$`(A1)`$   
The metric coefficients, stellar pressure $`p`$ and energy density $`ϵ`$ are obtained by integrating the standard set of equations with regular boundary conditions at the stellar centre:

$$\frac{dp}{dr}=\frac{(ϵ+p)(m+4\pi pr^3)}{r(r2m)},$$
$`(A2)`$   
$$\frac{dm}{dr}=4\pi r^2ϵ,$$
$`(A3)`$   
$$\frac{d\nu }{dr}=\frac{m+4\pi pr^3}{r(r2m)},$$
$`(A4)`$   
where $`m(r)`$ is defined from $`e^{2\lambda }=(12m/r)^1`$.

The equation of state is assumed to be polytropic, $`p=\kappa \rho ^{1+\frac{1}{N}}`$, where $`\kappa `$ and N are constants, and the rest-mass density $`\rho `$ is related to $`ϵ`$ by $`ϵ=\rho +Np`$.

### A.2 Equations for radial oscillations

For radial oscillations of spherical stars, only the Lagrangian displacement function $`\xi (r)`$ is needed to describe the physical perturbation (see Chapter 26 of Misner et al. 1973).

The equation of motion of the displacement ($`r^2e^\nu \xi \zeta `$) in the full perturbation theory is:

$$\zeta ^{\prime \prime }+A_f\zeta ^{}+B_f\zeta =0,$$
$`(A5)`$   
where $`A_f`$ and $`B_f`$ are defined as:

$$A_f\frac{p^{}}{p}\frac{2}{r}+\lambda ^{}+3\nu ^{},$$
$`(A6)`$   
and

$$B_f\frac{ϵ+p}{\mathrm{\Gamma }p}\left[(\nu ^{})^2+4\frac{\nu ^{}}{r}8\pi pe^{2\lambda }+\sigma ^2e^{2\lambda 2\nu }\right].$$
$`(A7)`$   
The prime after a variable refers its derivative with respect to $`r`$. Here $`\sigma `$ is the frequency and the adiabatic exponent $`\mathrm{\Gamma }`$ is defined by:

$$\mathrm{\Gamma }\frac{ϵ+p}{p}\frac{\mathrm{\Delta }p}{\mathrm{\Delta }ϵ},$$
$`(A8)`$   
where $`\mathrm{\Delta }`$ represents the Lagrangian perturbation. In general the adiabatic exponent need not coincide with $`1+1/N`$.

In the Cowling approximation, the equation of motion is instead ($`r^2e^\lambda \xi \eta `$):

$$\eta ^{\prime \prime }+A_c\eta ^{}+B_c\eta =0,$$
$`(A9)`$   
where $`A_c`$ and $`B_c`$ are defined as:

$$A_c\frac{p^{}}{p}\frac{2}{r}\lambda ^{}+\nu ^{},$$
$`(A10)`$   
and

$$B_c\frac{p^{}}{\mathrm{\Gamma }p}\left(\frac{2}{r}+\lambda ^{}\right)+\frac{ϵ+p}{\mathrm{\Gamma }p}(\nu ^{\prime \prime }+\sigma ^2e^{2\lambda 2\nu }).$$
$`(A11)`$   
These equations can be solved by the matching method: i.e., we have to search for $`\sigma ^2`$ which makes the Wronskian of the solutions, obtained by integrations from the stellar centre and from the surface, vanish at some matching point inside the star.

### A.3 Boundary conditions

The boundary condition for the equation of oscillations at the centre of the star is regularity of the variables. It requires $`\zeta ,\eta r^3`$ as $`r0`$.

At the surface of the star, we impose the boundary and regularity conditions which reduce to

$$\zeta ^{}+\frac{ϵ+p}{\mathrm{\Gamma }p^{}}[(\nu ^{})^2+4\nu ^{}r^1+\sigma ^2e^{2\lambda 2\nu }]\zeta =0,$$
$`(A12)`$   
or

$$\eta ^{}+\left[\frac{1}{\mathrm{\Gamma }}(2r^1+\lambda ^{})+\frac{ϵ+p}{\mathrm{\Gamma }p^{}}(\nu ^{\prime \prime }+\sigma ^2e^{2\lambda 2\nu })\right]\eta =0.$$
$`(A13)`$   
### A.4 Results

Keeping the polytropic index $`N`$ fixed, we vary the compactness (mass-to-radius ratio) of the model to construct an equilibrium sequence.

In Figs. A1 and A2 we show the typical sequences of the three lowest order modes. For all of the sequences, the eigenfrequency obtained by the Cowling approximation is larger than that from the full theory: i.e., the Cowling approximation overestimates the stability of the star. Before the turning point,<sup>2</sup><sup>2</sup>2Note that with this parametrization, the turning point does not correspond to the maximum mass configuration. Thus the zeroes of the fundamental radial modes are not the turning points. the two curves are nearly parallel. As expected, the relative error of the results obtained by the Cowling approximation becomes smaller for higher overtones.