Toward General-Relativistic Magnetohydrodynamics Simulations in Stationary Non-Vacuum Spacetimes

Accretion of magnetized gas on compact astrophysical objects such as black holes has been successfully modeled using general relativistic magnetohydrodynamic (GRMHD) simulations. These simulations have largely been performed in the Kerr metric, which describes the spacetime of a vacuum and stationary spinning black hole (BH) in general relativity (GR). The simulations have revealed important clues on the physics of accretion and jets near the BH event horizon, and have been used to interpret recent Event Horizon Telescope images of the supermassive BHs, M87$^*$ and Sgr A$^*$. GRMHD simulations require the spacetime metric in horizon-penetrating coordinates such that all metric coefficients are regular at the event horizon. The Kerr metric and its electrically charged spinning analog, the Kerr-Newman metric, are currently the only metrics available in such coordinates. We report here horizon-penetrating forms of a large class of stationary, axisymmetric, spinning metrics. These can be used to carry out GRMHD simulations of accretion on spinning, nonvacuum BHs and non-BHs within GR, as well as accretion on spinning objects described by non-GR metric theories of gravity.

The Event Horizon Telescope (EHT) collaboration recently demonstrated that it is now possible to "measure" the spacetime metric of astrophysical BHs.More concretely, we are able to measure the level of agreement of the spacetime geometry of the supermassive BHs, M87 * (The EHT Collaboration et al. 2019a,b;Psaltis et al. 2020) and Sgr A * (The EHT Collaboration et al. 2022a,b), with that of a Kerr BH (Kerr 1963), which is a vacuum, stationary solution of general relativity (GR; see, e.g., Wiltshire et al. 2009).Such investigations allow us, in principle, to test fundamental aspects of GR such as the no-hair hypothesis (see, e.g., Carter 1971), and whether astrophysical BHs truly have no surfaces and are devoid of matter (The EHT Collaboration et al. 2022b).
It has also been possible to test the level of agreement of EHT data with non-BH spacetime geometries (such as boson stars, wormholes, classical naked singularities, etc.) within GR, as well as both BH and non-BH spacetimes in non-GR theories (Völkel et al. 2021;Kocherlakota et al. 2021;The EHT Collaboration et al. 2022b;Vagnozzi et al. 2023).This program has seen intense activity, and current approaches use a combination of GRMHD simulations in Kerr BH spacetimes (see, e.g., The EHT Collaboration et al. 2019b, 2022b) as well as non-GR MHD ("xGRMHD") simulations in non-rotating spacetimes (Mizuno et al. 2018;Olivares et al. 2020;Fromm et al. 2021;Röder et al. 2023).In addition, a diverse range of instructive semi-analytic accretion flow models have been considered (see, e.g., Broderick et al. 2014;Shaikh et al. 2019;Shaikh & Joshi 2019;Narayan et al. 2019;Paul et al. 2020;Bauer et al. 2022;Özel et al. 2022;Younsi et al. 2023;Kocherlakota & Rezzolla 2022;Ayzenberg 2022;The EHT Collaboration et al. 2022b) to model the physics of accretion flows and their emission.
Nearly all the xGRMHD work to date has been restricted to spherically symmetric "static" spacetime metrics.This is because, apart from the spinning Kerr metric (Kerr 1963) and Kerr-Newman metric (Newman et al. 1965), the only metrics available in the horizon-penetrating "Kerr-Schild" coordinates required by current GRMHD codes (see, e.g., Font et al. 1998;McKinney & Gammie 2004;Sadowski et al. 2014;Porth et al. 2017;Liska et al. 2022 for a discussion of representative GRMHD codes) are static metrics.The one other exception is Nampalliwar et al. (2022), who use the Kerr-Schild form of the Johannsen (2013a) stationary metric, which is derived from the well-known Johannsen & Psaltis (2011) parametrized BH metric.
Astrophysical objects are expected to possess angular momentum.Hence their spacetimes will in general be axisymmetric (not spherically symmetric), and their metrics should be "stationary" (not static).For a more realistic confrontation of the underlying effective theory of gravity and fields in the vicinity of ultracompact objects to current observations, and in anticipation of future higher-quality data, it is necessary to be able to model accretion and emission processes in generic stationary spacetimes.However, other than the aforementioned Kerr, Kerr-Newman and Johannsen metrics, no other stationary spacetime metrics have been written down in the horizon-penetrating coordinates needed for computer simulations.The goal of the present paper is to eliminate this roadblock.
We introduce a general family of stationary and axisymmetric metrics (Azreg-Aïnou 2014a,b) which are derived from a corresponding family of "seed" spherically symmetric metrics.We then cast these stationary metrics in horizon-penetrating "Kerr-Schild" form.This is the main contribution of the present paper.Using the results described in this paper, we will report elsewhere the first high-resolution, 3D GRMHD simulations in spinning non-Kerr spacetimes (Chatterjee et al. in prep.a,i).
We begin in section 2 by introducing the general form of a (spherically-symmetric) static metric.This form includes the Schwarzschild metric as a special case but covers a wide range of other models as well.We pick out and discuss a few well-known examples of the latter.
In section 3 we introduce the Azreg-Aïnou (AA; Azreg-Aïnou 2014b) metric, which is a stationary, axisymmetric generalization of the static metric described in section 2. The AA metric was proposed as an ansatz, inspired by the Newman & Janis (1965) algorithm, to describe generic spinning spacetimes.This metric has several attractive features, which we touch on in Sec.1.1 and consider in greater detail in Sec. 3.
In section 4 we report the form of the AA metric in horizon-penetrating coordinates, as needed for most GRMHD codes, and in Sec. 5 we write down the 3 + 1 decomposition of the metric required for certain other codes.The spherical-polar Kerr-Schild coordinates we use to write the horizon-penetrating AA metric are adapted to the ingoing principal null congruence of spacetime.
Finally, in section 6 we summarize our findings and present a representative list of stationary and axisymmetric metrics that are obtained from the AA metric for specific choices of its metric functions.These metrics may be directly used in future GRMHD simulations.We also describe the steps needed to generate other stationary metrics in the future.Throughout we use geometrized units, G = c = 1, and the metric signature is (−, +, +, +).
Of necessity, this paper includes a fair amount of technical discussion about the AA metric and its various forms.A reader who is interested in utilizing various BH and non-BH metrics to explore with GRMHD simulations, can directly obtain the main results by reading only sections 2 and 6, along with tables 1 and 2.

Background on the Azreg-Aïnou Metric
While in general the description of an arbitrary stationary and axisymmetric spacetime requires ten independent metric functions, the AA metric remarkably makes use of only three functions.Furthermore, when considered as a generalization of a static metric, two of these functions are fixed by the original static metric functions.The one remaining free function is a conformal factor that can be fixed by solving the appropriate field equations.Moreover, as we discuss, there is a natu-rally attractive choice for this third function such that it is no longer a free function.Thus, the AA metric is both a powerful and surprisingly simple ansatz for spinning spacetimes.
These simplifications are partly a consequence of the AA metric describing only circular spacetimes, a particularly important subfamily of stationary and axisymmetric spacetimes which require at most five free functions (see, e.g., Sec.2.2 of Gourgoulhon 2010).The further reduction to just three metric functions (or indeed only two functions in the most natural from of the AA metric) is achieved via a novel modification in Azreg-Aïnou (2014b) to the original Newman-Janis algorithm, where the regularity of a coordinate transformation is used to eliminate the ad hoc "decomplexification step" that was required in Newman & Janis (1965) (cf.also Rajan 2016).
Circular spacetimes have been of central interest in explorations of non-Kerr spacetimes, and several popular parametrized metrics have been constructed to characterize and study their varied properties.The Johannsen-Psaltis (JP; Johannsen & Psaltis 2011; Johannsen 2013a) framework has had success, e.g., in establishing tests of the no-hair conjecture (Johannsen & Psaltis 2010) and of the post-Newtonian structure of astrophysical BH spacetimes (Psaltis et al. 2020).The Konoplya-Rezzolla-Zhidenko (KRZ;Rezzolla & Zhidenko 2014;Konoplya et al. 2016) framework has been used, e.g., to map well-known static metrics onto a single lowdimensional space with high accuracy to enable comparisons between spacetimes (Kocherlakota & Rezzolla 2020) and to test BH metrics with gravitationalwave as well as X-ray observations (Völkel & Barausse 2020;Cardenas-Avendano et al. 2020).
We note that while the JP and KRZ parametric metrics allow detailed explorations of black hole spacetimes, the AA ansatz metric can characterize also non-BH objects such as naked singularities, wormholes, boson stars, etc.Furthermore, it is clear to see from Johannsen (2013a) that the JP metric uses a fundamentally different ansatz compared to the AA metric to describe circular spacetimes, requiring four metric functions (see eq. 51 there), and motivated by imposing that geodesics be Lioville-separable (eq. 10 there), i.e., that they possess a Carter constant.This is similarly true for the KRZ metric, where five free functions are permitted (see eq. 7 there).Konoplya & Zhidenko (2021) explore a subclass of the KRZ metric that admits a Carter constant.As we show later, null geodesics in the AA metric always possess a Carter constant, and in the simplest form of the AA metric (which we favor) time-like geodesics also have a Carter constant.We note that sev-eral other parametrization frameworks have been constructed and used successfully to investigate observable effects due to modifications of the spacetime metric in the strong-field regime (see, e.g., Vigeland et al. 2011;Carson & Yagi 2020).For the status of parameterizations of non-circular spacetimes, we direct the reader to Delaporte et al. (2022).
In addition to its simplicity and appealing properties discussed above, our choice to use the class of AA metrics is further motivated by the knowledge that nearly all well-known solution metrics across various theories of gravity can be cast in this form.This metric, therefore, presents an excellent starting point for a broad forward-modeling study of the effect of the spacetime metric on various observables of interest.The Liouville-separability of the null geodesic equations makes computing various characteristic features of the spacetime such as the location of the photon shell or the shape of the shadow boundary (see, e.g., Shaikh 2019, Kocherlakota et al. 2021, Solanki et al. 2022) analytically tractable.Similarly, the separability of the timelike orbits facilitates a study of, e.g., equatorial Keplerian orbits, and allows us to go one step further with semi-analytic techniques.Details regarding geodesic orbits can be found in appendix C.

GENERAL CLASS OF STATIC METRICS
The line element ds 2 = ĝµν dx µ dx ν of an arbitrary static and spherically-symmetric metric ĝµν can be written in spherical-polar coordinates, x µ = (t, r, ϑ, ϕ), in the form where dΩ 2 2 = dϑ 2 + sin 2 ϑ dϕ 2 is the standard lineelement on a unit 2-sphere.This form of the metric involves three functions of the coordinate r: f (r), g(r), R(r).However, an arbitrary static metric can be defined using only two metric functions, i.e., any one of the metric functions, f , g, R, can be eliminated by an appropriate change of coordinates (see, e.g., Ch. 14 of Plebanski & Krasinski 2012).For example, it is possible to find a different radial coordinate r in which −ĝ tt ĝrr = g(r) = 1 (as we do in Table 1, and discuss further in Appendix A).However, the required coordinate transformation does not always lead to analytically convenient metric functions.We therefore keep an arbitrary g(r)> 0 in (1).
The function R(r) gives the curvature-radius or arealradius as a function of the radial coordinate r.In simple models, R(r) = r (see Table 1).When this is the case, the coordinate-radius r of a 2-sphere is its curvatureradius, κ G (r) = R (2) (r)/2 = 1/r 2 , as well as its arealradius, A(r) = 4πr 2 .Here we have used κ G , R (2) and A to denote the Gaussian curvature, the Ricci scalar, and the area of a 2-sphere, respectively.While it is Table 1.Metric functions of selected spherically-symmetric and static spacetimes.We have used the coordinate freedom to set −ĝttĝrr = g(r) = 1 in the general metric (1), hence we have only two metric functions, f (r) and R(r).The Schwarzschild and Reissner-Nordström (RN) spacetimes describe the vacuum and the electrovacuum black hole (BH) solutions of general relativity (GR), respectively.The Modified Hayward BH spacetime contains a regular BH and is generated by an anisotropic fluid.The GMGHS BH metric is a solution to the low energy theory of the heterotic string, and contains a scalar field (a dilaton, which vanishes asymptotically) as well as an electromagnetic field.The JNW naked singularity spacetime is generated by a scalar field.The JMN-1 naked singularity spacetime is an interior solution to the Schwarzschild spacetime and is generated by an anisotropic fluid.For easier comparison with the other five models, the JMN-1 metric is presented here with g(r) = 1 (see Appendix A for details), thus differing from the original form given in eq. ( 29) of Joshi et al. (2011).The parameter Q denotes an electromagnetic charge, L denotes a (de Sitter) length scale, R b ≥ 2.5M denotes a matching radius (σ = M/(R b − M ) and r b = (1 − σ)R b ), and ν denotes a scalar charge (r⋆ = 2M/(1 − ν) locates the curvature singularity).

Object
Theory Spacetime always possible to find a coordinate transformation to curvature coordinates in which R(r) = r, this again may only be tractable numerically for some models.Finding simple analytic forms for all the metric functions in such coordinates is not in general possible (see, e.g., Sec.IV D of Kocherlakota & Rezzolla 2020).Thus, for general convenience, we will leave R(r) free.
When the metric (1) describes a BH spacetime, its event horizon is a null, stationary surface.The normal to such a constant-r surface satisfies ĝµν ∂ µ r • ∂ ν r = 0, from which we find that the horizon is located at the outermost root of ĝrr = 0, i.e., at the outermost root of f (r).
In Table 1 we list the metric functions corresponding to six well-known and representative sphericallysymmetric static solutions in GR and in string theory.For uniformity, we have written all the metrics in coordinates such that g(r) = 1.The first two solutions listed are the canonical nonrotating BH solutions of GR, namely, the Schwarzschild BH and the charged Reissner-Nordström (RN) BH, which contain, respectively, a spacelike and a timelike curvature singularity (see, e.g., Poisson 2004).Several regular BH spacetime models have been proposed to mimic the desired effect of singularity-resolving physics (see, e.g., Bardeen 1968;Ayón-Beato & García 1998;Bronnikov 2001;Dymnikova 2004;Hayward 2006).The Modified Hayward BH (Zhou & Modesto 2023a; n = 4) spacetime included in Table 1 is generated by an anisotropic fluid with a radial equation of state, ω r = p r /ρ = −1, where ρ is the total energy density and p r is the radial pressure in the fluid rest-frame.The fluid is darkenergy-like, producing the pressure necessary to "avoid" a singularity.The charged Gibbons-Maeda-Garfinkle-Horowitz-Strominger (GMGHS) BH (Gibbons & Maeda 1988;Garfinkle et al. 1991) is different from the RN BH of GR because of the presence of a dilaton (a scalar field) that mediates the interaction between electromagnetism and gravity in the low-energy effective action of string theory.In particular, the central singularity of the GMGHS BH remains spacelike.
In addition to the above four BH solutions, Table 1 lists two other non-BH models.Over the years, non-BH solutions to various theories have been discussed in the literature and significant progress has been made to determine ways in which such objects can be distinguished from BHs using observations (see, e.g., Kocherlakota et al. 2021;The EHT Collaboration et al. 2022b).The Janis-Newman-Winicour naked singularity (JNW; Janis et al. 1968) spacetime is a static solution generated by a minimally-coupled, massless scalar field in GR.Indeed, this is also a solution to string theory when an electromagnetic field is absent (Virbhadra 1997) and is in addition a solution to the Brans-Dicke theory with the parameter ω = −1 (Kar 1997).Finally, the Joshi-Malafarina-Narayan-1 naked singularity (JMN-1; Joshi et al. 2011) is constructed in GR using an anisotropic fluid with vanishing radial pressure, p r = 0.

GENERALIZATION TO STATIONARY METRICS
Inspired by the Newman & Janis (1965) algorithm, Azreg-Aïnou (2014b) proposed a stationary and axisymmetric generalization of the spherically-symmetric metric (1) described in section 2. We will refer to this stationary metric as the AA metric and denote it by g µν (to distinguish it from ĝµν of the spherically symmetric metric).In Boyer-Lindquist (BL;Boyer & Lindquist 1967) coordinates, x µ = (t, r, ϑ, ϕ), 1 the AA metric is given by ds Here the parameter a is the spin of the central object and the "stationary" metric functions {F, ∆, Σ, Π} can be related to the "static" metric functions {f, g, R} in (1) as discussed below.This then leaves only the metric function X = X(r, ϑ) to be fixed by the field equations (Azreg-Aïnou 2014b).
While we focus in this paper on the AA metric, we note that different generalizations of the static metric (1) are possible, in principle, and can be achieved by modifying the complex coordinate transformation involved in the Newman & Janis (1965)-type solution-generating technique (see, e.g., Azreg-Aïnou 2014a).It is important to note, however, that the particular transformations used in Azreg-Aïnou (2014a,b) send the Schwarzschild, the Reissner-Nordström, and the GMGHS metrics, to their appropriate spinning generalizations, namely to the Kerr, the Kerr-Newman, and the Kerr-Sen metrics, respectively.We direct the reader to Erbin (2017) for a review of such solution-generating techniques, and to Sec. 7.1 of Wald (1984) for a general discussion of stationary, axisymmetric spacetimes, and on the construction of BL-like coordinates.
Expressions for the metric functions {F, ∆, Σ, Π} in the AA metric (2) are most conveniently written by first defining two auxiliary functions that are specific combinations of the original static metric functions, f (r), g(r), R(r), which can be written even more transparently as 1 Note that we do not distinguish between the labels for the radial coordinates used to write the non-spinning (1) and the spinning (2) metrics to avoid a proliferation of symbols.Furthermore, while the metric signature here (−, +, +, +) differs from that in Azreg-Aïnou (2014b), it can be checked that eq. 2 corresponds to eq. 16 there.
A and B are each functions of r alone, and in terms of them the functions {F, ∆, Σ, Π} in the AA metric (2) take the remarkably simple form Note that F and ∆ depend only on r, while Σ and Π are functions of r and ϑ.The spinning AA metric g µν (2) reduces to the nonspinning metric ĝµν (1) in the limit of vanishing spin (a = 0) for appropriate choices of X (see Azreg-Aïnou 2014a for a discussion on conformal vs. normal fluids).
We briefly illustrate now how the remaining metric function X can be fixed by considering, as an example, the case when the background matter is an anisotropic fluid in GR which flows around the spin axis.Since in its own rest-frame, {e µ (a) (a = 0 − 3)}, the fluid stressenergy-momentum tensor is diagonal, the Einstein equations imply that the Einstein tensor in this frame, G (a)(b) = G µν e µ (a) e ν (b) , must also be diagonal.However, it can be shown that two of the off-diagonal elements of the Einstein tensor, G (r)(ϑ) and G (t)(ϕ) , will in general be non-zero.Demanding that these two terms vanish yields one nonlinear partial differential equation (PDE) and one linear PDE, each of which involves the free metric function X, the auxiliary function A(r), and the fluid angular velocity Ω.These PDEs are given in eqs.15 and 18 of Azreg-Aïnou (2014a), respectively, where the fluid angular velocity was chosen to be Ω = a/(A+a2 ) (see e µ t in eq.16 there).It is interesting to note that this angular velocity does not correspond to the ZAMO angular velocity introduced above, Ω Z = 2aF/Π.It matches, however, the angular velocity of the principal null congruences of the spacetime introduced below (see eq. 8).For non-fluid matter models, additional equations of motion for the matter fields must be solved.For example, if the matter is a scalar field, the associated Klein-Gordon equations must additionally be solved.Similarly, if the matter is an electromagnetic field, the Maxwell equations have to be accounted for as well.We will not enter into a discussion of this topic but refer the reader to Erbin (2017).
For the purposes of the present paper, it suffices to note that the subclass of AA metrics with X = Σ is particularly interesting, as it corresponds to spacetimes that are asymptotically flat (see Appendix B for further details).Furthermore, it can be checked that this choice for the conformal factor X always solves the G (r)(ϑ) = 0 equation for axially-spinning matter.Finally, all of the BH solutions discussed in this paper have X = Σ, whether they are vacuum or contain scalar and/or electromagnetic and/or axion fields.Even though these spacetimes arise as "solution metrics" to a variety of field equations, X = Σ = A(r) + a 2 cos 2 ϑ consistently remains a valid choice.For this subclass of the AA metrics, the possibly divergent behavior of the Ricci R and the Kretschmann K scalars can be seen from the behavior of their denominators, which go roughly as ∼ Σ −3 and as ∼ Σ −6 respectively.Thus, when this metric describes a singular spacetime, the curvature singularity is located at r = r ⋆ such that Σ(r ⋆ , π/2) = 0, which corresponds to a ring singularity located in the equatorial plane.This last equation is equivalent to R(r ⋆ ) = 0.
When the general AA metric (2) describes a BH spacetime, as above, a horizon is present at every location r = r H > r ⋆ where g rr = ∆/Σ = 0. Equivalently, horizons are located at the real, positive (R(r H ) > 0) roots of ∆(r).Conversely, if no such roots exist, then the AA metric does not correspond to a BH spacetime (e.g., the a > M Kerr metric).The event horizon in particular is located at the largest such root.Since ∆ is a function of r alone, it is reassuring to find that the event horizon is indeed a round sphere in BL coordinates.Furthermore, the horizon angular velocity, Ω H := Ω Z (r H ), is simply given as Ω H = 2aF (r H )/(A(r H )+a 2 ) 2 = a/(A(r H )+a 2 ).Clearly, Ω H is independent of ϑ, i.e., an AA BH rotates like a rigid body, as a consequence of the weak rigidity theorem (see, e.g., Sec.8.4.4 of Straumann 2013).These desirable properties, which are well-known for the Kerr and Kerr-Newman metrics, are now seen to be true for the very wide class of AA stationary BH metrics.
It is interesting to note that the general AA metric (2) can also describe spacetimes which contain regions that admit closed timelike curves.In such regions, g ϕϕ < 0 (see, e.g., Sec.V.B of Johannsen 2013b), i.e., Π < 0 (e.g., the a > M Kerr metric).
We show in Appendix C that the null geodesic equations are Lioville-separable.Thus, the AA metric always admits a Carter constant for null geodesics.Furthermore, we also identify a subclass of AA metrics that additionally admit separable timelike geodesic equations.These require the conformal factor to take the form X(r, ϑ) = X r (r) + X ϑ (ϑ).For this subclass, a Carter constant exists for all geodesics.Concomitantly, arbitrary geodesics of all AA metrics with X = Σ possess Carter constants.
The rr-component of the AA metric (2) in BL coordinates diverges at the horizon (∆(r) = 0).The Ricci and Kretschmann scalars, however, reveal this to be merely a coordinate singularity, an artifact of the choice of coordinates.The singularity can be eliminated by a change of coordinates, and this is the topic we turn to next.

HORIZON-PENETRATING COORDINATES
In this section, we convert the AA metric, which is written in BL coordinates in equation 2, to a horizonpenetrating form in which no coordinate singularity is present at the horizon.This is the form that is most useful for GRMHD and xGRMHD simulations of astrophysical accretion flows.
The ingoing (−) and outgoing (+) principal null congruences (PNCs) of the AA spacetime consist of null geodesics x µ (λ) whose tangents ℓ µ ± = ẋµ = dx µ /dλ satisfy ℓ ϑ ± = θ = 0 and θ = 0 (see, e.g., Misner et al. 1973, Hioki & Miyamoto 2008 for further details).The 1−form fields (l ± ) µ associated with the PNCs in an AA spacetime can be expressed elegantly as, 2 where E is some constant.We now define the spherical-polar ingoing Kerr-Schild (siKS) coordinates, x μ = (τ, r, ϑ, φ), as those in which the ingoing PNC takes the form, To distinguish between the BL (x µ ) and the siKS (x μ) coordinate systems, we will use a bar for the indices of the latter.We are interested in the ingoing PNC since our coordinate system should penetrate the future horizon H + when one exists, and cover patches I and II of the Kruskal (Kruskal 1960) or Penrose-Carter (Penrose 1963;Carter 1966) diagrams.The Jacobian Λ μµ = ∂ µ x μ for the coordinate transformation from BL coordinates to siKS coordinates (i.e., dx μ = Λ μµ dx µ ) can be inferred from equations 9 and 10 to be, The time coordinate τ then is clearly a "tortoise" time coordinate (see, e.g., Blau 2023), and we recognize immediately that, in the limit of vanishing spin (a → 0), the siKS coordinates are analogous to the coordinates presented in Eddington (1924) and Finkelstein (1958).
Using the above coordinate transformation, we can write the AA metric in siKS coordinates as, where the asterisks denote components that are fixed by symmetry, g μν = g ν μ.Notice that the metric components now diverge only at the spacetime/curvature singularity, present at Σ = 0.The inverse AA metric in siKS coordinates is For completeness, we cast the general AA metric into its "classic" Kerr-Schild form in eq.B14 of Appendix B.
In siKS coordinates, the event horizon is still located at the outermost root of g 11 = g rr = 0, i.e., ∆(r) = 0 ( 14), but as we can see from eq. 13 the metric components g μν no longer diverge at ∆ = 0. Thus, the coordinate singularity at the horizon has been removed by changing from BL to siKS coordinates, showing that the latter are horizon-penetrating coordinates.Furthermore, a comparison of the Kerr metric in particular, in these coordinates, with that reported in Font et al. (1998), with k = 1 there, and in McKinney & Gammie (2004), shows that these are the horizon-penetrating coordinates canonically used in the context of GRMHD simulations.
Note that the original set of Kerr-Schild coordinates used to describe the Kerr spacetime (Kerr 1963;Kerr & Schild 2009) are adapted to the outgoing PNC in which l μ + ∝ δ μ r (see, e.g., Wiltshire et al. 2009;Azreg-Aïnou 2014b).The corresponding time coordinate u is the retarded time as well as a null affine parameter that parametrizes the outgoing PNC on future null infinity I + .Equivalently, the ingoing version of the original set of Kerr-Schild coordinates would be adapted to the ingoing PNC in which l μ − ∝ δ μ r with the corresponding time coordinate v being the advanced time as well as a null affine parameter that parametrizes the ingoing PNC on past null infinity I − .Clearly, the "original" ingoing Kerr-Schild system is quite different from the spherical-ingoing Kerr-Schild coordinates which we use in this section and which are matched to current GRMHD codes.From eq. 10, that the tangent to the ingoing PNC in our coordinates is given as 5. HORIZON-PENETRATING 3 + 1 FORM We can cast the AA metric when written in horizonpenetrating siKS coordinates into its associated 3 + 1 form (Arnowitt et al. 2008) by introducing a foliation into spacelike hypersurfaces Σ τ that are iso-surfaces of the scalar time function τ (see, e.g., Poisson 2004, Gourgoulhon 2007, Alcubierre 2008, Rezzolla & Zanotti 2013).The unit timelike normal to these hypersurfaces is where α := −1/g τ τ is called the lapse function, and from eq. 14 is given by The four-vector that is dual to the normal, i.e., n μ = g μν n ν = −αg μτ = [1/α, −α(2F/X), 0, 0], corresponds to the four-velocity of the local timelike Eulerian observer.The four-velocity of this observer in BL coordinates can be obtained as n µ = Λ µ μn μ, where Λ µ μ = Λ μµ −1 .It is straightforward to then check that the angular velocity of the Eulerian observer in BL coordinates Ω E := n ϕ /n t matches the angular velocity of the zero angular momentum observer introduced above Ω Z = −g tϕ /g ϕϕ = 2aF/Π exactly.However, the Eulerian observer has nonzero BL radial velocity, i.e., n r = −2αF/X(= n r ).With eq.16, the induced (Riemannian) metric γ on the spacelike hypersurfaces can now be introduced as The projection tensor, which is used to obtain the spatial components of any four-vector or tensor, is simply γ μν = g μ ᾱγ ᾱν = δ μν + n μn ν since γ μν n ν = 0.
Although the 1-form n μdx μ is collinear with dτ (see eq. 16), the dual four-vector n μ∂ μ is not, in general, collinear with ∂ τ (since −αg μτ ✚ ∝ δ μ τ ).Indeed, this noncollinearity is typically captured via a four-vector β μ, which is defined such that n μβ μ = 0: Since by construction β τ = 0, it is convenient to work instead with the spacelike three-vector (see, e.g., McKinney & Gammie 2004, Porth et al. 2017) where the index i takes values i = 1, 2, 3.This threevector measures the shift of the isolines of the three spatial coordinates with respect to the normal to the hypersurfaces and is therefore called the shift vector (see, e.g., Sec.1.3 of Gourgoulhon 2010).
The desired 3 + 1 form of the AA metric can now be written, using the lapse function α (eq.17), the shift vector β i (eq.20), and the components γ ij (eq.18) of the induced metric on the hypersurfaces Σ τ , as For completeness, we demonstrate how the temporal u τ and the spatial u i components of a vector u μ = (u τ , u i ) can be obtained using the unit normal and the projection tensor respectively.First, with Γ = −n μu μ, we have u τ = Γ/α.Second, introducing the purely spatial vectors (v 0 = 0) v i := (γ i μu μ)/Γ, we obtain , where v 2 = v i v i .When u corresponds to the fourvelocity of a fluid element, Γ has the significance of being the fluid Lorentz factor with respect to the Eulerian observer.An excellent summary of the 3+1 treatment relevant for GRMHD can be found in Gammie et al. (2003).

SUMMARY AND CONCLUSIONS
We have presented in equation ( 13), section 4, a horizon-penetrating Kerr-Schild form of the rather general stationary and axisymmetric Azreg-Aïnou (AA) metric.This form of the metric can be used as an input when performing general relativistic (GR) as well as non-GR (xGR) magnetohydrodynamics simulations of accretion flows onto compact objects in arbitrary metric theories of gravity.The AA metric can be used to generate a number of popular metrics that describe the stationary spacetimes in GR corresponding to electro-vacuum black holes (Kerr 1963;Newman et al. 1965), regular black holes such as the Kerr-Hayward models (Hayward 2006;Bambi & Modesto 2013;Zhou & Modesto 2023a,b), naked singularities such as the Kerr-Janis-Newman-Winicour spacetime (Solanki et al. 2022), etc.It can also be used to describe black hole solutions arising in alternative theories of gravity (e.g., Sen 1992).
To illustrate the ease of generating the Kerr-Schild forms of new stationary metrics, we now summarize the steps involved.We start with an arbitrary static and spherically-symmetric "seed" metric, e.g., the ones listed in Table 1, in standard spherical-polar coordinates, x µ = (t, r, ϑ, ϕ), As we noted in section 2, it is always possible to find a radial coordinate r in which g(r) = 1.Since this reduces the number of free functions, we used such coordinates for the static models listed in Table 1.However, the required transformation may not always lead to simple analytical forms for f (r) and R(r).When it does not, it is preferable to use the more general form ( 22).
The AA stationary generalization of the metric ( 22) in spherical-polar ingoing Kerr-Schild coordinates, x μ = (τ, r, ϑ, φ), takes the form, where the "stationary metric functions" {F, ∆, Σ, Π} are fixed by the "static metric functions" {f, g, R} in the seed metric ( 22) as, The parameter a describes the spin of the spacetime.
The only remaining freedom is the stationary metric function X(r, ϑ).
If the background matter in the stationary spacetime is a fluid flowing around the spin axis (ϑ = 0) in general relativity (GR), the unknown function X = X(r, ϑ) is fixed by solving two Einstein equations, G (r)(ϑ) = 0 and G (t)(ϕ) = 0, in the fluid rest frame (the energymomentum-stress components T (r)(ϑ) and T (t)(ϕ) of the fluid vanish in this frame).We will not enter into detailed implications of these equations here but note that X = Σ is always a solution of the first equation, G (r)(ϑ) = 0.The choice X = Σ is additionally appealing since spacetimes with this property are asymptotically-flat (provided of course that the original seed static metric is also asymptotically flat).Verifying that X = Σ also solves the second Einstein equation, G (t)(ϕ) = 0, in the fluid frame, and that the resulting energy-momentum-stress tensor is physically valid, i.e., satisfies various energy conditions, is an involved calculation and is beyond the scope of this paper.However, if one merely wants a novel stationary metric in Table 2. Metric functions for the stationary and axisymmetric generalizations of the static and spherically-symmetric spacetimes given in Table 1.The spinning generalizations of a nonspinning "seed" metric (1) is given in Boyer-Lindquist coordinates in eq. 2 and in horizon-penetrating spherical ingoing Kerr-Schild coordinates in eq.23.These involve five "stationary metric functions" {F, ∆, Σ, Π, X}, the first four of which are related to f and R via eq.24.We show below the spinning generalizations of the specific static solutions listed in Table 1.Since we used coordinates for the static metrics in which g(r) = 1, g does not appear in the expressions below.For brevity we only show the four metric functions 2F , ∆, Σ = X below.The fifth metric function is given by Π = (R 2 + a 2 ) 2 − ∆a 2 sin 2 ϑ.All the spacetime models listed here are asymptotically flat, and a denotes the usual spin parameter.A variety of physical scenarios are captured by the exemplary set listed here: the vacuum Kerr BH, the electromagnetically charged Kerr-Newman BH, the (modified) Kerr-Hayward spinning BH, the spinning JNW naked singularity, and the spinning JMN-1 naked singularity are solutions of general relativity, whereas the Kerr-Sen spacetime describes electromagnetically-charged BHs of string theory which are also charged under a dilaton field as well as an axion field.

Spacetime
which to carry out GRMHD simulations, and if one is not too concerned about necessarily satisfying energy conditions, then the metric ( 23) with X = Σ, and a suitable choice of the seed functions f (r), g(r), R(r), would be a reasonable choice.We should note that, for non-fluid models, additional equations of motion for the fields (e.g., the Klein-Gordon equations or the Maxwell equations) must be solved.The status and successes of such metric-solutiongenerating techniques in being able to yield solutions also to the nongravitational field equations has been discussed in Erbin (2017).
Table 2 presents a compilation of the metric functions of the stationary AA spacetimes corresponding to the six static "seed" models listed in Table 1.These are all models for which the choice X = Σ is valid.The stationary Kerr (1963) metric arises naturally as the AA generalization of the static Schwarzschild metric.It describes vacuum, spinning black holes (BHs) in GR and has been used extensively in GRMHD simulations.Similarly, the AA transformation automatically generates the Kerr-Newman metric (Newman et al. 1965), which describes electromagnetically-charged, spinning BHs in GR, as the spinning generalization of the non-spinning Reissner-Nordstrom BH.The Kerr-Hayward model proposed here (the third model in Table 2) is expected to describe spinning regular (i.e., no singularities) BHs in GR (see, e.g., Bambi & Modesto 2013;Zhou & Modesto 2023b).The spinning equivalent of the GMGHS BH solution in string theory is given by the Kerr-Sen metric (Sen 1992).The form of the Kerr-Sen metric shown in Table 2 can be put into the form presented in eq.20 of Xavier et al. (2020) by replacing r with r + Q 2 /M .It is worth noting that the spacetime of Kerr-Sen BHs contains not just electromagnetic fields but also a scalar field (dilaton) as well as an axion field.All of the BH solutions in Table 2 reduce to the Kerr BH in appropriate limits (Q → 0, L → 0).
The spinning Janis-Newman-Winicour (JNW) spacetime (the fifth model in Table 2) has been proposed in Solanki et al. (2022) and is expected to describe a spinning naked singularity spacetime containing scalar field matter.Finally, the spinning JMN-1 metric proposed here is expected to describe a spinning naked singularity spacetime containing an anisotropic fluid in GR.Careful explorations of these spacetimes, especially with a focus on the physicality of the underlying matter (e.g., how well energy conditions are satisfied), is left for future work.Vagnozzi et al. (2023) give an immense compilation of spacetime metrics including those corresponding to static and spherically-symmetric spacetimes.These metric functions describe alternative nonspinning black holes which are of considerable interest and can be studied using simulations.More importantly, each static model could potentially be used as a seed metric to obtain the horizon-penetrating form of the stationary metric that it may belong to via equations 23 and 24.The only open issue, which needs to be checked in each case, is whether the choice X = Σ is physically valid.
We conclude on an optimistic note, envisioning the widespread adoption of GRMHD simulations, facilitated by our easy-to-use metric formulation.A Kerr simulation library has recently been used successfully to infer the properties of plasma present in the close vicinity of Sgr A * , as well as the spacetime geometry of this supermassive BH (The EHT Collaboration et al. 2022a,b).While it might be a challenging endeavor in the near future to construct similar libraries of fully 3D simulations in the other stationary spacetimes described in this paper (because of the many additional spacetime parameters in these models such as Q, L, ν, σ), it should at least be possible to build extensive 2D GRMHD simulation libraries for these spacetimes.This will enable valuable exploration of magnetized relativistic gas dynamics in diverse spacetimes.
Since we are focused here on asymptotically-flat spacetimes (see, e.g., Adamo et al. 2009), let us note first that for the seed metric (1) to be asymptotically-flat, its metric functions must become R(r) = r, f (r) = 1, and g(r) = 1 in the limit r → ∞.Then, if we define ημν as, ) we will find that as r → ∞, the Riemann tensor associated with η vanishes identically for such asymptoticallyflat seeds.Thus, asymptotically, η is the flat Minkowski metric η, in disguise, i.e., lim r→∞ η = η.We stress that the Riemann tensor does not vanish in general at finite coordinate radii r the metric tensor η.Now, the AA metric in the siKS coordinates (13) can be expressed in terms of η and the tangent to the ingoing PNC ℓ − (10) everywhere as, We note further that ℓ − is null with respect to both g and η.We refer to the form of the metric above as the generalized (spherical ingoing) Kerr-Schild form of the AA metric (see Sec. 32.5 of Stephani et al. 2009 for related discussion).This can be seen to be conformally related to the "classic" Kerr-Schild form, as given in eq.1.1 of Kerr & Schild (2009).Therefore, for the AA metric itself to be asymptotically-flat, we require (a) lim r→∞ (2F/Σ) = 0, and (b) lim r→∞ (X/Σ) = 1.The first of these conditions is met due to the properties of the asymptotically-flat seed metric functions whereas the second is trivially met when X = Σ.Thus, asymptotically-flat seed metrics admit asymptoticallyflat stationary generalizations when X = Σ.Note that this latter condition (X = Σ) is a sufficient but not a necessary condition.It is instructive to compare the asymptotically Minkowski metric in our siKS coordinates (B13) with the analogous form of the Minkowski metric in the original (outgoing) Kerr-Schild coordinates for the Kerr metric, as given, e.g., in equations 1.7 and 1.13 of Wiltshire et al. (2009).

C. SEPARABILITY OF THE GEODESIC EQUATION FOR THE STATIONARY METRIC
The Lagrangian L describing a geodesic orbit x µ (λ) is given as 2L = u µ u µ , where u µ = dx µ /dλ = ẋµ is the four-velocity along the geodesic and λ is an affine parameter along it.Working now in Boyer-Lindquist coordinates (2), due to the two Killing symmetries of the spacetime that are generated by T = ∂ t and Φ = ∂ ϕ , we can find momenta, p µ = ∂ ẋµ L = ∂ u µ L = u µ , that are conserved, corresponding to the two cyclic variables.These are identified as being the energy E = −u µ T µ and the azimuthal angular momentum L = u µ Φ µ of the orbit respectively, 3 which can be used to obtain ṫ where we have introduced the first impact parameter, ξ := L/E, and the determinant of the tϕ−sector of the AA metric tensor in BL coordinates (2), det [g tϕ ] := g tt g ϕϕ − g 2 tϕ = −(X 2 /Σ 2 )∆ sin 2 ϑ.Now with eq.C16 we can write, It is easy to see then that the geodesic equation for null geodesics (2L = 0) is fully separable, X 2 θ2 E 2 = η 2 − (a sin ϑ − ξ csc ϑ) 2 =: Θ 0 (ϑ) , (C18) 2 − ∆η 2 =: R 0 (r) , (C19) where we have introduced a separation constant η 2 which is related to the Carter constant C through η 2 = C/E 2 .The Carter constant in turn can be used to demonstrate the existence of a Killing-Yano tensor and an associated hidden symmetry of the motion.The fundamental principal null congruences (PNCs) of the AA spacetime consist of null geodesics that satisfy θ = 0 and θ = 0 (see, e.g., Misner et al. 1973, Hioki & Miyamoto 2008; see also Sec. 2.3 of Adamo et al. 2009).This is equivalent to requiring that Θ 0 = 0 and ∂ ϑ Θ 0 = 0, which yields a solution η = 0 and ξ = a sin 2 ϑ.On a related note, the Newman & Penrose (1962) complex null tetrad adapted to the outgoing ( ṙ > 0) PNC for the AA metric can be found in Azreg-Aïnou (2014a).
If X(r, ϑ) is of the form X(r, ϑ) = X r (r) + X ϑ (ϑ), then the geodesic equation is separable even for nonnull orbits, For the special case when X = Σ these are given simply as X r (r) = A(r) and X ϑ (ϑ) = a 2 cos 2 ϑ.