On the other hand, advances in satellite instrumentation, e.g. the Rossi X-Ray Timing Explorer (RXTE), and the Advanced Satellite for Cosmology and Astrophysics (ASCA), are greatly stimulating, and are guiding theoretical research on accretion physics. The recent discovery of kHz quasi-periodic oscillations in low-mass X-ray binaries extends the frequency range over which these oscillations occur into timescales associated with the innermost regions of the accretion process (for a review see ). Moreover, in extragalactic sources spectroscopic evidence (broad iron emission lines) increasingly points to (rotating) black holes being the accreting central objects [215, 111, 42].
Additionally, most of the proposed theoretical models to explain -ray bursts involve a, possibly hyper-accreting, black hole at some point of the evolutionary paths : neutron star-neutron star and black hole-neutron star binaries, collapsars, black hole-white dwarf binaries, or common envelope evolution of compact binary systems. Developing the capability of performing accurate numerical simulations of time-dependent accretion flows in regions of strong gravitational fields, possibly dynamic during the hyper-accreting phase, is indeed of enormous interest.
Accretion theory is primarily based on the study of (viscous) stationary flows and their stability properties through linearized perturbations thereof. A well-known example is the solution consisting of isentropic constant angular momentum tori, unstable to a variety of non-axisymmetric global modes, discovered by Papaloizou and Pringle  (see  for a review of instabilities in astrophysical accretion disks). Establishing the nature of flow instabilities requires highly resolved and accurate time-dependent non-linear numerical investigations in strong gravitational fields. Such simulations have only been attempted, in general relativity, in a few cases and only for inviscid flows.
For a wide range of accretion problems, the Newtonian theory of gravity is adequate for the description of the background gravitational forces (see, e.g., ). The extensive experience with Newtonian astrophysics has shown that explorations of the relativistic regime benefit from the use of model potentials. Among those, we mention the Paczynski-Wiita pseudo-Newtonian potential for a Schwarzschild black hole , which gives approximations of general relativistic effects with accuracy of in regions from the black hole larger than the marginally stable radius, which corresponds to three times the Schwarzschild radius. Nevertheless, for comprehensive numerical work, a full (i.e. three-dimensional) formalism is required, able to cover also the maximally rotating hole. In rotating spacetimes the gravitational forces cannot be captured fully with scalar potential formalisms. Additionally, geometric regions such as the ergo-sphere would be very hard to model without a metric description. Whereas the bulk of emission occurs in regions with almost Newtonian fields, only the observable features attributed to the inner region may crucially depend on the nature of the spacetime.
In the following we present a summary of illustrative time-dependent accretion simulations in relativistic hydrodynamics. We concentrate on multidimensional simulations. In the limit of spherical accretion, exact stationary solutions exist for both Newtonian gravity  and relativistic gravity . Such solutions are commonly used as test-beds of time-dependent hydrodynamical codes, by analyzing whether stationarity is maintained during a numerical evolution [95, 125, 64, 184, 17].
Wilson's formulation has been extensively used in time-dependent numerical simulations of disk accretion. In  (see also ) Hawley and collaborators studied, in the test-fluid approximation and axisymmetry, the evolution and development of non-linear instabilities in pressure-supported accretion disks formed as a consequence of the spiraling infall of fluid with some amount of angular momentum. Their initial models were computed following the analytic theory of relativistic disks presented by Abramowicz et al. . The code used explicit second-order finite difference schemes with a variety of choices to integrate the transport terms of the equations (i.e. those involving changes in the state of the fluid at a specific volume in space). The code also used a staggered grid (with scalars located at the cell centers and vectors at the cell boundaries) for its suitability to difference the transport equations. Discontinuous solutions were stabilized with artificial viscosity terms.
With a three-dimensional extension of the axisymmetric code of Hawley, Smarr, and Wilson [94, 95], Hawley  studied the global hydrodynamic non-axisymmetric instabilities in thick, constant angular momentum accretion gas tori, orbiting around a Schwarzschild black hole. Such simulations showed that in radially wide, nearly constant angular momentum tori, the Papaloizu-Pringle instability saturates in a strong spiral pressure wave, not in turbulence. In addition, the simulations confirmed that accretion flows through the torus could reduce and even halt the growth of the global instability.
Igumenshchev and Beloborodov  have performed two-dimensional relativistic hydrodynamical simulations of inviscid transonic disk accretion onto a rotating (Kerr) black hole. The hydrodynamical equations follow Wilson's formulation but the code avoids the use of artificial viscosity. The advection terms are evaluated with an upwind algorithm which incorporates the PPM scheme  to compute the fluxes. Their numerical work confirms analytical expectations: (i) The structure of the innermost disc region strongly depends on the black hole spin, and (ii) the mass accretion rate goes as , being the energy gap at the cusp of the torus (i.e. , being the potential at the boundary of the torus) and the adiabatic index.
Yokosawa [237, 238], also using Wilson's formulation, studied the structure and dynamics of relativistic accretion disks and the transport of energy and angular momentum in magneto-hydrodynamical accretion onto a rotating black hole. In his code the hydrodynamic equations are solved using the Flux-Corrected-Transport (FCT) scheme  (a second-order flux-limiter method in smooth regions which avoids oscillations near discontinuities by reducing the magnitude of the numerical flux), and the magnetic induction equation is solved using the constrained transport method . The code contains a parametrized viscosity based on the -model . The simulations revealed different flow patterns, inside the marginally stable orbit, depending on the thickness h of the accretion disk. For thick disks with , being the radius of the event horizon, the flow pattern becomes turbulent.
Since those analytic studies numerical simulations by an increasing number of authors (see, e.g., [185, 25] and references therein) have extended the simplified analytic models and have helped to develop a thorough understanding of the hydrodynamic accretion scenario in its fully three-dimensional character. These investigations have revealed the formation of accretion disks and the appearance of non-trivial phenomena such as shock waves or flip-flop (tangential) instabilities.
Most of the existing numerical work has used Newtonian hydrodynamics to study the accretion onto non-relativistic stars. For compact accretors such as neutron stars or black holes the correct numerical modeling requires a general relativistic hydrodynamical description. Within the relativistic frozen star framework, wind accretion onto ``moving'' black holes was first studied in axisymmetry by Petrich et al. . In this work Wilson's formulation of the hydrodynamic equations was adopted. The integration algorithm was taken from  with the transport terms finite-differenced following the prescription given in . An artificial viscosity term of the form , with a being a constant, was added to the pressure terms of the equations.
More recently, an extensive survey of the morphology and dynamics of relativistic wind accretion past a Schwarzschild black hole was performed in [71, 70]. This investigation differs from  in both the use of a conservative formulation for the hydrodynamic equations (the Valencia formulation; see Section 2.1.3) and the use of advanced HRSC schemes. Axisymmetric computations were compared to , finding major differences in the shock location, opening angle, and accretion rates of mass and momentum. The reasons for the discrepancies may be diverse and related to the use of different formulations, numerical schemes and, possibly, to the grid resolution. In  canonical grid sizes were extremely coarse, of zones in r and respectively. The simulations presented in [71, 70] used much finer grids in every direction.
Non-axisymmetric two-dimensional studies, restricted to the equatorial plane of the black hole, were first performed in , motivated by the non-stationary patterns found in Newtonian simulations (see, e.g., ). The relativistic computations revealed that initially asymptotic uniform flows always accrete onto the hole in a stationary way which closely resembles the previous axisymmetric patterns.
Papadopoulos and Font  have recently presented a procedure which considerably simplifies the numerical integration of the general relativistic hydrodynamic equations near black holes. Their procedure is based on identifying classes of coordinates in which the black hole metric is free of coordinate singularities at the horizon (unlike the commonly adopted Boyer-Lindquist coordinates), independent of time, and admits a spacelike decomposition. With those coordinates the innermost radial boundary can be placed inside the horizon, allowing for a clean treatment of the entire (exterior) physical domain. In  Michel's (spherical) solution was re-derived using a particular coordinate system adapted to the black hole horizon, the Eddington-Finkelstein system. In Fig. 8 a representative sample of hydrodynamical quantities is plotted for this stationary solution. The solid lines correspond to the exact solution and the symbols correspond to the numerical solution. The solution is regular well inside the horizon at r =2 M . The steepness of the hydrodynamic quantities dominates the solution only near the real singularity.
In [73, 74] this approach was applied to the study of relativistic wind accretion onto rapidly rotating (Kerr) black holes. The effects of the black hole spin on the flow morphology were found to be confined to the inner regions of the black hole potential well. Within this region, the black hole angular momentum drags the flow, wrapping the shock structure around. An illustrative example is depicted in Fig. 9 . The left panel of this figure corresponds to a simulation employing the Kerr-Schild form of the Kerr metric, regular at the horizon. The right panel shows how the accretion pattern would look like, were the computation performed using the more common Boyer-Lindquist coordinates. The transformation induces a noticeable wrapping of the shock around the central hole. The shock would wrap infinitely many times before reaching the horizon. As a result, the computation in these coordinates would be much more challenging than in Kerr-Schild coordinates.
|Numerical Hydrodynamics in General Relativity
José A. Font
© Max-Planck-Gesellschaft. ISSN 1433-8351
Problems/Comments to firstname.lastname@example.org