A&A 373, 1043-1055 (2001)
DOI: 10.1051/0004-6361:20010673
J. M. Pittard - T. W. Hartquist - J. E. Dyson
Department of Physics & Astronomy, The University of Leeds, Woodhouse Lane, Leeds, LS2 9JT, UK
Received 6 March 2001 / Accepted 2 May 2001
Abstract
We present similarity solutions for adiabatic bubbles that are blown by
winds having time independent mechanical luminosities and that are each
mass-loaded by the hydrodynamic ablation of
distributed clumps. The mass loading is "switched-on'' at a specified
radius (with free-expansion of the wind interior to this point) and
injects mass
at a rate per unit volume proportional to
where
(1) if the flow is subsonic (supersonic) with respect to
the clumps.
In the limit of negligible mass loading a similarity solution found by
Dyson (1973) for expansion into a smooth ambient medium is
recovered. The presence of mass loading heats the flow, which
leads to a reduction in the Mach number of the supersonic
freely-expanding flow, and weaker jump conditions across
the inner shock. In solutions with large mass loading, it is possible for
the wind to connect directly to the contact discontinuity without first
passing through an inner shock, in agreement with previous hydrodynamic
simulations. In such circumstances, the flow may or
may not remain continuously supersonic with respect to the clumps.
For a solution that gives the mass of swept-up ambient gas to be less than
the sum of the masses of the wind and ablated material,
,
meaning that the exponent of the density profile
of the interclump medium must be at most slightly positive, with negative
values preferred. Maximum possible values for the ratio of ablated mass
to wind mass occur when mass loading
starts very close to the bubble center and when the flow is supersonic
with respect to the clumps over the entire bubble radius.
Whilst mass loading always reduces the temperature of the shocked wind,
it also tends to reduce the emissivity in the interior of the
bubble relative to its limb, whilst simultaneously increasing the central
temperature relative to the limb temperature. The maximum temperature in
the bubble often occurs near the onset of mass loading, and in some cases
can be many times greater than the post-inner-shock temperature.
Our solutions are potentially relevant to a wide range of astrophysical
objects, including stellar wind-blown bubbles, galactic winds, starburst
galaxy superwinds, and the impact of an AGN wind on its surrounding
environment. This work complements the earlier work of Pittard et al.
(2001) in which it was assumed that clumps were evaporated
through conductive energy transport.
Key words: hydrodynamics - shock waves - stars: mass-loss - ISM: bubbles - galaxies: active
Investigations into the self-similar nature of supernova remnant (McKee & Ostriker 1977; Chièze & Lazareff 1981; Dyson & Hartquist 1987) and stellar wind-blown bubble (Pittard et al.2001) evolution in tenuous media with embedded clumps have been performed. Detailed one-dimensional, time-dependent hydrodynamic models of specific WBBs associated with evolved stars have also been constructed (Arthur et al. 1993, 1994, 1996). In at least some isothermal mass-loaded winds, Arthur et al.(1994) and Williams et al. (1995, 1999) showed that a global shock does not occur in the region where the Mach number is large. Instead, a global shock occurs only after the wind has travelled far enough for mass loading to decelerate it sufficiently that its Mach number is less than a few. In some of these cases, a global shock within the mass loading region does not exist at all, and the wind continues until it encounters a termination shock caused by the interaction of the wind with an external medium.
Observations supporting mass loading in WBBs include the detection of blue-shifted absorption features of species with a range of ionization potentials in the spectrum of the central star of the Wolf-Rayet nebula RCW 58. Smith et al.(1984) argued that the observed velocity spread for the detected features is much larger than would be expected for K gas consisting only of stellar wind material decelerated through a terminal shock, and suggested that the velocity spread originates due to mass loading of the shocked wind through the entrainment of material from clumps. Crucially, one can obtain an estimate for the ratio of wind mass to entrained clump mass: the observed velocity spread implies a value of 40 or 50 to one. A detailed time-dependent, hydrodynamic model including nonequilibrium ionization structure is consistent with these conclusions (Arthur et al.1996). Spectroscopic data also provide evidence for transonic flows in the halo of core-halo planetary nebula (Meaburn et al.1991), again consistent with high mass loading rates. Finally, it is also suspected that the gradual deceleration of a wind by mass loading and the associated weakening of an inner shock may in fact contribute to the radio quietness of some WBBs, although this is a conjecture which currently cannot be proven (see Williams et al.1995, 1999). In Pittard et al.(2001, hereafter PDH), similarity solutions were derived for the structures and evolution of mass-loaded WBBs, under the assumption that conductive evaporation from embedded clumps was the dominant mass loading process, and that both evaporation from the cold swept-up shell and radiative losses were negligible. To obtain a similarity solution with these conditions, specific radial power-laws on the clump and interclump density distribution, and temporal power-laws on the wind mass-loss rate and terminal velocity, were required. Approximate similarity solutions for evaporatively mass-loaded WBBs with the assumptions of constant mass-loss rate and wind velocity, and an isobaric shocked wind region have previously been obtained by Weaver et al.(1977; evaporation from the cool swept-up shell) and Hanami & Sakashita (1987; mass-loading from clumps). A central assumption in both of these papers was that the shocked wind was approximately isobaric. However, PDH demonstrated that this was not necessarily a good assumption (e.g. see their Fig. 4). Indeed, imposing this condition is likely to set a limit to the amount of mass loading that can occur.
A central conclusion of PDH was that there exist maximum possible values for the ratio of evaporated mass to stellar wind mass, as a consequence of the evaporation rates dependence on temperature and the lowering of temperature by mass loading. In particular it was difficult to find ratios approaching what was observationally required. The work in this paper complements PDH by considering the case in which hydrodynamic ablation is the dominant mass addition process. As conductively driven evaporation has a very temperature sensitive rate, ablation is likely to regulate clump dispersal into lower temperature media.
Our solutions are again potentially relevant to such diverse objects as WBBs created by a faster wind interacting with a clumpy AGB superwind, by the wind of a young star interacting with surrounding molecular material, and the wind of an active galactic nucleus impacting its environment.
In this paper we present similarity solutions for WBBs mass-loaded by
hydrodynamic ablation. Our solution method is very similar to that used by
PDH, and we refer the reader to that paper for a discussion of
many of the assumptions involved. As a starting point, we also assume the
same qualitative shock structure in the WBB as shown in Fig. 1 of PDH:
that is a central wind source surrounded by a region of
unshocked wind, separated from an outer region of shocked
wind by an inner shock. The swept-up ambient medium is
assumed to radiate efficiently and collapse to a negligibly thin shell
coincident with the contact discontinuity separating the stellar and ambient
gas (cf. Dyson & de Vries 1972).
Figure 1: Results for , , and . In the limit of large ( i.e. negligible mass loading), the earlier results of Dyson (1973) are recovered. From the top the panels indicate density, velocity, internal energy, and Mach number (with respect to the stationary clumps) as a function of radius. Except for the Mach number, the panels are normalized to values of 1.0 at the contact discontinuity. | |
Open with DEXTER |
A central and fundamental difference between the conductive case considered in PDH and the ablative case examined here, however, is that here mass loading may also occur in the unshocked stellar wind. The effect is to heat this part of the flow, which leads to reduced Mach numbers, and weaker jump conditions across the inner shock. We shall see that for high mass loading rates, it is possible for the unshocked mass-loaded wind to connect directly with the contact discontinuity, without the presence of an inner shock.
In this work we assume that the flow can be treated as a single fluid. This requires that the ablated clump material merges with the global flow in the sense that its temperature, velocity, and density approach those of the surrounding tenuous material. Ablation by itself might require that the clump material mix microscopically with the wind to reach the density and temperature of the wind. However, though the mass-loss rate may be controlled by ablation, thermal conductivity will almost certainly be responsible for material, once it is stripped from a clump, reaching the physical state of the tenuous fast flowing medium. Thermal conduction accomplishes this phase transition without microscopic mixing, and acceleration to the global flow speed is effected by the response of stripped material to pressure gradients and viscous coupling, which may arise from a host of mechanisms including turbulence. Thus we envisage a two-stage process in which ablation controls the rate at which mass is stripped from the clump but conductivity becomes important for the merging into the global flow.
This receives support from the comparison of the conductively driven evaporation time of a spherical clump with and T=10^{4} K, embedded in a medium with the same pressure (which is typical of planetary nebulae, Wolf-Rayet wind-blown bubbles, and starburst superwinds) and T=10^{7} K (which is also typical of hot material in such regions) to the sound crossing time in the clump. The sound crossing time is somewhat (but not hugely) smaller than the conductively driven evaporation time for a large range of clump sizes, if the Cowie & McKee (1977) estimate of the mass-loss rate driven by saturated conduction is used. Consequently, ablation initiates mass-loss because it causes an increase in the surface area of a clump triggering more conductive heat transfer. Use of the analysis of Cowie & McKee (1977) and McKee & Cowie (1977) shows that for the assumptions given above, radiative losses do not hinder conductively driven evaporation unless the clump radius is greater than pc where is the number density ( ) of the global tenuous flow. Of course, a clump that does not cool radiatively after it is compressed by a shock will have a shorter sound crossing time and be ablated more rapidly. However, by assuming a clump temperature of 10^{4} K above we have established that the two-stage process is likely to be relevant in many environments if a clump radiatively cools after being shocked.
Observations of WBBs and PNe indicate that mass loading may not necessarily begin at zero distance from the wind source. For instance, in the young nebula Abell 30 the clumps appear almost all of the way down to the central star (Borkowski et al.1995), whilst in the more evolved Helix Nebula the clumps are all located closer to its edge (e.g. Meaburn et al.1996). We therefore include the radius at which mass-loading "switches on'' as a free parameter in our models, and assume that the wind undergoes free expansion interior to this radius. Hence each solution will consist of a region of supersonic wind with no mass loading and an adjacent region with mass loading.
One can imagine two possible causes for this minimum "mass-loading'' radius. In one scenario the clumps could have been ejected at low velocity from the central star at an earlier evolutionary stage. The ejection of clumps then abruptly stopped, so that at the time of observation they had travelled a finite distance from the central star. By this process a central region clear of clumps surrounded by a clumpy region can be generated. A second possibility is that clumps interior to the mass-loading radius have been completely destroyed by the wind. It seems reasonable to suppose that clumps located closest to the central star will be destroyed first, since they will have been subjected to the wind from the central star for the longest time. Then as the bubble or nebula evolves, clumps at ever increasing distance from the central star will be destroyed. Timescales for the destruction of clumps by ablation can be estimated from Hartquist et al.(1986) and Klein et al.(1994). Estimated destruction timescales vary from significantly less than to greater than the age of the bubble/PNe, in accord with the different spatial distribution of clumps in objects of differing age.
Regardless of which of the above scenarios is responsible for the existence of such a minimum mass-loading radius, this radius will physically increase with time. Our similarity solution requires that it increases in the same way as that of the contact discontinuity (i.e. , where is the radial dependence of the mass-loading). For most of the solutions presented in this paper, the minimum mass-loading radius scales with or close to t. Since on physical grounds we might expect it to scale as t, our solutions closely match this requirement.
In our solutions an inner shock may or may not be present - in the latter case the mass-loaded wind directly connects to the contact discontinuity, and the mass loading may be strong enough for the wind to become subsonic with respect to the clumps before the contact discontinuity is reached. If an inner shock is present, the postshock flow is by definition subsonic with respect to the shock, but may still be supersonic with respect to the clumps. In this case a number of different profiles for the Mach number are possible before the contact discontinuity is reached.
At the center of the bubble prior to the onset of mass loading,
we solve only the continuity and momentum equations, with the implicit
assumption that the thermal energy of the flow is negligible, whilst in
the mass loading regions we additionally solve the energy equation, and
include a source term for mass injection in the continuity equation. For
a
gas, the equations for the mass-loaded flow are:
(1) |
(2) |
(3) |
(4) |
(5) |
(6) |
(7) |
(13) |
(14) |
(15) |
(16) |
At the CD, we must also satisfy conservation of momentum. We take
(18) |
(20) |
M_{1} | |||||||||
0.1 | 0.271 | 0.46 | 0.014 | 0.469 | 0.059 | 102 | 1.454 | 376 | 4.79 |
0.2 | 0.268 | 0.46 | 0.013 | 0.469 | 0.059 | 101 | 1.371 | 390 | 4.79 |
0.33 | 0.266 | 0.46 | 0.013 | 0.469 | 0.059 | 99.7 | 1.289 | 408 | 5.00 |
0.5 | 0.268 | 0.46 | 0.013 | 0.469 | 0.059 | 99.2 | 1.223 | 426 | 5.56 |
0.9 | 0.264 | 0.46 | 0.013 | 0.469 | 0.059 | 99.0 | 1.152 | 450 | 11.9 |
0.02 | 0.82 | 0.17 | 0.63 | 0.18 | 0.022 | 2.95 | 6.50 | 0.51 | 1.60 |
0.1 | 0.76 | 0.22 | 0.53 | 0.22 | 0.028 | 2.95 | 4.48 | 0.80 | 1.45 |
0.2 | 0.73 | 0.24 | 0.49 | 0.24 | 0.030 | 2.90 | 3.84 | 0.92 | 1.34 |
0.5 | 0.66 | 0.25 | 0.45 | 0.27 | 0.033 | 2.67 | 2.55 | 1.23 | 1.12 |
0.9 | 0.67 | 0.24 | 0.39 | 0.33 | 0.041 | 2.85 | 1.75 | 2.17 | 2.53 |
(22) |
(23) |
(24) |
(26) |
The ratio of mass pick-up from the clumps to the mass swept-up from the interclump medium (which is related to the value of ) is obtained only after a particular solution has been found. For bubbles whose evolution is significantly altered by mass loading, we require that simultaneously and are both small.
The similarity equations were integrated with a fifth-order accurate adaptive step-size Bulirsch-Stoer method using polynomial extrapolation to infinitesimal step size. Once the CD was reached, the similarity variables were rescaled using the relationships defined in Eq. (21) so that . The mass, and kinetic and thermal energies of the bubble were calculated, as were the kinetic energy of the shell and the energy radiated from it. The correct normalization to satisfy global energy conservation was then obtained. Finally, for given values of and t, the similarity variables x, f, g and h may be scaled into the physical variables and .
We first checked our work against the solutions obtained by Dyson
(1973) for WBBs with no mass loading. For this comparison it was
required that
and that
was large. For
the radius of minimum mass loading increases as t.
Excellent agreement with Dyson (1973) was found
over a large range of
.
For large
(i.e. negligible
mass loading) the value of
has no effect on the resulting
solution. Figure 1 shows a sample solution.
Figure 2: Results for , , and various values of . For these parameters the mass loading is small but non-negligible, and hence demonstrates having a minor influence on the results. Results for five values of are shown: 0.1 (solid), 0.2 (dots), 0.33 (dashes), 0.5 (dot-dash), and 0.9 (dot-dot-dot-dash). The panels are the same as in Fig. 1. In Fig. 3 we show solutions where the value of has a much greater influence. | |
Open with DEXTER |
In Fig. 2 we show results for
,
,
with different values of
.
The mass loading in these results is
small, but not negligible, so that the precise value of
has
minor consequences for the
solution obtained. In Table 1 we tabulate important
parameters from these solutions. In particular, we find that values for
the ratio
and
the energy fractions are very similar. The fractional mass loading of
the bubble,
,
is somewhat more sensitive to the value of
,
and increases when the "turn-on'' radius for mass loading is decreased, as
one would expect. The preshock Mach number also reflects the degree of
preshock mass loading, again as expected. Interestingly, the profiles
of postshock Mach number are, however, very insensitive to the preshock
Mach number. With the assumption of constant ,
we note that since
varies for the above results it was important to see
whether this was a possible factor for some of the differences in these
solutions. Therefore we also investigated the effect of different values
of
whilst keeping
constant and varying
as
necessary. We find
that the resulting solutions are again fairly similar, and differ at about
the same level as the results in the top half of Table 1.
Therefore it seems that varying
has a direct effect on the
results, and not just an indirect influence through changes induced in
and/or .
In the following, therefore, we will vary
whilst keeping fixed.
Figure 3: Results for , , and various values of . For this value of , the mass loading is appreciable, and hence the solutions shown demonstrate having a major influence. Results for four values of are shown: 0.02 (solid), 0.1 (dots), 0.5 (dot-dash), and 0.9 (dot-dot-dot-dash). The panels again represent density, velocity, internal energy and Mach number. | |
Open with DEXTER |
At the other extreme, Fig. 3 shows results where the value
of
can have large consequences for the solutions obtained. This
occurs for values of
for which the resulting mass loading of the bubble
is important. Parameters from these solutions are tabulated in the lower
half of Table 1.
Although the ratio
remains relatively insensitive to the
value of
,
the fractional mass loading
,
the ratio of the
shell mass to the bubble mass
,
and the energy partition are
all significantly affected, as Table 1 shows. This is
expected given that most of the mass loading occurs near the minimum mass
loading radius for
.
Figure 4: Results for , . Here mass loading is so severe that the flow directly connects to the contact discontinuity. For the solution shown, the onset of mass loading occurs at . Once the flow starts mass loading, its velocity drops, its relative density increases, its kinetic energy drops, and its thermal energy increases. | |
Open with DEXTER |
An interesting consequence of mass loading the wind prior to the inner shock
is that if the mass loading is large, the Mach number of the flow can be
reduced so much that the flow directly connects to the contact discontinuity
without the presence of an inner shock. In Fig. 4
we show one solution where this happens. In this example,
the flow is continuously supersonic with respect to the clumps, and we
obtain
and
.
This perhaps indicates
that this morphology is most favoured at early stages of the evolution of
the bubble before a significant amount of ambient medium is swept-up. In
another example of "direct connection'', we find a solution where the Mach
number of the flow with respect to the clumps drops below unity (see
panel f of Fig. 6). In this case,
and
,
suggesting a much older bubble.
We also find that low values of
(
)
are apparently
needed to obtain directly-connected solutions.
Figure 5: Comparison of solutions with negligible (solid) and high (dots) mass loading rates, for a given value of (-2) and (0.606). The ratio of swept-up mass to bubble mass is 12.3 and 0.96 respectively, whilst the ratio of ablated mass to wind mass is 7.54 for the mass-loaded solution. The preshock Mach number for the mass-loaded solution is 1.006. In both cases, the similarity solution has been scaled such that the current wind parameters are appropriate for a WR star: and km s^{-1}. | |
Open with DEXTER |
In Fig. 5 we compare solution profiles for the
case where
and for negligible and high mass loading rates.
A
clump distribution is the most physically
plausible for many astrophysical objects. For this distribution
the minimum mass loading radius increases as t^{2/3}, which is not far
removed from a linear increase with t.
The solutions have been scaled from the similarity results by assuming
current values for the wind parameters appropriate for a WR star:
and
km s^{-1}. The mass-loss rate
increases with time as t^{2/3} whilst the wind velocity decreases as
t^{-1/3}, so the wind has evolved from a case more suitable to an O-star.
The assumed bubble age is
t = 10^{4} yr. The interclump medium varies
as
(i.e.
)
which is a flatter
distribution than for a constant velocity wind (r^{-2}). However,
since the wind mass-loss rate is increasing and the velocity decreasing,
an r^{-1/2} interclump density profile is physically realistic for
such a scenario.
In comparison with the negligibly mass-loaded solution, the increase in
density inside the bubble for the highly mass-loaded solution can be
clearly seen, along with a decrease in preshock velocity, preshock Mach
number, and postshock temperature. Interestingly, the Mach number of the
postshock flow is higher for the heavily mass-loaded solution.
Note also that the heavily mass-loaded solution extends to a larger bubble
radius. This is a consequence of the different preshock interclump densities:
for the heavily mass-loaded solution
,
whilst
for the negligibly mass-loaded solution.
Figure 6: The Mach number as a function of x for a number of different solutions. With respect to the clumps, the shocked region can be either entirely subsonic, entirely supersonic, or have one (or maybe more) sonic points. The panels show: a) shocked region with a monotonic decrease in Mach number ( , , , ); b) shocked region entirely supersonic with a monotonic increase in Mach number ( , , , ); c) postshock Mach number initially subsonic and decreasing, before levelling off and subsequently increasing through a sonic point ( , , , ); d) shocked region entirely subsonic with a decreasing Mach number followed by a modest increase close to the contact discontinuity ( , , , ); e) postshock flow continuously supersonic though with an initial decrease in Mach number ( , , , ); f) continuously supersonic flow (solid) directly connected to the contact discontinuity ( i.e. no inner shock) ( , , ) and a flow becoming subsonic (dots) directly connected to the contact discontinuity ( , , ). | |
Open with DEXTER |
In PDH it was shown that there is a
negative feed-back mechanism caused by the evaporation of mass from
embedded clumps, which set a maximum limit to the amount that a bubble
could be massloaded. Additionally, for the bubble mass to be greater than
the mass of the swept-up shell, a large radial dependence of the
mass loading was required (
). In this work we instead
find that small values of
are required to satisfy this condition.
Our simulations show that it is impossible to obtain a similarity solution
with
for
(e.g. for
,
). However, for smaller values of ,
solutions satisfying
can be found (e.g. for
(-3), we can obtain
(29)
,
where
(15)).
For a given value of
there appears to be a maximum
value for
.
It occurs when mass loading
starts very close to the bubble centre (i.e.
)
and
when the flow remains supersonic relative to the clumps over the entire
bubble radius.
Figure 7: Radial profiles of emissivity per unit volume and temperature, normalized to values of 1.0 at the limb. The upper panels show results where an inner shock is present, whilst the bottom panels show results where the mass loading region directly connects to the contact discontinuity. For the upper panels the results shown are for: , , , (solid); , , , (dots); , , , (dashes). For the lower panels the results are for: , , (solid); , , (dots). | |
Open with DEXTER |
We have also investigated the Mach number profile of the solutions that we obtain, of which some are already shown in Figs. 2 and 3. As for the evaporative mass loading in PDH, we note that a number of different forms are possible, and summarize these in Fig. 6.
Finally, radial profiles of the X-ray emissivity per unit volume and temperature have been calculated (Fig. 7). We assume that the emissivity , which is a good approximation over the temperature range (cf. Kahn 1976). The upper panels show results where an inner shock is present, and we will discuss these first. The solution with the solid line displays relatively constant values of and T, and is characterized by a high degree of mass loading ( ). Note that the highest temperature in the bubble does not occur postshock, but rather shortly after mass loading begins. The very high emissivity at the onset of mass loading is a consequence of the low temperature at this point and the assumed T^{-1/2} dependence of our emissivity, and is not physical. The solution with the dotted line has an interior dominated by free-expansion of the wind source. Mass loading only occurs over the final 20 per cent of the radius, and is relatively weak, and the maximum temperature in the bubble occurs immediately postshock. The solution with the dashed line has some characteristics of each of the other solutions. It is mass-loaded from small radius which again leads to the maximum bubble temperature occurring near the center. However, the mass loading is less severe ( ) and its radial dependence steeper ( ) than for the solution with the solid line, so that there is a large fall in density between the switch-on radius for mass loading and . This leads to a steep rise in emissivity from the inner shock to the bubble center in a similar fashion to the solution with the dotted line.
The bottom panels of Fig. 7 show results where there is no inner shock and where the mass loading region directly connects to the contact discontinuity. The Mach number profiles for these solutions are displayed in the bottom-right panel of Fig. 6. As can be seen, the temperature profiles are approximately flat over their respective mass loading regions, and are very similar in shape to each other. This contrasts with the behaviour of their emissivity profiles, where it is clear that the solution represented by the solid (dotted) line has generally decreasing (increasing) with radius.
A general result from Fig. 7, which was also found in the evaporatively mass-loaded bubbles of PDH, is that higher mass loading tends to decrease the central emissivity and increase the central temperature relative to the limb.
We first produced solutions with negligible mass loading and (which correspond to constant and ), which we compared with results obtained by Dyson (1973). Excellent agreement for the structure of the bubble was found. We also confirmed that for negligible mass loading, the value of had no effect on the resulting solutions, as desired.
We then investigated the changes in the structure of the bubble for different values of , , , and . The central conclusions are:
Acknowledgements
JMP would like to thank PPARC for the funding of a PDRA position, and Sam Falle and Rob Coker for helpful discussions. We would also like to thank an anonymous referee whose suggestions improved this paper. This work has made use of Nasa's Astrophysics Data System Abstract Service.