#
3D hydrodynamical models of point-symmetric planetary nebulae: the special case of H 1-67^{†}^{†}thanks: Based upon observations carried out at the Observatorio Astronómico Nacional on the Sierra San Pedro Mártir (OAN-SPM), Baja California, México.

###### Abstract

We present 3D hydrodynamical simulations of a precessing jet with a time-dependent ejection velocity or a time-dependent ejection density, interacting with a circumstellar medium given by a dense, anisotropic, and slow AGB wind, forming a torus. We explore a set of configurations with different values for the precession angle and number of ejections. The temporal evolution of these models is analised at times up to 1500 or 1800 yr. From our hydrodynamical models, we obtain position-velocity diagrams (PV diagrams) in the [N ii] 6583 line to be compared with high resolution observations of the planetary nebula H 1-67. From spectral data this object shows high-velocity jets and a point-symmetric morphology. With our synthetic PV diagrams we show that a precessing jet with a time-dependent ejection velocity or a time-dependent ejection density reproduce the point-symmetric morphological structure for this nebula if the precession cone angle is larger than 30. Our synthetic PV diagrams can be used to understand how the S-like morphology, also presented by other planetary nebulae, is formed. For H 1-67 we found a heliocentric velocity of -8.05 km s and height below the galactic plane of 451.6 pc.

###### keywords:

Planetary nebulae: general — hydrodynamics — methods: observational — methods: numerical — ISM: jets and outflows.^{†}

^{†}pagerange: 3D hydrodynamical models of point-symmetric planetary nebulae: the special case of H 1-67

^{†}

^{†}thanks: Based upon observations carried out at the Observatorio Astronómico Nacional on the Sierra San Pedro Mártir (OAN-SPM), Baja California, México.–A

^{†}

^{†}pubyear: 2002

## 1 Introduction

Planetary nebulae (PNe) are expanding regions of ionized gas, whose central stars had initial masses lying in the range of 0.8 - 8 M. In a traditional scenario, PNe are formed in a phase of mass loss of stars that are leaving the Asymptotic Giant Branch (AGB). These stars evolve towards high effective temperatures at about constant luminosity and ionize the ejected shell, afterwards they become a white dwarf (WD).

It is well known that PNe have a wide variety of shapes, showing shells, jets, torii, knots or internal structures and they are classified according to their morphology (Balick, 1987; Parker et al., 2006; Sahai et al., 2011). The main morphological groups are: round, elliptical, bipolar, and multipolar. Some PNe show a jet-like morphology which can be indicative of a PN evolving from a binary system (Soker & Livio, 1994; Bond, 2000; De Marco, 2009; Soker & Rappaport, 2000).

In this morphological scheme PNe catalogued as point-symmetric were introduced by Schwarz et al. (1993) and Stanghellini et al. (1993). These are defined
as PNe whose morphological components show point reflection symmetry with respect to the centre,
similar to a S-shape. Examples of point-symmetric PNe are IC 4634
(Guerrero et al., 2008), NGC 6309 (Vázquez
et al., 2008), and others. A very interesting case is the PN Fleming 1 (Fg 1) which shows a very extended S-shape where the filaments are moving at a high velocity of 75 km s (López et al., 1993).
See the catalogue by Sahai
et al. (2011) and the compilation by Guerrero (2000) for a list of point-symmetric PNe.

Point-symmetric PNe have been studied in previous works (Miranda &
Solf, 1992; López et al., 1993).
They are situated at scale height over the galactic plane of 310 pc (García-Segura
et al., 2002) which is smaller than the general scale height of disk PNe.
It has been suggested that the formation of this morphology is a consequence of a precessing collimated outflow (Guerrero, 2000) and a binary system (Livio &
Pringle, 1996).

Some models have been produced to understand this kind of morphology (García-Segura &
López, 2000; Rijkhorst
et al., 2004; Raga
et al., 1993, 2009).
Cliffe et al. (1995) calculated a hydrodynamical model for the formation of point-symmetric PNe and they concluded that the knots seen in these PNe can be interpreted as
resulting from the periodic ejection of dense material by a precessing jet.
3D numerical simulations like
the ones by Haro-Corzo et al. (2009) have shown that a jet ejected from a binary system can lead to morphological structures of point-symmetric PNe.

The ideas described above have been succesfully applied to, e.g., the well known point symmetric PN IC 4634 by Guerrero et al. (2008), who from very deep HST images and high-resolution spectroscopy, and by using the code Yguazú (described in §2.1), constructed a hydrodynamical model, which included a precessing fast collimated outflow interacting with nebular material. Due to they found low radial velocities of the S-shape structure ( 20 km s), these authors have assumed that the precession axis is on the plane of the sky.

Other well known point symmetric nebulae such as NGC 6309 and Fg 1 have been studied with deep images and kinematical models. It is worthy to say that the extended S-shaped Fg 1 resulted to possess a binary stellar system which is the most probable cause for its morphology (Boffin
et al., 2012). A close binary nucleus has also been found in the spectacular quadrupolar PN NGC 5189. This binary system is composed by a hot [WO 1] star and a massive WD (Manick
et al., 2015). The nebula contains two central toroidal structures (one IR and the other optical) associated to bipolar outflows and it shows multiple point-symmetric low-ionization knots or ansae and collimated outflows. Sabin
et al. (2012) found moderate expansion velocities in the three bubbles conforming the nebula, but the largest velocities found at both extremes of the nebula are about 40 km s and 50 km s. Other PNe with precessing collimated outflows are Hen 2-1, Hen 2-141, Hen 2-429, Hu 2-1, J320, M 1-61, NGC 5307, NGC 6881, PC19 among others.
For none of these objects reported as point-symmetric PN, except IC 4634, a specific hydrodynamical model has been constructed.

In this and in our previous work, our aim has been to analyse compact PNe presenting round or toroidal images but showing high velocity ejections larger than 70 km s observed in Position-Velocity diagrams. Due to their appearance, these objects have been classified as ‘compact’ or ‘elliptical’ PNe. The objects in our analysis have been selected from the San Pedro Mártir (SPM) Kinematical Catalogue of Galactic Planetary Nebulae (López et al., 2012).

By analysing position velocity (PV) diagrams of high resolution spectra of PNe that in images appear compact with no evident structure, we have found that some of them show high velocity ejections that seem to correspond to collimated material thrown through the poles of a torus. Such are the cases of M 1-32 and M 3-15 (Rechy-García et al., 2017). Other objects show point-symmetric features at high velocities.

One of the most interesting cases is the PN H1-67, presented in Figs. 1 and 2. The image and spectra were obtained with the Manchester Echelle Spectrograph (MES) attached to the 2.1-m telescope at Observatorio Astronómico Nacional, San Pedro Mártir (OAN-SPM), B.C., México, on the night 2013 May 18. A detailed analysis of these data is presented in §4. Here we want to show that this PN looks compact in direct images (see Fig. 1 and also the image presented by Reese & Zijlstra 2013), showing a relatively face-on “broken” torus with two bright condensations on the East and West sides, while the PV diagram (Fig. 2 top) obtained from the spectrum with the slit crossing the nebula through the central position, along the symmetry axis, at P.A. = 45, shows two central bright condensations corresponding to the broken torus and a S-shape ejection ending in two knots at about 100 km s. This velocity structure corresponds to a point-symmetric nebula, fact that is not noticed in the direct images. In addition, the PV diagram obtained with a slit crossing at P.A. = 0 shows the condensations corresponding to the broken torus and two knots at velocity of about 100 km s (Fig. 2 bottom). No S-shape is noticeable in this latter spectrum.

From the SPM Kinematical catalogue mentioned above, other compact nebulae such as M 2-36 and Wray 16-411, presented in §5, show similar S-shape velocity structures.

The main aim in the present study is to examine the formation of such point-symmetric morphologies. This can be carried out by analising the history of the injected gas, which shapes this kind of nebulae. Following Cliffe et al. (1995) we will study the effect of a periodic ejection of dense material in the form of a precessing jet, by means of hydrodynamical models. We calculate numerical simulations of a jet that changes its direction describing a precession cone in order to reproduce the S-like morphology and we include a time-dependent ejection velocity or a time-dependent ejection density, in order to explain the observed structures in these PNe. Furthermore, the jet interacts with a “clumpy” slow wind, equatorially concentrated emulating the presence of a broken dense torus.

This paper is organized as follows, we describe the hydrodynamical models and the initial settings in §2. In §3 the different morphologies obtained from the models are presented. In §4 we compare our results with observations of PN H 1-67 that presents S-like morphological features and compute specific models for this object. Other similar objects are presented in §5 and our conclusions are summarised in §6.

## 2 Numerical simulations

### 2.1 The hydrodynamical model and the Yguazú code

In order to reproduce the observational S-shape PV diagrams, we have carried out 3D hydrodynamical simulations taking into account the observed features of H1-67 described above: a broken dense torus and several clumps, at high and medium velocities, distributed along a filament with a S-like morphology.

From the H image it is evident that H 1-67 have a dense ‘broken” torus and that the high velocity material flows through the poles of the torus. In order to calculate a hydrodynamical model for this object, we have considered a circumstellar medium with a density distribution given by a dense and low-velocity AGB wind forming a torus. Subsequently, a second outflow is launched by the central star. Different types of second outflow were tested: an isotropic wind, an anisotropic wind (with anisotropical distributions in velocity or density), a cylindrical jet and a conical jet. The model that best explains the PV diagrams observed for H 1-67 consists of two components: a dense torus and a cilindrical jet.

Additionally, and in order to simulate a broken dense torus, we have built a clumpy wind by modulating the density by a fractal structure with a spectral index of 11/3, which generates a noisy structure and it is consistent with a turbulent interstellar medium (Esquivel et al., 2003). We imposed this clumpy density on the initial condition and this fluctuates with a 20 per cent of the mean density value.

Finally, to obtain the observed distribution of the clumps in the PV diagrams, a precession and a velocity ejection variability of the outflow were introduced. Alternatively, a model including a time-dependent ejection density of the outflow has been computed and it will be described ahead.

The numerical simulations were performed with the hydrodynamical code Yguazú (Raga et al., 2000). This is a 3D adaptive code which solves the gasdynamic equations employing the “flux vector splitting” scheme of van Leer (1982). Together with these equations, a system of equations for atomic/ionic species are also integrated in order to provide a cooling function. For these simulations, we have used the total abundances of N, O and S (relative to hydrogen), obtained by Escudero et al. (2004) for H1-67, which are similar to the ones presented by many PNe in the Galaxy.

### 2.2 Initial setup

Our numerical simulations use an adaptive
Cartesian grid. We computed several models varying the parameters of the outflow.
The sizes of the computational domain are of cm along
the x- and y-directions, and 210 cm along the z-axis for model 1, while for models 2, 3 and 4 (the latter one with a jet with density variations), the size along the y-axis increases to 210 cm (see Table 1).

At the beginning, following the work by Mellema
et al. (1991), to form the torus we have imposed on all computational domain the density distribution of an anisotropic AGB wind, which is given by:

(1) |

where is the distance to the centre of computational domain, is the polar angle measured with respect to the z-axis, and is the density for a reference radius , which is given by:

(2) |

being the mass loss rate and the terminal wind velocity, which was set to 20 km s. The angular dependence of the wind is described by the following equation:

(3) |

where is the parameter which determines the equator to pole density ratio, and gives the variation of the density from the equator ( plane or ) to the poles ( direction or ; direction or ).

A jet is imposed at the centre of the computational domain in a cylindrical volume with radius and a length (see Fig. 3). The jet direction describes a precession cone, with a semi-aperture angle and a period , where is a constant representing the factor of precession, and is a characteristic time given by:

(4) |

Furthermore, the jet velocity depends upon the time as:

(5) |

where is the mean jet velocity, is the amplitude of the velocity variation, is the time, and is the period of the velocity variation, which is given by , being a constant (the values of and used in our simulations are listed in Table 1). This jet has a constant density of 150 cm.

In the case of a jet with a time-dependent density ejection, the jet velocity is fixed to and the density is given by:

(6) |

where is the mean density of the ejection, is the amplitude of the density variation, is the time and is the period of the density variation. In this case, in Eq. 4 changes accordingly.

In all simulations the values of parameters and in Eq.(3) were chosen as 0.99 and 5.0, respectively, because these values are adequate to produce a highly- contrasted thin disk (Mellema et al., 1991). Such a disk (torus) constraints the stellar ejection helping to produce a bipolar nebula. In this torus the flow to external density ratio is much smaller than one in the equator ( = 90 in Eq.(3)) therefore the injected gas rapidly slows down and the ejection of the jet is prevented in this direction. Instead, this ratio is about 2 when = 0 (pole direction), therefore the jet escapes easily. We have also explored a model with a thick torus (with B=2) finding that in this case the jet is slowed, being necessary to employ larger integration times. However, this model is not successful in reproducing the “S” morphology observed in the PV diagrams.

The temperature for the AGB wind was fixed as 1000 K while its was of . The value of was chosen as 5 pixels of the grid at its maximum resolution, i.e., cm. For the jet with velocity variability we have set and as 140 km s and 0.5, respectively, while for the jet with density variability is 150 cm and is 0.7. The variability period is the same in both cases.

Initially we computed a model with a precessing cone with a semi amplitude angle equal to 15. Such a model did not produce a S-shape in the PV diagrams of the ejection and we discarded it. This model is not presented here, but helped us to understand that a larger is required. Afterwards we computed models with of 30 and 40 which show S-shape structures in the PV diagrams.

Different values of , and were used (see Table 1) to analise the PV diagrams produced. In the following we describe the properties of Model 1 with = 30 , =4 (the jet precess 4 times), =2 (two ejections per precession) and Model 2 with =40, =1.5 and =2.

Model | Computational domain | |||||
---|---|---|---|---|---|---|

, , and axis, (cm) | yr | yr | ||||

model 1 | (1.05, 1.05, 2.10) 10 | 4 | 2 | 30 | 230 | 1830 |

model 2 | (1.05, 2.10, 2.10) 10 | 1.5 | 2 | 40 | 690 | 2070 |

model 3 | (1.05, 2.10, 2.10) 10 | 1.5 | 4 | 40 | 345 | 2070 |

model 4 | (1.05, 2.10, 2.10) 10 | 1.5 | 4 | 40 | 304 | 1826 |

model with variable density. In column 6 is presented. |

## 3 Numerical results of hydrodynamical models

### 3.1 Density stratification evolution

Fig. 4 displays the temporal evolution of the electron density stratification at integration times of 500, 1000, and 1500 yr for models 1 and 2 mentioned above. These maps show density distribution projected on the xz-plane. The white arrows indicate the velocity field of the material in the jet.

The temporal evolution of the electron density maps on the yz-projection for these models, are shown in Fig. 5. As in Fig. 4, the maps for integration times of 500, 1000, and 1500 yr are displayed for model 1. For model 2, due to the precession angle is larger, = 40, we have enlarged the computational domain in the y-direction for the jet to fit into the computational domain. Therefore, we only show two evolution times for this model of 300 and 1500 yr.

For all density maps, both axes are given in units of 10 cm and the logarithmic color scale gives the number density in units of cm.

### 3.2 Synthetic PV diagrams

To compare with observed PV diagrams of point-symmetric planetary nebulae, from our hydrodynamic models we generate an array of synthetic PV diagrams at evolution time of 1500 yr, where we have considered different angles of inclination of the -axis of the model with respect to the line of sight (), as well as different counterclockwise rotations in the plane of the sky (). Besides we have considered both xz- and yz-projections.

#### 3.2.1 The xz-projection

Figs. 6 and 7 show the synthetic [N ii] 6583 PV diagrams (-projection) for models 1 y 2, respectively. Top, middle, and bottom panels of these figures display the PV diagrams for .

In Fig. 6, the PV diagrams for model 1 show two bright condensations with velocities close to 0 km s for , which correspond to the dense torus. As the angle increases these bright condensations appear slightly more separated in velocity (up to 20 km s for = 60). The structures associated with the ejections from the precessing jet, include a long filament at position 0 arcsec, which contains several ejections extending in velocity up to 200 km s and some faint knots, spatially separated from the centre by 2 or 3 arcsec
(see Fig. 6, bottom).
As the angle increases these knots show lower radial velocities, as a consequence of the projection effect (see Fig. 6, middle). In Fig. 6 top, where , the knots are overlapped with the emission associated with the dense torus.

The PV diagrams for model 2 are shown in Fig. 7. The behaviour of the bright condensations (dense torus) is similar to model 1 and also the filament shows the same behaviour as in model 1. For this model, the precession angle is 40, consequently the two knots, corresponding to the ejections, appear more spatially separated when (about 4 arc seconds from the centre, fig. 7 bottom). As increases these features show lower velocity and appear closer to the centre. When both features are overlapped with the torus, similar to model 1.

#### 3.2.2 The yz-proyection

The PV diagrams for models 1 and 2 obtained for the -projection are shown in Figs. 8 and 9, respectively. Both PV diagrams display similar morphologies. Several filaments are observed along the horizontal axis, which are produced by the variability of the jet ejections. All the filaments display a S-like shape. Model 1 has more filaments due it has a factor of precession of =4, while model 2 has =1.5. The filaments result more separated along the position axis for model 2 than for model 1. This is a consequence of a larger semi-aperture angle of precession (40 and 30, respectively).

## 4 The S-shape spectrum of H1-67: Observational Results

As it was said in the Introduction the PN H 1-67 presents S-shape spectra as obtained with the MES at OAN-SPM with the slit is located along its symmetry axis. Therefore this is a suitable object to compare its PV diagrams with our numerical simulations. Some characteristics of this and other similar objects are listed in Table 2.

H 1-67 is located towards the galactic bulge, showing a heliocentric radial velocity of 8.05 km s (this work). Its heliocentric distance is 5.88 kpc (Frew et al., 2016) and it is situated at a height of 451.6 pc (this work). It is ionised by a [WO 2] central star (Acker & Neiner, 2003). The total nebular abundance ratios for this PN are N/O = 0.56 and He/H = 0.125 (Escudero et al., 2004), therefore, this object is a nitrogen and helium rich PN (Peimbert Type I PN). Morphologically it has been classified as a bipolar nebula (Rees & Zijlstra, 2013). From our spectroscopic analysis it is found that the nebular morphology of H 1-67 is highly asymmetric, showing the appearance of broken face-on torus (Fig. 1) with low expansion velocity, and also knots and jets at high velocities of 100 km s.

In this work, we use high-resolution spectra to analyse the kinematical
behaviour of the gas in this [WO 2]PN and, by comparing it with the
hydrodynamical models presented above, we investigate if the morphological structure of
H 1-67 can be produced by the action a precessing jet
and a dense torus. We compare the
position-velocity diagrams obtained from the observations with those
predicted by the hydrodynamical models.

PN G | Name | Stellar | Dmean | z | N | O | |

Type | km s | kpc | pc | ||||

009.804.6 | H 1-67 | [WO 2] | 8.05 | 5.881.76 | 451.6 | 8.50 | 8.75 |

003.206.2 | M 2-36 | possible binary | +57.0 | 6.99 | – | 8.42 | 8.85 |

353.712.8 | Wray 16-411 | – | 50.0 | 6.141.75 | – | – | — |

In 12 + log /H (see Escudero et al. (2004) and Liu et al. (2001)) |

### 4.1 Long-Slit Spectroscopy

The Manchester Echelle Spectrometer (MES, Meaburn et al., 2003) is a long-slit echelle spectrometer, which uses narrow-band filters to isolate the orders containing the emission lines. The 90 Å bandwidth filter allowed us to obtain the H and [N ii] nebular emission lines. In this paper, we only show the behavior of the [N ii] 6583 line because H shows the same structure.

The Marconi 2 CCD was used as detector, which has a pixel size of 13.5 m. A binning of 22 pixel was applied in the observations and, using the secondary mirror f/7.5 which provides a plate scale of 13.2 /mm, the plate scale on the detector was 0.356 /pixel. The slit width used was 150 ( on the sky), providing a spectral resolution of 11 km s.

We obtained observations from two slit positions over H 1-67, one position with the slit at P.A. = 45 across the centre of the nebula corresponds to its symmetry axis.
The other position with the slit at P.A.= 0, lies on the East side. The slit positions are indicated
in Fig. 2 (left side) on H images obtained by us.
The integration time for the spectra was 900 s for both slit positions.
Immediately after every science observation,
exposures of a Th-Ar lamp were acquired for wavelength calibration. The internal precision
of the lamp calibrations
is better than 1.0 km s. The spectra were reduced using the standard
processes for the MES data, using IRAF^{1}^{1}1IRAF is distributed by the National Optical Astronomy Observatories, which is operated by the Association of Universities for Research in Astronomy, Inc., under contract to the National Science Foundation. tasks to correct for bias and calibrate in wavelength.

### 4.2 Position-Velocity diagrams of H 1-67

PV diagrams for this PN were obtained with the WIP
software (Morgan, 1995). In Fig. 2 right side, the PV diagrams
are shown for the slits located in the two positions, P.A. = 45 and P.A. = 0.
In the horizontal axis the heliocentric
radial velocity is represented and the vertical axis shows the spatial
direction. The emission line shown is [N ii] 6583.

In the PV diagram with the slit at P.A. = 45, two bright inclined condensations (at a heliocentric systemic velocity of 8.05 km s) are visible, they are labeled as B and C in Fig. 2 (top). Two knots at high velocity of 97.07 km s (relative to the systemic velocity) are seen, labeled as A and D in Fig. 2 (top). In the PV diagram with the slit at P.A. = 0, two bright condensations are visible (they are labeled as B’ and C’ in Fig. 2, bottom) and two knots at high velocity of 96.4 km s (relative to the systemic velocity) are seen, labeled as A’ and D’ in Fig. 2. In addition, comparing both PV diagrams, we observe that the spectra obtained with the slit at P.A. =45 shows a S-shape, differently from what is obtained with the slit at P.A. = 0. This suggests a bipolar, collimated jet whose ejection axis has changed following a precessing movement. The observations are not calibrated in flux, but we have added a scale in counts. In the observational PV diagrams, the ratio between the brightest (torus) and faintest (filament) regions is of the order of 20.

#### 4.2.1 Emission line profile

One-dimensional profile of the [N ii] emission-line for the two positions of the slit are displayed in Fig. 10. The wavelength is represented in the horizontal axis, and the vertical axis shows the intensity in counts. The features marked in Fig. 2 (for both positions) are shown in this figure. Both profiles are kinematically similar. We can distinguish two regions in this PN, an inner zone (B, C and B’, C’) and an outer zone (A, D and A’, D’). The inner zone that corresponds to the clumpy torus is expanding at about 18 km s while the outer zone corresponding to a precessing jet ending in knots, shows an expansion velocity of about 97 km s. Note that the radial heliocentric velocity in the knots marked as A,D and A’,D’ in both positions slits are similar. In Table 3 we present the kinematics for each feature.

P.A. = 45 | |||||
---|---|---|---|---|---|

Knot | A | D | B | C | |

km s | km s | km s | km s | ||

Vhelio | -103.7 | 90.40 | -26.3 | 10.6 | |

Vexp | 97.06 | 18.46 | |||

P.A. = 0 | |||||

Knot | A’ | D’ | B’ | C’ | |

km s | km s | km s | km s | ||

Vhelio | -103.4 | 89.3 | -28.6 | 7.3 | |

Vexp | 96.4 | 18.0 |

### 4.3 Computing models to reproduce H 1-67 observational results: Models 3 and 4

Analysing the synthetic PV diagrams obtained for the different models, we found that the S-shape is reached in the yz-projections when the precession angle is larger than 30. Whenever the angle is greater, the jet ejections open giving rise to a more defined S-shape. Also, the number of knots and filaments obtained depends on the number of precessions, , and the number of ejections in each turn, . Therefore, should be small to obtain less ejections in the PV diagrams (compare Fig. 8 and Fig. 9).

We notice that the observed PV diagrams shown in Fig. 2 for H 1-67, have some features similar to those obtained in the PV diagrams of model 2, for the yz-projection.
Due to this, we have carried out a 3rd model in order to reproduce better the main characteristics of the observed PV diagrams of H 1-67. In Model 3, which has =40 and =1.5, the ratio between the precession period and the velocity variability period, factor , was set to 4. Therefore model 3 is similar to model 2, but with more ejections per precession turn.

The electron density distribution map for this model in the xz-projection at evolution times of 300, 900, and 1500 yr, is presented in Fig. 11 (top) and the yz-projection is presented in Fig. 11 (bottom) at evolution times of 300 and 1500 yr. As expected, these maps are similar to the ones for model 2.

Fig. 12 shows the synthetic PV diagrams performed for the xz-projection considering angles =60, 40, 20 and =5, 10, 10. Fig. 13 shows the synthetic PV diagrams performed for the yz-projection considering angles =60, 40 and =18, 40.

For H 1-67, our observational PV diagram with the slit at 45 (Fig. 2, top) is well reproduced by the synthetic PV diagram with =60 and =40 (Fig. 13, middle). We reproduced the central condensations corresponding to the dense torus, separated almost 40 km s approximately. A clear point-symmetric morphology results from the presence of the precession jet. The synthetic PV diagrams show fast outflows up to velocities of 95 km s, with a S-shape that is similar to the observations.

We also obtained a good agreement with the PV diagram from observations along P.A. = 0 if the
jet axis was tilted at =60 and =18 (Fig. 13 top).

#### 4.3.1 Model 4, with variable ejections density

Instead of time-dependent ejection velocity for the jet, we have computed a model which considers a time-dependent ejection density for the jet, given by Eq. 6.
The parameters for this model are the same as those for Model 3 (the same , and ). In this case, the jet has constant velocity of 140 km s. The average density is 150 cm and the density variation is 70%.
The total integration time was of 1800 yr.

In Fig. 14 (top) the xz-projection of the electron density distribution map at 300, 900 and 1800 yr is shown while in the bottom figure, the yz-projection is presented at 300 and 1800 yr. These figures are similar to the ones presented for Model 3, although in this case the knots at high velocity appear denser and more defined.

The synthetic PV diagrams for this model are presented in Fig. 15 where three diagrams are showed with (, ) of (60,18), (60,40) and (40,40). The one corresponding to =60 and =18 reproduces well the observations with the slit at P.A.= 0. The second synthetic PV diagram for =60 and =40 reproduces well the observed PV diagram for H 1-67 for the slit at P.A.= 45. In both panels, the knots observed at 120 km s appear well defined in the model. The third synthetic PV diagram show =40 and =40, this PV diagram does not reproduce the observations.

#### 4.3.2 Synthetic Images

We have built synthetic images corresponding to the H, [N ii], [S ii] emissions from our model 4 for an evolution time of 1800 yr representing the image that could be observed with a very high spatial resolution like that of HST. Fig. 16 shows the synthetic images for the emission of [N ii] and [S ii]. H is not shown because presents the same structure. It should be noticed that projection does not include the rotation in the sky plane. The brightest emission corresponds to the torus and the filaments are observed above and below, showing a S-shape. These head of filament are moving at 130 km s, a velocity highly supersonic therefore they correspond to a shocked gas.

PN G | Name | Stellar | Outflow velocity | A | B | Spectral | Precession |
---|---|---|---|---|---|---|---|

type | km s | morphology | |||||

006.8+04.1 | M 3-15 | [WC 4] | 90 | 0.99 | 0.1 | Bipolar | No |

009.804.6 | H 1-67 | [WC 2-3] | 97 | 0.99 | 5 | point-symmetric | Yes |

006.8+04.1 | M 1-32 | [WC 4] pec | 180 | 0.9 | 5 | Bipolar | No |

model with variable density. In column 6 is presented. |

## 5 Other objects

Using the available spectra in The SPM Kinematic Catalogue of Galactic Planetary Nebulae (López et al., 2012), we have found two PNe, M 2-36 (PN G 003.206.2) and Wray 16-411 (PN G 353.712.8), with similar spectral characteristics to H 1-67. In their spectra they present a S-shape morphology with knots or filaments at high velocity (see Fig.17).

Some of their characteristics are listed in Table 2.
M 2-36 is located towards the galactic bulge, showing a radial velocity of 57 km s (Richer et al., 2017) and has an ADF^{2}^{2}2ADF is the Abundance Discrepancy Factor given by the ratio of ionic abundances derived from optical recombination lines and collisionaly excitation lines of 6.9 (Liu et al., 2001). Due to this large ADF, this PN is a candidate to have a binary central star (Corradi et al., 2015). Its heliocentric distance calculated by Zhang (1995) is 6.99
kpc. The total nebular abundance ratios for this PN are N/O = 0.95 and He/H = 0.135 (Liu et al., 2001), therefore, this object is a helium and nitrogen rich nebula, classified as a Peimbert Type I PN.

Wray 16-411 is located at a distance of 6.14 kpc (Frew et al., 2016). Its heliocentric radial velocity is 50 km s (Beaulieu et al., 1999). Weidmann et al. (2016) mention that this object could be the galaxy IRAS 18232â4031 due to its spiral morphology, but the spectral data in the SPM Catalogue (López et al., 2012) show clearly that it is a planetary nebula.

Our models could help us to reproduce these PNe with the same elements used here, a precessing jet with time-dependent velocity or density ejection, and a broken torus, but using other values for , and . This work is in process.

## 6 Discussion and Conclusions

In the literature there are several point-symmetric PN that present regular knotty ejections such as Fg 1 and Hen 2-90.
Sahai &
Nyman (2000) analyzed the case of Hen 2-90 which shows knotty filaments with many ejections along a straight line at both sides of a toroidal structure, in this case the jet seems to be not precessing. Sahai &
Nyman (2000) estimated an ejection time interval of the knots between 100-120 yrs for this nebula. Another case already mentioned is Fg 1, where the S-shape is very extended (the filaments span 2 arcmin to either side of the PN) and knotty. Its kinematical structure was analyzed by López et al. (1993). It is found that the knots were not emitted as regularly as in Hen 2-90. In this case the jet is precessing and from the figures presented by López et al. (1993) we estimated that the ejection time interval between knots is approximately larger than 1000 yrs.

From the above it is apparent that the existence of knotty jets is relatively common in point-symmetric PNe. In many cases the knots are emitted in time intervals of several hundreds of years and present masses between 10 and 10 M.

The knotty filaments presented by point-symmetric PNe are related to the mass-loss processes
in the central stars. It is clear that some AGB and post-AGB stars eject material not continuously but in regular (or irregular) mass-loss events. In the case of point-symmetric PNe, the ejection corresponds to collimated jets moving in different (opposite) directions, sometimes precessing and giving origin to S-shape structures.

In addition to these point-symmetric objects, there is a fraction of PNe and proto-PNe showing external structures in the form of rings and arcs around their central stars, occasionally regularly spaced. By analyzing a large number of images searching for external rings and arcs, Ramos-Larios et al. (2016) reported that about 8% of PNe present this type of structures. These authors investigated the spacing between rings and arcs and their number in 29 objects. A few of the PNe presenting rings and arcs, also present filaments and point-symmetric components like NGC 6543, a PN with a very complex morphology (Balick et al., 2001; Miranda & Solf, 1992).

Ramos-Larios et al. (2016) estimated that the averaged time-lapse between rings and arcs, in the PNe of their sample, is in the range from 500 to 1200 yr, which is similar to the time interval between ejections in the knotty point-symmetric PNe. It is possible that the mechanisms of ejection of both type of structures are similar. However the expansion velocities of the collimated outflows are in general large (several tenths of km s and in some cases larger than 1000 km s) while the expansion velocities of rings and arcs (basically unknown) are assumed to be the typical expansion velocity for an AGB envelope, of about 10 km s. Knotty jets and rings are repetitive phenomena occurring in PNe, and they present similar time scales but they correspond to different epochs of ejections, rings seem to be produced at the end of the AGB phase, while periodic jets seem to be produced when the PN is already formed. The morphology of both types of structures is also very different, ring and arcs being spherically symmetric, while collimated jets are filamentary structures. Therefore, except for the time scales, both phenomena seem unrelated.

Comparing the point-symmetric PNe mentioned above, we found that many of them have a disk or torus in the centre and a precessing jet flowing through the poles. Both components
are necessary to produce this kind of morphology. The structure of the knotty collimated outflows could be due to an ejection with time-dependent velocity or time-dependent density.

In this work we have carried out several 3D hydrodynamical models with the code Yguazú, in order to explore and understand the processes that give origin to S-shaped planetary nebulae. It is assumed that such a dynamical structure is a consequence of the effect of a highly collimated precessing jet, as it was mentioned above. Therefore our hydrodynamical simulations consist of high-velocity collimated outflows in a medium given by the noisy wind of an AGB, simulating the density distribution of a dense torus. This jet changes its direction forming a cone of precession and it could be treated as a time-dependent velocity ejection or a time-dependent density ejection.

Our main results from these models can be summarised as: A precession jet is required to obtain such a S-shape morphology. The semi-aperture angle of the precession cone should be equal or larger than 30 for obtaining this structure and the point-symmetric morphology is very evident when the angle is 40. In addition, the factor of precession () and the number of ejections () in each turn will determine the structure of the S-shaped filaments.

Two models are presented here, the first one with several turns () and ejections (), and a second one with smaller number of turns () and ejections (). The torus included in the models helps to produce a bipolar nebula, but has no effect on the S-shape. Only the precessing jet produces such a shape. The projection of this structure on the sky plane can be compared to real observations.

The PV diagrams computed for the models, in comparison with observed PV diagrams of the PN H 1-67, allowed us to produce a more refined model for this compact nebula. In this work we show that H 1-67 presents a spectroscopic point-symmetric morphology, even when its image does not show this morphology. For this nebula we obtained high resolution spectra, locating the slit along the symmetry axis of the nebula (P.A. = 45), finding that they present two bright condensations in the centre corresponding to a broken torus and a S-shape filament with knots at high velocity. Also, when the slit is located with a P.A = 0, this S-like morphology is not evident. From the models described above, we calculated new models varying different parameters trying to reproduce such a structure.

Synthetic PV diagrams were generated, considering different inclinations of the models with respect to the line of sight and on the plane of the sky. From comparing these synthetic PV with the observed PV diagrams of H 1-67, we found that our Model 3 with time-dependent velocity ejections, semi aperture angle of 40, , , and Model 4 with a time-dependent density ejections and the same parameters, are the ones that best reproduce the characteristics of the PV diagrams for this object, after an integration time of 1500 and 1800 yr respectively. According to our models, each ejection is launched every 345 yr (model 3) and 305 yr (model 4).

Both models, 3 and 4, show two bright condensations with velocities close to 0 km s, which are associated with the toroidal material. In the orientation = 60, = 40 the models show the observed S-shape ending in two knots at 100 and 120 km s respectively, these PV diagrams reproduce our observations of H 1-67 with the slit at 45. In the orientation = 60, = 18, two knots are seen at velocity of the order of 100 and 120 km s respectively, and no S-shape is observed. These PV diagrams reproduce our observations of H 1-67 with the slit at 0.

By using different sets of parameters, these hydrodynamical models could be applied to reproduce other compact PNe that present PV diagrams showing a S-shape structure, such as Wray 16-411 and M 2-36.

In summary our hydrodynamical model is a simple one. It consists of a bipolar jet moving inside an anisotropic AGB wind (with high density at the equator). Such a model helped us to reproduce the PV diagrams of the spectroscopic bipolar PNe M 1-32 and M 3-15 (Rechy-García et al., 2017). In this work we show that adding precession and time-dependent ejection velocity or density we can reproduce the point-symmetric PV diagram of H 1-67.

The PN H1-67 analyzed in this work presents a point-symmetric morphology seen from spectroscopic observations but not in an optical image due to the collimated outflow is projected onto the torus which is much brighter. Possibly the point-symmetric structure of H 1-67 could be observed in a very high-resolution image with such as we show in our synthetic image (Fig. 16). In this paper, the idea that point-symmetric PNe results from the periodic ejection of dense material by a precessing jet is corroborated.

Some PNe such as Hen 2-11 (Jones et al., 2014), Hen 2-155 (Jones et al., 2015), NGC 6326 and NGC 6778 (Miszalski et al., 2011) appear to possess filaments and jet-like structures around a binary stellar system. Therefore the point-symmetric structure that presents H 1-67 could be a consequence of a binary system. This should be investigated in the future.

## Acknowledgments

We thank the daytime and night support staff at the OAN-SPM, Gustavo Melgoza, Salvador Monroy, and Felipe Montalvo, for facilitating and helping to obtain our observations. Helpful comments by Dr. José Alberto López, Dr. Michael G. Richer, Dr. Christophe Morisset and Dr. Gloria Delgado-Inglada are deeply acknowledged. We thank the anonymous referee for her/his very useful suggestions and comments, which help us to improve the previous version of this manuscript. We thank Enrique Palacios (Cómputo-ICN) for mantaining the Linux servers where the hydrodynamical simulations were carried out. This work received finantial support from DGAPA-PAPIIT grants IN103117, IG100218 and CONACyT grant 241732. J.S.R.-G. acknowledges scholarship from CONACyT-México.

## References

- Acker & Neiner (2003) Acker A., Neiner C., 2003, AAP, 403, 659
- Balick (1987) Balick B., 1987, AJ, 94, 671
- Balick et al. (2001) Balick B., Wilson J., Hajian A. R., 2001, AJ, 121, 354
- Beaulieu et al. (1999) Beaulieu S. F., Dopita M. A., Freeman K. C., 1999, ApJ, 515, 610
- Boffin et al. (2012) Boffin H. M. J., Miszalski B., Rauch T., Jones D., Corradi R. L. M., Napiwotzki R., Day-Jones A. C., Köppen J., 2012, Science, 338, 773
- Bond (2000) Bond H. E., 2000, in Kastner J. H., Soker N., Rappaport S., eds, Asymmetrical Planetary Nebulae II: From Origins to Microstructures Vol. 199 of Astronomical Society of the Pacific Conference Series, . p. 115
- Cliffe et al. (1995) Cliffe J. A., Frank A., Livio M., Jones T. W., 1995, ApJl, 447, L49
- Corradi et al. (2015) Corradi R. L. M., García-Rojas J., Jones D., Rodríguez-Gil P., 2015, ApJ, 803, 99
- De Marco (2009) De Marco O., 2009, PASP, 121, 316
- Escudero et al. (2004) Escudero A. V., Costa R. D. D., Maciel W. J., 2004, AAP, 414, 211
- Esquivel et al. (2003) Esquivel A., Lazarian A., Pogosyan D., Cho J., 2003, MNRAS, 342, 325
- Frew et al. (2016) Frew D. J., Parker Q. A., Bojičić I. S., 2016, MNRAS, 455, 1459
- García-Segura et al. (2002) García-Segura G., Franco J., López J. A., Langer N., Rózyczka M., 2002 Rev. Mex. Â´ Astron. Astrofis., Vol. 12. p. 117
- García-Segura & López (2000) García-Segura G., López J. A., 2000, ApJ, 544, 336
- Guerrero (2000) Guerrero M. A., 2000, in Kastner J. H., Soker N., Rappaport S., eds, Asymmetrical Planetary Nebulae II: From Origins to Microstructures, ASP Conference Series, Vol. 199. p. 371
- Guerrero et al. (2008) Guerrero M. A., Miranda L. F., Riera A., Velázquez P. F., Olguín L., Vázquez R., Chu Y.-H., Raga A., Benítez G., 2008, ApJ, 683, 272
- Haro-Corzo et al. (2009) Haro-Corzo S. A. R., Velázquez P. F., Raga A. C., Riera A., Kajdic P., 2009, ApJl, 703, L18
- Jones et al. (2014) Jones D., Boffin H. M. J., Miszalski B., Wesson R., Corradi R. L. M., Tyndall A. A., 2014, AAP, 562, A89
- Jones et al. (2015) Jones D., Boffin H. M. J., Rodríguez-Gil P., Wesson R., Corradi R. L. M., Miszalski B., Mohamed S., 2015, AAP, 580, A19
- Liu et al. (2001) Liu X.-W., Luo S.-G., Barlow M. J., Danziger I. J., Storey P. J., 2001, MNRAS, 327, 141
- Livio & Pringle (1996) Livio M., Pringle J. E., 1996, ApJl, 465, L55
- López et al. (1993) López J. A., Meaburn J., Palmer J. W., 1993, ApJl, 415, L135
- López et al. (2012) López J. A., Richer M. G., García-Díaz M. T., Clark D. M., Meaburn J., Riesgo H., Steffen W., Lloyd M., 2012, RMxAA, 48, 3
- Manick et al. (2015) Manick R., Miszalski B., McBride V., 2015, MNRAS, 448, 1789
- Meaburn et al. (2003) Meaburn J., López J. A., Gutiérrez L., Quiróz F., Murillo J. M., Valdéz J., Pedrayez M., 2003, RMxAA, 39, 185
- Mellema et al. (1991) Mellema G., Eulderink F., Icke V., 1991, AAP, 252, 718
- Miranda & Solf (1992) Miranda L. F., Solf J., 1992, AAP, 260, 397
- Miszalski et al. (2011) Miszalski B., Jones D., Rodríguez-Gil P., Boffin H. M. J., Corradi R. L. M., Santander-García M., 2011, AAP, 531, A158
- Morgan (1995) Morgan J. A., 1995, in Shaw R. A., Payne H. E., Hayes J. J. E., eds, Astronomical Data Analysis Software and Systems IV Vol. 77 of Astronomical Society of the Pacific Conference Series, . p. 129
- Parker et al. (2006) Parker Q. A., Acker A., Frew D. J., Hartley M., Peyaud A. E. J., Ochsenbein F., Phillipps S., Russeil D., Beaulieu S. F., Cohen M., Köppen J., Miszalski B., Morgan D. H., Morris R. A. H., Pierce M. J., Vaughan A. E., 2006, MNRAS, 373, 79
- Raga et al. (1993) Raga A. C., Canto J., Biro S., 1993, MNRAS, 260, 163
- Raga et al. (2009) Raga A. C., Esquivel A., Velázquez P. F., Cantó J., Haro-Corzo S., Riera A., Rodríguez-González A., 2009, ApJl, 707, L6
- Raga et al. (2000) Raga A. C., Navarro-González R., Villagrán-Muniz M., 2000, RMxAA, 36, 67
- Ramos-Larios et al. (2016) Ramos-Larios G., Santamaría E., Guerrero M. A., Marquez-Lugo R. A., Sabin L., Toalá J. A., 2016, MNRAS, 462, 610
- Rechy-García et al. (2017) Rechy-García J. S., Velázquez P. F., Peña M., Raga A. C., 2017, MNRAS, 464, 2318
- Rees & Zijlstra (2013) Rees B., Zijlstra A. A., 2013, MNRAS, 435, 975
- Richer et al. (2017) Richer M. G., Suárez G., López J. A., García Díaz M. T., 2017, AJ, 153, 140
- Rijkhorst et al. (2004) Rijkhorst E.-J., Icke V., Mellema G., 2004, in Meixner M., Kastner J. H., Balick B., Soker N., eds, Asymmetrical Planetary Nebulae III: Winds, Structure and the Thunderbird, ASP Conference Series, Vol. 313. p. 472
- Sabin et al. (2012) Sabin L., Vázquez R., López J. A., García-Díaz M. T., Ramos-Larios G., 2012, RMxAA, 48, 165
- Sahai et al. (2011) Sahai R., Morris M. R., Villar G. G., 2011, AJ, 141, 134
- Sahai & Nyman (2000) Sahai R., Nyman L.-Å., 2000, ApJL, 538, L145
- Schwarz et al. (1993) Schwarz H. E., Corradi R. L. M., Stanghellini L., 1993, in Weinberger R., Acker A., eds, Planetary Nebulae, IAU Symp. 155. p. 214
- Soker & Livio (1994) Soker N., Livio M., 1994, ApJ, 421, 219
- Soker & Rappaport (2000) Soker N., Rappaport S., 2000, ApJ, 538, 241
- Stanghellini et al. (1993) Stanghellini L., Corradi R. L. M., Schwarz H. E., 1993, AAP, 279, 521
- Stanghellini & Haywood (2010) Stanghellini L., Haywood M., 2010, ApJ, 714, 1096
- van Leer (1982) van Leer B., 1982, in Krause E., ed., Numerical Methods in Fluid Dynamics Vol. 170 of Lecture Notes in Physics, Berlin Springer Verlag, . pp 507–512
- Vázquez et al. (2008) Vázquez R., Miranda L. F., Olguín L., Ayala S., Torrelles J. M., Contreras M. E., Guillén P. F., 2008, AAP, 481, 107
- Velázquez et al. (2012) Velázquez P. F., Raga A. C., Riera A., Steffen W., Esquivel A., Cantó J., Haro-Corzo S., 2012, MNRAS, 419, 3529
- Weidmann et al. (2016) Weidmann W. A., Schmidt E. O., Vena Valdarenas R. R., Ahumada J. A., Volpe M. G., Mudrik A., 2016, AAP, 592, A103
- Zhang (1995) Zhang C. Y., 1995, ApJS, 98, 659

## Appendix A Distance of H1-67 from GAIA

During the revision of this manuscript The GAIA Data Release 2 (GAIA DR2) have been published and there a parallax of 1.71 mas 0.50 mas is reported for H1-67.
This would indicate a heliocentric distance between 450 and 830 pc (585 pc ).
This distance is 10 times smaller than the 5880 pc reported by Frew
et al. (2016) which represents a parallax of 0.17 mas.
On the other hand (Stanghellini &
Haywood (2010) report
a distance of 7860 pc for this object, equivalent to a parallax of 0.13 mas).
The observed angular diameter of H 1-67 is 5.6 arcsec, equivalent to 0.15 pc at a distance of 5880 pc.
With the distance of 585 pc from GAIA DR2, the physical diameter of H 1-67 would 0.015 pc,
and the age, estimated as would be only 150 to 200 yr instead of 1500 to 2200 yr derived from the distance by Frew
et al. (2016). With these parameters H1-67 would be the youngest PN known in the Galaxy.

Other parameter that is largely affected with this change of distance is the ionised mass.
Calculating the ionised mass as M= it is obtained 0.00619 M with GAIA DR2 distance, and if the distance is 5.88 kpc (Frew) the ionized mass is 0.15 M.

The mass can be also calculated from the H flux with the expression

it is obtained 0.001 M for a distance of 585 pc, while the mass for a distance of 5.88 kpc
is of 0.1 M, which is much more reasonable. We have assumed I(H) = 1.510 erg cm s from Escudero
et al. (2004), multiplied by 3 to take into account that the slit was 2 arcsec wide
and the object diameter is about 6 arcsec.

All seems to indicate that the parallax reported by GAIA DR2 is incorrect for this object, because a too small radius, age, and ionised mass is derived with such a small distance. The parallax from GAIA DR2 could be incorrect for H1-67, possibly because it is an extended object with knots and filaments. The central star is too faint and possibly it was no detected by GAIA.