A Microscale Chemical Transport Model Simulation of an Ozone Episode in Detroit, Michigan

1. Introduction

The United States (U.S.) National Ambient Air Quality Standard (NAAQS) for ozone (O₃) is based on the maximum daily eight-hour average (8-h) O₃ concentration (MDA8 O₃) measured at a monitoring station, and is currently set by the U.S. Environmental Protection Agency (USEPA, Washington, DC, U.S.) at 70 parts per billion by volume (ppb) [1]. Attainment of the NAAQS is based on the three-year average of the annual 4th highest MDA8 O₃, referred to as the design value. Areas with design values that exceed the NAAQS at any station may be designated by the USEPA as in nonattainment of the O₃ standard within several categories of increasing seriousness from marginal to severe. In 2023, Southeast Michigan (SEMI) in the U.S. Great Lakes region, including the City of Detroit, was re-designated by the USEPA from moderate O₃ nonattainment to attainment based on a clean data determination for the three-year period of 2020-2022. This determination was partly based on an exceptional event demonstration by the State of Michigan attributing two days of high O₃ values in June 2022 to a Canadian wildfire.

SEMI’s attainment status could be impacted by a regulatory monitoring site located at East 7 Mile Rd in the City of Detroit. Fine-scale O₃ modeling for the Detroit area was performed with the Microscale Forward and Adjoint Chemical Transport (MicroFACT) model [2] to understand the impacts of local emission sources on ambient O₃ at East 7 Mile Rd station. The results may help to guide the selection of emission controls or other measures to attain the O₃ NAAQS in any future required State Implementation Plan.

This study examines an O₃ episode on 30 June 2022. During this episode, winds were predominantly from the southwest (the typical wind direction for MDA8 O₃ above 70 ppb in SEMI), with no forest fire smoke or lake breeze influence as determined from satellite imagery. Figure 1 shows twelve-hour (12-h) back trajectories arriving at or 500 m above the East 7 Mile Rd station at 14:00 Local Standard Time (LST), 30 June 2022 generated by the Hybrid Single-Particle Lagrangian Integrated Trajectory (HYSPLIT) model [3]. These trajectories indicate two things: 1) the potential for regional transport of O₃ and its precursors from industrialized areas not only in Michigan, but also in the neighboring states of Ohio and Indiana; and 2) the likelihood of elevated layers of enhanced O₃ due to the greater distance traveled over the industrialized U.S. Midwest by air parcels aloft relative to air parcels near the surface. MDA8 O₃ at East 7 Mile Rd station on 30 June 2022 was 76 ppb. While this concentration was the highest for the year in SEMI, it is not regulatorily significant because design values are based on the 4th highest MDA8 O₃. Nevertheless, understanding the physical mechanism behind the high O₃ measured during this episode may reveal insights useful to the design of future O₃ control strategies, if needed.

The EGLE monitoring network in SEMI includes the Allen Park station located to the southwest of East 7 Mile Rd (Figure 2), upwind of the latter during the 30 June 2022 episode. Between the two stations is the most heavily industrialized area in Michigan (Figure 3), including two natural gas-fired power plants (DG, DT), a refinery (MP), three automotive assembly facilities (FD, MA, JN), two steel mills (US, CC), a coking facility (ES), a lime plant (CL), a large wastewater treatment plant (WW), and a hazardous waste storage and treatment facility (UE). MDA8 O₃ at Allen Park station during the O₃ episode was 71 ppb, 5 ppb less than at East 7 Mile Rd., a concentration possibly influenced by long-range transport. The question addressed by this study is whether the O₃ gradient observed between the two stations on 30 June 2022 was due to the significant concentration of local emission sources downwind of Allen Park, and whether emission reductions at these sources would be beneficial as potential control strategies.

2. Methods

The MicroFACT model [2] used for this study is a sub-km horizontal resolution 3D Eulerian grid model with forward and inverse modes, the latter based on the 4D variational data assimilation technique [4,5]. It is driven by 3D building-sensitive winds derived from Los Alamos National Lab’s Quick Urban Industrial Complex (QUIC) model [6]. MicroFACT simulates pollutant transport by advection and turbulent diffusion, while chemistry is simulated with 116 gas-phase and 7 heterogeneous reactions. For this study, the recent parameterization of nitrous acid (HONO) formation based on smog chamber experiments by Song et al. [7] was added to the chemical mechanism to enhance radical production and HONO concentrations in line with recent field measurements [8]. The clear sky photolysis rates calculated using the method of Saunders et al. [9] were also enhanced for the effects of surface albedo on NO₂ and O₃ photolysis, assuming for simplicity that direct radiation dominates the former and diffuse radiation the latter [10].

MicroFACT can not only predict atmospheric concentrations of pollutants based on specified emissions but can also infer emissions at fine scale from ambient air concentrations. Olaguer previously used MicroFACT to perform inverse modeling of formaldehyde (HCHO) emissions from industrial and mobile sources [11] and methane emissions from natural gas pipeline leaks [12] based on mobile lab measurements during the 2021 Michigan-Ontario Ozone Source Experiment (MOOSE) [13] within the area shown in Figure 2.

For this study, MicroFACT was applied to a 24 km × 24 km, 400 m resolution horizontal grid as shown in Figure 3. The vertical grid had 20 layers from 0 to 1500 m above ground level, with parabolically increasing layer thickness with height (2 m at surface). The horizontal resolution of MicroFACT as implemented in this study is 10 times the typical highest resolution of regional air quality models used for regulatory applications in the U.S. Figure 3 shows contours of modeled O₃ concentration at 14:00 LST, 30 June 2022, the time of maximum hourly average O₃ simulated by MicroFACT at East 7 Mile Rd (see Table A1 in Appendix A). Note the significant horizontal gradients in O₃ apparent at fine scale that would be unresolved by coarser resolution models.

The QUIC model was used to generate hourly 3D wind fields consistent with flow around buildings at 20 m horizontal resolution based on background wind and other meteorological measurements at the Detroit Metropolitan Wayne County Airport (DTW) in Romulus, Michigan, roughly 10 km to the southwest of Allen Park. Airport measurements were run through the EPA AERMET pre-processor [14] to produce hourly values for wind speed and direction as well as turbulence parameters required by QUIC. The resulting 3D winds were then averaged and staggered with respect to the MicroFACT grid to ensure mass conservation (i.e., an Arakawa C-grid [15]). Figure 4 shows wind streamlines at 1 m above ground level generated by the QUIC model for 14:00 LST, 30 June 2022.

Hourly vertical temperature profiles in the domain interior were based on measurements at the Windsor International Airport (YQG) in Ontario, Canada located within the model grid. These measurements were processed with AERMET to produce boundary layer parameters required to infer atmospheric thermal structure. In the convectively unstable or neutral cases the temperature profile was based on the observed lapse rate. In the stable case, it was computed based on the Monin-Obukhov length, the surface heat flux, the roughness length, the friction velocity, and the inversion height [16].

U.S. annual point source emissions for 2022 and associated stack parameters were specified from Michigan’s State and Local Emissions Inventory System [17] with speciation based on the EPA Emissions Modeling Platform [18]. Point source emissions for Windsor, Ontario were supplied by Environment and Climate Change Canada (ECCC) for the available year 2017. Hourly mobile source emissions were based on 2022 runs of the USEPA Motor Vehicle Emission Simulator (MOVES) [19] by the Lake Michigan Air Directors Consortium (LADCO), while hourly biogenic emissions were based on LADCO 2022 runs of the Biogenic Emission Inventory System (BEIS) [20]. For 2022 nonpoint sources, annual combustion emissions were obtained from Wayne County, Michigan, while speciated Volatile Organic Compound (VOC) emissions from consumer and commercial products were obtained from the USEPA. Road lengths and population were used as surrogates to distribute emissions among the 400 m × 400 m grid cells. Note that point source emissions of HCHO were adjusted based on known HCHO:CO emission ratios (2-5%) according to the method described by Olaguer et al. [21], to be consistent with the findings of the MOOSE study regarding industrial primary HCHO [11,21].

Hourly species boundary conditions for 35 advected species and surrogates for lumped species (e.g., NO3X = NO₃+N₂O₅) were derived by combining measurements at Allen Park for surface O₃, CO, and total reactive nitrogen (NO_y) and simulated 1D vertical profiles at the same location from regional air quality models, specifically ECCC’s GEM-MACH [22] (2.5-km resolution for O₃, CO, HCHO, ethene, and most NO_y species) and a LADCO run of CAMx [23] (4-km resolution for most organic species and NO3X). While GEM-MACH had higher horizontal resolution than CAMx, the latter model’s lumping scheme for organic species was closer to that of MicroFACT, hence the necessity of separate sources for the boundary conditions.

Vertical ozone profiles derived from ozone sondes, ozone lidars, ground-based multi-axis differential optical absorption spectrometers, and satellites often show elevated layers of enhanced tropospheric ozone exceeding surface concentrations at altitudes below 1 km during pollution events in the planetary boundary layer [24,25,26,27]. Hourly O₃ profiles from the GEM-MACH model showed enhanced layers aloft in the late morning and afternoon, but with lower concentrations at the surface than measured at Allen Park (see Appendix B). The difference between the measured surface concentration and the GEM-MACH value was uniformly added throughout the vertical column to preserve the modeled layer structure (Table B1). Similar adjustments were made for CO and NO_y, with the latter distributed proportionally among relevant species (Table B2 and Table B3). Pollutants aloft were entrained to the surface as the planetary boundary layer expanded based on the vertical diffusivity parameterization of Delle Monache et al. [28] and boundary layer parameters derived from YQG airport data. Figure 5 shows the ozone boundary condition and the accompanying vertical diffusivity profile at 14:00 LST, 30 June 2022.

The simulation period was 5:00 to 19:00 LST with 20 s time steps. Initial conditions for species concentrations were set equal to the boundary conditions at 5:00 LST. The spin-up period consisted of the first four hours, during which the initial conditions were advected out of the domain (at a typical wind speed of ~3 m/s) and replaced by more realistic species concentrations consistent with local emissions and chemistry. Because of the low nocturnal mixing height, urban morphology was ignored in the wind field during the first hour to avoid computational instability due to taller structures in downtown Detroit. The ten-hour (10-h) period beyond spin-up from 9:00 to 19:00 LST coincided with measured O₃ concentrations at East 7 Mile Rd above 60 ppb (Table A2).

3. Results

The model inputs described above define the base case configuration. The simulation results for the base case were compared with measurements at the East 7 Mile Rd station for O₃, NO, NO₂, odd oxygen (O₃ + NO₂), and NO_y, as illustrated in Table 1 and Table 2 for the hours when observed O₃ exceeded 60 ppb. In these tables, the Normalized Mean Bias (NMB) was computed as the average of the difference between the modeled and measured hourly concentrations over the appropriate 8-h or 10-h period, divided by the corresponding measured 8-h or 10-h concentration. Comparisons between modeled concentrations and measurements for the individual hours are displayed in Table A1, Table A2 and Table A3.

MicroFACT successfully simulated 8-h O₃ for the periods when observed O₃ exceeded 60 ppb, although the simulated MDA8 O₃ occurred at 9:00 – 17:00 LST, whereas the observed MDA8 O₃ occurred at 11:00 – 19:00 LST. Reactive nitrogen species were over-predicted relative to observations, while the simulated odd oxygen (O_x) agreed well with observations (within 10%) despite the errors in NO_y species. The hourly comparison presented in Appendix A shows that the model underpredicted O₃ relative to observations in the hours just before sunset but compensated for this by overpredicting NO₂. During the real episode, the buildup of reactive nitrogen was prevented in the later hours of the day while O₃ remained at or above 75 ppb despite the weakening of photolytic conversion of NO₂ to O₃. The observed behavior of O₃ towards evening was perhaps due to greater persistence of O₃ entrainment from aloft in reality compared to the simulation.

The overprediction of NO_y was possibly due to errors in nitrogen oxide (NO_x = NO + NO₂) emissions or chemical loss mechanisms. In contrast, 2016 CMAQ model simulations of summertime NO_x by the USEPA for the continental U.S. have late morning to late afternoon NMB values roughly between -15% and -30% [29]. The discovery of much higher observed ratios of HONO:HNO₃ than simulated by current models during the MOOSE study [13] may indicate fundamental deficiencies in our current understanding of NO_y emissions and chemistry. Despite the inclusion of heterogeneous reactions on both aerosol and ground surfaces, a specified ratio of aerosol surface area to air volume of 0.0014 m² m⁻³ at the high end of measured values in the U.S. [30], and the addition of the HONO parameterization of Song et al. [7] in MicroFACT, the model did not generate HONO:HNO₃ ratios above 0.4, compared to mean ratios of 3 measured in Detroit during MOOSE.

Nitrous acid formation is closely tied to poorly understood heterogeneous chemical processes that may play an essential role in regulating the NO_y budget in urban airsheds. Akimoto and Tanimoto [31] noted frequent reports of disagreement between observations and corresponding air quality model simulations of HONO and HNO₃. They called for comprehensive field study measurements of NO_y species to improve the simulation of chemical processes in atmospheric models, with special emphasis on heterogeneous “renoxification” reactions involving HONO and HNO₃.

In addition to the base case simulation, a sensitivity experiment was conducted in which enhanced O₃ layers aloft were suppressed when surface O₃ at Allen Park exceeded 60 ppb. The measured hourly concentration at Allen Park in this case was set as the maximum value aloft in the O₃ boundary conditions. The result was that the modeled MDA8 O₃ at East 7 Mile Rd decreased by 15 ppb to 63 ppb, well below the value at Allen Park. This was because NO_x sources between Allen Park and East 7 Mile Rd titrated O₃ (which is recovered when NO₂ produced by NO titration of O₃ photolyzes further downwind), an effect diluted by the coarser resolution of regional air quality models. In the base case, O₃ entrained from aloft countered this titration by adding directly to surface O₃ and by increasing the production of radicals through the photolysis of subsiding ozone followed by H₂O + O(¹D) → 2OH. This is illustrated in Table 3 below, which shows differences in Ozone Production Efficiency (OPE) between the sensitivity experiment and the base case, as well as corresponding differences in O₃ and O_x averaged over 9:00 – 19:00 LST.

The OPE, defined by Chace et al. [32] as the slope of the linear regression between odd oxygen and reactive nitrogen reservoirs (NO_z = NO_y – NO_x), roughly indicates the amount of O₃ eventually generated per NO_x molecule through conversion of NO to NO₂ by atmospheric radicals before NO_x is terminated and stored in NO_z. Note how the OPE decreased by almost half between the base case and the sensitivity experiment. For comparison, the corresponding OPE computed from station measurements was 15.4, which may indicate that downward mixing of O₃ from aloft was even stronger than simulated by MicroFACT, especially late in the day.

In addition to the sensitivity experiment described above, additional experiments were conducted involving hypothetical emission reductions of NO_x (10%) and VOCs (10% and 50%) at stationary sources. These experiments are summarized in Table 4 below.

Clearly, NO_x titration results in an increase in O₃ due to NO_x reduction, so that local O₃ is radical-limited rather than NO_x-limited [33]. The weak impact of VOC reductions indicates that VOCs are not the major source of radicals compared to O₃ entrained from aloft. This does not rule out, however, that local VOC and NO_x controls may result in O₃ reductions further downwind of East 7 Mile Rd under southwesterly flow. It is also possible that emission reductions upwind of Allen Park but still within the SEMI region may decrease MDA8 O₃ at Allen Park and East 7 Mile Rd by mitigating O₃ transport aloft.

4. Discussion and Conclusions

Based on station measurements at Allen Park and East 7 Mile Rd, and corresponding MicroFACT and trajectory model simulation results for the 30 June 2022 O₃ episode in Detroit, Michigan, the following conclusions may be drawn:

• Local emission sources between Allen Park and East 7 Mile Rd tend to reduce MDA8 O₃ at the latter station below the corresponding concentration at Allen Park due to NO_x titration of O₃.
• Pollution layers aloft may counter O₃ titration when they are turbulently entrained to the surface, thus adding directly to ground-level O₃ and enhancing radical concentrations and O₃ production efficiency.
• Transport of O₃ around 500 m above ground level may significantly contribute to MDA8 O₃ above 70 ppb at East 7 Mile Rd during southwesterly wind flow that imports pollution from a wide region in the U.S. Midwest.

MicroFACT was able to demonstrate the above mechanism for the O₃ episode because of the model’s fine-scale resolution relative to coarser-resolution regional models, which dilute the NO_x titration effect and underestimate O₃ concentrations aloft by artificially diffusing plumes. The fine-scale simulation results show the importance of routine measurements of O₃ vertical profiles and boundary layer mixing height to enable evaluation of elevated pollution layers and their entrainment to the surface. Unfortunately, such measurements were unavailable for this study. (Note that a ceilometer is normally operated at the East 7 Mile Rd station. However, the instrument was offline on 30 June 2022.)

More explicit coupling between regional and fine-scale air quality models may improve our understanding of the physical mechanisms involved in high O₃ episodes in SEMI. Among the uncertainties that require further examination in future simulations of O₃ episodes are the details of NO_y chemistry, especially heterogeneous production and loss mechanisms that may limit NO_y concentrations before sunset.

While local VOC and NO_x controls in the modeled area may not significantly mitigate MDA8 O₃ at East 7 Mile Rd, SEMI at large may still be VOC-sensitive [34], so that a combination of VOC controls over a wider area and continued U.S. national reductions in NO_x emissions may help ensure O₃ attainment in SEMI by reducing either near-regional or longer-range O₃ transport in the lower planetary boundary layer.