A 3D Photoionization Model of the Extreme Planetary Nebula NGC 6302
We present a 3D photoionization model of the planetary nebula NGC 6302, one of the most complex and enigmatic objects of its kind. It’s highly bipolar geometry and dense massive disk, coupled with the very wide range of ions present, from neutral species up to Si, makes it one of the ultimate challenges to nebular photoionization modelling.
Our mocassin model is composed of an extremely dense geometrically thin circumstellar disk and a large pair of diffuse bipolar lobes, a combination which was necessary to reproduce the observed emission-line spectrum. The masses of these components, 2.2 M and 2.5 M respectively, gives a total nebular mass of 4.7 M, of which 1.8 M (39%) is ionized. Discrepancies between our model fit and the observations are attributed to complex density inhomogeneities in the nebula. The potential to resolve such discrepancies with more complex models is confirmed by exploring a range of models introducing small-scale structures. Compared to solar abundances helium is enhanced by 50%, carbon is slightly subsolar, oxygen is solar, and nitrogen is enhanced by a factor of 6. These all imply a significant 3rd dredge-up coupled with hot-bottom burning CN-cycle conversion of dredged-up carbon to nitrogen. Aluminium is also depleted by a factor of 100, consistent with depletion by dust grains.
The central star of NGC 6302 is partly obscured by the opaque circumstellar disk, which is seen almost edge-on, and as such its properties are not well constrained. However, emission from a number of high-ionization ‘coronal’ lines provides a strong constraint on the form of the high-energy ionizing flux. We model emission from the central star using a series of stellar model atmospheres, the properties of which are constrained from fits to the high-ionization nebular emission lines. Using a solar abundance stellar atmosphere we are unable to fit all of the observed line fluxes, but a substantially better fit was obtained using a 220,000 K hydrogen-deficient stellar atmosphere with and L. The H-deficient nature of the central star atmosphere suggests that it has undergone some sort of late thermal pulse, and fits to evolutionary tracks imply a central star mass of 0.73–0.82 M. Timescales for these evolutionary tracks suggest the object left the top of the asymptotic giant branch 2100 years ago, in good agreement with studies of the recent mass-loss event that formed one pair of the bipolar lobes. Based on the modelled nebular mass and central star mass we estimate the initial mass of the central star to be 5.5 M, in approximate agreement with that derived from evolutionary tracks.
keywords:planetary nebulae: individual: NGC 6302, ISM: abundances
Planetary nebulae (PNe) are one of the last evolutionary stages of the majority of low and intermediate mass stars (1-8 M). During the previous asymptotic giant branch (AGB) evolutionary phase the star ejects a large fraction of its mass during brief phases ( yrs) of high mass-loss (up to M yr). This mass-loss eventually exposes the central core of the star that radiates strongly in the UV and ionizes the surrounding nebula making it visible for study. Despite supposed spherically symmetric mass-loss while on the AGB, many PNe show point- or axi-symmetric shapes and complex morphologies that have been interpreted as evidence for the influence of a binary companion (e.g. De Marco, 2009). This and other outstanding issues in PNe research such as their small-scale structures (e.g. Gonçalves et al., 2001; Matsuura et al., 2009) and discrepancies between forbidden line and recombination line abundances for heavy elements (Liu et al., 2000; Tsamis et al., 2004; Wesson et al., 2005) suggest that out understanding of the late stages of stellar evolution is far from complete. Despite this, PNe are still regularly used as standard candles, metallicity indicators and tracers of stellar populations in distant galaxies (e.g. Buzzoni et al., 2006). If PNe are to be used in this way our understanding of their structure and evolution must be greatly improved, both through studies of large samples of such objects (e.g. Viironen et al., 2009; Sabin et al., 2010), their precursors the AGB stars (e.g. van der Veen & Habing, 1988; Zijlstra et al., 2001; Wright et al., 2008), and through individual case studies to understand the most complex examples (e.g. Gonçalves et al., 2006; Wareing et al., 2006).
NGC 6302 is one of the most complex and extreme examples of a planetary nebula. The nebula is broadly bipolar, with a second narrower pair of lobes visible further out (see Figure 1). It has a highly pinched waist that is typical of ‘butterfly’ bipolar PNe (Balick & Frank, 2002), as is a highly complex structure with many clumps and knots. The bipolar structure has been attributed to the presence of a very dense circumstellar disk that is seen almost edge-on and almost completely obscures the central star (, Matsuura et al., 2005; Szyszka et al., 2009), which has only recently been detected. The disk radiates strongly in the infrared (e.g. Lester & Dinerstein, 1984) with evidence for significant fractions of very large grains (Hoare et al., 1992) and crystalline silicates (Molster et al., 2001). NGC 6302 has been classified as a Type I PN, both according to the criteria or Peimbert & Torres-Peimbert (1983, N/O 0.5 and He/H 0.125) and the criteria of Kingsburgh & Barlow (1994) for a Milky-Way Type I PN, that N/O 0.8. The high helium and nitrogen abundances observed in the nebula imply a potentially massive progenitor, while observations of emission lines from very high ionization stages (e.g. Ashley & Hyland, 1988; Rowlands et al., 1994; Casassus et al., 2000) suggest an extremely hot central star.
Attempts to understand NGC 6302 have been hindered by its complexity and the inability to study its central star and therefore determine its evolutionary history. Aller et al. (1981) could find “no uniform density theoretical model that could reproduce the line intensities in NGC 6302 in any satisfactory fashion”, while Lame & Ferland (1991) could not produce a physically plausible model for the nebula without requiring shock-excitation, which Origlia et al. (1993) and Casassus et al. (2000) found strong evidence against. However, many of these models were forced to assume spherical symmetry due to the computational constraints of full 3-dimensional modeling, yet the inhomogeneous morphology of NGC 6302 clearly requires a 3-dimensional model if reliable conclusions are to be drawn.
In this paper we present a 3D photoionization model of NGC 6302 using the 3D Monte Carlo photoionization and radiative transfer code mocassin (Ercolano et al., 2003; Ercolano et al., 2005; Ercolano et al., 2008). The model presented here is a fully 3-dimensional model of NGC 6302 that includes both gas and dust, though we will only discuss the gas photoionization model in this paper. We will present the dust radiative transfer model in a future paper. In Section 2 we discuss the changes made to the mocassin code to enable the modeling of this object and outline the observations that our models are to reproduce. In Section 3 we describe the model used and the parameters used to build the model. In Section 4 we describe the modeling process and how different quantities varied according to changes in different parameters. In Section 5 we present the model results, including fits to the observed emission-line spectrum and the ionization and temperature structure of the nebula. Finally, in Section 6 we describe the resulting nebula properties and discuss the evolutionary history of the object.
2 The mocassin code
mocassin is a 3D photoionization and radiative transfer code that uses Monte Carlo methods to solve the radiative transfer in a photoionized nebula (Ercolano et al., 2003). The Monte Carlo approach allows it to fully treat all physical processes in 3D, self-consistently treating the diffuse component of the radiation field that 1D photoionization codes are unable to deal with. The code is highly versatile and can model any structure or morphology, introducing unlimited subgrids to resolve detail in areas of the model. It also allows multiple non-central ionizing sources and multiple gas chemistries. mocassin also has a complete radiative transfer dust model including the full treatment of all dust processes such as photoelectric heating and gas-grain collisions (Ercolano et al., 2005). The atomic database included opacity data from Verner et al. (1993) and Verner & Yakovlev (1995), energy levels, collisions strengths and transition probabilities from Version 5.2 of the CHIANTI111http://www.ukssdc.ac.uk/solar/chianti/ (Dere et al., 1997; Landi et al., 2006, and references therein) and the improved hydrogen and helium free-bound continuous emission data of Ercolano & Storey (2006). In this section we describe a number of changes made to the code as part of this work, as well as details of the computational facilities that the models were run on, and the observations that were to be reproduced by the models.
2.1 2-dimensional simulations in mocassin
mocassin is a 3-dimensional code, but many objects exhibit cylindrical symmetry along a rotation axis. Such objects may effectively be modeled in two dimensions, and , greatly reducing the computational and memory requirements of large and complicated models. We have developed and tested a 2-dimensional version of the mocassin code, which integrates the observations over from 0 to 360 (where is the angle between two directions in the plane). Since the nebula is symmetrical upon reflection in the plane, only positive values need be modeled, reducing the processing time further. This approximation still properly treats the diffuse component of the radiation field and thorough benchmarking has been performed that fully reproduces all of the mocassin benchmarks described in Ercolano et al. (2003). Although the models described here were run in 2D, we will often discuss the resulting 3D structure of the nebula, since the code reproduces the full 3D structure in the outputs. The final model was run in 3D to confirm the validity of the 2D modeling approach and no differences were found. The new routines are included in the latest public version of the code (Ercolano et al., 2008).
2.2 Computational details
The models presented in this work were run on the hiperspace facilities at University College London, including the Enigma SUN Microsystems V880 multiprocessor computer, the Keter SUN Microsystems V880 and V890 microprocessor computers and the Dell Legion Cluster, which uses Intel’s dual-core technology. Simulations were typically run with 8-16 processors and approximately 32 GB of memory. The necessary computational processing time varied depending on the complexity of the model and the required level of convergence of the Monte Carlo simulation, typically an iteration-to-iteration variation of the fraction of neutral hydrogen of 1% in 60% of the cells in the model.
2.3 Observations to reproduce
mocassin simulations produce several outputs that are used to compare the model results with observations. The outputs can also be adapted to the exact setup of the observations by simulating the slit or aperture through which the observations were made and the inclination angle of the nebula to the line of sight. As shown by Gonçalves et al. (2006) this is particularly important for long-slit observations of highly inhomogeneous nebulae such as NGC 7009 and the broad range of different measurements for NGC 6302 of the same emission line is testament to this. We have searched the literature for observations of NGC 6302 that span all wavelength regimes and that include lines from all the important elements and their ionization stages. The list in Table 1 was chosen based on spectral coverage and resolution. Where measurements overlap the most precise measurement was taken, or that which was accompanied by complementary information such as electron temperatures or densities derived from diagnostic line ratios.
|Region||Wavelength range||Telescope||Lines||Aperture / slit||Reference|
|Ultraviolet||1150-1975 Å||IUE SWP||7||10.3 23(oval)||Tsamis et al. (2003)|
|Ultraviolet||3426 Å||Mt. Lemmon 1.5-m||1||48(circular)||Rowlands et al. (1994)|
|Optical||3040-7400 Å||ESO 1.52-m||30||210 2||Tsamis et al. (2003)|
|Optical-IR||3300-8600 Å||Siding Springs 2.3-m||1||2(long-slit)||Groves et al. (2002)|
|Optical-IR||3500-20,000 Å||Cerro Tololo 92-cm||1||17 34||Danziger et al. (1973)|
|Near-infrared||2.7-4.8 m||UKIRT (echelle)||12||3 3||Casassus et al. (2000)|
|Mid-infrared||2.4-45 m||ISO SWS||14||14 20||Beintema & Pottasch (1999)|
|Mid-infrared||8-13 m||UKIRT N-band||5||4(circular)||Casassus et al. (2000)|
|Far-infrared||43-198 m||ISO LWS||7||40(circular)||Liu et al. (2001)|
The UV spectrum is mainly taken from the International Ultraviolet Explorer (IUE) spectrum presented by Tsamis et al. (2003) and complemented by observations of the [Ne v] 3426 Å line by Rowlands et al. (1994) for comparison with the [Ne v] lines in the mid-IR. The main optical emission-line spectrum was that of Tsamis et al. (2003), a deep optical spectrum made with a fixed slit aligned along the polar axis of the nebula that covers all aspects of the nebula. These observations include many of the diagnostic line ratios that will be important to diagnose the density and ionization structure within the nebula. This main line list is supplemented by lines from the far-red region using the published spectra of Danziger et al. (1973) and Groves et al. (2002). The infrared line list is split between the Infrared Space Observatory (ISO) observations using the Short Wavelength Spectrometer (SWS) presented by Beintema & Pottasch (1999) and the Long Wavelength Spectrometer (LWS) from Liu et al. (2001), and complemented by deep echelle spectroscopy from Casassus et al. (2000). These include many fine-structure lines in the mid- and far-infrared which have the advantage of only having a weak temperature dependence due to their low excitation energies.
In cases where the observed lines are known to be blends that are inseparable in the observations (e.g. the blend of [Ar iv] 4711 and [Ne iv] 4712), we accounted for the presence of the undesired component by calculating its intensity relative to other lines of the same species using the equib code.
Where possible all observations have been converted onto a scale relative to H for comparison with modeled line fluxes. For observations where this is not possible we compare observed fluxes, calculating model line fluxes using the distance to NGC 6302 of 1.17 kpc measured by Meaburn et al. (2008). We also compare the H flux for the entire nebula from our models with that measured by Tsamis et al. (2003) as F(H erg cm s. Using their derived reddening for NGC 6302 (c(H) we obtain a dereddened H flux of ergs cm s or L( ergs s, which can be compared with model outputs for the nebula.
2.4 Observational constraints on the central star luminosity
The central star of NGC 6302 was identified for the first time in deep Hubble Space Telescope images reported by Szyszka et al. (2009), who were able to perform photometric measurements of the central star in two filters and therefore make estimates of the surface temperature, extinction and stellar luminosity. However, at such high temperatures and extinctions, the exact solution is not well constrained and the authors were forced to use constraints on the age of the nebula, combined with evolutionary models, to fit the blackbody properties of the central star as (T,L) = (200,000 K, 2000 L). However, uncertainties on the photometry, evolutionary models, and extinction law led the authors to derive luminosities ranging from 860 – 24,000 L depending on these assumptions.
Another method for estimating the luminosity of the central ionizing source is based the sum of all emitted and observed radiation from the nebula. To perform this sum we have combined observations from a number of different observatories and wavelengths, the results of which are listed in Table 2. The IUE observations, which are for a arcsec aperture, were scaled up to the whole nebula using the intrinsic He ii 1640/4686 ratio (7 for K, cm, Case B, Osterbrock & Ferland, 2006) and the dereddened 4686 line flux from Tsamis et al. (2003). The UV – radio continuum emission is primarily free-free and bound-free emission, the former of which was estimated using the equations of Allen (1977), while the latter was calculated using the 2-photon continuum tool nebcont within dipso222http://star-www.rl.ac.uk/docs/sun50.htx/sun50.html. The infrared SED was taken from the ISO SWS and LWS spectra, with the ISO SWS spectrum was scaled up by a factor of 1.45 to match the LWS spectrum.
|Ultraviolet line fluxes||IUE||Feibelman (2001)||104|
|Optical line fluxes||ESO 1.52-m||Tsamis et al. (2003)||1324|
|Far-red line fluxes||Siding Springs 2.3-m||Groves et al. (2002)||73|
|Near-IR line fluxes||UKIRT||Casassus et al. (2000)||140|
|Free-free continuum||Theoretical||Allen (1977)||122|
|Infrared continuum||ISO||Molster et al. (2001)||3090|
The total luminosity of all components is 5700 L. This is likely to be a lower limit for the luminosity of the central star because the nebula and circumstellar disk are seen edge on and therefore a fraction of radiation from the central star may be escaping along the polar axis. Based on the geometry of the model nebula we have used, the opening angle of the circumstellar disk, and the reprocessing of radiation in the bipolar lobes, we estimate that the fraction of radiation escaping along the polar axis may be as much as a factor of two. We therefore estimate the luminosity of the central star to be L. This value is higher than previous estimates, which have either sampled only one component of the nebula (e.g. Matsuura et al., 2005, estimated a nebular luminosity of 3300 L based on infrared dust modeling) or have been estimated from evolutionary models (Szyszka et al., 2009, estimated 2000 L).
3 The model nebula
The final model described in this paper was built up over the course of this work starting with a simple pair of bipolar lobes and growing in complexity as both the models demanded and more observations were added to the list of observations to match. Every part of this model was necessary in some way to reproduce the observations described above. The principle components of the model are a large pair of bipolar lobes and a very dense disk-like distribution of circumstellar material. In this section we describe the different components of the NGC 6302 model and how this was constructed using the mocassin photoionization code.
Though observations show the presence of two pairs of bipolar lobes around NGC 6302 (e.g. Meaburn et al., 2005), we modeled only one pair for simplicity. Their similar orientations and the fact that most slits only cover the inner regions of the nebula suggest this will not significantly affect our results. The presence of large amounts of circumstellar material has been inferred from observations of a dark lane obscuring the central region of the nebula (Matsuura et al., 2005) and observations of molecular material in the core of NGC 6302 (Peretto et al., 2007). We attempted to model this material as both a circumstellar disk or a torus, the difference being that a torus has a scale height similar to its radial scale length, while a disk is intrinsically flat with a scale length significantly larger than its scale height. We found that the observations were better fit by a circumstellar disk (this will be discussed further in Section 6.1).
These structures were all modelled as cylindrically symmetric with the same axes as that of the bipolar lobes. This was a necessary simplification for this complicated nebula, but it should be noted that radio observations suggest that the disk is tilted with respect to the lobes, and most-likely warped (Gomez et al., 1989; Matsuura et al., 2005). All these components have a number of parameters that were initially estimated from previous observations of NGC 6302 but were allowed to vary to obtain the best fit during the modeling process. The final set of parameters used, and their relative uncertainties, are provided in Section 6.
3.1 The bipolar lobes
The bipolar lobes used in our model were developed from the model described by Dayal et al. (2000) for their model of MyCn 18. This geometrical model, shown in Figure 2, is based around a truncated paraboloid near the central star that develops into a cone-like structure further away from the star. The parameters indicated in Figure 2 are the latitude, , of transition from the paraboloid to the cone, the latitude, , of the end of the cone, the slope, , of the cone, the semi-height, , of the hourglass, and the waist-radius, , of the paraboloid. The shape of the walls of the hourglass can then be described using the equations
where the transition height, , is given by
and the radial distance from the hourglass center, , and the latitude, , of a point in the nebula are given by
The dimensions of the model nebula were estimated from images of NGC 6302 such as that in Figure 1 combined with an adopted inclination angle of (Meaburn et al., 2005) and distance of kpc (Meaburn et al., 2008). We estimated the height of the lobes to be cm, the latitudes of the parabola and cone as , and the slope of the cone as . These quantities were all constrained by images of the nebula and were not varied during the modeling process. Estimated measurement uncertainties are % for , % , and % for and . The waist radius, , was set to be the same value as the inner radius of the circumstellar disk and could be varied during the modeling process. Using these parameters the resulting nebula is quite wide (see Figure 3) as is seen in observations of NGC 6302. Finally, the density and density gradient in the lobes, as well as whether the lobes should be hollow (i.e. only edge-brightened) or filled was allowed to vary during the modeling process.
3.2 The circumstellar disk
To model the circumstellar disk we used the flared Keplerian disk model of Pascucci et al. (2004). In their model, the density distribution of the disk has the form
where is the radial distance and is the distance from the midplane. and are disk length and height parameters, which determine the scale length of the density distribution, and is a characteristic density for the disk. The radial density dependence of used by Pascucci et al. (2004) is also that found by Matsuura et al. (2005) from a best-fit to the spectral energy distribution (SED) in their radiative transfer models, though we allowed this to vary in our models. The flaring of the disk with distance from the central star is given by , utilising the flaring parameter , which controls the opening angle as a function of distance from the central star (see Figure 4). The disk extends from an inner radius, , to an outer radius, . All these parameters were allowed to vary in an attempt to find the best fitting model.
Initial parameters for the circumstellar disk were taken from the SED fitting of Matsuura et al. (2005) such as the height parameter, cm, and the inner radius, cm. For other parameters we used the initial estimates given by Pascucci et al. (2004). The outer radius was initially set as cm. The characteristic density, , was also allowed to vary during the modeling process.
3.3 The model grid
The bipolar lobes and circumstellar disk were mapped onto a 2D grid in the - plane, utilising the fact that the model nebula exhibits both cylindrical symmetry and is symmetric upon reflection in the plane (though, as noted above, the observed nebula is not cylindrically symmetric and is better described as point-symmetric). The number of cells in the grid and the size of the cells was varied during the modeling process to ensure that the ionization fronts and temperature and density gradients in the model were appropriately resolved such that the resulting emission line spectrum was a true representation of the model nebula. The final model consisted of a 2D grid of cells with the cell width varying from in the circumstellar disk, to in the outer regions of the bipolar lobes. This was then mapped onto a 3D grid of cells to confirm the validity of the 2D models and test the influence of density inhomogeneities in 3D (discussed below).
3.4 The central star
The central ionizing source of NGC 6302 was modeled using stellar model atmospheres produced using the Tübingen NLTE Model Atmosphere Package333Website: http://astro.uni-tuebingen.de/r̃auch/TMAP/TMAP.html. In the framework of the Virtual Observatory (VO, http://ww.ivoa.net), all spectral energy distributions from the TMAP model grids described here are available in VO compliant form from the VO service TheoSSA (http://vo.ari.uni-heidelberg.de/ssatr-0.01/TrSpectra.jsp?) provided by the German Astrophysical Virtual Observatory (GAVO, http://www.g-vo.org). (TMAP, Rauch & Deetjen, 2003; Werner et al., 2003). While it is sometimes true that the assumption of a blackbody is not worse than any other assumption in producing a photoionization model fit, at energies higher than 13.6 eV and 54.4 eV the differences between a blackbody flux and a stellar atmosphere can become so great that very large differences will result (e.g. Rauch, 1997, 2003; Armsdorfer et al., 2002). A grid of stellar atmosphere models were tested in our simulations, with effective temperatures ranging from 50,000 to 300,000 K and from 5-9. These were tested for both solar () and metal-poor (, ) elemental abundances (where the solar abundances used by TMAP are those from Asplund et al., 2005, and include all the elements up to nickel). Model atmosphere fluxes were also tried for ‘typical’ PG 1159 abundance ratios (He:C:N:O = 33:50:2:15 by mass), covering the same effective temperature and gravity ranges. Given the lack of knowledge about the central star of NGC 6302 we allowed all these properties to vary in our models.
3.5 Nebular abundances
No spatially resolved spectroscopic observations of NGC 6302 have been published, so there is no direct observational evidence for structures within the nebula with different elemental abundances. For this reason we used a homogeneous elemental abundance distribution in our models, despite mocassin’s ability to model regions with different chemistries. Our model used 14 elements, including all the principle elements that contribute to the heating and cooling of the nebula, as well as those responsible for the important density and temperature sensitive line ratios, and those required to reproduce the high-ionization coronal lines observed by Casassus et al. (2000). The initial abundances of He, C, N, O, Ne, S, Cl and Ar were taken from Tsamis et al. (2003), the abundances of Mg and Al were taken from the analysis of Casassus et al. (2000), while abundances for Na, Si and K were initially set to solar values (Lodders, 2003). All abundances were initially kept constant while structural parameters were being allowed to vary, but were then varied to obtain the best fit to the lines from each element (see final values in Table 5).
3.6 Dust modeling
NGC 6302 is known to have a very large dust component (e.g. Molster et al., 2001; Kemper et al., 2002) that will affect the radiative transport in the nebula. Any photoionization model of the nebula must therefore fully consider the influence of dust on the ionization structure of the gas. The models presented here therefore include a dust component that utilises mocassin’s full treatment of dust radiative transfer (e.g. Ercolano et al., 2005).The dust model will be fully discussed in a future paper, but we outline here the essential properties. Dust was included in the form of a number of dust species whose size distribution was adapted from the standard MRN model (Mathis et al., 1977). Optical constants were obtained from the literature and converted to absorption efficiencies using mocassin’s light scattering codes (that treat both spherical and elliptical dust grains). There is limited availability of UV and far-UV optical constants for common dust species in the literature, the most widely used being the amorphous silicates of Draine (2003), which were adopted in our model. Dust absorption within the nebula is significant in reducing the ionizing flux, with approximately half of all Lyman photons absorbed by the dust (in agreement with the observed fraction of luminosity coming from the IR SED, see Section 2.4).
4 The modeling process
In this section we describe the modeling process that was used to constrain the model properties by fitting the different observations. In particular we describe how different model parameters caused variations in the observable lines intensities or line ratios and how this led to the final model. These dependencies of lines or ratios on certain parameters were not absolute and changing other parameters could sometimes nullify the influence of certain parameters. This made the modeling process highly complicated and often required the parameter space around a model to be independently investigated at regular intervals in case new parameter dependencies had emerged due to changes in the model. Furthermore, in a model as complex as this there will inevitably be degeneracies in the models that can fit the observations, and some lines or ratios may remain not fitted. The objective here is to highlight which parts of the model are the most constrained by the observations and therefore represent the most important features of the model. The final model, while not a perfect fit, avoids many of the poor fits that certain parameters resulted in. We have divided this discussion, for the purposes of clarity, into the fundamental areas of the model (nebular structure, central star properties and nebula chemistry), but it should be noted that all the parameters were fitted simultaneously.
4.1 The nebular geometry and density
The geometry and density of the nebula, specifically its two main components, the circumstellar disk and the bipolar lobes, were refined by fitting the observed density- and temperature-sensitive line ratios. These line ratios respond to changes in the density and ionization structure of the nebula, but also highlight specific regions of the nebula according to their ionization potentials. The diagnostic line ratios can therefore be divided into those primarily responsive to changes in the low-density bipolar lobes, or the high-density circumstellar disk. These lists were not independent and many line ratios, particularly those sensitive to changes in the ionization structure of the nebula, were influenced by a large range of parameters. This was further complicated by the apertures and slits of the different observations of NGC 6302, each of which sampled a mixture of the two structural components. Despite this, many line ratios very highly responsive to specific parameters and by fitting these line ratios many structural parameters could be constrained. With these restrictions incorporated the variation of the remaining parameters was explored and the best-fitting model identified.
The circumstellar disk is modeled with seven parameters (, , , , , , ), though the characteristic density, , and the inner radius, , were the most influential in the resulting emission line spectrum and are therefore the most well constrained. This is because these parameters influence the region of the disk exposed to the full ionizing flux, while many of the other parameters affect regions of the disk beyond the ionization front. The line ratios that were responsive to changes in the circumstellar disk are a mixture of those that probe either the ionization structure, or regions of high temperature or density.
The most influential parameter was the characteristic density of the circumstellar disk, , which we fit to a value of cm, with an accuracy of approximately 20%. This parameter was influenced by nearly all the observed line ratios, though the strongest constraint came from the He ii / He i line ratio which was over-estimated by a factor 10 if cm, a disparity that could not be resolved by adjusting any other parameters. Many of the density-sensitive line ratios with high critical densities (i.e. those of [Cl iii], [Ar iv], and [K v] ) could also only be fit with a very dense circumstellar disk. Very little variation was observed in any of the line ratios when different density profiles were adopted, primarily because the circumstellar disk is highly optically thick. The initial choice of for the circumstellar disk radial density profile was made on the basis of the SED fitting of Matsuura et al. (2005), but or density profiles led to little difference.
The inner radius of the circumstellar disk, , was also very influential as it dictates the strength of the ionizing flux at the circumstellar disk. It was best constrained by a mixture of the [O ii] (3726 + 3729) / (7320 + 7330) and [Ar iv] density-sensitive line ratios. The former was over-estimated by a factor 5 if cm, and the latter increased significantly as was decreased. The final value of cm was chosen to balance these two ratios though it could be decreased if the [Ar iv] line ratio could be reduced by other methods. At a distance of 1.17 kpc this is equivalent to an angular radius of 0.68.
The thickness of the circumstellar disk (determined by ) influenced a large number of line ratios, most critically the [O ii] 3729 / 3726 ratio which became irretrievably over-estimated if the disk thickness was increased, and the [N ii] 5755 / (6584 + 6548) ratio that greatly increased if the thickness was significantly decreased. The best fit value of cm was constrained to an uncertainty of 50%.
The other disk parameters only influence material beyond the ionization front and therefore had a lot less influence on the modeled spectrum. The flaring parameter, , was only constrained to be by the [O iii] / [O ii] and [S ii] 4068 / (6731 + 6717) line ratios (which both increased significantly if the parameter was reduced), but no upper constraints were identified. It was therefore kept at its initial value, . The length parameter, , and the outer radius, , influenced a small number of line ratios, but none sufficiently enough to constrain them. The former parameter was therefore kept at its initial value of cm, while the latter parameter was adjusted based only on the results of radiative transfer dust-fitting (not discussed in this paper) and therefore we list its value as unconstrained here.
The circumstellar material was modeled here as a flat disk, but attempts were made to model it as a torus (i.e. a circumstellar distribution with a greater height and smaller length), as suggested by Matsuura et al. (2005) and Peretto et al. (2007). However, increasing the disk thickness to that suggested by Peretto et al. (2007) resulted in an unfittable [O ii] 3729 / 3726 ratio (as noted above). We therefore suggest that certainly in the inner regions where the ionized gas traces the structure of the circumstellar material, it is in the form of a thin disk. In the outer, neutral regions the circumstellar material may better represent a torus if traced by neutral material or dust.
The shape and dimensions of the bipolar lobes were adopted from images of the nebula such as that in Figure 1 with the distance to the nebula from Meaburn et al. (2008). These parameters were fixed throughout the modeling process and only the density of the lobes and the density profile were varied. The lobes were initially modeled with a low density of 2000 cm constrained by a number of line ratios that are particularly sensitive to density such as the [O ii] 3729 / 3726, [O iii] 4363 / 5007, [O iii] 52 / 88 m, and [S ii] (4068 + 4076) / (6731 + 6717) ratios (all sensitive to densities cm). However, none of these line ratios offered a particularly strong constraint on the final density, with the [O ii] 3729 / 3726 line ratio requiring densities in excess of 1000 cm, while the [O iii] 52 / 88 m line ratio increased if the density was significantly increased above 4000 cm. A final uncertainty of 50% is estimated for the final value of 2000 cm based on these constraints. No constraint could be found on the density profile of the bipolar lobes so they were kept at a constant density.
In response to an initially poor fit to a number of density-sensitive line ratios such as those of [Ar iv] and [K v], as well as the temperature-sensitive [N ii] line ratio, a third component was added to our model. This region, dubbed ‘the outflow’, is an inner region of the bipolar lobes with a higher density, , and with the same geometry but with an outer radius, , at which the density of the bipolar lobes dropped to the standard value. Introducing this component greatly reduced the [N ii], [Ar iv], and [K v], line ratios that had previously been significantly over-estimated. The final values of the two outflow parameters, cm and cm (11.5 at 1.17 kpc), are well constrained to 20%, due mainly to the first two line ratios, while the third offered less of a constraint. As an example of the modeling process Table 3 lists the line ratios most responsive to the introduction of, and changes to, the outflow component.
|Line ratio||Observed||Models with outflow density: (cm)|
|[N ii] 5755 / (6585+6546)||0.0290||0.133||0.0227||0.0231||0.0278|
|[O iii] (4959+5007) / [O ii] (3726+3729)||39.1||401||259||87.8||60.6|
|[S ii] 4068 / (6731 + 6717)||0.485||2.16||1.13||0.777||1.38|
|[Ar iv] 4741 / 4711||2.01||7.32||4.96||4.29||5.10|
|[K v] 4163 / 4123||1.07||8.38||4.61||4.32||4.12|
4.2 The central star and the IR coronal lines
Without any clear observations of the central star its parameters were particularly unconstrained at the beginning of the modeling process. While there had been many suggestions in the literature for a very high central star temperature (e.g. Ashley & Hyland, 1988; Casassus et al., 2000) we began the modeling process with a typical PN central star with a blackbody temperature of 90,000 K and solar abundances. While a large number of line ratios were responsive to changes in the central star temperature, the ionization structure diagnostics He ii / He i and [O iii] (4959 + 5007) / (3726 + 3729) were particularly responsive to it and would not come close to being fit unless the blackbody effective temperature was increased to at least 150,000 K.
At this temperature, with the peak of the SED in the far-UV, further changes produced little response in any line ratios. With the exception of the total H luminosity of the nebula, the stellar luminosity was also very unconstrained at these temperatures. The final stellar luminosity was therefore based on the best fit to the H line strength. For many PNe this may be a relatively unsatisfactory result, but given the high precision to which the distance to NGC 6302 is known, its absolute H luminosity and therefore the luminosity of its central source, is quite well constrained.
|(eV)||Solar composition||‘Typical’ H-def composition|
|[Mg v] 5.60 m||109||16.9||6.33||21.0||222||7.4||26.6||0.23|
|[Mg vii] 5.51 m||187||13.0||0.67||1.52||28.1||9.42||9.58||0.58|
|[Mg viii] 3.03 m||225||1.81||0.0||0.0||0.32||0.91||2.82||4.67|
|[Si vi] 1.96 m||167||23.5||14.1||21.7||429||31.6||73.4||1.34|
|[Si vii] 2.47 m||205||20.5||0.0||0.13||4.51||4.23||32.1||9.32|
|[Si ix] 3.93 m||303||0.030||0.0||0.0||0.0||0.002||1.04||139|
To constrain the stellar temperature further we used observations of a group of emission lines from high-ionization stage ions (the ‘infrared coronal lines’) that have been observed in the centre of NGC 6302 (e.g. Ashley & Hyland, 1988; Pottasch et al., 1996; Casassus et al., 2000). This group includes emission from ionization stages ranging up to Si and are unique to a very small number of PNe, with NGC 6302 showing the highest ionization species ever found in a PN. These were originally suspected to originate in high-velocity shocks in the nebula (e.g. Meaburn & Walsh, 1980; Lame & Ferland, 1991), but physical arguments against such levels of shock excitation (Oliva et al., 1996) and the lack of significantly broad wings on many of the lines (Casassus et al., 2000) argues for a photoionized origin.
The ionization potentials of many of these species lie in the extreme UV or X-ray region (see Table 4) and therefore place strong restrictions on the form of the high-energy ionizing flux from the central star. As discussed above, the ionizing flux from the central star was modeled using NLTE stellar model atmosphere fluxes for a range a range of central star abundances and surface temperatures. The effects of different temperatures and abundances on the fluxes of the six highest ionization stage ions is shown in Table 4. The effects of these changes on the lower ionization-stage nebula diagnostics was negligible.
Models using either a solar abundance or metal-poor stellar atmosphere were able to produce fits to many of the lines but were unable to reproduce the lines from the highest ionization stages, particularly the [Si ix] 3.93 m line. The strength of this line is underestimated by a factor even at K and an order-of-magnitude fit could not be obtained for any models with K. At such high stellar temperatures the lower ionization coronal lines (e.g. [Mg v] 5.60 m and [Si vi] 1.96 m) become overestimated by factors of several 100 or more. Models using hydrogen-deficient stellar atmospheres (Table 4) produced substantially better fits to these lines, with the best fit obtained using a 220,000 K stellar atmosphere. The fits to the magnesium lines presented in Table 4 are very good. The fits to the nebular silicon lines are slightly worse, with the 2.47 m line fitted well, but the 1.96 m line is overestimated by a factor of three and the 3.93 m line is overestimated by a factor of thirty.
Figure 5 shows the high-energy SEDs of 220,000 K solar-abundance and H-deficient model atmospheres. At the energies of the Mg and Si ionization potentials, there is a significantly higher flux in the H-deficient model atmospheres compared to the solar abundance model due to the lack of opacity from hydrogen. The absence of a number of low-abundance metals in these model atmospheres results in a reduced atmospheric opacity and therefore a higher flux at energies eV. Including them will likely reduce the photon flux at 330 eV and therefore resolve the problem of the over-predicted flux in the [Si ix] 3.93 m line in our models.
With the exception of a blackbody (which provides a better fit to the IR coronal lines, but is not an accurate representation of the ionizing flux, see e.g. Armsdorfer et al., 2002), the 220,000 K H-deficient stellar atmosphere provides the best fit to the observed IR coronal lines. Models with different surface gravities were tested but produced very minor changes in the model line fluxes. A value of produced the best fit, but is poorly constrained. The fit to the higher ionization stages of the silicon ions is not ideal, but it provides a closer fit than any of the solar composition models and is strong evidence that the central star of NGC 6302 is H-deficient.
4.3 The nebular chemistry
Once the main nebular and central star parameters were determined by fitting the density- and temperature-sensitive line ratios, attempts were made to achieve good fits to the strengths of individual emission lines by altering the chemical abundances of the nebula, which had originally been set to literature values for NGC 6302, primarily from Tsamis et al. (2003) and Casassus et al. (2000). The majority of abundances required small changes to fit the emission line fluxes, resulting in moderate variations in the nebula temperature structure which were then resolved by making other changes to the nebula structure. The sodium, magnesium, aluminium and silicon abundances were unchanged from their initial values either due to satisfactory fits, or insufficient numbers of lines to fully constrain their abundances. The final nebula abundances are listed in Table 5.
A full discussion of the dust modeling process will be reserved for a future paper, but as noted earlier dust was included in these models. This caused variations in a large number of line ratios due to significant changes in the nebular ionization structure.
5 Model results
In this section we present the best fitting model of NGC 6302 and compare the predicted emission-line fluxes from mocassin with those observed. The list of parameters used in the final model outlined in Section 3.1 of NGC 6302 are provided in Table 5. We also note how well constrained we estimate each parameter to be. The uncertainty on each parameter is estimated as the degree to which the parameter could vary without adversely affecting the model. For nebular abundances the uncertainty is based on the differences between observed and predicted line strengths for each element.
|Height of lobes,||cm||Estimated from images||10%|
|Transition latitude,||41.5||Estimated from images||25%|
|End of cone latitude,||49.6||Estimated from images||10%|
|Slope of cone,||66.0||Estimated from images||10%|
|Transition height,||cm||Estimated from images||25%|
|Waist radius,||cm||Model fit||20%|
|Lobe density||2000 cm||Model fit||50%|
|Density profile||constant||Model fit||-|
|Outflow density||20,000 cm||Model fit||20%|
|Outflow outer radius,||cm||Model fit||20%|
|Inner radius,||cm||Model fit||20%|
|Outer radius,||cm||Not well constrained||-|
|Length parameter,||cm||Model fit||50%|
|Height parameter,||cm||Model fit||20%|
|Flaring parameter,||1.125||Model fit||50%|
|Characteristic density,||80,000 cm||Model fit||10%|
|Central star properties:|
|Abundances||He:C:N:O = 33:50:2:15||Model fit||-|
|Effective temperature||220,000 K||Model fit||10%|
|Luminosity||14,300 L||Model fit||20%|
5.1 Fits to the emission line spectrum
|H/ erg cm s||776||852||T03||[Ne v] 3426||2.33||6.12||T03|
|H 4861||1.00||1.00||-||[Ne v] 14.3 m||634||1585||B99|
|He i 5876||0.175||0.144||T03||[Ne v] 24.3 m||308||434||B99|
|He i 4471||0.0592||0.0489||T03||[Ne vi] 7.64 m||315||1478||C00|
|He ii 4686||0.753||1.18||T03||[Na vi] 8.64 m||10.7||7.98||B99|
|[C ii] 157.7 m||13.7||22.7||L01||[Na vii] 4.69 m||4.00||2.13||C00|
|[N i] 5199||0.108||0.102||T03||[Mg v] 5.60 m||16.9||26.6||C00|
|[N ii] 5755||0.201||0.205||T03||[Mg vii] 5.51 m||13.0||9.58||C00|
|[N ii] 6548||1.74||2.19||T03||[Mg viii] 3.03 m||1.81||2.82||C00|
|[N ii] 6584||5.22||6.67||T03||[Al v] 2.88 m||0.0124||0.0742||C00|
|[N iii] 57.3 m||167||248||L01||[Al vi] 3.66 m||0.141||0.102||C00|
|N iv] 1486||5.34||5.86||T03||[Al viii] 3.69 m||0.16||0.0649||B99|
|N v 1240||13.82||19.0||T03||[Si ii] 34.6 m||17.4||325||B99|
|[O i] 5577||0.00232||0.00280||T03||[Si vi] 1.96 m||23.5||73.4||C00|
|[O i] 6300||0.243||0.338||T03||[Si vii] 2.47 m||20.5||32.1||C00|
|[O i] 63.2 m||284||135||L01||[Si ix] 3.93 m||0.030||1.74||C00|
|[O i] 145.5 m||10.5||8.27||L01||[S ii] 4068||0.145||0.218||T03|
|[O ii] 3726||0.300||0.158||T03||[S ii] 6716||0.101||0.0893||T03|
|[O ii] 3729||0.140||0.0517||T03||[S ii] 6731||0.198||0.191||T03|
|[O ii] 7320||0.0857||0.0866||T03||[S iii] 6312||0.0584||0.0328||T03|
|[O ii] 7330||0.0728||0.0713||T03||[S iii] 9533||0.717||0.666||B99|
|[O iii] 51.8 m||157||177||L01||[S iii] 18.7 m||57.7||55.3||B99|
|[O iii] 88.4 m||49||20.2||L01||[S iv] 10.5 m||44.0||18.9||C00|
|[O iii] 5007||12.9||13.8||T01||[Cl iii] 5517||0.00309||0.00391||T03|
|[O iii] 4959||4.29||4.62||T03||[Cl iii] 5537||0.00861||0.0127||T03|
|[O iii] 4363||0.396||0.316||T03||[Ar ii] 6.98 m||51.0||45.8||B99|
|[O iv] 25.9 m||322||472||B99||[Ar iii] 7751||0.0738||0.0474||G02|
|[Ne ii] 12.8 m||25.6||136||C00||[Ar iii] 7135||0.237||0.198||T03|
|[Ne iii] 15.6 m||377||186||B99||[Ar iv] 4741||0.209||0.191||T03|
|[Ne iii] 36.0 m||21.0||60.1||B99||[Ar iv] 4711||0.126||0.0445||T03|
|[Ne iii] 3967||0.276||0.338||T03||[Ar vi] 4.53 m||29.0||60.3||C00|
|[Ne iv] 2423||7.95||4.96||T03||[K v] 4123||0.00382||0.00150||T03|
|[Ne iv] 4724||0.0208||0.0177||T03||[K v] 4163||0.00407||0.00648||T03|
|[Ne iv] 4725||0.0182||0.0150||T03||[K vii] 3.19 m||0.279||0.323||C00|
Predicted emission line fluxes for our best-fitting model are given in Table 6, with strengths given relative to the intrinsic dereddened H flux where possible. The majority of the line strengths presented are in reasonable agreement with the observations, with fits to within 30%. In particular, nearly all the UV and optical line fluxes from Tsamis et al. (2003) are well reproduced. We obtain a reasonable fit to the N v 1240 Å resonance line which is often overestimated in photoionization modeling. Resonance lines must travel a much greater distance to escape from the nebula than do forbidden lines because of the large number of scatterings that they undergo and therefore can suffer greater attenuation due to dust. mocassin treats these lines separately and their absorption by dust is accounted for in this model (Ercolano et al., 2005).
The majority of infrared fine-structure line flux measurements are taken from spectra taken with the ISO SWS and LWS, presented by Beintema & Pottasch (1999) and Liu et al. (2001), respectively. The remaining infrared line flux measurements were presented by Casassus et al. (2000), from observations made with the UK Infrared Telescope (UKIRT). We note a significantly higher discrepancy between model and observations for line fluxes measured from ISO spectra than for line fluxes measured from the UKIRT spectra. This could be due to the combined effects of a misaligned pointing used for the ISO observations and the relatively small size of the ISO SWS aperture, as in the case of the NGC 3918 modeled by Ercolano et al. (2003). Beintema & Pottasch (1999) noted a significant discontinuity between short and long wavelength sections of the SWS exposures, which they attributed to misalignment of the pointing. Combined with the 3 pointing error on the arcsec ISO SWS aperture positioned over the central 10 diameter peak emission region suggests errors arising from such an offset could be high. As with all our predicted line fluxes, we have reproduced the relevant aperture or slit used for the observations and therefore any inaccuracies in the observational pointing would seriously affect a comparison of line fluxes, especially for an object like NGC 6302 with such strong density and temperature contrasts in its structure. The ISO LWS aperture is much larger, making any potential pointing error smaller, and we find a smaller discrepancy between our model results than for the SWS measurements, supporting this.
Many infrared lines from neutral species are predominantly emitted from the photo-dissociation regions (PDRs) of nebulae (e.g. Liu et al., 2001; Vasta et al., 2010) which are not treated by the current version of mocassin. Despite this Table 6 indicates that the [O i] 63 m and 145 m as well as the [C ii] 158 m lines show good agreement between the predicted and observed values, despite the fact that they originate from the neutral PDR. The current version of mocassin therefore appears to be matching these lines adequately. To complement this and for comparison with future studies we use our model to predict fluxes for two PDR lines without published fluxes, but which may be measured by future instruments. For the [C i] 369 and 609 m lines we derive fluxes for the entire nebula of and ergs cm s.
5.2 Fits to the density- and temperature-sensitive line ratios
|Line ratio||Observed||Model||Density (cm)||Reference|
|[C iii] 1907 / 1909||0.522||16000||B82|
|[O ii] 3729 / 3726||0.327||5000||G02|
|[O ii] (3726 + 3729) / (7320 + 7330)||1.33||5750||T03|
|[O iii] 52 / 88 m||8.76||1380||L01|
|[Ne v] 24.3 / 14.3 m||0.27||10000||B99|
|[S ii] 6731 / 6717||2.14||12900||T03|
|[Cl iii] 5537 / 5517||3.25||22450||T03|
|[Ar iv] 4741 / 4711||4.29||14900||T03|
|[K v] 4163 / 4122||4.32||10000||T03|
|Line ratio||Observed||Model||Temperature (K)||Reference|
|[N ii] 5755 / (6584 + 6548)||0.0231||14225||T03|
|[O i] 5577 / 6300||0.0083||10000||G02|
|[O iii] (4959 + 5007) / 4363||58.3||18400||T03|
|[O iii] 4363 / 5007||0.0229||17400||D73|
|[S ii] 4068 / (6731 + 6716)||0.777||10000||T03|
|[S iii] 6312 / 9533||0.0475||18300||B99|
|He i 5876 / 4471||2.94||T03|
|He ii 4686 / He i 5876||8.19||T03|
|[O iii] (4959 + 5007) / [O ii] (3726 + 3729)||87.8||T03|
Table 7 lists the model emission line ratios that are sensitive to changes in the density, temperature, or ionization structure of the nebula, and compares them to the observed line ratios of NGC 6302. Despite the wide range of electron densities and temperatures noted in the literature for NGC 6302 we have been able to reproduce the majority of emission line ratios using our model. The fits generally show good agreement to within 20%, with a small number only fit to within a factor of two. The most notably bad fit is the density-sensitive [K v] 4163 / 4122 line ratio, which is over-predicted by a factor of 4 and proved particularly difficult to fit.
The [Ar iv] and [K v] line ratios are both over-predicted by our model suggesting we have over-estimated the density of the circumstellar disk. However, many other line ratios constrained the density of the circumstellar disk and no suitable combination of parameters could be found using a lower-density circumstellar disk. Increasing the inner radius of the circumstellar disk offered a potential solution to fit these two line ratios, but increasing this parameter to cm resulted in extreme values of the [O ii] (3726 + 3729) / (7320 + 7330) line ratio. The introduction of an ‘outflow’ component (see Section 4.1) went a long way to resolving these discrepancies, but no simple model could be found to fit these lines accurately. Given this situation we were forced to leave these line ratios not fitted by our final model, but Section 5.2.1 discusses a possible solution to this problem. The poor fit to the He ii / He i line ratio is also worrisome, but again Section 5.2.1 presents a way to resolve this problem.
The majority of other line ratios are fit well and in some cases offer strong constraints on the geometry or density distribution of the nebula (see Section 4.1 for a discussion).
5.2.1 The effects of density inhomogeneities on the emission line spectrum
A number of the density- and temperature-sensitive line ratios discussed in Section 5.2, such as the He ii 4686 / He i 5876, [Ar iv 4741 / 4711, and [K v] 4163 / 4122 line ratios, are not well fit by our model. Exploring the parameter space in our model (see discussion in Section 4.1) revealed that almost every line ratio could be fit with some realistic combination of model parameters, but that no simple solution could be found that would fit all the line ratios. This does not necessarily suggest that a purely photoionization model solution does not exist for NGC 6302. A possible solution to this problem is that the morphology and density distribution of our model is an over-simplification of the true density distribution. The inhomogeneous appearance of the bipolar lobes of NGC 6302 (e.g. Figure 1) suggests that the density distribution is unlikely to be smooth, with many clumps, knots and small-scale structures visible (e.g. Matsuura et al., 2005, and Figure 1). Density enhancements can cause variations in density-sensitive line ratios as well as temperature-sensitive ratios though self-shielding. High densities can also induce collisional suppression of some collisionally-excited line, potentially invalidating the temperature sensitivity of some line ratios. The different densities determined from the [O iii] far-IR line ratio compared to the [Cl iii] and [Ar iv] line ratios, despite their similar ionization potentials (35.1eV compared to 23.8 and 40.7eV respectively), also suggests that moderate density inhomogeneities exist in the nebula and are responsible for quenching of the [O iii] far-IR lines by collisional de-excitation in the high-density regions (Liu et al., 2001).
Attempting to model this structure in detail would be impractical and largely unproductive compared to our attempts to characterise the overall structure of the nebula in Section 4.1. However, exploring the influence of such inhomogeneities will allow us to test the theory that they are responsible for the discrepancies between our model and the observations. To simulate these density inhomogeneities a number of ‘knots’ of various sizes and densities were simulated in our model using mocassin’s subgrid feature to resolve small-scale structures. For simplicity these knots were simulated as spherical over-densities in the structure of either the bipolar lobes or the circumstellar disk. Two parameters were used to model these structures, , a density factor representing the peak over-density of the knot, and , the radius of the knot. All the knots were modeled using a radial density profile. These parameters were allowed to vary over the ranges and , the latter of which was estimated from images of the nebula. The knots were distributed randomly across either the bipolar lobes, the circumstellar disk, or both, and the total number of knots was varied from 5–50.
The effect of these density enhancements on the nebula ionization structure depends on their position in the nebula and their density, but their main effect is to create shielded regions where lower ionization stages can exist. The most notable effects are therefore on the line ratios that probe the ionization structure of the nebula, particularly the He ii / He i line ratio which is significantly reduced if knots are introduced into the circumstellar disk. The other ionization-dependent line ratio, [O iii] (4959 + 5007) / [O ii] (3726 + 3729), was also found to decrease if knots were added to either the lobes or the disk, which had the effect of worsening the fit to this line ratio. However the magnitude of this change was smaller than that of the helium line ratio, suggesting it could be countered with other model changes.
Changes to density-sensitive line ratios were also observed when knots were added to the bipolar lobes, particularly the [Cl iii], [Ar iv], and [K v] line ratios that all probe high density regions, which were found to move closer to the observed ratios. The [O ii] and [S ii] line ratios, sensitive to lower-density regions, were less responsive. Adding knots to the circumstellar disk did not induce significant changes in these density-sensitive line ratios, most probably because either the upper critical densities of these line ratios are lower than the typical densities in the circumstellar disk or because the majority of the circumstellar disk is in a low ionization state.
In summary, introducing small-scale density inhomogeneities or ‘knots’ in the nebula structure could resolve some of the remaining discrepancies between our model and the observations. Knots in the circumstellar disk could resolve the over-predicted ionization-sensitive helium line ratio, while knots in the bipolar lobes could resolve differences in line ratios sensitive to high-densities. Despite this, a number of the remaining line ratios could not be matched within the parameter space explored (e.g. the low-density-sensitive [O ii] (3726 + 3729) / (7320 + 7330) line ratio).
5.3 The nebular ionization structure
|He / H||0.0645||0.0658|
|He / H||0.0477||0.0696|
|C / H||1.82(-5)||3.08(-5)|
|C / H||9.51(-5)||1.43(-5)|
|N / H||1.34(-4)||4.48(-5)|
|N / H||7.44(-5)||1.03(-4)|
|N / H||3.75(-5)||7.80(-5)|
|N / H||4.06(-5)||9.45(-5)|
|O / H||3.00(-4)||1.67(-4)|
|O / H||4.17(-5)||9.35(-5)|
|O / H||1.88(-5)||8.22(-5)|
|Ne / H||1.19(-5)||2.65(-5)|
|Ne / H||5.04(-5)||3.53(-5)|
|S / H||1.28(-6)||8.20(-7)|
|S / H||1.98(-6)||1.57(-6)|
|Cl / H||1.13(-7)||2.76(-8)|
|Ar / H||7.10(-7)||6.51(-7)|
|Ar / H||3.39(-7)||7.75(-7)|
A comparison of the relative ionic abundances in our model with those determined by Tsamis et al. (2003) is shown in Table 8. The ionization structure derived by those authors is in good agreement with our model, though there are a number of significant differences. The ionic fractions of O differ by a factor of four, while the fraction of Cl in our models is a factor of five larger than the empirical observational value. Fits to the [O iv] and [Cl iii] lines in Table 6 are generally good and don’t show major discrepancies. The likely explanation for this difference is that the model ionic abundances reported here are for the entire nebula, while those reported by Tsamis et al. (2003) only apply to the parts of the nebula caught under the slit.
Table 9 provides complete fractional ionic abundances for our NGC 6302 model for both the circumstellar disk and bipolar lobes. Figure 6 shows the ionic fractions for a number of elements along two sightlines that pass through the bipolar lobes (along the -axis) and the circumstellar disk (along the and axes). In the bipolar lobes nearly all ionization stages are seen somewhere, with a clear ionization sequence from the inner heavily ionized regions to the outer parts of the lobes where some material is neutral. Hydrogen is 73% ionized over the entire lobes, while helium is 51% singly ionized and 35% doubly ionized. In the circumstellar disk the ionization front occupies a thin region at the inner edge of the disk, almost one hundred times closer to the central star than in the bipolar lobes. Because of this the vast majority of the circumstellar disk is neutral with hydrogen 99% neutral in the disk. This is in agreement with the H images of Matsuura et al. (2005) that show very little emission from the disk of NGC 6302. Comparison between the [Ne v] images presented by Bohigas (1998) and our model Ne fractions also shows good agreement, with 13% of the neon in the bipolar lobes in the form of Ne, but less than 0.1% in the circumstellar disk.
Beyond the ionization front at a depth of cm the disk is optically thick and our model had difficulties converging due to the lack of photons. This is evident from the noise in the ionic fractions, predominantly of neutral species, beyond the ionization front. However, we found very little variation between Monte Carlo simulations in the total neutral and singly-ionized fractions of the majority of species, indicating that this noise did not affect our predicted emission line fluxes.
5.4 The nebular temperature structure
Mean electron temperatures weighted by ionic abundance have been calculated and are given in Table 10. The ionic temperatures for the same ions in different regions of the model show the influence of the density on the shielding of ions and therefore the temperature structure. In the optically thick circumstellar disk the neutral species all have temperatures 1000–3000 K, while those species in the optically thin bipolar lobes have temperatures around 7,000–10,000 K. The ionic temperatures increase towards higher ionization stages, but increase faster in the disk than in the lobes.
Because mocassin does not currently treat the neutral region fully, the temperatures of the neutral species likely represents the narrow transition region between ionized and neutral regions and not the highly optically thick and neutral region deep in the circumstellar disk. Beyond the ionization front, temperatures in the circumstellar disk drop to well below 1000 K and some cells have temperatures K. In these conditions complex molecules are likely to form, as observed by Peretto et al. (2007), but we are unable to model this with the current version of mocassin.
Electron temperature profiles are shown in Figure 6 as a function of radius along two axes. In the inner parts of the lobes the temperature reaches in excess of K (due to the lack of neutral H which would otherwise cause Ly- cooling) but drops to 10,000 K by a radial distance of cm. The temperatures in the disk do not reach such high levels, but do show the characteristic drops at the He, He, and H ionization fronts (the latter pair of which, at a radius of cm, have coincident ionization fronts as would be expected for such a hot central star). This range of temperatures is in full agreement with that found by many previous authors (e.g. Table 7), namely a nebula with a mean ion-weighted electron temperature of 19,400 K in the bipolar lobes, 3,800 K in the circumstellar disk and 12,000 K averaged over the entire nebula.
5.5 Comparison with observations of the central star
The central star of NGC 6302 was first detected in HST observations by Szyszka et al. (2009, and Figure 1), who attempted photometry in a number of narrow filters. Due to nebular contamination in many of the filters only two of these provide a reliable measurement of the stellar magnitude. Table 11 provides a comparison of their measured values with predictions from our model. To calculate the unreddened stellar magnitudes the K, model atmosphere was convolved with the relevant HST filter profiles and a spectrum of Vega. The resulting magnitudes are 6-10 magnitudes brighter than the observed photometry suggesting a considerable amount of extinction (e.g. Matsuura et al., 2005). Using an extinction law in the form given by Howarth (1983) we reddened the model photometry, varying the extinction to obtain agreement with the observed magnitudes. A very good fit was found using , fitting both magnitudes to within the observational uncertainties. This extinction is slightly larger than the found by Szyszka et al. (2009), but in agreement with other literature values (e.g. Matsuura et al., 2005, estimated , approximately ). The good agreement between observed and predicted photometric magnitudes and colours is strong verification of the central star luminosity and effective temperature we derived.
As illustrated by Figure 5 the high temperature and luminosity of the central star combined with its hydrogen-deficient nature results in a not-insignificant flux at UV and X-ray energies such that the central star of NGC 6302 should be a detectable X-ray source. This provides an observable test of both the high temperature and hydrogen deficiency of the central star. We have calculated the luminosity of the central star in the Chandra (0.5–8.0 keV) and XMM-Newton (0.3–4.5 keV) bands as erg s and erg s, respectively. X-ray absorption from neutral hydrogen in the circumstellar disk is calculated using the extinction of estimated above, and the conversion to X-ray absorbing cross-section from Ryter (1996). At a distance of 1.17 kpc this predicts absorbed fluxes of erg s cm for both X-ray bands (absorption due to neutral hydrogen is stronger at lower energies, off-setting the lower energy range of the XMM-Newton band for this soft X-ray source). Unfortunately the uncertainty on the predicted flux is relatively high due to uncertainties on the exact form of the stellar spectrum and the column of absorbing gas along the line of sight. However a H-rich central star of an equivalent effective temperature would be a significantly weaker X-ray source such that the detection of an X-ray point source in NGC 6302 would be sufficient to confirm the H-deficient nature of the central star. Since X-rays may also originate from shocked gas or collimated flows in the central cavity of the nebula (e.g. Kastner et al., 2008; Parkin et al., 2009), such observations must also distinguish between diffuse and point-source X-ray emission if the X-ray detection of the central star is to be observationally tested.
In this section we discuss the properties of the nebula and central star modeled here with particular focus on the parts of the model that are most strongly constrained by the observations. We discuss these in the context of previous studies of NGC 6302 and the general properties of planetary nebulae. Finally we suggest an evolutionary scenario for the formation and history of NGC 6302.
6.1 The nebula structure
NGC 6302 has long been known as a highly inhomogeneous nebula and our model matches this. Between the bipolar lobes and the inner edge of the circumstellar disk we model a density contrast of two orders of magnitude, from 2000 cm in the lobes to 300,000 cm in the circumstellar disk. This was necessary to reproduce many of the diagnostic emission line ratios, which were impossible to fit without including both the low and high density parts of the model. Such high densities are rare in PNe (e.g. Osterbrock & Ferland, 2006). In their study of PN densities Tsamis et al. (2003), NGC 6302 consistently displayed the highest densities from a range of diagnostic line ratios, with only NGC 5315 and IC 4191 showing comparably high densities from their sample.
|Volume of bipolar lobes||cm||20%|
|Mass of bipolar lobes||2.47 M||30%|
|Volume of circumstellar disk||cm||50%|
|Mass of circumstellar disk||2.22 M||50%|
|Total modeled nebular mass||4.69 M||40%|
|Total ionized mass||1.82 M||20%|
Derived total volumes and masses for our model nebula are listed in Table 12. The total ionized mass of our model is 1.82 M (predominantly located in the bipolar lobes), much higher than the typical ionized masses of Type I PNe (e.g. 0.54 M for six Magellanic Cloud Type I PNe, Barlow, 1987). The total mass of nebular material modeled (not including dust) is 4.7 M, also notably higher than the typical masses of PNe (Pottasch, 1996, studied 96 PNe and found masses in the range ), but not implausible. The final disk radius value of cm was constrained from fitting the infrared SED (not discussed here), but is unconstrained from photoionization modeling. The total disk mass is in good agreement with previous observations. Peretto et al. (2007) estimate the mass of the molecular region of NGC 6302 to be as much as M, in good agreement with our disk mass of 2.2 M, of which 97% (2.1 ) is neutral. However, Dinh-V-Trung et al. (2008) derive a significantly lower molecular mass of 0.12 M (scaled to a distance of 1.17 kpc), but note that this may be a lower limit if the CO temperature in the disk is higher than the 30 K they assume.
During the modeling process we attempted to model the circumstellar disk as a torus but found that the height of the circumstellar material was well constrained by the observations and that a flattened circumstellar disk was a better fit. Evidence for an extended disk-like morphology comes from Bohigas (1994) who observed differential extinction in the western lobe, which he suggested is viewed through the outer parts of the disk. Icke (2003) suggested that a warped disk structure is also necessary to form the multiple lobes seen in NGC 6302. A thin gas disk and a molecular torus could both be present in NGC 6302 if the torus of dust and molecular material were embedded within the outer regions of a sufficiently flared circumstellar disk. It is possible that any once torus-like structure would become thin and disk-like due to radiation pressure and stellar winds.
6.2 Nebular abundances
We fit our model of NGC 6302 with a helium abundance of 0.153, which is high even for a Type I PN (Kingsburgh & Barlow, 1994), indicative of a large amount of stellar processing in the progenitor star. According to the dredge-up models of Becker & Iben (1980), nebular material with this He abundance should correspond to a star with an initial mass M. The O/H abundance we derive for NGC 6302 () agrees within the uncertainty limits with the mean O/H ratio for Galactic disk PNe by Kingsburgh & Barlow (1994) and with the solar O/H ratios of Caffau et al. (2008) and Asplund et al. (2009). We do not confirm the depletion of oxygen by a factor of two that was found by Pottasch & Beintema (1999) and Tsamis et al. (2004) from empirical abundance analyses of collisionally excited UV, optical, and infrared lines. We attribute this difference to the neglect by empirical analyses of the large fraction of oxygen that resides in very high ionization stages in this unusually highly ionized nebula (see Table 9).
We find a C/O ratio of 0.43 for NGC 6302, indicating that the nebula is predominantly O-rich, in agreement with all other measurements of the C/O ratio in the literature (e.g. C/O, Aller et al. 1981, C/O, Pottasch & Beintema 1999, C/O, Casassus et al. 2000, and C/O, Tsamis et al. 2003), with the majority of values . Our derived C/H ratio of is 80% of the solar value listed by Asplund et al. (2009), indicating only a slight depletion of carbon, if any. However, the N/H ratio we derive for NGC 6302 () is a factor of 5.7 larger than the solar value given by Asplund et al. (2009), too large to have been produced by secondary conversion of initial carbon, while the combined (C+N+O)/H ratio for NGC 6302 of is 35% larger than the solar value of . We conclude that primary enrichment of nitrogen occurred in the precursor star of NGC 6302, via 3rd dredge-up enhancement of carbon, followed by hot bottom burning CN-cycle conversion of dredged-up carbon to nitrogen.
Several elements appear particularly depleted, such as magnesium (by a factor of two compared to the solar value) and aluminium (by a factor of one hundred). This is strong evidence for depletion onto dust grains, and given the evidence for large amounts of dust in NGC 6302 (e.g. Molster et al., 2001), is not unexpected.
6.3 The central star
The observation of emission lines from very high ionization stages and the nebula’s Type i status has long suggested that the central star of NGC 6302 was extremely hot. The hottest PN central stars directly measured to date have surface temperature of up to 180,000 K (Werner et al., 1997), which though high, are lower than all estimated temperatures for NGC 6302’s central star (e.g. Ashley & Hyland, 1988, estimated a temperature of 430,000 K).
Using a H-deficient stellar model atmosphere we find a best fit to the observed emission-line spectrum with central star properties = 220,000 K, , and L 14,300 L. This is the lowest temperature ever estimated for the central star of NGC 6302, due to our use of a H-deficient stellar atmosphere which has a larger flux at the energies necessary to significantly populate the highest ionization states from which emission lines are observed (see Figure 5). However, this luminosity is significantly larger than previous estimates, or our estimate of the total luminosity of 5700 L for NGC 6302, in Section 2.4, from integrating its dereddened observed emission across all wavelengths longwards of 1200 Å. That this is only 40% of the stellar luminosity needed by our bipolar model to match the observed line and continuum fluxes suggests that 60% of the stellar luminosity escapes through optically thin regions of the nebula.
The observed expansion of NGC 6302’s bipolar lobes (Meaburn et al., 2008; Szyszka et al., 2011) suggests it once had a fast wind that could affect the stellar spectrum. However strong mass-loss is usually only seen at lower effective temperatures and lower surface gravities, e.g. (e.g. Bohlin et al., 1982, find no evidence for stellar winds in the UV spectra of very hot PN central stars that Benedict et al. 2009 found to have surface gravities of ). As we will show in Section 6.4, evolutionary tracks for LTP and VLTP models imply a surface gravity of for the and we have derived, and this can also be derived from simple stellar relations by assuming a typical central star mass. Furthermore Casassus et al. (2000) could find no evidence for a hot wind or a wind blown cavity from observations of the IR coronal lines originating in the center of NGC 6302.
If the central star of NGC 6302 is H-deficient this implies that it has undergone some sort of late thermal pulse during its evolution, which has consumed the majority of the remaining surface hydrogen (e.g. Herwig et al., 1999). There are three possible scenarios for this: the ‘after final thermal pulse’ (AFTP) scenario, where an extra thermal pulse is experienced at the top of the AGB (Herwig, 2001); the ‘late thermal pulse’ (LTP) that occurs during the post-AGB evolution while H-burning is still active (Blöcker, 2001); or the ‘very late thermal pulse’ (VLTP) experienced by a hot white dwarf during its early cooling phase (Iben & Renzini, 1983). Only the VLTP scenario consumes the majority of the remaining hydrogen in the central star, while the other two scenarios suggest a dilution of the surface hydrogen to very low levels. We are currently unable to compare model atmospheres with H-deficient abundances with H-poor abundances, so we cannot use this method to separate these scenarios.
H-deficient central stars include the Wolf-Rayet central stars (WC, e.g. Hamann, 1997) and the hotter and rarer PG 1159 stars (Werner & Herwig, 2006). [WC]-type stars have been suggested to be the evolutionary precursors of PG 1159 stars (Górny & Tylenda, 2000), with the earlier stages characterised by strong stellar winds and interaction with the remaining nebulous material, and the latter stage exhibiting an often isolated central star at its hottest temperature. The high temperature and high surface gravity of the central star of NGC 6302 suggests it is similar to a PG 1159 star, and the presence of the remaining nebulous material can be explained if the star were particularly massive and has therefore evolved to this temperature before the nebula has dispersed.
There are many similarities between NGC 6302 and other PNe with H-deficient central stars that support this assertion. Zijlstra (2001) notes that nearly all PNe with [WC]-type central stars feature both O-rich and C-rich circumstellar material (i.e. dual-dust chemistry), as does NGC 6302 (Molster et al., 2001; Cohen, 2001). However, most H-deficient nebulae exhibit much stronger PAH (polycyclic aromatic hydrocarbon) emission than NGC 6302, though this could be due to the very low C/O ratio in NGC 6302 (e.g. Guzman-Ramirez et al., 2011, show that the production of PAHs is greatly reduced at low C/O ratios). Another similarity between NGC 6302 and PNe with H-poor or H-deficient central stars is that they both show highly inhomogeneous structures with many knots and filaments (see Figure 1 and Górny & Tylenda, 2000). Such knots may be associated with the complex turbulent velocities found by Acker et al. (2002) in the lobes of the majority of [WC]-type PNe. However the knots observed in NGC 6302 do not appear to be H-deficient (Szyszka et al., 2011), unlike those observed in other H-poor PNe (e.g. Abell 30 and Abell 78). This could indicate that the knots in NGC 6302 were ejected during an earlier, H-rich phase of mass-loss and have since been disrupted and accelerated by a later, H-poor mass-loss event.
6.4 Evolutionary tracks for late thermal pulse evolution
To further study the central star of NGC 6302 we have compared its properties determined from our models with those from evolutionary tracks for late thermal pulse models. This might allow us to determine the core mass of the central star, and therefore potentially its original mass, as well as its evolutionary state.
Figure 7 shows evolutionary tracks for the VLTP models of Blöcker (1995) and the LTP and VLTP models of Miller Bertolami & Althaus (2006). In the models of Blöcker (1995) the central star of NGC 6302 falls between the 0.61 and 0.84 M core mass tracks, at 0.82 M, corresponding to an initial mass of 4.8 M. In the Miller Bertolami & Althaus (2006) models the central star corresponds to a star with a core mass of 0.73 M, an initial mass of 3.7 M, and , close to our modeled value.
We can compare the timescales of these different evolutionary tracks with the known evolutionary history of NGC 6302, such as the last major mass-loss event 2200 yrs ago (Meaburn et al., 2008; Szyszka et al., 2011). The timescales for LTP events are considerably longer than those for VLTP events, which could suggest that the VLTP event is more likely given the short timescale since the last mass-loss event. Considering the VLTP model timescales, for the models of Blöcker (1995) for a 4.8 M initial mass star the position corresponds to an age of 2100 yrs since the star left the end of the AGB and 300 years since the VLTP event, suggesting that the eruptive event traced by Meaburn et al. (2008) was the last mass-loss event on the AGB before the central star evolved to hotter temperatures. For the models of Miller Bertolami & Althaus (2006) only ages since the VLTP event are provided, which, for a 3.7 M initial mass star correspond to 2000 yrs since the helium flash that caused this event, considerably longer than the equivalent evolutionary period in the Blöcker (1995) models.
If the timescales of Miller Bertolami & Althaus (2006) are correct this suggests that the eruptive event identified by Meaburn et al. (2008) is associated with a VLTP event. A VLTP mass-loss event however is not thought to be particularly large. van Hoof et al. (2007) estimated that the mass lost in the VLTP eruption of V334 Sgr might be as much as 0.01 M, at least an order of magnitude smaller than would be necessary for one pair of NGC 6302’s bipolar lobes. However, Althaus et al. (2008) found evidence that PG 1159 stars may have much thinner He-rich envelopes than previously thought and this might be due to considerable mass-loss during the VLTP event, though the helium layer is only 0.02 M thick and this puts an upper limit on what can be lost during such an event. Szyszka et al. (2011) note that the mass and momentum of the bipolar lobes is too high for them to have been ejected in a single mass-loss event and suggest that the recent mass-loss event may have accelerated gas that was previously ejected.
The large discrepancy in timescales between the models presented here for similar VLTP events implies a significant level of uncertainty in the current models that requires more accurate evolutionary tracks. Therefore at this point it is impossible to distinguish between the two sets of models as either could imply a timescale that agrees well with the evolutionary history of NGC 6302. Combining the derived core mass (0.74–0.82 M) with the nebular mass estimated from this work (4.7 M; Table 12) implies an initial stellar mass of at least 5.5 M. The initial stellar mass of 3.7–4.8 M from evolutionary tracks is not far from the 5.5 M that we derive.
6.4.1 Comparison with other known PG 1159 stars
A number of PG 1159 stars have measured effective temperatures and gravities (Werner & Herwig, 2006). Figure 8 shows these stars in the – plane (chosen because these quantities are better constrained than their luminosities), alongside the central star of NGC 6302 (using the value of fitted from Figure 7, which is better constrained than our modeled value), and the helium-deficient star H1504+65 (Werner et al., 2004). In addition we show evolutionary tracks for the Miller Bertolami & Althaus (2006) LTP/VLTP models. The position of the central star of NGC 6302 is in good agreement with that of many of the other PG 1159 stars. Is is clearly hotter than any of the known PG 1159 stars, and even hotter than the previously hottest known central star, H1504+65, though that source does appear to have had a higher initial mass and a previously higher temperature (it is now on the cooling track).
6.5 The evolutionary history of NGC 6302
From the evidence presented here we suggest that NGC 6302 was created by a star with an initial mass of at least 5 M, in agreement with both the chemical abundances and evolutionary tracks for the central star. The star may have had a lower-mass binary companion that diverted the AGB-phase mass-loss from isotropic to equatorial, forming the circumstellar disk (e.g. Bond & Livio, 1990; García-Arredondo & Frank, 2004).
The question of the age and stability of the circumstellar material is important for understanding the history of such a nebula. The presence of large dust grains (Hoare et al., 1992) and a high fraction of crystalline material (Molster et al., 2001) could indicate that the circumstellar disk is long-lived, providing a shielded and stable environment for these grains to build up. However, such dust characteristics could also be caused by high mass-loss rates and densities in the ejected material (Sylvester et al., 1999, found that the degree of crystallinity in OH/IR stars appears to scale with mass-loss rate). Peretto et al. (2007) measure an expansion velocity of 8 km s for the circumstellar torus and argue that it was ejected during a 5000 year mass-loss event that started 7500 years ago. Such a scenario requires mass-loss rates of M yr over this time period, potentially high enough to create the conditions for grain growth and crystallization in NGC 6302, but higher than current evolutionary models can sustain, even for shorter periods of time (e.g. Vassiliadis & Wood, 1993; Blöcker, 1995).
The kinematic age may be a lower limit if the disk has been accelerated by radiation pressure from the highly luminous central star. A simple calculation for the momentum imparted by radiation pressure on the circumstellar disk suggests that, based on the parameters derived here and the expansion velocity measured by Peretto et al. (2007), this would take 500,000 years at the current stellar luminosity (assuming that 10% of the stellar luminosity intercepts the circumstellar disk). This timescale is significantly longer than typical PNe lifetimes and as such unfeasible as the major acceleration mechanism for the circumstellar disk. A similar argument is presented by Peretto et al. (2007) who also conclude that alternative mechanisms are required and suggest interaction with a binary companion as a likely propulsion mechanism.
The circumstellar disk would have caused later mass-loss events to form the bipolar lobes. Matsuura et al. (2005) observed the circumstellar disk to be warped, which may be caused by an interaction with a binary companion, and according to the model of Icke (2003) would cause each pair of bipolar lobes to form with a different axis, as is observed. Hot-bottom burning at the base of the convective envelope (which occurs for stars more massive than 4 M, Iben & Renzini, 1983; Boothroyd et al., 1993), would convert carbon to nitrogen and prevent the star or the ejecta going through a prolonged C-rich phase.
At some point the star experienced some sort of late thermal pulse, ejecting a very small amount of material, causing the central star to become H-deficient and returning to the tip of the AGB. It is not currently possible to observationally distinguish between the different LTP scenarios, though the evolutionary timescales for the VLTP models of Miller Bertolami & Althaus (2006) show the strongest agreement with the known mass-loss history. Both sets of evolutionary tracks (above) suggest that the temperature of the star will continue to rise, potentially reaching a maximum temperature of 280,000 K, before entering the white dwarf cooling track again.
We have constructed a 3D photoionization model of the extreme planetary nebula NGC 6302 using the Monte Carlo photoionization and radiative transfer code mocassin. Our model reproduces the majority of emission line strengths and ratios, which place strong constraints on many of the stellar and nebular properties. The model consists of a very dense circumstellar disk where densities reach 300,000 cm at the inner edge, a large pair of bipolar lobes with a constant density of 2000 cm, and an intermediate component we have dubbed ‘the outflow’. This combination of components was required to match the majority of density- and temperature-sensitive line ratios from which a wide range of densities and temperatures have been observationally determined by previous authors. A number of line ratios remain not matched, which we suggest is due to complex density inhomogeneities throughout the nebula.
We derive a total nebular mass of 4.7 M, of which 1.8 M is ionized, almost entirely in the bipolar lobes. The C/O ratio for NGC 6302 is 0.43 indicating a predominantly O-rich nebula, in agreement with all other measurements in the literature. Carbon is found to be marginally under-abundant compared to the solar value, while nitrogen is significantly over-abundant. The combined (C+N+O)/H abundance is 35% larger than the solar value, all of which suggests that a 3rd dredge-up occurred in the precursor of NGC 6302, enriching the central star in carbon and nitrogen, followed by hot bottom burning CN-cycle conversion of carbon to nitrogen.
In modeling the central star we have incorporated NLTE model atmospheres to reproduce the ionizing flux. Fits to the optical emission-line spectrum imply an extremely high temperature for the central star, and we used the observed emission from high ionization-stage infrared ‘coronal’ lines to further constrain the form of the ionizing flux distribution. Using a solar abundance stellar atmosphere we were unable to fit all the observed line fluxes, overestimating many of the lower ionization stages while underestimating the higher stages. A substantially better fit was obtained using a hydrogen-deficient stellar atmosphere, implying that the central star of NGC 6302 is likely to be H-deficient.
Finally we compare the properties of the central star with evolutionary tracks for late thermal pulses that are capable of removing the majority of hydrogen from the central star. We find a good fit to a very late thermal pulse track for a star with an initial mass of 4–5 M. Timescales for this evolutionary model imply that the central star left the top of the AGB 2100 years ago, in good agreement with the age of one of the pairs of bipolar lobes determined by Meaburn et al. (2008) to be 2200 years. This mass is in reasonable agreement with the total modeled nebular mass plus central star mass of 5.5 M, taking into account contributions from the dust component. The mass is also in agreement with the chemical abundances in the nebula such as the high helium abundance, low C/O ratio and slightly enhanced (C+N+O)/H ratio.
A future paper will discuss the dust component of this model that used mocassin’s ability to fully model the radiative transfer of dust in an ionized nebula. The model includes the majority of observed amorphous and crystalline dust species, using multiple dust chemistry distributions, non-standard grain size distributions, and a new light scattering code for non-spherical dust grains.
The simulations presented in this work were performed in the DiRAC Facility jointly funded by STFC and the Large Facilities Captial Fund of BIS. The authors acknowledge support of the STFC funded Miracle Consortium (part of the DiRAC facility) in providing access to the UCL (University College London) Legion High Performance Computing Facility. The authors additionally acknowledge the support of UCL’s Research Computing team with the use of the Legion facility, particularly Jeremy Yates and Dugan Witherick. This work has also made us of the CHIANTI database, a collaborative project involving researchers at NRL (USA) RAL (UK), and the Universities of: Cambridge (UK), George Mason (USA), and Florence (Italy).
We are grateful to the referee, Albert Zijlstra, for a swift and helpful report that has improved the work presented here. We would also like to thank Pete Storey, Mikako Matsuura, Jonathan Rawlings, and Tim Gledhill for interesting discussions on this work. NJW acknowledges support from the STFC and an SAO Pre-doctoral Fellowship. BE acknowledges support from an STFC Advanced Fellowship. TR us supported by the German Aerospace Center (DLR) under grant 05 OR 0806.
- Acker et al. (2002) Acker A., Gesicki K., Grosdidier Y., Durand S., 2002, A&A, 384, 620
- Allen (1977) Allen K. W., 1977, Astrophysical quantities.. Cambridge University Press
- Aller et al. (1981) Aller L. H., Ross J. E., O’Mara B. J., Keyes C. D., 1981, MNRAS, 197, 95
- Althaus et al. (2008) Althaus L. G., Córsico A. H., Miller Bertolami M. M., García-Berro E., Kepler S. O., 2008, ApJL, 677, L35
- Armsdorfer et al. (2002) Armsdorfer B., Kimeswenger S., Rauch T., 2002, in Henney W. J., Franco J., Martos M., eds, Revista Mexicana de Astronomia y Astrofisica Conference Series Vol. 12, Effects of CSPN Models on PN Shell Modeling. pp 180–180
- Ashley & Hyland (1988) Ashley M. C. B., Hyland A. R., 1988, ApJ, 331, 532
- Asplund et al. (2005) Asplund M., Grevesse N., Sauval A. J., 2005, in Barnes III T. G., Bash F. N., eds, Cosmic Abundances as Records of Stellar Evolution and Nucleosynthesis Vol. 336 of ASP Conference Series, The Solar Chemical Composition.
- Asplund et al. (2009) Asplund M., Grevesse N., Sauval A. J., Scott P., 2009, ArA&A, 47, 481
- Balick & Frank (2002) Balick B., Frank A., 2002, ArA&A, 40, 439
- Barlow (1987) Barlow M. J., 1987, MNRAS, 227, 161
- Barral et al. (1982) Barral J. F., Canto J., Meaburn J., Walsh J. R., 1982, MNRAS, 199, 817
- Becker & Iben (1980) Becker S. A., Iben Jr. I., 1980, ApJ, 237, 111
- Beintema & Pottasch (1999) Beintema D. A., Pottasch S. R., 1999, A&A, 347, 942
- Benedict et al. (2009) Benedict G. F. et al., 2009, AJ, 138, 1969
- Blöcker (1995) Blöcker T., 1995, A&A, 299, 755
- Blöcker (2001) Blöcker T., 2001, ApSS, 275, 1
- Bohigas (1994) Bohigas J., 1994, A&A, 288, 617
- Bohigas (1998) Bohigas J., 1998, Revista Mexicana de Astronomia y Astrofisica, 34, 87
- Bohlin et al. (1982) Bohlin R. C., Harrington J. P., Stecher T. P., 1982, ApJ, 252, 635
- Bond & Livio (1990) Bond H. E., Livio M., 1990, ApJ, 355, 568
- Boothroyd et al. (1993) Boothroyd A. I., Sackmann I.-J., Ahern S. C., 1993, ApJ, 416, 762
- Buzzoni et al. (2006) Buzzoni A., Arnaboldi M., Corradi R. L. M., 2006, MNRAS, 368, 877
- Caffau et al. (2008) Caffau E., et al., 2008, A&A, 488, 1031
- Casassus et al. (2000) Casassus S., Roche P. F., Barlow M. J., 2000, MNRAS, 314, 657
- Cohen (2001) Cohen M., 2001, ApSS, 275, 103
- Danziger et al. (1973) Danziger I. J., Frogel J. A., Persson S. E., 1973, ApJL, 184, L29
- Dayal et al. (2000) Dayal A., Sahai R., Watson A. M., Trauger J. T., Burrows C. J., Stapelfeldt K. R., Gallagher III J. S., 2000, AJ, 119, 315
- De Marco (2009) De Marco O., 2009, PASP, 121, 316
- Dere et al. (1997) Dere K. P., Landi E., Mason H. E., Monsignori Fossi B. C., Young P. R., 1997, A&AS, 125, 149
- Dinh-V-Trung et al. (2008) Dinh-V-Trung Bujarrabal V., Castro-Carrizo A., Lim J., Kwok S., 2008, ApJ, 673, 934
- Draine (2003) Draine B. T., 2003, ApJ, 598, 1026
- Ercolano et al. (2005) Ercolano B., Barlow M. J., Storey P. J., 2005, MNRAS, 362, 1038
- Ercolano et al. (2003) Ercolano B., Barlow M. J., Storey P. J., Liu X.-W., 2003, MNRAS, 340, 1136
- Ercolano et al. (2003) Ercolano B., Morisset C., Barlow M. J., Storey P. J., Liu X.-W., 2003, MNRAS, 340, 1153
- Ercolano & Storey (2006) Ercolano B., Storey P. J., 2006, MNRAS, 372, 1875
- Ercolano et al. (2008) Ercolano B., Young P. R., Drake J. J., Raymond J. C., 2008, ApJS, 175, 534
- Feibelman (2001) Feibelman W. A., 2001, ApJ, 550, 785
- García-Arredondo & Frank (2004) García-Arredondo F., Frank A., 2004, ApJ, 600, 992
- Gomez et al. (1989) Gomez Y., Rodriguez L. F., Moran J. M., Garay G., 1989, ApJ, 345, 862
- Gonçalves et al. (2001) Gonçalves D. R., Corradi R. L. M., Mampaso A., 2001, ApJ, 547, 302
- Gonçalves et al. (2006) Gonçalves D. R., Ercolano B., Carnero A., Mampaso A., Corradi R. L. M., 2006, MNRAS, 365, 1039
- Górny & Tylenda (2000) Górny S. K., Tylenda R., 2000, A&A, 362, 1008
- Groves et al. (2002) Groves B., Dopita M. A., Williams R. E., Hua C.-T., 2002, Publications of the Astronomical Society of Australia, 19, 425
- Guzman-Ramirez et al. (2011) Guzman-Ramirez L., Zijlstra A. A., Níchuimín R., Gesicki K., Lagadec E., Millar T. J., Woods P. M., 2011, MNRAS, 414, 1667
- Hamann (1997) Hamann W.-R., 1997, in Habing H. J., Lamers H. J. G. L. M., eds, Planetary Nebulae Vol. 180 of IAU Symposium, Spectra of Wolf-Rayet type central stars and their analysis. pp 91
- Herwig (2001) Herwig F., 2001, ApSS, 275, 15
- Herwig et al. (1999) Herwig F., Blöcker T., Langer N., Driebe T., 1999, A&A, 349, L5
- Hoare et al. (1992) Hoare M. G., Roche P. F., Clegg R. E. S., 1992, MNRAS, 258, 257
- Howarth (1983) Howarth I. D., 1983, MNRAS, 203, 301
- Iben & Renzini (1983) Iben Jr. I., Renzini A., 1983, ArA&A, 21, 271
- Icke (2003) Icke V., 2003, A&A, 405, L11
- Kastner et al. (2008) Kastner J. H., Montez Jr. R., Balick B., De Marco O., 2008, ApJ, 672, 957
- Kemper et al. (2002) Kemper F., Molster F. J., Jäger C., Waters L. B. F. M., 2002, A&A, 394, 679
- Kingsburgh & Barlow (1994) Kingsburgh R. L., Barlow M. J., 1994, MNRAS, 271, 257
- Lame & Ferland (1991) Lame N. J., Ferland G. J., 1991, ApJ, 367, 208
- Landi et al. (2006) Landi E., Del Zanna G., Young P. R., Dere K. P., Mason H. E., Landini M., 2006, ApJS, 162, 261
- Lester & Dinerstein (1984) Lester D. F., Dinerstein H. L., 1984, ApJL, 281, L67
- Liu et al. (2001) Liu X.-W., et al., 2001, MNRAS, 323, 343
- Liu et al. (2000) Liu X.-W., Storey P. J., Barlow M. J., Danziger I. J., Cohen M., Bryce M., 2000, MNRAS, 312, 585
- Lodders (2003) Lodders K., 2003, ApJ, 591, 1220
- Loidl et al. (2001) Loidl R., Lançon A., Jørgensen U. G., 2001, A&A, 371, 1065
- Mathis et al. (1977) Mathis J. S., Rumpl W., Nordsieck K. H., 1977, ApJ, 217, 425
- Matsuura et al. (2009) Matsuura M., et al., 2009, ApJ, 700, 1067
- Matsuura et al. (2005) Matsuura M., Zijlstra A. A., Molster F. J., Waters L. B. F. M., Nomura H., Sahai R., Hoare M. G., 2005, MNRAS, 359, 383
- Meaburn et al. (2008) Meaburn J., Lloyd M., Vaytet N., López J. A., 2008, MNRAS
- Meaburn et al. (2005) Meaburn J., López J. A., Steffen W., Graham M. F., Holloway A. J., 2005, AJ, 130, 2303
- Meaburn & Walsh (1980) Meaburn J., Walsh J. R., 1980, MNRAS, 191, 5P
- Miller Bertolami & Althaus (2006) Miller Bertolami M. M., Althaus L. G., 2006, A&A, 454, 845
- Molster et al. (2001) Molster F. J., et al., 2001, A&A, 372, 165
- Oliva et al. (1996) Oliva E., Pasquali A., Reconditi M., 1996, A&A, 305, L21+
- Origlia et al. (1993) Origlia L., Moorwood A. F. M., Oliva E., 1993, A&A, 280, 536
- Osterbrock & Ferland (2006) Osterbrock D. E., Ferland G. J., 2006, Astrophysics of gaseous nebulae and active galactic nuclei, 2nd. ed, Sausalito, CA: University Science Books, 2006
- Parkin et al. (2009) Parkin E. R., Pittard J. M., Hoare M. G., Wright N. J., Drake J. J., 2009, MNRAS, 400, 629
- Pascucci et al. (2004) Pascucci I., et al., 2004, A&A, 417, 793
- Peimbert & Torres-Peimbert (1983) Peimbert M., Torres-Peimbert S., 1983, in D. R. Flower ed., Planetary Nebulae Vol. 103 of IAU Symposium, Type I planetary nebulae. pp 233–241
- Peretto et al. (2007) Peretto N., Fuller G., Zijlstra A., Patel N., 2007, A&A, 473, 207
- Pottasch (1996) Pottasch S. R., 1996, A&A, 307, 561
- Pottasch et al. (1996) Pottasch S. R., Beintema D., Dominguez-Rodriguez F. J., Schaeidt S., Valentijn E., Vandenbussche B., 1996, A&A, 315, L261
- Pottasch & Beintema (1999) Pottasch S. R., Beintema D. A., 1999, A&A, 347, 975
- Rauch (1997) Rauch T., 1997, A&A, 320, 237
- Rauch (2003) Rauch T., 2003, A&A, 403, 709
- Rauch & Deetjen (2003) Rauch T., Deetjen J. L., 2003, in I. Hubeny, D. Mihalas, & K. Werner ed., Stellar Atmosphere Modeling Vol. 288 of Astronomical Society of the Pacific Conference Series, Handling of Atomic Data. pp 103
- Rowlands et al. (1994) Rowlands N., Houck J. R., Herter T., 1994, ApJ, 427, 867
- Ryter (1996) Ryter C. E., 1996, ApSS, 236, 285
- Sabin et al. (2010) Sabin L., et al., 2010, PASA, 27, 166
- Sylvester et al. (1999) Sylvester R. J., Kemper F., Barlow M. J., de Jong T., Waters L. B. F. M., Tielens A. G. G. M., Omont A., 1999, A&A, 352, 587
- Szyszka et al. (2009) Szyszka C., Walsh J. R., Zijlstra A. A., Tsamis Y. G., 2009, ApJL, 707, L32
- Szyszka et al. (2011) Szyszka C., Zijlstra A. A., Walsh J. R., 2011, ArXiv e-prints
- Tsamis et al. (2003) Tsamis Y. G., Barlow M. J., Liu X.-W., Danziger I. J., Storey P. J., 2003, MNRAS, 345, 186
- Tsamis et al. (2004) Tsamis Y. G., Barlow M. J., Liu X.-W., Storey P. J., Danziger I. J., 2004, MNRAS, 353, 953
- van der Veen & Habing (1988) van der Veen W. E. C. J., Habing H. J., 1988, A&A, 194, 125
- van Hoof et al. (2007) van Hoof P. A. M., et al., 2007, A&A, 471, L9
- Vassiliadis & Wood (1993) Vassiliadis E., Wood P. R., 1993, ApJ, 413, 641
- Vasta et al. (2010) Vasta M., Barlow M. J., Viti S., Yates J. A., Bell T. A., 2010, MNRAS, 404, 1910
- Verner & Yakovlev (1995) Verner D. A., Yakovlev D. G., 1995, A&AS, 109, 125
- Verner et al. (1993) Verner D. A., Yakovlev D. G., Band I. M., Trzhaskovskaya M. B., 1993, Atomic Data and Nuclear Data Tables, 55, 233
- Viironen et al. (2009) Viironen K., et al., 2009, MNRAS, ?
- Wareing et al. (2006) Wareing C. J., et al., 2006, MNRAS, 366, 387
- Werner et al. (2003) Werner K., Deetjen J. L., Dreizler S., Nagel T., Rauch T., Schuh S. L., 2003, in Hubeny I., Mihalas D., Werner K., eds, Stellar Atmosphere Modeling Vol. 288 of Astronomical Society of the Pacific Conference Series, Model Photospheres with Accelerated Lambda Iteration. pp 31
- Werner et al. (1997) Werner K., Dreizler S., Heber U., Kappelmann N., Kruk J., Rauch T., Wolff B., 1997, in Schielicke R. E., ed., Reviews in Modern Astronomy Vol. 10, Ultraviolet Spectroscopy of Hot Compact Stars.. pp 219–252
- Werner & Herwig (2006) Werner K., Herwig F., 2006, PASP, 118, 183
- Werner et al. (2004) Werner K., Rauch T., Barstow M. A., Kruk J. W., 2004, A&A, 421, 1169
- Wesson et al. (2005) Wesson R., Liu X.-W., Barlow M. J., 2005, MNRAS, 362, 424
- Wright et al. (2008) Wright N. J., et al., 2008, MNRAS, 390, 929
- Zijlstra (2001) Zijlstra A. A., 2001, in Szczerba R., Górny S. K., eds, Astrophysics and Space Science Library Vol. 265, The Infrared [WC] Stars. pp 157
- Zijlstra et al. (2001) Zijlstra A. A., et al., 2001, MNRAS, 322, 280