2.1. Pure disk galaxies exist, and simulations show that they can naturally evolve to be barred galaxies

Pure disk galaxies exist. Kormendy et al. (Reference Kormendy, Drory, Bender and Cornell2010) obtained HST photometry of six nearby galaxies and measured their surface brightness profiles. They found an upperbound on the total (classical bulge + pseudobulge) stellar mass fraction of ~3% of the total disk mass. They also took an inventory of galaxies within 8 Mpc with $v_{\text{circ}} > 150\,\, \rm {km\, s}^{-1}$ , and found that 11 of 19 showed no evidence for a classical bulge, and four of the remaining eight may contain classical bulges contributing 5–12% of the stellar mass. It is the case, however surprisingly, that pure disk galaxies are common in the low-redshift universe, and thus it is plausible for the Milky Way to be one as well.

Simulations of disk galaxies consistently show that bars and subsequently buckling bars can occur spontaneously in disk galaxies (e.g. Combes & Sanders Reference Combes and Sanders1981; Raha et al. Reference Raha, Sellwood, James and Kahn1991; Athanassoula Reference Athanassoula2005; Martinez-Valpuesta, Shlosman, & Heller Reference Martinez-Valpuesta, Shlosman and Heller2006; Quillen et al. Reference Quillen, Minchev, Sharma, Qin and Di Matteo2014; Spinoso et al. Reference Spinoso, Bonoli and Dotti2017). The rate of evolution of the bar is sensitive to various variables, such as the ratio of the disk to halo mass (Combes & Sanders Reference Combes and Sanders1981). The bar may have a north–south asymmetry at first, but this rapidly dissipates (Raha et al. Reference Raha, Sellwood, James and Kahn1991). Overtime, the bars tend to grow in vertical extent, for example, via multiple, recurrent buckling episodes (Martinez-Valpuesta et al. Reference Martinez-Valpuesta, Shlosman and Heller2006, see Figure 3).

Separately, a bar does form in the gas-rich cosmological “zoom-in” simulation of Spinoso et al. (Reference Spinoso, Bonoli and Dotti2017). The bar's properties are a sensitive function of redshift and related to both gas transport and AGN feedback, suggesting the need for diligence in associating bar properties to Galactic archaeology.

The combination of these two facts, that pure disk galaxies are common in the local universe, and that simulations predict that pure disk galaxies can spontaneously evolve to have bars, render it a plausible model for the Milky Way’s bulge as well. The necessary initial conditions—a disk galaxy—are ubiquitous, and the necessary evolution is natural.

2.2. Radial velocity measurements from large spectroscopic surveys

The Milky Way is a barred galaxy (Dwek et al. Reference Dwek1995), however it has not been clear until recently how quantitatively dominant the bar is relative to the hypothetical classical bulge contribution. The advent of large spectroscopic surveys has allowed us to disentangle the bulge into plausible subcomponents.

Shen et al. (Reference Shen2010) fit a suite of N-body models to bulge radial velocity data from the BRAVA survey, specifically the mean radial velocity (and thus rotation curve) and radial velocity dispersion as a function of longitude and latitude. They recover a best-fit viewing angle, between the bar’s major axis and the line of sight between the Sun and the Galactic centre, of α_Bar = 20°, consistent with more recent measurements. More significantly, they show that models with significant classical bulges, with parameters defined to lie on the fundamental plane of ellipticals and bulges (Kormendy et al. Reference Kormendy, Fisher, Cornell and Bender2009), are not well fit by the data. They constrain the classical bulge mass of the Milky Way to be no more than 8% of the total mass of the disk.

Ness et al. (Reference Ness2013) analysed the data from the larger ARGOS survey, specifically 17 400 red giants with best-fit parameters including R _GC ⩽ 3.5 Kpc, which were spread over 30° of longitude and 20° of latitude. They find several features in the data consistent with N-body expectations for disk origin to the bulge. Bulge stars located between 5 and 10° from the plane ( $0.7\, \rm {Kpc} \lesssim |{\it z}| \lesssim 1.5\, \rm {Kpc}$ ) are cylindrically rotating, with the exception of the 5% of stars with [Fe/H] ⩽−1.0. They suggest that those stars are the inner extensions of the halo and the metal-weak thick disk. It is also the case that the thin disk includes very, very few stars with [Fe/H] ⩽−0.50 (Delgado Mena et al. Reference Delgado Mena2017), so it was always an implausible origin source for the metal-poor component of the bulge.

Di Matteo et al. (Reference Di Matteo2015) compared the ARGOS data to three different N-body models and found that the bulge could not be explained as having resulted from the evolution of a pure thin disk galaxy, and thus a joint origin with the thick disk is likely required as well. This may appear as requiring too many free parameters from models, but this is a parameter that exists in nature, the Milky Way is not a pure thin disk galaxy, it also has a thick disk. Di Matteo et al. (Reference Di Matteo2015) show that if the bulge was purely due to a thin disk buckling event, the curve of mean radial velocity versus longitude would decrease in amplitude, and the velocity dispersions would decrease, with increasing metallicity. The decrease in the amplitude of the rotation curve is seen for stars satisfying [Fe/H] ≳ −0.50, but not for lower metallicities. The velocity dispersions do decrease, but not as quickly as predicted from a pure thin disk model. Di Matteo et al. (Reference Di Matteo2015) suggest that the solution to this issue is that bulge stars with [Fe/H] ≲ −0.50 formed from the thick disk, rather than from the more metal-poor component of the thin disk.

Zasowski et al. (Reference Zasowski2016) investigated the skewness (third moment) and kurtosis (fourth moments) of the radial velocity distribution functions in APOGEE data towards the bulge. The correlation between skewness and mean velocity (first moment), a known diagnostic of bars, is only observed for the more metal-rich fraction of the stars, defined in that paper as stars with [Fe/H] ≳ −0.40. The data have a flat kurtotsis, Kurt (V) ≈0, consistent with that expected from models. The data have slightly higher velocity dispersion and slightly lower skewness than expected from N-body models of a simple single disk galaxy undergoing buckling.

The recent review of Di Matteo (Reference Di Matteo2016) goes over these lines of evidence in greater detail than possible here, and concludes that bulge stars with [Fe/H] ≳ −1.0 formed from buckling thick and thin disks.

2.3. The peanut/X-shaped spatial distribution of bulge stars

The sightline-dependent bifurcation of the apparent magnitude distribution of red clump stars (comprising ~99% of bulge horizontal branch stars, see Nataf et al. Reference Nataf2013) was first reported by Rattenbury et al. (Reference Rattenbury, Mao, Sumi and Smith2007). They speculated that it might be due to a distinct population lying at a different distance. Nataf et al. (Reference Nataf, Udalski, Gould, Fouqué and Stanek2010) and McWilliam & Zoccali (Reference McWilliam and Zoccali2010) showed that the double peak was spread across a large swath of the bulge, at large separations (|b| ≳ 5°) from the plane. McWilliam & Zoccali (Reference McWilliam and Zoccali2010) argued convincingly, that the feature had to be due to an X-shaped bulge (a very strong peanut shape), a feature of strong bars in N-body models.

Subsequent analyses confirmed the claim. Ness et al. (Reference Ness2012) found that the bright and faint red clumps had the same difference in mean radial velocity as expected from N-body models. Li & Shen (Reference Li and Shen2012) showed that the N-body model of Shen et al. (Reference Shen2010), already demonstrated to be a good match to radial velocity observations, also includes an X-shape. The relative brightness, and the dependence of the relative number counts on direction, was qualitatively similar to that seen in the data. Vásquez et al. (Reference Vásquez2013) showed that the proper motion distributions were consistent. Nataf et al. (Reference Nataf2015) showed that precision measurements of the brightness difference, relative number counts, and total number density of stars in the data from OGLE-III (Udalski et al. Reference Udalski, Szymanski, Soszynski and Poleski2008) were qualitatively matched by N-body models from Athanassoula (Reference Athanassoula2003) and Shen et al. (Reference Shen2010). An unambiguous image of the X-shaped bulge in integrated mid-infrared photometry of the Milky Way is shown in Figure 1, and the development of one from an N-body model is shown in Figure 2.

Figure 1. Figure 1 from Ness & Lang (Reference Ness and Lang2016), the X-shaped bulge is unambiguous in the integrated mid-infrared photometry of the Milky Way.

Figure 2. Figure 3 from Martinez-Valpuesta et al. (Reference Martinez-Valpuesta, Shlosman and Heller2006). An X-shaped bulge is a natural outcome of N-body models of disk galaxies undergoing the buckling instability.

Precision modelling of the photometric data by Wegg & Gerhard (Reference Wegg and Gerhard2013) has been followed by detailed dynamical modelling by Portail, Wegg, & Gerhard (Reference Portail, Wegg and Gerhard2015a) and Portail et al. (Reference Portail, Wegg, Gerhard and Martinez-Valpuesta2015b, Reference Portail, Gerhard, Wegg and Ness2017). A suite of N-body models coarsely consistent with Milky Way constraints was adjusted to be more consistent using the made-to-measure method (Syer & Tremaine Reference Syer and Tremaine1996), whereby weights of the particles in the model are shifted up or down to match the observational data. The end models always contain a significant X-shape, comprising 40–50% of the stellar mass of the bulge, but only dominant for stars with [Fe/H] ≳ −0.50. That is not a repetition of the results mentioned in the previous subsection, as Portail et al. (Reference Portail, Wegg and Gerhard2015a, Reference Portail, Wegg, Gerhard and Martinez-Valpuesta2015b, Reference Portail, Gerhard, Wegg and Ness2017) also constrained their N-body models to match the spectroscopic data from the APOGEE survey, which goes much closer to the plane (Ness et al. Reference Ness2016), and the global photometric parameters of the bulge which trace the number density and distance distribution function of stars (Wegg & Gerhard Reference Wegg and Gerhard2013).

Reiterating, this is an impressive series of observational tests for models to have passed. First, the alternative hypotheses (RGB bump, metallicity differences, etc.) do not work. Second, N-body model predictions qualitatively match the radial velocity offsets, the proper motion measurements, the mean brightness difference, and the dependence of the relative and total number of counts as a function of direction. This picture survives, and is in fact dominant for stars with [Fe/H] ≳ −0.50, when models are required to match observations from several different surveys.

2.4. The morphology of the long bar is another specific, precise prediction of these models

An issue of Galactic modelling in the last two decades is that of the long bar of the Milky Way (Benjamin et al. Reference Benjamin2005; Cabrera-Lavers et al. Reference Cabrera-Lavers2007, Reference Cabrera-Lavers, González-Fernández, Garzón, Hammersley and López-Corredoira2008). Evidence has accumulated from multiple investigations of a long, in-plane bar, with a half-length of ~4 Kpc, and an orientation angle of α_Bar = 45°. This represented a challenge to Galactic structure studies, as the orientation angle is far larger than the value found for the triaxial bulge towards the inner few Kpc, and the morphology is not consistent with the triaxial ellipsoid models used to fit for the bulge (Dwek et al. Reference Dwek1995; Stanek et al. Reference Stanek1997). It appears as a double-barred system where the two bars are of comparable length and are not aligned, which is an unstable configuration (Athanassoula Reference Athanassoula2012).

A suggestion to this resolution was independently proposed by Martinez-Valpuesta & Gerhard (Reference Martinez-Valpuesta and Gerhard2011) and Athanassoula (Reference Athanassoula2012). Their idea was the analytical triaxial ellipsoid models then widely used to model the bulge were only (somewhat) suitable for the peanut/X-shaped component, whereas bars in simulations often have long and thin extensions. This hypothesis made a specific prediction, that the difference in measured angle Δα_Bar, then around 20°, was due to observational errors. A smaller difference might remain if the ends of the bar develop interactions with spiral arms Martinez-Valpuesta & Gerhard (Reference Martinez-Valpuesta and Gerhard2011).

This prediction was confirmed by Wegg, Gerhard, & Portail (Reference Wegg, Gerhard and Portail2015), who combined photometry from Spitzer-GLIMPSE (Benjamin et al. Reference Benjamin2005), 2MASS (Skrutskie et al. Reference Skrutskie2006), VVV (Inoue & Saitoh Reference Inoue and Saitoh2012), and UKIDSS (Lucas et al. Reference Lucas2008) to make a global map of the bulge. The long bar had an angle of α_Bar = 28–33°, consistent with the values of α_Bar = 27 ± 2° and α_Bar ≈ 29° then recently measured by Wegg & Gerhard (Reference Wegg and Gerhard2013) and Cao et al. (Reference Cao, Mao, Nataf, Rattenbury and Gould2013), respectively, for the inner peanut/X-shaped component.

Thus, the issue of the long bar, which was previously a challenge to Galactic structure models, ended up being a triumph. The long bar is not only not surprising, but expected from simulations of buckling disk galaxies. Models predicted that the difference in orientation angle should shrink with better analysis, and it did.

2.5. The metallicity gradient can also be reproduced by buckling disk

The Galactic bulge has a vertical metallicity gradient, with stars further from the plane having a lower mean metallicity. Tiede, Frogel, & Terndrup (Reference Tiede, Frogel and Terndrup1995) used near-IR photometry and estimated a gradient in ∇[Fe/H] =−0.06 ± 0.03 dex deg⁻¹ or −0.43 ± 0.21 dex kpc⁻¹ between b = −3° and b = −12°. The metallicity gradient was independently and concurrently confirmed by Minniti et al. (Reference Minniti1995), and since confirmed numerous times with spectroscopic data (e.g. Zoccali et al. Reference Zoccali2008; Johnson et al. Reference Johnson, Rich, Fulbright, Valenti and McWilliam2011).

This was argued by Minniti et al. (Reference Minniti1995) to be evidence for a dissipative collapse, as models of spheroid formation predicted a radial metallicity gradient (Carlberg Reference Carlberg1984). However, Minniti et al. (Reference Minniti1995) acknowledged that ‘the alternative interpretation that the gradient itself is caused by the mixing of different components in the inner Galaxy cannot be ruled out’.

That second point, that the bulge has different components with both different metallicity distribution functions and scale heights, is very much the case, as should be clear from the literature evidence summarised in this review. Further, it also turns out that a metallicity gradient for the bulge can emerge from a pure disk galaxy. Martinez-Valpuesta & Gerhard (Reference Martinez-Valpuesta and Gerhard2013) evolved an N-body model of a pure disk galaxy which included an initial radial metallicity gradient. They chose [M/H](R) = +0.60 − 0.40(R) Kpc⁻¹, where ‘R’ is the galactocentric radius of particles at the start of the simulation. Given that particles from initial radii are scattered to different orbits, a metallicity gradient emerges which is qualitatively consistent with that observed for the whole bulge—there is impressive agreement with the global photometric metallicity maps from Gonzalez et al. (Reference Gonzalez2013), which are derived from the morphology of the red giant branch. That said, it is worth noting that the metallicity gradient required by Martinez-Valpuesta & Gerhard (Reference Martinez-Valpuesta and Gerhard2013) is extremely large.

This result is not altogether surprising as it has been known for a while that dynamical mixing processes do not erase gradients. For example, White (Reference White1980) demonstrated that simulations of mergers were expected to reduce, but not erase, gradients in metallicity. Sánchez et al. (Reference Sánchez2014) find that the metallicity gradients in interacting galaxies, measured in terms of oxygen abundance of H ii regions versus effective radius, is reduced by ~1/2 relative to non-interacting galaxies.

The metallicity gradient of the bulge, first largely argued to be evidence for dissipative collapse and thus evidence against a disk–galaxy origin, turns out not to be discriminating. A metallicity gradient can be matched by models of dissipative collapse, pure disk galaxies undergoing the buckling instability, a combination of the two, or as we will see, the clump-origin bulge scenario (Inoue & Saitoh Reference Inoue and Saitoh2012). The existence of a metallicity gradient, by itself, yields no significant constraints bulge formation.

2.6. Caveat against the bar scenario: Gas fractions and the challenge of initial conditions

The properties and evolution of bars in N-body simulations of pure exponential disks are a well-researched subject. However, real galaxies (such as the Milky Way) contain gas, and generally contained more gas in the past.

Athanassoula, Machado, & Rodionov (Reference Athanassoula, Machado and Rodionov2013) investigated the relative predicted properties of bars in N-body simulations with and without gas. They found that the gas-rich galaxies remain axisymmetric for longer. When they do develop bars, they do so at a slower rate, and end up much weaker. Their Figure 7 shows that bar strengths in gas-rich galaxies are predicted to end up lower even after the gas has been depleted. In contrast, as discussed in this section, the Milky Way has a very strong bar.

Indeed, the bar fraction is considerably lower at high redshift. Melvin et al. (Reference Melvin2014) find that it drops to 11% at z = 1, corresponding to a lookback time of 7.8 Gyr. Simmons et al. (Reference Simmons2014) also estimate a bar fraction of 11%, in the redshift range 0.50 < z < 2.0. By that time, most Milky Way bulge stars were already formed (Clarkson et al. Reference Clarkson2011; Bensby et al. Reference Bensby2017) and thus they already have a kinematic distribution.

The combination of these issues should give pause to the notion that a bar + buckling instability from a pure exponential disk is independently sufficient to explain the inner Milky Way’s dynamics. Separately from the issue of gas weakening bars, subsequent sections of this review will discuss how a minor classical bulge can actually strengthen the bar, and why the Milky Way is expected to be a former ‘clumpy galaxy’ which has implications for bulge formation.

3.1. A classical bulge is expected from theory

The most obvious advantages of this scenario are that it is predicted by straightforward theory of gravitational, dissipative collapse (Eggen et al. Reference Eggen, Lynden-Bell and Sandage1962), and hierarchical clustering in a cold dark matter universe (White & Rees Reference White and Rees1978; Kauffmann, White, & Guiderdoni Reference Kauffmann, White and Guiderdoni1993).

White & Rees (Reference White and Rees1978) proposed that most of the matter in the universe condensed early into small ‘dark’ objects. They suggested a model with Ω_m = 0.20 and dark matter making up 80% of the matter, impressively similar to the modern values (Planck Collaboration et al. 2014). Within this picture, the ‘pure disk galaxy’ that has particles and no gas is not a viable initial condition, as the universe and star formation begin with larger number of small haloes that coalesce and accrete additional small haloes.

As far as current, more up-to-date models are concerned, a classical bulge can be considered a requirement. Kobayashi & Nakasato (Reference Kobayashi and Nakasato2011) simulate 150 galaxies from various cosmologically motivated initial conditions. They formed disk structures in 48 of their galaxies, including five galaxies that had masses comparable to the Milky Way. Though some of their simulated galaxies had small classical bulges, none completely lacked a classical bulge. It is a concern that many local galaxies are found not to have classical bulges locally (Kormendy et al. Reference Kormendy, Drory, Bender and Cornell2010), but one can argue against that by saying that information is lost when looking at those galaxies in integrated light. That cannot be argued for the Milky Way.

3.2. The kinematics of metal-poor bulge stars

There has been much written in this review of the spectacular consistency between N-body models of buckling disks and the dynamics of bulge stars with [Fe/H] ≳ −1.0 or ≳ −0.50, depending on the study. These consistencies are not found for the metal-poor bulge stars. For example, the Galactic bar is either null (Dékány et al. Reference Dékány2013) or weak (Pietrukowicz et al. Reference Pietrukowicz2012, Reference Pietrukowicz2015) in the RR Lyrae stars, variable stars of standardisable distance that can be used to directly probe the metal-poor bulge. Dékány et al. (Reference Dékány2013) investigated 7 663 fundamental-mode RR Lyrae in I and K _s bands, and found that they do not trace a strong bar, but rather a more spheroidal and centrally concentrated distribution. No correlation is observed between dereddened distance and longitude, a correlation detected at high significance in the more metal-rich red clump stars which trace a bar (Stanek et al. Reference Stanek1997).

Pietrukowicz et al. (Reference Pietrukowicz2012) and Pietrukowicz et al. (Reference Pietrukowicz2015) obtained different results in their analysis of 16 836 and 27 258 RR Lyrae in the OGLE-III and OGLE-IV surveys, respectively. They do find a bar, but it is kinematically hotter than the bar measured in red clump stars, and does not show a peanut/X-shape at large separations from the plane.

Analysis of RR Lyrae kinematics further confirms these findings. Kunder et al. (Reference Kunder2016) analyse spectroscopic data for 947 RR Lyrae as part of the ongoing BRAVA-RR survey. These RR Lyrae, measured towards |l| ≲ 4° and −6° ≲ b ≲ −3° show higher velocity dispersions and weaker rotation than the metal-rich M-giants studied as part of the BRAVA survey. The velocity dispersion is ~15% higher, and the rotation is null or negligible. When the RR Lyrae are split into two metallicity bins with [Fe/H]=−0.75 marking the bifurcation point, no difference in rotation is measured, though the metal-poor RR Lyrae do have a higher velocity dispersion. The RR Lyrae radial velocity measurements are shown in Figure 3.

Figure 3. Left panel of Figure 2 from Kunder et al. (Reference Kunder2016). RR Lyrae stars show null or negligible Galactic rotation, as well as very high velocity dispersion, in contrast to the majority of bulge stars.

A further issue with the RR Lyrae is the very fact that they are RR Lyrae. The bulge RR Lyrae are measured to have a mean metallicity of [Fe/H]≈−1.0, whereas metallicities of [Fe/H]≈−1.60 are more typical of RR Lyrae in the globular clusters and in the halo (Dotter et al. Reference Dotter2010). More metal-rich stars have a higher turnoff mass at fixed age, and thus for these horizontal branch stars to be on the instability strip at a higher metallicity, they need to have a lower turnoff mass by other means, and thus likely a greater age. Lee (Reference Lee1992) estimated that the bulge RR Lyrae had to be Δt ~ 1.3 ± 0.30 Gyr older than halo stars, for stars formed with metallicities − 1.5 ≲ [Fe/H] ≲ -1.0. A bulge which is older than the halo is more consistent with a classical bulge then say, a buckling disk origin.

This argument also applies to bulge globular clusters. For example, NGC 6522, has a blue horizontal branch with a mean metallicity of [Fe/H]=−1.0 (Barbuy et al. Reference Barbuy2009), necessitating an extremely old age. This is a consistent pattern of bulge globular clusters (Barbuy et al. Reference Barbuy2006, Reference Barbuy2007). This information, combined with that of the RR Lyrae, strongly suggests that bulge stars of metallicity [Fe/H]=−1.0 are the oldest or among the oldest stellar populations in the Galaxy, and thus were around prior to there being a massive disk that could form a bar.

3.3. The theoretical interaction between classical bulges and buckling disks is consistent with observations

There is relatively sparse research on the theoretically predicted interaction of a central concentration of stars on a buckling disk, but what research is available now turns out to be consistent with the data. This is true both of the predictions of the classical bulge behaviour and the bar behaviour.

Athanassoula & Misiriotis (Reference Athanassoula and Misiriotis2002) studied three different N-body models of disk galaxies with variable initial central concentrations. The disk mass and disk-to-halo mass ratios are fixed. Two models have a non-centrally concentrated halo, one of those two also has a bulge. A third model has a centrally concentrated halo, but no bulge. Within this suite of models, the galaxy with a non-centrally concentrated halo and no bulge ends up forming the weakest bar, it does not develop cylindrical rotation, and its boxy shape does not evolve to an X-shape. Interestingly, the model with a centrally concentrated halo but no bulge ended up with the strongest bar. This is suggestive that greater central concentration was needed given the exceptionally strong bar of the Milky Way, though in and of itself, it is not sufficient as more models would be needed for a more convincing picture.

Conversely, Saha et al. (Reference Saha, Martinez-Valpuesta and Gerhard2012), Saha & Gerhard (Reference Saha and Gerhard2013), and Saha, Gerhard, & Martinez-Valpuesta (Reference Saha, Gerhard and Martinez-Valpuesta2016) have investigated how a joint classical bulge/buckling disk origin impacts the development of the bulge by means of N-body models. Saha et al. (Reference Saha, Martinez-Valpuesta and Gerhard2012) find that a small (7% of the total disk mass) classical bulge can pick up angular momentum from the larger rotating bar, and thus even develop into a triaxial object with cylindrical rotation. Saha & Gerhard (Reference Saha and Gerhard2013) show that the composite bulge always ends up rotating cylindrically, but may have deviations from cylindrical rotation at specific moments in its evolution. Curiously, they find that the final size of the composite bulges are reduced if the initial classical bulge has its own angular momentum. Saha et al. (Reference Saha, Gerhard and Martinez-Valpuesta2016) extend their prior results to show that even massive classical bulges can pick up as much specific angular momentum as low-mass classical bulges, but that the resulting rotation is non-cylindrical. All composite systems eventually form a boxy/peanut bulge.

Pérez-Villegas, Portail, & Gerhard (Reference Pérez-Villegas, Portail and Gerhard2017) model the evolution of a disk galaxy with a halo via an N-body model and find that the properties of the RR Lyrae are consistent with that of an inner halo, specifically their number density, their slow rotation, the lack of a peanut/X-shape in their spatial distribution, and their higher velocity dispersion. Within their model, only 12% of RR Lyrae end trapped on bar-like orbits, which can be considered as a prediction that the fraction will be very small once proper motions are available and once astronomers can compute orbits for these stars.

It is been established in the previous sections of this review that observations and analysis thereof have ruled out a predominantly classical bulge for the Milky Way bulge. However, the N-body models are also consistent with either a centrally concentrated inner halo, or an additional classical bulge in addition to that. Such a feature could help explain the Milky Way’s very strong bar as well as the behaviour of low-metallicity bulge stars. Further analysis and more data of metal-poor bulge stars are needed.

3.4. The bulge trends in the α-elements are not consistent with those of the thin and thick disks

The abundances of the α-elements relative to iron, [α/Fe], are a historic argument for a distinct formation to the bulge and the disk. That is because, the [α/Fe] abundance ratios trace the efficiency of star formation and possibly the initial mass function of a stellar population (Matteucci & Greggio Reference Matteucci and Greggio1986). A higher level of [α/Fe] can be due to a lower contribution to the chemical enrichment of interstellar gas from type Ia SNe (which take longer to form), and thus star formation would need to be faster.

The landmark study of McWilliam & Rich (Reference McWilliam and Rich1994) found that the trends for magnesium and titanium were enhanced by ≈0.30 dex relative to the solar neighbourhood trends over the full range of [Fe/H], whereas the trends for calcium and silicon were consistent with those of the disk. These different ratios suggested a different origin, and McWilliam & Rich (Reference McWilliam and Rich1994) said it may reflect a common enrichment process between bulges and ellipticals. Fulbright, McWilliam, & Rich (Reference Fulbright, McWilliam and Rich2007) compiled a more sophisticated analysis of a higher resolution, higher signal-to-noise sample. They found that the bulge has a magnesium trend elevated by ~0.30 dex relative to the disk, whereas the trends of oxygen, silicon, calcium, and titanium are slightly elevated relative to the disk. The bulge also has higher aluminum abundances at fixed iron abundance. Fulbright et al. (Reference Fulbright, McWilliam and Rich2007) conclude that the relative abundance offsets between the bulge and the disk are inconsistent with models where the bulge forms from the buckling of the disk. Separately, Fulbright et al. (Reference Fulbright, McWilliam and Rich2007) found that the metal-poor bulge stars also show higher mean abundances of silicon, calcium, and titanium than the halo. They concluded that the metal-poor bulge stars could not have formed from gas with the present-day halo composition.

The discussion shifted with the work of Alves-Brito et al. (Reference Alves-Brito, Meléndez, Asplund, Ramírez and Yong2010). That investigation used high resolution of optical spectra of 25 bulge giants and 55 comparison giants and analysed abundances in a homogeneous manner to minimise systematic offsets. Their results were that metal-poor bulge stars ([Fe/H] ≲ −0.50) have the same abundances as the thick disk, and more metal-rich stars have the same abundances as the thin disk.

The relative abundance offsets between the bulge and the disk have not converged to zero as more data have come in. Bensby et al. (Reference Bensby2017) compiled what is among the best datasets of bulge abundances, as they have high-resolution, high signal-to-noise abundances for 90 bulge stars located on the main-sequence turnoff and subgiant branch, analysed using the same methods as their comparison disk sample, which is composed of stars in the solar neighbourhood. Some elements trace the same abundance trends as the disk, but magnesium, titanium, and aluminum do not. They typically trace the upper end of the larger distributions spanned by the thick and thin disks. The abundance trend of [Mg/Fe] vs. [Fe/H] is shown in Figure 4.

Figure 4. From Figure 21 of Bensby et al. (Reference Bensby2017). The [Mg/Fe] vs. [Fe/H] abundance trend for bulge stars (coloured points) are shown superimposed on the disk trends (grey points). The magnesium abundances are elevated with respect to the disk abundances at all [Fe/H] values.

Relative abundance offsets between the bulge and the disk remain, which is a challenge to models of the Galaxy where the bulge is simply due to buckling thin and thick disks. As the abundance ratios are higher in magnesium in particular, suggesting a more rapid star formation, the suggestion is that the bulge formed faster than the disk, characteristic of an early, dissipative collapse.

4.1. This is what high-redshift galaxies actually look like

Guo et al. (Reference Guo2015) analysed 3 239 high-redshift, star-forming galaxies studied as part of the Cosmic Assembly Near-Infrared Deep Extragalactic Legacy Survey (CANDELS). They use a conservative definition for clumps—a clump has to contribute at least 8% of the UV light of a galaxy, which excludes smaller clumps. They also require clumps to be off-centre, which excludes central clumps and is thus limiting in our context as a central clump could obviously contribute to bulge formation. One of their mosaics of clumpy galaxies is shown in Figure 5.

Figure 5. Figure 4 from Guo et al. (Reference Guo2015) showing HST images of high-redshift galaxies. The norm is that of massive, star-forming clumps within disks, rather than simple and smooth exponential disks.

They find that galaxies with log(M _⋆./M_⊙) > 9.8 have a 55% probability of being clumpy at redshift z ~ 3, down to 15% at z ~ 0.5, which is largely due to the fact star formation declines with decreasing redshift. At all redshifts, the clump contribution to rest-frame UV light peaks at log(M _⋆./M_⊙) > 10.5—the current stellar mass of the Milky Way. Integrating over both clumpy and non-clumpy galaxies, they find that 4–10% of the star formation takes place within these massive clumps.

In other words, the clumpy galaxy is a very plausible assumption for the initial conditions of the Milky Way. Much has been said in this review that the pure disk galaxy can work, as Kormendy et al. (Reference Kormendy, Drory, Bender and Cornell2010) has pointed out. However, the smooth and massive exponential disks are widely seen in observations of the local universe. That clumpy galaxies are the norm for high-redshift observations suggests that most massive galaxies have passed through a clumpy phase.

4.2. Promising insights from simulations

The prevalence of clumps in high-redshift galaxies renders them a legitimate point of discussion for the origin of the Milky Way.

A way to test this is with comparison to the Milky Way, which Inoue & Saitoh (Reference Inoue and Saitoh2012) did. They used an N-body/SPH model to study the evolution of an isolated disk galaxy where the clumps migrate to the centre via dynamical friction and form a clump-origin bulge. The final bulge resembles what they call a pseudo bulge, and is referred to as a bar throughout this paper and most of the Milky Way literature. The surface density profile is nearly exponential, the final shape of the bar is boxy, and the rotation is significant. The resulting stars are old and metal rich, with a flat star-formation history in an interval of Δt ~ 2 Gyr followed by a rapid decline in star formation. They obtain a metallicity gradient that stretches across the full vertical extent of the bulge. All of these properties are qualitatively consistent with what is observed for the Milky Way bulge. However, these properties can be matched by other models, and further this is a comparison of a single clump-origin bulge model to those of the bulge. It is worthy of consideration, but it is far too premature to declare victory.

One concern is that of whether or not the clumps actually do migrate to the centre, as they do the simulations previously discussed in this review. This is largely a question for theory, as the migration duration is too long for the baseline of observations. In the simulations of Hopkins et al. (Reference Hopkins, Kereš, Murray, Quataert and Hernquist2012), the inclusion of their prescription for stellar feedback disrupted the clumps, and prevented them from migrating to the bulge where they can coalesce. This prediction is not reproduced by the simulations of Bournaud et al. (Reference Bournaud2014), who find that the ejection of stars in clumps due to stellar feedback is compensated by their accretion of gas from the gas-rich disks in which they are contained.

Mandelker et al. (Reference Mandelker2014) studied 770 snapshots of 20 simulated galaxies using adaptive-mesh refinement cosmological models of Galaxy evolution. The global number of clumps were consistent with those in observations, an important check of the violent-disk instability hypothesis of clump formation. They did not study the properties of the final resulting bulge, but did say that if the clumps can survive accretion onto the centre of the galaxy, they are expected to accrete gas from the surrounding interstellar medium (similarly to Bournaud et al. Reference Bournaud2014), and will thus show gradients in their mean properties with respect to separation from the centre of their galaxy, such as those measured by Förster Schreiber et al. (Reference Förster xsSchreiber2011) that clumps closer to the centres of their disks are redder, older, and more massive. Further, they found in their simulations was that a full 91% of galaxies develop a bulge clump. These are massive, typically equivalent to 40% of the disk mass, with 20% of the star formation, and gas fractions of less than 1%.

Mandelker et al. (Reference Mandelker2017) study 34 galaxies with more sophisticated prescriptions. Among their findings, they find that the inclusion of radiation pressure disrupts the smaller clumps, reducing their lifetimes to a few free-fall times, but that the more massive and dense clumps still nevertheless survive and migrate to the centre. The inclusion of radiation pressure reduces the number of long-lived clumps by 81%. Radiation pressure has little to no effect on the bulge clumps, with ~83% of simulated galaxies hosting a bulge clump.