Issue 
A&A
Volume 635, March 2020



Article Number  A143  
Number of page(s)  9  
Section  Astrophysical processes  
DOI  https://doi.org/10.1051/00046361/201937233  
Published online  24 March 2020 
Modeling the orbital motion of Sgr A*’s nearinfrared flares
^{1}
Max Planck Institute for Extraterrestrial Physics (MPE), Giessenbachstr.1, 85748 Garching, Germany
email: bauboeck@mpe.mpg.de
^{2}
LESIA, Observatoire de Paris, Université PSL, CNRS, Sorbonne Université, Univ. Paris Diderot, Sorbonne Paris Cité, 5 Place Jules Janssen, 92195 Meudon, France
^{3}
MaxPlanckInstitute for Astronomy, Königstuhl 17, 69117 Heidelberg, Germany
^{4}
1. Physikalisches Institut, Universität zu Köln, Zülpicher Str. 77, 50937 Köln, Germany
^{5}
Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France
^{6}
CENTRA and Universidade de Lisboa – Faculdade de Ciências, Campo Grande, 1749016 Lisboa, Portugal
^{7}
CENTRA and Universidade do Porto – Faculdade de Engenharia, 4200465 Porto, Portugal
^{8}
European Southern Observatory, KarlSchwarzschildStr. 2, 85748 Garching, Germany
^{9}
European Southern Observatory, Casilla 19001, Santiago 19, Chile
^{10}
Sterrewacht Leiden, Leiden University, Postbus 9513, 2300 RA Leiden, The Netherlands
^{11}
Departments of Physics and Astronomy, Le Conte Hall, University of California, Berkeley, CA 94720, USA
^{12}
School of Physics and Astronomy, Tel Aviv University, Tel Aviv 69978, Israel
^{13}
MaxPlanckInstitute for Radio Astronomy, Auf dem Hügel 69, 53121 Bonn, Germany
^{14}
Center for Computational Astrophysics, Flatiron Institute, 162 5th Ave., New York, NY 10010, USA
^{15}
JILA and Department of Astrophysical and Planetary Sciences, University of Colorado, Boulder, CO 80309, USA
^{16}
Department of Particle Physics and Astrophysics, Weizmann Institute of Science, Rehovot 7610001, Israel
Received:
2
December
2019
Accepted:
4
February
2020
Infrared observations of Sgr A* probe the region close to the event horizon of the black hole at the Galactic center. These observations can constrain the properties of lowluminosity accretion as well as that of the black hole itself. The GRAVITY instrument at the ESO VLTI has recently detected continuous circular relativistic motion during infrared flares which has been interpreted as orbital motion near the event horizon. Here we analyze the astrometric data from these flares, taking into account the effects of outofplane motion and orbital shear of material near the event horizon of the black hole. We have developed a new code to predict astrometric motion and flux variability from compact emission regions following particle orbits. Our code combines semianalytic calculations of timelike geodesics that allow for outofplane or elliptical motions with ray tracing of photon trajectories to compute timedependent images and light curves. We apply our code to the three flares observed with GRAVITY in 2018. We show that all flares are consistent with a hotspot orbiting at R ∼ 9 gravitational radii with an inclination of i ∼ 140°. The emitting region must be compact and less than ∼5 gravitational radii in diameter. We place a further limit on the outofplane motion during the flare.
Key words: black hole physics / Galaxy: center / accretion / accretion disks
GRAVITY has been developed by a collaboration of the Max Planck Institute for Extraterrestrial Physics, LESIA of Paris Observatory /CNRS/UPMC/Univ. Paris Diderot and IPAG of Université Grenoble Alpes/CNRS, the Max Planck Institute for Astronomy, the University of Cologne, the Centro de Astrofísica e Gravitação, and the European Southern Observatory.
© GRAVITY Collaboration et al. 2020
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Open Access funding provided by Max Planck Society.
1. Introduction
The black hole at the Galactic center provides a unique laboratory for exploring the physics of accretion and supermassive black holes. As it is the closest supermassive black hole, Sgr A* has the largest angular size of any black hole observable from Earth, with a Schwarzschild radius of R_{s} = 10.022 ± 0.020_{stat.} ± .032_{sys.} μas (Gravity Collaboration 2019). Additionally, the characteristic timescale of infalling matter near the black hole is on the order of minutes to hours (e.g., Baganoff et al. 2003; Meyer et al. 2009; Witzel et al. 2012, 2018; Hora et al. 2014), allowing variability and motion to be observed within a single night.
The GRAVITY instrument at the VLT Interferometer (Eisenhauer et al. 2008; Gravity Collaboration 2017) allows for highprecision astrometry of the Galactic center. This instrument has proven capable of an astrometric precision of 10 − 30 μas during bright flares (Gravity Collaboration 2018a). The observation of coherent circular motion of the centroid of emission during infrared flares lends new insight into the geometry and physical properties of matter in the accretion flow, as well as the spacetime near the black hole.
In the past, efforts to predict nearinfrared flares of the accretion onto Sgr A* have fallen into a few broad categories. The first of these involves detailed general relativistic magnetohydrodynamic simulations of a small volume around the black hole (e.g., Dolence et al. 2009; Dexter & Fragile 2013; Chan et al. 2015; Ball et al. 2016; Ressler et al. 2017). These simulations allow for selfconsistent modeling of the hydrodynamics of the accretion flow in the gravitational field of the black hole. However, the computational expense of these models preclude exploration of a large parameter space, making it difficult to constrain models based on observations. Simulations are usually able to produce only a limited number of orbits, making them sensitive to the initial conditions. Moreover, these models often fail to reproduce observed properties of Sgr A* such as the wavelengthdependent variability or the polarization fraction or angle during flares.
Yuan et al. (2003) developed semianalytic models of the spectral energy distribution and multiwavelength variability of Sgr A*. These models successfully match the global observed properties of the black hole, but they do not explain the mechanism causing flares nor do they predict the distribution of matter near the black hole.
Another class of recent models consist of ParticleinCell (PIC) simulations of a small region of the accretion flow (e.g., Sironi & Narayan 2015; Rowan et al. 2017, 2019). These models provide insight into the microphysical properties of the accretion flow, but they are computationally limited to very small regions. Predicting global observables, therefore, is not yet possible.
The last class of models consists of semianalytic calculations of hotspots in orbits around the central black hole (e.g., Broderick & Loeb 2005, 2006; Paumard et al. 2005; Meyer et al. 2006; Trippe et al. 2007; Hamaus et al. 2009; Zamaninasab et al. 2010). These models are computationally cheap enough to make it possible to explore a large parameter space of orbital configurations and black hole properties. However, they generally do not include any magnetohydrodynamic effects and, therefore, they cannot lend any insight into, for example, the mechanism that causes flares or the microphysical properties of the accretion flow.
In the past, semianalytic hotspot models have been used to perform numerical calculations of photon geodesics through the curved spacetime near the black hole, but they have been limited to particle orbits that can be calculated analytically. In practice, this has restricted calculations to a small family of orbits–generally, circular orbits in the equatorial plane of the black hole. However, there is no a priori reason to expect the transiently heated, relativistic electrons responsible for infrared flares to be confined to the midplane or to circular orbits.
In this paper, we describe the NERO code, which generalizes the hotspot model to allow for arbitrary particle orbits in the vicinity of a Kerr black hole. Models with noncircular, unbound, or outflowing orbits are similarly allowed in our simulations. We discuss the methods for calculating photon (Sect. 2) and particle (Sect. 3) geodesics to generate synthetic flare observables incorporating all relativistic effects. We use the resulting model to analyze the GRAVITY data from 2018 (Sect. 4). We consider both circular orbits as in past works (Gravity Collaboration 2018a) and obtain new constraints on outofplane motion (Sect. 5), as well as the size of the emission region due to shear (Sect. 6).
2. Ray tracing
The first step in calculating the appearance of the region near the black hole at the Galactic center is to find the paths followed by the photons as they travel towards a distant observer. We use the GEOKERR code (Dexter & Agol 2009) to calculate null geodesics between the region around the black hole and a distant image plane. We assume a Kerr metric characterized by the line element
where Σ = r^{2} + a^{2}cos^{2}θ, Δ = r^{2} − 2Mr + a^{2}, M is the mass of the black hole, and a = J/M^{2} is the dimensionless angular momentum. Here and throughout this paper, we use the convention G = c = 1.
There are two possible methods for finding particle trajectories in curved spacetimes. The more general approach is to directly solve the geodesic equations numerically. Alternatively, in certain astrophysically relevant metrics, the equations of motion can be reduced analytically to elliptical integrals (Cunningham & Bardeen 1973; Rauch & Blandford 1994). These integrals can then be calculated semianalytically to find the particle trajectories.
The GEOKERR code adopts the latter approach to solve the geodesic equations in the Kerr metric (Dexter & Agol 2009). Starting from an image plane at the location of a distant observer, the code solves the geodesic equations backwards in time to find the trajectories of the photons originating near the black hole that terminate on the image plane.
3. Relativistic orbits
The second component of the NERO code is the calculation of the relativistic orbits of test particles near the event horizon of the black hole. We are particularly interested in calculating a broad range of orbits without constraining particles to be in circular orbits or in the equatorial plane of the black hole.
Since we are interested in highly relativistic orbits close to the black hole, the usual Keplerian elements are a poor choice for characterizing orbits. Instead, we follow Hughes (2001) and identify each orbit by its energy E, z angular momentum L_{z}, and Carter constant Q in addition to choosing the spin a of the black hole. While the procedure we use is capable of calculating general orbits with arbitrary inclination and eccentricity, we confine our initial analysis to circular orbits. In this case, we use the orbital radius R as an orbital parameter along with L_{z}. We sample R between the innermost stable circular orbit and an arbitrary large value (we choose R = 14 for the analysis below).
In order to sample the range of allowed angular momenta, we first calculate the specific angular momentum of a prograde and retrograde orbit, respectively,
We can then create a family of orbits at a given radius by sampling L_{z} evenly between these bounds to obtain orbits with a range of inclinations to the spin vector of the black hole.
The condition of circularity on the orbits introduces the constraints,
where ℛ is the effective potential in the radial direction.
These constraints allow us to find the specific energy and specific Carter constant for each value of L_{z} and R. We can now find the energy and Carter constant of the orbit (Hughes 2001),
Once we have calculated the constants of motion, we check to make sure the orbit is stable. The criterion for stability is
Orbits in the region of the parameter space that does not satisfy Eq. (8) are unstable and therefore astrophysically irrelevant.
For each orbit we wish to model, we have now calculated the constants of motion E, L_{z}, and Q. The equations of motion for a test particle following a timelike geodesic through the Kerr metric in terms of these constants are given by (Carter 1968):
Here, τ is the proper time, Θ is the effective potential in the angular direction,
and T is
We calculate the solution to these equations of motion using the publiclyavailable YNOGKM code (Yang & Wang 2014). This code adapts the GEOKERR algorithm to semianalytically find particle trajectories in the KerrNewman metric for a given set of initial conditions. Since we are concerned with orbits about Sgr A*, we set the charge of the black hole to zero, reducing the metric to the Kerr solution given in Eq. (1).
The YNOGKM code takes as input the initial position and velocity of a test particle. As an initial position, we arbitrarily choose a point in the equatorial plane (θ = π/2) with a phase ϕ = π/2, with r set to the orbital radius R. The initial velocity in the BoyerLindquist frame can be calculated from the constants of motion:
where
We transform these velocities from BoyerLindquist coordinates into the locally nonrotating frame (LNRF, see Bardeen et al. 1972) for use by the YNOGKM code:
where
The result of the YNOGKM and GEOKERR codes is a set of position 4vectors corresponding to the trajectories of the test particle as well as photons propagating to the observer. The intersections of these trajectories correspond to a point in spacetime at which a photon is emitted by the test particle that reaches the observer. To simplify the task of finding these intersections, we assign the hotspot a Gaussian emissivity profile with a size s (set to a fiducial value of s = 0.5M):
where x is the distance from the center of the hotspot.
Given the positions and velocities of a particle as well as the constants of motion of the photon trajectory, we can calculate the Doppler shift of the intersection as
In terms of the metric elements defined above, we write this Doppler shift as
where is the usual relativistic factor,
With our assumption that the emission is optically thin, we can now find the flux observed at the image plane by integrating the emissivity scaled by the Doppler term backwards along the photon trajectories calculated from the GEOKERR code:
The time of arrival of the photons at the image plane is the sum of orbital time and the lighttravel time from the hotspot to the observer.
The end result of our simulation is a series of snapshots at increasing observer timesteps of the region immediately around the central black hole. For each image, we can then calculate a centroid position,
and a flux,
Figure 1 shows an example series of images of a hotspot passing behind the black hole. As the spot moves behind the black hole, the photon trajectories are strongly lensed in the direction of the observer, resulting in the corresponding peak in the lightcurve. Simultaneously, the centroid of the emission moves towards the location of the black hole.
Fig. 1. Sample image of a hotspot at different points in the orbit. Top row: image of the hotspot before, during, and after a strong lensing event. Topright plot: image averaged over the whole orbit. Lower left panel: lightcurve produced over several orbits, lower right panel: motion of the centroid of emission. This figure was produced with a hotspot in a circular orbit at R = 5M around a black hole with spin a = 0.9. The spin axis of the black hole is inclined i = 15° from the observer’s line of sight and the relative inclination between the spin axis and the orbital plane is θ = 50.54°. 
4. Orbit fitting
The GRAVITY instrument has found evidence of orbital motion around Sgr A*. Here we describe the process used to fit the observed motion with the NERO code to obtain constraints on the properties of the emitting region (Gravity Collaboration 2018a).
We construct a grid of model centroid motions across our parameter space. For the purposes of this work, we limit ourselves to fitting only for the orbital radius and inclination to the observer. Therefore, we calculate a range of models between R = 4 and R = 14 in increments of ΔR = 0.5 across all inclinations (1° < i < 180°, Δi = 5). Inclinations between 0° and 90° correspond to counterclockwise orbits on the sky, while inclinations above 90° correspond to clockwise orbits. At each point in this grid, we further vary the position angle of the orbit on the sky, the orbital phase at the start of the observation, and the x and yposition of the central mass. For each parameter combination, we calculate the χ^{2} difference between the model and the astrometric data. The location of the minimum χ^{2} value determines the best fit to the observation.
Figure 2 shows the bestfit models for each of the three flares observed during the summer of 2018. The reduced χ^{2} for each fit is 1.6, 4.5, and 2.3 for July 22, July 28, and May 27, respectively. In each case, the bestfit model lies mostly or entirely within the astrometric data points. This is due to our assumption that the orbiting hotspot follows a geodesic of motion, meaning that its orbital speed is determined by its distance from the black hole. At face value, the data from all three nights show motion that is faster than the expected orbital speed or equivalently, a radius that is larger than expected.
Fig. 2. Data and NERO fits to the flares from May 27, July 22, and July 28, 2018. Data points correspond to the centroid position observed by GRAVITY over the course of the flare. The solid lines show the orbits that best fits each of three flares. Colors indicate the time since the start of the flare. 
One effect that partially accounts for this tension is related to the scatter in the observed centroid points. Since the distance from the center of the orbit is a positive definite quantity, any scatter in the measurements will tend to increase the observed radius by a factor of . In the case of an orbit at R = 7M observed with 30 μas error bars, this effect implies a 10 μas bias in the observed radius. Other possible explanations are that the motion of the hotspot is not in the plane of the black hole but rather moving outward towards the observer (see Sect. 6 below) or that other forces, such as a pressure gradient or the magnetic field, have a significant effect on the motion of the hotspot.
We marginalize over the various parameters to find the 1, 2, and 3σ allowed regions, corresponding to the area containing 39%, 86%, and 98.8% confidence regions. Given the large values of the χ^{2} of the bestfit models, we scale the uncertainties of the model parameters such that the reduced χ^{2} of the bestfit model is 1. The first three panels of Fig. 3 show the constraints in the R − i plane for all three flares. The clockwise motion during the flares constrains the inclination of the orbit to be between ∼120° −175° for all but the May 27 flare. In the latter case, the bestfit inclination is also in this range (140°) but the allowed region extends over a larger fraction of the parameter space due to the poor astrometric precision of the data. Unless otherwise noted, we include a sin i prior on the inclination of the orbit to account for the geometric bias favoring edgeon orientations. We find that the astrometric fit favors models with moderate to high inclinations and radii between 7 and 12M. The allowed parameter spaces of all three flares overlap, indicating that a single model is consistent with all the flares.
Fig. 3. Confidence contours corresponding to 1, 2, and 3σ limits on the orbital radius and inclination for all three flares. The last panel shows the combined constraints, assuming the radius, inclination, and position angle are the same for all three flares. 
The contours shown in Fig. 3 are slightly less constraining than those found in Gravity Collaboration (2018a). This difference is due to the additional degrees of freedom allowed in the fit corresponding to the uncertainty of the position of the black hole. In previous fits, the center of the orbit was placed at the average of all the astrometric points over the course of the flare. Figure 4 illustrates the degeneracy between the xposition of the central mass and the radius of the orbit. This degeneracy indicates that the flare astrometry allows for models in which the flare represents only a small fraction of an orbit with a large radius, where the black hole is correspondingly offset from the center of emission. However, these models can be ruled out based on our knowledge of the position of the center of mass from the measured orbits of the Sstars. In particular, the orbit of the star S2 constrains the position of the center of mass to within ∼50 μas (Gravity Collaboration 2018b).
Fig. 4. Confidence contours showing the allowed values for the orbital radius and the xoffset of the center of the black hole for the July 22, 2018 flare. An offset of zero corresponds to the black hole lying at the mean location of emission averaged over the entire flare. 
The fits presented here only make use of the measurements of the centroid position of the emission throughout the flare. Two additional sources of information are available that can further constrain the parameters of the fit. The first is the absence of clear Doppler beaming during a portion of the flare. For orbits that are close to edgeon to the observer, the relativistic speed at which the hotspot moves should lead to significant brightening on the approaching side and dimming on the receding side. These effects are, of course, convolved with the intrinsic light curve of the flare, making it difficult to place precise constraints on the inclination. However, orbits with inclinations between 40° and 140° would lead to Doppler beaming comparable to the overall brightening of Sgr A* during the flare, making these scenarios unlikely (Gravity Collaboration 2018a).
Since we do not know the underlying mechanism for generating flares, it is impossible for us to construct light curves that can be compared directly to the data. However, the NERO code does calculate the Doppler boosting, lensing, and other relativistic effects of the motion of the hotspot on the observed flux. By combining these predictions with the observed lightcurves, we can find the intrinsic flare lightcurve as generated at the source. Figure 5 shows both the observed and inferred intrinsic lightcurves for each of the three flares discussed in this paper.
Fig. 5. Observed and inferred intrinsic lightcurves for the three GRAVITY flares. The flux is measured relative to the nearby star S2. Each dotted line corresponds to the lightcurve during the respective flare as observed by GRAVITY. The solid line shows the intrinsic lightcurve once all relativistic effects predicted from the bestfit orbital model have been remove. 
A further constraint arises from the observed change in the polarization angle during the course of the flare. The polarization angle changes on a timescale similar to the orbital timescale at the radii favored by the astrometric data. A detailed analysis of the constraints arising from the polarization will be addressed in a forthcoming paper (JiménezRosales et al., in prep.).
If we make the assumption that all of the flares share the same parameters, we can obtain a combined constraint for all of the flares. In order to calculate the combined constraint, we allow the orbital phase at the start of each observation as well as the x and yposition of the central mass to vary between observations. We hold the radius, inclination, and position angle of the model constant across all three observations. The last panel of Fig. 3 shows the constraint on the R − i parameter space from all three flares combined. The preferred radius is between 8 and 10M, and the inclination lies between i = 120° and i = 150°.
Figure 6 shows the bestfit region from the astrometry combined with the constraints from the lightcurve and the polarization. The contrast constraint is based on a conservative upper limit on the amount of Doppler boosting that is consistent with the observed lightcurve during the flare. The two hatched regions correspond to the radius range determined from the rate of change of the polarization angle and the behavior of the polarization angle in the Q−U plane. We note that this constraint is based on the assumption that the magnetic field is vertically oriented in the emitting region of the accretion flow (see Gravity Collaboration 2018a as well as JiménezRosales et al., in prep. for details).
Fig. 6. Combined constraints together with information from the polarization and the Doppler beaming. The blue astrometric constraints are as shown in Fig. 3. The green shaded region corresponds to models where the maximum contrast of the lightcurve is less than 2.5. The hatched polarization period constraint arises from the assumption that the change in polarization angle during the flare is caused by the motion of the hotspot through the magnetic field. The polarization loop constraint is based on a model with a vertical magnetic field and the requirement that the polarization angle matches the behavior observed during the flare. 
We also show the constraints on the orientation of the orbiting hotspots in Fig. 7. Interestingly, the allowed region of orbital orientations is broadly consistent with those of the stars in the clockwise stellar disk, as marked by the black cross in the figure (Gillessen et al. 2017, see also Bartko et al. 2009; Yelda et al. 2014). This clockwise disk consists of Ostars with significant stellar winds (Paumard et al. 2006). These winds have been conjectured as a possible source of matter to fuel the accretion onto the central black hole. If this is the case, the accretion flow is expected to roughly align with the stars in the clockwise disk (Ressler et al. 2018; Calderón et al. 2020). While our constraints cannot definitively confirm this picture, improved constraints from future observations may be able to pinpoint whether the accreting material originates in the winds from massive stars in the clockwise disk.
Fig. 7. Combined constraints on the orientation of the hotspot orbit from all three flares. The black marker indicates the approximate orientation of the clockwise disk of stars. 
5. Orbital shear
If the circular centroid motion during infrared flares is indeed due to the orbit of a hotspot, the hotspot is expected to shear over the course of its orbit due the different orbital velocity at different radii. Although the shape of the hotspot is unresolved, this shearing would manifest itself as an inward spiral of the centroid over the course of the flare. By incorporating the orbital shear into our model of the centroid motion and comparing this prediction to the observed astrometry, we can, in principle, obtain an upper limit on the size of the emitting region.
We allow for the possibility of shearing by modifying our model to account for the differential orbital velocity at different points in the emitting region. We define a parameter s corresponding to the radial extent of the hotspot. For a given orbital radius R, we calculate a series of particle trajectories with starting positions ranging from R − s to R + s. We then sum the emission from a Gaussian hotspot on each of these geodesics to find the total predicted lightcurve and centroid track.
We use the same procedure described above to fit the resulting shearedspot models to the data. In this case, we fix the inclination and restrict ourselves to orbital radii between 5M and 12.5M. We let the diameter of the emitting region vary between 1.0M and 6.5M. Given the poor astrometric precision of the flare observations on May 27 and July 28, we have restricted this analysis to the July 22 flare. Figure 8 shows contours of constant χ^{2} corresponding to 1, and 2σ limits on the orbital radius and spot extent. As the size of the spot increases, the quality of the fit decreases. Spots larger than 5M in diameter can be ruled out at the 1σ level, although larger spots are still within the 2σ limits from the current observations.
Fig. 8. 1 and 2σ constraints on the orbital radius and spot diameter for the July 22 flare. 
6. Outofplane motion
Although the detection of circular motion during flares would seem to favor the interpretation that flares originate in accreting material, it is, in principle, also consistent with material launched at the base of a jet. This material would retain some angular momentum from the accretion flow and therefore spiral outward from the jetlaunching region.
To investigate this possibility, we simulate the motion of material launched vertically from an accretion flow. In particular, we modify our existing models of orbital motion by adding a constant vertical component v_{z} to the velocity of the hotspot. We then renormalize the fourvelocity of the hotspot to conserve the relation u ⋅ u = −1. The result is a hotspot on a helical track originating in the equatorial plane. Since the volume of the region we simulate is limited for computational reasons, we have restricted ourselves to scenarios where v_{z} < 0.3c to prevent the hotspot from leaving the bounds of our simulation during a single orbit.
We again use the procedure described above to fit the parameters of the model to the July 22 flare. If the axis of a jet is significantly inclined with respect to the observer, the hotspot should display a systematic linear motion over time. The relatively circular motion measured by GRAVITY rules out this scenario. Therefore, we restrict the inclination to be greater than 160°. Given the faceon inclination inherently required by this scenario, we do not weight the probability distribution by sin i as above. We allow the vertical velocity to vary between 0.05c and 0.3c.
Figure 9 shows the allowed region of the parameter space for the July 22 flare. As expected, the combination of high inclination and large vertical velocity has been ruled out. For a sufficiently low inclination, however, the primary effect of the vertical velocity is to change the observed orbital period due to the propagation time of the photons and the decreasing distance to the hotspot. For faceon orbits, the observed centroid motion is consistent with a significant vertical velocity. In fact, adding a vertical velocity component marginally improves the quality of the fit, with the bestfit parameters being R = 9.0M, i = 175°, v_{z} = 0.15c.
Fig. 9. Constant χ^{2} contours corresponding the 1, and 2σ regions of the vertical velocity and inclination for the July 22 flare. We note that in this case we use a flat prior for the inclination of the orbit. 
7. Prospects for improved constraints
There are two possible avenues for improving the constraints on the flare properties from astrometric measurement. The first consists of continued observations of flares at the current level of precision. As shown in the last panel of Fig. 3, coadding several flares reduces the statistical error on the inferred parameters. Additional flare observations will further tighten the constraints. However, these improvements are contingent on the assumption that some parameters of the flares (i.e., the radius, inclination, and position angle) are consistent between multiple flares.
The second avenue for improved constraints stems from increasing the astrometric precision of measurements of individual flares. To simulate this scenario, we have generated simulated data sets at various levels of precision and used the NERO code to find the bestfit model and resulting constraints on the flare properties. Figure 10 shows the confidence contours for three such simulations. In each case, the underlying model is a hotspot orbiting at a radius of R = 8M and an inclination of i = 155°. Each simulated data set consists of ten data points chosen at time intervals corresponding to the observed times of the July 22nd flare. The x and ydata points are drawn from a normal distribution centered on the model value with a standard deviation corresponding to an error of 10, 20, and 30 μas, respectively.
Fig. 10. Constraints from simulated data with improved error bars. In each case, the underlying model is based on a hotspot orbiting at R = 8.0M with an inclination of i = 155°, as marked by the red cross in each panel. The astrometric error is 10 μas, 20 μas, and 30 μas, respectively. 
We then apply the same fitting procedure to our simulated data points as we use on the real data. The constraints from these fits are shown in Fig. 10. For small error bars, the radius and inclination are tightly constrained. We correctly recover the parameters of the underlying simulation, as marked in the figure. For a 30 μas error, which most closely corresponds to the current GRAVITY precision, the simulated data are consistent with radii between 8M and 14M. The inclination is only loosely constrained.
8. Summary and discussion
The NERO code combines the ability to numerically calculate relativistic orbits near a black hole with efficient raytracing between the vicinity of the black hole and a distant observer. The result is a hotspot model for fitting both astrometric and lightcurve data with the ability to account for the full general relativistic treatment of the spacetime near the black hole. It also takes into account a variety of physical scenarios, including hotspot orbits that are inclined to the spin axis, eccentric orbits, outofplane motion, and orbital shear from hotspots with a significant radial extent.
We apply this model to the three flares observed during the summer of 2018 with the GRAVITY instrument at the VLT. We find that all three flares are well fit by hotspots on circular orbits with radii between 7 and 10M and inclinations greater than ∼90°. A single model with R = 9M is consistent with all three observed flares at a range of inclinations.
While the observed centroid motions are consistent with circular orbital motion, the size of the observed motion on the sky seems to be systematically larger than what has been predicted by our models. Equivalently, the period of the motion is smaller than predicted by the observed size of the centroid track. One possible explanation for this discrepancy is that the hotspot is moving towards the observer at a significant fraction of the speed of light, causing the observed orbit to appear faster than it really is. We find that models with an additional vertical velocity component are consistent with the observations, but they are not significantly preferred over orbits in the equatorial plane.
We also investigate the effects of orbital shear on the predicted motion of the centroid. For a hotspot with a size comparable to its orbital radius, the differential velocity at different radii should lead to a centroid track that spirals inward over the course of the orbit. We find that smaller hotspots provide a better fit to the observations, but hotspots with diameters as large as 6M are allowed by the current astrometric precision.
Additional constraints on the viability of a model of orbital motion comes from the observed light curve and polarization of the flare. A hotspot orbiting at a high inclination to the observer is expected to show significant variation in brightness due to the difference in Doppler beaming between the approaching and receding parts of the orbit. However, this constraint is significantly complicated by the unknown emission mechanism of infrared flares, which leads to a large uncertainty in the expected light curve. In addition to astrometric motion during infrared flares, the polarization angle of the emission also varies on a timescale consistent with the orbital period of a hotspot. Using these observations to constrain the orbit requires a model of the magnetic field near the black hole and will be presented in a future work.
At the current level of precision, astrometric measurements during infrared flares provide strong evidence that matter is in orbital motion around the central black hole. The NERO code allows for relativistic orbits to be fit to the astrometric data in order to constrain the parameters of the orbits. The combination of multiple flare observations can reduce the uncertainties on the parameters of the orbiting material. Further observations of greater precision will be able to improve these constraints even further.
References
 Baganoff, F. K., Maeda, Y., Morris, M., et al. 2003, ApJ, 591, 891 [NASA ADS] [CrossRef] [Google Scholar]
 Ball, D., Özel, F., Psaltis, D., & Chan, C.K. 2016, ApJ, 826, 77 [NASA ADS] [CrossRef] [Google Scholar]
 Bardeen, J. M., Press, W. H., & Teukolsky, S. A. 1972, ApJ, 178, 347 [NASA ADS] [CrossRef] [Google Scholar]
 Bartko, H., Martins, F., Fritz, T. K., et al. 2009, ApJ, 697, 1741 [NASA ADS] [CrossRef] [Google Scholar]
 Broderick, A. E., & Loeb, A. 2005, MNRAS, 363, 353 [NASA ADS] [CrossRef] [Google Scholar]
 Broderick, A. E., & Loeb, A. 2006, ApJ, 636, L109 [NASA ADS] [CrossRef] [Google Scholar]
 Calderón, D., Cuadra, J., Schartmann, M., Burkert, A., & Russell, C. M. P. 2020, ApJ, 888, L2 [NASA ADS] [CrossRef] [Google Scholar]
 Carter, B. 1968, Phys. Rev., 174, 1559 [NASA ADS] [CrossRef] [Google Scholar]
 Chan, C.K., Psaltis, D., Özel, F., et al. 2015, ApJ, 812, 103 [NASA ADS] [CrossRef] [Google Scholar]
 Cunningham, C. T., & Bardeen, J. M. 1973, ApJ, 183, 237 [NASA ADS] [CrossRef] [Google Scholar]
 Dexter, J., & Agol, E. 2009, ApJ, 696, 1616 [NASA ADS] [CrossRef] [Google Scholar]
 Dexter, J., & Fragile, P. C. 2013, MNRAS, 432, 2252 [NASA ADS] [CrossRef] [Google Scholar]
 Dolence, J. C., Gammie, C. F., Mościbrodzka, M., & Leung, P. K. 2009, ApJS, 184, 387 [NASA ADS] [CrossRef] [Google Scholar]
 Eisenhauer, F., Perrin, G., Rabien, S., et al. 2008, in The Power of Optical/IR Interferometry: Recent Scientific Results and 2nd Generation, eds. A. Richichi, F. Delplancke, F. Paresce, & A. Chelli, 431 [CrossRef] [Google Scholar]
 Gillessen, S., Plewa, P. M., Eisenhauer, F., et al. 2017, ApJ, 837, 30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gravity Collaboration (Abuter, R., et al.) 2017, A&A, 602, A94 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gravity Collaboration (Abuter, R., et al.) 2018a, A&A, 618, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gravity Collaboration (Abuter, R., et al.) 2018b, A&A, 615, L15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gravity Collaboration (Abuter, R., et al.) 2019, A&A, 625, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hamaus, N., Paumard, T., Müller, T., et al. 2009, ApJ, 692, 902 [NASA ADS] [CrossRef] [Google Scholar]
 Hora, J. L., Witzel, G., Ashby, M. L. N., et al. 2014, ApJ, 793, 120 [Google Scholar]
 Hughes, S. A. 2001, Phys. Rev. D, 64, 064004 [NASA ADS] [CrossRef] [Google Scholar]
 Meyer, L., Eckart, A., Schödel, R., et al. 2006, A&A, 460, 15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Meyer, L., Do, T., Ghez, A., et al. 2009, ApJ, 694, L87 [NASA ADS] [CrossRef] [Google Scholar]
 Paumard, T., Perrin, G., Eckart, A., et al. 2005, Astron. Nachr., 326, 568 [NASA ADS] [Google Scholar]
 Paumard, T., Genzel, R., Martins, F., et al. 2006, ApJ, 643, 1011 [NASA ADS] [CrossRef] [Google Scholar]
 Rauch, K. P., & Blandford, R. D. 1994, ApJ, 421, 46 [NASA ADS] [CrossRef] [Google Scholar]
 Ressler, S. M., Tchekhovskoy, A., Quataert, E., & Gammie, C. F. 2017, MNRAS, 467, 3604 [NASA ADS] [CrossRef] [Google Scholar]
 Ressler, S. M., Quataert, E., & Stone, J. M. 2018, MNRAS, 478, 3544 [NASA ADS] [CrossRef] [Google Scholar]
 Rowan, M. E., Sironi, L., & Narayan, R. 2017, ApJ, 850, 29 [NASA ADS] [CrossRef] [Google Scholar]
 Rowan, M. E., Sironi, L., & Narayan, R. 2019, ApJ, 873, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Sironi, L., & Narayan, R. 2015, ApJ, 800, 88 [NASA ADS] [CrossRef] [Google Scholar]
 Trippe, S., Paumard, T., Ott, T., et al. 2007, MNRAS, 375, 764 [NASA ADS] [CrossRef] [Google Scholar]
 Witzel, G., Eckart, A., Bremer, M., et al. 2012, ApJS, 203, 18 [NASA ADS] [CrossRef] [Google Scholar]
 Witzel, G., Martinez, G., Hora, J., et al. 2018, ApJ, 863, 15 [NASA ADS] [CrossRef] [Google Scholar]
 Yang, X.L., & Wang, J.C. 2014, A&A, 561, A127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Yelda, S., Ghez, A. M., Lu, J. R., et al. 2014, ApJ, 783, 131 [NASA ADS] [CrossRef] [Google Scholar]
 Yuan, F., Quataert, E., & Narayan, R. 2003, ApJ, 598, 301 [NASA ADS] [CrossRef] [Google Scholar]
 Zamaninasab, M., Eckart, A., Witzel, G., et al. 2010, A&A, 510, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
All Figures
Fig. 1. Sample image of a hotspot at different points in the orbit. Top row: image of the hotspot before, during, and after a strong lensing event. Topright plot: image averaged over the whole orbit. Lower left panel: lightcurve produced over several orbits, lower right panel: motion of the centroid of emission. This figure was produced with a hotspot in a circular orbit at R = 5M around a black hole with spin a = 0.9. The spin axis of the black hole is inclined i = 15° from the observer’s line of sight and the relative inclination between the spin axis and the orbital plane is θ = 50.54°. 

In the text 
Fig. 2. Data and NERO fits to the flares from May 27, July 22, and July 28, 2018. Data points correspond to the centroid position observed by GRAVITY over the course of the flare. The solid lines show the orbits that best fits each of three flares. Colors indicate the time since the start of the flare. 

In the text 
Fig. 3. Confidence contours corresponding to 1, 2, and 3σ limits on the orbital radius and inclination for all three flares. The last panel shows the combined constraints, assuming the radius, inclination, and position angle are the same for all three flares. 

In the text 
Fig. 4. Confidence contours showing the allowed values for the orbital radius and the xoffset of the center of the black hole for the July 22, 2018 flare. An offset of zero corresponds to the black hole lying at the mean location of emission averaged over the entire flare. 

In the text 
Fig. 5. Observed and inferred intrinsic lightcurves for the three GRAVITY flares. The flux is measured relative to the nearby star S2. Each dotted line corresponds to the lightcurve during the respective flare as observed by GRAVITY. The solid line shows the intrinsic lightcurve once all relativistic effects predicted from the bestfit orbital model have been remove. 

In the text 
Fig. 6. Combined constraints together with information from the polarization and the Doppler beaming. The blue astrometric constraints are as shown in Fig. 3. The green shaded region corresponds to models where the maximum contrast of the lightcurve is less than 2.5. The hatched polarization period constraint arises from the assumption that the change in polarization angle during the flare is caused by the motion of the hotspot through the magnetic field. The polarization loop constraint is based on a model with a vertical magnetic field and the requirement that the polarization angle matches the behavior observed during the flare. 

In the text 
Fig. 7. Combined constraints on the orientation of the hotspot orbit from all three flares. The black marker indicates the approximate orientation of the clockwise disk of stars. 

In the text 
Fig. 8. 1 and 2σ constraints on the orbital radius and spot diameter for the July 22 flare. 

In the text 
Fig. 9. Constant χ^{2} contours corresponding the 1, and 2σ regions of the vertical velocity and inclination for the July 22 flare. We note that in this case we use a flat prior for the inclination of the orbit. 

In the text 
Fig. 10. Constraints from simulated data with improved error bars. In each case, the underlying model is based on a hotspot orbiting at R = 8.0M with an inclination of i = 155°, as marked by the red cross in each panel. The astrometric error is 10 μas, 20 μas, and 30 μas, respectively. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.