Impacts of physical parameterization on prediction of ethane concentrations for oil and gas emissions in WRF-Chem

Recent increases in natural gas (NG) production through hydraulic fracturing have called the climate benefit of switching from coal-fired to natural gas-fired power plants into question. Higher than expected levels of methane, non-methane hydrocarbons (NMHC), and NO_x have been observed in areas close to oil and NG operation facilities. Large uncertainties in the oil and NG operation emission inventories reduce the confidence level in the impact assessment of such activities on regional air quality and climate, as well as in the development of effective mitigation policies. In this work, we used ethane as the indicator of oil and NG emissions and explored the sensitivity of ethane to different physical parameterizations and simulation setups in the Weather Research and Forecasting with Chemistry (WRF-Chem) model using the US EPA National Emission Inventory (NEI-2011). We evaluated the impact of the following configurations and parameterizations on predicted ethane concentrations: planetary boundary layer (PBL) parameterizations, daily re-initialization of meteorological variables, meteorological initial and boundary conditions, and horizontal resolution. We assessed the uncertainties around oil and NG emissions using measurements from the FRAPPÉ and DISCOVER-AQ campaigns over the northern Front Range metropolitan area (NFRMA) in summer 2014. The sensitivity analysis shows up to 57.3 % variability in the normalized mean bias of the near-surface modeled ethane across the simulations, which highlights the important role of model configurations on the model performance and ultimately the assessment of emissions. Comparison between airborne measurements and the sensitivity simulations indicates that the model–measurement bias of ethane ranged from −14.9 to −8.2 ppb (NMB ranged from −80.5 % to −44 %) in regions close to oil and NG activities. Underprediction of ethane concentration in all sensitivity runs suggests an actual underestimation of the oil and NG emissions in the NEI-2011. An increase of oil and NG emissions in the simulations partially improved the model performance in capturing ethane and lumped alkanes (HC3) concentrations but did not impact the model performance in capturing benzene, toluene, and xylene; this is due to very low emission rates of the latter species from the oil and NG sector in NEI-2011.

1 Introduction

Recent advances in the unconventional natural gas (NG) production technology (hydraulic fracturing) have resulted in economical access to NG reserves in deep shale formations and a 36 % rise in NG production in the US from 2005 to 2014 (Lyon, 2015). The increase in the NG production, the decrease in the NG price, and the environmental advantages of NG-fired power plants over coal-fired power plants have made NG an important competitor for coal in the electricity generation sector. In 2015, NG and coal each had a 33 % share in the electricity generation in the US. It is predicted that NG's share in electricity generation will grow 1.5 % every year (Energy Information Administration of US Department of Energy, 2016; US Energy Information Administration, 2016). With the rapid increase in the unconventional oil and NG production, higher than expected levels of greenhouse gases, specifically methane, and air pollutants such as non-methane hydrocarbons (NMHC) and NO_x (from flaring or compressors, reboilers, pneumatic devices, trucks, and other equipment using fossil fuel) (Allen, 2016; Olaguer, 2012) have been observed in some places in the vicinity of oil and NG facilities. The high concentrations of these chemicals measured in many studies at different scales and in different regions suggest that official emission inventories (e.g., the Greenhouse Gas Inventory – GHGI and the Emission Database for Global Atmospheric Research – EDGAR) fail to capture the magnitude of emissions from unconventional extraction activities (Brandt et al., 2014). This underestimation by emission inventories has raised concerns regarding the climate implications of promoting NG as the “bridge fuel” (Alvarez et al., 2012; Howarth et al., 2011; Levi, 2013; McJeon et al., 2014), and its impacts on the air quality and public health (Halliday et al., 2016; McKenzie et al., 2012). Additionally, methane and NMHCs emitted from the oil and NG sector can degrade regional air quality and contribute to ozone formation on regional and global scales (Helmig et al., 2016). Outdated emission factors (EF), super-emitters in the production systems, and rapid growth in the production facilities are some of the reasons for the underestimation (Brandt et al., 2014; Lyon, 2015; Zavala-Araiza et al., 2015).

The Colorado northern Front Range (NFR), including the Denver metropolitan area, is located between the Rocky Mountains and the High Plains with a total population of about 4.8 million. In 2007, a large region of the NFR was declared in non-attainment of the National Ambient Air Quality Standard (NAAQS) for 8 h average ozone. Major sources of pollutants in this area are vehicle emissions, oil and NG operation, agriculture and feedlots, and power plants. In the past few years, oil and NG development has increased drastically in the NFR. NG production in Weld County has increased from 55.8×10⁶ m³ (1.97×10⁶ MCF – thousand cubic feet) to 181.8×10⁶ m³ (6.42×10⁶ MCF) from 2004 to 2016. The Wattenberg gas field in Weld County is close to populated regions and has the highest well density in the NFR with more than 25 000 active NG wells (Colorado Oil and Gas Conservation Commission, 2017). In the NFR, measured NMHCs are 18–77 times greater than the regional background as determined from the NOAA flask network (Thompson et al., 2014). High levels of NMHCs can cause health concerns at regional scales and can contribute significantly to the ozone pollution in the region (Cheadle et al., 2017; Gilman et al., 2013; McDuffie et al., 2016; Pétron et al., 2012; Pfister et al., 2017b; Thompson et al., 2014). Using box models constrained with observations, McDuffie et al. (2016) estimated that NFR oil and NG activities contribute ∼50 % to the regional volatile organic compound (VOC) OH reactivity and 20 % to the regional photochemical ozone production.

Mass balance approach methods have been widely used to estimate the emissions from oil and NG activities (Conley et al., 2016; Karion et al., 2015; Peischl et al., 2016; Pétron et al., 2012; Smith et al., 2015). This method cannot provide details on the spatial and temporal variability of emissions and has limitations in domains with complex atmospheric transport such as the NFR. High resolution three-dimensional atmospheric chemical transport models can better capture the variability in meteorology and chemistry in different domains. Paired with observations and using inverse modeling techniques, these models help evaluate the performance of emission inventories on high temporal and spatial scales (Barkley et al., 2017; Cui et al., 2014, 2017) and allow assessments of the impact of oil and NG activities on regional air quality. Ahmadov et al. (2015) used the Weather Research and Forecasting Model with Chemistry (WRF-Chem) to study high ozone episodes and emission reduction scenarios in the Uintah Basin. Their results show a strong underestimation of methane and VOC emissions in the National Emission Inventory 2011 (NEI-2011).

WRF-Chem provides users with different dynamical, physical, and chemical schemes (Grell et al., 2005; Skamarock et al., 2008). These choices can impact the performance of the model, specifically in regions with complex transport patterns (Saide et al., 2011). In order to assess the performance of emission inventories, it is critical to address the uncertainties derived from model configurations on simulated concentration fields. The goal of this study is to quantify the impact of WRF-Chem configurations on predicting the oil and NG emissions in the NFR. VOCs in the NFR have shown a clear source signature associated with oil and NG activities (Gilman et al., 2013; Pétron et al., 2014). Diverse air pollution sources and complex metrological patterns due to mountain–valley circulation, high elevation, and harsh terrain are some of the challenges for air quality modeling in this area. We use ethane, which has a simple chemical cycle and a lifetime of about 2 months, as a tracer for oil and NG (Helmig et al., 2016). The model and emission inventory performance are evaluated by comparing meteorological parameters as well as ethane and VOC concentrations to surface and airborne measurements. We explore the sensitivity of the modeled transport and ethane concentrations to different WRF-Chem physical parameterizations and setups. This work will be followed by the development of an inverse modeling technique to constrain the oil and NG emission rates by calculating optimal scaling factor for the emission inventory. Simulations discussed in this study will be used to calculate the variability of the optimal scaling factor. Moreover, to inform not only about the absolute magnitude in the ethane emissions but to further explore the feasibility to constrain other trace gas oil and NG emissions, we investigate CO and VOC emission estimates from the oil and NG sector and VOC ratios in the observations and in the model.

Figure 1 Terrain map of the WRF-Chem outer domain (d01) and inner domains (d02) and the location of observation sites. (a) The two nested domains designed for this study. (b) The zoomed in map of domain 2 with the location of several sites. Grey dots show the location of permitted wells (http://cogcc.state.co.us/, last acccess: 12 June 2018).

2 Method

2.1 Aircraft and ground-based observations

The National Science Foundation/National Center for Atmospheric Research (NSF/NCAR) Front Range Air Pollution and Photochemistry Éxperiment (FRAPPÉ) and the National Aeronautics and Space Administration (NASA) Deriving Information on Surface Conditions from Column and Vertically Resolved Observations Relevant to Air Quality (DISCOVER-AQ) campaigns were conducted in July and August 2014, in the NFR, Colorado. These two campaigns provide detailed and coherent airborne and ground-based measurements in this area, which can assist in the evaluation and improvement of chemical transport models and emission inventories. The NSF/NCAR C-130 collected extensive airborne measurements of various atmospheric constituents during the FRAPPÉ campaign. A total of 15 flights (∼80 flight hours) were conducted in the NFR with the goal of mapping the emissions and their transport and chemistry in this region. During the DISCOVER-AQ campaign, the NASA P-3B aircraft performed approximately 20 flights containing spiral ascents or descents over six key sites in the NFR to capture the vertical profiles of the atmospheric constituents and their diurnal variation. Ethane was measured on board of C-130 and P-3 aircraft. On C-130 aircraft, ethane was measured by the University of Colorado's CAMS instruments with a detection sensitivity of 15 ppt, the details of which are discussed in Richter et al. (2015). An Aerodyne Ethane-Mini spectrometer on P-3 was used to measure the ethane concentration (Yacovitch et al., 2014). Fried (2015) compared CAMS ethane measurements with sub-ppb precision with the Aerodyne measurements during wing tip comparisons and the agreement was within 9 %, corresponding to differences of less than 55 ppt.

The National Oceanic and Atmospheric Administration (NOAA), the Colorado Department of Public Health and Environment (CDPHE), and the National Park Services (NPS) operated numerous ground-level measurement sites during these two campaigns. In this work, we present ground-level measurements from the NOAA Boulder Atmospheric Observatory (BAO; 40.05^∘ N, 105.01^∘ W, 1584 m a.s.l. – above sea level), the NOAA Platteville site (PAO; 40.18^∘ N, −104.73^∘ W, 1523 m a.s.l.), Fort Collins (FC; 40.89^∘ N, −105.13^∘ W, 1572 m a.s.l.), NREL-Golden (Golden; 39.74^∘ N, −105.18^∘ W, 1833 m a.s.l.), and CDPHE wind measurements at Weld County tower (WC-Tower; 40.39^∘ N, −104.73^∘ W, 1483 m a.s.l.), Rocky Flats N (RF-N; 39.91^∘ N, −105.19^∘ W, 1803 m a.s.l.), Welch (39.64^∘ N, −105.14^∘ W, 1743 m a.s.l.), Chatfield (39.53^∘ N, −105.07^∘ W, 16 756 m a.s.l.), and Aurora-East (39.64^∘ N, −104.57^∘ W, 1802 m a.s.l.). BAO and PAO are located north of Denver and close to the Wattenberg gas field in Weld County (Fig. 1). Measurements of temperature, relative humidity, and wind speed and direction at 10, 100, and 300 m were recorded at BAO. Surface wind measurements from PAO (3 m) and WC-Tower (4 m) were used in this study. The planetary boundary layer (PBL) height was measured and calculated at PAO, FC, and Golden using micro-pulse lidar backscatter during the daytime (Compton et al., 2013).

2.2 WRF-Chem model

We used WRF-Chem 3.6.1 (Grell et al., 2005; Skamarock et al., 2008), a fully coupled online air quality and transport model, to investigate the sensitivity of modeled PBL, winds, temperature, relative humidity, and ethane concentrations to different physical parameterizations and configurations. Figure 1 illustrates the location of the two nested domains and the underlying terrain map. We used one-way nesting (i.e., the outer domain ran independently of the inner domain). The outer domain has a 12 km × 12 km horizontal resolution, and the inner domain has a 4 km × 4 km horizontal resolution. Both domains have 53 vertical levels with the domain top at 50 hPa (∼11 layers below 1 km). The outer domain is designed to capture the emission from the western US, and the inner domain includes Colorado and Utah. Sensitivity simulations start on 24 July and end on 18 August 2014. Table 1 shows a summary of the WRF-Chem configurations for this study, used in all sensitivity simulations. The Morrison double-moment scheme was selected as the microphysics option and Goddard shortwave (Chou and Suarez, 1999) and RRTMG longwave radiation schemes (Iacono et al., 2008) were used as shortwave and longwave radiation parameterizations, respectively. The Grell-Freitas convection scheme (Grell and Freitas, 2014) was used as convective parameterization for both the outer and inner domain. The inner domain falls into the “grey-scale” which means many of the assumptions used in convective parameterization will no longer be valid at this resolution. The Grell–Freitas convection scheme is a stochastic scale dependent convective parameterization based on the method proposed by Arakawa et al. (2011) and is designed for domains with a horizontal resolution of up to few kilometers. Comparisons between a simulation with resolved convection of the inner domain and a simulation using the Grell–Freitas convective parameterization in the inner domain showed similar performance in capturing transport (not shown). Thus, we used the Grell–Freitas convective scheme for both domains in all simulations to reduce the computation costs.

We selected the Regional Atmospheric Chemistry Mechanism chemistry using Earth System Research Laboratory (RACM-ESRL) (Stockwell et al., 1997) coupled to the Modal Aerosol Dynamics Model/Secondary Organic Aerosol Model (MADE/SORGAM). RACM includes 17 stable inorganics, 4 inorganic intermediates, 32 stable organic species, and 24 organic intermediates. RACM_ESRL (Kim et al., 2009) is an updated version of the RACM mechanism and includes 23 photolysis and 221 chemical reactions (Ahmadov et al., 2015). To reduce the computational costs, hydrocarbons with similar behavior are lumped together in the chemical mechanisms. For example, “HC3” in the RACM_ESRL mechanism includes alkanes such as propane, n-butane, isobutane, and acetylene (ethyne), and alcohols such as methanol and ethanol. “TOL” includes toluene and benzene. Ethane and methane are treated exclusively in the RACM_ESRL mechanism. More details regarding the reactions and lumping groups can be found in Stockwell et al. (1997). Chemical boundary conditions from Monitoring Atmospheric Composition and Climate reanalysis (MACC), available every 3 h, (Inness et al., 2013) and model outputs from RAQMS, available every 6 h, (Natarajan et al., 2012; Pierce et al., 2007) were used as chemical boundary and initial conditions in the simulations. The model outputs from these global models are specific to the simulation time (24 July to 18 August 2014) and are interpolated to the WRF-Chem domain and temporal resolution prior to starting the simulations. Ethane concentrations showed no strong sensitivity to the two different chemical initial and boundary conditions (i.e., RAQMS and MACC) and are not discussed further.

Table 1 Summary of basic WRF-Chem configuration.

WRF-Chem sensitivity tests

WRF-Chem provides users with a number of different dynamical, physical, and chemical schemes. Users can select schemes based on the physical properties of the domain of interest, goals of the study, and computational limitations. We evaluated the sensitivity of WRF-Chem to different physics options, such as the PBL parameterization, and configurations including the daily re-initialization of meteorological fields, different meteorological initial and boundary conditions, and varying horizontal resolution. Table 2 shows details on the sensitivity runs and lists the meteorological and chemical boundary conditions used for each run. The naming system for the simulations is based on the different settings (e.g., simulation 5-MnERi represents the simulation number “5”, PBL scheme “MYNN3”, meteorological initial and boundary condition “ERA-Interim”, chemical initial and boundary condition “RAQMS”, and daily re-initialization of meteorological fields “i”). The simulation IDs in Table 2 are used when discussing sensitivity tests in the paper.

Table 2 Summary of WRF-Chem configurations for sensitivity tests designed for this study. Sensitivity tests are divided by horizontal lines.

An accurate simulation of air pollution is dependent on a precise description of transport processes, meteorological conditions, and the PBL height (PBLH) (Cuchiara et al., 2014; Hu et al., 2010; Sobhani et al., 2018). Transport of pollutants within the domain depends on turbulent motions and vertical mixing within the PBL. WRF-Chem (3.6.1) has 11 different PBL schemes to address the closure problem in the simulation of turbulent motions. In general, PBL schemes can be classified into two main groups; local and non-local. A local PBL scheme estimates the turbulent fluxes of heat, momentum, and moisture from local mean and gradient flux values. In a non-local PBL scheme, non-local fluxes can influence fluxes in each grid; hence, these schemes are expected to better capture large-size eddies in the simulation (Stull, 1988). We tested one non-local and two local PBL schemes to understand the sensitivity of the model to PBL parameterization in a domain with high elevation and complex terrain. We used the Yonsei University (YSU) first-order scheme (Hong et al., 2006) as the non-local PBL scheme in the PBL1 (1-YFM) simulation. The local schemes used in the PBL2 (2-MjFM) and PBL3 (3-MnFM) simulations were Mellor–Yamada–Janjic (MYJ) 1.5 order (2.5 level) (Janjic, 2001; Janjic et al., 2000) and Mellor–Yamada–Nakanishi–Niino (MYNN3) 3rd level (Nakanishi and Niino, 2009).

WRF-Chem is a mesoscale model and requires initial and lateral boundary conditions from a larger-scale model. Usually, these initial and lateral boundary conditions are taken from the reanalysis products of larger-scale models optimized using assimilation techniques and observations. The choice of initial and lateral boundary condition products can impact the model performance (Angevine et al., 2012; Saide et al., 2011). We tested two different meteorological initial and boundary conditions, European Reanalysis (ERA-Interim) from the European Center for Medium-Range Weather Forecasts (ECMWF) in the Met5 (5-MnERi) simulation and NCEP's Global Forecast System (GFS) in the Met6 (6-MnFRi) simulation. ERA-Interim reanalysis is produced with 80 km by 80 km horizontal and 6 h temporal resolution (ECMWF, 2009), and NCEP FNL (final) operational global analysis is produced using GFS with 1 ^∘ × 1 ^∘ horizontal and 6 h temporal resolution (National Centers for Environmental Prediction, National Weather Service, and NOAA, 2000).

Simulations were performed for 24 days from 24 July to 18 August 2014. Initializing the meteorological fields in the simulation at the first time step with the larger-scale model values and running it for 24 days without any nudging will result in deviations from the larger-scale reanalysis products. Conversely, the lower resolution of the larger-scale models can lower the accuracy of WRF-Chem high-resolution simulations. To investigate this impact, we tested two different setups for WRF-Chem. In the Init4 (4-MnER) simulation, we initialized the meteorological fields at the first time step with larger-scale model values and ran the simulation freely for 24 days (“free run”). In the Init5 (5-MnERi) simulation, the meteorological fields were re-initialized every day at 18:00 UTC (12:00 LT – local time) and run for the next 30 h. The first 6 h of the simulation (18:00 to 00:00 UTC) were discarded to allow for the model to spin-up. In this setup, chemistry fields were recycled from previous cycles of simulations.

The sensitivity of the model to the horizontal resolution was examined by comparing the performance of the outer domain (12 km × 12 km) to the inner domain (4 km × 4 km) in the Hor5 (5-MnERi) simulation. In one-way nesting, the outer domain runs independently of the inner domain; thus, comparing the performance of the outer and inner domains is valid.

Table 3 Summary of model performance in capturing temperature (T) and relative humidity (RH) at BAO 100 m from 1 to 15 August 2014.

2.3 Emission inventory

NEI-2011 version 2 is a bottom-up emission inventory of US anthropogenic emissions. While we cannot expect the year 2011 inventory to fully represent the model year 2014, it was the only inventory available to the WRF-Chem user community at the time of this study. Emissions in this inventory are calculated based on fuel consumption, source activity, and emission factors reported by state, tribal, and local governing agencies (US Environmental Protection Agency, 2015). A processed version of NEI-2011 is available to the users, which includes emission of 76 species (50 speciated VOC compounds, 19 PM_2.5 aerosol species, and 7 primary species). NEI-2011 and emissions for only the oil and NG sector in the NEI-2011 were provided to us by Dr. Stuart McKeen (NOAA Earth Systems Laboratory, Boulder, CO). Table S1 in the Supplement includes details on the mapping table used to convert NEI-2011 species to the RACM and MADE/SORGAM chemical and aerosol mechanism. The separate oil and NG emission information was used to conduct an additional sensitivity simulation with perturbed oil and NG emission, which we used to study the sensitivity of modeled ethane concentrations as well as concentrations of VOCs and CO to the oil and NG emission sector. We used the Model of Emissions of Gases and Aerosols from Nature (MEGAN) for biogenic emission in all simulations (Guenther et al., 2012). Ethane does not have a significant biogenic source (Yacovitch et al., 2014); thus, we did not assess the impact of biogenic emissions in this study. Wildfire emissions were not included in the simulations, but this will have a negligible impact on the results as wildfires did not significantly influence the air quality in the NFR during the FRAPPÉ campaign (Valerino et al., 2017).

3 Results and discussion

We start with an evaluation of the overall performance of all simulations and later provide a detailed discussion on the different sets of sensitivity simulations. To evaluate the sensitivity of WRF-Chem to different physical parameterizations, we compared the simulated meteorological variables, such as temperature, relative humidity, wind fields, and PBLH, with measurements. The 27 and 28 July 2014 were reported as Denver cyclone episodes (Dingle et al., 2016; Valerino et al., 2017; Vu et al., 2016), and neither simulation captured the cyclone pattern and enhancements accurately on these 2 days. Thus, we only included the period from 1 to 15 August 2014 in our analysis to avoid skewing the results due to large model errors during the cyclone episodes. For quantitative comparison between the simulations, we used statistical measures including the correlation coefficient (R), the root mean square error (RMSE), the mean absolute error (MAE), the mean bias (MB), and the normalized mean bias (NMB). Definitions of these metrics can be found in the Supplement. We used the NMB as a proxy for model sensitivity to quantify the impact of model configuration on different variables. Variability of the NMB (calculated by subtracting the minimum NMB from the maximum NMB) in sensitivity tests can provide a range for uncertainties in the model cases independent of the model values.

Figure 2 Average diurnal cycle of temperature (a) , relative humidity (b), wind speed (c) and wind direction (d) for all tests and observations at BAO 100 m. Averages are calculated from 1 to 15 August 2014.

3.1 Evaluation of overall model performance

Table 3 includes the statistical measures for temperature and relative humidity in all the simulation tests at 100 m altitude at BAO. Figure 2 compares the diurnal cycles of measured temperature, relative humidity, wind speed, and wind direction at 100 m altitude at BAO with corresponding model values for all the simulation tests. While Fig. 2 provides an overview of all sensitivity tests, Fig. S2 in the Supplement separates each sensitivity test to provide a clearer test by test comparison. Similarly, Tables S2 to S5 include statistical measures and Figs. S1 and S3 show diurnal cycles of temperature, relative humidity, wind speed, and wind direction at BAO at 10 and 300 m. All model simulations capture the overall daily cycle in temperature and relative humidity well (Fig. 2 and Table 3). The variability across different sensitivity runs can be large, with modeled temperature varying by up to 6 ^∘C and the model–measurement NMB ranging from −3.9 % to 11.1 %. Relative humidity has larger variability among the simulations during nighttime compared to daytime. The NMB of relative humidity ranges from −29.7 % to 52.6 %.

Wind patterns vary significantly from daytime to nighttime. During the day, wind primarily blows from the east towards the Rocky Mountains with a slight southerly component. During the night, this pattern switches to predominantly westerly winds bringing cooler air to lower terrain. Wind measurements at the BAO at different altitudes (10, 100, 300 m) can help us better understand the wind pattern at higher model levels. Table 4 includes the mean and standard deviation of daytime and nighttime wind fields in the simulations and the observations at 100 m. Results for the 10 and 300 m levels at BAO from 1 August to 15 August 2014 are included in Tables S4 and S5, respectively. In addition to BAO, we investigated the wind sensitivity to physical parameterizations at two other sites that are close to oil and NG operations, WC-Tower and PAO (Fig. S5). At BAO, higher wind speeds were measured at higher elevations which is captured by the model. Overall, most simulations show skill in capturing diurnal cycles of wind speed and direction with better agreement with observations for daytime (Tables 4, S4, and S5). Overall, the model runs show fairly good performance in capturing temperature, relative humidity, and wind fields, especially for daytime. A higher sensitivity to the physical parameterization was observed for nighttime.

Ethane is predominantly emitted from oil and NG production sites (Helmig et al., 2016; Xiao et al., 2008) and is a valuable chemical tracer to study the transport patterns of oil and NG emissions. To evaluate the impact of vertical mixing intensity on the distribution of pollutants, we compared the vertical distribution of ethane concentrations between the simulations. Figure 3 shows the diurnal cycle of the averaged vertical cross section of ethane concentrations at PAO with the measured PBL height for each simulation.

Table 4 Summary of model performance in capturing wind speed and direction at BAO 100 m during the 1–15 August 2014 period.

Figure 3 Cross section of modeled ethane at PAO and the measured PBL height (black dots) averaged from 1 to 10 August 2014.

Figure 4 Vertical distribution of simulated and measured ethane in the NFR area separated by flight time. (a) C-130 AM 09:00 to 12:00 observation and the corresponding model values. (b) C-130 PM 12:00 to 18:00 observation and the corresponding model values. Measurement points were binned based on their elevation above the ground in 500 m intervals. The first bin includes all measurements below 1.5 km and the last bin includes all measurements above 3 km.

Complex local topography can cause localized transport patterns in the domain, which cannot be resolved at the model's 4 km × 4 km horizontal resolution. Pfister et al. (2017a) discuss the impacts of the complicated wind patterns in the NFR and the limitations of WRF-Chem simulations in capturing the transport during FRAPPÉ campaign in detail. To reduce the impact of localized influences on the sensitivity analysis we use airborne measurements which better represent the regional picture. The evaluation of modeled ethane concentrations with aircraft data provides information on the impact of different configurations on the transport of oil and NG emissions. Box and whisker plots of ethane concentrations at different elevations along the C-130 morning and afternoon flights are shown in Fig. 4. These plots limit the C-130 observation to the NFR region (east of −105.2 longitude) to reduce transport errors, and separate observations collected during the 09:00 to 12:00 (AM) and 12:00 to 18:00 (PM) flights to account for the diurnal changes in the PBLH. For this comparison, the hourly model output has been interpolated to the time and location of each 1 min average observation. Lower concentrations of ethane were measured during the PM flights compared to AM flights because of the higher PBLH and stronger vertical mixing in the afternoon (Fig. 3). Table 5 summarizes the mean and NMB of the ethane concentration for all simulations using ethane airborne measurements. In all simulations, the ethane concentrations are underpredicted by up to 3.3 ppb (NMB ranges between −63 % and −42 %) for the C-130 AM flights and up to 1.7 ppb (NMB ranges between −47.6 % and −29.5 %) for the C-130 PM flights. Overall, measured ethane concentrations, absolute biases, and absolute NMBs are higher for C-130 AM compared to C-130 PM. However, the differences between variability in NMBs for C-130 AM and C-130 PM are small, i.e., 21 % and 18.1 %.

Measurements from P-3 spirals focus on smaller regions and can capture the impact of local emissions. Figure 5 compares the vertical distribution of measured ethane concentrations against the corresponding model values (interpolated to the time and location of each 1 min average observation) for all the simulations at BAO and Platteville (PAO) spirals. Both sites are located close to oil and NG sources (Fig. 1), although urban emissions from Denver region can also reach BAO (Pfister et al., 2017a). Similar to C-130 observations, we illustrate the morning and afternoon data separately. Mean concentrations of up to 18.6 ppb (SD 2.8 ppb) were measured by P-3 aircraft, but these high values were not captured by the model and resulted in biases up to −14.9 ppb (NMB of −80.5 %) at PAO spirals and biases up to −7.16 ppb (NMB of −57.8 %) at BAO spirals. Similar to C-130 flights, higher measured ethane concentrations, absolute biases, and NMBs are observed for P-3 AM flights compared to PM flights. Higher absolute biases and larger variance at lower altitude in AM flights may be due to larger uncertainties in capturing the morning evolution of the PBL. Variability in NMBs across simulations are greater in the PM spirals (42.8 % at PAO and 57.3 % at BAO) compared to the AM spirals (36.5 % at PAO and 31.3 % at BAO).

While the model shows difficulty in representing the absolute magnitude in the ethane concentrations in all simulations at lower altitudes, most simulations capture the changes in the variance of ethane concentrations from lower to higher altitudes well – especially for the C-130 and P-3 BAO flights. The C-130 flights covered a larger region with varying flight patterns across the NFR; thus, less variability in the modeled ethane concentrations was observed in the C-130 flights compared to the P-3, which flew a repetitive pattern and the repeated spirals over the key surface locations that reflect a higher influence from localized emissions.

3.2 Sensitivity to planetary boundary layer parameterization

We evaluated the sensitivity of WRF-Chem meteorological fields and ethane concentrations to a non-local (YSU) and two local (MYJ and MYNN3) PBL schemes in the PBL1, PBL2, and PBL3 simulations, respectively. Table 2 includes details regarding the simulation configurations. The temperature at BAO changed little between the different PBL schemes and the model agrees with observations (Fig. 2). At all three altitudes, PBL1 had a small positive bias (errors less than 1 ^∘C) while PBL2 and PBL3 had a small negative bias (errors less than 1 ^∘C) (Tables 3 and S2). Relative humidity differed slightly between local and non-local PBL parameterizations. PBL1 captured relative humidity well, especially at lower altitudes (mean bias of 0.38 %, 1.47 %, and 4.93 % for 10, 100, and 300 m, respectively). PBL2 and PBL3 both overpredicted the relative humidity at all altitudes. The mean bias for PBL2 and PBL3 ranged from 11.12 % to 14.78 % and 6.61 % to 9.55 %, respectively.

Figure 5 Vertical distribution of ethane at the PAO (a, b) and BAO (c, d) sites measured during the P-3 spiral flights and the corresponding model values. Flights are separated by flight time. (a) and (c) show the P-3 AM flights that include 09:00 to 12:00 observations and the corresponding model values. (b) and (d) show P-3 PM that include 12:00 to 18:00 observations and model values. Measurement points were binned based on their elevation above the ground in 500 m intervals. The first bin includes all measurements below 1.5 km and the last bin includes all measurements above 3 km.

At all altitudes at BAO, PBL1 predicted higher wind speeds than observed in PBL2 and PBL3 (Figs. 2 and S1–S3). Wind direction does not vary significantly between PBL1, PBL2, and PBL3 at BAO tower and the model missed the southerly component of the afternoon winds. Figure S4 shows the 10 m average wind speed (from 1 to 11 August) in PBL1, PBL2, and PBL3 for daytime and nighttime and compares it with measurements. Higher daytime wind speed was predicted by PBL1 in the Colorado eastern plains, especially north of Denver and close to oil and NG operations. Figure S5 shows the averaged diurnal cycle of wind speed and wind direction at WC-Tower and PAO (sites close to oil and NG operations). At WC-Tower and PAO, PBL2 and PBL3 better captured the southerly component of afternoon winds compared to BAO.

Each PBL scheme in the WRF model uses different diagnostics to determine the PBLH. To have a consistent comparison of PBL height in the three simulations, we used the 1.5-theta-increase method to estimate the PBL height. In this method, the PBLH is the lowest altitude where the difference between minimum potential temperature and potential temperature is greater than 1.5 K (Hu et al., 2010; Nielsen-Gammon et al., 2008). Figure 6 shows the diurnal evolution of the PBLH as calculated using the 1.5-theta-increase method in the simulations. Observed PBLH at the PAO, FC, and Golden sites were retrieved from micro-pulse lidar backscatter profiles using covariance wavelet transform (CWT) (Compton et al., 2013). The PBLH in the PBL1 simulation is greater than PBL2, PBL3, and observations, and the bias is largest in the afternoon. Figure 3a–c shows that PBL1 distributed ethane higher into the atmosphere and that more dilution resulted in a lower ethane concentration within the PBL. Figure S6 shows up to 5 ppb higher surface ethane concentrations, on average, in simulations based on local PBL schemes (PBL2 and PBL3) compared with the simulation based on the non-local PBL scheme (PBL1).

Figure 6 Diurnal evolution of the PBL in the MYJ, MYNN3, and YSU schemes at the PAO (a), FC (b), and Golden (c) sites. PBLH was measured using micro-pulse lidar backscatter profiles during the daytime. Error bars represent the standard error.

The high bias in temperature, wind speed, and PBLH in PBL1, the non-local PBL scheme, suggests a strong vertical mixing that is more defined in the Colorado eastern plains and close to oil and NG activities. The local PBL schemes (i.e., PBL2 and PBL3) predict cooler and moister climates and a lower PBLH, which indicates less vertical mixing. This is consistent with previous works that compared local and non-local PBL schemes in the WRF model (Angevine et al., 2012; Hu et al., 2010).

The comparison between C-130 airborne measurements and modeled ethane concentrations across the NFR, as illustrated in Fig. 4, shows biases between −2.5 and −2.3 ppb for AM flights and between −1.7 and −1 ppb for PM flights. Lower NMB variability (4 %) was observed in the C-130 AM flight with the NMB ranging from −43.1 % to −47.1 % compared with C-130 PM that had a NMB variability of 18 % and NMB values ranging from −29.5 % to −47.6 %. Similarly to the C-130 comparison (Fig. 5) the simulations did not capture the high ethane values measured during P-3 BAO and P-3 PAO spirals. The sensitivity of modeled ethane profiles to the PBL scheme is larger in P-3 flights compared to C-130 flights, with a NMB variability of 14.1 % ranging from −58 % to −44 % for PAO AM flights and a NMB variability of 32.4 % ranging from −37.3 % to −69.7 % for the PAO PM flight. On average PBL1 predicted higher ethane concentrations during AM flights at lower altitudes compared with PBL2 and PBL3 (Fig. 3). Faster evolution of the morning PBL and stronger vertical mixing in PBL1 lofted pollutants (including ethane) higher into the atmosphere in the morning (Fig. 6). The rapid growth of the morning PBL in PBL1 resulted in higher concentration of ethane at higher altitudes (0.5 to 2 km) compared with PBL2 and PBL3.

3.3 Sensitivity to re-initialization

We investigated the impact of the daily initialization of meteorological fields on the model performance in capturing the transport of pollutants. For this, we conducted a sensitivity simulation (Init5) in which each daily cycle started at 18:00 UTC from ERA-Interim meteorological fields and ran for 30 h. In the comparison free-running simulation, Init4, we initialized the model at the first time step using the ERA-Interim model and ran the simulation from 24 July to 18 August 2014 freely. Physical configurations and meteorological and chemical initial and boundary conditions were kept the same for these two simulations (Table 2). Figure 2 shows an up to 3 ^∘C bias in the nighttime temperature in Init5, but good agreement with the measured temperature during the day. Init4 showed better skill in capturing the nighttime temperature compared to Init5, but predicted the lowest daytime temperature among all the simulations with a bias up to −3 ^∘C. On average, the NMB of the temperature at BAO 100 m is between 8.6% in Init5 and −6.0 % in Init4 (Table 3), which is the largest variability in the NMB for temperature across the simulations. Similar to the temperature, relative humidity showed a strong sensitivity to re-initialization. Init4 predicted the highest relative humidity, with a NMB of 39.2 % and Init5 predicted the lowest relative humidity, with a NMB of −26.5 % among the simulations at BAO 100 m (Table 3). Nighttime wind direction at BAO (Fig. 2), PAO, and WC-Tower (Fig. S2) had a strong southerly component in Init4 compared to Init5 and observations. In addition, Init4 predicted higher wind speeds compared with BAO measurements (Fig. 2) and Init5. Figure 7 shows higher wind speed on average at a 10 m altitude across the domain in Init4 compared with Init5 and measurements during both daytime and nighttime.

Figure 7 Measured (circles) and modeled (colored contours) wind speed at 10 m captured by Init4 (a, b) and Init5 (c, d) from 1 to 11 August 2014 and separated by daytime vs. nighttime.

When compared to C-130 AM ethane concentrations (Fig. 4), Init4 predicted the lowest ethane concentrations (a bias of −3.3 ppb and a NMB of −63 %) among all the simulations. This is likely due to the high bias in of wind speed in this simulation which resulted in lower concentrations of ethane (Fig. 7). The ethane bias is $\sim -\mathrm{2.5}$ ppb and the NMB is −47.9 % in Init5 during C-130 AM. Concentrations during the C-130 PM flights showed a weak sensitivity to re-initialization with the NMB ranging from −37.8 % (Init4) to −40.1 % (Init5). For the P-3 BAO and P-3 PAO spirals in both the AM and PM flights, Init4 had the lowest ethane values compared to all of the other simulations and compared to observations (Fig. 5). This resulted in the largest NMB variability across the simulations. During PAO AM, the NMB ranges between −80.5 % for Init4 and −53.2 % for Init5 (NMB variability of 27.3 %) and during PAO PM, the NMB ranges between −72.9 % for Init4 and −30.0 % for Init5 (NMB variability of 43.9 %).

3.4 Sensitivity to meteorological initial and boundary condition

We tested the performance of changing the meteorological initial and boundary conditions by comparing simulations using ERA-Interim (Met5) with simulations using NCEP-FNL (Met6). As was done for Met5, we initialized meteorological fields with the reanalysis fields every day allowing for a 6 h spin-up. To prepare meteorological initial and boundary conditions from global models, WRF interpolates these outputs to the designed domains. Figure S7 illustrates the differences between ERA-Interim and NCEP-FNL model outputs interpolated to the outer domain at the lowest model level and averaged from 1 to 15 August 2014. Overall, the wind speed predictions by these two global models are very similar with slightly (less than 1 m s⁻¹) higher prediction by NCEP-FNL. ERA-Interim and NCEP-FNL had larger discrepancies in temperature and relative humidity throughout the domain. Comparison with BAO observations (not shown) indicates similar performance for both models with somewhat lower temperature and higher relative humidity in ERA-Interim compared with NCEP-FNL. However, these discrepancies did not have a large impact on temperature and relative humidity in the WRF-Chem simulation. Figures 2 and S1–S3 indicate that the performance of the two simulations is comparable in capturing temperature and relative humidity with a better agreement with measurements during the day. Met5 had slightly higher temperature and lower relative humidity compared to Met6 and compared better to measurements especially during the night. This is because WRF-Chem only uses the global values as the initial and boundary values and resolves for atmospheric variables such as temperature and relative humidity in high resolution based on physical parameterizations set for the simulation.

The comparison of ethane measurements by the C-130 and P-3 aircraft with Met5 and Met6, shown in Figs. 4 and 5, respectively, also reflects an overall low sensitivity of the model performance to meteorological initial and boundary conditions for both the AM and PM flights. High sensitivity was observed during the P-3 PAO PM flight with an ethane NMB variability of 23.9 % where Met5 had a bias of −2.6 ppb (NMB of −30 %) and Met6 had a bias of −4.7 ppb (NMB of −53.9 %).

Table 5 Ethane mean, NMB, and NMB variability from C-130 and P-3 BAO and PAO airborne measurements below 2000 m and the corresponding model values.

3.5 Sensitivity to horizontal resolution

The two nested domains in the Hor5 simulation had a horizontal resolution of 12 km × 12 km (coarse) and 4 km × 4 km (fine). The one-way nesting method was used to prevent any feedback from the higher resolution inner domain on the outer domain. This means that while the outer domain provides the lateral boundary conditions to the inner domain, the higher-resolution fields from the inner domain do not alter the outer domain fields. To compare the impact of horizontal resolution, we compared the performance of the coarse domain with the fine domain in the same simulation (Hor5). Temperature and relative humidity did not show significant sensitivity to the horizontal resolution at BAO and PAO, and nor did surface winds at BAO (Fig. 2), PAO, or WC-Tower (Fig. S4). At altitudes of 100 and 300 m at BAO, the coarse domain predicted higher nighttime wind speeds compared to the fine domain and the measurements.

Averaged ethane concentrations along the C-130 flights (Fig. 4) do not vary significantly with horizontal resolution. However, higher differences are observed for the P-3 spirals. This might be due to the C-130 flights covering a larger area, which could average out the impact of horizontal resolution, whereas the P-3 spirals capture small-scale transport patterns in the domain more effectively. For the P-3 spirals (Fig. 5), the ethane NMB during BAO PM is +11.6 % for the fine domain and −21.6 % for the coarse domain. These values are −30 % and −55.1 % during the PAO PM flights for the fine and course domains, respectively.

Figure 8 Mean and mean bias ethane (a, b), CO (c, d), HC3 (e, f), and TOL (g, h) concentrations along the C-130 PM flights are limited to measurements below 2000 m a.g.l. and grids with more than four measurement points. The outline of Denver County (DEN) and the locations of BAO and PAO are marked on the underlying terrain map.

3.6 Oil and NG emission in the NFR

We assessed the performance of the model in capturing oil and NG emissions by focusing on ethane, which is mostly emitted from oil and NG emission sources, and on species with multiple emission sources such as CO and other VOCs. To investigate the contribution of oil and NG emissions to NFR air quality, we ran two additional simulations: in one, the emissions were based on the NEI-2011 as provided (base simulation or Em7), and in the other we doubled the oil and NG emissions (perturbed simulation or Em8).

Figure 8 shows the C-130 PM measurements and bias limited to altitudes below 2000 m and Fig. 9 displays scatterplots of the measured to modeled species concentrations limited to the NFR, below 2000 m, and with measured ethane greater than 2 ppb. Figure 8a and b illustrate high ethane concentrations in the vicinity of oil and NG facilities which were not captured by the model resulting in low biases. As can be expected, the simulated ethane concentrations show a high sensitivity to changes in the oil and NG emissions (Fig. S8). The highest sensitivity was observed for measurements taken over regions close to oil and NG sources, such as the P-3 PAO spirals. Ethane biases between Em7 and Em8 varied from −9.4 to −1 ppb (NMB from −50.8 % to −5.5 %) during P-3 PAO AM, and from −2.7 to +2.8 ppb (−31.2 % to +31.8 %) during P-3 PAO PM. Doubling oil and NG emissions lowered the absolute bias during the AM flights (NMB from −50.8 % to −5.5 %) but resulted in an overestimation of the ethane concentrations during the PM flights (NMB from −31.2 % to +31.8 %). One possible reason for the difference between the AM and PM biases might be the incorrect representation of the diurnal variation of ethane emission rates in NEI-2011. An inverse modeling technique, as will be subject of further studies, can be used to calculate optimal scaling factors for hourly ethane emissions with the goal to minimize the discrepancies between model and measurement.

Figure 9 Scatterplot of measured vs. corresponding model values of ethane (a), CO (b), HC3 (c), and TOL (d) along the C-130 PM flights which are limited to measurements in the NFR below 2000 m. Red diamonds represent the Em7 (base emissions) and blue circles represent Em8 (perturbed emissions). Red and blue lines show the best fit using the least square linear regression method for Em7 and Em8, respectively. Grey lines show the 1:1 lines.

CO is mostly emitted from combustion processes and is released from many different source sectors. Figure 8c shows CO enhancements over both Denver and oil and NG facilities. Biases along the C-130 flight tracks (Fig. 8d) show an overprediction of CO over Denver and west of Denver and an underprediction over the oil and NG facilities. The scatterplot in Fig. 9b reflects that an overall low bias in modeled CO that could be partly due to errors in capturing background CO. Doubling oil and NG emissions in Em8 only marginally increased the slope of the regression line indicating a low sensitivity of CO in the NFR to oil and NG emissions. This suggests that the source of the low bias in CO is likely related to other source categories and/or the model lateral boundary conditions.

In the RACM chemical mechanism, alkanes such as propane, n-butane, isobutane, and acetylene (ethyne), and alcohols such as methanol and ethanol are lumped under the “HC3” group (Stockwell et al., 1990). We compared the simulated HC3 concentrations with the sum of measured chemicals in the HC3 group during the C-130 flights. Similar to ethane, the highest values of HC3 were measured over oil and NG facilities (Fig. 8c). These enhancements were not captured in the model and resulted in low model biases (Fig. 7f). Comparison of measured HC3 with modeled values from Em7 and Em8 (Fig. 9c) confirms the low bias of HC3 and shows some increase in the slope of the regression line in Em8, albeit less pronounced compared to ethane.

Toluene and benzene are lumped together in the RACM chemistry under “TOL” (Stockwell et al., 1990). We compared simulated TOL with the sum of toluene and benzene concentrations observed during the C-130 flights. The transport sector is a strong source of toluene and benzene in the NFR as well as oil and NG activities. TOL enhancements were observed over oil and NG facilities and over Denver with higher values associated with oil and NG emissions (Fig. 8g). The model did not capture the enhancements in regions influenced by oil and NG emissions, but captured TOL values over Denver well (Fig. 7h). TOL showed very low sensitivity to perturbed oil and NG emissions as shown in Fig. 9d. TOL emissions from the oil and NG sector in the emission inventory used in this study (NEI-2011) were very low; thus, doubling oil and NG emissions did not increase TOL in the Em8 simulation. Similar to toluene and benzene, xylene enhancements were measured over oil and NG facilities and over Denver. The model underestimated xylene enhancements over oil and NG activities and overestimated these enhancements over Denver. Em8 with doubled oil and NG emissions showed very similar performance to Em7 which indicates low emission rates of xylene from the oil and NG sector in NEI-2011 (not shown).

Figure 10 illustrates the HC3 to TOL ratio measured along the C-130 PM flight, limited to the NFR region and altitudes below 2000 m, in addition to the corresponding model values. Figure 10a shows oil and NG influenced points with enhanced measured ethane (concentrations greater than 2 ppb). HC3 to TOL ratios in oil and NG influenced locations show inconsistency between measured (HC3 ∕ TOL = 68) and Em7 modeled ratios (HC3 ∕ TOL = 22) which was improved in the Em8 (HC3 ∕ TOL = 40.9). However, doubling oil and NG emission still resulted in underestimations of HC3, TOL, and their ratios in this region. Figure 10b shows urban influenced points with low measured ethane (concentrations less than 2 ppb). Modeled HC3 to TOL ratios (7.3 for Em7 and 8.9 for Em8) in the urban influenced locations did not show large sensitivity to oil and NG emissions and agreed well with the measurements (10.2). In both oil and NG and urban influenced regions models predicted lower than measured y intercepts which was not improved in Em8. Figure 9c also confirms the low bias (about −2 ppb) in background HC3 in the model. One reason for this offset could be the underestimation of the HC3 concentration in the lateral boundary condition fields or leakage from the NG distribution system which was not captured in the model.

The results suggest that HC3, toluene, benzene, and xylene from the oil and NG sector are significantly underestimated in NEI-2011. The low model bias for these species is more pronounced compared with the low model bias in ethane (Fig. 9). The inconsistency between these biases implies that NEI-2011 emission ratios might need to be changed and HC3, toluene, benzene, and that xylene oil and NG emissions would need to be increased by a larger factor than ethane.

4 Conclusion

We used WRF-Chem to understand the sensitivity of pollutant transport at a high horizontal resolution to different model configurations with a focus on oil and NG emissions. By conducting a range of different sensitivity simulations, we assessed the variability of meteorological variables such as temperature, relative humidity, and wind fields as well as of ethane concentrations (used as a tracer for the oil and NG sector) to different model configurations and parameterizations. The overall daily cycle of temperature and relative humidity was captured well in the simulations with NMB values ranging from −3.9 % to 11.1 % for temperature and from 29.7 % to 52.6 % for relative humidity. All simulations showed good skill in capturing daytime wind fields but showed higher biases for nighttime wind speeds.

Figure 10 Scatterplot of HC3 vs. TOL concentrations along the C-130 PM flights limited to measurements in the NFR below 2000 m altitude. (a) HC3 vs. TOL (when measured ethane is greater than 2 ppb) for measured and corresponding model values. (b) HC3 vs. TOL (when measured ethane is less than 2 ppb) for measured and corresponding model values. Grey circles represent measurements, red diamonds represent the Em7 (base emissions), and blue circles represent Em8 (perturbed emissions). Grey, red and blue lines show the best fit using least square linear regression method for observations, Em7, and Em8, respectively.

Table 5 summarizes the mean and NMB for ethane concentrations from C-130 and P-3 airborne measurements below 2000 m a.g.l. and the corresponding model values for all sensitivity tests. A significant underestimation of ethane in all simulations – especially in regions close to oil and NG activities – with biases up to −14.9 ppb (NMB up to −80.5 %) suggest that the emission inventory used (NEI-2011) underpredicts oil and NG emissions. The NMB variability (Table 5) was used as a proxy for variability in the model performance caused by model configurations. The NMB of the near-surface ethane concentration for aircraft flight patterns across sensitivity simulations varied by up to 57.3 % for P-3 BAO, by up to 42.8 % for P-3 PAO, and by up to 21.1 % for C-130 flights. The lower NMB variability during the C-130 flights may be due to the larger area coverage by this aircraft during the FRAPPÉ campaign and the irregular flight patterns. P-3 spirals, covering smaller regions within the domain during repetitive flight patterns, focused more on the local emissions and smaller-scale transport patterns and captured a larger ethane sensitivity to model configurations. The largest sensitivity occurred in the initialization test (comparing daily re-initialization with free-run simulation) with an ethane NMB variability of up to 57.3 %, followed by the horizontal resolution test (comparing horizontal resolution of 12 km × 12 km with 4 km × 4 km), and the PBL parameterization test (comparing local with non-local PBL schemes) with ethane NMB variabilities of up to 33.3 % and 32.4 %, respectively. To further investigate the performance of the model in capturing oil and NG emissions in the NFR we used a similar domain setup with 12 km × 12 km and 4 km × 4 km horizontal resolutions for outer and inner domains, respectively, daily re-initialization of meteorological variables with ERA-Interim model, and the MYNN3 PBL scheme.

We compared measured ethane, CO, lumped alkanes (HC3), lumped toluene and benzene (TOL), and xylene to corresponding modeled values and assessed the changes in the model performance when doubling oil and NG emissions. The model showed an underprediction of ethane with the original inventory and a strong sensitivity of ethane concentrations to oil and NG emissions. Doubling the oil and NG emissions resulted in an improvement during AM flights and an overestimation of ethane during the PM flights which suggests the possible incorrect representation of the diurnal variation of ethane emission rates in NEI-2011. The model tends to overestimate CO over the Denver region and underestimates CO over the oil and NG region. The low sensitivity of CO to oil and NG emissions indicates that CO in the region is predominantly emitted from sources other than oil and NG. Enhancements of HC3, TOL, and xylene over oil and NG facilities were not fully captured in the model and resulted in low biases. Doubling emissions from oil and NG emissions improved the model performance in capturing HC3, but still resulted in a low model bias. Although high values of TOL and xylene were measured over oil and NG facilities, the model did not capture these enhancements in either the simulations with base NEI-2011 emissions or in the simulations with the doubled oil and NG emissions. The inconsistency between the sensitivity of ethane, HC3, benzene, toluene, and xylene to the increase in oil and NG emissions and the mismatch between VOC ratios in the model and measurement suggest that oil and NG emission rates in NEI-2011 need to be scaled differently for these species. VOC ratios in the measurements can be used to update these ratios in the emissions inventory.

We recognize that using ethane as a tracer to assess the sensitivity of the WRF-Chem model to physical parameters can be limited by the biases in the emission inventory. Conducting WRF-Chem simulations using different physical parameterization and using NMB variability can help to reduce this limitation. The results presented reflect the challenges that one is faced with when attempting to improve emission inventories by contrasting measured with modeled concentrations, either through simple direct comparisons or more advanced methods, such as inverse modeling. Any uncertainties that arise from the model configuration will translate into the derived emission constraints, and it is important to be aware of the uncertainties resulting from different model setups. The WRF-Chem simulations and the knowledge gained from this study will be used to support inverse modeling studies aimed at improving estimates of emission from the oil and NG sector in the NFR.