Introduction

During the last decades, considerable efforts have been made to obtain exact analytical solutions suitable to modeling axisymmetric thin disks, within the framework of Newton’s and Einstein’s theories of gravity. Such models are of astrophysical interest because they can be used to model accretion disks, galaxies in thermodynamic equilibrium or galaxies with black holes centers (^{Bicák, Lynden-Bell and Pichon, 1993}, ^{Ledvinka, Zofka and Bicák, 1998}). Moreover, the addition of electromagnetic fields in those space-times allows studying neutron stars formation, white dwarfs, and quasars (^{Muñoz et. al., 2011}, ^{Alpar, 2001}).

Since the seminal works on exact solutions representing static thin disks carried out by Bonnor and Sackfield (^{Bonnor and Sackfield, 1968}) and Morgan and Morgan (^{Morgan and Morgan, 1969}), more realistic models have been proposed (^{Pichon and Lynden-Bell, 1996}, ^{González and Letelier, 2000}). The superposition of static and stationary thin disks with black holes at the center, has been considered by Lemos, Letelier and Semerak (^{Lemos and Letelier, 1993}, ^{Semerák, 2002}b, ²⁰⁰⁴). Vogt and Letelier studied the inclusion of electro- magnetic fields into thin disks made of dust (^{Vogt and Letelier, 2004}a) and of charged perfect fluid (^{Vogt and Letelier, 2004}b). Also interesting is the case of thick disks proposed by Gonzalez and Letelier (^{González and Letelier, 2004}), who extended the DCR method to include thick disks in their models (For the interested reader, a complete review of the state-of-the-art on relativistic disks was made by ^{Semerák, 2002}a).

Stability is an essential criterion to determine whether or not a model can be applied to describe an astrophysical system present in nature. In general, there are two approaches to study the stability of relativistic disks: The first option is based on analyzing the stability of particle orbits along geodesics (see e.g., (^{Letelier, 2003}) and (^{Vogt and Letelier, 2003})), while the second option consists in perturbing the energy-momentum tensor (see e.g., (^{Seguin, 1975})). From a theoretical point of view, the latter option is more rigorous, because in this case the collective behavior of the particles is taken into account. Working on this line, ^{Ujevic and Letelier (2004}) investigated the stability of three particular models for relativistic thin disks, performing a first-order perturbation analysis with variable coefficients. However, as a pathological result, the authors find that the total mass of the disk depends on the parameters of the DCR method, such that the boundary conditions are also dependent on these parameters.

With the aim to avoid the undesired dependencies between parameters and to observe the possible changes in the stability, in the present paper we redo the calculations made by ^{Ujevic and Letelier (2004)} using the Komar mass definition instead of the mass definition (along the paper we will call it total surface mass) introduced in ^{Vogt and Letelier (2003)}. The new results let us to choose the boundary condition for the cutoff radius such that it takes the maximum value allowing an appreciable and well-behaved perturbation on the disk.

Derivation of First-Order Perturbation Equations

Following ^{Ujevic and Letelier (2004)}, and for the sake of self-consistency of the current paper, in what follows we present the derivation of the first-order perturbation equations for relativistic thin disks. Let us start considering that the energy-momentum tensor for an isotropic fluid with a discoid shape and without heat flow can be written as

where δ denotes the Dirac delta function and

with σ the surface energy density, ^_pr and ^_pφ the radial and azimuthal pressure, respectively, and ^_Uµ , ^_Xµ , and ^_Yµ the nonzero components of the orthonormal tetrad. Assuming that the first-order perturbations in the Einstein field equations do not modify the background metric, the perturbed equation for the energy-momentum tensor reads as

Introducing the definition of energy-momentum tensor for a thin disk (1), in the perturbed equation (2), and integrating with respect to the coordinate z, we find

As a consecuence of the DCR method, the metric components should only depend on r and |z| (^{Vogt and Letelier, 2003}). Additionally, if we define the value of |z|, _z at z = 0 equals zero, the perturbed equation reduces to

Taking into account that the perturbed vectors must satisfy the orthonormality condition, and assuming that ^_δXφ = 0 (since the four-velocity and the thermodynamical variables do not depend on this quantity), we obtain

Due to the fact that the metric is static and axisymmetric, and the lack of z-dependences in ^_Tµν , all the coefficients depend only on the radial coordinate; therefore, the general perturbation can be chosen as

with k the wave number and w the angular frequency of the perturbation. We will focus on the term ^_δξµ (r), which has a significant role for the system stability.

Therefore, after substituting the perturbation (6) into equation (4), and replacing the conditions (5) in the resulting expression, the respective perturbed equations for t, r, and φ, read as

Finally, from the set of differential equations presented above, and the equation of state of a perfect fluid, δp = δσ( ^{_{p, r}} / ^{_{σ, r}} ), the differential equation for the perturbed energy density takes the form

whose coefficient A, B, and C, are given by

where ^_Ai , ^_Bi , and ^_Ci , are the factors multiplying the perturbed variables in equations (7), (8), and (9), respectively, with the indexes in numerical order according to their order of appearance.

Due to the cumbersome form of equation (10), it must be solved numerically. For this purpose, we shall impose Dirichlet boundary conditions, one at the center of the disk and the other one at the final boundary of the domain. However, it should be noted that the disk has an infinite radial extension, for this reason, it is necessary to introduce a cutoff on the radial coordinate.

Thermodynamic variables

In the previous section, we explicitly wrote down the perturbed equations for thin disks, nevertheless, such expressions are given in terms of the thermodynamic variables of the fluid, which requires finding explicit formulas for the surface energy density and the pressures on the disk. To this end, let us represent the matter distribution on a hypersurface Σ, defined by the function l(xα) = z, which divides the space-time into two regions: M ⁺ on top and M ⁻ at the bottom. Therefore, the normal vector to the hypersurface Σ is given by and the components of the metric tensor must satisfy

such that

where , should be understood as the metric tensors for the regions defined by z > 0 (M+) and z < 0 (M−), respectively.

By taking the limit z → 0, the discontinuities in the first derivatives of the metric tensor take the form

with bμν, the jump of the first derivative through Σ. Using the distributions method (^{Papapetrou, 1968}, ^{Poisson, 2004}), the metric can be written as

where Θ(l) is the usual Heaviside function. Performing the derivative of equation (14) we obtain

The last term in the right-hand side of equation (15) is singular, consequently the Christoffel symbols would not be defined as a distribution. An alternative way to avoid this difficulty, is to impose that the metric is continuous on the hypersurface, i.e., such that

whose derivative reads as

The corresponding Riemann curvature tensor is given by

where are the tensors defined in M^± and

with

From the last expression, it is clear that Riemann tensor, the Ricci tensor, and Ricci scalar on the hypersurface are

On the other hand, the energy-momentum tensor can be expressed as

where is the energy-momentum tensor associated with the hypersurface, and are the energy-momentum tensors associated to M±. Hence, we can write the Einstein field equations as follows

Given that the only source of gravitational field is a thin distribution of matter, the energy-momentum tensor satisfies

which, in quasi-cylindrical coordinates and considering the hypersurface Σ defined by the function , with normal vector , takes the form

From the above equation, and as noted by ^{Vogt and Letelier (2003)}, it can be shown that the non-zero components of the surface energy-momentum tensor are

Mass Definition

The mass concept in general relativity is not unique and there are several different definitions that are applicable under different circumstances. In the case of stationary spacetimes, the commonly accepted definitions are: the total volumetric (or surface) mass, Komar, ADM and Bondi-Sach masses. However, it is a well-known fact that for stationary spacetimes the ADM and Bondi-Sach masses are exactly alike (^{Poisson, 2004}). Besides, in the asymptotically flat case, it can be shown that the ADM and the Komar masses are equivalent (^{Jaramillo and Gourgoulhon, 2009}). Then, in the stationary, axisymmetric, asymptotically flat, vacuum solutions of Einstein’s equation, there are only two possibilities to choose from, the total volumetric (or surface) mass or the Komar mass. As commented in the introduction, the main difference between the present paper and the one by ^{Ujevic and Letelier (2004)}, is that they used the total surface mass, while along this paper, we use the Komar mass definition, which is independent of the methods used to build the particular solutions of Einstein’s equation.

A necessary condition for the definition of mass, is that must not involve any dependence with the specific choice of coordinates. This property is achieved for stationary and axially symmetric spacetimes, through the Komar formula Mk, which reads as

where Σ is spacelike hypersurface, is timelike vector normal to Σ, is timelike killing vector and is the determinat of metric associated to Σ. For a static and diagonal metric, the following relation holds

,

such that the Komar mass takes the form

The Komar mass definition differs from the total surface mass in the fact that the former one considers all the contributions to the energy-momentum tensor, while the surface mass, given by

takes only into account the surface energy density.

Case 1: Isotropic Schwarzschild Thin Disk

The isotropic Schwarzschild thin disk metric in quasi-cylindrical Weyl-Papapetrou coordinates was obtained by ^{Vogt and Letelier (2005)}, and can be written as

where m is a positive constant and R2 = r2 + (|z| + a)2.

Thus, the corresponding expressions for the pressure components and surface energy density (see e.g., ^{Vogt and Letelier (2003)}), are obtained by using equations (25, 26 y 27),

with

Let us define a general cutoff radius as the radial distance at which the matter within the thin disk formed up to such radius corresponds to the n% of the total matter of the infinite disk, i.e.

with MT the total mass of the infinite disk and Mrc the mass up to the cutoff radius. As noted in the previous section, the total mass value may depend on the considered definition, such that the cutoff radius could depend also on this choice. From (30), (31) and (32), and in accordance with Eq. (12) of ^{Vogt and Letelier (2003)}, we find

and

On the other hand, from the Komar mass definition (29) and using (31), (32) and (33), we get

and

Hence, the resulting expression for the cutoff radius when considering Eq. (30) is

while using Eq. (29), we get

In Fig. 1 we show a comparison between the cutoff radii and , using different values of the parameters m and a. Taking into account that values of n smaller than 0.9 are not used nor physically appropriate, we plot the difference in the range

Figure 1 Comparison between the cutoff radii rc and r*c using different values of the parameters m and a. 

In all that follows, except where especially noted, we define the cutoff radius rc using the 95% of the total mass of infinite disk,¹ and the remaining 5% of matter is distributed along the plane z = 0 from rc to infinity. Moreover, we assume that at r = 0, the perturbation equals 10% of the unperturbed energy density value and at r = rc the perturbation vanishes, i.e.

By setting m = 0.5 and a = 0.6, which are the set of parameters that best fit the expected energy density profile in a realistic model, we find . In Fig 2 we present the numerical solution to the differential equation (10) using the Komar mass definition (upper panel) and the total surface mass (lower panel). In this figure, we show the perturbed energy density profile for different values of w. It can be seen that the parameter w is proportional to the frequency of oscillations within the disk, i.e. the number of ring-like structures increases with increasing w. Also, it is worth noting that the oscillations quickly decay to zero regardless of the oscillation modes, however the oscillation amplitude is smaller for Eq. (29) (upper panel) than for Eq. (30) (upper panel).

Figure 2 Perturbed energy density profiles δσ˜ = √gzz δσ for different values of w using Eq. (29) (upper panel) and Eq. (30) (lower panel) in the isotropic Schwarzschild thin disk. 

On the other hand, in Fig. 3 we present the profile of the perturbed pressure. When fixing w and varying k, the amplitude for δp has an oscillatory behavior similar to the one of δσ. From Fig. 3 it can be seen that the oscillation amplitude is not only smaller for Eq. (29) (upper panel) than for Eq. (30) (upper panel), but also disappears in approximately the middle of the disk. It means that for the set of parameters here considered, the density and pressure in the disk are stable independently of the mass definition used.

Figure 3 Perturbed pressure profiles δp˜ = √gzz δp for the isotropic Schwarzschild thin disk using Eq. (29) (upper panel) and Eq. (30) (lower panel). 

A full analysis of the system requires also a comparison between the amplitude of the perturbed velocities with the escape velocity of the constituents in the disk. In Fig. 4 we show the profiles for the perturbed radial δUr velocities using Eq. (29) (upper panel) and Eq. (30) (lower panel), setting the respective cutoff radii for the 95.776% of the total mass of infinite disk². Both velocities exhibit an increase in the frequency of oscillations with increasing w; moreover, we see that the envelopes of the oscillating functions increase at the external radial boundary, however, the perturbed radial velocity grows faster in the lower panel than in the upper panel.

Figure 4 Profiles of the perturbed radial velocities for different oscillation modes in Schwarzschild isotropic thin disk using Eq. (29) (upper panel) and Eq. (30) (lower panel). The red dashed line denotes the value of the respective escape velocities. 

Concerning the escape velocity ^_νe , it is well known that it should be larger than the perturbed radial velocity, otherwise, the model will not have any astrophysical validity. The escape velocity required to overcome the force generated by a gravitational potential Φ can be calculated as ^_νe = √2Φ, which from the weak-field approximation is given as 2Φ = ^_gtt − 1. So, from (31), we find

For the set of parameters a = 0.6 and m = 0.5, the resulting escape velocities for Eq. (29) and Eq. (30) are ^_νe = 0.26 and ^_νe = 0.28, respectively. From Fig. 4 it is clear that using Eq. (29) the perturbed radial velocity is always less than the escape velocity (red dashed line), such that all particles remain inside the disk. Conversely, using Eq. (30) the perturbed radial velocity can be greater than the escape velocity (red dashed line), meaning that the particles may escape of the disk.

In Fig. 5 we show the profiles of the perturbed azimuthal velocities. In both cases the velocities are stable and some orders of magnitude smaller than the escape velocity. These results ratify that, under first order perturbations of the form (6) and using the Komar mass definition, the isotropic Schwarzschild thin disk has a stable behavior for a large set of parameters such that it can be used to describe astrophysical models.

Figure 5 Profiles of the perturbed azimuthal velocities for different oscillation modes in Schwarzschild isotropic thin disk using Eq. (29) (upper panel) and Eq. (30) (lower panel). 

Case 2: Chazy-Curzon Thin Disk

The Chazy-Curzon thin disk metric in Weyl coordinates is described as (^{Bicák, Lynden-Bell and Katz, 1993})

where the metric functions Φ y Λ are given by

with R ² = r ² + (|z| + a)² . From (42) and (43), and equations (25), (26), and (27), we obtain the expressions for the pressure and the surface energy density of the Chazy-Curzon thin disk,

where R₀, Λ₀ and Φ₀ are the functions evaluated at z = 0.

Substituting equations (43), (44) and (45), into the Komar mass equation (29), we obtain, accordingly, ^_Mk = m. A particular set of values that satisfy the energy conditions are a = 0.4 and m = 0.5, in which case ^_rc ≈ 8. In Fig. 6 we show the nu- merical solution to the differential equation (10) with boundary conditions at r = 0 and r = rc = 8. As can be noted, regardless of the values of w (upper panel) and k (lower panel), the profiles of the perturbed energy density grow rapidly to infinity. To rule out that the above effect is due to the chosen values of a, we perform an analysis varying δσ in the interval 0 ≤ m/a ≤ 1.3, observing the same tendency. The same procedure was performed using the total surface mass definition, obtaining the same behavior. These results show that the energy density, and hence the pressure and velocities, exhibit instabilities for this model and are independent of the mass definition used in the calculations.

Figure 6. Profiles of the perturbed energy density δσ˜ in the Chazy-Curzon thin disk using Eq. (29), for different values of w and k. The parameters have been set as a = 0.4 and m = 0.5. The insets show that for finite values of r, the perturbations tend to infinity. 

Case 3: Zipoy-Voorhees Thin Disk

The Zipoy-Voorhees thin disk metric in Weyl coordinates has the form (^{Ujevic and Letelier, 2004})

with metric functions Φ and Λ, given by

where Hence, by applying the same procedure used for the previous disks models, we obtain expressions for the surface energy density and the azimuthal pressure,

with _{^Ra 0}, _{^Rb 0}, Λ₀ and Φ₀ the respective functions evaluated at z = 0. Just like in the Chazy-Curzon model, the thin disk described by the Zipoy-Voorhees metric has no radial pressure. Replacing (47), (48) and (49) into equation (29), we can calculate the Komar mass for the Zipoy-Voorhees thin disk as ^_MT = m.

To be consistent with the previous models, the calculations were performed for a particular set of values that satisfy the energy conditions, m = 0.5, a = 1 and b = 2.15, such that the cutoff radius takes the value rc ≈ 31. In the same way as cases 1 and 2, in Fig. 7 we show the numerical solution to the differential equation (10) with boundary conditions at r = 0 and r = rc = 31.

Figure 7 Profiles of the perturbed energy density δσ˜ for different values of w (upper panel) and k (lower panel) in the Zipoy-Voorhees thin disk. The parameters have been set as a = 1, b = 2.15 and m = 0.5. The inset shows the tendency of the curves for larger scales. 

It can be seen from Fig. 7 that there exist instabilities for the Zipoy-Voorhees thin disk strongly amplified before they reach 10% of its cutoff radius, and regardless of the angular frecuency w (upper panel) or wave number k (lower panel) the perturbed energy density tends to infinity. The results for the Zipoy-Voorhees thin disk when using the total surface mass are not shown because they exhibit practically the same behavior than using the Komar mass. The independence of the results with the particular chosen values of a and b, was also analyzed observing the same tendency. Given the instability in the surface energy density, the other thermodynamic variables of the thin disk also are unstable.

Concluding Remarks

In this paper, we have reviewed the results obtained by ^{Ujevic and Letelier (2004)} on the stability of three particular models, using an alternative mass definition, i.e., the Komar mass. This formal concept of mass can be defined in any stationary spacetime, so it is applicable to all three particular thin disk models under consideration: the isotropic Schwarzschild, Chazy-Curzon, and Zipoy-Voorhees metrics. Assuming that the disks are made of a perfect fluid, we found the differential equation for the perturbed energy density and the nonzero components of the surface energy-momentum tensor. The differential equation was numerically solved, defining a cutoff radius ^_rc , such that the matter up to ^_rc is approximately 95% of the mass of the infinite disk and the remaining 5% of matter is distributed from ^_rc to infinity. Moreover, we assumed that at the center of the disk the perturbation equals 10% of the unperturbed energy density, and at the external cutoff radius the perturbation vanishes.

Once we derived the expressions for the pressure components, the surface energy density, and the metric functions for each particular case, we calculated the Komar mass. As a general result, we find that the mass parameter in each one of the metrics equals the Komar mass for the disk, while by using the total surface mass definition, as is the case of the reviewed paper, the total mass depends on the parameters of the DCR method. The use of the Komar definition lets us set the physical parameters that best fit the expected energy density profile in a realistic model and simultaneously satisfy the energy conditions. With this result, the cutoff radius only depends on the thermodynamic variables and the free parameters of each

As the main finding, we found that the cutoff radius is larger for the Komar mass definition than for the total surface mass. This result let us to increase the number of parameters that give place to stable Schwarzschild thin disk models. On the other hand, the Chazy-Curzon and Zipoy-Voorhees thin disk models are not stable under first- order perturbations, because the thermodynamic variables and fluid velocities tend to infinity for finite values of the radial coordinate. Such instabilities are a consequence of the lack of radial pressure and are not related to the definition of mass. We also have shown that the infinite radial extension of the disk can be the reason for the instability, as hypothesized by Ujevic et al., nevertheless, some instabilities can be artificially introduced into the model due to the use of a non-appropriate mass definition in the calculations.