SUPERNOVAE AND THEIR EXPANDING BLAST WAVES DURING THE EARLY EVOLUTION OF GALACTIC GLOBULAR CLUSTERS

In the past decade, huge observational evidence has deeply challenged the common paradigm of Galactic globular clusters (GGCs) as the prototype of single stellar populations. High-resolution spectroscopy (see Gratton et al. 2012 for a review) and high-accuracy photometry (see Piotto et al. 2012 and references therein) have actually shown that GGCs host samples of stars with distinct light-element chemical patterns, populating distinct evolutionary sequences evident in UV–optical color–magnitude diagrams (Sbordone et al. 2011; Piotto et al. 2015). In the bulk of GGCs, the distinct sub-populations in any given cluster—characterized by their peculiar chemical patterns—do not show any differences in their overall metallicity as commonly traced by the Fe abundance, and thus present no evidence of contamination by supernovae (SNe). The implication is that the homogeneously low [Fe/H] value in all their stars is probably original to their proto-cluster primordial cloud. Note, however, that given the mass of GGCs, ranging from a few times 10⁵ M_⊙ to several times 10⁶ M_⊙ and any reasonable initial mass function (IMF), the number of SNe II expected during their early 40 Myr of evolution, amounts from a few thousand to a few tens of thousands of events, and, even so, none of their products appear to have been trapped by the sub-populations. A large intrinsic difference in iron abundance ($[\mathrm{Fe}/{\rm{H}}]\geqslant 0.1\;\sim$ dex), where SNe have played a major role in the chemical enrichment history of the cluster, has been considered for a long time as a peculiarity of the most massive GGC ω Cen (see Marino et al. 2011 and references therein). This is now changing given the evidence of other massive GGCs with internal variations in their metallicity. Such is the case of M54 (Carretta et al. 2010), Terzan 5 (Ferraro et al. 2009), M2 (Yong et al. 2014), NGC 5824 (Da Costa et al. 2014), M22 (Marino et al. 2011), and NGC 5286 (Marino et al. 2015). A scenario(s) able to explain the events leading to star clusters hosting multiple stellar populations (with their intriguing photometric and spectroscopic peculiarities) is still a matter of debate (see Cassisi & Salaris 2013; Renzini 2013 and references therein). For instance, to understand why only some GCs were able to retain the ejecta of the SNe associated with a first stellar generation (1SG) and thus increase the metallicity of the following second (and in some cases also the third) generation is obviously of pivotal importance. At the same time, as discussed by Renzini (2013), the observed metallicity spread in clusters with "multiple [Fe/H] abundances" poses a stringent constraint on the efficiency with which these GCs were able to retained the SN ejecta: the fraction of the ejecta expelled by the whole population of SNe and retained by the most massive GCs is only of the order of a few percent.⁴

Here, we focus on issues that could provide important clues about the properties of star clusters during their early stages. The aim is to explain why only some GGCs retain a small fraction of the SN ejecta, while the bulk of GGCs did not retain them at all. For this, we explore the hydrodynamics of SN blast waves as they sweep the surrounding medium and look for the conditions that may inhibit their usually expected dispersing nature.

An explosion in an inhomogeneous atmosphere could lead to the acceleration of the leading shock into the density gradient. Numerical simulations (see Tenorio-Tagle & Bodenheimer 1988 and references therein) have considered an exponential or a Gaussian galactic HI disk surrounded by a large gaseous halo and a massive young stellar cluster at or near the galaxy plane, driving a supersonic wind. The latter first generates a strong shock into the surrounding medium that steady decelerates, in all directions, as it continuously sweeps more of the surrounding gas into an expanding shell. However, as shown by Koo & McKee (1992), if the shock and its shell reach one scale-height (a couple of hundred parsecs) moving still supersonically, then the shock accelerates. This is because as the violent expansion proceeds, the pressure driving the shock indeed decays, but if the density in the direction away from the galaxy plane falls even more rapidly, then the section of the shock moving in that direction, instead of decelerating, as in explosions in a constant density medium, it would accelerate (${V}_{\mathrm{shock}}\propto {(P/\rho )}^{0.5}$), and the more so, the smaller that ρ becomes as the expansion proceeds. The sudden acceleration initiates blowout: shock acceleration sets a Rayleigh–Taylor instability in the leading section of the shell and this soon causes its fragmentation. The multiple piercings on the shell allow for the venting, between shell fragments, of the hot ejecta out of the bubble. Once this happens, the shell stalls. Only a small fraction of the matter swept prior to blowout is accelerated away from the galaxy, while most of it (∼95%) remains in the perturbed fragmented shell at the height acquired prior to blowout (Mac Low et al. 1989).

If the efficiency of star formation in a first-generation GC has been about 50%, then all the stars including potential core-collapse SNe will evolve buried by the large amount of gas left over from star formation. Here, we assume that this conforms a centrally condensed cloud that spans across the cluster radius R_SC and presents a Gaussian density distribution: ${\rho }_{\mathrm{cl}}={\rho }_{0}\mathrm{exp}[{-(r/{R}_{{\rm{c}}})}^{2}],$ where ρ₀ and R_c are the central gas density and the cluster core radius and the whole distribution is surrounded by a low-density ambient medium with ${\rho }_{\mathrm{ISM}}\ll {\rho }_{0}.$ The mass of the cloud then is

$$\begin{eqnarray}{M}_{\mathrm{cl}} & = & 4\pi {\rho }_{0}{\displaystyle \int }_{0}^{{R}_{\mathrm{SC}}}\mathrm{exp}[{-(r/{R}_{{\rm{c}}})}^{2}]{r}^{2}{dr}\\ & = & 4\pi {\rho }_{0}{R}_{{\rm{c}}}^{3}\left[{\pi }^{1/2}\mathrm{erf}({R}_{\mathrm{SC}}/{R}_{{\rm{c}}})/4\right.\\ & & -\left.{R}_{\mathrm{SC}}/(2{R}_{{\rm{c}}})\mathrm{exp}[{-({R}_{\mathrm{SC}}/{R}_{{\rm{c}}})}^{2}]\right],\end{eqnarray} \tag{ 1 }$$

where $\mathrm{erf}({R}_{\mathrm{SC}}/{R}_{{\rm{c}}})$ is the error function, e.g., M_cl ≈ 2 × 10⁵ M_⊙, if one assumes, as in our reference model, the central gas number density n₀ = 10⁶ cm⁻³ and R_c = 1 pc. We also assume that most of the massive stars are binaries (as in de Mink et al. 2009). This is expected for a large stellar population restricted to a very small volume. Contrary to single stars, massive binaries deposit their H burning products with very low velocities, favoring mixing with the remaining gas without perturbing its overall density distribution. Thus, we envisage a strongly concentrated gas density distribution hardly affected by the Roche lobes around massive binaries. This allows us to look for the effects produced by a single SN. We ignored the short phases that follow mass transfer and that likely lead to excavated bubbles around exploding stars. Such structures, given the large cloud densities, would be small. However, they will extend the free-expansion phase of the SN ejecta, and upon the ejecta swept-up shell interaction, a sudden thermalization of the released kinetic energy would take place. Finally, and depending on the mass in the surrounding shell, the Sedov phase may be completely avoided, making the strongly radiative SNe enter its snowplow phase to then progress supersonically into our assumed cloud density distribution (see Tenorio-Tagle et al. 1991) similar to the initial condition in our calculations. Our approach differs from that of Krause et al. (2012, 2013), who consider a continuous central wind powered by single massive stars and sequential SNe and also by the energy delivered by accretion into black holes and neutron stars. Such an injection of energy led to the build up of superbubbles that in all their cases experienced blowout and thus led to no contamination of the gas left over from a 1SG, leaving without an explanation the Fe spread in the most massive GCs.

Here, we show that single SN blast waves evolving in a medium with a steep density gradient also naturally experience blowout, i.e., the sudden and continuous acceleration of the leading shock in the direction of lowest densities and the fragmentation of the leading section of the shell. The high pressure gas, the thermalized ejecta, rapidly escapes between shell fragments to closely follow the shock into the density gradient as this becomes more elongated. Given the size of the cluster (R_SC), the ejecta are then expelled into the surrounding medium soon after blowout. This leads to a rapid loss of energy and pressure from the volume swept by the SN shock, which strongly inhibits the lateral growth of the cavity. In fact, soon after blowout, the now larger pressure of the dense cloud matter surrounding the cavity would favor its collapse toward the elongated volume symmetry axis, leaving then no trace of the SN explosion. Explosions occurring within the cluster core, however, would require that their leading shocks reach the cluster core boundary and enter the density gradient while still progressing with a supersonic velocity for such an evolution to take place. In both cases, the SN products would be expelled from the cluster.

Only SN blast waves well contained near the cluster center while evolving into the dense background medium available for the formation of a 2SG would contribute to eventually enhancing its Fe abundance. This occurs if the SN blast wave does not supersonically reach the cluster core radius. In a dense proto-cluster cloud environment, the swept-up gas cools down and collapses into a cold, dense, and narrow shell very rapidly, making an early transition from the quasi-adiabatic Sedov to the snowplow evolution at ${t}_{0}=3{{kT}}_{{\rm{s}}}/(4{n}_{\mathrm{off}}{\rm{\Lambda }}({T}_{{\rm{s}}},Z)),$ where k is the Botzmann constant, ${\rm{\Lambda }}({T}_{{\rm{s}}},Z)$ is the cooling function, Z is the shocked gas metallicity, n_off is the gas number density at the explosion site, and T_s is the post-shock temperature: ${T}_{{\rm{s}}}=2{(\gamma -1)\eta {V}_{{\rm{s}}}^{2}/((\gamma +1)}^{2}k).$ Here, η = 14/23 m_H is the mean mass per particle in the shocked ionized plasma with 10 hydrogen atoms per helium atom, m_H is the proton mass, γ = 5/3 is the ratio of specific heats, and ${V}_{{\rm{s}}}=0.4{(\xi {E}_{0}/{\rho }_{{\rm{c}}})}^{1/5}{t}_{0}^{-3/5}$ is the velocity of the SN blast wave. Inside the cloud core radius R_c, where the density distribution is almost homogeneous, the velocity of the swept-up shell is (e.g., Pasko & Silich 1986)

$$\begin{eqnarray}&&{V}_{{\rm{s}}}=\displaystyle \frac{2}{7}{\left(\displaystyle \frac{147{E}_{{\rm{T}}0}{R}_{0}^{2}}{4\pi {\rho }_{{\rm{c}}}}\right)}^{1/7}{t}^{-5/7},\end{eqnarray} \tag{ 2 }$$

where ${R}_{0}={(\xi {E}_{0}/{\rho }_{{\rm{c}}})}^{1/5}{t}_{0}^{2/5}$ and ${E}_{{\rm{T}}0}=(\gamma +1){E}_{0}/(3\gamma -1)$ are the radius of the shock and the thermal energy driving the expansion at the end of the Sedov phase, E₀ = 10⁵¹ erg is the energy of the explosion, ρ_c is the central gas density in the proto-cluster cloud, and $\xi =75{(\gamma -1)(\gamma +1)}^{2}/(16\pi (3\gamma -1))$ (e.g., Bisnovatyi-Kogan & Silich 1995). The velocity drops to the sound speed value a_s, and the shell stalls at

$$\begin{eqnarray}&&{R}_{\mathrm{stall}}={\left(\displaystyle \frac{2}{7{a}_{{\rm{s}}}}\right)}^{2/5}{\left(\displaystyle \frac{147{E}_{{\rm{T}}0}{R}_{0}^{2}}{4\pi {\rho }_{{\rm{c}}}}\right)}^{1/5}.\end{eqnarray} \tag{ 3 }$$

The stalling radius R_stall is smaller than the core radius R_c if the central density exceeds the critical value:

$$\begin{eqnarray}&&{\rho }_{\mathrm{crit}}={\left(\displaystyle \frac{2}{7{a}_{{\rm{s}}}}\right)}^{10/7}{\left(\displaystyle \frac{147{E}_{{\rm{T}}0}}{4\pi }\right)}^{5/7}{\left(\displaystyle \frac{{\xi }^{2}{E}_{0}^{2}{t}_{0}^{4}}{{R}_{{\rm{c}}}^{25}}\right)}^{1/7}.\end{eqnarray} \tag{ 4 }$$

Equation (4) allows one to calculate the critical mass of the cloud M_crit. Only globular clusters with a gas mass available for a 2SG M > M_crit would be able to retain SN products and eventually enhance its Fe abundance. This may help to explain the observed Fe spread at the high-mass end of GGCs and the lack of it for the bulk of GCs. For a gas temperature T_cloud = 100 K and a cloud core radius between 0.5 pc ≤ R_c ≤ 2 pc, Equation (4) leads to a critical mass (2–5) × 10⁵ M_⊙.

To calculate the SN shock evolution for off-centered explosions, we used our 2D thin layer approximation code (Bisnovatyi-Kogan & Silich 1995). In all cases, the transition to the snowplow phase occurs when the shock radius is still small (0.045 pc $\leqslant \;{R}_{\mathrm{cool}}\leqslant 0.3$ pc), and thus for the calculations we adopted as initial conditions the shock radius and expansion velocity at the end of the Sedov evolution evaluated for the gas density ρ_off at the SN site.

Figure 1 shows our results for a proto-cluster cloud with a central gas density n₀ = 10⁶ cm⁻³ and R_c = 1 pc, with an SN explosion first at the center of the gas distribution and then at 1 pc and at 2 pc away from the center. The figure tracks the calculated shape of the SN blast wave at three different times (solid, dotted, and dashed lines), as well as the spherically symmetric Gaussian density distribution of the gas cloud, shown by thin concentric lines labeled with the logarithmic values of the local density. Also shown are the top pole shock velocity as a function of distance to the GC center. The latter plots display the strong shock acceleration and thus blowout as the shock enters the steep density distribution, whereas in the first case, the shock velocity goes to 0 km s⁻¹ (see Figure 1(b)) inside R_c, and thus the SN products are trapped to finally contaminate the remaining cloud. Other calculations assuming the explosion site to be ≤0.25 pc have lead to the same result: well-contained SN explosions, while explosions out of this radius sooner or later experience blowout. Comparing this volume with that of the cloud, the proportion of expected well-contained SNe is only 0.024%. More massive clouds with a larger trapping radius will lead to a larger proportion of SNe contaminating the cloud, as expected for the most massive GGCs.

Zoom In Zoom Out Reset image size

Note that as the SN blast wave sweeps the overtaken matter into a shell, it also unveils a large fraction of massive stars from the 1SG. These are able to photoionize the inner edge of the surrounding shell, causing a mayor increase on its pressure. Its temperature would be T_s ∼ 10⁴ K and its density is about 10⁶ cm⁻³, while in the cavity interior the gas would present a temperature around 10⁶ K and a density of the order of 10⁻² cm⁻³ and thus a strong inward champagne flow will take place (Tenorio-Tagle 1979; Franco et al. 1990). The collapse velocity of the cavity is ${v}_{\mathrm{collapse}}={c}_{{\rm{s}}}{({P}_{\mathrm{shell}}/{P}_{\mathrm{cavity}})}^{0.5},$ where c_s is the sound speed in the photoionized shell (∼10 km s⁻¹), which leads to v_collapse values larger than 10³ km s⁻¹. The collapse would be even faster for cavities left after blowout by off-centered explosions, as in all of these the exit of the hot gas will leave a negligible pressure in the cavity allowing the photoionized shell to rapidly restore the original density. The collapse of the cavity is shorter than the spacing between SNe, which for a 10⁶ M_⊙ cluster with a standard IMF occurs about every 3000 years. This implies that sequential SNe in clusters with a stellar mass ≤10⁶ M_⊙ will encounter a fully restored gas density distribution.

We have followed the consensus indicating that the peculiar light-element anti-correlations, the He enhancement ubiquitous to all GGCs, and particularly the Fe spread found in the most massive GC are signatures of 2SGs resultant from matter hugely polluted by the ejecta of a 1SG (see Renzini 2013 for a critical review on the various suggested 1SG polluters).

If one assumes a massive 1SG (M ∼ 10⁵–10⁶ M_⊙) and a Salpeter IMF, one expects during the first 40 Myr of evolution a core-collapse SN every ∼10⁴ to a few 10³ years, and their impact has been expected to be devastating, disrupting the mass left over from star formation. Here, however, we have shown that SNe exploding in a steep density gradient and those taking place within the densest central regions, but still being able to reach the cluster core radius (R_c) with supersonic velocities, lead to blast waves able to undergo blowout. The volume affected by the explosions becomes highly elongated as they progress into the gradient to reach the cluster edge within a few thousand to a few tens of thousand of years and then vent their energy and their metals into the medium surrounding the cluster. Such events lead then to no contamination of the gas left over from the formation of a 1SG. The rapid loss of pressure after discharging their metals into the surroundings inhibits the dispersal of the leftover cloud leaving then no trace within the clusters of the SN events. Note that blowout will also take place if the shock waves from other explosions find cavities similar to those shown in Figure 1, releasing their energy and metals within a thousand years out of the cluster volume. We have shown that only well-centered SN blast waves evolving into a dense medium would attain subsonic velocities before reaching the cluster core radius, and thus are likely to be retained by the cluster. Thus, if a 2SG then forms, this should present an Fe spread as observed in the most massive GGC. We have also shown that the cavities generated by SNe are likely to be refilled either after blowout or after the shell stalls, due to the large UV flux from all unveiled stars, which through photoionization provide the inner skin of the shell or cavity with a large pressure and cause the rapid restoration of the original density distribution. Thus, sequential SNe from massive clusters are likely to explode in very similar environments, and the majority of them will experience blowout.

The authors thank the referee for valuable comments. This study was supported by CONACYT—México, grants 167169 and 131913, and by the Spanish Ministry of Science and Innovation for the ESTALLIDOS collaboration (grants AYA2010-21887-C04-04 estallidos4 and AYA2013-47742-C4-2-P estallidos5). G.T.T. also acknowledges the Cátedra Severo Ochoa at the Instituto de Astrofi´sica de Canarias and the CONACYT grant 232876 for a Sabbatical leave. S.C. acknowledges financial support from PRIN-INAF 2014 (PI: S. Cassisi) and PRIN MIUR 2010-2011 (prot. 2010LY5N2T) and the friendly hospitality at the IAC. The authors appreciate the science and discussions among the participants of the ESTALLIDOS Star Formation Feedback Workshop (IAC, 2014 November) that triggered a good number of ideas.