First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole

Black holes are a fundamental prediction of the theory of general relativity (GR; Einstein 1915). A defining feature of black holes is their event horizon, a one-way causal boundary in spacetime from which not even light can escape (Schwarzschild 1916). The production of black holes is generic in GR (Penrose 1965), and more than a century after Schwarzschild, they remain at the heart of fundamental questions in unifying GR with quantum physics (Hawking 1976; Giddings 2017).

Black holes are common in astrophysics and are found over a wide range of masses. Evidence for stellar-mass black holes comes from X-ray (Webster & Murdin 1972; Remillard & McClintock 2006) and gravitational-wave measurements (Abbott et al. 2016). Supermassive black holes, with masses from millions to tens of billions of solar masses, are thought to exist in the centers of nearly all galaxies (Lynden-Bell 1969; Kormendy & Richstone 1995; Miyoshi et al. 1995), including in the Galactic center (Eckart & Genzel 1997; Ghez et al. 1998; Gravity Collaboration et al. 2018a) and in the nucleus of the nearby elliptical galaxy M87 (Gebhardt et al. 2011; Walsh et al. 2013).

Active galactic nuclei (AGNs) are central bright regions that can outshine the entire stellar population of their host galaxy. Some of these objects, quasars, are the most luminous steady sources in the universe (Schmidt 1963; Sanders et al. 1989) and are thought to be powered by supermassive black holes accreting matter at very high rates through a geometrically thin, optically thick accretion disk (Shakura & Sunyaev 1973; Sun & Malkan 1989). In contrast, most AGNs in the local universe, including the Galactic center and M87, are associated with supermassive black holes fed by hot, tenuous accretion flows with much lower accretion rates (Ichimaru 1977; Narayan & Yi 1995; Blandford & Begelman 1999; Yuan & Narayan 2014).

In many AGNs, collimated relativistic plasma jets (Bridle & Perley 1984; Zensus 1997) launched by the central black hole contribute to the observed emission. These jets may be powered either by magnetic fields threading the event horizon, extracting the rotational energy from the black hole (Blandford & Znajek 1977), or from the accretion flow (Blandford & Payne 1982). The near-horizon emission from low-luminosity active galactic nuclei (LLAGNs; Ho 1999) is produced by synchrotron radiation that peaks from the radio through the far-infrared. This emission may be produced either in the accretion flow (Narayan et al. 1995), the jet (Falcke et al. 1993), or both (Yuan et al. 2002).

When viewed from infinity, a nonrotating Schwarzschild (1916) black hole has a photon capture radius ${R}_{{\rm{c}}}=\sqrt{27}\,{r}_{{\rm{g}}}$, where ${r}_{{\rm{g}}}\equiv {GM}/{c}^{2}$ is the characteristic lengthscale of a black hole. The photon capture radius is larger than the Schwarzschild radius R_S that marks the event horizon of a nonrotating black hole, R_S ≡ 2 r_g. Photons approaching the black hole with an impact parameter b < R_c are captured and plunge into the black hole (Hilbert 1917); photons with b > R_c escape to infinity; photons with b = R_c are captured on an unstable circular orbit and produce what is commonly referred to as the lensed "photon ring." In the Kerr (1963) metric, which describes black holes with spin angular momentum, R_c changes with the ray's orientation relative to the angular-momentum vector, and the black hole's cross section is not necessarily circular (Bardeen 1973). This change is small (≲4%), but potentially detectable (Takahashi 2004; Johannsen & Psaltis 2010).

The simulations of Luminet (1979) showed that for a black hole embedded in a geometrically thin, optically thick accretion disk, the photon capture radius would appear to a distant observer as a thin emission ring inside a lensed image of the accretion disk. For accreting black holes embedded in a geometrically thick, optically thin emission region, as in LLAGNs, the combination of an event horizon and light bending leads to the appearance of a dark "shadow" together with a bright emission ring that should be detectable through very long baseline interferometery (VLBI) experiments (Falcke et al. 2000a). Its shape can appear as a "crescent" because of fast rotation and relativistic beaming (Falcke et al. 2000b; Bromley et al. 2001; Noble et al. 2007; Broderick & Loeb 2009; Kamruddin & Dexter 2013; Lu et al. 2014).

The observed projected diameter of the emission ring, which contains radiation primarily from the gravitationally lensed photon ring, is proportional to R_c and hence to the mass of the black hole, but also depends nontrivially on a number of factors: the observing resolution, the spin vector of the black hole and its inclination, as well as the size and structure of the emitting region. These factors are typically of order unity and can be calibrated using theoretical models.

Modern general-relativistic simulations of accretion flows and radiative transfer produce realistic images of black hole shadows and crescents for a wide range of near-horizon emission models (Broderick & Loeb 2006; Mościbrodzka et al. 2009; Dexter et al. 2012; Dibi et al. 2012; Chan et al. 2015; Mościbrodzka et al. 2016; Porth et al. 2017; Chael et al. 2018a; Ryan et al. 2018; Davelaar et al. 2019). These images can be used to test basic properties of black holes as predicted in GR (Johannsen & Psaltis 2010; Broderick et al. 2014; Psaltis et al. 2015), or in alternative theories of gravity (Grenzebach et al. 2014; Younsi et al. 2016; Mizuno et al. 2018). They can also be used to test alternatives to black holes (Bambi & Freese 2009; Vincent et al. 2016; Olivares et al. 2019).

VLBI at an observing wavelength of 1.3 mm (230 GHz) with Earth-diameter-scale baselines is required to resolve the shadows of the core of M87 (M87* hereafter) and of the Galactic center of Sagittarius A* (Sgr A*, Balick & Brown 1974), the two supermassive black holes with the largest apparent angular sizes (Johannsen et al. 2012). At 1.3 mm and shorter wavelengths, Earth-diameter VLBI baselines achieve an angular resolution sufficient to resolve the shadow of both sources, while the spectra of both sources become optically thin, thus revealing the structure of the innermost emission region. Early pathfinder experiments (Padin et al. 1990; Krichbaum et al. 1998) demonstrated the feasibility of VLBI techniques at ∼1.3 mm wavelengths. Over the following decade, a program to improve sensitivity of 1.3 mm-VLBI through development of broadband instrumentation led to the detection of event-horizon-scale structures in both Sgr A* and M87* (Doeleman et al. 2008, 2012). Building on these observations, the Event Horizon Telescope (EHT) collaboration was established to assemble a global VLBI array operating at a wavelength of 1.3 mm with the required angular resolution, sensitivity, and baseline coverage to image the shadows in M87* and Sgr A*.

In this paper, we present and discuss the first event-horizon-scale images of the supermassive black hole candidate M87* from an EHT VLBI campaign conducted in 2017 April at a wavelength of 1.3 mm. The accompanying papers give a more extensive description of the instrument (EHT Collaboration et al. 2019a, Paper II), data reduction (EHT Collaoration et al. 2019b, hereafter Paper III), imaging of the M87 shadow (EHT Collaboration et al. 2019c, hereafter Paper IV), theoretical models (EHT Collaboration et al. 2019d, hereafter Paper V), and the black hole mass estimate (EHT Collaboration et al. 2019e, hereafter Paper VI).

In Curtis (1918), Heber Curtis detected a linear feature in M87, later called a "jet" by Baade & Minkowski (1954). The jet is seen as a bright radio source, Virgo A or 3C 274 (Bolton et al. 1949; Kassim et al. 1993; Owen et al. 2000), which extends out to 65 kpc with an age estimated at about 40 Myr and a kinetic power of about 10⁴² to 10⁴⁵ erg s⁻¹ (de Gasperin et al. 2012; Broderick et al. 2015). It is also well studied in the optical (Biretta et al. 1999; Perlman et al. 2011), X-ray (Marshall et al. 2002), and gamma-ray bands (Abramowski et al. 2012). The upstream end of the jet is marked by a compact radio source (Cohen et al. 1969). Such compact radio sources are ubiquitous in LLAGNs (Wrobel & Heeschen 1984; Nagar et al. 2005) and are believed to be signatures of supermassive black holes.

The radio structures of the large-scale jet (Owen et al. 1989; de Gasperin et al. 2012) and of the core of M87 (Reid et al. 1989; Junor et al. 1999; Hada et al. 2016; Mertens et al. 2016; Kim et al. 2018b; Walker et al. 2018) have been resolved in great detail and at multiple wavelengths. Furthermore, the leveling-off of the "core-shift" effect (Blandford & Königl 1979), where the apparent position of the radio core shifts in the upstream jet direction with decreasing wavelength from increased transparency to synchrotron self-absorption, indicates that at a wavelength of 1.3 mm M87* is coincident with the supermassive black hole (Hada et al. 2011). The envelope of the jet limb maintains a quasi-parabolic shape over a wide range of distances from ∼10⁵ r_g down to ∼20 r_g (Asada & Nakamura 2012; Hada et al. 2013; Nakamura & Asada 2013; Nakamura et al. 2018; Walker et al. 2018).

VLBI observations at 1.3 mm have revealed a diameter of the emission region of ∼40 μas, which is comparable to the expected horizon-scale structure (Doeleman et al. 2012; Akiyama et al. 2015). These observations, however, were not able to image the black hole shadow due to limited baseline coverage.

Based on three recent stellar population measurements, we here adopt a distance to M87 of 16.8 ± 0.8 Mpc (Blakeslee et al. 2009; Bird et al. 2010; Cantiello et al. 2018, see Paper VI). Using this distance and the modeling of surface brightness and stellar velocity dispersion at optical wavelengths (Gebhardt & Thomas 2009; Gebhardt et al. 2011), we infer the mass of M87* to be $M={6.2}_{-0.6}^{+1.1}\times {10}^{9}$ M_⊙ (see Table 9 in Paper VI). On the other hand, mass measurements modeling the kinematic structure of the gas disk (Harms et al. 1994; Macchetto et al. 1997) yield $M={3.5}_{-0.3}^{+0.9}\times {10}^{9}$ M_⊙ (Walsh et al. 2013, Paper VI). These two mass estimates, from stellar and gas dynamics, predict a theoretical shadow diameter for a Schwarzschild black hole of ${37.6}_{-3.5}^{+6.2}\,\mu \mathrm{as}$ and ${21.3}_{-1.7}^{+5}\,\mu \mathrm{as}$, respectively.

The EHT (Paper II) is a VLBI experiment that directly measures "visibilities," or Fourier components, of the radio brightness distribution on the sky. As the Earth rotates, each telescope pair in the network samples many spatial frequencies. The array has a nominal angular resolution of λ/L, where λ is the observing wavelength and L is the maximum projected baseline length between telescopes in the array (Thompson et al. 2017). In this way, VLBI creates a virtual telescope that spans nearly the full diameter of the Earth.

To measure interferometric visibilities, the widely separated telescopes simultaneously sample and coherently record the radiation field from the source. Synchronization using the Global Positioning System typically achieves temporal alignment of these recordings within tens of nanoseconds. Each station is equipped with a hydrogen maser frequency standard. With the atmospheric conditions during our observations the coherent integration time was typically 10 s (see Figure 2 in Paper II). Use of hydrogen maser frequency standards at all EHT sites ensures coherence across the array over this timescale. After observations, recordings are staged at a central location, aligned in time, and signals from each telescope-pair are cross-correlated.

While VLBI is well established at centimeter and millimeter wavelengths (Boccardi et al. 2017; Thompson et al. 2017) and can be used to study the immediate environments of black holes (Krichbaum et al. 1993; Doeleman et al. 2001), the extension of VLBI to a wavelength of 1.3 mm has required long-term technical developments. Challenges at shorter wavelengths include increased noise in radio receiver electronics, higher atmospheric opacity, increased phase fluctuations caused by atmospheric turbulence, and decreased efficiency and size of radio telescopes in the millimeter and submillimeter observing bands. Started in 2009 (Doeleman et al. 2009a), the EHT began a program to address these challenges by increasing array sensitivity. Development and deployment of broadband VLBI systems (Whitney et al. 2013; Vertatschitsch et al. 2015) led to data recording rates that now exceed those of typical cm-VLBI arrays by more than an order of magnitude. Parallel efforts to support infrastructure upgrades at additional VLBI sites, including the Atacama Large Millimeter/submillimeter Array (ALMA; Matthews et al. 2018; Goddi et al. 2019) and the Atacama Pathfinder Experiment telescope (APEX) in Chile (Wagner et al. 2015), the Large Millimeter Telescope Alfonso Serrano (LMT) in Mexico (Ortiz-León et al. 2016), the IRAM 30 m telescope on Pico Veleta (PV) in Spain (Greve et al. 1995), the Submillimeter Telescope Observatory in Arizona (SMT; Baars et al. 1999), the James Clerk Maxwell Telescope (JCMT) and the Submillimeter Array (SMA) in Hawai'i (Doeleman et al. 2008; Primiani et al. 2016; Young et al. 2016), and the South Pole Telescope (SPT) in Antarctica (Kim et al. 2018a), extended the range of EHT baselines and coverage, and the overall collecting area of the array. These developments increased the sensitivity of the EHT by a factor of ∼30 over early experiments that confirmed horizon-scale structures in M87* and Sgr A* (Doeleman et al. 2008, 2012; Akiyama et al. 2015; Johnson et al. 2015; Fish et al. 2016; Lu et al. 2018).

For the observations at a wavelength of 1.3 mm presented here, the EHT collaboration fielded a global VLBI array of eight stations over six geographical locations. Baseline lengths ranged from 160 m to 10,700 km toward M87*, resulting in an array with a theoretical diffraction-limit resolution of ∼25 μas (see Figures 1 and 2, and Paper II).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We observed M87* on 2017 April 5, 6, 10, and 11 with the EHT. Weather was uniformly good to excellent with nightly median zenith atmospheric opacities at 230 GHz ranging from 0.03 to 0.28 over the different locations. The observations were scheduled as a series of scans of three to seven minutes in duration, with M87* scans interleaved with those on the quasar 3C 279. The number of scans obtained on M87* per night ranged from 7 (April 10) to 25 (April 6) as a result of different observing schedules. A description of the M87* observations, their correlation, calibration, and validated final data products is presented in Paper III and briefly summarized here.

At each station, the astronomical signal in both polarizations and two adjacent 2 GHz wide frequency bands centered at 227.1 and 229.1 GHz were converted to baseband using standard heterodyne techniques, then digitized and recorded at a total rate of 32 Gbps. Correlation of the data was carried out using a software correlator (Deller et al. 2007) at the MIT Haystack Observatory and at the Max-Planck-Institut für Radioastronomie, each handling one of the two frequency bands. Differences between the two independent correlators were shown to be negligible through the exchange of a few identical scans for cross comparison. At correlation, signals were aligned to a common time reference using an a priori Earth geometry and clock model.

A subsequent fringe-fitting step identified detections in correlated signal power while phase calibrating the data for residual delays and atmospheric effects. Using ALMA as a highly sensitive reference station enabled critical corrections for ionospheric and tropospheric distortions at the other sites. Fringe fitting was performed with three independent automated pipelines, each tailored to the specific characteristics of the EHT observations, such as the wide bandwidth, susceptibility to atmospheric turbulence, and array heterogeneity (Blackburn et al. 2019; Janssen et al. 2019, Paper III). The pipelines made use of standard software for the processing of radio-interferometric data (Greisen 2003; Whitney et al. 2004; McMullin et al. 2007, I. M. van Bemmel et al. 2019, in preparation).

Data from the fringe-fitting pipelines were scaled from correlation coefficients to a uniform physical flux density scale (in Jansky) by using an independent a priori estimate of the sensitivity of each telescope. The accuracies of the derived station sensitivities were estimated to be 5%–10% in amplitude, although certain uncharacterized losses (e.g., from poor pointing or focus) can exceed the error budget. By assuming total flux density values derived from ALMA interferometric data (Goddi et al. 2019) and utilizing array redundancy via network calibration (Paper III), we refined the absolute amplitude calibration of telescopes that are colocated and have redundant baselines, i.e., ALMA/APEX and JCMT/SMA.

The median scan-averaged signal-to-noise ratio for M87* was >10 on non-ALMA baselines and >100 on baselines to ALMA, leading to small statistical errors in visibility amplitude and phase. Comparisons between the three independent pipelines, the two polarizations, and the two frequency bands enabled estimation of systematic baseline errors of around 1° in visibility phase and 2% for visibility amplitudes. These small limiting errors remain after fitting station sensitivities and unknown station phases via self-calibration (Pearson & Readhead 1984) and affect interferometric closure quantities (Rogers et al. 1974; Readhead et al. 1980). Following data validation and pipeline comparisons, a single pipeline output was designated as the primary data set of the first EHT science data release and used for subsequent results, while the outputs of the other two pipelines offer supporting validation data sets.

The final calibrated complex visibilities V(u, v) correspond to the Fourier components of the brightness distribution on the sky at spatial frequency (u, v) determined by the projected baseline expressed in units of the observing wavelength (van Cittert 1934; Thompson et al. 2017). Figure 2 shows the (u, v) coverage and calibrated visibility amplitudes of M87* for April 11. The visibility amplitudes resemble those of a thin ring (i.e., a Bessel function J₀; see Figure 10.12 in Thompson et al. 2017). Such a ring model with diameter 46 μas has a first null at 3.4 Gλ, matching the minimum in observed flux density and is consistent with a reduced flux density on the longest Hawai'i–Spain baseline (JCMT/SMA-PV) near 8 Gλ. This particular ring model, shown with a dashed line in the bottom panel of Figure 2, is only illustrative and does not fit all features in the data. First, visibility amplitudes on the shortest VLBI baselines suggest that about half of the compact flux density seen on the ∼2 km ALMA–APEX baseline is resolved out by the interferometer beam (Paper IV). Second, differences in the depth of the first minimum as a function of orientation, as well as highly nonzero measured closure phases, indicate some degree of asymmetry in the source (Papers III, VI). Finally, the visibility amplitudes represent only half of the information available to us. We will next explore images and more complex geometrical models that can fit the measured visibility amplitudes and phases.

We reconstructed images from the calibrated EHT visibilities, which provide results that are independent of models (Paper IV). However, there are two major challenges in reconstructing images from EHT data. First, EHT baselines sample a limited range of spatial frequencies, corresponding to angular scales between 25 and 160 μas. Because the (u, v) plane is only sparsely sampled (Figure 2), the inverse problem is under-constrained. Second, the measured visibilities lack absolute phase calibration and can have large amplitude calibration uncertainties.

To address these challenges, imaging algorithms incorporate additional assumptions and constraints that are designed to produce images that are physically plausible (e.g., positive and compact) or conservative (e.g., smooth), while remaining consistent with the data. We explored two classes of algorithms for reconstructing images from EHT data. The first class of algorithms is the traditional CLEAN approach used in radio interferometry (e.g., Högbom 1974; Clark 1980). CLEAN is an inverse-modeling approach that deconvolves the interferometer point-spread function from the Fourier-transformed visibilities. When applying CLEAN, it is necessary to iteratively self-calibrate the data between rounds of imaging to solve for time-variable phase and amplitude errors in the data. The second class of algorithms is the so-called regularized maximum likelihood (RML; e.g., Narayan & Nityananda 1986; Wiaux et al. 2009; Thiébaut 2013). RML is a forward-modeling approach that searches for an image that is not only consistent with the observed data but also favors specified image properties (e.g., smoothness or compactness). As with CLEAN, RML methods typically iterate between imaging and self-calibration, although they can also be used to image directly on robust closure quantities immune to station-based calibration errors. RML methods have been extensively developed for the EHT (e.g., Honma et al. 2014; Bouman et al. 2016; Akiyama et al. 2017; Chael et al. 2018b; see also Paper IV).

Every imaging algorithm has a variety of free parameters that can significantly affect the final image. We adopted a two-stage imaging approach to control and evaluate biases in the reconstructions from our choices of these parameters. In the first stage, four teams worked independently to reconstruct the first EHT images of M87* using an early engineering data release. The teams worked without interaction to minimize shared bias, yet each produced an image with a similar prominent feature: a ring of diameter ∼38–44 μas with enhanced brightness to the south (see Figure 4 in Paper IV).

In the second imaging stage, we developed three imaging pipelines, each using a different software package and associated methodology. Each pipeline surveyed a range of imaging parameters, producing between ∼10³ and 10⁴ images from different parameter combinations. We determined a "Top-Set" of parameter combinations that both produced images of M87* that were consistent with the observed data and that reconstructed accurate images from synthetic data sets corresponding to four known geometric models (ring, crescent, filled disk, and asymmetric double source). For all pipelines, the Top-Set images showed an asymmetric ring with a diameter of ∼40 μas, with differences arising primarily in the effective angular resolutions achieved by different methods.

For each pipeline, we determined the single combination of fiducial imaging parameters out of the Top-Set that performed best across all the synthetic data sets and for each associated imaging methodology (see Figure 11 in Paper IV). Because the angular resolutions of the reconstructed images vary among the pipelines, we blurred each image with a circular Gaussian to a common, conservative angular resolution of 20 μas. The top part of Figure 3 shows an image of M87* on April 11 obtained by averaging the three pipelines' blurred fiducial images. The image is dominated by a ring with an asymmetric azimuthal profile that is oriented at a position angle ∼170° east of north. Although the measured position angle increases by ∼20° between the first two days and the last two days, the image features are broadly consistent across the different imaging methods and across all four observing days. This is shown in the bottom part of Figure 3, which reports the images on different days (see also Figure 15 in Paper IV). These results are also consistent with those obtained from visibility-domain fitting of geometric and general-relativistic magnetohydrodynamics (GRMHD) models (Paper VI).

Zoom In Zoom Out Reset image size

The appearance of M87* has been modeled successfully using GRMHD simulations, which describe a turbulent, hot, magnetized disk orbiting a Kerr black hole. They naturally produce a powerful jet and can explain the broadband spectral energy distribution observed in LLAGNs. At a wavelength of 1.3 mm, and as observed here, the simulations also predict a shadow and an asymmetric emission ring. The latter does not necessarily coincide with the innermost stable circular orbit, or ISCO, and is instead related to the lensed photon ring. To explore this scenario in great detail, we have built a library of synthetic images (Image Library) describing magnetized accretion flows onto black holes in GR¹⁴⁵ (Paper V). The images themselves are produced from a library of simulations (Simulation Library) collecting the results of four codes solving the equations of GRMHD (Gammie et al. 2003; Sa̧dowski et al. 2014; Porth et al. 2017; Liska et al. 2018). The elements of the Simulation Library have been coupled to three different general-relativistic ray-tracing and radiative-transfer codes (GRRT, Bronzwaer et al. 2018; Mościbrodzka & Gammie 2018; Z. Younsi et al. 2019, in preparation). We limit ourselves to providing here a brief description of the initial setups and the physical scenarios explored in the simulations; see Paper V for details on both the GRMHD and GRRT codes, which have been cross-validated for accuracy and consistency (Gold et al. 2019; Porth et al. 2019).

A typical GRMHD simulation in the library is characterized by two parameters: the dimensionless spin ${a}_{* }\equiv {Jc}/{{GM}}^{2}$, where J and M are, respectively, the spin angular momentum and mass of the black hole, and the net dimensionless magnetic flux over the event horizon $\phi \equiv {\rm{\Phi }}/{(\dot{M}{R}_{{\rm{g}}}^{2})}^{1/2}$, where Φ and $\dot{M}$ are the magnetic flux and mass flux (or accretion rate) across the horizon, respectively. Since the GRMHD simulations scale with the black hole mass, M is set only at the time of producing the synthetic images with the GRRT codes. The magnetic flux is generally nonzero because magnetic field is trapped in the black hole by the accretion flow and sustained by currents in the surrounding plasma.

These two parameters allow us to describe accretion disks that are either prograde (a_* ≥ 0) or retrograde (a_* < 0) with respect to the black hole spin axis, and whose accretion flows are either "SANE" (from "Standard and Normal Evolution," Narayan et al. 2012) with ϕ ∼ 1, or "MAD" (from "Magnetically Arrested Disk," Narayan et al. 2003) with ϕ ∼ 15.¹⁴⁶ In essence, SANE accretion flows are characterized by moderate dimensionless magnetic flux and result from initial magnetic fields that are smaller than those in MAD flows. Furthermore, the opening angles of the magnetic funnel in SANE flows are generically smaller than those in MAD flows. Varying a_* and ϕ, we have performed 43 high-resolution, three-dimensional and long-term simulations covering well the physical properties of magnetized accretion flows onto Kerr black holes.

All GRMHD simulations were initialized with a weakly magnetized torus orbiting around the black hole and driven into a turbulent state by instabilities, including the magnetorotational instability (Balbus & Hawley 1991), rapidly reaching a quasi-stationary state. Once a simulation was completed, the relevant flow properties at different times are collected to be employed for the further post-processing of the GRRT codes. The generation of synthetic images requires, besides the properties of the fluid (magnetic field, velocity field, and rest-mass density), also the emission and absorption coefficients, the inclination i (the angle between the accretion-flow angular-momentum vector and the line of sight), the position angle PA (the angle east of north, i.e., counter-clockwise on our images, of the projection on the sky of the accretion-flow angular momentum), the black hole mass M and distance D to the observer.

Because the photons at 1.3 mm wavelength observed by the EHT are believed to be produced by synchrotron emission, whose absorption and emission coefficients depend on the electron distribution function, we consider the plasma to be composed of electrons and ions that have the same temperature in the magnetically dominated regions of the flow (funnel), but have a substantially different temperature in the gas dominated regions (disk midplane). In particular, we consider the plasma to be composed of nonrelativistic ions with temperature T_i and relativistic electrons with temperature T_e. A simple prescription for the ratio of the temperatures of the two species can then be imposed in terms of a single parameter (Mościbrodzka et al. 2016), such that the bulk of the emission comes either from weakly magnetized (small R_high, Te ≃ Ti/R_high) or strongly magnetized (large R_high, Te ≃ Ti) regions. In SANE models, the disk (jet) is weakly (strongly) magnetized, so that low (high) R_high models produce most of the emission in the disk (jet). In MAD models, there are strongly magnetized regions everywhere and the emission is mostly from the disk midplane. While this prescription is not the only one possible, it has the advantage of being simple, sufficiently generic, and robust (see Paper V for a discussion of nonthermal particles and radiative cooling).

Since each GRMHD simulation can be used to describe several different physical scenarios by changing the prescribed electron distribution function, we have used the Simulation Library to generate more than 420 different physical scenarios. Each scenario is then used to generate hundreds of snapshots at different times in the simulation leading to more than 62,000 objects in the Image Library. From the images we have created model visibilities that correspond to the EHT observing schedule and compared them to the measured VLBI visibilities as detailed in Paper VI.

As an example, the top row of Figure 4 shows three GRMHD model snapshots from the Image Library with different spins and flow type, which fitted closure phases and amplitudes of the April 11 data best. For these models we produced simulated visibilities for the April 11 schedule and weather parameters and calibrated them with a synthetic data generation and calibration pipeline (Blecher et al. 2017; Janssen et al. 2019; Roelofs et al. 2019a). The simulated data were then imaged with the same pipeline used for the observed images. The similarities between the simulated images (bottom row of Figure 4) and the observed images (Figure 3) are remarkable.

Zoom In Zoom Out Reset image size

Overall, when combining all the information contained in both the Simulation Library and Image Library, the physical origin of the emission features of the image observed in M87* can be summarized as follows.

First, the observed image is consistent with the hypothesis that it is produced by a magnetized accretion flow orbiting within a few r_g of the event horizon of a Kerr black hole. The asymmetric ring is produced by a combination of strong gravitational lensing and relativistic beaming, while the central flux depression is the observational signature of the black hole shadow. Interestingly, all of the accretion models are consistent with the EHT image, except for the a_* = −0.94 MAD models, which fail to produce images that are sufficiently stable (i.e., the variance among snapshots is too large to be statistically consistent with the observed image).

Second, the north–south asymmetry in the emission ring is controlled by the black hole spin and can be used to deduce its orientation. In corotating disk models (where the angular momentum of the matter and of the black hole are aligned), the funnel wall, or jet sheath, rotates with the disk and the black hole; in counterrotating disk models, instead, the luminous funnel wall rotates with the black hole but against the disk. The north–south asymmetry is consistent with models in which the black hole spin is pointing away from Earth and inconsistent with models in which the spin points toward Earth.

Third, adopting an inclination of 17° between the approaching jet and the line of sight (Walker et al. 2018), the west orientation of the jet, and a corotating disk model, matter in the bottom part of the image is moving toward the observer (clockwise rotation as seen from Earth). This is consistent with the rotation of the ionized gas on scales of 20 pc, i.e., 7000 r_g (Ford et al. 1994; Walsh et al. 2013) and with the inferred sense of rotation from VLBI observations at 7 mm (Walker et al. 2018).

Finally, models with a_* = 0 are disfavored by the very conservative observational requirement that the jet power be P_jet > 10⁴² erg s⁻¹. Furthermore, in those models that produce a sufficiently powerful jet, it is powered by extraction of black hole spin energy through mechanisms akin to the Blandford–Znajek process.

In Paper VI, the black hole mass is derived from fitting to the visibility data of geometric and GRMHD models, as well as from measurements of the ring diameter in the image domain. Our measurements remain consistent across methodologies, algorithms, data representations, and observed data sets.

Motivated by the asymmetric emission ring structures seen in the reconstructed images (Section 5) and the similar emission structures seen in the images from GRMHD simulations (Section 6), we developed a family of geometric crescent models (see, e.g., Kamruddin & Dexter 2013) to compare directly to the visibility data. We used two distinct Bayesian-inference algorithms and demonstrate that such crescent models are statistically preferred over other comparably complex geometric models that we have explored. We find that the crescent models provide fits to the data that are statistically comparable to those of the reconstructed images presented in Section 5, allowing us to determine the basic parameters of the crescents. The best-fit models for the asymmetric emission ring have diameters of 43 ± 0.9 μas and fractional widths relative to the diameter of <0.5. The emission drops sharply interior to the asymmetric emission ring with the central depression having a brightness <10% of the average brightness in the ring.

The diameters of the geometric crescent models measure the characteristic sizes of the emitting regions that surround the shadows and not the sizes of the shadows themselves (see, e.g., Psaltis et al. 2015; Johannsen et al. 2016; Kuramochi et al. 2018, for potential biases).

We model the crescent angular diameter d in terms of the gravitational radius and distance, ${\theta }_{{\rm{g}}}\equiv {GM}/{c}^{2}D$, as d = αθ_g, where α is a function of spin, inclination, and R_high (α ≃ 9.6–10.4 corresponds to emission from the lensed photon ring only). We calibrate α by fitting the geometric crescent models to a large number of visibility data generated from the Image Library. We can also fit the model visibilities generated from the Image Library to the M87* data, which allows us to measure θ_g directly. However, such a procedure is complicated by the stochastic nature of the emission in the accretion flow (see, e.g., Kim et al. 2016). To account for this turbulent structure, we developed a formalism and multiple algorithms that estimate the statistics of the stochastic components using ensembles of images from individual GRMHD simulations. We find that the visibility data are not inconsistent with being a realization of many of the GRMHD simulations. We conclude that the recovered model parameters are consistent across algorithms.

Finally, we extract ring diameter, width, and shape directly from reconstructed images (see Section 5). The results are consistent with the parameter estimates from geometric crescent models. Following the same GRMHD calibration procedure, we infer values of θ_g and α for regularized maximum likelihood and CLEAN reconstructed images.

Combining results from all methods, we measure emission region diameters of 42 ± 3 μas, angular sizes of the gravitational radius θ_g = 3.8 ± 0.4 μas, and scaling factors in the range α = 10.7–11.5, with associated errors of ∼10%. For the distance of 16.8 ± 0.8 Mpc adopted here, the black hole mass is M = (6.5 ± 0.7) × 10⁹ M_⊙; the systematic error refers to the 68% confidence level and is much larger than the statistical error of 0.2 × 10⁹ M_⊙. Moreover, by tracing the peak of the emission in the ring we can determine the shape of the image and obtain a ratio between major and minor axis of the ring that is ≲4:3; this corresponds to a ≲10% deviation from circularity in terms of root-mean-square distance from an average radius.

Table 1 summarizes the measured parameters of the image features and the inferred black hole properties based on data from all bands and all days combined. The inferred black hole mass strongly favors the measurement based on stellar dynamics (Gebhardt et al. 2011). The size, asymmetry, brightness contrast, and circularity of the reconstructed images and geometric models, as well as the success of the GRMHD simulations in describing the interferometric data, are consistent with the EHT images of M87* being associated with strongly lensed emission from the vicinity of a Kerr black hole.

Table 1.  Parameters of M87*

Parameter	Estimate
Ring diameter a d	42 ± 3 μas
Ring width a	$\lt 20\,\mu \mathrm{as}$
Crescent contrast b	>10:1
Axial ratio a	<4:3
Orientation PA	150°–200° east of north
${\theta }_{{\rm{g}}}={GM}/{{Dc}}^{2}$ c	3.8 ± 0.4 μas
$\alpha =d/{\theta }_{{\rm{g}}}$ d	${11}_{-0.3}^{+0.5}$
Mc	(6.5 ± 0.7) × 109 M⊙
Parameter	Prior Estimate
D e	(16.8 ± 0.8) Mpc
M(stars) e	${6.2}_{-0.6}^{+1.1}\times {10}^{9}$ M⊙
M(gas) e	${3.5}_{-0.3}^{+0.9}\times {10}^{9}$ M⊙

Notes.

^aDerived from the image domain. ^bDerived from crescent model fitting. ^cThe mass and systematic errors are averages of the three methods (geometric models, GRMHD models, and image domain ring extraction). ^dThe exact value depends on the method used to extract d, which is reflected in the range given. ^eRederived from likelihood distributions (Paper VI).

Download table as:  ASCIITypeset image

A number of elements reinforce the robustness of our image and the conclusion that it is consistent with the shadow of a black hole as predicted by GR. First, our analysis has used multiple independent calibration and imaging techniques, as well as four independent data sets taken on four different days in two separate frequency bands. Second, the image structure matches previous predictions well (Dexter et al. 2012; Mościbrodzka et al. 2016) and is well reproduced by our extensive modeling effort presented in Section 6. Third, because our measurement of the black hole mass in M87* is not inconsistent with all of the prior mass measurements, this allows us to conclude that the null hypothesis of the Kerr metric (Psaltis et al. 2015; Johannsen et al. 2016), namely, the assumption that the black hole is described by the Kerr metric, has not been violated. Fourth, the observed emission ring reconstructed in our images is close to circular with an axial ratio ≲4:3; similarly, the time average images from our GRMHD simulations also show a circular shape. After associating to the shape of the shadow a deviation from the circularity—measured in terms of root-mean-square distance from an average radius in the image—that is ≲10%, we can set an initial limit of order four on relative deviations of the quadrupole moment from the Kerr value (Johannsen & Psaltis 2010). Stated differently, if Q is the quadrupole moment of a Kerr black hole and ΔQ the deviation as deduced from circularity, our measurement—and the fact that the inclination angle is assumed to be small—implies that ΔQ/Q ≲ 4 (ΔQ/Q = ε in Johannsen & Psaltis 2010).

Finally, when comparing the visibility amplitudes of M87* with 2009 and 2012 data (Doeleman et al. 2012; Akiyama et al. 2015), the overall radio core size at a wavelength of 1.3 mm has not changed appreciably, despite variability in total flux density. This stability is consistent with the expectation that the size of the shadow is a feature tied to the mass of the black hole and not to properties of a variable plasma flow.

It is also straightforward to reject some alternative astrophysical interpretations. For instance, the image is unlikely to be produced by a jet-feature as multi-epoch VLBI observations of the plasma jet in M87 (Walker et al. 2018) on scales outside the horizon do not show circular rings. The same is typically true for AGN jets in large VLBI surveys (Lister et al. 2018). Similarly, were the apparent ring a random alignment of emission blobs, they should also have moved away at relativistic speeds, i.e., at ∼5 μas day⁻¹ (Kim et al. 2018b), leading to measurable structural changes and sizes. GRMHD models of hollow jet cones could show under extreme conditions stable ring features (Pu et al. 2017), but this effect is included to a certain extent in our Simulation Library for models with R_high > 10. Finally, an Einstein ring formed by gravitational lensing of a bright region in the counter-jet would require a fine-tuned alignment and a size larger than that measured in 2012 and 2009.

At the same time, it is more difficult to rule out alternatives to black holes in GR, because a shadow can be produced by any compact object with a spacetime characterized by unstable circular photon orbits (Mizuno et al. 2018). Indeed, while the Kerr metric remains a solution in some alternative theories of gravity (Barausse & Sotiriou 2008; Psaltis et al. 2008), non-Kerr black hole solutions do exist in a variety of such modified theories (Berti et al. 2015). Furthermore, exotic alternatives to black holes, such as naked singularities (Shaikh et al. 2019), boson stars (Kaup 1968; Liebling & Palenzuela 2012), and gravastars (Mazur & Mottola 2004; Chirenti & Rezzolla 2007), are admissible solutions within GR and provide concrete, albeit contrived, models. Some of such exotic compact objects can already be shown to be incompatible with our observations given our maximum mass prior. For example, the shadows of naked singularities associated with Kerr spacetimes with $| {a}_{* }| \gt 1$ are substantially smaller and very asymmetric compared to those of Kerr black holes (Bambi & Freese 2009). Also, some commonly used types of wormholes (Bambi 2013) predict much smaller shadows than we have measured.

However, other compact-object candidates need to be analyzed with more care. Boson stars are an example of compact objects having circular photon orbits but without a surface or an event horizon. In such a spacetime, null geodesics are redirected outwards toward distant observers (Cunha et al. 2016), so that the shadow can in principle be filled with emission from lensed images of distant radio sources generating a complex mirror image of the sky. More importantly, accretion flows onto boson stars behave differently as they do not produce jets but stalled accretion tori that make them distinguishable from black holes (Olivares et al. 2019). Gravastars provide examples of compact objects having unstable photon orbits and a hard surface, but not an event horizon. In this case, while a single image of the accretion flow could in principle be very similar to that coming from a black hole, differences of the flow dynamics at the stellar surface (H. Olivares et al. 2019, in preparation), strong magnetic fields (Lobanov 2017), or excess radiation in the infrared (Broderick & Narayan 2006) would allow one to distinguish a gravastar from a black hole.

Altogether, the results derived here provide a new way to study compact-object spacetimes and are complementary to the detection of gravitational waves from coalescing stellar-mass black holes with LIGO/Virgo (Abbott et al. 2016). Our constraints on deviations from the Kerr geometry rely only on the validity of the equivalence principle and are agnostic about the underlying theory of gravity, but can be used to measure, with ever improved precision, the parameters of the background metric. On the other hand, current gravitational-wave observations of mergers probe the dynamics of the underlying theory, but cannot rely on the possibility of multiple and repeated measurements of the same source.

To underline the complementarity of gravitational-wave and electromagnetic observations of black holes, we note that a basic feature of black holes in GR is that their size scales linearly with mass. Combining our constraints on the supermassive black hole in M87 with those on the stellar-mass black holes detected via gravitational waves we can infer that this property holds over eight orders of magnitude. While this invariance is a basic prediction of GR, which is a scale-free theory, it need not be satisfied in other theories that introduce a scale to the gravitational field.

Finally, the radio core in M87 is quite typical for powerful radio jets in general. It falls on the fundamental plane of black hole activity for radio cores (Falcke et al. 2004), connecting via simple scaling laws the radio and X-ray properties of low-luminosity black hole candidates across vastly different mass and accretion rate scales. This suggests that they are powered by a scale-invariant common object. Therefore, establishing the black hole nature for M87* also supports the general paradigm that black holes are the power source for active galaxies.

We have assembled the EHT, a global VLBI array operating at a wavelength of 1.3 mm and imaged horizon-scale structures around the supermassive black hole candidate in M87. Using multiple independent calibration, imaging, and analysis methods, we find the image to be dominated by a ring structure of 42 ± 3 μas diameter that is brighter in the south. This structure has a central brightness depression with a contrast of >10:1, which we identify with the black hole shadow. Comparing the data with an extensive library of synthetic images obtained from GRMHD simulations covering different physical scenarios and plasma conditions reveals that the basic features of our image are relatively independent of the detailed astrophysical model. This allows us to derive an estimate for the black hole mass of M = (6.5 ± 0.7) × 10⁹ M_⊙. Based on our modeling and information on the inclination angle, we derive the sense of rotation of the black hole to be in the clockwise direction, i.e., the spin of the black hole points away from us. The brightness excess in the south part of the emission ring is explained as relativistic beaming of material rotating in the clockwise direction as seen by the observer, i.e., the bottom part of the emission region is moving toward the observer.

Future observations and further analysis will test the stability, shape, and depth of the shadow more accurately. One of its key features is that it should remain largely constant with time as the mass of M87* is not expected to change measurably on human timescales. Polarimetric analysis of the images, which we will report in the future, will provide information on the accretion rate via Faraday rotation (Bower et al. 2003; Marrone et al. 2007; Kuo et al. 2014; Mościbrodzka et al. 2017) and on the magnetic flux. Higher-resolution images can be achieved by going to a shorter wavelength, i.e., 0.8 mm (345 GHz), by adding more telescopes and, in a more distant future, with space-based interferometry (Kardashev et al. 2014; Fish et al. 2019; Palumbo et al. 2019; F. Roelofs et al. 2019b, in preparation).

Another primary EHT source, Sgr A*, has a precisely measured mass three orders of magnitude smaller than that of M87*, with dynamical timescales of minutes instead of days. Observing the shadow of Sgr A* will require accounting for this variability and mitigation of scattering effects caused by the interstellar medium (Johnson 2016; Lu et al. 2016; Bouman et al. 2018). Time dependent nonimaging analysis can be used to potentially track the motion of emitting matter near the black hole, as reported recently through interferometric observations in the near-infrared (Gravity Collaboration et al. 2018b). Such observations provide separate tests and probes of GR on yet another mass scale (Broderick & Loeb 2005; Doeleman et al. 2009b; Roelofs et al. 2017; Medeiros et al. 2017).

In conclusion, we have shown that direct studies of the event horizon shadow of supermassive black hole candidates are now possible via electromagnetic waves, thus transforming this elusive boundary from a mathematical concept to a physical entity that can be studied and tested via repeated astronomical observations.

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 Tecnología (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 (ERC) 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-1615796, AST-1614868, 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 (grant 2015-R1D1A1A01056807, the Global PhD Fellowship Grant: NRF-2015H1A2A1033752, and 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); the European Union's Horizon 2020 research and innovation programme under grant agreement No 730562 RadioNet; ALMA North America Development Fund; Chandra TM6-17006X.

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).

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 Astrófisica, Ó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 Hawaiian people.