Issue 
A&A
Volume 547, November 2012



Article Number  A29  
Number of page(s)  15  
Section  Interstellar and circumstellar matter  
DOI  https://doi.org/10.1051/00046361/201219515  
Published online  23 October 2012 
Discrepancies between the [O iii] and [S iii] temperatures in H ii regions
^{1} Département de Physique, de Génie Physique et d’Optique, Université Laval, Québec, QC, G1V 0A6, Canada
^{2} Instituto de Astronomía, Universidad Nacional Autónoma de México, D.F., Mexico
^{3} Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET), Argentina
^{4} Facultad de Ciencias Astronómicas y Geofísicas, Universidad Nacional de La Plata, Paseo del Bosque s/n, 1900 La Plata, Argentina
^{5} Departamento de Física Teórica, CXI, Universidad Autónoma de Madrid, 28049 Madrid, Spain
^{6} Research School of Astronomy and Astrophysics, Australian National University, Cotter Rd., Weston ACT 2611, Australia
^{7} Centro de Investigaciones de Astronomía, Av. Alberto Carnevalli, Mérida, Venezuela
^{8} Space Telescope Science Institute, Baltimore, Maryland 21218, USA
^{9} Instituto de Ciencias Nucleares, Universidad Nacional Autónoma de México, Ap. 70543, 04510 D.F., Mexico
Received: 30 April 2012
Accepted: 10 September 2012
Context. Analysis of published [O iii] and [S iii] temperatures measurements of emission line objects consisting of Hii galaxies, giant extragalactic Hii regions, Galactic Hii regions, and Hii regions from the Magellanic Clouds reveal that the [O iii] temperatures are higher than the corresponding values from [S iii] in most objects with gas metallicities in excess of 0.2 solar. For the coolest nebulae (the highest metallicities), the [O iii] temperature excess can reach ~3000 K.
Aims. We look for an explanation for these temperature differences and explore the parameter space of models with the aim of reproducing the observed trend of T_{O iii} > T_{S iii} in Hii regions with temperatures below 14 000 K.
Methods. Using standard photoionization models, we varied the ionization parameter, the hardness of the ionizing continuum, and the gas metallicities in order to characterize how models behave with respect to the observations. We introduced temperature inhomogeneities and varied their mean squared amplitude t^{2} . We explored the possibility of inhomogeneities in abundances by combining two models of widely different metallicity. We calculated models that consider the possibility of a nonMaxwellBoltzmann energy distribution (a κdistribution) for the electron energies. We also considered shock heating within the photoionized nebula.
Results. Simple photoionization calculations yield nearly equal [O iii] and [S iii] temperatures in the domain of interest. Hence these models fail to reproduce the [O iii] temperature excess. Models that consider temperature inhomogeneities, as measured by the mean squared amplitude t^{2} , also fail in the regime where T_{O iii} < 14 000 K. Three options remain that can reproduce the observed excess in T_{O iii} temperatures: (1) large metallicity inhomogeneities in the nebula; a (2) κdistribution for the electron energies; and (3) shock waves that propagate in the photoionized plasma at velocities ~60 km s^{1}.
Conclusions. The observed nebular temperatures are not reproduced by varying the input parameters in the pure photoionization case nor by assuming local temperature inhomogeneities. We find that (1) metallicity inhomogeneities of the nebular gas; (2) shock waves of velocities ≲ 60 km s^{1} propagating in a photoionized plasma; and (3) an electron energy distribution given by a κdistribution are successful in reproducing the observed excess in the [O iii] temperatures. However, shock models require proper 3D hydrodynamical simulations to become a fully developed alternative while models with metallicity inhomogeneities appear to fail in metalpoor nebulae, since they result in T_{rec}^{O++} ≳ T_{O iii} .
Key words: HII regions / ISM: lines and bands / shock waves / line: formation
© ESO, 2012
1. Introduction
The emission and absorbtion lines of starforming regions can inform us about the physical conditions of the gaseous media, such as abundances, temperature, and ionization degree. They also provide us with estimates of the ages, masses, and composition of the stellar population and the properties of the ionizing stellar clusters (see e.g. Hägele et al. 2011, hereafter H11; and references therein). The electron temperature is a property crucial for interpreting emission line intensities. When we investigated the [S iii] temperatures in the literature, the study of Hägele et al. (2006, hereafter H06), who compared the temperatures inferred from collisionally excited lines of [S iii] with those of [O iii], caught our attention. Reproducing these with photoionization models was set by some of us (RM and LB) as a potentially instructive exploratory exercise that turned out more challenging as we proceeded with our analysis. Our paper reports on this exploration and is a work in progress, mainly because the uncertainties in the data are substantial, as discussed in Sect. 2. In essence, in nebulae where measurement errors are the smallest, we find that the T_{O iii} temperatures exceed the values derived for T_{S iii} by up to ~3000 K, a trend that we were not able to reproduce using simple photoionization calculations. We consider it remarkable (yet with insight predictable) that the relationship between these two temperatures is quite insensitive to strong changes in the input parameters of photoionization models. We subsequently explored models with (1) metallicity inhomogeneities; (2) temperature inhomogeneities; (3) nonMaxwellBoltzmann energy distributions (hereafter a κdistribution), and completed our investigation with models that (4) combine photoionization with shock heating. We hope that the current work may motivate observers to undertake observational studies that are able to substantially reduce errors in the measurements of the nebular and auroral [S iii] lines.
Considering the extra difficulties in measuring the [S iii] lines in the infrared and the auroral counterpart at 6312 Å, what alternative species can be used to infer the electron temperature from an auroraltonebular line ratio? There are the two wellknown [N ii] or [O ii] temperature diagnostics, but they carry their own set of uncertainties, which renders an interpretation of these diagnostics ambiguous and more model dependent. For instance, the electron temperatures deduced from the [N ii] 5755 Å/(6548 Å + 6583 Å) ratio are frequently higher than those of the corresponding [O iii] lines (e.g. Tsamis et al. 2003). According to Stasinska (1980) and Copetti (2006), the increase in temperature outward as a function of nebular radius is likely caused by a combination of hardening of the radiation field with increasing optical depth, and stronger cooling from finestructure lines of [O iii] in the inner parts of the nebulae where O^{+2} dominates (see Stasinska 1980). This expectation is not fully confirmed, however, due to the significant contribution of recombination to the excitation of the auroral [N ii]λ5755 line (Rubin 1986; Tsamis et al. 2003) coupled with the potential presence of highdensity inclusions in nebulae (Viegas & Clegg 1994), all of which indicates that [N ii] fails to qualify as a suitable alternative. The [O ii] nebulartoauroral line ratio is even more susceptible to contribution from recombination and to the potential presence of highdensity inclusions in nebulae (Viegas & Clegg 1994), which may cause deceptively high temperatures to be derived from the [O ii] 3727 Å/(7320 Å + 7330 Å) ratio. Furthermore, this line ratio is more sensitive to variations in the extinction curve (R_{V}) or to secondorder effects from an inhomogeneous dust distribution, which standard reddening corrections cannot correct for. In comparison to N^{+} and O^{+}, the specie S^{+2} presents the clear advantage of strongly overlapping O^{+2} in the higher excitation region of a nebula. An appealing alternative might be the [Ar iii] auroraltonebular line ratio 5192 Å/(7136 Å + 7751 Å) (Keenan et al. 1988; Copetti 2006), but there are very little data available for these lines owing to their intrinsic weakness.
Studying physical properties of emission line objects at different redshifts with accurate and reliable methods is a fundamental issue when one aims to detect possible evolutionary effects such as systematic differences in the chemical composition of these objects. We discuss in Sect. 4 that the observation that the T_{O iii} temperatures are higher than T_{S iii} may be connected to gaseous nebula phenomena, which are not yet fully understood. One of these phenomena is the abundance discrepancy factor (ADF), which is the ratio of the abundance derived for a heavyelement ion from its optical recombination lines (ORLs) to the abundance derived for the same ion from its ultraviolet, optical, or infrared collisionally excited lines (CELs; Tsamis 2004). In Hii regions, the ADF is about 2 and has been the subject of various interpretations (Peimbert et al. 1995; Esteban et al. 2002; Tsamis & Péquignot 2005; PérezMontero & Díaz 2007; Tsamis et al. 2003, 2011; Nicholls et al. 2012; LópezSánchez et al. 2012, and references therein).
For the present study we carried out several models using the multipurpose photoionizationshock code mappings ie that allows us to explore five options: (1) simple photoionization models where we varied the main input parameters; (2) models that consider the effect of temperature inhomogeneities^{1}; (3) a combination of models of widely different metallicities; (4) models that consider κdistributions of electron energies (nonMaxwellBoltzmann energy distributions); and (5) models that combine shock excitation with photoionization. The article is structured as follows: Sect. 2 contains the data analysis, our modeling strategy is described in Sect. 3, and the comparison of models with the data is discussed in Sect. 4. Brief conclusions are presented in Sect. 5.
Fig. 1 Relation between T_{S iii} and T_{O iii} temperatures inferred from the [S iii] and [O iii] auroraltonebular line ratios for the objects studied by PM06, including the three Hii galaxies later observed by H06. The light gray dashed line depicts the locus of equal temperatures (T_{S iii} = T_{O iii}), while the solid line indicates the errorweighted leastsquares fit to the data (see Sect. 2). 
2. The data base
The analysis carried out in this paper relies on two temperaturesensitive line ratios that involve auroral and nebular line transitions of the same ion. These are [O iii] 4363 Å/(4959 Å + 5007 Å) and [S iii] 6312 Å/(9069 Å + 9532 Å). The inferred temperatures assuming the lowdensity limit hereafter are labeled T_{O iii} and T_{S iii}. Our sample corresponds to a compilation of data (99 objects) initially published by PérezMontero et al. (2006: hereafter PM06) and subsequently augmented by three Hii galaxies observed by H06. In total, the sample consists of the spectra of 42 Hii galaxies, 34 giant extragalactic Hii regions (GEHRs), 14 Galactic Hii regions, and 12 Hii regions from the Magellanic Clouds. Figure 1 shows the temperatures inferred from a total of 102 pairs of [O iii] and [S iii] line ratio measurements as well as the associated error bars. For clarity’s sake, the error bars are not drawn in the figures. Hereafter, we will use the terms Hiiregion or nebula as generic expressions representative of objects from the sample.
As is evident from the distribution, a correlation between the two temperatures is observed in Fig. 1. An errorweighted leastsquares fit of the data results in the following linear regression: (1)with a correlation coefficient of 0.85. The light gray dashed line depicts what would be the locus of equal temperature values in this plot and in subsequent figures. The correlation is markedly steeper than that inferred by Kennicutt et al. (2003) of T_{S iii} = 0.83 × T_{O iii} − 1700 K from a study of 20 Hii regions in the spiral galaxy M101 (these objects are included in our much larger sample). We note a systematic inequality between the two temperatures in Fig. 1. In particular, lower temperature objects T_{S iii} < 14 000 K tend to lie below the equal temperature dashed line. Towards the high temperature end (T_{S iii} > 14 000 K), however, the error bars and the overall dispersion of the different measurements become quite large. Due to this and to the absence of any clear trend in the data, our analysis will disregard the high temperature upper right quadrant in Fig. 1.
Given the distribution of the error bars and the increasing dispersion of the points toward higher temperatures, it is doubtful that a single linear regression provides an appropriate description of the current data set. We suggest that a bimodal distribution might be a better description of the data: 1) the lower left quadrant, which shows measurements that are correlated and possess smaller error bars; and 2) the upper right quadrant, where no clear trend is apparent and where the error bars are noticeably larger. Our study and the modeling presented below focus on the temperature behavior within this lower left quadrant and models the temperature departure from the dashed line in Fig. 1, that is, the domain T_{O iii} > T_{S iii} below T_{O iii} < 14 000 K. Even if the error bars had been underestimated by PM06, a proportionally higher noise amplitude would not erase the general trend that the [O iii] temperatures lie above those of [S iii]. Could this trend be the result of a bias affecting the weaker lines? Rola & Pelat (1994) established that the “measured” intensities in poor signaltonoise (S/N) line measurements are strongly biased toward overestimating their intrinsic values. The errors reported by PM06 for the vast majority of objects are markedly larger for T_{S iii} than for T_{O iii} . Because auroral lines are much weaker and the error in subtracting the underlying continuum larger, the errors pertaining to the weaker auroral lines are larger. The statistical bias studied by Rola & Pelat (1994) would in our case lead to an overestimate of [S iii]λ6312, hence of T_{S iii} , which is the opposite of the observed trend. Telluric absorption of the infrared [S iii] 9069 Å and 9532 Å lines or contamination of the 6312 Å line by atmospheric emission^{2} from nearby [O i] 6300 Å would both contribute to enhancing the inferred T_{S iii} values, in contrast to the observed trend. In our view, the spread in temperatures in the lower left quadrant is genuine, but we cannot rule out that unforeseen systematic errors are contributing to it.
Improving the accuracy of temperature measurements will be required before definitive conclusions can be drawn. The infrared [S iii] lines are intrinsically more difficult to measure. The risk of systematic errors caused by atmospheric telluric emission, which increases toward longer wavelengths near 9500 Å and beyond, remains a concern. As databases will grow, a significant improvement in the quality of the data should be the first priority. For instance, it might become preferable to discard data for which the 9532 Å/9069 Å ratio does not agree with the theoretical value of 2.46 (from Podobedova et al. 2009) within a reasonable S/N margin.
3. Model calculations
An abundance analysis of the complete sample has already been carried out by H06 and PM06, in which the authors took into account the measured nebular temperatures as well as the intensities of the forbidden lines with respect to the recombination lines. Such an abundance analysis will not be repeated in this paper. Rather, our aim will be to investigate whether one can model [O iii] temperatures that are higher than those derived from [S iii], as is observed in the lower left quadrant of Fig. 1. It follows that our study will focus on nebulae that tend to lie on the highmetallicity end of the data points. As an exercise, we can estimate the range of abundances that the objects in this quadrant cover. To do so, we can use the abundance calibration presented in Appendix A. It was derived from the data used by H06 and consists of a leastsquares fit of the observed [O iii] temperatures as a function of the deduced O/H abundance ratio. When adopting the value of T_{S iii} = T_{O iii} = 14 000 K as a marker that defines the lower left quadrant, this corresponds to Z ≈ 0.15 Z_{⊙} (see Appendix A). If we select instead the temperature at which the linear regression crosses the dashed line, which is at T_{O iii} = 12 400 K, the inferred metallicity is Z ≈ 0.25 Z_{⊙}. In summary, the nebulae whose temperatures we attempt to model have metallicities Z ≳ 0.2 Z_{⊙}.
In standard (steadystate) photoionization models, we disregard hydrodynamical effects and assume that equilibrium ionization and equilibrium temperature are reached at every point in the nebula. In that case, the temperature is determined locally by the condition that the energy is gained by photoionization and is balanced by the energy loss through line and continuum emission. As one increases the gas metallicity, the equilibrium temperature decreases. A systematic exploration of the [O iii] and [S iii] temperature behavior with respect to one another was initially carried out as we varied in turn the metallicity Z, the ionization parameter (U), and the spectral energy distribution (SED). It turned out to be a daunting challenge to have T_{O iii} significantly exceed T_{S iii} . With standard photoionization, we only obtained promising results when we considered combining ionized regions of greatly varying metallicities (Sect. 3.3). As an alternative to metallicity inhomogeneities, we subsequently initiated the exploration of three “extraneous” factors: (1) the possible contribution of shock excitation; (2) the presence of temperature inhomogeneities as introduced by Peimbert (1967; cf. PeñaGuerrero et al. 2012a, and references therein); and (3) the effects of a departure from a MaxwellBoltzman distribution for the electron energies (Nicholls et al. 2012: hereafter NDS12).
3.1. Photoionization calculations with mappings i
To compute nebular spectra, we used the multipurpose shock/photoionization code mappings ic described by Ferruit et al. (1997; see also Binette et al. 1985), in which the ionization and thermal structure of the nebula are computed in a standard fashion assuming local ionization and thermal equilibria throughout the nebula. The calculation proceeds at a progression of depths in the nebula until the entire ionization radiation has been absorbed. The updated version of the current code (mappings ie) now includes the option of considering a nonMaxwellBoltzmann distribution of the electron energies when calculating excitation rates, as described in Sect. 3.4.
3.1.1. Geometry
We adopted a simple slab geometry for the ionized gas, since for the lines of interest the integrated line spectrum does not appreciably differ from that calculated with spherical geometry (Shields 1974). This is the geometry preferred for the Orion nebula that Baldwin et al. (1991) modeled in detail. Hii regions are very chaotic in appearance, showing evidence of strong density inhomogeneities, with much of the emission coming from condensations or ionized “skins” of condensations. Modeling this level of complexity is beyond the current means of photoionization models and is undesirable, unless the observations were sufficiently detailed to constrain such models. Finally, the definition of the ionization parameter in the spherical case depends crucially on the Strömgren radius size (R_{s}), which depends on the volume filling factor of the gas and the radius R_{g} of the empty gas volume near the ionizing stars, both of which are poorly constrained if not downright illdefined.
We adopted the usual definition of the ionization parameter: (2)were c is the speed of light, h the Planck constant, and R_{g} the distance of the slab from the ionizing star, in other words, U is the ratio between the density of ionizing photons impinging on the slab (ϕ_{H}/c) and the slab H density n_{H} . With a slab geometry, we assume that the geometrical thickness R_{s} − R_{g} of the ionizationbounded slab is insignificant with respect to R_{g}.
In our grid of photoionization calculations, we varied U over a wide range from 0.001 to 0.046, although in practice a more limited range may suffice (e.g., Evans & Dopita 1985). Campbell (1988) found that most of the nebulae in her sample lie in the range − 1.9 ≳ lgU ≳ − 2.5. All photoionization calculations were ionizationbounded and we assumed the lowdensity regime, by adopting n_{H} = 10 cm^{3}. Some emission lines of low critical density, such as [S ii]λλ6716, 31 or [O ii]λλ3727, 29, would be affected by collisional deexcitation in those few cases where the nebular densities reach values as high as ~500 cm^{3}. Nevertheless, since the overall energy budget remains set by the majority of other lines with much higher critical density, the conclusions we reached about the T_{O iii} and T_{S iii} temperatures do not depend on the density for the overwhelming majority of Hii regions.
3.1.2. Reference set of abundances
We express the gas metallicity Z with respect to the solar abundance set of Asplund et al. (2006), where He/H = 0.085 and the other elements are C,N,O,Ne,Mg,Si,S,Ar,Ca,Fe: H = Z × (250,60,460,69,34,32,14,1.5,2,28) × 10^{6} Z_{⊙}, where Z is the metallicity in solar units. We computed nebular models of varying metallicities covering a range from 1% solar (Z = 0.01 Z_{⊙}) to more than twice solar (Z = 2.5 Z_{⊙}).
3.1.3. Stellar spectral energy distributions
We explored different stellar SEDs for the ionization source of the Hii regions. This did not turn out to be a sensitive parameter. In the models reported below, we compare photoionization models with an SED from a zeroage instantaneous stellar burst assuming a Chabrier (2003) initial mass function with M_{sup} = 100 M_{⊙}, and solar metallicity, obtained from the stellar population synthesis models by Bruzual & Charlot (2003). The stellar atmospheres used to derive the stellar ionizing SED are taken from the Kurucz spectra library (distributed by Lejeune et al. 1997) for stars with T_{eff} < 50 000 K, and the Rauch (1997) NLTE model atmospheres for the higher temperature stars (Magris et al. 2003, hereafter M03). We also explored the alternative case of a single stellar atmosphere of temperature T_{eff}, using the LTE atmosphere models of Hummer & Mihalas (1970, hereafter HM70).
3.2. Models that consider temperature inhomogeneities
It has been proposed that the differences between the temperatures from ORLs and CELs could be due to temperature inhomogeneities (e.g. Peimbert 1967). To evaluate how these inhomogeneities affected the line ratios, we adopted the temperature inhomogeneity characterization described in Binette & Luridiana (2000, hereafter BL00), who used different temperatures for the CELs in accordance with the energy separation of the atomic levels involved for each line emission calculation. BL00 did not assess the physical driver of these inhomogeneities. These authors simply assumed their existence and adopted a simple scheme in which the inhomogeneities corresponded to hotspots as a result of an unknown heating process, distinct from the photoelectric heating. These authors considered the global effects of temperature inhomogeneities on the emission lines, based on the statistical approximation introduced by Peimbert (1967), which quantifies the brightness increase or decrease of each line in terms of the departure of the temperature from the average temperature T_{o}. For each transition, Peimbert used an expansion in Taylor series of the temperature that takes the functional dependence with temperature of each line emissivity into account. Following Peimbert (1967), the mean nebular temperature T_{o} , is defined as (3)for a homogeneous metallicity nebula characterized by small temperature inhomogeneities, n_{e} is the electron density, T the electron temperature, and V the volume over which the integration is carried out. The rms amplitude t of the temperature inhomogeneities is defined as (4)We simplified the expression presented by Peimbert (1967), in which t^{2} depends on the density of the ionic species considered (n_{i}) while in the above equations we implicitly consider only ionized H (by setting n_{i} = n_{H + } = n_{e}).
In photoionization calculations, one determines the local equilibrium temperature at every point in the nebula, T_{eq} , that satisfies the condition that the cooling by radiative processes equals the heating due to the photoelectric effect. Assuming the existence of temperature inhomogeneities and adopting the scheme of BL00 that describes these as hotspots, we can define T_{eq} as the temperature floor above which putative temperature inhomogeneities take place. These inhomogeneities are taken into account only in the statistical sense, that is, through the use of appropriately “biased” temperatures that depend on the mean squared amplitude of the inhomogeneities and on the functional dependence of the (de)excitation rate on temperature. In this context, t^{2} also represents the amount of external energy that the hot spots contribute to the nebular energy budget (BL00). Computing the effects of the hotspots on the line intensities is a twostep process: first T_{o} is derived, and second the biased temperatures corresponding to each atomic transition are calculated. For the hotspot scenario that BL00 simulated, the following expression for T_{o} provided a satisfying parametric fit (5)A value of γ of ≃ − 15 was found to provide an adequate parametrization of the behavior of T_{o} .
To compute individual recombination line intensities in the presence of small inhomogeneities (which are proportional to T^{α}), the code uses a biased mean temperature T_{rec} given by (6)Because temperature inhomogeneities have a lesser impact on recombination lines than on collisionally excited lines, we adopted the simplification of considering a single constant value of the power α = −0.83 for all recombination processes.
For collisionally excited lines, when evaluating the excitation rates ( ∝ T^{βij} exp [ − ΔE_{ij}/kT] ) and deexcitation ( ∝ T^{βji}) of a given multilevel ion, each rate ij (population) or ji (depopulation) is calculated using T_{o} (instead of T_{eq} ), and then multiplied by the appropriate correction factor, either (7)in the case of excitation, or (8)in the case of deexcitation. k_{B} is the standard Boltzmann constant and β_{ji} = β_{ij}. These factors, adapted from the work of Peimbert et al. (1995), were applied to all collisionally excited transitions calculated by mappings ie. It typically results in a significant enhancement of collisional rates by an amount that strongly depends on the ratio ΔE_{ij}/kT_{o} for the selected transition ij.
3.3. Models that consider abundance inhomogeneities
We explored higher complexity photoionization calculations in the form of smallscale metallicity inhomogeneities within the ionized gas. One could imagine, for instance, enriched condensations embedded in a more diffuse gas of lower metallicity. A justification for such a scenario was provided by TenorioTagle (1996). It is based on the fact that the metalrich ejecta from type II supernovae (SNe) ought to follow a long excursion into the galactic environment before they are able to mix with the ISM. In this scenario, the violently ejected matter generates largescale superbubbles, which eventually cool nonuniformly, creating in the process multiple repressurizing shocks. This leads to the formation of small, dense cloudlets of metalrich gas, which eventually fall back toward the disk of the galaxy because of gravity. An efficient metallicity mixing with the surrounding unenriched ISM will only take place when these “metallic droplets” (or cloudlets) are caught within the nebulae associated to starforming regions. Stasinska et al. (2007) quantitatively investigated the hypothesis that the photoionization of these droplets can explain the ADFs observed in Hii regions. The authors derived bounds of 10^{13}–10^{15} cm on the droplet sizes to (1) explain that they have not yet been spatially resolved and (2) to ensure that they survive long enough their destruction by diffusion.
The possibility of metallicity inhomogeneities in Hii regions has been studied on two scales: (1) on a scale in which they have not yet been spatially resolved^{3} (Tsamis et al. 2003); and (2) on a scale associated with proplyds observed in Hii regions (Tsamis et al. 2011). An example of case (1) is given by the work of Tsamis & Péquignot (2005), who modeled the line spectrum of the central region of 30 Doradus in the LMC in detail. Their dual metallicity model consists of gas of approximately normal LMC composition (0.36 Z_{⊙}) coextensive with a hydrogenpoor component, in which most elements are enhanced by ≈ 0.9 dex. The relative weight of the hydrogenpoor component is relatively low, at an 8% level of the integrated Hβ flux. In addition to successfully fitting the ORLs from CNO elements (in addition to the UV/optical CELs), the model of Tsamis & Péquignot leads to greatly improved predictions for several other optical and IR lines as well. For planetary nebulae (PNe), the idea has been studied earlier (Liu et al. 2000; Ercolano et al. 2003; Tsamis et al. 2004; Bohigas 2009; Yuan et al. 2011), although the enriched gas component has a totally different origin than that envisaged by TenorioTagle (1996).
To study abundance inhomogeneities, we adopted a dual metallicity modeling approach in which we combined pairs of isobaric photoionization models that simply differ in metallicity. Although we did not consider all physical implications of embedded condensations, such as an intermingled diffuse field, thermal conduction, or a change in ionization parameter, our calculations suffice in evaluating how abundance inhomogeneities can affect the [S iii] and [O iii] temperatures. The technique we employed is the following: for a given pair of photoionization models, A and B, we multiplied all line fluxes of model A by a weight factor w_{i} and added them to the corresponding fluxes from model B after multiplying these by the complementary weight 1 − w_{i}. Subsequently, we calculated the line ratios of the resulting “averaged” model and derived physical quantities of interest, such as the “apparent” T_{S iii} and T_{O iii} temperatures. We built sequences of these models in which only w_{i} is varying along the sequence. We explored the case where models A and B differ only in their metallicity and the case where they differ in both U and Z. Since the results turned out to be qualitatively similar in our temperature diagram (provided the contrast in U was ≲ 10), we will only show cases where both models A and B share the same value in U . Our sequences can be defined by three parameters: the ionization parameter U_{A} , the metallicity Z_{A} of the (lower abundance) model A, and the metallicity contrast between models A and B given by θ = lg(Z_{B}/Z_{A}).
3.4. Models with a nonMaxwellBoltzmann energy distribution
Recently, NDS12 explored the possibility that electrons in Hii regions and PNe depart from a MaxwellBoltzmann equilibrium energy distribution. To represent the electron energy distributions, these authors adopted the parametrization given by a κdistribution, which is a generalized Lorentzian distribution initially introduced by Vasyliunas (1968). It is widely used in the studies of solar system plasmas (Livadiotis & McComas 2011; Livadiotis et al. 2011) where evidence of suprathermal electron energies abounds. NDS12 found that a κdistribution of electron energies resolved many of the difficultie encountered when attempting to reproduce the temperatures observed in nebulae. They concluded that mild departures from a MaxwellBoltzmann distribution is sufficient for reproducing the spread in temperature values found in diagnostic diagrams such as T_{rec} vs. ⟨T_{OII} + T_{NII}⟩ or T_{SII} vs. T_{OII} in Hii regions or T_{rec} vs. T_{O iii} in PNe. The κdistribution is a function of temperature and κ. In the limit as κ → ∞, the distribution becomes a MaxwellBoltzmann distribution. From examining data from Hii regions and PNe, NDS12 found that κ ≳ 10 is sufficient to account for most objects they compared their calculations with. These authors proposed various mechanisms that would result in a suprathermal distribution, such as (1) magnetic reconnection followed by the migration of highenergy electrons along field lines, and by the development of inertial Alfvén waves; (2) the ejection of highenergy electrons from the photoionization process itself (when the ultraviolet source is very hard, as in PNe or AGN); (3) or from photoionization of dust (Dopita & Sutherland 2000); and (4) by Xray ionization, resulting in highly energetic (~1 keV) innershell (Auger process) electrons (e.g., Aldrovandi & Gruenwald 1985).
To implement this new option in mappings ie, we proceeded as follows. In the case of collisionally excited lines of a given multilevel ion or atom (including Hi and He ii), each excitation rate ij is calculated first, using the collision strength Ω_{ij}, the excitation energy E_{ij}, and the local equilibrium temperature T = T_{eq}, as is customary (Osterbrock 1989). Next, each rate is multiplied by the corresponding correction factor (9)which is an expression taken from the work of NDS12 that corresponds to the ratio of the rate computed assuming a κdistribution over the rate that would be derived assuming a MaxwellBoltzmann distribution instead. To prevent numerical overflows, the natural logarithm of Eq. (9) is evaluated first. To compute the function Γ(z), we adopted the expression (Nemes 2010) (10)To derive the correction factors for the deexcitation rates ji, we used the above expression, but evaluated using ΔE_{ji} = 0. This coefficient’s value () is unique for a given combination of κ and T. It can be applied to thresholdfree processes when their cross section approximately behaves as ∝ E^{1}, as is the case with recombination rates. However, rather than applying correction factors to the recombination rates, we opted to use an effective recombination temperature T_{rec} that would result in the same rate coefficients as those obtained by applying the correction factor cf_{0}. Since recombination rates approximately vary as ∝ T^{ − 1/2}, it follows that the effective temperature to be used when computing recombination rates is given by T_{rec} ≃ T/(cf_{0})^{2}.
We also implemented the enhancement of collisional ionization that can result from κdistributions. There are two sets of correction factors: one applies to the excitationautoionization process (EA), and the other to direct collisional ionization process (DI). In the case of EA, we found from the work of Owocki & Scudder (1983) that Eq. (9) can be used to derive the correction factors when the excitation energy is replaced by the relevant autoionizing energies. These coefficients are then applied to the EA rates. In the case of DI, mappings ie uses the method proposed by Arnaud & Rothenflug (1985) to derive the ionizing rates, which is based on a fourparameter functional fit of the DI crosssections σ(E_{ij}) (their Eq. (1)), for which we lack an analytical prescription that would allow us to compute the correction factors when considering a κdistribution. We therefore opted for a numerical integration. For each electron shell j of ion i we numerically integrated σ(E_{ij}) assuming the selected κdistribution for the electron energies. We divided the result by the same numerical integration, assuming instead a MaxwellBoltzmann distribution. The resulting ratio defines the correction factor , which was then applied to the corresponding DI rate computed with the Arnaud & Rothenflug (1985) expression. For the calculations presented here, the enhancement of collisional ionization caused by a κdistribution was found to have no impact on the ionization structure of Hii regions.
3.5. Photoionization models with shock heating
To explore the possible impact of shocks in Hii regions, we used mappings ie to calculate steadystate planeparallel shock models, assuming a zero preshock magnetic field. This shock wave, by construction, is set to propagate in a nebula that is already photoionized by hot stars. The nebular gas constitutes the preshock zone with a temperature T_{pre} that is evaluated assuming photoionization equilibrium. We did not include the surrounding unshocked nebular emission in our calculations. Rather than assuming an arbitrary opacity ahead of the shock front, we adopted an unabsorbed ionizing SED. This means that we assumed that the shock wave occurs close to the inner nebula boundary, if homogeneous, or close to the inner boundary of the brightest filaments otherwise. Throughout the calculations, we assumed the lowdensity limit regime by adopting a density n_{pre} = 1.0 cm^{3}, which we term the preshock density. The basic input parameter of these models will be the postshock temperature T_{post} . Using the values of n_{pre} , T_{pre} and T_{post} , we computed the RankineHugoniot conditions to determine the shock velocity V_{s} and the postshock density n_{post} . Both quantities appear in Table 1. This density n_{post} is about four times that of n_{pre} when V_{s} corresponds to a high Mach number. Downstream of the shock front, as the postshock gas progressively cools, its density increases since the gas pressure is conserved throughout the postshock region. We adopted the previous zeroage instantaneous starburst model of M03 (Sect. 3.1.3) as stellar SED that is responsible for ionizing and heating the preshock zone. The ionization parameter is defined by Eq. (2) as before, except that n_{pre} substitutes n_{H} .
One caveat of our models is, however, that the putative shock wave is set to travel “inward”, toward the front face of the nebula, while one would favor the opposite, that is, we expect the shock wave to enter the nebula from the front face (that is, from the side exposed to the incoming stellar UV and stellar wind) and to propagate outward. The reason we are limited to the “inward” direction for the shock propagation has to do with the way the opacity is integrated, which is subsequently used to determine the radiative transfer of the ionizing flux. In mappings ie, this integration starts at the shock front and proceeds downstream in parallel with the temporal evolution of the shocked gas. The opacity is integrated along the flow. This way, the radiative transfer of the radiation emitted within the hot layers of the shock as well as that from the hot stars can be resolved. However, the opposite “outward” shock propagation might be preferable, for instance in scenarios where shocks are generated by a stellar wind or other stellar ejection mechanisms. In these scenarios where the ionizing photons enter downstream, from the backside of the shock, one needs to know from the start how the opacity behaves across the entire shocked layer. This requires an iterative method to ensure a proper radiative transfer determination. AldrovandiViegas & Contini (1985) developed such an iterative procedure for their code SUMA (Viegas & Contini 1997) and applied it to model the narrow line region (NLR) of active galactic nuclei (AGN). This option, however, is not available with the code mappings ie. To circumvent this limitation, we computed models for which the external ionizing flux is sufficiently intense so that the shocked layer’s depth remains a small fraction of the total ionized depth of the nebula in absence of shocks. When this condition is satisfied, the shocked layer can be expected to show a similar structure, whether it is propagating inward or outward, with respect to the direction of the dominant radiation field.
To determine the preionization conditions of the gas before it enters the shock, we assumed equilibrium ionization and equilibrium temperature. This is a straightforward assumption, since the preshock gas is immersed in a radiation field so intense that it is already highly ionized and at or near equilibrium temperature. The ion and electron temperatures were set to be equal in all shock models. Starting at the shock front, the adiabatic cooling (and recombination) of the plasma is followed only until such time as cooling by radiative processes has sufficiently decreased, so that it becomes equal again to the heating due to photoionization. Hence, the line intensities that we computed strictly apply to the shockheated layers (that is, the outofequilibrium layers beyond which steadystate photoionization eventually settles again, but this zone is not computed). We recall that the models do not include line emission from the photoionized preshock region either.
Fig. 2 Photoionization calculations with an SED corresponding to a zeroage instantaneous starburst model from M03 with M_{sup} = 100 M_{⊙}. A sequence of models of increasing metallicity Z is represented by the magenta solid line, for which U is constant at 0.01. The metallicities covered a range from Z = 0.01 Z_{⊙} to Z = 2.5 Z_{⊙}, where Z increases by 0.2 dex from model to model. Three sequences of models are superimposed in which U is increased (rather than Z), covering the range 0.001 to 0.046. Each sequence has a different metallicity, as follows: Z = 0.01 Z_{⊙} (cyan dotted line), Z = 0.25 Z_{⊙} (light green dotted line) and Z = 1.0 Z_{⊙} (blue dotted line). In this and subsequent figures, a light gray dashed line represents the locus of equal temperatures. 
To facilitate comparisons between photoionized plasma and steadystate shock waves, we will make use of column densities. For the assumed planeparallel geometry, the line emission measures are what models integrate to derive line luminosities. However, our interest in referring to integrated column densities is that the internal structure (i.e., the behavior with depth) of the ionization, the temperature, or the cumulative line ratios are the same when comparing homologous models and when the quantities of interest in these models are plotted as a function of column depth. In the familiar photoionization case, models with different densities are homologous, provided the ionization parameter U is the same, at least in the regime of sufficiently low densities where collisional deexcitation does not affect the thermal structure. In the case of shocks, constant V_{s} models (in the lowdensity regime) are homologous, and their structure remains unchanged as a function of column density, independently of the postshock density n_{post} . The reason is that although the internal energy of the shocked gas after the shock discontinuity scales as n_{post} , the local cooling rate scales as the density squared, which means that the cooling will proceed on a timescale that evolves as n_{H}^{1}, that is, decreasing with increasing postshock density. This keeps the behavior of the temperature or ionization invariant as a function of depth, when expressed as a column density. In the case of combined shockphotoionization models, models will be homologous if both U and V_{s} are kept constant when changing the input density. In the case of pure photoionization, the integrated column density of ionized H for an ionizationbounded model is given^{4} by , where α_{B} is the caseB recombination rate coefficient of hydrogen. For the ionization parameter value of U = 0.01 that we assumed in Sect. 4.5 for the combined shockphotoionization models, is 8.2 × 10^{20} cm^{2} in the pure photoionization case with Z = 1.0 Z_{⊙} (in absence of shocks). In the case of a steadystate shock wave, increases monotonically with shock velocity^{5}V_{s} without external photoionization. It is a manifestation of the longer elapsed time required for the shocked plasma to recombine and cool from T_{post} to ~100 K. In our calculations, for which the shock wave was immersed in a photoionized plasma, when integrated over the outofequilibrium layers only, is as low as 2.5 × 10^{19} cm^{2} for a shock velocity of 56 km s^{1}. Hence, the aforementioned condition of considering only shocked layers that are geometrically thin with respect to an ionizationbounded photoionized layer is clearly satisfied with the velocities considered in this work ( ≤ 65 km s^{1}).
4. Model results and discussion
In figures without the error bars, such as Fig. 2, the span in [O iii] temperatures to the right of the dashed line in the lower left quadrant stands out more clearly. At the lowtemperature end, T_{O iii} can exceed T_{S iii} by as much as ~3000 K. We now investigate which types of models or parameter values would cause the [O iii] temperatures to be appreciably higher than those from [S iii]. We begin with homogeneous photoionization calculations.
Fig. 3 Comparison of two sequences in metallicity that differ by the shape of their ionizing SED: a single temperature star of either T_{eff} = 55 000 K (light green line) or T_{eff} = 40 000 K (blue dotted line). The previous zeroage instantaneous starburst SED of M03 with M_{sup} = 100 M_{⊙} appears underneath (magenta line). Each sequence comprises 13 models with the metallicity increasing by a factor 0.2 dex from model to model, starting at Z = 0.01 Z_{⊙} and ending at 2.5 Z_{⊙}. The ionization parameter is U = 0.01 in all calculations. 
4.1. The pure homogeneous photoionization case
When the plasma temperature is set by the condition that equilibrium ionization and equilibrium temperature apply locally everywhere within the nebula, and that there are no inhomogeneities in either the local densities or the metallicities, one finds that varying the input parameters gives little leverage on the behavior of both temperatures. This is illustrated by the models presented in Fig. 2. The magenta solid line shows a sequence of 13 models (all with U = 0.01 and the zeroage instantaneous starburst SED of M03), in which only the metallicity increases along the sequence, by 0.2 dex from model to model, covering the range Z = 0.01 to 2.5 Z_{⊙}. Higher metallicity models result in cooler nebulae, while the lowest Z models occupy the upper right corner. What we infer from this plot is that below 14 000 K, there is no significant departure of the models from the equal temperature (dashed) line. Only beyond 16 000 K does a modest departure from the equal temperature line begin to occur. To complete the picture, we show three sequences, of six models each, starting at U = 0.001 and ending at 0.046, along which the ionization parameter increases by 0.33 dex from model to model. These three sequences are characterized by different metallicities Z of 0.01 Z_{⊙} (cyan dotted line), 0.25 Z_{⊙} (light green dotted line) and 1.0 Z_{⊙} (blue dotted line), respectively. They illustrate quantitatively that U is not a significant parameter, except at very low metallicities. In this case, the value of U matters because of the cooling provided by collisional excitation of Hi, which plays an increasingly important role in establishing the equilibrium temperature when metal cooling becomes progressively inoperative beyond Z ≲ 0.05 Z_{⊙}. Increasing U reduces the neutral fraction of hydrogen H^{0}/H, which diminishes the cooling that Hi collisional excitation can provide, resulting in a hotter nebula.
In Fig. 3, we compare the case of using an SED consisting of a single star atmosphere model (HM70) with T_{eff} = 55 000 K (light green solid line) with the case of using the previous M03 SED (magenta line). Both sequences are metallicity sequences covering the interval Z = 0.01 to 2.5 Z_{⊙}. We also explored stellar SEDs of lower T_{eff} or SEDs from M03 of lower M_{sup} and/or of an evolved population. To prevent cluttering in the figure with too many overlapping lines, these models are not plotted, except for the T_{eff} = 40 000 K HM70 stellar atmosphere SED (blue dotted line). Those models showed that reducing either M_{sup} or T_{eff} does not cause any improvement in reproducing the spread in the observed T_{S iii} and T_{O iii} temperatures of the lower left quadrant. The SEDs that we adopted for the calculations shown in Fig. 3 illustrate that a change in the continuum hardness is not effective in resolving the T_{S iii} –T_{O iii} temperature issue.
Fig. 4 Sequences of models in which t^{2} is increased from 0.004 to 0.16 by successive factors of 1.58 ( = 0.2 dex). All calculations were carried out with U = 0.01 and the zeroage instantaneous starburst model SED with M_{sup} = 100 M_{⊙} from M03. Each sequence corresponds to a different metallicity Z, which is indicated in the line color legend in Z_{⊙} units. 
4.2. Results of models with temperature inhomogeneities
Because mappings ie offers the possibility of simulating the effects of temperature inhomogeneities following the method described in Sect. 3.2, we explored various runs of models in which the mean squared amplitude of these inhomogeneities (t^{2}) is progressively increased^{6}. The t^{2} values that we explored range from 0.004 to 0.16, fully encompassing that found in Hii regions. For instance, in the sample of 28 nebulae used by PeñaGuerrero et al. (2012b, hereafter PG12) to recalibrate the Pagel method, 70% of the objects have t^{2} values within the interval 0.02 ≤ t^{2} ≤ 0.04, one object has a lower value, while the rest (25%) lie within the range 0.06–0.12. The results of our calculations are shown in Fig. 4. All calculation were carried out with U = 0.01 and the starburst SED from M03. Each t^{2} sequence comprises nine models, with t^{2} increasing by a factor 0.2 dex (=1.58) from model to model. The different t^{2}sequences differ only in their metallicity, whose value appears in Z_{⊙} units in the legend of Fig. 4. We found that the modifications to the procedure for calculating emissivities (with respect to pure photoionization) suffice to explain the growing divergence of lowmetallicity sequences from the equal temperature line as t^{2} is increased. Changes in ionization structure due to calculating the recombination rates at T_{rec} (see Eq. (6)), rather than at the equilibrium temperature, are minor.
Figure 4 shows that, for models within the lower left quadrant, increasing t^{2} increases the temperatures T_{S iii} and T_{O iii} by comparable amounts. As a result, the t^{2} sequences do not cover the spread in the data observed to the right of the equal temperature line in the lower left quadrant. Only within the upper right quadrant do the calculated temperatures diverge from the dashed line as t^{2} is increased. Hence, temperature inhomogeneities cannot account by themselves for the excess in T_{O iii} temperatures observed in the lower left quadrant.
Fig. 5 Sequences of models in which the relative weight w_{i} of two photoionization calculations is varied, as described in Sect. 3.3. For each sequence, the metallicity of its model A in solar units appears in the line color legend followed by the contrast θ in dex. An open circle indicates the position where the relative weight applied to model A and B is equal (w_{i} = 0.5). In all the above models, the ionization parameter is U = 0.01 and the SED corresponds to the zeroage instantaneous starburst model with M_{sup} = 100 M_{⊙}. 
4.3. Results of models with metallicity inhomogeneities
We now turn to the dual metallicity modeling approach described in Sect. 3.3, which has the advantage of not requiring an external contribution to the energy balance of the nebula. We report on sequences of models that have identical values in ionization parameter U_{A} = U_{B} = 0.01. Reducing this parameter to 0.001 did not affect the results in a noticeable way, except in cases where Z_{A} was ≪ 0.2 Z_{⊙}. We recall that along any sequence, only the relative weight w_{i} applied to models A and B varies, within the interval 0 ≤ w_{i} ≤ 1. In Fig. 5 we present examples of dual metallicity sequences. The θ = 0.7 dex sequence (cyan line) corresponds to a weighted average of a model A with Z_{A} = 0.2 Z_{⊙} with a model B with Z_{B} = 1.0 Z_{⊙}. The blue line sequence and the light green line sequence have the same Z_{A}, but a higher abundance contrast θ of 0.8 and 0.9 dex, respectively. The maroon dashed line sequence was calculated with a higher value Z_{A} = 0.32 Z_{⊙} and a metallicity contrast of 0.7 dex. We infer from Fig. 5 that sequences with a sufficiently high contrast θ can cover the observed domain of points to the right of the equal temperature (dashed) line. The locus of the various models also suggests that a metallicity contrast extending up to ≲ 1 dex is required to account for the nebulae with the largest temperature differences of T_{O iii} − T_{S iii}. A possible drawback might be that the modeling relies on some finetuning of w_{i} to reach the observed T_{O iii} value. To give an idea of how T_{S iii} and T_{O iii} varies with w_{i}, we use an open circle to indicate the position along each sequence where both models A and B have equal weight values of 0.5. This corresponds to the case where models A and B reprocess the same amount of ionizing photons (and emit about the same flux in recombination lines).
One interesting feature of models that consider abundance inhomogeneities is that the divergence from the equal temperature line becomes markedly smaller toward the hightemperature end. The reason is that beyond Z ≲ 0.1 Z_{⊙}, metals play a progressively vanishing role in affecting the plasma temperature. This behavior is illustrated by the seven sequences that are presented in Fig. 6, which extend over a much wider metallicity range compared to the previous figure. The abundance contrast θ is fixed at 1.0 dex and the metallicity of the model A of each sequence appears in Z_{⊙} units in the figure legend. Data of higher S/N would be essential to ascertain how nebulae behave in the lowmetallicity (hightemperature) quadrant of the diagram.
Fig. 6 Seven sequences of dual metallicity models covering a wider range in metallicities. The metallicity contrast between Z_{B} and Z_{A} is fixed at θ = 1.0 dex. Each sequence corresponds to a different metallicity Z_{A} value, which is indicated in the line color legend in Z_{⊙} units. An open circle along each sequence indicates the position where the relative weight applied to each model is equal (w_{i} = 0.5). In all the above sequences, the ionization parameter is U = 0.01, and the SED corresponds to the zeroage instantaneous starburst model with M_{sup} = 100 M_{⊙}. 
Fig. 7 Sequences of models in which the parameter κ has the values of 80, 40, 20, 10, and 5. The larger κ, the closer the model to standard models with a MaxwellBoltzmann electron energy distribution. These κ sequences differ in their metallicity Z (expressed in Z_{⊙} units in the line color legend). The dashed magenta line represents a metallicity sequence of standard photoionization models (κ = ∞) (same sequence as in Fig. 2). The light gray arrow illustrates the increase in the T_{S iii} and T_{O iii} temperatures when enhanced Hi excitation is left out from the calculations of the model with Z = 0.01 Z_{⊙} and κ = 5. In all calculations, we adopted U = 0.01 and an SED consisting of a zeroage instantaneous starburst model from M03 with M_{sup} = 100 M_{⊙}. 
4.4. Results of models with a κdistribution
After implementing the effects of a κdistribution on the line intensities and recombination processes as well as on the equilibrium temperature of a photoionized plasma, we proceeded to compute sequences of models where κ has the values of 80, 40, 20, 10, and 5 (which encompasses the range 6 − 20 that NDS12 explored). In Fig. 7, we present seven such sequences that differ only in their gas metallicity. The gas metallicity of each sequence appears in the figure legend, expressed in Z_{⊙} units. The dashed magenta line represents a metallicity sequence with a pure MaxwellBoltzmann distribution (the same Zsequence that was plotted in Fig. 2). All calculations assumed U = 0.01. Our results indicate that κdistributions are able to reproduce the spread in the observed values of T_{O iii} vs. T_{S iii} inside the lower left quadrant and therefore constitute an interesting alternative to the other options studied in this paper.
With κdistributions, the higher the excitation energy of a level, the larger the increase in the corresponding emission line intensity. This explains why [O iii]λ4363 grows considerably more than [S iii]λ6312. At a fixed κ value, the ratio between the respective T_{S iii} and T_{O iii} temperature increases is constant. Hence, independent of metallicities, models with the same κ lie at about the same distance from the equal temperature dashed line in Fig. 7. This is also the case for the other temperature diagrams presented by NDS12. Because mappings ie includes the enhancement in Hi excitation cooling due to the highenergy tail of κdistributions, T_{O iii} and T_{S iii} are seen to increase less and less in Fig. 7 when decreasing metallicity at a fixed κ value. A measure of this effect is represented by the arrow in Fig. 7, which represents the increase in the T_{S iii} and T_{O iii} temperatures when enhanced Hi excitation is left out from the calculations of the model with Z = 0.01 Z_{⊙} and κ = 5. Since collisional excitation of Hi is a strong function of the neutral fraction of hydrogen H^{0}/H (Luridiana et al. 2003), the importance of the effect described above at low metallicities is governed by the ionization parameter^{7}. If we were to increase U , it would result in a lower H^{0}/H fraction and this would reduce the importance of Hi excitation cooling present at very low metallicities and therefore result in a higher nebular temperature.
Fig. 8 Six sequences of combined shockphotoionization models along which the shock velocity V_{s} is increased in locked steps. The preshock conditions are set by photoionization equilibrium, using the parameter U = 0.01 and the M03 SED, and assuming thermal and ionization equilibria (see Sect. 3.5). Each sequence is characterized by a different metallicity Z value, which is expressed in Z_{⊙} units in the line color legend. For each model, the preshock temperature, the postshock density, the postshock temperature, and the resulting shock velocity are given in Table 1. The first two models of the Z = 0.5 Z_{⊙} and 1.0 Z_{⊙} sequences overlap too closely to be distinguished. 
4.5. Results from photoionization models with shock heating
A set of calculations was carried out in which we computed the line intensities of a planeparallel shock wave that propagates within a photoionized gas layer. To simplify the exercise, we adopted the zeroage instantaneous starburst SED and the ionization parameter value of U = 0.01 for the preshock photoionized gas in all calculations. More information about the models is presented in Sect. 3.5.
In Fig. 8 we present six sequences of combined shockphotoionization calculations, which differ by their metallicity Z (see the adopted values expressed in Z_{⊙} units in the figure legend). Preshock temperatures (T_{pre}) are given for each model in Table 1. Along each sequence, it is the postshock temperature (hence, the shock velocity) that is progressively increased from model to model. The exact values of the postshock temperature, T_{post} , and the associated shock velocities, V_{s} , are listed for each model in Table 1 (along with n_{post} ). The associated shock velocities typically cover the interval 20 to 60 km s^{1}. This velocity range is very similar to what is obtained in detailed kinematical studies of GEHRs and circumnuclear starforming regions (see e.g. Rozas et al. 2006a,b; Firpo et al. 2011; Hägele et al. 2010). The results depicted in Fig. 8 turn out to be quite encouraging. In the lower left quadrant, the models are surprisingly successful in covering the spread in object positions to the right of the equal temperature (dashed) line. As the value of the sequences’ metallicity decreases, the calculated temperatures move closer to the dashed line. A discussion on the behavior of T_{O iii} in combined shockphotoionization models can be found in Appendix B.
Shock heated photoionization models.
Despite the relative success of our calculations, it should be noted that they represent an oversimplification of a very complex problem. A subsequent analysis with 3D hydrodynamical simulations is necessary. As stated earlier, the line intensities and the temperatures calculated in the models shown in Fig. 8 correspond to the zone that is shockheated. That is, they include line emission from the shock discontinuity down to the point where cooling and heating (by the UV flux) become equal again, at which point the calculations stop^{8}. As indicated in Sect. 3.5, the total column density of an isobaric ionizationbounded photoionization model with U = 0.01 is 32 times that of the shocked layer of the 56 km s^{1} model from the Z = 1.0 Z_{⊙} sequence (last model of sequence No. 4 in Table 1). This is a concern, since it appears to suggest that one might need 32 successive shocks propagating across the whole nebula for the resulting emission lines to reflect the presence of shocks at the level shown in Fig. 8, a rather implausible proposition. Instead of moving inward, we assumed a possibly more realistic scenario in which the shocks are moving outward as a result of a stellar wind generated by the ionizing stars. In this case, the following scenario emerges. A high Mach number shock results in a compression factor of ~4 of the pressure (across the shock front) and, from that point on, the density increases because the pressure is conserved across the adiabatic flow, as the gas cools and the temperature decreases. At the point downstream where photoheating balances cooling again, the density has reached a value of about 4T_{post}/T_{pre}, where T_{pre} is an estimate of the local equilibrium temperature. This amounts to a compression factor of ~28. Hence the effective ionization factor near the backphotoionized zone is ~3.5 × 10^{4}, and even if this layer was kept photoionized, its [O iii]λ5007 emission contribution could not be significant. As the shock wave progresses outward within the nebula, the accumulating denser gas shell behind the shock will eventually quench the photon flux responsible for photoionizing the preshock zone ahead. The latter would quickly evolve into a lowexcitation fossil nebula (ahead of the shock) without significant emission in [O iii]λ5007. If we consider the back of the outwardmoving shock (facing the incoming ionizing radiation), we may expect the accumulating dense shell to break up into filaments as a result of instabilities.
Another complication that our simple planeparallel calculations does not consider is that Hii regions are characterized by a low volume filling factor (Osterbrock 1989). The fine spatial structure of Hii regions has so far not been fully resolved. Although progress has been significant, we have not yet reached the microstructure scale of nebular gas (e.g. O’Dell et al. 2003; MesaDelgado et al. 2008; MesaDelgado & Esteban 2010). It has been proposed that nebulae consist of gas condensations or filaments close to being ionizationbounded (e.g. Tsamis & Péquignot 2005). The shock scenario described above can in principle be transposed to the case of dense condensations. 3D hydrodynamical models, such as the one developed by Raga et al. (2008), would be required to simulate such a case and to determine the temperature of the diffuse gas in which the largescale shock surrounding the condensations (or the filaments) propagates.
To summarize, despite the caveats just discussed, we consider that shockphotoionization models have the potential of reproducing the observed spread in temperatures observed below the equal temperature line. These models show a tendency quite similar to that of abundance inhomogeneity sequences presented in Sect. 4.3, since T_{O iii} approaches T_{S iii} as the gas metallicities decrease (toward the upper right quadrant). This behavior is the opposite of that observed with temperature inhomogeneities as modeled with t^{2} parametrization (Sect. 4.2).
Fig. 9 Recombination temperature as a function of T_{O iii} for various models. Along each sequence, only the metallicity varies. Solid and dashed lines are used to represent . A dotted line of the same color depicts , but it is plotted only for cases where differs from by more than 3%. The ionization parameter is always U = 0.01. In the case of the 50 km s^{1} combined photoionizationshock model sequence (light green), Z assumes the values of 0.01 Z_{⊙}, 0.10 Z_{⊙}, 0.25 Z_{⊙} 0.49 Z_{⊙} 1.0 Z_{⊙}, and 1.6 Z_{⊙} (the higher the metallicity, lower is T_{rec} ). All other sequences cover the same metallicity range of 0.01 Z_{⊙} to 2.5 Z_{⊙}, with a fixed increment of 0.2 dex between successive models. Three sequences represent κdistributions (20, 10, or 5, see legend) and two sequences represent the case of temperature inhomogeneities (t^{2} = 0.05 or 0.10, see legend). The homogeneous pure photoionization sequence is represented by the magenta line (same as in Fig. 2). 
4.6. Recombination temperatures and the Balmer decrement
4.6.1. Behavior of the recombination temperatures
The abundance discrepancy factor can be defined as the ratio, or the logarithmic difference, between abundances of a same ion derived from ORLs and CELs. The ADF is a result of the large difference in the dependence of CELs and ORLs on temperature. To evaluate how the different models presented so far can account for the ADFs reported in the literature, we computed the recombination temperature characterizing each model. To reflect the fact that recombination line emissivities decrease with increasing temperature, we averaged the quantity for each selected ion X^{+i} using a weight across the nebular structure that is proportional to the product of this ion density with the electron density. We therefore define the average recombination temperature of species X^{+i} as . We derived in this way the three recombination temperatures , and , corresponding to the ions H^{+}, He^{+} and O^{++}, respectively. We adopted α = −0.83 (Eq. (6)) for H, He and − 0.87 for O. Since in any model and do not differ by more than 30 K, the temperature will not be discussed further.
To characterize the behavior of as a function of T_{O iii} , we present various sequences of models (solid lines) in Figs. 9 and 10. The metallicity is the only parameter that varies along any of the sequences. A dotted line is used when plotting , but it is only shown for cases where one model of the sequence has differing from by more than 3%. The behavior of in the case of κdistribution (with κ = 20, 10 and 5) is shown in Fig. 9 together with the case of t^{2} models (with t^{2} = 0.05 or 0.10). The Z sequence for the combined shockphotoionization case with V_{s} ≃ 50 km s^{1} is overlaid as well. Interestingly, the t^{2} sequences show a behavior similar to those of κdistributions. Since t^{2} = 0.10 is close to the highest value reported by PG12, and since the corresponding sequence (purple line) overlaps the κ = 10 sequence, we infer that values as low as κ ≃ 5 should not be required to reproduce even the largest ADFs observed in some Hii regions. NDS12 came to the same conclusion using different temperature diagnostics.
Fig. 10 Recombination temperatures as a function of T_{O iii} for models that consider metallicity inhomogeneities as discussed in Sect. 4.3. Along each sequence, only the metallicity varies (see legend for values of Z_{A} and θ). In no model does differs from by more than 3%. All four sequences from Fig. 5 with varying contrast θ are included as well as the seven sequences from Fig. 6 with θ = 1.0 (color–coding and nomenclature remain the same, see legend). An open circle indicates the position where w_{i} = 0.5. The ionization parameter is U = 0.01 throughout. 
All models representing metallicity inhomogeneities previously shown in Figs. 5 and 6 are displayed in Fig. 10 (see the figure line color legend for values of Z_{A} and θ). Wherever substantially exceeds T_{O iii} , the models imply a large ADF, which happens only with the higher metallicity sequences. A pertinent result from Fig. 10 is that for Z_{A} ≲ 0.05 sequences, the models diverge very little from the equal temperature gray dashed line. Moreover, the upward curvature^{9} above the equal temperature line (sequences Nos. 8–11 in Fig. 10) is opposite to that of the higher Z sequences, since in sequences Nos. 8–11 lie above T_{O iii} . These models become unsuitable to account for any ADF among metalpoor nebulae. If we increase the contrast θ much beyond the assumed value of 0.1 dex, we find that it is possible to reverse the sign of the curvature. For instance, for the navy sequence with Z_{A} ≲ 0.01, the metallicity contrast θ had to be increased to 2.4 dex for to become again lower than T_{O iii} . Even though obtaining is always possible with sufficiently high θ values, targeting metalpoor nebulae comes with the proviso that only low values of w_{i} are allowed. It appears to us that such a finetuning in the parameter space is an implausible requirement. Unlike the alternate mechanisms presented in previous Fig. 9, metallicity inhomogeneities cannot account for sizable ADFs in nebulae that have T_{O iii} ≥ 15 000 K, that is, in Hii regions that are metalpoor. Yet, in their recalibration of the Pagel method, PG12 reported five Hii regions beyond > 15 000 K. These nebulae are characterized by an equivalent t^{2} in the range 0.02–0.12, indicative of a temperature markedly lower than T_{O iii} .
Fig. 11 Calculations of the Balmer decrement as a function of metallicity for various models. The thick light gray line is a metallicity sequence representing the homogeneous pure photoionization case with U = 0.01 (same sequence as in Fig. 2). The (short vertical) black line shows a Usequence at a fixed metallicity Z = 0.01 Z_{⊙} (see cyan dotted line in Fig. 2). Three metallicity sequences representing the κdistribution case are shown. The corresponding κ value appears in the line color legend. The light green line represents a metallicity sequence for a 50 km s^{1} shock combined with photoionization. Two sequences (labeled 7 and 8 in the legend) are shown that represent models with metallicity inhomogeneities. The values of Z_{A} (in Z_{⊙} units) and θ (in dex) appear in the line color legend. An open circle indicates the position where w_{i} = 0.5. Yellow dashed lines represent the lowdensity recombination caseB Balmer decrement at three different temperatures (Osterbrock 1989). 
4.6.2. Enhanced Balmer decrement and dereddening
The reddening correction that is routinely applied to observed line intensities can be a concern, even if the value of the Balmer decrement adopted for determining A_{V} is apparently consistent with the subsequently derived [O iii] temperature. In particular, an excessive dereddening of the auroral and nebular line intensities would result in overestimates of the T_{O iii} and other CEL temperatures when these auroral lines occur at shorter wavelengths than the corresponding nebular lines^{10}. This inconsistency can have consequences for the case of κdistributions or large temperature inhomogeneities, since the intrinsic Balmer decrement is not set by T_{O iii} , but rather arises from a plasma at a considerably lower temperature, T_{rec} . For instance, the intrinsic Balmer decrement can approach values ≃3 when the abundances are above solar and T_{rec} ≪ T_{O iii}. A better understanding of temperatures in Hii regions may require that we relate the process of dereddening with temperature determination. To illustrate our concern, we carried out a comparison of the Balmer decrement from models that reproduce the spread in T_{O iii} . The results are shown^{11} in Fig. 11 where the calculated Balmer decrement is plotted as a function of nebular metallicity. Two metallicity inhomogeneity sequences with θ = 1 dex are plotted as a function of average metallicity defined as . These correspond to a model A with Z_{A} = 0.04 Z_{⊙} (purple line) and Z_{A} = 0.16 Z_{⊙} (dark gray line), respectively. We note that they do not depart appreciably from the homogeneous photoionization case. The combined photoionizationshock models show evidence at low metallicities of a slightly higher decrement. For κdistributions, the lower the adopted κ value, the higher the decrement toward low metallicities.
We recall that the lower left quadrant in our temperature diagram contains objects whose metallicities are likely to be higher than 0.2 Z_{⊙} (see Sect. 3). Yet, for our data set, the objects with the lowest T_{S iii} values can have corresponding T_{O iii} temperatures that exceed 10 000 K. At metallicities ~1–2 Z_{⊙}, which κdistribution models favor to reproduce the strongest excesses in T_{O iii} (red line in Fig. 7), the Balmer decrement predicted by these models in Fig. 11 is nearing or even exceeding 3. Hence, when comparing data with models that predict different behaviors of CEL temperatures, it would make sense to first apply a reddening correction consistent with the Balmer decrement that these same models predict.
5. Conclusions
An analysis of 102 pairs of [O iii] and [S iii] temperatures of Hii galaxies, giant extragalactic Hii regions, Galactic Hii regions, and Hii regions from the Magellanic Clouds was carried out. The [O iii] temperatures appear to be higher than the corresponding values derived from [S iii] in objects with T_{S iii} lower than 14 000 K. This trend is present in objects with metallicities Z ≳ 0.2. For the coolest nebulae (the highest metallicities) the [O iii] temperature excess can reach +3000 K. The results from grids of models that attempt to account for this excess can be summarized as follows:

1.
Simple photoionization calculations result in essentially equal [O iii] and [S iii] temperatures in the range of interest (T_{S iii} < 14 000 K). Varying either the ionization parameter U, the metallicity Z, or the hardness of the ionizing SED over a wide range does not alter this result.

2.
Including temperature inhomogeneities of large mean squared amplitude (t^{2}) does not result in values of T_{O iii} higher than T_{S iii} in the quadrant of reliable data (T_{S iii} < 14 000 K).

3.
Metallicity inhomogeneities can successfully reproduce the observed excess in T_{O iii} temperatures. By combining two photoionization models of widely differing metallicities (for instance a model with Z = 0.2 Z_{⊙} with another with Z = 1.26 Z_{⊙}), most of the observed spread in temperatures can be reproduced.

4.
Adopting a κdistribution for the electron energy distribution results in models that successfully reproduce the observed excess in T_{O iii} temperatures. Values of κ in the range ≳ 5 to 20 are favored.

5.
Shock waves that propagate in the photoionized gas can also account for the observed excess in T_{O iii} . The 1D models that we computed are simplistic and cannot be considered a fully autoconsistent description of shock waves propagating in a chaotic medium. 3D hydrodynamical calculations would be required that include proper radiative transfer in the context of outwardpropagating shocks.

6.
Both the shocked nebula model and the photoionization model with abundance inhomogeneities share the property of a diminishing excess in T_{O iii} with decreasing nebular metallicity, while the κdistributions result in a near constant T_{O iii} –T_{S iii} offset when κ is set to a fixed value. Through improved observations, it would be important to determine the behavior of the temperatures in the top right quadrant of the T_{S iii} vs. T_{O iii} diagram.

7.
We derived the recombination temperatures for the various models. Models that consider κdistributions or temperature inhomogeneities t^{2} result in recombination temperatures that are substantially lower than T_{O iii} over the whole metallicity range. On the other hand, models with metallicity inhomogeneities have higher than T_{O iii} in the upper right quadrant where T_{O iii} ≳ 15 000 K. Metallicity inhomogeneities are therefore only applicable to metallicityenriched nebulae and could not account for any sizable ADF in metalpoor nebulae (Z < 0.2 Z_{⊙}).
The alternative of “largescale” metallicity variations has also been studied by TorresPeimbert et al. (1990) and Kingdon & Ferland (1998). Interestingly, for the case of “smallscale” metallicity inhomogeneities, Zhang et al. (2007) showed that a small amount of cold gas is sufficient to account for the observed ORL intensities without the need of reproducing the high t^{2} values “empirically” inferred from Hii regions.
Dopita & Sutherland (1996) proposed the following relationship for a near solar metallicity plasma and in absence of magnetic field: .
Note that t^{2} as implemented in Sect. 3.2 is defined via a secondorder expansion of a Taylor series and as such is accurate only for low values of t^{2} < 0.1.
A metallicity sequence with t^{2} = 0.1 was calculated, but it is not shown for clarity in Fig. 11. The t^{2} sequence closely follows the thick light gray line.
Acknowledgments
L.B. acknowledges support from CONACyT grant CB128556. R.M. received financial support from project IN105511 from PAPIIT (UNAM). We are thankful to Jana Benda Klouda for providing us with the initial translation in English. We thank Diethild Starkmeth for proofreading the article. We acknowledge that many suggestions from an anonymous referee contributed to improve the usefulness and scope of the final version.
References
 Aldrovandi, S. M. V., & Contini, M. 1985, A&A, 149, 109 [NASA ADS] [Google Scholar]
 Aldrovandi, S. M. V., & Gruenwald, R. D. 1985, A&A, 147, 331 [NASA ADS] [Google Scholar]
 Arnaud, M., & Rothenflug, R. 1985, A&AS, 60, 425 [NASA ADS] [Google Scholar]
 Asplund, M., Grevesse, N., & Jacques Sauval, A. 2006, Nucl. Phys. A, 777, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Baldwin, J. A., Ferland, G. J., Martin, P. G., et al. 1991, ApJ, 374, 580 [NASA ADS] [CrossRef] [Google Scholar]
 Binette, L., & Luridiana, V. 2000, Rev. Mex. Astron. Astrofis., 36, 43 (BL00) [NASA ADS] [Google Scholar]
 Binette, L., Dopita, M. A., & Tuohy, I. R. 1985, ApJ, 297, 476 [NASA ADS] [CrossRef] [Google Scholar]
 Binette, L., Ferruit, P., Steffen, W., & Raga, A. C. 2003, Rev. Mex. Astron. Astrofis., 39, 55 [NASA ADS] [Google Scholar]
 Bohigas, J. 2009, Rev. Mex. Astron. Astrofis., 45, 107 [NASA ADS] [Google Scholar]
 Bruzual, A. G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
 Campbell, A. W. 1988, ApJ, 335, 644 [NASA ADS] [CrossRef] [Google Scholar]
 Chabrier, G. 2003, PASP, 115, 763 [NASA ADS] [CrossRef] [Google Scholar]
 Copetti, M. V. F. 2006, A&A, 453, 943 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Dopita, M. A., & Sutherland, R. S. 1996, ApJS, 102, 161 [NASA ADS] [CrossRef] [Google Scholar]
 Dopita, M. A., & Sutherland, R. S. 2000, ApJ, 539, 742 [NASA ADS] [CrossRef] [Google Scholar]
 Ercolano, B., Barlow, M. J., Storey, P. J., et al. 2003, MNRAS, 344, 1145 [NASA ADS] [CrossRef] [Google Scholar]
 Esteban, C., Peimbert, M., TorresPeimbert, S., & Rodríguez, M. 2002, ApJ, 581, 241 [NASA ADS] [CrossRef] [Google Scholar]
 Esteban, C., Bresolin, F., Peimbert, M., et al. 2009, ApJ, 700, 654 [NASA ADS] [CrossRef] [Google Scholar]
 Evans, I. N., & Dopita, M. A. 1985, ApJS, 58, 125 [NASA ADS] [CrossRef] [Google Scholar]
 Ferruit, P., Binette, L., Sutherland, R. S., & Pécontal, E. 1997, A&A, 322, 73 [NASA ADS] [Google Scholar]
 Firpo, V., Bosch, G., Hägele, G. F., Díaz, Á. I., & Morrell, N. 2011, MNRAS, 414, 3288 [NASA ADS] [CrossRef] [Google Scholar]
 Hägele, G. F., PérezMontero, E., Díaz, Á. I., Terlevich, E., & Terlevich, R. 2006, MNRAS, 372, 293 (H06) [NASA ADS] [CrossRef] [Google Scholar]
 Hägele, G. F., Díaz, Á. I., Terlevich, E., et al. 2008, MNRAS, 383, 209 (H08) [NASA ADS] [CrossRef] [Google Scholar]
 Hägele, G. F., Díaz, Á. I., Cardaci, M. V., Terlevich, E., & Terlevich, R. 2010, MNRAS, 402, 1005 [NASA ADS] [CrossRef] [Google Scholar]
 Hägele, G. F., GarcíaBenito, R., PérezMontero, E., et al. 2011, MNRAS, 414, 272 (H11) [NASA ADS] [CrossRef] [Google Scholar]
 Hummer, D. G., & Mihalas, D. M. 1970, MNRAS, 147, 339 (HM70) [NASA ADS] [CrossRef] [Google Scholar]
 Keenan, F. P., Kingston, A. E., & Johnson, C. T. 1988, A&A, 202, 253 [NASA ADS] [Google Scholar]
 Kennicutt, R. C., Jr.,Bresolin, F., & Garnett, D. R. 2003, ApJ, 591, 801 [CrossRef] [Google Scholar]
 Kingdon, J. B., & Ferland, G. J. 1998, ApJ, 506, 323 [NASA ADS] [CrossRef] [Google Scholar]
 Liu, X.W., Storey, P. J., Barlow, M. J., et al. 2000, MNRAS, 312, 585 [NASA ADS] [CrossRef] [Google Scholar]
 Lejeune, T., Cuisinier, F., & Buser, R. 1997, A&AS, 125, 229 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Livadiotis, G., & McComas, D. J. 2011, ApJ, 741, 88 [NASA ADS] [CrossRef] [Google Scholar]
 Livadiotis, G., McComas, D. J., Dayeh, M. A., et al. 2011, ApJ, 734, 1 [NASA ADS] [CrossRef] [Google Scholar]
 LópezSánchez, A. R., Dopita, M. A., Kewley, L. J., et al. 2012, MNRAS, in press [arXiv:1203.5021] [Google Scholar]
 Luridiana, V., Peimbert, A., Peimbert, M., & Cerviño, M. 2003, ApJ, 592, 846 [NASA ADS] [CrossRef] [Google Scholar]
 Magris, C. G., Binette, L., & Bruzual, A. G. 2003, ApJS, 149, 313 (M03) [NASA ADS] [CrossRef] [Google Scholar]
 MesaDelgado, A., & Esteban, C. 2010, MNRAS, 405, 2651 [NASA ADS] [Google Scholar]
 MesaDelgado, A., Esteban, C., & GarcíaRojas, J. 2008, ApJ, 675, 389 [NASA ADS] [CrossRef] [Google Scholar]
 Nicholls, D. C., Dopita, M. A., & Sutherland, R. S. 2012, ApJ, 752, 148 (NDS12) [Google Scholar]
 Nemes, G. 2010, Archiv der Mathematik, 95, 161 [CrossRef] [Google Scholar]
 O’Dell, C. R., Peimbert, M., & Peimbert, A. 2003, AJ, 125, 2590 [NASA ADS] [CrossRef] [Google Scholar]
 Osterbrock, D. E. 1989, in Astrophysics of Gaseous Nebulae and Active Galactic Nuclei (Mill Valley: University Science Books) [Google Scholar]
 Owocki, S. P., & Scudder, J. D. 1983, ApJ, 270, 758 [NASA ADS] [CrossRef] [Google Scholar]
 Peimbert, M. 1967, ApJ, 150, 825 [NASA ADS] [CrossRef] [Google Scholar]
 Peimbert, M., Luridiana, V., & TorresPeimbert, S. 1995, Rev. Mex. Astron. Astrofis., 31, 131 [NASA ADS] [Google Scholar]
 PeñaGuerrero, M. A., Peimbert, A., Peimbert, M., & Ruiz, M. T. 2012a, ApJ, 746, 115 [NASA ADS] [CrossRef] [Google Scholar]
 PeñaGuerrero, M. A., Peimbert, A., & Peimbert, M. 2012b, ApJ, 756, L14 (PG12) [NASA ADS] [CrossRef] [Google Scholar]
 PérezMontero, E., & Díaz, A. I. 2005, MNRAS, 361, 1063 [NASA ADS] [CrossRef] [Google Scholar]
 PérezMontero, E., & Díaz, Á. I. 2007, MNRAS, 377, 1195 [NASA ADS] [CrossRef] [Google Scholar]
 PérezMontero, E., Díaz, A. I., Vílchez, J. M., & Kehrig, C. 2006, A&A, 449, 193 (PM06) [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 PérezMontero, E., GarcíaBenito, R., Hägele, G. F., & Díaz, Á. I. 2010, MNRAS, 404, 2037 [NASA ADS] [Google Scholar]
 PérezMontero, E., Vílchez, J. M., Cedrés, B., et al. 2011, A&A, 532, A141 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Podobedova, L. I., Kelleher, D. E., & Wiese, W. L. 2009, J. Phys. Chem. Ref. Data, 38, 171 [NASA ADS] [CrossRef] [Google Scholar]
 Raga, A. C., Riera, A., Mellema, G., Esquivel, A., & Velázquez, P. F. 2008, A&A, 489, 1141 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rauch, T. 1997, A&A, 320, 237 [NASA ADS] [Google Scholar]
 Rola, C., & Pelat, D. 1994, A&A, 287, 676 [NASA ADS] [Google Scholar]
 Rozas, M., Richer, M. G., López, J. A., Relaño, M., & Beckman, J. E. 2006a, A&A, 455, 539 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rozas, M., Richer, M. G., López, J. A., Relaño, M., & Beckman, J. E. 2006b, A&A, 455, 549 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rubin, R. H. 1986, ApJ, 309, 334 [NASA ADS] [CrossRef] [Google Scholar]
 Shields, G. A. 1974, ApJ, 193, 335 [NASA ADS] [CrossRef] [Google Scholar]
 Stasińska, G. 1980, A&A, 84, 320 [NASA ADS] [Google Scholar]
 Stasińska, G., TenorioTagle, G., Rodríguez, M., & Henney, W. J. 2007, A&A, 471, 193 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 TenorioTagle, G. 1996, AJ, 111, 1641 [NASA ADS] [CrossRef] [Google Scholar]
 TorresPeimbert, S., Peimbert, M., & Peña, M. 1990, A&A, 233, 540 [NASA ADS] [Google Scholar]
 Tsamis, Y. G., & Péquignot, D. 2005, MNRAS, 364, 687 [NASA ADS] [CrossRef] [Google Scholar]
 Tsamis, Y. G., Barlow, M. J., Liu, X.W., Danziger, I. J., & Storey, P. J. 2003, MNRAS, 338, 687 [NASA ADS] [CrossRef] [Google Scholar]
 Tsamis, Y. G., Barlow, M. J., Liu, X.W., Storey, P. J., & Danziger, I. J. 2004, MNRAS, 353, 953 [NASA ADS] [CrossRef] [Google Scholar]
 Tsamis, Y. G., Walsh, J. R., Vílchez, J. M., & Péquignot, D. 2011, MNRAS, 412, 1367 [NASA ADS] [Google Scholar]
 Vasyliunas, V. M. 1968, J. Geophys. Res., 73, 2839 [NASA ADS] [CrossRef] [Google Scholar]
 Viegas, S. M., & Clegg, R. E. S. 1994, MNRAS, 271, 993 [NASA ADS] [CrossRef] [Google Scholar]
 Viegas, S. M., & Contini, M. 1997, IAU Colloq. 159: Emission Lines in Active Galaxies: New Methods and Techniques, 113, 365 [NASA ADS] [Google Scholar]
 Yuan, H.B., Liu, X.W., Péquignot, D., et al. 2011, MNRAS, 411, 1035 [NASA ADS] [CrossRef] [Google Scholar]
 Zhang, Y., Ercolano, B., & Liu, X.W. 2007, A&A, 464, 631 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
Appendix A: Sample abundance calibration
Figure A.1 shows the relationship between the electron temperatures inferred from the [O iii] lines and the total oxygen abundances 12 + lg(O/H). This sample was built from the literature and is a collection of objects for which measurements of the nebular and auroral emission lines of [O ii] and [O iii] were available. This sample comprises objects from PérezMontero & Díaz (2005) for which the determination of T_{O iii} and 12 + lg(O/H) was possible, using the direct method (see Hägele et al. 2008; H08), as well as the Hii galaxies from H06, H08, and H11. The solid line corresponds to a leastsquares fitting of the data, taking their individual errors into account:
where t_{e}([O iii]) is the [O iii] electron temperature in units of 10^{4} K. The O/H determination above is based on CELs only. It likely underestimates the true O/H abundance ratio by a factor ~2, which would likely have been inferred if the ADF had been determined from ORLs and had been taken into account in the above regression. If the electron energies exhibit a κ nonequilibrium distribution, this is expected to increase the oxygen abundance, as shown by NDS12. Likewise with models that consider temperature or abundance inhomogeneities.
Fig. A.1 Relation between t_{e}([O iii]) and total oxygen abundances (12 + lg(O/H)) for Hii galaxies from H06, H08 and H11 (solid circles), and Hii galaxies (open circles), GEHRs (upward triangles) and diffuse Hii regions in the Galaxy and the Magellanic Clouds (squares) from PérezMontero & Díaz (2005). The solid line is the actual fit to the data. The temperatures are in units of 10^{4} K. 
Appendix B: Behavior of combined shockphotoionization models when V_{s} is increased
In combined shockphotoionization models the span in temperature and in ionization stages of metals covers a much wider range than with photoionization alone. Furthermore, the temporal evolution of the plasma that characterizes shocks is another significant difference. To analyze why the T_{O iii} temperatures increase markedly with V_{s} , we compare two models from the cyan sequence with Z = 1.0 Z_{⊙} in Table 1. When one derives volumeaveraged quantities that are weighted by the density square of the ion O^{+2} or S^{+2}, it is found that the 45.3 km s^{1} and 56.1 km s^{1} models are characterized by a volumeaveraged ionic temperature ⟨T_{O++}⟩ of 8400 K and 9700 K, respectively, while for ⟨T_{S++}⟩ the corresponding values are 8070 K and 7980 K, respectively. Hence, not much separates the two models as far as average properties are concerned. In the 45.3 km s^{1} model, the mean doubly ionized fraction for oxygen ⟨n_{O++}/n_{O}⟩ is 0.164 while the equivalent fraction for sulphur ⟨n_{S++}/n_{S}⟩ is 0.889. In the 56.1 km s^{1} model, the mean doubly ionized fractions are 0.087 and 0.844, respectively. Photoionization is responsible for the small reverse effect of a reduction in ⟨n_{O++}/n_{O}⟩ because the mean density ⟨n_{H + }⟩ actually increases from 21.1 cm^{3} to 32.1 cm^{3}, (between the 45.3 km s^{1} and the 56.1 km s^{1} model, respectively), which has the effect of reducing the effective ionization parameter. We can conclude that combined photoionizationshock models in this velocity regime do not drastically alter the mean degree of ionization of oxygen or sulphur, or their mean temperatures where the doubly ionized ions reside. Indeed, both shock models start with the same ionization fraction at the shock front, of 0.93 for n_{O++}/n_{O} and 0.49 for n_{S++}/n_{S}, but they possess quite different postshock temperatures (35 000 K and 50 400 K, respectively).
To summarize, there are two main factors that account for the behavior of T_{O iii} in Fig. 8: (1) ⟨n_{O++}/n_{O}⟩ ≪ ⟨n_{S + + }/n_{S}⟩, which causes the much hotter region near the shock front to contribute proportionally more to the integrated line spectrum of the [O iii] lines than it does for the [S iii] lines, and (2) the strong dependance of the CELs on the factor ∝ T^{0.5}exp( − ΔE/kT) is responsible for favoring the copious emission of lines that have the highest energy transitions ΔE at the head of the shock where T ≲ T_{post}, as is the case for [O iii]λ4363. This is the main factor that accounts for the high 4363 Å/5007 Å line ratio, and hence for the higher values of T_{O iii} . [S iii]λ6312 on the other hand does not beneficiate much from either factor.
The decrease of both T_{O iii} and T_{S iii} at low velocities in sequences of very lowmetallicity gas (light green line and red line sequences in Fig. 8) reflects the increase in the mean density as the postshock temperature is successively increased, which causes a reduction of the effective ionization parameter and hence an increase in the neutral fraction H^{0}/H, which favors an increase in the cooling due to collisional excitation of Hi. This increased cooling in turn decreases ⟨T_{O++}⟩ and ⟨T_{S++}⟩. A similar mechanism is at work at very low metallicities in U sequences of homogeneous photoionization models (e.g. cyan dotted line in Fig. 2), as discussed in Sect. 4.1.
All Tables
All Figures
Fig. 1 Relation between T_{S iii} and T_{O iii} temperatures inferred from the [S iii] and [O iii] auroraltonebular line ratios for the objects studied by PM06, including the three Hii galaxies later observed by H06. The light gray dashed line depicts the locus of equal temperatures (T_{S iii} = T_{O iii}), while the solid line indicates the errorweighted leastsquares fit to the data (see Sect. 2). 

In the text 
Fig. 2 Photoionization calculations with an SED corresponding to a zeroage instantaneous starburst model from M03 with M_{sup} = 100 M_{⊙}. A sequence of models of increasing metallicity Z is represented by the magenta solid line, for which U is constant at 0.01. The metallicities covered a range from Z = 0.01 Z_{⊙} to Z = 2.5 Z_{⊙}, where Z increases by 0.2 dex from model to model. Three sequences of models are superimposed in which U is increased (rather than Z), covering the range 0.001 to 0.046. Each sequence has a different metallicity, as follows: Z = 0.01 Z_{⊙} (cyan dotted line), Z = 0.25 Z_{⊙} (light green dotted line) and Z = 1.0 Z_{⊙} (blue dotted line). In this and subsequent figures, a light gray dashed line represents the locus of equal temperatures. 

In the text 
Fig. 3 Comparison of two sequences in metallicity that differ by the shape of their ionizing SED: a single temperature star of either T_{eff} = 55 000 K (light green line) or T_{eff} = 40 000 K (blue dotted line). The previous zeroage instantaneous starburst SED of M03 with M_{sup} = 100 M_{⊙} appears underneath (magenta line). Each sequence comprises 13 models with the metallicity increasing by a factor 0.2 dex from model to model, starting at Z = 0.01 Z_{⊙} and ending at 2.5 Z_{⊙}. The ionization parameter is U = 0.01 in all calculations. 

In the text 
Fig. 4 Sequences of models in which t^{2} is increased from 0.004 to 0.16 by successive factors of 1.58 ( = 0.2 dex). All calculations were carried out with U = 0.01 and the zeroage instantaneous starburst model SED with M_{sup} = 100 M_{⊙} from M03. Each sequence corresponds to a different metallicity Z, which is indicated in the line color legend in Z_{⊙} units. 

In the text 
Fig. 5 Sequences of models in which the relative weight w_{i} of two photoionization calculations is varied, as described in Sect. 3.3. For each sequence, the metallicity of its model A in solar units appears in the line color legend followed by the contrast θ in dex. An open circle indicates the position where the relative weight applied to model A and B is equal (w_{i} = 0.5). In all the above models, the ionization parameter is U = 0.01 and the SED corresponds to the zeroage instantaneous starburst model with M_{sup} = 100 M_{⊙}. 

In the text 
Fig. 6 Seven sequences of dual metallicity models covering a wider range in metallicities. The metallicity contrast between Z_{B} and Z_{A} is fixed at θ = 1.0 dex. Each sequence corresponds to a different metallicity Z_{A} value, which is indicated in the line color legend in Z_{⊙} units. An open circle along each sequence indicates the position where the relative weight applied to each model is equal (w_{i} = 0.5). In all the above sequences, the ionization parameter is U = 0.01, and the SED corresponds to the zeroage instantaneous starburst model with M_{sup} = 100 M_{⊙}. 

In the text 
Fig. 7 Sequences of models in which the parameter κ has the values of 80, 40, 20, 10, and 5. The larger κ, the closer the model to standard models with a MaxwellBoltzmann electron energy distribution. These κ sequences differ in their metallicity Z (expressed in Z_{⊙} units in the line color legend). The dashed magenta line represents a metallicity sequence of standard photoionization models (κ = ∞) (same sequence as in Fig. 2). The light gray arrow illustrates the increase in the T_{S iii} and T_{O iii} temperatures when enhanced Hi excitation is left out from the calculations of the model with Z = 0.01 Z_{⊙} and κ = 5. In all calculations, we adopted U = 0.01 and an SED consisting of a zeroage instantaneous starburst model from M03 with M_{sup} = 100 M_{⊙}. 

In the text 
Fig. 8 Six sequences of combined shockphotoionization models along which the shock velocity V_{s} is increased in locked steps. The preshock conditions are set by photoionization equilibrium, using the parameter U = 0.01 and the M03 SED, and assuming thermal and ionization equilibria (see Sect. 3.5). Each sequence is characterized by a different metallicity Z value, which is expressed in Z_{⊙} units in the line color legend. For each model, the preshock temperature, the postshock density, the postshock temperature, and the resulting shock velocity are given in Table 1. The first two models of the Z = 0.5 Z_{⊙} and 1.0 Z_{⊙} sequences overlap too closely to be distinguished. 

In the text 
Fig. 9 Recombination temperature as a function of T_{O iii} for various models. Along each sequence, only the metallicity varies. Solid and dashed lines are used to represent . A dotted line of the same color depicts , but it is plotted only for cases where differs from by more than 3%. The ionization parameter is always U = 0.01. In the case of the 50 km s^{1} combined photoionizationshock model sequence (light green), Z assumes the values of 0.01 Z_{⊙}, 0.10 Z_{⊙}, 0.25 Z_{⊙} 0.49 Z_{⊙} 1.0 Z_{⊙}, and 1.6 Z_{⊙} (the higher the metallicity, lower is T_{rec} ). All other sequences cover the same metallicity range of 0.01 Z_{⊙} to 2.5 Z_{⊙}, with a fixed increment of 0.2 dex between successive models. Three sequences represent κdistributions (20, 10, or 5, see legend) and two sequences represent the case of temperature inhomogeneities (t^{2} = 0.05 or 0.10, see legend). The homogeneous pure photoionization sequence is represented by the magenta line (same as in Fig. 2). 

In the text 
Fig. 10 Recombination temperatures as a function of T_{O iii} for models that consider metallicity inhomogeneities as discussed in Sect. 4.3. Along each sequence, only the metallicity varies (see legend for values of Z_{A} and θ). In no model does differs from by more than 3%. All four sequences from Fig. 5 with varying contrast θ are included as well as the seven sequences from Fig. 6 with θ = 1.0 (color–coding and nomenclature remain the same, see legend). An open circle indicates the position where w_{i} = 0.5. The ionization parameter is U = 0.01 throughout. 

In the text 
Fig. 11 Calculations of the Balmer decrement as a function of metallicity for various models. The thick light gray line is a metallicity sequence representing the homogeneous pure photoionization case with U = 0.01 (same sequence as in Fig. 2). The (short vertical) black line shows a Usequence at a fixed metallicity Z = 0.01 Z_{⊙} (see cyan dotted line in Fig. 2). Three metallicity sequences representing the κdistribution case are shown. The corresponding κ value appears in the line color legend. The light green line represents a metallicity sequence for a 50 km s^{1} shock combined with photoionization. Two sequences (labeled 7 and 8 in the legend) are shown that represent models with metallicity inhomogeneities. The values of Z_{A} (in Z_{⊙} units) and θ (in dex) appear in the line color legend. An open circle indicates the position where w_{i} = 0.5. Yellow dashed lines represent the lowdensity recombination caseB Balmer decrement at three different temperatures (Osterbrock 1989). 

In the text 
Fig. A.1 Relation between t_{e}([O iii]) and total oxygen abundances (12 + lg(O/H)) for Hii galaxies from H06, H08 and H11 (solid circles), and Hii galaxies (open circles), GEHRs (upward triangles) and diffuse Hii regions in the Galaxy and the Magellanic Clouds (squares) from PérezMontero & Díaz (2005). The solid line is the actual fit to the data. The temperatures are in units of 10^{4} K. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.