{'fraction_non_alphanumeric': 0.04123817662138669, 'fraction_numerical': 0.03149954327135575, 'mean_word_length': 4.396326920019081, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Single-Star HII Regions as a Probe of Massive Star Spectral Energy Distributions

The shape of the OB-star spectral energy distribution is a critical component in many diagnostics of the ISM and galaxy properties. We use single-star HII regions from the LMC to quantitatively examine the ionizing SEDs from widely available CoStar, TLUSTY, and WM-basic atmosphere grids. We evaluate the stellar atmosphere models by matching the emission-line spectra that they predict from CLOUDY photoionization simulations with those observed from the nebulae. The atmosphere models are able to reproduce the observed optical nebular line ratios, except at the highest energy transitions>40 eV, assuming that the gas distribution is non-uniform. Overall we find that simulations using WM-basic produce the best agreement with the observed line ratios. The rate of ionizing photons produced by the model SEDs is consistent with the rate derived from the \Halpha\ luminosity for standard, log(g) = 4.0 models adopted from the atmosphere grids. However, there is a systematic offset between the rate of ionizing photons from different atmosphere models that is correlated with the relative hardness of the SEDs. In general WM-basic and TLUSTY atmosphere models predict similar effective temperatures, while CoStar predicts effective temperatures that are cooler by a few thousand degrees. We compare our effective temperatures, which depend on the nebular ionization balance, to conventional photospheric-based calibrations from the literature. We suggest that in the future, spectral type to effective temperature calibrations can be constructed from nebular data.


INTRODUCTION
Through the effects of radiative, mechanical and chemical feedback, massive stars play a critical role in shaping the galaxies in which they reside. These stars ionize their local interstellar medium (ISM) and form bright Hii regions that can be observed in distant galaxies, even where the individual stars are not resolved. Since the Hii region emission-line spectrum is sensitive to both the spectral energy distribution (SED) of the ionizing source and the properties of the nearby gas, they are used as diagnostics for galaxy properties, such as star formation rates and histories (e.g., Leitherer et al. 1999;Hunter & Massey 1990;Kennicutt et al. 2000), properties of the ionizing stellar population and the slope of the IMF (e.g., Baldwin et al. 1981;Kaler 1978;Rigby & Rieke 2004;Stasińska & Leitherer 1996;Copetti et al. 1986;Dufour 1975), and the chemical abundances and chemical evolution of galaxies (e.g., Kewley & Dopita 2002;Edmunds & Pagel 1984;Bresolin et al. 1999).
The shape of the stellar SED is particularly important for these diagnostics. For example, it determines the rate of ionizing photons (Q 0 ) emitted by massive stars, upon which the commonly used Hα star formation rate indicator depends (Kennicutt 1983). Other popular diagnostics use flux ratios of lines with different ionization potential as diagnostics for the effective temperature (T eff ) of the ionizing stars (e.g., Vílchez & Pagel 1988;Stoy 1933;Zanstra 1927). Thus, the shapes of massive star SEDs play key roles in deriving the physical properties of nebulae and stars from observed Hii region spectra.
Since massive stars radiate most of their flux at FUV wavelengths that are relatively inaccessible to observations, we are dependent on the predictions generated by stellar atmosphere models to describe the properties and SEDs of massive stars. To reproduce the stellar SED, the models need to take into account the effects of non-LTE conditions, stellar winds, and line-blanketing (e.g., Kudritzki & Hummer 1990;Schaerer & Schmutz 1994;Stasińska & Schaerer 1997;Lanz & Hubeny 2003;Pauldrach et al. 2001). Calculating these in detail is both challenging and computationally time-consuming. Therefore, the model atmospheres currently available incorporate these processes with different approximations, by balancing an exact treatment of the physics against faster computation. The different treatments affect the shape of the SED and the properties derived from them (e.g., Simón-Díaz & Stasińska 2008;Voges et al. 2008;Giveon et al. 2002). It is therefore crucial to know how well these different atmosphere models represent the true SED of these stars. Hii regions can be used as a test of these atmosphere models. As mentioned earlier, Hii region emission-line spectra strongly depend on the shape of the ionizing SED. Thus, a direct comparison between the emission lines from observed Hii regions and those predicted by photoionization simulations will reveal how well the different atmosphere models represent the SEDs of massive stars (e.g., Oey et al. 2000;Morisset et al. 2004;Esteban et al. 1993;Crowther et al. 1999;Stasińska & Schaerer 1997). For example, previous studies of both Galactic and Magellanic Clouds HII regions have shown that plane-parallel atmospheres are too soft and that line blanketing is a crucial process to include in the models (e.g., Martín-Hernández et al. 2002;Morisset et al. 2002;Giveon et al. 2002;Bresolin et al. 1999). This approach, however, is not always straightforward. To make this comparison, the metallicity and morphology of the gas in the nebula must be constrained. Furthermore, most studies use nebulae that contain many stars that are distributed throughout the ionized region. Properly accounting for these multiple ionizing sources is a major challenge to this approach and vastly limits the constraints that one can put on the SEDs (e.g., Ercolano et al. 2007). Our study circumvents this challenge by using single-star Hii regions in the nearby Large Magellanic Cloud (LMC). The Hii regions are spatially resolved and spherical, which makes the modeling more straightforward. We use photoionization simulations to evaluate how well atmosphere models represent the shape of the massive star SED.
1.1. Description of the model atmosphere codes Our goal in this work is to understand how well the OB atmosphere-model grids that are available to general users represent ionizing stars. For this reason, we use publicly available grids, rather than fitting each stellar spectrum in detail to obtain the appropriate model. We use the O-star grid presented in Smith et al. (2002), hereafter SNC02, in addition to the CoStar (Schaerer & de Koter 1997), TLUSTY (Lanz & Hubeny 2003), and WMbasic (Pauldrach et al. 2001) grids that are already available in the stellar atmosphere library of CLOUDY (Ferland et al. 1998). These massive-star atmosphere codes include important, but complex, physical processes such as non-LTE conditions and the effects of metal lines and winds on the transmitted spectrum. Including these processes in detail is computationally expensive. Therefore, research groups use various methods and algorithms to approximate some of these processes. The differences in the algorithms used result in non-negligible differences in the SEDs produced by different atmosphere codes (e.g., Simón-Díaz & Stasińska 2008).
All the atmosphere codes considered in this work solve for non-LTE radiation transfer. Typically this is accomplished by grouping lines of similar excitation energies together and applying the same non-LTE correction to the populations in that group (e.g., Hubeny & Lanz 1995;Pauldrach et al. 2001). In some atmosphere codes, such as CoStar and WM-basic, the non-LTE solution includes the effects of spherically expanding winds, while others, such as TLUSTY, assume a plane-parallel geometry.
Metal lines in the UV impact the emergent SED in two ways, line blocking and line blanketing. Line blocking refers to the absorption and scattering of the emergent flux due to the higher opacity in the line. Line blocking will increase the temperature in the deeper layers of the star because of the scattering of radiation back towards the star, an effect known as backwarming (Pauldrach et al. 2001). Line blanketing refers to the redistribution of energy to regions where the metal lines are not so densely packed (Pauldrach et al. 2001). The main method used to include line blocking and blanketing is the opacity sampling method. For this method, the opacities are evaluated for a grid of frequency points. The approximation approaches the exact solution as the code increases the number of frequency points it samples. CoStar incorporates the opacity sampling using ∼ 50Åwide bands in the Monte Carlo radiative transfer solution (Schaerer & Schmutz 1994). Of the three codes considered here, CoStar has the most approximate treatment.
The TLUSTY atmosphere code includes 180,000-200,000 frequency points in the opacity sampling grid (Lanz & Hubeny 2003). The WM-basic atmosphere code is solved in two parts. First, the radiative transfer is solved with a fast approximate treatment that samples ∼ 1, 000 frequency points in the Lyman continuum. This fast solution is repeated iteratively and generates starting values for the final solution. The final solution consists of fewer iterations that solve the radiative transfer with an exact treatment of the line-blanketing (Pauldrach et al. 2001).
In addition to line-blanketing by metal lines, stellar winds change the shape of the ionizing SED at high energies (e.g., Sellmaier et al. 1996). This is particularly important near the 54.4eV ionization potential of Heii (Gabler et al. 1989). Both CoStar and WM-basic calculate the non-LTE solution for expanding atmospheres. In the WM-basic atmosphere code, the standard opacity sampling is modified to take into account line shifts due to the expanding winds. These line shifts effectively increase the frequency range that can be blocked by a given line (Pauldrach et al. 2001).
While SNC02 is generated with the WM-basic code, they use a different set of assumptions for the stellar parameters in calculating the SEDs. Specifically, they determine the appropriate mass-loss rate and terminal velocities a priori using the empirical relations from Prinja et al. (1990), Lamers et al. (1995) and Kudritzki & Puls (2000). To ensure the final atmospheres have these parameters, they manually adjust the radiative acceleration by changing the force multipliers used in the code (Smith et al. 2002). In contrast, the WM-basic grid included in CLOUDY obtains the the mass-loss rate and wind terminal velocity from the solution of the atmosphere code itself (Pauldrach et al. 1998). We note that the SNC02 grid was previously implemented in STAR-BURST99 (Leitherer et al. 1999).

Observations
To evaluate the atmosphere grids described above, we first construct a sample of single-star Hii regions. We use the narrow-band emission-line images from the Magellanic Clouds Emission Line Survey (MCELS; Smith et al. 2005) to select small, spherical Hii regions that are likely ionized by a single star. Table 1 lists the properties of these Hii regions. The first five columns give the MCELS (Pellegrini et al. 2012), DEM (Davies et al. 1976) and stellar designations, and the (J2000) positions. The next three columns list observed properties of the ionizing stars. The ninth column has the observed radius of the Hii region in parsecs as measured from the Hα MCELS image. Figure 1 shows both Hα and three-color composites in [Oiii], [Sii], and Hα of the Hii regions in our sample.
We used the Inamori-Magellan Areal Camera & Spectrograph (IMACS) on the Magellan Baade Telescope at Las Campanas Observatory to obtain both long slit spectra and Bessell B and V images. Our data were collected on the nights of 2008 Jan 29-31 using the f/4 configuration. In this setup, IMACS has an eight-chip mosaic CCD that has a total of 8,000×8,000 15 µm pixels, each of which corresponds to 0.11 . The seeing was good over the observing run, resulting in a final spatial resolution  Davies et al. (1976) b ID from Henize (1956) c Unless otherwise noted, ID from OGLE-III (Udalski et al. 2008) d ID from Massey (2002) e R.A. and Dec correspond to that of the ionizing star. f Measurement error on V is 0.05 mag and the systematic error is discussed in the text. g SpT is inferred but not directly observed, see §3.4. h Magnitude from the OGLE-III (Udalski et al. 2008) of ∼ 1.6 .
The long slit observations consist of 3 × 1200s exposures using a 0.7 slit. We use the 600 l/mm grating, which has a spectral resolution R ∼ 2730 at Hα. The data are binned by 2 and 4 in the dispersion and spatial directions, respectively, which results in ∆ ∼ 0.76Å per pixel. For the first night of the observing run the spectral coverage is 3700-5900Å, while for the nights of 2008 Jan 30-31 the wavelength coverage is 3700 -6740Å due to a different grating tilt.
We use standard IRAF 3 procedures for the data reduction. The spectra are extracted with the IRAF task apextract using separate apertures for the nebula and star. The extraction includes a local background subtraction. We extract the nebular spectra in two or more apertures. These apertures are selected to exclude any stars along the slit, and they otherwise span the entire nebula. We sum those spectra to obtain the total emission from the nebula along the slit. Figures 2 and 3 show the stellar and nebular spectra for the objects in our sample. The intensities are scaled to an arbitrary value for presentation, and the gaps between the chips are assigned a value equal to the continuum level. We flux calibrate the spectra using standard stars LTT 3218, LTT 1788, LTT 2415, LTT 2754and EG 21 (Hamuy et al. 1994. We measure the emission line fluxes from the nebular spectra ( Figure 3) with the IRAF task splot assuming gaussian line profiles. We use the reddening equation: to find the reddening coefficient, c(Hβ).
Here, I(λ)/I(Hβ) and I 0 (λ)/I 0 (Hβ) are the observed and intrinsic Balmer ratios, respectively, and we use the reddening law of Cardelli et al. (1989) to determine f λ and f Hβ . For the reddening solution we set the ratio of total to selective extinction, R V = 3.45, which is appropriate for the LMC (Gordon et al. 2003). The dereddened line strengths, relative to Hβ, and c(Hβ) values for each nebula are listed Table 2. In Table 2, we also include the Hβ flux that we measure in the slit. Note that the slit width is 0.7 relative to the 0.5-2 sizes of the nebulae. Therefore, the Hβ flux in Table 2 is not representative of the Hβ flux in the nebula. We calculate the emission-line measurement error as σ 2 cont N + σ 2 cont EW/∆ (Gonzalez-Delgado et al. 1994). Here, σ cont is the rms in the continuum near the line, N is width of the measured line in pixels, EW is the equivalent width of the line and ∆ is the dispersion of the spectra inÅ per pixel (Gonzalez-Delgado et al. 1994). The reddening error and the flux calibration error together contribute < 10% error for all our sample, and ∼ 5% for most of our sample. These sources of error are combined in quadrature to obtain the flux error listed in Table 2.
We assign spectral types (SpT) to the stars from the rectified stellar spectra ( Figure 2) following the criteria of Walborn & Fitzpatrick (1990). The spectral type assignments are based on the independent spectral typing by four individuals and are accurate to half a spectral type. The earliest SpT is O5.5, and the latest is B0.5. Most of the stars are luminosity class V, but our sample also contains two B supergiants (Table 1).
In addition, we collect single exposures of 20s and 10s in filters B and V, respectively. We obtain photometry from these images using the IRAF task apphot. Ten of our targets are in common with stars in the OGLE-III survey (Udalski et al. 2008). We compare our measured V magnitudes to the OGLE-III magnitudes and find a mean difference of -0.07 mag with a standard deviation of 0.15 mag. We also have two nebulae, MCELS L 43 and L 52, whose stars are in common with the Massey (2002) UBVR survey of the Magellanic Clouds, and we find agreement in V to 0.01 and 0.09 mag respectively. As a further check, we compare the Massey (2002) and OGLE III V magnitudes for ∼ 50 stars with V between   13 and 14 mag. The mean difference (V Massey −V OGLEIII ) and standard deviation are -0.05 and 0.11 mag, respectively. From this we see that our values are bracketed by the published literature values. We obtain absolute V magnitudes, M V , after correcting our observed V for the measured extinction and assuming an LMC distance modulus of 18.48 (Westerlund 1997).
2.2. OB companions One important consideration is the possible contamination from OB companions in our sample. The binary fraction for massive stars is thought to be 40-70% (e.g., Sana et al. 2009Sana et al. , 2011 and possibly as high as 90% (Kiminki & Kobulnicky 2012). Any additional OB stars in the nebula will both contribute to the ionizing photon budget and affect the shape of the ionizing SED. For OB stars, Q 0 changes by 0.1 to 0.6 dex for each full step in SpT (Smith et al. 2002). Thus, the impact of a companion on Q 0 is maximized for equal mass stars and drops rapidly with later SpT.
The available data put some constraints on possible companion stars. For the objects in our sample, the nearest resolved stars have magnitudes that are > 2 mag fainter than our target stars, thus ruling out contributions to Q 0 from companions at large distances from the central star. However, we can still have companions on scales smaller than our spatial resolution, 1.6 which corresponds to ∼ 0.5 pc. In fact, two of our objects, MCELS L 346 and MCELS L 390, are confirmed eclipsing binaries, and we discuss them in detail in §3.4.
For the rest of the sample, the observed magnitudes will include contributions from both stars, if a companion is present. Thus, the luminosities derived from our observed magnitudes should be the sum of the luminosities from all components and will be representative of the system. Furthermore, for a binary or companion, one would expect to have spectral type diagnostics representing a mixture of the two spectral types, except in the case of an equal mass binary. We note that none of our stellar spectra show evidence for composite SpT. If any of our objects are equal mass binaries, Q 0 from the atmosphere models would be around half that of Q 0 from L(Hα). We discuss Q 0 in more detail in §3.3, but we note here that we do not find support for equal mass binaries from the Q 0 comparison. Thus, while we cannot conclusively rule out binaries from the rest of our sample, the data available suggest that any binaries present will not significantly affect our results.

Method
To evaluate the stellar atmosphere models, we compare the Hii regions in our sample to the predictions of photoionization simulations. As discussed in the Introduction, the line emission from Hii regions depends primarily on the metallicity (Z ), the ionization parameter (U ), and the SED of the ionizing star. When Z and U are constrained, differences between the predictions from photoionization simulations and the observed Hii region spectrum can be directly linked to differences in the shape of the model SED and that of the actual ionizing star.
Photoionization calculations are performed with version 08.00 of CLOUDY (Ferland et al. 1998). To set up the CLOUDY simulations, we need to match the nebular abundance and ionization parameter to the observed nebulae. In MCELS L 28 and L 351, we detect the auroral [Oiii] λ4363 line. We use the ratio [Oiii] λ4959, 5007/[Oiii] λ4363 as input for the IRAF task temden to obtain the electron temperature, T e . With T e in hand, we derive the oxygen abundance (Table 2) using the IRAF task ionic. For this calculation we assume T e ([Oii]) = T e ([Oiii]). To set the abundance of S, C, Ne, Ar, Si and Fe, we use the relations from McGaugh (1991), which relate the elemental abundances to that of oxygen. For the rest of our objects, [Oiii] λ4363 falls in the gap between CCD chips. In these cases, we adopt the mean LMC abundances measured by Garnett (1999). For all our objects, except for MCELS L 346, L 351, and DEM L 283b, we calculate the nitrogen abundance from the ratio [Nii]/[Oii], as described in Pérez- Montero & Díaz (2005). Those values are shown in Table 2. We note that except for MCELS L 345, the mean difference between the calculated log(N/O) and the Garnett (1999) values is 0.17 dex with a standard deviation of 0.11 dex. This is comparable to the uncertainty in both the relation (Pérez-Montero & Díaz 2005) and in the measured abundances above.
The ionization parameter is defined by where Q 0 is the rate of H-ionizing photons, c is the speed of light, n e is the electron density, R S is the radius of the Strömgren sphere and is the filling factor of the gas. Essentially, U describes the ionizing photon density relative to the gas density. The ionization parameter will depend on both the type of star producing the ionizing radiation and the gas distribution around the star. As can be seen in Figure 1, the nebulae are Strömgren spheres, which means we can assume simple spherical geometry for the simulated nebulae. The uncertain nebular parameters that control U are those that describe the radial density profile of the nebula: the initial (R inner ) and outer (R nebula ) radii of the cloud, hydrogen density (n H ), and . We use a combination of emission-line diagnostics and the MCELS Hα images to set these values. To constrain the density distribution, we examine the MCELS images and the Hα and [Oiii] λ5007 line profiles along the slit. The [Oiii] λ5007 and Hα profiles indicate that R inner is between 0.1×R nebula and 0.5×R nebula . We measure R nebula from the Hα MCELS images of the nebulae and use those values to set the inner and outer radii in the CLOUDY simulations. To obtain n H , we assume n H = n e , since hydrogen is fully ionized in a typical Hii region. The [Oii] λλ3726, 3729 and [Sii] λλ6716, 6731 densitysensitive doublets both indicate that nearly all our objects are in the low-density regime (n e < 100 cm −3 ) below which these diagnostics are no longer sensitive to the density. We also estimate the density using the Hα emission measure based on the Hα photometry, and find n e 30 cm −3 , which is consistent with our limit from the [Sii] and [Oii] diagnostics. Based on these constraints, we use n H ≤ 100 cm −3 as an upper limit for the densities in all our simulations. Additionally, the Hii regions in our sample are optically thick (Pellegrini et al. 2012). Therefore, we set n H and such that the ionizing photons are absorbed within the observed nebular size (Table 1). Finally all our simulated nebulae include 5 km s −1 turbulence.
In addition to the properties discussed above, dust will affect the nebula. In particular, photoelectric heating from dust can contribute as much as 30% of the total heating (van Hoof et al. 2004). Therefore, we include both graphite and silicate dust in our simulations. We adopt LMC gas distributions and a dust-go-gas ratio, (2001), which were introduced into CLOUDY in Pellegrini et al. (2011). Finally, since our observations are based on long slit spectroscopy, we implement the CLOUDY command slit, which predicts the line ratios that would be observed through a slit across the simulated nebula.
We run CLOUDY simulations with the same set of nebular parameters for each stellar atmosphere model. In addition to specifying which atmosphere grid to use for the source, we must also assign the luminosity and T eff of the star. We calculate stellar luminosity from M V using where M bol, is the solar bolometric magnitude and BC is the the bolometric correction from Martins et al. (2005), except for the case of MCELS L 345, the B0.5 Iab star, where we use the bolometric correction from Crowther et al. (2006). The SpT-T eff conversion is model-dependent and we leave the T eff as a free parameter in our simulations. We note that since the BC is temperature dependent, the input stellar luminosity changes with T eff to maintain consistency. For all our dwarf stars, we assume log(g) = 4.0. We determine the best T eff for each star to be the one that produces the correct balance of ionizing flux at both high and low energies. We accomplish this task by plotting the predicted emission-line flux as a function of ionization potential. An example is shown in Figure 4, in which we plot the results of six simulations that differ in T eff . As expected, as T eff decreases from 43,000 K to 35,000 K, the predictions for [Oii] and [Nii] increase, while the predictions for [Oiii] and [Neiii] decrease. Thus, the overall slope of the points decreases with decreasing T eff . We define the optimal model at the turnover T eff , where the slope is flat. For the models shown in Figure 4, this would be around T eff =39,000 K (yellow circles). We also prioritize models that match the observed Hei lines in our spectrum. The ratio of the Hei to Hi recombination lines reflects the level of He ionization in the nebula. As long as He remains partially ionized, this ratio will be primarily dependent on the SED of the ionizing source (Kennicutt et al. 2000). We note that we do not detect Heii λ4686 in any of our spectra. Therefore, we can set an upper limit on T eff by requiring a non-detection of that line in the simulations as well. Figure 1 shows that our HII regions, to first approximation, are simple Strömgren spheres. We use the method described in §2.3 and the models described in §1.1 to first generate a grid of model Hii regions with uniform densities. Figures 5 and 6 show the results using this prescription for MCELS L 28 and L 43, which are ionized by an O5.5 V and an O8 V star, respectively. These figures are representative of the results from the rest of the sample. The left panels of the figures show [Oiii] λ5007 vs.

Uniform Density Models
[Oii] λ3727, which reflects variations of U. In nebulae with high U there will be more [Oiii] relative to [Oii], and the points will occupy the top left of the plot. As U decreases, [Oii] becomes stronger relative to [Oiii], and the points will move toward the lower right. The right panels show [Oii] λ3727 vs. Hei λ5876 . This panel illustrates how well [Oii] is predicted at the appropriate T eff , which is traced by Hei (Kennicutt et al. 2000). The solid and dotted lines show how the line ratios change when the R inner and T eff are changed, respectively.
It is immediately apparent that these models do not reproduce the observed emission-line ratios. In the simulations using CoStar and WM-basic atmosphere models, the simulation tracks generally run below and to the left of the observed values ( Figure 5; left panel). This indicates that not enough [Oii] is produced when [Oiii] is well predicted and suggests that U is too high. Furthermore, from the righthand plots, we see that simply changing the T eff is not a viable solution; the predicted [Oii] flux does not match the observations for the entire range in T eff . This is the case for all our objects except for MCELS L 43 ( Figure 6) and MCELS L 344 (not shown). In contrast, the simulated nebulae ionized by TLUSTY atmospheres have lower ionization parameters. The softer TLUSTY SED is able to reproduce the observed [Oiii] and [Oii] lines in MCELS L 28 ( Figure 5), MCELS L 351, and MCELS L 394. For the rest of the objects, the TLUSTY SED is too soft, and the simulations have an ionization parameter that is too low (Figure 6).
We rule out the adopted metallicity as the source of the discrepancy between the simulations and observations, as follows. The metallicity of the gas strongly affects the observed line ratios. As the metallicity decreases, there are fewer metals to provide cooling, and the nebula is hotter. Thus, the under-prediction of the low ions seen above could indicate that the adopted metallicities are too low. To explore how this might affect our results, we generate a CLOUDY grid for MCELS L 28 in which we change the metallicity from log ( We find that changing the metallicity will not reconcile the discrepancy in the low ions. Although increasing the metallicity does decrease the predicted [Oiii]/[Oii] ratio for a given T eff (Figure 7(a)), it does so by decreasing [Oiii] (Figure 7(b)). Meanwhile, [Oii] remains under-predicted for the whole range of metallicities (Figure 7(c)). Therefore, even though changing the metallicity might bring the overall [Oiii]/[Oii] closer to the observed, it cannot account for the discrepancies between the simulations and the observations for the individual lines. We note that decreasing log(O/H) by 0.1 dex, from the observed to the mean LMC value, will change the T eff required to reproduce the observed emission-line fluxes by 500 K. This is comparable to the 250-500 K uncertainty in determining T eff using the emission-line ratios. In most of the spectra, the chip gap falls on or slightly red-ward of Hγ. We assign this region with a value equal to the mean continuum level. The best-fitting model is defined as the one with the flattest overall slope, with most points lying between the red lines. In this case, T eff =39,000 K (yellow circles) has the best fit.
Since we do not have measured metallicities for most of the objects in our sample, this 500 K is included in the error on the best fitting T eff , discussed below.
Another possible cause of the general under-prediction of the low ions is the nebular structure. In this case, there is not enough dense gas in our simulations at large radii to receive a diluted radiation field, which is needed to reproduce the emission from lines with low ionization potential. This could be because R inner is set too close to the star, or because the nebula has a clumpy distribution of denser gas. An increase of R inner increases the [Oii] and [Nii] flux, while it simultaneously suppresses the [Oiii] and [Neiii] flux. At face value, this would move the predictions closer to the observations. However, even if we vary the inner radius from 0.1 to 0.5 R nebula , we still are not able to simultaneously reproduce the two sets of lines ( Figure 5). Furthermore, the inner radius is constrained to be 0.25 R nebula in all our objects, except DEM L 283b, based on the spatially resolved [Oiii] λ5007 line profiles and the MCELS images. Thus, to within the ∼0.5 pc spatial resolution of our spectral data, these Hii regions do not appear to be hollow shells in projection. Therefore, neither the metallicity nor R inner can explain the discrepancy in these uniform density models. In the next section, we explore the effects of a clumpy gas distribution.

Density Fluctuations
In Hii regions, non-uniform gas and temperature distributions are often invoked to explain the observed properties. In fact, direct imaging of Hii regions reveals complex morphology, which includes density gradients and/or clumps and filaments with high density contrast. Non-uniform gas density has long been invoked to explain observed parameters (e.g., Osterbrock & Flather 1959). One example is the long-standing problem reconciling the discrepant values for T e that result from different di-  agnostics. This has been of particular concern because the uncertainty in T e translates into significant uncertainty in the derived abundances. One suggestion is that this problem arises from temperature fluctuations in the ionized gas, which are described by the t 2 parameter (Peimbert 1967). These fluctuations could arise from small clumps of dense gas within the nebulae (e.g., Viegas & Clegg 1994;Williams 1992;Liu & Danziger 1993) or from fluctuations in chemical abundance rather than density  (Rubin 1989;Kingdon & Ferland 1995;Giammanco et al. 2004;Tsamis & Péquignot 2005).
We therefore explore how changing the density structure affects line emission by considering a clumpy medium. For our analysis here, we assume that dense clumps of gas are uniformly spread throughout the nebula, with a vacuum between them (Osterbrock & Flather 1959). In reality, the clumps will be interspersed with lower density, diffuse gas. However, the emission measure is proportional to n 2 e . Therefore, most of the observed flux is coming from the dense clumps, and we can treat the inter-clump space as though it were a vacuum. In the simulations, these clumps are treated in a statistical sense; both the volume emissivity of the gas and the optical depth along the line of sight are decreased by in each zone calculated by the code (Ferland et al. 1998).
In our final simulation grid, T eff and are varied, with the density chosen as discussed above. A decrease in the filling factor of the gas necessitates an increase in the density of that gas to maintain the nebular radius (Equation 2). In this scenario, dense gas is present at small radii to experience a strong radiation field, as required for the lines with high ionization potential. At the same time, more dense gas is present at larger radii from the star that receives a diluted radiation field to reproduce ions with low ionization potential. Our constraints on the sizes of the simulated nebulae, combined with our upper limit on the density, effectively provide us with lower limits for n H and . We cannot decrease the density past the value that leads to correctly sized nebulae when = 1.0, and we cannot decrease below the value that corresponds to n H = 100 cm −3 . We note that the filling factor set by the density limit, ∼ 0.10 in most of our objects, is consistent with the 10 −1 −10 −3 commonly used in the literature (Kennicutt 1984;Giammanco et al. 2004;Hunt & Hirashita 2009). Furthermore, the nebulae in our sample exhibit substructure (Figure 1) that is consistent with this level of density fluctuations. Table 3 shows the parameters of the best fitting model for each nebula. These parameters include R inner and log(g) for each nebula and star, as well as the T eff , , and n H that correspond to the best simulation for each atmosphere model. Figure 8 shows the results of these simulations. In this Figure, we plot the ratio of the predicted emission-line flux to the observed, F pre /F obs , as a function of the ionization potential of the emission line. Each model atmosphere corresponds to a different color and symbol. The error bars represent the observational error, which includes the line measurement, reddening correction and flux calibration errors. As discussed in §2.3, these simulations are selected by finding the combination of free parameters ( , T eff ) that reproduces the observed Hei emission-lines, and yields the flattest slope with most points lying within ±20% of unity, denoted by the red horizontal lines. The 20% takes into account the observational variance, which is estimated from MCELS L 28, L 32, L 35 and L 344, for which we have multiple nights of data and/or multiple slit positions (Figures 8(a), 8(i), 8(j), and 8(e)). The additive error on our selected T eff and for a given R inner is ∼ 500 − 1000 K and 0.1, respectively. The error for T eff includes the contribution based on the uncertainty in our metallicity assumptions, as discussed in §3.1. While they are included in Figure 8 for completeness, we find that the [Oi] λ6300 and [Sii] λ6716, 6731 emission lines are not good diagnostics for distinguishing between the SEDs. These lines are primarily dependent on nebular conditions, and less so on the particular atmosphere chosen. This can be seen in Figure 9,   tion potential of this line is 10.60 eV, which is below that of H. In this regime, the differences between atmosphere models are small, and so the line will be more sensitive to changes in the gas density and morphology.
[Oi] λ6300, on the other hand, is produced via charge-exchange and, in the conditions of a typical nebula, depends critically on the density of H + (Osterbrock & Ferland 2006).

Rate of Ionizing Photons
One important comparison is between the rate of ionizing photons produced by the observed nebulae, Q 0,Hα , and that of the best fitting atmosphere model, Q 0,SED . To calculate Q 0,Hα , we obtain Hα photometry from the MCELS, Hα emission-line image ( Table 1). The rate of ionizing photons calculated from this photometry is listed in column 2 of Table 4. We note that the MCELS narrowband images are not continuum subtracted. As a result, the error on Q 0,Hα ∼ 20% (Pellegrini et al. 2012). Columns 3-6 show the ratio of Q 0,SED and the measured Q 0,Hα . It is immediately apparent that the Q 0,SED for dwarf stars is offset from Q 0,Hα by +5, -12, -12, and -37 % on average for TLUSTY, WM-basic, SNC02, and CoStar, respectively. This excludes MCELS L394, for which Q 0,SED is 3-5 times higher than Q 0,Hα (see §3.4), and the two supergiants. From this diagnostic, TLUSTY, best represents the ionizing source. However, both WM- Figure 9.
[Oi] λ6300 and [Sii] λ6716 vs Hei λ6876 for MCELS L 28. Dotted lines are models with constant R inner =0.10, 0.25, and 0.50 R nebula and solid lines have constant T eff = 35, 37, 39 kK, except WM-basic, which also includes T eff = 40 kK. The values of R inner increase going upward and T eff increases from left to right. basic and SNC02 also do well, particularly within the 20% error on the Hα photometry. We note that the offset does depend on atmosphere model, but we find no trend with spectral type or A V .
One important factor that affects this comparison is the assumed gravity. In our simulations, we use log(g)= 4.0 for all luminosity class V stars. However, these stars may have log(g) in the range from 3.8 to 4.1 (e.g., Herrero et al. 1999;Massey et al. 2005;Martins et al. 2005). Decreasing log(g) by just 0.1 dex can increase Q 0,SED by as much as a factor of 2. This change more than compensates for any discrepancy shown in Table 4, and supports that Q 0,SED is consistent with Q 0,Hα , within the uncertainty.
We note that T eff will affect Q 0,SED . Both the bolometric correction and the shape of the SED play a role in determining Q 0 , and both properties depend on T eff . As noted in §3.2, the error on T eff is between 500 and 1000 K. If T eff is too low by 1000 K, Q 0,SED would increase by 10-15%. However, a change in T eff of 1000 K significantly affects the ionization structure of the nebula, which worsens many of the fits shown in Figure 8. Thus, uncertainty in T eff does not significantly contribute to any offset in Q 0,SED /Q 0,Hα .
In addition to the general offsets in Q 0,SED /Q 0,Hα , we also find differences between model grids. In general, TLUSTY has the highest Q 0,SED /Q 0,Hα , with most values near unity, and CoStar has the lowest, with many values near 0.6. WM-basic and SNC02 fall between them. The trend between atmosphere models reflects the important role that the relative hardness of the SED plays in the temperature selection process, which in turn affects Q 0 . For a given effective temperature, different atmosphere models will produce different Q 0 (e.g., Voges et al. 2008;Simón-Díaz & Stasińska 2008). On top of that, we are sensitive to the shape of the ionizing SED because our atmosphere-model selection criteria depend on matching the emission-line ratios. This means that for a softer atmosphere model (e.g., TLUSTY), we select simulations with higher effective temperatures than we would for simulations using harder atmosphere models (e.g., CoStar). This effect leads to the trend seen above: the best-fitting models for TLUSTY have higher Q 0 than those of CoStar. Therefore, while the systematic offset in Q 0 is most likely caused by uncertainty in log(g), it may also suggest that the WM-basic, SNC02, and CoStar atmosphere models are generally harder than the observed stars. However, we note that this is opposite from what we find based on the [Neiii] line, as we discuss in §4.1.

Individual Objects
In this section we briefly comment on the individual objects in our sample. Unless otherwise noted, we assume the mean LMC abundances from Garnett (1999). Table  2 lists the Helium abundances, which we calculate from the Hei lines using the relations from Benjamin et al. (2002). As a reference, Table 3 shows the parameters of our best fitting models.
MCELS L 28 is ionized by OGLE-III 8225 (Udalski et al. 2008), which is the earliest spectral type in our sample, an O5.5 V star (Figures 1(a) and 8(a)). In this nebula, we detect the auroral [Oiii] λ4363 emission line, from which we derive log(O/H) = −3.51. The simula-  Figure 8(a) clearly shows this is not the case in the = 0.10 models. If we had selected the = 1 model as our final best fit, the corresponding T eff would be 39,000 K. This is 3,500 K cooler than the adopted T eff =42,500 K selected with the lower filling factor. This T eff difference highlights the importance of using as many emission-line diagnostics as available to determine the fits. The best fitting T eff are 39500, 42500, 41500, and 41500 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively. MCELS L 52 is ionized by M2002 19696 (Massey 2002), an O6.5 V star (Figures 1(f) and 8(b)). We find excellent agreement for lines with ionization potential (IP) below 40 eV, and all the atmosphere models produce consistent results. For [Neiii], whose IP = 40.96, this agreement does not hold. The simulation using CoStar, the hardest atmosphere, over-predicts [Neiii], and TLUSTY under-predicts it. Interestingly, we find somewhat different nebular parameters produce the best fitting simulations between CoStar and the other atmosphere models. For CoStar, a slightly higher filling factor and lower density, 0.25 and 25 cm −3 , respectively, work best. Whereas, = 0.10 and n H = 60 cm −3 produce better fits with the other atmosphere models. We note that this difference is close to our error in determining of 0.10, and, thus, might not be significant. The best fitting T eff are 36000, 38500, 38750, and 38500 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 351 is ionized by OGLE-III 43846, which is also an O6.5 V star and, with a radius of 11 pc, it is the largest nebula in our sample (Figures 1(b) and 8(c)). We detect [Oiii] λ4363 in this object, from which we calculate T e and log(O/H) (Table 2). However, our [Oiii] λ4363 detection is at a very poor signal-to-noise, and our result is consistent with the mean LMC metallicity. Therefore, we use the mean LMC metallicity in the CLOUDY simulations. As for MCELS L 52, we find excellent agreement between simulations and observations for emission lines with IP < 40 eV. The best fitting T eff are 36500, 39500, 39000, and 39000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively. DEM L 283b is ionized by OGLE-III 44989, an O6.5 V((f)) star (Figures 1(c) and 8(d)). Since we do not have IMACS photometry for this object, we use the OGLE-III V -band magnitude to calculate the stellar luminosity. From the [Oiii] λ5007 line profile, we assign R inner = 0.25 R nebula . Overall, we find good fits for all atmosphere models for lines with IP < 40 eV and = 0.10. The best fitting T eff are 36500, 39000, 39000, and 39000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 344 is ionized by OGLE-III 45830, an O6.5-7 V((f)) star (Figures 1(d) and 8(e)). We note that the nebular structure diverges from a Strömgren sphere to the east in Figure 1(d). However, our slit position runs N-S and covers the regular portion of the nebula. Overall we find excellent fits for models where IP < 40 eV. The measured [Neiii] flux shows a 45% difference between the spectra from 2008 Jan 29 and 31, which results in a large scatter for [Neiii] in Figure 8(e). This apparent discrepancy is caused by the large observational uncertainty. Interestingly, for Jan 31, the simulation using CoStar is the only one to match the observations, and the others under-predict [Neiii] by > 50%. The best fitting T eff are 37000, 39500, 39000, and 39000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 43 is ionized by M2002 17251, an O8 V star (Figures 1(e) and 8(f)). As noted in §3.1, MCELS L 43 is one of the few objects for which [Oiii]/[Oii] is well predicted by simulations using = 1 with atmosphere models from WMbasic, SNC02 and CoStar. However, the right panel of Figure 6 shows that these models cannot simultaneously match the Hei lines, which are a more reliable T eff indicator. At the temperatures required to match Hei, [Oiii] λ5007 is over predicted by a factor of 2-3. The fit improves in all emission lines when we decrease the filling factor to 0.10 (Figure 8(f)). The best fitting T eff are 37000, 39000, 39000, and 39000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 390 is ionized by OGLE-III 28307, an O9 V star (Figures 1(g) and 8(g)). For this object, the simulation using CoStar reproduces the relative amounts of [Oiii] and [Neiii], while the others do not. This is in contrast to the rest of our sample, for which the CoStar models are too hard. This star is an eclipsing binary, according to the OGLE-III catalog of eclipsing binaries (Graczyk et al. 2011). However, the second star only contributes 10% of the I band flux. Futhermore, Table 4 shows agreement between Q 0,SED and Q 0,Hα . Therefore, the second star is not significant. The best fitting T eff are 35000, 35000, 37000, and 37000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 394 is ionized by OGLE-III 28239, a star whose spectral type is between an O9 and a B1.5 V, based on the noisy stellar spectrum (Figure 1(g) and 8(h)). However, we assign this star the earlier SpT in this range because the observed [Neiii] rules out stars with SpT of B0 V or later. Even with = 0.10, it is still challenging to reproduce the [Sii], [Oii] and [Nii] lines. As seen in Figure  8(h) the predictions lie at the edge of our 20% tolerance limit. Additionally, all the atmosphere models, except for CoStar, do not produce enough flux at the energies needed to reproduce the observed [Neiii] emission-line (Figure 8(h)). The best fitting T eff are 37000, 39000, 39000, and 39000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
Unlike the other objects, Q 0,Hα in MCELS L 394 is a factor of a few less than Q 0,SED (Table 4). This could indicate that the nebula is optically thin. If the nebula were optically thin, we expect to observe less [Oii] and [Sii] than if the nebula were optically thick. Based on our modeling, however, we see the opposite. We observe too much [Oii] and [Sii] relative to the optically thick model. Therefore, while the Hα photometry may suggest an optically thin nebula, it is clear that this does not adequately explain all the observations. MCELS L 32 is ionized by OGLE-III 8229, a B0 V star (Figure 1(a) and 8(i)). [Neiii] is not detected in this nebula, and thus provides a firm upper limit for the T eff of the ionizing star. Both of the grids generated with the WM-basic code include atmosphere models down to T eff = 30,000 K, which is adequate for this B0 star. However, the CoStar grid in CLOUDY does not have dwarf B star models. The coolest CoStar dwarf, T eff =34,000 K, is still a bit too hard to match the observed emission lines. We optimized the T eff to reproduce the Hei lines, with the result that [Oiii] is under-predicted. The best fitting T eff are 34000, 33000, 34000, and 34000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.
MCELS L 35 is ionized by OGLE-III 8203, one of our two supergiants, a B0 Ib star (Figure 1(a) and 8(j)). In addition to varying the T eff and , we vary the surface gravity between log(g) = 2.5 − 3.5, as appropriate for supergiants (Crowther et al. 2006;Trundle et al. 2004). TLUSTY is the only atmosphere code which includes this whole range of gravities and temperatures. The lowest gravity in the available grids for WM-basic and SNC02 are 3.0 and 2.95, respectively. In the case of CoStar, the lowest available gravity for T eff >26,000 K is log(g) = 2.86. For TLUSTY, we find that simulations using log(g) = 2.6 − 3.0 all yield good fits, although at different temperatures. For example, for log(g) = 2.6, the best fitting T eff = 27,250 K, whereas for log(g) = 3.0, T eff = 28,500 K. For WM-basic, SNC02 and CoStar we use log(g) = 3.5 and find T eff = 31,000 K for all three atmospheres.
MCELS L 345 is ionized by OGLE-III 50093, the other supergiant in our sample, a B0.5 Iab star ( Figure   1(d) and 8(k)). Based on the [Nii]/[Oii] ratio, this object has a much higher nitrogen abundance than the rest of our sample, log(N/O) = −0.77 compared to the mean LMC value of log(N/O) = −1.5. This higher abundance is confirmed by the simulations because those that use the mean value are unable to reproduce the [Nii]/[Oii] line ratio.
[Oiii] λ5007 is not detected in this object, which places an upper limit of T eff =26,000 K. As discussed for MCELS L 32, we are only able to test the TLUSTY atmosphere models for this object. Note, the points along the bottom of Figure 8(k) are lines for which we only have upper limits on the line strength, and the simulations are consistent with the data. We find the best fitting T eff falls in the range of 22,000 -24,000 K.
MCELS L 346 is ionized by an eclipsing binary, MA-CHO 81.9725.16 (Alcock et al. 1997). The spectral type determined from our observed data is B1-1.5 V ( Figure  2). However, we do detect [Neiii] in this object, despite the late SpT. The strengths of the [Neiii] and [Oiii] lines suggest that this nebula is ionized by a star of much earlier spectral type. Based on the [Oiii]/Hβ and [Neiii]/Hβ line ratios compared to the rest of our sample, the ionizing star in this nebula likely has a SpT between O8 V and O9.5 V. This analysis is further supported by the measured log(Q 0,Hα )=48.48, which falls into the [47.99-48.65] range spanned by the O8 -9 V stars in our sample. Considering the short period of this binary, 1.2 days, it is possible that we observed the nebula when the earlier type star was eclipsed by a later type companion. Since the early type star will dominate the ionizing SED, we model this object with a single atmosphere model.
We note that one side of this nebula has an irregular morphology, while the other is more Strömgren-like (Figure 1(h)). Therefore, we only use the emission-line ratios from the portion of the slit that lies along the regular half of the nebula. The best fitting models are shown in Figure 8(l). These models have R inner = 0.25 R nebula , = 0.10, and n H = 100 cm −3 . We find excellent fits for WM-basic, TLUSTY, and SNC02 except at [Neiii] λ3869, while CoStar over-predicts both [Oiii] and [Neiii]. The best fitting T eff are 36,000, 37,500, 38,000, and 37,000 K for CoStar, TLUSTY, WM-basic, and SNC02, respectively.

Comparison of Different Atmosphere Models
In the previous section we discussed the simulations for each object individually. However, the question at hand is: how well do the atmosphere models reproduce the ionizing population, in general? In this section, we discuss the trends we find across our sample. To look at the results collectively, we plot F pre /F obs for individual emission lines as a function of the spectral type of the ionizing star. These plots are shown in Figure 10.
With the exception of [Neiii] λ3869 (Figure 10(d)), Hii region simulations using CoStar are consistent with the other atmosphere models and the observations ( Figure  10). For high ionization potential lines, such as [Neiii] shown in Figure 10(d), simulations with the CoStar atmospheres over-produce the line emission. Other studies have similarly found that the CoStar SED is too hard at high energies (e.g., Oey et al. 2000;Morisset et al. 2004;Mokiem et al. 2004). Additionally, we find that a small range of T eff reproduces the observed nebular emission for most of our range in SpT (Figure 12). This implies that the SEDs generated by CoStar are less sensitive to changes in the atmospheres between stars of different SpT.
As with CoStar, the TLUSTY models reproduce most of the emission lines up to the ionization potential of Ne +2 . At the IP of Ne +2 , 41 eV, the modeled Hii region spectrum under-predicts the line emission, and it does so to a greater degree than the other atmosphere models. This suggests that the TLUSTY atmosphere is not producing enough flux at high energies and reflects differences in the incident SEDs, as shown in Figure 11. Figure 11 compares the SEDs produced by different stellar atmosphere codes for two temperatures, T eff =35,000 K and T eff =41,000 K. At both T eff , the high energy slope of the TLUSTY SED is steeper than the other atmosphere models. Similarly, Simón-Díaz & Stasińska (2008) show that for ionizing SEDs with 30, 000 ≤ T eff ≤ 40, 000 K, TLUSTY produces less flux at high energies than model atmospheres with expanding winds. Using IR emission lines, Morisset et al. (2004) also find that plane-parallel model atmospheres, TLUSTY included, are too soft between 27 eV and 41 eV to reproduce the observed nebular emission. This confirms that the treatment of the expanding wind plays an important role in reproducing the flux of ionizing photons with higher energies. Gabler et al. (1989) showed that self-consistent treatment of the wind was particularly important for the continuum near the Heii ionizing edge at 54.4 eV. It was also shown by Sellmaier et al. (1996) that the [Neiii] problem, in which photoionization models consistently under-predicted [Neiii] λ3869 (e.g., Simpson et al. 1995Simpson et al. , 1986Mathis et al. 1985;Rubin et al. 1991), is not an issue in simulations using non-LTE atmospheres that have a wind extension. We note that while there is a general trend for [Neiii] to be under predicted, the scatter is large and the degree to which [Neiii] is under predicted in our TLUSTY simulations is much less than that seen for earlier generations of stellar atmosphere models.
At the opposite energy regime, simulations using TLUSTY yield more [Oii] emission than the other atmosphere models, and typically over-predict [Oii] 3727Å by 15-30% (Figure 10(a)). While TLUSTY is softer (Simón-Díaz & Stasińska 2008), it is surprising that this trend is not seen in the other lines with similar ionization potential, such as [Nii] λ6584 (Figure 10(b)). One explanation for this could be the relative shapes of the SEDs above the ionization potentials of N +2 (29.6 eV) and O +2 (35.1 eV), since the flux in that energy range controls the ionization from O + to O +2 and N + to N +2 . For energies greater than 35 eV, the luminosity of the TLUSTY SED drops relative to other models, at least for T eff < 40, 000 K (Simón-Díaz & Stasińska 2008). Thus, the slight excess of [Oii] may result from the softer SED.
Simulations using atmospheres generated with the WM-basic code are the closest to the observed values for most emission lines (Figure 10a-c). However, Figure  10(d) shows that the predictions for [Neiii] λ3869 have more scatter than what is seen for the emission lines of lower IP. The ratio of predicted and observed flux for [Neiii] ranges from 0.7 to 1.7, with most points underpredicted by 20-30%. As expected, the atmosphere mod-els presented in the SNC02 grid are consistent with the WM-basic grid implemented in CLOUDY, although we find that they occasionally require a slightly different T eff . This is likely due to the different stellar properties (stellar radius, mass loss rate and wind terminal velocity) used in generating the atmosphere model grid, as discussed in §1.1. We note that the temperature differences, when they occur, are in the range of 250-750 K. Although the predictions from the two grids are close, the WM-basic grid does slightly better than the SNC02 grid in reproducing [Oiii] for late SpT stars.

Comparison between SpT and T eff
The calibration between spectral type and effective temperature is an important consequence of the atmosphere modeling. Ordinarily, the calibration is based on careful fitting of photospheric absorption lines in model SEDs to reproduce high-resolution stellar spectra. Therefore, the SpT-T eff calibration depends on the stellar atmosphere models used in the fitting. This bias is especially apparent when one compares the calibrations based on atmospheres codes without line blanketing (e.g., Vacca et al. 1996) to calibrations based on atmospheres with line blanketing (e.g., Martins et al. 2005;Massey et al. 2005). The latter calibrations assign effective temperatures that are several thousand degrees cooler than the previous calibrations.
In contrast to this method, the effective temperatures from Table 3 depend on reproducing the ionization balance in the nebula that the star ionizes. The points in Figure 12 show the effective temperatures that we generate as a function of the SpTs from Table 1. For comparison, the lines show calibrations from the literature for dwarf stars. The typical error on the literature calibrations is ∼ 1000 K (Martins et al. 2005). Figure 12 does not show MCELS L 346 since the SpT of the ionizing star is uncertain. MCELS L 394, for which the stellar spectrum is noisy, is placed at SpT = O9.5 with a large error bar. Based on the nebular spectrum, however, this object is probably at the upper end of the shown error bar.
We find that, in general, WM-basic, SNC02, and TLUSTY predict similar T eff , while CoStar predicts cooler T eff . This trend is consistent with our findings in §3.3; the CoStar SEDs are harder than the other atmosphere models, and therefore require a lower T eff to reproduce the observed emission-line spectra. For the earliest stars, the best fitting T eff selected for TLUSTY is slightly hotter than that of WM-basic and SNC02. At late SpT we see a trend toward the opposite. Since the treatment of the winds is the biggest difference between the TLUSTY and WM-basic atmospheres, this may reflect changes in the wind properties toward late SpT stars.
We note that the two supergiants in our sample have T eff between 24 kK and 31 kK. This large range spans, and is consistent with, the expected T eff for Magellanic class I stars 25-32 kK (e.g., Mokiem et al. 2007;Trundle et al. 2004).
The effective temperatures presented in Figure 12 are consistent with the relations in the literature, within the errors. However, we note that the points appear to follow a shallower slope. This is primarily caused by the stars of later spectral type. The O5.5 -O7 V stars follow the LMC calibration of Massey et al. (2005)  line), but the later spectral types prefer the hotter temperatures of the Vacca et al. (1996) calibration (dotted line).
One contributing factor to the T eff trend is the assumed log(g). We calculate an expected variation in log(g) for dwarf stars across our spectral range of 0.05 dex, based on the modeled stellar parameters from Schaerer et al. (1993). If there is an error in our spectroscopically determined luminosity class, it will affect the best fitting T eff selected. Specifically, if any of the stars classified as dwarfs are actually giants or supergiants, the T eff will be 1000-3000 K cooler than the T eff chosen for a dwarf star. The difference between the log(g) appropriate for dwarfs and supergiants increases for later SpT. Therefore, changes to the SED caused by uncertainty in log(g) will be strongest for stars of late spectral type (e.g., Martins et al. 2005). We find 0-1000 K decrease in T eff for MCELS L 28, an O5.5 V, going from log(g)=4.0 to log(g)=3.0, but a 2000-3000 K decrease in T eff for MCELS L 20, an O8 V.
Another consideration is the metallicity for which these SpT-T eff calibrations are valid, since lower metallicity stars will have higher T eff (Massey et al. 2005;Mokiem et al. 2004). All of the literature calibrations, except for the thick-dashed line, are calibrations that use solar metallicity stars. Thus, we expect most of these calibrations to lie at cooler T eff than the points in Figure  12. While this is the case for the calibrations based on line blanketed models, which we see by comparing the thick-dashed line to the thin-dashed and solid lines, we note the scatter between the points is larger than the difference in T eff between the LMC and the MW calibrations. Furthermore, the scatter is comparable to the difference between calibrations based on atmospheres with and without line blanketing.
The conventional approach to obtaining the SpT-T eff calibration is based on photospheric lines and is sensitive to the SED at the photosphere. Our method, however, depends on the ionization balance of the nebula and is also sensitive to changes in the SED due to layers outside the photosphere. Thus, SpT-T eff relations derived from nebular data have the potential to probe factors in the upper layers of the atmosphere models that do not match the conditions in the ionizing stars. Our sample is small, and as seen earlier, there are a variety of uncertainties, limiting its usefulness in this regard, but with larger samples and better data, this approach promises an alternative and important way to refine SpT-T eff calibrations.

CONCLUSIONS
Diagnostics from Hii region emission-line spectra are used extensively to determine the physical conditions of the ISM of galaxies, as well as understand their stellar populations and their chemical evolution. The shape of the ionizing SED is one of the primary properties that determines the Hii region emission-line spectrum. However, the intrinsic shapes of massive star SEDs are uncertain, and we depend on atmosphere models to describe them. It is therefore crucial to understand and quantify how well massive star atmosphere models represent the ionizing sources.
In this study, we compare a sample of single-star Hii regions from the nearby LMC to photoionization simulations ionized by atmosphere models from widelyavailable CoStar, TLUSTY, and WM-basic atmosphere grids. We select our sample of small, Strömgren spherelike Hii regions from narrowband MCELS images. We obtain spectra and photometry of this sample using IMACS on the Magellan Baade Telescope at Las Campanas Observatory. We assign spectral types to the ionizing stars in our sample and find a range from O5.5 V to B0 V with 10 of our 12 Hii regions ionized by dwarf stars. While we choose Hii regions that are likely to be ionized by a single star, contamination by OB companions cannot be entirely ruled out. Two of the objects in our sample, MCELS L 346 and MCELS L 390, are confirmed eclipsing binaries. MCELS L 346 contains a nearly equal mass binary, and is modeled with two ionizing stars. The rest of our sample show no obvious signs of OB companions, such as composite spectral types or large discrepancies between Q 0,Hα and Q 0,SED , and as discussed in §2.2 the binary status will not affect our results.
We evaluate the stellar atmosphere models by matching the emission-line spectra predicted by CLOUDY photoionization simulations with that observed in the nebulae. Our nebular simulations show that uniform Strömgren spheres cannot reproduce the observations. In these cases, [Sii], [Oii] and [Nii] are under-predicted, which suggests that the ionization parameter is too high in these simulations. We explore the effects of changing R inner and Z, but find that neither parameter can resolve this issue. Additionally, we find that [Sii] and [Oi] depend strongly on the nebular parameters, but only weakly on the particular SED used in the simulation. In contrast to the uniform Strömgren spheres, simulations that assume a clumpy medium reproduce the observed data. The clumpy medium is described by specifying a filling factor of dense gas and assuming that the remaining volume is a vacuum. In general, we find that a filling factor of ∼ 0.10 produces simulations that agree with observations.
We compare the rate of ionizing photons predicted by the best fitting atmosphere, assuming a constant log(g), with the rate calculated from L(Hα). Within the errors due to uncertainty in log(g), which can be as high as a factor of two, we find that the two values are consistent. However, we do find a systematic offset between the rate of ionizing photons predicted by different atmosphere models. TLUSTY produces the highest Q 0 , while WM-basic and SNC02 have slightly lower rates. In contrast, the Q 0 from CoStar is a factor of two lower than the other atmosphere models. The hardness of the CoStar SED can explain this offset. We select the best-fitting atmosphere models by matching the ionization balance in the nebula. Since the CoStar models have harder SEDs, the best-fitting CoStar atmospheres are a few thousand degrees cooler and produce fewer ionizing photons as a result.
We evaluate the stellar atmosphere models by comparing the predicted to observed emission line ratios. For emission lines with IP ≤ 36 eV, we find that simulations ionized by WM-basic and SNC02 atmospheres, reproduce the observations. Simulations using TLUSTY are also consistent, with the exception of [Oii] λ3727. For this line, TLUSTY over-predicts the amount of [Oii] (a) T eff =41,000 K (b) T eff =35,000 K Figure 11. SEDs for stars with T eff = 41,000 K and log(L bol ) = 37.09 (top) and T eff = 35,000 K and log(L bol ) = 37.78 (bottom) for each stellar atmosphere model. Figure 12. A comparison of the T eff from this paper to different SpT-T eff calibrations. Symbols represent temperatures derived from this work, while the lines are SpT-T eff calibrations for dwarf stars in the literature. Dash-dotted, solid, thin-and thick-dashed, and dotted correspond to calibrations from Conti (1973), Martins et al. (2005), Massey et al. (2005), and Vacca et al. (1996) respectively. We note that the thick dashed line is a calibration for the LMC while the rest are derived from Galactic stars. The hollow symbols correspond to the luminosity class I stars.
λ3727 in the nebula. Surprisingly, despite having the most approximate treatment of line blanketing, simulations using CoStar reproduce the emission lines up to [O III], at 35.1 eV, but they do so at lower T eff than the other model atmospheres.
At energies > 36 eV, the predictions from different atmosphere models diverge from each other and the observations. We find significant scatter between the predictions and observations of [Neiii] λ3869 from all atmosphere models. CoStar, the hardest atmosphere overpredicts the [Neiii] line by more than a factor of two in most of our objects. In contrast, the other atmosphere models range from over-to under-predicting [Neiii] by 70%. Simulations ionized by TLUSTY, with its softer SED and plane-parallel geometry, systematically produce less [Neiii] than simulations using WM-basic and SNC02.
Finally, we compare the SpT-T eff calibrations in the literature with the best fitting T eff from the simulations in this study. Our results for stars with SpT earlier than O8 V fall along the LMC calibration of Massey et al. (2005), while stars later than O8 V seem to be hotter than the calibrations predict. The SpT-T eff calibrations from the literature are based on photospheric lines, while the effective temperatures obtained from this work depend the ionization balance of the nebulae. The comparison of these two methods has the potential to reflect differences between the outer atmospheres of the models and the actual ionizing stars. However, there are a variety of uncertainties and our sample is small, which limits its usefulness in this regard. Further study with a larger sample and detailed stellar modeling is needed to determine if the deviations at late SpT are significant.
This work uses single-star HII regions to test the predictions of stellar atmosphere models. By using singlestar Strömgren spheres, we dramatically reduce the free parameters involved in evaluating atmosphere models with observed nebulae. Future directions that could improve on the groundwork presented here include leveraging a larger range in T eff by including stars earlier than O5.5 V, using a detailed fit to high resolution stellar spectra to generate the input model SEDs, and evaluating model atmospheres at different metallicities.