Triaxial orbit-based modelling of the Milky Way Nuclear Star Cluster

Triaxial orbit-based modelling of the Milky Way Nuclear Star Cluster


We construct triaxial dynamical models for the Milky Way nuclear star cluster using Schwarzschild’s orbit superposition technique. We fit the stellar kinematic maps presented in Feldmeier et al. (2014). The models are used to constrain the supermassive black hole mass M, dynamical mass-to-light ratio , and the intrinsic shape of the cluster. Our best-fitting model has M = () , = (0.90) /, and a compression of the cluster along the line-of-sight. Our results are in agreement with the direct measurement of the supermassive black hole mass using the motion of stars on Keplerian orbits. The mass-to-light ratio is consistent with stellar population studies of other galaxies in the mid-infrared. It is possible that we underestimate M and overestimate the cluster’s triaxiality due to observational effects. The spatially semi-resolved kinematic data and extinction within the nuclear star cluster bias the observations to the near side of the cluster, and may appear as a compression of the nuclear star cluster along the line-of-sight. We derive a total dynamical mass for the Milky Way nuclear star cluster of M = (2.1)  10 within a sphere with radius = 2 = 8.4 pc. The best-fitting model is tangentially anisotropic in the central = 0.5-2 pc of the nuclear star cluster, but close to isotropic at larger radii. Our triaxial models are able to recover complex kinematic substructures in the velocity map.

Galaxy: center; kinematics and dynamics.

1 Introduction

The Milky Way nuclear star cluster is the ideal object to study the dynamics of a stellar system around a supermassive black hole. At a distance of 8 kpc it is close enough to resolve the individual stars, and measure discrete velocities in three dimensions. Modelling the stellar kinematics can constrain the mass distribution of the star cluster, and reveal the presence of a central dark massive object. In the special case of our own Galaxy, it it possible to observe Keplerian orbits of stars around a dark, point-mass-like object in the Galactic centre. These observations constrain this dark object to be a supermassive black hole with a mass of (4.1 0.6)  (Ghez et al., 2008), (4.3 0.39)  (Gillessen et al., 2009), or (4.02 0.20)  (Boehle et al., 2016). Unfortunately, similar high-resolution observations are not yet possible in other galaxies.

Already in the 1970s the requirement of a central supermassive black hole in the Galactic centre was discussed to explain observational data (e.g. Oort, 1977). Several studies used stellar radial velocities to constrain the mass distribution in the Galactic centre (e.g. Rieke & Rieke, 1988; McGinn et al., 1989; Sellgren et al., 1990; Haller et al., 1996; Genzel et al., 1996). Also stellar proper motions were used to study the Galactic centre mass distribution (Schödel et al., 2009). Several studies combined radial velocity and proper motion data (Trippe et al., 2008; Do et al., 2013; Fritz et al., 2016). The mass distribution was derived using the spherical Jeans (1922) equations or the projected mass estimators of Bahcall & Tremaine (1981) for spherical systems. These studies found that a central dark mass of 25  10 is required to explain the observations.

Together with the increase of observational data, also the modelling became more advanced. Trippe et al. (2008) included the rotation of the nuclear star cluster in the modelling, although the rotation velocity of their data was too high (Schödel et al., 2009; Feldmeier et al., 2014). Feldmeier et al. (2014) and Chatzopoulos et al. (2015a) studied the Milky Way nuclear star cluster using axisymmetric Jeans models. Chatzopoulos et al. (2015a) showed the advantages of axisymmetric models over spherical Jeans models, which cannot explain the observed asymmetry of the velocity dispersion of proper motions parallel and perpendicular to the Galactic plane. The nuclear star cluster appears to be flattened in its light distribution (Schödel et al., 2014) as well as in the kinematics (Chatzopoulos et al., 2015a). Most studies showed that the nuclear star cluster kinematics is in agreement with isotropy (Schödel et al., 2009; Do et al., 2013; Chatzopoulos et al., 2015a), although the uncertainties are quite large (e.g. Do et al. 2013). All these models assumed a constant mass-to-light ratio for the light distribution of the cluster.

In this study we relax the assumption of axisymmetry and use triaxial orbit-based Schwarzschild (1979) models. Orbit-based models make no assumptions on the velocity anisotropy of the stellar motions, as Jeans models do. Further, the higher moments of the kinematics can also be included (Rix et al., 1997), which is important to break the degeneracy of mass and anisotropy in dynamical models.

Orbit-based models are commonly used to analyse line-of-sight velocity data of other galaxies (e.g. van der Marel et al., 1998; Gebhardt et al., 2000; Valluri et al., 2005; van den Bosch et al., 2008), and are an excellent tool to detect and measure the masses of supermassive black holes and dark matter halos. For extragalactic systems, the data are usually obtained from integrated light observations. Each data point contains the accumulated kinematics of many stars, weighted by their respective brightness. However, modelling the dynamics of integrated light data may be prone to systematic uncertainties, and bias the results of the central black hole mass. Therefore, it is interesting to test dynamical models on systems for which we know the central black hole mass from other independent measurements. The Milky Way nuclear star cluster is a good object for this kind of test. Also megamaser disc galaxies are useful to validate stellar dynamical black hole measurements. Black hole mass measurements from megamasers are very precise with uncertainties of only about 10 per cent. However, there is currently only one megamaser disc galaxy with a stellar dynamical black hole mass measurement (van den Bosch et al., 2016), NGC 4258. Different dynamical studies found either a 15 per cent lower or a 25 per cent higher black hole mass than the maser measurement (Siopis et al., 2009; Drehmer et al., 2015).

We use the triaxial orbit-based code by van den Bosch et al. (2008) to model the light distribution and line-of-sight kinematics of the Milky Way nuclear star cluster. We use the spectroscopic data cube constructed by Feldmeier et al. (2014) for the kinematic data, and derive a surface brightness distribution using Spitzer 4.5  and NACO band images. We assume a galactocentric distance of 8 kpc (Malkin, 2012) and a position angle 3140 East of North (J2000.0 coordinates, Reid & Brunthaler 2004) with respect to the Galactic plane. This paper is organised as follows: We describe the kinematic and photometric data in Section 2. The dynamical models are introduced in Section 3. Section 4 discusses the results, and Section 5 summarizes the main conclusions.

Figure 1: Kinematic data (top row) and respective uncertainties (bottom row). The columns denote velocity , velocity dispersion , Gauss-Hermite moments , and . White pixels are due to excluded bright stars.

2 Description of the Data

2.1 Kinematic data

The line-of-sight velocity distribution (LOSVD) provides constraints on the dynamical structure of stellar systems. To extract this information, we used the near-infrared band spectroscopic data cube of Feldmeier et al. (2014), which has a pixel scale of 2.22 arcsec pixel. We used the data cube that was cleaned of foreground stars and bright stars. The cleaned data cube contains the light of the old red giant star population.

We fitted the LOSVD as in Feldmeier et al. (2014), i.e. on the stellar CO absorption lines (2.29022.365 ) with the IDL routine pPXF (Cappellari & Emsellem, 2004) and the high resolution spectra of Wallace & Hinkle (1996) as template stars. We applied the same spatial binning as Feldmeier et al. (2014), resulting in 175 spatial bins. Feldmeier et al. (2014) fitted only the velocity and velocity dispersion . However, we fitted in addition also higher moments of the LOSVD, in particular the Gauss-Hermite parameters and . We added noise to each of the 175 integrated light spectra in 100 Monte Carlo simulation runs and obtained a distribution for each moment of the LOSVD. The mean and standard deviation of the Monte Carlo distribution are taken as measurement and 1 uncertainty of the kinematics.

Since the Milky Way nuclear star cluster is at a distance of only 8 kpc, the spectroscopic observations are spatially semi-resolved. Bright stars can be resolved individually, and contribute a large fraction of the flux. For that reason we used the cleaned data cube of Feldmeier et al. (2014), where bright stars were excluded. However, the kinematic maps still show stochastic shot noise. As a consequence, the difference of the kinematics in adjacent bins can be higher than their uncertainties, which causes problems when we model the kinematics. The stochastic noise can be mistaken for signal, and this means the best fit will be achieved by modelling the shot noise. To prevent this, we increased our kinematic uncertainties such that the difference of the measurement in two adjacent bins () is less than the sum of their uncertainties (). We did this for the velocity , velocity dispersion , , and , and find that it is required for about 68 per cent of the kinematic data uncertainties. Additionally, we point-symmetrised the kinematics using the procedure of van den Bosch & de Zeeuw (2010). The median uncertainties or , , , and are 24.6 km s, 18.4 km s, 0.15, and 0.17. Our resulting kinematic maps are consistent with the maps of Feldmeier et al. (2014). We find rotation in the velocity map of approximately 50 km s and an increase in the velocity dispersion from about 65 km s towards =135 km s at the centre. The kinematic maps are shown on the top row of Fig. 1, the uncertainties are shown on the bottom row.

2.2 Imaging data and surface brightness distribution

The light distribution of the nuclear star cluster traces the stellar density. We require the two-dimensional light distribution of the red giant stars, which are our kinematic tracers. The extinction is high at optical wavelengths in the Galactic centre ( 30 mag, Scoville et al., 2003; Gao et al., 2013), therefore we used near- and mid-infrared images.

For the central 40.4 arcsec  40.4 arcsec (1.6 pc  1.6 pc) we used the high-resolution NACO -band mosaic of Schödel et al. (2009), which has a spatial scale of 0.027 arcsec pixel. We preferred the band over the band in order to avoid light from gas emission lines in the band (Br  and He I, Paumard et al., 2004). Our kinematic tracers are cool late-type stars, but there are also more than 100 hot, young stars located in the centre of the cluster, within a projected radius  = 0.5 pc (12.8 arcsec, Paumard et al. 2006). We masked out the young stars from the image with a 15 pixel radius. For the bright red supergiant IRS 7 we used a larger mask with a 30 pixel radius. Beyond the central 0.5 pc, the nuclear star cluster light is dominated by cool stars, and the contribution of young stars is negligible (Feldmeier-Krause et al., 2015).

For the large-scale light distribution, we used Spitzer IRAC images (Stolovy et al., 2006). These images were corrected for dust extinction and PAH emission by Schödel et al. (2014). We used the extinction and emission corrected 4.5  image to measure the light distribution. The image was smoothed to a scale of 5 arcsec pixel, and extends over 270 pc  200 pc. We excluded a central circle with  = 0.6 pc (15.4 arcsec) to avoid contribution from ionised gas emission and young stars. In addition we masked out the young Quintuplet star cluster (Figer et al., 1999), and the dark 20-km s-cloud M-0.13-0.08 (García-Marín et al., 2011).

We used the MGE_FIT_SECTORS package (Cappellari, 2002) to derive the surface brightness distribution. The Multi-Gaussian-Expansion model (Emsellem et al., 1994) has the advantage that it can be deprojected analytically. We measured the photometry of the two images along the major axis and the minor axis. We assumed that the cluster’s major axis is aligned along the Galactic plane, as found by Schödel et al. (2014), and constant. The centre is the position of Sgr A*, which is the radio source associated with the Galactic centre supermassive black hole. We fitted a scale factor to match the photometry of the two images in the region where they overlap (1627.8 arcsec). Then we measured the photometry on each image along 12 angular sectors, and converted the NACO photometry to the Spitzer flux. Assuming four-fold symmetry, the measurements of four quadrants are averaged on elliptical annuli with constant ellipticity. Using the photometric measurements of the two images, we fitted a set of two-dimensional Gaussian functions, taking the point-spread-function (PSF) of the NACO image into account.

A comparison with the surface brightness profile of Fritz et al. (2016, their Fig. 2) showed that our profile is steeper in the central 30 arcsec. A possible reason is the small overlap region of the Spitzer and NACO images, and that the Spitzer flux could be too high at the centre. Maybe the PAH emission correction of the Spitzer image was too low. The mid-infrared dust emission is significant out to almost one arcmin. Fritz et al. (2016) used NACO - and -band images in the central  = 20 arcsec. Out to 1 000 arcsec (39 pc) they used Hubble Space Telescope WFC3 data (M127 and M153 filters) and public VISTA Variables in the Via Lactea Survey images ( and bands, Saito et al. 2012). We lowered the intensities of the central Gaussians by scaling our averaged profile to the one-dimensional flux density profile of Fritz et al. (2016). As a result the amplitudes of the inner Gaussians become smaller, but the outer Gaussians (>100 arcsec 4 pc) are nearly unchanged. We list the components of the Multi-Gaussian Expansion in Table 1 and plot the surface brightness profile in Fig. 2 (upper panel). We also show the projected axial ratio as a function of radius in the lower panel of Fig. 2. Out to the central 1 pc, is increasing from 0.4 to 0.7. Schödel et al. (2014) found a mean axial ratio of 0.710.02 for the nuclear star cluster. This is in agreement with our maximum value of . However, decreases at larger radii, as the contribution from the nuclear stellar disc becomes more important and the light distribution therefore flatter.

Table 1: The Multi-Gaussian-Expansion (MGE) fit parameters for the 4.5  Spitzer/IRAC dust extinction and PAH emission corrected image in combination with the NACO -band mosaic scaled to Spitzer flux. is the peak surface brightness used in the dynamical modelling, is the standard deviation, and is the axial ratio of the Gaussian components. is the peak surface brightness before scaling to Fritz et al. (2016).
Figure 2: Upper Panel: Surface brightness profile derived from a dust extinction and PAH emission corrected Spitzer/IRAC 4.5 m image and NACO -band mosaic for the centre, scaled to the measurements of Fritz et al. (2016, blue crosses). The black full line denotes the MGE fit along the major axis, and the red dashed line along the minor axis. Lower panel: Projected axial ratio as a function of .

We note that there are three main differences with the surface brightness distribution derived by Feldmeier et al. (2014): (1) We used an -band instead of a -band NACO image to avoid ionised gas emission; (2) We masked young stars in the NACO image to match the distribution of stars used as kinematic tracers; and (3) We scaled the central photometry to the flux density data of Fritz et al. (2016) to avoid a possible overestimation of the central flux when scaled to the Spitzer image. All three changes influence only the central part of the surface brightness distribution, as ionised gas emission and light from young stars are only important in the central parsec.

The surface brightness profile is deprojected to obtain the three-dimensional spatial stellar light distribution. In general, the deprojection of a surface brightness profile is non-unique (e.g. Rybicki, 1987; Franx, 1988), and our deprojection is only one possible solution. The MGE deprojection produces smooth intrinsic densities, which are in agreement with the photometric observations (Cappellari, 2002).

3 Dynamical Models of the Milky Way nuclear star cluster

3.1 Schwarzschild’s method

Orbit-based models or Schwarzschild models are a useful tool to model the dynamics of stellar systems by orbit superposition. The first step of Schwarzschild’s method is to integrate the equations of motion for a representative library of stellar orbits in a gravitational potential . Then the observables for each orbit are computed, considering projection, PSF convolution and pixel binning. The next step is to find orbital weights to combine the orbits such that they reproduce the observed data. Schwarzschild models are a powerful tool to recover the intrinsic kinematical structure and the underlying gravitational potential (Schwarzschild, 1979; van de Ven et al., 2008; van den Bosch & van de Ven, 2009). We refer the reader for further details to van den Bosch et al. (2008) for implementation and van de Ven et al. (2008) for verification of the triaxial Schwarzschild code.

Mass model

We calculated orbits in the combined gravitational potential of a supermassive black hole and the star cluster , inferred from the imaging data. As we run triaxial models, there are three intrinsic shape parameters, , , and , for the cluster. The shape parameters characterise the axial ratios for the long, intermediate and short axes , , and . They are defined as , , and , where is the length of the longest axis projected on the sky. Thus, represents the compression of due to projection on the sky. Each set of axial ratios refers to a set of viewing angles (, see also van den Bosch et al. 2008). The surface brightness distribution is deprojected given the intrinsic shape parameters , and multiplied with the dynamical mass-to-light ratio to get the intrinsic stellar mass density . From Poisson’s equation one calculates the gravitational potential. We did this for different values of the black hole mass M, dynamical mass-to-light ratio , and different shape parameters. In total our model has five free parameters, M, , , , and .

Besides the considered stellar population and the supermassive black hole, there are other components within the nuclear star cluster, which we neglected. We measured a dynamical mass-to-light ratio, which combines the stellar mass-to-light ratio with other components. These components are the young stars, ionised gas, neutral gas, and dark matter. The young stars are at a distance of about 0.5 pc from the supermassive black hole. The lower limit of the total mass of young stars is 12 000 . However, the total enclosed extended mass in the same region is 10 (Oh et al., 2009; Feldmeier et al., 2014), and the mass of the supermassive black hole is 4  10. The mass of the young stars is therefore probably negligible. The hot ionised gas has a mass of only a few 100  (Ferrière, 2012), and cannot influence the stellar dynamics significantly. The neutral gas in the circum-nuclear disc may contribute more mass, estimates range from 10 (Etxaluze et al., 2011; Requena-Torres et al., 2012) to 10 (Christopher et al., 2005), though this is probably the upper limit (Genzel et al., 2010). The circum-nuclear disc extends over a distance of about 1 pc to more than 5 pc from the centre. At 5 pc, the total enclosed mass is 10 (McGinn et al., 1989; Feldmeier et al., 2014). We decided to neglect the mass distribution of the circum-nuclear disc in our dynamical models, since it is very uncertain, and makes up only 0.1 to 10 per cent of the enclosed mass. The contribution of dark matter to the nuclear star cluster mass is also neglected. Linden (2014) showed that the fraction of dark matter in the central 100 pc of the Milky Way is about 6.6 per cent, assuming the traditional dark matter profile of Navarro et al. (1996).

Orbit library

The orbit library should be as general as possible and representative for the gravitational potential. We assumed that the orbits are regular and that three integrals of motion, , , and , are conserved. The orbit families consist of box orbits, which can cross the centre and have an average angular momentum of zero, and three types of tube orbits, which avoid the centre. The tube orbits are divided in short-axis-tube orbits, which have non-zero mean angular momentum around the short axis, outer and inner long-axis-tube orbits, which have non-zero mean angular momentum around the long axis. The orbit grid should sample the entire phase space. It has to be dense enough to suppress discreteness noise, but integration has to be done in a reasonable amount of computing time.

Figure 3: Comparison of the observed stellar surface brightness and kinematics (top row) and the best-fitting Schwarzschild model. The columns denote surface brightness, velocity , velocity dispersion , Gauss-Hermite moments , and .
Figure 4: The surface brightness and kinematics of the best-fitting model (red filled diamond symbols) and the data (open black squares and error bars) as a function of radius. Data points in the Galactic East are plotted at negative radii, Galactic West at positive radii. From top to bottom: surface brightness, velocity , velocity dispersion , and higher Gauss-Hermite moments and .

We followed van den Bosch et al. (2008) and sampled the orbit energy using a logarithmic grid in radius. Each energy is linked to the radius by calculating the potential at = (, 0, 0). We sampled = 35 energies calculated from in logarithmic steps ranging from = to = , i.e. 3.16 arcsec to 4.4 degree or 0.12 pc to 616.5 pc. We note that the outer radius is about 3.5 times the outermost Gaussian of the MGE fit. We tested lower values of the inner radius but found consistent results. For each energy, the starting point of an orbit was selected from a linear grid over 14 values each. For details on the orbit sampling we refer to van den Bosch et al. (2008). In total, we have = 6860 orbits. Each orbit was integrated over 200 periods, and sampled on 100 000 points per orbit. For each orbit we stored the intrinsic and projected properties. The projected orbits are stored in a grid, with PSF convolution and pixel size of the observed data taken into account. The velocities were stored in 183 bins between 7.4  and +7.4 . These numbers guarantee a proper sampling of the observed velocity profiles (Cretton et al., 1999).

Solving the orbital weight distribution

The model has to fit the kinematic data, the intrinsic and the projected mass distribution. The fit was done by finding a linear combination of the orbits, and solving for orbital weights . Each orbital weight corresponds to a mass on the respective orbit , and the weights are therefore non-negative. We used the non-negative least-squares (NNLS) algorithm of Lawson & Hanson (1974) which was also used by Rix et al. (1997), van der Marel et al. (1998), and Cretton et al. (1999). One of the fitting constraints is to make sure that the model is self-consistent. It is required that the orbit superposition reproduces the intrinsic and projected aperture masses within two per cent, which is the typical accuracy of the observed surface brightness (van den Bosch et al., 2008).

3.2 Constraining the input parameters

We ran 4899 models with different parameter combinations of M, , . The black hole mass M was sampled in logarithmic steps of 0.2 from 5.5 to 7.5, starting with 6.3 (i.e. M 2). The mass-to-light ratio was linearly sampled between 0.1 and 2.0 with steps of 0.04, with a starting value of 0.6 (in units of ). The starting model had () = (0.29, 0.84, 0.99). We sampled different combinations of () with a step size of (0.01, 0.02, 0.01), with the boundaries 0.05 < <0.29, 0.40 < <0.99, and 0.70 < 0.99. When fitting the orbital weights, we already applied the self-consistency criteria to make sure that the photometry is fitted with a high accuracy of at least two per cent for each model. The best fit of the five parameters M, , , however, was found by calculating the of the kinematic moments. After the starting model we computed different combinations of (M, , ) in the grid. We used a analysis to find regions with good fits and started new models, with a finer sampling at low . This was done in several iterations, and it was not necessary to compute each point in the (M, , ) grid.

We measured four kinematic moments in 175 spatial bins of the spectroscopic map in Sect. 2.1. In each bin we measured the velocity , velocity dispersion , and the higher Gauss-Hermite moments , and . Thus, the number of observables is the number of spatial bins times the number of kinematic moments, in our case N = 175 4 = 700. We minimised the function


where and denote the kinematic measurements, the respective measurement uncertainties, , , , and the model kinematics, to find the best-fitting model.

3.3 Modelling results

The best-fitting model

Figure 5: Illustration of the fitted parameter space, was computed with Eq. 1. Each symbol denotes a model, the coloured symbols are models with < = 37.3, black diamonds are models with > . The 1, 2, and 3 colours corresponding to = 5.9, 11.3 and 18.2, are denoted. The black cross denotes the best-fitting model.
Figure 6: The values plotted against the five free parameters (M, , , , ). was computed with Eq. 1. The best-fitting model is denoted as blue asterisk, the 1, 2, and 3 confidence levels, corresponding to = 5.9, 11.3, and 18.2, are denoted as blue lines. The red line denotes = 37.3.

Our best-fitting parameters are M = , = 0.90, = 0.28, = 0.64, = 0.99. This corresponds to best-fitting viewing angles = 80°, = 79°, = 91°. We show the surface brightness map and the symmetrised kinematic maps in Fig. 3. The upper row are the data, the lower row are the maps of the best-fitting model. The misalignment of the kinematic rotation axis with respect to the photometry, and the perpendicular rotating substructure at 20 arcsec (0.8 pc) found by Feldmeier et al. (2014) are well reproduced in the model velocity map. We also show the data and the best-fitting model in Fig. 4, the panels denote surface brightness, , , and . The surface brightness map is reproduced within one per cent. The highest and lowest values of , , , and the lowest values of in the data have higher absolute values than in the best-fitting model, but are consistent within their uncertainties. The best fit has = 290. With M = 5 fitted parameters and N = 4175=700 observational constraints, this means = 0.42. That is less than one is partially due to the large uncertainties of the kinematics, and the fact that the kinematic measurements are correlated.

We illustrate the distribution of for the 4899 models in Fig. 5. We plot each combination of parameters. Red colours denote low , bluer, smaller symbols denote high . The black cross denotes the best-fitting model. The observed projected flattening of the surface brightness profile constrains the viewing angles and thus also the intrinsic shape parameters. In particular, a flat 1 means that the stellar system is observed along one of the principal planes (van den Bosch et al., 2008). We obtained as lowest value = 0.30, this limits the possible projections and puts an upper limit on the parameter . Likewise, the value of = 0.99 is the boundary value of the grid. The values of and denote the respective minimum values of the radially varying intrinsic axial ratios , and over the entire radial range of the photometry, i.e. the nuclear stellar disc and the embedded nuclear star cluster. The upper left panel of Fig. 5 shows that for each value of , the best-fitting is approximately 0.90. A similar behaviour is found with and . There is only a slight increase of the best-fitting with higher values of . At the same time, the best-fitting values of M do not show a strong dependence on or (second row) in the allowed parameter ranges. The intrinsic shape parameters do not influence our best fit for M, as this measurement is mostly made from the inner bins and the outer bins contribute little. The outer bins, however, contribute to the intrinsic shape fit. The supermassive black hole mass and the dynamical mass-to-light ratio are correlated. For higher values of , a lower M fits the data.

We show how depends on the different parameters in Fig. 6. The best-fitting model, which has the lowest , is marked as blue asterisk symbol. The blue lines denote the 1, 2, and 3 confidence limits, corresponding to = 5.9, 11.3, and 18.2. The red line illustrates the standard deviation of itself, i.e. = 37.3, where N = 700, and M = 5. This value was used as confidence limit by van den Bosch & van de Ven (2009). In Table 2 we list the 1 and 3 uncertainties.

parameter best fit 1 3 unit
Table 2: The best-fitting model results and the 1 and 3 uncertainties, corresponding to = 5.9 and 18.2.
Figure 7: Enclosed total mass in spherical shells as a function of radius, in units of and in logarithmic scaling. The black line denotes the enclosed mass with = 0.90 and M = , the grey shaded contours the 3 uncertainty. The horizontal line denotes a supermassive black hole with the mass M = 4 . The vertical, dotted line denotes the outer edge of the kinematic data, the vertical, solid line the effective radius. We also plot the results for the enclosed mass from previous studies. We scaled the masses to a distance = 8.0 kpc, if the study assumed a different Galactocentric distance: McGinn et al. (1989, diamonds, assumed = 8.5 kpc), Lindqvist et al. (1992, upward triangles, = 8.5 kpc), Deguchi et al. (2004, squares), Trippe et al. (2008, x-symbol), Oh et al. (2009, leftfacing triangles), Schödel et al. (2014, asterisk), Feldmeier et al. (2014, blue dashed line), Chatzopoulos et al. (2015a, rightfacing triangle, = 8.3 kpc), and Fritz et al. (2016, downward triangle, = 8.2 kpc).
Figure 8: Intrinsic axial ratios as a function of radius . The blue dot-dashed line denotes , the black solid line , the green dashed line the triaxiality , the shaded regions the 1 uncertainties, the vertical, solid line = 4.2 pc.

Mass profile and intrinsic shape

We show the enclosed total mass as a function of the projected radius in Fig. 7, grey shaded contours are the 3 uncertainties. The mass was computed within spherical shells. We also plot the results of various other studies. Most studies assumed a spherical shape of the cluster, Feldmeier et al. (2014) and Chatzopoulos et al. (2015a) assumed axisymmetry. Some of the studies used also different Galactocentric distances, so we scaled the masses using = 8.0 kpc. Our results are in agreement with several studies in the central 100 arcsec. However, we obtain a lower enclosed mass than Trippe et al. (2008) and Oh et al. (2009), who used the Jeans equation for a spherical system to obtain the enclosed mass. At larger radii 400 arcsec (15.5 pc) beyond the reach of our kinematic data, we obtained a higher mass than Lindqvist et al. (1992), but we are in agreement at =750–1600 arcsec. Their data extend to larger radii, but their assumption of spherical symmetry does no longer hold at such large radii.

The enclosed total mass within a sphere with the radius = 8.4 pc, i.e. about two times the effective radius of the nuclear star cluster, is M = (2.10.7)  10, here we give the 3 uncertainty.

The black hole influences the stellar kinematics only at the centre of the nuclear star cluster. Out to = 53 arcsec (2 pc), the best-fitting mass of the black hole (M = ) is higher than the enclosed stellar mass of our best-fitting model. Merritt (2004) defined the radius of influence of a black hole as the radius where the enclosed stellar mass equals two times the black hole mass. With this definition and a black hole mass of 4 , we obtain = (104) arcsec, i.e. approximately (4.0) pc. This result is in agreement with Alexander (2005), who found = 3 pc. The kinematic measurements at larger radii have little influence on the black hole mass measurement, but are important to constrain the orbital structure and dynamical mass-to-light ratio.

The shape of the nuclear star cluster is illustrated in Fig. 8. We show the intrinsic axial ratios and as a function of radius . The axial ratio is low in the centre (=0.3), increases to =0.8 at 35 arcsec, and then decreases to =0.6 at =150 arcsec. The best-fitting shape parameter, =0.28, is approximately the central value. The axial ratio is also low in the centre (0.65), increases to =0.9 at =40 arcsec, and then decreases to =0.75 at =150 arcsec. Though, the uncertainty of is rather high, and the decrease not significant. Also the best-fitting axial ratio, =0.64, is close to the central value. We plot the triaxiality , it varies between 0.45 and 0.7. The triaxiality increases from =70 arcsec to the outer radius =150 arcsec. This may be a signature of the increasing influence of the triaxial Galactic bulge at larger radii (Tsatsi et al., 2017). However, given the large uncertainty of , this increase is not significant. A constant or a decreasing are not excluded by the large uncertainties.

Figure 9: Top: Anisotropy as a function of radius . Negative values denote tangential anisotropy, positive values radial anisotropy. Bottom: Fraction of stellar mass per orbit type as a function of average orbit radius . The green, solid line denotes long-axis-tube orbits; the blue, dot-dashed line short-axis-tube orbits; the red, dashed line box orbits; the vertical, solid line = 4.2 pc.
Figure 10: Orbit density with angular momentum (top), i.e. rotation around the short axis, and (bottom), i.e. rotation around the long axis as a function of average orbit radius . Dark, blue colour indicates higher orbit density. The vertical line denotes = 4.2 pc.

Internal dynamics

The best-fitting model has tangential anisotropy in the centre of the cluster. The value of the anisotropy = is negative, where is the tangential velocity dispersion and is the radial velocity dispersion. We show the anisotropy as a function of radius in Fig. 9, top panel. We plot the mean anisotropy of the models within the 1 uncertainty limit. The uncertainty of is about 0.1. The plot extends to the outer edge of the kinematic data at 150 arcsec. The cluster kinematics becomes nearly isotropic at radii >70 arcsec.

We show the angular momentum distribution of the orbits in Fig. 10. The colours denote the density of orbits passing radius with mean angular momentum (top panel) or (bottom panel). The plot of denotes rotation about the short z-axis. Orbits with 0 are contributed by short-axis-tube orbits, while long-axis-tube orbits have = 0. On the other hand, denotes rotation about the long x-axis (bottom panel), and orbits with 0 are contributed by long-axis-tube orbits. Short-axis-tube orbits have = 0. Long-axis-tube orbits are most important in the central 20–60 arcsec and at larger radii 80 arcsec. Short-axis-tube orbits, which contribute in total more mass than long-axis-tube orbits, are most important at = 60–140 arcsec. We illustrate the distribution of the stellar mass on the different orbit types also in Fig. 9 (bottom panel) as a function of radius. Most stars (>50 per cent) are on short-axis-tube orbits, i.e. they orbit the minor axis. Long-axis-tube orbits contribute about 40 per cent to the luminosity and thus stellar mass in the central 30 arcsec. They produce the perpendicular rotating substructure at  20 arcsec (0.8 pc) found by Feldmeier et al. (2014). At larger radii, long-axis-tube orbits contribute only about 30 per cent to the stellar mass. Box orbits contribute little mass in the centre (<10 per cent), but their relative mass increases towards larger radii. At = 150 arcsec (5.8 pc), they contribute 20 per cent.

4 Discussion

4.1 Difference of the resulting black hole mass

The currently best black hole mass estimate is (4.02 0.20) , derived by Boehle et al. (2016) from Keplerian stellar orbits around the supermassive black hole. Using axisymmetric Jeans models and the same spectroscopic data as our study, Feldmeier et al. (2014) found a lower value of M = (1.7) . Our best fit using triaxial Schwarzschild models is () . This measurement is consistent with the direct measurements of Boehle et al. (2016) within the 1 uncertainty limit. The result is also in agreement with the lower black hole mass of Feldmeier et al. (2014). We derived a 3 lower limit for the black hole of 0.7 , and an upper limit of 5.4 . We briefly discuss the model degeneracies, possible reasons for the different black hole mass measurements, and why our results are closer to the direct measurement than the black hole mass derived by Feldmeier et al. (2014).

Model degeneracies

Some model parameters seem to be correlated. This becomes clear when looking at Fig. 5. The best-fitting value of apparently increases with increasing dynamical mass-to-light ratio (second column of the first row). However, the value of has little effect on M, as can be seen in the second column of the second row in Fig. 5 . With a higher value of , the best-fitting M increases slightly. At larger , the -contours of M and broaden. This means that for a more oblate axisymmetric cluster with closer to one, M and are not as well constrained as with smaller values of . The black hole mass M=4.1  can be obtained with a higher value of =0.82 and =1.22, this model is within the 1 uncertainties.

The dynamical mass-to-light ratio is inversely correlated with the black hole mass (fourth column of the first row in Fig. 5). The higher , i.e. the more massive the cluster, the less massive is the black hole. This degeneracy is often obtained in dynamical models. Valluri et al. (2004) found that the degeneracy of M depends on how well the black hole’s sphere of influence is resolved. The measurement of is better constrained when the data extend to larger radii, provided that is constant over the entire field. We have several kinematic data bins within the radius of influence of the supermassive black hole, and our data extend to one effective radius. This may not be sufficient to put strong constraints on . To get agreement with the measurement of (4.02 0.20)  (Boehle et al., 2016), we would require a lower value of 0.75, which is within the 1 uncertainties.

Influence of the surface brightness profile

The shape of the surface brightness profile is important to estimate the mass of the supermassive black hole. The surface brightness profile has to represent the density of the kinematic tracer. We excluded young stars and ionised gas from the surface brightness profile, as these components contribute little mass compared to the cool, old stars we used as kinematic tracers. Excluding these components results in a lower surface brightness and stellar mass in the centre compared to Feldmeier et al. (2014). The stellar mass we obtain at = 48 arcsec (1.8 pc) is 1.3  less. Our black hole mass is therefore higher, and closer to the direct measurement of M 4 . We ran the same axisymmetric models (Cappellari, 2008), using the same kinematic data as Feldmeier et al. (2014), but our surface brightness distribution from Table 1. The best fit is obtained with M = () , = 0.89 , and a constant anisotropy of = 0.3. This result is in agreement with the triaxial Schwarzschild models, and confirms that the surface brightness profile has a strong influence on the results of the black hole mass and dynamical mass-to-light ratio.

The deprojection of the surface brightness profile is non-unique and only one possible solution. We applied the MGE method, which produces smooth intrinsic densities (Cappellari, 2002). We assumed that the deprojected density is smooth, as this is what we observe in other galaxies, where the clusters are observed from various viewing angles.

van den Bosch (1997) studied how much central density can be hidden in an oblate axisymmetric galaxy without effects on the projected surface brightness. They found that the percentage of hidden density in the centre of a Staeckel potential is 0 per cent for inclination =90, and 10 per cent for 80 of the total galaxy mass. The percentage of hidden density can increase in steeper cusps. However, the effect on the dynamics is still negligible. For our models, a hidden density is equivalent to a spatially varying mass-to-light ratio.

Spatially varying mass-to-light ratio

We assumed a constant dynamical mass-to-light ratio for the Schwarzschild models. We obtained = 0.90 (1 uncertainty). The dynamical mass-to-light ratio combines the stellar mass-to-light ratio with other components, it is sensitive to the presence of gas or dark matter.

Our best-fitting value of = 0.9 is consistent with stellar-population studies. Norris et al. (2014) found =0.8–1.2 at 4.6  for stellar populations with ages >7.0 Gyr and metallicities [] ranging from -1.0 dex to +0.3 dex. For populations with lower metallicity []=–2.18 dex and 13 Gyr age, can be even 1.5. However, for younger stellar populations 5 Gyr, Norris et al. (2014) obtained lower values 0.6. Our measurement of is averaged over the entire field of the kinematic data. We cannot exclude that the stellar age or metallicity changes over the range of the kinematic data. Stellar population studies of the red giant population were so far confined to the central 1–2 pc (e.g. Pfuhl et al., 2011; Do et al., 2015; Feldmeier-Krause et al., 2017). The stars in this region are mostly older than 5 Gyr and rather metal-rich. Our knowledge of the stellar population at the outer region of our field is based on only a few bright stars (e.g. Blum et al., 2003; Feldmeier et al., 2014). But these stars are brighter and probably younger than our kinematic tracer stars. However, the mass-to-light ratio for old stars in the mid-infrared varies modestly with age and metallicity in comparison to the optical mass-to-light ratio (Meidt et al., 2014). Therefore we do not expect a change of by more than 0.4 within the cluster. Should vary with radius, our mass profile (Fig. 7) would have a different shape. For example, if was lower in the centre than outside, this would increase M, and there would be less mass in the stellar component.

However, the stellar mass-to-light ratio may also increase towards the central =0.5 pc, as massive stellar remnants may migrate to the centre. The mass and distribution of dark stellar remnants, i.e. stellar-mass black holes and neutron stars, in the central parsec of the nuclear star cluster is uncertain. For a top-heavy initial mass function, there could be >1  in dark remnants (Morris, 1993), though Löckmann et al. (2010) found a lower mass of about 1  10 for a canonical initial mass function.

In our models we neglected the mass of molecular gas in the circum-nuclear disc. The molecular gas may contribute 10. The gas disc extends from 17 pc along the Galactic plane, but only to 3 pc along the minor axis (Ferrière, 2012). Thus, the molecular gas is located in the central part of our spectroscopic field, but absent in the North. If the gas contributes significantly to the cluster mass, our assumption of spatially constant would be violated, and would be higher than for a stellar component alone. When we assume the maximum gas mass of , the value of a constant decreases to about 0.85, which is within our 1 uncertainty limit.

The spatial distribution of dark matter in the Galactic centre is uncertain. A classical cuspy Navarro et al. (1996) dark matter profile results in a dark matter fraction of about 6.6 per cent in the central 100 pc (Linden, 2014). However, black hole accretion, dark matter annihilation, and scattering alter the shape of the dark matter distribution in the Galactic centre. Vasiliev & Zelnikov (2008) found that these effects produce a shallower dark matter profile in the central 2 pc than further out. The dark matter mass inferred from the classical cusp is reduced by up to 50 per cent in the central 2 pc. The contribution of dark matter to the nuclear star cluster mass should therefore be negligible. Although the dark matter distribution may be different from the luminous baryonic matter, and the dynamical mass-to-light ratio for that reason not spatially constant, the effect on the cluster mass distribution should be only minor.

4.2 Triaxial cluster shape

Our best-fitting model has axial ratios of  = = 0.28 , = = 0.64 , and = = 0.99 . These axial ratios correspond to viewing angles = 80°, = 79°, and = 91°. The angle denotes the polar viewing angle, the azimuthal viewing angle, and is the misalignment angle between photometric major axis and the projected intrinsic long axis (van den Bosch et al., 2008; van den Bosch & van de Ven, 2009). For the best-fitting model the angle between the cluster’s major axis and the line-of-sight is about 79°. The cluster’s shape is illustrated in Fig. 8, the triaxiality paramter varies between 0.45 and 0.7, with an average at 0.6. An oblate axisymmetric system has = 0, a prolate axisymmetric system has = 1.

Also the Milky Way’s bulge is triaxial, the axial ratios are = 0.26 and = 0.63 (Wegg & Gerhard, 2013). The shape was derived from the density of red clump stars in the central 800 pc of the bulge. The Milky Way bulge is much larger than the nuclear star cluster, and extends out to about 2.5 kpc. At the outer edge of our data, at =150 arcsec, we obtain for the intrinsic shape parameters =0.60 and =0.75. Both and are higher than found in the bulge. However, they are also decreasing, though the decrease of is not significant. The bulge has a peanut or X-shape (Nataf et al., 2010; McWilliam & Zoccali, 2010). The angle between the bulge major axis and the line-of-sight to the Galactic centre is about 27° (Rattenbury et al., 2007; Wegg & Gerhard, 2013), while we obtained 79° for the nuclear star cluster. There are also indications for another bar within the Galactic bulge from star count data in the inner 140 pc (Alard, 2001) or even 560 pc (Nishiyama et al., 2005). Rodriguez-Fernandez & Combes (2008) derived an angle 60-75° for a thick triaxial nuclear bar with axial ratios 0.55 and 0.75.

One possible scenario for nuclear star cluster formation is that massive star clusters (10–10) formed in the galactic disc, migrated to the galaxy’s centre and merged (Neumayer et al., 2011; Guillard et al., 2016). Simulations of multiple star cluster mergers and of star cluster accretion on a nuclear stellar component can produce triaxial nuclear star clusters (Bekki et al., 2004; Hartmann et al., 2011; Perets & Mastrobuono-Battisti, 2014). However, so far no observational study was able to constrain the triaxial shape of nuclear star clusters in general. Hartmann et al. (2011) constrained the shape of two nuclear star clusters and found agreement with an axisymmetric shape. For the Milky Way nuclear star cluster, may be as high as 0.94 within the 3 uncertainties. Thus, a nearly axisymmetric shape is consistent with our data.

4.3 Caveats and Considerations

Regime of semi-resolved populations

We used integrated light spectroscopy to measure the stellar kinematics. This is the common approach for extragalactic systems, which have a distance of several Mpc. The measured kinematics are weighted by the respective luminosities of different stars. As the Milky Way nuclear star cluster is only 8 kpc distant, we are in the regime of semi-resolved populations. The brightest stars can be resolved individually, and these stars contribute a large fraction of the flux. In consequence, individual spatial bins can be dominated by a single star. Instead of measuring the spectrum of an ensemble of stars, one measures a spectrum in which a large percentage of the flux is contributed by one single star. This causes shot noise, and high differences between neighbouring spatial bins. We accounted for this problem by excluding the brightest stars from the spectroscopic map. This method helps to significantly reduce the intrinsic scatter of the velocity dispersion (Lützgendorf et al., 2011; Bianchini et al., 2015). We further increased the kinematic uncertainties such that the data in two neighbouring bins have consistent values within their uncertainties. This helps to prevent that the models fit only stochastic shot noise. Due to the large kinematic uncertainties, the intrinsic shape parameters , , , and the dynamical mass-to-light ratio are not very well constrained, and have large error margins.

At a distance of only 8 kpc, also the relative distances of the stars become more important. A star located on the near side of the nuclear star cluster, at a distance d = 7.9 kpc, contributes 1.05 times more flux than a star with the same absolute magnitude at the far side of the cluster, at d = 8.1 kpc. In an extragalactic system, the distance of a star at the near side and the distance of a star at the far side with respect to the observer are approximately the same, as the system is farther away. For a galaxy at d = 5 Mpc, a relative difference of 200 pc changes the flux only by a factor 1.00008. Even foreground stars that belong to the outer parts of the stellar system contribute roughly the same flux as a star with the same magnitude that is located in the galactic nucleus.

Interstellar extinction

Another observational complication is interstellar extinction in the Galactic centre, which varies on arcsecond scales (Schödel et al., 2010). In particular, the field of view of the kinematic data contains the so-called 20-km s-cloud (M-0.13-0.08, e.g. García-Marín et al., 2011) in the Galactic southwest. It lies at a projected distance of about 70 arcsec (3 pc) from the centre, and probably about 5 pc in front of Sgr A* (Ferrière, 2012). This cloud blocks the light from stars of the nuclear star cluster. We cannot access the kinematics of stars behind this cloud. There is also interstellar dust within a projected distance of 20 arcsec (0.8 pc) from the centre, i.e. within the radius of influence of the black hole. This dust causes extinction within the nuclear star cluster by up to 0.8 mag (Chatzopoulos et al., 2015b). As a consequence, the two effects of dimming by distance and by extinction add up and stars that lie on the far side of the nuclear star cluster appear even more faint than the stars on the near side.


Both the semi-resolved stellar population and the inter-cluster extinction cause that our observations are biased to the near side of the nuclear star cluster. As a consequence, we measured a lower limit of the velocity dispersion. Feldmeier et al. (2014) found that the velocity dispersion in the projected radial range 6 arcsec < < 20 arcsec is smaller compared to the velocity dispersion computed from proper motion data of Schödel et al. (2009), which is based on resolved stars. For resolved stars, the velocity dispersion is not weighted by the flux of the stars. An underestimated velocity dispersion means that the black hole mass measurement is biased to lower values.

This observational bias also influences the measurements of , and . In particular, the cluster may appear compressed along the line-of-sight, and thus the value of = 0.64 may be too low. As a consequence, = 0.90 would be underestimated (see second column of the first row in Fig. 5).

Influence of figure rotation

The Galaxy rotates, and with it the nuclear star cluster. In a non-axisymmetric, rotating system, centrifugal and Coriolis forces play a role. However, figure rotation and the resulting forces were not included in our triaxial models. Figure rotation influences the stellar orbits. The prograde and retrograde tube orbits no longer fill the same volumes, while the box orbits acquire net mean angular momentum (e.g. Heisler et al., 1982; Schwarzschild, 1982; Sellwood & Wilkinson, 1993; Skokos et al., 2002). As a result, orbit-based, tumbling, triaxial models are computationally expensive. Other than an early attempt by Zhao (1996), no such models have been constructed that take into account kinematic data. It is difficult to predict how our results would change in a rotating model. The inferred orbital structure will be affected (depending on the tumbling speed of the nuclear star cluster), but our results on the mass distribution are likely to be fairly robust, as the assumption of a constant mass-to-light ratio is probably more important.

5 Summary and Outlook

We constructed for the first time triaxial orbit-based Schwarzschild models of the Milky Way nuclear star cluster. We used the spectroscopic integrated light maps of Feldmeier et al. (2014) to measure the cluster kinematics of the central 60 pc of the Milky Way. As photometry we used Spitzer 4.5 and NACO -band images, and measured a two-dimensional surface brightness distribution. We excluded young stars, avoided gas emission and dark clouds in the photometric data. Our triaxial models were based on the code by van den Bosch et al. (2008). Our best-fitting model contains a black hole of mass M = () , a dynamical mass-to-light ratio of = (0.90) /, and shape parameters = 0.28, = 0.64, and = 0.99. Our black hole mass measurement is in agreement with the direct measurement of (4.02 0.20)  by Boehle et al. (2016). We obtain a total cluster mass M = (2.10.7)  10within a spherical shell with radius = 2  = 8.4 pc. The best-fitting model is tangentially anisotropic in the central = 0.5-2 pc of the nuclear star cluster, but close to isotropic at larger radii. The model is able to recover the long-axis rotation in the central = 0.8 pc found by Feldmeier et al. (2014), and the misalignment of the kinematic rotation axis from the photometric minor axis.

There are several possible ways to extend the dynamical models in the future. One way is to include a component for the neutral gas disc inside the nuclear star cluster. If the gas mass is close to the upper limit of 10, the dynamical mass-to-light ratio would probably decrease slightly, and in return would slightly increase the black hole mass. Modelling a spatially varying mass-to-light ratio may provide a better representation of the cluster’s intrinsic properties. Further, proper motions can be included in combination with discrete line-of-sight velocities, as shown by van de Ven et al. (2006) and van den Bosch et al. (2006) for axisymmetric Schwarzschild models. Watkins et al. (2013) extended axisymmetric Jeans models and implemented discrete kinematic data without binning. Using discrete data means that the stars are not weighted by their luminosities. This prevents the previously discussed bias towards the near side of the cluster.


GvdV acknowledges partial support from Sonderforschungsbereich SFB 881 “The Milky Way System” (subproject A7 and A8) funded by the German Research Foundation. RS acknowledges funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/2007-2013) / ERC grant agreement n. 614922. We thank Remco van den Bosch for helpful discussions about the project. We finally thank the anonymous reviewer for useful comments and suggestions.


  1. thanks: Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile (289.B-5010(A)).
  2. pubyear: 2016
  3. pagerange: Triaxial orbit-based modelling of the Milky Way Nuclear Star Cluster thanks: Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile (289.B-5010(A)).Triaxial orbit-based modelling of the Milky Way Nuclear Star Cluster thanks: Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile (289.B-5010(A)).
  4. thanks: Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile (289.B-5010(A)).
  5. thanks: Based on observations collected at the European Organisation for Astronomical Research in the Southern Hemisphere, Chile (289.B-5010(A)).


  1. Alard C., 2001, \aap, 379, L44
  2. Alexander T., 2005, \physrep, 419, 65
  3. Bahcall J. N., Tremaine S., 1981, \apj, 244, 805
  4. Bekki K., Couch W. J., Drinkwater M. J., Shioya Y., 2004, \apjl, 610, L13
  5. Bianchini P., Norris M. A., van de Ven G., Schinnerer E., 2015, \mnras, 453, 365
  6. Blum R. D., Ramírez S. V., Sellgren K., Olsen K., 2003, \apj, 597, 323
  7. Boehle A., Ghez A. M., Schödel R., Meyer L., Yelda S., Albers S., Martinez G. D., Becklin E. E. et al, 2016, \apj, 830, 17
  8. Cappellari M., 2002, \mnras, 333, 400
  9. —, 2008, \mnras, 390, 71
  10. Cappellari M., Emsellem E., 2004, \pasp, 116, 138
  11. Chatzopoulos S., Fritz T. K., Gerhard O., Gillessen S., Wegg C., Genzel R., Pfuhl O., 2015a, \mnras, 447, 948
  12. Chatzopoulos S., Gerhard O., Fritz T. K., Wegg C., Gillessen S., Pfuhl O., Eisenhauer F., 2015b, \mnras, 453, 939
  13. Christopher M. H., Scoville N. Z., Stolovy S. R., Yun M. S., 2005, \apj, 622, 346
  14. Cretton N., de Zeeuw P. T., van der Marel R. P., Rix H.-W., 1999, \apjs, 124, 383
  15. Deguchi S., Imai H., Fujii T., Glass I. S., Ita Y., Izumiura H., Kameya O., Miyazaki A. et al, 2004, \pasj, 56, 261
  16. Do T., Kerzendorf W., Winsor N., Støstad M., Morris M. R., Lu J. R., Ghez A. M., 2015, \apj, 809, 143
  17. Do T., Martinez G. D., Yelda S., Ghez A., Bullock J., Kaplinghat M., Lu J. R., Peter A. H. G. et al, 2013, \apjl, 779, L6
  18. Drehmer D. A., Storchi-Bergmann T., Ferrari F., Cappellari M., Riffel R. A., 2015, \mnras, 450, 128
  19. Emsellem E., Monnet G., Bacon R., 1994, \aap, 285, 723
  20. Etxaluze M., Smith H. A., Tolls V., Stark A. A., González-Alfonso E., 2011, \aj, 142, 134
  21. Feldmeier A., Neumayer N., Seth A., Schödel R., Lützgendorf N., de Zeeuw P. T., Kissler-Patig M., Nishiyama S. et al, 2014, \aap, 570, A2
  22. Feldmeier-Krause A., Kerzendorf W., Neumayer N., Schödel R., Nogueras-Lara F., Do T., de Zeeuw P. T., Kuntschner H., 2017, \mnras, 464, 194
  23. Feldmeier-Krause A., Neumayer N., Schödel R., Seth A., Hilker M., de Zeeuw P. T., Kuntschner H., Walcher C. J. et al, 2015, \aap, 584, A2
  24. Ferrière K., 2012, \aap, 540, A50
  25. Figer D. F., McLean I. S., Morris M., 1999, \apj, 514, 202
  26. Franx M., 1988, \mnras, 231, 285
  27. Fritz T. K., Chatzopoulos S., Gerhard O., Gillessen S., Genzel R., Pfuhl O., Tacchella S., Eisenhauer F. et al, 2016, \apj, 821, 44
  28. Gao J., Li A., Jiang B. W., 2013, Earth, Planets, and Space, 65, 1127
  29. García-Marín M., Eckart A., Weiss A., Witzel G., Bremer M., Zamaninasab M., Morris M. R., Schödel R. et al, 2011, \apj, 738, 158
  30. Gebhardt K., Richstone D., Kormendy J., Lauer T. R., Ajhar E. A., Bender R., Dressler A., Faber S. M. et al, 2000, \aj, 119, 1157
  31. Genzel R., Eisenhauer F., Gillessen S., 2010, Reviews of Modern Physics, 82, 3121
  32. Genzel R., Thatte N., Krabbe A., Kroker H., Tacconi-Garman L. E., 1996, \apj, 472, 153
  33. Ghez A. M., Salim S., Weinberg N. N., Lu J. R., Do T., Dunn J. K., Matthews K., Morris M. R. et al, 2008, \apj, 689, 1044
  34. Gillessen S., Eisenhauer F., Trippe S., Alexander T., Genzel R., Martins F., Ott T., 2009, \apj, 692, 1075
  35. Guillard N., Emsellem E., Renaud F., 2016, \mnras, 461, 3620
  36. Haller J. W., Rieke M. J., Rieke G. H., Tamblyn P., Close L., Melia F., 1996, \apj, 456, 194
  37. Hartmann M., Debattista V. P., Seth A., Cappellari M., Quinn T. R., 2011, \mnras, 418, 2697
  38. Heisler J., Merritt D., Schwarzschild M., 1982, \apj, 258, 490
  39. Jeans J. H., 1922, \mnras, 82, 122
  40. Lawson C. L., Hanson R. J., 1974, Solving least squares problems
  41. Linden T., 2014, in IAU Symposium, Vol. 303, The Galactic Center: Feeding and Feedback in a Normal Galactic Nucleus, Sjouwerman L. O., Lang C. C., Ott J., eds., pp. 403–413
  42. Lindqvist M., Habing H. J., Winnberg A., 1992, \aap, 259, 118
  43. Löckmann U., Baumgardt H., Kroupa P., 2010, \mnras, 402, 519
  44. Lützgendorf N., Kissler-Patig M., Noyola E., Jalali B., de Zeeuw P. T., Gebhardt K., Baumgardt H., 2011, \aap, 533, A36
  45. Malkin Z., 2012, ArXiv e-prints
  46. McGinn M. T., Sellgren K., Becklin E. E., Hall D. N. B., 1989, \apj, 338, 824
  47. McWilliam A., Zoccali M., 2010, \apj, 724, 1491
  48. Meidt S. E., Schinnerer E., van de Ven G., Zaritsky D., Peletier R., Knapen J. H., Sheth K., Regan M. et al, 2014, \apj, 788, 144
  49. Merritt D., 2004, Coevolution of Black Holes and Galaxies, 263
  50. Morris M., 1993, \apj, 408, 496
  51. Nataf D. M., Udalski A., Gould A., Fouqué P., Stanek K. Z., 2010, \apjl, 721, L28
  52. Navarro J. F., Frenk C. S., White S. D. M., 1996, \apj, 462, 563
  53. Neumayer N., Walcher C. J., Andersen D., Sánchez S. F., Böker T., Rix H.-W., 2011, \mnras, 413, 1875
  54. Nishiyama S., Nagata T., Baba D., Haba Y., Kadowaki R., Kato D., Kurita M., Nagashima C. et al, 2005, \apjl, 621, L105
  55. Norris M. A., Meidt S., Van de Ven G., Schinnerer E., Groves B., Querejeta M., 2014, \apj, 797, 55
  56. Oh S., Kim S. S., Figer D. F., 2009, Journal of Korean Astronomical Society, 42, 17
  57. Oort J. H., 1977, \araa, 15, 295
  58. Paumard T., Genzel R., Martins F., Nayakshin S., Beloborodov A. M., Levin Y., Trippe S., Eisenhauer F. et al, 2006, \apj, 643, 1011
  59. Paumard T., Maillard J.-P., Morris M., 2004, \aap, 426, 81
  60. Perets H. B., Mastrobuono-Battisti A., 2014, \apjl, 784, L44
  61. Pfuhl O., Fritz T. K., Zilka M., Maness H., Eisenhauer F., Genzel R., Gillessen S., Ott T. et al, 2011, \apj, 741, 108
  62. Rattenbury N. J., Mao S., Sumi T., Smith M. C., 2007, \mnras, 378, 1064
  63. Reid M. J., Brunthaler A., 2004, \apj, 616, 872
  64. Requena-Torres M. A., Güsten R., Weiß A., Harris A. I., Martín-Pintado J., Stutzki J., Klein B., Heyminck S. et al, 2012, \aap, 542, L21
  65. Rieke G. H., Rieke M. J., 1988, \apjl, 330, L33
  66. Rix H.-W., de Zeeuw P. T., Cretton N., van der Marel R. P., Carollo C. M., 1997, \apj, 488, 702
  67. Rodriguez-Fernandez N. J., Combes F., 2008, \aap, 489, 115
  68. Rybicki G. B., 1987, in IAU Symposium, Vol. 127, Structure and Dynamics of Elliptical Galaxies, de Zeeuw P. T., ed., p. 397
  69. Saito R. K., Hempel M., Minniti D., Lucas P. W., Rejkuba M., Toledo I., Gonzalez O. A., Alonso-García J. et al, 2012, \aap, 537, A107
  70. Schödel R., Feldmeier A., Kunneriath D., Stolovy S., Neumayer N., Amaro-Seoane P., Nishiyama S., 2014, \aap, 566, A47
  71. Schödel R., Merritt D., Eckart A., 2009, \aap, 502, 91
  72. Schödel R., Najarro F., Muzic K., Eckart A., 2010, \aap, 511, A18
  73. Schwarzschild M., 1979, \apj, 232, 236
  74. Schwarzschild M., 1982, \apj, 263, 599
  75. Scoville N. Z., Stolovy S. R., Rieke M., Christopher M., Yusef-Zadeh F., 2003, \apj, 594, 294
  76. Sellgren K., McGinn M. T., Becklin E. E., Hall D. N., 1990, \apj, 359, 112
  77. Sellwood J. A., Wilkinson A., 1993, Reports on Progress in Physics, 56, 173
  78. Siopis C., Gebhardt K., Lauer T. R., Kormendy J., Pinkney J., Richstone D., Faber S. M., Tremaine S. et al, 2009, \apj, 693, 946
  79. Skokos C., Patsis P. A., Athanassoula E., 2002, \mnras, 333, 847
  80. Stolovy S., Ramirez S., Arendt R. G., Cotera A., Yusef-Zadeh F., Law C., Gezari D., Sellgren K. et al, 2006, Journal of Physics Conference Series, 54, 176
  81. Trippe S., Gillessen S., Gerhard O. E., Bartko H., Fritz T. K., Maness H. L., Eisenhauer F., Martins F. et al, 2008, \aap, 492, 419
  82. Tsatsi A., Mastrobuono-Battisti A., van de Ven G., Perets H. B., Bianchini P., Neumayer N., 2017, \mnras, 464, 3720
  83. Valluri M., Ferrarese L., Merritt D., Joseph C. L., 2005, \apj, 628, 137
  84. Valluri M., Merritt D., Emsellem E., 2004, \apj, 602, 66
  85. van de Ven G., de Zeeuw P. T., van den Bosch R. C. E., 2008, \mnras, 385, 614
  86. van de Ven G., van den Bosch R. C. E., Verolme E. K., de Zeeuw P. T., 2006, \aap, 445, 513
  87. van den Bosch F. C., 1997, \mnras, 287, 543
  88. van den Bosch R., de Zeeuw T., Gebhardt K., Noyola E., van de Ven G., 2006, \apj, 641, 852
  89. van den Bosch R. C. E., de Zeeuw P. T., 2010, \mnras, 401, 1770
  90. van den Bosch R. C. E., Greene J. E., Braatz J. A., Constantin A., Kuo C.-Y., 2016, \apj, 819, 11
  91. van den Bosch R. C. E., van de Ven G., 2009, \mnras, 398, 1117
  92. van den Bosch R. C. E., van de Ven G., Verolme E. K., Cappellari M., de Zeeuw P. T., 2008, \mnras, 385, 647
  93. van der Marel R. P., Cretton N., de Zeeuw P. T., Rix H.-W., 1998, \apj, 493, 613
  94. Vasiliev E., Zelnikov M., 2008, \prd, 78, 083506
  95. Wallace L., Hinkle K., 1996, \apjs, 107, 312
  96. Watkins L. L., van de Ven G., den Brok M., van den Bosch R. C. E., 2013, \mnras, 436, 2598
  97. Wegg C., Gerhard O., 2013, \mnras, 435, 1874
  98. Zhao H., 1996, \mnras, 283, 149
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minumum 40 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description