Using the motion of S2 to constrain vector clouds around SgrA ∗

The dark compact object at the centre of the Milky Way is well established to be a supermassive black hole with mass 𝑀 • ∼ 4 . 3 · 10 6 𝑀 ⊙ , but the nature of its environment is still under debate. In this work, we used astrometric and spectroscopic measurements of the motion of the star S2, one of the closest stars to the massive black hole, to determine an upper limit on an extended mass composed of a massive vector field around Sagittarius A*. For a vector with effective mass 10 − 19 eV ≲ m s ≲ 10 − 18 eV, our Markov Chain Monte Carlo analysis shows no evidence for such a cloud, placing an upper bound 𝑀 cloud ≲ 0 . 1% 𝑀 • at 3 𝜎 confidence level. We show that dynamical friction exerted by the medium on S2 motion plays no role in the analysis performed in this and previous works, and can be neglected thus.


INTRODUCTION
Since the star S2 has been discovered orbiting the Galactic Center (GC) (Schödel et al. 2002;Ghez et al. 2003;Gillessen et al. 2009Gillessen et al. , 2017)), its orbital motion has been largely and extensively used to constrain the properties of the supermassive black hole (SMBH) Sagittarius A * (SgrA * ) and the environment around it.S2 is part of the so-called S-cluster, which currently counts up to tens of detected stars (Sabha et al. 2012;Habibi et al. 2017;GRAVITY Collaboration et al. 2022).
The astrometric and spectroscopic data collected by two independent groups showed that the dynamics of S-stars is entirely dominated by the presence of a compact source with  • ∼ 4.3 • 10 6  ⊙ at a distance of  0 ∼ 8.3 kpc.There is overwhelming evidence that the compact source is a SMBH (Schödel et al. 2002;Ghez et al. 2008;Genzel et al. 2010;Genzel 2021;GRAVITY Collaboration 2019b, 2022).Very strong arguments that the central dark mass is indeed an SMBH come from the measurement of the Schwarzschild precession in the orbit of S2 (GRAVITY Collaboration 2020), from the observations of near-IR flares in correspondence with the innermost circular orbit of the SMBH (GRAVITY Collaboration 2018; Abuter et al. 2023) and by the image released by the Event Horizon Telescope collaboration, which is compatible with the expected image of a Kerr BH (Akiyama et al. 2022).
The physics of horizons is so puzzling that any further evidence for their existence is welcome and provides important information on the scales at which new physics sets in.Currently, it is challenging to use orbits of S-stars around the GC to test the nature of the compact source itself and to distinguish it from other possible models, such as boson stars, dark matter (DM) cores or wormholes, which have similar features to BHs (Amaro-Seoane et al. 2010;Grould et al. 2017a;Boshkayev & Malafarina 2019;Della Monica & de Martino 2022;De Laurentis et al. 2022).Note, however, that the optical appearance of hot spots (or stars) close to the accretion zone of SgrA * , may differ significantly should an horizon be absent (Rosa et al. 2022).
Equally important is the nature of the environment around SMBHs, in particular around SgrA * .Dark matter (DM) is expected to cluster at the center of galaxies leading to "overdensities" (Gondolo & Silk 1999;Sadeghian et al. 2013), which might leave an imprint in the motion of stars.S-stars are currently the main observational tool we have to look into this inner region of our Galaxy and thus they must be exploited to gain as much information as possible from their motion.For this and other reasons, the possibility of an extended mass distribution around SgrA * have been studied (Lacroix 2018;Bar et al. 2019;Heißel et al. 2022;GRAVITY Collaboration 2022;Foschi et al. 2023).Specifically, GRAVITY Collaboration (2022) derived an upper limit of  ∼ 4000  ⊙ ∼ 0.1% • for a density distribution described by a Plummer profile with length-scale  0 = 0.3 ′′ .
A special, and interesting, model for dark matter concerns ultralight bosons.These arise in a variety of scenarios, for instance the "string axiverse" (Arvanitaki et al. 2010;Arvanitaki & Dubovsky 2011;Marsh 2016) or as a hidden U(1) gauge boson, a generic feature of extensions of the Standard Model (Goodsell et al. 2009;Jaeckel & Ringwald 2010).In fact, such fields can exist and grow even if they are only a minute component of DM, as they are amplified via a mechanism known as BH superradiance (Brito et al. 2015b).In this process, the light boson extracts rotational energy away from the spinning BH, depositing it in a "bosonic cloud", which can acquire a sizeable fraction of the BH mass.For a fundamental boson of mass   the key parameter controlling the superradiant growth and energy extraction is the mass coupling  =  •   .
In a recent work (Foschi et al. 2023), we investigated the possibility that a massive scalar field clusters around SgA * in the form of a cloud (GRAVITY Collaboration 2019a).We showed that for the range of (dimensionless) mass couplings, 0.01 ≲  ≲ 0.045 (which corresponds to a mass of the scalar field of 6 • 10 −19 eV ≲ m s ≲ 3 • 10 −18 eV) we are able to constrain the mass of the cloud to be  cloud ≲ 0.1%  • , recovering the upper bound found in GRAVITY Collaboration (2022).
Here, we focus on a similar system: a massive vector cloud.As scalar fields, massive vector fields can form bound states around Kerr BHs, giving rise to stationary clouds.At the linear level and using the small coupling approximation, it has been shown that the superradiant instability is triggered on a timescale   ∝  −7 for vector clouds when compared to the scalar case of   ∝  −9 (Pani et al. 2012;Brito et al. 2015b;Cardoso et al. 2018;Endlich & Penco 2017).Hence vector clouds grow much faster than their scalar counterparts and the field's mass   needed to make them grow in a timescale smaller than the cosmic age is much smaller, making them more likely to be observed.
In this work we will use the astrometric and spectroscopic data of star S2 collected at the Very Large Telescope (VLT) to constrain the fractional mass of a possible vector cloud around SgrA * .

SETUP
In this work, we consider a massive vector field   described by the Lagrangian and   satisfies the Proca equation of motion     =  2   .If the Compton wavelength of the vector field is much larger than the Schwarzschild radius   =  • , the bound states of the field oscillate with frequency   ≃  and can be written as (Baryakhtar et al. 2017) In the limit  ≫   , the Proca equation becomes a Schröedinger-like equation, and the Ψ 0 component can be expressed in terms of Ψ  .
Since the radial part of the potential is spherically symmetric, Ψ  can be decomposed as where the  ℓ,   (, ) are the so-called pure-orbital vector spherical harmonics (Thorne 1980;Santos et al. 2020).
The fundamental mode of the field, which is also the mode that grows fastest due to superradiant mechanisms (Baryakhtar et al. 2017) is given by ℓ = 0,  =  = 1 and  = 0.At leading order in  we can neglect  0 and consider only the spatial components of the field, which can be written as (Chen et al. 2023) From this profile, we can compute the energy-momentum tensor (Herdeiro et al. 2016) and take the Newtonian limit, i.e. neglecting all the spatial derivatives and assuming a real field, obtaining: which coincides with the expression in Chen et al. (2023).
As done in Foschi et al. (2023), we can integrate the energy density in Eq. ( 5) to relate the amplitude of the field Ψ 0 with the mass of the vector cloud: From the energy density in Eq. ( 5) we can get the potential generated by the cloud solving Poisson's equation: ∇ 2   = 4 and using the spherical harmonic decomposition of Poisson & Will (2012) to get: where we have defined Λ =  cloud / • .

Effects of the cloud on S2 orbit with osculating elements
We start our analysis of the effects of vector cloud on S2 motion using the method of osculating elements that can be found in Poisson & Will (2012).The basic idea is to treat the effect of the vector cloud as a perturbation of the Newtonian acceleration, assuming that the Keplerian description of the orbit is still approximately true.In this way, we are able to express the equations of motion in terms of the Keplerian elements (, , , , Ω, M 0 ) (eccentricity, semi-major axis, inclination, argument of the periastron, longitude of the ascending node and mean anomaly at epoch, respectively), which would be constant in a pure Newtonian setup, and see how the perturbing force modifies them.In order to do so, we introduce a vectorial basis adapted to the orbital motion of the binary system BH-S2: where n = r/, e  = h/ℎ with h := r × v and  is orthogonal to both n and e  .We also assume that the mass of the star is negligible compared to the BH mass  • .
The perturbing force can be decomposed as: The variation of the orbital elements in terms of the perturbing force components is given in Kopeikin et al. (2011); Poisson & Will (2012) and we report it for completeness in Appendix A.
Once the variation in time of the orbital elements is known, one can compute the secular change of the orbital element   over a complete orbit using: where and (11)

Effect of the vector cloud alone
Due to the spherical symmetry of the energy distribution in Eq. ( 5), the only non-zero component of f  is the radial one: while S  = W  = 0.

Inclusion of the 1PN correction
Since the Schwarzschild precession has been detected at 8 confidence level by the GRAVITY collaboration (GRAVITY Collaboration 2020, 2022), it is interesting to see how the previous results change if we include the first Post Newtonian (PN) correction to the equations of motion.
This corresponds to having a total acceleration where with  =  r,  =  r,   θ,   sin  φ and  = ||.
The decomposition of the acceleration in Eq. ( 14) into the basis (n, , e  ) has been done in Poisson & Will (2012) and here we report the result: and W 1PN = 0.In order to express everything in terms of the orbital elements, we need to use the expressions for ,  and  reported in Sec.10.1.3 of Poisson & Will (2012).
In this second case we set Λ = 10 −3 , which corresponds to the current upper limit obtained by the GRAVITY collaboration for the fractional mass of an extended mass distribution around SgrA * (GRAV-ITY Collaboration 2022; Foschi et al. 2023).

Data
The set of available data  is the same as in Foschi et al. (2023).

Fitting approach
The next step is to obtain a best-fit value for the fractional mass Λ for different coupling  values.The procedure followed in this work is exactly the same as the one reported in Foschi et al. (2023).Specifically, we solve the equations of motion in Eq. ( 13) using the initial conditions reported in Appendix B. The solutions of this set of equations are given in the BH reference frame and must be projected into the observer reference frame using the three Euler angles Ω, , .
Following Grould et al. (2017b) we can define a new reference frame { ′ ,  ′ ,  obs } such that  ′ = DEC,  ′ = R.A. are the collected astrometric data,  obs points towards the BH and   obs corresponds to the radial velocity (see Appendix C for details about how to perform the rotation of the reference frame).
Moreover, it is true that S2 motion happens mostly in a Newtonian regime, i.e. with  ≪ 1, but near the periastron, it reaches a total space velocity  ∼ 10 −2 .In this region, relativistic effects become important and can not be neglected.For this reason, we correct the radial velocity coming from Eq. ( 13), including both the relativistic Doppler shift and the gravitational redshift (Abuter et al. 2018).
Finally, we also consider the so-called Rømer delay, which is the difference between the observational dates and the actual emission dates of the signal due to the finite speed of light.Details about how to include Rømer delay and relativistic effects are reported in Appendix D.
For any given value of , we fit for the following set of parameters, The additional parameters { 0 ,  0 ,   0 ,   0 ,   0 } characterise the NACO/SINFONI data reference frame with respect to Sgr A* (Plewa et al. 2015).We refer the reader to Appendix E for more details about the MCMC implementation.

Variation of the orbital elements
In Figure 1 we show the variation of the orbital elements Δ  /Λ due to the presence of the vector cloud for different values of the coupling , as described in Sec.2.1.The secular change is negligible for both the eccentricity  and the semi-major axis .
The change in the mean anomaly at epoch M 0 is instead proportional to , increasing monotonically.M 0 is directly related to the time of pericenter passage   : a larger mean anomaly at the epoch corresponds to a later pericenter passage.
The only meaningful change in the orbital elements is found in Δ, which quantifies the precession effect on the orbit, with  the argument of pericenter.First of all, we observe that Δ < 0 always.This is a consequence of the fact that the presence of an extended mass within the orbit of S2 would produce a retrograde precession of the orbit (Heißel et al. 2022).
Unsurprisingly, its maximum variation is found in the range 0.003 ≲  ≲ 0.03 .
Indeed, as in the case of scalar clouds (Foschi et al. 2023), this behaviour is expected if we compute the effective peak position of the energy distribution in Eq. ( 5), which, for the values of  reported in Eq. ( 18), corresponds to 5 • 10 2  • ≲  peak ≲ 5 • 10 4  • , i.e. it roughly matches the orbital range of S2 (3 . This result shows that the maximum variation in  is found when the star crosses regions of higher (vector) density, while its orbit remains basically unaffected if the cloud is located away from its apoastron or too close to the central BH mass.
In Figure 2 we show the variation of the orbital elements when the 1PN correction is included in the equations of motion, as described in Sec.2.1.2.Opposite to the previous case, here, the variation of the argument of the pericenter Δ can be either positive or negative, according to the value of .Indeed now the retrograde precession induced by the vector cloud is compensated by the (prograde) Schwarzschild precession due to the 1PN correction in the equations of motion, and its maximum value corresponds to Δ ≃ −1.8 ′ , which is smaller than the previous case with Λ = 10 −3 (Δ ≃ −6 ′ ).

Limit on the fractional mass Λ
Before running the MCMC algorithm we determine the initial guesses for the parameters listed in Eq. ( 17).We performed a simple  2 minimization using the Python package lmfit.minimize(Newville et al. 2016) with Levenberg-Marquardt method.In Figure 3 we report the best-fit values of Λ with relative 1 uncertainties, and we compare the range of  with the effective peak position of the cloud in Eq. ( 19).The smallest uncertainties for Λ are found roughly in the range of Eq. ( 18), which is slightly different from the scalar cloud case (Foschi et al. 2023) and in agreement with the orbital variation reported in Figure 2.
After performing the MCMC analysis, we look for the maximum likelihood estimator (MLE) Λ, which in this case corresponds to the value that maximises the posterior density distribution reported in Figure 4, as a consequence of using flat priors and a Gaussian likelihood.
In Table 1 we report the values of Λ with relative 1 uncertainties When the posterior distribution is found to be non-normal and peaked at zero, we estimated the 1 (3) confidence interval looking for that value of Λ such that roughly the 68% (99%) of (Λ|) lies below that value.When  ≳ 0.3, the distribution of Λ start to be flat, with a sudden drop around Λ ≃ 10 −2 .One can show that for flat distributions in an interval [, ], the mean is given by ( − )/2 while the variance is ( − ) 2 /12 (Bailer-Jones 2017).We report those values in Table 1.However, what is important to notice in these cases is that for  ≳ 0.03 ( peak ≲ 550  • ), it is not possible to determine a unique value for Λ that best fits the data, confirming the expectation from the  2 minimisation.
When  is in the range of Eq. ( 18) the posterior distributions of Λ are Gaussian whose means and standard deviations are reported in Table 1.For all cases considered in this range, Λ ∼ 10 −3 with 1 uncertainties roughly of the same order of magnitude.This makes all the Λ values derived from the MCMC analysis compatible with zero within the 3 confidence level.In addition to this, the associated Bayes factors always have log  < 2. This result, according to the literature (Kass & Raftery 1995), shows no statistical evidence in favour of the BH plus vector cloud model with respect to the nonperturbative case where no cloud is present.Hence we derive an upper limit of Λ ≲ 10 −3 at 3 confidence level.This upper bound imposes a limit on the superradiant growth, that in general would lead to transfer up to ∼ O (10)% of the BH mass into the vector cloud (Brito et al. 2015a;East & Pretorius 2017;Herdeiro et al. 2022).Here we showed that for a field's effective mass of   ∼ 10 −19 −10 −18 eV, the mass of the cloud around SgrA * can not exceed the limit  cloud ≲ 0.1% • .For a BH spinning with / ∼ 0.5 (an indicative value), the growth timescale of the cloud can vary between 10 5 − 10 10 yrs, exact values depend on the effective mass   .This estimate is below the age of the Universe ( age ∼ 10 10 yrs), making the superradiant process and our constraints relevant.In Appendix  F we report the corner plots of two illustrative cases ( = 0.01,  = 0.001) to show the correlations between parameters.

Inclusion of environmental effects
All the above results are obtained neglecting the backreaction effects of the matter on the motion of S2.Indeed, the presence of a matter distribution induces a gravitational drag force on the body moving in it, with the consequence that part of the material is dragged along the motion producing dynamical friction force on the main body (Chandrasekhar 1983;Ostriker 1999).It has been shown that dynamical friction induced by ultralight bosons may play a significant role in the strong regime (Traykova et al. 2021;Vicente & Cardoso 2022).
Here we investigated whether dynamical friction affects S2 motion too.
In a Newtonian setup, including the dynamical friction force means adding the following two components to the equations of motion (Macedo et al. 2013): where  2 =  2 +  2  2 , since we have assumed that the motion of S2 happens on the equatorial plane ( = /2) of the central SMBH.
The term  DF has been derived in Ostriker (1999) for a perturber in linear motion and it reads: where  is the density of the matter distribution in Eq. ( 5),   is the mass of the star S2 that we take to be   = 14  ⊙ and   is the speed of sound in the medium which constitutes the environment.Kim & Kim (2007) showed that Eq. ( 21) correctly reproduces the results obtained for circular orbits if one substitutes  → 2 ().
Despite the orbit of S2 is far from being circular, we are going to use Eqs.(20) in a first approximation.
We tested four different values of the speed of sound   for both the supersonic (  = 10 −6 ,   = 10 −3 ) and the subsonic (  = 0.1,   = 0.03) regimes, for different values of .We set Λ = 10 −3 , since this corresponds to the maximum allowed value of the fractional mass, but results scale linearly with it.
We found that results are independent on   and that the maximum difference in both the astrometry and the radial velocity with respect to the case where no dynamical friction is implemented is always negligible.
In Figure 5 we report the absolute difference in DEC, R.A. and radial velocity in the supersonic case with   = 10 −3 .Overall, the effect of dynamical friction is at most 10 −5 mas in the astrometry and ≈ 10 −3 km/s in the radial velocity, and in both cases is reached around the periastron passages.Overall, it remains well below the current (and future) instrument precision and can be neglected.
We performed the same analysis for the scalar cloud model implemented in Foschi et al. (2023) and the Plummer density profile tested in GRAVITY Collaboration (2022) too.In both cases, we found similar results to Figure 5 and hence we conclude that dynamical friction effects can be safely neglected.
Along the same line, one can try to compute the effect that regular gas around SgrA * has on S2 orbit.In Gillessen et al. (2018), the A. and radial velocity between the case where dynamical friction is implemented in the supersonic case with   = 10 −3 and the case where no dynamical friction is present.We set Λ = 10 −3 , but results scale linearly with Λ.The difference is maximum around the periastron passages and minimum at the apoastron (black dotted line).Overall, they remain far below the current instrument threshold, whatever the value of .
authors detected a drag force acting on the gas cloud G2 orbiting around SgrA * and they derived an estimate for the number density of the ambient.Here we used their same formulation for the drag force, meaning where  = 1,  is the relative velocity between the medium and the star, that, following Gillessen et al. (2018), is assumed to be equal to the velocity of the star itself and   parametrizes the strength of the drag force and it is related to the normalized number density of the gas ambient.In Gillessen et al. (2018) they derived   ∼ 10 −3 , which is the value used in this work as well.In this case no vector cloud is present (Λ = 0) and only the force contribution due to the presence of gas is considered.The maximum difference induced by the drag force exerted by the gas ambient on the astrometry and the radial velocity of S2 is of order ∼ 10 −6 mas and ∼ 10 −3 km/s, respectively.Hence, also the contribution due to regular gas around SgrA * has a negligible effect on S2.We also note that the difference induced by the presence of gas is comparable with the effect produced by dynamical friction.Hence, even with the development of future instruments and the advent of GRAVITY+, it will still be hard to disentangle the two effects.

CONCLUSIONS
In this paper we investigated the possibility that a vector cloud of superradiant origin clusters around the SMBH SgrA * , extending the analysis on scalar clouds performed in Foschi et al. (2023).Specifically, we considered a massive vector field, which gives rise to a spherically symmetric cloud and in Sec.3.1 we investigated the imprints of such a cloud in S2's orbital elements.The MCMC analysis in Sec.3.2 confirmed the current upper bound for the fractional mass of Λ ≲ 0.1% • , recovering previous results on extended masses (GRAVITY Collaboration 2022;Foschi et al. 2023).Despite the range of field's masses that can be tested with S2 motion is roughly the same in both the scalar and vector cloud case (10 −18 eV ≲ m s ≲ 10 −19 eV), in the latter those values can effectively engage a superradiant instability in a timescale shorter than the cosmic age.This strongly constrains the mass of a possible superradiant cloud at the GC, improving the theoretical bound that can lead to have masses up to two order of magnitude larger (Brito et al. 2015a;East & Pretorius 2017;Herdeiro et al. 2022).
Moreover, the effect of the environment on S2 orbit was also investigated for the first time.We considered both the dynamical friction exerted by the medium on the star, and the effect of ambient gas around SgrA * .In both cases, the effect on the astrometry and the radial velocity are negligible.This analysis was also extended to the scalar cloud case considered in Foschi et al. (2023) and to the Plummer profile of GRAVITY Collaboration (2022), showing that even in those cases both effect can be neglected.However, since the difference in the astrometry and the radial velocity induced by those effects is of the same order of magnitude, it will be difficult to separate them even with the advent of future instrumentation.minimum found by minimize we skipped the burning-in phase and we used the last 80% of the chains to compute the mean and standard deviation of the posterior distributions.The convergence of the MCMC analysis is assured by means of the auto-correlation time   , i.e. we ran  iterations such that  ≫ 50   .

APPENDIX F: CORNER PLOTS
Here we report the corner plots for two representative values of  ( = 0.01 and  = 0.001), to show the behaviour of the parameters when the cloud is located in and outside S2's orbital range.The strong correlation between Λ and the periastron passage   when  = 0.01 can be understood following the argument of Heißel et al. (2022): the presence of an extended mass will induce a retrograde precession in the orbit that will result in a positive shift of the periastron passage time, needed to compensate the (negative) shift in the initial true anomaly.Indeed, when considering the Schwarzschild precession, which instead induces a prograde precession (hence a positive initial shift in the true anomaly),   will undergo a negative shift, as can be seen from the strong anti-correlation between  SP and   reported in GRAVITY Collaboration (2020).
This paper has been typeset from a T E X/L A T E X file prepared by the author.

Figure 1 .
Figure 1.Variation of the orbital elements Δ  /Λ over an entire orbit for different values of the coupling constant  when only the vector cloud is present.The maximum variation in Δ/Λ is roughly found in the range 0.003 ≲  ≲ 0.03.

Figure 2 .
Figure 2. Variation of the orbital elements Δ  over an entire orbit for different values of the coupling constant  when one includes the Schwarzschild precession in the equation for the osculating elements.Here Λ = 10 −3 .The maximum variation is still found in 0.003 ≲  ≲ 0.03.

Figure 3 .Figure 4 .
Figure 3. Best-fit values for Λ and relative 1 uncertainties as function of the coupling  obtained minimizing the  2 .The grey dashed line represents the effective peak position of the vector cloud given by Eq. (19), while the orange band gives the orbital range of S2.

Figure 5 .
Figure5.Absolute difference in DEC, R.A. and radial velocity between the case where dynamical friction is implemented in the supersonic case with   = 10 −3 and the case where no dynamical friction is present.We set Λ = 10 −3 , but results scale linearly with Λ.The difference is maximum around the periastron passages and minimum at the apoastron (black dotted line).Overall, they remain far below the current instrument threshold, whatever the value of .

Table 1 .
Maximum Likelihood Estimator Λ with associated 1 error and Bayes factors log 10  for different values of .The measurements for each  are not independent (the same orbit was used to derive them) and therefore cannot be combined to derive a more stringent upper limit.For non-normal distributions we report Λ 1 and Λ 2 defined such that  (Λ  < Λ 1 |) ≈ 68% and  (Λ  < Λ 2 |) ≈ 99% of  (Λ  |).
Heavens et al. (2017)ue of the Bayes factor log .The latter is obtained computing the marginal likelihoods by making use of the Python package MCEvidence developed inHeavens et al. (2017)and it is defined as  = ( |  )/(| 0 ), where   represents the BH plus vector cloud model while  0 corresponds to the non perturbative one.

Table E2 .
Plewa et al. (2015)d in the MCMC analysis.Initial guesses Θ 0  coincide with the best-fit parameters found by minimize.andrepresent the mean and the standard deviation of the distributions, respectively, and they come fromPlewa et al. (2015).