%
%\documentclass[12pt]{report}
\documentclass[12pt,preprint]{aastex}  %% double-column, single-spaced document
%\documentclass[12pt,article]{emulateapj}  %% double-column, single-spaced document
%\documentclass[12pt]{article}
%\documentstyle{report}
%\usepackage{graphics,graphicx,sectsty}
\usepackage{graphics,graphicx}


\voffset -0.5in
%\hoffset 0.5in
\setlength{\textheight}{9.0in}
\setlength{\textwidth}{6.5in}
\setlength{\topmargin}{0.0cm}
%\addtolength{\oddsidemargin}{-1.5cm}
%\addtolength{\evensidemargin}{-1.5cm}
\newcommand{\mso}{$M_{\odot}$ }
\newcommand{\rso}{$R_{\odot}$ }
\newcommand{\msoyr}{$M_{\odot}/{\rm year}$}
\newcommand{\be}{\begin{equation}}
\newcommand{\ee}{\end{equation}}
\newcommand{\etal}{{\em et al.}\ }
\newcommand{\cmm}{{cm$^{-1}$}\ }
\newcommand{\asec}{\mbox{\rlap{$^{\prime\prime}$}.\hbox to 2pt{}}}
%\pagestyle{empty} %suppresses page numbers 

%\pagestyle{headings}
%\pagenumbering{arabic}
%\begin{spacing}{1.66}
%\include{chap1}
%\end{spacing}
%\begin{spacing}{1.24}
%\include{bibliography}
%\end{spacing}

%\documentclass[manuscript]{aastex} % one-column, double-spaced document
%\documentclass[12pt,preprint]{aastex}
%\documentclass[letterpaper,12pt]{article}
%
%
%\setlength{\textwidth}{210mm}
%\setlength{\textheight}{9.0in}
%\setlength{\oddsidemargin}{25mm} \setlength{\evensidemargin}{25mm}
%\setlength{\topmargin}{-13mm}
%\setlength{\unitlength}{1mm}
%\def\gap{\rule[-5mm]{0mm}{1mm}}
%\newsavebox{\yesno}
%\newcounter{rxpage}
%\renewcommand{\thepage}{{page} }

\setcounter{secnumdepth}{4}
%\def\subsubsection{\@startsection{subsubsection}{3}{0pt}{-3.25ex plus
%-1ex minus -.2ex}{1.5ex plus .2ex}{\normalsize\it}}

%\def\subsubsection{\@startsection{subsubsection}{3}{0pt}{0pt}{0pt}{\normalsize\it}}

\renewcommand{\thesubsection}{\arabic{section}.\arabic{subsection}}
%\renewcommand{\thesubsubsection}{\arabic{section}.\arabic{subsection}.\arabic{subsubsection} }

%\renewcommand{\thepage}{\roman{page}}
\renewcommand{\thepage}{\arabic{page}}
\renewcommand{\thesection}{\arabic{section}}


\begin{document}


%\documentstyle[12pt,aaspp4,includegraphics[width=5.0in,angle=270]]{article}
%\tighten
%\begin{document}
%
%
\title{ NONAXISYMMETRIC INSTABILITIES OF SELF-GRAVITATING DISKS. $I$ TOROIDS }
%
%\author{ James N. Imamura and Katherine Hadley }
%
%\affil{Institute of Theoretical Science and Department of Physics,
%University of Oregon, Eugene, OR  97403}
%
\begin{abstract} 


We study the onset and development of nonaxisymmetric 
gravito-rotation instabilities in polytropic 
toroids with power law angular velocity distributions.
Nonaxismmetric modes with azimuthal dependence 
$e^{im\phi}$, where $m$ is the azimuthal mode number and 
$\phi$ is the aziumthal angle are considered. Earlier studies of 
self-gravitating annuli have identified three principal nonaxisymmetric 
disk modes, $I$ modes, $P$ modes, and $J$ modes. $I$ and $J$ 
modes appear only in self-gravitating disk systems.
For polytropic self-gravitating toroids with 
power law angular velocity distributions,
we find that instability sets in through
$m$ = 2 $I$ modes at $T/|W|$ $\sim$ 0.16-0.18 where
$T$ is the rotational kinetic
energy and $W$ is the gravitational energy. 
Instability in the $m$ = 2 $I$ mode
peaks in strength around $T/|W|$ = 0.22-0.23 and then 
weakens. For toroids with polytropic index $n$ = 3/2 and power 
law angular velocity distributions with indices 
$q$ = 3/2 and 2, the $m$ = 2 $I$ mode 
becomes stable at $T/|W|$ $\sim$ 0.26. Just before the $m$ = 2 
mode stabilizes, an $m$ = 1 $I$ mode instability sets in 
and dominates the $m$ = 2 mode 
for $T/|W|$ $\buildrel>\over\sim$ 
0.25. However, around when the $m$ = 1 mode starts to 
dominate the $m$ = 2 $I$ mode, $J$ modes appear.
$J$ modes with $m$ = 2, 3, and 4 become unstable
for $T/|W|$ $\sim$ 0.25-0.26, 0.23-0.25, and 0.25-0.26, 
respectively. The $m$ $\ge$ 3 $J$ modes quickly dominate the
$m$ = 1 $I$ mode and control the evolution of 
toroids with $T/|W|$ $\buildrel>\over\sim$ 0.26. Neither
$I$ mode nor $J$ mode instabilities produce prompt fission 
in the nonlinear regime.

\end{abstract}


\keywords{Hydrodynamics: Stability; Stars: Rotation}



\section{ INTRODUCTION }



Self-gravitating toroids and disks arise in a wide range 
of astrophysical environments playing crucial roles in the 
evolution of systems on galactic to stellar scales. The 
nature and stability properties of 
self-gravitating toroids and disks are not completely 
understood. Extensive work on infinitesimally 
thin accretion disks and self-gravitating annuli 
have been performed. Much less attention has
been directed toward thick self-gravitating toroids
and disks. Cold (thin) disks which, essentially, are
two-dimensional (2D), exhibit a rich spectrum of 
modes ({\it e.g.,} Lin \& Papaloizou 1995). It is
not clear how much of the 2D work carries over to
hot (thick) disks (toroids). The extra degree-of-freedom 
may introduce new types of modes as well as affect the 
properties of existing modes. Here, we begin an investigation of 
the nonaxisymmetric instabilities of thick, three-dimensional 
(3D) self-gravitating disks and toroids. In this paper, we consider 
self-gravitating toroids. In later works we study star/disk systems 
where the central star is treated as a point mass and where it is
resolved and we model the coupling of internal modes excited 
in the star with nonaxisymmetric modes excited in the disk. 


Self-gravitating toroids have been the subject of 
a handful of studies over the years (Wong 1974; 
Eriguchi \& Hachisu 1984; Hachisu \& Eriguchi 1985;
Hachisu, Tohline, \& Eriguchi
1987, 1988; Tohline \& Hachisu 1990; Woodward, Tohline 
\& Hachisu 1994). The earliest works which examined the properties
of incompressible fluids (cf. Wong 1974; Eriguchi \& Hachisu
1984 and references therein), showed that two ring
sequences branched off the Maclaurin sequence, the
one-ring and two-ring sequences (Hachisu \& Eriguchi 1985).
Similarly compressible toroids with polytropic equations of state
and angular momentum distributions of Maclaurin spheroids (the
angular momentum distribution of a uniform density, uniformly rotating
sphere) branch off the Maclaurin-like
sequence at $T/|W|$ = 0.41 (Hachisu, Tohline, \& Eriguchi 1987, 1988).
In both the compressible and incompressible cases, toroids appear
well after $T/|W|$ = 0.27 is reached, the bar mode stability limit
for rotating polytropes. Maclaurin-like toroids are thus likely to
be strongly dynamically unstable and not steps in secular
evolution of collapsing interstellar clouds. They may,
however, represent metastable states reached during
inviscid axisymmetric dynamical collapses (Tohline 1980).
Not as cleanly connected to the star formation problem but 
perhaps more useful from a theoretical standpoint are 
polytropic toroids with power law angular velocity distributions,
$\Omega(\varpi)$ $\propto$ $\varpi^{-q}$ ({\it e.g.}, Tohline 2002). 
Such toroids may arise after instabilities and viscosity
have acted and they may, in fact, be more likely 
to exist in nature because for a broad range of $q$, 
stable toroids exist down to zero
angular momentum unlike toroids with Maclaurin
spheroid angular momentum distributions where
equilibria exist only for high angular momentum.


Self-gravitating polytropic toroids with power law
angular velocity distributions $\Omega(\varpi)$ have been
investigated by Hachisu, Tohline, \& Eriguchi (1987,1988), Tohline
\& Hachisu (1990), and Woodward, Tohline, \& Hachisu (1994).
Tohline \& Hachisu (1990) and Woodward, Tohline, \& Hachisu (1994)
showed that barlike dynamic instabilities set in as early as
$T/|W|$ = 0.17 in toroids, while Hachisu, Tohline, \& Eriguchi (1987,1988) 
argued that instability could set in as early as $T/|W|$ = 0.142
based on time scale arguments. Both thresholds
are well below $T/|W|$ = 0.27. In this paper, we perform an
extensive study of dynamic nonaxisymmetric instabilities 
in self-gravitating polytropic toroids. 
We determine the general properties of the fastest 
growing modes for unstable equilibrium models and compare
and contrast the nonaxisymmetric modes that arise in toroidal
structures to those that arise in polytropic stars 
({\it cf.}, Toman {\it et al}. 1998, Imamura {\it et al}. 2000), 
and to the nonaxisymmetric modes isolated in incompressible 
narrow toroids (ICTs) ({\it e.g.,} Andalib, Tohline, \&
Christodoulou 1997). 
Parameter space for incompressible narrow toroids (ICTs) with 
constant specific angular momentum has been explored, {\it e.g.,} 
see Andalib, Tohline, \& Christodoulou (1997). The stability properties 
of toroids have been less well-studied. TH and WTH considered
the stability properties of even $m$ modes for seven toroids 
with $n$ = 3/2 for $q$ = 3/2, 7/4, and 2. We perform an
extensive exploration of parameter space for self-gravitating toroids 
and consider both even and odd $m$ modes.
For selected toroids, we follow the development of instability
in the linear regime,  
through the quasi-linear regime,
and into the fully 
nonlinear regime. 



The remainder of our paper is organized as follows. In 
$\S$2, our numerical methods are presented. In $\S$3, 
our numerical results are compared to other works. In
$\S$4, our conclusions are summarized.



\section{ EQUILIBRIUM TOROID MODELS }



The inviscid hydrodynamics equations are
\begin{eqnarray}
0 & = & \frac {\partial\rho} {\partial t} + \nabla\cdot(\rho {\bf v}), \\
0 & = & \rho\left( \frac {\partial} {\partial t} + 
{\bf v}\cdot\nabla\right) {\bf v}
+ \nabla P +\rho\nabla \Phi_g, \\
0 & = & \left( \frac {\partial} {\partial t} + {\bf v}\cdot\nabla\right)
\left( P\rho^{-\gamma} \right) -(\gamma-1)\left({\cal H}-{\cal L}\right)
\end{eqnarray}
where $t$ is time, $\rho$ is the density, ${\bf v}$ is the velocity, $P$ is
the pressure, $\Phi_g$ is the gravitational potential, $\gamma$ is the
adiabatic index, and ${\cal H}$ and ${\cal L}$ are the local heating 
and cooling rates. In steady-state, we drop the time derivatives 
and assume axisymmetry. For an adiabatic gas 
${\cal H}$ = ${\cal L}$ and we replace the energy equation with 
the polytropic pressure-density 
relationship, $P(\rho)$ = $K\rho^\gamma$ where
$K$ is the polytropic constant, 
$\gamma$ = $1+1/n$, and $n$ is the polytropic index.
The velocity field is defined with power law angular velocity
distribution $\Omega(\varpi)$ $=$ $\Omega_{\circ}(\varpi/R_{+})^{-q}$,
where $\Omega_{\circ}$ is a normalization constant,
$R_+$ is the outer radius of the disk, $\varpi$ is the cylindrical
radial coordinate and $q$ is a constant.
We consider toroids with $n$ = 3/2 and $q$ = 1, 3/2, and 2.
The equations are solved on 
cylindrical grids with uniformly sized cells and 
typical dimension $n_{\varpi}\times n_z$
= $256\times256$. We use a Hachisu-type (Hachisu 1986) 
self-consistent field code (Pickett, Durisen, \& Davis 1996).
Families of toroids are defined by $(n,q)$-pairs. 
Family members 
are parameterized by the ratio of the radius of the inner edge of the 
toroid $R_-$ to the radius of the outer edge of the toroid $R_+$.
There is a one-to-one correspondence
between the ratio $R_-/R_+$ and $T/|W|$ (and the angular momentum $J$)
for toroids with power law $\Omega(\varpi)$ 
distributions. For large $R_-/R_+$, toroids are more
symmetric with cross-sections approaching circular in shape.
Properties of the equilibrium toroids are
summarized in Table 1. In the far right column is listed 
$\epsilon$ = $(R_{max}-R_-)/R_{max}$, a parameter which
measures the width of the toroids. For toroids with small 
$T/|W|$ (small $R_-/R_+$), $\epsilon$ $\rightarrow$ 1. For toroids
with large $T/|W|$ (large $R_-/R_+$), 
$\epsilon$ $\rightarrow$ 0. 
The ICT approximations are most accurate when 
$\epsilon$ $\ll$ 1. Unless
otherwise noted, all quantities are presented in polytrope units, 
the unit system where $G$ = $K$ = $M$ = 1. More extensive 
compilations of the properties of polytropic toroids may be 
found in Hachisu, Tohline, \& Eriguchi (1987, 1988).



Equilibrium toroids display several characteristic radii,
radii important to understanding the nature of their eigenmodes 
({\it e.g.,} Papaloizou \& Savonije 1991). Among these are: 
(a) the location of the density maximum of the 
toroid $R_{max}$; (b) the extremum of 
the vortensity distribution $R_{\lambda}$, where
the vortensity ${\bf \lambda}$ is defined as the mass 
column density $\Sigma$ divided by the local vorticity.
For rotation on cylinders and power law $\Omega(\varpi)$,
\begin{equation}
{\bf \lambda} =
\left(\frac{\Sigma} {(2-q)\Omega}\right)\;\hat{{\bf z}}.
\end{equation}
(c) The boundaries of the region where 
the Toomre Q-parameter,
\begin{equation}
{\rm Q} =  \frac {c_s\kappa}{\pi{\rm G}\Sigma}
\end{equation}
drops below one.
Here, $c_s$ is the sound speed, $\kappa$ is the epicyclic
frequency (see Equation [10]), and $\Sigma$ is the surface density. $Q$
reaches maxima at the inner and outer edges of the toroid. 
with minimum 
near the density maximum. Local axisymmetric modes may be unstable 
where $Q$ $\buildrel<\over\sim$ 1.



\section{ NONAXISYMMETRIC INSTABILITIES IN TOROIDS }



\subsection{ Nonaxisymmetric Disk Modes}



Self-gravitating annuli in star/disk systems
exhibit various nonaxisymmetric modes, 
the self-gravitating $I$ and $J$ modes, 
and the $P$ modes ({\it e.g.,} see Andalib, Tohline, \& Christodoulou 1997).
(i) $P$ modes were discovered in studies of non-self-gravitating disks 
by Papaloizou \& Pringle (1982). Subsequent works showed that 
they persist into the self-gravitating 
regime but weaken greatly in strength with increasing 
self-gravity. They are driven by 
edge waves which couple across a forbidden region around 
corotation, where corotation falls at or near the location 
of the density maximum $\rho_{max}$. We
do not find $P$ modes in rapidly rotating toroids.
(ii) $I$ modes arise in self-gravitating disks and are
driven by the merger of two waves which carry equal but 
opposite amounts of angular momentum. One wave is driven by 
self-gravity and represents azimuthal compressions, the
other is due to free epicyclic motions which
have high frequency in the rotating frame. Corotation for $I$ modes falls
well outside the location of $\rho_{max}$.
(iii) $J$ modes appear only in self-gravitating disks. They are 
driven by the merger of two waves 
both of which are driven by self-gravity. $J$ 
modes are related to Jeans instabilities.
Corotation for $J$ modes falls near the location of $\rho_{max}$. 
Instability in $I$ modes sets in well before instability in $J$ modes 
in ICTs for toroids with increasing angular momentum and given 
star mass to disk mass ratio, $M_*/M_d$.



\subsection{ Numerical Modeling of Nonaxisymmetric Instabilities}



We find the fastest growing mode with given $m$ in the linear
regime by solving an
initial value problem (IVP) where the evolution of
random, low amplitude, Eulerian perturbations 
are followed away from equilibrium in unstable
toroidal configurations.
The IVP evolution equations are derived
by linearizing the inviscid hydrodynamics equations 
(Equs. [1]-[3]) using Eulerian
perturbations of the form
\begin{eqnarray}
\rho &=& \rho_o+\rho_1(\varpi,z,t){\rm exp}(im\phi),\\
v_{\varpi} &=& v_{\varpi,1}(\varpi,z,t){\rm exp}(im\phi),\\
v_{\phi} &=& \Omega\varpi+v_{\phi,1}(\varpi,z,t){\rm exp}(im\phi),\\
v_{z} &=& v_{z,1}(\varpi,z,t){\rm exp}(im\phi),
\end{eqnarray}
where $\rho_o$ and $\Omega$ are the equilibrium density 
and angular velocity, respectively. 
The linearized hydrodynamic equations are then
\begin{eqnarray}
\frac{\partial\rho_1}{\partial t}&=&-im\Omega\rho_1-
\rho_o\frac{v_{\varpi,1}}{\varpi}- v_{\varpi,1}\frac{\partial
\rho_o}{\partial\varpi}- v_{z,1}\frac{\partial\rho_o}{\partial 
z}-\rho_o(\frac{\partial
 v_{\varpi,1}}{\partial\varpi} 
 +\frac{im}{\varpi} v_{\phi,1}
 +\frac{\partial v_{z,1}}{\partial z}) \\
\frac{\partial v_{\varpi,1}}{\partial t}&=&-im\Omega v_{\varpi,1}
+ 2\Omega v_{\phi,1}
-\frac{\gamma P_o}{\rho_o^2}\frac{\partial
\rho_1}{\partial\varpi}-(\gamma-2)\frac{\rho_1}{\rho_o^2}
\frac{\partial P_o}
{\partial\varpi} 
-\frac{\partial\Phi_1}{\partial\varpi},\\
\frac{\partial v_{\phi,1}}{\partial t}&=&-im\Omega v_{\phi,1}
-\frac{v_{\varpi,1}}{\varpi}\frac{\partial(\Omega\varpi^2)}{\partial\varpi}
-\frac{im}{\varpi}\frac{\gamma P_o}{\rho_o^2}\rho_1-\frac{im}{\varpi}\Phi_1,
\end{eqnarray}
and
\begin{eqnarray}
\frac{\partial v_{z,1}}{\partial t}&=&-im\Omega v_{z,1}-\frac{\gamma P_o}{\rho_o^2}
\frac{\partial\rho_1}{\partial z}-(\gamma-2)\frac{\rho_1}{\rho_o^2}
\frac{\partial P_o}
{\partial z}-\frac{\partial\Phi_1} {\partial z}.
\end{eqnarray}
In the above, we suppressed the $``(\varpi,z,t)$'' from the 
notation describing the perturbations, and we let $P_o$ = 
$K\rho_o^{\gamma}$ where $\gamma$ = 1+1/$n$. The linearized equations are 
complex. We solve them as pairs of real equations by substituting 
\begin{eqnarray}
\rho_1&=&\eta_R+i\eta_I,\\
v_{\varpi,1}&=&T_R+iT_I,\\
v_{\phi,1}&=&V_R+iV_I,\\
v_{z,1}&=&W_R+iW_I,\\
\Phi_1&=&\Phi_R+i\Phi_I
\end{eqnarray}
into the linearized equations. We find 
\begin{eqnarray}
\frac{\partial\eta_R} {\partial t} &=&m\Omega\eta_I-T_R
\frac{\rho_o} {\varpi} (1+\frac{\varpi}{\rho_o}
\frac{\partial\rho_o} {\partial\varpi})
-W_R\frac{\partial\rho_o} 
{\partial z}-\rho_o(\frac{\partial T_R} {\partial\varpi}
-\frac{m} {\varpi}
V_I+\frac{\partial W_R} {\partial z}),\\
\frac{\partial\eta_I} {\partial t} &=&-m\Omega\eta_R-T_I\frac{\rho_o}{\varpi}(1+
\frac{\varpi}{\rho_o}\frac{\partial\rho_o}{\partial\varpi}
-W_I\frac{\partial\rho_o} {\partial z}-\rho_o(\frac{\partial T_I}{\partial\varpi}
+\frac{m}{\varpi}
V_R+\frac{\partial W_I}{\partial z}),\\
\frac{\partial T_R}{\partial t}&=&m\Omega T_I
+ 2\Omega V_R-\gamma \frac{P_o}{\rho_o^2}\frac{\partial
\eta_R}{\partial\varpi}
-(\gamma-2)\frac{\eta_R}{\rho_o^2}\frac{\partial P_o}
{\partial\varpi}-\frac{\partial\Phi_R}{\partial\varpi},\\
\frac{\partial T_I}{\partial t}&=&-m\Omega T_R
+ 2\Omega V_I-\gamma \frac{P_o}{\rho_o^2}\frac{\partial
\eta_I}{\partial\varpi}
-(\gamma-2)\frac{\eta_I}{\rho_o^2}\frac{\partial P_o}
{\partial\varpi}-\frac{\partial\Phi_I}{\partial\varpi},\\
\frac{\partial V_R}{\partial t}&=&m\Omega V_I
-\frac{T_R}{\varpi}\frac{\partial(\Omega\varpi^2)}{\partial\varpi}
+\frac{m}{\varpi}[\frac{\gamma P_o}{\rho_o^2}\eta_I+\Phi_I],\\
\frac{\partial V_I}{\partial t}&=&-m\Omega V_R
-\frac{T_I}{\varpi}\frac{\partial(\Omega\varpi^2)}{\partial\varpi}
-\frac{m}{\varpi}[\frac{\gamma P_o}{\rho_o^2}\eta_R+\Phi_R],\\
\frac{\partial W_R}{\partial t}&=&m\Omega W_I-\frac{\gamma P_o}{\rho_o^2}
\frac{\partial\eta_R}{\partial z}-(\gamma-2)\frac{\eta_R}{\rho_o^2}
\frac{\partial P_o}
{\partial z}-\frac{\partial\Phi_R}{\partial z},
\end{eqnarray}
and
\begin{eqnarray}
\frac{\partial W_I}{\partial t}&=&-m\Omega W_R-\frac{\gamma P_o}{\rho_o^2}
\frac{\partial\eta_I}{\partial z}-(\gamma-2)\frac{\eta_I}{\rho_o^2}
\frac{\partial P_o}
{\partial z}-\frac{\partial\Phi_I}{\partial z}.
\end{eqnarray}              

The Poisson equation for the gravitational perturbations is 
\begin{equation}
\nabla^2\Phi_1 = 4\pi G\rho_1,
\end{equation}
and so the $\Phi_R$ and $\Phi_I$ may be calculated from $\eta_R$ and 
$\eta_I$ individually.  We use a Poisson equation solver taken from the 
numerical hydrodynamics code of J. Tohline (Tohline 1980) to
calculate $\Phi_R$ and $\Phi_I$.


The following boundary conditions are imposed:  (1) the perturbations
are forced to be symmetric about the equatorial plane; (2) the 
perturbed radial and azimuthal velocities are set to 0 on the 
rotation axis; (3) the density perturbation is unconstrained  
on the surface of the unperturbed polytrope; and (4) the 
velocity perturbations are set to 0
on the surface of the unperturbed polytrope. 
The continuity of the gravitational potential at the surface 
of the polytrope is guaranteed by imposing boundary conditions on 
the potential on the surface of the cylinder defined 
by the maximum $\varpi$ and $z$ of the computational grid. The
potential on the surface of the cylinder is found by solving the 
integral form for the perturbed gravitational potential. 


The linearized evolution equations 
are solved as follows (Toman {\it et al.} 1998). 
We discretize space forming square cells in 
($\varpi$, z)-space.
The $\rho$'s are defined at cell centers and the $v$'s 
at cell vertices. The equations are solved at cell centers
by approximating the spatial derivatives with finite differences 
while leaving the time derivatives continuous. This leads to
a set of coupled total differential equations in the time, $t$.  
The system of coupled total differential equations is 
advanced in time using a Runge-Kutta algorithm.
The simulations are 
performed on the same grids as used for the equilibrium 
calculations. 


Exponentially growing global modes
are identified as follows. We plot the phase and
the logarithm of the amplitude of the density perturbation 
at three locations in the equatorial plane of the toroid. 
If the perturbation settles into the same pattern
of growth at the three radii, we 
start to monitor the full structure of the 
disturbance to see if the toroid has fallen into global growth. 
If the growth is uniform, then the growth rate 
for the mode is determined by a least squares fit to the exponential
growth portion of the evolution and the frequency of the mode is
found by a least squares fit to the linear portion of the phase
evolution. 


The early nonlinear evolution of unstable barlike modes in polytropes
is dominated by angular momentum transport which arises from 
a gravitational self-interaction torque
between a central bar and {\it ejected} spiral arms 
(Imamura, Picket, \& Durisen 2001).
Imamura, Pickett, \& Durisen explicitly demonstrated
that this interaction torque limited the
growth of bar-like modes in polytropes using 
a quasi-linear (QL) formulation for the torque. 
The QL self-interaction torque on
the annular cylinder between $\varpi$ 
and $\varpi+{\rm d}\varpi$ in the toroid is given by
\begin{equation} 
\Upsilon = m \pi \int|\delta\rho||\delta\Phi|\;
{\rm sin}(\phi_{\rho}-\phi_{\Phi}){\rm dz}
\end{equation} 
where $\delta\rho$ and $\delta\Phi$ are the perturbed
density and gravitational potential, and
$\phi_{\rho}$ and $\phi_{\Phi}$ are the phases of the
density and gravitational perturbations (Imamura, Pickett,
\& Durisen 2001). The QL torque determines the size
and properties of the bars which form and
how the angular momentum 
$J$ is redistributed by the nonlinear development of the
bar mode instability in polytropes,
features which allowed Imamura, Durisen, \& Pickett (2001) 
to predict the
initial development of the $m$ = 2 barlike mode in the nonlinear
regime based only on the linear eigenfunctions.
We apply the QL theory
to toroids.


The nonlinear development of nonaxisymmetric modes in selected 
toroids is followed using a version of the second-order,
three-dimensional (3D) hydrodynamics computer code
described in Pickett, Durisen, \& Davis (1996, hereafter PDD). 
The numerical hydrodynamics code is explicit and
second-order in space and time. The PDD code is extended to include
an energy equation using artifical viscosity to handle shocks. 
In our version of the code, the hydrodynamics equations 
(Equations [1]-[3]) are solved in conservative form,
\begin{equation}
\frac {\partial\rho}{\partial t} = -\nabla\cdot{\rho{\bf v}},
\end{equation}
\begin{equation}
\frac {\partial\rho{\bf v}}{\partial t} = -\nabla\cdot{\rho{\bf vv}} -\nabla 
(P+P_Q) -\rho\nabla\Phi_g,
\end{equation}
and
\begin{equation}
\frac {\partial\epsilon^{1/\gamma}}{\partial t} 
= -\nabla\cdot{\epsilon^{1/\gamma}{\bf v}}
+ \frac {\epsilon^{1/\gamma-1}}{\gamma}\Gamma_Q
\end{equation}
where, $P_Q$ and $\Gamma_Q$ are the
artificial viscosity (AV) dissipative terms inserted
to handle shock formation (see Norman \& Hawley 1985). 
$P_Q$ and $\Gamma_Q$ are given by  
\begin{equation}
P_Q = \rho \left(Q_{\varpi\varpi}^2+Q_{\phi\phi}^2+Q_{zz}^2\right),
\end{equation}
and
\begin{equation}
\Gamma_Q = - \rho \left( Q_{\varpi\varpi} \frac {\partial v_{\varpi}} 
{\partial\varpi}
+Q_{\phi\phi} \frac{1}{\varpi} \frac {\partial v_{\phi}} {\partial\phi}
+Q_{zzi}\frac{\partial v_{z}}{\partial z}\right),
\end{equation}
where 
\begin{equation}
Q_{jj} = \cases {\begin{array}{ll}
C_{Q}(\Delta v_j)^2 & \Delta v_j \leq 0 \cr
0 & \Delta v_j  > 0 \cr,
\end{array}}
\end{equation}
$C_{Q}$ is a constant of order unity and 
$\Delta v_j$ are set to zero in rarefactions so that 
artificial viscosity acts only to dissipate energy.
The hydrodynamic equations are solved on a uniform 
cylindrical grid with second-order transport in space and time.
The numerical grid is again chosen to 
match that used in the equilibrium calculations and the
equilibrium models described in $\S2.1$ can be input directly 
into the nonlinear hydrodynamics computer code. Boundary
conditions on the grid are outflow at the largest $r$                  
and $z$, and reflection through the equatorial plane. The
gravitational potential $\Phi_g$ is found by a direct solution of 
the Poisson equation. Boundary conditions at large distance are  
found by an expansion of $\Phi_g$ on the boundary of the grid in 
spherical harmonics (up to $\ell$ = $m$ = 10). As the bar mode
develops and the model 
expands, the grid is enlarged radially and vertically when
necessary. The maximum grid size used is (384,128,128). 


\subsection{ Energy and Work Integrals }


The perturbed energy equation is given by 
\begin{equation}
\frac{\partial\left<{\cal E}\right>}{\partial t} =
\sigma_R + \sigma_{\Pi} + \sigma_{\Phi},
\end{equation}
where $\left<{\cal E}\right>$ is the mode energy, the 
sum of the kinetic energy and acoustic energy of the 
disturbance
\begin{equation}
\left<{\cal E}\right> =
0.5\;\rho_{\circ} 
\left<\delta v\cdot\delta v\right>
+ 0.5\gamma(P_{\circ}/\rho_{\circ}^2)\left<\delta\rho\right>^2.
\end{equation}
Here, quantities enclosed in brackets, 
$\left<{\rm X}\right>$, indicate averages of the quantity X over $\phi$, and
$\sigma_R$, $\sigma_{\Phi}$, and $\sigma_{\Pi}$ are defined as follows.
$\sigma_R$ is the Reynolds stress given by
\begin{equation}
\sigma_{R} =
-\rho_{\circ}\varpi(\frac{\partial\Omega}{\partial\varpi})\left<\delta 
v_{\varpi}\delta v_{\phi}\right>.
\end{equation} 
The Reynolds stress $\sigma_R$ is a measure of the energy 
in the shear flow of the equilibrium model transfered to the 
disturbance. The stress $\sigma_{\Pi}$
measures the redistribution of energy with in the toroid 
by acoustic wave fluxes carried by the disturbance. $\sigma_{\Pi}$
is given by
\begin{equation}
\sigma_{\Pi} =
-\nabla\cdot\left<\delta{\rm P}\delta{\bf v}\right>.
\end{equation}
$\sigma_{\Phi}$ is a measure of the work performed 
on the disturbance by gravitational forces. $\sigma_{\Phi}$ is given by
\begin{equation}
\sigma_{\Phi} =
-\rho_{\circ}\left<\delta 
{\bf v}\cdot\nabla(\delta\Phi_d)\right>.
\end{equation}
Although the work 
integrals are second-order, they are expressed using only 
the linear eigenfunctions. The energy equation shows regions where 
instability is driven and is damped in the toroid. For normal
modes, the growth time in terms of
work integrls is given as
\begin{equation}
\tau = 2\frac{\int\left<{\cal E}\right>d^3x} {\int(\sigma_R+\sigma_{\Phi})d^3x}.
\end{equation}
The acoustic wave flux term $\sigma_{\Pi}$ does not appear in 
the expression for $\tau$
because it goes to zero at the surface and so 
integrates to zero over the volume of the toroid.


\section{ STABILITY RESULTS }


\subsection{ Linear Regime }


\subsubsection{ Eigenvalues }


The instability regions in ($T/|W|$,$R_-/R_+$)-space are mapped in this
section. To facilitate comparison of our results with previous workers, 
the oscillation frequencies $\omega_{m}$ and growth 
rates $\tau_{m}^{-1}$ are recast in the dimensionless forms 
\begin{equation}
y_1(m)= \frac {\left|\omega_{m}\right|} {\Omega_{max}}-m,\;\;{\rm and}\;\; 
y_2(m)= \frac {1} {\tau_{m}\Omega_{max}}
\end{equation}
(Kojima 1986). If $y_1(m)$ $<$ 0, corotation falls oustide the density 
maximum and if $y_1(m)$ $>$ 0, corotation falls inside the density 
maximum. The eigenvalues for the
$m$ = 2, 3, and 4 modes for $n$ = 3/2 and $q$ = 1, 3/2, and 2 toroids 
and the eigenvalues for the $m$ = 1 mode for $n$ = 3/2, $q$ = 3/2
toroids are presented in Figure 1.


Take the ($n,q$) = (3/2,3/2) toroids as representative of 
toroid behavior. Instability sets in at $T/|W|$ $\approx$ 0.163
through an $m$ = 2 mode.
The growth rate, $y_2(2)$, increases strongly 
after threshold is passed. Its rate of increase
slows beyond $T/|W| \sim$ 0.18, peaking around $T/|W| \approx$
0.22-0.23. After $T/|W|$ = 0.22-0.23, the
growth rate declines until by $T/|W|$
$\approx$ 0.25-0.26, the $m$ = 2 mode stabilizes.
The mode's oscillation frequencies are roughly 
constant, $y_1(2)$ $\sim$ -1.1 near threshold, increasing slowly to
-0.91 by $T/|W|$ $\sim$ 0.26. Corotation for the 
$m$ = 2 modes falls well 
outside $\rho_{max}$, $R_{co}$ $\sim$ 1.5-1.6 $R_{max}$. 
$R_{co}$ approaches $R_{+}$ as the $m$ = 2 mode stabilizes.
Although instability sets in through and is dominated by the 
$m$ = 2 barlike mode, one-armed spirals ($m$ = 1 modes) 
are also unstable in toroids over this $T/|W|$ range. They do not 
dominate the $m$ = 2 mode but are clearly unstable for $T/|W|$ 
$\buildrel>\over\sim$ 0.225 where their growth times are $\tau_1$ 
$\buildrel<\over\sim$ 1.5 $\tau_{max}$. They may become unstable earlier
but because their observed growth times are long, it is not clear 
if this early growth is physical or numerical in origin. We find $\tau_1$ 
$\sim$ 10 $\tau_{max}$ for the $T/|W|$ $\sim$ 0.192 toroid.
The $m$ = 1 modes all have
fast oscillation frequencies, $y_1(1)$ $\sim$ 1 with 
corotation between the inner edge of
the disk and the density maximum of the toroid. They may be
the inner edge $I$ modes discussed in Christodoulou (1994). The $m$ = 1 
modes grow faster than $m$ = 2
modes in toroids with $T/|W|$ $\buildrel>\over\sim$ 
0.25. However, $m$ = 1 modes are not expected to play large
roles in the evolution of toroids. Around the time they start to
dominate the $m$ = 2 mode, $J$ modes appear.
The $J$ mode instability thresholds are
$T/|W|$ $\sim$  0.25-0.26,
0.23-0.25 and 0.25-0.26 
for $m$ = 2, 3, and 4, respectively.
The $m$ = 3 $J$ mode 
appears first and dominates the $m$ = 1 mode by $T/|W|$
$\sim$ 0.26. The $m$ = 3 $J$ mode gives way to the $m$ = 4 $J$ mode 
for $T/|W|$ $>$ 0.28. 
The $J$ modes all have 
faster pattern frequencies than do the 
previously discussed $m$ = 2 modes. For
the $m$ = 2 $J$ mode, the oscillation 
frequency is roughly twice that of the lower $T/|W|$ $m$ = 
2 modes. All of the $J$ modes 
have $y_1(m)$ $\sim$ -0.2 to 0.2 so that 
their corotation radii fall 
close to the location of their density maxima, 
$R_{co}/R_{max}$ $\approx$ 0.95-1.1. 


The general depedence of the stability limits on $T/|W|$ and
the location of corotation suggests that the lower $T/|W|$ modes
are $I$ modes, although there is some uncertainty about this
identification for the $m$ = 2 modes in the $T/|W|$ $<$ 0.19 toroids
(as we note in the next section).
The general stability properties of these $I$ and $J$ modes are not 
qualitatively affected by $q$. The effect of $q$ is to 
change the quantitative locations of the instability
thresholds in $T/|W|$-space. We find that instability
sets in earliest for $q$ = 3/2 toroids, then for $q$ = 2
toroids, and lastly for $q$ = 1 toroids (see Figure 1).
For both $I$ and $J$ modes, the spread in $T/|W|$ for the
instability thresholds is small, $<$ 0.02 for $q$ = 1 to 2. 


The onset of instability in toroids at low $T/|W|$ ($T/|W|$ $<$ 0.27,
the commonly quoted stability limit for barlike modes in polytropes)
was predicted by Hachisu, Tohline, \& Eriguchi (1987)
and verified by Hachisu \& Tohline (1990). Hachisu, Tohline, 
\& Eriguchi argued that gravitationally driven instabilities 
would set in when the ratio of the local free-fall time, 
$\tau$ = $\sqrt{3\pi/32G\rho}$, and 
the sound crossing time $\tau_s$ = $2\pi\varpi/c_s$, where $c_s$
is the local sound speed, fell below 0.08-0.11. Based on this
argument, they predicted that toroids would be unstable for $T/|W|$ 
$\buildrel>\over\sim$ 0.184, 0.151, and 0.142 for 
$q$ = 1, 3/2, and 2. The prediction
for $q$ = 1 is close to that which we found; $q$ = 1 toroids 
are unstable for $T/|W|$ $\buildrel>\over\sim$ 0.18. 
The other predictions are reasonable but less accurate. 
We find limits of $T/|W|$ $\sim$ 0.163 and 
0.17 for $q$ = 3/2 and 2 toroids. Further, 
Hachisu, Tohline, \& Eriguchi argued that the 
effects of the Coriolis force would shift the
stability limits to higher $T/|W|$. We do not see this 
effect. Coriolis forces weaken with increasing $q$, but 
we see the largest differences from the predictions of 
Hachisu, Tohline, \& Eriguchi (1987) for the $q$ = 2 toroids, 
the toroids with the largest $q$ we considered.


\subsubsection{ Eigenfunctions }


The geometric properties of the eigenmodes are 
discussed in this section with particular attention paid to the 
differences in appearance of $I$ and $J$ modes.  
We first look at how the
$m$ = 2 mode eigenfunctions 
change with increasing $T/|W|$ for $(n,q)$ = (3/2,3/2) toroids.
The amplitude of the relative density eigenfunction
$|\delta\rho/\rho_{\circ}|$ in the 
equatorial plane, the azimuthally averaged gravitational self-interaction
torques $\Upsilon$, the constant phase locus for 
the density eigenfunction $\delta\rho$
in the equatorial plane, and the perturbed density and velocity 
fields in the equatorial plane in the frame of the bar 
for a series of $m$ = 2 barlike modes 
are presented in Figure 2. In Figures 3 and 4, we show the perturbed
energies and stresses, the angular momentum, and self-interaction torque
$\gamma$ for the $m$ = 2 modes shown in Figure 2.
In addition to the relative density eigenfunction, 
we also show the eigenfunction 
\begin{equation}
{\cal W} = 
\gamma\left[\frac{P_{\circ}}{\rho_{\circ}^2}\right]\delta\rho+\delta\Phi_g
\end{equation}
({\it e.g.,} see Papaloizou \& Pringle 1982). Both $|\delta\rho|/\rho_{\circ}$
and ${\cal W}$ are normalized to their peak values for presentation in the
figures.
${\cal W}$ is a useful representation near the instability threshold as 
can be seen from the dispersion relation for an eigenmode with time 
and azimuthal dependence given by e$^{i(m\phi-\omega t)}$,
\begin{equation}
\delta\rho = 
\frac {\partial}{\partial\varpi}
\left[\frac{\varpi\rho_{\circ}}{D}(
{\sigma}\frac{\partial}{\partial\varpi}
-\frac{2m\Omega}{\varpi}){\cal W}\right]
+\frac{\partial}{\partial z}
\left[\frac{\varpi\rho_{\circ}}{\sigma}
\frac{\partial}{\partial z}{\cal W}\right]
+\left[\frac{m\rho_{\circ}\kappa^2}{2\Omega D}\frac{\partial}{\partial\varpi}
-\frac{m^2\rho_{\circ}\sigma}{\varpi D}
+\frac{\varpi\rho_{\circ}\sigma}{c_s^2}\right]{\cal W}.
\end{equation}
Here, $\sigma = \omega+m\Omega$, $D = \sigma^2-\kappa^2$, and 
$\kappa$, the epicyclic frequency, is given by
\begin{equation}
\kappa^2 =
2(2-q)\Omega^2
\end{equation}
for power law $\Omega(\varpi)$. Corotation falls where the real part of 
$\sigma$, $\sigma_R$ is equal to zero. The 
inner and outer Lindblad resonances fall where $\sigma_R$ $\pm$ 
$\kappa$ = 0. For power law $\Omega(\varpi)$, the Lindblad resonances 
and the location of corotation are related as
\begin{equation}
\frac {R_{LR}} {R_{co}} = \left(\frac{m}{m\pm\sqrt{4-2q}}\right)^{-1/q}.
\end{equation}
For toroids with constant specific angular
momentum ($q$ = 2), $\kappa^2$ = 0, $D$ = 
$\sigma^2$, and corotation and the Lindblad resonances are degenerate.
We show these and other critical radii for the $m$ = 2
modes of $(n,q)$ = (3/2,3/2) toroids in Figure 5.
Detailed properties of the $m$ = 2 modes are given in Table 2.


First, consider the $m$ = 2 modes 
over the range 0.16 $<$ $T/|W|$ $<$ 0.26.
Near $T/|W|$ = 0.16, the relative density 
eigenfunction $|\delta\rho/\rho_{\circ}|$ shows a 
phase shift of $\delta$ $\sim$ $\pi/3$ at the location where 
$|\delta\rho/\rho_{\circ}|$ exhibits a sharp minimum, $R_{\phi}$,
not where ${\cal W}$ passes through 
minimum, $R_{\cal W}$. The perturbed toroid
appears as a coherent bar which undergoes 
a sharp phase advance at $R_{\phi}$
outside of which the figure of the fluid
settles into another coherent bar which extends to its outer edge. 
The phase swing at $R_{\phi}$ falls inside $R_{max}$, $R_{co}$,
$R_{\Upsilon}$, $R_{\lambda}$ and both Lindblad resonances. The 
inner Lindblad resonance $R_{ILR}$ does, however, fall near $R_{max}$.
The $m$ = 2 mode changes character near $T/|W|$ = 0.19, the
model where $R_{OLR}$ first moves outside the outer edge of the toroid.
For the eigenfunctions in $T/|W|$ $\buildrel>\over\sim$ 0.19 toroids,
the phase change switches from
leading to trailing with figures now nicely
modeled as coherent central bars with rapid $\delta$
$\sim$ -$\pi/2$ phase lags at $R_{\phi}$ followed
by nearly coherent bars to their outer edges. 
For these toroids, $R_{\phi}$ continues to fall
inside $R_{\Upsilon}$, $R_{co}$, and $R_{\lambda}$, and
$R_{ILR}$ continues to fall near $R_{max}$. Suggestively, 
the phase change at $R_{\phi}$ is roughly where 
the Toomre Q-parameter first drops below 1 as one moves 
outward from the inner edge of the toroid. Corotation 
does not coincide with any particular location for these toroids and
approaches the outer edge of the toroid for $T/|W|$ $\sim$
0.26, where the $m$ = 2 mode stabilizes. For these toroids, the 
evanescent wave region, the area where $Q$ $<$ 1 is very broad 
and does not always tightly bracket corotation as needed for 
scenarios which rely on the {\it tunneling} 
of edge waves which carry negative and positive energies
through corotation as the amplification mechanism for instability
({\it e.g.}, Nakagawa \& Sekiya 1994).


The eigenfunction ${\cal W}$
varies smoothly through $R_{\phi}$ with a minimum at
$R_{\cal W}$. The minimum location, $R_{\cal W}$,
falls near where $\Upsilon$ changes sign which, for 
toroids near threshold is just inside $R_{co}$.
The location of the vortensity extremum,
$R_{\lambda}$, does not coincide with either $R_{co}$ or 
$R_{ILR}$ and the $m$ = 2 modes are apparently 
not related to the vortensity modes of self-gravitating thin 
disks studied by Papaloizou \& Savonije (1991). A necessary 
condition for the vortensity modes is that $R_{co}$ = $R_{\lambda}$ 
at threshold to make corotation singlarities ineffective.
As $T/|W|$ increases and the toroid moves beyond
threshold the sharp minimum in ${\cal W}$ smooths out but
still appears to be associated with $R_{\Upsilon}$.
The minimum smooths because resonances are
weakened when instability becomes stronger;
by $T/|W|$ = 0.19, $\tau_2/P_2$ $\sim$ 0.25. The sharp 
minimum in $|\delta\rho/\rho_{\circ}|$, and hence the rapid phase 
change persists beyond threshold and is apparently not related to 
such resonances.


Representative $m$ = 1 eigenfunctions for 
$T/|W|$ = 0.192 and 0.253 toroids are shown in Figure 6.
The appearance of ${\cal W}$ and
$\delta\rho/\rho_{\circ}$ are similar in shape to the 
$m$ = 2 modes seen in 
$T/|W|$ $>$ 0.19 toroids. 
The  $\delta\rho/\rho_{\circ}$ show sharp minima outside 
the location of $\rho_{max}$, $R_{\phi}$ = 1.17 $R_{max}$,
where the eigenmodes undergo rapid phase swings
of $\delta$ $\sim$ -$\pi$. The 0.192 ${\cal W}$ shows a 
sharp minimum at 0.66 $R_{max}$, just outside corotation 
$R_{co}$ = 0.63 $R_{max}$. The 
0.253 ${\cal W}$ shows a shallower minimum at
0.68 $R_{max}$ which again falls near, but slightly outside cortation,
$R_{co}$ = 0.65 $R_{max}$. Based on the appearance of the mode and 
location of these and other key radii, {\it e.g.,} the 
vortensity extrema fall near 
$R_{max}$, it appears that the 
$m$ = 1 modes are related to $I$ modes (see 
Christodoulou 1994) rather than to the {\it eccentric} $m$ = 1 
modes of Adams, Ruden, \& Shu (1989) aa expected because the
{\it eccentric} mode of Adams, Ruden, \& Shu relies on a
central star moving out-of-phase with the motion of the disk 
center-of-mass. 


$J$ modes emerge around $T/|W|$ $\sim$ 0.25 and the barlike 
eigenfunctions again show distinct changes. The eigenfunctions 
still exhibit phase lags at $R_{\phi}$ of
$\delta$ $\sim$ -$\pi/2$, but the transitions are more rounded than 
for the lower $T/|W|$ modes and the perturbations lead to objects 
composed of coherent central bars which trail smoothly into
pairs of spiral arms.
The dimensionless oscillation frequencies of $J$ modes are 
also faster than the lower $T/|W|$ modes and
corotation falls closer to $R_{max}$. 
Further, the location of the phase swing,
$R_{\phi}$, coincides with $R_{\Upsilon}$, $R_{co}$, 
$R_{{\cal W}}$, and $R_{\lambda}$
for these high $T/|W|$ $I$ modes and
Lindblad resonances re-appear. The Lindblad resonances are roughly
symmetrically placed about $R_{co}$ for the $T/|W|$ = 0.273 toroid.
As $T/|W|$ increases, the Lindblad resonances move away
from $R_{co}$ and both fall outside the toroid for $T/|W|$ $>$ 
0.30. The Q-barriers, the locations where Q passes through 1, now 
tightly bracket corotation. The evanescent wave
zone around $R_{co}$ is narrow for these $J$ modes 
increasing the probability of wave {\it tunneling} through 
corotation and amplification may occur through 
{\it super}-reflection of waves at corotation. Further, 
after examination of the azimuthally-averaged angular momentum
of the $m$ = 2 modes for these toroids, we see that the
angular momentum carried by the eigenmodes
changes sign at $R_{\phi}$ (see Figure 4). The change in sign
of the angular momentum is consistent with the suggestion
that instability is driven by the coupling of waves across 
the corotation singularity (Nakagawa \& Sekiya 1992). This behavior
is not seen in the lower $T/|W|$ $I$ modes. Further, for the $T/|W|$ $<$
0.19 $m$ = 2 modes, the angular momentum perturbation is positive at 
both the inner and outer regions of the toroid; it is negative only 
where the self-interaction torque $\gamma$ changes sign.


\subsubsection{ $\left<{\cal E}\right>$ and the Work Integrals }


The azimuthally averaged perturbed energy, $\left<{\cal E}\right>$,
composed of the perturbed
kinetic energy ${\cal E}_k$ and the perturbed
enthalpy ${\cal E}_h$, and the
stresses $\sigma_R$, $\sigma_{\Phi}$, and $\sigma_{\Pi}$ 
for the $m$ = 1 and 2 modes shown in Figure 2 are
presented in Figure 3. For toroids with $T/|W|$ $<$ 0.192,
the kinetic energy growth rate 
$2\gamma{\cal E}_k$ and the acoustic wave stress 
$\sigma_{\Pi}$ peak at the same radius while the 
the enthalphy growth rate 
$2\gamma{\cal E}_h$ and the Reynolds stres  
$\sigma_{R}$ peak at the same radius.
The $2\gamma{\cal E}_k$ peak falls closer to the inner edge of 
the toroid than does that of 
$2\gamma{\cal E}_h$.  The enthalphy growth rate $2\gamma{\cal E}_h$
shows a second smaller peak outside $R_{max}$ which grows in size and 
importance as $T/|W|$ increases. For these low $T/|W|$ toroids, 
the strongest driving occurs near the inner edge of the toroid.
Notice that the gravitational stress $\sigma_{\Phi}$
peaks near $R_{max}$ and that $\sigma_{\Pi}$ 
is negative where $\sigma_{\Phi}$
is the largest. The region where $\sigma_{\Pi}$ is $<$ 0 is  
where energy is transferred from the perturbation to 
the background flow
consistent with the fact that the central region of the 
toroid has $Q$ $<$ 1 and is an evanescent wave propagation region.


As $T/|W|$ increases beyond $/T/|W|$ $\approx$ 0.19, the secondary 
peak in $2\gamma{\cal E}_k$ increases in size relative to the main inner
peak, and $\sigma_R$ and $\sigma_{\Phi}$ increase in importance
compared to $\sigma_{\Pi}$. The $\sigma_{\Pi}$ and $\sigma_R$ 
continue to peak close to the inner edge of the toroid, however, 
the peak in $\sigma_{\Pi}$ moves inside that of 
$\sigma_R$.  It is also now more obvious that 
$\sigma_{\Phi}$ peaks where $\left( {\cal E}_k+{\cal E}_h\right)$ 
reaches a local maximum, which falls 
near $\rho_{max}$. The large region where $\sigma_{\Pi}$ $<$ 0 is 
again the region of the toroid where $Q$ $<$ 1. This effect
is much more pronounced in the $T/|W|$ $>$ 0.19 toroids than in 
the $T/|W|$ $<$ 0.19 toroids. 


For the $T/|W|$ $>$ 0.25-0.26 toroids, where $J$ modes dominate, 
the enthalpy term $2\gamma{\cal E}_h$ shows two strong maxima.
Near threshold, the maximum near the inner edge of the toroid dominates
but as $T/|W|$ increases the two maxima become more comparable 
in size. This is consistent with the fact that the
$Q$ $<$ 1 region narrows around density maximum for 
$T/|W|$ $>$ 0.25-0.26. The gravitational stress
$\sigma_{\Phi}$ peaks at the minimum in $2\gamma{\cal E}_h$,
where $\sigma_{\Pi}$ is the most negative. This is also near the
density maximum. This suggests that gravity puts energy into the 
disturbance at the maximum density and drives instability. The Reynolds stress 
still peaks where $2\gamma{\cal E}_K$ peaks, but its amplitude is much smaller
in comparison to $\sigma_{\Phi}$ than for the low $T/|W|$ toroids.


We show ${\cal R}$ = $\int\sigma_{\Phi}{\rm d}^3x/\int\sigma_{R}{\rm d}^3x$ 
in Figure 7. ${\cal R}$ $\sim$ 2 for toroids with $T/|W|$ $<$
0.25-0.26 toroids. For toroids with larger $T/|W|$, 
${\cal R}$ is larger, ${\cal R}$ $>$ 4 for toroids
with $T/|W|$ $>$ 0.26. This behavior is consistent
with the suggestion that the high $T/|W|$ modes are 
the gravity driven $J$ modes, and the low
$T/|W|$ modes are the $I$ modes, modes
{\it intermediate} in nature between the shear-driven
P modes and the gravity driven $J$ modes.


\subsection{ Nonlinear Regime }


We run simulations of $(n,q)$ = (3/2,3/2) toroids 
with $T/|W|$ = 0.192, 0.253, and 0.272 into the 
nonlinear regime to study the long-term evolution of 
toroids dominated by $I$ and $J$ modes. The $T/|W|$  
= 0.192 toroid is linearly unstable to barlike $m$ = 2 modes 
and likely stable to $m$ = 1 perturbations.
We cannot unambiguously assert that $m$ = 1 modes
are stable, however, because of numerical uncertainty.
$m$ = 1 modes do become more unstable with increasing $T/|W|$. 
For the $T/|W|$ = 0.253 toroid, the $m$ = 1 mode
dominates the $m$ = 2 $I$ mode in the linear regime. We show
the $m$ = 1, 2, and 3 linear density eigenmodes for the $T/|w|$
= 0.253 toroid in Figure xxx.
For the $T/|W|$ = 0.272 toroid, the $m$ = 1 mode is still 
unstable but the $m$ =  3 $J$ mode dominates it in the
linear regime. We show the $m$ = 1, 2, 3, and 4 density 
eigenmodes in Figure xxx for the $T/|W|$ = 0.272 toroid. 
We describe each nonlinear simulation in 
the following sections.


\subsubsection{ $T/|W|$ = 0.192 Toroid }


The low $m$ Fourier amplitudes ${\cal A}_m$,
cos$(\phi_{\rho})$ for the $m$ = 2 mode, and
$T/|W|$ are shown in
Figure xxx, where the ${\cal A}_m$ are
the integrals of the Fourier amplitudes over the volume 
of the toroid (see Pickett, Durisen, \& Davis 1996),
and $\phi_2$ is the phase for the $m$ = 2 density disturbance. 
All low-$m$ nonaxisymmetric modes participate in an early 
{\it settling} phase during which their amplitudes grow to 
$\sim$ $10^{-5}$ where they level off. During this phase,
the toroid initially expands inward and outward, 
but quickly re-settles in an equilibrium nearly identical
to that of the original model. This initial transient phase is
probably caused by slight adjustments made to the equilibrium model
because of meridional motions. The equilibrium assumes rotation on
cylinders. We find that meridional motions arise right after we
start the simulation. The kinetic energy in the meridional motions 
is $<$ 0.01 \% than that contained in the rotation and thermal motions,
small not zero.


The $m$ = 2 mode organizes soon after the initial settling phase
and starts to grow
by $t$ $\sim$ $0.5 \tau_{max}$. The $m$ = 2 mode settles into
steady global growth by $t$ = 2-3 $\tau_{max}$.
The eigenvalues for the growing barlike mode are nearly
identical to those of the linear barlike mode.
We find $P_{2}$ = 2.0 $\tau_{max}$ and
$\tau_{2}$ $\sim$ 0.57 $\tau_{max}$ for the 
nonlinear simulation. The eigenvalues for the
linear simulation are $P_{2}$ = 2.09 $\tau_{max}$ and
$\tau_{2}$ = 0.515 $\tau_{max}$. All other modes are found
to be stable.\footnotemark
\footnotetext
{
The $m$ = 1 mode is at best marginally unstable. It has fast
pattern period =  0.49 $\tau_{max}$ and long growth time $\tau_1$ 
$\sim$ 10 $\tau_{max}$.
}
The density eigenfunction $|\delta\rho/\rho_{\circ}|$ 
strongly resembles
the linear eigenfunction (compare Figures xxx and xxx) at
early times. $|\delta\rho/\rho_{\circ}|$ 
shows maxima at the inner and outer edges of the toroid 
with a sharp minimum at $R_{\phi}$ = 0.7 $R_{max}$ 
which marks the location of the $\pi/2$ phase shift in 
the figure of the bar. The linear eigenfunction shows the 
sharp minimum at $R_{\phi}$ = 0.74 $R_{max}$.
Corotation falls at $\sim$ 1.6 $R_{max}$;
$R_{co}$ = 1.63 $R_{max}$ for the linear barlike mode. 


During the earliest part of the nonlinear simulation, 
the mass and angular momentum do not consistently
flow inward or outward as the mode organizes. They alternately 
flow across $R$ $\sim$ 1.3 $R_{max}$, roughly 
the radius where $\Upsilon$ changes sign for the linear
eigenfunction, $R_{\Upsilon}$ =  1.37 $R_{max}$. 
In Figure xxx, we show
the mass $\Delta m_c$ and angular momentum $\Delta J$ contained 
in cylindrical annuli between $\varpi$ and $\varpi+{\rm d}\varpi$, 
and the radial mass flux $\rho U$, and the radial angular
momentum flux $|{\bf {A}}|U$ as functions of $\varpi$ where ${\bf U}$ 
is the the radial velocity and ${\bf A}$ is the angular momentum 
density. Then, as
the barlike mode grows in amplitude, the mass and angular
momentum flows settle in steady patterns in which they
flow inward and outward across 0.87 $R_{max}$, near the radius 
where the rapid phase swing occurs and 
the angular momentum of the linear eigenfunction
changes sign. The torque, defined as the 
change in the local angular momentum
$\Upsilon_{nl}$,
has a form similar to that of the 
linear gravitationl self-interaction torque 
$\Upsilon$. As the bar
grows in amplitude, the innermost
and outermost regions of the toroid gain angular momentum 
while the central region sheds angular momentum. The 
outer sign change occurs near 1.37 $R_{max}$. 
The gravitational self-interaction torque as discussed in 
$\S3.2$, appears to control the
angular momentum redistribution in the nonlinear simulation. 


Nonlinear effects start to alter the toroid structure
when the bar reaches amplitude ${\cal A}_2$ 
$\sim$ 0.1 at around $t$ = 4.1 $\tau_{max}$. 
At this time even $m$ $>$ 2 disturbances
which have the same pattern 
periods as the $m$ = 2 mode does, start to grow.
The $m$ $>$ 2 disturbances are most likely the 
result of a nonsinusoidal $m$ = 2 waveform. As the bar 
strengthens, the inner radius of the toroid moves inward
filling in the center of the toroid. 
Once this happens, the object, strictly speaking,
is no longer a toroid. Not surprisingly, 
we find that the nonlinear eigenfunction differs
strongly from the linear eigenfunction when this occurs
(see Figures xxx and xxx).
The nonlinear eigenfunction now generally increases in 
amplitude from the rotation axis outward,
peaking at the equatorial radius. However, a sharp minimum still remains
in the central region of the toroid. The minimum falls inside the radius
where the minimum occurs in the linear eigenfunction and
a pronounced off-center density maximum remains. Most
of the mass is still contained 
between radii 0.44 $R_{max}$ and 2.4 $R_{max}$, close to 
the original inner and outer radii of the toroid.
The $m$ = 2 disturbance saturates at amplitude
${\cal A}_2$ $\sim$ 0.96 around $t$ $\sim$ 5.6 $\tau_{max}$.
Near saturation, a small disk with much less than 1 \% of the 
mass and angular momentum of the toroid is ejected.


Imamura {\it et al.} (2001) showed that saturation of barlike 
$f$ modes in polytropes was caused by the gravitational 
self-interaction torque $\Upsilon$ described in $\S3.2$. 
The self-interaction torque led to the ejection of significant 
amounts of mass from the
polytrope, the remaining mass then losing angular momentum to
the ejected disk until it reached a configuration marginally 
stable to the barlike instability. If this is also true for toroids, 
then we would expect that the mass inside $R_{\Upsilon}$ would
shed $J$ until it reached the $J$ of a toroid (or other oblate 
figure) marginally 
stable to the barlike instabiity. Given this and the QL torque
given by Equation (6), Imamura {\it et al.} would argue that 
the barlike mode would grow until the angular momentum transport
time
\begin{equation}
\tau_j = \left| \frac{j_b-j_{lim}}{\Upsilon} \right|,
\end{equation}
where $j_b$ is the angular momentum contained in the region 
inside $R_{\Upsilon}$ and $j_{lim}$ is the angular 
momentum of the marginally stable axisymmetric object
with $M_b$.
Using the QL torque and the peak ${\cal A}_2$ from the 
nonlinear simulation, the bar mass 
$M_b$ = 0.629 and angular momentum $j_b$ = 0.395 from the 
quasi-linear simulation, 
and the limiting angular momentum for the marginally stable
toroid from the linear analysis $j_{lim}$ = 0.332, we find
$\tau_j$ = 1.2 $\tau_{max}$. The bar mode linear growth time
$\sim$ 0.57 $\tau_{max}$ and the self-interaction torque
based on the linear eigenfunction is too weak to saturate the
growth of the barlike mode if the target state is the 
marginally stable toroid. Because saturation 
occurs after the toroid fills-in, however, there are distinct 
differences between the nonlinear and linear eigenfunctions as 
noted earlier. Near saturation, although the nonlinear
$\Upsilon_{nl}$ (see Figure xxx) resembles 
the linear $\Upsilon$ (see Figure xxx), there are quantitative
differences. If we use 
the nonlinear torque $\Upsilon_{nl}$, we do find
$\tau_j$ = 1.2 $\tau_{max}$ similar to the prediction of the QL 
theory. The gravitational self-interaction torque, 
at least initially, tries to drive the toroid to the marginally 
stable axisymmetric toroid state.


After saturation, 
$\Upsilon_{nl}$ changes character. $R_{\Upsilon}$ remains at 
1.37 $R_{max}$ but the region inside $R_{\Upsilon}$
now gains angular momentum and the region outside $R_{\Upsilon}$ 
loses angular momentum. The bar amplitude
reaches its first minimum, ${\cal A}_2$ $\sim$ 0.47 around
6.9 $\tau_{max}$, after which it rebounds and 
$\Upsilon_{nl}$ reverts to its original form where the region 
inside $R_{\Upsilon}$ $\sim$ 1.37 $R_{max}$ loses angular 
momentum and the region outside $R_{\Upsilon}$ gains angular 
momentum. The barlike mode reaches its second maximum, 
${\cal A}_2$ $\sim$ 0.95, around $t$ $\sim$ 7.5 $\tau_{max}$ 
where $\Upsilon_{nl}$ again pivots about 1.37 $R_{max}$ and 
another small disk is ejected.
The toroid continues to follow this slow damping, 
oscillatory pattern until $t$ $\sim$ 9 $\tau_{max}$ 
when the damping
rate increases, $T/|W|$ starts to 
increase, and the evolution of the central region changes. The evolution
switches from one which resembles  a damping toroidal 
instability to one characterized as a secular change to a
new equilibrium state, a filled-center polytrope equilibrium. 
The central region starts to rotate in the retrograde sense with a 
strong increase in the outward mass flux
and an outward pulse of negative angular momentum. The
mass and angular momentum redistribution permanently alter the 
structure and specific angular momentum distribution of the central
object. Around $t$ $\sim$ 13 $\tau_{max}$, 
$T/|W|$ stops growing and the simulation settles to a 
configuration approaching a distorted, stable, steady star/disk system. 
We stop the simulation after $\sim$ 19 $\tau_{max}$.


The final density structure and velocity field in the frame of the 
bar, and the constant density phase locus for the bar in the equatorial 
plane are shown in Figure xxx. The final specific angular momentum 
distribution and $\Omega$ are compared to the specific 
angular momentum distribution of a Maclaurin spheroid 
and the original $\Omega$ in Figure xxx. The final 
state has ${\cal A}_2$ $\sim$ 0.12 and $T/|W|$ $\sim$ 0.193, well 
beyond the instability threshold for $(n,q)$ = (3/2,3/2) toroids
but below the usually quoted filled-center polytrope threshold
of $T/|W|$ = 0.27.
The central object has a core in near uniform rotation with an 
envelope with power law angular velocity distribution, 
$\Omega$ $\propto$ $\varpi^{-1.1}$. The central object has an
off-center density maximum. The uniformly rotating core contains 
$\sim$ 10-12 \% of the mass with the transition to the uniform 
velocity region falling near $\varpi$ $\sim$ 0.5 $R_{max}$. 
The central object is surrounded by a low mass disk. The transition
region between the central object and the disk where the angular 
velocity increases contains 0.75 \% of the mass and $\sim$ 
1.62 \% of the angular momentum. The disk contains 0.65 \% of the 
mass and $\sim$ 2.7 \% of the angular momentum with 
$\Omega$ $\propto$ $\varpi^{-1.75}$.  The barlike instability 
does not lead to a prompt disruption of the toroid
rather, it leads to the formation of a stable star/disk system.
However, the final state has off-center density maxima with
nonaxisymmetric structure.


\subsubsection{ $T/|W|$ = 0.253 Toroid }


The low $m$ Fourier amplitudes ${\cal A}_m$ and
$T/|W|$ for the evolving toroid are shown in Figure xxx.
The ${\cal A}_m$ quickly grow to 
amplitudes of 
$\sim$ $10^{-5}$ where they then level off. During this time the toroid
structure adjusts but quickly re-settles into an equilibrium nearly 
identical to that of the original model. The
$m$ = 1, 2, and 3 modes then organize and start to grow
after $t$ $\sim$ 0.5 $\tau_{max}$. The modes
settle into steady growth by $t$ $\sim$ 1 $\tau_{max}$. We find
$P_{1}$ = 0.525 $\tau_{max}$, $P_{2}$ = 1.88 $\tau_{max}$, and 
$P_{3}$ = 1.36 $\tau_{max}$ with growth times
$\tau_{1}$ = 0.38 $\tau_{max}$, $\tau_{2}$ $\sim$ 0.45 $\tau_{max}$,
and $\tau_{3}$ = 0.52 $\tau_{max}$ for the nonlinear simulation.
The linear eigenvalues are $P_1$ = 0.52 $\tau_{max}$, $P_2$ =
1.84 $\tau_{max}$, and $P_3$ = 1.19 $\tau_{max}$, and
$\tau_1$ = 0.40 $\tau_{max}$, $\tau_2$ = 
0.73 $\tau_{max}$, and $\tau_3$ = 0.71 $\tau_{max}$. The nonlinear
oscillation frequencies are very close to the linear frequencies,
however, the $m$ = 2 and 3 growth rates are more similar to those 
found for slightly lower $T/|W|$ toroids.
This occurs because the $T/|W|$ = 0.253 toroid is 
in the region where the $m$ = 2 $I$ mode stabilizes and the 
$m$ = 3 $J$ mode first becomes unstable. Consequently, small 
changes in the equilibrium toroid could lead to relatively large 
changes in the growth rates for the $m$ = 2 and 3 modes.

 
At early times, the $m$ = 1 and 2 mode nonlinear eigenfunctions are 
similar in shape to the $m$ = 1 and 2 linear eigenfunctions. We show the 
$m$ = 1 and 2 $|\delta\rho/\rho_{\circ}|$ and their constant phase loci 
in the equatorial plane in Figure xxx. Corotation for the $m$
= 1 and 2 modes falls near $\sim$ 0.8 $R_{max}$.


The $m$ = 1 mode with its faster growth rate eventually dominates 
the $m$ = 2 mode; it overtakes the $m$ = 2 mode around
$t$ $\sim$ 3.3 $\tau_{max}$. The $m$ = 2 mode saturates at $t$ $\sim$
4.1 $\tau_{max}$ when ${\cal A}_2$ $\sim$ 0.076. The $m$ = 1 mode continues
to grow saturating at $t$ $\sim$ 4.5 $\tau_{max}$ when ${\cal A}_1$ $\sim$ 0.xx. 
Nonlinear effects set in before the $m$ = 1 saturation as the toroid starts to 
spread when ${\cal A}_1$ $\sim$ 0.1. After saturation the toroid continues to 
evolve eventually settling into a nearly asixymmetric quasi-steady state 
with ${\cal A}_1$ $\sim$ 0.xx.



\subsubsection{ $T/|W|$ = 0.272 Toroid }


The (3/2,3/2), $T/|W|$ = 0.272 toroid is representative of 
models dominated by $J$ modes rather than $I$ modes.
The ${\cal A}_m$ and $T/|W|$ for 
the simulation are shown in Figure xxx. There 
is an initial {\it settling} phase in which growth in all
modes is seen during 
which the toroid adjusts slightly but returns quickly 
to its original equilibrium structure. The
$m$ = 3 and 4 modes are the first to organize and 
grow. Both modes appear after $\sim$ 0.75 $\tau_{max}$ has
elapsed, settling into steady growth by $t$ = 1 $\tau_{max}$.
The $m$ = 4 mode dominates early, but
the $m$ = 3 mode overtakes it at $t$ = 1.8 $\tau_{max}$. 
In addition, growth is also seen 
in $m$ = 1
starting at $t$ $\sim$ 1.5 $\tau_{max}$, 
\footnotemark
\footnotetext{
In the linear regime, $P_1$ = 0.507 $\tau_{max}$ and $\tau_1$ = 
0.293 $\tau_{max}$. Although $m$ = 1 growth is faster in the linear
regime, the $m$ = 3 and 4 modes organize more quickly in the nonlinear
simulation and dominate the evolution.
}
In this early phase of the nonlinear simulation, 
$P_3$ = 1.10 $\tau_{max}$
and $P_4$ = 1.08 $\tau_{max}$, and $\tau_3$ = 
0.20 $\tau_{max}$, and $\tau_4$ = 0.26 $\tau_{max}$.
The nonlinear eigenvalues agree well with the linear 
results where the $m$ = 3 and 4 modes have
$P_2$ = $1.11\tau_{max}$ 1.10 and
$P_3$ = $1.19\tau_{max}$ 1.10, and
$\tau_{2}$ = $0.205\tau_{max}$ = 0.205 and
$\tau_{3}$ = $0.248\tau_{max}$. The $m$ = 3 and 4
modes are nearly commensurate.
The {\it nonlinear} pattern periods and growth times 
match the linear ones to within 5 \% - 10 \%.
The $m$ = 3 eigenfunction in the equatorial
plane at time $t$ = 1.98 $\tau_{max}$ is shown in 
Figure xxx. At this time ${\cal A}_3$ $\sim$ 0.1. The
disturbance strongly resembles the 
linear eigenfunction (igure xxx) with
maxima at the inner and outer edges of the toroid and
a rounded minimum near the density maximum, $\varpi$ 
$\sim$ 0.65 $R_{max}$. Corotation falls nearly at the 
density maximum, $R_{co}$ $\sim$ 1.1 $R_{max}$ as
characteristic of $J$ modes.


The angular momentum redistriubtion function quickly settles into a 
steady form. $\Upsilon_{nl}$ is similar to the linear $\Upsilon_{l}$ 
in that the outermost region of the toroid gains angular
momentum from the inner region,  however, the innermost
region also gains angular momentum in the nonlinear toroid.
The gravitational self-interaction torque plays a role 
in the evolution of the $m$ = 3 mode as it approaches saturation,
but it may not be the dominant effect as in the $m$ = 2 $I$ mode 
simulation discussed in $\S4.2.1$. 
Nonlinear effects, such as growth in high order modes appear
between $t$ = 1.5 and 2 $\tau_{max}$ when the $m$ = 3 
amplitude has reached ${\cal A}_3$ $\sim$ 0.1. 
Growth in $m$ = 3 saturates around $t$ $\sim$ 2.7 $\tau_{max}$
when ${\cal A}_3$ $\approx$ 0.416 and $T/|W|$ attains 
its first minimum at $T/|W|$ $\sim$ 0.26. By
saturation the toroid has already undergone significant
evolution in that the center of the toroid has filled-in
and the outer radius of the toroid has increased in size 
by more than 50 \%. The $m$ = 3 eigenfunction no longer 
resembles the linear eigenfunction (see Figures xxx and xxx);
the $m$ = 3 disturbance looks more 
like that of a filled-center, star-like polytrope.


The toroid continues to evolve with $T/|W|$ oscillating in size 
but, on average, declining while the outer part of the disk 
continues to expand. The $m$ = 2 mode continues to grow after
the $m$ = 3 mode saturates and starts to decline. 
Around $t$ = 3.3 $\tau_{max}$, the $m$ = 2 mode becomes larger
than the $m$ = 3 mode. The equatorial radius of the toroid (the 
ejected disk's edge) now extends to a distance around a factor 
of four greater than that of the original toroid. The toroid 
then appears to be settling in a new filled-center equilibrium
where an $m$ = 2 disturbance dominates. The $T/|W|$ reaches a 
second maximum at $T/|W|$ = 0.292. $m$ = 2 continues to grow a little 
after this point and reaches maximum ${\cal A}_2$ $\sim$ 0.265 around 
0.05 $\tau{max}$ after the peak in $T/|W|$. The toroid then continues
to {\it oscillate} as $T/|W|$ approaches 0.xxx and ${\cal A}_2$ approaches
0.xxx . During this time, the $m$ = 1 mode is comparable to or
larger than the $m$ = 2 and 3 modes, reaching maximum amplitude 
${\cal A}_1$ $\sim$ 0.163.
The simulation is stopped at $t$ $\sim$ 6 $\tau_{max}$ 
because of poor energy conservation; the toroid has lost 10 \%
of its energy by 6 $\tau_{max}$. 


\section{ DISCUSSION }


\subsection{ Comparison to Previous Works }


Hachisu \& Tohline (1990) and Woodward, Tohline, \& Hachisu (1994) performed
a handful of calculations directly comparable to ours in their study of the 
properties of nonaxisymmetric gravito-rotation modes in
self-gravitating toroids and star/disk systems. Tohline {\it et al.} 
performed simulations of self-gravitating toroids for 
$n$ = 3/2 and $q$ = 3/2, 7/4, or 2 using a fully nonlinear,
three-dimensional, numerical hydrodynamics computer code.
They did not perform an extensive study of
parameter space, however. They modeled only even $m$ modes 
for seven toroid models.
Our results, where they overlap with the linear $I$ mode results given 
in Woodward, Tohline, \& Hachisu (1994), agree very well.\footnotemark
\footnotetext{
The Woodward, Tohline, \& Hachisu calculations were performed using a 
nonlinear hydrodynamics computer code. So, strictly speaking, they were
never in the linear regime. Linear in the context of the nonlinear 
simulations of WTH means the portion of the evolution where exponential
growth was observed. 
}
The $\omega_2$ are within 5 \% and the $\tau_2$ within 5-10 \%
for the $I$ modes with $T/|W|$ $<$ 0.26. The eigenvalues agree less well 
with the one $J$ mode simulation presented in Woodward, Tohline, \& 
Hachisu (1994). The toroid model has $T/|W|$ = 0.27 for which, because 
they modeled only even $m$ modes, Woodward, Tohline, \& Hachisu
missed the fastest growing $J$ mode, the $m$ = 3 mode.


A more extensive study of self-gravitating toroids was performed by 
Andalib, Tohline, \& Christodoulou (1997) who
studied $P$, $I$, and $J$ 
modes in narrow, incompressible, self-gravitating annuli with
constant specific angular momentum and circular meridonal cross-sections. 
The set of  assumptions forms the 
{\it ICT} approximation. Although not directly comparable to our work 
because of these simplifying assumptions, {\it ICT} studies offer an useful 
framework for a discussion of our
study of the stability properties of self-gravitating
toroids. Andalib, Tohline, \& Christodoulou (1997) present
stability conditions using $T/|W|$ as the stability parameter
(see their Figures 6 and 7 for $I$ modes and Figures 10 and 11 for
$J$ modes). (i) Andalib, Tohline, and Christodoulou (1997) found
for $I$ modes that $m$ = 1 modes
were the most unstable, instability setting in for 
$T/|W|$ $>$ 0.175858. $m$ = 2 modes are less unstable,
instability setting in for $T/|W|$ $>$ 0.2620.
(ii) For the $J$ modes, Andalib, Tohline, and Christodoulou found 
that the $m$ = 4 and 5 modes were the most unstable with instability
setting in for $T/|W|$ $>$ 0.3018. Other $m$ $>$ 2
$J$ modes became unstable only for higher $T/|W|$. Our
results are in qualitative agreement with those of 
Andalib, Tohline, \& Christodoulou (1997) despite the fact that the 
toroids used in our study 
populate a regime in parameter space where the ICT assumptions are 
suspect. Toroids with $T/|W|$ $<$ 0.27 have $\epsilon$ $>$
0.3-0.4. For the $I$ modes, we find that $m$ = 2 modes are the most
unstable with instability setting in for $T/|W|$ 
$>$ 0.16-0.18, and that $m$ = 1 $I$ modes are more stable with  
instability setting in at $T/|W|$ 
$\sim$ 0.22. For $J$ modes, we find 
that $m$ = 2, 3, and 4 $J$ modes are unstable
for $T/|W|$ $>$ 0.23-0.26. We thus find that instability sets in at lower
$T/|W|$ for both $I$ and $J$ modes, and that instability sets in at 
higher $m$ for $I$ modes and at lower
$m$ for $J$ modes than found by Andalib, Tohline, \& Christodoulou. 
These differences are probably not physical, however. 
Near instability thresholds, toroids are wide and 
the ICT theory breaks down.  The breakdown 
will affect low $m$ modes more strongly than 
it will high $m$ modes. This occurs because 
low $m$ modes likely have longer wavelengths as shown by studies of 
uniformly rotating incompressible fluids which find that 
nonaxisymmetric modes with $l$ = $|m|$ are the
most unstable for a given model. Here, $l$ (and $m$) are the 
angular quantum numbers for the associated Legendre polynomials.
Consequently, low $m$ modes will have long wavelengths in the
meridional plane as well as in
azimuthal direction. The accuracy of their modeling then 
relies on how well the surface of the toroid is resolved.
The accuracy of the high $m$ modes, with their shorter wavelengths,
relies more on the accurate modeling of the interior 
volume of the toroid. If results from the two
calculations are to differ then, they will
most likely differ for low $m$ modes near threshold
which is indeed where we find the discrepancies. Overall, however, we  
find that the results of our exact calculations of the properties of 
the $I$ and $J$ modes in self-gravitating toroids and those from 
ICT studies, although not in quantitative agreement, exhibit
good qualitative agreement so that
ICTs are useful tools for the study of the nonaxisymmetric 
modes of disks.


The location of the instability thresholds for toroids does not 
depend as cleanly on other stability parameters used in ICT studies 
of star/disk systems. For example,  
(i) we find that the location of $I$ mode instability thresholds
vary strongly with $\epsilon$, the width of the toroid
(see Table 1). We find that
$\epsilon$ $\approx$ 0.91, 0.64, and 0.55 at the $m$ = 2 instability
thresholds for $q$ = 1, 3/2, and 2 toroids. (ii) We find 
that the locations of the
instability thresholds in terms of the self-gravity parameter, $p$,
although consistent with the results presented in Andalib, Tohline, 
\& Chrisotodoulou, show significant differences. Here,
$p$, the self-gravity parameter, is given by 
$p^2$ = $4\pi G\rho_{max}/\Omega_{max}^2$. In terms of $p$, 
we find that the
$I$ mode instability thresholds fall at $p$ = 2.7, 2.7,
and 3.5 for toroids with $q$ = 1, 3/2, and 2,
and that $J$ modes 
dominate $I$ modes when $p$ $>$ 
4.6, 5.2, and 5.7 for toroids with $q$ = 1, 3/2, and 2.
Andalib, Tohline, \& Chrisotdoulou (1997)
found that $I$ modes are unstable for $p$ $>$ 2.97, in 
good agreement with our results, but that 
$J$ modes are unstable for $p$ $>$ 5.205 but only dominate 
$I$ modes when $p$ $>$ 7.526. The predictions of Andalib,
Tohline, \& Christodoulou (1997) are most useful when
given in $T/|W|$-space for toroids, reasonably good when $p$
is used, and least useful in terms of $\epsilon$. 


Toroids are susceptible to $I$ and $J$ modes in the linear regime,
depending on the $T/|W|$ of the toroid but, {\it How does instability 
ultimately manifest itself?}  We follow the evolution
of unstable toroids into the nonlinear regime and find that 
toroids do not undergo prompt fission as a result of 
the development of $I$ mode instabilities 
({\it cf.}, Tohline \& Hachisu 1992, Woodward, Tohline, \& 
Hachisu 1994) or $J$ mode instabilities. The amplitudes never reach
values greater than ${\cal A}_m$ $\sim$ 0.6. The ultimate outcome of 
$I$ mode and $J$ mode instabilities appears to be quasi-stable, damping
slightly nonaxisymmetric figures, barlike
in the unstable $I$ mode regime or a higher order figure in the 
unstable $J$ mode regime. Further evolution of the toroids must 
probably then be drive by secular dissipative processes, such as radiative 
cooling or viscosity.



\subsection{ Protostellar Star/Disk Systems }



Low mass star formation occurs in cold, overdense regions in Giant 
Molecular Clouds (GMCs).
Cold cloud cores contain 1 to 2 $M_{\odot}$
and have large specific angular momenta, specific 
angular momenta much larger than found in any star
in the Galaxy ({\it e.g.,} see Tohline 2002). Cloud cores 
are typically evolved in the sense that they are centrally 
condensed. If the clouds are only margninally Jeans unstable 
then formation of star/disk systems is the likely outcome
of the collapse of thse cloud cores rather than a single 
star ({\it e.g.,} see Tohline 2002). 
For uniform density cloud cores which, initially,  
are in uniform  rotation, a small equilibrium object first forms 
from material with low specific angular momentum, material which
was initially close to the rotation axis of the cloud. 
Later, material with larger specific
angular momentum flows inward but is unable to accrete onto the
central object. The material forms a circumstellar
disk which slowly increases in mass and size with time. 
The central core generally contains 
less than a few percent of the total cloud mass under typical 
conditions (Adams \& Shu 1989). The system evolves as the circumstellar disk 
continues to accrete, growing in mass as the central core more or less 
remains at its original mass. The formation of stellar-sized
objects is then determined by how the mass contained in the
disk finds its way onto the central core. A likely transfer 
mechanism is the development of gravito-rotation driven 
nonaxismmetric instabilities in the disk. Such instabilities may arise
as the disk mass accretes and grows. Studies of star/disk systems
which do not include the effects of accretion
find that instability sets in when $M_*/M_d$ falls 
below $\sim$ 1 ({\it e.g.,} ref.).  


The work of Andalib, Tohline, \& Christodoulou (1997) suggests that toroids 
represent the limit of star/disk systems where the mass of the central
star goes to zero and do not form a set of singular disk solutions. 
If this is indeed correct, then results from our exact simulations of
instabilities in star/disk systems offer insight on the stability of
protostellar systems at late times when massive disks have formed around
small, centrally condensed cores, systems where the central mass compared
to the disk mass $M_*/M_d$ $\ll$ 1, such as are thought to arise under
typical conditions in star forming clouds (Shu \& Adams 198x).
We are currently investigating this question. 
Results for a simulation of a $(n,q)$ = (3/2,3/2) toroid
with $M_*/M_d$ = 0.05 and $R_-/R_+$ =
0.173 ($T/|W|$ = 0.219, $R_+$ = 4.973, $j$ = 0.7660, $\rho_{max}$ =
0.01233, $R_{max}$ = 2.3596, $\tau_{max}$ = 49.434 to be fully discussed
in a later paper (Hadley, Hickok, \& Imamura 200x) are now presented. 
For the barlike $m$ = 2 mode, $(\omega_2,\tau_2^{-1})$ 
=  -1,0354,2.3826(0.3792), $R_{\phi}$ = 1.6547, $R_{\lambda}$ = 3.16
with eigenfunction structure as shown in Figure xxx.
These results for the barlike mode are similar to 
those found for the (3/2,3/2) toroid with
$R_-/R_+$ $=$ 0.173 ($T/|W|$ $=$ 0.192) 
indicating that for this $M_*/M_d$ ratio, barlike modes in 
toroids behave like those in star/disk systems. In this sense, 
toroids are the limiting solution for star/disk
systems where the ratio of the stellar mass to 
the disk mass $M_*/M_d$ $\rightarrow$ 0. The
$m$ = 1 mode complicates this discussion as it shows 
different behavior for the system with $M_*/M_d$ $=$ 0.05
compared to the toroid. The toroid's
$m$  = 1 mode has $(\omega_1,\tau_1^{-1})$ = 
(-2.044,0.1098[9.52]) and so, at best, is only weakly
unstable. For the star/disk system with a low mass
central star, the $m$ = 1 mode dominates
the barlike mode with its $(\omega_1,\tau_1^{-1})$ = 
(-0.871,3.065[0.326]).  The $m$ = 1 mode
clearly grows much faster with corotation near $\rho_{max}$. 
Its eigenfunction, however, is similar in form to that of the $T/|W|$ 
= 0.19 toroid. So, although the nonaxisymmetric instabilities we 
find for toroids are apparently the limiting 
form for those of low mass star/disk systems, there are 
differences. The evolution of systems without central 
stars (toroids) are dominated by barlike instabilities while 
systems with even tiny central stars, $M_*/M_d$ = 0.05,
are dominated by $m$ = 1 modes. This and other differences 
are currently under investigation 
(Hadley, Imamura, \& Hickok 2009).


Toroids with the specific angular momentum distribution of a 
Maclaurin spheroid may form through the collapse of 
nearly uniform density rotating cores whose mass is 
significantly larger than the local Jeans mass 
({\it e.g.}, see Tohline 2002). Hachisu, Tohline,
\& Eriguchi (1987) investigated this process in 
(F,$T/|W|$)-space, where F is a parameter formed from 
$G, K,$ J, and $M$,
\begin{equation}
{\rm F} = [K^nG^{3-2n}M^{10-4n}J^{2n-6}]^{\frac{1}{n+1}}.
\end{equation}
F can be given in terms of 
the ratio of the thermal energy and gravitational
energy of the cold condensed cloud core, 
$\alpha_{\circ}$, and the ratio of the rotational kinetic 
energy and gravitational
energy of the cloud, $\beta_{\circ}$ = $T/|W|$,
\begin{equation}
{\rm F}  = 
1.049 \alpha_{\circ}^{0.6}\beta_{\circ}^{-0.143}.
\end{equation}
Clouds whose masses are greater than the local Jeans mass 
have $\alpha_{\circ}$ $\ll$ 1/2. Such cold, 
clouds are expected to collapse homologously with constant
F. In a homologous collapse,  $\beta$ = $T/|W|$ increases.
Along a sequence of collapsing spheroids, the barlike mode 
bifurcates at $T/|W|$ $\approx$ 0.27 
({\it e.g.}, Tassoul 1978), and ring-like sequences bifurcate 
at $T/|W|$ $\approx$ 0.41. Numerical simulations found, however,
that if spheroids were able to bypass the barlike instability, 
they became unstable to ring formation at $T/|W|$
$\approx$ 0.39, not $T/|W|$ = 0.41. This seeming discrepancy 
arises because the $T/|W|$ for the ring-like sequence is 
not a single-valued function of the angular momentum $j$. Two distinct 
equilibria are possible except when $j$ $\approx$ 2.86 and $T/|W|$ 
$\approx$ 0.39, where the two branches merge (see Figure 1 in 
Hachisu, Tohline, \& Eriguchi 1987).
Hachisu, Tohline, \& Eriguchi (1987) argued that collapsing clouds 
would evolve from low-to-high $T/|W|$ and eventually
settle in the lowest available $T/|W|$ equilibrium.
For small F (large $j$), they argued that condensed core 
collapses would approach Maclaurin-like toroidal equilibria rather 
than spheroidal equilibria. We calculated and then
determined the stability properties of selected Maclaurin-like
toroids. Our results are summarized in Tables 3 \& 4 where we give 
properties of the equilibrium toroids and the results of linear
stability analyses of selected equilibrium models. We find that
$m$ = 1, 2, 3, and 4 modes are all linearly
unstable in the models tested. This is not surprising given that the
Macluarin toroids tested all have $T/|W|$ $>$ 0.37.


\section{ CONCLUSIONS }


We studied the dynamic nonaxisymmetric
instabilities of thick, self-gravitating, inviscid polytropic
toroids using linear, quasi-linear, and nonlinear 
techniques. Self-gravitating toroids are susceptible to a broad range  
of gravito-rotation nonaxisymmetric instabilities. 
For polytropic toroids with power law angular velocity distributions
and nonaxismmetric modes with azimuthal dependence 
$e^{im\phi}$, where $m$ is the azimuthal mode number and 
$\phi$ is the aziumthal angle, we find that
instability sets in through
$m$ = 2 $I$ modes at $T/|W|$ $\sim$ 0.16-0.18 where
$T$ is the rotational kinetic
energy and $W$ is the gravitational energy. 
Instability in the $m$ = 2 $I$ mode
peaks in strength around $T/|W|$ = 0.22-0.23 and then 
weakens. For toroids with polytropic index $n$ = 3/2 and power 
law angular velocity distributions with indices 
$q$ = 3/2 and 2, the $m$ = 2 $I$ mode 
becomes stable at $T/|W|$ $\sim$ 0.26. Just before the $m$ = 2 
mode stabilizes, an $m$ = 1 $I$ mode instability sets in 
and is expected to dominate the $m$ = 2 mode 
for $T/|W|$ $\buildrel>\over\sim$ 
0.25. However, around the $T/|W|$ where the $m$ = 1 mode starts to 
dominate the $m$ = 2 $I$ mode, $J$ modes appear.
$J$ modes with $m$ = 2, 3, and 4 become unstable
for $T/|W|$ $\sim$ 0.25-0.26, 0.23-0.25, and 0.25-0.26, 
respectively, $m$ $\ge$ 3 $J$ modes quickly dominate the
$m$ = 1 $I$ mode and are expected to control the evolution of 
toroids with $T/|W|$ $\buildrel>\over\sim$ 0.26. 
Toroids are thus susceptible to $I$ and $J$ modes in the linear regime,
depending on the $T/|W|$ of the toroid however, we find 
that toroids do not undergo prompt fission as a result of 
the development of either $I$ or $J$ mode instabilities 
into the nonlinear regime. The ultimate outcome of $I$
and $J$ mode instabilities appears to be quasi-stable, damping
nonaxisymmetric figures with off-center maxima, either barlike
in form for the $I$ modes or tri-cornered or
square figures for the $J$ modes. Consequently, the further evolution
of such figures driven by cooling may lead to the onset of elliptical
instabilities and eventually to fission (Lebovitz \& Lfischitz 1996,
Ou \& Tohline 2006).


\newpage


\begin{center}
REFERENCES
\end{center}

\setlength{\parindent}{0em}

Adams, F.C., Ruden S.P. \& Shu F.H. 1989, Ap. J., 347, 959

Durisen, R.H., Gingold, R.A., Tohline, J.E., \& Boss, A.P. 1986
Ap. J., 305, 281

Durisen, R.H. \& Imamura, J.N. 1981, Ap. J., 243, 612

Durisen, R.H. \& Tohline, J.E. 1985, in Protostars and Planets II, eds.
  D.C. Black and M.S. Matthews (Tucson: U. Arizona)

Eriguchi, Y. \& Hachisu, I. 1984, Pub. Astr. Soc. J., 36, 491

Friedman, J.L. \& Schutz, B.F. 1978a, Ap. J., 221, 937

\underline{\hspace{1in}}. 1978b, Ap. J., 222, 281

Hachisu, I. 1986 Ap. J. S., 61, 479

Hachisu, I. \& Eriguchi, Y. 1985, As. Sp. Sci., 99, 71

Hachisu, I., Tohline, J., \& Eriguchi, Y. 1987, Ap. J., 323, 592

\underline{\hspace{1.0in}}. 1988, Ap. J. Suppl., 66. 315

Imamura, J.N., Friedman, J.L., \& Durisen, R.H. 1985,
Ap. J., 294, 474

Imamura, J.N., Toman, J., Durisen, R. H., Pickett, B., K., \&
Yang, S. 1995, Ap. J., 444, 363

Kojima, Y. 1986, Progr. Theor. Phys., 75, 251

Kojima, Y. 1989, M. N. R. A. S., 236, 589

Lin, D.N.C. \& Papaloizou, J.C.B. 1995, Ann. Rev. Astr. As., 33, 505

Lynden-Bell, D. \& Katz, J. 1981, Proc. Roy. Soc. London A, 378, 202

Pickett, B. K., Durisen, R. H., \& Davis, G.A. 1996, 
Ap. J., 458, 714

Pickett, B. K., Durisen, R. H., \& Link, R. 1997, Icarus, 126, 243

Tohline, J.E. 1980, Ap. J., 235, 866

\underline{\hspace{1.0in}}. 2002, Ann. Rev. Astr. As., ,

Tohline, J.E., \& Hachisu, I. 1990, Ap. J., 361, 394

Toman, J., Imamura, J.N., Pickett, B.K., \& Durisen, R.H. 1998,
Ap. J., 497, 

Wong, C.-Y. 1974, Ap. J., 190, 675

Woodward, J.W., Tohline, J.E.,\& Hachisu, I. 1994, Ap. J., 420, 247

\clearpage

\begin{deluxetable}{cccccccccc}
\tablecaption{ Equilibrium Models } 
\tablehead{
\colhead{$(n,q)$} &
\colhead{$\displaystyle\frac{T}{|W|}$} &
\colhead{$\rho_{max}$} &
\colhead{$\tau_{max}$} &
\colhead{$R_{max}$} &
\colhead{$R_{\lambda}$} &
\colhead{$j$} &
\colhead{$\displaystyle\frac{R_-}{R_+}$} &
\colhead{$R_+$} &
\colhead{$\epsilon$ }
}
\tablewidth{0pt}
\startdata
(3/2,1) &0.171  &0.01967 &33.66 &1.434 &1.434& 0.5532 & 0.024  & 4.261 &0.929\\
        &0.179  &0.017532&36.667&1.5625&1.583& 0.5872 & 0.0315 & 4.4428&0.910\\
        &0.192  &0.0146  &42.38 &1.806 &1.806& 0.644  & 0.04724& 4.73  &0.876 \\
        &0.215  &0.01014 &57.36 &2.395 &2.395& 0.76647& 0.09016& 5.4555&0.795\\
        &0.231  &0.00765 &72.369&2.940 &2.942& 0.87293& 0.134  & 6.0950&0.722\\
        &0.253  &0.004777&107.76&4.1068&4.107& 1.0811 & 0.220  & 7.4244&0.602 \\
        &0.260  &0.004077&123.56&4.5915&4.592& 1.1607 & 0.252  & 7.9606&0.563 \\
        &0.265  &0.003627&136.44&4.9714&4.97 & 1.2229 & 0.276  & 8.3902&0.534 \\
        &0.271  &0.003108&155.52&5.5141&5.514& 1.3098 & 0.307  & 9.0069&0.499 \\
        &0.280  &0.002372&196.69&6.6297&6.630& 1.4760 & 0.362  &10.238 &0.441 \\
        &0.302  &0.001171&368.44&10.737&10.74& 2.0094 & 0.500  &13.6225&0.366\\
        &0.321  &0.000612&656.18&16.552&16.55& 2.6593 & 0.611  &20.8035&0.232\\
        &0.364  &0.0002128&1673.3&33.61&33.61& 4.2596 & 0.754  &38.412 &0.138 \\
\\
(3/2,3/2)&0.163 &0.01967 &33.66 &1.702 &1.922& 0.525  & 0.126  & 4.138 &0.636 \\
        &0.168  &0.01747 &41.298&1.801 &2.019& 0.5480 & 0.134  & 4.256 &0.683 \\
        &0.171  &0.0169  &42.29 &1.844 &2.048& 0.5596 & 0.138  & 4.316 &0.677 \\
        &0.180  &0.01468 &47.14 &2.0369&2.252& 0.60587& 0.154  & 4.5583&0.655 \\
        &0.192  &0.0124  &53.79 &2.291 &2.838& 0.66354& 0.173  & 4.871 &0.632 \\
        &0.216  &0.00849 &72.78 &2.977 &3.48 & 0.81046& 0.223  & 5.712 &0.572\\
        &0.230  &0.00650 &90.34 &3.566 &3.735& 0.9251 & 0.260  & 6.411 &0.533\\
        &0.253  &0.00401 &133.9 &4.91  &5.342& 1.164  & 0.333  & 7.982 &0.459\\
        &0.263  &0.003175&163.14&5.739 &6.282& 1.2956 & 0.370  & 8.9077&0.426\\
        &0.273  &0.002485&200.39&6.733 &6.891& 1.4475 & 0.409  & 10.031&0.391 \\
        &0.280  &0.00203 &236.95&7.666 &8.059& 1.583  & 0.442  & 11.08 &0.361\\
        &0.300  &0.001132&391.48&11.240&11.36& 2.0472 & 0.535  & 14.99 &0.287\\
        &0.321  &0.000586&697.64&17.356&17.36& 2.7252 & 0.634  & 21.452&0.215\\
%        &0.369  &0.000222&1607.8&32.871&32.87& 4.231  & 0.756  & 37.40 &0.137\\
\\
(3/2,2) &0.172  &0.01504  &50.051 &2.1217&\nodata &0.5651&0.213 &4.4723 &0.551\\
        &0.181  &0.013159 &54.942 &2.3153&\nodata &0.6131&0.228 &4.72366&0.516\\
        &0.192  &0.0111   &63.248 &2.6125&\nodata &0.6780&0.249 &5.0757 &0.516\\
        &0.216  &0.007424 &86.513 &4.0674&\nodata &0.8444&0.299 &6.0355 &0.556\\
        &0.231  &0.005658 &107.03 &5.0674&\nodata &0.9712&0.335 &6.8193 &0.549\\
        &0.253  &0.003351 &155.32 &5.49  &\nodata &1.2192&0.398 &8.477  &0.385\\
        &0.264  &0.002751 &191.78 &6.476 &\nodata &1.3742&0.433 &9.5919 &0.359\\
        &0.276  &0.002078 &241.41 &7.7476&\nodata &1.5632&0.472 &11.0246&0.328\\
        &0.300  &0.0010295&439.20 &12.213&\nodata &2.1340&0.571 &15.8681&0.258\\
        &0.320  &0.0005638&733.96 &17.994&\nodata &2.7719&0.651 &22.0660&0.202\\
\enddata
\end{deluxetable}
\clearpage


%\begin{deluxetable}{cccccc}
%\tablecaption{ Eigenvalues: ($y_1(m),y_2(m)$) }
%\tablehead{
%\colhead{$(n,q)$} &
%\colhead{$\displaystyle\frac{T}{|W|}$} &
%\colhead{ $m$ = 2 } &
%\colhead{ $m$ = 3 } &
%\colhead{ $m$ = 4 } &
%\colhead{ jmax}
%}
%\tablewidth{0pt}
%\startdata
%(3/2,1)   & 0.171 & \nodata      & \nodata     & \nodata    &   256 \\
%          & 0.179 &-0.940,0.0216 & \nodata     & \nodata    &   256 \\
%          & 0.192 &-1.027,0.127  & \nodata     & \nodata    &   256 \\
%          & 0.215 &-0.996,0.378  & \nodata     & \nodata    &   256 \\
%          & 0.231 &-0.996,0.439  &-0.779,0.0493& \nodata    &   256 \\
%          & 0.253 &-0.945,0.362  &-0.982,0.127 & \nodata    &   256 \\
%          & 0.260 &-0.905,0.285  &-0.499,0.2566& \nodata    &   256 \\
%          & 0.265 &-0.781,0.130  &-0.497,0.455 &-0.429,0.211&   256 \\
%          & 0.271 &-0.534,0.157  &-0.4385,0.611&0.4803,0.3947&  256 \\
%          & 0.280 &-0.347,0.296  &-0.3421,0.8417&-0.4381,0.72648& 256 \\
%          & 0.302 &+0.0266,0.600 &-0.1147,1.3314&-0.19815,1.4849& 256 \\
%          & 0.321 &+0.1589,0.902 &-0.0055,1.6844&-0.0713,2.0306 & 256 \\
%          & 0.364 &+0.0260,1.379 &3.0566,2.102  &+0.03358,2.624 & 256 \\
%(3/2,3/2) & 0.163 &-0.9124,0.0624& \nodata     & \nodata    &   256 \\
%          & 0.171 &-0.9834,0.128 & \nodata     & \nodata    &   256 \\
%          & 0.180 &-1.019,0.230  & \nodata     & \nodata    &   256 \\
%          & 0.192 &-1.038,0.321  & \nodata     & \nodata    &   256 \\
%          & 0.216 &-1.035,0.403  & \nodata     & \nodata    &   256 \\
%          & 0.230 &-1.002,0.390  & \nodata     & \nodata    &   256 \\
%          & 0.253 &-0.914,0.217  &-0.483,0.224 & \nodata    &   256 \\
%          &       &-0.898,0.179  &             & \nodata    &   512 \\
%          & 0.263 & \nodata      &-0.395,0.549 &-0.380,0.324&   256 \\
%          & 0.273 &-0.240,0.291  &-0.286,0.778 &-0.374,0.643&   256 \\
%          &       &-0.           &-0.315,0.800 &-0.         &   512 \\
%          & 0.280 &-0.170,0.387  &-0.252,0.932 &-0.338,0.874&   256 \\
%          &       &-0.235,0.503  &- &-0.&   512 \\
%          & 0.300 &-0.048,0.722  &-0.1061,1.35 &-0.1761,1.50&   256 \\
%          & 0.321 &+0.1463,0.9666&+0.0226,1.712&-0.0372,2.0643& 256 \\
%          & 0.369 &+0.117,1.389  &+0.1001,2.076&+0.0623,2.5737 &   256 \\
%(3/2,2)   & 0.172 &-1.045,0.0798 & \nodata     & \nodata    &   256 \\
%          & 0.181 &-1.0707,0.219 & \nodata     & \nodata    &   256 \\
%          & 0.192 &-1.067,0.305  & \nodata     & \nodata    &   256 \\
%          & 0.216 &-1.040,0.394  & \nodata     & \nodata    &   256 \\
%          & 0.231 &-1.009,0.383  & \nodata     & \nodata    &   256 \\
%          & 0.253 &-0.907,0.171  &-0.4794,0.240& \nodata    &   256 \\
%          & 0.264 &  \nodata     &-0.3757,0.602&-0.386,0.388&   256 \\
%          & 0.276 &-0.154,0.376  &-0.272,0.863 &-0.3552,0.7644& 256 \\
%          & 0.300 &+0.03958,0.74678&-0.0507,1.3347&-1266,1.551& 256 \\
%          & 0.320 &+0.193,0.945  &+0.0319,1.6990&-0.0313,2.053& 256 \\
%\enddata 
%\end{deluxetable}
%\clearpage

\begin{deluxetable}{ccccccccccc}
\tablecaption{ Barlike ($m$ = 2) Mode Properties }
\tablehead{
\colhead{$(n,q)$} &
\colhead{$\displaystyle\frac{T}{|W|}$} &
\colhead{$\left(\displaystyle\frac{\omega_2}{\Omega_{max}},
\displaystyle\frac{\tau_2}{\tau_{max}}\right)$} &
\colhead{$\displaystyle\frac{R_{\phi}}{R_{max}} $} &
\colhead{$M_{\Upsilon}$} &
\colhead{$\displaystyle\frac{j_{\Upsilon}}{j} $} &
\colhead{$\displaystyle\frac{R_{\Upsilon}}{R_{max}} $} &
\colhead{$\displaystyle\frac{R_{ILR}}{R_{max}}$} &
\colhead{$\displaystyle\frac{R_{co}}{R_{max}} $} &
\colhead{$\displaystyle\frac{R_{OLR}}{R_{max}}$} 
}
\tablewidth{0pt}
\startdata
(3/2,1)   & 0.179 &-1.072,7.309   & 0.419 &0.556 & 0.654 & 1.800 & 0.551
& 1.88 &\nodata\\
          & 0.192 &-0.973,1.255   & 0.420 &0.629 & 0.515 & 1.506 & 0.602
& 2.054& \nodata\\
          & 0.215 &-1.005,0.434   & 0.605 &0.650 & 0.545 & 1.286 & 0.583
& 1.99 &\nodata\\
          & 0.231 &-1.004,0.363   & 0.755 &0.664 & 0.564 & 1.31  & 0.588
& 2.00 &\nodata \\
          & 0.253 &-1.055,0.440   & 0.922 &0.593 & 0.504 & 1.17  & 0.491
& $>R_{+}$ &\nodata\\
          & 0.265 &-1.219,1.228   & 0.980 &0.533 & 0.452 & 1.092 & 0.479
& 1.637 &\nodata\\
          & 0.271 &-1.466,1.013   & 1.019 &0.523 & 0.447 & 1.070 & \nodata
& 1.367 &\nodata\\
          & 0.280 &-1.653,0.537   & 1.024 &0.515 & 0.449 & 1.054 & \nodata
& 1.276 &\nodata\\
          & 0.302 &-2.02656,0.2654& 1.032 &0.506 & 0.458 & 1.0245& \nodata
& 0.9872&\nodata\\
          & 0.321 &-2.159,0.176   & 1.030 &0.493 & 0.469 & 1.0210& \nodata 
&0.924 &\nodata \\
          & 0.364 &-2.024,0.115   & 0.942 &0.493 & 0.473 & 0.939 & \nodata 
&0.932 &\nodata \\
(3/2,3/2) & 0.163 &-1.0826,2.551 & 0.656 &0.435 & 0.411 & 1.404 & 0.942 
& 1.498 & 1.962\\
          & 0.171 &-1.017,1.240  & 0.677 &0.514 & 0.483 & 1.40  & 0.988
& 1.57 & 2.057\\
          & 0.180 &-0.980,0.692  & 0.692 &0.582 & 0.546 & 1.394 & 1.015
& 1.613& 2.135\\
          & 0.192 &-0.962,0.497  & 0.732 &0.629 & 0.595 & 1.37  & 1.026 
& 1.63 &\nodata\\
          & 0.216 &-0.954,0.400  & 0.842 &0.663 & 0.616 & 1.27  & 1.0316  
& 1.64 &\nodata\\
          & 0.231 &-0.997,0.397  & 0.887 &0.630 & 0.588 & 1.19  & 1.000
& 1.59 &\nodata\\ 
          & 0.253 &-1.086,0.730  & 0.968 &0.557 & 0.517 & 1.10  & 0.950
& 1.51 &\nodata\\
          & 0.273 &-1.76,0.546   & 1.035 &0.505 & 0.474 & 1.04  & 0.689
& 1.095& 1.435\\
          & 0.280 &-1.833,0.439  & 1.034 &0.501 & 0.472 & 1.03  & 0.663
& 1.054& 1.381\\
          & 0.300 &-1.951,0.220  & 1.026 &0.505 & 0.483 & 1.0231& \nodata
& 1.014&\nodata\\
          & 0.321 &-2.1463,0.2131& 1.019 &0.485 & 0.469 & 1.0083& \nodata
& 0.9507&\nodata\\
%          & 0.369 &-2.117,0.1146 & 1.009 &0.481 & 0.470 & 1.0024& \nodata 
%&0.9644& \nodata\\
(3/2,2)   & 0.172 &-0.956,2.021  & 0.768 &0.578 & 0.579 & 1.381 &\nodata
& 1.45 &\nodata\\
          & 0.181 &-0.9286,0.7292& 0.791 &0.632 & 0.633 & 1.38  &\nodata
& 1.464 &\nodata\\
          & 0.192 &-0.889,0.705  & 0.824 &0.616 & 0.616 & 1.35  &\nodata
& 1.50 &\nodata\\
          & 0.216 &-0.960,0.404  & 0.889 &0.654 & 0.6533& 1.208 &\nodata
& 1.446&\nodata\\
          & 0.231 &-0.9912,0.4161& 0.927 &0.595 & 0.595 & 1.138 &\nodata
& 1.424&\nodata\\
          & 0.253 &-1.0927,0.931 & 0.976 &0.515 & 0.515 & 1.063 &\nodata
& 1.35 &\nodata\\
          & 0.276 &-1.8462,0.423 & 1.034 &0.499 & 0.500 & 1.034 &\nodata
& 1.039&\nodata\\
          & 0.300 &-2.03958,0.213& 1.020 &0.487 & 0.488 & 1.0153&\nodata
& 0.9907&\nodata\\
          & 0.320 &-2.193  ,0.169& 1.019 &0.493 & 0.495 & 1.011 &\nodata
& 0.956 &\nodata
\enddata
\end{deluxetable}

\clearpage

%
%\begin{deluxetable}{ccccccc}
%\tablecaption{ ${\it Critical}$ Radii: $(n,q)$=(3/2,3/2) Toroids }
%\tablehead{
%\colhead{$\displaystyle\frac{T}{|W|}$} &
%\colhead{$\displaystyle\frac{R_{\phi}}{R_{max}} $} &
%\colhead{$\displaystyle\frac{R_Q}{R_{max}} $} &
%\colhead{$\displaystyle\frac{R_{\Upsilon}}{R_{max}} $} &
%\colhead{$\displaystyle\frac{R_{ILR}}{R_{max}}$} &
%\colhead{$\displaystyle\frac{R_{co}}{R_{max}} $} &
%\colhead{$\displaystyle\frac{R_{OLR}}{R_{max}}$} 
%}
%\tablewidth{0pt}
%\startdata
%0.163 & 0.656 & 0.764,1.915 & 1.404 & 0.942 
%& 1.498 & 1.962\\
%0.171 & 0.677 & 0.781,1.838 & 1.40  & 0.988
%& 1.57 & 2.057\\
%0.180 & 0.692 & 0.795,1.753 & 1.394 & 1.015
%& 1.613& 2.135\\
%0.192 & 0.732 & 0.8075,1.663 & 1.37  & 1.026 
%& 1.63 &\nodata\\
%0.216 & 0.842 & 0.846,1.5015 & 1.27  & 1.0316 
%& 1.64 &\nodata\\
%0.231 & 0.887 & 0.855,1.402 & 1.19  & 1.000
%& 1.59 &\nodata\\
%0.253 & 0.968 & 0.896,1.295 & 1.10  & 0.950
%& 1.51 &\nodata\\
%0.273 & 1.035 & 0.906,1.206 & 1.04  & 0.689
%& 1.095& 1.435\\
%0.280 & 1.034 & 0.9105,1.1767 & 1.03  & 0.663
%& 1.054& 1.381\\
%0.300 & 1.026 & 0.934,1.121 & 1.0231& \nodata
%& 1.014&\nodata\\
%0.321 & 1.019 & 0.9507,1.083 & 1.0083& \nodata
%& 0.951&\nodata\\
%0.369 & 1.009 & 0.967,1.0465 & 1.0024& \nodata
%&0.9644& \nodata\\
%\enddata
%\end{deluxetable}
%
%\clearpage
%

\begin{deluxetable}{cccccccc}
\tablecaption{Equilibrium Models: Maclaurin Sequences } 
\tablehead{
\colhead{$n$} &
\colhead{$T/|W|$} &
\colhead{$\rho_{max}$} &
\colhead{$\Omega_{max}$} &
\colhead{$R_{max}$} &
\colhead{$J$} &
\colhead{$R_+$} &
\colhead{$R_-/R_+$}
}
\tablewidth{0pt}
\startdata
3/2 &0.4066 & 2.25$\times10^{-4}$& 0.00816&
                        10.6 & 3.01&   55.2& 0.008 \nl
3/2 &0.4062 & 2.27$\times10^{-4}$& 0.00859&
                        11.0 & 3.00&   54.9& 0.016 \nl
3/2 &0.4027 & 2.47$\times10^{-4}$& 0.00831&
                        11.0  & 2.94&  52.8& 0.048 \nl
3/2 &0.3978 & 2.70$\times10^{-4}$& 0.00875&
                        11.4  & 2.88&  50.7& 0.080 \nl
3/2 &0.3929 & 2.86$\times10^{-4}$& 0.00759&
                        12.3  & 2.85&  49.8& 0.112 \nl
3/2 &0.3884 & 2.90$\times10^{-4}$& 0.00711&
                        13.6  & 2.86&  50.0& 0.144 \nl
3/2 &0.3842 & 2.84$\times10^{-4}$& 0.00653 &
                        15.3  & 2.90&  51.7& 0.176 \nl
3/2 &0.3802 & 2.67$\times10^{-4}$& 0.00523 &
                        17.1  & 2.99&  54.7& 0.208 \nl
3/2 &0.3765 & 2.41$\times10^{-4}$& 0.00467 &
                        20.1  & 3.14&  59.8& 0.240 \nl
3/2 &0.3729 & 2.04$\times10^{-4}$& 0.00353 &
                        24.2  & 3.37&  68.7& 0.272 \nl
3/2 &0.3695 & 1.53$\times10^{-4}$& 0.00261 &
                        32.2  & 3.81&  87.5& 0.304 \nl
\nl
\enddata
\end{deluxetable}


\begin{deluxetable}{ccccccc}
\tablecaption{Eigenvalues: Maclaurin Sequences }
\tablehead{
\colhead{$n$} &
\colhead{R$_-$/R$_+$} &
\colhead{$T/|W|$} &
\colhead{$m=1$} &
\colhead{$m=2$} &
\colhead{$m=3$} &
\colhead{$m=4$} 
}
\tablewidth{0pt}
\startdata
3/2,m1& 0.008&0.4066
               & -0.00953,0.00360
               & -0.0105,0.00503 
               & -0.0209,0.00574 
               & \nodata \nl
3/2,m4& 0.080&0.3978
               & -0.0104,0.00330 
               & -0.0104,0.00521
               & -0.0212,0.00564 
               & \nodata \nl
3/2,m8& 0.208&0.3820
               & -0.0172,0.00320 
               & -0.00641,0.00487
               & -0.017251,0.003249
               & \nodata \nl
\nl
\enddata
\end{deluxetable}
\clearpage

%
\begin{figure}
\includegraphics[width=6.0in]{eigenvaluesm2-4.ps}
\caption{
Low $m$ eigenvalues for toroids with $q$ = 1, 3/2, and 2. 
The left column shows $y_1$, 
the growth rates and the right column shows $y_2$, the oscillation
frequencies. The squares with a dot are for $m$ = 1, 
"+" are for $m$ = 2, the "$\times$" are for $m$ = 3, and 
the "$*$" are for $m$ = 4. Legend: a) q = 1; b) q = 3/2; c) q = 2.
}
\end{figure}

%## bar modes

\begin{figure}
\includegraphics[width=5.0in]{combine_FIG1}
\caption{
The $m$ = 2 mode eigenfunctions and constant phase loci in 
the equatorial plane for the 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.163, 0.192, 0.253, and 0.273 toroids
from top-to-bottom.
}
\end{figure}

%## bar modes work 

\begin{figure}
\includegraphics[width=5.0in]{combine_FIG2}
\caption{
The perturbed energies, $\left<{\cal E}\right>_k$ and
$\left<{\cal E}\right>_hk$ and, the streses, $\sigma_R$, 
$\sigma_{\Pi}$, and $\sigma_{\Phi}$ for the
$m$ = 2 modes of the (n,q) = (3/2,3/2), $T/|W|$ = 
0.163, 0.192, 0.253, and 0.273 toroids,
from top-to-bottom.
}
\end{figure}

%## bar modes angular momentum and self-interaction torque 

\begin{figure}
\includegraphics[width=5.0in]{combine_FIG3}
\caption{
The self-interaction gravitational torque $\gamma$ and
the azimuthally averaged angular momentum of the perturbed 
modes $\delta J$ 
for the (n,q) = (3/2,3/2), $T/|W|$ =
0.163, 0.192, 0.253, and 0.273 toroids,
from top-to-bottom.
}
\end{figure}

%#   critical radii

\begin{figure}
\includegraphics[width=5.0in]{critical_radii.ps}
\caption{
Critical radii for the $m$ = 2 barlike modes for the $(n,q)$ = 
(3/2,3/2) toroid sequence.
}
\end{figure}

%# One-Armed Spiral eigenfunctions and work

\begin{figure}
\includegraphics[width=5.0in,angle=270]{combine_m1_rho_FIGa.ps}
\caption{
The $m$ = 1 mode eigenfunction and properties for the 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.253 toroid.
}
\end{figure}

\begin{figure}
\includegraphics[width=5.0in]{combine_m1_work_FIGb.ps}
\caption{
The $m$ = 1 mode eigenfunction properties for the 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.253 toroid.
}
\end{figure}

\clearpage

%##  end of m = 1 mode eigenfunctions and work

%##  Ratio of integrated gravity and Reynolds stress

\begin{figure}
\includegraphics[width=5.0in,angle=270]{work_integrals.ps}
\caption{
The ratio ${\cal R}$ = $(\int\sigma_{{\Phi}}d^3x)/(\int\sigma_{R}d^3x)$
for the $m$ = 2 mode for the $(n,q)$ = (3/2,3/2) toroid sequence.
}
\end{figure}

%## Ratio of integrated gravity and Reynolds stress

%## Nonlinear evolutions

%## T/|W| = 0.192

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_192_eigen_amp.ps}
\caption{The ${\cal A}_m$ for $m$ = 1 to 8, cos$(\phi_{\rho})$ for
a point in the equatorial plane, and $T/|W|$ for the 
$(n,q)$ = (3/2,3/2) $T/|W|$ = 0.19 toroid.
}
\end{figure} 

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_192_nl_eigen_rho.ps}
\caption{The $m$ = 2
eigenfunction $|\delta\rho|$ in the equatorial plane 
of the $(n,q)$ = (3/2,3/2), $T/|W|$ = 0.192 toroid at
$t$ = 3.8 $\tau_{max}$ and xxx $\tau_{max}$.
}
\end{figure} 
\clearpage

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_192_nl_eigen_phase.ps}
\caption{The constant phase locus for the $m$ = 2 mode 
in the equatorial plane of
the $\delta\rho$ for the
$(n,q)$ = (3/2,3/2), $T/|W|$ = 0.192 toroid at
$t$ = 3.8 $\tau_{max}$ and xxx $\tau_{max}$. The circles mark the 
inner and outer radii of the equilibrium toroid.
}
\end{figure} 
\clearpage

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_192_eigen_omega.ps}
\caption{
The specific angular momentum distribution $h(m_c)$ and angular velocity
distribution $\Omega$ for the final model ($t$ $\sim$ 19 $\tau_{max}$) in the 
$(n,q)$ = (3/2,3/2), $T/|W|$ = 0.19 toroid nonlinear simulation.  In the
left panel, the blue line has slope $\varpi^{-1.05}$ and the green line
has slope $\varpi^{-1.73}$. In the right panel, the green line shows the
specific angular momentum distribution of a Maclaurin spheroid while 
the red line is the specific angular momentum distribution of the final
model.
}
\end{figure} 

%## T/|W| = 0.192

\clearpage

%## T/|W| = 0.253

\begin{figure}
\includegraphics[width=5.0in]{combine_b253_FIG5.ps}
\caption{
The $m$ = 1, 2, and 3 mode eigenfunctions for the 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.253 toroid.
}
\end{figure}


\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15b253_eigen_rho.ps}
\caption{
$m$ = 1 and 2 relative density eigenfunctions in the equatorial
plane for the $(n,q)$ = 
(3/2,3/2) toroid with $T/|W|$ = 0.253 at time $t$ = 3.3 $\tau_{max}$.
The central circle in the phase plots shows the location of the 
maximum density $\rho_{max}$.
}
\end{figure}


\begin{figure}
\includegraphics[width=5.0in]{combine_b273_FIG4.ps}
\caption{
The $m$ = 1, 2, 3, and 4 mode eigenfunctions for the 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.273 toroid.
}
\end{figure}

%## linear eigenfunctions

\clearpage

%## nonlinear simulations results

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_272_nl_eigen_rho.ps}
\caption{The $m$ = 3 
eigenfunction $|\delta\rho|$ in the equatorial plane 
of the $(n,q)$ = (3/2,3/2), $T/|W|$ = 0.272 toroid at
$t$ = 1.98  $\tau_{max}$ and 2.84 $\tau_{max}$.
}
\end{figure} 
\clearpage

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_272_nl_eigen_phase.ps}
\caption{The constant phase locus for the $m$ = 3 mode 
in the equatorial plane of
$\delta\rho$ for the
$(n,q)$ = (3/2,3/2), $T/|W|$ = 0.272 toroid at
$t$ = 1.98 $\tau_{max}$ and 2.84 $\tau_{max}$. The circles mark the inner
and outer radii of the equilibrium toroid.
}
\end{figure} 
\clearpage


%## nonlinear simulation results

\clearpage

%## start M_*/M_d = 0.05

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_m0p05_b219_m1_torque.ps}
\caption{
The $m$ = 1 mode eigenfunction and properties for the 
star/disk system with $M_*/M_d$ = 0.05 and disk with 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.219 toroid.
}
\end{figure}

\begin{figure}
\includegraphics[width=5.0in,angle=270]{n15q15_m0p05_b219_m2_torque.ps}
\caption{
The $m$ = 2 mode eigenfunction and properties for the 
star/disk system with $M_*/M_d$ = 0.05 and disk with 
($n,q$) = (3/2,3/2), $T/|W|$ = 0.219 toroid.
}
\end{figure}

%## end M_*/M_d = 0.05

\end{document}