# The cometary cavity created by an aligned streaming environment/collimated outflow interaction

###### Abstract

We present a “thin shell” model of the interaction of a biconical outflow and a streaming environment (aligned with the direction of the flow), as well as numerical (axisymmetric) simulations of such an interaction. A similar situation, although in a more complex setup, takes place at the head of the cometary structure of Mira. Thus, for most of the numerical simulations we explore parameters consistent with the observed bipolar outflow from Mira B. For these parameters, the interaction is non-radiative, so that a rather broad jet/streaming environment interaction region is formed. In spite of this, a reasonable agreement between the thin-shell analytic model and the numerical simulations is obtained.

D. López-Cámara, A. Esquivel, A. C. Raga, P. F. Velázquez, A. Rodríguez-González: Instituto de Ciencias Nucleares, Universidad Nacional Autónoma de México, Ap. 70-543, 04510 D.F., México () J. Cantó: Instituto de Astronomía, Universidad Nacional Autónoma de México, Ap. 70-264, 04510 D.F., México \listofauthorsD. López-Cámara, A. Esquivel, J., Cantó, A. C. Raga, P. F. Velázquez, A. Rodríguez-González \indexauthorLópez-Cámara, D. \indexauthorEsquivel, A. \indexauthorCantó, J. \indexauthorRaga, A. C. \indexauthorVelázquez, P. F. \indexauthorRodríguez-González, A. \resumen Presentamos un modelo analítico de “capa delgada” de la interacción de un chorro bicónico y un medio ambiente en movimiento (alineado con la dirección del chorro), así como simulaciones numéricas (axisimétricas) de dicha interacción. Una situación similar, aunque en un escenario más complejo, ocurre en la cabeza de la estructura cometaria de Mira. Por esta razón, en la mayoría de las simulaciones numéricas exploramos parámetros que son consistentes con las observaciones del flujo bipolar de Mira B. Para estos parámetros la interacción es no-radiativa, lo que resulta en una zona de interacción jet/medio ambiente bastante ancha. A pesar de esto, encontramos que el modelo analítico de capa delgada describe de manera satisfactoria la morfología básica del flujo. \addkeywordhydrodynamics \addkeywordISM: dynamics \addkeywordISM: Jets and outflows \addkeywordStars: AGB and post-AGB \addkeywordStars: Mass loss

## 0.1 Introduction

Prompted by the discovery of a striking cometary tail attached to Mira by Martin et al. (2007), a series of theoretical and observational studies of this object have been made. Theoretical works include Wareing et al. (2007); Raga et al. (2008); Raga & Cantó (2008); Esquivel et al. (2010), all of which are devoted to study the interaction between the post-AGB wind (from Mira A) and a streaming interstellar medium (ISM). This streaming medium corresponds to the motion of Mira through the galactic plane. Observations of the cometary tail have been done in several wavebands, including 21cm (Matthews et al., 2008) and the IR (Ueta, 2008).

Optical (H) observations by Meaburn et al. (2009) revealed a fast bipolar outflow, which probably arises from the less luminous component, Mira B. This outflow is embedded in the head of the cometary structure. It has radial velocities of km s (with respect to Mira), is relatively broad, and its projection onto the plane of the sky is approximately aligned with the direction of Mira’s proper motion.

The mass loss rate from Mira A has estimated values from to Myr (Knapp & Morris, 1985; Bowers & Knapp, 1988). Part of this mass loss is accreted onto a disk around Mira B (and eventually onto Mira B itself). Ireland et al. (2007) estimated a Myr accretion rate. Part of this accretion rate would then be redirected into the bipolar outflow from Mira B.

In the Mira A/B system, we therefore appear to have a bipolar outflow from Mira B which first interacts with the wind from Mira A and at larger distances emerges from the wind into the streaming ISM region. However, in the present paper we consider the more simple problem of an aligned (bi-conical) bipolar outflow/streaming environment interaction (i. e., neglecting the presence of the wind from Mira A). This is a first step to modeling the high velocity outflow observed by Meaburn et al. (2009).

We present an analytic, stationary, thin shell model of the interaction of a bi-conical outflow and an aligned, streaming environment in Section 0.2. Time-dependent, axisymmetric numerical simulations of this problem are described in Section 0.3, and their results presented in Section 0.4. With the simulations we show the characteristics of the flow in greater detail, and evaluate the level of agreement with the analytic, thin shell solution. Finally, the implications of this model for the interpretation of Mira’s outflow are discussed in Section 0.5.

## 0.2 The analytic model

Let us consider a source that ejects a bipolar, conical outflow (with a half-opening angle , velocity and mass loss rate for each outflow lobe), immersed in a uniform, streaming environment of density , which travels parallel to the outflow axis at a velocity . The situation is shown in the schematic diagram of Figure 1. We assume that as the material from the outflow encounters the streaming environment a thin shell of well mixed material is formed (the thick solid line shown in Fig 1). This layer has a bow shock shape , where is the spherical radius and the angle measured from the -axis, see Figure 1).

At a position along the bow shock, the material flowing along the thin shell has a mass rate:

(1) |

where the first term on the right represents the mass fed into the shell by the conical jet, and the second term is the mass from the streaming environment (with being the cylindrical radius, see Figure 1).

The radial momentum rate of the material in the shell is:

(2) |

with a contribution only from the jet. The axial momentum (along the -direction) flowing along the shell is:

(3) |

In equations (2) and (3), and are the - and -components of the velocity of the (well mixed) material flowing along the thin shell. One can also write an expression for the angular momentum rate of the material flowing along the shell:

(4) |

where is the velocity in the -direction (of the material flowing along the thin shell), which is given by:

(5) |

Now, we consider a conical, bipolar outflow with a mass loss rate given by:

(6) |

With this form for , it is possible to carry out the integrals in equations (1-3), allowing to obtain the shape of the bow shock:

(7) |

(8) |

(9) |

where

(10) |

is the on-axis standoff distance between the bow shock and the jet source.

## 0.3 The numerical simulations

We have computed seven different axisymmetric gasdynamic simulations of aligned bi-conical jet/streaming environment interactions with the walicxe code, which is described in detail by Esquivel et al. (2010). The code integrates the gasdynamic equations with an approximate Riemann solver (in this case we used a hybrid HLL-HLLC algorithm, see Harten et al. 1983; Toro et al. 1994; Esquivel et al. 2010). A hydrogen ionization rate equation is solved along with the gasdynamic equations in order to include the radiative losses through a parametrized cooling function, that depends on the density, temperature and hydrogen ionization fraction (Raga & Reipurth, 2004).

The walicxe code has a block based adaptive grid, with blocks of a constant number of cells (in our case, of cells), which are refined by successive factors of 2. Most of our simulations were run at a medium resolution (labeled “mr” in Table 1). These consist of two root blocks aligned in the axial direction, which are allowed to have 8 levels of refinement. The resolution at the finest level would correspond to cells (axial and radial, respectively) in a uniform grid. In addition, we ran low and high resolution versions of model 1 (labeled “lr” and “hr” in Table 1). The low resolution run has 7 levels of refinement ( cells in a uniform grid), while the high resolution run has 9 levels of refinement ( cells in a uniform grid).

Models 1-6 have a computational box of the same size, (radial and axial). Therefore the “mr” runs have a maximum resolution of cm, and the “lr” and “hr” runs have resolutions of cm and cm, respectively. In model 7 we were more interested in resolving the structure at the head of the shock than in the long cometary tail. To achieve a higher resolution ( cm) with the same number of computational cells, in these models we have reduced the size of the computational box by a factor of .

The simulations were performed in a cylindrical computational grid with an axial extent (with the jet source placed at ) and extending radially from the symmetry axis (, where a reflecting boundary condition is applied) out to . An inflow condition is applied in the positive grid boundary, where a streaming environment is injected with a constant velocity in the direction. Outflow boundaries were applied at the negative grid boundary and at . In models 1-3 and , in models 4-6 and , in model 7 and .

The parameters of the simulations are summarized in Table 1. The jets are imposed in two cones with a half-opening angle and an initial outer radius equal to cm for models 1-6, and cm for model 7. The outflows have a velocity (see Table 1) and a temperature K. The density inside the source follows a profile set by the mass loss rate . The environment has an initial temperature K, a density and a velocity . The values for the parameters of models 1-6 are roughly consistent with the jet of Mira B while model 7 (with a denser streaming environment) was chosen to produce a more radiative interaction region. Both the jet and the environment are initially neutral, except for a seed electron density assumed to come from the presence of singly ionized carbon.

In Table 1 we have also included the corresponding value of the bow-shock standoff distance (see equation 10), and the cooling distances associated with the double shock structure at . The cooling distance corresponds to the forward shock (with the streaming ambient medium), while corresponds to the cooling behind the reverse shock. In order to calculate these cooling distances, we have used the fit:

(11) |

to the cooling distances to K given in the tabulation of self-consistent preionization, plane-parallel shock models of Hartigan et al. (1987). In equation (11), is the pre-shock atom+ion number density, and is the shock velocity.

Model | tan | tan | resolution | ||||||
---|---|---|---|---|---|---|---|---|---|

[] | [km s] | [ s] | [km s] | [cm] | [ cm] | ||||

1 | 10 | 200 | 1.00 10 | 125 | 0.044 | 2.97 | 402.4 | 11.9 | lr,mr,hr |

2 | 20 | 200 | 3.9710 | 125 | 0.044 | 2.97 | 194.9 | 5.7 | mr |

3 | 30 | 200 | 8.8210 | 125 | 0.044 | 2.97 | 122.9 | 3.6 | mr |

4 | 10 | 200 | 1.00 10 | 125 | 0.044 | 9.39 | 127.2 | 3.7 | mr |

5 | 20 | 200 | 3.9710 | 125 | 0.044 | 9.39 | 61.6 | 1.8 | mr |

6 | 30 | 200 | 8.8210 | 125 | 0.044 | 9.39 | 38.9 | 1.1 | mr |

7 | 10 | 40 | 1.0010 | 100 | 1.0 | 3.48 | 0.28 | 1.3 | mr |

## 0.4 Results

In Figures 3 - 5, we show the time-evolution (, 3, and 5 yr, respectively) obtained from model 1 (see Table 1). For each time, we show the corresponding mesh configuration, the density (in cm) and the temperature (in K) (top, middle, and bottom panels, respectively). Thanks to the adaptive grid of the code we used less than 10% (for yr) and 20% (for yr) of the finest grid in the domain, thus reducing notably the required computer time.

From the time-sequence shown in Figures 3 - 5, we see that the flow develops a cometary-shaped structure, due to the presence of the side-streaming environmental wind (which flows parallel to the -axis, from right to left). The head of the cometary structure approaches the steady state star/bow shock stagnation distance (see §2) at a time yr.

In Figure 6, we show single time frames ( yr) of the density and temperature stratifications (top and bottom parts, respectively, of each of the three panels) obtained from models 1 through 3 (see Table 1). In each plot we include the corresponding “thin shell” analytic solution (see Section 0.2).

Even though models 1 through 3 differ in and (see Table 1), they result in basically the same flow morphology. Model 3 (which has the largest mass loss rate and opening angle) has the broadest flow. Independently of the and values, all three models coincide with the theoretical stagnation point ( cm), and they also have tails with the same basic morphology. It is clear that the numerical simulations produce broad interaction regions, which are a direct result of the fact that the post-shock cooling distances (, and , see Table 1) are not small compared to the width of the flows. We find that the thin layer analytic solution lies along a line dividing the computed flows into two regions: an internal low-density, hot region (corresponding to the jet cocoon) and a high-density external, cold region (corresponding to the shocked, streaming ambient medium).

The fact that the cooling distances are large compared to the cross section of the jet indicates that the flows are non-radiative. In this case it is remarkable that the analytical solution does successfully divide the hot/low-density from the cold/high-density regions. This is not necessarily expected, as the analytical model assumes that the material is well mixed in a thin layer, which would arise more naturally in radiative flows (where cooling occurs in a small region).

The density and temperature stratifications (at yr) for models 4-6 are shown in Figure 7. These models, with mass loss rates one order of magnitude larger than those from the first three models ( Ms, for details see Table 1) have smaller cooling distances than models 1-3, however they are still in the non-radiative regime. The three models roughly match the analytically derived stagnation point ( cm), and produce broader flows for larger values of the opening angle .

In order to check how well the numerical simulations reproduced the theoretical stagnation point in the thin shell analytical solution from Section 0.2, we show a zoom of the density structures of the jet heads for all models in Figures 8-9. We find that a steady state is reached (with small discrepancies, more evident in models 4-7, Fig. 9), with a shock structure approximately centered on the stagnation radius (see Table 1). One can also see that the place at which the solution curves sharply coincides nicely with the end of the dense wall of the jet.

In order to illustrate the effect of the numerical resolution, we computed two extra simulations with the setup of model 1. The resulting density stratifications for an integration time yr, are shown in Figure 10.

Even though more complex structures are obtained for increasing resolutions, the basic morphology of the flow (stagnation point, cocoon size, density profile, etc.) is the same in all cases. Taking the number of grid points across the jet as an estimate of the Reynolds number of the simulation (Re), in these models we only reach Re (in our hr run). To avoid Reynolds number dependency, the simulations should have resolutions larger by a factor of at least two orders of magnitude greater, which is presently unattainable. However, from the results shown in Figure 8, it is clear that the mr resolution appears to produce the correct basic morphology for the flow. We must note that this result is similar to the one found by Wareing et al. (2007). In this, the authors simulated the interaction of an AGB star moving through the ISM (with also Re ). The key result was that independently of the Reynolds number, the ISM ram pressure stripped material from the head of the AGBÕs bow shock. Such stripping as well as the presence of Rayleigh-Taylor instabilities produced turbulent clump structures that eventually moved down stream to the tail of the bow shock. This mechanism even though in smaller scale, is also present in our simulations (see Figures 3 - 5).

As we have mentioned earlier, the analytic model assumes that the material of the environment and the jet mix well in a thin shell. Such a situation is more easily met in radiative flows, characterized by a thin cooling layer. However, our simulations with parameters appropriate for the jet of Mira B are basically non-radiative. In order to check what happens with a more radiative flow we ran model 7. In this model, the reverse shock is indeed radiative () while the cooling distance of the forward shock is of the order of the cross section of the jet () (see Table 1). This model has a 10 times smaller computational domain, and therefore a 10 times larger resolution than the mr models (see Section 3). In Figure 11, we overlap the thin layer analytic solution on top of the density and temperature stratifications obtained from model 7 at various times ( km s, and km s). For this model the stagnation point is cm, and the cooling distances are , and . As can be seen from the time-sequence, a steady state is not reached in this model. The model develops a turbulent structure which is not stable, and the distance from the jet source to the bow-shock oscillates around .

## 0.5 Discusion

Motivated by the recently observed bipolar jet from Mira (Meaburn et al., 2009), we study the interaction of a steady bi-conical outflow with a streaming environment. We use this scenario to perform 2D axisymmetrical simulations and an analytical model of such an interaction.

We have studied the particular case in which the bi-conical jets are aligned with the streaming environment (as shown in Figure 1). This of course is an idealization of the situation which could possibly be found in Mira, which is unlikely to have a perfect alignment between the outflow axis and the direction of Mira’s motion. It allows, however, simple axisymmetrical analytical and numerical solutions which are useful as a first exploration of the more complex, less symmetric problem of the jet in Mira B.

The proposed “thin shell” analytical solution, based on the momentum equilibrium between the outflow material from the conical jet and the ambient medium, traces the shape of the shock (see Equations 7-9). The model assumes that the interaction occurs in a very thin layer in which the material from the jet and the ISM mixes well. This situation is likely to occur in radiative flows, there the size of the cooling distances are small compared with the flow dimensions. This, however is not the case of the resulting interaction of the outflow launched from Mira B and its surrounding environment.

The models reach a stationary state with a fixed stand-off distance between the jet head and the outflow source. This configuration might be relevant for the collimated outflow from Mira, in which the forward directed jet head has a very low proper motion with respect to Mira (see Meaburn et al. (2009)).

We find that the non-radiative flow configurations (that result from the parameters deduced from the bipolar outflow from Mira) have broad interaction regions which approximately coincide with the predictions from the analytic model (in which a thin shell interaction is assumed). Also, the analytical solutions lie close to the contact discontinuity dividing the material in the jet cocoon from the shocked environment region.

We have tested the convergence of the numerical simulations for increasing resolutions. We find that at different resolutions we obtain similar large-scale flow structures, but with different structures at smaller scales (in which the Kelvin-Helmholtz instabilities associated with shear layers are active). Better convergence is expected only for resolutions higher by at least two orders of magnitude, in which the simulations would start to approach the “high Reynolds number regime” appropriate for the real, astrophysical flow.

We end by noting again that the present model of a bipolar jet/streaming ISM interaction is not completely consistent with the outflow from the Mira system (Meaburn et al., 2009). This is because in the binary Mira AB system, the mass loss rate from Mira A’s uncollimated wind (with 10s, Ireland et al. (2007)) clearly dominates the interaction with the ambient medium. The mass loss from the bi-conical jet (10s, inferred from Martin et al. (2007)), will simply modify the morphology created by the interaction from the isotropic wind and the ambient medium, except in the region in which the leading jet from Mira B emerges from the Mira A wind bubble. This is the region to which our model would in principle apply.

Future models of the outflow from Mira should include both the wind from the primary star and the collimated outflow from the secondary. This more complex flow will clearly be the topic of future papers on Mira’s outflow/ISM interaction.

## Acknowledgments

DLC work has been possible thanks to a DGAPA postdoctoral grant from the UNAM. We also acknowledge support from CONACyT grants 61547, 101256, and 101975. We thank Enrique Palacios and Martín Cruz for supporting the servers in which the calculations of this paper were carried out.

## References

- Bowers & Knapp (1988) Bowers, P. F., & Knapp, G. R. 1988, ApJ, 332, 299
- Esquivel et al. (2010) Esquivel, A., et al. 2010, ApJ, 725, 1466
- Harten et al. (1983) Harten, A., Lax, P.D., van Leer, B. 1983, SIAM review, 25, 357
- Hartigan et al. (1987) Hartigan, P., Raymond, J. C., Hartmann, L. W. 1987, ApJ, 316, 323
- Ireland et al. (2007) Ireland, M. J., et al. 2007, ApJ, 662, 651
- Knapp & Morris (1985) Knapp, G. R., & Morris, M. 1985, ApJ, 292, 640
- Martin et al. (2007) Martin, D. C., et al. 2007, Nature, 448, 780
- Matthews et al. (2008) Matthews, L. D., Libert, Y., Gérard, E., Le Bertre, T., & Reid, M. J. 2008, ApJ, 684, 603
- Meaburn et al. (2009) Meaburn, J., López, J. A., Boumis, P., Lloyd, M., & Redman, M. P. 2009, A&A, 500, 827
- Raga & Cantó (2008) Raga, A. C., & Cantó, J. 2008, ApJL, 685, L141
- Raga et al. (2008) Raga, A. C., Cantó, J., De Colle, F., Esquivel, A., Kajdic, P., Rodríguez-González, A., & Velázquez, P. F. 2008, ApJL, 680, L45
- Raga & Reipurth (2004) Raga, A.C., & Reipurth, B. 2004, RevMexAA, 40, 15
- Toro et al. (1994) Toro, E.F., Spruce, M.,& Speares, W. 1994, Shock Waves, 4, 25
- Ueta (2008) Ueta, T. 2008, ApJL, 687, L33
- Wareing et al. (2007) Wareing, C. J., Zijlstra, A. A., & O’Brien, T. J. 2007, ApJL, 660, L129
- Wareing et al. (2007) Wareing, C. J., Zijlstra, A. A., O’Brien, T. J., & Seibert, M. 2007, ApJL, 670, L125