SHAPING THE RED RECTANGLE PROTO-PLANETARY NEBULA BY A PRECESSING JET

It has been suggested that bipolar or axisymmetric planetary nebulae (PNe) or proto-planetary nebulae (PPNe) can be the result of a binary central system (Bond et al. 1978; Livio et al. 1979; Soker & Livio 1994). This suggestion is supported by the increasing number of binary systems found inside bipolar PNe. De Marco et al. (2008) carried out a survey showing that at least 10%–15% of PNe have binary systems with short orbital periods (<3 days). Although this fraction seems low in comparison with the large fraction of PNe that apparently have been shaped by binary interactions, it has been considered a lower limit because this survey cannot detect large orbital periods.

Bipolar PNe and PPNe can be formed by two jets (or collimated fast winds (CFW)) ejected by one of the components of a binary system (e.g., Morris 1987; Soker & Rappaport 2000). Sahai & Trauger (1998) showed that these jets can actually play a primary role in forming bipolar PPNe (see also Balick & Frank 2002). This binary system could be composed of a white dwarf (WD; the secondary) that accretes material from an asymptotic giant branch (AGB) or post-AGB star (the primary), forming an accretion disk. Two jets or CFWs are then launched almost perpendicular to the accretion disk plane (Frank & Blackman 2004; Frank et al. 2007).

Jets (continuous, time-dependent, or precessing) have been successfully used to explain the morphology of several PNe and PPNe such as PPN CRL 618 (Lee & Sahai 2003), PPN Hen 3-1475 (Velázquez et al. 2004; Riera et al. 2004), PN K 3–35 (Velázquez et al. 2007), and PN IC 4634 (Guerrero et al. 2008). Recently, Raga et al. (2009) studied how mirror and point-symmetric morphologies can be obtained from a precessing jet from a binary system. Haro-Corzo et al. (2009) analyzed the influence of the orbital motion on the line emission of these objects, finding that it could display mirror or point-symmetric distributions, depending on the orientation with respect to the observer.

The Red Rectangle PPN belongs to the bipolar PPNe group because it exhibits a global biconical or "hourglass" shape (Cohen et al. 2004). The binary system HD 44179, which has an orbital period of 322 days, a semi-major axis a less than 1 AU, and a low eccentricity of 0.34 (Men'shchikov et al. 2002), is located at the center of this nebula. This PPN has been called unique because of its biconical morphology and the presence of the so-called ladder rungs. These are linear features almost perpendicular to the nebula symmetry axis, which join some eddy-like structures observed close to the edge of the biconical structure. Soker (2005) presented an analytical model showing that the Red Rectangle morphology is not unique because it can be explained by the same mechanism that is used to model other PNe or PPNe. He considers a CFW blown by the companion star in a binary system. This CFW is emitted intermittently and several double-ring features are formed (Soker 2002) when the CFW material collides with the AGB wind shell. These features can be associated with the "ladder rungs."

In this work, we follow the idea given by Soker (2005), in which an outflow launched by the companion of an AGB star is responsible for the final shape of the Red Rectangle PPN. However, in contrast with Soker (2005), we use a three-dimensional hydrodynamical simulation and propose an alternative mechanism: a continuous but precessing jet (CFW) evolving into an AGB stellar wind. The precession period is a multiple of the orbital period (see Section 2) and produces a jet with a time-dependent structure. We also analyze whether the orbital motion produces important effects on the final PPN shape.

In Section 2 we present the adopted scenario for the orbital and precession motion of the binary system of the Red Rectangle. In Section 3 we describe the initial conditions of the numerical simulation and in Section 4 we discuss our main results of the resulting images. Finally, we summarize our results in Section 5.

The binary system within this PPN (HD 44179) is composed of a post-AGB star (the primary) and a companion which probably is a WD star. Both stars follow elliptical orbits around their center of mass. There are several determinations of the orbital parameters, obtained from spectroscopic studies, which give periods τ_o in the [298, 318] day range and eccentricities e between 0.38 and 0.45 (Van Winckel et al. 1995; Waelkens et al. 1996; Waters et al. 1998). Men'shchikov et al. (2002) recomputed the orbital parameters from original observations and obtained τ_o = 322 days, e = 0.37, and a₁ sin i = 0.34 AU, where a₁ is the semi-major axis of the primary's orbit and i the angle between the orbital axis and the line of sight. They also estimated stellar masses of 0.57 and 0.35  M_☉ for the primary and the secondary components, respectively. In this work we have adopted the orbital parameters given by Men'shchikov et al. (2002).

The general scenario is a WD star that accretes material from the primary, forming an accretion disk, and a collimated outflow or jet (which is launched with a velocity v_j). The jet axis precesses around the z-direction (perpendicular to the orbital plane), describing a cone with a semi-aperture angle α. The precession is retrograde with respect to the orbital motion (Terquem et al. 1999) and its period can be estimated, as a first approach, considering the Figure 1 of Raga et al. (2009). Employing this figure and the masses of this binary system, we have that M₂/M₁ ≃ 0.6 and the resulting precession period τ_p must be more than 10 times the orbital period. In order to model the general "hourglass" morphology of the Red Rectangle nebula, we consider that the orbital plane lies on the x–y plane, and perpendicular to this plane we consider a retrograde precessing jet, which is launched by the companion star. The accretion process is not simulated.

Images of the Red Rectangle PPN show an angular extension of ∼20 arcsec, which is equivalent to 14300 AU or 2.1 × 10¹⁷ cm (Cohen et al. 2004), if a distance to this object of 710 pc is considered (Men'shchikov et al. 2002). At the edge of the biconical structures, several more or less periodic features often referred to as "eddies" (Cohen et al. 2004) are observed. These structures appear to be joined by parallel lines, the so-called ladder rungs, giving the nebula a resemblance to a ladder. These "ladder rungs" are separated by angular distances between 0.9 and 1.5 arcsec, equivalent to 9.5 × 10¹⁵ cm and 1.6 × 10¹⁶ cm, respectively (right panel of Figure 2 of Cohen et al. 2004). These features could be a result of the precession and/or orbital motion of the jet source. The effects of precession are apparently at a scale (Raga et al. 2009; Masciadri & Raga 2002)

$$\begin{equation} D_p= v_j \tau _p, \end{equation} \tag{ 1 }$$

where v_j is the velocity of the collimated outflow from the companion star. On the other hand, the effects of the orbital motion have an associated spatial scale

$$\begin{equation} D_o= v_j \tau _o, \end{equation} \tag{ 2 }$$

where τ_o is the orbital period.

D_o and D_p correspond to the lengths of spiral loops associated with the orbital and precession motion, respectively. From Equation (2) and considering the previous estimates of the distance between successive "ladder rungs," we obtain velocities for the outflow between 2800 km s⁻¹ and 5000 km s⁻¹, which seem to be high for typical outflows from PN or PPN (of ∼100 km s⁻¹). Doing the same estimate with Equation (1) and setting τ_p = 20 τ_o, one finds values for v_j in the [160–290]  km s⁻¹ range. These estimates show that the large-scale morphology of the Red Rectangle could be produced by the precession motion alone. For the case of the AGB wind, it can also be shown that the structure produced by the orbital motion occurs only in scales smaller than those in the Red Rectangle nebula (Cantó et al. 1999). For these reasons we will neglect the effect of the orbital motion in the model of the Red Rectangle nebula. The binary nature of the central source therefore has the sole effect of producing a precession of the ejection direction.

The effect of the precession can be modeled as a variation in the ejection velocity of the jet. In Cartesian coordinates (with the precession axis aligned with z):

$$\begin{equation} v_x=v_j \sin \alpha \cos (\omega _p t), \end{equation} \tag{ 3 }$$

$$\begin{equation} v_y=-v_j \sin \alpha \sin (\omega _p t), \end{equation} \tag{ 4 }$$

$$\begin{equation} v_z=v_j \cos \alpha, \end{equation} \tag{ 5 }$$

where ω_p is the precession angular velocity and t is the time. For the case in which the orbital motion would also be important, it would be necessary to add the x- and y-component of the orbital speed, v_ox(t) and v_oy(t), to the right-hand side of Equations (3) and (4), respectively.

We carried out the three-dimensional numerical simulation with the yguazú-a hydrodynamical code (Raga et al. 2000), which integrates the gas-dynamic equations with a second-order accurate scheme (in time and space) using the "flux-vector splitting" method of van Leer (1982) on a binary adaptive grid. A rate equation for atomic hydrogen is integrated together with the gas-dynamic equations to include the radiative losses through a parametrized cooling function that depends on the density, temperature, and hydrogen ionization fraction (Raga & Reipurth 2004). The time-dependent cooling function for a gas parcel that cools at a constant (atom+ion) number density n = 1 cm⁻³ is plotted in Figure 1.

Zoom In Zoom Out Reset image size

To model the Red Rectangle nebula, we employed a computational domain of (1.1, 1.1, 2.2) × 10¹⁷ cm along the x-, y-, and z-axis, respectively. Five levels of refinement were allowed in the adaptive Cartesian grid, achieving a resolution of ∼4.3 × 10¹⁴ cm at the finest level, corresponding to (256 × 256 × 512) pixels in an uniform grid.

The jet (assumed to be produced by the companion star) evolves into a circumstellar medium (CSM) previously shaped by the dense and slow wind from the AGB primary star. As discussed above, we disregard the influence of the orbital motion on the AGB wind, and for the initial time we impose an isotropic analytical solution of an isotropic constant velocity AGB wind that fills the entire computational domain. This wind has the following density distribution law:

$$\begin{equation} \rho _w=\frac{\dot{M}_w}{4 \pi \ r^2 \ v_w}, \end{equation} \tag{ 6 }$$

where r is the distance from the primary star, $\dot{M}_w$ is the mass-loss rate of the AGB wind, and v_w is the terminal wind velocity. We have assumed $\dot{M}_w=2\times 10^{-6}{\,{M_{\odot }}\ \rm { yr^{-1}}}$, v_w = 15  km s⁻¹, and a constant temperature T_w = 100 K. After several time-integration steps, the temperature acquires a decreasing profile with distance r.

The jet/counter-jet system is injected at the center of the computational domain inside a cylindrical volume with radius r_j and length l_j, both equal to 3 × 10¹⁵ cm (equivalent to 6 pixels in the finest grid). The radius of the jet is larger than the real one at the base of simulated flow, which is a requirement introduced by the limited resolution of the simulations. In other astrophysical flows, such as in HH jets, the initial ejection radius of the jets might be as small as a few AU, but at distances of ∼10¹⁶ cm from the source, their radii expand to 200–500 AU. The jet simulations of these objects are initialized with such a radius (100 AU) rather than the approximately 1–3 AU initial radius of a real HH jet. In the case of the Red Rectangle, the real initial radius of the jet would be a fraction of an AU (the order of the size of the accretion disk around the companion star of HD 44179; Witt et al. 2009), but the width deduced from the width of the observed "rungs" would be approximately 02 (equivalent to 140 AU at a distance of 710 pc). This represents an expansion factor for the width of the jet beam similar to the one seen in HH jets.

The axis of the jet precesses describing a cone with a semi-aperture angle α = 30° (chosen so as to obtain a morphology resembling the Red Rectangle) and a precession period τ_p = 17.6 yr (20 times the orbital period). With the resolution of our numerical simulation it is not possible to resolve spatially the WD orbit, but we have taken care that the jet precession period is well resolved in time (the time step was limited such that the precession is discretized at ∼1° intervals). The jet velocity v_j and its number density n_j were set to 300  km s⁻¹ and 5 × 10⁴ cm⁻³, respectively. With these values, the corresponding mass injection rate of the jet is ${\dot{M}}_j = 1.4\times 10^{-6} \,M_{\odot }\ \rm {yr^{-1}}$, which is one order of magnitude lower than the accretion mass rate estimated by Witt et al. (2009) for the secondary star of HD 44179. In 200 yr (the total integration time of our simulation), both jets inject a total mass of 5.8 × 10⁻⁴ M_☉ to the surrounding CSM. At t = 0 yr, the jet velocity was (0, 150, 260)  km s⁻¹, along the x-, y-, and z-directions. Both the jet and the environment are atomic and initially neutral, except for a seed electron density coming from singly ionized carbon. The jet temperature was set to 10³ K.

3.1. Simulating Synthetic Scattered Light Intensity Maps

In order to compare our hydrodynamical results with observations of the Red Rectangle nebula, we have generated synthetic scattered light intensity maps, which are obtained by using

$$\begin{equation} I_{\lambda }=\int F_{\lambda }\ n\ \frac{\sigma _{{\rm sc}}}{4\pi }\ dl= B_{\lambda }(T_*)\ R^2_*\ \frac{\sigma _{{\rm sc}}}{4}\ \int \ \frac{n}{r^2}\ dl, \end{equation} \tag{ 7 }$$

where B_λ is the Planck's function at wavelength λ = 6130 Å (at the center of the F622W HST filter) and stellar temperature T_*, R_* is the stellar radius, σ_sc is the scattering cross section by dust, n is the number density obtained from the three-dimensional numerical simulation, r is the distance from the source located at the center of the computational domain, and dl indicates that the integration is carried out along the line of sight. In Equation (7), σ_sc can be calculated from

$$\begin{equation} \sigma _{{\rm sc}}= \frac{n_g}{n} \ \pi \ r^2_g \ Q_{{\rm sc}}, \end{equation} \tag{ 8 }$$

n_g and r_g being the dust number density and grain radius, respectively, and Q_sc an efficiency factor. Considering T_* = 7700 K and R_* = 3.1 × 10¹² cm (Witt et al. 2009), n_g/n = 3 × 10⁻¹⁵ and r_g = 4 × 10⁻⁵ cm (Höfner 2008; corresponding to a mass gas-to-dust ratio of ≈830 for grains with a density of 3.3 g cm⁻³), and Q_sc = 0.8, Equation (7) results in I_λ = A∫(n/r²)dl, with A = 1.83 × 10⁸ erg s⁻¹ cm² sr⁻¹ Å⁻¹.

Synthetic maps were obtained with the SHAPE code (Steffen et al. 2011)⁶ from our numerical results (see Section 3.1).

Figure 2 displays the synthetic maps obtained for an integration time of 200 yr, considering that the line of sight is tilted 10° (Men'shchikov et al. 2002) with respect to the y-axis (xz-projection, left panel) or with respect to the x-axis (yz-projection, right panel). Both projections show a clear biconical morphology, with a size along the symmetry axis (z-axis) of 2 × 10¹⁷ cm. This value corresponds to an angular size of ≃ 187, which is similar to that of the F622W PC image obtained by Cohen et al. (2004). This biconical structure is observed in our synthetic maps as a double X shape which corresponds to the walls of the cocoon produced by the jet. As the jet evolves around the precession axis both the jet and the cocoon describe a helical path, quite evident in Figure 2. Due to the small inclination with respect to the plane of the sky, this helical path looks almost horizontal in some parts, which can be associated with the "ladder rungs" observed in the Red Rectangle PPN images (Cohen et al. 2004; Soker 2005). The distance between successive "rungs" is similar to the observed one. Furthermore, at the edges of the conical structures, the helical structure appears more curved and resembles the eddy-like features identified by Cohen et al. (2004).

Zoom In Zoom Out Reset image size

The inner wall of the biconical structure can be identified with the inner X region of the observations, although in the simulated maps these do not intersect at the center. The outer wall of the cocoon shows a single-hourglass structure, while observations revealed several parabolic features. This discrepancy may be due to the fact that we have not considered an AGB wind with several episodes of emission, which would generate a structure with multiple concentric shells, where the jet propagates.

There are several V-shaped structures that seem to be pointing to the central source. They are clearly seen in the xz-projection of Figure 2, while in the yz-projection they can be confused with the eddy-like features. These features are a direct result of the initialization of the simulation, and they are produced when the jet is suddenly "turned on" at t = 0. These structures are not observed in the Red Rectangle images, because of a more gradual "turning on" of the jet is likely to have occurred in the real case.

Figure 3 shows that in a central region, of approximately 3'' in diameter, there is an emission peak, which corresponds to the base of the jet (the bright jet region within 10¹⁶ cm from the source, present in both projections of Figure 2). This peak represents an increase of a factor of approximately 200 in flux with respect to the periphery of the nebula. Comparing with the bottom panel of Figure 4 of Cohen et al. (2004), it is seen that the inner region (of 3'' diameter) of the Red Rectangle has a peak with a flux increase of a factor 1000. Therefore, if this nebula had a jet-like structure with the intensity predicted from our model, it would not be evident in the Hubble Space Telescope (HST) images against the bright central component due to the circumbinary disk around the HD 44179 (Men'shchikov et al. 2002; Witt et al. 2009), which is not included in our simulations. Notwithstanding, a qualitative and quantitative agreement is achieved between the observed flux density profile (see Figure 4 of Cohen et al. 2004) and our simulated profile, for the regions beyond the innermost 3''. After converting the synthetic profile given in units of erg cm⁻² s⁻¹ Å⁻¹ sr⁻¹ to erg cm⁻² s⁻¹ Å⁻¹ pixel⁻¹, considering the HST resolution of 00456 pixel⁻¹, our simulated maps give a flux a factor of 3.3 lower than the observed one. Therefore, we obtain a flux of the same order of magnitude as the observed flux in the outer region of the nebula. A local minimum is observed at r ≃ 0 in the synthetic emission profile. This minimum is produced because we have disregarded the emission coming from the cylinder where the jet condition is imposed at each time-integration step.

Zoom In Zoom Out Reset image size

In Figure 4, we explore how the general morphology changes at different angles i between the precession axis of the biconical structure and the line of sight. The biconical or hourglass morphology appears at angles i  ⩾  60°. At i = 30°, the general shape becomes in two partial shells that looks similar for instance to the PPN He 2–104 (Sahai & Trauger 1998). The nebular morphology becomes a single shell with some ripples when i = 0° which could resembles the morphology of PPNe He 2–131 and He 2–138 (Sahai & Trauger 1998). Our results support the idea that several additional morphologies of PPNe or PNe may be produced with the same mechanism, as was previously suggested by Soker (2005). In this case, the mechanism is a precessing jet.

Zoom In Zoom Out Reset image size

We performed a three-dimensional numerical simulation of a precessing jet, moving into a CSM previously filled by an AGB stellar wind, with the idea of reproducing the observed large-scale morphology of the Red Rectangle PPN.

The scenario presented in this work is based on the fact that a binary system is located inside of this PN (HD 44179). The companion star receives material from the primary one, forming an accretion disk. A jet is assumed to be launched from the companion star, whose axis precesses in a cone with semi-aperture angle α and period τ_p. This setup is similar to the one used by Soker (2005), with the difference that in that work an intermittent outflow was considered instead of a precessing jet.

Reasonable agreement with observations is achieved with a precession angle α = 30°, a period τ_p = 17.6 yr (20 times the orbital period), and with a jet velocity of 300 km s⁻¹. Synthetic scattered light maps obtained from our numerical results, at an integration time of 200 yr, display a global biconical morphology with an angular size of 187 (considering a distance to this object of 710 pc) which is comparable with the observed one. An internal helical structure is seen, which can be associated with the observed "ladder rungs." In addition, our models also reproduce the qualitative and quantitative behavior of the observed brightness profile and the contrast in the intensity of emission from regions close to the nebula center and its periphery. The orbital speed (of the order of 40 km s⁻¹ at periastron) has a negligible effect on the global morphology of the PPN.

The "Hourglass" morphology can be observed as an incomplete double, or a single shell, depending on the point of view of the observer. This demonstrates that a single mechanism, in this case a precessing jet, can produce several of the morphologies observed in PN and PPN, in agreement with the previous results of Soker (2005).

In summary, a precessing jet model can reproduce the size and general morphology of the biconical structure, the so-called ladder rungs and eddies, and the order of magnitude of the observed flux of Red Rectangle PPN. The orbital motion has no direct effect on the global PN morphology. This result supports previous conclusions by Soker (2005), who claimed that an ad hoc model is not necessary to explain the peculiar morphology of this PPN.

The authors thank the referee for useful comments. Authors also thank Hans Van Winckel for his very useful suggestions and discussions. We thank Enrique Palacios for maintaining the Linux servers where the simulation was performed. P.F.V., A.C.R., A.E., and J.C. acknowledge support from CONACyT grants 61547, 101356, and 101975, and UNAM DGAPA grant PAPIIT IN119709. W.S. acknowledges support from UNAM DGAPA grant PAPIIT IN100410. The work of A.R. was supported by the MICINN grants AYA2008-06189-C03 and AYA2008-04211-C02-01 (co-funded with FEDER funds).