General Relativistic Magnetohydrodynamic Code Research Articles

Hybrid GRMHD and force-free simulations of black hole accretion

ABSTRACT We present a new approach for stably evolving general relativistic magnetohydrodynamic (GRMHD) simulations in regions where the magnetization $\sigma =b^2/\rho c^2$ becomes large. GRMHD codes typically struggle to evolve plasma above $\sigma \approx 100$ in simulations of black hole accretion. To ensure stability, GRMHD codes will inject mass density artificially to the simulation as necessary to keep the magnetization below a ceiling value $\sigma _{\rm max}$. We propose an alternative approach where the simulation transitions to solving the equations of general relativistic force-free electrodynamics (GRFFE) above a magnetization $\sigma _{\rm trans}$. We augment the GRFFE equations in the highly magnetized region with approximate equations to evolve the decoupled field-parallel velocity and plasma energy density. Our hybrid scheme is explicit and easily added to the framework of standard-volume GRMHD codes. We present a variety of tests of our method, implemented in the GRMHD code koral, and we show results from a 3D hybrid GRMHD + GRFFE simulation of a magnetically arrested disc (MAD) around a spinning black hole. Our hybrid MAD simulation closely matches the average properties of a standard GRMHD MAD simulation with the same initial conditions in low magnetization regions, but it achieves a magnetization $\sigma \approx 10^6$ in the evacuated jet funnel. We present simulated horizon-scale images of both simulations at 230 GHz with the black hole mass and accretion rate matched to M87*. Images from the hybrid simulation are less affected by the choice of magnetization cut-off $\sigma _{\rm cut}$ imposed in radiative transfer than images from the standard GRMHD simulation.

Viscous torque in turbulent magnetized active galactic nucleus accretion disks and its effects on the gravitational waves of extreme mass ratio inspirals

The merger of supermassive black holes produces millihertz gravitational waves (GWs), which are potentially detectable by the future Laser Interferometer Space Antenna (LISA). Such binary systems are usually embedded in an accretion disk environment at the center of an active galactic nucleus (AGN). Recent studies suggest the plasma environment imposes measurable imprints on the GW signal if the mass ratio of the binary is around q ∼ 10−4 − 10−3. The effect of the gaseous environment on the GW signal is strongly dependent on the disk’s parameters; therefore, it is believed that future low-frequency GW detections will provide us with precious information about the physics of AGN accretion disks. We investigated this effect by measuring the viscous torque via modeling of the evolution of magnetized tori around the primary massive black hole. Using the general relativistic magnetohydrodynamic HARM-COOL code, we performed 2D and 3D simulations of weakly magnetized, thin accretion disks, with a possible truncation and transition to advection-dominated accretion flow. We studied the angular momentum transport and turbulence generated by the magnetorotational instability. We quantified the disk’s effective alpha viscosity and its evolution over time. We applied our numerical results to quantify the relativistic viscous torque on a hypothetical low-mass secondary black hole via a 1D analytical approach, and we estimated the GW phase shift due to the gas environment.

Bridging Scales in Black Hole Accretion and Feedback: Magnetized Bondi Accretion in 3D GRMHD

Fueling and feedback couple supermassive black holes (SMBHs) to their host galaxies across many orders of magnitude in spatial and temporal scales, making this problem notoriously challenging to simulate. We use a multi-zone computational method based on the general relativistic magnetohydrodynamic (GRMHD) code KHARMA that allows us to span 7 orders of magnitude in spatial scale, to simulate accretion onto a non-spinning SMBH from an external medium with a Bondi radius of R B ≈ 2 × 105 GM •/c 2, where M • is the SMBH mass. For the classic idealized Bondi problem, spherical gas accretion without magnetic fields, our simulation results agree very well with the general relativistic analytic solution. Meanwhile, when the accreting gas is magnetized, the SMBH magnetosphere becomes saturated with a strong magnetic field. The density profile varies as ∼r −1 rather than r −3/2 and the accretion rate is consequently suppressed by over 2 orders of magnitude below the Bondi rate . We find continuous energy feedback from the accretion flow to the external medium at a level of . Energy transport across these widely disparate scales occurs via turbulent convection triggered by magnetic field reconnection near the SMBH. Thus, strong magnetic fields that accumulate on horizon scales transform the flow dynamics far from the SMBH and naturally explain observed extremely low accretion rates compared to the Bondi rate, as well as at least part of the energy feedback.

Study accretion and ejection using a new GPU-accelerated GRMHD code

AbstractWe study disks and jets in various accretion states (SANE and MAD) using novel, GPU-accelerated general-relativistic magneto-hydrodynamic (GR-MHD) code which we developed, based on HARM. This code, written in CUDA-c and uses OpenMP to parallelize multi-GPU setups, allows high resolution simulations of accretion disks and the formation and structure of jets without the need of multi-node supercomputer infrastructure. A 2563 simulation is well within the reach of an Nvidia DGX-V100 server, with the computation being a factor about 100 times faster if only the CPU was used.We use this code to examine several disk structures, wind and jet properties in the MAD and SANE states. In the MAD state, we find that the magnetic flux threading the horizon mostly depends on the spin of the BH. This implies that the jet structure and power are strong functions of the spin, with non-spinning BHs have the widest jets.

GRaM-X: a new GPU-accelerated dynamical spacetime GRMHD code for Exascale computing with the Einstein Toolkit

We present GRaM-X (General Relativistic accelerated Magnetohydrodynamics on AMReX), a new GPU-accelerated dynamical-spacetime general relativistic magnetohydrodynamics (GRMHD) code which extends the GRMHD capability of Einstein Toolkit to GPU-based exascale systems. GRaM-X supports 3D adaptive mesh refinement (AMR) on GPUs via a new AMR driver for the Einstein Toolkit called CarpetX which in turn leverages AMReX, an AMR library developed for use by the United States DOE’s Exascale Computing Project. We use the Z4c formalism to evolve the Einstein equations and the Valencia formulation to evolve the equations of GRMHD. GRaM-X supports both analytic as well as tabulated equations of state. We implement TVD and WENO reconstruction methods as well as the HLLE Riemann solver. We test the accuracy of the code using a range of tests on static spacetime, e.g. 1D magnetohydrodynamics shocktubes, the 2D magnetic rotor and a cylindrical explosion, as well as on dynamical spacetimes, i.e. the oscillations of a 3D Tolman-Oppenheimer-Volkhof star. We find excellent agreement with analytic results and results of other codes reported in literature. We also perform scaling tests and find that GRaM-X shows a weak scaling efficiency of ∼40%–50% on 2304 nodes (13824 NVIDIA V100 GPUs) with respect to single-node performance on OLCF’s supercomputer Summit.

Nozzle Shocks, Disk Tearing, and Streamers Drive Rapid Accretion in 3D GRMHD Simulations of Warped Thin Disks

The angular momentum of gas feeding a black hole (BH) may be misaligned with respect to the BH spin, resulting in a tilted accretion disk. Rotation of the BH drags the surrounding spacetime, manifesting as Lense–Thirring torques that lead to disk precession and warping. We study these processes by simulating a thin (H/r = 0.02), highly tilted () accretion disk around a rapidly rotating (a = 0.9375) BH at extremely high resolutions, which we performed using the general-relativistic magnetohydrodynamic code H-AMR. The disk becomes significantly warped and continuously tears into two individually precessing subdisks. We find that mass accretion rates far exceed the standard α-viscosity expectations. We identify two novel dissipation mechanisms specific to warped disks that are the main drivers of accretion, distinct from the local turbulent stresses that are usually thought to drive accretion. In particular, we identify extreme scale height oscillations that occur twice an orbit throughout our disk. When the scale height compresses, “nozzle” shocks form, dissipating orbital energy and driving accretion. Separate from this phenomenon, there is also extreme dissipation at the location of the tear. This leads to the formation of low-angular momentum “streamers” that rain down onto the inner subdisk, shocking it. The addition of low-angular momentum gas to the inner subdisk causes it to rapidly accrete, even when it is transiently aligned with the BH spin and thus unwarped. These mechanisms, if general, significantly modify the standard accretion paradigm. Additionally, they may drive structural changes on much shorter timescales than expected in α-disks, potentially explaining some of the extreme variability observed in active galactic nuclei.

Quantifying uncertainties in general relativistic magnetohydrodynamic codes

In this paper, we show that similar open-source codes for general relativistic magnetohydrodynamic (GRMHD) produce different results for key features of binary neutron star mergers. First, we present a new open-source version of the publicly available IllinoisGRMHD code that provides support for realistic, finite temperature equations of state. After stringent tests of our upgraded code, we perform a code comparison between GRHydro, IllinoisGRMHD, Spritz, and WhiskyTHC, which implement the same physics, but slightly different computational methods. The benefit of the comparison is that all codes are embedded in the EinsteinToolkit suite, hence their only difference is algorithmic. We find similar convergence properties, fluid dynamics, and gravitational waves, but different merger times, remnant lifetimes, and gravitational wave phases. Such differences must be resolved before the post-merger dynamics modeled with such simulations can be reliably used to infer the properties of nuclear matter especially in the era of precision gravitational wave astronomy.

Studying Postmerger Outflows from Magnetized-neutrino-cooled Accretion Disks

Neutrino-cooled accretion flow around a spinning black hole, produced by a compact binary merger, is a promising scenario for jet formation and launching magnetically driven outflows. Based on GW170817 gravitational wave detection by LIGO and Virgo observatories, followed by electromagnetic counterparts, this model can explain the central engine of the short-duration gamma-ray bursts and kilonova radiations. Using the open-source general relativistic magnetohydrodynamic HARM-COOL code, we evolved several 2D magnetized accretion disk–black hole models with a realistic equation of state in the fixed curved spacetime background. We applied the particle tracer technique to measure the properties of the outflows. The disk and black hole's initial parameters are chosen in a way to represent different possible postmerger scenarios of the merging compact objects. Our simulations show a strong correlation between the black hole's spin and ejected mass. Generally, mergers producing massive disks and rapidly spinning black holes launch stronger outflows. We observed our models generate winds with moderate velocity (v/c ∼ 0.1–0.2) and a broad range of electron fractions. We use these results to estimate the luminosity and light curves of possible radioactively powered transients emitted by such systems. We found that the luminosity peaks within the range of 1040–1042 erg s−1, which agrees with previous studies for disk wind outflows.

CuHARM: A New GPU-accelerated GRMHD Code and Its Application to ADAF Disks

We introduce a new GPU-accelerated general relativistic magnetohydrodynamic code based on HARM, which we call cuHARM. The code is written in CUDA-C and uses OpenMP to parallelize multi-GPU setups. Our code allows us to run high-resolution simulations of accretion disks and the formation and structure of jets without a need for multinode supercomputer infrastructure. A 2563 simulation is well within the reach of an Nvidia DGX-V100 server, with the computation being about 10 times faster than if only the CPU is used. We use this code to examine several disk structures all in the SANE state. We find that (i) increasing the magnetic field in the SANE state does not affect the mass accretion rate; (ii) simultaneously increasing the disk size and the magnetic field, while keeping the ratio of energies fixed, leads to the destruction of the jet once the magnetic flux through the horizon decreases below a certain limit (this demonstrates that the existence of a jet is a function of the initial mass distribution, and not of the initial intensity of the magnetic field, since the magnetorotational instability dictates the evolution of the magnetic field); and (iii) the structure of the jet is a weak function of the adiabatic index of the gas, with relativistic gas tending to have a wider jet.

H-AMR: A New GPU-accelerated GRMHD Code for Exascale Computing with 3D Adaptive Mesh Refinement and Local Adaptive Time Stepping

General relativistic magnetohydrodynamic (GRMHD) simulations have revolutionized our understanding of black hole accretion. Here, we present a GPU-accelerated GRMHD code H-AMR with multifaceted optimizations that, collectively, accelerate computation by 2–5 orders of magnitude for a wide range of applications. First, it introduces a spherical grid with 3D adaptive mesh refinement that operates in each of the three dimensions independently. This allows us to circumvent the Courant condition near the polar singularity, which otherwise cripples high-resolution computational performance. Second, we demonstrate that local adaptive time stepping on a logarithmic spherical-polar grid accelerates computation by a factor of ≲10 compared to traditional hierarchical time-stepping approaches. Jointly, these unique features lead to an effective speed of ∼109 zone cycles per second per node on 5400 NVIDIA V100 GPUs (i.e., 900 nodes of the OLCF Summit supercomputer). We illustrate H-AMR's computational performance by presenting the first GRMHD simulation of a tilted thin accretion disk threaded by a toroidal magnetic field around a rapidly spinning black hole. With an effective resolution of 13,440 × 4608 × 8092 cells and a total of ≲22 billion cells and ∼0.65 × 108 time steps, it is among the largest astrophysical simulations ever performed. We find that frame dragging by the black hole tears up the disk into two independently precessing subdisks. The innermost subdisk rotation axis intermittently aligns with the black hole spin, demonstrating for the first time that such long-sought alignment is possible in the absence of large-scale poloidal magnetic fields.

Handing off the outcome of binary neutron star mergers for accurate and long-term postmerger simulations

We perform binary neutron star (BNS) merger simulations in full dynamical general relativity with IllinoisGRMHD, on a Cartesian grid with adaptive-mesh refinement. After the remnant black hole has become nearly stationary, the evolution of the surrounding accretion disk on Cartesian grids over long timescales (1s) is suboptimal, as Cartesian coordinates over-resolve the angular coordinates at large distances, and the accreting plasma flows obliquely across coordinate lines dissipating angular momentum artificially from the disk. To address this, we present the Handoff, a set of computational tools that enables the transfer of general relativistic magnetohydrodynamic (GRMHD) and spacetime data from IllinoisGRMHD to HARM3D, a GRMHD code that specializes in modeling black hole accretion disks in static spacetimes over long timescales, making use of general coordinate systems with spherical topology. We demonstrate that the Handoff allows for a smooth and reliable transition of GRMHD fields and spacetime data, enabling us to efficiently and reliably evolve BNS dynamics well beyond merger. We also discuss future plans, which involve incorporating advanced equations of state and neutrino physics into BNS simulations using the \handoff approach.

Testing Relativistic Reflection Models with GRMHD Simulations of Accreting Black Holes

X-ray reflection spectroscopy is currently one of the leading techniques for studying the inner part of accretion disks around black holes, measuring black hole spins, and even testing fundamental physics in strong gravitational fields. However, the accuracy of these measurements depends on the reflection models employed for the spectral analysis, which are sometimes questioned. In this work, we use a general relativistic magnetohydrodynamic code to generate a thin accretion disk in Kerr spacetime and ray-tracing techniques to calculate its relativistically broadened reflection spectrum. We simulate NuSTAR observations and we test the capability of current reflection models to recover the correct input parameters. Our study shows that we can measure the correct input parameters in the case of high inclination angle sources, while we find some minor discrepancies when the inclination angle of the disk is low.

Implementing a new recovery scheme for primitive variables in the general relativistic magnetohydrodynamic code Spritz

General relativistic magnetohydrodynamic (GRMHD) simulations represent a fundamental tool to probe various underlying mechanisms at play during binary neutron star (BNS) and neutron star (NS) - black hole (BH) mergers. Contemporary flux-conservative GRMHD codes numerically evolve a set of conservative equations based on `conserved' variables which then need to be converted back into the fundamental (`primitive') variables. The corresponding conservative-to-primitive variable recovery procedure, based on root-finding algorithms, constitutes one of the core elements of such GRMHD codes. Recently, a new robust, accurate and efficient recovery scheme called RePrimAnd was introduced, which has demonstrated the ability to always converge to a unique solution. The scheme provides fine-grained error policies to handle invalid states caused by evolution errors, and also provides analytical bounds for the error of all primitive variables. In this work, we describe the technical aspects of implementing the RePrimAnd scheme into the GRMHD code Spritz. To check our implementation as well as to assess the various features of the scheme, we perform a number of GRMHD tests in three dimensions. Our tests, which include critical cases such as a NS collapse to a BH as well as the early evolution (~50 ms) of a Fishbone-Moncrief BH-accrection disk system, show that RePrimAnd is able to support magnetized, low density environments with magnetic-to-fluid pressure ratios as high as 10^4, in situations where the previously used recovery scheme fails.

Spritz: general relativistic magnetohydrodynamics with neutrinos

We here present a new version of the publicly available general relativistic magnetohydrodynamic (GRMHD) code Spritz, which now includes an approximate neutrino leakage scheme able to handle neutrino cooling and heating. The leakage scheme is based on the publicly available ZelmaniLeak code, with a few modifications in order to properly work with Spritz. We discuss the involved equations, physical assumptions, and implemented numerical methods, along with a large battery of general relativistic tests performed with and without magnetic fields. Our tests demonstrate the correct implementation of the neutrino leakage scheme, paving the way for further improvements of our neutrino treatment and the first application to magnetized binary neutron star mergers. We also discuss the implementation in the Spritz code of high-order methods for a more accurate evolution of hydrodynamical quantities.

A Numerical Scheme for General Relativistic Radiation Magnetohydrodynamics Based on Solving a Grid-based Boltzmann Equation

Abstract We develop a general relativistic radiation magnetohydrodynamics (GR-RMHD) code inazuma in which the time-dependent radiation transfer equation (frequency-integrated Boltzmann equation) is solved in curved spacetime. The Eddington tensor is derived from the specific intensity, and we solve the zeroth and first moment equations in order to update the radiation fields. Therefore, our code can solve the radiation field around relativistic compact objects more appropriately than an approximation method like the M1 closure scheme. The numerical scheme of magnetohydrodynamics is the same as that of our previous code. In some test calculations for propagating radiation and radiation hydrodynamics in flat spacetime, our code shows similar results to our previous work. Radiation propagation in curved spacetime is also properly solved for. We also show the radiation transport from the super-Eddington accretion disk around the black hole. The disk structure, such as the density, velocity, and temperature, is fixed by the model obtained using the GR-RMHD simulation with the M1 method. We found that the difference between our scheme and the M1 method appears in the optically thin outflow region around the rotation axis while the radiation field is almost the same in the optically thick disk region.

Observational signatures of disc and jet misalignment in images of accreting black holes

ABSTRACT Black hole (BH) accretion is one of nature’s most efficient energy extraction processes. When gas falls in, a significant fraction of its gravitational binding energy is either converted into radiation or flows outwards in the form of BH-driven jets and disc-driven winds. Recently, the Event Horizon Telescope (EHT), an Earth-sized submillimetre radio interferometer, captured the first images of M87’s BH. These images were analysed and interpreted using general relativistic magnetohydrodynamics (GRMHD) models of accretion discs with rotation axes aligned with the BH spin axis. However, since infalling gas is often insensitive to the BH spin direction, misalignment between accretion disc and BH spin may be a common occurrence in nature. In this work, we use the general relativistic radiative transfer code bhoss to calculate the first synthetic radio images of (highly) tilted disc/jet models generated by our GPU-accelerated GRMHD code h-amr . While the tilt does not have a noticeable effect on the system dynamics beyond a few tens of gravitational radii from the BH, the warping of the disc and jet can imprint observable signatures in EHT images on smaller scales. Comparing the images from our GRMHD models to the 43 and 230 GHz EHT images of M87, we find that M87 may feature a tilted disc/jet system. Further, tilted discs and jets display significant time variability in the 230 GHz flux that can be further tested by longer-duration EHT observations of M87.

Spritz: a new fully general-relativistic magnetohydrodynamic code

The new era of multimessenger astrophysics requires the capability of studying different aspects of the evolution of compact objects. In particular, the merger of neutron star binaries is a strong source of gravitational waves and electromagnetic radiation, from radio to γ-rays, as demonstrated by the detection of GW170817 and its electromagnetic counterparts. In order to understand the physical mechanisms involved in such systems, it is necessary to employ fully general relativistic magnetohydrodynamic (GRMHD) simulations able to include the effects of a composition and temperature dependent equation of state describing neutron star matter as well as neutrino emission and reabsorption. Here, we present our new code named Spritz that solves the GRMHD equations in 3D Cartesian coordinates and on a dynamical spacetime. The code can support tabulated equations of state, taking into account finite temperature effects and allowing for the inclusion of neutrino radiation. In this first paper, we present the general features of the code and a series of tests performed in special and general relativity to assess the robustness of the basic GRMHD algorithms implemented. Among these tests, we also present the first comparison between a non-staggered and a staggered formulation of the vector potential evolution, which is used to guarantee the divergence-less character of the magnetic field. With respect to other publicly available GRMHD codes, Spritz combines the robust approach of a staggered formulation of the vector potential together with the use of an equation of state driver (EOS_Omni) that can allow the code to use finite temperature equations of state. A next version of the code will fully test the EOS_Omni driver by coupling it with a neutrino leakage scheme.

Two-moment scheme for general-relativistic radiation hydrodynamics: a systematic description and new applications

ABSTRACT We provide a systematic description of the steps necessary – and of the potential pitfalls to be encountered – when implementing a two-moment scheme within an implicit–explicit (IMEX) scheme to include radiative-transfer contributions in numerical simulations of general-relativistic (magneto-)hydrodynamics (GRMHD). We make use of the M1 closure, which provides an exact solution for the optically thin and thick limits, and an interpolation between these limits. Special attention is paid to the efficient solution of the emerging set of implicit conservation equations. In particular, we present an efficient method for solving these equations via the inversion of a 4 × 4-matrix within an IMEX scheme. While this method relies on a few approximations, it offers a very good compromise between accuracy and computational efficiency. After a large number of tests in special relativity, we couple our new radiation code, frac, with the GRMHD code bhac to investigate the radiative Michel solution, namely, the problem of spherical accretion on to a black hole in the presence of a radiative field. By performing the most extensive exploration of the parameter space for this problem, we find that the accretion’s efficiency can be expressed in terms of physical quantities such as temperature, T, luminosity, L, and black hole mass, M, via the expression $\varepsilon =(L/L_{\rm Edd})/(\dot{M}/\dot{M}_{\rm Edd})= 7.41\times 10^{-7}\left(T/10^6\, \mathrm{K}\right)^{0.22} \left(L/L_\odot \right)^{0.48} \left(M/M_\odot \right)^{0.48}$, where LEdd and $\dot{M}_{\mathrm{Edd}}$ are the Eddington luminosity and accretion rate, respectively. Finally, we also consider the accretion problem away from spherical symmetry, finding that the solution is stable under perturbations in the radiation field.

Numerical relativity in spherical coordinates: A new dynamical spacetime and general relativistic MHD evolution framework for the Einstein Toolkit

We present SphericalNR, a new framework for the publicly available Einstein Toolkit that numerically solves the Einstein field equations coupled to the equations of general relativistic magnetohydrodynamics (GRMHD) in a 3+1 split of spacetime in spherical coordinates without symmetry assumptions. The spacetime evolution is performed using reference-metric versions of either the Baumgarte-Shapiro-Shibata-Nakamura equations or the fully covariant and conformal Z4 system with constraint damping. We have developed a reference-metric version of the Valencia formulation of GRMHD with a vector potential method, guaranteeing the absence of magnetic monopoles during the evolution. In our framework, every dynamical field (both spacetime and matter) is evolved using its components in an orthonormal basis with respect to the spherical reference-metric. Furthermore, all geometric information about the spherical coordinate system is encoded in source terms appearing in the evolution equations. This allows for the straightforward extension of Cartesian high-resolution shock-capturing finite volume codes to use spherical coordinates with our framework. To this end, we have adapted GRHydro, a Cartesian finite volume GRMHD code already available in the Einstein Toolkit, to use spherical coordinates. We present the full evolution equations of the framework, as well as details of its implementation in the Einstein Toolkit. We validate SphericalNR by demonstrating it passes a variety of challenging code tests in static and dynamical spacetimes.

Large-scale poloidal magnetic field dynamo leads to powerful jets in GRMHD simulations of black hole accretion with toroidal field

ABSTRACT Accreting black holes (BHs) launch relativistic collimated jets, across many decades in luminosity and mass, suggesting the jet launching mechanism is universal, robust, and scale-free. Theoretical models and general relativistic magnetohydrodynamic (GRMHD) simulations indicate that the key jet-making ingredient is large-scale poloidal magnetic flux. However, its origin is uncertain, and it is unknown if it can be generated in situ or dragged inward from the ambient medium. Here, we use the GPU-accelerated GRMHD code h-amr to study global 3D BH accretion at unusually high resolutions more typical of local shearing box simulations. We demonstrate that turbulence in a radially extended accretion disc can generate large-scale poloidal magnetic flux in situ, even when starting from a purely toroidal magnetic field. The flux accumulates around the BH till it becomes dynamically important, leads to a magnetically arrested disc (MAD), and launches relativistic jets that are more powerful than the accretion flow. The jet power exceeds that of previous GRMHD toroidal field simulations by a factor of 10 000. The jets do not show significant kink or pinch instabilities, accelerate to γ ∼ 10 over three decades in distance, and follow a collimation profile similar to the observed M87 jet.