Representing sub-grid scale variations in nitrogen deposition associated with land use in a global Earth System Model : implications for present and future nitrogen deposition fluxes over North America

Reactive nitrogen (N) emissions have increased over the last 150 years as a result of greater fossil fuel combustion and food production. The resulting increase in N deposition can alter the function of ecosystems, but characterizing its ecological impacts remains challenging, in part because of uncertainties in model-based estimates of N dry deposition. Here, we leverage the tiled structure of the land component (LM3) of the Geophysical Fluid Dynamics Laboratory (GFDL) 5 Earth System Model to represent the impact of physical, hydrological, and ecological heterogeneities on the surface removal of chemical tracers. We show that this framework can be used to estimate N deposition at more ecologically-relevant scales (e.g., natural vegetation, water bodies) than from the coarse-resolution global chemistry–climate model (GFDL-AM3). Focusing on North America, we show that the faster removal of N over forested ecosystems relative to cropland and pasture implies 10 that coarse resolution estimates of N deposition from global models systematically underestimate N deposition to natural vegetation by 10 to 30% in the Central and Eastern US. Neglecting the subgrid scale heterogeneity of dry deposition velocities also results in an underestimate (overestimate) of the amount of reduced (oxidized) nitrogen deposited to water bodies. Overall, changes in land cover associated with human activities are found to slow down the removal of N from the atmosphere, 15 causing a reduction in the dry oxidized, dry reduced, and total N deposition over the contiguous US of 8%, 26%, and 6%, respectively. We also find that the reduction in the overall rate of removal of N associated with land-use change tends to increase N deposition on the remaining natural vegetation and facilitate N export to Canada. We show that subgrid scale differences in the surface removal of oxidized and reduced nitrogen imply that near-term (2010–2050) changes in oxidized (-47%) 20 1 Atmos. Chem. Phys. Discuss., https://doi.org/10.5194/acp-2018-572 Manuscript under review for journal Atmos. Chem. Phys. Discussion started: 19 July 2018 c © Author(s) 2018. CC BY 4.0 License.

The goal of this study is to develop a framework to diagnose ecosystem-specific N dry deposition fluxes within a global chemistry climate model on decadal to centennial timescales.First we describe the coupling of the Geophysical Fluid Dynamics Laboratory (GFDL) land model (LM3) to the GFDL atmospheric chemistry-climate model (AM3) to represent the impact of natural (e.g., vegetation type, soil, and canopy wetness) and man-made (e.g., deforestation, cropping) heterogeneities on dry deposition.We then show that the tiled structure of LM3 can be leveraged to derive N deposition on a more ecologically relevant scale (e.g., deposition on water bodies or natural vegetation).Finally, we discuss how this framework can be used to better represent the impact of land-use change and future trends in N emissions on N deposition.

Model description
We use an updated version of the GFDL-AM3 (Donner et al., 2011;Naik et al., 2013;Paulot et al., 2016) to simulate atmospheric dynamics and chemistry.Except for the treatment of dry deposition, the model configuration is identical to the one recently described by Paulot et al. (2016) and Paulot et al. (2017), including updates to wet deposition and the chemistry of sulfate and nitrate.The horizontal resolution of the model is 200 km with 48 vertical levels.
In AM3, the surface removal of chemical tracers is calculated using a prescribed monthly climatology of dry deposition velocities (Naik et al., 2013;Paulot et al., 2016).The lack of a dynamic representation of dry deposition reduces the ability of the model to capture the impact of past and future variability in environmental conditions (e.g., drought, climate change; Wu et al., 2016) and land-use change on atmospheric chemistry.We note that these limitations are not specific to AM3 but affect all chemical transport models that do not include a comprehensive land model (Ellis et al., 2013;Ran et al., 2017).
Here, we describe the development of a new model, in which dry deposition of gaseous and aerosol species is calculated within the dynamic vegetation model LM3 (Shevliakova et al., 2009;Milly et al., 2014).The combined model will be referred to as AM3-LM3-DD hereafter.
LM3 is a comprehensive climate land model that includes detailed representations of vegetation dynamics and hydrology and is designed to be run over decadal to centennial timescales under both historical and future conditions.LM3 can be run both coupled with AM3 and in stand-alone mode with prescribed meteorological fields (Milly et al., 2014).
In LM3, the heterogeneity of the land surface and vegetation is represented using a sub-grid mosaic of tiles (Shevliakova et al., 2009;Malyshev et al., 2015) as illustrated in Fig. 1.Each tile has distinct energy and moisture balances for a vegetation-snow-soil column, biophysical properties, and exchanges of radiant and turbulent fluxes with the over-Figure 1.Schematic representation of the resistance scheme used to represent the dry deposition of gaseous tracers for each tile.R a , R b,i , R ac , R m , R s , and R sf,i are the aerodynamic resistance, laminar resistance, canopy aerodynamic resistance, mesophyll resistance, stomatal resistance, and surface resistance, respectively.The g, s, and v indexes (i) refer to ground, stem, and vegetation.Note that for clarity, deposition on soil and vegetation that are covered by snow or liquid water are not shown.lying atmosphere.LM3 predicts physical, biogeochemical, and ecological characteristics for each sub-grid land surface tile from the top of the vegetation canopy to the bottom of the soil column, including leaves and canopy temperature, canopy-air specific humidity, stomatal conductance, snow cover and depth, runoff, vertical distribution of soil moisture, ice, and temperature.The land-use history is prescribed from the Hurtt et al. (2011) reconstruction for each grid cell in terms of annual transition rates among four distinct land-use types: undisturbed (hereafter referred to as natural), crops, pastures, and secondary vegetation.Secondary vegetation is defined in LM3 as the vegetation recovering after land-use and land-cover changes and not currently managed.This includes all abandoned agricultural land as well as the land where wood was harvested at least once in prior years.The model keeps track of different recovery states by creating a secondary vegetation tile every time a disturbance occurs and simulating the subsequent vegetation regrowth in the tile.To avoid unrestricted growth of the number of tiles, the number of secondary vegetation tiles is limited to 10 per grid cell in the configuration of LM3 used here.When more than 10 secondary vegetation titles exist in a grid cell, secondary vegetation tiles with similar properties are merged (Shevliakova et al., 2009), while preserving water, energy, and carbon balances.Land properties that affect the surface removal of chemical tracers, such as snow cover, canopy liquid water and snow mass, surface and canopy temperature, leaf area index (LAI), stomatal conductance, and vegetation height are all prognostic (Shevliakova et al., 2009).Vegetation carbon is partitioned into five pools: leaves, fine roots, sapwood, heartwood, and labile storage.The model simulates changes in vegetation and soil carbon pools, as well as the carbon exchange among these pools and the atmosphere.The sizes of the pools are modified daily depending on the carbon uptake according to a set of allocation rules.Additionally, the model simulates changes in the vegetation carbon pools due to phenological processes, natural mortality, and fire.LAI is determined by vegetation leaf biomass and specific leaf area, prescribed for each vegetation type.Each vegetated tile has a unique vegetation type (C 3 grass, C 4 grass, temperate decidu-ous, coniferous, or tropical vegetation), which is determined based on biogeographical rules that take into account environmental conditions as well as vegetation biomass in each tile (Shevliakova et al., 2009).The fraction of the canopy covered by liquid water (f l ) and snow (f s ) are estimated from the intercepted canopy liquid water mass (w l ) and snow mass (w s ) following Bonan (1996): where W l,max = 0.02 kg m −2 and W s,max = 0.2 kg m −2 are the maximum liquid water and snow holding capacities, respectively.If both snow and liquid water are present simultaneously, water and snow are assumed to be distributed independently of each other.
The representation of management practices is important in determining the impact of land-use change on dry deposition, as it affects the vegetation type, and the seasonality of the vegetation cover.In LM3, crop harvesting and pasture grazing are performed annually at the end of the calendar year (Malyshev et al., 2015).Previous work has shown that this treatment contributes to an underestimate of the impact of management on land cover (Malyshev et al., 2015).To address these biases, we make the following modifications.For pasture, we assume that 10 % of leaf biomass is removed daily by grazing, provided LAI exceeds 2 to avoid overgrazing.This higher grazing frequency and intensity are needed to avoid the excessive growth of vegetation biomass on pasture in the tropics and midlatitudes, a problem which was noted in previous versions of LM3 (Malyshev et al., 2015) leading to misclassification of pasture vegetation cover as forests (Malyshev et al., 2015).LM3 does not estimate the cropping schedule (e.g., Bondeau et al., 2007), so we specify planting and harvesting dates from the global monthly irrigated and rainfed crop area climatology (Portmann et al., 2010).The impact of management practices on the timing and magnitude of agricultural emissions (e.g., Paulot et al., 2014) is not accounted for in AM3-LM3-DD.
The tiled structure of LM3 is especially useful to diagnose fluxes to areas, such as natural vegetation or water bodies, which are generally not well represented by the average properties of the grid box, in which they are located, because of their small geographical extent (Fig. S1 in the Supplement).
Briefly, the aerodynamic resistance (R a ) to the exchange of tracers between the canopy and the atmosphere is determined using the Monin-Obukhov similarity theory.Within the canopy, the aerodynamic resistances to the ground (R ac,g ) and to the vegetation (R ac,v ) are independent of the chemical tracer and taken from Erisman (1994) and Choudhury and Monteith (1988), respectively.
where SAI, h, V, and u are the stem area index (unitless), the height of the vegetation (in meters), the normalized wind (m s −1 ) at the top of the canopy, and the friction velocity (m s −1 ), respectively.Note that unlike Erisman (1994), we include SAI in the calculation of R ac,g , which tends to reduce deposition to the ground in winter.
We focus next on the representation of the dry deposition of gases, which is much faster than that of fine particles (Zhang et al., 2002).
Following Jensen and Hummelshøj (1995) and Jensen and Hummelshøj (1997), the canopy laminar resistance (R b,v ) is defined as where lw is the characteristic obstacle length of the canopy (in m, Table S1), ν is the kinematic viscosity, and D X is the diffusivity of species X.Following Hicks et al. (1987), the stem laminar resistance is where P r is the Prandtl number, Sc(X) is the Schmidt number, i.e., the ratio of the kinematic to the mass diffusivity (Sc ∝ D −1 X ), and κ is the von Kármán constant (κ = 0.4).Similarly, the ground surface laminar resistance is where u g is the friction velocity near the ground (Loubet et al., 2006).The mesophyll resistance is expressed following Wesely (1989): The stomatal resistance (R s (X)) is calculated as where M(X) is the molecular weight of species X and R s (H 2 O) is the stomatal resistance for water vapor, calculated according to the Leuning model (Leuning, 1995;Milly et al., 2014).This model accounts for the impact of water stress and CO 2 concentration, which have been shown to modulate the response of surface ozone to drought (Huang et al., 2016) and CO 2 increase (Sanderson et al., 2007).Cuticle (v), stem (s), and ground (g) resistances for species X are parameterized based on SO 2 and O 3 : where R sf,i (SO 2 ) and R sf,i (O 3 ) are tabulated resistances (Table S1) for each surface type, α(X) and β(X) are weighting factors (Table S2) estimated using the solubility (for α) and reactivity (for β) of X (Wesely, 1989;Zhang et al., 2002), s(T ) is a temperature adjustment factor (Zhang et al., 2003), and γ (X) is a co-deposition adjustment, which reflects changes in R sf,i (X) associated with surface acidity (Erisman et al., 1994;Massad et al., 2010;Neirynck et al., 2011;Wu et al., 2016).Here, we use the parameterizations of Massad et al. (2010) for NH 3 and Simpson et al. (2003) for SO 2 : To avoid unrealistic oscillations in v d (NH 3 ) and v d (SO 2 ), we estimate the acid ratio (r SN ) using the ratio of the 24 h integrated total dry deposition of acids to the dry deposition of ammonia and ammonium, rather than using the ratio of their surface concentrations (Massad et al., 2010;Simpson et al., 2003).
The bidirectional exchange of ammonia is not represented in AM3-LM3-DD (Massad et al., 2010;Flechard et al., 2013).This reflects in part uncertainties in the emission potential of vegetation and the lack of detailed treatment of agricultural activities in LM3 (Riddick et al., 2016).We thus expect AM3-LM3-DD to overestimate NH 3 dry deposition in source regions (Zhu et al., 2015;Sutton et al., 2007).

Experimental design
We perform two sets of global simulations representative of present-day (circa 2010) and future (2050) conditions.For present-day conditions, AM3-LM3-DD is run from 2007 to 2010 using 2007 as a spin-up.The model is forced with observed sea surface temperatures and sea ice cover, and land use from the Representative Concentration Pathway 8.5 scenario (RCP8.5;Riahi et al., 2011).Anthropogenic emissions are from the Hemispheric Transport of Air Pollution 2 (HTAPv2; Janssens-Maenhout et al., 2015).Natural emissions are based on Naik et al. (2013), except for isoprene emissions, which are calculated interactively using the Model of Emissions of Gases and Aerosols from Nature (MEGAN; Guenther et al., 2006).This simulation will be referred to as R2010 hereafter.An additional sensitivity experiment is performed (R2010_no_lu) in which natural vegetation is assumed to cover all vegetated tiles (i.e., no human land use).In both experiments, horizontal winds are nudged to those from the National Centers for Environmental Prediction reanalysis (Kalnay et al., 1996) to minimize meteorological variability between R2010 and R2010_no_lu.
For 2050, we use the vegetation, sea surface temperatures, and sea ice cover simulated by the GFDL-CM3 model under the RCP8.5 scenario in 2050 (Levy et al., 2013).RCP8.5 anthropogenic emissions for 2050 are used (Lamarque et al., 2011) except for NH 3 , where we use the spatial distribution and seasonality of HTAPv2 emissions following Paulot et al. (2016).The model is run for 10 years with land use fixed to year 2050, and we use the average of the last 9 years to minimize the impact of internal variwww.atmos-chem-phys.net/18/17963/2018/Atmos.Chem.Phys., 18, 17963-17978, 2018 ability.This simulation will be referred to as R2050 hereafter.We perform two additional sensitivity experiments to characterize how land-use change (R2050_2010lu) and climate (R2050_2010climate) contribute to the change in deposition velocity between R2010 and R2050.The different model runs are summarized in Table 1.
3 Results and discussion

Evaluation of simulated v d against observations
The resistance approach for calculating dry deposition velocities implemented in AM3-LM3-DD is similar to that used in most chemical transport models.However differences in implementations can result in large differences between simulated deposition velocities (Wu et al., 2018).To illustrate these differences, Fig. 2 shows the sensitivity of v d (SO 2 ) and v d (NH 3 ) to temperature, wetness, and surface acidity in three global models: MOZART (Emmons et al., 2010), GEOS-Chem (Wang et al., 1998), and AM3-LM3-DD.Under dry conditions, GEOS-Chem and AM3-LM3-DD produce identical results for v d (SO 2 ), with the temperature dependence driven by that of the stomatal conductance.At low and high temperatures, v d (NH 3 ) is faster in AM3-LM3-DD than GEOS-Chem, which reflects small differences in the assumed surface pH (6.35 and 6.6, respectively).In contrast, MOZART assumes a surface pH = 5 and accounts for changes in the effective solubility of SO 2 and NH 3 with temperature, similar to Nguyen et al. (2015).The increase in solubility with decreasing temperature results in faster v d (X) at cold temperatures in MOZART, while the lower pH increases v d (NH 3 ) and decreases v d (SO 2 ).The impact of surface wetness on v d (X) is only considered in MOZART and AM3-LM3 DD.In MOZART the presence of dew more than doubles v d (SO 2 ) but reduces v d (NH 3 ) below 25 • C. In contrast, both v d (NH 3 ) and v d (SO 2 ) increase in AM3-LM3-DD when the canopy is wet, which is supported by observations (Erisman et al., 1994(Erisman et al., , 1999;;Massad et al., 2010).AM3-LM3-DD also accounts for the modulation of R sf,v (SO 2 ) and R sf,v (NH 3 ) by the acidity of the surface.Our results suggest that when α SN = 2 , i.e., when the deposition of acids is twice as large as the deposition of bases, the impact of codeposition can be greater than that of canopy wetness.Our comparison suggests that the implementation of the Wesely scheme in MOZART, AM3-LM3 DD, and GEOS-Chem produce similar v d (SO 2 ) and v d (NH 3 ) (within 50 %) under dry conditions and for temperatures close to 20 • C.However, differences in the sensitivity of v d (SO 2 ) and v d (NH 3 ) to environmental conditions (temperature, wetness, acidity) can result in large differences (> 2).Such differences highlight the need for detailed evaluation of v d (X) across a wide range of conditions and chemical species (Wu et al., 2018).
We first evaluate the simulated present-day (R2010) v d (SO 2 ) against a compilation of field-based v d (SO 2 ) observations (Table S3).We sample the simulated monthly v d (SO 2 ) at the location of the measurements in the tile that best represents the type of vegetation reported in the observations.When observations are available, we further distinguish between daytime and nighttime as well as wet and dry conditions.For daytime and nighttime observations, we sample the model from 08:00 LT to 17:00 LT and 22:00 LT to 04:00 LT, respectively.For wet conditions, we sample the model when the canopy wetness is greater than 10 %. Figure 3 shows observed and simulated v d (SO 2 ) grouped among the four types of vegetation simulated by LM3 (deciduous, coniferous, tropical, and grass).
Simulated deposition velocities generally fall within a factor of 2 of the observations, with better agreement during the day than at night, when the model is biased high.This uncertainty range is similar to that reported by Wu et al. (2018) in different dry deposition models.More specifically, AM3-LM3-DD qualitatively captures the range of deposition velocities over forested ecosystems, including the slower deposition of SO 2 in winter than in summer and under dry than under wet conditions in deciduous forests, and the fast removal of SO 2 over coniferous forests.However, the model fails to capture the elevated v d (SO 2 ) (> 1 cm s −1 ) reported by several studies over grasslands.This may reflect uncertainties in the representation of ammonia emissions (e.g., no sub-grid heterogeneity), which could result in an underestimate of SO 2 -NH 3 co-deposition over crops or fertilized grasslands (Nemitz et al., 2001;Flechard et al., 2013).
Figure 4 shows the observed deposition velocities for HNO 3 , a range of organic nitrates (ISOPN, MVKN, PROPNN) derived from isoprene photooxidation (Paulot et al., 2009), HCN, and H 2 O 2 .We refer the reader to Nguyen et al. (2015) for information regarding the site and Caltech observations.We compare these observations with the simulated deposition velocities at this site decomposed into its stomatal, cuticle (wet and dry), stem, and ground components.
To facilitate the comparison between simulated and observed deposition velocities, we use meteorological fields (wind speed, temperature, precipitation, and downward radiation) from the Modern-Era Retrospective analysis for Research and Applications (MERRA) (Rienecker et al., 2011) to drive a stand-alone version of LM3-DD.This provides a more accurate representation of the site conditions than using meteorological fields simulated by AM3.
The compounds measured by Nguyen et al. (2015) have different chemical properties, allowing us to evaluate the representation of different deposition pathways in AM3-LM3-DD.In particular, HNO 3 and H 2 O 2 have negligible cuticu-  (Hicks et al., 1987), and R s (Wesely, 1989) for all models.Solar irradiation increases linearly from 0 to 800 W m −2 with temperature.We neglect deposition to the ground and stems.Codep refers to the decrease in R sf,v (SO 2 ) and R sf,v (NH 3 ) associated with base and acid deposition, respectively.For illustrative purposes, the ratio of acid to base deposition is set to 0.5 for SO 2 and 2 for NH 3 .The lifetimes of SO 2 and NH 3 are estimated assuming a boundary layer height of 900 m.GEOS-Chem and AM3-LM3-DD produce identical results for SO 2 under dry conditions.lar resistance (R surf,v 0) (Nguyen et al., 2015), such that v d (X) [R a +R b,v (X)] −1 (ground deposition is negligible).Figure 4 shows that LM3-DD captures both v d (H 2 O 2 ) and v d (HNO 3 ) well, including the faster deposition of H 2 O 2 relative to HNO 3 , consistent with the dependence of R b on 1/D X ∝ √ MW (X) (Eq.5).In contrast, the low solubility and low reactivity at the leaf surface of HCN produces a large non-stomatal resistance (R sf,v > 1 s m −1 , Nguyen et al., 2015), such that v d (HCN) R s (HCN) −1 .A comparison of observed and modeled v d (HCN) suggests that the Leuning model captures the stomatal conductance well at this site.Since R a , R b,v , and R s are well represented over the measurement period, we use observations of v d (ISOPN), v d (MVKN), and v d (PROPNN) at this site to estimate α and β for these organic nitrates (Eq.10).We find that α = 7 and β = 1 provide a reasonable fit for all organic nitrates.These parameters imply that the deposition of isoprene-derived organic nitrates is primarily controlled by dry cuticles with small contributions of stomata and stems.We note that these parameters imply a much greater solubility and reactivity of organic nitrogen than in other models (e.g., α = 0, β = 0.5 in AURAMS; Zhang et al., 2002).While we use these parameters globally, such large differences warrant further investigations, as the deposition of organic nitrogen may account for over 25 % of the overall N deposition but remains rarely measured (Jickells et al., 2013).
Finally, we note that the comparison against SOAS observations points to a significant high bias in simulated nighttime deposition velocity.During this time period, the deposition is dominated by wet cuticles, which reflects the formation of dew in LM3.Since this bias is found for all species including those with little surface resistance (H 2 O 2 and HNO 3 ), it is likely to be associated with an underestimate of the stability of the nocturnal boundary layer.

Impact of land heterogeneities on present-day N deposition
Figure 5 shows the simulated dry deposition of oxidized N (dominated by HNO 3 ) and reduced N (dominated by NH 3 ) as well as the total N deposition (wet+dry) in North America.As noted in previous studies (Zhang et al., 2012;Lamarque et al., 2013), the overall pattern of N deposition mirrors the underlying distribution of NH 3 and NO emissions, with high deposition in the Northeast and a greater contribution of reduced nitrogen to N deposition in the US Midwest and North Carolina than elsewhere in the eastern US.S3.
The grid-cell average dry deposition represents the areaweighted sum of the deposition fluxes to the tiles that comprise each grid cell.Figure 5 (middle column) shows that N deposition over natural vegetation is generally greater than the grid-cell average, which is consistent with faster deposition velocities over forests relative to grasslands (Finkelstein, 2001;Hicks, 2006; and Fig. S1).Overall, the simulated total N deposition to natural ecosystems exceeds the grid-box average deposition by 10 % to 30 % over most of the eastern and central US.This enhancement is largest in regions where land-use change has caused a large decrease in vegetation height and LAI (e.g., in the US Midwest and Northeast, Fig. S2) and smallest in regions with little agricultural activity (e.g., most of Canada) or where managed vegetation differs little in height and LAI from natural vegetation (e.g., in the western US, Fig. S2). Figure 5 (middle column) also shows that the dry deposition of NH x exhibits a greater enhancement over natural vegetation than the dry deposition of NO y , consistent with the greater sensitivity of v d (NH 3 ) than v d (HNO 3 ) to surface properties (Fig. S3).The enhancements of the dry deposition of NH x over natural vegetation is likely to be underestimated in AM3-LM3-DD as the surface bidirectional exchange of NH 3 tends to reduce its deposition in source regions.
Figure 5 (right column) also shows that water bodies receive more reduced N but less oxidized N through dry deposition than the grid-box average.These differences can be attributed to the large effective solubility of NH 3 in freshwater, which results in lower R sf,g (NH 3 ) than over vegetated surfaces (R sf,g (HNO 3 ) is low over all surfaces).Our model suggests that v d (HNO 3 ) is generally slower over water bodies than over vegetated surfaces because of the lower roughness    height of water bodies (see Fig. S3).The westward increase in the ratio of NH 3 to NO emissions thus results in water bodies receiving less N than the average grid cell in the eastern US and Canada but more in the central and western US.

Impact of anthropogenic land-use change on present-day N deposition
Figure 6 shows the change in dry NO y , dry NH x , and total N deposition associated with anthropogenic land-use change, which is estimated by comparing R2010 and R2010_no_lu.We find that anthropogenic land-use change reduces dry NO y , dry NH x , and total N deposition over the contiguous US by 8 %, 26 %, and 6 %, respectively.The reduction in N deposition associated with anthropogenic land-use change is largest in the central and eastern US, where deforestation has caused a large reduction in LAI and vegetation height (Fig. S2).
While anthropogenic land use is estimated to reduce the overall N deposition in the contiguous US, we find that it tends to increase the surface concentration of reactive nitrogen species, which leads to greater N deposition on the remaining natural vegetation.Figure 6 shows that land use has important implications for N deposition at national parks, which are best represented by natural vegetation tiles.For instance, we find that anthropogenic land-use change is as-sociated with a 14 % reduction in the overall N deposition in the region of Shenandoah National Park, but an increase of 9 % on natural vegetation in the same grid box.The slower removal of N near source regions also facilitates N export to remote regions, such as eastern Canada, where N deposition (primarily through wet deposition) increases by more than 10 %.This suggests that anthropogenic land-use change in North America has contributed to the increase of N deposition to natural ecosystems both near source regions and in remote receptor regions.

Implications for future N deposition
Figure 7 shows the simulated difference between N deposition in 2008-2010 (R2010) and 2050 (R2050).This difference reflects changes in anthropogenic emissions as well as changes in climate and land properties induced by climate and land-use change.Total N deposition is projected to increase by 9 % over the contiguous US.Most of the increase is driven by greater deposition in the Midwest and western US associated with higher NH 3 emissions (+40 %).In contrast, N deposition is projected to decrease in the eastern US following the decrease of NO emissions (−47 %, mostly in the eastern US).
We find a small increase (< 10 %) in the deposition velocity of HNO 3 over most of the US between R2010 and R2050 (Fig. S4).This is attributed to a reduction in the land fraction devoted to agriculture in RCP8.5 between 2010 and 2050 (Davies-Barnard et al., 2014), which results in taller vegetation and higher LAI.The impact of this change in land use between 2010 and 2050 is larger for v d (NH 3 ), which increases by more > 10 % over most of the Midwest and eastern US.However, in the eastern and Midwest US, this increase is more than compensated by a reduction in acid deposition, which results in an overall decrease of v d (NH 3 ) of 10 % to 20 % over most of the eastern US.This highlights the need to better characterize the impact of the co-deposition of acids and ammonia on the removal of ammonia to improve projection of future N deposition.
Figure 7 also shows that trends in N deposition simulated for all land types tend to be amplified over natural vegetation, because of the faster deposition velocities as discussed earlier.In contrast, water bodies are projected to experience an increase in N deposition over most of the US, including in regions which experience an overall decrease in N deposition.This contrast is driven by the faster removal of NH 3 over water relative to managed vegetation, which results in greater sensitivity to changes in the emissions of reduced N.The different responses of N deposition on natural tiles and water tiles are important for projections of N deposition in national parks, where N deposition to both vegetation and water bodies is of concern.For instance, the changes in N deposition to natural vegetation from 2010 to 2050 at Voyageurs and Shenandoah national parks are 30 % greater than simulated in the grid box where they are located, while N deposition to water bodies in the Shenandoah region is projected to increase by 16 %, even though overall N deposition for the grid decreases by 18 % in this region.

Conclusions
Our study highlights the importance of accounting for surface heterogeneities and anthropogenic land use in modulating the magnitude and trend of N deposition.Here, we leverage the tiled structure of the GFDL land model to efficiently represent the sub-grid scale heterogeneity of surface properties and their evolution in a changing climate.We have shown that the shift of N emissions from oxidized to reduced N in North America will exacerbate the sensitivity of N deposition to small-scale heterogeneities, which highlights the need to improve the representation of non-stomatal surface resistances (R sf,v , R sf,s , and R sf,g ) including their modulation by canopy wetness and acidity (Flechard et al., 2013;Wentworth et al., 2016;Wu et al., 2018).
Our approach is best suited to long timescales (decadal to centennial) and is complementary to ongoing efforts to improve the representation of present-day N deposition using a combination of high-resolution models and observations (Schwede and Lear, 2014).Future work will aim at coupling the representation of dry deposition presented here to the N cycle in the GFDL land model (Gerber et al., 2010), which will enable us to represent the bidirectional exchange of NH 3 (Nemitz et al., 2001;Flechard et al., 2013;Bash et al., 2013) and improve our understanding of the impact of N deposition on ecosystems and carbon cycling (Magnani et al., 2007;Janssens et al., 2010;Fleischer et al., 2013Fleischer et al., , 2015)).
Data availability.Model outputs are available upon request to Fabien Paulot.Instruction to run the AM3 model are available at https://www.gfdl.noaa.gov/am3/(last access: 14 December 2018).
Author contributions.FP designed the research, developed and implemented the dry deposition model with help from SM, and performed the simulations and analysis.TN and JDC collected the observations at SOAS and estimated the deposition velocities.FP wrote the paper with inputs from SM, TN, JDC, ES, and LWH.
Competing interests.The authors declare that they have no conflict of interest.
Disclaimer.The statements, findings, conclusions, and recommendations are those of the authors and do not necessarily reflect the views of NOAA.

Figure 2 .
Figure 2. Simulated deposition velocity of NH 3 and SO 2 over a coniferous forest (LAI = 5, u = 0.5 m s −1 , RH = 80 %) at different canopy wetness and temperatures.To facilitate comparison across models, we use the same R a = 20 s m −1 , R b(Hicks et al., 1987), and R s(Wesely, 1989) for all models.Solar irradiation increases linearly from 0 to 800 W m −2 with temperature.We neglect deposition to the ground and stems.Codep refers to the decrease in R sf,v (SO 2 ) and R sf,v (NH 3 ) associated with base and acid deposition, respectively.For illustrative purposes, the ratio of acid to base deposition is set to 0.5 for SO 2 and 2 for NH 3 .The lifetimes of SO 2 and NH 3 are estimated assuming a boundary layer height of 900 m.GEOS-Chem and AM3-LM3-DD produce identical results for SO 2 under dry conditions.

Figure 3 .
Figure 3. Observed and simulated deposition velocities of SO 2 for different vegetation types.The symbol shape indicates the canopy status: wet (upward pointing triangle), dry (downward point triangle), and circle (average).The symbol fill indicates the time period: filled (day), half-filled (day+night), and empty (night).The monthly diurnal cycle of deposition velocities simulated by AM3-LM3-DD (R2010 simulation) is sampled at each observation site in the tile that best represents the observed ecosystem accounting for the month, time of day, and canopy wetness status when the observations were collected.References for the different sites are given in TableS3.

Figure 4 .
Figure 4. Observed (circles ± standard deviation) and simulated (bars) dry deposition velocities for several nitrogen-containing species and hydrogen peroxide over the Talladega National Forest (southeastern US) in June 2013 (5 days; Nguyen et al., 2015).The bar colors indicate the contribution of the different surfaces to the overall surface removal of the chemical tracer.

Figure 5 .
Figure 5. Simulated reactive nitrogen deposition (left column) from dry oxidized nitrogen deposition (top row), dry reduced nitrogen deposition (middle row), and total nitrogen deposition (bottom row) over the 2008-2010 period.The ratio between the deposition on selected land types and the grid cell average deposition is shown in the middle (for natural vegetation) and right columns (water bodies).

Figure 6 .
Figure 6.Relative change in the 2008-2010 average land deposition of dry oxidized nitrogen (a), dry reduced nitrogen (b), and total nitrogen (c) associated with anthropogenic land-use change.The relative change is shown as (with land use − without land use) / with land use.From top right to bottom right, the percentages indicate the change in N deposition at Banff National Park (cross), Voyageurs National Park (circle), and Shenandoah National Park (star) at the grid-box level and on natural vegetation, a better proxy for these parks.

Figure 7 .
Figure 7. Simulated change in reactive nitrogen deposition from 2010 to 2050 in the RCP8.5 scenario at the grid-box level (a), on natural tiles (b), and on water bodies (c).From top right to bottom right, the percentages indicate the change in N deposition at Banff National Park (cross), Voyageurs National Park (circle), and Shenandoah National Park (star) for each land type.The fractional change in N deposition over the contiguous US is indicated in the inset (bottom left).