Draft version April 11, 2019
Typeset using L
A
T
E
X default style in AASTeX62
The Event Horizon General Relativistic Magnetohydrodynamic Code Comparison Project
Oliver Porth,
1, 2
Koushik Chatterjee,
1
Ramesh Narayan,
3, 4
Charles F. Gammie,
5
Yosuke Mizuno,
2
Peter Anninos,
6
John G. Baker,
7, 8
Matteo Bugli,
9
Chi-kwan Chan,
10, 11
Jordy Davelaar,
12
Luca Del Zanna,
13, 14
Zachariah B. Etienne,
15, 16
P. Chris Fragile,
17
Bernard J. Kelly,
18, 19, 7
Matthew Liska,
1
Sera Markoff,
20
Jonathan C. McKinney,
21
Bhupendra Mishra,
22
Scott C. Noble,
23, 7
H
́
ector Olivares,
2
Ben Prather,
24
Luciano Rezzolla,
2
Benjamin R. Ryan,
25, 26
James M. Stone,
27
Niccol
`
o Tomei,
13, 14
Christopher J. White,
28
and
Ziri Younsi
29, 2
—
Kazunori Akiyama,
30, 31, 32, 4
Antxon Alberdi,
33
Walter Alef,
34
Keiichi Asada,
35
Rebecca Azulay,
36, 37, 34
Anne-Kathrin Baczko,
34
David Ball,
10
Mislav Balokovi
́
c,
3, 4
John Barrett,
31
Dan Bintley,
38
Lindy Blackburn,
3, 4
Wilfred Boland,
39
Katherine L. Bouman,
3, 4, 40
Geoffrey C. Bower,
41
Michael Bremer,
42
Christiaan D. Brinkerink,
12
Roger Brissenden,
3, 4
Silke Britzen,
34
Avery E. Broderick,
43, 44, 45
Dominique Broguiere,
42
Thomas Bronzwaer,
12
Do-Young Byun,
46, 47
John E. Carlstrom,
48, 49, 50, 51
Andrew Chael,
3, 4
Shami Chatterjee,
52
Ming-Tang Chen,
41
Yongjun Chen (
陈
永
军
),
53, 54
Ilje Cho,
46, 47
Pierre Christian,
10, 3
John E. Conway,
55
James M. Cordes,
52
Geoffrey, B. Crew,
31
Yuzhu Cui,
56, 57
Mariafelicia De Laurentis,
58, 2, 59
Roger Deane,
60
Jessica Dempsey,
38
Gregory Desvignes,
34
Jason Dexter,
61
Sheperd S. Doeleman,
3, 4
Ralph P. Eatough,
34
Heino Falcke,
12
Vincent L. Fish,
31
Ed Fomalont,
30
Raquel Fraga-Encinas,
12
Bill Freeman,
62, 63
Per Friberg,
38
Christian M. Fromm,
2
Jos
́
e L. G
́
omez,
33
Peter Galison,
64, 65, 4
Charles F. Gammie,
5
Roberto Garc
́
ıa,
66
Olivier Gentaz,
66
Boris Georgiev,
44
Ciriaco Goddi,
12, 67
Roman Gold,
2
Minfeng Gu (
顾
敏
峰
),
53, 68
Mark Gurwell,
3
Kazuhiro Hada,
56, 57
Michael H. Hecht,
31
Ronald Hesper,
69
Luis C. Ho (
何
子
山
),
70, 71
Paul Ho,
35
Mareki Honma,
56, 57
Chih-Wei L. Huang,
35
Lei Huang (
黄
磊
),
53, 68
David H. Hughes,
72
Shiro Ikeda,
73, 32, 74, 75
Makoto Inoue,
35
Sara Issaoun,
12
David J. James,
3, 4
Buell T. Jannuzi,
10
Michael Janssen,
12
Britton Jeter,
44
Wu Jiang (
江
悟
),
53
Michael D. Johnson,
3, 4
Svetlana Jorstad,
76, 77
Taehyun Jung,
46, 47
Mansour Karami,
78, 79
Ramesh Karuppusamy,
34
Tomohisa Kawashima,
32
Garrett K. Keating,
3
Mark Kettenis,
80
Jae-Young Kim,
34
Junhan Kim,
10
Jongsoo Kim,
46
Motoki Kino,
32, 81
Jun Yi Koay,
35
Patrick, M. Koch,
35
Shoko Koyama,
35
Michael Kramer,
34
Carsten Kramer,
42
Thomas P. Krichbaum,
34
Cheng-Yu Kuo,
82
Tod R. Lauer,
83
Sang-Sung Lee,
46
Yan-Rong Li (
李
彦
荣
),
84
Zhiyuan Li (
李
志
远
),
85, 86
Michael Lindqvist,
55
Kuo Liu,
34
Elisabetta Liuzzo,
87
Wen-Ping Lo,
88, 35
Andrei P. Lobanov,
34
Laurent Loinard,
89, 90
Colin Lonsdale,
31
Ru-Sen Lu (
路
如
森
),
53, 34
Nicholas R. MacDonald,
34
Jirong Mao (
毛
基
荣
),
91, 92, 93
Daniel P. Marrone,
10
Alan P. Marscher,
94
Iv
́
an Mart
́
ı-Vidal,
55, 95
Satoki Matsushita,
35
Lynn D. Matthews,
31
Lia Medeiros,
10, 96
Karl M. Menten,
34
Izumi Mizuno,
38
James M. Moran,
3, 4
Kotaro Moriyama,
56
Monika Moscibrodzka,
12
Cornelia Mu
̈
ller,
12, 34
Hiroshi Nagai,
97, 98
Neil M. Nagar,
99
Masanori Nakamura,
35
Gopal Narayanan,
100
Iniyan Natarajan,
101
Roberto Neri,
42
Chunchong Ni,
44
Aristeidis Noutsos,
34
Hiroki Okino,
56, 102
Gisela N. Ortiz-Le
́
on,
34
Tomoaki Oyama,
56
Feryal
̈
Ozel,
10
Daniel C. M. Palumbo,
3, 4
Nimesh Patel,
3
Ue-Li Pen,
103, 104, 105, 106
Dominic W. Pesce,
3, 4
Vincent Pi
́
etu,
42
Richard Plambeck,
107
Aleksandar PopStefanija,
100
Jorge A. Preciado-L
́
opez,
43
Dimitrios Psaltis,
10
Hung-Yi Pu,
43
Venkatessh Ramakrishnan,
99
Ramprasad Rao,
41
Mark G. Rawlings,
38
Alexander W. Raymond,
3, 4
Bart Ripperda,
2
Freek Roelofs,
12
Alan Rogers,
31
Eduardo Ros,
34
Mel Rose,
10
Arash Roshanineshat,
10
Helge Rottmann,
34
Alan L. Roy,
34
Chet Ruszczyk,
31
Kazi L.J. Rygl,
87
Salvador S
́
anchez,
108
David S
́
anchez-Arguelles,
109, 72
Mahito Sasada,
110, 56
Tuomas Savolainen,
111, 112, 34
F. Peter Schloerb,
100
Karl-Friedrich Schuster,
42
Lijing Shao,
71, 34
Zhiqiang Shen (
沈
志强
),
53, 54
Des Small,
80
Bong Won Sohn,
46, 47, 113
Jason SooHoo,
31
Fumie Tazaki,
56
Paul Tiede,
43, 44
Remo P.J. Tilanus,
67, 12, 114
Michael Titus,
31
Kenji Toma,
115, 116
Pablo Torne,
108, 34
Tyler Trent,
10
Sascha Trippe,
117
Shuichiro Tsuda,
118
Ilse van Bemmel,
80
Huib Jan van Langevelde,
80, 119
Daniel R. van Rossum,
12
Jan Wagner,
34
John Wardle,
120
Jonathan Weintroub,
3, 4
Norbert Wex,
34
Robert Wharton,
34
Maciek Wielgus,
3, 4
George N. Wong,
24
Qingwen Wu (
吴
庆
文
),
121
Ken Young,
3
Andr
́
e Young,
12
Feng Yuan (
袁
峰
),
53, 68, 122
Ye-Fei Yuan (
袁
业
飞
),
123
J. Anton Zensus,
34
Guangyao Zhao,
46
Shan-Shan Zhao,
12, 85
and Ziyan Zhu
65
(The Event Horizon Telescope Collaboration)
1
Anton Pannekoek Institute for Astronomy, University of Amsterdam, Science Park 904, 1098 XH, Amsterdam, The Netherlands
2
Institut f ̈ur Theoretische Physik, Goethe-Universit ̈at Frankfurt, Max-von-Laue-Straße 1, D-60438 Frankfurt am Main, Germany
Corresponding author: Oliver Porth
o.porth@uva.nl
arXiv:1904.04923v1 [astro-ph.HE] 9 Apr 2019
2
GRMHD community and the EHTC
3
Center for Astrophysics
|
Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA
4
Black Hole Initiative at Harvard University, 20 Garden Street, Cambridge, MA 02138, USA
5
Department of Astronomy; Department of Physics; University of Illinois, Urbana, IL 61801 USA
6
Lawrence Livermore National Laboratory, Livermore, CA 94550, USA
7
Gravitational Astrophysics Laboratory, NASA Goddard Space Flight Center, Greenbelt, MD 20771, USA
8
Joint Space-Science Institute, University of Maryland, College Park, MD 20742, USA
9
IRFU/D ́epartement d’Astrophysique, Laboratoire AIM, CEA/DRF-CNRS-Universit ́e Paris Diderot, CEA-Saclay F-91191, France
10
Steward Observatory and Department of Astronomy, University of Arizona, 933 N. Cherry Ave., Tucson, AZ 85721, USA
11
Data Science Institute, University of Arizona, 1230 N. Cherry Ave., Tucson, AZ 85721, USA
12
Department of Astrophysics, Institute for Mathematics, Astrophysics and Particle Physics (IMAPP), Radboud University, P.O. Box
9010, 6500 GL Nijmegen, The Netherlands
13
Dipartimento di Fisica e Astronomia, Universit`a di Firenze e INFN – Sez. di Firenze, via G. Sansone 1, I-50019 Sesto F.no, Italy
14
INAF, Osservatorio Astrofisico di Arcetri, Largo E. Fermi 5, I-50125 Firenze, Italy
15
Department of Mathematics, West Virginia University, Morgantown, WV 26506, USA
16
Center for Gravitational Waves and Cosmology, West Virginia University, Chestnut Ridge Research Building, Morgantown, WV
26505, USA
17
Department of Physics & Astronomy, College of Charleston, 66 George St., Charleston, SC 29424, USA
18
Department of Physics, University of Maryland, Baltimore County, Baltimore, MD 21250, USA
19
CRESST, NASA Goddard Space Flight Center, Greenbelt, MD 20771, USA
20
Anton Pannekoek Institute for Astronomy & GRAPPA, University of Amsterdam, Postbus 94249, 1090GE Amsterdam, The
Netherlands
21
H2O.ai, 2307 Leghorn St., Mountain View, CA 94043
22
JILA, University of Colorado and National Institute of Standards and Technology, 440 UCB, Boulder, CO 80309-0440, USA
23
Department of Physics and Engineering Physics, The University of Tulsa, Tulsa, OK 74104, USA
24
Department of Physics, University of Illinois, 1110 West Green St, Urbana, IL 61801, USA
25
CCS-2, Los Alamos National Laboratory, P.O. Box 1663, Los Alamos, NM 87545, US
26
Center for Theoretical Astrophysics, Los Alamos National Laboratory, Los Alamos, NM, 87545, USA
27
Department of Astrophysical Sciences, Princeton University, Princeton, NJ 08544, USA
28
Kavli Institute for Theoretical Physics, University of California Santa Barbara, Kohn Hall, Santa Barbara, CA 93107, USA
29
Mullard Space Science Laboratory, University College London, Holmbury St. Mary, Dorking, Surrey, RH5 6NT, United Kingdom
30
National Radio Astronomy Observatory, 520 Edgemont Rd, Charlottesville, VA 22903, USA
31
Massachusetts Institute of Technology Haystack Observatory, 99 Millstone Road, Westford, MA 01886, USA
32
National Astronomical Observatory of Japan, 2-21-1 Osawa, Mitaka, Tokyo 181-8588, Japan
33
Instituto de Astrof ́ısica de Andaluc ́ıa - CSIC, Glorieta de la Astronom ́ıa s/n, E-18008 Granada, Spain
34
Max-Planck-Institut f ̈ur Radioastronomie, Auf dem H ̈ugel 69, D-53121 Bonn, Germany
35
Institute of Astronomy and Astrophysics, Academia Sinica, 11F of Astronomy-Mathematics Building, AS/NTU No. 1, Sec. 4,
Roosevelt Rd, Taipei 10617, Taiwan, R.O.C.
36
Departament d’Astronomia i Astrof ́ısica, Universitat de Val`encia, C. Dr. Moliner 50, E-46100 Burjassot, Val`encia, Spain
37
Observatori Astron`omic, Universitat de Val`encia, C. Catedr ́atico Jos ́e Beltr ́an 2, E-46980 Paterna, Val`encia, Spain
38
East Asian Observatory, 660 N. A’ohoku Pl., Hilo, HI 96720, USA
39
Nederlandse Onderzoekschool voor Astronomie (NOVA), PO Box 9513, 2300 RA Leiden, the Netherlands, Niels Bohrweg 2, 2333 CA
Leiden, the Netherlands
40
California Institute of Technology, 1200 East California Boulevard, Pasadena, CA 91125, USA
41
Institute of Astronomy and Astrophysics, Academia Sinica, 645 N. A’ohoku Place, Hilo, HI 96720, USA
42
Institut de Radioastronomie Millim ́etrique, 300 rue de la Piscine, 38406 Saint Martin d’H`eres, France
43
Perimeter Institute for Theoretical Physics, 31 Caroline Street North, Waterloo, ON, N2L 2Y5, Canada
44
Department of Physics and Astronomy, University of Waterloo, 200 University Avenue West, Waterloo, ON, N2L 3G1, Canada
45
Waterloo Centre for Astrophysics, University of Waterloo, Waterloo, ON N2L 3G1 Canada
46
Korea Astronomy and Space Science Institute, Daedeok-daero 776, Yuseong-gu, Daejeon 34055, Republic of Korea
47
University of Science and Technology, Gajeong-ro 217, Yuseong-gu, Daejeon 34113, Republic of Korea
48
Kavli Institute for Cosmological Physics, University of Chicago, Chicago, IL 60637, USA
49
Department of Astronomy and Astrophysics, University of Chicago, Chicago, IL 60637, USA
50
Department of Physics, University of Chicago, Chicago, IL 60637, USA
51
Enrico Fermi Institute, University of Chicago, Chicago, IL 60637, USA
52
Cornell Center for Astrophysics and Planetary Science, Cornell University, Ithaca, NY 14853, USA
53
Shanghai Astronomical Observatory, Chinese Academy of Sciences, 80 Nandan Road, Shanghai 200030, China
54
Key Laboratory of Radio Astronomy, Chinese Academy of Sciences, Nanjing 210008, China
Event Horizon Code Comparison
3
55
Department of Space, Earth and Environment, Chalmers University of Technology, Onsala Space Observatory, SE-439 92 Onsala,
Sweden
56
Mizusawa VLBI Observatory, National Astronomical Observatory of Japan, 2-12 Hoshigaoka, Mizusawa, Oshu, Iwate 023-0861, Japan
57
Department of Astronomical Science, The Graduate University for Advanced Studies (SOKENDAI), 2-21-1 Osawa, Mitaka, Tokyo
181-8588, Japan
58
Dipartimento di Fisica ”E. Pancini”, Universit ́a di Napoli ”Federico II”, Compl. Univ. di Monte S. Angelo, Edificio G, Via Cinthia,
I-80126, Napoli, Italy
59
INFN Sez. di Napoli, Compl. Univ. di Monte S. Angelo, Edificio G, Via Cinthia, I-80126, Napoli, Italy
60
Department of Physics, University of Pretoria, Lynnwood Road, Hatfield, Pretoria 0083, South Africa; Centre for Radio Astronomy
Techniques and Technologies, Department of Physics and Electronics, Rhodes University, Grahamstown 6140, South Africa
61
Max-Planck-Institut f ̈ur Extraterrestrische Physik, Giessenbachstr. 1, D-85748 Garching, Germany
62
Department of Electrical Engineering and Computer Science, Massachusetts Institute of Technology, 32-D476, 77 Massachussetts Ave.,
Cambridge, MA 02142, USA
63
Google Research, 355 Main St., Cambridge, MA 02142, USA
64
Department of History of Science, Harvard University, Cambridge, MA 02138, USA
65
Department of Physics, Harvard University, Cambridge, MA 02138, USA
66
Institut de Radioastronomie Millim ́etrique, 300 rue de la Piscine, 38406 Saint Martin d’H`eres, France
67
Leiden Observatory - Allegro, Leiden University, P.O. Box 9513, 2300 RA Leiden, The Netherlands
68
Key Laboratory for Research in Galaxies and Cosmology, Chinese Academy of Sciences, Shanghai 200030, China
69
NOVA Sub-mm Instrumentation Group, Kapteyn Astronomical Institute, University of Groningen, Landleven 12, 9747 AD Groningen,
The Netherlands
70
Department of Astronomy, School of Physics, Peking University, Beijing 100871, China
71
Kavli Institute for Astronomy and Astrophysics, Peking University, Beijing 100871, China
72
Instituto Nacional de Astrof ́ısica,
́
Optica y Electr ́onica. Apartado Postal 51 y 216, 72000. Puebla Pue., M ́exico
73
The Institute of Statistical Mathematics, 10-3 Midori-cho, Tachikawa, Tokyo, 190-8562, Japan
74
Department of Statistical Science, The Graduate University for Advanced Studies (SOKENDAI), 10-3 Midori-cho, Tachikawa, Tokyo
190-8562, Japan
75
Kavli Institute for the Physics and Mathematics of the Universe, The University of Tokyo, 5-1-5 Kashiwanoha, Kashiwa, 277-8583,
Japan
76
Institute for Astrophysical Research, Boston University, 725 Commonwealth Ave., Boston, MA 02215
77
Astronomical Institute, St.Petersburg University, Universitetskij pr., 28, Petrodvorets,198504 St.Petersburg, Russia
78
Perimeter Institute for Theoretical Physics, 31 Caroline Street North, Waterloo, Ontario N2L 2Y5, Canada
79
University of Waterloo, 200 University Avenue West, Waterloo, Ontario N2L 3G1, Canada
80
Joint Institute for VLBI ERIC (JIVE), Oude Hoogeveensedijk 4, 7991 PD Dwingeloo, The Netherlands
81
Kogakuin University of Technology & Engineering, Academic Support Center, 2665-1 Nakano, Hachioji, Tokyo 192-0015, Japan
82
Physics Department, National Sun Yat-Sen University, No. 70, Lien-Hai Rd, Kaosiung City 80424, Taiwan, R.O.C
83
National Optical Astronomy Observatory, 950 North Cherry Ave., Tucson, AZ 85719, USA
84
Key Laboratory for Particle Astrophysics, Institute of High Energy Physics, Chinese Academy of Sciences, 19B Yuquan Road,
Shijingshan District, Beijing, China
85
School of Astronomy and Space Science, Nanjing University, Nanjing 210023, China
86
Key Laboratory of Modern Astronomy and Astrophysics, Nanjing University, Nanjing 210023, China
87
Italian ALMA Regional Centre, INAF-Istituto di Radioastronomia, Via P. Gobetti 101, 40129 Bologna, Italy
88
Department of Physics, National Taiwan University, No.1, Sect.4, Roosevelt Rd., Taipei 10617, Taiwan, R.O.C
89
Instituto de Radioastronom ́ıa y Astrof ́ısica, Universidad Nacional Aut ́onoma de M ́exico, Morelia 58089, M ́exico
90
Instituto de Astronom
́
Ia, Universidad Nacional Aut ́onoma de M ́exico, CdMx 04510, M ́exico
91
Yunnan Observatories, Chinese Academy of Sciences, 650011 Kunming, Yunnan Province, China
92
Center for Astronomical Mega-Science, Chinese Academy of Sciences, 20A Datun Road, Chaoyang District, Beijing, 100012, China
93
Key Laboratory for the Structure and Evolution of Celestial Objects, Chinese Academy of Sciences, 650011 Kunming, China
94
Institute for Astrophysical Research, Boston University, 725 Commonwealth Ave., Boston, MA 02215, USA
95
Centro Astron ́omico de Yebes (IGN), Apartado 148, 19180 Yebes, Spain
96
Department of Physics, Broida Hall, University of California Santa Barbara, Santa Barbara, CA 93106, USA
97
National Astronomical Observatory of Japan, Osawa 2-21-1, Mitaka, Tokyo 181-8588, Japan
98
Department of Astronomical Science, The Graduate University for Advanced Studies (SOKENDAI), Osawa 2-21-1, Mitaka, Tokyo
181-8588, Japan
99
Astronomy Department, Universidad de Concepci ́on, Casilla 160-C, Concepci ́on, Chile
100
Department of Astronomy, University of Massachusetts, 01003, Amherst, MA, USA
101
Centre for Radio Astronomy Techniques and Technologies, Department of Physics and Electronics, Rhodes University, Grahamstown
6140, South Africa
4
GRMHD community and the EHTC
102
Department of Astronomy, Graduate School of Science, The University of Tokyo, 7-3-1 Hongo, Bunkyo-ku, Tokyo 113-0033, Japan
103
Canadian Institute for Theoretical Astrophysics, University of Toronto, 60 St. George Street, Toronto, ON M5S 3H8, Canada
104
Dunlap Institute for Astronomy and Astrophysics, University of Toronto, 50 St. George Street, Toronto, ON M5S 3H4, Canada
105
Canadian Institute for Advanced Research, 180 Dundas St West, Toronto, ON M5G 1Z8, Canada
106
Perimeter Institute for Theoretical Physics, 31 Caroline Street North, Waterloo, ON N2L 2Y5, Canada
107
Radio Astronomy Laboratory, University of California, Berkeley, CA 94720
108
Instituto de Radioastronom ́ıa Milim ́etrica, IRAM, Avenida Divina Pastora 7, Local 20, 18012, Granada, Spain
109
Consejo Nacional de Ciencia y Tecnolog ́ıa, Av. Insurgentes Sur 1582, 03940, Ciudad de M ́exico, M ́exico
110
Hiroshima Astrophysical Science Center, Hiroshima University, 1-3-1 Kagamiyama, Higashi-Hiroshima, Hiroshima 739-8526, Japan
111
Aalto University Department of Electronics and Nanoengineering, PL 15500, 00076 Aalto, Finland
112
Aalto University Mets ̈ahovi Radio Observatory, Mets ̈ahovintie 114, 02540 Kylm ̈al ̈a, Finland
113
Department of Astronomy, Yonsei University, Yonsei-ro 50, Seodaemun-gu, 03722 Seoul, Republic of Korea
114
Netherlands Organisation for Scientific Research (NWO), Postbus 93138, 2509 AC Den Haag , The Netherlands
115
Frontier Research Institute for Interdisciplinary Sciences, Tohoku University, Sendai 980-8578, Japan
116
Astronomical Institute, Tohoku University, Sendai 980-8578, Japan
117
Department of Physics and Astronomy, Seoul National University, Gwanak-gu, Seoul 08826, Republic of Korea
118
Mizusawa VLBI Observatory, National Astronomical Observatory of Japan, Hoshigaoka 2-12, Mizusawa-ku, Oshu-shi, Iwate 023-0861,
Japan
119
Leiden Observatory, Leiden University, Postbus 2300, 9513 RA Leiden, The Netherlands
120
Physics Department, Brandeis University, 415 South Street, Waltham , MA 02453
121
School of Physics, Huazhong University of Science and Technology, Wuhan, Hubei, 430074, China
122
School of Astronomy and Space Sciences, University of Chinese Academy of Sciences, No. 19A Yuquan Road, Beijing 100049, China
123
Astronomy Department, University of Science and Technology of China, Hefei 230026, China
(Received April 9, 2019)
Submitted to ApJS
ABSTRACT
Recent developments in compact object astrophysics, especially the discovery of merging neutron
stars by LIGO, the imaging of the black hole in M87 by the Event Horizon Telescope (EHT) and
high precision astrometry of the Galactic Center at close to the event horizon scale by the GRAVITY
experiment motivate the development of numerical source models that solve the equations of general
relativistic magnetohydrodynamics (GRMHD). Here we compare GRMHD solutions for the evolution
of a magnetized accretion flow where turbulence is promoted by the magnetorotational instability
from a set of nine GRMHD codes:
Athena++
,
BHAC
,
Cosmos++
,
ECHO
,
H-AMR
,
iharm3D
,
HARM-Noble
,
IllinoisGRMHD
and
KORAL
. Agreement between the codes improves as resolution increases, as measured
by a consistently applied, specially developed set of code performance metrics. We conclude that the
community of GRMHD codes is mature, capable, and consistent on these test problems.
Keywords:
numerical methods - magnetohydrodynamics (MHD) - magnetic fields - general relativity
- black hole physics
1.
INTRODUCTION
Fully general relativistic models of astrophysical sources are in high demand, not only since the discovery of gravita-
tional waves emitted from merging stellar mass black holes (Abbott et al. 2016). The need for an accurate description
of the interplay between strong gravity, matter and electromagnetic fields is further highlighted by the recent detection
of electromagnetic counterpart radiation to the coalescence of a neutron star binary (Abbott et al. 2017). Our own
effort is motivated by the Event Horizon Telescope (EHT) project, which allows direct imaging of hot, luminous plasma
near a black hole event horizon (Doeleman et al. 2008; Goddi et al. 2017; EHT Collaboration 2019a). The main targets
of the EHT are the black hole at the center of the Milky Way (also known by the name of the associated radio source,
Sgr A*, e.g. Lu et al. (2018)) and the black hole at the center of the galaxy M87 with associated central radio source
M87* (Doeleman et al. 2012; Akiyama et al. 2015). In order to extract information on the dynamics of the plasma
that lies within a few
GM/c
2
of the event horizon (
M
≡
black hole mass) as well as information about the black hole’s
Event Horizon Code Comparison
5
gravitational field, it is necessary to develop models of the accretion flow, associated winds and relativistic jets, and
the emission properties in each of the components.
Earlier semi-analytic works (Narayan & Yi 1995; Narayan et al. 1998; Yuan et al. 2002) have provided with the general
parameter regime of the galactic center by exploiting spectral information. For example, Mahadevan & Quataert (1997)
demonstrated that the electrons and ions are only weakly collisionally coupled and unlikely in thermal equilibrium.
Also key parameters like the accretion rate are typically estimated based on simple one-dimensional models (Marrone
et al. 2007). They have solidified the notion that the accretion rate in Sgr A* is far below the Eddington limit
̇
M
Edd
=
L
Edd
/
(0
.
1
c
2
)
'
2
M/
(10
8
M
)
M
yr
−
1
where
L
Edd
= 4
πGMc/σ
T
is the Eddington luminosity (with
σ
T
being
the Thomson electron cross section). New observational capabilities like mm- and IR- interferometry, as provided by
the EHT and GRAVITY (Gravity Collaboration et al. 2018) collaborations now allow to go much closer to the source
which requires a description of general relativistic and dynamical (hence time-dependent) effects.
The most common approach to dynamical relativistic source modeling uses the ideal general relativistic magneto-
hydrodynamic (GRMHD) approximation. It is worth reviewing the nature and quality of the two approximations
inherent in the GRMHD model. First, the plasma is treated as a fluid rather than a collisionless plasma. Second, the
exchange of energy between the plasma and the radiation field is neglected.
The primary EHT sources Sgr A* and M87* fall in the class of low-luminosity active galactive nuclei (AGN) and
accrete with
̇
M/
̇
M
Edd
.
10
−
6
(Marrone et al. 2007) and
̇
M/
̇
M
Edd
.
10
−
5
(Kuo et al. 2014) far below the Eddington
limit. In both cases the accretion flow is believed to form an optically thin disk that is geometrically thick and therefore
has temperature comparable to the virial temperature (see Yuan & Narayan 2014 for a review). The plasma is at
sufficiently high temperature and low density that it is collisionless: ions and electrons can travel
GM/c
2
along
magnetic field lines before being significantly deflected by Coulomb scattering, while the effective mean free path
perpendicular to field lines is the gyroradius, which is typically
GM/c
2
. A rigorous description of the accreting
plasma would thus naively require integrating the Boltzmann equation at far greater expense than integrating the fluid
equations. Full Boltzmann treatments of accretion flows are so far limited to the study of localized regions within the
source (e.g. Hoshino 2015; Kunz et al. 2016). Global models that incorporate nonideal effects using PIC-inspired closure
models suggest, however, that the effects of thermal conduction and pressure anisotropy (viscosity) are small (Chandra
et al. 2015, 2017; Foucart et al. 2017), and thus that one would not do too badly with an ideal fluid prescription.
For Sgr A*, radiative cooling is negligible (Dibi et al. 2012). For M87 radiative cooling is likely important (e.g.
Mo ́scibrodzka et al. 2011; Ryan et al. 2018; Chael et al. 2018). Cooling through the synchrotron process and via
inverse Compton scattering primarily affects the electrons, which are weakly coupled to the ions and therefore need
not be in thermal equilibrium with them. To properly treat the radiation field for the non-local process of Compton
scattering requires solving the Boltzmann equation for photons (the radiative transport equation) in full (e.g. Ryan
et al. 2015) or in truncated form with “closure”. A commonly employed closure is to assume the existence of a frame in
which the radiation field can be considered isotropic, yielding the “M1” closure (Levermore 1984) for which a general
relativistic derivation is shown for example in S ̧adowski et al. (2013). As expected, the computational demands imposed
by the additional “radiation fluid” are considerable. It may however be possible to approximate the effects of cooling
by using a suitable model to assign an energy density (or temperature) to the electrons (Mo ́scibrodzka et al. 2016b).
Again an ideal fluid description, which automatically satisfies energy, momentum, and particle number conservation
laws is not a bad place to start.
It is possible to write the GRMHD equations in conservation form. This enables one to evolve the GRMHD
equations using techniques developed to evolve other conservation laws such as those describing nonrelativistic fluids
and magnetized fluids. Over the last decades, a number of GRMHD codes have been developed, most using conservation
form, and applied to a large variety of astrophysical scenarios (Hawley et al. 1984; Koide et al. 1999; De Villiers &
Hawley 2003; Gammie et al. 2003; Baiotti et al. 2005; Duez et al. 2005; Anninos et al. 2005; Ant ́on et al. 2006; Mizuno
et al. 2006; Del Zanna et al. 2007; Giacomazzo & Rezzolla 2007; Tchekhovskoy et al. 2011; Radice & Rezzolla 2013;
Radice et al. 2014; S ̧adowski et al. 2014; McKinney et al. 2014; Etienne et al. 2015; White et al. 2016; Zanotti &
Dumbser 2015; Meliani et al. 2016; Liska et al. 2018a).
Despite the conceptual simplicity of the MHD equations, the non-linear properties which allow for shocks and
turbulence render their treatment difficult. This is particularly true for the case study considered here: in state-of-the-
art simulations of black hole accretion, angular momentum transport is provided by Maxwell- and Reynolds- stresses
of the orbiting plasma. MHD turbulence is seeded by the magnetorotational instability (MRI) in the differentially
rotating disk (Balbus & Hawley 1991, 1998) and gives rise to chaotic behavior which hinders strict convergence of
6
GRMHD community and the EHTC
the solutions. Nonetheless, it can be demonstrated that certain global properties of the solutions exhibit signs of
convergence.
Another challenge is posed by the “funnel” region near the polar axis where low angular momentum material will be
swallowed up by the black hole (e.g. McKinney 2006). The strong magnetic fields that permeate the black hole (held
in place by the equatorial accretion flow) are able to extract energy in the form of Poynting flux from a rotating black
hole, giving rise to a relativistic “jet” (Blandford & Znajek 1977; Takahashi et al. 1990). The ensuing near-vacuum
and magnetic dominance are traditionally difficult to handle for fluid-type simulations, but analytic calculations (e.g.
Pu et al. 2015) and novel kinetic approaches (Parfrey et al. 2018) can be used to validate the flow in this region.
Due to their robustness when dealing e.g. with super-sonic jet outflows, current production codes typically employ
high-resolution shock-capturing schemes in a finite-volume or finite-difference discretization (Font 2008; Rezzolla &
Zanotti 2013; Mart ́ı & M ̈uller 2015). However, new schemes, for example based on discontinuous Galerkin finite-element
methods, are continuously being developed (Anninos et al. 2017; Fambri et al. 2018).
Given the widespread and increasing application of GRMHD simulations, it is critical for the community to evaluate
the underlying systematics due to different numerical treatments and demonstrate the general robustness of the results.
Furthermore, at the time of writing, all results on the turbulent black hole accretion problem under consideration are
obtained without explicitly resolving dissipative scales in the system (called the implicit large eddy simulation (ILES)
technique). Hence differences are in fact expected to prevail even for the highest resolutions achievable. Quantifying
how large the systematic differences are is one of the main objectives in this first comprehensive code comparison of
an accreting black hole scenario relevant to the EHT science case. This work has directly informed the generation
of the simulation library utilized in the modeling of the EHT 2017 observations (EHT Collaboration 2019b). We use
independent production codes that differ widely in their algorithms and implementation. In particular, the codes that
are being compared are:
Athena++
,
BHAC
,
Cosmos++
,
ECHO
,
H-AMR
,
iharm3D
,
HARM-Noble
,
IllinoisGRMHD
and
KORAL
.
These codes are described further in section 3 below.
The structure of this paper is as follows: in section 2.1 we introduce the GRMHD system of equations, define the
notation used throughout this paper and briefly discuss the astrophysical problem. Code descriptions with references
are given in section 3. The problem setup is described in section 4, where code-specific choices are also discussed.
Results are presented in section 5 and we close with a discussion in section 6.
2.
ASTROPHYSICAL PROBLEM
Let us first give a brief overview of the problem under investigation and the main characteristics of the accretion
flow.
2.1.
GRMHD system of equations
For clarity and consistency of notation, we here give a brief summary of the ideal GRMHD equations in a coordinate
basis (
t,x
i
) with four-metric (
g
μν
) and metric determinant
g
. As customary, Greek indices run through [0
,
1
,
2
,
3] while
Roman indices span [1
,
2
,
3]. The equations are solved in (geometric) code units (i.e. setting the gravitational constant
and speed of light to unity
G
=
c
= 1) where, compared to the standard Gauss cgs system, the factor 1
/
√
4
π
is further
absorbed in the definition of the magnetic field. Hence in the following, we will report times and radii simply in units
of the mass of the central object M.
The equations describe particle number conservation;
∂
t
(
√
−
gρu
t
) =
−
∂
i
(
√
−
gρu
i
)
(1)
where
ρ
is the rest-mass density and
u
μ
is the four-velocity; conservation of energy-momentum:
∂
t
(
√
−
g T
t
ν
)
=
−
∂
i
(
√
−
g T
i
ν
)
+
√
−
g T
κ
λ
Γ
λ
νκ
,
(2)
where Γ
λ
νκ
is the metric connection; the definition of the stress-energy tensor for ideal MHD:
T
μν
MHD
= (
ρ
+
u
+
p
+
b
2
)
u
μ
u
ν
+
(
p
+
1
2
b
2
)
g
μν
−
b
μ
b
ν
(3)
where
u
is the fluid internal energy,
p
the fluid pressure, and
b
μ
is the magnetic field four-vector; the definition of
b
μ
in terms of magnetic field variables
B
i
, which are commonly used as “primitive” or dependent variables:
b
t
=
B
i
u
μ
g
iμ
,
(4)
Event Horizon Code Comparison
7
b
i
= (
B
i
+
b
t
u
i
)
/u
t
;
(5)
and the evolution equation for
B
i
, which follows from the source-free Maxwell equations:
∂
t
(
√
−
gB
i
) =
−
∂
j
(
√
−
g
(
b
j
u
i
−
b
i
u
j
))
.
(6)
The system is subject to the additional no-monopoles constraint
1
√
−
g
∂
i
(
√
−
g B
i
) = 0
,
(7)
which follows from the time component of the homogeneous Maxwell equations. For closure, we adopt the equation
of state of an ideal gas. This takes the form
p
= (ˆ
γ
−
1)
u
, where ˆ
γ
is the adiabatic index. More in-depth discussions
of the ideal GRMHD system of equations can be found in the various publications of the authors, e.g. Gammie et al.
(2003); Anninos et al. (2005); Del Zanna et al. (2007); White et al. (2016); Porth et al. (2017).
To establish a common notation for use in this paper, we note the following definitions: the magnetic field strength
as measured in the fluid-frame is given by
B
:=
√
b
α
b
α
. This leads to the definition of the magnetization
σ
:=
B
2
/ρ
and the plasma-
β β
:= 2
p/B
2
. In addition, we denote with Γ the Lorentz factor with respect to the normal observer
frame.
2.2.
The magnetised black hole accretion problem
We will now discuss the most important features of the problem at hand and introduce the jargon that has developed
over the years. A schematic overview with key aspects of the accretion flow is given in Figure 1.
At very low Eddington rate
̇
M
∼
10
−
6
̇
M
Edd
, the radiative cooling timescale becomes longer than the accretion
timescale. In such radiatively inefficient accretion flows (RIAF), dynamics and radiation emission effectively decouple.
For the primary EHT targets, Sgr A* and M87*, this is a reasonable first approximation and hence purely non-radiative
GRMHD simulations neglecting cooling can be used to model the data. For a RIAF, the protons assume temperatures
close to the virial one which leads to an extremely “puffed-up” appearance of the tenuous accretion
disk
.
In the polar regions of the black hole, plasma is either sucked in or expelled in an outflow, leaving behind a highly
magnetized region called the
funnel
. The magnetic field of the funnel is held in place by the dynamic and static pressure
of the disk. Since in ideal MHD, plasma cannot move across magnetic field lines (due to the frozen-in condition), there
is no way to re-supply the funnel with material from the accretion disk and hence the funnel would be completely
devoid of matter if no pairs were created locally. In state-of-the-art GRMHD calculations, this is the region where
numerical
floor
models inject a finite amount of matter to keep the hydrodynamic evolution viable.
The general morphology is separated into the components of
i)
the
disk
which contains the bound matter
ii)
the
evacuated
funnel
extending from the polar caps of the black hole and the
iii) jet sheath
which is the remaining
outflowing matter. In Figure 1, the regions are indicated by commonly used discriminators in a representative simu-
lation snapshot: the blue contour shows the bound/unbound transition defined via the geometric Bernoulli parameter
u
t
=
−
1
1
, the red contour demarcates the funnel-boundary
σ
= 1 and the green contour the equipartition
β
= 1
which is close to the bound/unbound line along the disk boundary (consistent with McKinney & Gammie (2004)). In
McKinney & Gammie (2004) also a disk-
corona
was introduced for the material with
β
∈
[1
,
3], however as this choice
is arbitrary, there is no compelling reason to label the corona as separate entity in the RIAF scenario.
Since plasma is evacuated within the funnel, it has been suggested that unscreened electric fields in the charge
starved region can lead to particle acceleration which might fill the magnetosphere via a pair cascade (e.g. Blandford
& Znajek 1977; Beskin et al. 1992; Hirotani & Okamoto 1998; Levinson & Rieger 2011; Broderick & Tchekhovskoy
2015). The most promising alternative mechanism to fill the funnel region is by pair creation via
γγ
collisions of seed
photons from the accretion flow itself (e.g. Stepney & Guilbert 1983; Phinney 1995). Neither of these processes is
included in current state of the art GRMHD simulations, however the efficiency of pair formation via
γγ
collisions can
be evaluated in post-processing as demonstrated by Mo ́scibrodzka et al. (2011).
1
There are various ways to define the unbound material in the literature. For example, McKinney & Gammie (2004) used the
geometric
Bernoulli parameter
u
t
. The
hydrodynamic
Bernoulli parameter used for example by Mo ́scibrodzka et al. (2016a) is given by
−
hu
t
where
h
=
c
2
+
p
+
u
denotes the specific enthalpy. Chael et al. (2018); Narayan et al. (2012) included also
magnetic
and
radiative
contributions.
Certainly the geometric and hydrodynamic prescriptions underestimate the amount of outflowing material.
8
GRMHD community and the EHTC
0
10
20
30
r
KS
sin
KS
[M]
10
5
0
5
10
15
20
25
30
r
KS
cos
KS
[M]
log
10
7
6
5
4
3
2
1
0
0
5
10
15
20
25
30
r
KS
sin
KS
[M]
10
5
0
5
10
15
20
25
30
log
10
4
3
2
1
0
1
2
0
5
10
15
20
r
KS
sin
KS
[M]
10
5
0
5
10
15
20
25
30
Funnel
Jet Sheath
Funnel Wall
Disk
Figure 1.
Views of the radiatively inefficient turbulent black hole accretion problem at
t
KS
= 10 000 M against the Kerr-Schild
coordinates (subscript KS).
Left:
logarithmic rest-frame density (hue) and rendering of the magnetic field structure using line-
integral convolution (luminance), showing ordered field in the funnel region and turbulence in the disk.
Center:
the logarithm of
the magnetization with colored contours indicating characteristics of the flow. The magnetized funnel is demarcated by
σ
= 1,
(red), the disk is indicated by
β
= 1 (green) and the geometric Bernoulli criterion (
u
t
=
−
1) is given as blue solid line in the
region outside of the funnel.
Right:
schematic of the main components. In these plots, the black hole horizon is the black disk
and the ergosphere is shown as black contour. The snapshot was obtained from a simulation with
BHAC
.
Turning back to the morphology of the RIAF accretion, Figure 1, one can see that between evacuated funnel
demarcated by the
funnel wall
(red) and bound disk material (blue), there is a strip of outflowing material often also
referred to as the
jet sheath
(Dexter et al. 2012; Mo ́scibrodzka & Falcke 2013; Mo ́scibrodzka et al. 2016a; Davelaar
et al. 2018). As argued by Hawley & Krolik (2006), this flow emerges as plasma from the disk is driven against the
centrifugal barrier by magnetic and thermal pressure (which coined the alternative term
funnel wall jet
for this region).
In current GRMHD based radiation models as utilized e.g. in EHT Collaboration (2019b), as the density in the funnel
region is dominated by the artificial floor model, the funnel is typically excised from the radiation transport. The
denser region outside the funnel wall remains which naturally leads to a limb-brightened structure of the observed
M87 “jet” at radio frequencies (e.g. Mo ́scibrodzka et al. 2016a; Chael et al. 2018; Davelaar et al. 2019 in prep.). In the
mm-band (EHT Collaboration 2019a), the horizon scale emission originates either from the body of the disk or from
the region close to the funnel wall, depending on the assumptions on the electron temperatures (EHT Collaboration
2019b).
In RIAF accretion, a special role is played by the horizon penetrating magnetic flux Φ
BH
: normalized by the accretion
rate
φ
:= Φ
BH
/
√
̇
M
, it was shown that a maximum for the magnetic flux
φ
max
≈
15 (in our system of units) exists
which depends only mildly on black hole spin, but somewhat on the disk scale height (with taller disks being able to
hold more magnetic flux, Tchekhovskoy et al. 2012). Once the magnetic flux reaches
φ
max
, accretion is brought to a
near-stop by the accumulation of magnetic field near the black hole (Tchekhovskoy et al. 2011; McKinney et al. 2012)
leading to a fundamentally different dynamic of the accretion flow and maximal energy extraction via the Blandford &
Znajek (1977) process. This state is commonly referred to as Magnetically Arrested Disk (MAD, Bisnovatyi-Kogan &
Ruzmaikin 1976; Narayan et al. 2003) to contrast with the Standard and Normal Evolution (SANE) where accretion
is largely unaffected by the black hole magnetosphere (here
φ
∼
few). While the MAD case is certainly of great
scientific interest, in this initial code comparison we focus on the SANE case for two reasons:
i)
the SANE case is
already extensively discussed in the literature and hence provides the natural starting point
ii)
the MAD dynamics
poses additional numerical challenges (and remedies) which render it ill-suited to establish a baseline agreement of
GMRHD accretion simulations.
3.
CODE DESCRIPTIONS
Event Horizon Code Comparison
9
In this section, we give a brief, alphabetically ordered overview of the codes participating in this study, with notes
on development history and target applications. Links to public release versions are provided, if applicable.
3.1.
Athena++
Athena++
is a general-purpose finite-volume astrophysical fluid dynamics framework, based on a complete rewrite of
Athena
(Stone et al. 2008). It allows for adaptive mesh refinement in numerous coordinate systems, with additional
physics added in a modular fashion. It evolves magnetic fields via the staggered-mesh constrained transport algorithm
of Gardiner & Stone (2005) based on the ideas of Evans & Hawley (1988), exactly maintaining the divergence-
free constraint. The code can use a number of different time integration and spatial reconstruction algorithms.
Athena++
can run GRMHD simulations in arbitrary stationary spacetimes using a number of different Riemann solvers
(White et al. 2016). Code verification is described in White et al. (2016) and a public release can be obtained from
https://github.com/PrincetonUniversity/athena-public-version
3.2.
BHAC
The BlackHoleAccretionCode (
BHAC
) first presented by Porth et al. (2017) is a multidimensional GRMHD module
for the MPI-AMRVAC framework (Keppens et al. 2012; Porth et al. 2014; Xia et al. 2018).
BHAC
has been designed
to solve the equations of general-relativistic magnetohydrodynamics in arbitrary spacetimes/coordinates and exploits
adaptive mesh refinement techniques with an oct-tree block-based approach. The algorithm is centred on second order
finite volume methods and various schemes for the treatment of the magnetic field update have been implemented,
on ordinary and staggered grids. More details on the various
∇·
B
preserving schemes and their implementation in
BHAC
can be found in Olivares et al. (2018). Originally designed to study black hole accretion in ideal MHD,
BHAC
has
been extended to incorporate nuclear equations of state, neutrino leakage, charged and purely geodetic test particles
(Bacchini et al. 2018, 2019) and non-black hole fully numerical metrics. In addition, a non-ideal resistive GRMHD
module is under development (e.g. Ripperda et al. 2019). Code verification is described in Porth et al. (2017).
3.3.
Cosmos++
Cosmos++
(Anninos et al. 2005; Fragile et al. 2012, 2014) is a parallel, multidimensional, fully covariant, modern
object-oriented (C++) radiation hydrodynamics and MHD code for both Newtonian and general relativistic astrophys-
ical and cosmological applications. Cosmos++ utilizes unstructured meshes with adaptive (
h
-) refinement (Anninos
et al. 2005), moving-mesh (
r
-refinement) (Anninos et al. 2012), and adaptive order (
p
-refinement) (Anninos et al. 2017)
capabilities, enabling it to evolve fluid systems over a wide range of spatial scales with targeted precision. It includes
numerous hydrodynamics solvers (conservative and non-conservative), magnetic fields (ideal and non-ideal), radiative
cooling and transport, geodesic transport, generic tracer fields, and full Navier-Stokes viscosity (Fragile et al. 2018).
For this work, we utilize the High Resolution Shock Capturing scheme with staggered magnetic fields and Constrained
Transport as described in Fragile et al. (2012). Code verification is described in Anninos et al. (2005).
3.4.
ECHO
The origin of the
Eulerian Conservative High-Order
(
ECHO
) code dates back to the year 2000 (Londrillo & Del Zanna
2000; Londrillo & Del Zanna 2004), when it was first proposed a shock-capturing scheme for classical MHD based on
high-order finite-differences reconstruction routines, one-wave or two-waves Riemann solvers, and a rigorous enforce-
ment of the solenoidal constraint for staggered electromagnetic field components (the
Upwind Constraint Transport
,
UCT). The GRMHD version of
ECHO
used in the present paper is described in Del Zanna et al. (2007) and preserves
the same basic characteristics. Important extensions of the code were later presented for dynamical spacetimes (Buc-
ciantini & Del Zanna 2011) and non-ideal Ohm equations (Bucciantini & Del Zanna 2013; Del Zanna et al. 2016; Del
Zanna & Bucciantini 2018). Specific recipes for the simulation of accretion tori around Kerr black holes can be found
in Bugli et al. (2014, 2018). Further references and applications may be found at
www.astro.unifi.it/echo
. Code
verification is described in Del Zanna et al. (2007).
3.5.
H-AMR
H-AMR
is a 3D GRMHD code which builds upon
HARM
(Gammie et al. 2003; Noble et al. 2006) and the public
code
HARM-PI
(https://github.com/atchekho/harmpi) and has been extensively rewritten to increase the code’s speed
and add new features (Liska et al. 2018a; Chatterjee et al. 2019).
H-AMR
makes use of GPU acceleration in a na-
tively developed hybrid CUDA-OpenMP-MPI framework with adaptive mesh refinement (AMR) and locally adaptive