Metallicity Properties of the Galactic Bulge Stars Near and Far: Expectations from the Auriga Simulation

In Fig. 2, we divide the stars in AU23 into different Galactic latitude slices, from $b=-9^{\circ}$ to $b=9^{\circ}$ with $1^{\circ}$ height for each slice. We further superpose the opposite latitudinal slices ($\pm b$) to increase the number statistics. In each panel of Fig. 2 ([Fe/H] color-coded), the contour represents the number density distribution of stars younger than 11 Gyr to highlight the position of the double peaks, useful to identify the stars belonging to the X-shape (similar to observations). In order to better visualize the density contours, larger bins are used in the slices at $|b|>6^{\circ}$. Including all the stars will produce a smooth distribution. In the slices at $|b|\simeq 5^{\circ}$, there are two clear over-density regions, corresponding to the double peaks in the apparent magnitude distributions of the RC stars as in observations (McWilliam & Zoccali, 2010; Nataf et al., 2010) and previous simulation works (e.g. Li & Shen 2012; Ness et al. 2012).

Close to the mid-plane, such over-density regions gradually approach to each other and merger together in slices $|b|\leq 3^{\circ}$ (P1 and P2). In slices $|b|\geq 4^{\circ}$ the double peaks are visible only in $|l|\lesssim 2.5^{\circ}$ for $|b|=4^{\circ}$ and $|l|\lesssim 3.5^{\circ}$ for $|b|=5^{\circ}$, and different peaks dominate different sides (P3).

From the metallicity distribution in Fig. 2, we can see that in slices $|b|\geq 4^{\circ}$, the two over-density regions are more metal-rich than the region between them. In slices $|b|\leq 3^{\circ}$, the stellar metallicity is higher close to GC and along the major axis of the bar. The slices closer to the mid-plane are more metal-rich. In slices $|b|>4^{\circ}$, the near side (bright RCs) seems slightly more metal-rich than the far-side (faint RCs).

According to the stellar distance from the GC, we can divide the stars into two groups: 4 kpc$<d<$8 kpc as the near side, and 8 kpc$<d<$12 kpc as the far side. By analyzing the metallicity properties, we can get a qualitative picture of the metallicity difference between the near and far sides of the bulge.

Fig. 3 shows the [Fe/H] distribution in the $l-b$ diagram of the two populations. As shown in the left two panels, the most metal-rich region of the two populations are located at different Galactic longitude sides. The right panel shows the [Fe/H] difference of the two populations, with $\Delta$[Fe/H]=[Fe/H]_near $-$ [Fe/H]_far (the difference between previous two panels). In the high latitude region ($|b|>6^{\circ}$) the near side is more metal-rich than the farther one at all longitudes. This is due to the negative vertical metallicity gradient of the MW; along certain line-of-sight towards the Galactic bulge region, the fainter/farther peak is located at larger vertical height above the disk plane than brighter/closer one. In addition, in the right panel we can see a horizontal [Fe/H] gradient in the intermediate latitude region. Close to the mid-plane ($|b|\leq 6^{\circ}$) at negative longitude side ($l<0^{\circ}$), the metallicity difference even becomes negative. This is unexpected, since the vertical metallicity gradient of the disk should indicate the near side to be more metal-rich than the far side.

To better understand the variation of $\Delta$[Fe/H] across the bulge region, we show the faced-on and edge-on maps of AU23 color-coded with the mean [Fe/H] in the left panel of Fig. 4. The metallicity is larger along the major axis of the bar than that along the minor axis in the face-on map. The more metal-rich (younger) stars are more concentrated around the bar major axis with higher metallicity, while the metal-poor (older) stars show rounder distribution with smaller ellipticity, as shown in Fig. A1 of Fragkoudi et al. 2020. It can help to understand the phenomena in Fig. 3. The right panels in Fig. 4 show the latitude-distance distributions of stars in two vertical slices at $l=\pm 3^{\circ}$. The black solid line marks the position of the GC distance at 8 kpc. For the sight line at $l=3^{\circ}$, close to the mid-plane ($|b|\leq 6^{\circ}$), the near side is dominated by the stars in the major axis of the bar (metal-rich), while the far side is populated by the stars in the minor axis, which is more metal-poor. On the contrary, for the sight line at $l=-3^{\circ}$ and $|b|\leq 6^{\circ}$, the near side is mainly populated by stars in the minor axis region (metal-poor) and the far side is mainly dominated by the major axis region of the bar (metal-rich). This difference further enhances the horizontal metallicity gradient close to the mid-plane. Interestingly, in the left panels of Fig. 4 the most metal-rich region of the bar seems to coincide with the central boxy core identified in Li & Shen (2015) for numerical buckled bars.

In the region away from the mid-plane, the metallicities of both the near and far sides decrease, with the far side showing a steeper latitudinal gradient. Therefore, at both positive and negative longitude regions, the near side gradually becomes more metal-rich than the far side. This is demonstrated in the lower panel of an edge-on projection of the model with two lines of sight at $b=\pm 5^{\circ}$. The vertical metallicity gradient is more important to result in a positive metallicity difference between the near and far sides. At such latitude range we can still find a longitudinal gradient of the metallicity difference. The negative longitude side generally shows smaller (or even negative) metallicity difference. This difference arises from the fact that at negative longitude region, the near side is mainly contributed by the more metal-poor minor axis of the bar than the more metal-rich far side around the major axis. The latitudinal metallicity gradient decreases faster in the far side, to gradually result in a positive metallicity difference at larger latitude regions.

The results in Section 3.2 hint that the bright (near) and faint (far) peaks of the X-shape may show similar patterns of the metallicity difference across the Galactic bulge region. To quantify such phenomenon more precisely, we transform the stellar spatial distribution into distance histograms for different $(l,\,\,b)$ fields, similar to Fig. 2 in Li & Shen (2012). To better visualize the double peaks in the distance distribution, we mainly focus on stars younger than 11 Gyr (as shown in Fig. 2) to determine the peak position. Note that the double peaks can be identified in the region $|l|<3^{\circ}$ and $3^{\circ}\lesssim|b|\lesssim 8^{\circ}$. Then we apply a Gaussian mixture model provided by scikit-learn (Pedregosa et al., 2012), which implements an expectation-maximization algorithm, to identify the location of the double peaks. More details are shown in the Appendix. In the following analysis, we pick out stars (including all age populations) that lie within 0.25 kpc distance from the peaks to study their metallicity.

The left column of Fig. 5 shows the mean metallicity and metallicity differences as a function of the longitude in different latitude slices from $|b|=4^{\circ}$ to $8^{\circ}$ of the bright and faint double peaks. For lower Galactic latitude slices ($|b|=4^{\circ},\,\,5^{\circ},\,\,6^{\circ}$), with $l$ increasing, [Fe/H] of the brighter (fainter) peak gradually increases (decreases). $\Delta$[Fe/H] gradually increases from negative (at $l<0^{\circ}$) to positive values (at $l>0^{\circ}$) due to the metallicity difference between the major and minor axis of the bar demonstrated in Fig. 4. In slices $|b|\geq 6^{\circ}$ almost all the metallicity difference values are positive due to the vertical metallicity gradient.

The right column of Fig. 5 shows the similar result but in different Galactic longitude slices ($l\,=\,-3^{\circ},\,-1^{\circ},\,0^{\circ},\,1^{\circ},\,3^{\circ}$). At each slice, the brighter peak always has a shallower latitudinal metallicity gradient than the fainter one, since at the same latitude the fainter/farther stars are located at larger vertical heights than the brighter/closer stars. In the $l=-3^{\circ}$ slice, $\Delta$[Fe/H] gradually increases from negative to positive values with increasing $|b|$. For $l\geq 0^{\circ}$ slices, at different latitudes, $\Delta$[Fe/H] values are all positive. This has been demonstrated in Fig. 4: at $l\geq 0^{\circ}$, the near side is mainly contributed by the major axis part of the bar, which is more metal-rich than the minor axis part to produce a positive metallicity difference.

Therefore, according to the simulation, we can see that the near side is not always more metal-rich than the far side. The near side is even more metal-poor than the far side at $-3^{\circ}<l<-1^{\circ}$ and $|b|<6^{\circ}$. As we have explained in Section 3.2, it is the metallicity difference between the minor axis (metal-poor) and major axis (metal-rich) of the bar to result in negative $\Delta$[Fe/H] values in such region. In the other fields, the near side is generally more metal-rich, which is mainly due to the vertical metallicity gradient in the bulge/bar structure. Although this model may not truly reflect the chemical properties of the Galactic bulge, it is still valuable to shed light on the spatial variation of the metallicity differences between the near and far sides across the bulge region.

The metallicity difference between the double peaks in the simulation is consistent with the previous observations (Joo et al., 2017; Lee et al., 2019). We note that the observed value of metallicity difference between the doubel RCs in Lim et al. (2021) is $0.149\pm 0.036$ dex at $(l,\,\,b)=(-1^{\circ},\,\,-8.5^{\circ})$ that is comparable to the same window in AU23, where $\Delta$[Fe/H]=[Fe/H]${}_{\rm{near}}-$[Fe/H]${}_{\rm{far}}=0.157$ dex. Our results confirm that the buckled bar can naturally explain the metallicity difference between the two RC peaks. A multi-population classical bulge component, i.e., the MCB scenario, is not necessarily needed to explain the observations.

Most importantly, we find a gradient of the metallicity difference along the longitudinal direction at $|b|\leq 6^{\circ}$ in the simulation shown in Figs. 3 and 5. Such a horizontal gradient is mainly caused by two factors. The first one is that the metallicity is higher close to the major axis of the bar, and lower around the minor axis. The tilted bar angle is the second factor. However, such a phenomenon is not expected in the MCB scenario; a prominent spheroidal classical bulge tends to produce similar $\Delta$[Fe/H] values between the double RCs at the opposite $\pm l$ regions. Particularly at $|b|=5^{\circ}$ region, we find a monotonically increasing pattern of $\Delta$[Fe/H] with $l$. This trend should not be expected in the MCB scenario. In addition, along the latitudinal direction, we note that the vertical metallicity gradient in brighter peaks is shallower than fainter ones. In the MCB scenario, the bright and faint peaks should have no significant difference in the vertical metallicity gradient.

The metallicity difference map between the near and far sides in Fig. 3 shows similar pattern to the profiles in Fig. 5. We apply linear fitting to the result in Fig. 3 to get a more quantitative estimation of the longitudinal (horizontal) and latitudinal (vertical) gradients of $\Delta$[Fe/H] with error estimated by Monte Carlo simulation. The results are listed in Table 1. We can see that both directions generally show positive gradients except for the latitudinal gradients at positive $l$. The longitudinal gradient is highest at smaller latitude values ($|b|=3^{\circ}$) and decreases towards higher latitudes. The latitudinal gradients are much higher at the negative longitudinal side than the positive longitudinal side with slightly negative gradients.

Our results show that a self-consistent buckled bar in AU23 can account for the metallicity difference between the bright and faint RCs in the Galactic bulge with a similar value to the observation at the bulge field $(l,\,\,b)=(-1^{\circ},-8.5^{\circ})$. Compared to the MCB scenario, a boxy/peanut bulge can naturally explain the structure, kinematics, and metallicity properties of the bulge stars. On the other hand, the MCB scenario has to include an unbuckled but very thick bar (extending to $|b|\sim 8^{\circ}$) to produce the Galactic longitude dependence of RCs at high $|b|$ region (Joo et al., 2017). Such a bar configuration is actually quite unlikely to be found in the MW or other disk galaxies. Future observations of the longitudinal/latitudinal metallicity difference gradients will be the smoking gun evidence to discriminate the MCB scenario.

The bar tilting angle ($\alpha$) plays an important role in generating the metallicity difference gradient along the longtitudinal directions as shown in Figs. 3 and 4. Different bar tilting angle should result in different $\Delta$[Fe/H] pattern. From the simulation, we choose different tilting angles of the bar to produce the metallicity difference maps in the bulge region. The longitudinal gradient of $\Delta$[Fe/H] could be used to constrain the bar tilting angle.

Fig. 6 shows the $\Delta$[Fe/H] maps for different bar tilting angle. In the $0^{\circ}$ case, $\Delta$[Fe/H] values are generally positive across the bulge region at $|b|>2^{\circ}$. There is no horizontal gradient in $\Delta$[Fe/H]. For intermediate bar angles, there is a negative $\Delta$[Fe/H] region close to the mid-plane and a horizontal $\Delta$[Fe/H] gradient at $3^{\circ}<|b|<9^{\circ}$. The gradient $d\Delta$[Fe/H]/$d\,l$ seems sensitive to the bar tilting angle. Fig. 7 shows the gradient $d\Delta$[Fe/H]/$d\,l$ at different $|b|$ slices for different $\alpha$. The gradient curves at $|b|<6^{\circ}$ show a sinusoidal pattern, with the steepest gradient at $\alpha\sim 45^{\circ}$. For higher $|b|$ slices ($|b|\geq 7^{\circ}$), the gradient is much shallower and varies mildly at different bar tilting angles.

With the gradient curves in Fig. 7, future observations of the bulge stars could help to constrain the bar tilting angle. It is suggested to measure the gradients in multiple latitudinal slices lower than $6^{\circ}$ since the other curves at higher latitudes in Fig. 7 are relatively flat with low amplitudes ($<0.005~{}{\rm dex/deg}$). If the gradient in only one latitudinal slice is used, even a small error will result in a large uncertainty of the bar tilting angle. As listed in Table 1, the gradient at $|b|=5^{\circ}$ is 0.012 dex/deg with the error of 0.002, corresponding to $\alpha=15^{\circ}\sim 45^{\circ}$. Therefore, to better constrain the bar angle, the error of the gradients should be small, and measurements at multiple latitudinal slices are also needed.