pense-bête de bruno sanchiz

Accueil > Physique > Trous Noirs :l’image de M87 > M87-5 Origine Physique de l’anneau asymétrique

M87-5 Origine Physique de l’anneau asymétrique

Publié le 15 avril 2019, dernière mise-à-jour le 1er mai 2019, 2 visites, visites totales.

Origine Physique de l’anneau asymétrique

Le télescope Event Horizon (EHT) a cartographié à 1,3 mm la source radio compacte au centre de la galaxie elliptique M87 avec une résolution angulaire sans précédent. Nous examinons ici les implications physiques de l’anneau asymétrique observé dans les données de l’EHT de 2017.
À cette fin, nous construisons une vaste bibliothèque de modèles basés sur des simulations magnétohydrodynamiques relativistes générales (GRMHD) et des images de synthèse produites par des lancés de rayons de relativité générale.
Nous comparons les visibilités observées avec cette bibliothèque et confirmons que l’anneau asymétrique est cohérent avec les prévisions antérieures de forte lentille gravitationnelle de l’émission synchrotron d’un plasma chaud en orbite autour de l’horizon des événements de trou noir.
Le rayon de l’anneau et l’asymétrie de l’anneau dépendent respectivement de la masse et du spin du trou noir, et les deux devraient donc être stables lorsqu’ils seront observés lors de futures campagnes EHT. Dans l’ensemble, l’image observée correspond aux prévisions concernant l’ombre d’un trou noir de Kerr en rotation, telles que prédites par la relativité générale. Si le spin du trou noir et le jet à grande échelle de M87 sont alignés, alors le vecteur de spin du trou noir est dirigé loin de la Terre. Les modèles de notre bibliothèque de trous noirs pas en rotation ne concordent pas avec les observations car ils ne produisent pas de jets suffisamment puissants. Dans le même temps, dans ces modèles qui produisent un jet suffisamment puissant, ce dernier est alimenté par extraction de l’énergie de spin du trou noir par le biais de mécanismes similaires au processus de Blandford-Znajek. Nous considérons brièvement des alternatives au trou noir pour l’objet compact central.
L’analyse des données de polarisation EHT existantes et des données prises simultanément à d’autres longueurs d’onde permettra bientôt de nouveaux tests des modèles GRMHD, de même que les futures campagnes EHT à 230 et 345 GHz.

introduction

En 1918, la galaxie Messier 87 (M87) fut observée par Curtis et fut trouvée avoir ’’’un curieux rayon droit… apparemment lié au noyau par une fine ligne de matière’’’ (Curtis 1918, p. 31). Le rayon de Curtis est maintenant connu pour être un jet, s’étendant de l’échelle sous-parsec à plusieurs kpc, et peut être observé à travers le spectre électromagnétique, de la radio jusqu’aux rayons γ. Les observations d’interférométrie à base très longues (VLBI) qui zooment sur le noyau, sondant des échelles angulaires progressivement plus petites à des fréquences progressivement plus élevées jusqu’à 86 GHz avec la matrice ?? Global mm-VLBI (GMVA ; par exemple, Hada et al. 2016 ; Boccardi et al. 2017 ; Kim et al. 2018 ; Walker et al 2018) ont révélé que le jet émerge d’un noyau central. Les modèles de la distribution de vitesse stellaire impliquent une masse pour le noyau central $M \ approx 6,2 \ times 10 ^ 9 \, M _ \ odot $$ à une distance de $16,9 ​​\, \ mathrm Mpc $$ ( Gebhardt et al., 2011) ; les modèles de raies d’émission de gaz ionisé d’échelle d’une seconde d’arc impliquent une masse qui est inférieure d’un facteur deux environ (Walsh et al. 2013).

Le modèle conventionnel de l’objet central dans M87 est un trou noir entouré par un flux d’accrétion de disque optiquement fin et géométriquement épais (par exemple, Ichimaru 1977 ; Rees et al. 1982 ; Narayan & Yi 1994, 1995 ; Reynolds et al. 1996).
Le pouvoir radiatif du flux d’accrétion dérive finalement de l’énergie de liaison gravitationnelle du plasma entrant. Il n’existe pas de modèle consensuel pour le lancement de jets, mais les deux scénarios principaux sont que le jet est un flux dominé magnétiquement qui est finalement alimenté en exploitant l’énergie de rotation du trou noir (Blandford & Znajek 1977) et que le jet est un vent collimaté magnétiquement depuis le disque d’accrétion environnant (Blandford et Payne, 1982 ; Lynden-Bell, 2006).

Les observations VLBI de M87 aux fréquences $ \gtrsim 230 \, \mathrm GHz $$ avec le Télescope Horizon des Evénements (EHT) peuvent résoudre des échelles angulaires de dizaines de $ \ mu \ mathrm as$$, comparables à l’échelle de l’horizon des événements (Doeleman et al. 2012 ; Akiyama et al. 2015 ; EHT Collaboration et al. 2019a, 2019b, 2019c, ci-après désigné : documents I, II et III).
Ils ont donc le pouvoir de sonder la nature de l’objet central et de tester des modèles de lancement de jet. En outre, les observations EHT peuvent limiter les paramètres physiques clés du système, notamment la masse et le spin du trou noir, le taux d’accrétion et le flux magnétique piégé par l’accrétion de plasma dans le trou noir.

Dans cette lettre, nous adoptons l’hypothèse de travail selon laquelle l’objet central est un trou noir décrit par la métrique de Kerr, avec une masse M et un spin sans dimension $a_ * , -1 \ lt a _ * \ lt 1 $$.
Ici, $* \ equiv Jc / GM ^ 2 $$, où J, G et c sont respectivement le moment angulaire du trou noir, la constante de gravitation et la vitesse de la lumière. Dans notre convention, $* \ lt 0 $$ implique que le moment angulaire du flux d’accrétion et celui du trou noir sont antialignés. À l’aide de modèles magnétohydrodynamiques relativistes généraux (GRMHD) pour le flux d’accrétion et d’images synthétiques de ces simulations produites par des calculs de transfert radiatif relativiste général, nous vérifions si les résultats de la campagne d’observation EHT de 2017 (ci-après EHT2017) sont compatibles avec l’hypothèse du trou noir. .

Cette lettre est organisée comme suit. Dans la section 2, nous passons en revue les principales caractéristiques des observations et fournissons des estimations par ordre de grandeur des conditions physiques de la source. Dans la section 3, nous décrivons les modèles numériques. Dans la section 4, nous décrivons notre procédure de comparaison des modèles avec les données de manière à prendre en compte la variabilité du modèle. Dans la section 5, nous montrons qu’un grand nombre de modèles ne peuvent pas être rejetés sur la base des seules données EHT. Dans la section 6, nous combinons les données EHT avec d’autres contraintes sur l’efficacité radiative, la luminosité des rayons X et la puissance des jets, et montrons que cette dernière contrainte élimine tous les modèles $a_ * = 0 $$. Dans la section 7, nous discutons des limites de nos modèles et des solutions de rechange aux modèles de trou noir de Kerr. Dans la section 8, nous résumons nos résultats et expliquons comment une analyse plus poussée des données EHT existantes, des futures données EHT et des observations complémentaires sur plusieurs longueurs d’onde accentueront les contraintes sur les modèles.

Review and Estimates

Dans la Collaboration EHT et al. (2019d ; ci-après Paper IV), nous présentons des images générées à partir des données EHT2017 (pour plus de détails sur le tableau, la campagne d’observation de 2017, la corrélation et l’étalonnage, voir les documents II et III). Une image représentative est reproduite dans le panneau de gauche de la figure 1.

Figure 1. Left panel : an EHT2017 image of M87 from Paper IV of this series (see their Figure 15). Middle panel : a simulated image based on a GRMHD model. Right panel : the model image convolved with a $20\,\mu \mathrmas$ FWHM Gaussian beam. Although the most evident features of the model and data are similar, fine features in the model are not resolved by EHT.

Quatre caractéristiques de l’image sur le panneau de gauche de la figure 1 jouent un rôle important dans notre analyse : (1) la géométrie en anneau, (2) la température de brillance de pointe, (3) la densité de flux totale et (4) la asymétrie de la bague. Nous considérons maintenant chacun à son tour.

(1) La source compacte présente un anneau lumineux avec une zone sombre centrale sans composants étendus significatifs. Cela présente une similitude remarquable avec la structure longtemps prédite pour l’émission optiquement mince d’un plasma chaud entourant un trou noir (Falcke et al. 2000). Le trou central entouré d’un anneau brillant résulte d’une forte lentille gravitationnelle (par exemple, Hilbert 1917 ; von Laue 1921 ; Bardeen 1973 ; Luminet 1979). Le « cercle de photons » correspond aux lignes de vision qui passent près des orbites de photons (instables) (voir Teo 2003), s’attardent près de l’orbite des photons et ont donc un long trajet à travers le plasma émetteur. Ces lignes de visée apparaîtront relativement brillantes si le plasma émetteur est optiquement mince. La dépression du flux central est ce qu’on appelle le"trou noir"("trou noir") (Falcke et al. 2000) et correspond aux lignes de mire qui se terminent à l’horizon des événements. L’ombre pourrait être vue contrairement à l’émission environnante provenant du flux d’accrétion ou du contre-jet à lentilles dans M87 (Broderick & Loeb 2009).

L’anneau de photons est presque circulaire pour toutes les rotations des trous noirs et toutes les inclinaisons de l’axe de rotation des trous noirs par rapport à la ligne de mire (par exemple, Johannsen & Psaltis 2010). Pour un _ * = 0 $ trou noir de masse $ et de distance D, le rayon angulaire de l’anneau du photon dans le ciel est
Équation (1)

où nous avons réduit à la masse la plus probable de Gebhardt et al. (2011) et une distance de 6,9 ​​\, \ mathrm Mpc $ (voir aussi EHT Collaboration et al. 2019e, (ci-après Article VI ; Blakeslee et al. 2009 ; Bird et al. 2010 ; Cantiello et al. 2018). Le rayon angulaire de l’anneau du photon pour les autres inclinaisons et les valeurs de _ * $ diffère d’au plus 13% de l’équation (1), et la majeure partie de cette variation se produit à - | a _ * | \ ll 1 $ ( par exemple, Takahashi 2004 ; Younsi et al 2016). À l’évidence, le rayon angulaire de l’anneau de photons observé est approximativement égal à $ \ sim 20 \, \ mu \ mathrm as $ (Figure 1 et papier IV), prédiction du modèle de trou noir donné dans l’équation (1).

(2) The observed peak brightness temperature of the ring in Figure 1 is $T_b,pk\sim 6\times 10^9\,\rmK$, which is consistent with past EHT mm-VLBI measurements at 230 GHz (Doeleman et al. 2012 ; Akiyama et al. 2015), and GMVA 3 mm-VLBI measurements of the core region (Kim et al. 2018). Expressed in electron rest-mass (me) units, $\rm\Theta

_b,pk\equiv k_\rmB

T_b,pk/(m_ec^2)\simeq 1$, where $k_\rmB

$ is Boltzmann’s constant. The true peak brightness temperature of the source is higher if the ring is unresolved by EHT, as is the case for the model image in the center panel of Figure 1.

L’émission de 1,3 mm de M87 illustrée à la figure 1 devrait être générée par le processus synchrotron (voir Yuan & Narayan 2014 et les références qui y figurent) et dépend donc de la fonction de distribution d’électrons (eDF). Si le plasma émetteur a un eDF thermique, il est caractérisé par une température électronique
$T_e\geqslant T_b$, or $\rm\Theta

_e\equiv k_\rmB

T_e/(m_ec^2)\gt 1$, because $\rm\Theta

_e\gt \rm\Theta

_b,pk$ if the ring is unresolved or optically thin.

La température de luminosité observée est-elle conforme à ce que l’on pourrait attendre de modèles phénoménologiques de la source ? Les modèles d’écoulement d’accrétion inefficaces sur le plan radiatif de M87 (Reynolds et al. 1996 ; Di Matteo et al. 2003) produisent une émission en mm dans un beignet de plasma géométriquement épais autour du trou noir. Le plasma émetteur est sans collision : la diffusion de Coulomb est faible à ces faibles densités et hautes températures. Par conséquent, les températures des électrons et des ions ne doivent pas nécessairement être les mêmes (par exemple, Spitzer 1962). Dans les modèles de flux d’accrétion inefficaces par radiation, la température des ions est légèrement inférieure à la température viriale des ions,

Equation (2)

where $r_\rmg

\equiv GM/c^2$ is the gravitational radius, r is the Boyer–Lindquist or Kerr–Schild radius, and mp is the proton mass. Most models have an electron temperature $T_e\lt T_i$ because of electron cooling and preferential heating of the ions by turbulent dissipation (e.g., Yuan & Narayan 2014 ; Mościbrodzka et al. 2016). If the emission arises at $\sim 5\,r_\rmg

$, then $\rm\Theta

_e\simeq 37(T_e/T_i)$, which is then consistent with the observed $\rm\Theta

_b,pk$ if the source is unresolved or optically thin.

(3)La densité de flux totale dans l’image à.3 \, \ mathrm mm $ est $ \ simeq 0.5 $ Jy. Avec quelques hypothèses, nous pouvons utiliser ceci pour estimer la densité de nombre d’électrons et l’intensité du champ magnétique B dans la source. Nous adoptons un modèle simple, sphérique, à une zone pour la source avec rayon
$r\simeq 5\,r_\rmg

$, pressure $n_ikT_i+n_ekT_e=\beta _\rmp

B^2/(8\pi )$ with $\beta _\rmp

\equiv p_\mathrmgas/p_\mathrmmag\sim 1$, $T_i\simeq 3T_e$, and temperature $\theta _e\simeq 10\theta _b,pk$, ce qui est conforme à la discussion dans (2) ci-dessus. En fixant ne = ni (c’est-à-dire en supposant un plasma d’hydrogène totalement ionisé), les valeurs de B et ne nécessaires pour produire la densité de flux observée peuvent être trouvées en résolvant une équation non linéaire (en supposant un angle moyen entre le champ et la ligne de visée, 60 °). La solution peut être approchée comme une loi de puissance :

Equation (3)

Equation (4)

en supposant que = 6.2 \ times 10 ^ 9 \, M _ \ odot $ et = 16.9 \, \ mathrm Mpc $, et en utilisant l’émissivité thermique approximative de Leung et al. (2011). Ensuite, la profondeur optique du synchrotron à.3 \, \ mathrm mm $ est 0.2. On peut maintenant estimer un taux d’accrétion à partir de (3) en utilisant

Equation (5)

assuming spherical symmetry. The Eddington accretion rate is
Equation (6)

where $L_\mathrmEdd\equiv 4\pi GMcm_p/\sigma _T$ is the Eddington luminosity ($\sigma _T$ is the Thomson cross section). Setting the efficiency $\epsilon =0.1$ and $M=6.2\times 10^9\,M_\odot $, $\dotM_\mathrmEdd=137\,M_\odot \,\mathrmyr^-1$, and therefore $\dotM/\dotM_\mathrmEdd\,\sim 2.0\times 10^-5$.

Cette estimation est similaire mais légèrement supérieure à la limite supérieure déduite des propriétés de polarisation linéaire à 230 GHz de M87 (Kuo et al. 2014).

(4)L’anneau est plus brillant dans le sud que le nord. Cela peut s’expliquer par une combinaison de mouvement dans la source et de faisceau Doppler. Comme exemple simple, nous considérons un anneau lumineux optiquement mince tournant avec la vitesse v et un vecteur de moment angulaire incliné à un angle de vue i
> 0° to the line of sight.Ensuite, le côté approchant de l’anneau est dopé par Doppler et le côté reculant est atténué en Doppler, produisant un contraste de luminosité de surface d’un ordre égal à un si v est relativiste. L’approche du jet à grande échelle dans M87 est orientée ouest-nord-ouest (angle de position $ \ mathrm PA \ environ 288 ^ \ circ ; $ dans le papier VI, on l’appelle
$\mathrmPA_\mathrmFJ$), ou à droite et légèrement en haut dans l’image. Walker et al. (2018) ont estimé que l’angle entre le jet en approche et la ligne de mire était de 17 °. Si l’émission est produite par un anneau tournant avec un vecteur moment cinétique orienté le long de l’axe du jet, le plasma au sud s’approche de la Terre et le plasma au nord recule. Cela implique une circulation dans le sens des aiguilles d’une montre du plasma dans la source, tel que projeté sur le plan du ciel. Ce sens de rotation est cohérent avec le sens de rotation dans le gaz ionisé à des échelles d’arcsecondes (Harms et al. 1994 ; Walsh et al. 2013). Notez que l’asymétrie de l’anneau est compatible avec l’asymétrie déduite d’observations à 43 GHz du rapport de luminosité entre les côtés nord et sud du jet et du contre-jet (Walker et al. 2018).

Toutes ces estimations donnent une image de la source qui correspond remarquablement aux attentes du modèle de trou noir et aux modèles GRMHD existants (par exemple, Dexter et al. 2012 ; Mościbrodzka et al. 2016). Ils suggèrent même une sensation de rotation du gaz près du trou noir. Une comparaison quantitative avec les modèles GRMHD peut en révéler davantage.

models

Conformément à la discussion à la section 2, nous adoptons maintenant l’hypothèse de travail selon laquelle M87 contient un flux d’accrétion turbulent et magnétisé entourant un trou noir de Kerr. Pour tester quantitativement cette hypothèse par rapport aux données EHT2017, nous avons généré une bibliothèque de simulation de modèles GRMHD idéaux 3D dépendant du temps. Pour générer efficacement cette bibliothèque coûteuse en calculs et avec des contrôles indépendants des résultats, nous avons utilisé plusieurs codes différents qui ont évolué pour correspondre aux conditions initiales à l’aide des équations de GRMHD idéal. Les codes utilisés incluent BHAC (Porth et al. 2017), H-AMR (Liska et al. 2018 ; K. Chatterjee et al. 2019, en préparation), iharm (Gammie et al. 2003) et KORAL (Sa̧dowski et al. 2013b, 2014). On trouvera une comparaison de ces codes et d’autres codes GRMHD dans O. Porth et al. 2019 (en préparation), qui montre que les différences entre les intégrations d’un modèle d’accrétion standard avec des codes différents sont moins importantes que les fluctuations dans les simulations individuelles

À partir de la bibliothèque de simulation, nous avons généré une grande bibliothèque d’images synthétiques. Des instantanés des évolutions de GRMHD ont été produits en utilisant les schémas de traçage de rayons relativistes généraux (GRRT) ipole (Mościbrodzka & Gammie 2018), RAPTOR (Bronzwaer et al. 2018) ou BHOSS (Z. Younsi et al. 2019b, en préparation). On peut trouver une comparaison de ces codes et d’autres codes GRRT dans Gold et al. (2019), qui montre que les différences entre les codes sont faibles.

Dans les modèles GRMHD, le gros de l’émission de 1,3 mm est produit à l’intérieur du $ \ lesssim 10 \, r _

r

du trou noir, où les modèles peuvent atteindre un état statistiquement stable. Il est donc possible de calculer des modèles radiatifs prédictifs pour cette composante compacte de la source sans représenter avec précision le flux d’accrétion à tous les rayons.

Nous notons que les modèles actuels de pointe pour M87 sont des modèles GRMHD de rayonnement qui incluent la rétroaction radiative et la thermodynamique électron-ion (Ryan et al. 2018 ; Chael et al. 2019). Ces modèles sont trop coûteux en calcul pour une vaste étude de l’espace des paramètres, de sorte que, dans la présente Lettre, nous ne considérons que les modèles GRMHD non radiatifs avec un traitement paramétré de la thermodynamique des électrons.

3.1. Simulation Library
Toutes les simulations GRMHD sont initialisées avec un tore de plasma faiblement magnétisé en orbite autour du trou équatorial du trou noir (par exemple, De Villiers et al. 2003 ; Gammie et al. 2003 ; McKinney & Blandford 2009 ; Porth et al. 2017). Nous ne considérons pas les modèles inclinés, dans lesquels le moment angulaire du flux d’accrétion est mal aligné avec la rotation du trou noir. Les limites de cette approche sont examinées à la section 7.

The initial torus is driven to a turbulent state by instabilities, including the magnetorotational instability (see e.g., Balbus & Hawley 1991). In all cases the outcome contains a moderately magnetized midplane with orbital frequency comparable to the Keplerian orbital frequency, a corona with gas-to-magnetic-pressure ratio $\beta _\rmp

\equiv p_\mathrmgas/p_\mathrmmag\sim 1$, and a strongly magnetized region over both poles of the black hole with $B^2/\rho c^2\gg 1$. We refer to the strongly magnetized region as the funnel, and the boundary between the funnel and the corona as the funnel wall (De Villiers et al. 2005 ; Hawley & Krolik 2006). All models in the library are evolved from t = 0 to $t=10^4\,r_\rmg

c^-1$.

Le résultat de la simulation dépend de l’intensité et de la géométrie initiales du champ magnétique, dans la mesure où celles-ci affectent le flux magnétique traversant le disque, comme indiqué ci-dessous. Une fois la simulation lancée, le disque passe dans un état turbulent et perd la mémoire de la plupart des détails des conditions initiales. Cet état de turbulence décontracté se trouve à l’intérieur d’un rayon caractéristique qui se développe au cours de la simulation. Pour être sûr de ne représenter que les régions assouplies, nous établissons des instantanés pour les comparer avec les données de \ times 10 ^ 3 \ le t / r _ \ rm g

c ^ - 1 \ le 10 ^ 4 $.

GRMHD models have two key physical parameters. The first is the black hole spin $a_* $, $-1\lt a_* \lt 1$. The second parameter is the absolute magnetic flux $\rm\Phi

_\mathrmBH$ crossing one hemisphere of the event horizon (see Tchekhovskoy et al. 2011 ; O. Porth et al. 2019, in preparation for a definition). It is convenient to recast $\rm\Phi

_\mathrmBH$ in dimensionless form $\phi \equiv \rm\Phi

_\mathrmBH\left(\dotMr_\rmg

^2c\right)^-1/2$.110

Le flux magnétique phgr est non nul car le champ magnétique est envoyé dans l’horizon des événements par le flux d’accrétion et maintenu par les courants dans le plasma environnant. À $ \ phi \ gt \ phi _ \ max \ sim 15 $, 111 simulations numériques montrent que le flux magnétique accumulé explose, écarte le flux d’accrétion et s’échappe (Tchekhovskoy et al. 2011 ; McKinney et al. 2011). 2012). Les modèles avec $ \ phi \ sim 1 $ sont habituellement appelés modèles d’évolution standard et normale (SANE ; Narayan et al. 2012 ; Sa̧dowski et al (2013a)) ; les modèles avec $ \ phi \ sim \ phi _ \ max $ sont classiquement appelés modèles de disques à arrêts magnétiques (MAD ; Igumenshchev et al. 2003 ; Narayan et al. 2003).

La bibliothèque de simulation contient des modèles SANE avec _ * = - 0,94 $, −0,5, 0, 0,5, 0,75, 0,88, 0,94, 0,97, et 0,98 et des modèles MAD avec _ * = - 0,94 $, −0,5, 0, 0,5, 0,75 et 0,94. La bibliothèque de simulation occupe 23 To d’espace disque et contient un total de 43 simulations GRMHD, certaines répétées à plusieurs résolutions avec plusieurs codes, avec des résultats cohérents (O. Porth et al. 2019, en préparation).

3.2. Image Library Generation

Pour produire des images modèles à partir des simulations afin de les comparer aux observations EHT, nous utilisons GRRT pour générer un grand nombre d’images de synthèse et de produits de données VLBI dérivés. Pour créer les images de synthèse, nous devons spécifier les éléments suivants : (1) le champ magnétique, le champ de vitesse et la densité en fonction de la position et de la durée ; (2) les coefficients d’émission et d’absorption en fonction de la position et du temps ; et (3) l’angle d’inclinaison entre le vecteur moment cinétique du flux d’accrétion et la ligne de visée i, l’angle de position $ \ mathrm PA $, la masse du trou noir $ et la distance D à l’observateur. Dans ce qui suit, nous discutons chaque entrée. Le lecteur qui n’est intéressé que par une description de haut niveau de la bibliothèque d’images peut passer à la section 3.3.

(1) GRMHD models provide the absolute velocity field of the plasma flow. Nonradiative GRMHD evolutions are invariant, however, under a rescaling of the density by a factor $\mathscrM$. In particular, they are invariant under $\rho \to \mathscrM\rho $, field strength $B\to \mathscrM

^1/2B$, and internal energy $u\to \mathscrMu$ (the Alfvén speed $B/\rho ^1/2$ and sound speed $\propto \sqrtu/\rho $ are invariant). That is, there is no intrinsic mass scale in a nonradiative model as long as the mass of the accretion flow is negligible in comparison to $M$.112 We use this freedom to adjust $\mathscrM$ so that the average image from a GRMHD model has a 1.3 mm flux density ≈0.5 Jy (see Paper IV). Once $\mathscrM$ is set, the density, internal energy, and magnetic field are fully specified.

The mass unit $\mathscrM$ determines $\dotM$. In our ensemble of models $\dotM$ ranges from $2\times 10^-7\dotM_\mathrmEdd$ to $4\times 10^-4\dotM_\mathrmEdd$. Accretion rates vary by model category. The mean accretion rate for MAD models is $\sim 10^-6\dotM_\mathrmEdd$. For SANE models with $a_* \gt 0$ it is $\sim 5\times 10^-5\dotM_\mathrmEdd ;$ and for $a_* \lt 0$ it is $\sim 2\times 10^-4\dotM_\mathrmEdd$.

(2) Les distributions d’énergie spectrale radioélectrique observées (SED) et les caractéristiques de polarisation de la source indiquent clairement que l’émission de 1,3 mm correspond à un rayonnement synchrotron, typique des noyaux galactiques actifs (AGN). Les coefficients d’absorption et d’émission du synchrotron dépendent du eDF. Dans ce qui suit, nous adoptons un modèle thermique relativiste pour le format eDF (distribution Maxwell-Jüttner ; Jüttner 1911 ; Rezzolla & Zanotti 2013). Nous discutons des limites de cette approche à la section 7.

Tous nos modèles de M87 ont une densité et une température suffisamment basses et élevées pour que le plasma soit sans collision (voir Ryan et al. 2018, pour une discussion du couplage Coulomb en M87). Par conséquent, Te n’est probablement pas égal à la température de L’ion Ti, qui est fournie par les simulations. Nous définissons Te en utilisant la densité de GRMHD ρ, la densité d’énergie interne u, et le plasma
$\beta _\rmp

$ using a simple model :
Equation (7)

où nous avons supposé que le plasma est composé d’hydrogène, les ions sont non relativiste, et les électrons sont relativistes. Here $R\equiv T_i/T_e$ and
Equation (8)

This prescription has one parameter, $R_\mathrmhigh$, and sets $T_e\simeq T_i$ in low $\beta _\rmp

$ regions and $T_e\simeq T_i/R_\mathrmhigh$ in the midplane of the disk. It is adapted from Mościbrodzka et al. (2016) and motivated by models for electron heating in a turbulent, collisionless plasma that preferentially heats the ions for $\beta _\rmp

\gtrsim 1$ (e.g., Howes 2010 ; Kawazura et al. 2018).

(3)Nous devons spécifier l’inclinaison de l’observateur i, l’orientation de l’observateur à travers l’angle de position $\mathrmPA$, la masse du trou noir$, et la distance D à la source. Les contraintes Non-EHT sur i, $\mathrmPA$, et $ sont considérées ci-dessous ; nous avons généré des images à =12^\circ ,17^\circ ,22^\circ ,158^\circ ,163^\circ $, et 168° et quelques-unes à i = 148°. L’angle de position (PA) peut être modifié en tournant simplement l’image. Toutes les caractéristiques des modèles que nous avons examinés, y compris $\dotm$, sont insensibles aux petites variations de I. La morphologie de l’image dépend de si je suis supérieur ou inférieur à 90°, comme nous le montrerons ci-dessous.\n

Les images modèles sont générées avec un champ de vision de 60 \ times 160 \, \ mu \ mathrm as $ et \ mu \ mathrm as $ pixels, ce qui est petit par rapport au $ \ sim 20 \, \ mu \ mathrm as $ résolution nominale de EHT2017. Notre analyse est insensible aux changements de champ de vision et d’échelle de pixel.

Pour $, nous utilisons la valeur la plus probable du travail sur la ligne d’absorption stellaire,.2 \ times 10 ^ 9 M _ \ odot $ (Gebhardt et al. 2011). Pour la distance D, nous utilisons 6,9 \, \, \ mathrm Mpc $, ce qui est très proche de celui utilisé dans le papier VI. Le rapport /(c^2BuchD)3.62\,\mu \ mathrm as $ (ci-après M / D) détermine l’échelle angulaire des images. Pour certains modèles, nous avons également généré des images avec = 3.5 \ times 10 ^ 9 \, M _ \ odot $ pour vérifier que les résultats de l’analyse ne sont pas prédéterminés par la masse du trou noir entrée.

3.3. Image Library Summary

La bibliothèque d’images contient environ 60 000 images. Nous générons des images de 100 à 500 fichiers de sortie distincts à partir de chacun des modèles GRMHD à chacun des _ \ mathrm high = 1,10,20,40,80 $ et 160. En comparant les données, nous ajustons le $ \ mathrm PA $ par rotation, ainsi que le flux total et l’échelle angulaire de l’image en redimensionnant simplement les images à partir des paramètres standard de la bibliothèque d’images (voir la figure 29 du document VI). Les tests indiquent que les comparaisons avec les données sont insensibles à la procédure de redimensionnement sauf si le facteur de mise à l’échelle angulaire ou le facteur de mise à l’échelle du flux est grand.113

The comparisons with the data are also insensitive to image resolution.114

Les figures 2 et 3 présentent un ensemble représentatif d’images de la bibliothèque d’images présentant une moyenne temporelle. Ces figures montrent clairement que la modification des paramètres _ * $, phgr et _ \ mathrm high $ peut changer la largeur et l’asymétrie de l’anneau à photons et introduisent des structures supplémentaires extérieures et intérieures à l’anneau à photons.

Figure 2. Time-averaged 1.3 mm images generated by five SANE GRMHD simulations with varying spin ($a_* =-0.94$ to $a_* =+0.97$ from left to right) and $R_\mathrmhigh$ ($R_\mathrmhigh=1$ to $R_\mathrmhigh=160$ from top to bottom ; increasing $R_\mathrmhigh$ corresponds to decreasing electron temperature). The colormap is linear. All models are imaged at i = 163°. The jet that is approaching Earth is on the right (west) in all the images. The black hole spin vector projected onto the plane of the sky is marked with an arrow and aligned in the east–west direction. When the arrow is pointing left the black hole rotates in a clockwise direction, and when the arrow is pointing right the black hole rotates in a counterclockwise direction. The field of view for each model image is $80\,\mu \mathrmas$ (half of that used for the image libraries) with resolution equal to $1\mu \mathrmas$/pixel (20 times finer than the nominal resolution of EHT2017, and the same employed in the library images).

Figure 3. Same as in Figure 2 but for selected MAD models

The location of the emitting plasma is shown in Figure 4, which shows a map of time- and azimuth-averaged emission regions for four representative $a_* \gt 0$ models. For SANE models, if $R_\mathrmhigh$ is low (high), emission is concentrated more in the disk (funnel wall), and the bright section of the ring is dominated by the disk (funnel wall).115 Appendix B shows images generated by considering emission only from particular regions of the flow, and the results are consistent with Figure 4.

Figure 4. Binned location of the point of origin for all photons that make up an image, summed over azimuth, and averaged over all snapshots from the simulation. The colormap is linear. The event horizon is indicated by the solid white semicircle and the black hole spin axis is along the figure vertical axis. This set of four images shows MAD and SANE models with $R_\mathrmhigh=10$ and 160, all with $a_* =0.94$. The region between the dashed curves is the locus of existence of (unstable) photon orbits (Teo 2003). The green cross marks the location of the innermost stable circular orbit (ISCO) in the equatorial plane. In these images the line of sight (marked by an arrow) is located below the midplane and makes a 163° angle with the disk angular momentum, which coincides with the spin axis of the black hole.

Figures 2 and 3 show that for both MAD and SANE models the bright section of the ring, which is generated by Doppler beaming, shifts from the top for negative spin, to a nearly symmetric ring at $a_* =0$, to the bottom for $a_* \gt 0$ (except the SANE $R_\mathrmhigh=1$ case, where the bright section is always at the bottom when i > 90°). That is, the location of the peak flux in the ring is controlled by the black hole spin : it always lies roughly 90 degrees counterclockwise from the projection of the spin vector on the sky. Some of the ring emission originates in the funnel wall at $r\lesssim 8\,r_\rmg

$. The rotation of plasma in the funnel wall is in the same sense as plasma in the funnel, which is controlled by the dragging of magnetic field lines by the black hole. The funnel wall thus rotates opposite to the accretion flow if $a_* \lt 0$. This effect will be studied further in a later publication (Wong et al. 2019). The resulting relationships between disk angular momentum, black hole angular momentum, and observed ring asymmetry are illustrated in Figure 5.

Figure 5. Illustration of the effect of black hole and disk angular momentum on ring asymmetry. The asymmetry is produced primarily by Doppler beaming : the bright region corresponds to the approaching side. In GRMHD models that fit the data comparatively well, the asymmetry arises in emission generated in the funnel wall. The sense of rotation of both the jet and funnel wall are controlled by the black hole spin. If the black hole spin axis is aligned with the large-scale jet, which points to the right, then the asymmetry implies that the black hole spin is pointing away from Earth (rotation of the black hole is clockwise as viewed from Earth). The blue ribbon arrow shows the sense of disk rotation, and the black ribbon arrow shows black hole spin. Inclination i is defined as the angle between the disk angular momentum vector and the line of sight.

The time-averaged MAD images are almost independent of $R_\mathrmhigh$ and depend mainly on $a_* $. In MAD models much of the emission arises in regions with $\beta _\rmp

\sim 1$, where $R_\mathrmhigh$ has little influence over the electron temperature, so the insensitivity to $R_\mathrmhigh$ is natural (see Figure 4). In SANE models emission arises at $\beta _\rmp

\sim 10$, so the time-averaged SANE images, by contrast, depend strongly on $R_\mathrmhigh$. In low $R_\mathrmhigh$ SANE models, extended emission outside the photon ring, arising near the equatorial plane, is evident at $R_\mathrmhigh=1$. In large $R_\mathrmhigh$ SANE models the inner ring emission arises from the funnel wall, and once again the image looks like a thin ring (see Figure 4).

Figure 6 and the accompanying animation show the evolution of the images, visibility amplitudes, and closure phases over a $5000\,r_\rmg

c^-1\approx 5\,\mathrmyr$ interval in a single simulation for M87. It is evident from the animation that turbulence in the simulations produces large fluctuations in the images, which imply changes in visibility amplitudes and closure phases that are large compared to measurement errors. The fluctuations are central to our procedure for comparing models with the data, described briefly below and in detail in Paper VI.

Figure 6. Single frame from the accompanying animation. This shows the visibility amplitudes (top), closure phases plotted by Euclidean distance in 6D space (middle), and associated model images at full resolution (lower left) and convolved with the EHT2017 beam (lower right). Data from 2017 April 6 high-band are also shown in the top two plots. The video shows frames 1 through 100 and has a duration of 10 s.

The timescale between frames in the animation is $50\,r_\rmg

c^-1\simeq 18$ days, which is long compared to EHT2017 observing campaign. The images are highly correlated on timescales less than the innermost stable circular orbit (ISCO) orbital period, which for $a_* =0$ is $\simeq 15\ r_\rmg

c^-1\simeq 5$ days, i.e., comparable to the duration of the EHT2017 campaign. If drawn from one of our models, we would expect the EHT2017 data to look like a single snapshot (Figures 6) rather than their time averages (Figures 2 and 3).
4. Procedure for Comparison of Models with Data

As described above, each model in the Simulation Library has two dimensionless parameters : black hole spin $a_* $ and magnetic flux phgr. Imaging the model from each simulation adds five new parameters : $R_\mathrmhigh$, i, $\mathrmPA$, $M$, and D, which we set to $16.9\,\mathrmMpc$. After fixing these parameters we draw snapshots from the time evolution at a cadence of 10 to $50\,r_\rmg

c^-1$. We then compare these snapshots to the data.

The simplest comparison computes the $\chi _\nu ^2$ (reduced chi square) distance between the data and a snapshot. In the course of computing $\chi _\nu ^2$ we vary the image scale M/D, flux density Fν, position angle $\mathrmPA$, and the gain at each VLBI station in order to give each image every opportunity to fit the data. The best-fit parameters $(M/D,F_\nu ,\mathrmPA)$ for each snapshot are found by two pipelines independently : the Themis pipeline using a Markov chain Monte Carlo method (A. E. Broderick et al. 2019a, in preparation), and the GENA pipeline using an evolutionary algorithm for multidimensional minimization (Fromm et al. 2019a ; C. Fromm et al. 2019b, in preparation ; see also Section 4 of Paper VI for details). The best-fit parameters contain information about the source and we use the distribution of best-fit parameters to test the model by asking whether or not they are consistent with existing measurements of M/D and estimates of the jet $\mathrmPA$ on larger scales.

The $\chi _\nu ^2$ comparison alone does not provide a sharp test of the models. Fluctuations in the underlying GRMHD model, combined with the high signal-to-noise ratio for EHT2017 data, imply that individual snapshots are highly unlikely to provide a formally acceptable fit with $\chi _\nu ^2\simeq 1$. This is borne out in practice with the minimum $\chi _\nu ^2=1.79$ over the entire set of the more than 60,000 individual images in the Image Library. Nevertheless, it is possible to test if the $\chi _\nu ^2$ from the fit to the data is consistent with the underlying model, using ’’’Average Image Scoring’’’ with Themis (Themis-AIS), as described in detail in Appendix F of Paper VI). Themis-AIS measures a $\chi _\nu ^2$ distance (on the space of visibility amplitudes and closure phases) between a trial image and the data. In practice we use the average of the images from a given model as the trial image (hence Themis-AIS), but other choices are possible. We compute the $\chi _\nu ^2$ distance between the trial image and synthetic data produced from each snapshot. The model can then be tested by asking whether the data’s $\chi _\nu ^2$ is likely to have been drawn from the model’s distribution of $\chi _\nu ^2$. In particular, we can assign a probability p that the data is drawn from a specific model’s distribution.

In this Letter we focus on comparisons with a single data set, the 2017 April 6 high-band data (Paper III). The eight EHT2017 data sets, spanning four days with two bands on each day, are highly correlated. Assessing what correlation is expected in the models is a complicated task that we defer to later publications. The 2017 April 6 data set has the largest number of scans, 284 detections in 25 scans (see Paper III) and is therefore expected to be the most constraining.116
5. Model Constraints : EHT2017 Alone

The resolved ring-like structure obtained from the EHT2017 data provides an estimate of M/D (discussed in detail in Paper VI) and the jet $\mathrmPA$ from the immediate environment of the central black hole. As a first test of the models we can ask whether or not these are consistent with what is known from other mass measurements and from the orientation of the large-scale jet.

Figure 7 shows the distributions of best-fit values of M/D for a subset of the models for which spectra and jet power estimates are available (see below). The three lines show the M/D distribution for all snapshots (dotted lines), the best-fit 10% of snapshots (dashed lines), and the best-fit 1% of snapshots (solid lines) within each model. Evidently, as better fits are required, the distribution narrows and peaks close to $M/D\sim 3.6\,\mu \mathrmas$ with a width of about $0.5\mu \mathrmas$.

Figure 7. Distribution of M/D obtained by fitting Image Library snapshots to the 2017 April 6 data, in $\mu \mathrmas$, measured independently using the (left panel) Themis and (right panel) GENA pipelines with qualitatively similar results. Smooth lines were drawn with a Gaussian kernel density estimator. The three lines show the best-fit 1% within each model (solid) ; the best-fit 10% within each model (dashed) ; and all model images (dotted). The vertical lines show $M/D=2.04$ (dashed) and $3.62\,\mu \mathrmas$ (solid), corresponding to M = 3.5 and $6.2\times 10^9\,M_\odot $. The distribution uses a subset of models for which spectra and jet power estimates are available (see Section 6). Only images with $a_* \gt 0$, i > 90° and $a_* \lt 0$, i < 90° (see also the left panel of Figure 5) are considered

The distribution of M/D for the best-fit $\lt 10 \% $ of snapshots is qualitatively similar if we include only MAD or SANE models, only models produced by individual codes (BHAC, H-AMR, iharm, or KORAL), or only individual spins. As the thrust of this Letter is to test the models, we simply note that Figure 7 indicates that the models are broadly consistent with earlier mass estimates (see Paper VI for a detailed discussion). This did not have to be the case : the ring radius could have been significantly larger than $3.6\,\mu \mathrmas$.

We can go somewhat further and ask if any of the individual models favor large or small masses. Figure 8 shows the distributions of best-fit values of M/D for each model (different $a_* $, $R_\mathrmhigh$, and magnetic flux). Most individual models favor M/D close to $3.6\,\mu \mathrmas$. The exceptions are $a_* \leqslant 0$ SANE models with $R_\mathrmhigh=1$, which produce the bump in the M/D distribution near $2\mu \mathrmas$. In these models, the emission is produced at comparatively large radius in the disk (see Figure 2) because the inner edge of the disk (the ISCO) is at a large radius in a counter-rotating disk around a black hole with $| a_* | \sim 1$. For these models, the fitting procedure identifies EHT2017’s ring with this outer ring, which forces the photon ring, and therefore M/D, to be small. As we will show later, these models can be rejected because they produce weak jets that are inconsistent with existing jet power estimates (see Section 6.3).

Figure 8. Distributions of M/D and black hole mass with $D=16.9\,\mathrmMpc$ reconstructed from the best-fit 10% of images for MAD (left panel) and SANE (right panel) models (i = 17° for $a_* \le 0$ and 163° for $a_* \gt 0$) with different $R_\mathrmhigh$ and $a_* $, from the Themis (dark red, left), and GENA (dark green, right) pipelines. The white dot and vertical black bar correspond, respectively, to the median and region between the 25th and 75th percentiles for both pipelines combined. The blue and pink horizontal bands show the range of M/D and mass at $D=16.9\,\mathrmMpc$ estimated from the gas dynamical model (Walsh et al. 2013) and stellar dynamical model (Gebhardt et al. 2011), respectively. Constraints on the models based on average image scoring (Themis-AIS) are discussed in Section 5. Constraints based on radiative efficiency, X-ray luminosity, and jet power are discussed in Section 6.

Figure 8 also shows that M/D increases with $a_* $ for SANE models. This is due to the appearance of a secondary inner ring inside the main photon ring. The former is associated with emission produced along the wall of the approaching jet. Because the emission is produced in front of the black hole, lensing is weak and it appears at small angular scale. The inner ring is absent in MAD models (see Figure 3), where the bulk of the emission comes from the midplane at all values of $R_\mathrmhigh$(Figure 4).

We now ask whether or not the PA of the jet is consistent with the orientation of the jet measured at other wavelengths. On large ( mas) scales the extended jet component has a PA of approximately 288° (e.g., Walker et al. 2018). On smaller ($\sim 100\,\mu \mathrmas$) scales the apparent opening angle of the jet is large (e.g., Kim et al. 2018) and the PA is therefore more difficult to measure. Also notice that the jet PA may be time dependent (e.g., Hada et al. 2016 ; Walker et al. 2018). In our model images the jet is relatively dim at 1.3 mm, and is not easily seen with a linear colormap. The model jet axis is, nonetheless, well defined : jets emerge perpendicular to the disk.

Figure 9 shows the distribution of best-fit PA over the same sample of snapshots from the Image Library used in Figure 7. We divide the snapshots into two groups. The first group has the black hole spin pointed away from Earth (i > 90° and $a_* \gt 0$, or i < 90° and $a_* \lt 0$). The spin-away model PA distributions are shown in the top two panels. The second group has the black hole spin pointed toward Earth (i > 90 and $a_* \lt 0$ or i > 90° and $a_* \lt 0$). These spin-toward model PA distributions are shown in the bottom two panels. The large-scale jet orientation lies on the shoulder of the spin-away distribution (the distribution can be approximated as a Gaussian with, for Themis (GENA) mean 209 (203)° and $\sigma _\mathrmPA\,=54\,(55)^\circ ;$ the large-scale jet PA lies $1.5\sigma _\mathrmPA$ from the mean) and is therefore consistent with the spin-away models. On the other hand, the large-scale jet orientation lies off the shoulder of the spin-toward distribution and is inconsistent with the spin-toward models. Evidently models in which the black hole spin is pointing away from Earth are strongly favored.

Figure 9. Top : distribution of best-fit PA (in degree) scored by the Themis (left) and GENA (right) pipelines for models with black hole spin vector pointing away from Earth (i > 90° for $a_* \gt 0$ or i < 90° for $a_* \lt 0$). Bottom : images with black hole spin vector pointing toward Earth (i < 90° for $a_* \gt 0$ or i > 90° for $a_* \lt 0$). Smooth lines were drawn with a wrapped Gaussian kernel density estimator. The three lines show (1) all images in the sample (dotted line) ; (2) the best-fit 10% of images within each model (dashed line) ; and (3) the best-fit 1% of images in each model (solid line). For reference, the vertical line shows the position angle $\mathrmPA\sim 288^\circ $ of the large-scale (mas) jet Walker et al. (2018), with the gray area from (288 – 10)° to (288 + 10)° indicating the observed PA variation.

The width of the spin-away and spin-toward distributions arises naturally in the models from brightness fluctuations in the ring. The distributions are relatively insensitive if split into MAD and SANE categories, although for MAD the averaged PA is $\langle \mathrmPA\rangle =219^\circ $, $\sigma _\mathrmPA=46^\circ $, while for SANE $\langle \mathrmPA\rangle =195^\circ $ and $\sigma _\mathrmPA=58^\circ $. The $a_* =0$ and $a_* \gt 0$ models have similar distributions. Again, EHT2017 data strongly favor one sense of black hole spin : either $| a_* | $ is small, or the spin vector is pointed away from Earth. If the fluctuations are such that the fitted PA for each epoch of observations is drawn from a Gaussian with $\sigma _\mathrmPA\simeq 55^\circ $, then a second epoch will be able to identify the true orientation with accuracy $\sigma _\mathrmPA/\sqrt2\simeq 40^\circ $ and the Nth epoch with accuracy $\sigma _\mathrmPA/\sqrtN$. If the fitted PA were drawn from a Gaussian of width $\sigma _\mathrmPA=54^\circ $ about $\mathrmPA=288^\circ $, as would be expected in a model in which the large-scale jet is aligned normal to the disk, then future epochs have a >90% chance of seeing the peak brightness counterclockwise from its position in EHT2017.

Finally, we can test the models by asking if they are consistent with the data according to Themis-AIS, as introduced in Section 4. Themis-AIS produces a probability p that the $\chi _\nu ^2$ distance between the data and the average of the model images is drawn from the same distribution as the $\chi _\nu ^2$ distance between synthetic data created from the model images, and the average of the model images. Table 1 takes these p values and categorizes them by magnetic flux and by spin, aggregating (averaging) results from different codes, $R_\mathrmhigh$, and i. Evidently, most of the models are formally consistent with the data by this test.

aThe Average Image Scoring (Themis-AIS) is introduced in Section 4. bflux : net magnetic flux on the black hole (MAD or SANE). c $a_* $ : dimensionless black hole spin. d $\langle p\rangle $ : mean of the p value for the aggregated models. e $N_\mathrmmodel$ : number of aggregated models. f $\mathrmMIN(p)$ : minimum p value among the aggregated models. g $\mathrmMAX(p)$ : maximum p value among the aggregated models.

One group of models, however, is rejected by Themis-AIS : MAD models with $a_* =-0.94$. On average this group has p = 0.01, and all models within this group have $p\leqslant 0.04$. Snapshots from MAD models with $a_* =-0.94$ exhibit the highest morphological variability in our ensemble in the sense that the emission breaks up into transient bright clumps. These models are rejected by Themis-AIS because none of the snapshots are as similar to the average image as the data. In other words, it is unlikely that EHT2017 would have captured an $a_* =-0.94$ MAD model in a configuration as unperturbed as the data seem to be.

The remainder of the model categories contain at least some models that are consistent with the data according to the average image scoring test. That is, most models are variable and the associated snapshots lie far from the average image. These snapshots are formally inconsistent with the data, but their distance from the average image is consistent with what is expected from the models. Given the uncertainties in the model—and our lack of knowledge of the source prior to EHT2017—it is remarkable that so many of the models are acceptable. This is likely because the source structure is dominated by the photon ring, which is produced by gravitational lensing, and is therefore relatively insensitive to the details of the accretion flow and jet physics. We can further narrow the range of acceptable models, however, using additional constraints.
6. Model Constraints : EHT2017 Combined with Other Constraints

We can apply three additional arguments to further constrain the source model. (1) The model must be close to radiative equilibrium. (2) The model must be consistent with the observed broadband SED ; in particular, it must not overproduce X-rays. (3) The model must produce a sufficiently powerful jet to match the measurements of the jet kinetic energy at large scales. Our discussions in this Section are based on simulation data that is provided in full detail in Appendix A.
6.1. Radiative Equilibrium

The model must be close to radiative equilibrium. The GRMHD models in the Simulation Library do not include radiative cooling, nor do they include a detailed prescription for particle energization. In nature the accretion flow and jet are expected to be cooled and heated by a combination of synchrotron and Compton cooling, turbulent dissipation, and Coulomb heating, which transfers energy from the hot ions to the cooler electrons. In our suite of simulations the parameter $R_\mathrmhigh$ can be thought of as a proxy for the sum of these processes. In a fully self-consistent treatment, some models would rapidly cool and settle to a lower electron temperature (see Mościbrodzka et al. 2011 ; Ryan et al. 2018 ; Chael et al. 2019). We crudely test for this by calculating the radiative efficiency $\epsilon \equiv L_\mathrmbol/(\dotMc^2)$, where $L_\mathrmbol$ is the bolometric luminosity. If it is larger than the radiative efficiency of a thin, radiatively efficient disk,117 which depends only on $a_* $ (Novikov & Thorne 1973), then we reject the model as physically inconsistent.

We calculate $L_\mathrmbol$ with the Monte Carlo code grmonty (Dolence et al. 2009), which incorporates synchrotron emission, absorption, Compton scattering at all orders, and bremsstrahlung. It assumes the same thermal eDF used in generating the Image Library. We calculate $L_\mathrmbol$ for 20% of the snapshots to minimize computational cost. We then average over snapshots to find $\langle L_\mathrmbol\rangle $. The mass accretion rate $\dotM$ is likewise computed for each snapshot and averaged over time. We reject models with epsilon that is larger than the classical thin disk model. (Table 3 in Appendix A lists epsilon for a large set of models.) All but two of the radiatively inconsistent models are MADs with $a_* \geqslant 0$ and $R_\mathrmhigh=1$. Eliminating all MAD models with $a_* \geqslant 0$ and $R_\mathrmhigh=1$ does not change any of our earlier conclusions.
6.2. X-Ray Constraints

As part of the EHT2017 campaign, we simultaneously observed M87 with the Chandra X-ray observatory and the Nuclear Spectroscopic Telescope Array (NuSTAR). The best fit to simultaneous Chandra and NuSTAR observations on 2017 April 12 and 14 implies a $2\mbox10\,\mathrmkeV$ luminosity of $L_

\rmX

_\mathrmobs

\,=4.4\pm 0.1\times 10^40\,\,\mathrmerg\,\rms

^-1$. We used the SEDs generated from the simulations while calculating $L_\mathrmbol$ to reject models that consistently overproduce X-rays ; specifically, we reject models with $\mathrmlogL_

\rmX

_\mathrmobs

\lt \mathrmlog\langle L_\rmX

\rangle -2\sigma (\mathrmlogL_\rmX

)$. We do not reject underluminous models because the X-rays could in principle be produced by direct synchrotron emission from nonthermal electrons or by other unresolved sources. Notice that $L_\rmX

$ is highly variable in all models so that the X-ray observations currently reject only a few models. Table 3 in Appendix A shows $\langle L_\rmX

\rangle $ as well as upper and lower limits for a set of models that is distributed uniformly across the parameter space.

In our models the X-ray flux is produced by inverse Compton scattering of synchrotron photons. The X-ray flux is an increasing function of $\tau _TT_e^2$ where τT is a characteristic Thomson optical depth ($\tau _T\sim 10^-5$), and the characteristic amplification factor for photon energies is $\propto T_e^2$ because the X-ray band is dominated by singly scattered photons interacting with relativistic electrons (we include all scattering orders in the Monte Carlo calculation). Increasing $R_\mathrmhigh$ at fixed $F_\nu (230\,\ \mathrmGHz)$ tends to increase $\dotM$ and therefore τT and decrease Te. The increase in Te dominates in our ensemble of models, and so models with small $R_\mathrmhigh$ have larger $L_\rmX

$, while models with large $R_\mathrmhigh$ have smaller $L_\rmX

$. The effect is not strictly monotonic, however, because of noise in our sampling process and the highly variable nature of the X-ray emission.

The overluminous models are mostly SANE models with $R_\mathrmhigh\leqslant 20$. The model with the highest $\langle L_\rmX

\rangle =4.2\,\times 10^42\,\,\mathrmerg\,\rms

^-1$ is a SANE, $a_* =0$, $R_\mathrmhigh=10$ model. The corresponding model with $R_\mathrmhigh=1$ has $\langle L_\rmX

\rangle =2.1\,\times 10^41\,\,\mathrmerg\,\rms

^-1$, and the difference between these two indicates the level of variability and the sensitivity of the average to the brightest snapshot. The upshot of application of the $L_\rmX

$ constraints is that $L_\rmX

$ is sensitive to $R_\mathrmhigh$. Very low values of $R_\mathrmhigh$ are disfavored. $L_\rmX

$ thus most directly constrains the electron temperature model.
6.3. Jet Power

Estimates of M87’s jet power ($P_\mathrmjet$) have been reviewed in Reynolds et al. (1996), Li et al. (2009), de Gasperin et al. (2012), Broderick et al. (2015), and Prieto et al. (2016). The estimates range from 1042 to $10^45\,\,\mathrmerg\,\rms

^-1$. This wide range is a consequence of both physical uncertainties in the models used to estimate $P_\mathrmjet$ and the wide range in length and timescales probed by the observations. Some estimates may sample a different epoch and thus provide little information on the state of the central engine during EHT2017. Nevertheless, observations of HST-1 yield $P_\mathrmjet\sim 10^44\ \,\mathrmerg\,\rms

^-1$ (e.g., Stawarz et al. 2006). HST-1 is within $\sim 70\,\mathrmpc$ of the central engine and, taking account of relativistic time foreshortening, may be sampling the central engine $P_\mathrmjet$ over the last few decades. Furthermore, the 1.3 mm light curve of M87 as observed by SMA shows $\lesssim 50 \% $ variability over decade timescales (Bower et al. 2015). Based on these considerations it seems reasonable to adopt a very conservative lower limit on jet power $\equiv P_\mathrmjet,\min =10^42\ \,\mathrmerg\,\rms

^-1$.

To apply this constraint we must define and measure $P_\mathrmjet$ in our models. Our procedure is discussed in detail in Appendix A. In brief, we measure the total energy flux in outflowing regions over the polar caps of the black hole in which the energy per unit rest mass exceeds $2.2\,c^2$, which corresponds to βγ = 1, where $\beta \equiv v/c$ and γ is Lorentz factor. The effect of changing this cutoff is also discussed in Appendix A. Because the cutoff is somewhat arbitrary, we also calculate $P_\mathrmout$ by including the energy flux in all outflowing regions over the polar caps of the black hole ; that is, it includes the energy flux in any wide-angle, low-velocity wind. $P_\mathrmout$ represents a maximal definition of jet power. Table 3 in Appendix A shows $P_\mathrmjet$ as well as a total outflow power $P_\mathrmout$.

The constraint $P_\mathrmjet\gt P_\mathrmjet,\min =10^42\,\mathrmerg\,\rms

^-1$ rejects all $a_* =0$ models. This conclusion is not sensitive to the definition of $P_\mathrmjet$ : all $a_* =0$ models also have total outflow power $P_\mathrmout\,\lt 10^42\,\mathrmerg\,\rms

^-1$. The most powerful $a_* =0$ model is a MAD model with $R_\mathrmhigh=160$, which has $P_\mathrmout=3.7\times 10^41\,\mathrmerg\,\rms

^-1$ and $P_\mathrmjet$ consistent with 0. We conclude that our $a_* =0$ models are ruled out.

Can the $a_* =0$ models be saved by changing the eDF ? Probably not. There is no evidence from the GRMHD simulations that these models are capable of producing a relativistic outflow with $\beta \gamma \gt 1$. Suppose, however, that we are willing to identify the nonrelativistic outflow, whose power is measured by $P_\mathrmout$, with the jet. Can $P_\mathrmout$ be raised to meet our conservative threshold on jet power ? Here the answer is yes, in principle, and this can be done by changing the eDF. The eDF and $P_\mathrmout$ are coupled because $P_\mathrmout$ is determined by $\dotM$, and $\dotM$ is adjusted to produce the observed compact mm flux. The relationship between $\dotM$ and mm flux depends upon the eDF. If the eDF is altered to produce mm photons less efficiently (for example, by lowering Te in a thermal model), then $\dotM$ and therefore $P_\mathrmout$ increase. A typical nonthermal eDF, by contrast, is likely to produce mm photons with greater efficiency by shifting electrons out of the thermal core and into a nonthermal tail. It will therefore lower $\dotM$ and thus $P_\mathrmout$. A thermal eDF with lower Te could have higher $P_\mathrmout$, as is evident in the large $R_\mathrmhigh$ SANE models in Table 3. There are observational and theoretical lower limits on Te, however, including a lower limit provided by the observed brightness temeprature. As Te declines, ne and B increase and that has implications for source linear polarization (Mościbrodzka et al. 2017 ; Jiménez-Rosales & Dexter 2018), which will be explored in future work. As Te declines and ne and ni increase there is also an increase in energy transfer from ions to electrons by Coulomb coupling, and this sets a floor on Te.

The requirement that $P_\mathrmjet\gt P_\mathrmjet,\min $ eliminates many models other than the $a_* =0$ models. All SANE models with $| a_* | =0.5$ fail to produce jets with the required minimum power. Indeed, they also fail the less restrictive condition $P_\mathrmout\gt P_\mathrmjet,\min $, so this conclusion is insensitive to the definition of the jet. We conclude that among the SANE models, only high-spin models survive.

At this point it is worth revisiting the SANE, $R_\mathrmhigh=1$, $a_* =-0.94$ model that favored a low black hole mass in Section 5. These models are not rejected by a naive application of the $P_\mathrmjet\gt P_\mathrmjet,\min $ criterion, but they are marginal. Notice, however, that we needed to assume a mass in applying the this criterion. We have consistently assumed $M=6.2\times 10^9\,M_\odot $. If we use the $M\sim 3\times 10^9\,M_\odot $ implied by the best-fit M/D, then $\dotM$ drops by a factor of two, $P_\mathrmjet$ drops below the threshold and the model is rejected.

The lower limit on jet power $P_\mathrmjet,\min =10^42\,\mathrmerg\,\rms

^-1$ is conservative and the true jet power is likely higher. If we increased $P_\mathrmjet,\min $ to $3\times 10^42\,\mathrmerg\,\rms

^-1$, the only surviving models would have $| a_* | =0.94$ and $R_\mathrmhigh\geqslant 10$. This conclusion is also not sensitive to the definition of the jet power : applying the same cut to $P_\mathrmout$ adds only a single model with $| a_* | \lt 0.94$, the $R_\mathrmhigh=160$, $a_* =0.5$ MAD model. The remainder have $a_* =0.94$. Interestingly, the most powerful jets in our ensemble of models are produced by SANE, $a_* =-0.94$, $R_\mathrmhigh=160$ models, with $P_\mathrmjet\simeq 10^43\,\,\mathrmerg\,\rms

^-1$.

Estimates for $P_\mathrmjet$ extend to $10^45\,\mathrmerg\,\rms

^-1$, but in our ensemble of models the maximum $P_\mathrmjet\sim 10^43\,\,\mathrmerg\,\rms

^-1$. Possible explanations include : (1) $P_\mathrmjet$ is variable and the estimates probe the central engine power at earlier epochs (discussed above) ; (2) the $P_\mathrmjet$ estimates are too large ; or (3) the models are in error. How might our models be modified to produce a larger $P_\mathrmjet$ ? For a given magnetic field configuration the jet power scales with $\dotMc^2$. To increase $P_\mathrmjet$, then, one must reduce the mm flux per accreted nucleon so that at fixed mm flux density $\dotM$ increases.118 Lowering Te in a thermal model is unlikely to work because lower Te implies higher synchrotron optical depth, which increases the ring width. We have done a limited series of experiments that suggest that even a modest decrease in Te would produce a broad ring that is inconsistent with EHT2017 (Paper VI). What is required, then, is a nonthermal (or multitemperature) model with a large population of cold electrons that are invisible at mm wavelength (for a thermal subpopulation, $\rm\Theta

_e,\mathrmcold\lt 1$), and a population of higher-energy electrons that produces the observed mm flux (see Falcke & Biermann 1995). We have not considered such models here, but we note that they are in tension with current ideas about dissipation of turbulence because they require efficient suppression of electron heating.

The $P_\mathrmjet$ in our models is dominated by Poynting flux in the force-free region around the axis (the ’’’funnel’’’), as in the Blandford & Znajek (1977) force-free magnetosphere model. The energy flux is concentrated along the walls of the funnel.119 Tchekhovskoy et al. (2011) provided an expression for the energy flux in the funnel, the so-called Blandford–Znajek power $P_\mathrmBZ$, which becomes, in our units,
Equation (9)

where $f(a_* )\approx a_* ^2\left(1+\sqrt1-a_* ^2\right)^-2$ (a good approximation for $a_* \lt 0.95$) and $\dotM_\mathrmEdd=137\,M_\odot \,\mathrmyr^-1$ for $M=6.2\times 10^9\,M_\odot $. This expression was developed for models with a thin disk in the equatorial plane. $P_\mathrmBZ$ is lower for models where the force-free region is excluded by a thicker disk around the equatorial plane. Clearly $P_\mathrmBZ$ is comparable to observational estimates of $P_\mathrmjet$.

In our models (see Table 3) $P_\mathrmjet$ follows the above scaling relation but with a smaller coefficient. The ratio of coefficients is model dependent and varies from 0.15 to 0.83. This is likely because the force-free region is restricted to a cone around the poles of the black hole, and the width of the cone varies by model. Indeed, the coefficient is larger for MAD than for SANE models, which is consistent with this idea because MAD models have a wide funnel and SANE models have a narrow funnel. This also suggests that future comparison of synthetic 43 and 86 GHz images from our models with lower-frequency VLBI data may further constrain the magnetic flux on the black hole.

The connection between the Poynting flux in the funnel and black hole spin has been discussed for some time in the simulation literature, beginning with McKinney & Gammie (2004 ; see also McKinney 2006 ; McKinney & Narayan 2007). The structure of the funnel magnetic field can be time-averaged and shown to match the analytic solution of Blandford & Znajek (1977). Furthermore, the energy flux density can be time-averaged and traced back to the event horizon. Is the energy contained in black hole spin sufficient to drive the observed jet over the jet lifetime ? The spindown timescale is $\tau =(M-M_\mathrmirr)c^2/P_\mathrmjet$, where $M_\mathrmirr\equiv M\left(\left(1+\sqrt1-a_* ^2\right)/2\right)^1/2$ is the irreducible mass of the black hole. For the $a_* =0.94$ MAD model with $R_\mathrmhigh=160$, $\tau =7.3\times 10^12\,\mathrmyr$, which is long compared to a Hubble time ($\sim 10^10$ yr). Indeed, the spindown time for all models is long compared to the Hubble time.

We conclude that for models that have sufficiently powerful jets and are consistent with EHT2017, $P_\mathrmjet$ is driven by extraction of black hole spin energy through the Blandford–Znajek process.
6.4. Constraint Summary

We have applied constraints from AIS, a radiative self-consistency constraint, a constraint on maximum X-ray luminosity, and a constraint on minimum jet power. Which models survive ? Here we consider only models for which we have calculated $L_\rmX

$ and $L_\mathrmbol$. Table 2 summarizes the results. Here we consider only i = 163° (for $a_* \geqslant 0$) and i = 17° (for $a_* \lt 0$). The first three columns give the model parameters. The next four columns show the result of application of each constraint : Themis-AIS (here broken out by individual model rather than groups of models), radiative efficiency ($\epsilon \lt \epsilon _\mathrmthin\mathrmdisk$), $L_\rmX

$, and $P_\mathrmjet$.

aflux : net magnetic flux on the black hole (MAD, SANE). b $a_* $ : dimensionless black hole spin. c $R_\mathrmhigh$ : electron temperature parameter. See Equation (8). dAverage Image Scoring (Themis-AIS), models are rejected if $\langle p\rangle \leqslant 0.01$. See Section 4 and Table 1. eepsilon : radiative efficiency, models are rejected if epsilon is larger than the corresponding thin disk efficiency. See Section 6.1. f $L_\rmX

$ : X-ray luminosity ; models are rejected if $\langle L_\rmX

\rangle 10^-2\sigma \gt 4.4\times 10^40\,\,\mathrmerg\,\rms

^-1$. See Section 6.2. g $P_\mathrmjet$ : jet power, models are rejected if $P_\mathrmjet\leqslant 10^42\,\,\mathrmerg\,\rms

^-1$. See Section 6.3.

The final column gives the logical AND of the previous four columns, and allows a model to pass only if it passes all tests. Evidently most of the SANE models fail, with the exception of some $a_* =-0.94$ models and a few $a_* =0.94$ models with large $R_\mathrmhigh$. A much larger fraction of the MAD models pass, although $a_* =0$ models all fail because of inadequate jet power. MAD models with small $R_\mathrmhigh$ also fail. It is the jet power constraint that rejects the largest number of models.
7. Discussion

We have interpreted the EHT2017 data using a limited library of models with attendant limitations. Many of the limitations stem from the GRMHD model, which treats the plasma as an ideal fluid governed by equations that encode conservation laws for particle number, momentum, and energy. The eDF, in particular, is described by a number density and temperature, rather than a full distribution function, and the electron temperature Te is assumed to be a function of the local ion temperature and plasma $\beta _\rmp

$. Furthermore, all models assume a Kerr black hole spacetime, but there are alternatives. Here we consider some of the model limitations and possible extensions, including to models beyond general relativity.
7.1. Radiative Effects

Post-processed GRMHD simulations that are consistent with EHT data and the flux density of 1.3 mm emission in M87 can yield unphysically large radiative efficiencies (see Section 6). This implies that the radiative cooling timescale is comparable to or less than the advection timescale. As a consequence, including radiative cooling in simulations may be necessary to recover self-consistent models (see Mościbrodzka et al. 2011 ; Dibi et al. 2012). In our models we use a single parameter, $R_\mathrmhigh$, to adjust Te and account for all effects that might influence the electron energy density. How good is this approximation ?

The importance of radiative cooling can be assessed using newly developed, state-of-the-art general relativistic radiation GRMHD (’’’radiation GRMHD’’’) codes. Sa̧dowski et al. (2013b ; see also Sa̧dowski et al. 2014, 2017 ; McKinney et al. 2014) applied the M1 closure (Levermore 1984), which treats the radiation as a relativistic fluid. Ryan et al. (2015) introduced a Monte Carlo radiation GRMHD method, allowing for full frequency-dependent radiation transport. Models for turbulent dissipation into the electrons and ions, as well as heating and cooling physics that sets the temperature ratio Ti/Te, have been added to GRMHD and radiative GRMHD codes and used in simulations of Sgr A* (Ressler et al. 2015, 2017 ; Chael et al. 2018) and M87 (Ryan et al. 2018 ; Chael et al. 2019). While the radiative cooling and Coulomb coupling physics in these simulations is well understood, the particle heating process, especially the relative heating rates of ions and electrons, remains uncertain.

Radiation GRMHD models are computationally expensive per run and do not have the same scaling freedom as the GRMHD models, so they need to be repeatedly re-run with different initial conditions until they produce the correct 1.3 mm flux density. It is therefore impractical to survey the parameter space using radiation GRMHD. It is possible, however, to check individual GRMHD models against existing radiation GRMHD models of M87 (Ryan et al. 2018 ; Chael et al. 2019).

The SANE radiation GRMHD models of Ryan et al. (2018) with $a_* =0.94$ and $M=6\times 10^9\,M_\odot $ can be compared to GRMHD SANE $a_* =0.94$ models at various values of $R_\mathrmhigh$. The radiative models have $\dotM/\dotM_\mathrmEdd=5.2\times 10^-6$ and $P_\mathrmjet\,=5.1\times 10^41\,\mathrmerg\,\rms

^-1$. The GRMHD models in this work have, for $1\leqslant R_\mathrmhigh\leqslant 160$, $0.36\times 10^-6\leqslant \dotM/\dotM_\mathrmEdd\leqslant 20\times 10^-6$, and $0.22\leqslant P_\mathrmjet/(10^41\,\,\mathrmerg\,\rms

^-1)\leqslant 12$ (Table 3). Evidently the mass accretion rates and jet powers in the GRMHD models span a wide range that depends on $R_\mathrmhigh$, but when we choose $R_\mathrmhigh$ = 10 − 20 they are similar to what is found in the radiative GRMHD model when using the turbulent electron heating model (Howes 2010).

Figure 1. Left panel : an EHT2017 image of M87 from Paper IV of this series (see their Figure 15). Middle panel : a simulated image based on a GRMHD model. Right panel : the model image convolved with a $20\,\mu \mathrmas$ FWHM Gaussian beam. Although the most evident features of the model and data are similar, fine features in the model are not resolved by EHT.

Download figure :
Standard image High-resolution image Export PowerPoint slide

Four features of the image in the left panel of Figure 1 play an important role in our analysis : (1) the ring-like geometry, (2) the peak brightness temperature, (3) the total flux density, and (4) the asymmetry of the ring. We now consider each in turn.

(1) The compact source shows a bright ring with a central dark area without significant extended components. This bears a remarkable similarity to the long-predicted structure for optically thin emission from a hot plasma surrounding a black hole (Falcke et al. 2000). The central hole surrounded by a bright ring arises because of strong gravitational lensing (e.g., Hilbert 1917 ; von Laue 1921 ; Bardeen 1973 ; Luminet 1979). The so-called ’’’photon ring’’’ corresponds to lines of sight that pass close to (unstable) photon orbits (see Teo 2003), linger near the photon orbit, and therefore have a long path length through the emitting plasma. These lines of sight will appear comparatively bright if the emitting plasma is optically thin. The central flux depression is the so-called black hole ’’’shadow’’’ (Falcke et al. 2000), and corresponds to lines of sight that terminate on the event horizon. The shadow could be seen in contrast to surrounding emission from the accretion flow or lensed counter-jet in M87 (Broderick & Loeb 2009).

The photon ring is nearly circular for all black hole spins and all inclinations of the black hole spin axis to the line of sight (e.g., Johannsen & Psaltis 2010). For an $a_* =0$ black hole of mass $M$ and distance D, the photon ring angular radius on the sky is
Equation (1)

where we have scaled to the most likely mass from Gebhardt et al. (2011) and a distance of $16.9\,\mathrmMpc$ (see also EHT Collaboration et al. 2019e, (hereafter Paper VI ; Blakeslee et al. 2009 ; Bird et al. 2010 ; Cantiello et al. 2018). The photon ring angular radius for other inclinations and values of $a_* $ differs by at most 13% from Equation (1), and most of this variation occurs at $1-| a_* | \ll 1$ (e.g., Takahashi 2004 ; Younsi et al. 2016). Evidently the angular radius of the observed photon ring is approximately $\sim 20\,\mu \mathrmas$ (Figure 1 and Paper IV), which is close to the prediction of the black hole model given in Equation (1).

(2) The observed peak brightness temperature of the ring in Figure 1 is $T_b,pk\sim 6\times 10^9\,\rmK$, which is consistent with past EHT mm-VLBI measurements at 230 GHz (Doeleman et al. 2012 ; Akiyama et al. 2015), and GMVA 3 mm-VLBI measurements of the core region (Kim et al. 2018). Expressed in electron rest-mass (me) units, $\rm\Theta

_b,pk\equiv k_\rmB

T_b,pk/(m_ec^2)\simeq 1$, where $k_\rmB

$ is Boltzmann’s constant. The true peak brightness temperature of the source is higher if the ring is unresolved by EHT, as is the case for the model image in the center panel of Figure 1.

The 1.3 mm emission from M87 shown in Figure 1 is expected to be generated by the synchrotron process (see Yuan & Narayan 2014, and references therein) and thus depends on the electron distribution function (eDF). If the emitting plasma has a thermal eDF, then it is characterized by an electron temperature $T_e\geqslant T_b$, or $\rm\Theta

_e\equiv k_\rmB

T_e/(m_ec^2)\gt 1$, because $\rm\Theta

_e\gt \rm\Theta

_b,pk$ if the ring is unresolved or optically thin.

Is the observed brightness temperature consistent with what one would expect from phenomenological models of the source ? Radiatively inefficient accretion flow models of M87 (Reynolds et al. 1996 ; Di Matteo et al. 2003) produce mm emission in a geometrically thick donut of plasma around the black hole. The emitting plasma is collisionless : Coulomb scattering is weak at these low densities and high temperatures. Therefore, the electron and ion temperatures need not be the same (e.g., Spitzer 1962). In radiatively inefficient accretion flow models, the ion temperature is slightly less than the ion virial temperature,
Equation (2)

where $r_\rmg

\equiv GM/c^2$ is the gravitational radius, r is the Boyer–Lindquist or Kerr–Schild radius, and mp is the proton mass. Most models have an electron temperature $T_e\lt T_i$ because of electron cooling and preferential heating of the ions by turbulent dissipation (e.g., Yuan & Narayan 2014 ; Mościbrodzka et al. 2016). If the emission arises at $\sim 5\,r_\rmg

$, then $\rm\Theta

_e\simeq 37(T_e/T_i)$, which is then consistent with the observed $\rm\Theta

_b,pk$ if the source is unresolved or optically thin.

(3) The total flux density in the image at $1.3\,\mathrmmm$ is $\simeq 0.5$ Jy. With a few assumptions we can use this to estimate the electron number density ne and magnetic field strength B in the source. We adopt a simple, spherical, one-zone model for the source with radius $r\simeq 5\,r_\rmg

$, pressure $n_ikT_i+n_ekT_e=\beta _\rmp

B^2/(8\pi )$ with $\beta _\rmp

\equiv p_\mathrmgas/p_\mathrmmag\sim 1$, $T_i\simeq 3T_e$, and temperature $\theta _e\simeq 10\theta _b,pk$, which is consistent with the discussion in (2) above. Setting ne = ni (i.e., assuming a fully ionized hydrogen plasma), the values of B and ne required to produce the observed flux density can be found by solving a nonlinear equation (assuming an average angle between the field and line of sight, 60°). The solution can be approximated as a power law :
Equation (3)

Equation (4)

assuming that $M=6.2\times 10^9\,M_\odot $ and $D=16.9\,\mathrmMpc$, and using the approximate thermal emissivity of Leung et al. (2011). Then the synchrotron optical depth at $1.3\,\mathrmmm$ is 0.2. One can now estimate an accretion rate from (3) using
Equation (5)

assuming spherical symmetry. The Eddington accretion rate is
Equation (6)

where $L_\mathrmEdd\equiv 4\pi GMcm_p/\sigma _T$ is the Eddington luminosity ($\sigma _T$ is the Thomson cross section). Setting the efficiency $\epsilon =0.1$ and $M=6.2\times 10^9\,M_\odot $, $\dotM_\mathrmEdd=137\,M_\odot \,\mathrmyr^-1$, and therefore $\dotM/\dotM_\mathrmEdd\,\sim 2.0\times 10^-5$.

This estimate is similar to but slightly larger than the upper limit inferred from the 230 GHz linear polarization properties of M87 (Kuo et al. 2014).

(4) The ring is brighter in the south than the north. This can be explained by a combination of motion in the source and Doppler beaming. As a simple example we consider a luminous, optically thin ring rotating with speed v and an angular momentum vector inclined at a viewing angle i > 0° to the line of sight. Then the approaching side of the ring is Doppler boosted, and the receding side is Doppler dimmed, producing a surface brightness contrast of order unity if v is relativistic. The approaching side of the large-scale jet in M87 is oriented west–northwest (position angle $\mathrmPA\approx 288^\circ ;$ in Paper VI this is called $\mathrmPA_\mathrmFJ$), or to the right and slightly up in the image. Walker et al. (2018) estimated that the angle between the approaching jet and the line of sight is 17°. If the emission is produced by a rotating ring with an angular momentum vector oriented along the jet axis, then the plasma in the south is approaching Earth and the plasma in the north is receding. This implies a clockwise circulation of the plasma in the source, as projected onto the plane of the sky. This sense of rotation is consistent with the sense of rotation in ionized gas at arcsecond scales (Harms et al. 1994 ; Walsh et al. 2013). Notice that the asymmetry of the ring is consistent with the asymmetry inferred from 43 GHz observations of the brightness ratio between the north and south sides of the jet and counter-jet (Walker et al. 2018).

All of these estimates present a picture of the source that is remarkably consistent with the expectations of the black hole model and with existing GRMHD models (e.g., Dexter et al. 2012 ; Mościbrodzka et al. 2016). They even suggest a sense of rotation of gas close to the black hole. A quantitative comparison with GRMHD models can reveal more.
3. Models

Consistent with the discussion in Section 2, we now adopt the working hypothesis that M87 contains a turbulent, magnetized accretion flow surrounding a Kerr black hole. To test this hypothesis quantitatively against the EHT2017 data we have generated a Simulation Library of 3D time-dependent ideal GRMHD models. To generate this computationally expensive library efficiently and with independent checks on the results, we used several different codes that evolved matching initial conditions using the equations of ideal GRMHD. The codes used include BHAC (Porth et al. 2017), H-AMR (Liska et al. 2018 ; K. Chatterjee et al. 2019, in preparation), iharm (Gammie et al. 2003), and KORAL (Sa̧dowski et al. 2013b, 2014). A comparison of these and other GRMHD codes can be found in O. Porth et al. 2019 (in preparation), which shows that the differences between integrations of a standard accretion model with different codes is smaller than the fluctuations in individual simulations.

From the Simulation Library we have generated a large Image Library of synthetic images. Snapshots of the GRMHD evolutions were produced using the general relativistic ray-tracing (GRRT) schemes ipole (Mościbrodzka & Gammie 2018), RAPTOR (Bronzwaer et al. 2018), or BHOSS (Z. Younsi et al. 2019b, in preparation). A comparison of these and other GRRT codes can be found in Gold et al. (2019), which shows that the differences between codes is small.

In the GRMHD models the bulk of the 1.3 mm emission is produced within $\lesssim 10\,r_\rmg

$ of the black hole, where the models can reach a statistically steady state. It is therefore possible to compute predictive radiative models for this compact component of the source without accurately representing the accretion flow at all radii.

We note that the current state-of-the-art models for M87 are radiation GRMHD models that include radiative feedback and electron-ion thermodynamics (Ryan et al. 2018 ; Chael et al. 2019). These models are too computationally expensive for a wide survey of parameter space, so that in this Letter we consider only nonradiative GRMHD models with a parameterized treatment of the electron thermodynamics.
3.1. Simulation Library

All GRMHD simulations are initialized with a weakly magnetized torus of plasma orbiting in the equatorial plane of the black hole (e.g., De Villiers et al. 2003 ; Gammie et al. 2003 ; McKinney & Blandford 2009 ; Porth et al. 2017). We do not consider tilted models, in which the accretion flow angular momentum is misaligned with the black hole spin. The limitations of this approach are discussed in Section 7.

The initial torus is driven to a turbulent state by instabilities, including the magnetorotational instability (see e.g., Balbus & Hawley 1991). In all cases the outcome contains a moderately magnetized midplane with orbital frequency comparable to the Keplerian orbital frequency, a corona with gas-to-magnetic-pressure ratio $\beta _\rmp

\equiv p_\mathrmgas/p_\mathrmmag\sim 1$, and a strongly magnetized region over both poles of the black hole with $B^2/\rho c^2\gg 1$. We refer to the strongly magnetized region as the funnel, and the boundary between the funnel and the corona as the funnel wall (De Villiers et al. 2005 ; Hawley & Krolik 2006). All models in the library are evolved from t = 0 to $t=10^4\,r_\rmg

c^-1$.

The simulation outcome depends on the initial magnetic field strength and geometry insofar as these affect the magnetic flux through the disk, as discussed below. Once the simulation is initiated the disk transitions to a turbulent state and loses memory of most of the details of the initial conditions. This relaxed turbulent state is found inside a characteristic radius that grows over the course of the simulation. To be confident that we are imaging only those regions that have relaxed, we draw snapshots for comparison with the data from $5\times 10^3\le t/r_\rmg

c^-1\le 10^4$.

GRMHD models have two key physical parameters. The first is the black hole spin $a_* $, $-1\lt a_* \lt 1$. The second parameter is the absolute magnetic flux $\rm\Phi

_\mathrmBH$ crossing one hemisphere of the event horizon (see Tchekhovskoy et al. 2011 ; O. Porth et al. 2019, in preparation for a definition). It is convenient to recast $\rm\Phi

_\mathrmBH$ in dimensionless form $\phi \equiv \rm\Phi

_\mathrmBH\left(\dotMr_\rmg

^2c\right)^-1/2$.110

The magnetic flux phgr is nonzero because magnetic field is advected into the event horizon by the accretion flow and sustained by currents in the surrounding plasma. At $\phi \gt \phi _\max \sim 15$,111 numerical simulations show that the accumulated magnetic flux erupts, pushes aside the accretion flow, and escapes (Tchekhovskoy et al. 2011 ; McKinney et al. 2012). Models with $\phi \sim 1$ are conventionally referred to as Standard and Normal Evolution (SANE ; Narayan et al. 2012 ; Sa̧dowski et al (2013a)) models ; models with $\phi \sim \phi _\max $ are conventionally referred to as Magnetically Arrested Disk (MAD ; Igumenshchev et al. 2003 ; Narayan et al. 2003) models.

The Simulation Library contains SANE models with $a_* =-0.94$, −0.5, 0, 0.5, 0.75, 0.88, 0.94, 0.97, and 0.98, and MAD models with $a_* =-0.94$, −0.5, 0, 0.5, 0.75, and 0.94. The Simulation Library occupies 23 TB of disk space and contains a total of 43 GRMHD simulations, with some repeated at multiple resolutions with multiple codes, with consistent results (O. Porth et al. 2019, in preparation).
3.2. Image Library Generation

To produce model images from the simulations for comparison with EHT observations we use GRRT to generate a large number of synthetic images and derived VLBI data products. To make the synthetic images we need to specify the following : (1) the magnetic field, velocity field, and density as a function of position and time ; (2) the emission and absorption coefficients as a function of position and time ; and (3) the inclination angle between the accretion flow angular momentum vector and the line of sight i, the position angle $\mathrmPA$, the black hole mass $M$, and the distance D to the observer. In the following we discuss each input in turn. The reader who is only interested in a high-level description of the Image Library may skip ahead to Section 3.3.

(1) GRMHD models provide the absolute velocity field of the plasma flow. Nonradiative GRMHD evolutions are invariant, however, under a rescaling of the density by a factor $\mathscrM$. In particular, they are invariant under $\rho \to \mathscrM\rho $, field strength $B\to \mathscrM

^1/2B$, and internal energy $u\to \mathscrMu$ (the Alfvén speed $B/\rho ^1/2$ and sound speed $\propto \sqrtu/\rho $ are invariant). That is, there is no intrinsic mass scale in a nonradiative model as long as the mass of the accretion flow is negligible in comparison to $M$.112 We use this freedom to adjust $\mathscrM$ so that the average image from a GRMHD model has a 1.3 mm flux density ≈0.5 Jy (see Paper IV). Once $\mathscrM$ is set, the density, internal energy, and magnetic field are fully specified.

The mass unit $\mathscrM$ determines $\dotM$. In our ensemble of models $\dotM$ ranges from $2\times 10^-7\dotM_\mathrmEdd$ to $4\times 10^-4\dotM_\mathrmEdd$. Accretion rates vary by model category. The mean accretion rate for MAD models is $\sim 10^-6\dotM_\mathrmEdd$. For SANE models with $a_* \gt 0$ it is $\sim 5\times 10^-5\dotM_\mathrmEdd ;$ and for $a_* \lt 0$ it is $\sim 2\times 10^-4\dotM_\mathrmEdd$.

(2) The observed radio spectral energy distributions (SEDs) and the polarization characteristics of the source make clear that the 1.3 mm emission is synchrotron radiation, as is typical for active galactic nuclei (AGNs). Synchrotron absorption and emission coefficients depend on the eDF. In what follows, we adopt a relativistic, thermal model for the eDF (a Maxwell-Jüttner distribution ; Jüttner 1911 ; Rezzolla & Zanotti 2013). We discuss the limitations of this approach in Section 7.

All of our models of M87 are in a sufficiently low-density, high-temperature regime that the plasma is collisionless (see Ryan et al. 2018, for a discussion of Coulomb coupling in M87). Therefore, Te likely does not equal the ion temperature Ti, which is provided by the simulations. We set Te using the GRMHD density ρ, internal energy density u, and plasma $\beta _\rmp

$ using a simple model :
Equation (7)

where we have assumed that the plasma is composed of hydrogen, the ions are nonrelativistic, and the electrons are relativistic. Here $R\equiv T_i/T_e$ and
Equation (8)

This prescription has one parameter, $R_\mathrmhigh$, and sets $T_e\simeq T_i$ in low $\beta _\rmp

$ regions and $T_e\simeq T_i/R_\mathrmhigh$ in the midplane of the disk. It is adapted from Mościbrodzka et al. (2016) and motivated by models for electron heating in a turbulent, collisionless plasma that preferentially heats the ions for $\beta _\rmp

\gtrsim 1$ (e.g., Howes 2010 ; Kawazura et al. 2018).

(3) We must specify the observer inclination i, the orientation of the observer through the position angle $\mathrmPA$, the black hole mass $M$, and the distance D to the source. Non-EHT constraints on i, $\mathrmPA$, and $M$ are considered below ; we have generated images at $i=12^\circ ,17^\circ ,22^\circ ,158^\circ ,163^\circ $, and 168° and a few at i = 148°. The position angle (PA) can be changed by simply rotating the image. All features of the models that we have examined, including $\dotM$, are insensitive to small changes in i. The image morphology does depend on whether i is greater than or less than 90°, as we will show below.

The model images are generated with a $160\times 160\,\mu \mathrmas$ field of view and $1\mu \mathrmas$ pixels, which are small compared to the $\sim 20\,\mu \mathrmas$ nominal resolution of EHT2017. Our analysis is insensitive to changes in the field of view and the pixel scale.

For $M$ we use the most likely value from the stellar absorption-line work, $6.2\times 10^9M_\odot $ (Gebhardt et al. 2011). For the distance D we use $16.9\,\,\mathrmMpc$, which is very close to that employed in Paper VI. The ratio $GM/(c^2D)=3.62\,\mu \mathrmas$ (hereafter M/D) determines the angular scale of the images. For some models we have also generated images with $M=3.5\times 10^9\,M_\odot $ to check that the analysis results are not predetermined by the input black hole mass.
3.3. Image Library Summary

The Image Library contains of order 60,000 images. We generate images from 100 to 500 distinct output files from each of the GRMHD models at each of $R_\mathrmhigh=1,10,20,40,80$, and 160. In comparing to the data we adjust the $\mathrmPA$ by rotation and the total flux and angular scale of the image by simply rescaling images from the standard parameters in the Image Library (see Figure 29 in Paper VI). Tests indicate that comparisons with the data are insensitive to the rescaling procedure unless the angular scaling factor or flux scaling factor is large.113

The comparisons with the data are also insensitive to image resolution.114

A representative set of time-averaged images from the Image Library are shown in Figures 2 and 3. From these figures it is clear that varying the parameters $a_* $, phgr, and $R_\mathrmhigh$ can change the width and asymmetry of the photon ring and introduce additional structures exterior and interior to the photon ring.
Figure 2.

Figure 2. Time-averaged 1.3 mm images generated by five SANE GRMHD simulations with varying spin ($a_* =-0.94$ to $a_* =+0.97$ from left to right) and $R_\mathrmhigh$ ($R_\mathrmhigh=1$ to $R_\mathrmhigh=160$ from top to bottom ; increasing $R_\mathrmhigh$ corresponds to decreasing electron temperature). The colormap is linear. All models are imaged at i = 163°. The jet that is approaching Earth is on the right (west) in all the images. The black hole spin vector projected onto the plane of the sky is marked with an arrow and aligned in the east–west direction. When the arrow is pointing left the black hole rotates in a clockwise direction, and when the arrow is pointing right the black hole rotates in a counterclockwise direction. The field of view for each model image is $80\,\mu \mathrmas$ (half of that used for the image libraries) with resolution equal to $1\mu \mathrmas$/pixel (20 times finer than the nominal resolution of EHT2017, and the same employed in the library images).

Download figure :
Standard image High-resolution image Export PowerPoint slide
Figure 3.

Figure 3. Same as in Figure 2 but for selected MAD models.

Download figure :
Standard image High-resolution image Export PowerPoint slide

The location of the emitting plasma is shown in Figure 4, which shows a map of time- and azimuth-averaged emission regions for four representative $a_* \gt 0$ models. For SANE models, if $R_\mathrmhigh$ is low (high), emission is concentrated more in the disk (funnel wall), and the bright section of the ring is dominated by the disk (funnel wall).115 Appendix B shows images generated by considering emission only from particular regions of the flow, and the results are consistent with Figure 4.
Figure 4.

Figure 4. Binned location of the point of origin for all photons that make up an image, summed over azimuth, and averaged over all snapshots from the simulation. The colormap is linear. The event horizon is indicated by the solid white semicircle and the black hole spin axis is along the figure vertical axis. This set of four images shows MAD and SANE models with $R_\mathrmhigh=10$ and 160, all with $a_* =0.94$. The region between the dashed curves is the locus of existence of (unstable) photon orbits (Teo 2003). The green cross marks the location of the innermost stable circular orbit (ISCO) in the equatorial plane. In these images the line of sight (marked by an arrow) is located below the midplane and makes a 163° angle with the disk angular momentum, which coincides with the spin axis of the black hole.

Download figure :
Standard image High-resolution image Export PowerPoint slide

Figures 2 and 3 show that for both MAD and SANE models the bright section of the ring, which is generated by Doppler beaming, shifts from the top for negative spin, to a nearly symmetric ring at $a_* =0$, to the bottom for $a_* \gt 0$ (except the SANE $R_\mathrmhigh=1$ case, where the bright section is always at the bottom when i > 90°). That is, the location of the peak flux in the ring is controlled by the black hole spin : it always lies roughly 90 degrees counterclockwise from the projection of the spin vector on the sky. Some of the ring emission originates in the funnel wall at $r\lesssim 8\,r_\rmg

$. The rotation of plasma in the funnel wall is in the same sense as plasma in the funnel, which is controlled by the dragging of magnetic field lines by the black hole. The funnel wall thus rotates opposite to the accretion flow if $a_* \lt 0$. This effect will be studied further in a later publication (Wong et al. 2019). The resulting relationships between disk angular momentum, black hole angular momentum, and observed ring asymmetry are illustrated in Figure 5.
Figure 5.

Figure 5. Illustration of the effect of black hole and disk angular momentum on ring asymmetry. The asymmetry is produced primarily by Doppler beaming : the bright region corresponds to the approaching side. In GRMHD models that fit the data comparatively well, the asymmetry arises in emission generated in the funnel wall. The sense of rotation of both the jet and funnel wall are controlled by the black hole spin. If the black hole spin axis is aligned with the large-scale jet, which points to the right, then the asymmetry implies that the black hole spin is pointing away from Earth (rotation of the black hole is clockwise as viewed from Earth). The blue ribbon arrow shows the sense of disk rotation, and the black ribbon arrow shows black hole spin. Inclination i is defined as the angle between the disk angular momentum vector and the line of sight.

Download figure :
Standard image High-resolution image Export PowerPoint slide

The time-averaged MAD images are almost independent of $R_\mathrmhigh$ and depend mainly on $a_* $. In MAD models much of the emission arises in regions with $\beta _\rmp

\sim 1$, where $R_\mathrmhigh$ has little influence over the electron temperature, so the insensitivity to $R_\mathrmhigh$ is natural (see Figure 4). In SANE models emission arises at $\beta _\rmp

\sim 10$, so the time-averaged SANE images, by contrast, depend strongly on $R_\mathrmhigh$. In low $R_\mathrmhigh$ SANE models, extended emission outside the photon ring, arising near the equatorial plane, is evident at $R_\mathrmhigh=1$. In large $R_\mathrmhigh$ SANE models the inner ring emission arises from the funnel wall, and once again the image looks like a thin ring (see Figure 4).

Figure 6 and the accompanying animation show the evolution of the images, visibility amplitudes, and closure phases over a $5000\,r_\rmg

c^-1\approx 5\,\mathrmyr$ interval in a single simulation for M87. It is evident from the animation that turbulence in the simulations produces large fluctuations in the images, which imply changes in visibility amplitudes and closure phases that are large compared to measurement errors. The fluctuations are central to our procedure for comparing models with the data, described briefly below and in detail in Paper VI.

Figure 6. Single frame from the accompanying animation. This shows the visibility amplitudes (top), closure phases plotted by Euclidean distance in 6D space (middle), and associated model images at full resolution (lower left) and convolved with the EHT2017 beam (lower right). Data from 2017 April 6 high-band are also shown in the top two plots. The video shows frames 1 through 100 and has a duration of 10 s.

(An animation of this figure is available.)

Download figure :
Video Standard image High-resolution image

The timescale between frames in the animation is $50\,r_\rmg

c^-1\simeq 18$ days, which is long compared to EHT2017 observing campaign. The images are highly correlated on timescales less than the innermost stable circular orbit (ISCO) orbital period, which for $a_* =0$ is $\simeq 15\ r_\rmg

c^-1\simeq 5$ days, i.e., comparable to the duration of the EHT2017 campaign. If drawn from one of our models, we would expect the EHT2017 data to look like a single snapshot (Figures 6) rather than their time averages (Figures 2 and 3).
4. Procedure for Comparison of Models with Data

As described above, each model in the Simulation Library has two dimensionless parameters : black hole spin $a_* $ and magnetic flux phgr. Imaging the model from each simulation adds five new parameters : $R_\mathrmhigh$, i, $\mathrmPA$, $M$, and D, which we set to $16.9\,\mathrmMpc$. After fixing these parameters we draw snapshots from the time evolution at a cadence of 10 to $50\,r_\rmg

c^-1$. We then compare these snapshots to the data.

The simplest comparison computes the $\chi _\nu ^2$ (reduced chi square) distance between the data and a snapshot. In the course of computing $\chi _\nu ^2$ we vary the image scale M/D, flux density Fν, position angle $\mathrmPA$, and the gain at each VLBI station in order to give each image every opportunity to fit the data. The best-fit parameters $(M/D,F_\nu ,\mathrmPA)$ for each snapshot are found by two pipelines independently : the Themis pipeline using a Markov chain Monte Carlo method (A. E. Broderick et al. 2019a, in preparation), and the GENA pipeline using an evolutionary algorithm for multidimensional minimization (Fromm et al. 2019a ; C. Fromm et al. 2019b, in preparation ; see also Section 4 of Paper VI for details). The best-fit parameters contain information about the source and we use the distribution of best-fit parameters to test the model by asking whether or not they are consistent with existing measurements of M/D and estimates of the jet $\mathrmPA$ on larger scales.

The $\chi _\nu ^2$ comparison alone does not provide a sharp test of the models. Fluctuations in the underlying GRMHD model, combined with the high signal-to-noise ratio for EHT2017 data, imply that individual snapshots are highly unlikely to provide a formally acceptable fit with $\chi _\nu ^2\simeq 1$. This is borne out in practice with the minimum $\chi _\nu ^2=1.79$ over the entire set of the more than 60,000 individual images in the Image Library. Nevertheless, it is possible to test if the $\chi _\nu ^2$ from the fit to the data is consistent with the underlying model, using ’’’Average Image Scoring’’’ with Themis (Themis-AIS), as described in detail in Appendix F of Paper VI). Themis-AIS measures a $\chi _\nu ^2$ distance (on the space of visibility amplitudes and closure phases) between a trial image and the data. In practice we use the average of the images from a given model as the trial image (hence Themis-AIS), but other choices are possible. We compute the $\chi _\nu ^2$ distance between the trial image and synthetic data produced from each snapshot. The model can then be tested by asking whether the data’s $\chi _\nu ^2$ is likely to have been drawn from the model’s distribution of $\chi _\nu ^2$. In particular, we can assign a probability p that the data is drawn from a specific model’s distribution.

In this Letter we focus on comparisons with a single data set, the 2017 April 6 high-band data (Paper III). The eight EHT2017 data sets, spanning four days with two bands on each day, are highly correlated. Assessing what correlation is expected in the models is a complicated task that we defer to later publications. The 2017 April 6 data set has the largest number of scans, 284 detections in 25 scans (see Paper III) and is therefore expected to be the most constraining.116
5. Model Constraints : EHT2017 Alone

The resolved ring-like structure obtained from the EHT2017 data provides an estimate of M/D (discussed in detail in Paper VI) and the jet $\mathrmPA$ from the immediate environment of the central black hole. As a first test of the models we can ask whether or not these are consistent with what is known from other mass measurements and from the orientation of the large-scale jet.

Figure 7 shows the distributions of best-fit values of M/D for a subset of the models for which spectra and jet power estimates are available (see below). The three lines show the M/D distribution for all snapshots (dotted lines), the best-fit 10% of snapshots (dashed lines), and the best-fit 1% of snapshots (solid lines) within each model. Evidently, as better fits are required, the distribution narrows and peaks close to $M/D\sim 3.6\,\mu \mathrmas$ with a width of about $0.5\mu \mathrmas$.
Figure 7.

Figure 7. Distribution of M/D obtained by fitting Image Library snapshots to the 2017 April 6 data, in $\mu \mathrmas$, measured independently using the (left panel) Themis and (right panel) GENA pipelines with qualitatively similar results. Smooth lines were drawn with a Gaussian kernel density estimator. The three lines show the best-fit 1% within each model (solid) ; the best-fit 10% within each model (dashed) ; and all model images (dotted). The vertical lines show $M/D=2.04$ (dashed) and $3.62\,\mu \mathrmas$ (solid), corresponding to M = 3.5 and $6.2\times 10^9\,M_\odot $. The distribution uses a subset of models for which spectra and jet power estimates are available (see Section 6). Only images with $a_* \gt 0$, i > 90° and $a_* \lt 0$, i < 90° (see also the left panel of Figure 5) are considered.

Download figure :
Standard image High-resolution image Export PowerPoint slide

The distribution of M/D for the best-fit $\lt 10 \% $ of snapshots is qualitatively similar if we include only MAD or SANE models, only models produced by individual codes (BHAC, H-AMR, iharm, or KORAL), or only individual spins. As the thrust of this Letter is to test the models, we simply note that Figure 7 indicates that the models are broadly consistent with earlier mass estimates (see Paper VI for a detailed discussion). This did not have to be the case : the ring radius could have been significantly larger than $3.6\,\mu \mathrmas$.

We can go somewhat further and ask if any of the individual models favor large or small masses. Figure 8 shows the distributions of best-fit values of M/D for each model (different $a_* $, $R_\mathrmhigh$, and magnetic flux). Most individual models favor M/D close to $3.6\,\mu \mathrmas$. The exceptions are $a_* \leqslant 0$ SANE models with $R_\mathrmhigh=1$, which produce the bump in the M/D distribution near $2\mu \mathrmas$. In these models, the emission is produced at comparatively large radius in the disk (see Figure 2) because the inner edge of the disk (the ISCO) is at a large radius in a counter-rotating disk around a black hole with $| a_* | \sim 1$. For these models, the fitting procedure identifies EHT2017’s ring with this outer ring, which forces the photon ring, and therefore M/D, to be small. As we will show later, these models can be rejected because they produce weak jets that are inconsistent with existing jet power estimates (see Section 6.3).
Figure 8.

Figure 8. Distributions of M/D and black hole mass with $D=16.9\,\mathrmMpc$ reconstructed from the best-fit 10% of images for MAD (left panel) and SANE (right panel) models (i = 17° for $a_* \le 0$ and 163° for $a_* \gt 0$) with different $R_\mathrmhigh$ and $a_* $, from the Themis (dark red, left), and GENA (dark green, right) pipelines. The white dot and vertical black bar correspond, respectively, to the median and region between the 25th and 75th percentiles for both pipelines combined. The blue and pink horizontal bands show the range of M/D and mass at $D=16.9\,\mathrmMpc$ estimated from the gas dynamical model (Walsh et al. 2013) and stellar dynamical model (Gebhardt et al. 2011), respectively. Constraints on the models based on average image scoring (Themis-AIS) are discussed in Section 5. Constraints based on radiative efficiency, X-ray luminosity, and jet power are discussed in Section 6.

Download figure :
Standard image High-resolution image Export PowerPoint slide

Figure 8 also shows that M/D increases with $a_* $ for SANE models. This is due to the appearance of a secondary inner ring inside the main photon ring. The former is associated with emission produced along the wall of the approaching jet. Because the emission is produced in front of the black hole, lensing is weak and it appears at small angular scale. The inner ring is absent in MAD models (see Figure 3), where the bulk of the emission comes from the midplane at all values of $R_\mathrmhigh$(Figure 4).

We now ask whether or not the PA of the jet is consistent with the orientation of the jet measured at other wavelengths. On large ( mas) scales the extended jet component has a PA of approximately 288° (e.g., Walker et al. 2018). On smaller ($\sim 100\,\mu \mathrmas$) scales the apparent opening angle of the jet is large (e.g., Kim et al. 2018) and the PA is therefore more difficult to measure. Also notice that the jet PA may be time dependent (e.g., Hada et al. 2016 ; Walker et al. 2018). In our model images the jet is relatively dim at 1.3 mm, and is not easily seen with a linear colormap. The model jet axis is, nonetheless, well defined : jets emerge perpendicular to the disk.

Figure 9 shows the distribution of best-fit PA over the same sample of snapshots from the Image Library used in Figure 7. We divide the snapshots into two groups. The first group has the black hole spin pointed away from Earth (i > 90° and $a_* \gt 0$, or i < 90° and $a_* \lt 0$). The spin-away model PA distributions are shown in the top two panels. The second group has the black hole spin pointed toward Earth (i > 90 and $a_* \lt 0$ or i > 90° and $a_* \lt 0$). These spin-toward model PA distributions are shown in the bottom two panels. The large-scale jet orientation lies on the shoulder of the spin-away distribution (the distribution can be approximated as a Gaussian with, for Themis (GENA) mean 209 (203)° and $\sigma _\mathrmPA\,=54\,(55)^\circ ;$ the large-scale jet PA lies $1.5\sigma _\mathrmPA$ from the mean) and is therefore consistent with the spin-away models. On the other hand, the large-scale jet orientation lies off the shoulder of the spin-toward distribution and is inconsistent with the spin-toward models. Evidently models in which the black hole spin is pointing away from Earth are strongly favored.
Figure 9.

Figure 9. Top : distribution of best-fit PA (in degree) scored by the Themis (left) and GENA (right) pipelines for models with black hole spin vector pointing away from Earth (i > 90° for $a_* \gt 0$ or i < 90° for $a_* \lt 0$). Bottom : images with black hole spin vector pointing toward Earth (i < 90° for $a_* \gt 0$ or i > 90° for $a_* \lt 0$). Smooth lines were drawn with a wrapped Gaussian kernel density estimator. The three lines show (1) all images in the sample (dotted line) ; (2) the best-fit 10% of images within each model (dashed line) ; and (3) the best-fit 1% of images in each model (solid line). For reference, the vertical line shows the position angle $\mathrmPA\sim 288^\circ $ of the large-scale (mas) jet Walker et al. (2018), with the gray area from (288 – 10)° to (288 + 10)° indicating the observed PA variation.

Download figure :
Standard image High-resolution image Export PowerPoint slide

The width of the spin-away and spin-toward distributions arises naturally in the models from brightness fluctuations in the ring. The distributions are relatively insensitive if split into MAD and SANE categories, although for MAD the averaged PA is $\langle \mathrmPA\rangle =219^\circ $, $\sigma _\mathrmPA=46^\circ $, while for SANE $\langle \mathrmPA\rangle =195^\circ $ and $\sigma _\mathrmPA=58^\circ $. The $a_* =0$ and $a_* \gt 0$ models have similar distributions. Again, EHT2017 data strongly favor one sense of black hole spin : either $| a_* | $ is small, or the spin vector is pointed away from Earth. If the fluctuations are such that the fitted PA for each epoch of observations is drawn from a Gaussian with $\sigma _\mathrmPA\simeq 55^\circ $, then a second epoch will be able to identify the true orientation with accuracy $\sigma _\mathrmPA/\sqrt2\simeq 40^\circ $ and the Nth epoch with accuracy $\sigma _\mathrmPA/\sqrtN$. If the fitted PA were drawn from a Gaussian of width $\sigma _\mathrmPA=54^\circ $ about $\mathrmPA=288^\circ $, as would be expected in a model in which the large-scale jet is aligned normal to the disk, then future epochs have a >90% chance of seeing the peak brightness counterclockwise from its position in EHT2017.

Finally, we can test the models by asking if they are consistent with the data according to Themis-AIS, as introduced in Section 4. Themis-AIS produces a probability p that the $\chi _\nu ^2$ distance between the data and the average of the model images is drawn from the same distribution as the $\chi _\nu ^2$ distance between synthetic data created from the model images, and the average of the model images. Table 1 takes these p values and categorizes them by magnetic flux and by spin, aggregating (averaging) results from different codes, $R_\mathrmhigh$, and i. Evidently, most of the models are formally consistent with the data by this test.

Table 1. Average Image Scoringa Summary
Fluxb $a_* $ c $\langle p\rangle $ d $N_\mathrmmodel$ e $\mathrmMIN(p)$ f $\mathrmMAX(p)$ g
SANE −0.94 0.33 24 0.01 0.88
SANE −0.5 0.19 24 0.01 0.73
SANE 0 0.23 24 0.01 0.92
SANE 0.5 0.51 30 0.02 0.97
SANE 0.75 0.74 6 0.48 0.98
SANE 0.88 0.65 6 0.26 0.94
SANE 0.94 0.49 24 0.01 0.92
SANE 0.97 0.12 6 0.06 0.40
MAD −0.94 0.01 18 0.01 0.04
MAD −0.5 0.75 18 0.34 0.98
MAD 0 0.22 18 0.01 0.62
MAD 0.5 0.17 18 0.02 0.54
MAD 0.75 0.28 18 0.01 0.72
MAD 0.94 0.21 18 0.02 0.50

Notes.
aThe Average Image Scoring (Themis-AIS) is introduced in Section 4. bflux : net magnetic flux on the black hole (MAD or SANE). c $a_* $ : dimensionless black hole spin. d $\langle p\rangle $ : mean of the p value for the aggregated models. e $N_\mathrmmodel$ : number of aggregated models. f $\mathrmMIN(p)$ : minimum p value among the aggregated models. g $\mathrmMAX(p)$ : maximum p value among the aggregated models.

Download table as : ASCIITypeset image

One group of models, however, is rejected by Themis-AIS : MAD models with $a_* =-0.94$. On average this group has p = 0.01, and all models within this group have $p\leqslant 0.04$. Snapshots from MAD models with $a_* =-0.94$ exhibit the highest morphological variability in our ensemble in the sense that the emission breaks up into transient bright clumps. These models are rejected by Themis-AIS because none of the snapshots are as similar to the average image as the data. In other words, it is unlikely that EHT2017 would have captured an $a_* =-0.94$ MAD model in a configuration as unperturbed as the data seem to be.

The remainder of the model categories contain at least some models that are consistent with the data according to the average image scoring test. That is, most models are variable and the associated snapshots lie far from the average image. These snapshots are formally inconsistent with the data, but their distance from the average image is consistent with what is expected from the models. Given the uncertainties in the model—and our lack of knowledge of the source prior to EHT2017—it is remarkable that so many of the models are acceptable. This is likely because the source structure is dominated by the photon ring, which is produced by gravitational lensing, and is therefore relatively insensitive to the details of the accretion flow and jet physics. We can further narrow the range of acceptable models, however, using additional constraints.
6. Model Constraints : EHT2017 Combined with Other Constraints

We can apply three additional arguments to further constrain the source model. (1) The model must be close to radiative equilibrium. (2) The model must be consistent with the observed broadband SED ; in particular, it must not overproduce X-rays. (3) The model must produce a sufficiently powerful jet to match the measurements of the jet kinetic energy at large scales. Our discussions in this Section are based on simulation data that is provided in full detail in Appendix A.
6.1. Radiative Equilibrium

The model must be close to radiative equilibrium. The GRMHD models in the Simulation Library do not include radiative cooling, nor do they include a detailed prescription for particle energization. In nature the accretion flow and jet are expected to be cooled and heated by a combination of synchrotron and Compton cooling, turbulent dissipation, and Coulomb heating, which transfers energy from the hot ions to the cooler electrons. In our suite of simulations the parameter $R_\mathrmhigh$ can be thought of as a proxy for the sum of these processes. In a fully self-consistent treatment, some models would rapidly cool and settle to a lower electron temperature (see Mościbrodzka et al. 2011 ; Ryan et al. 2018 ; Chael et al. 2019). We crudely test for this by calculating the radiative efficiency $\epsilon \equiv L_\mathrmbol/(\dotMc^2)$, where $L_\mathrmbol$ is the bolometric luminosity. If it is larger than the radiative efficiency of a thin, radiatively efficient disk,117 which depends only on $a_* $ (Novikov & Thorne 1973), then we reject the model as physically inconsistent.

We calculate $L_\mathrmbol$ with the Monte Carlo code grmonty (Dolence et al. 2009), which incorporates synchrotron emission, absorption, Compton scattering at all orders, and bremsstrahlung. It assumes the same thermal eDF used in generating the Image Library. We calculate $L_\mathrmbol$ for 20% of the snapshots to minimize computational cost. We then average over snapshots to find $\langle L_\mathrmbol\rangle $. The mass accretion rate $\dotM$ is likewise computed for each snapshot and averaged over time. We reject models with epsilon that is larger than the classical thin disk model. (Table 3 in Appendix A lists epsilon for a large set of models.) All but two of the radiatively inconsistent models are MADs with $a_* \geqslant 0$ and $R_\mathrmhigh=1$. Eliminating all MAD models with $a_* \geqslant 0$ and $R_\mathrmhigh=1$ does not change any of our earlier conclusions.
6.2. X-Ray Constraints

As part of the EHT2017 campaign, we simultaneously observed M87 with the Chandra X-ray observatory and the Nuclear Spectroscopic Telescope Array (NuSTAR). The best fit to simultaneous Chandra and NuSTAR observations on 2017 April 12 and 14 implies a $2\mbox10\,\mathrmkeV$ luminosity of $L_

\rmX

_\mathrmobs

\,=4.4\pm 0.1\times 10^40\,\,\mathrmerg\,\rms

^-1$. We used the SEDs generated from the simulations while calculating $L_\mathrmbol$ to reject models that consistently overproduce X-rays ; specifically, we reject models with $\mathrmlogL_

\rmX

_\mathrmobs

\lt \mathrmlog\langle L_\rmX

\rangle -2\sigma (\mathrmlogL_\rmX

)$. We do not reject underluminous models because the X-rays could in principle be produced by direct synchrotron emission from nonthermal electrons or by other unresolved sources. Notice that $L_\rmX

$ is highly variable in all models so that the X-ray observations currently reject only a few models. Table 3 in Appendix A shows $\langle L_\rmX

\rangle $ as well as upper and lower limits for a set of models that is distributed uniformly across the parameter space.

In our models the X-ray flux is produced by inverse Compton scattering of synchrotron photons. The X-ray flux is an increasing function of $\tau _TT_e^2$ where τT is a characteristic Thomson optical depth ($\tau _T\sim 10^-5$), and the characteristic amplification factor for photon energies is $\propto T_e^2$ because the X-ray band is dominated by singly scattered photons interacting with relativistic electrons (we include all scattering orders in the Monte Carlo calculation). Increasing $R_\mathrmhigh$ at fixed $F_\nu (230\,\ \mathrmGHz)$ tends to increase $\dotM$ and therefore τT and decrease Te. The increase in Te dominates in our ensemble of models, and so models with small $R_\mathrmhigh$ have larger $L_\rmX

$, while models with large $R_\mathrmhigh$ have smaller $L_\rmX

$. The effect is not strictly monotonic, however, because of noise in our sampling process and the highly variable nature of the X-ray emission.

The overluminous models are mostly SANE models with $R_\mathrmhigh\leqslant 20$. The model with the highest $\langle L_\rmX

\rangle =4.2\,\times 10^42\,\,\mathrmerg\,\rms

^-1$ is a SANE, $a_* =0$, $R_\mathrmhigh=10$ model. The corresponding model with $R_\mathrmhigh=1$ has $\langle L_\rmX

\rangle =2.1\,\times 10^41\,\,\mathrmerg\,\rms

^-1$, and the difference between these two indicates the level of variability and the sensitivity of the average to the brightest snapshot. The upshot of application of the $L_\rmX

$ constraints is that $L_\rmX

$ is sensitive to $R_\mathrmhigh$. Very low values of $R_\mathrmhigh$ are disfavored. $L_\rmX

$ thus most directly constrains the electron temperature model.
6.3. Jet Power

Estimates of M87’s jet power ($P_\mathrmjet$) have been reviewed in Reynolds et al. (1996), Li et al. (2009), de Gasperin et al. (2012), Broderick et al. (2015), and Prieto et al. (2016). The estimates range from 1042 to $10^45\,\,\mathrmerg\,\rms

^-1$. This wide range is a consequence of both physical uncertainties in the models used to estimate $P_\mathrmjet$ and the wide range in length and timescales probed by the observations. Some estimates may sample a different epoch and thus provide little information on the state of the central engine during EHT2017. Nevertheless, observations of HST-1 yield $P_\mathrmjet\sim 10^44\ \,\mathrmerg\,\rms

^-1$ (e.g., Stawarz et al. 2006). HST-1 is within $\sim 70\,\mathrmpc$ of the central engine and, taking account of relativistic time foreshortening, may be sampling the central engine $P_\mathrmjet$ over the last few decades. Furthermore, the 1.3 mm light curve of M87 as observed by SMA shows $\lesssim 50 \% $ variability over decade timescales (Bower et al. 2015). Based on these considerations it seems reasonable to adopt a very conservative lower limit on jet power $\equiv P_\mathrmjet,\min =10^42\ \,\mathrmerg\,\rms

^-1$.

To apply this constraint we must define and measure $P_\mathrmjet$ in our models. Our procedure is discussed in detail in Appendix A. In brief, we measure the total energy flux in outflowing regions over the polar caps of the black hole in which the energy per unit rest mass exceeds $2.2\,c^2$, which corresponds to βγ = 1, where $\beta \equiv v/c$ and γ is Lorentz factor. The effect of changing this cutoff is also discussed in Appendix A. Because the cutoff is somewhat arbitrary, we also calculate $P_\mathrmout$ by including the energy flux in all outflowing regions over the polar caps of the black hole ; that is, it includes the energy flux in any wide-angle, low-velocity wind. $P_\mathrmout$ represents a maximal definition of jet power. Table 3 in Appendix A shows $P_\mathrmjet$ as well as a total outflow power $P_\mathrmout$.

The constraint $P_\mathrmjet\gt P_\mathrmjet,\min =10^42\,\mathrmerg\,\rms

^-1$ rejects all $a_* =0$ models. This conclusion is not sensitive to the definition of $P_\mathrmjet$ : all $a_* =0$ models also have total outflow power $P_\mathrmout\,\lt 10^42\,\mathrmerg\,\rms

^-1$. The most powerful $a_* =0$ model is a MAD model with $R_\mathrmhigh=160$, which has $P_\mathrmout=3.7\times 10^41\,\mathrmerg\,\rms

^-1$ and $P_\mathrmjet$ consistent with 0. We conclude that our $a_* =0$ models are ruled out.

Can the $a_* =0$ models be saved by changing the eDF ? Probably not. There is no evidence from the GRMHD simulations that these models are capable of producing a relativistic outflow with $\beta \gamma \gt 1$. Suppose, however, that we are willing to identify the nonrelativistic outflow, whose power is measured by $P_\mathrmout$, with the jet. Can $P_\mathrmout$ be raised to meet our conservative threshold on jet power ? Here the answer is yes, in principle, and this can be done by changing the eDF. The eDF and $P_\mathrmout$ are coupled because $P_\mathrmout$ is determined by $\dotM$, and $\dotM$ is adjusted to produce the observed compact mm flux. The relationship between $\dotM$ and mm flux depends upon the eDF. If the eDF is altered to produce mm photons less efficiently (for example, by lowering Te in a thermal model), then $\dotM$ and therefore $P_\mathrmout$ increase. A typical nonthermal eDF, by contrast, is likely to produce mm photons with greater efficiency by shifting electrons out of the thermal core and into a nonthermal tail. It will therefore lower $\dotM$ and thus $P_\mathrmout$. A thermal eDF with lower Te could have higher $P_\mathrmout$, as is evident in the large $R_\mathrmhigh$ SANE models in Table 3. There are observational and theoretical lower limits on Te, however, including a lower limit provided by the observed brightness temeprature. As Te declines, ne and B increase and that has implications for source linear polarization (Mościbrodzka et al. 2017 ; Jiménez-Rosales & Dexter 2018), which will be explored in future work. As Te declines and ne and ni increase there is also an increase in energy transfer from ions to electrons by Coulomb coupling, and this sets a floor on Te.

The requirement that $P_\mathrmjet\gt P_\mathrmjet,\min $ eliminates many models other than the $a_* =0$ models. All SANE models with $| a_* | =0.5$ fail to produce jets with the required minimum power. Indeed, they also fail the less restrictive condition $P_\mathrmout\gt P_\mathrmjet,\min $, so this conclusion is insensitive to the definition of the jet. We conclude that among the SANE models, only high-spin models survive.

At this point it is worth revisiting the SANE, $R_\mathrmhigh=1$, $a_* =-0.94$ model that favored a low black hole mass in Section 5. These models are not rejected by a naive application of the $P_\mathrmjet\gt P_\mathrmjet,\min $ criterion, but they are marginal. Notice, however, that we needed to assume a mass in applying the this criterion. We have consistently assumed $M=6.2\times 10^9\,M_\odot $. If we use the $M\sim 3\times 10^9\,M_\odot $ implied by the best-fit M/D, then $\dotM$ drops by a factor of two, $P_\mathrmjet$ drops below the threshold and the model is rejected.

The lower limit on jet power $P_\mathrmjet,\min =10^42\,\mathrmerg\,\rms

^-1$ is conservative and the true jet power is likely higher. If we increased $P_\mathrmjet,\min $ to $3\times 10^42\,\mathrmerg\,\rms

^-1$, the only surviving models would have $| a_* | =0.94$ and $R_\mathrmhigh\geqslant 10$. This conclusion is also not sensitive to the definition of the jet power : applying the same cut to $P_\mathrmout$ adds only a single model with $| a_* | \lt 0.94$, the $R_\mathrmhigh=160$, $a_* =0.5$ MAD model. The remainder have $a_* =0.94$. Interestingly, the most powerful jets in our ensemble of models are produced by SANE, $a_* =-0.94$, $R_\mathrmhigh=160$ models, with $P_\mathrmjet\simeq 10^43\,\,\mathrmerg\,\rms

^-1$.

Estimates for $P_\mathrmjet$ extend to $10^45\,\mathrmerg\,\rms

^-1$, but in our ensemble of models the maximum $P_\mathrmjet\sim 10^43\,\,\mathrmerg\,\rms

^-1$. Possible explanations include : (1) $P_\mathrmjet$ is variable and the estimates probe the central engine power at earlier epochs (discussed above) ; (2) the $P_\mathrmjet$ estimates are too large ; or (3) the models are in error. How might our models be modified to produce a larger $P_\mathrmjet$ ? For a given magnetic field configuration the jet power scales with $\dotMc^2$. To increase $P_\mathrmjet$, then, one must reduce the mm flux per accreted nucleon so that at fixed mm flux density $\dotM$ increases.118 Lowering Te in a thermal model is unlikely to work because lower Te implies higher synchrotron optical depth, which increases the ring width. We have done a limited series of experiments that suggest that even a modest decrease in Te would produce a broad ring that is inconsistent with EHT2017 (Paper VI). What is required, then, is a nonthermal (or multitemperature) model with a large population of cold electrons that are invisible at mm wavelength (for a thermal subpopulation, $\rm\Theta

_e,\mathrmcold\lt 1$), and a population of higher-energy electrons that produces the observed mm flux (see Falcke & Biermann 1995). We have not considered such models here, but we note that they are in tension with current ideas about dissipation of turbulence because they require efficient suppression of electron heating.

The $P_\mathrmjet$ in our models is dominated by Poynting flux in the force-free region around the axis (the ’’’funnel’’’), as in the Blandford & Znajek (1977) force-free magnetosphere model. The energy flux is concentrated along the walls of the funnel.119 Tchekhovskoy et al. (2011) provided an expression for the energy flux in the funnel, the so-called Blandford–Znajek power $P_\mathrmBZ$, which becomes, in our units,
Equation (9)

where $f(a_* )\approx a_* ^2\left(1+\sqrt1-a_* ^2\right)^-2$ (a good approximation for $a_* \lt 0.95$) and $\dotM_\mathrmEdd=137\,M_\odot \,\mathrmyr^-1$ for $M=6.2\times 10^9\,M_\odot $. This expression was developed for models with a thin disk in the equatorial plane. $P_\mathrmBZ$ is lower for models where the force-free region is excluded by a thicker disk around the equatorial plane. Clearly $P_\mathrmBZ$ is comparable to observational estimates of $P_\mathrmjet$.

In our models (see Table 3) $P_\mathrmjet$ follows the above scaling relation but with a smaller coefficient. The ratio of coefficients is model dependent and varies from 0.15 to 0.83. This is likely because the force-free region is restricted to a cone around the poles of the black hole, and the width of the cone varies by model. Indeed, the coefficient is larger for MAD than for SANE models, which is consistent with this idea because MAD models have a wide funnel and SANE models have a narrow funnel. This also suggests that future comparison of synthetic 43 and 86 GHz images from our models with lower-frequency VLBI data may further constrain the magnetic flux on the black hole.

The connection between the Poynting flux in the funnel and black hole spin has been discussed for some time in the simulation literature, beginning with McKinney & Gammie (2004 ; see also McKinney 2006 ; McKinney & Narayan 2007). The structure of the funnel magnetic field can be time-averaged and shown to match the analytic solution of Blandford & Znajek (1977). Furthermore, the energy flux density can be time-averaged and traced back to the event horizon. Is the energy contained in black hole spin sufficient to drive the observed jet over the jet lifetime ? The spindown timescale is $\tau =(M-M_\mathrmirr)c^2/P_\mathrmjet$, where $M_\mathrmirr\equiv M\left(\left(1+\sqrt1-a_* ^2\right)/2\right)^1/2$ is the irreducible mass of the black hole. For the $a_* =0.94$ MAD model with $R_\mathrmhigh=160$, $\tau =7.3\times 10^12\,\mathrmyr$, which is long compared to a Hubble time ($\sim 10^10$ yr). Indeed, the spindown time for all models is long compared to the Hubble time.

We conclude that for models that have sufficiently powerful jets and are consistent with EHT2017, $P_\mathrmjet$ is driven by extraction of black hole spin energy through the Blandford–Znajek process.
6.4. Constraint Summary

We have applied constraints from AIS, a radiative self-consistency constraint, a constraint on maximum X-ray luminosity, and a constraint on minimum jet power. Which models survive ? Here we consider only models for which we have calculated $L_\rmX

$ and $L_\mathrmbol$. Table 2 summarizes the results. Here we consider only i = 163° (for $a_* \geqslant 0$) and i = 17° (for $a_* \lt 0$). The first three columns give the model parameters. The next four columns show the result of application of each constraint : Themis-AIS (here broken out by individual model rather than groups of models), radiative efficiency ($\epsilon \lt \epsilon _\mathrmthin\mathrmdisk$), $L_\rmX

$, and $P_\mathrmjet$.

Table 2. Rejection Table
Fluxa $a_* $ b $R_\mathrmhigh$ c AISd epsilone $L_\rmX

$ f $P_\mathrmjet$ g
SANE −0.94 1 Fail Pass Pass Pass Fail
SANE −0.94 10 Pass Pass Pass Pass Pass
SANE −0.94 20 Pass Pass Pass Pass Pass
SANE −0.94 40 Pass Pass Pass Pass Pass
SANE −0.94 80 Pass Pass Pass Pass Pass
SANE −0.94 160 Fail Pass Pass Pass Fail
SANE −0.5 1 Pass Pass Fail Fail Fail
SANE −0.5 10 Pass Pass Fail Fail Fail
SANE −0.5 20 Pass Pass Pass Fail Fail
SANE −0.5 40 Pass Pass Pass Fail Fail
SANE −0.5 80 Fail Pass Pass Fail Fail
SANE −0.5 160 Pass Pass Pass Fail Fail
SANE 0 1 Pass Pass Pass Fail Fail
SANE 0 10 Pass Pass Pass Fail Fail
SANE 0 20 Pass Pass Fail Fail Fail
SANE 0 40 Pass Pass Pass Fail Fail
SANE 0 80 Pass Pass Pass Fail Fail
SANE 0 160 Pass Pass Pass Fail Fail
SANE +0.5 1 Pass Pass Pass Fail Fail
SANE +0.5 10 Pass Pass Pass Fail Fail
SANE +0.5 20 Pass Pass Pass Fail Fail
SANE +0.5 40 Pass Pass Pass Fail Fail
SANE +0.5 80 Pass Pass Pass Fail Fail
SANE +0.5 160 Pass Pass Pass Fail Fail
SANE +0.94 1 Pass Fail Pass Fail Fail
SANE +0.94 10 Pass Fail Pass Fail Fail
SANE +0.94 20 Pass Pass Pass Fail Fail
SANE +0.94 40 Pass Pass Pass Fail Fail
SANE +0.94 80 Pass Pass Pass Pass Pass
SANE +0.94 160 Pass Pass Pass Pass Pass
MAD −0.94 1 Fail Fail Pass Pass Fail
MAD −0.94 10 Fail Pass Pass Pass Fail
MAD −0.94 20 Fail Pass Pass Pass Fail
MAD −0.94 40 Fail Pass Pass Pass Fail
MAD −0.94 80 Fail Pass Pass Pass Fail
MAD −0.94 160 Fail Pass Pass Pass Fail
MAD −0.5 1 Pass Fail Pass Fail Fail
MAD −0.5 10 Pass Pass Pass Fail Fail
MAD −0.5 20 Pass Pass Pass Pass Pass
MAD −0.5 40 Pass Pass Pass Pass Pass
MAD −0.5 80 Pass Pass Pass Pass Pass
MAD −0.5 160 Pass Pass Pass Pass Pass
MAD 0 1 Pass Fail Pass Fail Fail
MAD 0 10 Pass Pass Pass Fail Fail
MAD 0 20 Pass Pass Pass Fail Fail
MAD 0 40 Pass Pass Pass Fail Fail
MAD 0 80 Pass Pass Pass Fail Fail
MAD 0 160 Pass Pass Pass Fail Fail
MAD +0.5 1 Pass Fail Pass Fail Fail
MAD +0.5 10 Pass Pass Pass Pass Pass
MAD +0.5 20 Pass Pass Pass Pass Pass
MAD +0.5 40 Pass Pass Pass Pass Pass
MAD +0.5 80 Pass Pass Pass Pass Pass
MAD +0.5 160 Pass Pass Pass Pass Pass
MAD +0.94 1 Pass Fail Fail Pass Fail
MAD +0.94 10 Pass Fail Pass Pass Fail
MAD +0.94 20 Pass Pass Pass Pass Pass
MAD +0.94 40 Pass Pass Pass Pass Pass
MAD +0.94 80 Pass Pass Pass Pass Pass
MAD +0.94 160 Pass Pass Pass Pass Pass

Notes.
aflux : net magnetic flux on the black hole (MAD, SANE). b $a_* $ : dimensionless black hole spin. c $R_\mathrmhigh$ : electron temperature parameter. See Equation (8). dAverage Image Scoring (Themis-AIS), models are rejected if $\langle p\rangle \leqslant 0.01$. See Section 4 and Table 1. eepsilon : radiative efficiency, models are rejected if epsilon is larger than the corresponding thin disk efficiency. See Section 6.1. f $L_\rmX

$ : X-ray luminosity ; models are rejected if $\langle L_\rmX

\rangle 10^-2\sigma \gt 4.4\times 10^40\,\,\mathrmerg\,\rms

^-1$. See Section 6.2. g $P_\mathrmjet$ : jet power, models are rejected if $P_\mathrmjet\leqslant 10^42\,\,\mathrmerg\,\rms

^-1$. See Section 6.3.

Download table as : ASCIITypeset image

The final column gives the logical AND of the previous four columns, and allows a model to pass only if it passes all tests. Evidently most of the SANE models fail, with the exception of some $a_* =-0.94$ models and a few $a_* =0.94$ models with large $R_\mathrmhigh$. A much larger fraction of the MAD models pass, although $a_* =0$ models all fail because of inadequate jet power. MAD models with small $R_\mathrmhigh$ also fail. It is the jet power constraint that rejects the largest number of models.
7. Discussion

We have interpreted the EHT2017 data using a limited library of models with attendant limitations. Many of the limitations stem from the GRMHD model, which treats the plasma as an ideal fluid governed by equations that encode conservation laws for particle number, momentum, and energy. The eDF, in particular, is described by a number density and temperature, rather than a full distribution function, and the electron temperature Te is assumed to be a function of the local ion temperature and plasma $\beta _\rmp

$. Furthermore, all models assume a Kerr black hole spacetime, but there are alternatives. Here we consider some of the model limitations and possible extensions, including to models beyond general relativity.
7.1. Radiative Effects

Post-processed GRMHD simulations that are consistent with EHT data and the flux density of 1.3 mm emission in M87 can yield unphysically large radiative efficiencies (see Section 6). This implies that the radiative cooling timescale is comparable to or less than the advection timescale. As a consequence, including radiative cooling in simulations may be necessary to recover self-consistent models (see Mościbrodzka et al. 2011 ; Dibi et al. 2012). In our models we use a single parameter, $R_\mathrmhigh$, to adjust Te and account for all effects that might influence the electron energy density. How good is this approximation ?

The importance of radiative cooling can be assessed using newly developed, state-of-the-art general relativistic radiation GRMHD (’’’radiation GRMHD’’’) codes. Sa̧dowski et al. (2013b ; see also Sa̧dowski et al. 2014, 2017 ; McKinney et al. 2014) applied the M1 closure (Levermore 1984), which treats the radiation as a relativistic fluid. Ryan et al. (2015) introduced a Monte Carlo radiation GRMHD method, allowing for full frequency-dependent radiation transport. Models for turbulent dissipation into the electrons and ions, as well as heating and cooling physics that sets the temperature ratio Ti/Te, have been added to GRMHD and radiative GRMHD codes and used in simulations of Sgr A* (Ressler et al. 2015, 2017 ; Chael et al. 2018) and M87 (Ryan et al. 2018 ; Chael et al. 2019). While the radiative cooling and Coulomb coupling physics in these simulations is well understood, the particle heating process, especially the relative heating rates of ions and electrons, remains uncertain.

Radiation GRMHD models are computationally expensive per run and do not have the same scaling freedom as the GRMHD models, so they need to be repeatedly re-run with different initial conditions until they produce the correct 1.3 mm flux density. It is therefore impractical to survey the parameter space using radiation GRMHD. It is possible, however, to check individual GRMHD models against existing radiation GRMHD models of M87 (Ryan et al. 2018 ; Chael et al. 2019).

The SANE radiation GRMHD models of Ryan et al. (2018) with $a_* =0.94$ and $M=6\times 10^9\,M_\odot $ can be compared to GRMHD SANE $a_* =0.94$ models at various values of $R_\mathrmhigh$. The radiative models have $\dotM/\dotM_\mathrmEdd=5.2\times 10^-6$ and $P_\mathrmjet\,=5.1\times 10^41\,\mathrmerg\,\rms

^-1$. The GRMHD models in this work have, for $1\leqslant R_\mathrmhigh\leqslant 160$, $0.36\times 10^-6\leqslant \dotM/\dotM_\mathrmEdd\leqslant 20\times 10^-6$, and $0.22\leqslant P_\mathrmjet/(10^41\,\,\mathrmerg\,\rms

^-1)\leqslant 12$ (Table 3). Evidently the mass accretion rates and jet powers in the GRMHD models span a wide range that depends on $R_\mathrmhigh$, but when we choose $R_\mathrmhigh$ = 10 − 20 they are similar to what is found in the radiative GRMHD model when using the turbulent electron heating model (Howes 2010).

Table 3. Model Table
Flux phgr Spin $R_\mathrmhigh$ $L_\mathrmbol/\left(\dotM\,c^2\right)$ $L_\rmX

$ (cgs) $P_\mathrmjet$ (cgs) $P_\mathrmjet/\left(\dotM\,c^2\right)$ $P_\mathrmout$ (cgs) $P_\mathrmout/\left(\dotM\,c^2\right)$ $P_\mathrmjet,\mathrmem/P_\mathrmjet$ $\dotM/\dotM_\mathrmEdd$
SANE 1.02 −0.94 1 1.27 × 10−2 $3.18_\gt 0.20^\lt 49.55\times 10^41$ 1.16 × 1042 5.34 × 10−3 1.19 × 1042 5.48 × 10−3 0.84 2.77 × 10−5
SANE 1.02 −0.94 10 1.6 × 10−3 $9.62_\gt 1.44^\lt 64.42\times 10^40$ 4.94 × 1042 5.34 × 10−3 5.07 × 1042 5.48 × 10−3 0.84 1.19 × 10−4
SANE 1.02 −0.94 20 6.09 × 10−4 $3.26_\gt 0.90^\lt 11.86\times 10^40$ 5.8 × 1042 5.34 × 10−3 5.96 × 1042 5.48 × 10−3 0.84 1.39 × 10−4
SANE 1.02 −0.94 40 2.45 × 10−4 $8.89_\gt 1.56^\lt 50.53\times 10^39$ 7.02 × 1042 5.34 × 10−3 7.21 × 1042 5.48 × 10−3 0.84 1.69 × 10−4
SANE 1.02 −0.94 80 1.33 × 10−4 $2.65_\gt 0.39^\lt 18.26\times 10^39$ 8.89 × 1042 5.34 × 10−3 9.13 × 1042 5.48 × 10−3 0.84 2.13 × 10−4
SANE 1.02 −0.94 160 7.12 × 10−5 $6.36_\gt 0.73^\lt 55.27\times 10^38$ 1.2 × 1043 5.34 × 10−3 1.23 × 1043 5.48 × 10−3 0.84 2.87 × 10−4
SANE 1.11 −0.5 1 1.62 × 10−2 $1.97_\gt 0.98^\lt 3.94\times 10^41$ 2.62 × 1040 1.86 × 10−4 3.84 × 1040 2.72 × 10−4 0.88 1.81 × 10−5
SANE 1.11 −0.5 10 2.17 × 10−3 $1.94_\gt 0.69^\lt 5.40\times 10^41$ 1.95 × 1041 1.86 × 10−4 2.85 × 1041 2.72 × 10−4 0.88 1.34 × 10−4
SANE 1.11 −0.5 20 6.69 × 10−4 $3.72_\gt 1.80^\lt 7.72\times 10^40$ 2.26 × 1041 1.86 × 10−4 3.31 × 1041 2.72 × 10−4 0.88 1.56 × 10−4
SANE 1.11 −0.5 40 2.47 × 10−4 $9.44_\gt 6.67^\lt 13.37\times 10^39$ 2.62 × 1041 1.86 × 10−4 3.83 × 1041 2.72 × 10−4 0.88 1.81 × 10−4
SANE 1.11 −0.5 80 1.26 × 10−4 $1.23_\gt 0.33^\lt 4.58\times 10^39$ 3.2 × 1041 1.86 × 10−4 4.68 × 1041 2.72 × 10−4 0.88 2.21 × 10−4
SANE 1.11 −0.5 160 7.86 × 10−5 $3.72_\gt 0.83^\lt 16.68\times 10^38$ 4.21 × 1041 1.86 × 10−4 6.16 × 1041 2.72 × 10−4 0.88 2.9 × 10−4
SANE 0.99 0 1 3.17 × 10−2 $2.08_\gt 0.02^\lt 194.22\times 10^41$ 2.24 × 1036 4.4 × 10−8 5.22 × 1039 1.03 × 10−4 1.01 6.5 × 10−6
SANE 0.99 0 10 1.88 × 10−2 $4.2_\gt 0.04^\lt 425.40\times 10^42$ 4.38 × 1037 4.4 × 10−8 1.02 × 1041 1.03 × 10−4 1.01 1.27 × 10−4
SANE 0.99 0 20 5.83 × 10−3 $1.57_\gt 0.06^\lt 39.69\times 10^42$ 8.02 × 1037 4.4 × 10−8 1.87 × 1041 1.03 × 10−4 1.01 2.33 × 10−4
SANE 0.99 0 40 7.8 × 10−4 $8.92_\gt 1.92^\lt 41.45\times 10^40$ 9.16 × 1037 4.4 × 10−8 2.14 × 1041 1.03 × 10−4 1.01 2.66 × 10−4
SANE 0.99 0 80 1.69 × 10−4 $2.5_\gt 0.33^\lt 19.17\times 10^39$ 1.03 × 1038 4.4 × 10−8 2.41 × 1041 1.03 × 10−4 1.01 3 × 10−4
SANE 0.99 0 160 1.08 × 10−4 $3.44_\gt 0.89^\lt 13.32\times 10^38$ 1.23 × 1038 4.4 × 10−8 2.87 × 1041 1.03 × 10−4 1.01 3.57 × 10−4
SANE 1.10 0.5 1 4.97 × 10−2 $5.5_\gt 0.88^\lt 34.41\times 10^40$ 2.57 × 1039 1.63 × 10−4 9.19 × 1039 5.86 × 10−4 0.88 2.01 × 10−6
SANE 1.10 0.5 10 5.98 × 10−3 $4.73_\gt 0.25^\lt 88.59\times 10^40$ 1.91 × 1040 1.64 × 10−4 6.84 × 1040 5.86 × 10−4 0.88 1.5 × 10−5
SANE 1.10 0.5 20 3.33 × 10−3 $3.83_\gt 0.30^\lt 49.18\times 10^40$ 4.09 × 1040 1.64 × 10−4 1.47 × 1041 5.86 × 10−4 0.88 3.2 × 10−5
SANE 1.10 0.5 40 1.74 × 10−3 $2.52_\gt 0.28^\lt 22.73\times 10^40$ 8.02 × 1040 1.64 × 10−4 2.87 × 1041 5.86 × 10−4 0.88 6.28 × 10−5
SANE 1.10 0.5 80 6.95 × 10−4 $7.84_\gt 0.67^\lt 91.92\times 10^39$ 1.27 × 1041 1.64 × 10−4 4.55 × 1041 5.86 × 10−4 0.88 9.95 × 10−5
SANE 1.10 0.5 160 2.78 × 10−4 $1.37_\gt 0.08^\lt 22.85\times 10^39$ 1.69 × 1041 1.63 × 10−4 6.06 × 1041 5.86 × 10−4 0.88 1.33 × 10−4
SANE 1.64 0.94 1 1.4 $2.38_\gt 0.02^\lt 359.03\times 10^41$ 2.2 × 1040 7.76 × 10−3 3.38 × 1040 1.19 × 10−2 0.82 3.63 × 10−7
SANE 1.64 0.94 10 2.7 × 10−1 $2.79_\gt 0.02^\lt 508.99\times 10^41$ 1.4 × 1041 7.76 × 10−3 2.15 × 1041 1.19 × 10−2 0.82 2.31 × 10−6
SANE 1.64 0.94 20 1.74 × 10−1 $5.75_\gt 0.02^\lt 1685.98\times 10^41$ 3.22 × 1041 7.76 × 10−3 4.94 × 1041 1.19 × 10−2 0.82 5.31 × 10−6
SANE 1.64 0.94 40 7.2 × 10−2 $4.71_\gt 0.01^\lt 2490.36\times 10^41$ 5.97 × 1041 7.76 × 10−3 9.17 × 1041 1.19 × 10−2 0.82 9.84 × 10−6
SANE 1.64 0.94 80 2.38 × 10−2 $1.42_\gt 0.00^\lt 860.83\times 10^41$ 8.87 × 1041 7.76 × 10−3 1.36 × 1042 1.19 × 10−2 0.82 1.46 × 10−5
SANE 1.64 0.94 160 8.45 × 10−3 $3.22_\gt 0.01^\lt 1687.88\times 10^40$ 1.23 × 1042 7.76 × 10−3 1.89 × 1042 1.19 × 10−2 0.82 2.03 × 10−5
MAD 8.04 −0.94 1 7.61 × 10−1 $2.12_\gt 0.25^\lt 17.74\times 10^41$ 1.36 × 1042 2.09 × 10−1 1.6 × 1042 2.46 × 10−1 0.75 8.32 × 10−7
MAD 8.04 −0.94 10 7.54 × 10−2 $5.76_\gt 0.49^\lt 68.06\times 10^40$ 1.97 × 1042 2.09 × 10−1 2.32 × 1042 2.46 × 10−1 0.75 1.21 × 10−6
MAD 8.04 −0.94 20 3.76 × 10−2 $2.27_\gt 0.18^\lt 29.09\times 10^40$ 2.38 × 1042 2.09 × 10−1 2.8 × 1042 2.46 × 10−1 0.75 1.46 × 10−6
MAD 8.04 −0.94 40 2.07 × 10−2 $6.18_\gt 0.49^\lt 77.36\times 10^39$ 3 × 1042 2.09 × 10−1 3.54 × 1042 2.46 × 10−1 0.75 1.84 × 10−6
MAD 8.04 −0.94 80 1.17 × 10−2 $1.32_\gt 0.07^\lt 26.36\times 10^39$ 3.99 × 1042 2.09 × 10−1 4.71 × 1042 2.46 × 10−1 0.75 2.45 × 10−6
MAD 8.04 −0.94 160 6.52 × 10−3 $2.57_\gt 0.14^\lt 46.76\times 10^38$ 5.7 × 1042 2.09 × 10−1 6.73 × 1042 2.46 × 10−1 0.75 3.5 × 10−6
MAD 12.25 −0.5 1 2.96 × 10−1 $1.39_\gt 0.17^\lt 11.56\times 10^41$ 3.43 × 1041 4.91 × 10−2 6.04 × 1041 8.64 × 10−2 0.82 8.95 × 10−7
MAD 12.25 −0.5 10 4.53 × 10−2 $2.43_\gt 0.30^\lt 19.86\times 10^40$ 5.31 × 1041 4.92 × 10−2 9.33 × 1041 8.64 × 10−2 0.82 1.38 × 10−6
MAD 12.25 −0.5 20 2.67 × 10−2 $8.18_\gt 0.86^\lt 77.51\times 10^39$ 6.45 × 1041 4.92 × 10−2 1.13 × 1042 8.64 × 10−2 0.82 1.68 × 10−6
MAD 12.25 −0.5 40 1.69 × 10−2 $2.17_\gt 0.21^\lt 22.33\times 10^39$ 8.07 × 1041 4.92 × 10−2 1.42 × 1042 8.64 × 10−2 0.82 2.1 × 10−6
MAD 12.25 −0.5 80 1.07 × 10−2 $4.87_\gt 0.47^\lt 50.76\times 10^38$ 1.05 × 1042 4.92 × 10−2 1.85 × 1042 8.64 × 10−2 0.82 2.74 × 10−6
MAD 12.25 −0.5 160 6.43 × 10−3 $1.09_\gt 0.17^\lt 7.06\times 10^38$ 1.46 × 1042 4.92 × 10−2 2.57 × 1042 8.64 × 10−2 0.82 3.81 × 10−6
MAD 15.44 0 1 2.67 × 10−1 $1.22_\gt 0.10^\lt 14.60\times 10^41$ 0.0 0.0 8.39 × 1040 1.51 × 10−2 0.00 7.12 × 10−7
MAD 15.44 0 10 4.53 × 10−2 $1.86_\gt 0.11^\lt 31.55\times 10^40$ 0.0 0.0 1.39 × 1041 1.51 × 10−2 0.00 1.18 × 10−6
MAD 15.44 0 20 2.81 × 10−2 $5.98_\gt 0.35^\lt 101.81\times 10^39$ 0.0 0.0 1.71 × 1041 1.51 × 10−2 0.00 1.46 × 10−6
MAD 15.44 0 40 1.85 × 10−2 $1.63_\gt 0.10^\lt 27.75\times 10^39$ 0.0 0.0 2.15 × 1041 1.51 × 10−2 0.00 1.82 × 10−6
MAD 15.44 0 80 1.21 × 10−2 $3.51_\gt 0.20^\lt 61.34\times 10^38$ 0.0 0.0 2.77 × 1041 1.51 × 10−2 0.00 2.35 × 10−6
MAD 15.44 0 160 7.63 × 10−3 $8.06_\gt 0.81^\lt 80.62\times 10^37$ 0.0 0.0 3.73 × 1041 1.51 × 10−2 0.00 3.17 × 10−6
MAD 15.95 0.5 1 5.45 × 10−1 $1.57_\gt 0.21^\lt 11.98\times 10^41$ 4.64 × 1041 1.16 × 10−1 6.74 × 1041 1.69 × 10−1 0.85 5.11 × 10−7
MAD 15.95 0.5 10 9.45 × 10−2 $2.71_\gt 0.20^\lt 36.30\times 10^40$ 8.07 × 1041 1.16 × 10−1 1.17 × 1042 1.69 × 10−1 0.85 8.89 × 10−7
MAD 15.95 0.5 20 5.54 × 10−2 $9.67_\gt 0.74^\lt 126.69\times 10^39$ 1.02 × 1042 1.16 × 10−1 1.49 × 1042 1.69 × 10−1 0.85 1.13 × 10−6
MAD 15.95 0.5 40 3.5 × 10−2 $3.3_\gt 0.28^\lt 39.01\times 10^39$ 1.32 × 1042 1.16 × 10−1 1.92 × 1042 1.69 × 10−1 0.85 1.45 × 10−6
MAD 15.95 0.5 80 2.22 × 10−2 $8_\gt 0.70^\lt 91.84\times 10^38$ 1.74 × 1042 1.16 × 10−1 2.52 × 1042 1.69 × 10−1 0.85 1.92 × 10−6
MAD 15.95 0.5 160 1.35 × 10−2 $1.79_\gt 0.38^\lt 8.44\times 10^38$ 2.38 × 1042 1.16 × 10−1 3.46 × 1042 1.69 × 10−1 0.85 2.62 × 10−6
MAD 12.78 0.94 1 3.65 $5.19_\gt 0.62^\lt 43.60\times 10^41$ 1.97 × 1042 8.23 × 10−1 2.29 × 1042 9.55 × 10−1 0.80 3.07 × 10−7
MAD 12.78 0.94 10 3.68 × 10−1 $1.3_\gt 0.13^\lt 13.22\times 10^41$ 3.04 × 1042 8.23 × 10−1 3.52 × 1042 9.55 × 10−1 0.80 4.73 × 10−7
MAD 12.78 0.94 20 1.79 × 10−1 $5_\gt 0.44^\lt 56.22\times 10^40$ 3.73 × 1042 8.23 × 10−1 4.33 × 1042 9.55 × 10−1 0.80 5.81 × 10−7
MAD 12.78 0.94 40 9.43 × 10−2 $1.54_\gt 0.11^\lt 22.13\times 10^40$ 4.74 × 1042 8.23 × 10−1 5.5 × 1042 9.55 × 10−1 0.80 7.38 × 10−7
MAD 12.78 0.94 80 5.19 × 10−2 $3.74_\gt 0.17^\lt 80.85\times 10^39$ 6.26 × 1042 8.23 × 10−1 7.27 × 1042 9.55 × 10−1 0.80 9.75 × 10−7
MAD 12.78 0.94 160 2.82 × 10−2 $6.97_\gt 0.26^\lt 186.48\times 10^38$ 8.75 × 1042 8.23 × 10−1 1.02 × 1043 9.55 × 10−1 0.80 1.36 × 10−6

Download table as : ASCIITypeset images : 1 2

We have also directly compared the Te distribution in the emitting region, and the radiation GRMHD model is quite close to the $R_\mathrmhigh=10$ model. The resulting images are qualitatively similar, with an asymmetric photon ring that is brighter in the south and a weak inner ring associated with the funnel wall emission as in Figure 2. The radiation GRMHD SANE model, like all our nonradiative GRMHD SANE models (except the $R_\mathrmhigh=160$ model), would be ruled out by the condition $P_\mathrmjet\gt 10^42\,\mathrmerg\,\rms

^-1$.

The MAD radiation GRMHD models of Chael et al. (2019) with $a_* =0.94$ and $M=6.2\times 10^9M_\odot $ can be compared to GRMHD MAD $a_* =0.94$ models at various values of $R_\mathrmhigh$. Chael et al. (2019) uses two dissipation models : the Howes (2010, hereafter H10) model of heating from a Landau-damped turbulent cascade, and the Rowan et al. (2017, hereafter R17) model of heating based on simulations of transrelativistic magnetic reconnection. The (H10, R17) models have $\dotM/\dotM_\mathrmEdd=(3.6,2.3)\times 10^-6$ and $P_\mathrmjet=(6.6,13)\,\times 10^42\,\mathrmerg\,\rms

^-1$. The GRMHD models have, for $1\leqslant R_\mathrmhigh\,\leqslant 160$, $0.13\times 10^-6\leqslant \dotM/\dotM_\mathrmEdd\leqslant 1.4\times 10^-6$ and $2.3\leqslant P_\mathrmjet/(10^42\,\,\mathrmerg\,\rms

^-1)\leqslant 8.8$ (Table 3). In the radiation GRMHD MAD models $\dotM$ lies in the middle of the range spanned by the nonradiative GRMHD models, and jet power lies at the upper end of the range spanned by the nonradiative GRMHD models. The Te distributions in the radiative and nonradiative MAD models differ : the mode of the radiation GRMHD model Te distribution is about a factor of 3 below the mode of the Te distribution in the $R_\mathrmhigh=20$ GRMHD model, and the GRMHD model has many more zones at $\rm\Theta

_e\sim 100$ that contribute to the final image than the radiation GRMHD models. This difference is a consequence of the $R_\mathrmhigh$ model for Te : in MAD models almost all the emission emerges at $\beta _\rmp

\lesssim 1$, so $R_\mathrmhigh$, which changes Te in the $\beta _\rmp

\gt 1$ region, offers little control over Te in the emission region. Nevertheless, the jet power and accretion rates are similar in the radiative and nonradiative MAD models, and the time-averaged radiative and nonradiative images are qualitatively indistinguishable. This suggests that the image is determined mainly by the spacetime geometry and is insensitive to the details of the plasma evolution.

This review of radiative effects is encouraging but incomplete : it only considers a limited selection of models and a narrow set of observational constraints. Future studies of time dependence and polarization are likely to sharpen the contrast between radiative and nonradiative models.
7.2. Nonthermal Electrons

Throughout this Letter we have considered only a thermal eDF. While a thermal eDF can account for the observed emission at mm wavelengths in M87 (e.g., Mościbrodzka et al. 2016 ; Prieto et al. 2016 ; Ryan et al. 2018 ; Chael et al. 2019), eDFs that include a nonthermal tail can also explain the observed SED (Broderick & Loeb 2009 ; Yu et al. 2011 ; Dexter et al. 2012 ; Li et al. 2016 ; Davelaar et al. 2018 ; J. Davelaar et al. 2019, in preparation).

The role of nonthermal electrons (and positrons) in producing the observed compact emission is not a settled question, and cannot be settled in this first investigation of EHT2017 models, but there are constraints. The number density, mean velocity, and energy density of the eDF are fixed or limited by the GRMHD models. In addition, the eDF cannot on average sustain features that would be erased by kinetic instabilities on timescales short compared to $r_\rmg

c^-1$. Some nonthermal eDFs increase $F_\nu /\dotM$ in comparison to a thermal eDF, implying lower values of $\dotM$ than quoted above (Ball et al. 2018 ; Davelaar et al. 2018 ; J. Davelaar et al. 2019, in preparation). These lower values of $\dotM$ can slightly change the source morphology, e.g., by decreasing the visibility of the approaching jet (e.g., Dexter et al. 2012).

One can evaluate the influence of nonthermal eDFs in several ways. For example, it is possible to study simplified, phenomenological models. Emission features due to the cooling of nonthermal electrons may then reveal how and where the nonthermal electrons are produced (Pu et al. 2017). Emission features created by the injection of nonthermal electrons within GRMHD models of the jet and their subsequent cooling will be studied separately (T. Kawashima et al. 2019, in preparation). The effect of nonthermal eDFs can also be studied by post-processing of ideal GRMHD models if one assumes that the electrons have a fixed, parameterized form such as a power-law distribution (Dexter et al. 2012) or a κ-distribution (Davelaar et al. 2018 ; J. Davelaar et al. 2019, in preparation). These parameterized models produce SEDs that agree with radio to near-infrared data, but they are approximations to the underlying physics and do not resolve the microscopic processes that accelerate particles. One can also include dissipative processes explicitly in the GRMHD models, including scalar resistivity (Palenzuela et al. 2009 ; Dionysopoulou et al. 2013 ; Del Zanna et al. 2016 ; Qian et al. 2017 ; Ripperda et al. 2019), heat fluxes and viscosities (pressure anisotropies ; Chandra et al. 2015 ; Ressler et al. 2015 ; Foucart et al. 2017), and particle acceleration (e.g., Chael et al. 2017). Ultimately special and general relativistic particle-in-cell codes (Watson & Nishikawa 2010 ; Chen et al. 2018 ; Levinson & Cerutti 2018 ; Parfrey et al. 2019) will enable direct investigations of kinetic processes.
7.3. Other Models and Analysis Limitations

We have used a number of other approximations in generating our models. Among the most serious ones are as follows.

(1) Fast Light Approximation. A GRMHD simulation produces a set of dump files containing the model state at a single global (Kerr–Schild) coordinate time. Because the dynamical time is only slightly longer than the light-crossing time, in principle one needs to trace rays through a range of coordinate times, i.e., by interpolation between multiple closely spaced dump files. In practice this is difficult because a high cadence of output files is required, limiting the speed of the GRMHD simulations and requiring prohibitively large data storage. In addition, the cost of ray tracing through multiple output files is high. Because of this, we adopt the commonly used fast light approximation in which GRMHD variables are read from a single dump file and held steady during the ray tracing. Including light-travel time delays produces minor changes to the small-scale image structure and to light curves (e.g., Dexter et al. 2010 ; Bronzwaer et al. 2018 ; Z. Younsi et al. 2019b, in preparation), although it is essential for the study of variability on the light-crossing timescale.

(2) Untilted Disks. We have assumed that the disk angular momentum vector and black hole spin vector are (anti-) aligned. There is no reason for the angular momentum vector of the accretion flow on large scales to align with the black hole spin vector, and there is abundant evidence for misaligned disks in AGNs (e.g., Miyoshi et al. 1995). How might disk tilt affect our results ?

Tilting the disk by as little as 15° is enough to set up a standing, two-armed spiral shock close to the ISCO (Fragile & Blaes 2008). This shock directly affects the morphology of mm wavelength images, especially at low inclination, in models of Sgr A* (Dexter & Fragile 2013, especially Figure 5), producing an obvious two-armed spiral pattern on the sky. If this structure were also present in images of tilted models of M87, then it is possible that even a modest tilt could be ruled out.

If a modest tilt is present in M87 it is unlikely to affect our conclusion regarding the sign of black hole spin. That conclusion depends on emission from funnel wall plasma in counter-rotating ($a_* \lt 0$) disks. The funnel wall plasma is loaded onto funnel plasma field lines by local instabilities at the wall and then rotates with the funnel and therefore the black hole (Wong et al. 2019). The funnel wall is already unsteady, fluctuating by tens of degrees in azimuth and in time, so a modest tilt seems unlikely to dramatically alter the funnel wall structure.

Is there observational evidence for tilt in M87 ? In numerical studies of tilted disks the jet emerges perpendicular to the disk (Liska et al. 2018), and tilted disks are expected to precess. One might then expect that a tilted source would produce a jet that exhibits periodic variations, or periodic changes in jet direction with distance from the source, as seen in other sources. There is little evidence of this in M87 (see Park et al. 2019 for a discussion of possible misaligned structure in the jet). Indeed, Walker et al. (2018) saw at most small displacements of the jet with time and distance from the source at mas scales. In sum, there is therefore little observational motivation for considering tilted disk models.

Tilted disk models of M87 are an interesting area for future study. It is possible that the inner disk may align with the black hole via a thick-disk variant of the Bardeen & Petterson (1975) effect. Existing tilted thick-disk GRMHD simulations (e.g., Fragile et al. 2007 ; McKinney et al. 2013 ; Shiokawa 2013 ; Liska et al. 2018) show some evidence for alignment and precession (McKinney et al. 2013 ; Shiokawa 2013 ; Liska et al. 2018), but understanding of the precession and alignment timescales is incomplete. It will be challenging to extend the Image Library to include a survey of tilted disk models, however, because with tilted disks there are two new parameters : the two angles that describe the orientation of the outer disk with respect to the black hole spin vector and the line of sight.

(3) Pair Production. In some models of M87 the mm emission is dominated by electron-positron pairs within the funnel, even close to the horizon scale (see Beskin et al. 1992 ; Levinson & Rieger 2011 ; Mościbrodzka et al. 2011 ; Broderick & Tchekhovskoy 2015 ; Hirotani & Pu 2016). The pairs are produced from the background radiation field or from a pair-cascade process following particle acceleration by unscreened electric fields, which we cannot evaluate using ideal GRMHD models. We leave it to future work to assess whether or not these models can plausibly suppress emission from the disk and funnel wall, and simultaneously produce a sufficiently powerful jet.

(4) Numerical Treatment of Low-density Regions. Virtually all MHD simulations, including ours, use a ’’’floor’’’ procedure that resets the density if it falls below a minimum value. If this is not done, then truncation error accumulates dramatically in the low-density regions and the solution is corrupted. If the volume where floors are activated contains only a small fraction of the simulation mass, momentum, and energy, then most aspects of the solution are unaffected by this procedure (e.g., McKinney & Gammie 2004).

In regions where the floors are activated the temperature of the plasma is no longer reliable. This is why we cut off emission from regions with $B^2/\rho \gt 1$, where floors are commonly activated. In models where floors are only activated in the funnel (e.g., most SANE models), the resulting images are insensitive to the choice of cutoff $B^2/\rho $. In MAD models the regions of low and high density are mixed because lightly loaded magnetic field lines that are trapped in the hole bubbles outward through the disk. In this case emission at $\nu \,\gt 230$ GHz can be sensitive to the choice of cutoff $B^2/\rho $ Chael et al. (2019). The sense of the effect is that greater cutoff $B^2/\rho $ implies more emission at high frequency. Our use of a cutoff $B^2/\rho =1$ is therefore likely to underestimate mm emission and therefore overestimate $\dotM$ and $P_\mathrmjet$. Accurate treatment of the dynamics and thermodynamics of low-density regions and especially sharp boundaries between low- and high-density regions is a fundamental numerical problem in black hole accretion flow modeling that merits further attention.
7.4. Alternatives to Kerr Black Holes

Although our working hypothesis has been that M87 contains a Kerr black hole, it is interesting to consider whether or not the data is also consistent with alternative models for the central object. These alternatives can be grouped into three main categories : (i) black holes within general relativity that include additional fields ; (ii) black hole solutions from alternative theories of gravity or incorporating quantum effects ; (iii) black hole ’’’mimickers,’’’ i.e., compact objects, both within general relativity or in alternative theories, whose properties could be fine-tuned to resemble those of black holes.

The first category includes, for example, black holes in Einstein–Maxwell–dilaton-axion gravity (e.g., García et al. 1995 ; Mizuno et al. 2018), black holes with electromagnetic or Newman–Unti-Tamburino (NUT) charges (e.g., Grenzebach et al. 2014), regular black holes in nonlinear electrodynamics (e.g., Abdujabbarov et al. 2016), black hole metrics affected by a cosmological constant (e.g., Dymnikova 1992) or a dark matter halo (e.g., Hou et al. 2018), and black holes with scalar wigs (e.g., Barranco et al. 2017) or hair (e.g., Herdeiro & Radu 2014). While the shadows of this class of compact objects are expected to be similar to Kerr and therefore cannot be ruled out immediately by current observations (Mizuno et al. 2018), the most extreme examples of black holes surrounded by massive scalar field configurations should produce additional lobes in the shadow or disconnected dark regions (Cunha et al. 2015). As these features are not found in the EHT2017 image, these alternatives are not viable models for M87.

The second category comprises black hole solutions with classical modifications to general relativity, as well as effects coming from approaches to quantum gravity (see, e.g., Moffat 2015 ; Dastan et al. 2016 ; Younsi et al. 2016 ; Amir et al. 2018 ; Eiroa & Sendra 2018 ; Giddings & Psaltis 2018). These alternatives have shadows that are qualitatively very similar to those of Kerr black holes and are not distinguishable with present EHT capabilities. However, higher-frequency observations, together with the degree of polarization of the emitted radiation or the variability of the accretion flow, can be used to assess their viability.

Finally, the third category comprises compact objects such as spherically symmetric naked singularities (e.g., Joshi et al. 2014), superspinars (Kerr with $| a_* | \gt 1$, which are axisymmetric spacetime with naked singularities), and regular horizonless objects, either with or without a surface. Examples of regular surfaceless objects are : boson stars (Kaup 1968), traversable wormholes, and clumps of self-interacting dark matter (Saxton et al. 2016), while examples of black hole mimickers with a surface are gravastars (Mazur & Mottola 2004) and collapsed polymers (Brustein & Medved 2017), to cite only a few. Because the exotic genesis of these black hole mimickers is essentially unknown, their physical properties are essentially unconstrained, thus making the distinction from black holes rather challenging (see, however, Chirenti & Rezzolla 2007, 2016). Nevertheless, some conclusions can drawn already. For instance, the shadow of a superspinar is very different from that of a black hole (Bambi & Freese 2009), and the EHT2017 observations rule out any superspinar model for M87. Similarly, for certain parameter ranges, the shadows of spherically symmetric naked singularities have been found to consist of a filled disk with no dark region120 in the center (Shaikh et al. 2019) ; clearly, this class of models is ruled out. In the same vein, because the shadows of wormholes can exhibit large deviations from those of black holes (see, e.g., Bambi 2013 ; Nedkova et al. 2013 ; Shaikh 2018), a large portion of the corresponding space of parameters can be constrained with the present observations.

A comparison of EHT2017 data with the boson star model, as a representative horizonless and surfaceless black hole mimicker, and a gravastar model as a representative horizonless black hole mimicker, will be presented in Olivares et al. (2019a). Both models produce images with ring-like features similar to those observed by EHT2017, which are consistent with the results of Broderick & Narayan (2006), who also consider black hole alternatives with a surface. The boson star generically requires masses that are substantially different from that expected for M87 (H. Olivares et al. 2019b, in preparation), while the gravastar has accretion variability that is considerably different from that onto a black hole.

In summary, because each of the many exotic alternatives to Kerr black holes can span an enormous space of parameters that is only poorly constrained, the comparisons carried out here must be considered preliminary. Nevertheless, they show that the EHT2017 observations are not consistent with several of the alternatives to Kerr black holes, and that some of those models that produce similar images show rather different dynamics in the accretion flow and in its variability. Future observations and more detailed theoretical modeling, combined with multiwavelength campaigns and polarimetric measurements, will further constrain alternatives to Kerr black holes.
8. Conclusion

In this Letter we have made a first attempt at understanding the physical implications of a single, high-quality EHT data set for M87. We have compared the data to a library of mock images produced from GRMHD simulations by GRRT calculations. The library covers a parameter space that is substantially larger than earlier model surveys. The results of this comparison are consistent with the hypothesis that the compact 1.3 mm emission in M87 arises within a few $r_\rmg

$ of a Kerr black hole, and that the ring-like structure of the image is generated by strong gravitational lensing and Doppler beaming. The models predict that the asymmetry of the image depends on the sense of black hole spin. If this interpretation is accurate, then the spin vector of the black hole in M87 points away from Earth (the black hole spins clockwise on the sky). The models also predict that there is a strong energy flux directed away from the poles of the black hole, and that this energy flux is electromagnetically dominated. If the models are correct, then the central engine for the M87 jet is powered by the electromagnetic extraction of free energy associated with black hole spin via the Blandford–Znajek process.

In our models, M87’s compact mm emission is generated by the synchrotron mechanism. Our ability to make physical inferences based on the models is therefore intimately tied to the quality of our understanding of the eDF. We have used a thermal model with a single free parameter that adjusts the ratio of ion to electron temperature in regions with plasma $\beta _\rmp

\gt 1$ (i.e., regions where magnetic pressure is less than gas pressure). This simple model does not span the range of possible plasma behavior. The theory of high temperature, collisionless plasmas must be better understood if this core physical uncertainty of sub-Eddington black hole accretion is to be eliminated. At present our understanding is inadequate, and alternative eDF models occupy a large, difficult-to-explore parameter space with the potential to surprise. Despite these uncertainties, many of the models produce images with similar morphology that is consistent with EHT2017 data. This suggests that the image shape is controlled mainly by gravitational lensing and the spacetime geometry, rather than details of the plasma physics.

Although the EHT2017 images are consistent with the vast majority of our models, parts of the parameter space can be rejected on physical grounds or by comparison with contemporaneous data at other wavelengths. We reject some models because, even though all models are variable, some models are too variable to be consistent with the data. We can also reject models based on a radiative efficiency cut (the models are not self-consistent and would cool quickly if radiative effects were included), an X-ray luminosity cut using contemporaneous Chandra and NuSTAR data, and on a jet-power cut. The requirement that the jet power exceed a conservative lower limit of $10^42\,\mathrmerg\,\rms

^-1$ turns out to eliminate many models, including all models with $a_* =0$.

We have examined the astrophysical implications of only a subset of EHT2017 data ; much remains to be done, and there are significant opportunities for further constraining the models. EHT2017 data includes tracks from four separate days of observing ; each day is $2.8\,r_\rmg

c^-1$ (see Paper IV). This timescale is short compared to the decorrelation timescale of simulated images, which is $\sim 50\,r_\rmg

c^-1$, and smaller than the light-crossing time of the source plasma. Analysis techniques that use short-timescale variations in the data will need to be developed and are likely to recover new, more stringent constraints on the model from the EHT2017 data set. EHT2017 took polarized data as well. Our simulations already predict full polarization maps, albeit for our simple eDF model. Comparison of model polarization maps of the source with EHT2017 data are likely to sharply limit the space of allowed models (Mościbrodzka et al. 2017). Finally, in this Letter the only multiwavelength companion data that we consider are X-ray observations. Simultaneous data are available at many other wavelengths, from the radio to the gamma-rays, and is likely to further limit the range of acceptable models and guide the implementation of predictive electron physics models.

In this Letter we have focused on the time-dependent ideal GRMHD model. Physically motivated, semi-analytic models including nonthermal emission have not been applied yet and will be discussed in future papers (A. E. Broderick et al. 2019b, in preparation ; T. Kawashima et al. 2019, in preparation ; H.-Y. Pu et al. 2019, in preparation).

We have also not yet considered how the physical properties of the jet are constrained by lower-frequency VLBI observations, which constrain jet kinematics (Mertens et al. 2016 ; Britzen et al. 2017 ; Hada et al. 2017 ; Kim et al. 2018 ; Walker et al. 2018), the jet width profile (Asada & Nakamura 2012 ; Hada et al. 2013 ; Nakamura et al. 2018), the total jet power at kilo-parsec scale (Owen et al. 2000 ; Stawarz et al. 2006), the jet power (e.g., Kino et al. 2014, 2015), the core shift (Hada et al. 2011), and the symmetric limb-brightening structure (Takahashi et al. 2018 ; Kim et al. 2018). The jet width profile is potentially very interesting because it depends on the magnetic flux phgr : the jet internal magnetic pressure $\propto \phi ^2$. We therefore expect (and see in our numerical simulations ; see Figure 4) that MAD jets are wider at the base than SANE jets. Future theoretical work will help connect the ring-like structure seen in EHT2017 to the large-scale jet (M. Nakamura et al. 2019, in preparation).

A second epoch of observations ($\gtrsim 50\,r_\rmg

c^-1\sim 2$ weeks after EHT2017, when the models suggest that source structure will decorrelate) will increase the power of the average image analysis to reject models. The EHT2017 data were able to reject one entire category of models with confidence : high magnetic flux (MAD), retrograde, high-spin models. Other categories of models, such as the low magnetic flux, high-spin models, are assigned comparatively low probabilities by the average image scoring scheme. Data taken later, more than a decorrelation time after EHT2017 (model decorrelation times are of order two weeks), will provide an independent realization of the source. The probabilities attached to individual models by average image scoring will then multiply. For example, a model with probability 0.05 that is assigned probability 0.05 in comparison to a second epoch of observation would then have probability $0.05^2=2.5\times 10^-3$, and would be strongly disfavored by the average image scoring criterion (see Section 4).

Future EHT 345 GHz campaigns (Paper II) will provide excellent constraints, particularly on the width of the ring. The optical depth on every line of sight through the source is expected to decrease (the drop is model and location dependent). In our models this makes the ring narrower, better defined, easier to measure accurately from VLBI data, and less dependent on details of the source plasma model.

Certain features of the model are geometric and should be present in future EHT observations. The photon ring is a persistent feature of the model related to the mass and distance to the black hole. It should be present in the next EHT campaign unless there is a dramatic change in $\dotM$, which would be evident in the SED. The asymmetry in the photon ring is also a persistent feature of the model because, we have argued, it is controlled by the black hole spin. The asymmetry should therefore remain in the southern half of the ring for the next EHT campaign, unless there is a dramatic tilt of the inner accretion flow. If the small-scale and large-scale jet are aligned, then EHT2017 saw the brightest region at unusually small PA, and future campaigns are likely (but not certain) to see the peak brightness shift further to the west. Future 230 GHz EHT campaigns (Paper II) will thus sharply test the GRMHD source models.

Together with complementary studies that are presently targeting either the supermassive black hole candidate at the Galactic Center (Eckart & Genzel 1997 ; Ghez et al. 1998 ; Gravity Collaboration et al. 2018a, 2018b) or stellar-mass binary black holes whose gravitational-wave emission is recorded by the LIGO and Virgo detectors (Abbott et al. 2016), the results provided here are consistent with the existence of astrophysical black holes. More importantly, they clearly indicate that their phenomenology, despite being observed on mass scales that differ by eight orders of magnitude, follows very closely the one predicted by general relativity. This demonstrates the complementarity of experiments studying black holes on all scales, promising much imp roved tests of gravity in its most extreme regimes.

The authors of this Letter thank the following organizations and programs : the Academy of Finland (projects 274477, 284495, 312496) ; the Advanced European Network of E-infrastructures for Astronomy with the SKA (AENEAS) project, supported by the European Commission Framework Programme Horizon 2020 Research and Innovation action under grant agreement 731016 ; the Alexander von Humboldt Stiftung ; the Black Hole Initiative at Harvard University, through a grant (60477) from the John Templeton Foundation ; the China Scholarship Council ; Comisión Nacional de Investigación Científica y Tecnológica (CONICYT, Chile, via PIA ACT172033, Fondecyt 1171506, BASAL AFB-170002, ALMA-conicyt 31140007) ; Consejo Nacional de Ciencia y Tecnológica (CONACYT, Mexico, projects 104497, 275201, 279006, 281692) ; the Delaney Family via the Delaney Family John A. Wheeler Chair at Perimeter Institute ; Dirección General de Asuntos del Personal Académico—Universidad Nacional Autónoma de México (DGAPA—UNAM, project IN112417) ; the European Research Council Synergy Grant ’’’BlackHoleCam : Imaging the Event Horizon of Black Holes’’’ (grant 610058) ; the Generalitat Valenciana postdoctoral grant APOSTD/2018/177 ; the Gordon and Betty Moore Foundation (grants GBMF-3561, GBMF-5278) ; the Istituto Nazionale di Fisica Nucleare (INFN) sezione di Napoli, iniziative specifiche TEONGRAV ; the International Max Planck Research School for Astronomy and Astrophysics at the Universities of Bonn and Cologne ; the Jansky Fellowship program of the National Radio Astronomy Observatory (NRAO) ; the Japanese Government (Monbukagakusho : MEXT) Scholarship ; the Japan Society for the Promotion of Science (JSPS) Grant-in-Aid for JSPS Research Fellowship (JP17J08829) ; JSPS Overseas Research Fellowships ; the Key Research Program of Frontier Sciences, Chinese Academy of Sciences (CAS, grants QYZDJ-SSW-SLH057, QYZDJ-SSW-SYS008) ; the Leverhulme Trust Early Career Research Fellowship ; the Max-Planck-Gesellschaft (MPG) ; the Max Planck Partner Group of the MPG and the CAS ; the MEXT/JSPS KAKENHI (grants 18KK0090, JP18K13594, JP18K03656, JP18H03721, 18K03709, 18H01245, 25120007) ; the MIT International Science and Technology Initiatives (MISTI) Funds ; the Ministry of Science and Technology (MOST) of Taiwan (105-2112-M-001-025-MY3, 106-2112-M-001-011, 106-2119-M-001-027, 107-2119-M-001-017, 107-2119-M-001-020, and 107-2119-M-110-005) ; the National Aeronautics and Space Administration (NASA, Fermi Guest Investigator grant 80NSSC17K0649) ; the National Institute of Natural Sciences (NINS) of Japan ; the National Key Research and Development Program of China (grant 2016YFA0400704, 2016YFA0400702) ; the National Science Foundation (NSF, grants AST-0096454, AST-0352953, AST-0521233, AST-0705062, AST-0905844, AST-0922984, AST-1126433, AST-1140030, DGE-1144085, AST-1207704, AST-1207730, AST-1207752, MRI-1228509, OPP-1248097, AST-1310896, AST-1312651, AST-1337663, AST-1440254, AST-1555365, AST-1715061, AST-1614868, AST-1615796, AST-1716327, OISE-1743747, AST-1816420) ; the Natural Science Foundation of China (grants 11573051, 11633006, 11650110427, 10625314, 11721303, 11725312, 11873028, 11873073, U1531245, 11473010) ; the Natural Sciences and Engineering Research Council of Canada (NSERC, including a Discovery Grant and the NSERC Alexander Graham Bell Canada Graduate Scholarships-Doctoral Program) ; the National Youth Thousand Talents Program of China ; the National Research Foundation of Korea (the Global PhD Fellowship Grant : grants NRF-2015H1A2A1033752, 2015-R1D1A1A01056807, the Korea Research Fellowship Program : NRF-2015H1D3A1066561) ; the Netherlands Organization for Scientific Research (NWO) VICI award (grant 639.043.513) and Spinoza Prize SPI 78-409 ; the New Scientific Frontiers with Precision Radio Interferometry Fellowship awarded by the South African Radio Astronomy Observatory (SARAO), which is a facility of the National Research Foundation (NRF), an agency of the Department of Science and Technology (DST) of South Africa ; the Onsala Space Observatory (OSO) national infrastructure, for the provisioning of its facilities/observational support (OSO receives funding through the Swedish Research Council under grant 2017-00648) the Perimeter Institute for Theoretical Physics (research at Perimeter Institute is supported by the Government of Canada through the Department of Innovation, Science and Economic Development Canada and by the Province of Ontario through the Ministry of Economic Development, Job Creation and Trade) ; the Russian Science Foundation (grant 17-12-01029) ; the Spanish Ministerio de Economía y Competitividad (grants AYA2015-63939-C2-1-P, AYA2016-80889-P) ; the State Agency for Research of the Spanish MCIU through the ’’’Center of Excellence Severo Ochoa’’’ award for the Instituto de Astrofísica de Andalucía (SEV-2017-0709) ; the Toray Science Foundation ; the US Department of Energy (USDOE) through the Los Alamos National Laboratory (operated by Triad National Security, LLC, for the National Nuclear Security Administration of the USDOE (Contract 89233218CNA000001)) ; the Italian Ministero dell’Istruzione Università e Ricerca through the grant Progetti Premiali 2012-iALMA (CUP C52I13000140001) ; ALMA North America Development Fund ; Chandra TM6-17006X and DD7-18089X ; a Sprows Family VURF Fellowship ; a NSERC Discovery Grant ; the FRQNT Nouveaux Chercheurs program ; CIFAR ; the NINS program of Promoting Research by Networking among Institutions(Grant Number 01421701) ; MEXT as a priority issue (Elucidation of the fundamental laws and evolution of the universe) to be tackled by using post-K Computer and JICFuS ; part of this work used XC50 at the Center for Computational Astrophysics, National Astronomical Observatory of Japan.This work used the Extreme Science and Engineering Discovery Environment (XSEDE), supported by NSF grant ACI-1548562, and CyVerse, supported by NSF grants DBI-0735191, DBI-1265383, and DBI-1743442. XSEDE Stampede2 resource at TACC was allocated through TG-AST170024 and TG-AST080026N. XSEDE JetStream resource at PTI and TACC was allocated through AST170028. The simulations were performed in part on the SuperMUC cluster at the LRZ in Garching, on the LOEWE cluster in CSC in Frankfurt, and on the HazelHen cluster at the HLRS in Stuttgart. This research was enabled in part by support provided by Compute Ontario (http://computeontario.ca), Calcul Quebec (http://www.calculquebec.ca) and Compute Canada (http://www.computecanada.ca). Results in this Letter are based in part on observations made by the Chandra X-ray Observatory (Observation IDs 20034, 20035, 19457, 19458, 00352, 02707, 03717) and the Nuclear Spectroscopic Telescope Array (NuSTAR ; Observation IDs 90202052002, 90202052004, 60201016002, 60466002002). The authors thank Belinda Wilkes, Fiona Harrison, Pat Slane, Joshua Wing, Karl Forster, and the Chandra and NuSTAR scheduling, data processing, and archive teams for making these challenging simultaneous observations possible. We thank the staff at the participating observatories, correlation centers, and institutions for their enthusiastic support. This Letter makes use of the following ALMA data : ADS/JAO.ALMA#2016.1.01154.V. ALMA is a partnership of the European Southern Observatory (ESO ; Europe, representing its member states), NSF, and National Institutes of Natural Sciences of Japan, together with National Research Council (Canada), Ministry of Science and Technology (MOST ; Taiwan), Academia Sinica Institute of Astronomy and Astrophysics (ASIAA ; Taiwan), and Korea Astronomy and Space Science Institute (KASI ; Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, Associated Universities, Inc. (AUI)/NRAO, and the National Astronomical Observatory of Japan (NAOJ). The NRAO is a facility of the NSF operated under cooperative agreement by AUI. APEX is a collaboration between the Max-Planck-Institut für Radioastronomie (Germany), ESO, and the Onsala Space Observatory (Sweden). The SMA is a joint project between the SAO and ASIAA and is funded by the Smithsonian Institution and the Academia Sinica. The JCMT is operated by the East Asian Observatory on behalf of the NAOJ, ASIAA, and KASI, as well as the Ministry of Finance of China, Chinese Academy of Sciences, and the National Key R&D Program (No. 2017YFA0402700) of China. Additional funding support for the JCMT is provided by the Science and Technologies Facility Council (UK) and participating universities in the UK and Canada. The LMT project is a joint effort of the Instituto Nacional de Astrofísica, Óptica y Electrónica (Mexico) and the University of Massachusetts at Amherst (USA). The IRAM 30-m telescope on Pico Veleta, Spain is operated by IRAM and supported by CNRS (Centre National de la Recherche Scientifique, France), MPG (Max-Planck-Gesellschaft, Germany) and IGN (Instituto Geográfico Nacional, Spain). The SMT is operated by the Arizona Radio Observatory, a part of the Steward Observatory of the University of Arizona, with financial support of operations from the State of Arizona and financial support for instrumentation development from the NSF. Partial SPT support is provided by the NSF Physics Frontier Center award (PHY-0114422) to the Kavli Institute of Cosmological Physics at the University of Chicago (USA), the Kavli Foundation, and the GBMF (GBMF-947). The SPT hydrogen maser was provided on loan from the GLT, courtesy of ASIAA. The SPT is supported by the National Science Foundation through grant PLR- 1248097. Partial support is also provided by the NSF Physics Frontier Center grant PHY-1125897 to the Kavli Institute of Cosmological Physics at the University of Chicago, the Kavli Foundation and the Gordon and Betty Moore Foundation grant GBMF-947.The EHTC has received generous donations of FPGA chips from Xilinx Inc., under the Xilinx University Program. The EHTC has benefited from technology shared under open-source license by the Collaboration for Astronomy Signal Processing and Electronics Research (CASPER). The EHT project is grateful to T4Science and Microsemi for their assistance with Hydrogen Masers. This research has made use of NASA’s Astrophysics Data System. We gratefully acknowledge the support provided by the extended staff of the ALMA, both from the inception of the ALMA Phasing Project through the observational campaigns of 2017 and 2018. We would like to thank A. Deller and W. Brisken for EHT-specific support with the use of DiFX. We acknowledge the significance that Maunakea, where the SMA and JCMT EHT stations are located, has for the indigenous Hawaiʻian people.
Appendix A : Table of Simulation Results

Below we provide a table of simulation results for models with a standard inclination of 17° between the approaching jet and the line of sight. In the notation of this Letter this corresponds to i = 17° for $a_* \lt 0$ or i = 163° for $a_* \geqslant 0$. The table shows models for which we were able to calculate $L_\mathrmbol$ and $L_\rmX

$. When $M$ is needed to calculate, e.g., $P_\mathrmjet$, we assume $M=6.2\times 10^9\,M_\odot $.

The first, third, and fourth columns in the table identify the model parameters : SANE or MAD based on dimensionless flux, $a_* $, and $R_\mathrmhigh$. Once these parameters are specified, an average value of $\dotM$ for the model, which is shown in last column, can be found from the requirement that the average flux density of 1.3 mm emission is 0.5 Jy (see Paper IV). This $\dotM$ is shown in units of the Eddington accretion rate $\dotM_\mathrmEdd=137M_\odot \,\mathrmyr^-1$. The measured average dimensionless magnetic flux phgr is shown in the second column. Notice that phgr is determined solely from the GRMHD simulation and is independent of the mass scaling $\mathscrM$ and the mass $M$ used to fix the flux density. It is also independent of the electron thermodynamics ($R_\mathrmhigh$).

The fifth column shows the radiative efficiency, which is the bolometric luminosity $L_\mathrmbol$ over $\dotMc^2$. Here $L_\mathrmbol$ was found from a relativistic Monte Carlo radiative transport model that includes synchrotron emission, Compton scattering (all orders), and bremsstrahlung. The Monte Carlo calculation makes no approximations in treating the Compton scattering (see Dolence et al. 2009). Bremsstrahlung is negligible in all models.

The sixth column shows predicted X-ray luminosity $L_\rmX

$ in the $2\mbox10\,\mathrmkeV$ band. This was calculated using the same relativistic Monte Carlo radiative transport model as for $L_\mathrmbol$. There are three numbers in this column : the average $\langle L_\rmX

\rangle $ (left) of the 20 sample spectra used in the calculation, and a maximum and minimum value. The maximum and minimum are obtained by taking the standard deviation $\sigma (\mathrmlog_10L_\rmX

)$ and setting the maximum (minimum) to $10^+2\sigma \langle L_\rmX

\rangle $ ($10^-2\sigma \langle L_\rmX

\rangle $).

The seventh column shows the jet power
Equation (10)

The integral is evaluated at $r=40\,r_\rmg

$ for SANE models and $r=100\,r_\rmg

$ for MAD models. These radii were chosen because they are close to the outer boundary of the computational domain. Here $\rm\Delta t$ is the duration of the time-average, $-T^r_t$ is a component of the stress-energy tensor representing outward radial energy flux, g is the determinant of the (covariant) metric, ρ is the rest-mass density, and ur is the radial component of the four-velocity. Here we use Kerr–Schild $t,r,\theta ,\phi $ for clarity ; in practice, the integral is evaluated in simulation coordinates. The quantity in parentheses is the outward energy flux with the rest-mass energy flux subtracted off. The θ integral is done after time averaging and azimuthal integration over the region where
Equation (11)

Here βγ would be the radial four-velocity as $r\to \infty $ if the flow were steady and all internal magnetic and internal energy were converted to kinetic energy. In Table 3 we use $(\beta \gamma )_\mathrmcut^2=1$ to define the jet. This is equivalent to restricting the jet to regions where the total energy per unit rest-mass (including the rest-mass energy) exceeds $\sqrt5c^2\simeq 2.2c^2$.

The ninth column shows the total outflow power $P_\mathrmout$, defined using the same integral as in Equation (10), but with the θ integral carried out over the entire region around the poles where there is steady outflow (and $\theta \lt 1$, although the result is insensitive to this condition). $P_\mathrmout$ thus includes both the narrow, fast, relativistic jet and any wide-angle, slow, or nonrelativistic outflow. It is the maximal $P_\mathrmjet$ under any definition of jet power.

Finally, the tenth column shows the ratio of the electromagnetic to total energy flux in the jet. In most cases this number is close to 1 ; i.e., the jet is Poynting dominated. This measurement is sensitive to the numerical treatment of low-density regions in the jet where the jet can be artificially loaded with plasma by numerical ’’’floors’’’ in the GRMHD evolution. More accurate treatment of the funnel would raise values in this column.

Our choice of $(\beta \gamma )_\mathrmcut^2$, and therefore $P_\mathrmjet$, is somewhat arbitrary. To probe the sensitivity of $P_\mathrmjet$ to $(\beta \gamma )_\mathrmcut^2$, Figure 10 shows the ratio $P_\mathrmjet$/$P_\mathrmout$ (which is determined by the GRMHD model and is thus independent of the electron thermodynamics, i.e., $R_\mathrmhigh$) as a function o

Figure 10. Ratio $P_\mathrmjet$/$P_\mathrmout$ as a function of the outflow velocity cutoff parameter $\beta \gamma _\mathrmcut$. Evidently, as the cut is decreased, so that the maximum asymptotic speed of the jet flow is decreased, an increasing fraction of $P_\mathrmout$ is classified as $P_\mathrmjet$. Our nominal cutoff is $\beta \gamma =1$, which corresponds to $\beta \equiv v^r/c=1/\sqrt2$. Using this definition, $P_\mathrmjet$ for $a_* =0$ models is small because the energy flux in the relativistic outflow is small.

The eighth and tenth columns show the jet and outflow efficiency. This is determined by the GRMHD evolution, i.e., it is independent of electron thermodynamics ($R_\mathrmhigh$). It is $\gt 0.1$ only for MAD models with $a_* \geqslant 0.5$.

The eleventh column shows the fraction of $P_\mathrmjet$ in Poynting flux. This fraction is large for all models, and meaningless for the $a_* =0$ models, which have $P_\mathrmjet$ that is so small that it is difficult to measure accurately.

The problem of defining $P_\mathrmjet$ and $P_\mathrmout$ has been discussed extensively in the literature (e.g., Narayan et al. 2012 ; Yuan et al. 2015 ; Mościbrodzka et al. 2016), where alternative definitions of unbound regions and of the jet have been used, some based on a fluid Bernoulli parameter $B_e\,\equiv -u_t(\rho +u+p)/\rho -1$, while others use μ (the ratio of energy flux to rest mass flux), which is directly related to our βγ.
Appendix B : Image Decomposition into Midplane, Nearside, and Farside Components

In Section 3.3 we presented representative images from the Image Library spanning a broad range of values in both $a_* $ and $R_\mathrmhigh$. It was noted that for SANEs with low values of $R_\mathrmhigh$ the emission is concentrated more in the midplane, whereas for larger values of $R_\mathrmhigh$ this emission is concentrated in the funnel wall. In particular, Figure 4 presented temporal- and azimuthal-averaged images of the point of origin of photons comprising images from $a_* =0.94$ MAD and SANE simulations with $R_\mathrmhigh=10$ and 160.

Figure 11 presents the decomposition of the four images in Figure 4 into components that we refer to as : midplane (material within 32fdg7 of the midplane), nearside (material within 1 radian, or 57fdg3, of the polar axis nearest to the observer), and farside (material within 1 radian of the polar axis furthest from the observer). From inspection of the first three models (rows) in Figure 11, the ratio of nearside to farside flux in the simulations is small (compared to the midplane) and of order unity and the midplane emission is dominant, as in Figure 4.

Figure 11. Decomposition of time-averaged 1.3 mm images from Figure 4 into midplane, nearside, and farside components (MAD and SANE models with $a_* =0.94$). Each model (row) of the figure corresponds to a simulation in Figure 4. The percentage of the total image flux from each component is indicated in the bottom right of each panel. The color scale is logarithmic and spans three decades in total flux with respect to the total image from each model, chosen in order to emphasize both nearside and farside components, which are nearly invisible when shown in a linear scale. The field of view is $80\,\mu \mathrmas$.

However, for the SANE, $R_\mathrmhigh=160$ model the farside emission contributes a flux that exceeds that produced from the midplane, and is significantly brighter than the nearside emission. This is in agreement with what is seen in the bottom-right panel of Figures 4, and can be understood to arise from the SANE model possessing an optically thin disk and bright funnel wall in the $R_\mathrmhigh=160$ case, compared to SANE, $R_\mathrmhigh=10$, as also seen in Figures 2 and 3. Due to the reduced opacity along the line of sight in this case, mm photons can pass through both the intervening nearside material and the midplane without significant attenuation, before reaching the photospheric boundary in the farside component (where τ 1), where they originate. The image decomposition and its application to M87’s image structure will be explored further in Z. Younsi et al. 2019a (in preparation).

M87-6 L’ombre et la masse du trou noir central
Footnotes

110

phgr is determined by the outcome of the simulation and cannot be trivially predicted from the initial conditions, but by repeated experiment it is possible to manipulate the size of the initial torus and strength and geometry of the initial field to produce a target phgr.

111

In Heaviside units, where a factor of $\sqrt4\pi $ is absorbed into the definition of B, $\phi _\max \simeq 15$. In the Gaussian units used in some earlier papers, $\phi _\max \simeq 50$.

112

For a black hole accreting at the Eddington rate, the ratio of the accreting mass onto a black hole mass is $\sim 10^-22(M/M_\odot ) ;$ in our models mass accretion rate is far below the Eddington rate.

113

In particular the distribution of best-fit M/D, which is defined in Section 4, have mean and standard deviation of $M/D=3.552\pm 0.605\,\mu \mathrmas$ when the images are made with an input $M/D=3.62\,\mu \mathrmas$, and $3.564\pm 0.537\,\mu \mathrmas$ when the images are made with an input $M/D=2.01\,\mu \mathrmas$. We have also checked images made with an input 1.3 mm flux ranging from 0.1 to 1.5 Jy and find relative changes in M/D and $\mathrmPA$ of less than 1%.

114

In particular, doubling the image resolution changes the mean best-fit M/D by 7 nano-arcsec, and the best-fit $\mathrmPA$ by 0fdg3.

115

In GRMHD models the jet core is effectively empty and the density is set by numerical ’’’floors.’’’ In our radiative transfer calculations emission from regions with $B^2/\rho \gt 1$ is explicitly set to zero.

116

Paper I and Paper IV focus instead on the April 11 data set.

117

The thin disk radiative efficiency is 0.038 for $a_* =-1$, 0.057 for $a_* =0$, and 0.42 for $a_* =1$. See Equations (2.12) and (2.21) of Bardeen et al. (1972) ; the efficiency is $1-E/\mu _p$, where $\mu _p$ is the rest mass of the particle. The rejected model list is identical if instead one simply rejects all models with $\epsilon \gt 0.2$.

118

The compact mm flux density could be a factor of 2 larger than our assumed 0.5 Jy. That would raise $P_\mathrmjet$ by slightly less than a factor of 2.

119

The total energy flux inside a cone of opening angle $\theta _0$ is proportional to $\sin ^4\theta _0$ in the Blandford & Znajek (1977) monopole model if the field strength is fixed, and $\sin ^2\theta _0$ if the magnetic flux is fixed.

120

The width of the ring, the central flux depression, and a quantitative discussion of the black hole shadow can be found in Paper VI.

[bruno sanchiz]