Extension of a gaseous dry deposition algorithm to oxidized volatile organic compounds and hydrogen cyanide for application in chemistry transport models

The dry deposition process refers to flux loss of an atmospheric pollutant due to uptake of the pollutant by the Earth's surfaces, including vegetation, underlying soil, and any other surface types. In chemistry transport models (CTMs), the dry deposition flux of a chemical species is typically calculated as the product of its surface layer concentration and its dry deposition velocity (V_d); the latter is a variable that needs to be highly empirically parameterized due to too many meteorological, biological, and chemical factors affecting this process. The gaseous dry deposition scheme of Zhang et al. (2003) parameterizes V_d for 31 inorganic and organic gaseous species. The present study extends the scheme of Zhang et al. (2003) to include an additional 12 oxidized volatile organic compounds (oVOCs) and hydrogen cyanide (HCN), while keeping the original model structure and formulas, to meet the demand of CTMs with increasing complexity. Model parameters for these additional chemical species are empirically chosen based on their physicochemical properties, namely the effective Henry's law constants and oxidizing capacities. Modeled V_d values are compared against field flux measurements over a mixed forest in the southeastern US during June 2013. The model captures the basic features of the diel cycles of the observed V_d. Modeled V_d values are comparable to the measurements for most of the oVOCs at night. However, modeled V_d values are mostly around 1 cm s⁻¹ during daytime, which is much smaller than the observed daytime maxima of 2–5 cm s⁻¹. Analysis of the individual resistance terms and uptake pathways suggests that flux divergence due to fast atmospheric chemical reactions near the canopy was likely the main cause of the large model–measurement discrepancies during daytime. The extended dry deposition scheme likely provides conservative V_d values for many oVOCs. While higher V_d values and bidirectional fluxes can be simulated by coupling key atmospheric chemical processes into the dry deposition scheme, we suggest that more experimental evidence of high oVOC V_d values at additional sites is required to confirm the broader applicability of the high values studied here. The underlying processes leading to high measured oVOC V_d values require further investigation.

1 Introduction

Atmospheric pollutants impact human health and can also cause detrimental effects on sensitive ecosystems (Wright et al., 2018). Quantifying atmospheric deposition for atmospheric pollutants is needed to estimate their lifetimes in air and deposition rates to ecosystems. In the mass continuity equation of a chemistry transport model (CTM), atmospheric deposition is calculated separately for dry and wet deposition fluxes. Dry deposition refers to the removal process through which pollutants are taken up by the Earth's surface, and this process, while being quite slow, is a continuous process happening all the time, even during precipitation. In contrast, wet deposition is fast but episodic, and pollutants need to first be incorporated into hydrometeors before being delivered to the surface via precipitation. The amount of dry deposition of a pollutant of interest is typically calculated as the product of its ambient concentration and its dry deposition velocity (V_d), with V_d being calculated using empirically developed dry deposition schemes (Wesely and Hicks, 2000). In most V_d formulations, turbulent and diffusion effects are parameterized as aerodynamic and quasi-laminar resistance, respectively, above and sometimes also inside the canopy. Uptake effects by canopies, underlying soils, and any other surface types are parameterized as canopy (or surface) resistance, which includes several flux pathways such as to stomatal, cuticle, and soil. All of these flux pathways can be simultaneously affected by meteorological, biological, and chemical factors, most of which cannot be explicitly considered and are thus highly empirically parameterized in existing dry deposition schemes, which are known to have large uncertainties even for the most commonly studied chemical species such as O₃, SO₂, and more commonly measured nitrogen species with relatively rich flux datasets (Flechard et al., 2011; Wu et al., 2012, 2018)

Existing dry deposition schemes have thus far considered a small number of oxidized volatile organic compounds (oVOCs). Due to the lack of field flux data for oVOCs, V_d of these species is typically parameterized based on physicochemical properties, taking SO₂ and O₃ as references (Wesely, 1989; Zhang et al., 2003). In these existing schemes, V_d values of most oVOCs are of a similar order of magnitude as or slightly smaller than that of V_d of O₃. However, higher daytime V_d values for certain oVOCs than predicted by these schemes were recently reported by two studies (Karl et al., 2010; Nguyen et al., 2015). In one study Karl et al. (2010) found that V_d values of oVOCs calculated using existing schemes are about a factor of 2 lower than those based on canopy-level concentration gradient measurements over six forest and shrubland sites. V_d in their study was calculated from an inverse Lagrangian transport model with concentration gradient data as model input. The ratios of magnitudes between V_d(oVOCs) and V_d(O₃) in the study of Karl et al. (2010) are similar to those of Zhang et al. (2003) in that V_d(oVOCs) are slightly smaller than V_d(O₃) in both cases. However, the typical daytime V_d(O₃) over vegetated canopies is around 1 cm s⁻¹ in the literature from numerous studies (see summary in Silva and Heald, 2018), and the value in Karl et al. (2010) is much higher (e.g., up to 2.4 cm s⁻¹ at canopy top). One hypothesis explaining both high V_d(O₃) and high V_d(oVOCs) would be the reaction of O₃ with oVOCs, which depends on the chemical structure of the oVOC, but data required for validating this hypothesis are still lacking. We thus suspect that the very high V_d(oVOCs) values presented in Karl et al. (2010) were likely caused by atmospheric chemical processes not typically considered in the dry deposition process. High V_d(oVOCs) values were also observed over a temperate mixed forest in the southeastern US in a more recent short-term study (Nguyen et al., 2015), which again were suspected to be caused by atmospheric chemical reactions near vegetation surface. The flux measurements themselves also contain uncertainty. For example, Wu et al. (2015) showed that different measurement methods (e.g., flux gradient versus eddy correlation) resulted in very different daytime V_d(O₃) over the same forest canopy.

Hydrogen cyanide (HCN) is one of the most abundant cyanides present in the atmosphere (Singh et al., 2003) and is considered a biomass burning marker (Bunkan et al., 2013), but few existing studies have considered its dry deposition, which is critical to estimating the total sinks and atmospheric lifetimes of cyanides.

To meet the demands of modeling a large number of organic compounds in CTMs (Kelly et al., 2019; Moussa et al., 2016; Paulot et al., 2018; Pye et al., 2015; Xie et al., 2013), existing or newly developed air–surface exchange and dry deposition schemes need to be expanded to include additional oVOCs. At this stage with very limited knowledge of oVOC V_d, air–surface exchange models based on various theoretical and/or measurement approaches should be developed so that these models can be made available to the scientific community where such models are urgently needed as well as for future evaluation and improvement should more flux measurements become available. For example, Nguyen et al. (2015) modified the Wesely (1989) scheme to fit the flux data. A more sophisticated model, with a bottom-up approach, was adopted in Nizzetto and Perlinger (2012) to handle air–canopy exchange of semivolatile organic compounds.

The original dry deposition scheme of Zhang et al. (2003) includes 9 inorganic species and 22 organic species. Most of these 22 organic species are oVOCs formed from oxidation of nonmethane hydrocarbons. To take advantage of the recent flux dataset of a large number of oVOCs and HCN collected over a temperate forest (Nguyen et al., 2015), the present study extends the Zhang et al. (2003) scheme by including 12 additional oVOC species and HCN while keeping the same original model structure and theory. These additional oVOCs include hydroxymethyl hydroperoxide, peroxyacetic acid, organic hydroxy nitrates, and other multifunctional species that are mainly formed from the oxidation of biogenic VOCs (e.g., isoprene and monoterpenes). Model parameters for these newly included species are theoretically constrained based on the effective Henry's law constants and oxidizing capacities of the individual species as well as by considering the measured V_d values. Such an approach provides a top-down determination of V_d through comparison with measured (bottom-up) fluxes. Model–measurement comparison is conducted for V_d as well as resistance components and uptake pathways, results from which identify the major causes of model–measurement discrepancies. This study provides a computer code that is potentially useful for CTMs handling these oVOCs.

2 Methodology

2.1 Brief description of the V_d formulation

In the scheme of Zhang et al. (2003), V_d is calculated as follows:

$$\begin{array}{}\text{(1)}& V_{\mathrm{d}}\left(z\right)={\left(R_{\mathrm{a}}\left(z\right)+R_{\mathrm{b}}+R_{\mathrm{c}}\right)}^{-\mathrm{1}},\end{array}$$

where R_a is the aerodynamic resistance, R_b the quasi-laminar sub-layer resistance, R_c the surface resistance, and z the reference height above the vegetation. R_c is parameterized as

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{\mathrm{1}}{R_{\mathrm{c}}}}={\displaystyle \frac{\mathrm{1}-W_{\mathrm{st}}}{R_{\mathrm{s}}+R_{\mathrm{m}}}}+{\displaystyle \frac{\mathrm{1}}{R_{\mathrm{ns}}}},\text{(3)}& {\displaystyle \frac{\mathrm{1}}{R_{\mathrm{ns}}}}={\displaystyle \frac{\mathrm{1}}{R_{\mathrm{ac}}+R_{\mathrm{g}}}}+{\displaystyle \frac{\mathrm{1}}{R_{\mathrm{cut}}}},\end{array}$$

where R_s is the canopy stomatal resistance, R_m the mesophyll resistance, R_ns the non-stomatal resistance including resistance for uptake by leaf cuticles (R_cut) and by soil or ground litter (R_g), R_ac in-canopy aerodynamic resistance, and W_st the fraction of stomatal blocking under wet conditions.

R_s is calculated as follows:

$$\begin{array}{}\text{(4)}& {\displaystyle \frac{\mathrm{1}}{R_{\mathrm{s},\mathrm{i}}}}=G_{\mathrm{s}}\left(\mathrm{PAR}\right)f\left(T\right)f\left(D\right)f\left(\mathrm{\Psi }\right){\displaystyle \frac{D_{\mathrm{i}}}{D_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}}}.\end{array}$$

Here G_s(PAR) is the unstressed canopy stomatal conductance for water vapor, a function of photosynthetically active radiation (PAR). The dimensionless functions f(T), f(D), and f(ψ) range from 0 to 1, representing the fractional degree of stomatal closure caused by the stress from temperature, water vapor pressure deficit, and leaf water potential, respectively. $D_{{\mathrm{H}}_{\mathrm{2}}\mathrm{O}}$ and D_i are the molecular diffusivities for water vapor and the gas of interest, respectively.

R_cut and R_g for any chemical species are scaled to those of SO₂ and O₃ with two species (i)-dependent scaling parameters α(i) and β(i):

$$\begin{array}{}\text{(5)}& {\displaystyle \frac{\mathrm{1}}{R_{\mathrm{cut}/\mathrm{g}}\left(i\right)}}={\displaystyle \frac{\mathit{\alpha }\left(i\right)}{R_{\mathrm{cut}/\mathrm{g}}\left({\mathrm{SO}}_{\mathrm{2}}\right)}}+{\displaystyle \frac{\mathit{\beta }\left(i\right)}{R_{\mathrm{cut}/\mathrm{g}}\left({\mathrm{O}}_{\mathrm{3}}\right)}}.\end{array}$$

Details of the R_s-related formulas were described in Zhang et al. (2002), R_ns-related formulas in Zhang et al. (2003), and R_a and R_b formulas in Wu et al. (2018).

2.2 Extension of the scheme to additional chemical species

Dry deposition of a gaseous compound to most canopy types is mainly through non-stomatal uptake during nighttime and through both non-stomatal and stomatal uptake during daytime. The non-stomatal uptake depends on water solubility and reactivity of the species, which can be quantified by its effective Henry's law constant (H^*) and oxidizing capacity, respectively (Wesely, 1989; Zhang et al., 2002).

In the Supplement, Table S1 lists H^* values and Table S2 lists the oxidizing capacities for oVOCs and HCN considered in the present study. As shown in Eq. (5) above, two model parameters (α and β) are needed for every chemical species to calculate the non-stomatal uptake, with α being dependent on H^* and β dependent on oxidizing capacity. Initial α values were first given based on the relative magnitudes of H^* of all the chemical species and that of SO₂. Considering that the majority of the chemical species are very reactive, a value of 1.0 was used for β for most species and smaller values for a few less reactive species. α and β values were then adjusted based on the agreement of nighttime V_d between modeled values and measured fluxes obtained from a forest site in the southeastern US during summer (Nguyen et al., 2015). When adjusting α and β values, two rules were first applied: (1) the trends in α (or β) values between different chemical species should be consistent with the trends of their log (H^*) (or oxidizing capacity) (see Fig. S1 in the Supplement for the finalized α versus log (H^*)), and (2) modeled mean and median nighttime V_d should be mostly within a factor of 2.0 of the measured values (see discussion in Sect. 3.2 below). Only after these two rules were satisfied were the possible maximum α and β values chosen to reduce the gap between the modeled and measured daytime V_d, knowing that model-predicted V_d values were mostly lower than the measured ones. The finalized α and β values for the additional 12 oVOCs and HCN are listed in Table 1.

Table 1 List of model parameters needed in the scheme of Zhang et al. (2003) to simulate the dry deposition velocity of additional oVOC species and HCN: α and β are scaling parameters for non-stomatal resistance, and R_m is mesophyll resistance.

^a Treated as one group of compounds in the field measurements due to instrument limitations and have the same parameter values in the model. ^b The β value for HCOOH in Zhang et al. (2003) is 0.0, and here it is given as 0.2 to be consistent with other oVOC species (which would make no difference since the α value of 2 would dominate the non-stomatal resistance).

Model parameters chosen for the additional oVOCs and HCN can produce the magnitude of nighttime V_d for nearly all the chemical species, but they inevitably underpredicted daytime V_d for several oVOC species with very high measured daytime V_d values. We designed the model parameters this way due to the following considerations: (1) some of the chemical processes causing flux loss at the surfaces may be treated separately in the mass continuity equation in chemical transport models, (2) some of the oVOCs may also experience bidirectional air–surface exchange, and (3) more flux measurements are needed to confirm if the very high daytime flux for certain oVOCs is a universal phenomenon, noting that the existing data used here were from a short period of several days and over only one surface type.

Besides α and β, another chemical-species-dependent parameter that needs to be arbitrarily chosen is R_m. R_m for HCN was set to 100 s m⁻¹ based on its effective Henry's law constants and oxidizing capacities. Karl et al. (2010) found that enzymatic conversion can be an efficient pathway for the immobilization of oVOCs (e.g., methacrolein and methyl vinyl ketone, acetaldehyde, methacrolein) within the leaf interior, besides dissolution and oxidation, which suggests that the magnitude of R_m for oVOCs is minimal. Thus, the R_m for the oVOCs was set to 0 s m⁻¹ (Table 1).

2.3 Field flux data

The fluxes of 16 atmospheric compounds (including 13 oVOC species, HCN, hydrogen peroxide – H₂O₂, and nitric acid – HNO₃) were measured using the eddy covariance (EC) technique at the Centreville (“CTR”) Southeastern Aerosol Research and Characterization Study (SEARCH) site (hereinafter referred to as CTR). Note that formic acid (HCOOH) is the only overlapping oVOC species between the original Zhang et al. (2003) scheme and the flux measurement dataset. The CTR site (Brent, Alabama; 32.90^∘ N, 87.25^∘ W) is surrounded by a grassy field to the south and a temperate mixed forest that is part of the Talladega National Forest in all the other directions. The forest canopy is comprised of needleleaf coniferous (shortleaf, longleaf, and loblolly pine; ∼ 40 %) and broadleaf deciduous (primarily oak, sweetgum, and hickory; ∼ 60 %) tree species. The canopy height near the tower is on average ∼ 10 m with a leaf area index (LAI) of ∼ 4.7 m² m⁻². A 20 m metal walk-up tower is used as the main structure supporting instruments that measured the eddy covariance fluxes and related meteorological variables. The sonic anemometer and the gas inlet were mounted at a height of about 22 m, facing north toward the forest. Mixing ratios of gas-phase compounds were measured with negative-ion chemical ionization mass spectrometry (CIMS) at 8 Hz or faster. A database of half-hourly V_d for 16 atmospheric compounds covering 5 non-continuous days in June 2013 was obtained at the site. During these periods, the predominant winds were northerly, which is ideal to sample air from the forest (Fig. S2), and the requirement for energy balance closure was met (see Nguyen et al., 2015). At CTR, it was typically humid (RH 50 %–80 %) and warm (28–30 ^∘C) in the daytime during the experiment (Fig. S3). A comprehensive description of the V_d dataset, data processing protocols, instrumental methods, uncertainty analysis, and site characterizations can be found in Nguyen et al. (2015).

3 Results and discussion

3.1 Comparison of modeled resistance components

3.1.1 Atmospheric resistances (R_a and R_b)

For very reactive and soluble substances such as HNO₃ and H₂O₂, R_c is often assumed to be close to 0 (Hall and Claiborn, 1997; Meyers et al., 1989; Valverde-Canossa et al., 2006; Wesely and Hicks, 2000). The analysis of the measurement data showed that the daytime averaged V_d for HNO₃ and H₂O₂ fit the rate of deposition well without surface resistance ($V_{\mathrm{d}}=\mathrm{1}/[R_{\mathrm{a}}+R_{\mathrm{b}}]$) (Nguyen et al., 2015), which supports the assumption of near-zero R_c for HNO₃ and H₂O₂ over the mixed deciduous–coniferous CTR site under a humid environment. Therefore, the measured V_d of HNO₃ and H₂O₂ can be used to evaluate the modeled atmospheric resistances for those species (the sum of R_a and R_b). R_a represents the resistance for turbulent transport between the reference height and the surface and is not chemical-compound-specific. R_b quantifies the resistance for the mass transfer across the thin layer of air in contact with surface elements and is a function of the molecular diffusivity of a specific compound (Wesely and Hicks, 1977). In theory, the differences in R_b between any two gaseous species are only determined by differences in their molecular diffusivity at any given turbulent condition.

Figure 1 compares the modeled average diel variations of V_d for HNO₃ and H₂O₂ against observations. The measured V_d values for HNO₃ and H₂O₂ peaked around noon at about 4 and 6 cm s⁻¹, respectively, and were less than 1 cm s⁻¹ during the night. The model reproduced the diel pattern and captured the peak V_d values at noon well. During the early nighttime (hours 19–23), the modeled V_d values for HNO₃ and H₂O₂ were on the order of 1 cm s⁻¹, which is much higher than the measurements (< 0.2 cm s⁻¹). During the night, R_a dominates atmospheric resistance as it is usually much larger than R_b in magnitude. This discrepancy between the measurement and the model during the early night could be due to the stability correction functions used in the R_a calculation (the equations can be found in the article by Wu et al., 2018), which is subject to large uncertainties under nocturnal stable conditions (Högström, 1988). The measurements indicated that H₂O₂ deposited slightly faster than HNO₃, and the model reproduces this well, as shown in Fig. 1. Modeled R_b for H₂O₂ is always smaller than that for HNO₃ due to the smaller molecular weight and the larger molecular diffusivity. Overall, the model was in good agreement with the measurements regarding V_d for HNO₃ and H₂O₂, implying that the parameterization for atmospheric resistances (R_a and R_b) was reasonable for the site during the study period.

Figure 1 Comparison of the observed and modeled average diel variations of dry deposition velocities (V_d) for HNO₃ and H₂O₂. The shaded area indicates the standard deviation of the observations. The model assumes that surface resistances (R_c) for HNO₃ and H₂O₂ are zero.

3.1.2 Stomatal resistance (R_s)

Over vegetated areas, gas molecules can exit and enter the leaf through the stomata by molecular diffusion, similar to the leaf–air exchange of water vapor and CO₂. In dry deposition models, R_s for water vapor is estimated using evapotranspiration stomatal submodels, an approach that is also popular in the land surface and climate communities. R_s is extended to any gas species using the ratio of molecular diffusivity of the species of interest to that of water vapor (Pleim and Ran, 2011; Wesely and Hicks, 2000). Figure 2 compares the modeled canopy stomatal conductance ($G_{\mathrm{s}}=\mathrm{1}/R_{\mathrm{s}}$) for water vapor against the observation-based estimates. The observation-based G_s was estimated by using the inversion of the Penman–Monteith (P-M) equation (Monteith and Unsworth, 1990), which calculates R_s for water vapor by using measured water vapor fluxes and related meteorological data (e.g., humidity, temperature). The evaporation from soil water and liquid water on the vegetation surfaces is usually a minor contribution to the total water vapor flux observed above a forest canopy during daytime in summer. It was assumed that 85 % of the water vapor flux originated from transpiration in this study, following that used in the study of Turnipseed et al. (2006) at Duke Forest, North Carolina. Note that a value of 90 % was used by Clifton et al. (2017) at Harvard Forest, Massachusetts. The uncertainty of the calculated R_s related to the uncertainty in water vapor flux portion (on an order of 10 %) is much smaller than the differences between the modeled and the observation-based stomatal conductance (by a factor of 2) as discussed below.

Figure 2 Comparison of observation-based and modeled average diel variations of stomatal conductance (G_s) for water vapor. The shaded area indicates the standard deviation of the observation-based G_s(H₂O) estimated by the P-M method. “Mod-R_s,min” refers to a model sensitivity test in which R_s,min was reduced by 40 %.

As shown in Fig. 2, the model reproduced the basic diel pattern in G_s (i.e., highest values between 08:00 and 11:00), but the peak value is only about half of the observation-based values. The Jarvis-type stomatal submodel (Jarvis, 1976) is known for its linear dependence on the prescribed minimum stomatal resistance (R_s,min), a term that is subject to large uncertainties (Kumar et al., 2011; Wu et al., 2018, 2011). A series of tests conducted by iteratively adjusting the R_s,min values showed the modeled G_s to be in better agreement with observations if R_s,min was decreased by 40 % (Fig. 2). Modeled G_s with the adjusted R_s,min was in good agreement with the observation-based values most of the time, though the modeled values were slightly smaller than the observation-based estimates around noon. Analysis of the R_s parameterization indicates that this discrepancy was related to the stress function for water vapor pressure deficit (VPD) used in the Jarvis-type stomatal submodel, which may overpredict the stress on stomatal opening due to high VPD around noon.

3.1.3 Non-stomatal resistance (R_ns)

To assess if the non-stomatal resistance (R_ns) parametrization (Eq. 3) is reasonable, modeled $\mathrm{1}/R_{\mathrm{ns}}$ (defined as G_ns) values are compared with the non-stomatal portion of the flux, the inverse of which is termed the residual conductance (G_residual). G_residual includes all processes influencing deposition aside from R_a, R_b, R_m, and R_s, calculated as $[V_{\mathrm{d}}^{-\mathrm{1}}-(R_{\mathrm{a}}+R_{\mathrm{b}}){]}^{-\mathrm{1}}-(R_{\mathrm{s}}+R_{\mathrm{m}}{)}^{-\mathrm{1}}$. Here V_d is from the observations, R_a and R_b are calculated by the model driven by the observed meteorology, R_s is the observation-based estimate by the P-M method adjusted by the molecular diffusivity of each gas (similar to Eq. 4), and R_m is listed in Table 1. The uncertainties in individual resistance terms of Zhang et al. (2003) and several other dry deposition schemes have been thoroughly assessed by Wu et al. (2018), from which we believe G_residual estimated using the above formula is meaningful although with large uncertainties. The estimated G_residual can provide useful information on the $\mathrm{flux}/V_{\mathrm{d}}$ resulting from processes such as deposition to the leaf cuticle and ground (i.e., non-stomatal) or chemical loss due to reactions within and near the canopy that lead to flux divergence.

Figure 3 compares the observation-based G_residual for each oVOC species or HCN against the corresponding modeled non-stomatal conductance (G_ns) under different conditions. The mean and median values are presented in Table S3. During the nighttime when the canopy surface was dry (no dew), the G_residual for oVOC species ranged from 0.08 to 0.18 cm s⁻¹, and the modeled G_ns was comparable in magnitude. When the surface was wet from dew formation on leaves and needles, the oVOC species showed an increase in G_residual by 55 %–440 % compared to the nighttime dry surface. The model captured the increases in non-stomatal uptake when the surface become wet with dew, although it may underestimate (e.g., HDC₄, INP, HCN) or overestimate (e.g., PAA, DHC₄, HCOOH) the wetness effects. During the daytime of the study period, no precipitation was recorded at the CTR site (Fig. S3) and the canopy surface was dry. The mean G_residual for oVOCs ranged from 0.5 to 8.7 cm s⁻¹ during the daytime, which is much higher than the modeled G_ns for most species (0.2–1 cm s⁻¹). Figure S4 presents the diel variations of G_residual and G_ns, and it shows that the modeled G_ns had smaller diel variations than those of G_residual; large differences in magnitude can be seen during the daytime. The modeled G_ns showed a peak during the early morning (around 07:00), which may be due to the enhanced non-stomatal uptake by dew-wetted surfaces.

Figure 3 Box plot of the observation-based residual conductance (G_residual) and the modeled non-stomatal conductance (G_ns) during the nighttime dry period (ND, n=88), nighttime wet period (NW, n=40), and daytime dry period (DD, n=85). In each box, the central mark is the median, and the edges of the box are the 25th and 75th percentiles. The filled dots represent the arithmetical mean of data between the 25th and 75th percentiles. Daytime is 09:00–17:00 (local time) and nighttime is 20:00–06:00 (local time). The wet surface conditions were determined in the model driven by the observations of relative humidity, precipitation rate, friction velocity, and temperature.

3.2 Evaluation of modeled deposition velocities

Figure 4 shows a model–measurement comparison of diel V_d of the oVOCs and HCN, and Table 2 presents the statistical results of the comparison. As described in Sect. 2, the assigned α and β values should first produce reasonable nighttime V_d. Modeled nighttime mean V_d values were very close to measurements for the majority of the chemical species, although the differences were somewhat larger for the median values (Table 2). Three species (HAC, HPALD, PROPNN) still had 50 % lower modeled than measured nighttime mean V_d but had slightly higher modeled than measured nighttime median V_d. In contrast, modeled daytime mean V_d values were more than 50 % lower than the measured values for four species (HMHP, PAA, HPALD, ISOPOOH/IEPOX) and were also significantly lower for several other species. Only three species (MTNP, HCN, HCOOH) had comparable modeled and measured V_d for both daytime and nighttime. One species (DHC₄) had slightly lower modeled than measured daytime mean or median V_d, but with an opposite trend for nighttime V_d.

Table 2 Statistical results of the observed and modeled dry deposition velocity (V_d) for oVOCs and HCN (cm s⁻¹)^*.

^* Note: N is the number of samples; R is the correlation coefficient between observation (Obs) and model simulation (Mod); “Mod-R_s,min” refers to a sensitivity test in which R_s,min was reduced by 40 %; daytime is 09:00–17:00 (local time) and nighttime is 20:00–06:00 (local time). Median values are provided in parentheses, following arithmetic mean values.

Figure 4 Comparison of averaged diel cycles of observed and modeled dry deposition velocities (V_d) of oVOCs and HCN. The shaded area indicates the standard deviation of the observations. “Mod-R_s,min” refers to a sensitivity test in which R_s,min was reduced by 40 %.

The model reproduced the basic features of the diurnal pattern of the observations, showing the highest values during the day and the lowest values at night. Correlation coefficients between the measurement and the model ranged from 0.52 to 0.77. At night, the measured V_d for the oVOCs remained relatively low, typically ranging from 0.1–0.5 cm s⁻¹, and the model produced the same magnitudes for most of the species. During the daytime, the model can only capture the magnitudes of the measured V_d for a few species (e.g., HCN, HCOOH, MTNP, DHC₄), the peak V_d values of which were less than 1.5 cm s⁻¹. For the other species, the measured peak V_d values were in the range of 2 to 5 cm s⁻¹, while the modeled results were below 1 cm s⁻¹. As shown in Sect. 3.1.2, the modeled G_s was likely underestimated when compared to the simultaneous measurements of water vapor flux. Adjusting G_s higher by 67 % (through reducing R_s,min by 40 %) can only increase the modeled V_d of the oVOCs by 10 %–40 % during the daytime (see the sensitivity test in Fig. 4), and the peak values were still mostly below 1 cm s⁻¹. Figure 5 shows that the model captured the differences in measured V_d for the oVOCs to some extent. The model–measurement agreements were good for species with the measured mean V_d below 0.5 cm s⁻¹, above which the discrepancy increased. For the measurements, the mean values were significantly larger than the median values, especially for the fast-deposited species, indicating that the distribution of the measured V_d values skewed to the right (high values). The model has better agreement with the measurements by comparing the median versus mean values.

Figure 5 Box plot of observed and modeled hourly dry deposition velocities (V_d) of oVOCs and HCN. In each box, the central mark is the median, the edges of the box are the 25th and 75th percentiles, and the whiskers extend to the 10th and 90th percentiles. The filled dots represent the arithmetical mean of all the data. “Mod-R_s,min” refers to a sensitivity test in which R_s,min was reduced by 40 %.

3.3 Fast chemical reactions as potential causes of the daytime model–measurement discrepancies

At night when stomata are mostly closed and atmospheric chemical reactions are largely inhibited, the measured fluxes above the canopy should better represent non-stomatal surface uptake. In the presence of sunlight, fast chemical reactions between the inlet and canopy could make a significant or even dominant contribution to the measured fluxes of reactive species (Cape et al., 2009; Farmer and Cohen, 2008; Wolfe et al., 2011). The impact of fast chemical reactions on surface fluxes should be different for different chemical species. To verify this hypothesis, two chemical species (HAC and PAA) having similar molecular weights (74 and 76 Da, respectively) but very different daytime fluxes were compared (Fig. 6). Their similar molecular diffusivities (controlled by molecular weight) suggest that they should be transferred through the quasi-laminar sub-layer and taken up through leaf stomata at similar rates, resulting in similar resistance components of R_b and R_s. Note that R_a is universal to any trace gases and R_m is assumed to be negligible. Thus, the differences between their V_d should be caused by their different non-stomatal sinks. At night, V_d values were similar between HAC and PAA (median values: 0.04 cm s⁻¹) over dry surfaces. When the surfaces were wet due to dew formation, V_d for both HAC and PAA increased (median values: 0.30–0.48 cm s⁻¹). In contrast, V_d(PAA) was much higher than V_d(HAC) during daytime, suggesting additional or larger sinks for PAA compared to HAC. The reactivity parameters listed in Table S2 in the Supplement also suggest that PAA is more reactive than HAC. Thus, fast chemical processing and subsequent flux divergence above the canopy likely caused the large discrepancies between the measured and modeled V_d for the reactive oVOC compounds during the daytime.

Figure 6 Scatter plot of the measured dry deposition velocities (V_d) for hydroxyacetone (HAC) and peroxyacetic acid (PAA) during nighttime (20:00–06:00, local time) and daytime (09:00–17:00, local time). The shaded (white) cycles correspond to wet (dry) surface conditions.

Chemical processes can indeed cause flux divergence or convergence at the surface, which has been supported by growing evidence from field measurements (e.g., Farmer and Cohen, 2008; Min et al., 2014; Wolfe et al. 2009). For example, Wolfe et al. (2009) suggested that the differences in loss rate between the inlet and canopy may be an important contributor to the measured net flux of peroxyacetyl nitrate, irrespective of turbulent timescales. Photochemical OH production is reduced within canopies, which in turn slows down the oxidation of volatile organic compounds and the photolysis of organic nitrates. The oVOCs measured at the CTR site are mainly produced from the oxidation of isoprene and monoterpenes (Nguyen et al., 2015). Most of the oVOCs are quite chemically reactive and can undergo fast oxidation (e.g., multifunctional carbonyls), decomposition (e.g., HMHP), or photolysis (e.g., organic nitrates) (Müller et al., 2014; Nguyen et al., 2015). Vertical gradients in the chemical production and loss rates below the inlet can exhibit chemical flux divergence, which contributes to the net flux above the canopy. Quantifying the effects of chemical processing on the net flux would require a multi-layer model with resolved emission, deposition, turbulent diffusion, and chemical processes throughout the canopy, which is recommended for future studies (e.g., Ashworth et al., 2015; Bryan et al., 2012; Stroud et al., 2005; Wolfe and Thornton, 2011; Zhou et al., 2017).

Quantifying V_d as the ratio of flux to concentration at one measurement height only ($V_{\mathrm{d}}=F/C_{zr}$), rather than as the ratio of flux to the concentration difference at the measurement height and the surface ($V_{\mathrm{d}}=F/[C_{zr}-C_{\mathrm{0}}]$), although commonly employed in analyzing eddy covariance flux measurements, is a simplification. It is valid for (1) matter that disappears nearly completely by reactions at the surface and (2) unstable or neutral conditions. Most chemical species considered here may satisfy the first condition. With regards to the second condition, our analysis is based on the assumption that, under stable conditions at nighttime, concentrations observed at the measurement height change in relation to the fluxes measured at this height. However, no relation between measured concentration and flux is typically observed due to the presence of a shallow stable boundary layer, connection between the stable free atmosphere and stable boundary layer by internal gravity waves, ground inversions, and low-level jets, leading to intermittent turbulence at the measurement height containing a gravity wave signal and non-steady-state conditions (Foken, 2017). Future efforts to model oVOC and HCN deposition velocities above forest canopies should be based on neutral or unstable boundary layer flux measurements only or, for example, on a modified Bowen ratio flux measurement in which concentrations are measured at two heights in the constant flux layer. Such an approach can provide a means to compute a measured deposition velocity of a surface-reactive substance as proportional to the ratio between the measured flux and the measured concentration difference.

4 Summary and recommendations

The number of chemical species simulated in chemical transport models (CTMs) has been increasing with increasing computer power. Among these, oVOCs and HCN are important groups of atmospheric pollutants for which dry deposition processes need to be treated as accurately as possible so that their inputs to ecosystems (noting that some oVOCs are organic nitrogen) and their roles in other atmospheric chemistry processes (e.g., formation of ozone and secondary organic aerosols) can be assessed. Earlier dry deposition schemes considered very few oVOCs and need to be extended for more species. Dry deposition of HCN was assumed to be negligible in some CTMs (e.g., Moussa et al., 2016). The present study first generated the effective Henry's law constant and oxidizing capacity, the two key physical and chemical properties that are considered to control the dry deposition process (Wesely and Hicks, 2000), for 12 oVOC species and HCN. Two scaling factors for the non-stomatal resistance and one for the mesophyll resistance were applied to individual oVOCs and HCN for calculating their respective V_d.

The modeled nighttime V_d agrees well with the measured data for most of the oVOCs, suggesting that the current non-stomatal parameterization scheme is a reasonable approach. The stomatal conductance for water vapor, with adjusted (reduced) R_s,min, also agrees well with measured values. However, the modeled peak V_d values during daytime are only a fraction (0.2–0.5) of the measured values for some of the oVOCs, suggesting that fast atmospheric chemical processes likely contributed to the total measured fluxes. In practice, these additional fluxes during daytime can be modeled as non-stomatal uptake, and better model–measurement agreement can be obtained by adjusting the non-stomatal parameterization scheme (e.g., Müller et al., 2018; Paulot et al., 2018). However, using this approach will produce unreasonably high values for the solubility parameter and overpredict V_d during nighttime if the same non-stomatal formulas are used for both day and nighttime (as is the case in the existing schemes). More importantly, the high measured V_d values have only been observed at relatively few sites during very short periods (Karl et al., 2010; Nguyen et al., 2015). More evidence is needed to parameterize V_d for oVOCs to different land use categories over entire seasons. Until then, the conservative estimates of V_d such as modeled in this study are still recommended for use in CTMs. The model parameters chosen for V_d of these oVOCs provide the best-known representation of their respective physicochemical properties, and the modeled V_d values fall within the range of the low-end values of the available measurements.

Future field studies should focus on conducting flux measurements of oVOC compounds with the highest uncertainties, such as those that are most chemically reactive in the atmosphere or most rapidly taken up by wet surfaces. Additional measurements are also needed in different ecosystems to inform the representativeness of the high oVOC V_d reported by Nguyen et al. (2015) and Karl et al. (2010). Furthermore, concurrent chemical measurements of oxidants such as O₃ and radicals are needed to quantify flux divergence due to fast within- and near-canopy chemical reactions. Future dry deposition schemes should include additional biochemical processes and species-dependent parameters for non-stomatal uptake, including enzymatic reactions (Karl et al., 2010), the octanol–air partitioning coefficients to account for cavity formation and polar intermolecular interactions with leaf surfaces and reservoirs (Nizzetto and Perlinger, 2012), and the enhancement and reduction effects due to soil and leaf moisture. Chemical processes within the canopy airspace could also be coupled with emission and deposition schemes to realistically simulate chemical fate and transport, including bidirectional fluxes of the reactive compounds discussed here and less reactive compounds such as methanol. Such an approach would require specification of chemical conditions within and near the canopy as well as in-canopy radiation and airflow. While more computationally intensive, the results presented here reinforce the need for such advanced models to explicitly resolve the non-stomatal processes contributing to the net atmosphere–biosphere exchange of reactive compounds. Above all, intercomparison studies should first be conducted for existing models that can handle oVOC dry deposition processes to quantify the magnitudes of uncertainties in the simulated V_d as well as the associated ambient concentrations and deposition fluxes.