text:
# Hydromagnetic Stability of a Slim Disk in a Stationary Geometry

## 1 Introduction

The process of accretion onto compact objects has been long recognized as the primary mechanism to power the most luminous events in space. In the traditional picture of Shakura & Sunyaev (1973) and Novikov & Thorne (1973), entropy is generated and radiated locally from the free energy available in a shear flow with a Keplerian angular velocity profile. Two salient oversimplifications of this framework have been the focus of intense research and progress in the last decade: energy advection by the flow and free-energy tapping and angular momentum transport through magnetohydrodynamical processes.

The magnetorotational instability or MRI (Velikov 1959, Chandrasekhar 1961, Balbus & Hawley 1991), justifies the long-sought mechanism for efficient, turbulent transport of angular momentum that enables accretion disks to operate with astrophysically interesting mass accretion rates (Pringle 1981). The importance of this process cannot be overstated: By catalyzing accretion into gravitational wells, the MRI enables a plethora of astrophysical phenomena to occur, from protostar formation inside molecular clouds to jet launching in quasars. The MRI also holds the key to understand the extraction of free energy from the differential shear flow of otherwise hydrodynamically stable disks (Balbus et al. 1999, Godon & Livio 2000).

On the observational front, the wealth of high-quality data from spectral and timing devices aboard space-borne high-energy observatories has turned out the most compelling evidence yet of accretion onto black holes. The discoveries of pairs of high frequency quasi-periodic oscillations in RXTE X-ray timing data from microquasars GRO J1655-40 and GRO 1915+105 (Strohmayer 2001a,b) have brought the spotlight to hydrodynamical models of adiabatic global excitations of the inner disk, a.k.a. diskoseismology models (Perez et al. 1997) or relativistic precession models (Stella et al. 1999). Interestingly, the role of magnetic fields has been largely ignored in spite of clear evidence that QPO’s, being non-thermal, hard X-ray phenomena, likely do not originate in the accretion disk proper but rather on a magnetically active accretion disk corona (R. Blandford, Priv.Comm.). Likewise, the recent report (Wilms et al. 2001) of the detection of a very broad Fe K$`\alpha `$ feature on the XMM-EPIC spectrum of MCG-6-30-15, has made a very strong case for the inadequacy of standard models of energy deposition in accretion disks. The proposed exits to this paradox–extraction of black hole spin energy (Blandford & Znajek 1977), or non-zero torque at the marginally stable circular orbit radius, $`r_{\mathrm{ms}}`$ (Agol & Krolik 2000)– both rely on magnetic coupling between a standard disk and the flow inside $`r_{\mathrm{ms}}`$.

On a more exotic front, theoretical progress in our understanding of accretion processes at the most extreme imaginable conditions–stellar-mass black holes hyper-accreting at twelve orders of magnitude above the Eddington limit–requires attention to a detailed physical account of highly relativistic accretion flows. Aside from the potential to explain gamma-ray burst phenomenology, such studies are chiefly relevant to assess the likelihood of “failed” supernovae as gravitational wave sources (Fryer et al. 2001). Indeed, when neutrino trapping occurs at $`\dot{M}\text{ }>\text{ }1M_{}`$ sec-1 (Popham et al. 1998), the associated dynamical stress will mimic the effects of radiation stress in standard disks where clumpy accretion ensues (Turner et al. 2001). If the mass fraction in the clumps is large, prolific gravitational waves will be emitted from the mass quadrupole moment associated with the bulk motion of large mass over-densities. Such a scenario will also lead to excitations of the black hole’s geometry which, at high values of the spin parameter $`a`$, can produce highly characteristic, monochromatic black hole ringing as the geometry settles towards a quiescent Kerr state (Araya-Góchez 2003). Remarkably, the expectation of a large mass fraction in the clumps is reasonable and justifiable by the physical picture of near-hole accretion presented herein.

An outstanding issue yet to be addressed in light of recent theoretical progress is our view of black hole accretion inside the marginally stable orbit, the putative inner boundary of standard accretion disk theory. In particular, very little is concretely known about the inertial effects of strong gravity on the relevant MHD processes. Previous work has either assumed pure hydrodynamical flow (and energetically negligible energy release) or, alternatively, laminar flow under ideal MHD conditions (Krolik 1999, Gammie 1999). Krolik (1999) has brought out an interesting point: Under mere flux freezing conditions the assumption of ballistic orbits in the plunging region is never self-consistent; when the radial velocity component is significant, the magnetic field energy density becomes comparable to the rest mass energy density of the matter.

In this paper, we address the issue of stability of the magnetic field (co-moving frame) in a stationary, axially symmetric background geometry. Curiously, the two key developments in accretion disk theory of the last few years may have come of age to properly address the problem at hand: Inside $`r_{\mathrm{ms}}`$ the accretion flow will be mildly advective, with a slightly sub-Keplerian angular velocity profile and possibly supported in part by the radial pressure gradient of a hot MHD fluid with significant relativistic enthalpy (see Popham & Gammie (1998) solutions for moderate values of $`\alpha `$ and advected fraction $`f`$). In this spirit, we argue in §5 that the natural evolution of the MRI inside the marginally stable orbit is at least consistent with this view.

The (magnetohydro)dynamics of black hole accretion comprises two important aspects that have received relatively little attention: the effects of radiation pressure (see, however, Blaes & Socrates and Turner et al. 2001), and the effects of strong gravity (see footnote 7 of Gammie & Popham 1998). We will address the former problem in a future paper (Araya-Góchez & Vishniac 2002) while concentrating on general relativity in this one. As a background, §2 looks at the Lagrangian displacement vector field formulation of the MRI concentrating on inertial and compressibility effects. In section 3, we develop a fully covariant theory of the instability. The intention is to build a theoretical framework from first principles in order to avoid missing any subtleties associated to the full incorporation of gravitational effects (e.g., reference is made to the Cowling approximation and to the fixing of the gauge associated with the component of the Lagrangian diplacement along the fluid’s four-velocity). The elastic reponse of the field is computed by noting that the surface of invariance of the Faraday tensor attributes mathematically identical variational properties to the two four-vectors that span it: the magnetic field four-vector and the fluid’s four-velocity. We then make the minimal modifications to the relativistic fluid equations that allow for the inclusion of a coherent magnetic field and undertake a local stability analysis of this field in the medium of a slim disk around a rotating black hole, while suppressing compression. The role of compressibility in a photon gas is then briefly assessed.

## 2 A Lagrangian Formulation of the MRI in Compressible Media

The MRI is essentially a local instability. In the frame of the fluid, the interplay of inertial “forces” with the elastic coupling of fluid elements creates an unstable situation to the redistribution of specific angular momentum, $`\mathrm{}`$. Without the elastic coupling provided by the bending of field lines, such inertial forces–namely, the shear(tide) and the coriolis terms–induce radial epicyclic motions while preserving specific angular momentum in collisionless fluids (e.g., stars in the Galaxy). This is related to the Rayleigh criterium for stability of a differentially rotating fluid: $`r^3`$ d$`{}_{r}{}^{}\mathrm{}_{}^{2}=\kappa ^20`$, where $`\kappa `$ is the frequency of epicyclic motions.

In the weak field limit, one may construct a dispersion relation quite independently of the specific magnetic field topology: Highly sub-thermal fields, $`v_{_{\mathrm{Alf}}}/c_s1`$, guarantee that the instability is truly local222 When the field is non-negligible, $`k_{}^1`$ may approach the disk’s pressure scale height and in the case of supra-thermal toroidal fields, non-axisymmetric modes have fastest growing wave numbers that may approach the inverse radial scale-length (Foglizzo & Tagger 1995). , occurring at large values of $`k_{}\mathrm{\Omega }/v_{_{\mathrm{Alf}}}(𝐤\mathrm{𝟏}_{\widehat{𝐁}})`$. In this simplified approach, the global disk structure is ignored (no curvature nor radial structure) and the response of the field amounts to nothing more than providing a restoring force to displacements from equilibrium (Balbus & Hawley 1992). Indeed, in the horizontal regime of Lagrangian displacement two orthogonal field topologies yield nearly identical mathematical dispersion relations for wavemodes: axisymmetric perturbations of a meridional field and non-axisymmetric perturbations of a toroidal field. The former case corresponds to the “classical” Balbus-Hawley instability and it’s physical relevance is free of controversy. The relative importance of the latter analyses is a more subtle issue.

A somewhat technical point–well discussed in the review by Balbus & Hawley (1998)–is the non-locality induced by shear on wave-modes with Eulerian coordinate phase-dependencies. For $`k_\phi 0`$ modes, shear evolves the radial component of wavenumbers according to $`k_r(t)=k_{0r}[\mathrm{d}_{\mathrm{ln}r}\mathrm{\Omega }]k_\phi t`$, which means that modes that could be “unstable” are only so, transiently. The maximum instantaneous growth rate occurs when $`k_r0`$ and matches that of the local axisymmetric modes. In the end, this issue turns out to be more academic than practical but it stresses the importance of treating the instability locally, in co-moving coordinates. The down side is that this greatly complicates global approaches that rely on eigenmode solutions in Eulerian coordinates extrinsic to the fluid. On the other hand, in a local approach azimuthal wavenumbers are no longer discrete (Ogilvie & Pringle 1996) and consequently, neither are the co-moving frequencies (see below).

A related issue concerns the relevance of non-axisymmetric mode analyses when the magnetic field is not purely toroidal. Balbus & Hawley (1998) argue that the strict ordering of wavenumber components (and narrow phase space) necessitated to achieve fastest growth: $`k_\phi k_rk_\theta `$, ensue in violent poloidal Alfvénic couplings that promptly take over the dynamics. Non-axisymmetric modes, however, are important for at least two key reasons: a- the ordering is not so restrictive when the fields are not weak (as needed to explain $`\alpha `$ values of a few tenths), and b- compressive, non-axisymmetric modes are fundamental to examine energy deposition when radiation stress becomes significant (Araya-Góchez & Vishniac 2002). Moreover, because the dispersion relations relate simply (at least in the horizontal regime of fastest growth), it is rather lucrative to examine both cases at once.

Aiming to formulate a fully covariant relativistic theory of the MRI in §3, this section conducts the same in three dimensions. The linear stability analysis is carried out in terms of the Lagrangian displacement vector field, $`𝝃`$. Foglizzo (1995) has stressed the usefulness of this approach to account for the polarization of compressive MHD modes. A simple meridional stratification profile sets the physical scale-length of the problem: d$`{}_{z}{}^{}\mathrm{ln}\rho =^1`$, with gas, radiation (and possibly magnetic) pressures tracking the unperturbed density profile $`\rho \mathrm{\Omega }=p_{\mathrm{r}_+\mathrm{g}}+p__𝐁`$. The problem naturally splits in two parts: computation of the inertial/geometric terms §2.1, and computation of the body forces from gas, radiation and electromagnetic stresses §2.2. We avoid going into the rotating frame from the onset in order to preserve a transparent connection to a “universal” standard of rest frame (to be associated with Boyer-Lindquist coordinates).

### 2.1 inertial terms

Inertial accelerations are geometrically imprint in the connection terms for the covariant derivatives of the Eulerian velocity components. For spherical coordinate motion $`(\dot{r},\dot{\phi },\dot{\theta })(V^r,V^\phi ,V^\theta )`$, the only non-trivial connections are $`\mathrm{\Gamma }_{\phi \phi }^r\mathrm{\Gamma }_{r\phi }^\phi `$. Denoting the Lagrangian time derivative by $`\mathrm{d}_t_t+𝐕\mathbf{}`$, the three components of Euler’s equation read

$`\mathrm{d}_tV^r=_tV^r+V^jV_{,j}^r+(r)V^\phi V^\phi `$ $`=`$ $`g^{rj}𝐟_j`$ (1)
$`\mathrm{d}_tV^\phi =_tV^\phi +V^jV_{,j}^\phi +(2/r)V^rV^\phi `$ $`=`$ $`g^{\phi j}𝐟_j`$ (3)
$`\mathrm{and}\mathrm{d}_tV^\theta =_tV^\theta +V^jV_{,j}^\theta `$ $`=`$ $`g^{\theta j}𝐟_j`$ (5)
$`\mathrm{where}𝐟{\displaystyle \frac{1}{\rho }}\mathbf{}p+{\displaystyle \frac{1}{4\pi \rho }}𝐉\times 𝐁,`$ (7)
$`g^{ij}`$ is the flat-space metric for spherical coordinates, and $`\mathbf{}`$ will stand for the covariant derivative hereon.

Assuming an equilibrium from purely azimuthal (but differential) bulk motion V= $`\mathrm{\Omega }\mathrm{𝟏}_\phi `$, an Eulerian perturbation of such a state $`(V^r,V^\phi ,V^z)(v^r,\mathrm{\Omega }+v^\phi ,v^z)`$, leads to the usual equations associated with a rotating frame and its coriolis and centrifugal terms. For coordinate motion, Euler equations for the perturbations of the fluid read

$`(_t+\mathrm{\Omega }_\phi )v^r2r\mathrm{\Omega }v^\phi `$ $`=`$ $`g^{rj}\delta 𝐟_j`$ (8)
$`(_t+\mathrm{\Omega }_\phi )v^\phi +(2{\displaystyle \frac{\mathrm{\Omega }}{r}}+\mathrm{\Omega }_{,r})v^r`$ $`=`$ $`g^{\phi j}\delta 𝐟_j`$ (10)
$`(_t+\mathrm{\Omega }_\phi )v^\theta `$ $`=`$ $`g^{\theta j}\delta 𝐟_j`$ (12)
where $`\delta 𝐟`$ stands for the Eulerian perturbation of the sum of specific body forces. The standard form of these equations, e.g., for non-coordinate motion (see, Chandrasekhar 1961), may be obtained from Eqs  above by “dimensionalizing” $`V^\phi `$, (i.e., in the second Eq, multiplying by $`r`$ and completing the differential while recalling that the covariant derivative and the metric commute $`[\mathbf{},g_{ij}]=\mathrm{}`$).

Next, one switches dynamical variables from the Euler velocity perturbation, $`𝐯`$, to the Lagrangian displacement, $`𝝃`$, using the first order relation between Lagrangian and Eulerian variations, $`\stackrel{~}{\mathrm{\Delta }}=\delta +𝝃\mathbf{}`$, whilst denoting333 The tilde indicates that this form of Lagrangian displacement–which is generally non-unique–has had its gauge “fixed” in accordance to the non-relativistic regime. Mathematically, this amounts to a choice of Universal time direction, $`\mathrm{𝟏}_t`$ (e.g., unaffected by the fluid’s motion), while adopting the gauge fixing condition $`\xi \mathrm{𝟏}_t\mathrm{}`$.  $`\stackrel{~}{\mathrm{\Delta }}𝐕\mathrm{d}_t𝝃`$ and $`\delta 𝐕𝐯`$ (see, e.g., Chandrasekhar & Lebovitz 1964, Lynden-Bell & Ostriker 1967)

$$𝐯=\{_t+𝐕\mathbf{}\}𝝃(𝝃\mathbf{})𝐕i\sigma 𝝃\xi ^r\mathrm{\Omega }_{,r}\mathrm{𝟏}_\phi .$$
(13)
The algebraic relation follows from the assumption of differential rotation and from writing exp$`i(\omega t+m\phi +k_zz)`$ dependencies for $`𝝃`$. Note that the connection coefficients in Eq  cancel one another and that $`\sigma \omega +m\mathrm{\Omega }`$ denotes the co-moving frequency of the perturbations.

These geometrical equations have their more traditional equivalents in the so-called shearing sheet approximation where a co-moving, “locally Cartesian frame” $`(\widehat{r},\widehat{\phi },\widehat{\theta })(x,y,z)`$, is used along with the linearized shear velocity field, $`𝐕(x)=[\mathrm{d}_{\mathrm{ln}r}\mathrm{\Omega }]x\mathrm{𝟏}_y`$, to treat the problem locally while introducing the coriolis terms by hand. Defining the Cartesian derivative operator $`\stackrel{}{}`$, then the equivalent to Eq  is $`(_t+𝐕\stackrel{}{})𝝃=𝐯+𝝃\stackrel{}{}𝐕`$, which has Galilean invariance in the sense that Lagrangian time derivatives produce co-moving frequency factors in the dispersion relation: $`(_t+𝐕\stackrel{}{})𝝃d_t𝝃i\sigma 𝝃`$.

The equations of motion (EoM) for the (coordinate) Lagrangian displacement are

$`(\sigma ^2+2\mathrm{\Omega }r\mathrm{\Omega }_{,r})\xi ^r2r\mathrm{\Omega }i\sigma \xi ^\phi `$ $`=`$ $`\delta 𝐟^r`$ (14)
$`\sigma ^2\xi ^\phi +2{\displaystyle \frac{\mathrm{\Omega }}{r}}i\sigma \xi ^r`$ $`=`$ $`\delta 𝐟^\phi `$ (16)
$`\sigma ^2\xi ^\theta `$ $`=`$ $`\delta 𝐟^\theta .`$ (18)
Note that the Eulerian shear term, $`\mathrm{𝟏}_\phi `$, becomes the tide term, $`\mathrm{𝟏}_r`$, in terms of $`𝝃`$.

Let us re-cast these equations in a more compact form

$$\ddot{\xi }^i+2\mathrm{\Gamma }_{jk}^iV^j\dot{\xi }^k2\mathrm{\Gamma }_{jk}^iV^j(v\dot{\xi })^k=\delta 𝐟^i$$
(19)
where each over-dot stands for a factor of $`i\sigma `$ (from a Lagrangian time derivative). In the shearing sheet approximation, these equations correspond to the Hill equations for non-coordinate motion (Chandrasekhar 1961, Balbus & Hawley 1992) $`\ddot{𝝃}+2𝛀\times \dot{𝝃}+2r\mathrm{\Omega }\mathrm{\Omega }_{,r}\xi _x\mathrm{𝟏}_x=\delta 𝐟`$.

### 2.2 compressibility

The Lagrangian perturbation of mass density and the Eulerian perturbation of the field follow from mass and magnetic flux conservations (recall the non-relativistic relation $`\stackrel{~}{\mathrm{\Delta }}=\delta +𝝃\mathbf{}`$)

$$\frac{\mathrm{\Delta }\rho }{\rho }=\mathbf{}𝝃,\delta 𝐁=\mathbf{}\times (𝝃\times 𝐁).$$
(20)
The latter equation includes possible gradients of the background field $`𝐁\mathrm{}`$; however, in the spirit of examining the instability as a local phenomenon the global structure of the field is ignored herein. The Lorentz force variation (co-moving frame) may then be written as

$`\delta \left({\displaystyle \frac{1}{4\pi \rho }}𝐉\times 𝐁\right)`$ $`=`$ $`v_{_{\mathrm{Alf}}}^2\times \left\{\mathbf{}(\mathbf{}𝝃)+_{\widehat{𝐁}}^2𝝃_{\widehat{𝐁}}\mathbf{}(\mathrm{𝟏}_{\widehat{𝐁}}𝝃)\mathrm{𝟏}_{\widehat{𝐁}}_{\widehat{𝐁}}(\mathbf{}𝝃)\right\}`$ (23)
$`\stackrel{i𝐤}{}v_{_{\mathrm{Alf}}}^2\times \{(k_i\xi _ik_{\widehat{𝐁}}\xi _{\widehat{𝐁}})𝐤+k_{\widehat{𝐁}}^2𝝃\mathrm{𝟏}_{\widehat{𝐁}}k_{\widehat{𝐁}}k_i\xi _i\}`$
where the scalar operator $`_{\widehat{𝐁}}\mathrm{𝟏}_{\widehat{𝐁}}\mathbf{}`$, and $`\mathrm{𝟏}_{\widehat{𝐁}}`$ is a unit vector in the direction of the unperturbed field. Note that the term $`(k_i\xi _ik_{\widehat{𝐁}}\xi _{\widehat{𝐁}})k_{}\xi _{}`$ may be interpreted as a restoring force due to the compression of field lines (distinct from line bending, Foglizzo & Tagger 1995).

The Lagrangian variation of the specific pressure gradient contains two terms (Lynden-Bell & Ostriker 1967): one $`\mathrm{\Delta }\rho ^1`$ and another $`\mathrm{\Delta }\mathbf{}p_{\mathrm{r}_+\mathrm{g}}`$. In terms of the displacement vector, the first term is proportional to the equilibrium value of $`\mathbf{}p_{\mathrm{r}_+\mathrm{g}}`$ which is negligible444 the variation of the mass density is also negligible when the focus is on the effects of radiative heat conduction: Loss of pressure support out of compressive modes involves only the pressure term (Araya-Góchez & Vishniac 2002).  in the local treatment (proportional to a radial gradient). For the same reason, the Eulerian and Lagrangian variations of the pressure “force” are identical.

The thermodynamic pressure term is then given by

$$\delta \left(\frac{1}{\rho }\mathbf{}p_{\mathrm{r}_+\mathrm{g}}\right)=\mathrm{\Gamma }\frac{p_{\mathrm{r}_+\mathrm{g}}}{\rho }\mathbf{}(\mathbf{}𝝃)\stackrel{i𝐤}{}c_s^2𝐤k_i\xi _i,$$
(24)
where, for heterogeneous media, $`\mathrm{\Gamma }d_{[\mathrm{ln}\rho ]}\mathrm{ln}p`$ represents a generalized adiabatic index (see, e.g., Chandrasekhar 1939, Mihalas & Mihalas 1984).

Putting the above equations together, one gets the EoM in Fourier-space

$$\ddot{\xi }^i+2\mathrm{\Gamma }_{jk}^iV^j\dot{\xi }^k2\mathrm{\Gamma }_{jk}^iV^j(v\dot{\xi })^k=g^{ij}[\{c_s^2𝐤+v_{_{\mathrm{Alf}}}^2(𝐤k_{\widehat{𝐁}}\mathrm{𝟏}_{\widehat{𝐁}})\}(𝐤𝝃)+v_{_{\mathrm{Alf}}}^2(k_{\widehat{𝐁}}^2𝝃k_{\widehat{𝐁}}\xi _{\widehat{𝐁}}𝐤)]_j.$$
(25)
which agrees with the matrix de-composition of Foglizzo & Tagger (1995) in the case of a purely toroidal field embedded in a gas with adiabatic index $`\mathrm{\Gamma }=1`$.

Reckoning of fluid compressibility has complicated somewhat the equations of motion. Yet, these generally unwieldy equations simplify greatly in the regime of fastest growth (a.k.a. the horizontal regime) and for two ideal field topologies of interest. When the field is meridional, the fastest growth modes have $`\xi _{\widehat{𝐁}}\mathrm{}`$, and $`𝐤k_{\widehat{𝐁}}\mathrm{𝟏}_{\widehat{𝐁}}`$ thus yielding a simple isotropic elastic response $`v_{_{\mathrm{Alf}}}^2(ik_{\widehat{𝐁}})^2𝝃`$.

Alternatively, when the field is purely toroidal, the meridional component of Eq  yields an anisotropy constraint: $`v_{_{\mathrm{Alf}}}^2(k_{}\xi _{})=c_s^2(𝐤𝝃)`$ (Foglizzo & Tagger 1995), which allows for a straightforward solution in this regime. Defining $`\mathrm{\Lambda }`$ through

$$1\mathrm{\Lambda }=\frac{\mathbf{}𝝃}{ik_{}\xi _{}}=\frac{2\mathrm{\Theta }}{\mathrm{\Gamma }+2\mathrm{\Theta }}$$
(26)
where $`\mathrm{\Theta }p_{_{𝐁_\phi }}/p_{\mathrm{r}_+\mathrm{g}}`$, the dispersion relation out of Eq  reads

$$\widehat{\sigma }^4\{(\mathrm{\Lambda }+1)\widehat{q}_{_{\widehat{𝐁}}}^2+\widehat{\chi }^2\}\widehat{\sigma }^2+\mathrm{\Lambda }\widehat{q}_{_{\widehat{𝐁}}}^2\{\widehat{q}_{_{\widehat{𝐁}}}^2+4\widehat{\mathrm{A}}\}=\mathrm{},$$
(27)
where all frequencies are normalized to the rotation rate, $`\widehat{\mathrm{A}}\frac{1}{2}\mathrm{d}_{\mathrm{ln}r}\mathrm{ln}\mathrm{\Omega }`$ is the Oort A “constant”, $`\widehat{\chi }^24(1+\widehat{\mathrm{A}})`$ is the squared of the epicyclic frequency, and $`\widehat{q}_{_{\widehat{𝐁}}}(𝐤𝐯_{\mathrm{Alf}})/\mathrm{\Omega }`$ is a frequency related to the component of the wave vector along the field (in velocity units).

The non-axisymmetric modes of fastest growth conform with (Araya-Góchez & Vishniac 2002)

$$\widehat{q}_{_{\widehat{𝐁}}}^2=2\widehat{\mathrm{A}}+(\frac{1+\mathrm{\Lambda }}{2\mathrm{\Lambda }})\times \left\{\frac{2\mathrm{\Lambda }\widehat{\mathrm{A}}^2}{𝒟}\right\},\mathrm{where}𝒟1+(\frac{1\mathrm{\Lambda }}{2})\widehat{\mathrm{A}}+\sqrt{1+(1\mathrm{\Lambda })\widehat{\mathrm{A}}},$$
(28)
while the expression in the curly brackets identifies with the negative root of the dispersion relation.

Note that the compressibility of non-axisymmetric modes is imprint on the deviations of $`\mathrm{\Lambda }`$ from unity. From Eq  one reads that the degree of compression of these modes gets stronger with the (toroidal) field strength and, naturally, with a softer equation of state. In the companion paper (Araya-Góchez & Vishniac 2002), we find that when radiation pressure begins to dominate the disk dynamics, an “ultra soft” effective (adiabatic) index accentuates the effects of mode compressibility. On the other hand, setting $`\mathrm{\Lambda }1`$ and re-orienting the field vertically produces the standard (incompressible) dispersion form for the Balbus-Hawley instability of a meridional field in the horizontal regime.

## 3 General Relativistic Effects in the Cowling Limit

In contrast to the Newtonian case, the formulation of a covariant theory of accretion disk oscillations requires more than mere application of the Lagrangian rate of change operator dt (or its relativistic counterpart dτ) to the Eulerian velocity perturbation. This is insufficient to carry out a normal mode analysis because of the freedom associated with the choice of coordinates. It is much more useful and proper to free the eigenmodes from the coordinate representation; treating them rather as intrinsic to the physical system. It is here that a Lagrangian construction comes in handy.

A working covariant definition of the Lagrangian displacement is that of a vector field that moves a fluid element’s world line from its unperturbed position in spacetime to its perturbed one. The fundamental relation between the Lagrangian (following the fluid’s world line) and Eulerian (taken at a fixed coordinate point) variational operators is

$$\mathrm{\Delta }=\delta +_\xi $$
(29)
where $``$ stands for the Lie derivative.

An elemental use of this relation involves particle number conservation (Schutz & Sorkin 1977): with the use of a number flux density $`𝒩^\nu n\sqrt{g}U^\nu `$, such a law reads $`\mathrm{\Delta }𝒩\mathrm{}`$, where $`g=det|g^{\mu \nu }|`$ and $`U^\nu `$ is the four-velocity of the fluid. In the Cowling approximation, $`\delta g\mathrm{}`$, and in the absence of co-moving sources (or sinks) of particles, one gets for the variations of the four-velocity of an ideal fluid:

$$\mathrm{\Delta }U0=\delta U+_\xi U$$
(30)
which demonstrates that Eulerian perturbations of the four velocity, $`\delta Uu`$, obey $`u^\nu =_\xi U^\nu `$.

The connection to the Newtonian limit is recuperated upon identifying $`cU\mathbf{}`$ with the convective (or material) rate of change $`cU\mathbf{}\stackrel{c\mathrm{}}{}(_t+𝐕\mathbf{})`$ so that, with $`\stackrel{~}{\mathrm{\Delta }}𝐕\mathrm{d}_t𝝃`$, one has

$$\stackrel{~}{\mathrm{\Delta }}=\delta +𝝃\mathbf{}$$
(31)
(see, e.g., Chandrasekhar 1964, Lynden-Bell & Ostriker 1967, and compare $`\delta 𝐕𝐯`$ with Eq ).

The fluid particles that constitute a thin accretion disk (with negligible radial inflow) embedded in a Kerr spacetime geometry have unperturbed four velocity $`U^\nu =\gamma (\mathrm{𝟏}_t+\mathrm{\Omega }\mathrm{𝟏}_\phi )`$ where $`\mathrm{𝟏}_t=(1,0,0,0)`$ and $`\mathrm{𝟏}_\phi =(0,0,1,0)`$ are the Killing vector fields of the stationary, axisymmetric geometry and where $`\gamma `$ is the “redshift” factor of the fluid elements at fixed radius $`\gamma =U^t=\mathrm{d}_\tau t`$. In terms of the Lagrangian displacement vector field, each Lie derivative with respect to one of the Killing vector fields of the geometry “brings down” a wavenumber co-factor in the dispersion relation (modulo spatial gradients of the four velocity). This leads to an algebraic relation between $`\xi `$ and $`u\delta U`$ in the case of a differentially rotating fluid:

$`u`$ $`=`$ $`_{\gamma (\mathrm{𝟏}_t+\mathrm{\Omega }\mathrm{𝟏}_\phi )}\xi `$ (32)
$`=`$ $`\dot{\xi }\gamma \xi ^r\mathrm{\Omega }_{,r}\mathrm{𝟏}_\phi U\xi \mathrm{ln}\gamma .`$ (34)
Here $`\dot{\xi }i\sigma \gamma \xi ^\nu `$, with $`\sigma =\omega +m\mathrm{\Omega }`$ the co-moving frequency of the perturbation as measured at asymptotic infinity (see Ipser & Lindblom 1992).

The relativistic generalization of Eq  is found upon projecting $`\xi `$ on the 3-surface perpendicular to the (unperturbed) four velocity. With $`h^{\alpha \beta }U^\alpha U^\beta +g^{\alpha \beta }`$ the projection operator, one has

$$h_\beta ^\alpha u^\beta \widehat{u}^\alpha =i\sigma \gamma \widehat{\xi }^\alpha \gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi ^\alpha $$
(35)
while fixing the gauge freedom associated with the component of the Lagrangian displacement perpendicular to space-like hypersurfaces (Schutz & Sorkin 1977), i.e. along the local “time” direction. Note that the requirement of unit normal for the perturbed velocity, $`U+u`$, under the Cowling approximation fixes the gauge accordingly: $`2U^\alpha u_\alpha =U^\alpha U^\beta \delta g_{\alpha \beta }0`$.

Dynamical conservation laws for an ideal fluid in the presence of a large-scale electromagnetic field are written succinctly through the Einstein-Maxwell equation (Cowling approximation)

$$T_{;\nu }^{\mu \nu }F^{\mu \nu }J_\nu =0,$$

where the first term stands for the matter stress and the second equals the Maxwell stress. The notation is standard fare: $`F`$ is the electromagnetic field tensor and $`J=neU`$ is the four-current. Ideal MHD makes things easy by stating that the electric field in the co-moving frame vanishes everywhere. Since the latter is the contraction of the field tensor with the four-velocity, it follows that $`F_{\mu \nu }J^\nu =0`$, and the four-acceleration from the Maxwell stress vanishes as well (but not its perturbation).

We shall concern ourselves with the material stress first.

Denoting the relativistic enthalpy of the fluid by $`\varrho \rho +\epsilon +p`$, it is straightforward to show that for a non-dissipative, ideal fluid such that $`T^{\mu \nu }=\varrho U^\mu U^\nu +pg^{\mu \nu }`$,

$$T_{;\nu }^{\mu \nu }=\varrho \mathrm{d}_\tau U^\mu +g^{\mu \nu }p_{,\nu }\mathrm{where}\mathrm{d}_\tau U$$
(36)
stands for the generalization of the convective rate of change, i.e. the Lagrangian proper-time derivative. The projection of this equation along the four-velocity states energy conservation while the perpendicular components express conservation of momentum.

The specific Eulerian perturbation of the four-acceleration (normalized to the enthalpy) looks like $`uU+Uu`$, and one can use Eq  to switch the dynamical variable in favor of the projected $`\xi `$’s555 note that with this form of the stress-energy tensor, $`h_\beta ^\alpha \mathrm{d}_\tau U^\beta \mathrm{d}_\tau U^\alpha `$, i.e. the four-acceleration automatically lies in proper space-like hypersurfaces. :

$`\widehat{u}U+h(U\widehat{u}^\alpha )`$ $`=`$ $`h(U\dot{\xi }\dot{\xi }U)+2\dot{\xi }U`$ (37)
$``$ $`\gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi UhU(\gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi )`$ (39)
where the silly hats on the $`\xi `$’s (signifying projected components) were dropped.

The first term on the r.h.s. of Eq  may be readily identified with (the projection of) the Lie derivative of $`\dot{\xi }`$ along $`U`$. Defining $`q\mathrm{𝟏}_t+\mathrm{\Omega }\mathrm{𝟏}_\phi `$ (so that $`U=\gamma q`$), one computes

$`_U\dot{\xi }`$ $`=`$ $`\gamma _q\dot{\xi }q\dot{\xi }\gamma `$ (40)
$`=`$ $`\ddot{\xi }\gamma \dot{\xi }^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi +\dot{\xi }U\mathrm{ln}\gamma U\dot{\xi }\mathrm{ln}\gamma `$ (42)
where $`\ddot{\xi }(i\sigma \gamma )^2\xi `$. Note that the last term disappears upon (re)projection onto proper space-like hypersurfaces.

The second term on the r.h.s. of Eq  is easily evaluated, $`\dot{\xi }U^\alpha =\gamma \dot{\xi }^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi ^_\alpha +\mathrm{\Gamma }_{\mu \nu }^\alpha U^\mu \dot{\xi }^\nu `$, and the third simply involves a projected affine connection. Evaluation of the (non-projected) last term yields four pieces:

$`U(\gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi )`$ $`=`$ $`\xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi U\gamma +\gamma \mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi U\xi ^r`$
$`+`$ $`\gamma \xi ^r\widehat{\mathrm{𝟏}}_\phi U\mathrm{\Omega }_{,r}+\gamma \xi ^r\mathrm{\Omega }_{,r}U\widehat{\mathrm{𝟏}}_\phi .`$
Under the premise of negligible radial motion, in the second piece above $`U\xi ^r\dot{\xi }^r+𝒪(U^r)`$, while the third is $`𝒪(U^r)`$. Likewise, pairing of all terms proportional to the logarithmic gradient of the redshift factor yield same order (negligible) corrections $`(\dot{\xi }\gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi )U\mathrm{ln}\gamma 𝒪(U^r)`$. Moreover, the last term above involves the same connection coefficient as the third term on the r.h.s. of Eq .

When all this is said and done, one gets for the specific Eulerian perturbation of the four-acceleration:

$`\widehat{u}U+h(U\widehat{u})\stackrel{𝒪(U^1)}{}\ddot{\xi }^\nu +2\mathrm{\Gamma }_{\alpha \beta }^\nu U^\alpha \dot{\xi }^\beta 2\mathrm{\Gamma }_{\alpha \beta }^\nu U^\alpha (\widehat{u}\dot{\xi })^\beta .`$ (45)
The resemblance with Eq  is remarkable but not accidental.

The shear(tidal) term is embodied by the third term on the r.h.s.: $`\widehat{u}\dot{\xi }=\gamma \xi ^r\mathrm{\Omega }_{,r}\widehat{\mathrm{𝟏}}_\phi `$. Note the non-trivially hatted unit vector $`\widehat{\mathrm{𝟏}}_\phi ^\nu =h_\mu ^\nu \mathrm{𝟏}_\phi ^\mu =\mathrm{𝟏}_\phi ^\nu +U^\nu U_\phi `$. We evaluate this term first using the standard form of the Kerr metric in the equatorial plane (Boyer-Lindquist coordinates):

$$ds^2=\frac{𝒟}{𝒜}dt^2+r^2𝒜(d\phi \omega dt)^2+\frac{1}{𝒟}dr^2,$$
(46)
with $`\omega 2a/𝒜r^3`$ the rate of frame dragging by the hole and where the metric functions of the radial BLF coordinate are written as relativistic corrections (e.g. Novikov and Thorne 1973):

$$𝒜1+a^2/r^2+2a^2/r^3,\mathrm{and}𝒟12/r+a^2/r^2,$$

in normalized geometrical units ($`c=G=M_{_{\mathrm{bh}}}=1`$).

In expanded form, the projection of the Killing vector associated with the azimuthal symmetry is

$$\widehat{\mathrm{𝟏}}_\phi =[1+\stackrel{~}{\gamma }^2\stackrel{~}{r}v^{\stackrel{~}{\phi }}\mathrm{\Omega }]\mathrm{𝟏}_\phi +\stackrel{~}{\gamma }^2\stackrel{~}{r}v^{\stackrel{~}{\phi }}\mathrm{𝟏}_t,$$

where $`\stackrel{~}{\gamma }=\gamma \sqrt{𝒟/𝒜}`$ is the redshift factor relative to “locally non-rotating observers” (Bardeen, Press & Teukolsky 1972) and $`\stackrel{~}{r}r𝒜/\sqrt{𝒟}`$ is the radius of gyration for the physical velocity in that frame, $`v^{\stackrel{~}{\phi }}=\stackrel{~}{r}(\mathrm{\Omega }\omega )`$.

Evaluation of the tidal term is a bit lengthy but straightforward

$`2\mathrm{\Gamma }_{\alpha \beta }^rU^\alpha (\widehat{u}\dot{\xi })^\beta `$ $`=`$ $`{\displaystyle \frac{4}{r^3}}\left\{\frac{1}{2}\gamma ^2𝒟\mathrm{d}_{\mathrm{ln}r}\mathrm{\Omega }\right\}\left[\left({\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}}a\right)+\stackrel{~}{\gamma }^2\stackrel{~}{r}v^{\stackrel{~}{\phi }}\left(1{\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }}}_+\right)\left(1{\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }}}_{}\right)\right]\xi ^r`$ (47)
where $`\mathrm{\Omega }_\pm =\pm (r^{3/2}\pm a)^1`$ refer to prograde and retrograde circular orbits and where the expression in the curly brackets equals (minus) the shear of the congruence of circular, equatorial geodesics (Novikov and Thorne 1973).

Next, to evaluate the coriolis terms one finds $`\xi ^t`$ from the gauge fixing condition $`\xi U0`$. Accordingly one finds

$`2\mathrm{\Gamma }_{\alpha \beta }^rU^\alpha \dot{\xi }^\beta `$ $`=`$ $`2\gamma {\displaystyle \frac{𝒟}{r^2}}\left({\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}}a+(1a\mathrm{\Omega }){\displaystyle \frac{r^2𝒜(\mathrm{\Omega }\omega )}{1\frac{2}{r}(1a\mathrm{\Omega })}}\right)\dot{\xi }^\phi ,`$ (48)
$`2\mathrm{\Gamma }_{\alpha \beta }^\phi U^\alpha \dot{\xi }^\beta `$ $`=`$ $`2\gamma {\displaystyle \frac{1}{r^4𝒟}}\left({\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}}a+2r^2\mathrm{\Omega }\right)\dot{\xi }^r`$ (50)
$`\mathrm{and}`$ (51)
$`2\mathrm{\Gamma }_{\alpha \beta }^tU^\alpha \dot{\xi }^\beta `$ $`=`$ $`2\gamma {\displaystyle \frac{1}{r^4𝒟}}\left(r^2+a^2a\mathrm{\Omega }(3r^2+a^2)\right)\dot{\xi }^r.`$ (52)
Eq  for the (Eulerian) perturbation of the four-acceleration, $`a^\mu =\widehat{u}U^\mu +h(U\widehat{u}^\mu )`$, was derived for the components of the Lagrangian displacement in a coordinate frame that is fixed with respect to distant stars, i.e. in the Boyer-Lindquist “frame”. However, because the instability is local (at least for weak fields in thin disks) one needs to transform the components of Eq  for manipulation in terms of the local tetrad carried by co-moving observers. This simply involves (matrix) multiplication by the basis vectors of such a tetrad (e.g., Novikov & Thorne 1973). In our notation, the relevant basis vectors are $`e_\alpha ^{\widehat{r}}=1/\sqrt{𝒟}(0,1,0,0)`$ and $`e_\alpha ^{\widehat{\phi }}=\gamma r\sqrt{𝒟}(\mathrm{\Omega },0,0,1)`$ (note that transformation to the local tetrad yields equations for non-coordinate motion, i.e. equivalent to motion in a local Cartesian basis).

Transformation of the $`r`$-component is trivial (since one needs to transform both the acceleration and the displacement vector in the basic EoM below, the radial scale, $`1/\sqrt{𝒟}`$, has no net effect):

$`\sqrt{𝒟}a^{\widehat{r}}`$ $`=`$ $`\ddot{\xi }^r+2\gamma {\displaystyle \frac{𝒟}{r^2}}\left({\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}}a+(1a\mathrm{\Omega }){\displaystyle \frac{r^2𝒜(\mathrm{\Omega }\omega )}{1\frac{2}{r}(1a\mathrm{\Omega })}}\right)\dot{\xi }^\phi `$ (53)
$``$ $`{\displaystyle \frac{4}{r^3}}\left\{\frac{1}{2}\gamma ^2𝒟\mathrm{d}_{\mathrm{ln}r}\mathrm{\Omega }\right\}\left[\left({\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}}a\right)+\stackrel{~}{\gamma }^2\stackrel{~}{r}v^{\stackrel{~}{\phi }}\left(1{\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }}}_+\right)\left(1{\displaystyle \frac{\mathrm{\Omega }}{\mathrm{\Omega }}}_{}\right)\right]\xi ^r.`$ (55)
On the other hand, using the local azimuthal base vector, Eqs  and the gauge fixing condition $`\xi U=0`$, computation of the local $`\phi `$-component, $`\mathrm{\Omega }a^t+a^\phi `$, is a bit more involved (again, the radial scale factors out in the EoM and does not affect the dispersion relation)

$$\frac{1}{\gamma r\sqrt{𝒟}}a^{\widehat{\phi }}=\left(\frac{1\frac{2}{r}(12a\mathrm{\Omega })r^2𝒜\mathrm{\Omega }^2}{1\frac{2}{r}(1a\mathrm{\Omega })}\right)\ddot{\xi }^\phi 2\gamma \frac{1}{r^4𝒟}\left(\frac{\mathrm{\Omega }}{\mathrm{\Omega }_+\mathrm{\Omega }_{}}a+(1a\mathrm{\Omega })\mathrm{\Omega }(3r^2+a^2)\right)\dot{\xi }^r.$$
(56)
Let us take a look at the Maxwell stress next.

The fundamental premise of ideal MHD may be stated rather succinctly: in the rest frame of the fluid currents will flow uninhibited to (instantaneously) cancel any hint of an electric field. A relativistic generalization of ideal MHD may be achieved by a similar covariant (albeit imperfect) postulate: $`E_{\mathrm{r}.\mathrm{f}.}=FU\mathrm{}`$, e.g., the Faraday field tensor is “purely magnetic” in the fluid frame.

Such postulate brings a few mathematical consequences (Phinney 1983):
i- the second electromagnetic invariant vanishes everywhere,

$$\frac{1}{4}F_{\mu \nu }^{\mu \nu }\mathrm{}$$
(57)
ii- the Faraday field tensor is Lie transported along the worldlines of the fluid,

$$_UF\mathrm{}$$
(58)
iii- the other zero eigenvector of the field tensor is the (space-like) magnetic field $`BU`$ with $`^{\mu \nu }ϵ^{\mu \nu \alpha \beta }F_{\alpha \beta }`$ the dual to the field tensor,

$$F(U)\mathrm{}$$
(59)
and iv- the field tensor is also invariant when transported along the magnetic field four-vector:

$$__UF\mathrm{}$$
(60)
Note further that, since the four-velocity and the four-magnetic field are orthogonal, $`UB=0`$, properties ii and iv above define a 2-surface of invariance for the Faraday tensor

$$_{_{aU+bB}}F\mathrm{}$$
(61)
where a and b are arbitrary real numbers.

Aside from the intrinsic (physical) difference in their spacetime orientation, the mathematical similarities between $`U`$ and $`B`$ are uncanny.

Let us go back to the Newtonian case for a moment. From our definition of the Lagrangian variation of the three-velocity: $`\mathrm{d}_t𝝃\stackrel{~}{\mathrm{\Delta }}𝐕`$, one finds the equation governing the Lagrangian change of three-velocity: $`\stackrel{~}{\mathrm{\Delta }}𝐕=(𝐕\mathbf{})𝝃`$, (e.g. Eq ). As noted in the footnote, the difference between $`\stackrel{~}{\mathrm{\Delta }}`$ and $`\mathrm{\Delta }`$ is related to the choice of gauge for $`\xi ^\alpha `$. The induction equation of non-relativistic MHD yields a virtually identical relation for the magnetic field variation (in a frame where the fluid was originally at rest) which is spoiled by fluid compressibility: $`\stackrel{~}{\mathrm{\Delta }}𝐁=(𝐁\mathbf{})𝝃𝐁\mathbf{}𝝃`$. Nevertheless, making use of the continuity equation and weighting the field by the inverse of the fluid’s mass density $`\stackrel{~}{𝐁}𝐁/\rho `$ cleans up its connection to the displacement vector field $`\stackrel{~}{\mathrm{\Delta }}\stackrel{~}{𝐁}=\stackrel{~}{𝐁}\mathbf{}𝝃`$. If only conservative forces, $`Uf0`$, act on the fluid, it can be shown that use of energy conservation in lieu of mass conservation simply swaps the fluid’s rest mass density by the relativistic enthalpy (a world scalar), above. Thus, we choose to work below with a specific measure of the magnetic four-vector weighted by the inverse of the fluid’s relativistic enthalpy $`\stackrel{~}{B}\frac{1}{\varrho }(U)`$. Such combination of observables (and its perturbation) occurs naturally in the problem at hand.

Applying the Lagrangian variational operator, c.f. Eq , on $`\stackrel{~}{B}`$ under the constraints from ideal MHD noted above, Eq , yields (contrast this with Eq)

$$\mathrm{\Delta }\stackrel{~}{B}\mathrm{}=\delta \stackrel{~}{B}+_\xi \stackrel{~}{B}.$$
(62)
This equation states a manifestly covariant expression for the Eulerian perturbation of the (enthalpy-weighted) Faraday tensor under ideal MHD constraints: $`\stackrel{~}{b}=_{\stackrel{~}{B}}\xi `$. Demanding that the total magnetic field four-vector be orthogonal to the (unperturbed) four-velocity is equivalent to projecting its Eulerian perturbation into proper spacelike hypersurfaces: $`\stackrel{~}{b}^\mu \widehat{b}^\mu =h_\nu ^\mu _{\stackrel{~}{B}}\xi ^\nu `$. Again, we suppress the hats below while tacitly imposing the condition $`\xi U\mathrm{}`$ throughout.

The Eulerian perturbation of the specific measure of the Lorentz force, $`\delta (\stackrel{~}{T}_{_{\mathrm{e}.\mathrm{m}.}})_{;\nu }^{\mu \nu }=\delta F^{\mu \nu }\stackrel{~}{J}_\nu +F^{\mu \nu }\delta \stackrel{~}{J}_\nu `$, may now be written in terms of the Lagrangian displacement but the general expressions are not particularly illuminating. Evaluated in the frame where the fluid was originally at rest, the Eulerian perturbations of the field tensor and of the four-current depend linearly on the components of $`b`$: $`\delta F=F(b)`$ and $`4\pi \delta \stackrel{~}{J}=\mathrm{d}\delta \stackrel{~}{F}(b)`$.

We proceed by assuming negligible gradients of the background specific field ($`\stackrel{~}{B}\mathrm{}`$):

$$\stackrel{~}{b}=\stackrel{~}{B}\xi \xi \stackrel{~}{B}\stackrel{~}{B}(ik)\xi ,$$
(63)
and, consistent with this assumption, we also ignore the $`\delta F\stackrel{~}{J}`$ term in the perturbation of the Maxwell stress (i.e. gradients of the background field tensor ($`J`$) gentler than those of the perturbations).

The simple “linear poking” of the field tensor may now be written in a manifestly covariant manner

$$F^{\mu \nu }\delta \stackrel{~}{J}_\nu =\frac{1}{\varrho }(Bik)^2\xi ^\mu .$$
(64)
Naturally, evaluation of the elastic response of the field is straightforward in the rest frame of the fluid where one has $`F^{\mu \nu }\delta \stackrel{~}{J}_\nu (v_{_{\mathrm{Alf}}}k_{\widehat{𝐁}})^2\xi ^\mu `$. The only difference with the non-relativistic analog is that the Alfvén speed is now weighted by the relativistic enthalpy of the fluid $`\varrho v_{_{\mathrm{Alf}}}^2\frac{1}{2}F^{\mu \nu }F_{\mu \nu }`$.

We are now all geared up to put together the pieces of the puzzle. In terms of the Lagrangian displacement vector field, the r.h.s.’s of Eqs \[55 & 56\] are to be balanced by the elastic response of the field tensor to the poking by $`\xi `$, c.f. Eq  (note that the radial scales of the transformation into the rest frame of the fluid cancel one another). This balance is locally equivalent to $`a^\mu =q_{\widehat{𝐁}}^2\xi ^\mu `$, i.e. the covariant components of the fluid’s acceleration respond to a force proportional to the displacement vector (with the unnormalized “spring constant” $`q_{\widehat{𝐁}}=(v_{_{\mathrm{Alf}}}k_{\widehat{𝐁}})`$ provided by the field). By construction, both of these vectors are orthogonal to $`U`$ and collinear. Furthermore, since $`\ddot{\xi }^\mu (i\gamma \sigma )^2\xi ^\mu `$ and $`\gamma \sigma `$ is a world scalar to be identified with the true co-moving frequency (as measured by an observer riding along with the fluid), it follows that $`\ddot{\xi }^{\widehat{\phi }}(i\gamma \sigma )^2\xi ^{\widehat{\phi }}`$. With these relations and the aforementioned equations for the tidal and coriolis terms, one arrives to lengthy component equations for $`\xi ^r`$ and $`\xi ^\phi `$, for general $`\mathrm{\Omega }U^\phi /U^t`$ and negligible radial flow.

In the case of circular geodesic flow, the equations simplify beautifully (horizontal regime)

$`\ddot{\xi }^r2\gamma {\displaystyle \frac{𝒟}{r^{1/2}}}\mathrm{\Omega }_\pm \left({\displaystyle \frac{r^33r^2\pm 2ar^{3/2}}{r^{3/2}\pm a2r^{1/2}}}\right)\dot{\xi }^\phi {\displaystyle \frac{4}{r^{3/2}}}\left\{\frac{3}{4}\gamma ^2𝒟r^{3/2}\mathrm{\Omega }_\pm ^2\right\}\xi ^r`$ $`=`$ $`q_{\widehat{𝐁}}^2\xi ^r`$ (65)
$`\ddot{\xi }^\phi +2\gamma {\displaystyle \frac{1}{r^{5/2}𝒟}}\mathrm{\Omega }_\pm \left(r^{3/2}\pm a2r^{1/2}\right)\dot{\xi }^r`$ $`=`$ $`q_{\widehat{𝐁}}^2\xi ^\phi .`$ (67)
These immediately yield the sought after dispersion relation near a rotating hole

$$(\gamma \sigma )^4\left[4\gamma ^2\mathrm{\Omega }_\pm ^2\left(𝒞_\pm \frac{3}{4}𝒟\right)+2q_{\widehat{𝐁}}^2\right](\gamma \sigma )^2+q_{\widehat{𝐁}}^2\left[q_{\widehat{𝐁}}^24\left\{\frac{3}{4}\gamma ^2𝒟\mathrm{\Omega }_\pm ^2\right\}\right]=\mathrm{}$$
(68)
where $`𝒞_\pm 1\frac{3}{r}\pm \frac{2a}{r^{3/2}}`$ corresponds to the $`𝒞`$ function of Novikov & Thorne (1973) for prograde orbits.

Factoring out the extrinsic666 As defined, $`\mathrm{\Omega }U^\phi /U^t`$ reflects motion as observed in the Boyer-Lindquist frame, i.e., in a frame extrinsic to the fluid. It follows that the timescale associated with $`\mathrm{\Omega }^1`$ does not reflect a proper dynamical timescale.  dynamical frequency, $`\mathrm{\Omega }_\pm `$, one arrives to the normalized dispersion relation (with $`\gamma \sigma \mathrm{\Omega }_\pm \widehat{\sigma }`$)

$$\widehat{\sigma }^4\left[\widehat{q}_{_{\widehat{𝐁}}}^2+\widehat{\chi }_\pm ^2\right]\widehat{\sigma }^2+\widehat{q}_{_{\widehat{𝐁}}}^2\left[\widehat{q}_{_{\widehat{𝐁}}}^2+4\widehat{\mathrm{A}}\right]=\mathrm{}$$
(69)
where

$$\widehat{\mathrm{A}}\left\{\frac{3}{4}\gamma ^2𝒟\right\}\mathrm{and}\widehat{\chi }_\pm ^2=4\gamma ^2\left(𝒞_\pm \frac{3}{4}𝒟\right)$$

denote the normalized shear parameter and (co-moving) epicycle frequency (note that $`\frac{1}{\gamma }\widehat{\chi }`$ corresponds to the well known result of epicycle frequency as measured at asymptotic infinity). One thus sees that with the proper generalizations of the epicycle frequency and shear parameter, the local dispersion relation is identical with the Newtonian case in the limit of no fluid compression and $`\mathrm{𝟏}_{\widehat{𝐁}}𝝃\mathrm{}`$ (i.e., the “classical” Balbus-Hawley instability, Eq).

Using the relation $`\gamma ^2=(1\pm a/r^{3/2})^2𝒞_\pm ^1`$ for cold, circular, geodesic flow (Novikov and Thorne 1973), one finds the fastest growing modes to conform with

$$\widehat{q}_{_{\widehat{𝐁}}}^2=1\frac{1}{16}\widehat{\chi }^4=1\left(1\pm \frac{a}{r^{\frac{3}{2}}}\right)^4\left\{1\frac{3}{4}\frac{𝒟}{𝒞_\pm }\right\}^2$$
(70)
which remains finite and close to the Newtonian value of $`\frac{15}{16}`$ for all radii outside the ISCO (and for any value of the rotation parameter).

To attach meaning to the polynomial functions that appear naturally in the dispersion relation for the magnetorotational instability, recall the range of radii that define particle dynamics in the Kerr geometry (Bardeen et al. 1972):

1. The marginally stable circular orbit (a.k.a. the ISCO), $`r_{\mathrm{ms}}`$, corresponds to the root of $`\widehat{\chi }_\pm =0`$.
2. The radius of the circular photon orbit, $`r_{\mathrm{ph}}`$, is where $`𝒞_\pm =0`$.
3. The event horizon, $`r_+`$, happens at the outer root of $`𝒟=0`$.

One therefore has the following ordering of radii for any value of the rotation parameter $`a`$: $`r_{\mathrm{ms}}>r_{\mathrm{ph}}>r_+`$. As remarked by Bardeen et al. ’72, when $`a=1`$, the proper radial distance between these radii is non-zero in spite of “coinciding with the horizon”, i.e. in spite of laying at the same Boyer-Lindquist radial coordinate.

Inspection of Eq  now shows that $`\widehat{q}_{_{\widehat{𝐁}}}0^+`$ as $`rr_{\mathrm{ph}}^+`$ so the most unstable MRI modes go to large scale just outside the photon orbit. Moreover, utilizing that expression for $`\widehat{q}_{_{\widehat{𝐁}}}`$ in the unstable root of the dispersion relation, one finds the growth rate (or frequency!) to be given by

$$\widehat{\sigma }^2=\left\{\frac{3}{4}\frac{𝒟}{𝒞_\pm }\right\}^2\left[\left(1\pm \frac{a}{r^{\frac{3}{2}}}\right)^4\frac{8}{3}\frac{𝒞_\pm }{𝒟}\left(\pm 2\frac{a}{r^{\frac{3}{2}}}+5\frac{a^2}{r^3}\pm 4\frac{a^3}{r^{\frac{9}{2}}}+\frac{a^4}{r^6}\right)+\left(\frac{4}{3}\frac{𝒞_\pm }{𝒟}\right)^2\left(4\frac{a^2}{r^3}\pm 4\frac{a^3}{r^{\frac{9}{2}}}+\frac{a^4}{r^6}\right)\right].$$
(71)
For a non-rotating hole,

$$\widehat{q}_{_{\widehat{𝐁}}}0\mathrm{@}r=r_{\mathrm{ph}}(1+\frac{1}{5}),\mathrm{and}\mathrm{the}\mathrm{local}\mathrm{growth}\mathrm{rate}\widehat{\sigma }=\frac{3}{4}\frac{𝒟}{𝒞_\pm }2$$

while for a rotating hole, the MRI quenching radii (for fastest growing modes) also occur just outside the circular photon orbit and may be readily extracted from the above relations. In Figs 1 and 2, we plot the general relativistic modifications the fastest growing linear wavemodes, wavenumbers and growth rates respectively, as functions of radius and for different values of spin parameter $`a`$.

To go beyond this point, one would need to address global effects arising, for instance, from field curvature terms (see, e.g., Curry & Pudritz 1995, Ogilvie & Pringle 1996) and from the non-negligible radial velocity profile. Further investigation of the nature of the global instability is beyond the scope of this paper.

## 4 Discussion

The MRI–in its most simple, local, incompressible variant–is found to operate virtually unabated down to the marginally stable orbit for massive particles. This radius is nearly coincident with the putative inner boundary of standard, thin accretion disks in the Kerr geometry. A vanishing epicycle frequency at $`r_{\mathrm{ms}}`$ means that the fastest growing wavenumbers tend to be of a bit smaller scale, $`\widehat{q}_{_{\widehat{𝐁}}}^2:\frac{15}{16}1𝒪(ar^{3/2})`$, while growing faster than classically, $`i\widehat{\sigma }:\frac{3}{4}1+𝒪(ar^{3/2})`$. The effects of strong gravity become truly significant only in a regime where circular, cold, geodesic flow is unstable (i.e., where $`\widehat{\chi }_\pm ^2<0`$).

Recall that particle trajectories with $`U^\phi /U^t=\mathrm{\Omega }_\pm `$ exist inside $`r_{\mathrm{ms}}`$ and all the way down to $`r_{\mathrm{ph}}`$ but, in the presence of turbulent velocity fluctuations, body forces such as a radial pressure gradient would be required to confine the flow to such circular orbits. Although very little is concretely known about the accretion flow inside $`r_{\mathrm{ms}}`$, two rather robust remarks may be ascertained: The flow inside $`r_{\mathrm{ms}}`$ cannot be supported centrifugally and it must therefore deviate from a standard thin disk. In addition, depending on the timescale for infall, the flow may not have time to cool significantly and advection of entropy will become progressively more important as $`r_+`$ is approached. A robust prediction of this paper is the expectation that free energy tapping from the differential shear flow goes on in the region immediately below $`r_{\mathrm{ms}}`$.

One may envisage the situation inside the ISCO to evolve from a Mildly Advective Accretion Flow (MAAF) to a fully Advection Dominated Accretion Flow (ADAF) as the photon orbit is approached. In fully or partly advective accretion flows, such as those modeled by Popham & Gammie (1998), the angular velocity profile “peaks” precisely at $`r_{\mathrm{ph}}`$ and quickly drops therein to match the angular velocity at $`r_+`$. More importantly, when cooling by advection of entropy is moderately important–say, for advection fractions $`f`$ a few percent–the angular velocity profile departs very slowly from circular geodesic flow, $`U^\phi /U^t\mathrm{\Omega }_\pm `$ down to a region below the marginally bound orbit. The transition from nearly Keplerian to plunging orbits can be clearly seen in one of the very few global slim disk models where the cooling fraction is calculated explicitly: the 1D models of Popham et al. (1998, albeit in the exotic scenario of a hyper-accreting black hole). In these models the radial velocity component is non-negligible when compared with the local sound speed (the sonic point generally occurs below $`r_{\mathrm{ms}}`$, even near $`r_{\mathrm{mb}}`$ for low values of $`\alpha `$), but $`v^{\stackrel{~}{r}}`$ is generally smaller than $`v^{\stackrel{~}{\phi }}`$ down to the region below $`r_{\mathrm{mb}}`$. (Note that the radial speed in the corotating frame, e.g. Gammie & Popham’s (1998) $`V`$, is related to the same in the locally non-rotating frame by $`v^{\stackrel{~}{r}}=\gamma _\phi ^1V`$.)

The major limitation of the work presented herein is the presumption of negligible radial flow which greatly simplifies matters from the onset (see Eq ). At this point, it is unclear how much the results will change when full consideration is made for the radial inflow. Since the changes could be qualitatively significant–recent reports negate the reversal of the centrifugal force when the radial speed overwhelms the azimuthal component (Mukhopadhyay & Prasanna 2001, Prasanna 2001)–this point should be the subject of close scrutiny in a future paper. Meanwhile, the adoption of an angular velocity profile corresponding to circular equatorial geodesic orbits seems a reasonable rough approximation in view of the above observations of advective flows. In this spirit, we argue below that the natural evolution of the MRI inside the marginally stable orbit is at least consistent with this assumption.

Assume, in quasi-linear fashion, that the time- and length-scales provided by the linear dispersion relation reflect the growth and size of the dominant turbulent eddies to within factors of order unity to a few. Provided that $`v^{\stackrel{~}{r}}\text{ }<\text{ }v^{\stackrel{~}{\phi }}`$, simple linear growth/non-linear decay arguments (e.g., Araya-Góchez 1999) can be used to predict a predominantly toroidal field topology: The MRI constantly promotes radial/azimuthal field growth from “horizontal” velocity fluctuations, $`\xi ^{\widehat{r}}\xi ^{\widehat{\phi }}`$, while the coherent, background azimuthal shear flow converts this field into toroidal field at twice the rate of radial field generation. In the rest frame of the fluid, the tapping of free-energy associated with the shear flow becomes very rapid as the flow turns relativistic. Indeed, a co-moving observer measures the shear parameter, $`2\mathrm{A}/\mathrm{\Omega }_\pm `$ to be $`\frac{3}{2}\gamma ^2𝒟\frac{3}{2}(1\pm a/r^{3/2})^2𝒟/𝒞_\pm `$, higher than the “Keplerian” frequency associated with the global dynamical timescale as seen at large distances (the redshift factor comes in because we chose to measure the angular frequency in terms of Boyer-Lindquist coordinates).

The ratio $`𝒟/𝒞_\pm `$ represents a gauge of the relative strength of two inertial terms, shear and coriolis. Setting aside the issue of radial flow for a moment, our dispersion relation suggests that as material approaches the region just outside of the photon orbit where $`𝒞_\pm `$ vanishes, the slow branch of the dispersion relation (i.e. the MRI) is stabilized by the predominance of shear over the coriolis terms. Recall that the location of the circular photon orbit is the place where the centrifugal force reverses its direction: Inside $`r_{\mathrm{ph}}`$, increasing the velocity of a test particle pulls it in further (see, e.g., Abramowicz & Prasanna 1990 and references therein). The limit of $`\widehat{q}_{_{\widehat{𝐁}}}0^+`$ means that what was essentially a local instability becomes a global phenomenon. Although such a regime is formally outside the scope of the local analysis, one can anticipate a few rather interesting qualitative consequences.

At first glance, the dispersion relation Eq  shows the appearance of an interchange, radially buoyant mode (T. Foglizzo priv. comm., Araya-Góchez 1999). More likely, this would simply imply the need for a steep radial stratification profile. Indeed, if the coherence lengthscale of the field were to reach the comoving length associated with the radial scale, $`𝒜r/\sqrt{𝒟}`$, the disk could make a transition from centrifugally driven to magnetically driven: MRI modulated dynamics guarantee that the Alfvén speed associated with the toroidal field at this large scale would be comparable to the orbital speed. Moreover, the field generated at large scales is less prompt to decay through reconnection and also more buoyant. This has very important consequences for the energy fraction going into–and persisting in–electromagnetic channels.

The radial velocity profile will very likely change the expected outcome once the radial velocity becomes supersonic or super-Alfvenic, but some of the qualitative features of the this analysis may carry over when the full problem is solved, analytically or otherwise. If so, in this part of the so-called “plunging region” of the flow, the turbulent eddies will tend to grow larger while the field direction will tend to track the surfaces of null angular velocity gradients (no longer purely toroidal). The implied field topology is that of large-scale horizontal field domains.

### 4.1 effects of radiation stress and neutrino trapping

Precise assessment of the dynamical role of radiation in the general relativistic regime is hampered by the breakdown of one key assumption made to simplify the “linear poking” on the Faraday field tensor: Use of the enthalpy weighted specific four magnetic field in Eq . On the other hand, one expects a photon gas–semi-contained by a neutral plasma through Compton scattering–to comprise a rather funny MHD fluid where the magnetic field is truly frozen only to the co-moving volume associated with the mass density but for which pressure perturbations do not behave adiabatically. It follows that when the fluid transitions into radiation pressure domination, compressive modes (e.g. toroidal field, non-axisymmetric modes) may lose pressure support in an unfavorable range of wavenumber phase-space (Agol & Krolik 1998). One can prove that the MRI falls squarely in such radiative heat conduction damping regime (Blaes & Socrates 2001). Araya-Góchez & Vishniac (2002) show that the behavior of the energy equation is in some (algebraic) sense “quasi-adiabatic” for exponentially growing, non-propagating modes. Mathematically, this means that a real, analytical, slow-varying function of the scale of the perturbations, $`\stackrel{~}{\mathrm{\Gamma }}(i\stackrel{~}{𝐤}^2/\stackrel{~}{k}^0)`$, can be used to treat the energy equation in quasi-adiabatic fashion. Radiative heat conduction isotropizes the modes and, to zeroth order, one can use such quasi-adiabatic index in Eq  to anticipate that the effects of radiative heat conduction out of compressive toroidal modes is to increment the threshold of shear parameter where $`\widehat{q}_{_{\widehat{𝐁}}}0^+`$ from $`2`$ to $`\widehat{\mathrm{A}}1+2𝒟/(1+\mathrm{\Lambda })`$. Nevertheless, since $`\widehat{\mathrm{A}}𝒟/𝒞_\pm `$ and $`𝒞_\pm 0`$ @ $`r_{\mathrm{ph}}`$, the increase in shear threshold in this setting is rather inconsequential.

Note further that the qualitative nature of energy deposition in radiation pressure dominated fluids is insensitive to the details of the (global) cooling but it is explicitly sensitive to the optical thickness of the relevant eddies. Thus, upon the onset of neutrino trapping in the neutrino cooling regime of hyper-accreting black holes, one may reasonably expect MRI modulated dynamics at $`p_\nu \text{ }>\text{ }p_{\mathrm{r}_+\mathrm{g}}`$ (gas and radiation are tightly coupled) to resemble the standard disk case when $`p_{\mathrm{rad}}\text{ }>\text{ }p_{\mathrm{gas}}`$. Turner et al. 2001 report that the non-linear outcome of the MRI in this setting is a porous medium with drastic density contrasts as to cheat the Eddington limit at high accretion rates. Under nearly constant total pressure and temperature, the non-linear regime shows that density enhancements anti-correlate with azimuthal field domains (just as expected from the linear theory) and that turbulent eddies live for about a dynamical timescale while mass clumps are destroyed through collisions or by running through a localized region of shear.

Since the turbulent eddies in the disk are largely instabilities of the toroidal field (at moderate values of the field), large-scale horizontal field domains near the marginally bound orbit would naturally force the baryonic component of the accretion flow into spatially segregated, massive clumps that occur near the nodes of non-axisymmetric (toroidal) MRI eddies (Araya-Góchez & Vishniac 2002). This expectation motivates the picture of massive clumpy accretion suggested in the introduction.

## 5 Ending Notes

In summary, this work shows that the MRI is virtually unaffected by strong gravity outside the innermost stable circular orbit. Secondly, it indicates that the instability becomes non-local inside this region. Indeed, the MRI may leave behind a large scale, ordered field as the fluid heads in towards the circular photon orbit (with an orientation that tracks surfaces of null angular velocity gradients). Assuming incompressibility and the angular velocity profile of circular geodesic flow, the fastest growing modes die off while going to large scales at a radius just inside the marginally bound orbit. Accountability of compressibility as required to address the effects of radiation stress will bring the critical MRI quenching radius in, slightly closer to the photon orbit. Radiation stress, when significant, will diminish the growth rate while increasing the threshold of shear parameter to quench the MRI.

Radial inflow will affect the global field topology but the details depend on poorly understood fluid trajectories in a region where cold, circular geodesic flow is unstable. As was pointed out by Krolik (1999), the standard assumption of ballistic orbits is never self-consistent for ideal MHD accretion inside $`r_{\mathrm{ms}}`$. Indeed, when magnetic turbulence is the culprit of angular momentum transport in the disk, the magnetic field energy density must become comparable to the rest-mass energy density of the fluid in the plunging region. Yet, unlike Krolik’s suggestion, we do not believe that linear Alfvén waves could efficiently transport energy from inside $`r_{\mathrm{ms}}`$; the magnetic field there is still highly unstable, and the range of stability of such waves is limited by inertial forces.

On the other hand, in this paper we demonstrate that energy deposition and angular momentum transport through the MRI go on virtually unscathed in the region just below $`r_{\mathrm{ms}}`$. An important note is the promptitude of this process at near Eddington rates since the MRI directly feeds the photon bath through compressive damping of the modes (Araya-Góchez & Vishniac 2002). Energy deposition into the radiation field thus occurs on the MRI timescale! On the other hand, near $`r_{\mathrm{ph}}`$ the flow will inevitably transition into advective cooling. Assessing the magnetic field dynamics in the region $`r_{\mathrm{ms}}>r\text{ }>\text{ }r_{\mathrm{mb}}`$ is essential to predict the efficiency of accretion, and to address some large scale effects such as jet launching and disk-hole coupling.

At highly super-Eddington accretion rates (such as those expected in the prompt stages of hyper-accreting black hole formation), the fluid may possess non-trivial amounts of internal energy per unit rest mass of baryons. For such a hot MHD fluid, $`r_{\mathrm{ms}}`$ does not represent a significant boundary to the disk/flow and such may occur rather closer to $`r_{\mathrm{ph}}`$. This stresses the importance of addressing MHD processes in the region above the circular photon orbit. Along these lines, we have motivated the provocative conjecture that copious gravitational wave losses ensue through black hole ringing when a hyper-accreting black hole enters the accretion regime where neutrino trapping occurs. This argument, which combines linear regime phenomenology with the latest numerical results from accretion in radiation-stress dominated environs, leads to a picture of near-hole accretion where large-scale horizontal field domains channel the flow into massive clumps that “thump” the hole.

Lastly, note that a strong, toroidal field topology is ripe ground for MHD instabilities that promote poloidal field generation such as the Parker and radial interchange instabilities (in the vertical and horizontal regime respectively). This instabilities could provide a physical justification for desirable field topologies invoked in jet launching and the Blandford-Znajek processes.

## Acknowledgments

It is a great pleasure to acknowledge stimulating and instructive discussions with Lee Lindblom, Roger Blandford, Ethan Vishniac, Thierry Foglizzo and Sterl Phinney. I am also very thankful to an anonymous referee for useful comments which helped improve the manuscript. I am greatly indebted with Caltech’s Theoretical Astrophysics and Relativity Group for their hospitality.