EQUILIBRIUM SEQUENCES AND GRAVITATIONAL INSTABILITY OF ROTATING ISOTHERMAL RINGS

In this section, we explore equilibrium sequences of rotating isothermal bodies in the presence of both self-gravity and external gravity. These bodies can take a spheroid-like or ring-like configuration when the total angular momentum is small or large. We assume that equilibrium bodies are rotating rigidly at angular frequency Ω₀ about its symmetry axis that is aligned in the $\hat{z}$-direction. The equation of steady equilibrium then reads

$$\begin{eqnarray}&&{c}_{s}^{2}{\rm{\nabla }}\mathrm{ln}\rho +{\rm{\nabla }}{{\rm{\Phi }}}_{{\rm{eff}}}=0,\end{eqnarray} \tag{ 1 }$$

where ρ is the density, ${c}_{s}$ is the isothermal speed of sound, and Φ_eff is the effective potential defined by

$$\begin{eqnarray}&&{{\rm{\Phi }}}_{{\rm{eff}}}={{\rm{\Phi }}}_{e}+{{\rm{\Phi }}}_{s}-\displaystyle \frac{1}{2}{{\rm{\Omega }}}_{0}^{2}{R}^{2},\end{eqnarray} \tag{ 2 }$$

with R being the cylindrical radial distance from the rotation axis. In Equation (2), Φ_e represents the external gravitational potential, while Φ_s is the self-gravitational potential satisfying the Poisson equation

$$\begin{eqnarray}&&{{\rm{\nabla }}}^{2}{{\rm{\Phi }}}_{s}=4\pi G\rho .\end{eqnarray} \tag{ 3 }$$

For nuclear rings in barred galaxies, Φ_e is provided mainly by a dark halo as well as a stellar disk and a bulge. The last term in Equation (2) is the centrifugal potential.

We assume ${\rm{\nabla }}{{\rm{\Phi }}}_{e}={{\rm{\Omega }}}_{e}^{2}{\boldsymbol{R}}$, so that the external gravity alone can make a body rotate at constant angular frequency Ω_e. For nuclear rings, this approximation is valid if the rings are geometrically thin. Equation (1) is then integrated to yield

$$\begin{eqnarray}&&{c}_{s}^{2}\mathrm{ln}\rho +{{\rm{\Phi }}}_{s}-\displaystyle \frac{1}{2}{{\rm{\Omega }}}_{s}^{2}{R}^{2}=C,\end{eqnarray} \tag{ 4 }$$

where C is constant, and

$$\begin{eqnarray}&&{{\rm{\Omega }}}_{s}^{2}\equiv {{\rm{\Omega }}}_{0}^{2}-{{\rm{\Omega }}}_{e}^{2}.\end{eqnarray} \tag{ 5 }$$

Note that Ω_s corresponds to the equilibrium angular frequency of an isolated self-gravitating ring in the absence of the external gravity. Our aim is to obtain ρ satisfying Equations (3) and (4) simultaneously.

To obtain equilibrium structure of slowly rotating stars, Ostriker & Mark (1968) introduced an efficient SCF method that solves Equations (3) and (4) alternatively and iteratively. In the SCF method, one first takes a trial distribution for ρ and solves the Poisson equation to find Φ_s, which in turn yields new ρ from Equation (4). Calculations are repeated until the trial and new density distributions agree within a tolerance. Hachisu (1986a) extended the original SCF method to make it suitable for rapidly rotating polytropes. Here, we closely follow Hachisu's method to determine isothermal equilibria.

Following Hachisu (1986a), we let ${\rho }_{c}$ and ${R}_{{\rm{A}}}$ denote the maximum density and the maximum radial extent in the equatorial plane of an equilibrium object, respectively. We introduce the following dimensionless variables: $\widehat{\rho }\equiv \rho /{\rho }_{c},\,$ $\,\widehat{R}\equiv R/{R}_{{\rm{A}}},\,$ $\,{\widehat{{\rm{\Phi }}}}_{s}\equiv {{\rm{\Phi }}}_{s}/({{GR}}_{{\rm{A}}}^{2}{\rho }_{c}),$ and ${\widehat{{\rm{\Omega }}}}_{s}\equiv {{\rm{\Omega }}}_{s}/{(G{\rho }_{c})}^{1/2}.$ Then, Equation (4) reduces to

$$\begin{eqnarray}&&\alpha \mathrm{ln}\widehat{\rho }=\widehat{C}-{\widehat{{\rm{\Phi }}}}_{s}+\displaystyle \frac{1}{2}{\widehat{{\rm{\Omega }}}}_{s}^{2}{\widehat{R}}^{2},\end{eqnarray} \tag{ 6 }$$

where $\widehat{C}$ is a dimensionless constant and

$$\begin{eqnarray}&&\alpha \equiv {c}_{s}^{2}/({{GR}}_{{\rm{A}}}^{2}{\rho }_{c})\end{eqnarray} \tag{ 7 }$$

measures the relative importance of the thermal to gravitational potential energies.

In the SCF method, it is crucial to solve the Poisson equation accurately and efficiently. For spheroid-like configurations, it is customary to employ a multipole expansion technique on spherical polar coordinates. On the other hand, it is more efficient to utilize toroidal coordinates for ring-like configurations, especially for slender rings. In Appendix A, we describe the methods to find Φ for given ρ both in spherical and toroidal coordinates.

In the case of a polytropic equation of state, an object in equilibrium achieves vanishing density at a finite radius, and is thus self-truncated. On the other hand, an isothermal object in steady equilibrium would extend to infinite distance without an external medium. In reality, gas clouds or rings are usually in pressure equilibrium with their surrounding hot rarefied medium that provides a confining external pressure P_ext. Fixing P_ext is equivalent to choosing $\widehat{C}$ in Equation (6), or to placing the boundaries where $\rho ={P}_{{\rm{ext}}}/{c}_{s}^{2}$.

Following Hachisu (1986a), we let (positive) ${R}_{{\rm{B}}}$ denote the radial distance of the boundary along the z-axis for spheroid-like configurations. For ring-like configurations, ${R}_{{\rm{B}}}$ takes the negative of the radial distance to the inner boundary in the equatorial plane. Then, Equation (6) requires

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{s}^{2}=\left\{\begin{array}{l}2[{\widehat{{\rm{\Phi }}}}_{s}(1,\pi /2)-{\widehat{{\rm{\Phi }}}}_{s}({\widehat{R}}_{{\rm{B}}},0)],\\ \quad {\rm{for}}\,{\rm{spheroid \mbox{-} like}}\,{\rm{bodies}},\\ 2[{\widehat{{\rm{\Phi }}}}_{s}(1,\pi /2)-{\widehat{{\rm{\Phi }}}}_{s}(-{\widehat{R}}_{{\rm{B}}},\pi /2)]/(1-{\widehat{R}}_{{\rm{B}}}^{2}),\\ \quad {\rm{for}}\,{\rm{ring \mbox{-} like}}\,{\rm{bodies}},\end{array}\right.\end{eqnarray} \tag{ 8 }$$

where

$$\begin{eqnarray}&&{\widehat{R}}_{{\rm{B}}}\equiv {R}_{{\rm{B}}}/{R}_{{\rm{A}}}.\end{eqnarray} \tag{ 9 }$$

with $0\lt {\widehat{R}}_{{\rm{B}}}\leqslant 1$ for spheroid-like bodies and $-1\lt {\widehat{R}}_{{\rm{B}}}\lt 0$ for ring-like bodies. Since ${\widehat{{\rm{\Omega }}}}_{s}$ depends on ${\widehat{R}}_{{\rm{B}}}$ through Equation (8), an isothermal equilibrium can be completely specified by two parameters: α and ${\widehat{R}}_{{\rm{B}}}$.

Because ${R}_{{\rm{A}}}$ is defined to be the maximum radial extent in the meridional plane, the existence of an equilibrium demands ${\widehat{{\rm{\Phi }}}}_{s}-{\widehat{{\rm{\Omega }}}}_{s}^{2}{\widehat{R}}^{2}/2$ in Equation (6) to be an increasing function of $\widehat{R}$ near $\widehat{R}=1$: the boundary should otherwise retreat to a smaller radius where the thermal pressure is equal to ${P}_{{\rm{ext}}}$. Because the potential minimum occurs inside $\widehat{R}=1$, this requires that the equilibrium should be sufficiently self-gravitating and/or have small enough ${\widehat{{\rm{\Omega }}}}_{s}$. As will be presented in Section 3.2, an isothermal equilibrium turns out to be nonexistent for fairly small $| {R}_{{\rm{B}}}|$ because self-gravity is not strong enough or the angular frequency is too large to form gravitationally bound objects.

In Appendix A.3, we compare the results based on the potential expansions in the spherical and toroidal coordinates for ring-like equilibria, and show that the two methods agree with each other when ${\widehat{R}}_{{\rm{B}}}\gtrsim -0.86$, while the multipole expansion in the spherical coordinates overestimates ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ for smaller ${\widehat{R}}_{{\rm{B}}}$. When we present the results in Section 3, therefore, we employ the multipole expansion with l_max = 10 in the toroidal coordinates for flattened equilibria with ${\widehat{R}}_{{\rm{B}}}\leqslant -0.8$, while adopting the spherical multipole expansion with l_max = 50 for any other equilibria.

As a domain of computation, we consider a meridional cross-section of an equilibrium body, and divide it into N_r × N_a cells. Here, N_r and N_a refer, respectively, to the mesh numbers over $0\leqslant \widehat{r}\leqslant 1.2$ in the r-direction and $0\leqslant \theta \leqslant \pi /2$ in the θ-direction of the spherical coordinates, or over 2.5 ≤ σ ≤ 9.0 in the σ-direction and over 0 ≤ τ ≤ π in the τ-direction of the toroidal coordinates.

Initially, we take $\widehat{\rho }=1$ when $\widehat{r}\leqslant 1$ and $\widehat{\rho }=0$ otherwise for spheroid-like configurations, and $\widehat{\rho }=1$ when ${\widehat{R}}_{{\rm{B}}}\leqslant \widehat{r}\sin \theta \leqslant 1$ and $\widehat{\rho }=0$ otherwise for ring-like configurations in spherical coordinates. When we use the toroidal coordinates, we set the focal length equal to $\widehat{a}\equiv a/{R}_{{\rm{A}}}=(1-{\widehat{R}}_{{\rm{B}}})/2$, and take $\widehat{\rho }=1$ in the regions with ${(x-a)}^{2}\,+\,{z}^{2}\leqslant {a}^{2}$ and $\widehat{\rho }=0$ otherwise. We then calculate ρ_l from Equation (56) or (63), and Φ using Equation (57) or (62) based on Gaussian and Newton–Cotes quadratures (e.g., Press et al. 1988). Next, we calculate ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ from Equation (8) and then update $\widehat{\rho }$ from Equation (6). In each iteration on the toroidal mesh, $\widehat{a}$ is set to move to the location of the maximum density. We repeat the calculations using the updated density until the relative difference in ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ from two successive iterations is smaller than 10⁻⁶. We employ N_r × N_a = 1024 × 512 cells, and it typically takes less than 20 iterations to obtain a converged solution.

Once we find an equilibrium configuration ρ(r, θ) in the spherical coordinates, it is straightforward to calculate its volume V, mass M, angular momentum J, kinetic energy T, and gravitational potential energy W via

$$\begin{eqnarray}&&M=2\pi \int \rho {r}^{2}\sin \theta {drd}\theta ,\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}&&J=2\pi \int \rho {\rm{\Omega }}{r}^{4}{\sin }^{3}\theta {drd}\theta ,\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&T=\pi \int \rho {{\rm{\Omega }}}^{2}{r}^{4}{\sin }^{3}\theta {drd}\theta ,\end{eqnarray} \tag{ 12 }$$

and

$$\begin{eqnarray}&&W=\pi \int \rho {{\rm{\Phi }}}_{s}{r}^{2}\sin \theta {drd}\theta .\end{eqnarray} \tag{ 13 }$$

Note that T = JΩ/2 for rigidly rotating bodies. These quantities will be evaluated and used to draw the equilibrium sequences of isothermal rings in Section 3.

Eriguchi & Sugimoto (1981) found that an incompressible body in axial symmetry takes a form of a Maclaurin spheroid when rotating slowly, and bifurcates into a one-ring sequence when the total angular momentum exceeds a critical value (see also Chandrasekhar 1967; Bardeen 1971; Hachisu 1986a). As a test of our SCF method, we first apply it to construct equilibrium configurations in the incompressible limit, which is attained by taking α ≫ 1. By comparing our results with Eriguchi & Sugimoto (1981) for the Maclaurin spheroids and incompressible rings, we can check the accuracy of our SCF method.

Figure 1 plots the boundaries in the meridional plane of axially symmetric incompressible bodies in (a) spheroid-like and (b) ring-like equilibrium for some selected values of ${\widehat{R}}_{{\rm{B}}}$. For all cases, α = 10⁵ is taken. For $0.158\lesssim {\widehat{R}}_{{\rm{B}}}\leqslant 1$, an equilibrium is exactly a Maclaurin spheroid with an ellipticity $e={(1-{\widehat{R}}_{{\rm{B}}}^{2})}^{1/2}$. The cases with $0\lt {\widehat{R}}_{{\rm{B}}}\lesssim 0.158$ result in concave "hamburgers" that are somewhat more flared at intermediate $\widehat{R}$ (∼0.6–0.8) than at the symmetry axis (e.g., Eriguchi & Sugimoto 1981; Hachisu 1986a). When ${\widehat{R}}_{{\rm{B}}}\lt 0$, on the other hand, equilibrium bodies take a form of rotating rings (or tori), with a larger $| {\widehat{R}}_{{\rm{B}}}|$ corresponding to a more slender ring. Table 1 lists the dimensionless quantities ${\widehat{{\rm{\Omega }}}}_{s}$, $\widehat{M}=M/({R}_{{\rm{A}}}^{3}{\rho }_{c})$, $\widehat{J}=J/({G}^{1/2}{R}_{{\rm{A}}}^{5}{\rho }_{c}^{3/2})$, $\widehat{T}=T/({{GR}}_{{\rm{A}}}^{5}{\rho }_{c}^{2})$, and $\widehat{W}=W/({{GR}}_{{\rm{A}}}^{5}{\rho }_{c}^{2})$ for incompressible bodies.

Zoom In Zoom Out Reset image size

Table 1.  Properties of Axially Symmetric Incompressible Bodies in Steady Equilibrium

${\widehat{R}}_{{\rm{B}}}$	${\widehat{{\rm{\Omega }}}}_{s}^{2}$	$\widehat{M}$	$\widehat{J}$	$\widehat{T}$	$-\widehat{W}$
1.0	0.000E+0	4.189E+0	0.000E+0	0.000E+0	1.054E+1
0.9	3.263E−1	3.767E+0	8.599E−1	2.456E−1	8.811E+0
0.8	6.307E−1	3.349E+0	1.063E+0	4.221E−1	7.220E+0
0.7	9.082E−1	2.927E+0	1.115E+0	5.311E−1	5.729E+0
0.6	1.140E+0	2.510E+0	1.071E+0	5.718E−1	4.385E+0
0.5	1.314E+0	2.093E+0	9.581E−1	5.491E−1	3.179E+0
0.4	1.405E+0	1.671E+0	7.910E−1	4.688E−1	2.121E+0
0.3	1.380E+0	1.254E+0	5.880E−1	3.454E−1	1.253E+0
0.2	1.192E+0	8.389E−1	3.660E−1	1.998E−1	5.902E−1
0.1	1.014E+0	7.885E−1	3.609E−1	1.817E−1	5.012E−1
0.0	1.008E+0	8.504E−1	3.998E−1	2.007E−1	5.726E−1
−0.1	1.019E+0	8.796E−1	4.162E−1	2.100E−1	6.103E−1
−0.2	9.845E−1	9.403E−1	4.528E−1	2.246E−1	6.844E−1
−0.3	8.743E−1	9.397E−1	4.537E−1	2.121E−1	6.703E−1
−0.4	7.189E−1	8.615E−1	4.076E−1	1.728E−1	5.558E−1
−0.5	5.459E−1	7.182E−1	3.235E−1	1.195E−1	3.845E−1
−0.6	3.794E−1	5.366E−1	2.214E−1	6.818E−2	2.162E−1
−0.7	2.321E−1	3.439E−1	1.224E−1	2.947E−2	9.076E−2
−0.8	1.130E−1	1.696E−1	4.661E−2	7.835E−3	2.311E−2
−0.9	3.204E−2	4.654E−2	7.537E−3	6.745E−4	1.911E−3

Note. The number behind E indicates the exponent of the power of 10.

Download table as:  ASCIITypeset image

Figure 2 plots (a) the square of the normalized angular velocity ${\omega }_{s}^{2}\equiv {{\rm{\Omega }}}_{s}^{2}/(4\pi G\langle \rho \rangle )$ and (b) the ratio of the kinetic to gravitational potential energy $t\equiv T/| W|$ as functions of the square of the normalized angular momentum ${j}^{2}\equiv {J}^{2}/(4\pi {{GM}}^{10/3}\langle \rho {\rangle }^{-1/3})$. Here, $\langle \rho \rangle$ denotes the volume-averaged density. The black and blue solid lines are the spheroid-like and ring-like equilibria, respectively, that we obtain by taking α = 10⁵. The black dashed lines plot the theoretical predictions

$$\begin{eqnarray}&&{\omega }_{s}^{2}=\displaystyle \frac{{(1-{e}^{2})}^{1/2}}{2{e}^{2}}\left[(3-2{e}^{2})\displaystyle \frac{{\sin }^{-1}e}{e}-3{(1-{e}^{2})}^{1/2}\right],\end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray}&&t=\displaystyle \frac{3}{2{e}^{2}}-1-\displaystyle \frac{3{(1-{e}^{2})}^{1/2}}{2e{\sin }^{-1}e},\end{eqnarray} \tag{ 15 }$$

and

$$\begin{eqnarray}&&{j}^{2}=\displaystyle \frac{4{\omega }_{s}^{2}}{25}{\left(\displaystyle \frac{3}{4\pi }\right)}^{4/3}{(1-{e}^{2})}^{-2/3},\end{eqnarray} \tag{ 16 }$$

of the Maclaurin spheroid sequence (e.g., Equations (3.234), (3.236), and (3.239) of Lang 1999), which are in excellent agreement with our results at the small-j part. The filled stars at j² = 2.233 × 10⁻² (occurring at ${\widehat{R}}_{{\rm{B}}}=0.158$) mark the bifurcation point where the Maclaurin sequence branches out into the concave hamburgers and then into the one-ring sequence after the filled circles at j² = 2.183 × 10⁻² (or ${\widehat{R}}_{{\rm{B}}}=0$). The insets zoom in the regions with 0.020 ≤j² ≤ 0.024 around the bifurcation point for clarity. For comparison, we also plot the ${\omega }_{s}^{2}$–j² and $T/| W|$–j² relationships adopted from Table 1 of Eriguchi & Sugimoto (1981) as the red dotted lines, which are slightly (∼1%–2%) different from our results (see also Hachisu 1986a). These discrepancies are presumably due to the insufficient resolution used by these authors in solving the Poisson equation.³

Zoom In Zoom Out Reset image size

We now present density distributions of isothermal objects in axisymmetric equilibrium. We first visit non-rotating isothermal spheres truncated by an external pressure. We then explore how rotation changes equilibrium structures.

3.2.1. Bonnor–Ebert Spheres

Consider non-rotating, self-gravitating isothermal spheres with Ω_s = 0, namely Bonner–Ebert spheres, in hydrostatic equilibrium. Equations (3) and (4) are then combined to yield

$$\begin{eqnarray}&&\displaystyle \frac{1}{{\xi }^{2}}\displaystyle \frac{d}{d\xi }\left({\xi }^{2}\displaystyle \frac{d\psi }{d\xi }\right)=\exp (-\psi ),\end{eqnarray} \tag{ 17 }$$

where $\psi =({{\rm{\Phi }}}_{s}-C)/{c}_{s}^{2}$ and $\xi ={(4\pi G{\rho }_{c}/{c}_{s}^{2})}^{1/2}\,r\ =(4\pi /{\alpha )}^{1/2}\,\widehat{r}$ are the dimensionless potential and radius, respectively. Equation (17) can be solved subject to the proper boundary conditions, ψ = dψ/dξ = 0 at ξ = 0, to give the density distribution $\rho ={\rho }_{c}\exp (-\psi )$, which shall be compared with the results of our SCF method.

As a second test, we apply our SCF method to obtain density distributions of isothermal spheres by setting ${\widehat{R}}_{{\rm{B}}}=1$. Figure 3(a) plots the resulting density profiles for α = 1, 0.1, and 0.01 as functions of $\widehat{r}$: all the cases are truncated at $\widehat{r}=1$. Note that ρ varies more steeply for smaller α in order to compensate for a smaller sound speed in balancing self-gravity. Note also that when drawn against ξ, as shown in Figure 3(b), all the curves lie well (within 1%) on the inner parts of the solution of Equation (17) plotted as a dashed line, confirming the accuracy of our SCF method. A smaller α corresponds to a larger truncation radius in ξ.

Zoom In Zoom Out Reset image size

3.2.2. Rigidly Rotating Isothermal Equilibria

By taking ${\widehat{R}}_{{\rm{B}}}$ less than unity, we obtain equilibrium density of isothermal objects in rigid rotation. Figure 4 presents the resulting density distributions on the meridional plane for such equilibria with differing ${\widehat{R}}_{{\rm{B}}}$. The left, middle, and right columns are for α = 1, 0.1, and 0.01, respectively. Figure 5 plots the exemplary density profiles along the $\widehat{R}$-axis (solid lines) and the $\widehat{z}$-axis (dotted lines) for the spheroid-like configurations with ${\widehat{R}}_{{\rm{B}}}=0.6$ and the ring-like configurations with ${\widehat{R}}_{{\rm{B}}}=-0.4$. Clearly, an equilibrium is more centrally concentrated for smaller α. The vertical extent of an equilibrium body is smaller than the horizontal extent for both spheroid-like and ring-like objects. Unlike Bonner–Ebert spheres whose density distributions are independent of α when expressed in terms of ξ, we find that the density profiles of rotating isothermal objects with different α along the $\widehat{z}$- or $\widehat{R}$-axis cannot be expressed by a single function of ${(4\pi /\alpha )}^{1/2}\widehat{z}$ or ${(4\pi /\alpha )}^{1/2}\widehat{R}$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 6 plots the variations of the square of the angular velocity ${\widehat{{\rm{\Omega }}}}_{s}^{2}$, the total mass $\widehat{M}$, the averaged density $\langle \widehat{\rho }\rangle$, the kinetic energy $\widehat{T}$, and the gravitational potential energy $\widehat{W}$ as functions of ${\widehat{R}}_{{\rm{B}}}$ for isothermal equilibria with α = 1, 0.1, and 0.01, in comparison with the incompressible cases. Both ${\widehat{{\rm{\Omega }}}}_{s}$ and $\widehat{M}$ increase as $| {\widehat{R}}_{{\rm{B}}}|$ decreases from unity. For spheroid-like configurations, this is because an equilibrium body becomes more flattened and occupies a smaller volume as it rotates faster. On the other hand, ring-like configurations attain a larger volume and mass with decreasing $| {\widehat{R}}_{{\rm{B}}}|$, and thus requires larger ${\widehat{{\rm{\Omega }}}}_{s}$ to balance self-gravity. Note, however, that $\widehat{M}$ is not a monotonically decreasing function of ${\widehat{R}}_{{\rm{B}}}$ for ring-like configurations due to complicated dependence of their volume on ${\widehat{R}}_{{\rm{B}}}$ (e.g., Figure 1(b)). For α = 0.1, for example, $\widehat{M}$ increases as ${\widehat{R}}_{{\rm{B}}}$ moves away from −1, is maximized at ${\widehat{R}}_{{\rm{B}}}=-0.30$, and starts to decreases afterwards. Obviously, $\langle \widehat{\rho }\rangle =1$ for incompressible configurations due to uniform density. Overall, $\langle \widehat{\rho }\rangle$ increases with decreasing ${\widehat{R}}_{{\rm{B}}}$ and tends to unity as ${\widehat{R}}_{{\rm{B}}}$ approaches −1. The dependency of $\widehat{T}$ and $\widehat{W}$ on ${\widehat{R}}_{{\rm{B}}}$ is closely related to that of ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ and $\widehat{M}$, respectively. For $\alpha \lesssim 0.1$, $\widehat{M}$ and $\widehat{W}$ are insensitive to ${\widehat{R}}_{{\rm{B}}}\gtrsim 0.6$ since the density in the outer parts is very small. All the quantities are smaller with smaller α due to a stronger density concentration: these values are also listed in Table 2 for some selected values of ${\widehat{R}}_{{\rm{B}}}$.

Zoom In Zoom Out Reset image size

Table 2.  Properties of Axially Symmetric Isothermal Bodies in Steady Equilibrium

${\widehat{R}}_{{\rm{B}}}$	${\widehat{{\rm{\Omega }}}}_{s}^{2}$	$\widehat{M}$	$\langle \widehat{\rho }\rangle$	$\widehat{J}$	$\widehat{T}$	$-\widehat{W}$
			α = 1
1.0	0.000E+0	1.806E+0	4.310E−1	0.000E+0	0.000E+0	2.140E+0
0.9	1.861E−1	1.735E+0	4.611E−1	2.566E−1	5.536E−2	2.028E+0
0.8	3.805E−1	1.657E+0	4.964E−1	3.551E−1	1.095E−1	1.900E+0
0.7	5.826E−1	1.564E+0	5.392E−1	4.209E−1	1.606E−1	1.746E+0
0.6	7.810E−1	1.455E+0	5.901E−1	4.599E−1	2.032E−1	1.561E+0
0.5	9.653E−1	1.322E+0	6.512E−1	4.709E−1	2.313E−1	1.335E+0
0.4	1.111E+0	1.152E+0	7.243E−1	4.453E−1	2.347E−1	1.054E+0
0.3	1.171E+0	9.247E−1	8.065E−1	3.672E−1	1.987E−1	7.133E−1
0.1	9.806E−1	7.472E−1	9.465E−1	3.316E−1	1.642E−1	4.531E−1
0.0	9.681E−1	8.052E−1	9.350E−1	3.662E−1	1.802E−1	5.162E−1
−0.1	9.728E−1	8.281E−1	9.305E−1	3.781E−1	1.865E−1	5.440E−1
−0.2	9.272E−1	8.744E−1	9.202E−1	4.044E−1	1.947E−1	5.950E−1
−0.3	8.176E−1	8.696E−1	9.173E−1	4.030E−1	1.822E−1	5.767E−1
−0.4	6.726E−1	7.999E−1	9.232E−1	3.644E−1	1.494E−1	4.811E−1
−0.5	5.156E−1	6.745E−1	9.360E−1	2.946E−1	1.058E−1	3.402E−1
−0.6	3.633E−1	5.119E−1	9.523E−1	2.064E−1	6.222E−2	1.971E−1
−0.7	2.257E−1	3.338E−1	9.700E−1	1.171E−1	2.782E−2	8.563E−2
−0.8	1.115E−1	1.671E−1	9.854E−1	4.562E−2	7.617E−3	2.246E−2
−0.9	3.193E−2	4.636E−2	9.961E−1	7.494E−3	6.695E−4	1.897E−3
			α = 0.1
1.0	0.000E+0	2.490E−1	5.941E−2	0.000E+0	0.000E+0	5.219E−2
0.9	3.749E−2	2.480E−1	6.607E−2	1.141E−2	1.104E−3	5.258E−2
0.8	8.080E−2	2.467E−1	7.483E−2	1.696E−2	2.410E−3	5.292E−2
0.7	1.319E−1	2.440E−1	8.718E−2	2.177E−2	3.954E−3	5.295E−2
0.6	1.912E−1	2.385E−1	1.053E−1	2.573E−2	5.625E−3	5.212E−2
−0.2	5.916E−1	4.694E−1	4.925E−1	1.546E−1	5.944E−2	1.831E−1
−0.3	5.060E−1	4.825E−1	4.910E−1	1.645E−1	5.852E−2	1.861E−1
−0.4	4.242E−1	4.684E−1	5.224E−1	1.632E−1	5.316E−2	1.709E−1
−0.5	3.439E−1	4.281E−1	5.802E−1	1.499E−1	4.394E−2	1.407E−1
−0.6	2.632E−1	3.596E−1	6.605E−1	1.224E−1	3.141E−2	9.898E−2
−0.7	1.809E−1	2.632E−1	7.617E−1	8.246E−2	1.754E−2	5.379E−2
−0.8	9.929E−2	1.478E−1	8.706E−1	3.804E−2	5.993E−3	1.764E−2
−0.9	3.093E−2	4.481E−2	9.627E−1	7.130E−3	6.270E−4	1.774E−3
			α = 0.01
1.0	0.000E+0	2.018E−2	4.831E−3	0.000E+0	0.000E+0	4.390E−4
0.9	3.318E−3	2.006E−2	5.349E−3	2.255E−4	6.494E−6	4.400E−4
0.8	7.235E−3	1.984E−2	6.058E−3	3.313E−4	1.409E−5	4.391E−4
0.7	1.199E−2	1.946E−2	7.109E−3	4.155E−4	2.275E−5	4.344E−4
0.6	1.779E−2	1.862E−2	8.857E−3	4.593E−4	3.064E−5	4.188E−4
−0.4	1.084E−1	7.682E−2	9.693E−2	1.105E−2	1.819E−3	5.609E−3
−0.5	8.838E−2	8.399E−2	1.134E−1	1.363E−2	2.027E−3	6.225E−3
−0.6	7.534E−2	8.618E−2	1.551E−1	1.509E−2	2.071E−3	6.293E−3
−0.7	6.356E−2	8.293E−2	2.362E−1	1.517E−2	1.913E−3	5.714E−3
−0.8	4.852E−2	6.856E−2	4.028E−1	1.229E−2	1.354E−3	3.937E−3
−0.9	2.391E−2	3.357E−2	7.206E−1	4.664E−3	3.583E−4	1.005E−3

Note. The number behind E indicates the exponent of the power of 10.

Download table as:  ASCIITypeset image

The dependencies of ${\widehat{{\rm{\Omega }}}}_{s}$ and $\widehat{M}$ on ${\widehat{R}}_{{\rm{B}}}$ make an equilibrium sequence with fixed α ( ≤ 1) cease to exist for an intermediate range of ${\widehat{R}}_{{\rm{B}}}$: $0.13\lt {\widehat{R}}_{{\rm{B}}}\lt 0.27$ for α = 1, $-0.14\lt {\widehat{R}}_{{\rm{B}}}\lt 0.51$ for α = 0.1, and $-0.40\lt {\widehat{R}}_{{\rm{B}}}\lt 0.59$ for α = 0.01, with the corresponding boundaries marked by filled circles in Figure 6. This is unlike the incompressible bodies for which any value of $0\leqslant | {\widehat{R}}_{{\rm{B}}}| \leqslant 1$ readily yields an equilibrium. As mentioned in Section 2.2, the presence of a steady equilibrium requires large enough $| {\widehat{{\rm{\Phi }}}}_{s}|$ and/or small enough ${\widehat{{\rm{\Omega }}}}_{s}$ for ρ to be a decreasing function of $\widehat{R}$ near $\widehat{R}=1$ (see Equation (6)). The absence of an isothermal equilibrium for intermediate ${\widehat{R}}_{{\rm{B}}}$ results from the fact that the self-gravitational potential is too weak to overcome the centrifugal potential in establishing gravitationally bound objects.

Figure 7 plots the square of the normalized angular velocity ${\omega }_{s}^{2}$ as well as the energy ratio $t=T/| W|$ as functions of the square of the normalized angular momentum j² for isothermal equilibria with α = 1, 0.1, and 0.01. The incompressible cases with α = 10⁵ are compared as dashed lines. Comparison of Figure 7 with Figures 10 and 11 of Hachisu (1986a) reveals that the ${\omega }_{s}^{2}$–j² and t–j² relationships of isothermal objects with α ∼ 0.01–1 are very close to those of polytropes with index n ∼ 0.1–1.5. Hachisu (1986a) found that spheroid-like polytropic equilibria can be possible only in the region ${\omega }_{s}^{2}+5{j}^{2}\lt 0.185$ when j² < 0.02. Our results suggest that spheroid-like isothermal equilibria can exist in the shaded region in Figure 7(a), which is bounded by

$$\begin{eqnarray}&&{\omega }_{s}^{2}+4.12{j}^{2}=0.172,\end{eqnarray} \tag{ 18 }$$

and the incompressible ${\omega }_{s}^{2}$–j² relation.

Zoom In Zoom Out Reset image size

Here, we focus on the properties of slender isothermal rings whose minor axis is much shorter than the major axis. Density distributions of such rings can be obtained by taking ${\widehat{R}}_{{\rm{B}}}$ close to −1 in our SCF method. Using a perturbation analysis, Ostriker (1964b) derived approximate expressions for the density and angular frequency of both polytropic and isothermal rings in steady equilibrium. Our objective in this subsection is to compare the results of our SCF method with Ostriker (1964b).

For a ring-like configuration with ${\widehat{R}}_{{\rm{B}}}\lt 0$, we define its major axis R₀ and minor axis η₀ as

$$\begin{eqnarray}&&{R}_{0}={ \mathcal V }/(2\pi { \mathcal A }),{\rm{}}\,{\rm{and}}\,{\rm{}}{\eta }_{0}={({ \mathcal A }/\pi )}^{1/2},\end{eqnarray} \tag{ 19 }$$

where ${ \mathcal V }$ and ${ \mathcal A }$ refers to the total volume and the meridional cross-section, respectively, occupied by the equilibrium body. We plot in Figure 8(a) as solid lines ${\widehat{R}}_{0}={R}_{0}/{R}_{{\rm{A}}}$ and ${\widehat{\eta }}_{0}={\eta }_{0}/{R}_{{\rm{A}}}$ together with ${\widehat{R}}_{0}+{\widehat{\eta }}_{0}$ resulting from the SCF method as functions of ${\widehat{R}}_{{\rm{B}}}$ for α = 10⁵, 1, 0.1, and 0.01. Note that ${\widehat{\eta }}_{0}\simeq (1+{\widehat{R}}_{{\rm{B}}})/2$ and ${\widehat{R}}_{0}\simeq (1-{\widehat{R}}_{{\rm{B}}})/2$, resulting in ${\widehat{R}}_{0}+{\widehat{\eta }}_{0}\simeq 1$, as expected, for ${\widehat{R}}_{{\rm{B}}}\lesssim -0.6$. This indicates that R₀ and η₀ defined in Equation (19) describe the real major and minor axes of slender rings reasonably well.

Zoom In Zoom Out Reset image size

Using a perturbation analysis, Ostriker (1964b) showed that an incompressible, slender ring in the absence of external gravity should obey

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{s}^{2}=\displaystyle \frac{M}{2\pi {R}_{0}^{3}{\rho }_{c}}\left[\mathrm{ln}\left(\displaystyle \frac{8{R}_{0}}{{\eta }_{0}}\right)-\displaystyle \frac{5}{4}\right].\end{eqnarray} \tag{ 20 }$$

For an isothermal ring of infinite extent (without pressure truncation), he also derived

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{s}^{2}=\displaystyle \frac{{M}_{\infty }}{2\pi {R}_{0}^{3}{\rho }_{c}}\left[\mathrm{ln}\left(\displaystyle \frac{8{R}_{0}}{{\eta }_{1/2}}\right)-2\right],\end{eqnarray} \tag{ 21 }$$

which is valid when ${\eta }_{1/2}/{R}_{0}\ll 1$. Here, ${M}_{\infty }$ denotes the total mass and ${\eta }_{1/2}={({M}_{\infty }/2{\pi }^{2}{R}_{0}{\rho }_{c})}^{1/2}$ is the half-mass radius. Figure 8(b) plots Equation (20) as dotted lines, which are in good agreement with the results of our SCF method, shown as solid lines, for α ≳ 1 at ${\widehat{R}}_{{\rm{B}}}\lesssim -0.4$ and even for α as small as 0.01 at ${\widehat{R}}_{{\rm{B}}}\lesssim -0.8$. Figure 8(b) also plots Equation (21) as dashed lines after taking $M={M}_{\infty }$, which matches well our SCF results only for α = 0.01. The discrepancy between Equation (21) and our results is in fact expected since isothermal rings considered in the present paper are truncated by external pressure. Since rings with α = 0.01 are highly concentrated, however, Equation (21) can still be a good approximation for truncated slender rings.

As the bottom panels of Figure 4 show, the density distributions of slender rings at the meridional plane appear almost circularly symmetric with respect to the point (R, z) = (R₀, 0). Ostriker (1964b) showed that the meridional density distribution is, to the zeroth order in η_1/2, given by

$$\begin{eqnarray}&&{\rho }_{{\rm{sr}}}=\displaystyle \frac{{\rho }_{c}}{{[1+{(\eta /H)}^{2}/8]}^{2}},\end{eqnarray} \tag{ 22 }$$

where η denotes the distance from the density maximum and

$$\begin{eqnarray}&&H\equiv {c}_{s}/{(4\pi G{\rho }_{c})}^{1/2}\end{eqnarray} \tag{ 23 }$$

is the characteristic ring thickness. Note that Equation (22) is also the solution for non-rotating isothermal cylinders of infinite length along its symmetry axis (e.g., Ostriker 1964a; Nagasawa 1987; Inutsuka & Miyama 1992). Figure 9(a) compares Equation (22) (black) with the density profiles from the SCF method for slender rings with ${\widehat{R}}_{{\rm{B}}}=-0.8$ (or with ${\widehat{\eta }}_{0}=0.1$) along the radial (blue) and vertical (red) directions from the density maximum. The cases with α = 1, 0.1, and 0.01 are shown as dashed, dotted, and solid lines, respectively. The relative errors, ${\rho }_{{\rm{sr}}}/\rho -1$, given in Figure 9(b) are only a few percents in the most of the dense regions, demonstrating that Equation (22) is a good approximation to the true density distributions of slender isothermal rings.

Zoom In Zoom Out Reset image size

The basic hydrodynamic equations governing evolution of isothermal gas read

$$\begin{eqnarray}&&\displaystyle \frac{\partial \rho }{\partial t}+{\rm{\nabla }}\cdot (\rho {\boldsymbol{v}})=0,\end{eqnarray} \tag{ 25 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\boldsymbol{v}}}{\partial t}+{\boldsymbol{v}}\cdot {\rm{\nabla }}{\boldsymbol{v}}=-\displaystyle \frac{{c}_{s}^{2}}{\rho }{\rm{\nabla }}\rho -{{\rm{\Omega }}}_{e}^{2}{\boldsymbol{R}}-{\rm{\nabla }}{{\rm{\Phi }}}_{s},\end{eqnarray} \tag{ 26 }$$

together with Equation (3).

To analyze a linear stability of a ring, it is convenient to introduce the new curvilinear coordinates $(\eta ,\lambda ,\phi )$, as depicted in Figure 10. The new coordinates are related to the Cartesian coordinate system (x, y, z) through

$$\begin{eqnarray}&&\left(\begin{array}{c}x\\ y\\ z\end{array}\right)=\left(\begin{array}{c}({R}_{0}+\eta \cos \lambda )\cos \phi \\ ({R}_{0}+\eta \cos \lambda )\sin \phi \\ \eta \sin \lambda \end{array}\right).\end{eqnarray} \tag{ 27 }$$

The coordinate η is the distance from a reference circle of radius R₀ located in the horizontal plane, while λ is the polar angle measured from the horizontal plane; ϕ is the usual cylindrical azimuthal angle. The new coordinates are orthogonal, and reduce to the usual spherical polar coordinates $(\eta ,\pi /2-\lambda ,\phi )$ in the limit of ${R}_{0}/\eta \ll 1$.

Zoom In Zoom Out Reset image size

Appendix B derives gas dynamical equations for slender rings in the new curvilinear coordinates. Fluid variables are in general three-dimensional, and finding dispersion relations of three-dimensional perturbations applied to a ring is a daunting task. However, gravitationally unstable modes we seek in the present work are dominated by the velocity components in the η- and ϕ-directions, without much involving gas motions in the λ-direction. For simplicity, therefore, we take v_λ = 0 and assume that all physical quantities are independent of λ. We also take $\cos \lambda =1$, which allows us to fully capture the rotational effect of the ring material around the symmetry axis. We note however that our method, by construction, is unable to handle shear that may exist in the rotation of the ring. In Section 4.5, we will use direct numerical simulations to show that shear does not significantly affect gravitational instabilities of slender rings.

Under these circumstances, Equations (73)–(77) can be simplified to

$$\begin{eqnarray}&&\displaystyle \frac{\partial \rho }{\partial t}+\displaystyle \frac{1}{\eta }\displaystyle \frac{\partial (\eta \rho {v}_{\eta })}{\partial \eta }+\displaystyle \frac{1}{{R}_{0}}\displaystyle \frac{\partial (\rho {v}_{\phi })}{\partial \phi }=0,\end{eqnarray} \tag{ 28 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {v}_{\eta }}{\partial t}+\left({v}_{\eta }\displaystyle \frac{\partial }{\partial \eta }+\displaystyle \frac{{v}_{\phi }}{{R}_{0}}\displaystyle \frac{\partial }{\partial \phi }\right){v}_{\eta }-\displaystyle \frac{{v}_{\phi }^{2}}{{R}_{0}}\\ &&\quad =-\,\displaystyle \frac{{c}_{s}^{2}}{\rho }\displaystyle \frac{\partial \rho }{\partial \eta }-{{\rm{\Omega }}}_{e}^{2}{R}_{0}-\displaystyle \frac{\partial {{\rm{\Phi }}}_{s}}{\partial \eta }\end{eqnarray} \tag{ 29 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {v}_{\phi }}{\partial t}+\left({v}_{\eta }\displaystyle \frac{\partial }{\partial \eta }+\displaystyle \frac{{v}_{\phi }}{{R}_{0}}\displaystyle \frac{\partial }{\partial \phi }\right){v}_{\phi }+\displaystyle \frac{{v}_{\eta }{v}_{\phi }}{{R}_{0}}=-\displaystyle \frac{{c}_{s}^{2}}{{R}_{0}\rho }\displaystyle \frac{\partial \rho }{\partial \phi }-\displaystyle \frac{1}{{R}_{0}}\displaystyle \frac{\partial {{\rm{\Phi }}}_{s}}{\partial \phi }\end{eqnarray} \tag{ 30 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{{\rm{\Phi }}}_{s}}{\partial {\eta }^{2}}+\displaystyle \frac{1}{\eta }\displaystyle \frac{\partial {{\rm{\Phi }}}_{s}}{\partial \eta }+\displaystyle \frac{1}{{R}_{0}^{2}}\displaystyle \frac{{\partial }^{2}{{\rm{\Phi }}}_{s}}{\partial {\phi }^{2}}=4\pi G\rho .\end{eqnarray} \tag{ 31 }$$

The initial density distribution (i.e., Equation (24)) satisfies Equations (29) and (31) as long as ${{\rm{\Omega }}}_{s}^{2}\ll {{\rm{\Omega }}}_{0}^{2}$, a condition easily met for nuclear rings.

We now consider small-amplitude perturbations applied to the initial equilibrium. Denoting the background quantities and the perturbations using the subscripts "0" and "1," respectively, we linearize Equations (28)–(31) to obtain

$$\begin{eqnarray}&&\left(\displaystyle \frac{\partial }{\partial t}+{{\rm{\Omega }}}_{0}\displaystyle \frac{\partial }{\partial \phi }\right){\rho }_{1}+\displaystyle \frac{1}{\eta }\displaystyle \frac{\partial (\eta {\rho }_{0}{v}_{\eta 1})}{\partial \eta }+\displaystyle \frac{{\rho }_{0}}{{R}_{0}}\displaystyle \frac{\partial {v}_{\phi 1}}{\partial \phi }=0,\end{eqnarray} \tag{ 32 }$$

$$\begin{eqnarray}&&\left(\displaystyle \frac{\partial }{\partial t}+{{\rm{\Omega }}}_{0}\displaystyle \frac{\partial }{\partial \phi }\right){v}_{\eta 1}-2{{\rm{\Omega }}}_{0}{v}_{\phi 1}=-\displaystyle \frac{\partial {\chi }_{1}}{\partial \eta },\end{eqnarray} \tag{ 33 }$$

$$\begin{eqnarray}&&\left(\displaystyle \frac{\partial }{\partial t}+{{\rm{\Omega }}}_{0}\displaystyle \frac{\partial }{\partial \phi }\right){v}_{\phi 1}+2{{\rm{\Omega }}}_{0}{v}_{\eta 1}=-\displaystyle \frac{1}{{R}_{0}}\displaystyle \frac{\partial {\chi }_{1}}{\partial \phi },\end{eqnarray} \tag{ 34 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{{\rm{\Phi }}}_{s1}}{\partial {\eta }^{2}}+\displaystyle \frac{1}{\eta }\displaystyle \frac{\partial {{\rm{\Phi }}}_{s1}}{\partial \eta }+\displaystyle \frac{1}{{R}_{0}^{2}}\displaystyle \frac{{\partial }^{2}{{\rm{\Phi }}}_{s1}}{\partial {\phi }^{2}}=4\pi G{\rho }_{1},\end{eqnarray} \tag{ 35 }$$

where

$$\begin{eqnarray}&&{\chi }_{1}\equiv {c}_{s}^{2}\displaystyle \frac{{\rho }_{1}}{{\rho }_{0}}+{{\rm{\Phi }}}_{s1}.\end{eqnarray} \tag{ 36 }$$

We assume that any perturbation, q₁, varies in space and time as

$$\begin{eqnarray}&&{q}_{1}(\eta ,\phi ,t)={q}_{1}(\eta )\exp ({im}\phi -i\omega t),\end{eqnarray} \tag{ 37 }$$

with ω and m being the frequency and azimuthal mode number of the perturbations, respectively. Substituting Equation (37) into Equations (32)–(35), one obtains

$$\begin{eqnarray}&&\displaystyle \frac{{{dv}}_{\eta 1}}{d\eta }+\displaystyle \frac{d\mathrm{ln}(\eta {\rho }_{0})}{d\eta }{v}_{\eta 1}=i{\omega }_{D}\displaystyle \frac{{\rho }_{1}}{{\rho }_{0}}-i\displaystyle \frac{m}{{R}_{0}}{v}_{\phi 1},\end{eqnarray} \tag{ 38 }$$

$$\begin{eqnarray}&&{v}_{\eta 1}=\displaystyle \frac{i}{{\omega }_{D}^{2}-4{{\rm{\Omega }}}_{0}^{2}}\left(\displaystyle \frac{2m{{\rm{\Omega }}}_{0}}{{R}_{0}}{\chi }_{1}-{\omega }_{D}\displaystyle \frac{d{\chi }_{1}}{d\eta }\right),\end{eqnarray} \tag{ 39 }$$

$$\begin{eqnarray}&&{v}_{\phi 1}=\displaystyle \frac{1}{{\omega }_{D}^{2}-4{{\rm{\Omega }}}_{0}^{2}}\left(\displaystyle \frac{m{\omega }_{D}}{{R}_{0}}{\chi }_{1}-2{{\rm{\Omega }}}_{0}\displaystyle \frac{d{\chi }_{1}}{d\eta }\right),\end{eqnarray} \tag{ 40 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}{{\rm{\Phi }}}_{s1}}{d{\eta }^{2}}+\displaystyle \frac{1}{\eta }\displaystyle \frac{d{{\rm{\Phi }}}_{s1}}{d\eta }-\displaystyle \frac{{m}^{2}}{{R}_{0}^{2}}{{\rm{\Phi }}}_{s1}=4\pi G{\rho }_{1},\end{eqnarray} \tag{ 41 }$$

where ${\omega }_{D}\equiv \omega -m{{\rm{\Omega }}}_{0}$ is the Doppler-shifted frequency. Eliminating v_η1 and v_ϕ1 in favor of χ₁ from Equations (38)–(40), one finally obtains

$$\begin{eqnarray}&&\displaystyle \frac{{d}^{2}{\chi }_{1}}{d{\eta }^{2}}+\displaystyle \frac{d\mathrm{ln}(\eta {\rho }_{0})}{d\eta }\displaystyle \frac{d{\chi }_{1}}{d\eta }-\left[\displaystyle \frac{2{{\rm{\Omega }}}_{0}}{{\omega }_{D}}\displaystyle \frac{m}{{R}_{0}}\displaystyle \frac{d\mathrm{ln}(\eta {\rho }_{0})}{d\eta }+\displaystyle \frac{{m}^{2}}{{R}_{0}^{2}}\right]{\chi }_{1}\\ &&\quad =-\,({\omega }_{D}^{2}-4{{\rm{\Omega }}}_{0}^{2})\displaystyle \frac{{\rho }_{1}}{{\rho }_{0}}.\end{eqnarray} \tag{ 42 }$$

Equations (41) and (42) constitute a set of coupled equations that can be solved simultaneously for eigenvalue ω subject to proper boundary conditions.

Since Equations (41) and (42) are second-order differential equations, we need to have five constraints in order to determine ω unambiguously. Since this is a linear problem, we are free to choose the amplitude of one variable arbitrarily. Two conditions come from the inner boundary by the requirements

$$\begin{eqnarray}&&{\left.\displaystyle \frac{d{\chi }_{1}}{d\eta }\right|}_{\eta =0}={\left.\displaystyle \frac{d{{\rm{\Phi }}}_{s1}}{d\eta }\right|}_{\eta =0}=0,\end{eqnarray} \tag{ 43 }$$

for regular solutions at η = 0.

The remaining two conditions can be obtained from the outer boundary. Perturbations given in Equation (37) also disturb the ring surface to

$$\begin{eqnarray}&&\eta ={\eta }_{0}+{\eta }_{1}\exp ({im}\phi -i\omega t),\end{eqnarray} \tag{ 44 }$$

with amplitude η₁, implying that

$$\begin{eqnarray}&&{v}_{\eta 1}({\eta }_{0})=-i\omega {\eta }_{1}.\end{eqnarray} \tag{ 45 }$$

The pressure equilibrium at the disturbed surface, ${c}_{s}^{2}{\rho }_{0}({\eta }_{0})+{P}_{1}({\eta }_{0})={P}_{{\rm{ext}}}$, with P₁ representing the perturbed pressure (e.g., Nagasawa 1987), requires

$$\begin{eqnarray}&&{\left.\displaystyle \frac{d{\rho }_{0}}{d\eta }\right|}_{{\eta }_{0}}{\eta }_{1}+{\rho }_{1}({\eta }_{0})=0,\end{eqnarray} \tag{ 46 }$$

which is a third boundary condition.

To derive a fourth boundary condition, we follow Goldreich & Lynden-Bell (1965) to assume that the perturbed mass near the outer boundary is restricted to a thin annulus such that ${\rho }_{1}={\rho }_{0}{\eta }_{1}\delta (\eta -{\eta }_{0})$ (see also, Nagasawa 1987; Kim et al. 2012a). Integrating Equation (41) from η = η₀ to η = η₀ + η₁, one obtains

$$\begin{eqnarray}&&\displaystyle \frac{d{{\rm{\Phi }}}_{s1}^{+}}{d\eta }-\displaystyle \frac{d{{\rm{\Phi }}}_{s1}^{-}}{d\eta }=4\pi G{\rho }_{0}{\eta }_{1}\ \ \ \mathrm{at}\ \ \ \eta ={\eta }_{0},\end{eqnarray} \tag{ 47 }$$

where the superscripts "+" and "−" indicate the potentials evaluated just outside and inside the ring surface, respectively. Assuming that the region outside the ring is filled with an extremely hot, tenuous gas, ${{\rm{\Phi }}}_{s1}^{+}$ should satisfy Equation (35) with ρ₁ = 0. The regular solution at infinity can be expressed as

$$\begin{eqnarray}&&{{\rm{\Phi }}}_{s1}^{+}={{AK}}_{0}\left(\displaystyle \frac{m}{{R}_{0}}\eta \right),\ \ \ \mathrm{for}\ \ \ \eta /{\eta }_{0}\geqslant 1,\end{eqnarray} \tag{ 48 }$$

where A is a constant and K_n is the second-kind modified Bessel function of order n. The condition that the gravitational potential should be continuous across the surface gives $A\ ={{\rm{\Phi }}}_{s1}({\eta }_{0})/{K}_{0}(m{\eta }_{0}/{R}_{0})$. Plugging Equation (48) into Equation (47) gives

$$\begin{eqnarray}&&\displaystyle \frac{d{{\rm{\Phi }}}_{s1}}{d\eta }+4\pi G{\rho }_{0}{\eta }_{1}=-\displaystyle \frac{m}{{R}_{0}}\displaystyle \frac{{K}_{1}(m{\eta }_{0}/{R}_{0})}{{K}_{0}(m{\eta }_{0}/{R}_{0})}{{\rm{\Phi }}}_{s1}\ \ \ \mathrm{at}\ \ \ \eta ={\eta }_{0}.\end{eqnarray} \tag{ 49 }$$

Equations (43), (46), and (49) are our complete set of the boundary conditions.

By writing Equations (41) and (42) into a dimensionless form, one can see that the problem of finding the dimensionless eigenvalue $\widehat{\omega }\equiv \omega /{(G{\rho }_{c})}^{1/2}$ is completely specified by four dimensionless parameters, η₀$/$R₀, ${R}_{0}/H$, ${\widehat{{\rm{\Omega }}}}_{0}\equiv {{\rm{\Omega }}}_{0}/{(G{\rho }_{c})}^{1/2}$, and m. Nuclear rings typically have ${\eta }_{0}\sim 0.1\,{\rm{kpc}}$, ${R}_{0}\sim 1\,{\rm{kpc}}$, $M\sim 4\times {10}^{8}\,{M}_{\odot }$, and Ω₀ ∼ 100–200 ${\text{km s}}^{-1}$ kpc⁻¹. The corresponding mean density is $\langle \widehat{\rho }\rangle =M/(2{\pi }^{2}{\eta }_{0}^{2}{R}_{0})\sim 1\ \times {10}^{-22}\,{\text{g cm}}^{-3}$, and the characteristic ring thickness is

$$\begin{eqnarray}&&H=35{\rm{pc}}\ \left(\displaystyle \frac{{c}_{s}}{10\,{\text{km s}}^{-1}}\right){\left(\displaystyle \frac{{\rho }_{c}}{1\times {10}^{-22}{{\rm{g}}{\rm{cm}}}^{-3}}\right)}^{-1/2},\end{eqnarray} \tag{ 50 }$$

We therefore choose η₀/R₀ = 0.1 and R₀/H = 30 as our standard set of parameters, but also vary R₀/H and ${\widehat{{\rm{\Omega }}}}_{0}$ to explore the effects of ring thickness and rotation on the ring stability. For slender rings with η₀/R₀ = 0.1, R₀/H is related to α through $\alpha \simeq 10.4{({R}_{0}/H)}^{-2}$.

As a normalization condition, we take Re(χ₁) = Im(χ₁) = 1 at the outer boundary. For given m and ${\widehat{{\rm{\Omega }}}}_{0}$, we first choose two trial values for $\widehat{\omega }$ and Φ_s1 at η = η₀, and calculate $d{\chi }_{1}/d\eta$ and dΦ_s1/dη at the outer boundary from Equations (46) and (49). Next, we integrate Equations (41) and (42) from η = η₀ to η = 0 and check if the two conditions in Equation (43) are satisfied. If not, we vary Φ_s1(η₀) and $\widehat{\omega }$ one by one and repeat the calculations until the inner boundary conditions are fulfilled within tolerance.

Figures 11(a)–(c) plots the imaginary parts of eigenfrequencies for gravitationally unstable modes for isothermal rings with R₀/H = 10, 30, and 60 (or α = 1.04 × 10⁻¹, 1.15 × 10⁻², and 2.88 × 10⁻³) as functions of the azimuthal mode number m and the rotational angular frequency ${\widehat{{\rm{\Omega }}}}_{0}$. For all cases, the ring thickness is fixed to η₀/R₀ = 0.1. When ${\widehat{{\rm{\Omega }}}}_{0}=0$, the maximum growth rates are Im $({\widehat{\omega }}_{{\rm{\max }}})=0.88$, 1.01, and 1.19, which are achieved at m_max = 6, 10, and 18, for the rings with R₀/H = 10, 30, and 60, respectively. Note that the dispersion relations with ${\widehat{{\rm{\Omega }}}}_{0}=0$ are identical to those of axisymmetric modes in an infinite isothermal cylinder, presented by Nagasawa (1987), as long as m/R₀ is replaced with the vertical wavenumber k. Without rotation, $\widehat{\omega }$ is always an imaginary number, corresponding to pure instability. When ${\widehat{{\rm{\Omega }}}}_{0}\ne 0$, however, eigenfrequencies are complex numbers with the real parts almost linearly proportional to m, corresponding to overstability, as exemplified in Figure 11(d) for R₀/H = 30. This is a generic property of any instability occurring in a non-static background medium (e.g., Matsumoto et al. 1994; Kim et al. 2014, 2015).

Zoom In Zoom Out Reset image size

It is apparent that a ring with larger R₀/H (or smaller α) is more unstable owing to a smaller sound speed and/or a larger ring mass. Overall, rotation tends to stabilize the instability, reducing both the growth rate and the unstable range of m, although their dependence on ${\widehat{{\rm{\Omega }}}}_{0}$ is not simple. When R₀/H = 10 (or 30), the reduction of the growth rate due to rotation is larger at larger (or smaller) m, making m_max shifted to a smaller (or larger) value as ${\widehat{{\rm{\Omega }}}}_{0}$ increases. In the case of R₀/H = 60, on the other hand, rotation simply reduces the growth rate, without much affecting the unstable range of m. The instability is completely suppressed when ${\widehat{{\rm{\Omega }}}}_{0}\gtrsim 0.81$, 0.64, and 3.70 for R₀/H = 10, 30, and 60, respectively. Figure 12 plots the critical angular frequency ${\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}$ against α as the solid line, with the upper-right region corresponding to the stable regime. Note that ${\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}$ is almost constant at ∼0.7 for α ≳ 0.01 and increases rapidly as α decreases.

Zoom In Zoom Out Reset image size

It is interesting to compare our results for ${\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}$ with those from other simple estimates. The Toomre stability parameter ${Q}_{T}={\kappa }_{0}{c}_{s}/(\pi G{{\rm{\Sigma }}}_{0})$ has usually been invoked to judge whether a flattened system under consideration is gravitationally stable or not. For a ring, it is unclear how to choose Σ₀ since the ring surface density varies with R. If we simply take ${{\rm{\Sigma }}}_{0}=2{\rho }_{c}H$, about a half of the maximum surface density ${{\rm{\Sigma }}}_{{\rm{\max }}}\ =2{\int }_{0}^{\infty }{\rho }_{{\rm{sr}}}d\eta ={2}^{1/2}\pi {\rho }_{c}H$ across the ring center, the Toomre condition for marginal stability Q_T = 1 corresponds to ${\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}={(\pi /4)}^{1/2}\approx 0.89$, independent of α, which is plotted as the horizontal dotted line in Figure 12. The critical angular frequency from the Toomre condition is close to the results of our stability analysis for α ≳ 0.01, but deviates considerably for smaller α. Strictly speaking, the Toomre condition is valid only for thin disks that are infinitesimally thin in the vertical direction but infinitely extended in the horizontal direction. Even the thin-disk gravity underestimates self-gravity of highly concentrated rings with α ≲ 0.01.

Elmegreen (1994) presented a local dispersion relation for gravitational instability of nuclear rings by treating them as being thin and locally cylindrical. In Appendix C, we solve Equations (41) and (42) for local waves that vary very rapidly in the azimuthal direction (i.e., $m/{R}_{0}\gg | d\mathrm{ln}(\eta {\rho }_{0})/d\eta |$) but remain constant in the η-direction (i.e., $d{\chi }_{1}/d\eta =0$). The resulting dispersion relation (Equation (84)) is the same as the one given in Elmegreen (1994) (in the absence of magnetic fields and gas accretion). The critical angular frequency for local waves is then given by

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}^{2}=\mathop{\max }\limits_{m}\left\{\pi \left[1-\displaystyle \frac{m{\eta }_{0}}{{R}_{0}}{K}_{1}\left(\displaystyle \frac{m{\eta }_{0}}{{R}_{0}}\right)\right]-\displaystyle \frac{\alpha }{4}{m}^{2}\right\},\end{eqnarray} \tag{ 51 }$$

which is plotted in Figure 12 as the blue dashed line for η₀/R₀ = 0.1. Although ${\widehat{{\rm{\Omega }}}}_{0,\mathrm{crit}}$ from Equation (51) is similar to the results of the full analysis for α ∼ 3 × 10⁻², it underestimates the latter considerably for α ≲ 5 × 10⁻³. This is because rings with smaller α are increasingly more strongly concentrated that the approximations of constant ρ₀ and χ₁ over η become invalid.⁴ For $\alpha \propto {(H/{R}_{0})}^{2}\to 0$, rings can be approximated as strongly concentrated cylinders for which motions along the symmetry axis do not affect their gravitational instability much. Very large values of ${\widehat{{\rm{\Omega }}}}_{0}$ are required to stabilize such small-α rings.

We thus conclude that both the Toomre condition and the local results cannot adequately describe the critical angular frequencies of nuclear rings, especially when α is very small. In Section 5.2, we will apply the stability curve derived from the full stability analysis to observed rings in real galaxies.

To check the results of our linear stability analysis as well as to explore the effect of differential rotation, we conduct direct numerical simulations for gravitational instability of a slender isothermal ring using the GIZMO code (Hopkins 2015, 2016). GIZMO is a second-order accurate magnetohydrodynamics code based on a Lagrangian, mesh-free, Godunov-type method, and conserves mass, momentum, and energy almost exactly. In our calculations, the basic equations of hydrodynamics are solved by the Meshless Finite-Mass Method known to conserve angular momentum very well.

To obtain an initial particle distribution, we first use the SCF method to construct the equilibrium density distribution of the rigidly rotating ring with R₀/H = 30 and η₀/R₀ = 0.1 in the absence of external gravity. The rotation frequency of this ring is found to be ${\widehat{{\rm{\Omega }}}}_{s}=0.22$. The initial particle positions are then sampled by a rejection technique that uses Halton's quasi-random sequences over usual random numbers in order to reduce Poisson noises (e.g., Press et al. 1988). We then impose a radially varying external gravity ${\widehat{{\rm{\Omega }}}}_{e}^{2}(R)\widehat{{\boldsymbol{R}}}$ to boost the angular velocity of the ring particles according to

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{0}={({\widehat{{\rm{\Omega }}}}_{s}^{2}+{\widehat{{\rm{\Omega }}}}_{e}^{2})}^{1/2}=0.30{\left(\displaystyle \frac{R}{{R}_{0}}\right)}^{-q}.\end{eqnarray} \tag{ 52 }$$

Here, q is the rate of shear in the ring rotation, such that q = 0 and 1 corresponds to rigid-body and flat rotation, respectively. We vary q from 0 to 1.5 to study the effect of shear on the growth of gravitational instability.

To represent a hot tenuous external medium, we also distribute low-mass particles/cells outside the ring but inside a cylindrical volume with radius R/R₀ = 2 and height h/R₀ = 2. We let the external medium follow the rotation law given in Equation (52), and adjust its density distribution to balance the centrifugal force and gravity of the ring. We ensure a pressure equilibrium between the ring and the external medium that has 100 times lower density than the ring at the contact surfaces. The number of particles/cells for the ring and external medium is 5 × 10⁵ each.

At the very initial phase of evolution ($\tau \equiv t{(G{\rho }_{c})}^{1/2}\lesssim 0.25$), we iron out any residual Poisson sampling noises by introducing an artificial damping force $f\propto -{v}_{\eta }\widehat{\eta }$. The system thus evolves from a smooth, steady equilibrium state, without undergoing violent expansion or contraction, and gradually picks up gravitationally unstable modes. Figure 13 plots snapshots of the projected density ${\rm{\Sigma }}=\int \rho {dz}$ onto the equatorial plane at τ = 0, 5, 7, and 8.8 for a rigidly rotating model with q = 0. Defining the line density as ${ \mathcal L }\equiv \int {\rm{\Sigma }}{dR}$, we calculate the amplitude ${{ \mathcal L }}_{m}$ of each azimuthal mode m via Fourier transform at each time. Figure 14 displays ${{ \mathcal L }}_{m}$ as functions of τ for even-m modes in the q = 0 model. Note that the modes with m = 8–12 dominate during the linear growth phase (4 ≲ τ ≲ 6), resulting in 11 clumps at the end of the run. The growth rates of these modes are all similar at ∼(0.80–0.82) (Gρ_c)^1/2, as indicated as the dashed-line segment, consistent with the linear results shown in Figure 11(b). This indirectly confirms that the assumptions made in our stability analysis are quite reasonable.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The left panels of Figure 15 plot the snapshots of surface density in the equatorial plane at the end of the runs (at τ = 8.8, 9.0, and 9.1) for models with q = 0.5, 1.0, and 1.5 from top to bottom, while the right panels give the temporal evolution of ${{ \mathcal L }}_{m}$ for even-m modes. In the q = 0.5 model, the m = 8 and 10 modes dominate almost equally, while the models with q = 1.0 and 1.5 are dominated by the m = 10 and m = 8 mode, respectively. Note that the growth rates of the dominant modes in all models are very close to $0.81{(G{\rho }_{c})}^{1/2}$, marked by the dashed-line segment in each panel. The number of clumps produced as a result of gravitational instability is 10 or 11, insensitive to q, demonstrating that shear does not affect the character of gravitational instability of slender rings.

Zoom In Zoom Out Reset image size

Nuclear rings at centers of barred galaxies exhibit strong star formation activities. They are thought to undergo gravitational instability when sufficiently massive. To study their equilibrium properties and stability to gravitational perturbations, we approximate nuclear rings as isothermal objects. We modify the SCF method of Hachisu (1986a) to make it suitable for an isothermal equation of state, and construct equilibrium sequences of rigidly rotating, self-gravitating, isothermal bodies. A steady equilibrium is uniquely specified by two dimensionless parameters: α and ${\widehat{R}}_{{\rm{B}}}$ (see Equations (7) and (9)). The former is the measure of the thermal energy relative to gravitational potential energy of an equilibrium body, while the latter corresponds to the ellipticity for spheroid-like configurations or the thickness for ring-like configurations. We take a convention that ${\widehat{R}}_{{\rm{B}}}$ is positive (or negative) for spheroid-like (or ring-like) objects.

To test our SCF method, we first apply it to the case of rotating incompressible bodies, and confirm that our method is able to reproduce the Maclaurin spheroid sequence when $0.158\leqslant {\widehat{R}}_{{\rm{B}}}\leqslant 1$. With improved resolution, our method gives more accurate results than those obtained by Eriguchi & Sugimoto (1981) and Hachisu (1986a) for the concave hamburger sequence with $0\lesssim {\widehat{R}}_{{\rm{B}}}\lt 0.158$. Our method also successfully reproduces isothermal Bonnor–Ebert spheres, with larger α corresponding to a higher degree of central density concentration.

We then use our SCF method to obtain the density distributions of rotating isothermal equilibria on the meridional plane, as illustrated in Figure 4. We calculate the dependence on ${\widehat{R}}_{{\rm{B}}}$ of various dimensionless quantities such as the rotational angular frequency ${\widehat{{\rm{\Omega }}}}_{s}$, the total mass $\widehat{M}$, the mean density $\langle \widehat{\rho }\rangle$, the total kinetic energy $\widehat{T}$, and the gravitational potential energy $\widehat{W}$. These values are tabulated in Table 2 and given graphically in Figure 6. We find that an equilibrium density profile is more centrally concentrated for smaller α. Unlike the incompressible bodies, not all values of ${\widehat{R}}_{{\rm{B}}}$ result in an isothermal equilibrium configuration. Spheroid-like equilibria exist only for ${\widehat{R}}_{{\rm{B,1}}}\leqslant {\widehat{R}}_{{\rm{B}}}\leqslant 1$, while ring-like (or hamburger-like) configurations are possible only for $-1\lt {\widehat{R}}_{{\rm{B}}}\lt {\widehat{R}}_{{\rm{B,2}}}$: otherwise, the centrifugal potential is too large to form gravitationally bound objects. The critical ${\widehat{R}}_{{\rm{B}}}$ values are found to be ${\widehat{R}}_{{\rm{B,1}}}=0.27$, 0.51, and 0.59, and ${\widehat{R}}_{{\rm{B,2}}}=0.13$, −0.14, and −0.40 for α = 1, 0.1, and 0.01, respectively.

In general, ${\widehat{{\rm{\Omega }}}}_{s}$ is a decreasing function of $| {\widehat{R}}_{{\rm{B}}}|$. This is naturally expected for spheroid-like configurations since faster rotation leads to a more flattened equilibrium. As ${\widehat{R}}_{{\rm{B}}}$ approaches −1, on the other hand, ring-like configurations becomes less massive and thus requires weaker centrifugal force to balance self-gravity. Due to stronger central concentration, ${\widehat{{\rm{\Omega }}}}_{s}$, $\widehat{M}$, $\langle \widehat{\rho }\rangle$, $\widehat{T}$, and $| \widehat{W}|$ all become smaller as α decreases. For $\alpha \lt 0.1$, $\widehat{M}$ and $\widehat{W}$ are insensitive to ${\widehat{R}}_{{\rm{B}}}\gtrsim 0.6$ for spheroid-like equilibria since the density in the outer parts becomes vanishingly small. For a given value of the normalized angular momentum j, the normalized angular frequency ω_s becomes smaller with decreasing α, although the energy ratio $T/| W|$ is insensitive to α.

The density distribution of finite slender rings obtained by our SCF method for ${\widehat{R}}_{{\rm{B}}}\lesssim -0.6$ is found to be well approximated by Equation (22), which is also the solution for static, isothermal cylinders of infinite extent. This indicates that the rotation as well as geometrical curvature effect are insignificant in determining an equilibrium for rings with the major axis R₀ much longer than the minor axis η₀. The equilibrium angular frequency for isothermal slender rings with α ≳ 0.1 is well described by Equation (20) applicable to truncated incompressible rings (Ostriker 1964b).

To explore gravitational instability of nuclear rings, we calculate the growth rates of non-axisymmetric modes with azimuthal mode number m by assuming that the rings are slender with η₀/R₀ = 0.1, and that perturbations are independent of the polar angle λ in the meridional plane. In the absence of rotation, the resulting dispersion relations are the same as those of axisymmetric modes for an infinite isothermal cylinder studied by Nagasawa (1987) if m/R₀ is taken equal to the wavenumber in the direction along the cylinder (see Figure 11). Only large-scale modes can be gravitationally unstable, and the unstable range of m as well as the maximum growth rate increase with decreasing α.

Rotation tends to stabilize the gravitational instability, reducing both the growth rates and the unstable ranges of m. The instability is completely suppressed when ${\widehat{{\rm{\Omega }}}}_{0}$ exceeds the critical value that is relatively constant at ∼0.7 for α ≳ 0.01 and increases rapidly with decreasing α (see Figure 12). The simple estimates of the critical angular frequencies from the Toomre condition as well as the local dispersion relation are smaller than the results of our full stability analysis for α ≲ 5 × 10⁻³ due to underestimation of self-gravity at the ring centers. Shear turns out to be unimportant for the gravitational instability of rings as long as they are slender.

Mazzuca et al. (2008) analyzed photometric data of a sample of nuclear rings to estimate the ages of Hα-emitting star clusters and found that about a half of their sample contains an age gradient of star clusters along the azimuthal direction. Since nuclear rings with age gradient are thought to be gravitationally stable and form stars preferentially at the contact points, it is interesting to apply the results of our linear stability analysis to the observed rings to tell whether they are really stable.

Table 3 lists the properties of 19 observed nuclear rings in galaxies with a noticeable bar, compiled from the literature where the information on the presence/absence of an age gradient is available.⁵ Column (1) lists each galaxy name. Columns (2) and (3) give the semimajor axis and ellipticity of nuclear rings adopted from Comerón et al. (2010). Column (4) lists the rotational velocity v_rot adopted from Mazzuca et al. (2008) or from the references given in Column (9). Column (5) lists the total gas mass M_g in the ring from Sheth et al. (2005) or the references in Column (9) only for the galaxies with available data; we otherwise take M_g = 4 × 10⁸ ${M}_{\odot }$ as a reference value. Column (6) indicates the presence or absence of an azimuthal age gradient of star clusters adopted from Mazzuca et al. (2008), Allard et al. (2006), Sandstrom et al. (2010), and Brandl et al. (2012): a question mark is used when it is difficult to characterize the age distribution. Columns (7) and (8) give α and ${\widehat{{\rm{\Omega }}}}_{0}$ calculated by

$$\begin{eqnarray}&&\alpha =0.01{\left(\displaystyle \frac{{c}_{s}}{10{\text{km s}}^{-1}}\right)}^{2}{\left(\displaystyle \frac{{M}_{g}}{4\times {10}^{8}{M}_{\odot }}\right)}^{-1}\left(\displaystyle \frac{{R}_{0}}{1\,{\rm{kpc}}}\right),\end{eqnarray} \tag{ 53 }$$

and

$$\begin{eqnarray}&&{\widehat{{\rm{\Omega }}}}_{0}=2.1\left(\displaystyle \frac{{v}_{{\rm{rot}}}}{200\,{\text{km s}}^{-1}}\right){\left(\displaystyle \frac{{M}_{g}}{4\times {10}^{8}{M}_{\odot }}\right)}^{-0.5}{\left(\displaystyle \frac{{R}_{0}}{1{\rm{kpc}}}\right)}^{0.5},\end{eqnarray} \tag{ 54 }$$

after taking ${R}_{0}={R}_{1}{(1-{e}^{2})}^{1/4}$ corresponding to the geometric means of the major and minor axes of eccentric rings, ${c}_{s}=10\,{\text{km s}}^{-1}$, and η₀ = R₀/10. We replace ρ_c with $\langle \rho \rangle$ since the ring central density is difficult to constrain observationally.

Table 3.  Properties of Observed Nuclear Rings

	R1		vrot	Mg
Galaxy	(kpc)	e	(km s−1)	(107 ${M}_{\odot }$)	Age Grad.	α	${\widehat{{\rm{\Omega }}}}_{0}$	Ref.
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)
NGC 473	1.69	0.06	125	40	Yes	1.69E−2	1.71	⋯
NGC 613	0.40	0.26	115	40	?	3.93E−3	0.76	⋯
NGC 1097	0.97	0.32	220	140	No	2.70E−3	1.20	(a), (b)
NGC 1300	0.40	0.15	155	40	No	3.98E−3	1.03	⋯
NGC 1343	1.97	0.30	80	40	Yes	1.92E−2	1.17	⋯
NGC 1530	1.20	0.80	180	40	Yes	9.30E−3	1.82	⋯
NGC 2903	0.16	0.32	60	35	?	1.78E−3	0.27	(c)
NGC 3351	0.15	0.11	120	31	?	1.89E−3	0.55	(c)
NGC 4303	0.35	0.11	90	42	No	3.32E−3	0.54	⋯
NGC 4314	0.56	0.31	160	21	Yes	1.04E−2	1.71	(d), (e)
NGC 4321	0.87	0.32	170	51	Yes	6.64E−3	1.45	⋯
NGC 5248	0.65	0.20	150	42	Yes	6.13E−3	1.23	⋯
NGC 5728	1.10	0.23	180	40	Yes	1.09E−2	1.97	⋯
NGC 5905	0.39	0.14	150	40	?	3.88E−3	0.98	⋯
NGC 5953	1.00	0.43	150	40	?	9.50E−3	1.54	⋯
NGC 6951	0.56	0.17	160	40	Yes	5.56E−3	1.25	⋯
NGC 7552	0.34	0.15	150	40	No	3.38E−3	0.92	(f)
NGC 7716	1.20	0.04	150	40	No	1.20E−2	1.72	⋯
NGC IC14	0.68	0.36	204	40	Yes	6.57E−3	1.74	⋯

Note. Columns (2) and (3) give the semimajor axis and ellipticity of nuclear rings adopted from Comerón et al. (2010). Column (4) is the rotational velocity adopted from Mazzuca et al. (2008) or from the references given in Column (9). Column (5) is the total gas mass in the ring from Sheth et al. (2005) or from references in Column (9); we take M_g = 4 × 10⁸ ${M}_{\odot }$ if no information is available. Column (6) cites the age distribution: "Yes" and "No" for the presence and absence of an azimuthal age gradient, respectively, and "?" for uncertain cases, adopted from Allard et al. (2006), Mazzuca et al. (2008), Sandstrom et al. (2010), and Brandl et al. (2012). Columns (7) and (8) give α and ${\widehat{{\rm{\Omega }}}}_{s}$ calculated by Equations (53) and (54). Column (9) is the references for v_rot or M: (a) Onishi et al. (2015), (b) Hsieh et al. (2011), (c) Mazzuca et al. (2011), (d) Garcia-Barreto et al. (1991), (e) Benedict et al. (1996), (f) Brandl et al. (2012).

Download table as:  ASCIITypeset image

In Figure 16, we plot ${\widehat{{\rm{\Omega }}}}_{0}$ against α for the rings listed in Table 3 using various symbols, with numbers indicating galaxy names. Overall, rings with larger α tend to have larger ${\widehat{{\rm{\Omega }}}}_{0}$. Blue circles represent rings with an azimuthal age gradient, while red diamonds are for those with no age gradient. Rings for which the age distribution cannot be judged are marked by star symbols. It is apparent that all rings with an azimuthal age gradient are located at the stable regime, while all rings with no age gradient, except NGC 7716, correspond to unstable configurations. These results are consistent with two modes of star formation proposed by Böker et al. (2008), such that rings sufficiently massive or rotating sufficiently slowly form stars in the popcorn style caused by gravitational instability, and thus do not show an apparent age gradient. On the other hand, star formation in stable rings may occur preferentially at the contact points to exhibit an azimuthal age gradient of star clusters like pearls on a string.

Zoom In Zoom Out Reset image size

The ring models we have considered so far ignored the effects of magnetic fields that are pervasive in galaxies (e.g., Beck et al. 1996; Fletcher et al. 2011). In spiral galaxies, the presence of toroidal magnetic fields is known to play a destabilizing role in forming giant clouds inside spiral arms where tension forces from bent field lines resist the stabilizing effect of the Coriolis force (e.g., Balbus 1988; Kim & Ostriker 2001, 2002, 2006). In addition, magnetic fields are likely to reduce the degree of central density concentration by exerting pressure forces. It will be interesting to study how magnetic fields embedded in nuclear rings change the critical angular frequencies for gravitational instability compared to those of unmagnetized rings.

For spheroid-like configurations, it is convenient to solve Equation (3) in spherical polar coordinates (r, θ, ϕ). When an equilibrium is axially symmetric, we may expand the density using the Legendre function ${P}_{l}(\cos \theta )$ as

$$\begin{eqnarray}&&\rho (r,\theta )=\displaystyle \sum _{l=0}^{\infty }{\rho }_{l}(r){P}_{l}(\cos \theta ),\end{eqnarray} \tag{ 55 }$$

with the coefficients ρ_l given by

$$\begin{eqnarray}&&{\rho }_{l}(r)=\displaystyle \frac{2l+1}{2}{\int }_{0}^{\pi }\rho (r,\theta ){P}_{l}(\cos \theta )\sin \theta d\theta .\end{eqnarray} \tag{ 56 }$$

Then, Equation (3) has the series solution of the form

$$\begin{eqnarray}&&{\rm{\Phi }}(r,\theta )=-4\pi G\displaystyle \sum _{l=0}^{{l}_{{\rm{\max }}}}\displaystyle \frac{{P}_{l}(\cos \theta )}{2l+1}\\ &&\quad \times \left[\displaystyle \frac{1}{{r}^{l+1}}{\displaystyle \int }_{0}^{r}{{dr}}^{\prime }{r}^{\prime l+1}{\rho }_{l}({r}^{\prime })+{r}^{l}{\displaystyle \int }_{r}^{\infty }\displaystyle \frac{{{dr}}^{\prime }}{{r}^{\prime l-1}}{\rho }_{l}({r}^{\prime })\right],\end{eqnarray} \tag{ 57 }$$

(see, e.g., Equation (2.95) of Binney & Tremaine 2008; Sirotkin & Kim 2009). The summation in Equation (57) is truncated typically at l_max = 50 in our computation.

The multipole expansion using Legendre functions described above becomes inefficient for highly flattened systems such as slender rings/tori because it requires a high cutoff number l_max for accurate potential evaluation. To resolve such equilibria, the SCF method in the spherical coordinates also needs a very large computation grid. If we instead employ toroidal coordinates, flattened equilibria are well resolved with a relatively small grid and a very small number of multipole terms.

The toroidal coordinates (τ, σ, ϕ) are defined by

$$\begin{eqnarray}&&\left(\begin{array}{c}x\\ y\\ z\end{array}\right)=\left(\begin{array}{c}a\sinh \sigma \cos \phi \\ a\sinh \sigma \sin \phi \\ a\sin \tau \end{array}\right)/f(\sigma ,\tau ),\end{eqnarray} \tag{ 58 }$$

where a is a constant (focal length) and

$$\begin{eqnarray}&&f(\sigma ,\tau )\equiv \cosh \sigma -\cos \tau .\end{eqnarray} \tag{ 59 }$$

Note that $\sigma \in [0,\infty )$ is a (dimensionless) radial distance, $\tau \in [0,2\pi )$ is a poloidal angle, and $\phi \in [0,2\pi )$ is equal to the usual cylindrical azimuthal angle. The constant-σ surface is a torus with a circular cross-section: it becomes a focal ring with radius a when $\sigma =\infty$.

In this coordinate system, the Green's function for a Laplacian operator can be expanded based on the half-integer degree Legendre functions ${P}_{l-1/2}^{m}$ and ${Q}_{l-1/2}^{m}$ as

$$\begin{eqnarray}&&\displaystyle \frac{1}{| {\boldsymbol{x}}-{{\boldsymbol{x}}}^{\prime }| }=\displaystyle \frac{{({{ff}}^{\prime })}^{1/2}}{a\pi }\displaystyle \sum _{l,m=0}^{\infty }{\epsilon }_{l}{\epsilon }_{m}{(-i)}^{m}\displaystyle \frac{{\rm{\Gamma }}\left(l-m+\tfrac{1}{2}\right)}{{\rm{\Gamma }}\left(l+m+\tfrac{1}{2}\right)}\\ &&\quad \times \cos [l(\tau -{\tau }^{\prime })]\cos [m(\phi -{\phi }^{\prime })]\\ \quad \times \left\{\begin{array}{ll}{P}_{l-1/2}^{m}(\cosh \sigma ){Q}_{l-1/2}^{m}(\cosh {\sigma }^{\prime }),{\rm{}}\,{\rm{for}}\,{\rm{}} & {\sigma }^{\prime }\gt \sigma ,\\ {P}_{l-1/2}^{m}(\cosh {\sigma }^{\prime }){Q}_{l-1/2}^{m}(\cosh \sigma ),{\rm{}}\,{\rm{for}}\,{\rm{}} & {\sigma }^{\prime }\lt \sigma ,\end{array}\right.\end{eqnarray} \tag{ 60 }$$

where ${f}^{\prime }\equiv f({\sigma }^{\prime },{\tau }^{\prime })$, and _k is equal to unity for k = 1 and two otherwise (e.g., Morse & Feshbach 1953; Cohl & Tohline 1999; Cohl et al. 2000).⁶

The gravitational potential is then given by

$$\begin{eqnarray}{\rm{\Phi }}(\sigma ,\tau ) & = & -G{\displaystyle \int }_{V}\displaystyle \frac{\rho ({{\boldsymbol{x}}}^{\prime })}{| {\boldsymbol{x}}-{{\boldsymbol{x}}}^{\prime }| }{d}^{3}{x}^{\prime }\\ & = & -G{\displaystyle \int }_{0}^{2\pi }d{\phi }^{\prime }{\displaystyle \int }_{0}^{2\pi }d{\tau }^{\prime }{\displaystyle \int }_{0}^{\infty }d{\sigma }^{\prime }\displaystyle \frac{\rho ({{\boldsymbol{x}}}^{\prime })}{| {\boldsymbol{x}}-{{\boldsymbol{x}}}^{\prime }| }\displaystyle \frac{{a}^{3}\sinh {\sigma }^{\prime }}{f{({\sigma }^{\prime },{\tau }^{\prime })}^{3}}.\end{eqnarray} \tag{ 61 }$$

For an axisymmetric mass distribution ρ(σ, τ), Equation (61) reduces to

$$\begin{eqnarray}\displaystyle \frac{{\rm{\Phi }}(\sigma ,\tau )}{4{a}^{2}{Gf}{(\sigma ,\tau )}^{1/2}} & = & -\displaystyle \sum _{l=0}^{{l}_{\max }}{\epsilon }_{l}{P}_{l-1/2}(\cosh \sigma )\cos (l\tau )\\ & & \times {\displaystyle \int }_{\infty }^{\sigma }d{\sigma }^{\prime }\sinh {\sigma }^{\prime }{Q}_{l-1/2}(\cosh {\sigma }^{\prime }){\rho }_{l}({\sigma }^{\prime })\\ & & -\displaystyle \sum _{l=0}^{{l}_{\max }}{\epsilon }_{l}{Q}_{l-1/2}(\cosh \sigma )\cos (l\tau )\\ & & \times {\displaystyle \int }_{0}^{\sigma }d{\sigma }^{\prime }\sinh \sigma ^{\prime} {P}_{l-1/2}(\cosh {\sigma }^{\prime }){\rho }_{l}({\sigma }^{\prime }),\end{eqnarray} \tag{ 62 }$$

where

$$\begin{eqnarray}&&{\rho }_{l}({\sigma }^{\prime })={\int }_{0}^{\pi }\displaystyle \frac{\rho ({\sigma }^{\prime },{\tau }^{\prime })\cos (l{\tau }^{\prime })}{f{({\sigma }^{\prime },{\tau }^{\prime })}^{5/2}}d{\tau }^{\prime },\end{eqnarray} \tag{ 63 }$$

and we assume a reflection symmetry of ρ(σ, τ) with respect to the τ = 0 plane.

To compare the rate of convergence against l_max as well as accuracy between the results from the two different multipole expansions, we calculate ring-like equilibria with ${\widehat{R}}_{{\rm{B}}}\leqslant -0.8$ on both spherical and toroidal meshes. Figure 17(a) plots ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ as a function of l_max for an equilibrium with ${\widehat{R}}_{{\rm{B}}}=-0.8$ and α = 0.01. The blue and red lines correspond to the results with the spherical and toroidal multipole expansions, respectively. The inset zooms in the results with the toroidal mesh, showing that the SCF method converges to ${\widehat{{\rm{\Omega }}}}_{0}^{2}=0.048$ very rapidly when the potential is expanded in the toroidal coordinates. On the other hand, the SCF method with the spherical multipole expansions converges quite slowly, requiring l_max ≳ 42 to obtain the solution within 1% of the value from the toroidal multipole expansion shown as the dotted line.

Zoom In Zoom Out Reset image size

Figure 17(b) compares ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ obtained by the spherical mutipole expansion with l_max = 50 and the toroidal mutipole expansion with l_max = 10 for ring-like equilibria with $-0.80\geqslant {\widehat{R}}_{{\rm{B}}}\geqslant -0.99$ and α = 10⁻² (solid lines) and α = 10⁵ (dotted lines). The relative errors given in Figure 17(c) show that the two methods agree within 1% for ${\widehat{R}}_{{\rm{B}}}\gtrsim -0.86$, and the spherical expansions overestimate ${\widehat{{\rm{\Omega }}}}_{s}^{2}$ for smaller ${\widehat{R}}_{{\rm{B}}}$. This is mainly because such flattened equilibria are not well resolved in the spherical mesh. For this reason, we employ the toroidal potential expansion in constructing equilibria with ${\widehat{R}}_{{\rm{B}}}\leqslant -0.8$ presented in Section 3.