Hydromagnetic Stability of a Slim Disk in a Stationary Geometry

The magnetorotational instability originates from the elastic coupling of ﬂuid elements in orbit around a gravitational well. Since inertial accelerations play a fundamental dynamical role in the process, one may expect substantial modiﬁcations by strong gravity in the case of accretion on to a black hole. In this paper, we develop a fully covariant, Lagrangian displacement vector ﬁeld formalism with the aim of addressing these issues for a disk embedded in a stationary geometry with negligible radial ﬂow. This construction enables a transparent connection between particle dynamics and the ensuing dispersion relation for MHD wave modes. The MRI–in its incompressible variant–is found to operate virtually unabated down to the marginally stable orbit; the putative inner boundary of standard accretion disk theory. To get a qualitative feel for the dynamical evolution of the ﬂow below r ms , we assume a mildly advective accretion ﬂow such that the angular velocity proﬁle departs slowly from circular geodesic ﬂow. This exercise suggests the turbulent eddies will occur at spatial scales approaching the radial distance while tracking the surfaces of null angular velocity gradients. The implied ﬁeld topology, namely large-scale horizontal ﬁeld domains, should yield strong mass segregation at the displacement nodes of the non-linear modes when radiation stress dominates the local disk structure (an expectation supported by quasi-linear arguments and by the non-linear behavior of the MRI in a non-relativistic setting). Under this circumstance, baryon-poor ﬂux in horizontal ﬁeld domains will be subject to radial buoyancy and to the Parker instability, thereby promoting the growth of poloidal ﬁeld.


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 highenergy 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α 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, rms (Agol & Krolik 2000)-both rely on magnetic coupling between a standard disk and the flow inside rms.
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Ṁ > ∼1M⊙ sec −1 , 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 rms 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  solutions for moderate values of α 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 ). 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.

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, ℓ. 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 drℓ 2 = κ 2 ≥ 0, where κ 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 Alf /cs ≪ 1, guarantee that the instability is truly local 1 , occurring at large values of k ≃ Ω/v Alf (≡ k · 1B). 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ϕ = 0 modes, shear evolves the radial component of wavenumbers according to kr(t) = k0r −[d ln r Ω] kϕt, which means that modes that could be "unstable" are only so, transiently. The maximum instantaneous growth rate occurs when kr → 0 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ϕ ≪ kr ≪ k θ , ensue in violent poloidal Alfvénic couplings that promptly take over the dynamics. Nonaxisymmetric 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 α values of a few tenths), and bcompressive, 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: dz ln ρ = H −1 , with gas, radiation (and possibly magnetic) pressures tracking the unperturbed density profile ρHΩ = pr + g + p B . 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).

inertial terms
Inertial accelerations are geometrically imprint in the connection terms for the covariant derivatives of the Eulerian velocity components. For spherical coordinate motion (ṙ,φ,θ) −→ (V r , V ϕ , V θ ), the only non-trivial connections are Γ r ϕϕ ∧ Γ ϕ rϕ . Denoting the Lagrangian time derivative by dt ≡ ∂t + V · ∇, the three components of Euler's equation read g ij is the flat-space metric for spherical coordinates, and ∇ will stand for the covariant derivative hereon. Assuming an equilibrium from purely azimuthal (but differential) bulk motion V= Ω1ϕ, an Eulerian perturbation of such a state (V r , V ϕ , V z ) −→ (v r , Ω + v ϕ , 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 where δf 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 [1] above by "dimensionalizing" V ϕ , (i.e., in the second Eq, multiplying by r and completing the differential while recalling that the covariant derivative and the metric commute [∇, gij] = ∅). Next, one switches dynamical variables from the Euler velocity perturbation, v, to the Lagrangian displacement, ξ, using the first order relation between Lagrangian and Eulerian variations,∆ = δ + ξ · ∇, whilst denoting 2∆ V ≡ dtξ and δV ≡ v (see, e.g., Chandrasekhar & Lebovitz 1964, Lynden-Bell & Ostriker 1967 (3) The algebraic relation follows from the assumption of differential rotation and from writing exp i(ωt + mϕ + kzz) dependencies for ξ. Note that the connection coefficients in Eq [3] cancel one another and that σ .
= ω + mΩ 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" (r,φ,θ) → (x, y, z), is used along with the linearized shear velocity field, V(x) = [d ln r Ω] x1y, to treat the problem locally while introducing the coriolis terms by hand. Defining the Cartesian derivative operator ∂, then the equivalent to Eq [3] is (∂t + V · ∂)ξ = v + ξ · ∂ V, which has Galilean invariance in the sense that Lagrangian time derivatives produce co-moving frequency factors in the dispersion relation: The equations of motion (EoM) for the (coordinate) Lagrangian displacement are Note that the Eulerian shear term, ∝ 1ϕ, becomes the tide term, ∝ 1r, in terms of ξ.
Let us re-cast these equations in a more compact form where each over-dot stands for a factor of iσ (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 ξ + 2 Ω ×ξ + 2rΩΩ,r ξx 1x = δf .

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∆ = δ + ξ · ∇) The latter equation includes possible gradients of the background field ∂B = ∅; 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 where the scalar operator ∇B ≡ 1B · ∇, and 1B is a unit vector in the direction of the unperturbed field. Note that the term (kiξi − kBξB) . = k ⊥ ξ ⊥ 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 ∝ ∆ρ −1 and another ∝ ∆∇pr + g. In terms of the displacement vector, the first term is proportional to the equilibrium value of ∇pr + g which is negligible 3 in the local treatment (proportional to a radial gradient). For the same reason, the Eulerian and Lagrangian variations of the pressure "force" are identical.
Putting the above equations together, one gets the EoM in Fourier-spacë 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 Γ = 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 ξB . = ∅, and k ≃ kB1B thus yielding a simple isotropic elastic response ∝ −v 2 Alf (ikB) 2 ξ. Alternatively, when the field is purely toroidal, the meridional component of Eq [9] yields an anisotropy constraint: (Foglizzo & Tagger 1995), which allows for a straightforward solution in this regime. Defining Λ through where Θ ≡ p Bϕ /pr + g, the dispersion relation out of Eq [9] readŝ where all frequencies are normalized to the rotation rate,Â ≡ 1 2 d ln r ln Ω is the Oort A "constant",χ 2 ≡ 4(1 +Â) is the squared of the epicyclic frequency, andqB ≡ (k · v Alf )/Ω 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) 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 Λ from unity. From Eq [10] 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 Λ . = 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.

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 where L 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 N ν ≡ n √ −g U ν , such a law reads ∆N . = ∅, where g = det |g µν | and U ν is the four-velocity of the fluid. In the Cowling approximation, δg . = ∅, 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: which demonstrates that Eulerian perturbations of the four velocity, δU ≡ u, obey u ν = −L ξ U ν . The connection to the Newtonian limit is recuperated upon identifying c U · ∇ with the convective (or material) rate of change c U · ∇ c→∞ −→ (∂t + V · ∇) so that, with∆V ≡ dtξ, one has (see, e.g., Chandrasekhar 1964, Lynden-Bell & Ostriker 1967, and compare δV ≡ v with Eq [3]). The fluid particles that constitute a thin accretion disk (with negligible radial inflow) embedded in a Kerr spacetime geometry have unperturbed four velocity U ν = γ(1t + Ω1ϕ) where 1t = (1, 0, 0, 0) and 1ϕ = (0, 0, 1, 0) are the Killing vector fields of the stationary, axisymmetric geometry and where γ is the "redshift" factor of the fluid elements at fixed radius γ = U t = dτ 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 ξ and u ≡ δU in the case of a differentially rotating fluid: Hereξ ≡ iσγξ ν , with σ = ω + mΩ the co-moving frequency of the perturbation as measured at asymptotic infinity (see Ipser & Lindblom 1992). The relativistic generalization of Eq [3] is found upon projecting ξ on the 3-surface perpendicular to the (unperturbed) four velocity. With h αβ ≡ U α U β + g αβ the projection operator, one has 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 α uα = − U α U β δg αβ . = 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) 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 = n e U 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µν J ν = 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 ̺ ≡ ρ + ε + p, it is straightforward to show that for a non-dissipative, ideal fluid such that T µν = ̺U µ U ν + pg µν , 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.
Under the premise of negligible radial motion, in the second piece above U·∇ξ r ≃ξ r + O(U r ), while the third is O(U r ). Likewise, pairing of all terms proportional to the logarithmic gradient of the redshift factor yield same order (negligible) corrections (ξ − γξ r Ω,r1ϕ)U·∇ ln γ ≃ O(U r ). Moreover, the last term above involves the same connection coefficient as the third term on the r.h.s. of Eq [19]. When all this is said and done, one gets for the specific Eulerian perturbation of the four-acceleration: The resemblance with Eq [5] is remarkable but not accidental. The shear(tidal) term is embodied by the third term on the r.h.s.:û −ξ = −γξ r Ω,r1ϕ. Note the non-trivially hatted unit 4 note that with this form of the stress-energy tensor, h α β dτ U β . = dτ U α , i.e. the four-acceleration automatically lies in proper space-like hypersurfaces.
c 2002 RAS, MNRAS 000, 000-000 vector1 ν ϕ = h ν µ 1 µ ϕ = 1 ν ϕ + U ν Uϕ. We evaluate this term first using the standard form of the Kerr metric in the equatorial plane (Boyer-Lindquist coordinates): with ω ≡ 2a/Ar 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): A ≡ 1 + a 2 /r 2 + 2a 2 /r 3 , and D ≡ 1 − 2/r + a 2 /r 2 , in normalized geometrical units (c = G = M bh = 1). In expanded form, the projection of the Killing vector associated with the azimuthal symmetry iŝ 1ϕ = [1 +γ 2r vφ Ω]1ϕ +γ 2r vφ1t, whereγ = γ D/A is the redshift factor relative to "locally non-rotating observers" (Bardeen, Press & Teukolsky 1972) and r ≡ rA/ √ D is the radius of gyration for the physical velocity in that frame, vφ =r(Ω − ω). Evaluation of the tidal term is a bit lengthy but straightforward where Ω± = ±(r 3/2 ± 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 ξ t from the gauge fixing condition ξ · U . = 0. Accordingly one finds Eq [21] for the (Eulerian) perturbation of the four-acceleration, a µ =û·∇ U µ +h (U·∇û µ ), 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 [21] 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 er α = 1/ √ D(0, 1, 0, 0) and eφ α = γr √ D(−Ω, 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/ √ D, has no net effect): On the other hand, using the local azimuthal base vector, Eqs [24] and the gauge fixing condition ξ · U = 0, computation of the local ϕ-component, ∝ −Ωa t + a ϕ , is a bit more involved (again, the radial scale factors out in the EoM and does not affect the dispersion relation) 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 r.f. = F · U . = ∅, 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, ii-the Faraday field tensor is Lie transported along the worldlines of the fluid, iii-the other zero eigenvector of the field tensor is the (space-like) magnetic field B ≡ F · U with F µν ≡ ǫ µναβ F αβ the dual to the field tensor, and iv-the field tensor is also invariant when transported along the magnetic field four-vector: Note further that, since the four-velocity and the four-magnetic field are orthogonal, U · B = 0, properties ii and iv above define a 2-surface of invariance for the Faraday tensor 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: dtξ ≡∆V, one finds the equation governing the Lagrangian change of three-velocity:∆V = (V·∇) ξ, (e.g. Eq [15]). As noted in the footnote, the difference between∆ and ∆ is related to the choice of gauge for ξ α . The induction equation of nonrelativistic 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:∆B = (B·∇) ξ − B ∇·ξ. Nevertheless, making use of the continuity equation and weighting the field by the inverse of the fluid's mass densityB ≡ B/ρ cleans up its connection to the displacement vector field∆B =B·∇ξ. If only conservative forces, U · f . = 0, 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 enthalpyB ≡ 1 ̺ (F ·U ). Such combination of observables (and its perturbation) occurs naturally in the problem at hand. Applying the Lagrangian variational operator, c.f. Eq [13], onB under the constraints from ideal MHD noted above, Eq [31], yields (contrast this with Eq[14]) This equation states a manifestly covariant expression for the Eulerian perturbation of the (enthalpy-weighted) Faraday tensor under ideal MHD constraints:b = LBξ. 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:b µ →b µ = h µ ν LBξ ν . Again, we suppress the hats below while tacitly imposing the condition ξ · U . = ∅ throughout. The Eulerian perturbation of the specific measure of the Lorentz force, δ(T e.m. ) µν ;ν = δF µνJ ν + F µν δJν, 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: δF = F (b) and 4π δJ = d · δF (b).
We proceed by assuming negligible gradients of the background specific field (∇B . and, consistent with this assumption, we also ignore the δF ·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 Naturally, evaluation of the elastic response of the field is straightforward in the rest frame of the fluid where one has F µν δJν . = −(v Alf kB) 2 ξ µ . The only difference with the non-relativistic analog is that the Alfvén speed is now weighted by the relativistic enthalpy of the fluid ̺v 2 Alf ≡ 1 2 F µν Fµν . 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 [25 & 26] are to be balanced by the elastic response of the field tensor to the poking by ξ, c.f. Eq [34] (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 µ = −q 2 B ξ µ , i.e. the covariant components of the fluid's acceleration respond to a force proportional to the displacement vector (with the unnormalized "spring constant" qB = (v Alf kB) provided by the field). By construction, both of these vectors are orthogonal to U and collinear. Furthermore, sinceξ µ ≡ (iγσ) 2 ξ µ and γσ 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ξφ ≡ (iγσ) 2 ξφ. With these relations and the aforementioned equations for the tidal and coriolis terms, one arrives to lengthy component equations for ξ r and ξ ϕ , for general Ω ≡ U ϕ /U t and negligible radial flow.
In the case of circular geodesic flow, the equations simplify beautifully (horizontal regime) These immediately yield the sought after dispersion relation near a rotating hole where C± ≡ 1 − 3 r ± 2a r 3/2 corresponds to the C function of Novikov & Thorne (1973) for prograde orbits. Factoring out the extrinsic 5 dynamical frequency, Ω±, one arrives to the normalized dispersion relation (with γσ ≡ Ω±σ) wherê A ≡ − 3 4 γ 2 D andχ 2 ± = 4γ 2 C± − 3 4 D denote the normalized shear parameter and (co-moving) epicycle frequency (note that 1 γχ 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 1B · ξ . = ∅ (i.e., the "classical" Eq[11]). Using the relation γ 2 = (1 ± a/r 3/2 ) 2 C −1 ± for cold, circular, geodesic flow (Novikov and Thorne 1973), one finds the fastest growing modes to conform witĥ which remains finite and close to the Newtonian value of 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): (i) The marginally stable circular orbit (a.k.a. the ISCO), rms, corresponds to the root ofχ± = 0. (ii) The radius of the circular photon orbit, r ph , is where C± = 0. (iii) The event horizon, r+, happens at the outer root of D = 0.
One therefore has the following ordering of radii for any value of the rotation parameter a: rms > r 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 [38] now shows thatqB → 0 + as r → r + ph so the most unstable MRI modes go to large scale just outside the photon orbit. Moreover, utilizing that expression forqB in the unstable root of the dispersion relation, one finds the growth rate (or frequency!) to be given by For a non-rotating hole, qB → 0 @ r = r ph (1 + 1 5 ), and the local growth rateσ = 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.

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 rms means that the fastest growing wavenumbers tend to be of a bit smaller scale,q 2 B : 15 16 → 1 − O(ar −3/2 ), while growing faster than classically, iσ : 3 4 → 1 + O(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χ 2 ± < 0). Recall that particle trajectories with U ϕ /U t = Ω± exist inside rms and all the way down to r 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 rms, two rather robust remarks may be ascertained: The flow inside rms 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 rms.
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 , the angular velocity profile "peaks" precisely at r ph and quickly drops therein to match the angular velocity at r+. More importantly, when cooling by advection of entropy is moderately importantsay, for advection fractions f ≃ a few percent-the angular velocity profile departs very slowly from circular geodesic flow, U ϕ /U t ≃ Ω± 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 rms, even near r mb for low values of α), but vr is generally smaller than vφ down to the region below r mb . (Note that the radial speed in the corotating frame, e.g.  V , is related to the same in the locally non-rotating frame by vr = γ −1 ϕ V .) The major limitation of the work presented herein is the presumption of negligible radial flow which greatly simplifies matters from the onset (see Eq [16]). 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 vr < ∼ vφ, 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, ξr ≃ −ξφ, 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, 2A/Ω± to be 3 2 γ 2 D ≃ 3 2 (1 ± a/r 3/2 ) 2 D/C±, 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 D/C± 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 C± 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 ph , increasing the velocity of a test particle pulls it in further (see, e.g., Abramowicz & Prasanna 1990 and references therein). The limit ofqB → 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 [37] 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, Ar/ √ D, 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.

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 [34]. 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,Γ(ik 2 /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 [12] to anticipate that the effects of radiative heat conduction out of compressive toroidal modes is to increment the threshold of shear parameter whereqB → 0 + from −2 to −Â → 1+2D/(1+Λ). Nevertheless, sinceÂ ∝ D/C± and C± → 0 @ r 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ν > ∼pr + g (gas and radiation are tightly coupled) to resemble the standard disk case when p rad > ∼pgas . 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), largescale 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.

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 rms. 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 rms; 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 rms. 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 ph the flow will inevitably transition into advective cooling. Assessing the magnetic field dynamics in the region rms > r > ∼ r 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, rms does not represent a significant boundary to the disk/flow and such may occur rather closer to r 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.