General Relativistic Force-Free Electrodynamics: A New Code and Applications to Black Hole Magnetospheres

0.45 MB·English

by  Jonathan C. McKinney

Mon.Not.R.Astron.Soc.000,1–14(2005) Printed5February2008 (MNLATEXstylefilev2.2)
General Relativistic Force-Free Electrodynamics: A New
Code and Applications to Black Hole Magnetospheres
⋆
Jonathan C. McKinney
Institute forTheory and Computation, Harvard-Smithsonian Centerfor Astrophysics, 60 Garden Street, MS51, Cambridge, MA 02138, USA
Accepted 2006January18.Received2005Dec8;inoriginalform2005Dec8
6 ABSTRACT
0
0 The force-free limit of magnetohydrodynamics (MHD) is often a reasonable ap-
2 proximation to model black hole and neutron star magnetospheres. We describe a

n general relativistic force-free (GRFFE) formulation that allows general relativistic
a magnetohydrodynamic (GRMHD) codes to directly evolve the GRFFE equations of
J motion.Established,accurate,andwell-testedconservativeGRMHDcodescansimply

9 addanewinversionpieceofcodetotheirexistingcode,whilecontinuingtouseallthe
1 already-developed facilities present in their GRMHD code. We show how to enforce
theE·B=0constraintandenergyconservation,andweintroduceasimplifiedgeneral

1 modelofthedissipationoftheelectricfieldtoenforcetheB2−E2 >0constraint.We
v
alsointroduceasimplifiedyetgeneralmethodtoresolvecurrentsheets,withoutmuch
0
reconnection, over many dynamical times. This formulation is incorporated into an
1
existing GRMHD code (HARM), which is demonstrated to give accurate and robust
4
GRFFE results for Minkowski and black hole space-times.
1
0
Key words: black hole physics, galaxies: jets, gamma rays: bursts, MHD, stars:
6
neutron, magnetic field, outflows, methods: numerical
0
/
h
p
-
o 1 INTRODUCTION metric modes or stability against reconnection in current
r sheets if present.
t General relativistic force-free electrodynamics (GRFFE)
s Onesimpletechniquetoseekstationarysolutionsandto
a is the low-inertia limit of general relativistic mag-
studytime-dependentstabilityistodirectlynumericallyin-
: netohydrodynamics (GRMHD) (see, e.g., Komissarov
v tegratetheGRFFEequationsofmotion.BoththeGRMHD
2002b, 2004a). Neutron star and black hole mag-
i and GRFFE equations of motion can be written as a set
X netospheres exhibit regions of space that are very
of conservation laws that can be directly integrated. Con-
r nearly force-free, and self-consistent moderately realis-
a servative numerical GRMHD methods, such as of HARM
tic quasi-analytic solutions exist that describe the ideal
(Gammie et al. 2003a), use so-called “primitive” quantities
MHD or force-free environment of such systems (see,
(P) to define so-called “conserved” quantities (U), fluxes
e.g., Blandford & Znajek 1977; Contopoulos et al. 1999;
(F),andsource(S)terms.Thetemporalintegrationisdeter-
Goodwin et al. 2004; Beskin & Nokhrina 2005; Gruzinov
minedbysolvingthesetof“conservation”equations,which
2005). It is generally difficult to solve the GRFFE equa-
can bewritten as
tions to find even stationary solutions, except with sim-
plified assumptions, that apply to astrophysical systems.
∂Up(P) ∂Fpi(P)
For example, despite recent progress in studies of neu- = +Sp(P), (1)
tron star magnetospheres, no self-consistent analytic solu- ∂t − ∂xi
tion considers general relativistic effects. Also, the solution
ofBlandford & Znajek(1977)istheonlyself-consistent an- where p labels the conservation equation and i is the spa-
alytic force-free solution that has a realistic Poynting jet. tial index. The set of source terms (Sp) accounts for the
Also, analytic solutions typically assume stationarity and connection coefficients or other sources of mass-energy-
axisymmetry and so rarely address the global or local sta- momentum. Thus conservation is explicitly true as long as
bilityofthesolutionsagainsttime-dependentornonaxisym- the source terms vanish. For an axisymmetric, stationary
metric and as written in HARM, energy and angular mo-
mentum are explicitly conserved to machine error because
the source terms vanish. HARMhas been successfully used
⋆ E-mail:[email protected] to study black hole accretion flows, winds, and Poynting
(cid:13)c 2005RAS
2 Jonathan C. McKinney
jets (Gammie et al. 2003a; Gammie, Shapiro, & McKinney ods directly evolve the electric and magnetic fields (primi-
2004;McKinney & Gammie 2004; McKinney 2005a,b,c). tivequantitiesare E,B ,thefield-evolutionapproach)and
{ }
While GRFFE is mathematically just the no-inertia have to explicitly check that the velocity is less than the
limit of GRMHD, numerical truncation errors limit theuse speed of light (i.e. magnetic energy density is greater than
of GRMHD codes in evolving systems with regions where electric energy density, B2 E2 >0) since the electric and
−
the magnetic energy density greatly exceeds the rest-mass magnetic field are evolved with no constraint on their evo-
energy (see, e.g., McKinney & Gammie 2004; Komissarov lution (Komissarov 2002b; Krasnopolsky et al. 2005). The
2004b, 2005a,b). Conversely, the GRFFE equations of mo- method we describe evolves thedrift velocity and magnetic
tion do not describe the rest-mass motion along field lines field(primitivequantitiesare v,B ,(velocity-evolutionap-
{ }
or the thermal energy of the particles. So long as the iner- proach),explicitlybreaksdownwhenthatconstraintisvio-
tia of particles remains negligible, theGRFFE equations of lated and this ensures a physicaland causal evolution.
motion properly describe the magnetic field geometry and We show that this velocity-evolution approach can
motion. be directly used by GRMHD codes, and explicitly dis-
A GRFFE code can be used to studyneutron star and cuss how this is done for conservative numerical meth-
black hole magnetospheres. For example, the origin of the ods. Many GRMHD codes have recently been devel-
collimation and stability of astrophysical jets remains un- oped (Koide et al. 1999; Komissarov 1999; Gammie et al.
explained. A GRFFE code can be used to examine the ori- 2003a; DeVilliers & Hawley 2003; Shibata& Sekiguchi
gin of the collimation and the stability of any model for a 2005; Anninoset al. 2005; Antonet al. 2005). Our dis-
Poynting-dominated jet. cussion applies to any GRMHD method, but focuses on
Existing jet solutions may not be stable and self- conservative-basedcodes.WeshowthataseparateGRFFE
collimation of isolated jets may not work (Spruit 1996; code does not have to be developed, and many of the tools
Okamoto1999,2000,2003).Collimationbyapredominately and methods used to make the GRMHD code perform well
toroidal field (hoop stress) may lead to a pinch/kink in- carry directly over tothe GRFFE version.
stability, whereas poloidal collimation works only if the 2showshowaconservativeGRMHDcodecanbeused
§
jet is surrounded by an extended disk wind (Spruit 1996; toevolvetheGRFFEequationsofmotion.Simplifiedmodels
Begelman & Li 1994; Begelman 1998). The solution to the of dissipation in GRFFE are discussed in order to handle
collimation problem may be that the disk wind collimates regions heretheelectric field dominates themagnetic field.
the Poynting jet. By using both GRMHD and GRFFE nu- 3 includes a series of standard tests of the GRFFE
§
merical models, the origin of collimation, acceleration, and code.
stability of jets can be examined. Thus it is more efficient 4discussesmodelsofblackholemagnetospheres.This
§
tohaveasinglecodebeabletoperform bothGRMHDand sectiondemonstratestheusefulnessoftheGRFFEequations
GRFFE studies. of motion and themodel of dissipation.
Even a black hole system accreting a thick disk 5 summarizes theresults of thepaper.
§
has a force-free magnetosphere (McKinney & Gammie TheGRMHDnotationfollowsMisner et al.(1973)and
2004). Such a magnetosphere quantitatively agrees with we use Heaviside-Lorentz units unless otherwise specified,
the solution by Blandford & Znajek (1977) for low which is like Gaussian units without the 4π. For exam-
black hole spins and agrees qualitatively for high spins ple, the 4-velocity components are uµ (contravariant) or
(McKinney & Gammie 2004; Komissarov 2004b, 2005a). uµ (covariant). For a black hole with angular momentum
Gamma-ray bursts (GRBs), active galactic nuclei (AGN), J =jGM2/c,j =a/M isthedimensionlessKerrparameter
and x-ray binaries probably exhibit both a Poynting- with 16j 61. The contravariant metric components are
dominated jet and a jet from the disk, as suggested by gµν an−dcovariantcomponentsaregµν.Thecomovingenergy
GRMHD numerical models (McKinney 2005b,c). In partic- densityisb2/2,wherebµ isthecomovingmagneticfield.See
ular,ofalljet mechanisms,theBlandford-Znajek drivenjet Gammie et al.(2003a);McKinney & Gammie(2004)forde-
is the only one that can clearly produce an ultrarelativis- tails.
tic jet (McKinney 2005a). Thus is it important to study
theBlandford-Znajekprocessindetail.AGRFFEcodecan
help determine the stability of the Poynting-dominated jet
generated by the Blandford-Znajek effect. For example, a
2 FORCE-FREE ELECTRODYNAMICS
GRFFE code was used to studythe puremonopole version
of the Blandford-Znajek solution, which is found to be sta-
This section shows how a conservative-based algorithm de-
ble (Komissarov 2001, 2002a). However, general Poynting-
signed for solving the ideal GRMHD equations of motion
dominated jets produced by the Blandford-Znajek effect
based upon the velocity and a magnetic field can be used
may be violently unstable to pinch and kink instabilities
tosolvetheGRFFEequationsof motion. Inparticular, the
(see, e.g., Li 2000, but see also Tomimatsu et al. 2001). A
electric and magnetic field are used to derive the velocity
GRFFEcodeprovesvaluableinstabilityanalysesbyavoid-
of the frame in which the electric field vanishes. This al-
ing the difficulty of analytically crossing the light cylinder
lowsaconservativeidealGRMHDnumericalcodetoevolve
(see, e.g., Tomimatsu et al. 2001), which manifests itself as
the force-free equations of motion by simply modifying the
a singularity in theGrad-Shafranov equation for stationary
inversion from conserved quantities to primitive quantities.
solutions.
Theforce-free electrodynamics equationsofmotion are
ThegoalofthispaperistoformulatetheGRFFEequa-
the8 Maxwell equations:
tions of motion such that a conservative GRMHDcode can
be used with small modifications. Present numerical meth- Fµν = Jν, (2)
µ
∇ −
(cid:13)c 2005RAS,MNRAS000,1–14
General RelativisticForce-FreeElectrodynamics:A New CodeandApplicationsto BlackHole Magnetospheres 3
where Fµν are the components of the Faraday tensor and HARMusestheprimitivequantitiestodefineso-called“con-
Jµ are thecomponents of the current density,and served” quantities (U) and the fluxes (F), which are both
∗ µν closed-form operations. These conserved quantities can be
F =0, (3)
∇µ evolved forward in time, but then an inversion to primitive
where ∗Fµν 1ǫµναβF are the components of the dual quantities is required to easily define the fluxes and other
Faraday(or≡Ma2xwell) teαnβsor,Fµν 1ǫµναβ∗F ,ǫαβδγ required quantitiesto determine thenext update.
≡−2 αβ ≡ NoclosedformsolutionappearstoexistinGRMHD,so
( 1/√ g)[αβγδ], and ǫ √ g[αβγδ], where brack-
αβδγ
− − ≡ − mostinversionmethodsrelyoniterativeproceduressuchas
etsdenotetheunitlengthcompletelyantisymmetrictensor.
Newton’s method. This approach can be used in force-free
These8equationsdefine6evolution equationsand2differ-
electrodynamics, but may not always work due to known
ential constraints.
problems(suchaspoorlyconditionedorsingularJacobians)
By the antisymmetry and duality of the Faraday and
withNewton’smethod.However,aclosed-formsolutiondoes
Maxwell tensors, theybe expressed in component form as
exist for force-free electrodynamics, as described below.
Fαβ fηαEβ fηβEα hBγηδǫαβγδ (4) Theonly conserved quantitiesof relevance in force-free
≡ − −
electrodynamics are the lab-frame momentums Tt (the en-
and µ
ergy evolution equation giving the µ = t term is actually
∗Fαβ hηαBβ+hηβBα fEγηδǫαβγδ, (5) redundant) and the lab-frame magnetic field Bµ (only 3
≡− −
componentsare independent).If onecould obtain Eµ,then
where f and h are arbitrary independent constants
onecouldreconstructtheFaradayfromequation(4)orany
that we set to be f = h = 1 consistent with the
other quantities from Eµ and Bµ. It is straight-forward to
conventions in Misner et al. (1973), and see also, e.g.,
show that if EαB =0, then
Baumgarte & Shapiro (2003). Here Eα η Fαβ, Bα α
β
ηβ∗Fβα, and ηµ is an arbitrary 4-vector≡. This gives tha≡t Eα=ǫαβγδBβTγξηξηδ/(B2η4), (10)
Fµν∗Fµν = 4EµBµ( η2) and FµνFµν = 2(B2 E2)( η2). asshownbyasubstitutionofequation(6)andequation(4)
For a time-like ηµ,−the sign conventions are a−s in sp−ecial into this expression. Only the spatial components of η Tµ
µ ν
relativity. andBµareneededifonechoosesaspecialformofηµ.Notice
Theelectromagneticstressenergytensorisconstructed that EαB = 0 is explicitly true, therefore the degeneracy
α
from theFaraday tensor and must be quadraticin the field condition of E B = 0 can be preserved to machine error
·
strengths, symmetric, and is divergenceless in vacuum by regardlessofthetruncationerrorinTorB.Alsonoticethat
Maxwell’s equations. This uniquetensor is equation (10) projects out a component perpendicular to
1 time,thespatialfield,andthemomentums.Forafixedfield
Tµν =FµλFνλ− 4gµνFλκFλκ, (6) Bi,only 2componentsoftheelectricfield areindependent.
Onemaychoosetohaveη2 = 1suchthatηµ onlyde-
wheregµν arethecomponentsofthemetric.Thedualitybe- −
pendsonthemetric.Onechoiceisη = α,0,0,0 ,where
tween the Faraday and Maxwell tensors and thedefinitions α 1/√ gtt. Then ηµ =(1/α) 1, µβi {,−where βi} α2gti.
of Eµ and Bµ givethat ≡ − { − } ≡
Thisdefinesazeroangularmomentum(ZAMO)frame.This
B2+E2 choice of ηµ makes it possible to only require Tt and Bi to
Tµν = (ηµην+Pµν) (7) i
2 obtain Eα in equation (10).
(cid:18) (cid:19)
η2 (BµBν+EµEν) Another interesting choice for ηµ is to have Eν ≡
− − η Fµν = 0. In the ideal GRMHD equations of motion, for
µ
η(cid:0)αEβ(cid:1)Bκ ηµǫναβκ+ηνǫµαβκ , a fluid velocity uµ, this choice corresponds to the electric
−
fieldinthecomovingframeuµ beingeν u Fµν =0.Then
(cid:16) (cid:17) µ
where the projection operator Pµν = ( η2)gµν + ηµην. η2 = 1 since the fluid velocity is tim≡e-like. In force-free
For example, the energy density in frame−ηµ is ηµηνTµν = electro−dynamics, there is no unique4-velocity that satisfies
η4(E2 +B2)/2, which for a time-like ηµ is the same ex- eν = 0, but one such frame is constructed below that is
pression as in special relativity. If the electromagnetic field uniquely always time-like with a minimum Lorentz factor
is the only source of stress-energy, then equations (2) are with respect to theframe with 4-velocity ηµ.
equivalenttotheenergy-momentumconservationequations Asshownnext,anytwoframeswith4-velocitiesuµ and
Tµν =J Fµν =0, (8) ηµ can be easily related to determine a 4-velocity of the
∇µ µ frame in which eν = 0 for the Faraday defined in terms
where only 2 of the 4 components of equation (8) are inde- of ηµ. Thus the metric and Faraday alone determine a 4-
pendentbecause equations(8) implies velocity that allows oneto use theideal GRMHD Faraday,
Fµν∗Fµν =4EµBµ(−η2)=0. (9) Fαβ ≡−bγuδǫαβγδ (11)
For more details see Komissarov (2002b, 2004a). and Maxwell,
∗Fαβ uαbβ+uβbα, (12)
≡−
2.1 Inversion where bν u ∗Fµν. This formulation and sign conventions
µ
≡
Conservative numerical GRMHD methods, such as of arethesameasusedinHARM.InHARM,anew4-velocity
HARM, operate primarily on so-called “primitive” quan- is introduced that is unique by being related to a physical
tities (P): fluid density, fluid internal energy, coordinate observer for any space-time and has well-behaved interpo-
fluid velocity, and the lab-frame coordinate magnetic field. lated values,
(cid:13)c 2005RAS,MNRAS000,1–14
4 Jonathan C. McKinney
u˜i ui γηi, (13) Thesecondterminequation(17)representsthepurelyspa-
≡ −
tial “E B drift.” Note that there is no evolutionary con-
where γ = −uαηα and so ut = γ/α. This additional term straint ×on Tµν that forces B2 E2 > 0, and when this
representsthespatialdrift oftheZAMOframe definedear- −
is violated the force-free model is no longer physical. The
lier. Onecan show that γ =(1+q2)1/2 with q2 ≡giju˜iu˜j. value of vφ coincides with the “field rotation frequency”
To obtain the 4-velocity,notice that
Ω =F /F =F /F for stationary axisymmetric flows
F tr rφ tθ θφ
0=eα= (v Eβ)ηα+Eα+v B η ǫαβγδ, (14) for which vr =vθ =0.
− β β γ δ Now the inversion from Tt,Bi Ei,Bi and then
where vβ uβ/(uαηα). The most general form of the 4- vi,Bi completely define{seiquati}on→s({11)and}(12)for a
≡ − →{ }
velocity that satisfies theabove is general relativistic force-free electrodynamics evolution us-
ing a conservativeGRMHD code.
ǫ ησEξBη η
vβ =G − βσξηη2B2 +H ηβ2 Obviously this formulation can be also used to study
(cid:18) − (cid:19) (cid:18)− (cid:19) specialrelativisticmodelsaswell.Noticethatinspecialrel-
E B ativitythatthederivedvelocityexpressionreducestotheso-
+K ηβ2E2!+L ηβ2B2!, (15) called “electromagnetic 3-velocity”v=E×B/B2 =S/B2,
− − where S=E B is thePoyntingflux.
where G,pH, K, and L arpe functions to be determined. Asusedi×nHARM,thisformulationpreserves B=0
∇·
Eachterm representsoneoffourorthogonaldirections, and and E B=0 to round-off error for both theGRMHD and
·
when multiplied by arbitrary functions represent the most GRFFE equations of motion. Notice that this differs from
general solution. The only nontrivial term is the antisym- otherformulationsthatonlypreserve B=0andE B=0
∇· ·
metric product between the last term in equation (14) and to truncation error (Komissarov 2002b, 2004a)
the first term in equation (15). Substitution of this vβ The fact that the field-evolution approach does not
into equation (14) gives that G = 1. Since vβηβ = 1, break down might be considered an advantage when seek-
then H = 1. With uβ = γvβ, then u2 = 1 gives t−hat ing stationary solutions. In such a case the evolution may
−
γ = η2B2/((1 K2)B2 E2), which simply defines have regions that only transiently have B2 E2 < 0 and
γ = −puα−ηα as the L−orentz fa−ctor between the two frames. the integration can pass smoothly through t−his region into
All terms proportional to K are orthogonal to each other a physical solution space. A related advantage of the field-
and to Eα, and so in general K =0. Hence, uαEα =0. To evolutionapproachisthatRunge-Kuttatemporalevolution
determineL,noticethat bα in∗Fµν in equation(12) can be can recover the correct temporal order of accuracy. That
written as is,withoutcharacteristicinterpolation,HARMusesRunge-
PαBβ Kutta to time step to achieve higher order temporal accu-
bα= β , (16) racy. Runge-Kutta is only first order accurate for the first
γ
substep, but after all substeps are completed, the method
wherePα=δα+uαu istheprojectiontensor.Thusbα has is accurate to arbitrary order. The velocity-evolution ap-
β β β
terms proportional to Bα and uα. So extra terms added to proach can yield an unphysical result for a substep and be
uα proportional to Bα vanish due to the antisymmetry of unabletocontinueortreattheresultasaviolationofforce-
the Maxwell and so do not contribute to the stress-energy freeelectrodynamics,whilethefield-evolutionapproachcan
tensorortheequationsofmotion.Thus,thefunction Lpa- avoid such first order errors and recover to arbitrary order
rameterizes the arbitrary velocity component along a field accuracy.However,wehavefoundthevelocity-evolutionap-
line,andwechooseL=0tominimizetheLorentzfactor of proachtobesufficient.Thisissueoffallingoutsidethelight
the frame. Hence, uαB = 0 and so bα = Bα/γ. With this cone is the same issue one encounters when evolving the
α
choiceofL,ifB2 E2>0,thentheframewithuµ isalways GRMHD equations of motion, and in that case the author
−
time-like. Thus theframe definedby the4-velocity knows of no method that does not require the velocity at
somepoint duringtheintegration, so theissueis treated as
uα = B2 ηα ǫαβγδηβEγBδ (17) a generic one that simply requires a more accurate integra-
rB2−E2!(cid:18) − B2 (cid:19) tion.
istime-likeforanyforce-freeelectrodynamicsolution.Thus,
byconstruction,wehaveshownthatinforce-freeelectrody- 2.2 Currents
namics that it is possible to boost into a time-like frame
Thissectionshowsthatthecurrentscanbecomputedwith-
where the electric field vanishes and thusthe Poynting flux
out time derivatives, which is numerically convenient to
vanishes (see also, e.g. Komissarov 2002b). This also shows
avoidstoringdataatprevioustimes.InidealMHDorforce-
thatforce-freeelectrodynamicsisacausallimit ofGRMHD
free electrodynamics there are many dependent ways of
as long as B2 E2 > 0. This 4-velocity also represents a
− equally describing the same physics. Researchers often in-
unique covariant definition of the “field-line velocity,” and
voke, such as in discussions of current closure, the current
also describes the field-linevelocity even in ideal MHD.
and the magnetic field as a more intuitive set of quantities
A GRMHD code may only need the coordinate lab-
than the electromagnetic fields. The current could be com-
frame3-velocity.Sinceut =γηt,thenfortheearlierdefined
putedfrom
ZAMOframeηµ,thecoordinatelab-frame3-velocityisgiven
by Jα Fαβ, (19)
≡ ;β
ui [ijk]E B whichapparentlyrequirestimederivatives.Inforce-freeelec-
vi = βi+α2 j k (18)
≡ ut − √ gB2 trodynamics, however, the current (like the 4-velocity) sits
−
(cid:13)c 2005RAS,MNRAS000,1–14
General RelativisticForce-FreeElectrodynamics:A New CodeandApplicationsto BlackHole Magnetospheres 5
in the null space of Fαβ, i.e. JαFβ =0. Hence, JαE =0. then explicitly treats the cell interface discontinuity appro-
α α
Thesameanalysis as insection 2.1 leadstothesameresult priatelyandreducestoatrivialformforsmoothflows.How-
for Jα except the function L is no longer arbitrary, where ever,whilethehydrodynamicequationsofmotion allow ar-
here bitraryleftandrightinitialstates,theelectromagnetic field
must obey general well-known jump conditions as a mani-
J Bβ η2
L= β − , (20) festation of the Bianchi identities and the antisymmetry of
(cid:18) ρq,η (cid:19) r B2 ! theFaradaytensor(see,e.g.,chpt.15ofMisner et al.1973).
where ρ Jαη , which only actually requires spatial Alternativelystated,electrodynamicscanbewrittenasaset
q,η α
≡ ofHamilton-Jacobi equationsinvolvingthevectorpotential
derivativesifη = α,0,0,0 .Pluggingequations(4)into
µ
equations (19) and{c−ontracting}J with Bβ gives and the electric field. Arbitrary interpolation of arbitrary
β
quantitiesdoesnotgenerallyenforcethesejumpconditions.
JβBβ =Bβ Eαηβ;α ηαEβ;α (Bγηδǫβαγδ);α , (21) Thus, these jump conditions must be self-consistently en-
− −
forced byinterpolatingappropriatequantitiesintheappro-
(cid:16) (cid:17)
wherethefirstandlasttermsonlyrequirespatialderivatives priate way.
with thechosen ηµ. Using HARM directly operates on the magnetic field rather
∗ αβ than the magnetic vector potential, so the scheme must
F =0, (22)
;β enforce the electrodynamic jump conditions on the elec-
and contracting with E , the second term in equation (21) tric and magnetic fields. For schemes that only interpolate
α
can bewritten as in space and not time, one only requires the spatial jump
conditions. The jump conditions across a one-dimensional
−BβηαEβ;α=BαEβηβ;α−Eβ;αEγηδǫβαγδ. (23) space-like surface for a single lab-frame time (t) are ob-
Now, for the chosen η , each term in equation (21) only tained by integrating Maxwell’s equations. For simplicity
requires spatial derivatµives. So the current can be written, define √ gE˜α tβFαβ and √ gB˜α tβ∗Fβα. Also, con-
− ≡ − ≡
only actually requiringspatial derivatives, as sider three arbitrary orthogonal space-like vectors A,B,C
that describe a space-like volume and the time-like (often
Jα =ρ ηα ǫαβγδηβEγBδ +Bα JβBβ , (24) but not always Killing) vector tα = 1,0,0,0 . For the ho-
q,η − B2 B2 mogeneous Maxwell equations, { }
(cid:18) (cid:19) (cid:18) (cid:19)
where finally ∗ µν
(√ gF )
JβBβ = EαBβ(ηβ;α+ηα;β) 0=Z −√−g ,νdΣµ, (29)
+ (Bα;βBγ−Eα;βEγ) ηδǫαβγδ , (25) where ∗Fij = [ijk]E˜k. First, let dΣµ = ǫµαβγtαdAβdBγ be
where the asymmetry between E(cid:16)α and B(cid:17)α just expresses a 1-form 2-volume for a single lab-frame time (see, e.g.,
the asymmetry in Maxwell’s equations for Jα. One can use Lichnerowicz 1967, 1976; Anile 1989). Then equation (29)
gives for the two spatial parallel components that
equation (17) to write
Jα =uα ρq,η +Bα JβBβ . (26) ∆C √−gE˜AdA[CAB]=0, (30)
γ B2 Z
(cid:18) (cid:19) (cid:18) (cid:19)
wheretheC-directionischosenasperpendiculartothesur-
That is, there is a perpendicular drift current (J ) and a
⊥ face and never is a sum implied for [CAB], which only
field-aligned current (J ), i.e.
k gives the 3-signature for arbitrary, but fixed, A,B,C cor-
Jα =Jα+Jα. (27) responding to any 3 spatial directions. Choosing instead
⊥ k
dΣ =ǫ dAαdBβdCγ,onehasforthecontravariantfield
Equation (26) is Ohm’s law in dissipationless GRFFE. In µ µαβγ
that
the special relativistic regime this reduces to a lab-frame
current density of ∆ √ gBCdAdB=0. (31)
C
−
E B( E)+B(B ( B) E ( E)) Z
J= × ∇· · ∇× − · ∇× , (28)
B2 Likewise, for theinhomogeneous Maxwell equations,
where η Jα =ρ = αJt = E is the charge density (√ gFµν)
in the frαame moqv,iηng w−ith 4-ve−lo∇cit·y ηµ. As in general, the JµdΣµ = −√ g ,νdΣµ, (32)
Z Z −
special relativistic equation obviously has no time deriva-
tives. where Fij =[ijk]B˜k. With dΣµ =ǫµαβγtαdAβdBγ, then
∆ √ gB˜ dA=[CAB] √ gJBdAdC, (33)
C A
2.3 Jump Conditions − −
Z Z
Conservative schemes are often based upon 1D piece-wise for a possible surface current JB δ(C)KB, where upper
≡
constantGodunovmethodsthatachievehigherthanfirstor- (lower) A,B,C denotes the contravariant (covariant) com-
der accuracy by relying on a one-dimensional interpolation ponents parallel to the surface. Choosing instead dΣµ =
from, e.g., grid cell centers to cell interfaces. The interpo- ǫµαβγdAαdBβdCγ for space-like orthogonal A, B, and C,
lated states areassumed tobeapproximately constant over one hasfor thecontravariant field that
thecellface,otherwiseageneralizedRiemannproblemmust
∆ √ gECdAdB= η Jµ√ gdAdBdC, (34)
be solved (see, e.g. Toro 1999). The 1D Godunov scheme C − µ −
Z Z
(cid:13)c 2005RAS,MNRAS000,1–14
6 Jonathan C. McKinney
forapossiblesurfacechargeρ =η Jµ δ(C)σ ,where enforce energy conservation and compare the nonconserva-
q,η µ q,η
≡
σ isthesurfacechargeintheframemovingwith4-velocity tiveintegration withaconservativeonetogaugetheactual
q,η
ηµ and the C-direction is taken to beacross the surface. effect of thetruncation error.
Foraninfinitesimalsurface√ gdAdB,lengths√ gdA One requires a solution for Tt as a function of an ar-
− − i
and √ gdB, or for point values of the fields, these expres- bitrary set of three components of Tt. Then an arbitrary
− µ
sions reduceto choice can be made to move the truncation error from en-
ergy to a momentum. For axisymmetric, stationary space-
∆ E˜ [CAB]=0, (35)
C A times the natural choice of momenta are the radial and θ
∆CBC =0, (36) momenta. This relationship between Tµt that only depends
otherwise on Bµ can be obtained from
and
E2+B2
∆CB˜A=[CAB]KB, (37) Tνµηµην =η4 2 , (39)
∆ EC =σ , (38) where equation (10) also gives that
C q,η
awnhticvhs.acpoavratrifarnotmcotmhepocnonenctesr,nawndithE˜,√B˜−vgs,.tEhe,Bcoinstridaevnatrii-- E2= −ηαBηβ2ηT6γαTδβ!(cid:18)gγδ− ηηγ2ηδ − BBγB2δ(cid:19). (40)
cal to the Cartesian special relativistic expressions for the
jumpconditions.Noticethatnodistinctionismadebetween Equation(40)doesnotactuallydependonTtt forourchoice
boundandfreechargesorcurrents,sotheseexpressionsare ofηµ.Equation (39)givesTtt intermsofTit andBi.Notice
generally true. thatcontractingwithonlyην showsthatηνTνµ justreduces
the 1T-hDeceeqlulaintitoenrsfa(c3e5s.)Fthorrocuegnhte(r3e8d)smchuestmbese tphreesecrovnetdinuat- tthoaηtαoBnβlyTβαtw=o0arfeorinνde6=pet.ndTehnitsfsohrowasfitxheadtBofi.thHeotwherveeer,Tfito’rs
itiescanbeenforcedbyusingthesameinterpolationstencil anarbitrarymagneticfield,theevolutionofallthreeTit’sis
for the left/right interpolations to the cell interface from required to avoid singular expressions.
theleft/rightcellcenters.Thisisthemethodchosenforthe To replace any particular Tit with Ttt, one solves the
GRFFE version of HARM. These constraints on the fields quadraticequation(39)forthatspatialcomponentinterms
formanimplicitconstraintonthedriftvelocityandthemag- ofTttandtheremainingspatialcomponents.Thisnewsetof
netic field. Notice that equations (37) and (38) enforce no effectiveTit’scanthenbepluggedintoequation(10).Notice
specificconstraintunlessthereisanenforcedsurfacecharge that in general the quadratic equation gives two solutions.
or surface current. Also notice that these constraints must This degeneracy is introduced by using the energy, which
also be preserved in ideal GRMHD. lacks directional information. In practice it is sufficient to
Interpolatingtheelectricandmagneticfieldcanleadto use the solution closest to the one obtained originally from
unphysicalinterfacestatesifB2 E2 <0isunconstrainedat integrationofonlythespatialparts(Tit).Thisprocedurecan
the interface. To avoid this, one−can interpolate γE˜ and γ beusedtokeepenergyandangularmomentumconservedto
α
separately andthenreconstructE˜ foragiveninterpolated machine error, unlike in prior methods (Komissarov 2002b;
α
γ at the interface. This guarantees that the interface state Krasnopolsky et al.2005).Caremustbetakenfornumerical
is physical and reasonably similar to thecenter states. methods with a large truncation error in the energy, since
An alternative scheme can be designed with a stag- theevolutionofthespatialandtemporalstress-energycom-
gered grid that automatically enforces these constraints ponents may be disparate. The origin of this larger trunca-
(Del Zannaet al.2003).However,bothmethodsinvolvethe tion error is thelarger nonlinearity of the energy compared
samenumberofinterpolations sinceintheircasetheymust tothemomentums.Lackofenergy-momentumconservation
interpolate the interface field components to the center be- can also beduetodissipative processes, asdescribed in the
fore performing the inversion from conserved to primitive next section.
quantities. They must also use a stencil that guarantees no
discontinuities at their cell center. Thus both methods are
equivalent. 2.5 Dissipation in Force-Free Electrodynamics
There is no evolutionary constraint in dissipationless force-
free electrodynamics that preserves B2 E2 > 0, whose
2.4 Energy Conservation −
violation is taken as evidence that the plasma being mod-
The formulation above for the inversion in force-free elec- elled would have a nonnegligible inertial back-reaction on
trodynamics only explicitly requires Tt and not Tt, which the electric field. This typically occurs in current sheets or
i t
is a similar feature to other methods (Komissarov 2002b; inregionswheretheinertiawouldrestrictthefieldtohavea
Krasnopolsky et al.2005).Sincethetruncationerrorineach velocityofv<c.Aphysicalsystemrespondsbydissipating
Tt is independent, the conserved quantity associated with theelectric field into other forms of energy.
µ
energy conservation (Tt) is generally inconsistent with Tt. Thedissipation oftheelectricfieldisdeterminedbyan
t i
Hence, energy is only conserved to truncation error. As for Ohm’s law. For magnetospheres with an ample supply of
any numerical scheme that only conserves energy to trun- charges (i.e. not “charge-starved”) the Ohm’s law in force-
cation error, this error can be used to gauge the reliability freeelectrodynamicsiswell-approximatedbyalargeconduc-
of the integration. However, in steady-state problems this tivity along the magnetic field and a vanishing conductiv-
truncationerrormaybesecularandbuild-upandleadtoun- ity perpendicular to the magnetic field (Komissarov 2004a,
realisticsolutions.Itisfruitfultoformulatetheinversionto 2005b). This reduces Ohm’s law to the condition E B=0
·
(cid:13)c 2005RAS,MNRAS000,1–14
General RelativisticForce-FreeElectrodynamics:A New CodeandApplicationsto BlackHole Magnetospheres 7
andtheonly perpendicularcurrentisthedrift currentwith drifts.Thisisusefulinstudyingsystemsforwhichtheeffect
velocity given by equation (18). of reconnection is uncertain and the current sheet may be
To model this dissipation, Komissarov (2005b) use the stable.Theaboveprescriptioncanbegeneralized tosetany
prescriptionthatifB2 E2<0,thentheyintroducealarge arbitrary drift speed into the current sheet. Thus, physical
−
cross-field conductivity. They also have to modify the drift models of the current sheet and reconnection speed can be
velocitytokeepv<c.Theproblemwiththisprescriptionis includedinforce-freemodelsaslongasinertiaplaysnoother
that there is no dissipation until v >c. Indeed, in currents role than in thesheet.
sheets the rate of dissipation is related to the drift velocity At present this is only implemented for a priori known
and is allowed to be v c, leading to the fastest possible locations of thecurrentsheet,although ageneral algorithm
→
reconnection rate. should be similar to reconnection-capturing methods (see,
Lyubarsky (2005) study the relativistic Sweet-Parker e.g., Stone& Pringle 2001). In the tests below that have a
sheet and Petschek configuration and determine that the current sheet, there are 4 numerical grid zones in a current
rate of reconnection is much less than the speed of light, sheet region that are forced to obey the above condition.
contrary to previous estimates (see, e.g. Lyutikov 2003; Thisguaranteesthatthestencil,usedbythereconstruction
Lyutikov& Uzdensky 2003). This suggests that dissipation and other dissipative procedures in HARM, does not cou-
should limit thedrift of field into a current sheet. plequantities diffusively across the sheet. This also ensures
First,inertiallossesduetothedissipationoftheelectric that upon convergence testing that this numerical scheme
field in a relativistic flow, such as beyond thelight cylinder to avoid reconnection in current sheets plays no role in the
ofarotatingcompactobject,canbe“captured”bylimiting results.
the Lorentz factor of the drift velocity. This is a similar This formulation ensures that B2 E2 > 0 and that
−
approach taken in GRMHD numerical models in order to current sheets can be forced to be mostly stable, unlike in
avoid significant numerical errors at large Lorentz factors. Komissarov (2002b); Krasnopolsky et al. (2005) and is im-
A simple prescription is to limit γ uαη such that 1 6 provedcompared toKomissarov (2005b).
α
γ2 6 γ2 . From 1 = uαu , this≡m−eans replacing B2 in
equatiomna(x18) (and−the similαar B2 in the second term of
equation 17) with 2.6 Quasi-GRMHD and Stationary Force-Free
( η2)E2B2 WhilethispaperdescribesaGRFFEformulation, sincethe
P2 = −1 1 (41) samecodealsohasaGRMHDformulation,onecanimagine
s − γm2ax hybrid schemes that integrate the decoupled equations of
onlywhen 1>γ2 >γm2ax. Thisgives alimited 3-velocity of motion in the stiff regime where b2/ρ0 ≫ 1, where ρ0 is
therest-massdensity.AllGRMHDnumericalschemeshave
vi = βi+α2[ijk]EjBk, (42) difficulties with this regime, so a decoupled evolution may
− √ gP2 prove useful for studying the first nonzero order effect of
−
inertiainaforce-freefield.Thishybridmethodwillbeused
that has a continuous transition between standard force-
in future work, but themethod is outlined here.
free electrodynamics and dissipative electrodynamics. The
In the force-free limit, one may still retain the evolu-
4-velocitygivenbyequation(17)is,byconstruction,limited
tion of the rest-mass and internal energy density with no
to always be time-like and have γ 6 γ . The energy-
max
backreactionontothefield(seealsoMestel & Shibata1994;
momentum lost in such a limiting procedure is assumed to
Contopoulos 1995; Contopoulos et al. 1999). The method
begainedbyinertialmassintheformof,e.g.,(rapidlylost)
is to evolve the full GRMHD equations of motion, but to
thermalorfield-aligned kineticenergy that hasnoeffecton
determine the field-perpendicular velocity components and
thefield.
fieldfromtheforce-freeequationsalone.Inthiscaseonecan
Second,currentsheetsdissipateduetoadvectionoffield
readilyobtainthefield-alignedvelocitycomponentfromthe
intothesheetandsubsequentreconnection andcancelation
GRMHDinversion.Thatis,inidealGRMHDingeneral,the
of the magnetic field. Even under exactly symmetric con-
fluidvelocitymaybebrokenintoafield-perpendicularanda
ditions, numerical round-off error quickly leads to random
field-alignedvelocityorequivalentlyintoanelectromagnetic
velocity drifts across the current sheet. This magnetic field
and a matter velocity,
advectioncanbelimitedoravoidedbylimitingthedriftve-
locity perpendicular to the current sheet. For example, the uα =uα +uα =uα +uα . (44)
FL ⊥ k EM MA
3-velocitygivenbyequation(42)canbefurthermodifiedto
have a small or vanishing component perpendicular to the Alsoofinterestisthefluidmotioninastationaryforce-
current sheet. That is, if ni is the spatial normal vector of free field, which allows a study of the Lorentz factor along
the current sheet plane at a particular time-slice, then we a force-free field line (Mestel & Shibata 1994; Contopoulos
can set 1995; Contopoulos et al. 1999). For an stationary, axisym-
metric model the energy (Bernoulli) equation alone deter-
njvigij =0 (43) minesthefield-alignedvelocityandthefrozen-inconditions,
withinsomeinfinitesimalregionboundingthecurrentsheet. urFL = uθFL = uφFL−ΩFutFL, (45)
ThischangesOhm’slawinequation(27)tohaveavanishing Br Bθ Bφ
conductivityperpendiculartofieldlines,i.e.thespatialpart
ofJαvanishessuchthatthe“E Bdrift”vanishesalongthe apply to thefluid velocity and by uEαMBα=0,
⊥ ×
nsiogrnmifiacladnitrercetcioonnnoefctthioencubryreanvtosidhienegt.aTnhoimseaxlopuliscintluymaevroiicdasl urEM = uθEM = uφEM −ΩFutEM, (46)
Br Bθ Bφ
(cid:13)c 2005RAS,MNRAS000,1–14
8 Jonathan C. McKinney
apply to the electromagnetic velocity. These equations can is used for all tests except an additional current sheet test
besolved for theelectromagnetic 3-velocity to obtain with parabolic spatial interpolation.
Notice that apart from the initial condition given in
tµ+Ω φµ
vi =Ω φi BiB F (47) Komissarov (2002b,2004a),forthefastwaveonerequiresa
EM F − µ(cid:18) B2 (cid:19) relation determined earlier in their paper: Ez =C µfBy,
wheretµandφµarethetimeandφKillingvectorsandB = where C is constrained such that B2 E2 <0 (C−=1 was
µ −
gµνBν. That is, for stationary, axisymmetric solutions the chosen) and µf =+1 was chosen by them.
magnetic field and Ω alone determinetheelectromagnetic Figure 1 shows the suite of tests demonstrated in
F
velocity, unlikein general as given by equation (18). Komissarov (2002b).Forthetop3panels,thesolid linede-
notestheanalyticsolution attheinitialandfinaltime.The
diamonds denote the numerical solution at the final time.
The top panel shows B2 for the fast wave. The head of the
wave is more resolved than in Komissarov (2002b), but the
3 ALGORITHM TESTS
tail of the wave is slightly less resolved. The second and
A parameter space test of the analytic GRFFE inversion third panels from the top show B2 and B3 for the degen-
described above was performed to evaluate the range of al- erate Alfv´en wave problem. The code does well to capture
lowed Lorentz factors the inversion is capable of handling. both fast and Alfv´en waves. The next two panels show B2
The procedure is to start with a known vi and Bi for a and E1 for thethree-waveproblem. The waves are resolved
specific point in space-time, determine Tt, and then invert similarlyasinKomissarov(2002b).Thebottompanelshows
µ
to get vi. Double precision is used for all quantities. The b2/B2 =(B2 E2)/B2 forthesmoothed B2 E2 0test
− − →
typicalfailure modeisthatthevelocityisdeterminedtobe problem. Overall, the GRFFE version of HARM performs
space-like, and this is nearly coincident with a significant comparably to the code by Komissarov (2002b) based on
increase in theerrorin theinversion.Ofinterest isthegen- more complicated Riemann solvers.
eral maximum value of the Lorentz factor below which no Figure 2 shows the suite of tests demonstrated in
failures occur. Komissarov (2004a). The top panel shows Bz for a station-
In Minkowski space-time tests, the maximum Lorentz ary Alfv´en wave. This wave is much more resolved than in
factorisut 1010 beforeroundofferrorcausestheinversion Komissarov (2004a) and indicates that the code’s effective
to suddenly∼have significant error. Below this ut, the rela- resistive diffusion coefficient is quite low. The next panel
tiveerrorissimilartomachineerror.Asanextremetest,the shows Ez and By for a current sheet problem as described
space-time point is chosen in Kerr-Schild coordinates with inKomissarov(2004a)withB0 =0.5.Thefeaturesarewell-
a Kerr spin parameter of a/M = 0.9375 for a point on the resolved. The two bottom panels show the second current
horizon at θ = π/4. For grid-aligned flows the maximum sheetproblemwithB0 =2.Theupperofthetwoisusingthe
Lorentz factor is ut 1010 as before. However, for arbi- MC limiter, while the lower of the two is using a parabolic
trary flow directions,∼the inversion can fail for ut & 2200. interpolation (Colella 1984), which gives similar results to
Most astrophysical flows of interest have ut <2000, so this Komissarov (2004a). In the fast wave region the MC lim-
is sufficient for our purposes. Otherwise a smaller machine iter gives a Lorentz factor of Γ = 1.90, while the parabolic
precision should be used. method gives Γ = 2000, the largest allowed Lorentz factor.
The GRFFE formulation is coupled to HARM to test Foreithermethod,theregionatx1=0.5withinthecurrent
the formulation’s ability to handle typical force-free prob- sheet reaches B2 E2 0, but the Lorentz factor limiter
− ∼
lems. The GRFFE formulation is used in HARM to per- (herelimited toγmax =2000) keepsthecodestabledespite
formtheMinkowskispace-timetestcalculationsthatarede- thepresenceofthecurrentsheet.Nocomplicateddissipation
scribed by Komissarov (2002b, 2004a). These tests include: model had to be included to achieve such a result.
1) a fast wave ; 2) a degenerate Alfv´en wave ; 3) a three-
wave problem ; 4) a problem that evolves to B2 E2 <0 ;
−
5) a standing Alfv´en wave ; and 6) two current sheet prob-
4 PHYSICAL MODELS
lems.SinceourGRFFEformulationassumesE B=0,their
·
non-degenerate Alfv´en wave test is not considered. This section considers astrophysical models for which the
Their notation of the “wave frame” fields E′ and B′ GRFFE approximation is a reasonable one. The GRFFE
are equivalentto theHARMnotation for eµ=0 and bµ for formulation is used in HARM.
E′ = 0, and otherwise the wave frame, E′, and B′ can be Komissarov (2004a) study the Wald (Wald 1974) and
processed through the inversion routine to setup an initial split-monopole BZ77 (Blandford & Znajek 1977) solutions
lab-frame Bi and vi from their E′, B′, and speed µ. The for slowly and moderately rapidly rotating black holes. Of
otherproblemstheysetuponlyrequirethelab-frameE and particular interest is whether the Wald and split-monopole
B, and again our GRFFE inversion can be used to obtain solutions can be better represented compared to as shown
the initial lab-frame Bi and vi as long as the problem is in figures 3 and 4 of Komissarov (2004a). We consider the
degenerate. Note that our sign conventions for E and B monopoleandsplit-monopoleforslowlyrotatingblackholes
agree with theirs. We estimate that they used about 200 in this paper.
grid zones for their tests and so use the same number for The other models are considered in a separate paper
our tests. We use the same final time, box size, and plot (McKinney2005e),whichdemonstratesourGRFFEformu-
labels for easy comparisons. We shift the x= 0 position to lation’s ability to handle the current sheet in the actual
besimilarforeachsetoftestsinKomissarov(2002b,2004a). split-monopole solution of Blandford & Znajek (1977) and
ACourantfactorof0.9andthemonotonizedcentrallimiter theWaldproblemwitharapidlyrotatingblackholespin.In
(cid:13)c 2005RAS,MNRAS000,1–14
General RelativisticForce-FreeElectrodynamics:A New CodeandApplicationsto BlackHole Magnetospheres 9
Figure1.SuiteoftestsdescribedbyKomissarov(2002b).Solidcurvesrepresenttheinitialandfinalanalyticsolutionforthetopthree
panels.Diamondsrepresentthenumericalsolutionatthefinaltime,exceptfortheplotofb2 thatshowsbothtimes.Testsareasfollows
from top to bottom: 1) B2 for fast wave ; 2) B2 and B3 for degenerate Alfv´en wave ; 3) B2 and E1 for three-wave problem ; and 4)
initialandfinalsmoothed b2=B2−E2→0problem.
thatpaperfortheWaldproblem,wewereabletoreachsim- In this paper we focus on slowly rotating black hole
ilarsolutionsKomissarov (2005a)whousedanMHDmodel magnetospheres that require general relativity and so
to avoid significant reconnection in the current sheet that full GRFFE. The interesting quantities that are plotted
developed in their force-free models Komissarov (2004a). throughout the following sections are Ω F /F , which
F tθ θφ
≡
is also Ω = F /F for a stationary, axisymmetric flow.
F tr rφ
Neutron star magnetospheres are studied in a sepa-
This quantityoften appears asa ratio totheblack hole an-
ratepaper(McKinney2005d),whichdemonstratesthatthe
gular velocity of ΩH a/(2r+). The radial and θ magnetic
GRFFE code can be used to study pulsar magnetospheres ≡
field strengths for these plots is defined as in Komissarov
eveninthepresenceofacurrentsheet.Thatpapershowswe (2004a), with Bi ∗Fit. Also interesting is the conserved
are able to avoid the problems encountered by Komissarov ≡
(2005b)withforce-freeandthecurrentsheet.Wefoundsim-
ilar results to, but more accurate than, their ideal MHD
results.
(cid:13)c 2005RAS,MNRAS000,1–14
10 Jonathan C. McKinney
Figure 2. Suite of tests described by Komissarov (2004a). Solid curves represent the initial and final analytic solution except in the
bottom panel. Pluses and diamonds represent the numerical solution at the final time. Tests are as follows from top to bottom: 1) Bz
forstationaryAlfv´enwave;2)Ez andBy forcurrentsheetwithB0=0.5;3,4)Ez,b2,andBy forcurrentsheetwithB0=2.Thirdis
with2ndordermethodandfourthiswith4thordermethod.
toroidal flux of √ gB˜ B ∗F = √ gFrθ.1 This of a purely rotational velocity Ω = Ω leads to a null tra-
φ φ tφ F
is because the elec−tromag≡netic e≡nergy flux is−Fi = Ti = jectory, i.e. from uαu = 1,
E − t α −
BiΩ B andtheelectromagneticangularmomentumflux
−is Fi =F Fφi/Ω . Also of interest is the magnetic vector po- gtt+2Ωgφt+Ω2gφφ =0. (48)
L E F
tential (A ), whose contours are plotted and represent flow
φ
surfaces. For axisymmetric, stationary flows these surfaces 4.1 Black Hole Magnetosphere: Blandford-Znajek
define surfaces of constant ΩF, Bφ, FEi/Bi, and FLi/Bi in Monopole and Split-Monopole
dissipationlessforce-freeelectrodynamics.Finally,alsoofin-
terest are the light surfaces, which are defined for the case Komissarov (2002b) demonstrated the stability of the pure
monopole solution of Blandford & Znajek (1977) by study-
ing one hemisphere. Our GRFFE formulation generates
1 Noticetheslightchange,forsimplicity,innotationforBφ and quantitatively similar results. Krasnopolsky et al. (2005)
Bi fromthispointonward. also use a method based on HARM and study the depen-
(cid:13)c 2005RAS,MNRAS000,1–14