Recent Developments Concerning the Gravitational Potential of the Milky Way. III. A Closer Look at the RAR Model

I am primarily an extragalactic astronomer – someone who studies galaxies outside our own. Our home Galaxy is a subject in its own right. Naturally, I became curious how the Milky Way appeared in the light of the systematic behaviors we have learned from external galaxies. I first wrote a paper about it in 2008; in the process I realized that I could use the RAR to infer the distribution of stellar mass from the terminal velocities observed in interstellar gas. That’s not necessary in external galaxies, where we can measure the light distribution, but we don’t get a view of the whole Galaxy from our location within it. Still, it wasn’t my field, so it wasn’t until 2015/16 that I did the exercise in detail. Shortly after that, the folks who study the supermassive black hole at the center of the Galaxy provided a very precise constraint on the distance there. That was the one big systematic uncertainty in my own work up to that point, but I had guessed well enough, so it didn’t make a big change. Still, I updated the model to the new distance in 2018, and provided its details on my model page so anyone could use it. Then Gaia data started to pour in, which was overwhelming, but I found I really didn’t need to do any updating: the second data release indicated a declining rotation curve at exactly the rate the model predicted: -1.7 km/s/kpc. So far so good.

I call it the RAR model because it only involves the radial force. All I did was assume that the Milky Way was a typical spiral galaxy that followed the RAR, and ask what the mass distribution of the stars needed to be to match the observed terminal velocities. This is a purely empirical exercise that should work regardless of the underlying cause of the RAR, be it MOND or something else. Of course, MOND is the only theory that explicitly predicted the RAR ahead of time, but we’ve gone to great lengths to establish that the RAR is present empirically whether we know about MOND or not. If we accept that the cause of the RAR is MOND, which is the natural interpretation, then MOND over-predicts the vertical motions by a bit. That may be an important clue, either into how MOND works (it doesn’t necessarily follow the most naive assumption) or how something else might cause the observed MONDian phenomenology, or it could just be another systematic uncertainty of the sort that always plagues astronomy. Here I will focus on the RAR model, highlighting specific radial ranges where the details of the RAR model provide insight that can’t be obtained in other ways.

The RAR Milky Way model was fit to the terminal velocity data (in grey) over the radial range 3 < R < 8 kpc. Everything outside of that range is a prediction. It is not a prediction limited to that skinny blue line, as I have to extrapolate the mass distribution of the Milky Way to arbitrarily large radii. If there is a gradient in the mass-to-light ratio, or even if I guess a little wrong in the extrapolation, it’ll go off at some point. It shouldn’t be far off, as V(R) is mostly fixed by the enclosed mass. Mostly. If there is something else out there, it’ll be higher (like the cyan line including an estimate of the coronal gas in the plot that goes out to 130 kpc). If there is a bit less than the extrapolation, it’ll be lower.

From 8 to 19 kpc, the Gaia data as realized by Zhao et al. fall bang on the model. They evince exactly the slowly declining rotation curve that was predicted. That’s pretty good for an extrapolation from R < 8 kpc. I’m not aware of any other model that did this well in advance of the observation. Indeed, I can’t think of a way to even make a prediction with a dark matter model. I’ve tried this – a lot – and it is as easy to come up with a model whose rotation curve is rising as one that is falling. There’s nothing in the dark matter paradigm that is predictive at this level of detail.

Beyond R > 19 kpc, the match of the model and Zhou et al. realization of the data is not perfect. It is still pretty damn good by astronomical standards, and better than the Keplerian dotted line. Cosmologists would be wetting themselves with excitement if they could come this close to predicting anything. Heck, they’re known to do that even when they’re obviously wrong*.

If the difference between the outermost data and the blue line is correct, then all it means is that we have to tweak the model to have a bit less mass than assumed in the extrapolation. I call it a tweak because it would be exactly that: a small change to an assumption I was obliged to make in order to do the calculation. I could have assumed something else, and almost did: there is discussion in the literature that the disk of the Milky Way is truncated at 20 kpc. I considered using a mass model with such a feature, but one can’t make it a sharp edge as that introduces numerical artifacts when solving the Poisson equation numerically, as this procedure depends on derivatives that blow up when they encounter sharp features. Presumably the physical truncation isn’t unphysically sharp anyway, rather being a transition to a steeper exponential decline as we sometimes see in other galaxies. However, despite indications of such an effect, there wasn’t enough data to constrain it in a way useful for my model. So rather than introduce a bunch of extra, unconstrained freedom into the model, I made a straight extrapolation from what I had all the way to infinity in the full knowledge that this had to be wrong at some level. Perhaps we’ve found that level.

That said, I’m happy with the agreement of the data with the model as is. The data become very sparse where there is even a hint of disagreement. Where there are thousands of stars per bin in the well-fit portion of the rotation curve, there are only tens per bin outside 20 kpc. When the numbers get that small, one has to start to worry that there are not enough independent samples of phase space. A sizeable fraction of those tens of stars could be part of the same stellar stream, which would bias the results to that particular unrepresentative orbit. I don’t know if that’s the case, which is the point: it is just one of the many potential systematic uncertainties that are not represented in the formal error bars. Missing those last five points by two sigma is as likely to be an indication that the error bars have been underestimated as it is to be an indication that the model is inadequate. Trying to account for this sort of thing is why the error bars of Jiao et al. are so much bigger than the formal uncertainties in the three realization papers.

That’s the outer regions. The place where the RAR model disagrees the most with the Gaia data is from 5 < R < 8 kpc, which is in the range where it was fit! So what’s going on there?

Again, the data disagree with the data. The stellar data from Gaia disagree with the terminal velocity data from interstellar gas at high significance. The RAR model was fit to the latter, so it must per force disagree with the former. It is tempting to dismiss one or the other as wrong, but do they really disagree?

In order to build the model depicted above, I chose to split the difference between the first and fourth quadrant terminal velocity data. I fit them separately in McGaugh (2016) where I made the additional point that the apparent difference between the two quadrants is what we expect from an m=2 mode – i.e., a galaxy with spiral arms. That means these velocities are not exactly circular as commonly assumed, and as I must per force assume to build the model. So I split the difference above in the full knowledge that this is not the exact circular velocity curve of the Galaxy, it’s just the best I can do at present. This is another example of the systematic uncertainties we encounter: the difference between the first and fourth quadrant is real and is telling us that the galaxy is not azimuthally symmetric – as anyone can tell by looking at any spiral galaxy, but is a detail we’d like to ignore so we can talk about disk+dark matter halo models in the convenient limit of axisymmetry.

Though not perfect – no model is – the RAR model Milky Way is a lot better than models that ignore spiral structure entirely, which is basically all of them. The standard procedure assumes an exponential disk and some form of dark matter halo. Allowance is usually made for a central bulge component, but it is relatively rare to bother to include the interstellar gas, much less consider deviations from a pure exponential disk. Having adopted the approximation of an exponential disk, one inevitably get a smooth rotation curve like the dashed line below:

The common assumption of exponential disk precludes the possibility of fitting the bumps and wiggles observed in the terminal velocities. These occur because of deviations from a pure exponential profile caused by features like spiral arms. By making this assumption, the variations in mass due to spiral arms is artificially smoothed over. They are not there by assumption, and there is no way to recover them in a dark matter fit that doesn’t know about the RAR.

Depending on what one is trying to accomplish, an exponential model may suffice. The Bovy & Rix model shown above is perfectly reasonable for what they were trying to do, which involved the vertical motions of stars, not the bumps and wiggles in the rotation curve. I would say that the result they obtain is in reasonable agreement with the rotation curve, given what they were doing and in full knowledge that we can’t expect to hit every error bar of every datum of every sort. But for the benefit of the chi-square enthusiasts who are concerned about missing a few data points at large radii, the reduced chi-squared of the Bovy & Rix model is 14.35 while that of the RAR model is 0.6. A good fit is around 1, so the RAR model is a good fit while the smooth exponential is terrible – as one can see by eye in the residual inset: the smooth exponential model gets the overall amplitude about right, but hits none of the data. That’s the starting point for every dark matter model that assumes an exponential disk; even if they do a marginally better job of fitting the alleged Keplerian downturn, they’re still a lot worse if we consider the terminal velocity data, the details of which are usually ignored.

If instead we pay attention the details of the terminal velocity data, we discover that the broad features seen there in are pretty much what we expect for the kinematic signatures of photometrically known spiral arms. That is, the mass density variations inferred by fitting the RAR correspond to spiral arms that are independently known from star counts. We’ve discussed this before.

If we accept that the bumps and wiggles in the terminal velocities are tracers of bumps and wiggles in the stellar mass profiles, as seen in external galaxies, then we can return to examining the apparent discrepancy between them and the stellar rotation curve from Gaia. The latter follow from an application of the Jeans equation, which helps us sort out the circular motion from the mildly eccentric orbits of many stars. It includes a term that depends on the gradient of the density profile of the stars that trace the gravitational potential. If we assume an exponential disk, then that term is easily calculated. It is slowly and smoothly varying, and has little impact on the outcome. One can explore variations of the assumed scale length of the disk, and these likewise have little impact, leading us to infer that we don’t need to worry about it. The trouble with this inference is that it is predicated on the assumption of a smooth exponential disk. We are implicitly assuming that there are no bumps and wiggles.

The bumps and wiggles are explicitly part of the RAR model. Consequently, the gradient term in the Jeans equation has a modest but important impact on the result. Applying it to the Gaia data, I get the black points:

The velocities of the Gaia data in the range illustrated all go up. This systematic effect reconciles the apparent discrepancy between the stellar and gas rotation curves. The red points are highly discrepant from the gray points, but the black points are not. All it took was to drop the assumption of a smooth exponential profile and calculate the density gradient numerically from the data. This difference has a more pronounced impact on rotation curve fits than any of the differences between the various realizations of the Gaia DR3 data – hence my cavalier attitude towards their error bars. Those are not the important uncertainties.

Indeed, I caution that we still don’t know what the effective circular velocity of the potential is. I’ve made my best guess by splitting the difference between the first and fourth quadrant terminal velocity data, but I’ve surely not got it perfectly right. One might view the difference between the quadrants as the level at which the perfect quantity is practically unknowable. I don’t think it is quite that bad, but I hope I have at least given the reader some flavor for some of the hidden systematic uncertainties that we struggle with in astronomy.

It gets worse! At small radii, there is good reason to be wary of the extent to which terminal velocities represent circular motion. Our Galaxy hosts a strong bar, as artistically depicted here:

Bars are a rich topic in their own right. They are supported by non-circular orbits that maintain their pattern. Consequently, one does not expect gas in the region where the bar is to be on circular orbits. It is not entirely clear how long the bar in our Galaxy is, but it is at least 3 kpc – which is why I have not attempted to fit data interior to that. I do, however, have to account for the mass in that region. So I built a model based on the observed light distribution. It’s a nifty bit of math to work out the equivalent circular velocity corresponding to a triaxial bar structure, so having done it once I’ve not been keen to do it again. This fixes the shape of the rotation curve in the inner region, though the amplitude may shift up and down with the mass-to-light ratio of the stars, which dominate the gravitational potential at small radii. This deserves its own close up:

Here is another place where the terminal velocities disagree with the stellar data. This time, it is because the terminal velocities do not trace circular motion. If we assume they do, then we get what is depicted above, and for many years, that was thought to be the Galactic rotation curve, complete with a pronounced classical bulge. Many decades later, we know the center of the Galaxy is not dominated by a bulge but rather a bar, with concominant non-circular motions – motions that have been observed in the stars and carefully used to reconstruct the equivalent circular velocity curve by Portail et al. (2017). This is exactly what we need to compare to the RAR model.

Note that 2008, when the bar model was constructed, predates 2017 (or the 2016 appearance of the preprint). While it would have been fair to tweak the model as the data improved, this did not prove necessary. The RAR model effectively predicted the inner rotation curve a priori. That’s a considerably more impressive feat than getting the outer slope right, but the model manages both sans effort.

No dark matter model can make an equivalent boast. Indeed, it is not obvious how to do this at all; usually people just make a crude assumption with some convenient approximation like the Hernquist potential and call it a day without bothering to fit the inner data. The obvious prediction for a dark matter model overshoots the inner rotation curve, as there is no room for the cusp predicted in cold dark matter halos – stars dominate the central potential. One can of course invoke feedback to fix this, but it is a post hoc kludge rather than a prediction, and one that isn’t supposed to apply in galaxies as massive as the Milky Way. Unless it needs to, of course.

So, lets’s see – the RAR model Milky Way reconciles the tension between stellar and interstellar velocity data, indicates density bumps that are in the right location to correspond to actual spiral arms, matches the effective circular velocity curve determined for stars in the Galactic bar, correctly predicted the slope of the rotation curve outside the solar circle out to at least 19 kpc, and is consistent with the bulk of the data at much larger radii. That’s a pretty successful model. Some realizations of the Gaia DR3 data are a bit lower than predicted, but others are not. Hopefully our knowledge of the outer rotation curve will continue to improve. Maybe the day will come when the data have improved to the point where the model needs to be tweaked a little bit, but it is not this day.

---

*To give one example, the BICEP II experiment infamously claimed in March of 2014 to have detected the Inflationary signal of primordial gravitational waves in their polarization data. They held a huge press conference to announce the result in clear anticipation of earning a Nobel prize. They did this before releasing the science paper, much less hearing back from a referee. When they did release the science paper, it was immediately obvious on inspection that they had incorrectly estimated the dust foreground. Their signal was just that – excess foreground emission. I could see that in a quick glance at the relevant figure as soon as the paper was made available. Literally – I picked it up, scanned through it, saw the relevant figure, and could immediately spot where they had gone wrong. And yet this huge group of scientists all signed their name to the submitted paper and hyped it as the cosmic “discovery of the century”. Pfft.