Updates to the Noah Land Surface Model in WRF-CMAQ to Improve Simulated Meteorology, Air Quality, and Deposition.

Regional, state, and local environmental regulatory agencies often use Eulerian models to investigate the potential impacts on pollutant deposition and air quality from changes in land use, anthropogenic and natural emissions, and climate. The Noah land surface model (LSM) in the Weather Research and Forecasting (WRF) model is widely used with the Community Multiscale Air Quality (CMAQ) model for such investigations, but there are many inconsistencies that need to be changed so that they are consistent with dry deposition and emission processes. In this work, the Noah LSM in WRFv3.8.1 is improved in its linkage to CMAQv5.2 by adding important parameters to the WRF/Noah output, updating the WRF soil and vegetation reference tables that influence CMAQ wet and dry photochemical deposition processes, and decreasing WRF/Noah's top soil layer depth to be consistent with CMAQ processes (e.g., windblown dust and bidirectional ammonia exchange). The modified WRF/Noah-CMAQ system (both off-line and coupled) impacts meteorological predictions of 2-m temperature (T2; increases and decreases), 2-m mixing ratio (Q2; decreases), and 10-m wind speed (WSPD10; decreases) in the United States. These changes are mostly driven by leaf area index values and aerodynamic roughness lengths updated in the vegetation tables based on satellite data, with additional impacts from soil tables updated based on recent soil data. Improvements in the consistency in the treatment of land surface processes between CMAQ and WRF resulted in improvements in both estimated meteorological (e.g., T2, WSPD10, and latent heat fluxes) and chemical (e.g., ozone, sulfur dioxide, and windblown dust) model estimates.

Introduction

The Noah land surface model (LSM) is one of the most robust and well‐established models of its kind for meteorological and climate modeling. The Noah LSM is the result of long‐term, overlapping and continued development in land surface modeling initiated by the four agencies that comprise its name: National Center for Atmospheric Research (NCAR), Oregon State University, the U.S. Air Force, and National Centers for Environmental Prediction's (NCEP's) Office of Hydrology. The Noah LSM has been widely developed, applied, and evaluated for more than three decades (e.g., Chen & Dudhia, 2001; Chen et al., 1996, 1997, 2007; Ek et al., 2003; Li et al., 2013; Mitchell et al., 2004; Niu et al., 2011; Pan & Mahrt, 1987; Yin et al., 2015). The Noah LSM includes four soil temperature and soil moisture layers, as well as a vegetation canopy model, snow prediction, evapotranspiration, and soil drainage and runoff. It may be used either as a stand‐alone model or embedded within a parent atmospheric model to simulate the radiative and heat fluxes between the ground and the atmosphere. This coupling allows the land surface to simultaneously respond to short‐term weather events and long‐term climate changes, while the land heat and water storage anomalies can impact climate predictability. The Noah LSM is currently used in operational NCEP's weather and climate prediction models (i.e., Global Forecast System and North American Model), at the Air Force Weather Agency (Zheng et al., 2014), and across the weather and climate modeling research communities in both retrospective and future scenario cases (e.g., Gao et al., 2012; Salathé et al., 2010; Trail et al., 2013; Wang & Kotamarthi, 2015). The Noah LSM has also been used in both short‐term weather and air quality forecasting in National Oceanic and Atmospheric Administration‐Environmental Protection Agency (NOAA‐EPA)'s National Air Quality Forecasting Capability (NAQFC) since 2007 and 2009 for ozone (O3) and fine particulate matter (PM2.5) predictions, respectively (Lee et al., 2017; Mathur et al., 2008; Stajner et al., 2012).

In the Weather Research and Forecasting (WRF) model (Powers et al., 2017; Skamarock & Klemp, 2008), there are several options to extend the Noah LSM to tailor it for additional snow and ice processes, to improve representations of urban areas, to increase subgrid representation of various land use (LU) classifications, and to add complexity in simulating hydrological processes.

The WRF model with the Community Multiscale Air Quality Model (CMAQ; hereafter WRF‐CMAQ) is a robust modeling system that allows for chemistry and physics of pollutant transport and fate to be solved either sequentially (i.e., off‐line) as described by Byun and Schere (2006) or simultaneously (i.e., coupled) using the framework of Wong et al. (2012). In an off‐line WRF‐CMAQ framework, the WRF model is run and processed by the Meteorology‐Chemistry Interface Processor (MCIP) (Otte & Pleim, 2010), and then CMAQ is run using the MCIP output. The coupled WRF‐CMAQ framework exchanges information between WRF and CMAQ throughout the WRF simulation and can optionally include two‐way CMAQ aerosol feedbacks on WRF meteorology (see Figure 1 in Wong et al., 2012, and discussion therein). CMAQv5.0 and v5.1 were extensively evaluated (Appel et al., 2013, 2017), and CMAQv5.2 (Environmental Protection Agency, EPA, 2017) included major updates to photochemistry, photolysis rates, aerosol processes, aqueous and heterogeneous chemistry, lightning interactions, and transport processes. The Noah LSM option cannot accurately be used within the coupled WRF‐CMAQ model system, as it contains missing information that is pivotal to calculating atmospheric deposition and returns erroneous air quality (e.g., O3) concentrations. While the Noah LSM is available in the off‐line WRF‐CMAQ system, the necessary information is parameterized and adds additional uncertainty in air quality predictions. Thus, the current study seeks to improve the WRF‐CMAQ model from advanced treatments of land surface processes and land‐atmosphere interactions in WRF/Noah, while more firmly establishing both the off‐line and coupled WRF‐CMAQ connection to better support the large number of applications that use the Noah LSM in WRF, thus expanding the models use in retrospective, forecasting, and climate applications.

While many WRF‐CMAQ meteorological‐air quality studies use the Pleim‐Xiu (PX) LSM in WRF (Noilhan & Planton, 1989; Pleim & Xiu, 1995; Xiu & Pleim, 2001), the WRF/PX‐CMAQ model is not suitable for examining future climate/air quality applications as it relies on constraining the surface fluxes using observed near‐surface moisture and temperature (Pleim & Gilliam, 2009; Pleim & Xiu, 2003). Thus, the WRF/PX‐CMAQ configuration is best suited for retrospective applications (Gilliam & Pleim, 2010). On the other hand, the off‐line WRF/Noah and CMAQ models (hereafter WRF/Noah‐CMAQ) have been applied to investigate future climate‐air quality changes over the United States (e.g., Campbell et al., 2018a, 2018b; Fann et al., 2015; Gao et al., 2013; Gonzalez‐Abraham et al., 2015; Spero et al., 2016; Zhang et al., 2016). The WRF/Noah and CMAQ models are also important components of an integrated energy system with agricultural, meteorological, and air/water quality models to investigate climate, emissions, and land use (LU) change impacts on future pollutant deposition and nutrient cycles (i.e., a one‐biosphere approach), which then can also include economic drivers, human health, and ecosystem endpoints (Cooter et al., 2016).

In this work we modify both the off‐line and coupled (without aerosol feedbacks) versions of WRF/Noah version 3.8.1 (WRFv3.8.1) and CMAQ version 5.2 (CMAQv5.2) (EPA, 2017) to increase its utility for air quality and ecosystem deposition studies. This work demonstrates that the modified off‐line and coupled WRF/Noah‐CMAQ model system improves overall model performance compared to the unmodified off‐line WRF/Noah‐CMAQ system for some meteorological and chemical variables. The following sections present the WRF/Noah‐CMAQ model configurations and simulation design used to test the modified off‐line and coupled systems (section 2), the motivation, details, and impact assessments for the individual WRF/Noah‐CMAQ modifications (section 3), the evaluation of the modified WRF/Noah‐CMAQ air quality predictions (section 4), and conclusions (section 5). A list of the commonly used acronyms in this paper (Text S1), as well as additional figures (Figures S1–S11) and tables (Tables S1–S5) used to supplement the main analyses, are also included in the supporting information.

Methodology and Model Setup

The model configuration and simulation design are motivated by the individual modifications to the off‐line and coupled WRFv3.8.1‐CMAQv5.2 model system (hereafter referred to as WRF/Noah‐CMAQ), which include using stomatal (R ) and aerodynamic resistance (R ) from WRF/Noah output directly into CMAQ, updates to the WRF/Noah vegetation and soil input tables and their influence in CMAQ, and tests of modified WRF/Noah soil layer depths to better align with those used in CMAQ processes.

The effects of using a tiled (i.e., mosaic) LU approach in WRF/Noah (Li et al., 2013), and its impacts on CMAQ predictions, are also compared to the dominant LU approach. In short, the mosaic approach accounts for subgrid‐scale land cover characterization by way of a user‐selected number of tiles (N) within each model grid cell, which contrasts with the dominant LU approach that only accounts for the single most abundant LU tile. Here the N mosaic tiles is set to eight from independent tests (not shown) of this model configuration/domain (discussed below) that indicate N = 8 results in about 97% of all model grid cells having 99% of their LU categories represented. The normalized LU fraction for the N = 8 tiles is then used in the area‐weighted average of the Noah mosaicked diagnostic variables for each model grid cell (i.e., the eight LU types occupy 100% of grid cell). Hence, in this work, the WRF/Noahv3.8.1 mosaic code was also modified to provide LU fraction weighted averages for R and leaf area index (LAI), which are important inputs to CMAQ. Although there is no Noah mosaicked R variable, it is still impacted by mosaicked roughness length calculations (more details in section 3).

As part of this work, the canopy moisture content (CMC) in WRF/Noah‐CMAQ was also corrected because WRFv3.8.1 contains an important unit conversion error. While this change is benign for most WRF/Noah applications, it can be significant in WRF/Noah‐CMAQ for dry deposition of gas species that are highly soluble and depend strongly on surface wetness (e.g., sulfur dioxide; SO2). Consequently, this error in WRFv3.8.1 estimates canopy moisture to be 3 orders of magnitude too small and results in substantially underpredicted dry deposition of soluble gaseous species in CMAQ. In all subsequent analyses, the modified off‐line and coupled WRF/Noah‐CMAQ model reflects the CMC × 103 correction, which has been added to the upcoming release of WRFv4.

Finally, when a more realistic LAI and vegetation coverage is used in WRF, CMAQ exhibits a high O3 bias in sparsely vegetated areas that has been attributed to the O3 deposition to soil mechanism (Ran et al., 2016). To ameliorate this bias, the resistance of O3 deposition to soil from Mészáros et al. (2009) was adopted. This model more realistically scales the resistance between 200 and 500 s/m by the predicted soil moisture content as opposed to a constant resistance of 667 s/m with the unmodified CMAQ. Overall, the spatial impact of the Mészáros et al. (2009) parameterization is to almost linearly increase dry deposition of O3 to soil as soil moisture decreases (not shown), which can help reduce high O3 biases. This change has also been applied to all modified off‐line and coupled WRF/Noah‐CMAQ simulations in this work.

Model Configuration

The modified off‐line and coupled (one‐way) WRF/Noah‐CMAQ model systems (see section 1 for discussion of off‐line and coupled systems) are tested in an application over North America (centered over the contiguous United States; CONUS) at a horizontal resolution of 12 × 12 km (see Figure 1 for domain boundaries) spanning a simulation period of 21 May to 31 August 2011. Table 1 provides the WRF‐CMAQ domain specifications, model configurations, meteorological and chemical initial/boundary conditions, anthropogenic and biogenic emissions, and the pertinent references.

Figure 1

Average (21–31 May 2011) incremental changes due to the modified SOILPARM table in WRF/Noah‐CMAQ (runs A4 and A3) for (a) soil_H2O, (b) LH_flx, (c) Q2, (d) TSLB, (e) SH_flx, and (f) T2. SOILPARM = soil hydraulic parameter table; WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; soil_H2O = soil water content; LH_flx = latent heat flux; Q2= 2‐m mixing ratio; TSLB = topmost soil layer temperature; SH_flx = sensible heat flux; T2 = 2‐m temperature.

Table 1

Off‐Line and Coupled WRFv3.8.1‐CMAQv5.2 Model Domain, Configurations, Inputs and References

Model attribute	Configuration	Reference (s)
Simulation period	21 May to 31 August 2011 (10‐day spin‐up)	n/a
Domain	Continental United States; Center = 40°N, 97°W	n/a
Horizontal resolution	12 km	n/a
Vertical resolution	35 layers from surface to 50 hPa	n/a
Land use/Cover data	IGBP‐Modified MODIS 20‐category	https://modis.gsfc.nasa.gov/; Friedl et al. (2002, 2010)
Meteorological ICs and BCs	0.141° × ~0.141° ECMWFa‐ Operational Model Analysis ‐ 2011	https://rda.ucar.edu/datasets/ds113.0/
Chemical ICs and BCs	Goddard Earth Observing System model with Chemistry (GEOS‐Chem)	http://acmg.seas.harvard.edu/geos/
Anthropogenic Emissions	U.S. EPA 2011 National Emissions inventory	http://www.epa.gov/ttnchie1/eiinformation.html
Biogenic emissions	Biogenic Emissions Inventory System version 3.14	Schwede et al. (2005) and Vukovich and Pierce (2002)
Microphysics	Morrison 2‐Moment	Morrison et al. (2009)
PBL physics scheme	Yonsei University (YSU)	Hong (2010) and Hong et al. (2006)
Cumulus physics parameterization	Kain‐Fritsch	Kain (2004)
Shortwave and longwave radiation	Rapid Radiative Transfer Model for GCMs (RRTMG)	Clough et al. (2005) and Iacono et al. (2008)
Land surface model	Modified Noah (both dominant and mosaic LU)	Chen and Dudhia (2001), Ek et al. (2003), Li et al. (2013), Tewari et al. (2004), and this work
Surface layer	Revised MM5 Monin‐Obukhov (MO)	Grell et al. (1994), Jimenez et al. (2012), and Monin and Obukhov (1954)
Gas‐phase chemistry	CB05‐TUCL	Sarwar and Bhave (2007), Sarwar et al. (2006), Sarwar et al. (2008), Sarwar, Appel, et al. (2011), Whitten et al. (2010), and Yarwood et al. (2005)
Aqueous‐phase chemistry	CMAQ AQCHem updates	Alexander et al. (2009), Martin and Good (1991), Sarwar, Fahey, et al. (2011), and Sarwar et al. (2013)
Aerosol module/size	AERO6/3 modes	Appel et al. (2013)
Other model attributes	‐FDDA using analysis nudging	Deng et al. (2007) and Stauffer and Seaman (1994)
	‐Subgrid cumulus feedback	Alapaty et al. (2012) and Herwehe et al. (2014)
	‐In‐line photolysis	Binkowski et al. (2007)
	‐In‐line bi‐directional NH3 exchange	Bash et al. (2010, 2013), Cooter et al. (2012), and Massad et al. (2010)
	‐In‐line lightning NO emissions	Allen et al. (2012)
	‐ Physics‐based wind‐blown dust emissions scheme	Foroutan et al. (2017)
	‐Ocean/surf‐zone mask for sea‐salt emissions	Kelly et al. (2010)

Note. EPA = Environmental Protection Agency; WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; IGBP = International Geosphere‐Biosphere Programme; MODIS = Moderate Resolution Imaging Spectroradiometer.

European Centre for Medium‐Range Weather Forecasts, 2011‐Research Data Archive at the National Center for Atmospheric Research. https://doi.org/10.5065/D6ZG6Q9F. Accessed‐10 Feb 2017.

Simulation Design

To investigate the individual and combined impacts of the modifications to WRF/Noah‐CMAQ discussed above, the simulation design consists of 10 different simulation runs, which are divided into two major categories of short‐impact assessment and summer 2011 evaluation runs. Table 2 summarizes the 10 simulations, their duration, details, abbreviation, and respective analyses.

Table 2

Off‐Line and Coupled WRFv3.8.1‐CMAQv5.2 Simulation Design

WRF‐CMAQ simulation sets		Descriptions
Run #	Run details	Abbreviation	Analyses
Short Impact Assessments‐Set A = 21–31 May 2011 (5‐day spin‐up only for evaluation)
1	Offline WRF‐CMAQ with Unmodified Noah	A1	Unmodified/Default
2	Coupled WRF‐CMAQ with modified Noah R s and R a	A2	Impacts of coupled WRF/Noah‐CMAQ modifications
3	Same as A2 but with modified Noah Mosaic LU approach	A3	Impacts of coupled WRF/Noah mosaic modifications
4	Same as A3 but with updated Soil Hydraulic Table	A4	Impacts of soil hydraulic table updates
5	Same as A4 but with updated Vegetation Table	A5	Impacts of vegetation table updates
6	Same as A5 but with modified soil depth layers	A6	Impacts of “thin” (1 cm) soil layer 1 depth
7	Same as A6 but with in‐line wind‐blown dust module	A7	Test thin soil layer 1 depth on wind‐blown dust
Summer Evaluation Simulations‐Set B = 21 May to 31 August 2011 (10‐day spin‐up)
8	Offline WRF‐CMAQ with Default Noah Mosaic	B1	Unmodified/Default
9	Same as B1, but with all modifications in A2 ‐ A6	B2	Impacts of all mods in Runs A2 – A6 for off‐line WRF‐CMAQ
10	Same as B2, but for coupled WRF‐CMAQ	B3	Impacts of all mods in Runs A2 – A6 for coupled WRF‐CMAQ

Note. Table includes the abbreviations used in the text and additional details regarding their respective analyses. All modified runs A2–A6 and B2 and B3 contain the CMC × 103 correction and the Mészáros et al. (2009) parameterization for dry deposition of O3 to soil. WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality.

The incremental runs (A1–A7) are for 21–31 May 2011 and are used to demonstrate the impacts on meteorological, atmospheric dry deposition, and chemical variables for each WRF/Noah‐CMAQ modification, as well as a test of the modified WRF/Noah soil layer depths with the physics‐based windblown (WB) dust module in CMAQ (Foroutan et al., 2017). For the impact assessments (runs A1–A7), we analyze over the entire 21–31 May 2011 period, and there is no spin‐up period. The short simulation time and lack of spin‐up for the impact assessments were chosen due to the large computational expense of running numerous 12‐km coupled WRF‐CMAQ simulations over CONUS to investigate the incremental impacts of all model changes. Otherwise for model evaluation of the impact assessments, a relatively short 5‐day spin‐up period for runs A1–A7 is included (valid for 26–31 May 2011). The full suite of modifications is then evaluated for both the off‐line and coupled WRF/Noah‐CMAQ system (B1–B3) over the 2011 meteorological summer season (1 June to 31 August), using a 10‐day model spin‐up of 21–31 May. While a full year evaluation considering all seasons is ideal, under the scope of this study the summer season was chosen due to the importance of photochemical O3 production, enhanced heat and moisture fluxes from the surface, and the significant stomatal pathway of O3 and SO2 during the summer in CONUS. The results from the 2011 summer runs (section 4) also in part corroborate the choice of smaller spin‐up times for the impact assessment runs. Figure S1 in the supporting information shows the spatial plots of the unmodified off‐line WRF/Noah‐CMAQ (Table 2; A1) averaged between 21 and 31 May 2011, for all the meteorological, deposition, and chemical variables that undergo impact assessments due to the WRF/Noah‐CMAQ modifications presented in section 3.

We note that because the coupled WRF/Noah‐CMAQ does not work without using the modified R and R from WRF/Noah directly into CMAQ, these modifications (both with dominant and mosaic LU) are made first (A2 and A3). However, the incremental impact assessments in section 3 start by describing the updated soil (A4) and vegetation impacts (A5), as these have important WRF surface/meteorological impacts that influence R and R (and CMAQ deposition) in the evaluation of the combined off‐line and coupled WRF/Noah‐CMAQ modifications discussed in section 4 (B1–B3).

Observations and Evaluation Protocol

The impact assessment and summer 2011 simulations described in Table 2 are evaluated against observational networks that include measurements of meteorological states (temperature, moisture, and wind speed), soil conditions (moisture and temperature), surface fluxes (latent and sensible heat), and air quality concentrations (ozone, sulfur dioxide, ammonia, and fine and coarse particulate matter). The networks include the following: (1) Aerometric Information Retrieval System‐Air Quality Subsystem (AIRS‐AQS; https://www.epa.gov/aqs), (2) Clean Air Status and Trends Network (CASTNET, http://www.epa.gov/castnet/), (3) the Interagency Monitoring of Protected Visual Environments (IMPROVE; http://vista.cira.colostate.edu/improve/), (4) the Ammonia Monitoring Network (AMoN), and (5) the Ameriflux network of eddy covariance flux towers (discussed in section 3.3). The evaluation protocol uses average spatial analyses, time series, and common statistical evaluations that include root‐mean‐square error (RMSE), mean bias (MB), normalized mean bias (NMB), mean error (ME), normalized mean error (NME), and Pearson's correlation coefficient (R). Figure S2 provides the site locations for each network used for the different analyses that are area averaged in section 3. The evaluations against the network sites are all based on a distance‐weighting matching between the closest model grid cell for each observational site available.

Updates to WRF/Noah‐CMAQ

Updates to the Soil and Vegetation Parameters

WRF/Noah uses the soil hydraulic (SOILPARM) and vegetation (VEGPARM) tables to define parameters that are functions of the soil texture and LU category, respectively, and are influential drivers of the land surface characteristics (e.g., R and R ) and fluxes, and the prognostic WRF meteorological and CMAQ chemical fields. In this work, we specifically update the National Resources Conservation Service (NRCS) State Soil Geographic (STATSGO, i.e., “STAS”; https://www.nrcs.usda.gov) soil hydraulic parameters in the SOILPARM table, and NASA's MODIS‐based (https://modis.gsfc.nasa.gov) LU characteristics in the VEGPARM table. Figure S3 shows spatial plots of the dominant 17 (out of 19) STAS soil and 20 International Geosphere‐Biosphere Programme (IGBP)‐MODIS (Friedl et al., 2002, 2010) LU categories found in the model domain, which are updated for their characteristics in the SOILPARM and VEGPARM tables, respectively.

SOILPARM Table

The soil hydraulic parameters (i.e., SOILPARM table) help describe the water‐holding characteristics of soil and have been used in mostly agricultural areas (Saxton et al., 1986). These values affect water and energy fluxes in WRF/Noah and, consequently, the meteorological and air quality simulations with WRF/Noah‐CMAQ. Although increasing amounts of soil data have become available in the National Cooperative Soil Survey, National Cooperative Soil Characterization Database (USDA‐NRCS), the SOILPARM table used in WRF/Noah has not been updated for more than 25 years (Kishné et al., 2017). The default SOILPARM table is updated here following Kishné et al., which is based on 6,749 USDA‐NRCS soil samples of percent clay, silt, and sand; bulk density; and soil moisture content (SMC) at oven dry, air dry, and 33 kPa matric potential within a study area covering Texas and adjacent areas. For other parameters not measured, the values in Kishné et al. are based on pedotransfer functions (PTFs) from Cosby et al. (1984) and Chen and Dudhia (2001). Kishné et al. updated eight of the ten soil parameters in the SOILPARM table, and only table parameters “F11” (soil thermal diffusivity/conductivity coefficient) and “QTZ” (soil quartz content) were unchanged due to lack of supporting documentation and because of missing soil data, respectively (Table 3). In this study, the SOILPARM parameters for sand, loam, silt, and clay soil texture categories (Table 3; Categories 1–12) are used from Kishné et al. (2017), and in addition, we also updated SOILPARM values for soils that were excluded from Kishné et al., such as organic material, bedrock, land ice, playa, lava, and white sand horizons (Table 3; Categories 13–19 but excluding water). The parameters updated (i.e., filled‐in) for these soil types include the SMC for ceased evaporation for the top soil layer (DRYSMC), field capacity (REFSMC), and the wilting point (WILTSMC) using similar PTFs (Campbell & Norman, 1998; Cosby et al., 1984). As soil hydraulic parameters are important to processes simulated in CMAQ, such as WB dust and bidirectional ammonia (NH3) exchange (BIDI‐NH3), these values are also updated within CMAQ to be consistent with those in Table 3. While the Texas measurements used to update the hydraulic properties in the SOILPARM table in Kishné et al. cover a large range of U.S. soil types (Figure S3a), their extrapolation to other regions remains an approximation, which is compounded by the inability to represent the subgrid spatial variability of percent soil textural fractions.

Table 3

Modified WRF SOILPARM Table Comparison

Category	Type	BB	DRYSMC	MAXSMC	REFSMC	SATPSI	SATDK	SATDW	WLTSMC
	New, Default	Model calibration parameter	Dry soil moisture threshold	Porosity	Field Capacity	Sat. soil matric potential	Sat. soil conductivity	Sat. soil diffusivity	Wilting point soil moisture
1	Sand	3.36 , 2.79	0.004 , 0.01	0.402 , 0.339	0.086 , 0.236	0.047 , 0.069	2.40E−5 , 4.66E−5	9.36E − 6 , 6.08E−7	0.024 , 0.01
2	Loamy sand	4.06 , 4.26	0.01 , 0.028	0.396 , 0.421	0.142 , 0.383	0.063 , 0.036	1.72E − 5 , 1.41E−5	1.11E − 5 , 5.14E−6	0.057 , 0.028
3	Sandy loam	4.85 , 4.74	0.016 , 0.047	0.413 , 0.434	0.213 , 0.383	0.105 , 0.141	1.01E − 5 , 5.23E−6	1.18E − 5 , 8.05E−6	0.081 , 0.047
4	Silt loam	5.72 , 5.33	0.023 , 0.084	0.456 , 0.476	0.303 , 0.36	0.399 , 0.759	2.20E − 6 , 2.81E−6	9.99E − 6 , 2.39E−5	0.123 , 0.084
5	Silt	4.18 , 5.33	0.01 , 0.084	0.438 , 0.476	0.346 , 0.383	0.568 , 0.759	1.30E − 6 , 2.81E−6	7.78E − 6 , 2.39E−5	0.064 , 0.084
6	Loam	6.01 , 5.25	0.022 , 0.066	0.44 , 0.439	0.274 , 0.329	0.219 , 0.355	4.12E − 6 , 3.38E−6	1.20E − 5 , 1.43E−5	0.128 , 0.066
7	Sandy clay loam	7.03 , 6.77	0.029 , 0.067	0.416 , 0.404	0.288 , 0.314	0.139 , 0.135	7.08E − 6 , 4.45E−6	1.58E − 5 , 9.90E−6	0.168 , 0.067
8	Silty clay loam	8.49 , 8.72	0.039 , 0.12	0.457 , 0.464	0.35 , 0.387	0.517 , 0.617	1.48E − 6 , 2.03E−6	1.41E − 5 , 2.37E−5	0.212 , 0.12
9	Clay loam	8.2 , 8.17	0.036 , 0.103	0.449 , 0.465	0.335 , 0.382	0.288 , 0.263	3.02E − 6 , 2.45E−6	1.52E − 5 , 1.13E−5	0.196 , 0.103
10	Sandy clay	8.98 , 10.73	0.037 , 0.1	0.425 , 0.406	0.355 , 0.338	0.171 , 0.098	5.13E − 6 , 7.22E−6	1.93E − 5 , 1.87E−5	0.239 , 0.1
11	Silty Clay	10.24 , 10.39	0.052 , 0.126	0.467 , 0.468	0.392 , 0.404	0.588 , 0.324	1.27E − 6 , 1.34E−6	1.61E − 5 , 9.64E−6	0.264 , 0.126
12	Clay	11.56 , 11.55	0.058 , 0.138	0.506 , 0.468	0.428 , 0.412	0.483 , 0.468	1.66E − 6 , 9.74E−7	1.74E − 5 , 1.12E−5	0.285 , 0.138
13	Organic material	5.25, 5.25	0.003 , 0.066	0.439, 0.439	0.286 , 0.329	0.355, 0.355	3.38E−6, 3.38E−6	1.43E−5, 1.43E−5	0.118 , 0.066
14	Water	0, 0	0, 0	1, 1	0, 0	0, 0	0, 0	0, 0	0, 0
15	Bedrock	2.79, 2.79	0.001 , 0.006	0.2, 0.2	0.05 , 0.17	0.069, 0.069	1.41E−4, 1.41E−4	1.36E−4, 1.36E−4	0.009 , 0.006
16	Other (land‐ice)	4.26, 4.26	0.01 , 0.028	0.421, 0.421	0.145 , 0.283	0.036, 0.036	1.41E−5, 1.41E−5	5.14E−6, 5.14E−6	0.049 , 0.028
17	Playa	11.55, 11.55	0.147 , 0.03	0.468, 0.468	0.395 , 0.454	0.468, 0.468	9.74E−7, 9.74E−7	1.12E−5, 1.12E−5	0.264 , 0.03
18	Lava	2.79, 2.79	0.001 , 0.006	0.2, 0.2	0.05 , 0.17	0.069, 0.069	1.41E−4, 1.41E−4	1.36E−4, 1.36E−4	0.009 , 0.006
19	White sand	2.79, 2.79	0.001 , 0.01	0.339, 0.339	0.084 , 0.236	0.069, 0.069	4.66E−5, 4.66E−5	6.08E−7, 6.08E−7	0.015 , 0.01

Note. Table includes both unmodified (italic) and modified (bold) WRF/Noah soil hydraulic parameters that are part of Weather Research and Forecasting's SOILPARM table. Increases in the modified table are colored red while the decreases are colored blue, and no changes are indicated in black. Modified values for sand and clay soil texture categories (1–12) are adapted from Kishné et al. (2017). All other categories are updated based on descriptions provided in this work. SOILPARM = soil hydraulic parameter table.

An important aspect of the modified SOILPARM tables is that decreased DRYSMC for most categories (Table 3) reflects that soils could dry via evaporation to a lower threshold, whereas the DRYSMC is identical to WILTSMC in the unmodified tables. Therefore, the previous values found in the table are unrealistic as soil moisture should reach the wilting point before further decreasing to where evaporation from soil ceases.

The changes to the saturated air entry water potential (SATPSI) alter the partitioning of precipitation between runoff and infiltration (e.g., reduced values for silt loam and silt), whereas saturated hydraulic conductivity (SATDK) and diffusivity (SATDW) differences affect the water filling and pressure movements within different soils, respectively. Here modified values have larger movement in sandy soils, for example. The effect of reducing potential (REFSMC–WILTSMC) and total available water (MAXSMC–WILTSMC) by 35–76% (Kishné et al., 2017) in the modified table is exacerbated by increasing WILTSMC in all soil categories except for silt, so that the modified tables reduce the potential soil water available for plants to transpire in all soil texture categories (except silt) compared to the unmodified tables (Figure 5b in Kishné et al., 2017). Overall, drier soils and reduced plant transpiration (from less potential and total available water) from most of the modified SOILPARM categories will alter the WRF/Noah surface moisture fluxes.

Average spatial analysis of the impact assessment of these changes in the modified WRF/Noah‐CMAQ (Table 2; run A4) shows that the reduced MAXSMC, REFSMC, and DRYSMC for many soil categories (Table 3) reduces the threshold for soil evaporation and results in many areas of decreased soil water content (soil_H2O %), i.e., drier soils, and increased latent heat flux (LH_flx) and diagnosed 2‐m water vapor mixing ratio (Q2) (Figure 1). Other regions of increased soil_H2O % (e.g., northern Nebraska, northern Michigan, Florida), are dominated by a very fine “sand” classification (Figure S3a), which has a larger MAXSMC in the updated SOILPARM table. The enhanced soil evaporation affects the partitioning between sensible heat flux (SH_flx) and LH_flx (i.e., decreased Bowen Ratio; B = SH_flx/LH_flx), and thus, there are also widespread decreases in topmost soil layer temperature (TSLB), SH_flx, and diagnosed 2‐m temperature (T2) due to the modified SOILPARM table.

We note that there are small impacts of the SOILPARM changes on the CMAQ deposition and chemical variables analyzed in the short‐term, 10‐day assessment run, and thus, their results are not shown. Table S1 provides a summary of the model performance impacts for the updated SOILPARM table, as well as for the unmodified off‐line WRF/Noah‐CMAQ and each additional modification made to the coupled WRF/Noah‐CMAQ system. Incorporation of the modified SOILPARM table both slightly improves (e.g., SH_flx) and degrades (e.g., LH_flx) model performance (e.g., NMB, RMSE, and R) for the meteorological variables impacted (Figure 1). Finally, as with all other analyses that follow, there are small impacts on precipitation for the impact assessment runs, and results are not included here.

VEGPARM Table

Vegetation data (i.e., VEGPARM table) is used for retrospective WRF simulations, and is useful to capture changes in vegetation parameters when simulations include LU change or future projections; however, these data can be better informed from satellite observations (Bousetta et al., 2015). To assess the existing tabular data in the unmodified Noah VEGPARM table, Noah LAI estimates were evaluated against MODIS monthly mean LAI data distributed with WRFv3.8.1 and then the VEGPARM table was repopulated based on MODIS satellite observations. MODIS LAI estimates have shown biases, particularly overestimates in sparsely vegetated areas and during leaf off conditions, but MODIS LAI also contains the most complete LAI coverage (Garrigues et al., 2008) and is a data set distributed and supported by the WRF modeling system. LAI is estimated within Noah based on equation (1) but constrained using minimum and maximum values found in the VEGPARM table: Where F veg is the green vegetation coverage fraction distributed with the WRF Preprocessing System (AVHRR estimates in this case), F LU is the fraction of the grid cell that is composed of LU category i, n is the number of LU categories in the grid cell, and LAI min and LAI max are the minimum and maximum LAIs associated with each LU category i and listed in VEGPARM. This assumes that the maximum LAI coincides with the maximum vegetation fraction. When using the unmodified WRF VEGPARM table associated with MODIS LU data in the WRFv3.8.1/Noah model, the calculated LAI is nearly twice as high as the satellite LAI observations throughout most of the study area (Figure S4).

Here this overestimation is rectified, first, by calculating the satellite LAI at the minimum and maximum vegetation fractions, and respectively, from the multiyear monthly LAI and vegetation coverage satellite data distributed with WRF. Subsequently, we estimate the appropriate maximum and minimum LAI values using multiple least squares regressions. The LAI at the maximum and minimum vegetation coverage can be expressed in equations (2) and (3) (following linear regression), respectively:

Here and can be estimated from the satellite‐derived LAI and F veg, and then the VEGPARM LAI min and LAI max tabular data for each LU category can be estimated using least squares linear regression model similar to Bousetta et al. (2015) where LAI min and LAI max are the regression coefficients. The model was evaluated using a tenfold cross validation where the LAI arrays were randomly sampled to create 11 unique arrays, the linear regression model was run for 10 of those arrays, and the VEGPARM LAI min and LAI max data were the averages of the iterations where the regression coefficient was significant, using a p < 0.05 criterion. Except for a few LU categories, all 10 of the regressions' coefficients were typically significant at p < 0.001 (Tables S2 and S3). These revised VEGPARM values reduced the maximum and minimum calculated VEGPARM LAI differences when compared to MODIS values from 108.9% to 9.3% and 351.8% to 0.0%, respectively, compared to the observations (n = 13,316) not used in the regressions. As a result, the January and July LAI spatially averaged biases when compared to satellite values were reduced from 382.2% to 25.0% and from 92.6% to 2.4% and the normalized mean error was reduced from 100.0% to 20.6% and from 381.6% to 45.6%, respectively (Figure S4). Note that this regression was done for a CONUS domain and the MODIS LU data.

The maximum aerodynamic roughness length (z ) in the updated VEGPARM table, , was estimated from 1‐km resolution remotely sensed lidar canopy heights from the Geoscience Laser Altimeter System (GLAS) from 2005 as documented in Simard et al. (2011). The canopy heights were logarithmically averaged (equal to the difference divided by the logarithm of their quotient) over the 12‐km grid cell where the fractional forest land use types were greater than zero in this application. The average height determined using the same methodology as the LAI was used to estimate the canopy height (z topv) in the VEGPARM table. was assumed to be 7% of the z topv in the range of observations reported by Nakai et al. (2008). The spatially explicit data was then mapped to the VEGPARM table using a tenfold cross validation similar to the LAI updates. The minimum z () was then estimated from using the ratio of to from the original VEGPARM table, as leaf off conditions could not be assessed from the lidar data.

Table 4 shows a comparison of the modified VEGPARM table's LAI min and LAI max, and , and z topv used in this work compared to the unmodified table (all other VEGPARM parameters unchanged). Results would vary depending on the LU data and model domain selected. The reduced LAI in the VEGPARM table (Table 4 and Figure S4) results in a reduction in canopy transpiration and LH_flx, increased SH_flx, and shifts in B that perturb the surface energy balance in the Noah LSM (Figure 2). Consequently, there are less clouds and more downward shortwave radiation at the surface, especially in the central CONUS regions, and the surface energy balance in the Noah LSM requires that there be some compensating increases in the TSLB, SH_flx (i.e., increased B), and diagnosed T2 (neglecting wind speed changes; Figures 2d–2f). In contrast to SOILPARM updates, there are minimal short‐term impacts of the VEGPARM updates on soil_H2O % changes (not shown). The increased T2 due to the VEGPARM updates also leads to widespread deepening of the planetary boundary layer height (PBLH; not shown). There are widespread decreases in 10‐m wind speed (WSPD10; Figure 2a), which are a consequence of increased z associated with some of the LU categories in the updated VEGPARM table (Table 4). The decrease in LAI changes the partitioning of total evaporation between bare soil and canopy transpiration that dominates the LH_flx decrease (Figure 2b), while the decrease in wind speed impacts the coupling between the skin layer and first model layer (i.e., the surface exchange coefficient) that impacts the SH_flx increase (Chen & Dudhia, 2001). The result of these changes is widespread T2 increases (warmer) and Q2 decreases (drier; Figures 2c and 2f).

Table 4

Modified WRF VEGPARM Table Comparison

Category	Type	LAI min	LAI max	zomin	zomax	z topv
	New, Default	Minimum LAI	Maximum LAI	Minimum aerodynamic roughness length	Maximum aerodynamic roughness length	Canopy top height
1	Evergreen Needleleaf Forest	1.00 , 5.00	3.00 , 6.40	1.41 , 0.50	1.41 , 0.50	20.1, 17.0
2	Evergreen Broadleaf Forest	3.85 , 3.08	6.24 , 6.48	1.03 , 0.50	1.03 , 0.50	14.7, 35.0
3	Deciduous Needleleaf Forest	0.01 , 1.00	2.00 , 5.16	0.50 , 0.50	0. 50, 0.50	5.0 , 14.0
4	Deciduous Broadleaf Forest	0.60 , 1.85	5.68 , 3.31	0.53 , 0.50	1.89 , 0.50	27.0 , 20.0
5	Mixed Forests	0.79 , 2.80	3.90 , 5.50	0.87 , 0.20	1.60 , 0.50	22.9 , 18.0
6	Closed Shrublands	0.10 , 0.50	1.00 , 3.66	0. 01, 0.01	0.19 , 0.05	1.90 , 0.50
7	Open Shrublands	0.14 , 0.60	0.50 , 2.60	0. 01, 0.01	0.18 , 0.06	1.80 , 0.50
8	Woody Savannas	1.44 , 0.50	3.39 , 3.66	0. 01, 0.01	0. 05, 0.05	0.50, 0.50
9	Savannas	0.88 , 0.50	2.69 , 3.66	0.05 , 0.15	0.05 , 0.15	0.50, 0.50
10	Grasslands	0.17 , 0.52	0.62 , 2.90	0.06 , 0.10	0.06 , 0.12	0.57 , 0.50
11	Permanent Wetlands	0.10 , 1.75	1.50 , 5.72	0.01 , 0.30	0.01 , 0.30	0.10 , 0.00
12	Croplands	0.16 , 1.56	2.00 , 5.68	0.16 , 0.05	0.16 , 0.15	0.50, 0.50
13	Urban and Built‐up	0.45 , 1.00	1.75 , 1.00	1.34 , 0.50	1.34 , 0.50	0.00, 0.00
14	Cropland/Natural Vegetation Mosaic	0.38 , 2.29	2.26 , 4.29	0.15 , 0.05	0.15 , 0.14	1.48 , 0.50
15	Snow and Ice	0.04 , 0.01	0.44 , 0.01	0.001, 0.001	0.001, 0.001	0.00, 0.00
16	Barren or Sparsely Vegetated	0.01 , 0.10	0.01 , 0.75	0.16 , 0.01	0.16 , 0.01	0.02, 0.02
17	Water	0.01, 0.01	0.01, 0.01	0.0001, 0.0001	0.0001, 0.0001	0.00, 0.00
18	Wooded Tundra	0.02 , 0.41	0.71 , 3.35	0.06 , 0.30	0.06 , 0.30	0.64 , 10.0
19	Mixed Tundra	0.05 , 0.41	1.56 , 3.35	0.42 , 0.15	0.42 , 0.15	4.24 , 5.00
20	Barren Tundra	0.41, 0.41	3.35, 3.35	0. 05, 0.05	0. 10, 0.10	0.02, 0.02

Note. Table includes comparison of the default (italic) and updated (bold) Weather Research and Forecasting/Noah VEGPARM LAI min, LAI max, , , and z topv for this model domain configuration and IGBP MODIS‐20 Noah LU data set (Table 1). Increases in the updated table are in red text, while the decreases are in blue text. VEGPARM = vegetation parameter table; IGBP = International Geosphere‐Biosphere Programme; MODIS = Moderate Resolution Imaging Spectroradiometer.

Figure 2

Average (21–31 May 2011) incremental changes due to the modified VEGPARM table in WRF/Noah‐CMAQ (runs A5 and A4) for (a) WSPD10, (b) LH_flx, (c) Q2, (d) TSLB, (e) SH_flx, and (f) T2. VEGPARM = vegetation parameter table; WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; WSPD10 = 10‐m wind speed; LH_flx = latent heat flux; Q2= 2‐m mixing ratio; TSLB = topmost soil layer temperature; SH_flx = sensible heat flux; T2 = 2‐m temperature.

Incorporating more realistic satellite‐based vegetation characteristics in the updated VEGPARM table shows either consistency or improvement in NMB, RMSE, and R statistical performances of T2, LH_flx, and WSPD10 (Table S1). The most notable improvement of these variables is WSPD10, where WRF overpredictions are reduced due to increased z in most categories for the updated VEGPARM table (Table 4), and consequently there is a decrease in NMB of 18%, decrease in RMSE of 0.3 m/s, and an increase in R from 0.61 to 0.67 (Table S1). The largest reductions in bias and error for T2 and WSPD10 occur in the observation ranges >30 °C and <8 m/s, respectively (observational binned analyses not shown). There are other notable meteorological performance improvements for LH_flx due to reduced LAI in the updated VEGPARM table, which impact the heat flux partitioning (i.e., decreased LH_flx and increased SH_flx), and consequently lead to a decrease in NMB and reduced RMSE (Table S1). This change, however, results in an increase in NMB and increased RMSE for SH_flx. Thus, the introduction of vegetation more closely related to MODIS satellite observations has a significant effect on the thermodynamic state (i.e., tuning) of the WRF/Noah model up to this point, while a more exhaustive investigation (outside the scope here) on other model factors contributing to the LH_flx and SH_flx dichotomy is recommended. We note that a recent set of four papers based on the Clouds Above the United States and Errors at the Surface (CAUSES) project identified that the simulation of the evaporative fraction (EF = LH_flx/[SH_flx+LH_flx]) at the surface and simulation of deep convective clouds are main contributors to numerical weather prediction models' systematic warm bias in the central CONUS (see Steiner, 2018, and references contained within).

Gases such as O3, SO2, and nitric acid (HNO3) are governed by unidirectional deposition, and they are strongly dependent on surface characteristics. Thus, there are also impacts of the VEGPARM table updates on modified WRF/Noah‐CMAQ dry deposition (DDEP) and mixing ratios of O3, SO2, and HNO3 (Figure 3). The reduced LAI in the updated VEGPARM table acts to decrease the DDEP_O3 (domain relative change ~ −27%) and DDEP_SO2 in the central and western CONUS, with some increases in the eastern CONUS (Figures 3a and 3b). The decrease in DDEP_SO2 is strongly mitigated by the CMC × 103 correction (Figure 3b), which increases DDEP of soluble gases (especially SO2), as it was applied to all modified WRF/Noah‐CMAQ simulations (Table 2). As DDEP_HNO3 is limited by R , the reduced low‐level wind speeds from the VEGPARM table updates also reduces the R in the modified coupled WRF/Noah‐CMAQ model used in this simulation (see sections 2.2. and 3.2), and thus results in widespread increases in DDEP_HNO3 (Figure 3c).

Figure 3

Average (21–31 May 2011) incremental changes due to the modified VEGPARM table in WRF/Noah‐CMAQ (runs A5 and A4) for (a) DDEP_O3, (b) DDEP_SO2, (c) DDEP_HNO3, (d) DDEP_NH3, (e) O3, (f) SO2, (g) HNO3, and (h) NH3. VEGPARM = vegetation table; WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; DDEP = dry deposition.

These modifications induce changes in O3 and SO2 mixing ratios with a range of about −5 to +14 ppb and −1 to +3 ppb, respectively (Figures 3e and 3f). There are also offsetting effects of DDEP_O3 increases and O3 mixing ratio decreases due to the updated Mészáros et al. (2009) parameterization for dry deposition of O3 to soil in the relative dry regions of CONUS (e.g., southwest CONUS). The changes in HNO3 mixing ratio are not strongly tied to the unidirectional DDEP_HNO3 changes; rather HNO3 is driven by changes in its precursor gases (i.e., nitrogen dioxide, NO2, and hydroxyl radical, OH), which are most strongly impacted by the indirect effects of the updated VEGPARM table on such gas concentrations (not shown). The resulting change in HNO3 mixing ratio has a domain‐wide range of about −0.2 to +0.6 ppb (Figure 3g). The variation of HNO3 may also be related to NH3 concentration in regions (e.g., parts of central and western United States, including California) that are characterized by high NH3 and particulate nitrate levels in winter/early spring (see Figure 5 in Campbell et al. (2015)), as typically decreases (increases) in NH3 will result in decreases (increases) particulate nitrate and increases (decreases) in HNO3.

The BIDI‐NH3 exchange option is used in the modified WRF/Noah‐CMAQ (Table 1), and thus the NH3 flux is governed by bidirectional exchange and atmosphere‐canopy compensation points, which are a function of R and R (and LAI) that are recalculated (defined) in CMAQ (Bash et al., 2010, 2013; Cooter et al., 2012). When the compensation point is larger than the ambient concentration the flux is evasion and when it is lower the flux is deposition. The compensation point is generally most sensitive to the soil emission potential in agricultural LU and vegetation emission potentials in nonagricultural LU. Thus, the changes in DDEP_NH3 due to the modifications in this work are somewhat different than changes in DDEP_O3 and DDEP_SO2. In regions of significant DDEP_NH3 (and enhanced NH3 concentrations) during the spring season there are widespread increases in DDEP_NH3 due to the CMC × 103 correction. Dependent on the CONUS region considered, however, the VEGPARM table updates can also either increase (as in the southwest and northeast CONUS) or decrease (as in the central CONUS) the DDEP_NH3 due to the indirect impacts of VEGPARM‐induced meteorological/soil parameter changes on BIDI‐NH3 exchange (Figure 3d).

There are widespread decreases in NH3 due to the VEGPARM table updates (Figure 3h), which are due to the parameterization of the NH3 flux as a gradient‐dependent process between atmospheric NH3 and estimated NH3 in soil and vegetation media. The changes to WRF primarily alter the rate at which NH3 is transferred across these gradients. For example, in an area with high soil NH4 + content and NH3 compensation points, a reduction in the rate of transfer across this gradient would reduce the evasive flux and consequential ambient NH3 mixing ratio. This is what happened with the VEGPARM table update and the reduction in the ambient NH3 concentrations in the Midwest.

WRF/Noah Stomatal and Aerodynamic Resistances

The unmodified WRF/Noah does not make available either the bulk R or R that could be used as conductance (g  = 1/R or g  = 1/R ) to compute the dry deposition velocity (V ) for gas species in CMAQ using the electrical resistance analog model “M3Dry” (Pleim et al., 2001). In the unmodified off‐line WRF/Noah‐CMAQ, R or R are instead parameterized in MCIP using static LU type look‐up tables based on the PX LSM and M3Dry, which are inconsistent with the set of vegetation parameters used to calculate WRF/Noah fluxes. These parameterizations in MCIP therefore contain additional uncertainty, especially for chemical species that have a significant stomatal deposition pathway (e.g., O3 and SO2). Here WRF/Noah is modified so that its time‐varying R and R are made available from WRF and then passed directly to CMAQ to ensure consistent treatment of meteorological (heat and moisture; evapotranspiration) and chemical fluxes as in the PX LSM (Pleim & Ran, 2011), while simplifying the CMAQ dry deposition calculations for WRF/Noah. The formulation of R in WRF/Noah is based on the empirical Jarvis model (Chen et al., 1996; Jacquemin & Noilhan, 1990; Jarvis, 1976; Noilhan & Planton, 1989), which is similarly parameterized in the off‐line MCIP framework discussed above.

In a spatial comparison of the unmodified off‐line and modified coupled WRF/Noah‐CMAQ R (with a set of runs each using dominant and mosaic/tiled LU; runs A2 and A3 in Table 2), there are notable spatial differences especially in the west‐northwest and east‐southeast parts of CONUS, and in parts of Canada (Figures 4a and 4b). The differences in R stem from the fact that the unmodified (off‐line) WRF/Noah‐CMAQ calculations of R in MCIP are typically performed using hourly WRF output without considering seasonality, as the minimum R (), which is a LU‐dependent parameter in the Jarvis model, is only based on the assumed maximum LAI for each LU category (Figure 4a). Furthermore, the effects of soil moisture on the Jarvis R calculation in the unmodified system are also based on a LU category look‐up table that uses a single loam soil texture, which excludes the effects of precipitation on soil moisture. In the modified coupled WRF/Noah‐CMAQ (Figure 4b), however, R is calculated at every 1‐min WRF integration time step when both the vegetation fraction and LAI are greater than zero and depends on prognostic WRF/Noah variables such as solar radiation, air temperature and atmospheric water vapor deficit at the lowest model level, and soil moisture fields that include effects of precipitation that are pivotal to the Noah LSM water balance. Both the unmodified and modified WRF/Noah‐CMAQ use initial table values for from LU‐based tables, but differences arise from the modulation of by the different environmental stress functions in the Jarvis model at different update frequencies, that is, 60‐min WRF output in the unmodified off‐line system versus 1‐min time step in the modified coupled WRF/Noah‐CMAQ system.

Figure 4

The 26–31 May 2011 average spatial comparisons of R for the (a) unmodified WRF/Noah‐CMAQ (run A1), (b) difference in modified‐unmodified WRF/Noah‐CMAQ (runs A2 and A1), and (c) difference in modified‐unmodified WRF/Noah‐CMAQ with mosaic (runs A3–A1). Panel (d) shows comparisons of the domain‐wide averaged diurnal pattern of g for the unmodified (run A1: black), modified (run A2: red), and modified with mosaic WRF/Noah‐CMAQ (run A3: green). WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality.

The overall enhanced spatial variability (and lower values) for R in the relatively dry grasslands, open shrublands, and savannas within the western CONUS (e.g., eastern Montana, eastern Wyoming, and Washington; see Figure S3 for LU classification) are a result of lower for these LU categories in the modified WRF/Noah‐CMAQ that are taken from the VEGPARM table, compared to the taken from the tables used in the off‐line MCIP parametrization (Table S4). Increasing soil moisture due to precipitation during this period (not accounted for in the unmodified WRF/Noah‐CMAQ) enhances the differences in R in the relatively drier western CONUS region. For example, there are impacts from precipitation that occurred between 26 and 31 May 2011 just west of the Cascade Range in Washington and western Oregon (not shown), where the R values are considerably lower for the modified WRF/Noah‐CMAQ. The larger modified R values that are found in the south and southeastern CONUS (e.g., Florida), are dominated by a mixture of sand and sandy loam in the top soil layer (Figure S3a), and were typical of lower soil moisture during this period (Figure S1e) compared to areas in northern CONUS, which are dominated by silt and loam in the top soil layer (Figure S3a) and are typical of higher soil moisture (Figure S1e). These differences increase R in the southern and southeastern CONUS, resulting in more local variability in R for the modified WRF/Noah‐CMAQ.

The modified R for Noah mosaic LU representation (Li et al., 2013), which calculates the grid cell R (also following the Jarvis model) based on a weighted average of the fractional LU coverages, features additional spatial heterogeneity, as well as more urban LU representation within grid cells that had only rural LU types with the dominant treatment (Mallard & Spero, 2018). In this case, the mosaic leads to more grid cell representation of LU categories (e.g., croplands) that have lower  compared to that for forests (Table S4), which are prevalent in some Mid‐Atlantic states (e.g., West Virginia) and parts of eastern CONUS (<95°W) under a dominant LU representation (Figures S3b).

There is a more pronounced diurnal pattern and peak g for the modified WRF/Noah‐CMAQ (both dominant and mosaic LU) following sunrise, which is strongly impacted by the lower for many of the widespread LU categories across CONUS in the modified system (Table S4), and because the dynamic photosynthetic active radiation (PAR) is updated every model time step in the Jarvis calculation of R , compared to an hourly basis in the off‐line MCIP parameterization of R

The unmodified off‐line WRF/Noah‐CMAQ model also parameterizes R using hourly WRF output in MCIP following Monin‐Obukhov (MO) similarity theory (with atmospheric stability corrections) prior to its use in CMAQ. R is not included in the unmodified coupled WRF/Noah‐CMAQ model. Thus, to improve its physical consistency between WRF/Noah and CMAQ, R is now directly calculated in each model time step in the modified coupled WRF/Noah‐CMAQ system following Garland (1977) as follows: where WSPD is the magnitude of the horizontal wind speed in the lowest model layer; u and v are the eastward and northward components of the wind, respectively; w* is the convective velocity scale (i.e., stability) correction (over land: Beljaars, 1995; over water: Wyngaard, 1988); V is the subgrid‐scale velocity correction (Mahrt & Sun, 1995); and u the friction velocity. An additional stability correction term on the right‐hand side of equation (4) in the original form of R in Garland (1977; not shown) is zero in neutral and stable conditions (more details found in section 2.3 of Nemitz et al., 2009). In unstable conditions when the stability correction term becomes larger than zero, its impact is captured in equation (4) by using the WRF predicted w* in WSPD. Thus, the R and aerodynamic conductance (g  = 1/R ) in modified WRF/Noah‐CMAQ (Figure 5) is computed using a simple closed‐form parameterization for R that closely resembles MO similarity theory, is supported by the surface layer schemes in WRF, and has minimal numerical noise.

Figure 5

Same as in Figure 4 but for aerodynamic resistance (R ) and conductance (g ).

Because both the unmodified off‐line MCIP and modified coupled WRF/Noah‐CMAQ (equation (4)) use MO similarity theory, R is comparable in both simulations (Figures 5a–5c), with smaller and more systematic differences than R . The areas of widespread larger R values in the west and in parts of the east CONUS, are driven by slight differences in its R formulation in equation (4), including the addition of the w* and V sg corrections to WSPD, which increase WSPD and R . There is also an impact from higher or lower 1‐min instantaneous (i.e., model time step) WSPD values in equation (4) that are calculated by the modified coupled WRF/Noah‐CMAQ model, compared to hourly WSPD values that are used in the unmodified off‐line WRF/Noah‐CMAQ. This modification allows for a stronger diurnal effect of transitioning between relatively laminar (higher WSPD) and turbulent flow (lower WSPD), thus leading to both spatially and temporally localized increases and decreases in R that are somewhat hidden in the domain‐wide average, especially when analyzing the daily trends in the east and west CONUS, where there is evidence of day‐to‐day and spatial variability (Figure S5). There is a consistently lower domain‐averaged aerodynamic conductance (g  = 1/R ) diurnal profile for the modified and modified w/mosaic compared to the default WRF/Noah‐CMAQ (Figure 5d), and only very slight differences in the spatially averaged dominant and mosaic LU diurnal profiles of g .

Using modified coupled WRF/Noah‐CMAQ R and R (dominant LU) leads to changes in atmospheric deposition of gases compared to the unmodified off‐line MCIP‐CMAQ system (Figure 6). The decreases in R over much of the CONUS lead to increases in DDEP_O3 and DDEP_SO2, especially over eastern CONUS (Figures 6a and 6b; <95°W). These increases are also impacted by the CMC × 103 correction. Implementation of the mosaic R updates have only slight impacts on further increasing the DDEP_O3 and DDEP_SO2 in eastern CONUS (not shown). As noted for the VEGPARM updates (section 3.1.2), DDEP_HNO3 is not dependent on surface resistance (i.e., chemical, physical, or biological characteristics), and the moderate (domain relative changes of <5%) DDEP_HNO3 decreases (Figure 6c) are instead determined by the changes in R (Figure 5b).

Figure 6

Average (21–31 May 2011) incremental changes due to the modified R (dominant land use) and R in WRF/Noah‐CMAQ (runs A2 and A1) for (a) DDEP_O3, (b) DDEP_SO2, (c) DDEP_HNO3, (d) DDEP_NH3, (e) O3, (f) SO2, (g) HNO3, and (h) NH3. WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; DDEP = dry deposition.

The spatial distribution of DDEP_O3 and DDEP_SO2 are qualitatively similar to combined changes in their respective O3 and SO2 mixing ratios (Figures 6e and 6f). The spatial changes in HNO3 mixing ratio (Figure 6g), however, are again not well correlated with the widespread increases in DDEP_HNO3. Rather they are governed by widespread increases in NO2 due to reduced DDEP_NO2 in central CONUS (not shown), decreases in OH due to significant decreases in O3 in eastern CONUS (<95°W; see Figure 6e), and impacts of changes in NH3 and particle nitrate formation on resulting HNO3 concentrations (see discussion of Figures 3d and 3e in comparison).

There are widespread increases (and localized decreases) in both DDEP_NH3 (Figure 6d) and NH3 (Figure 6h) that are not due to the R and R updates, as these variables are re‐calculated in the BIDI‐NH3 model in CMAQ. We note that the R recalculation for BIDI‐NH3 model results in similar values as the Noah estimate, and the R recalculation is also similar, but with higher  for agricultural crops. These differences are thus a result of an interplay of the CMC × 103 correction over areas with considerable agriculture (i.e., croplands; Figure S3b) and NH3 evasion in this simulation. In these regions, there is an increase in the rate of transport across the soil‐vegetation‐air gradient, which results in an increase in DDEP_NH3 and NH3 in this region despite a decrease in the other modeled species.

Modified Soil Layer Depths

Soil moisture and temperature in WRF/Noah are vital inputs to CMAQ. The topmost soil layer depth in Noah (10 cm) is not shallow enough to be consistent with CMAQ processes for WB dust (Pleim & Xiu, 1995) and BIDI‐NH3 that maintains the soil ammonium mass balance in the top soil layer to represent surface fertilizer application (Bash et al., 2013; Cooter et al., 2012). When running the unmodified off‐line WRF/Noah‐CMAQ, the soil moisture and temperature are adjusted in CMAQ from 10‐ to 1‐cm depth based on Darmenova et al. (2009), which creates additional uncertainty in model predictions. In this work, we modified the WRF structure of the default 4 Noah soil layer depths from 10, 40, 100, and 200 cm to depths of 1, 10, 100, and 200 cm, so that the modified WRF/Noah can provide prognostic soil temperature and moisture from the top two soil layers (1 and 10 cm) directly to CMAQ. In this case, the soil layer depths are adjusted in the WRF preprocessor module that controls both the depths in the WRF input data, as well as their depths in WRF/Noah.

Comparing an impact assessment of the “thinner” topmost soil layer in the modified WRF/Noah‐CMAQ (Table 2; run A6) to observations from the topmost soil layer depth at 40 AmeriFlux sites (https://ameriflux.lbl.gov/; Baldocchi et al., 2001; Schaefer et al., 2012) shows that it leads to a larger range (i.e., stronger diurnal surface signal) in TSLB reaching warmer temperatures, and also lower soil_H2O% with drier soils compared to the unmodified 10‐cm depth. This leads to an overall worse WRF model performance for the thinner soil layer, with higher NMB and RSME, and lower R for TSLB, soil_H2O %, and SH_flx compared to Ameriflux (Figure S6). We note that the 40 Ameriflux sites have variable topmost soil layer depth measurements (~2–10 cm) at different sites (see Table 1 in Qiu et al., 2016 for a comparison across 36 Ameriflux sites), whereas the values from WRF/Noah are extracted at fixed depths. In an additional comparison over the summer 2011, however, there is good agreement in the AmeriFlux soil_H2O % and modified WRF/Noah results at 1 cm, with a bias very close to zero (not shown).

The warmer TSLB due to a thinner topmost soil layer depth (at 1 cm) also has an impact on BIDI‐NH3, which leads to increases in DDEP_NH3, and NH3, especially near agricultural regions on the central Mid‐West United States (not shown).

Using the thinner 1‐cm layer also allows for a unique opportunity to test an updated CMAQ WB dust module (Foroutan et al., 2017) in the modified WRF/Noah‐CMAQ model (Table 2: run A7), which predicts a topmost layer soil moisture at 1‐cm depth that is more appropriate for dust generation (compared to 10 cm for the default Noah LSM; see Darmenova et al., 2009, for discussion).

Figure 7 shows the unmodified 10 cm model (A1), modified 1 cm (A6) without WB dust, and modified 1‐cm model with WB dust (A7) runs during a 29 May 2011 dust event in southwest CONUS (see Foroutan et al., 2017) compared against 12 AQS coarse (PM10) and fine mode PM (PM2.5) sites within the state of New Mexico for hourly time series and biases, and against six IMPROVE PM2.5 sites for daily average stacked (i.e., percent composition) bar chart comparisons (see Figures S2d and S2e for site locations).

Figure 7

The 29 May 2011 time series and biases of unmodified (run A1; red), modified with no windblown dust (run A6; blue), and modified with windblown dust (1 cm) (run A7; green) hourly (a, b) PM10 and (c, d) PM2.5 against averaged observations over 12 Air Quality Subsystem sites (gray) within the state of New Mexico. Panel (e) compares 29 May 2011 24‐hr average PM2.5 composition (i.e., stacked bar plots) from runs A6 and A7 against an average of the six IMPROVE sites in New Mexico. IMPROVE = Interagency Monitoring of Protected Visual Environments; RMSE = root‐mean‐square error; WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality.

Incorporation of the new dust module (Foroutan et al., 2017) in the modified WRF/Noah‐CMAQ system leads to an overall good model response of enhanced PM10 and PM2.5 concentrations, and lower bias/error (compared to no dust in A1 and A6) during the 29 May 2011 dust event that led to average observed PM10 concentrations above >250 μg m3 (Figures 7a and S7b; ~17:00 LT), while incorporating dust leads to 92% (85%) and 92% (100%) of the 13 and 12 AQS sites having decreased MB (ME) for PM10 and PM2.5, respectively (not shown).

The increase in percent soil composition for the WB dust run also leads to a better agreement with the average of the six IMPROVE sites. While outside the scope in this paper, a more consistent comparison of dust performance in modified WRF/Noah (both at topmost soil layer depths of 10 vs. 1 cm) is recommended; however, we note that overall the modified WRF/Noah‐CMAQ with new dust module captures the 29 May 2011 dust event well.

Summer 2011 Evaluation of the Modified WRF/Noah‐CMAQ

In this section we evaluate the impacts and performance changes for the suite of modifications (described in section 3) to the off‐line (Table 2; run B2) and coupled WRF/Noah‐CMAQ (B3) model against the unmodified off‐line WRF/Noah‐CMAQ (B1) model for 2011 summer (June, July, and August; JJA). There are widespread increases in O3 over the west and central United States, with decreases over the east for O3, and widespread decreases in SO2 across CONUS, which are most concentrated in Ohio Valley and northeastern United States (Figures 8a and 8d). The changes lead to reductions in MB and ME for O3 and SO2 in CONUS (Table S5), especially over the eastern CONUS and Mid‐Atlantic states (e.g., Virginia, Pennsylvania, New York), where there is reduced MB and ME found at a predominant number (>60–80% in some cases) of AQS sites for the modified off‐line and coupled model.

Figure 8

Average summer (June–August) 2011 spatial difference (modified‐unmodified WRF/Noah‐CMAQ) plots of hourly O3, SO2, and NH3 and corresponding MB (mean bias) and ME (mean error) difference plots against hourly AQS O3 and SO2 observations (ppb) and weekly AMoN observations (μg m−3). Relative percent changes are shown in the difference plots and are based on Figures S7–S9 in the supporting information. WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; AQS = Air Quality Subsystem; AMoN = Ammonia Monitoring Network; CONUS = continental United States.

The results are very similar for the modified off‐line and coupled WRF/Noah‐CMAQ models, runs B2 and B3, respectively (see Table S5 and the spatial plots of absolute values and relative percent in Figures S7–S11), which further corroborates the short‐term impact assessments that compared the modified coupled against unmodified off‐line WRF/Noah‐CMAQ for air quality predictions (runs A1–A6). Thus, hereafter the results for either runs B2–B1 or B3–B1 will be uniformly referred to as modified‐unmodified WRF/Noah‐CMAQ.

In the analyses that follow, we use the results from the short‐term impact assessments to infer the relative effects of the incremental changes on the combined modifications in the 2011 summer runs. This method is inherently limited by some nonuniform behavior of each modification and its impact from late spring to summer, where some individual impacts may either be further enhanced or dampened. However, we note that given the nature of the modifications, there is confidence that the overall direction of change is consistent across the short‐term impact and full summer 2011 runs. This is further substantiated and discussed in the following analyses.

The enhanced DDEP_O3 and DDEP_SO2 in the modified model, which are due to mainly lower R and the CMC × 103 correction (discussion in section 3), reduces the SO2 and O3 mixing ratio, and demonstrates that these changes can reduce typical O3 overpredictions by similarly configured CMAQ simulations in eastern CONUS (e.g., Appel et al., 2017; Canty et al., 2015; Herwehe et al., 2011), which are likely impacted by an inaccurate dry deposition loss pathway. In the central and western CONUS, however, the large reductions in LAI in the modified VEGPARM table (Table 4 and Figure S4) lead to widespread reductions in DDEP_O3 and DDEP_SO2 that are impactful enough on O3 to result in increased concentrations and further O3 overpredictions for the modified WRF/Noah‐CMAQ, especially over the agricultural regions in the central and upper Midwest (Figure 8a). This may be due to either too coarse MODIS AVHRR resolution and/or underestimated in situ measured values of LAI in agricultural regions with spatially heterogeneous land use, which negatively affected the modified VEGPARM LAI values used here (Tian et al., 2002). The enhancements in SO2 for the modified model (B3–B1) occur in locations where the VEGPARM table changes overcompensate for the CMC × 103 correction factor, which can then lead to decreases in DDEP_SO2 that become more prevalent in the summer season. We note that the O3 MB and ME decrease further for the modified model in the west‐southwest CONUS (Figures 8b‐c), attributable to the updated DDEP_O3 to soil parameterization (Mészáros et al., 2009) under more realistic vegetation fractions (Ran et al., 2016).

There are improved diurnal O3 and SO2 patterns for the modified compared to unmodified WRF/Noah‐CMAQ in the daytime hours over CONUS for O3 (especially in the east) and for most hours for SO2 (Figure 9). Also, the CONUS‐wide averaged statistics show a majority (> 50%) of AQS sites have lower MB or ME for the modified model (Table S5).

Figure 9

Average summer (June‐July‐August) 2011 diurnal (a, b) O3 and (c, d) SO2 for unmodified (run B1; blue) and modified (run B3; red) WRF/Noah‐CMAQ simulations against AQS observations (gray) taken over CONUS (left) and east CONUS (right). WRF = Weather Research and Forecasting; CMAQ = Community Multiscale Air Quality; AQS = Air Quality Subsystem; CONUS = continental United States.

There are strong increases in average NH3 mixing ratio for the modified WRF/Noah‐CMAQ model over predominantly agricultural regions in the central and eastern CONUS (Figure 8g), where moderate to high fertilizer application occurs during the summer months (e.g., the U.S. Corn/Soybean Belts). This leads to increased MB and ME at most AMoN sites in the eastern and central United States. These changes are driven by reduced LAI in the modified VEGPARM table and a thinner topmost soil layer that induce temperature increases (see T2 max relative percent change ~ +24%; Figures S10d and S10e), which consequently alters the BIDI‐NH3 exchange (Bash et al., 2013) in the modified WRF/Noah‐CMAQ model. In much of the western CONUS that consists of minimally fertilized or unfertilized vegetation, there are widespread decreases in temperature (Figures S10d and S10e), increases in DDEP_NH3, (Figure S11) and decreases in NH3 for modified WRF/Noah‐CMAQ compared to the unmodified model (Figure 8g).

The interactions between the BIDI‐NH3 exchange model and the modifications to WRF/Noah‐CMAQ tend to exacerbate NH3 overpredictions, which consequently leads to overall performance degradation for weekly NH3 concentrations compared to the unmodified model for the AMoN network (Table S5). We note, however, that there is a low sample size for the AMoN comparison during this period (2011 JJA total ~ 100 sample pairs), and that the domain‐averaged modified‐unmodified WRF/Noah‐CMAQ NH3 differences (mean relative percent difference ~5%; Figure S9d and S9e) are comparable to the reported precision of AMoN NH3 measurements that were compared against collocated high‐resolution (in time) ground‐based samplers (e.g., mean relative percent difference ~9%) in Puchalski et al. (2015). While beyond the scope of this work, a more comprehensive NH3 analysis/evaluation for the modified WRF/Noah‐CMAQ system against a suite of satellite and ground‐based observations is recommended.

Conclusions

In this work, the WRF/Noah‐CMAQ modeling system has been improved for both off‐line and coupled applications. The individual and combined changes and impacts were analyzed, and model performance was evaluated for key meteorological and chemical (i.e., air quality) variables. The major changes to the modified model are to update the Noah soil (SOILPARM; following Kishné et al., 2017) and vegetation (VEGPARM) tables in WRF, and directly use the WRF/Noah stomatal (R ) and aerodynamic resistances (R ) in CMAQ. Other modifications include correcting a long‐standing error in the WRF/Noah canopy moisture content, revising the CMAQ dry deposition of O3 to soil (Mészáros et al., 2009), and revising the soil layer depths in WRF/Noah to better align with the soil depths in CMAQ.

Extensive impact assessment (i.e., sensitivity) tests and analyses illustrate that the major changes to T2, Q2, and WSPD10 in the modified WRF/Noah result from reducing LAI in the updated VEGPARM table (based on MODIS satellite data), updates to the SOILPARM table, and using a thinner top soil layer in the southwest CONUS. The modified WRF/Noah is most impactful on WSPD10, where the reduced LAI and enhanced roughness lengths for most LU categories in the updated VEGPARM table reduce NMB by up to 18%. Furthermore, T2 and LH_flx are either comparably simulated or improved by using the updated VEGPARM table and the mosaic land use representation (compared to dominant land use), but there is degradation in SH_flx performance.

In the modified off‐line and coupled WRF/Noah‐CMAQ simulations, the modified R and R and updated VEGPARM table drive the overall changes in dry deposition and mixing ratios of O3, SO2, and HNO3, and model performance is improved across much of the United States for O3 and SO2 compared to the unmodified off‐line WRF/Noah‐CMAQ system. In the eastern United States, the mean bias (error) for the modified WRF/Noah‐CMAQ system is reduced at 60% (53%) and 62% (69%) of the AQS sites for O3 and SO2, respectively. There are also large differences in the spatial patterns of NH3, which are due to the interactions of the modified system and the more complex nature of the BIDI‐NH3 exchange in CMAQ. The changes in NH3 evasion and deposition were driven entirely by the changes in the soil/meteorology conditions as the underlying assumptions of the mosaic/tiled deposition and emission processes in CMAQ with bidirectional NH3 exchange had not been modified to be completely consistent with the modified WRF/Noah. Supporting a more consistent parameterization of the NH3 bidirectional tiled emissions and deposition with WRF land surface models is an area of research currently under development. Although beyond the scope of this study, other limitations in this work include inherent model biases from other sources (e.g., model boundary conditions), as well as neglecting aerosol feedbacks in the coupled WRF‐CMAQ model. While preliminary, reducing the topmost soil layer depth to 1 cm to align with CMAQ shows a definitive impact on BIDI‐NH3 exchange, while also allowing the ability of the modified WRF/Noah‐CMAQ system to well simulate the timing and peak particle concentration during a windblown dust event.

In summary, the WRF/Noah‐CMAQ system has been updated to improve its physical consistency, include data tables that reflect more up‐to‐date data sets and prior studies, and increase its utility for applications to both retrospective and future air quality‐climate studies. To our knowledge, the modifications to the R and R introduced in this paper are the first to enable the use of the coupled WRF/Noah‐CMAQ model with scientific consistency for such applications. Furthermore, including aerosol feedbacks in the modified coupled WRF/Noah‐CMAQ would be a very interesting development/application study for a full year simulation, where the results of such an analysis would be unprecedented for this modeling system.

Supporting Information S1

Click here for additional data file.