gmd-14-6309-2021.pdf

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

the Creative Commons Attribution 4.0 License.

FORest Canopy Atmosphere Transfer (FORCAsT) 2.0: model

updates and evaluation with observations at a mixed forest site

Dandan Wei

, Hariprasad D. Alwe

, Dylan B. Millet

, Brandon Bottorff

, Michelle Lew

, Philip S. Stevens

Joshua D. Shutter

, Joshua L. Cox

, Frank N. Keutsch

, Qianwen Shi

, Sarah C. Kavassalis

, Jennifer G. Murphy

Krystal T. Vasquez

, Hannah M. Allen

, Eric Praske

, John D. Crounse

, Paul O. Wennberg

, Paul B. Shepson

10,11

Alexander A. T. Bui

, Henry W. Wallace

, Robert J. Griffin

, Nathaniel W. May

, Megan Connor

Jonathan H. Slade

, Kerri A. Pratt

, Ezra C. Wood

, Mathew Rollings

15,a

, Benjamin L. Deming

15,b

Daniel C. Anderson

14,c

, and Allison L. Steiner

Department of Climate and Space Sciences and Engineering, University of Michigan, Ann Arbor, MI, USA

Department of Soil, Water and Climate, University of Minnesota, Twin Cities, St. Paul, MN, USA

Department of Chemistry, Indiana University, Bloomington, IN, USA

School of Public and Environmental Affairs, Indiana University, Bloomington, IN, USA

Department of Chemistry and Chemical Biology, Harvard University, Cambridge, MA, USA

Department of Chemistry, University of Toronto, Toronto, Ontario, Canada

Division of Chemistry and Chemical Engineering, California Institute of Technology, Pasadena, CA, USA

Division of Geological and Planetary Sciences, California Institute of Technology, Pasadena, CA, USA

Divisions of Engineering and Applied Science and Geological and Planetary Science,

California Institute of Technology, Pasadena, CA, USA

Department of Chemistry, Purdue University, West Lafayette, IN, USA

School of Marine and Atmospheric Sciences, Stony Brook University, Stony Brook, NY, USA

Department of Civil and Environmental Engineering, Rice University, Houston, TX, USA

Department of Chemistry, University of Michigan, Ann Arbor, MI, USA

Department of Chemistry, Drexel University, Philadelphia, PA, USA

Department of Chemistry, University of Massachusetts, Amherst, MA, USA

now at: Department of Chemistry, University of California, Berkeley, CA, USA

now at: Department of Chemistry and Cooperative Institute for Research in Environmental Sciences,

University of Colorado, Boulder, CO, USA

now at: Universities Space Research Association, Columbia, MD and NASA Goddard Space Flight Center,

Greenbelt, MD, USA

Correspondence:

Dandan Wei (dandanwe@umich.edu)

Received: 31 March 2021 – Discussion started: 23 April 2021

Revised: 27 July 2021 – Accepted: 20 September 2021 – Published: 21 October 2021

Abstract.

The FORCAsT (FORest Canopy Atmosphere

Transfer) model version 1.0 is updated to FORCAsT 2.0 by

implementing five major changes, including (1) a change

to the operator splitting, separating chemistry from emis-

sion and dry deposition, which reduces the run time of the

gas-phase chemistry by 70 % and produces a more real-

istic in-canopy profile for isoprene; (2) a modification of

the eddy diffusivity parameterization to produce greater and

more realistic vertical mixing in the boundary layer, which

ameliorates the unrealistic simulated end-of-day peaks in

isoprene under well-mixed conditions and improves day-

time air temperature; (3) updates to dry deposition veloci-

ties with available measurements; (4) implementation of the

Reduced Caltech Isoprene Mechanism (RCIM) to reflect the

current knowledge of isoprene oxidation; and (5) extension

of the aerosol module to include isoprene-derived secondary

Published by Copernicus Publications on behalf of the European Geosciences Union.

6310

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

organic aerosol (iSOA) formation. Along with the opera-

tor splitting, modified vertical mixing, and dry deposition,

RCIM improves the estimation of first-generation isoprene

oxidation products (methyl vinyl ketone and methacrolein)

and some second-generation products (such as isoprene

epoxydiols). Inclusion of isoprene in the aerosol module in

FORCAsT 2.0 leads to a 7 % mass yield of iSOA. The most

important iSOA precursors are IEPOX and tetrafunctionals,

which together account for

86 % of total iSOA. The iSOA

formed from organic nitrates is more important in the canopy,

accounting for 11 % of the total iSOA. The tetrafunctionals

compose up to 23 % of the total iSOA formation, highlight-

ing the importance of the fate (i.e., dry deposition and gas-

phase chemistry) of later-generation isoprene oxidation prod-

ucts in estimating iSOA formation.

1 Introduction

Forests cover 30 % of the land surface and play an impor-

tant role in the Earth system through exchanges of energy,

water, carbon dioxide, and reactive chemical species with

the atmosphere (Bonan, 2008). Forest canopies emit large

amounts of volatile organic compounds (VOCs) into the

atmosphere (Guenther et al., 2006) that drive atmospheric

chemistry (e.g., Chameides et al., 1992; Taraborrelli et al.,

2012) and are precursors to climate-relevant species such as

ozone (e.g., Wolfe et al., 2011) and particulate matter (e.g.,

Palm et al., 2016). In addition, forest canopies serve as a

major sink of VOCs through dry deposition (e.g., Nguyen

et al., 2015). These bidirectional exchanges and their in-

fluences on atmospheric chemistry are complicated by the

three-dimensional structure of forest canopies, which cre-

ates turbulent flows significantly different from the overlying

atmospheric boundary layer (e.g., Gao et al., 1993; Patton

et al., 2001) and affects the vertical transport and the chem-

istry of trace gases (e.g., Kaser et al., 2015).

The complexity and interplay of these chemical and phys-

ical processes challenge our understanding of forest-driven

climate impacts on local, regional, and global scales (e.g.,

Spracklen et al., 2008; Vanwalleghem and Meentemeyer,

2009). Improving our understanding of the chemical and

physical processes governing forest–atmosphere interactions

at a local scale is helpful to generalize the net impact of

the terrestrial biosphere on chemistry and climate at broader

scales. Canopy–chemistry models that explicitly represent

the physical, chemical, and biological processes of an indi-

vidual forest canopy are useful tools to investigate the chem-

ically relevant interactions between forests and the atmo-

sphere at local scales. These canopy–chemistry models cal-

culate the environmental variables that drive emissions, dry

deposition, turbulent mixing, and chemical reactions verti-

cally throughout a canopy at very fine resolutions (e.g., on

the order of meters), while atmospheric chemical transport

models approximate canopy processes through parameteri-

zations and operate on the 10–200 km scale. Several forest

canopy–chemistry models (Stroud, 2005; Forkel et al., 2006;

Boy et al., 2011; Wolfe and Thornton, 2011) have been de-

veloped to study chemically relevant forest–atmosphere ex-

changes with a focus on the gas-phase chemical processes.

The FORCAsT model version 1.0 (FORest Canopy Atmo-

sphere Transfer) (Ashworth et al., 2015) is one of the few

canopy models currently capable of simulating the forma-

tion of secondary organic aerosol (SOA) from biogenic VOC

oxidation.

Over the past decade, the understanding of isoprene chem-

istry under a wide range of NO

conditions and their im-

pact on atmospheric particles has greatly expanded. Specif-

ically, understanding of isoprene oxidation under low-NO

conditions has improved (Wennberg et al., 2018), and proper

representation of isoprene oxidation and isoprene-derived

SOA formation in canopy–chemistry models is now recog-

nized to be important for a more accurate understanding of

forest–atmosphere exchange. Isoprene-hydroxy-peroxy radi-

cals (ISOPOO), produced by addition of a hydroxyl radical

(OH) across one of the double bonds followed by the rapid

addition of molecular oxygen (O

), react with nitric oxide

(NO), hydroperoxy radicals (HO

), and other peroxy radi-

cals or undergo unimolecular isomerization. Historically, the

dominant fate of ISOPOO was thought to be reaction with

NO, as mechanisms were developed for urban locations and

the NO loss pathway dominates in polluted regions. Under

low-NO

conditions common in forested regions, unimolec-

ular chemistry and reaction with HO

are also important.

The first-generation product of the reaction of ISOPOO with

, hydroxy hydroperoxide (ISOPOOH), is an important

SOA precursor following their oxidation to epoxydiol prod-

ucts (Surratt et al., 2006; Paulot et al., 2009). The other novel

low-NO

pathway recently elucidated is ISOPOO isomer-

ization, which can sustain elevated OH concentrations un-

der low-NO conditions (Peeters and Müller, 2010; Crounse

et al., 2011; Teng et al., 2017; Møller et al., 2019). These

different branches of ISOPOO pathways produce a different

ensemble of oxygenated compounds with low volatility and

are thus crucial for accurate prediction of the environmental

and climate impacts of isoprene chemistry. In addition to the

gas-phase fate of isoprene, field studies found evidence of C5

compounds in ambient particles (Claeys, 2004; Kleindienst

et al., 2007), and the modeling of isoprene-derived SOA has

been significantly advanced in the past decade (e.g., Comper-

nolle et al., 2009; Marais et al., 2016; Gaston et al., 2014).

Wennberg et al. (2018) compile a comprehensive isoprene

mechanism incorporating the current knowledge of isoprene

chemistry and include the necessary isoprene SOA precur-

sors. The Reduced Caltech Isoprene Mechanism (RCIM), a

reduced version with the same product yields of known com-

pounds and minimal simplifications beyond lumping of iso-

meric compounds and removal of minor (

2 % yield) path-

ways, is developed concurrently with the explicit mechanism

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6311

and is more suitable for implementing in canopy–chemistry

models.

Here, we develop FORCAsT 1.0 to FORCAsT 2.0 to in-

corporate important updates to gas-phase isoprene chem-

istry, isoprene-derived SOA formation, and two physical

components of the model. Specifically, major updates in-

clude (i) separating the integration of chemistry from the

emission and dry deposition (also called operator splitting)

to provide more realistic representations of vertical gradi-

ents in the forest canopy and to make the chemical mod-

ule more flexible with future chemical mechanism updates,

(ii) improving the vertical mixing parameterization in the

boundary layer, (iii) updating the dry deposition velocities

for chemical species with available measurements (Nguyen

et al., 2015), (iv) implementing the RCIM to reflect the cur-

rent understanding of isoprene fate under low-NO

condi-

tions (Wennberg et al., 2018), and (v) extending the MPMPO

(Model to Predict the Multiphase Partitioning of Organics)

aerosol module (Griffin et al., 2005; Ashworth et al., 2015) to

include isoprene-derived SOA formation. We evaluate FOR-

CAsT 2.0’s performance against FORCAsT 1.0 (Ashworth

et al., 2015) and the observations from the AMOS (Atmo-

spheric Measurements of Oxidants in Summer) field cam-

paign conducted at the University of Michigan Biological

Station (UMBS) during the summer of 2016.

2 Model description

The FORCAsT model, based on the CACHE canopy model

by Forkel et al. (2006), is a one-dimensional model that cou-

ples atmospheric chemistry (gas-phase and gas–particle par-

titioning) and canopy processes. The vertical resolution of

FORCAsT can be configured with a minimum of 20 and a

maximum of 60 vertical layers. In this study, the number of

model levels is set to 40, with 20 layers within the canopy, 8

layers representing the boundary layer (

∼

1 km), and the re-

maining layers extending to the lower troposphere (

∼

4 km).

In addition to the aboveground layers, the model includes

15 soil layers for computing soil heat and moisture storage,

as well as exchange with the atmosphere and root extraction

(Forkel et al., 2006).

The FORCAsT model parameterizes the processes of ra-

diative transfer, chemical species emission, advection, depo-

sition, vertical exchange, and chemistry and then integrates

the energy and mass balance equations at each vertical layer

in the canopy. The parameterization of the radiative trans-

fer, emission, advection, and turbulent mixing in the canopy

layers has been described extensively in Bryan et al. (2012)

and Ashworth et al. (2015) and remains unchanged in the up-

dated version 2.0. Here we describe the major updates to the

model, including the operator splitting, gas-phase chemistry,

gas–particle partitioning, and some aspects of the dry depo-

sition and turbulent mixing.

2.1 Operator splitting

Processes such as emission, turbulent mixing, dry deposition,

and chemical reactions occur at the same time in the atmo-

sphere. However, in numerical models, these processes (i.e.,

operators) are split and integrated over time and/or space in

sequence, which is commonly referred to as operator split-

ting. It is generally faster to integrate the operators sepa-

rately than to compute the solution when the operators are

treated together (Lapointe et al., 2020). This computational

efficiency comes at the cost of an error introduced by the

splitting. In the context of atmospheric chemistry modeling,

model accuracy is affected by the order in which the opera-

tors are applied (Santillana et al., 2016) and by the integra-

tion time steps of the operators (also called operator duration)

(Philip et al., 2016). In prior studies of chemical transport

models, the operator duration causes greater differences in

concentrations of reactive emitted species such as nitrogen

oxides (up to 5 times; Philip et al., 2016) than the order of

operators (up to 10 %; Santillana et al., 2016).

In FORCAsT 1.0, the order of operation is from vertical

transport to chemistry, with emissions and dry deposition in-

tegrated within the chemistry solver (Table 1). The chemistry

solver typically dominates the computational cost of the sim-

ulations (Lapointe et al., 2020). To increase computational

efficiency and to allow for flexible chemical mechanisms in

the future, we separate the chemical solver from emission

and dry deposition and integrate the latter two operators with

the vertical transport (Table 1). The impacts of this opera-

tor splitting on FORCAsT 2.0’s performance are discussed

in Sect. 3.1.

2.2 Turbulent mixing

In the surface layer (roughly 10 % of the boundary layer), the

eddy diffusivity

is commonly defined as a simple function

of height

of the form

κu

∗

, where

is the von Kar-

man constant and

∗

is the friction velocity (Stull, 1988).

In the rest of the boundary layer, defining

is not as clear

for numerical models. Several approaches have been used in

the literature to define the

profile, such as linearly de-

creasing from the surface layer to the top of the boundary

layer (Estoque, 1963), exponentially decreasing with height

from the surface layer, and finding an interpolating poly-

nomial passing through prescribed points with predefined

slopes (O’Brien, 1970). The general approach for approxi-

mating

in the boundary layer has used a power-law de-

pendence on

z/z

, where

is the top of the boundary layer

and scale parameters are derived from similarity theory or

empirically. A commonly used shape function is of the form

(Troen and Mahrt, 1986)

κu

∗

(

−

)

(1)

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6312

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

Table 1.

The order of operation in the FORCAsT model.

∂C

∂t

denotes the time evolution of concentrations of chemical species

FORCAsT 1.0

FORCAsT 2.0

∂C

∂t

transport

∂C

∂t

transport

emission

dry deposition

⇓

∂C

∂t

chemistry

emission

dry deposition

∂C

∂t

chemistry

The exponent

2 has been commonly used in the litera-

ture (Nissanka et al., 2018), while values between 2 and 3

agree with different observed profiles.

The eddy diffusivity (

) parameterization based on

mixing-length theory in FORCAsT 1.0 is described as fol-

lows, with greater detail presented in Forkel et al. (2006) and

Bryan et al. (2012):

∂u

∂z

(2)

κz

(3)

∂u

∂z

∗

κz

(4)

where

is the mixing length,

is the mean wind speed,

is the von Karman constant (0.41),

is the height above the

ground, and

∗

is the friction velocity. The

∂u

∂z

is derived from

the common logarithmic expression for the boundary layer

(Eq. 4). Combing Eqs. (2), (3), and (4) and rearranging the

terms leads to

κu

∗

(

κz

)

(5)

where

is a function of the height

as shown below.











z < h

max

(

−

), h

−

(6)

Here,

is the canopy height,

is the height of the bound-

ary layer,

is the geostrophic wind at the top of the bound-

ary layer set to 17 m s

−

, and

is the Coriolis parameter.

This parameterization yields a relatively small value of

(i.e.,

10 m

−

) in the boundary layer and thus implicitly

produces a low boundary layer height (

∼

250–300 m around

14:00 local time) in FORCAsT 1.0. In addition, this mixing

parameterization in FORCAsT 1.0 produces an unrealistic

end-of-day peak in isoprene under well-mixed conditions.

To produce more realistic

in the boundary layer and

boundary layer height, we adopt the parameterization de-

scribed in Eq. (1) to calculate the

in the boundary layer.

In the present study, we use a cubic power of height

(i.e.,

2) for

in the boundary layer:

new

{

κu

∗

φ(z/L)

z < z

sfc

aκu

∗

−

)

φ(z/L)

, z

sfc

(7)

where

sfc

is the surface layer height (here assumed to be

10 % of the boundary layer height

), and the constant

−

sfc

)

23 is used to ensure a continuous transition

new

from the surface layer to the boundary layer. The

is calculated as a function of the sensible heat flux at the top

of the canopy (Stull, 1988):

√

(

)

∫

′

)

cpy

(8)

where

is a standard entrainment parameter (0.2), and

the lapse rate in the free atmosphere (0.0065 K m

−

). The

integral

∫

′

)

cpy

is approximated to

′

)

cpy

, where

is the elapsed time. The calculation for

starts when the

heat flux

′

)

cpy

first becomes positive in the morning and

continues 1 h longer after

′

)

cpy

becomes negative. The

nighttime

is set to 200 m.

The stability function

φ(z/L)

is given in Eq. (9) (Nissanka

et al., 2018).

φ(z/L)







(

−

z/L)

−

, z/L <

z/L,

z/L >

z/L

(9)

Overall, the new parameterization produces a larger

, a

more realistic boundary layer height, and thus a more real-

istic air temperature (see details in Sect. 3.2).

2.3 Dry deposition

The dry deposition velocity to the canopy foliage (

) is cal-

culated using the resistance model (Meyers and Baldocchi,

1988; Wesely, 1989) in FORCAsT 1.0, with details provided

in Bryan et al. (2012) and Ashworth et al. (2015). Recent

work by Nguyen et al. (2015) reports dry deposition veloci-

ties based on measured fluxes and concentrations for 16 at-

mospheric compounds above a southeastern United States

forest, suggesting that the parameterization in the FORCAsT

1.0 underestimates dry deposition velocities for these oxy-

genates (Table 2). We adopt this newer parameterization in

FORCAsT 2.0, wherein the major revisions to the resistance

model include (1) the addition of the aerodynamic resis-

tance

, (2) the formulation of the molecular diffusion

and (3) the addition of temperature dependence to the mes-

ophyll resistance

and the cuticular resistance

. Details

of the revised parameterization can be found in Nguyen et al.

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6313

Table 2.

Estimates of the dry deposition velocities (cm s

−

) to

the canopy foliage for relevant molecules at 12:00 local time. The

measurement-based deposition velocities for a temperate forest are

provided as a reference (Nguyen et al., 2015).

Species

Revised

Old

Measurement-

resistance

based

model

5.2

0.2

HNO

4.6

7.7

HCOOH

1.6

3.3

HAC

0.7

0.1

HMHP

4.7

1.0

IHN

1.3

0.0

HPALD

3.2

−

ISOPOOH/IEPOX

3.4

1.0

PROPNN

2.0

INP

1.3

NA – not available.

(2015). The revised resistance model is highly sensitive to

Henry’s law coefficient (Nguyen et al., 2015), which is un-

known for many species. Here, we apply the revised resis-

tance model to species whose Henry’s law coefficients are

available and measurement-based

exists to validate the

estimates, including for H

(hydrogen peroxide), HNO

(nitric acid), HCOOH (formic acid), HAC (hydroxyacetone),

PROPNN (propanone nitrate), HMHP (hydroxymethyl hy-

droperoxide), INP (isoprene nitrooxy hydroperoxide), IHNs

(isoprene hydroxy nitrates), HPALD (hydroperoxy aldehy-

des), IEPOX (epoxydiols), and ISOPOOH (isoprene hy-

droxy hydroperoxides) (Table 2). The revised deposition ve-

locity parameterization is also applied to all the isomers

for INP, IHN, ISOPOOH, IEPOX, and HPALD in RCIM.

Note that the revised dry deposition scheme yields a much

smaller

(0.35 cm s

−

) for methyl vinyl ketone (MVK) and

methacrolein (MACR) than the observation-based estimate

(up to 2.4 cm s

−

; Karl et al., 2010). Therefore, we prescribe

the observation-based 2.4 cm s

−

for MVK and MACR in

FORCAsT. Because the resistance model assumes that each

canopy layer is one “big leaf” perpendicular to the sunlight,

the estimates of dry deposition velocities for all compounds

are then scaled by the leaf area distribution in the canopy.

Species other than those listed above use the old parame-

terization in FORCAsT 1.0. The comparison of dry depo-

sition velocities between the old and new parameterizations

is shown in Table 2, with the new parameterization generally

increasing the deposition velocity with the exception of nitric

and formic acids. As input data and measurements become

available for other species, the revised parameterization can

be evaluated against new observations and applied to other

species.

2.4 Isoprene gas-phase mechanism: the Reduced

Caltech Isoprene Mechanism (RCIM)

We replace the chemical reactions for isoprene oxidation in

the Caltech Atmospheric Chemistry Mechanism (CACM) in

FORCAsT 1.0 with RCIM. The original CACM was devel-

oped for application to urban conditions (CACM0.0) (Griffin

et al., 2002, 2005), and the version incorporated into FOR-

CAsT 1.0 was updated for low-NO

conditions (Ashworth

et al., 2015), hereinafter referred to as CACM1.0 in the fol-

lowing sections. RCIM is version 4.3 of the “reduced” mech-

anism in the Wennberg et al. (2018) mechanism repository

(Bates and Wennberg, 2017), including the essential chem-

istry required to accurately simulate isoprene oxidation under

remote conditions in the atmosphere. Compiled concurrently

with the full explicit mechanism, RCIM groups isomers with

similar reaction rates and products and lumps minor path-

ways (

2 % branching ratio) into the major channels to min-

imize the number of species and reactions while retaining an

accurate description of the oxidative fate of those grouped

species and lumped pathways (Wennberg et al., 2018). The

reduced mechanism includes 119 species and 221 reactions,

in contrast to the 385 species and 810 reactions in the explicit

isoprene mechanism in Wennberg et al. (2018) and the 113

reactions for isoprene in CACM1.0.

RCIM treats the initial system of allylic and peroxy radi-

cals formed following the addition of OH to isoprene dynam-

ically. Older mechanisms, including CACM1.0 in FORCAsT

1.0, implicitly used fixed distributions of isoprene-hydroxy-

peroxy radicals (ISOPOO) derived from experiments per-

formed under high-NO conditions. Addition of O

to allylic

radicals under ambient conditions is in fact a reversible pro-

cess, resulting in a dynamic system with interconversion be-

tween the six major ISOPOO isomers (two subgroups of

three defined by a common OH position) (Peeters et al.,

2009; Teng et al., 2017). Wennberg et al. (2018) represent

the reversibility of O

addition in the explicit mechanism by

including 10 species (i.e., four allylic radicals and six major

ISOPOOs; see Fig. 3 in Wennberg et al., 2018) and 69 re-

actions. The reduced mechanism, RCIM, retains this novel

treatment of the ISOPOO system (Wennberg et al., 2018), al-

though it simplifies the 10-species radical system to two ma-

jor ISOPOO isomers, i.e.,

(

−

)

−

ISOPOO and

(

−

)

−

ISOPOO.

RCIM includes important updates to the formation and

fates of isoprene hydroxyl nitrates (IHNs) through pressure-

and temperature-dependent parameterizations of the branch-

ing ratios and a new structure–activity relationship for cal-

culating the formation of nitrates from multifunctional per-

oxy radicals without measured yields (Wennberg et al.,

2018). The dynamic representation of the ISOPOO isomers

also contributes to higher production of IHN than previ-

ous mechanisms when included in global models (Bates

and Jacob, 2019). In addition, RCIM includes 12 dis-

tinct C5 tetrafunctional compounds with unique combina-

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6314

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

tions of functional groups. They are dihydroxy hydroper-

oxy nitrate, carbonyl hydroxy nitrooxy hydroperoxide, di-

hydroxy carbonyl nitrate, hydroxy hydroperoxy dialdehyde,

hydroxy hydroperoxy dinitrate, dihydroxy dinitrate, carbonyl

hydroperoxy-diol, carbonyl hydroxy dinitrate, dihydroxy hy-

droperxy epoxide, dihydroxy dihydroperoxide, hydroxy ni-

trooxy dihydroperoxide, and hydroxy nitrooxy hydroper-

oxy epoxide. Each C5 tetrafunctional compound represents

a variety of isomers. Chamber experiments suggest that

multifunctionals contribute to isoprene-derived SOA (iSOA)

(Ng et al., 2008; Schwantes et al., 2019). Additional as-

pects of RCIM relative to CACM1.0 include decreased C5-

hydroperoxy-aldehyde (HPALD) yields following the 1,6-

H shifts of the Z

−

peroxy radicals in Teng et al.

(2017) and additional intramolecular H shifts, including

rapid peroxy–hydroperoxy shifts, resulting in higher OH re-

cycling under low-NO conditions (Wennberg et al., 2018).

These mechanism changes are manifested by changes in gas-

phase isoprene oxidation products (Sect. 3.3).

2.5 Isoprene-derived secondary organic aerosol

FORCAsT simulates the partitioning of condensable species

into the particle phase using the Model to Predict the Mul-

tiphase Partitioning of Organics (MPMPO; Griffin et al.,

2005). In FORCAsT 1.0, 99 out of the 300 prognostic species

in CACM1.0 are treated as condensable and lumped into

12 surrogate species according to their structures, sources

(biogenic or anthropogenic), volatilities, and dissociative ca-

pabilities (Ashworth et al., 2015). Specifically, the 12 sur-

rogate groups include 6 anthropogenic aromatic groups, 4

monoterpene-derived biogenic surrogates, and 1 group com-

posed of non-volatile dimers of multifunctional acids from

monoterpene oxidation (i.e., phthalic acid). Additionally,

one surrogate group based on keto-propanoic acid and ox-

alic acid (formed from oxidation of methyl vinyl ketone

and methacrolein) is considered condensable in CACM1.0-

MPMPO. However, explicit formation of SOA from isoprene

is missing in FORCAsT 1.0. We incorporate six new surro-

gate groups in MPMPO to account for the isoprene-derived

SOA (iSOA) in FORCAsT 2.0 (Table 3).

As the iSOA precursors are small organic molecules (num-

ber of carbon atoms

≤

5) with numerous functional groups,

they are expected to be highly hydrophilic. Under humid con-

ditions representative of the summertime PROPHET (Pro-

gram for Research on Oxidants: Photo-chemistry, Emis-

sions, Transport) boundary layer, aqueous aerosol provides

a medium for reactive uptake (Surratt et al., 2009), and thus

iSOA is likely aqueous (Couvidat and Seigneur, 2011; Er-

vens et al., 2011). Marais et al. (2016) proposed a mecha-

nism for irreversible reactive uptake of iSOA precursors by

preexisting aqueous aerosols, dependent on the Henry’s law

coefficients, that is widely used in chemical transport mod-

els. Bates and Jacob (2019) estimate a global iSOA yield by

mass of 25 % using this mechanism coupled with RCIM. On

the other hand, dry chamber experiments (relative humid-

ity

10 %) suggest up to 11 % of iSOA yield (Kroll et al.,

2006), suggesting a partitioning between the gas phase and a

non-aqueous organic phase. Thornton et al. (2020) estimate

an upper bound of the non-aqueous iSOA yield of 3 % under

atmospheric conditions using a volatility-driven gas–particle

partitioning implemented in a box model.

The MPMPO in FORCAsT only considers organic aerosol

species and assumes the partitioning of the gases into two

particulate phases: a purely organic aerosol and an aque-

ous aerosol with associated molecular and ionic components.

Equilibrium between the gas and aerosols is assumed for or-

ganic species. The equilibrium organic aerosol-phase mass

concentration of an individual species

(μg m

−

air), is

given by the following relationship (Griffin et al., 2005):

L,i

(10)

where

air μg

−

) is the partitioning coefficient

that describes the phase distribution of the condensing or-

ganic species (Pankow, 1994),

is its corresponding gas-

phase concentration (μg m

−

air), and

is the total con-

centration (μg m

−

air) of organic aerosol mass available to

act as the partitioning medium.

is the ideal gas constant

−

atm mol

−

is temperature (K),

is the average molecular weight (g mol

−

) of the absorbing

organics including both primary and secondary organic com-

pounds,

L,i

is the pure component vapor pressure (atm) of

species

, and

is the activity coefficient of species

in the

organic phase. The factor of 10

converts grams (g) to micro-

grams (μg).

The aqueous-phase concentration of species

(μg m

−

air), is given by (Griffin et al., 2005):

LWC

(11)

where

LWC

the

aerosol

liquid

water

content

(μg H

O m

−

air),

is the Henry’s law coefficient

(

μg μg

−

)/(

μg m

−

air

)

), and

is the activity coef-

ficient of organic species

in the aqueous phase normalized

by that at infinite dilution. The aqueous-phase equilibrium is

constrained by dissociation of the dissolved organic species,

for which pH is needed to calculate the concentration of the

charged ions (see Eqs. 13 and 14 in Ashworth et al., 2015).

The LWC in Eq. (11) is calculated offline as the sum

of aerosol liquid water associated with inorganic species

(LWC

inorg

) and organics (LWC

org

). The LWC

inorg

and pH are

calculated by the Extended AIM Aerosol Thermodynamics

Model model II (E-AIM, http://www.aim.env.uea.ac.uk/aim/

model2/model2d.php, last access: 14 February 2021) (Clegg

et al., 1998) using measurements of NH

(g) and PM

sul-

fate, nitrate, and ammonium (Table 4), following the method

of Murphy et al. (2017). The LWC

org

is calculated according

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6315

Table 3.

Properties of the six new iSOA surrogate groups in FORCAsT 2.0. The Henry’s law constants below are from Sander (2015) and

Safieddine et al. (2017) at 298 K. Temperature dependence of Henry’s law constants is included in FORCAsT.

Surrogate species

Representative

Henry’s law coefficient

Surrogate group

molar mass

(at 298 K; M atm

−

)

group name

(g mol

−

)

GLYX (C

), MGLY (C

)

72.0

GLYX

IEPOXt, IEPOXc, IEPOXD (C

), ICHE (C

118.0

IEPOX

HMML (C

)

ISOP1N4OH, ISOP1OH4N

147.0

IHN

ISOP1OH2N, ISOP3N4OH (C

)

INPB, INPD (C

)

163.0

INP

MVK3OOH4N, MACR2OOH3N (C

)

149.0

MVK3N4OH, MVK3OH4N (C

)

MACR2OH3N, MACR2N3OH (C

)

C5 tetrafunctional compounds

208.0

Tetra

to Petters and Kreidenweis (2007):

LWC

org

−

(12)

where

org

is the organic mass concentration from high-

resolution time-of-flight aerosol mass spectrometer (HR-

ToF-AMS) measurements,

org

is the density of organics

(1.4 g cm

−

; Cerully et al., 2015),

is the density of wa-

ter, RH is the relative humidity, and

org

is the hygroscopic-

ity growth for organics taken as 0.08 in this study. The

org

is derived from the HR-ToF-AMS data using the method de-

scribed in Cerully et al. (2015):

org

inorg

(13)

where

org

and

inorg

are the volume fractions of organic and

inorganic species calculated using AMS mass fraction data,

and

inorg

is the hygroscopicity growth determined from the

speciated inorganic concentrations and

inorg

for individual

inorganic species from Padró et al. (2010).

is the total hy-

groscopicity growth and has been shown to be insensitive to

location and organic fraction (Padró et al., 2010). Here we

use a value of 0.23 for

based on measurements from pre-

vious studies (Padró et al., 2010; Cerully et al., 2015). We

calculate LWC as the sum of LWC

inorg

and LWC

org

based on

surface observations and vertically scale LWC by the relative

humidity to extrapolate to other model heights.

The thermodynamic model UNIFAC (UNIversal Func-

tional group Activity Coefficient) is used to calculate the ac-

tivity coefficients

in Eqs. (10) and (11) (Fredenslund et al.,

1975). The standard UNIFAC parameters (e.g., alkane group)

are found in Hansen et al. (1991), Balslev and Abildskov

(2002), and Wittig et al. (2003). The UNIFAC parameters for

the functional groups of nitrate and hydroperoxide are taken

from Compernolle et al. (2009). The missing UNIFAC pa-

rameters (e.g., the unknown interaction parameters) for some

functional groups (such as nitrate and hydroperoxide) are set

to zero, which introduces uncertainties in estimating the ac-

tivity coefficient for the new SOA surrogate groups.

3 Model evaluation with observations

The performance of FORCAsT 2.0 is evaluated against the

observations obtained at the PROPHET (Program for Re-

search on Oxidants: Photo-chemistry, Emissions, Transport)

tower at UMBS during the 2016 AMOS field campaign. Full

details of the PROPHET tower and the UMBS site can be

found in Millet et al. (2018). Measurement details and refer-

ences are listed in Table 4. The results presented in this study

are based on a 2 d model simulation for the 2 sunny days of

22–23 July 2016; the first day (22 July) is a well-mixed day

and the second day (23 July) is relatively stagnant based on

micrometeorological analysis (Wei et al., 2020). The model

spin-up time is 24 h. The 2 d period is relatively hot for the

site, with a mean high of 30.0

◦

C and a mean daily temper-

ature of 23.9

◦

C compared to the monthly averages of 25.4

and 21.0

◦

C, respectively. The canopy structure in the model

is identical to that used in Bryan et al. (2012) and Ashworth

et al. (2015). The input data are based on the measurements

during AMOS 2016, with the model driven by observed PAR

(photosynthetically active radiation), standard deviation of

the vertical velocity (

), friction velocity (

∗

), and the cal-

culated aerosol liquid water content. The initial concentra-

tions for the chemical species are taken from the measure-

ments when available, including ozone (O

), nitric oxide

(NO), nitrogen dioxide (NO

), formaldehyde (CH

O), and

methyl vinyl ketone and methacrolein (MVK

MACR). The

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6316

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

isoprene emission factor is increased from prior FORCAsT

1.0 studies (Bryan et al., 2012; Ashworth et al., 2015), as the

last measurement campaign (i.e., CABINEX 2009) occurred

during a relatively cool period and utilized an emission factor

roughly a factor of 2–3 lower than previously observed at the

site (Pressley, 2005; Unger, 2013). The 2 d period evaluated

in this study is warmer, and increasing the isoprene emission

factor is consistent with the effects of prior day temperature

(Guenther et al., 2006). The gas-phase chemical mechanism

in FORCAsT includes explicit treatment of two monoter-

pene surrogate species:

-pinene (APIN) and

-limonene

(DLMN). Light dependence of monoterpene emissions is in-

cluded in the 2 d simulation by changing the synthesis emis-

sion factors from 0 to 0.4 nmol m

−

for the two monoter-

pene surrogates, as this has been observed by Ortega et al.

(2007) at the site and this emission change improved the

agreement with the measured monoterpene concentrations

for a more realistic evaluation of the monoterpene-derived

SOA. In the following sections, we describe the impact of

model changes (Sect. 2.0) on the physical and chemical en-

vironment in the forest canopy during this 2 d period.

3.1 Operator splitting

As described in Sect. 2.1, the chemistry is separated from the

operators of emission and dry deposition in FORCAsT 2.0

(Table 1). The splitting of chemistry from the emission and

dry deposition leads to weaker gradients in the vertical pro-

files for emitted species such as isoprene (Fig. 1). The old

order of operators simulates higher concentrations between

0.7 and 1.0

z/h

where the emission occurs (Fig. 1a–c), likely

due to the higher production rates of isoprene resulting from

emission (on the order of 10

−

–10

−

ppbv s

−

depending on

the leaf area distribution) compared to the reaction rates (on

the order of 10

−

ppbv s

−

for isoprene

OH) in the same

solver, resulting in the vertical concentration profiles resem-

bling the emission profiles. The new order of operators sim-

ulates more realistic in-canopy concentration gradients for

the emitted species such as isoprene. Specifically, the verti-

cal profile of isoprene before sunrise is better captured by the

new order of operation with an RMSE (root mean square er-

ror) less than 0.1 (Fig. 1a, d). The midday in-canopy profiles

are also improved by the new order, while simulated concen-

trations are higher above the canopy than observed (Fig. 1b,

e). The overall RMSE for the midday case (Fig. 1e) suggests

that the results for the new order agree slightly better with

the observations. The differences between the two orders are

up to 16 % around the height of 0.85

z/h

in the midday case

and are in good agreement with Santillana et al. (2016), who

report a 10 % difference. In the evening hours, the emission

rates and chemistry quickly decrease due to the reduced ra-

diation, leading to more fluctuations in the observed profile

that are challenging for models to capture (Fig. 1c).

In addition, the time for running the model with the gas-

phase chemistry alone (i.e., without running the aerosol mod-

ule) is reduced from 10 to 3 min. In chemical models, the

computational cost is mainly due to the chemical solver,

which has a very small internal time step (about 1.2 s in our

case). However, the transport solver has a time step of 1 min;

therefore, moving the emission and dry deposition from the

chemical solver to the transport solver reduces the run time

for the gas-phase chemistry by 70 %, making the model more

computationally efficient for adding increasingly complex

chemical mechanisms in the future. In summary, the differ-

ences of the modeled gradients in isoprene concentrations

between the two operator orders are relatively small (up to

16 %), with the new order having a more realistic in-canopy

profile and a higher computational efficiency.

The impacts of the operator splitting on the gas concen-

trations are correlated with the chemical lifetime (Fig. 2).

The OH and nitric oxide (NO) concentrations increase by

160 %–180 % and 130 % with the new order of operator, re-

spectively, while CO (carbon monoxide) only differs by 13 %

(Fig. 2). The increased NO, CH

O (formaldehyde), and O

using the new order improve agreement with the observations

(Fig. 4). In addition, the vertical gradients in concentrations

between in- and above-canopy decrease with increasing life-

time, with 20 % difference for OH and almost zero for CO

(Fig. 2). Overall, the results here draw attention to the influ-

ences of the operator splitting on reactive trace gases such

as OH and NO, which are critically important for accurately

predicting the gas-phase chemistry and aerosol formation. In

summary, the advantages of the new operator order from the

modeling perspective include (i) generally better agreement

with the observations for the critical species NO and OH,

(ii) improved in-canopy gradients for emitted species such

as isoprene, and (iii) higher computational efficiency. There-

fore, we implement the new order of operations in FORCAsT

2.0.

3.2 Revised vertical mixing

The revised vertical mixing parameterization produces a

larger eddy diffusivity (

) and a more realistic boundary

layer height (Fig. 3a). The parameterization of

in the

boundary layer remains challenging, as

is a derived pa-

rameter that is analogous to diffusion yet is not completely

accurate for boundary layer turbulent mixing. Despite its lim-

itations, it is still a useful approximation when a full solution

for turbulence is computationally expensive to implement.

Kumar and Sharan (2012) compiled estimated values for

in the boundary layer based on previous studies, suggesting

that the magnitude of the

ranges from 60 to 200 m

−

under weakly unstable conditions (

z/L

=−

2). In FORCAsT

1.0,

peaks at 16 m

−

at 150 m above the ground during

the daytime (Fig. 3a). This leads to weak mixing, resulting

in an unrealistic end-of-day peak in isoprene concentrations

around sunset for the well-mixed day of 22 July (Fig. 3c).

The revised parameterization produces a larger

that falls

within the lower range of previously reported values (Ku-

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6317

Table 4.

2016 UMBS PROPHET observations utilized in this study for model evaluation, including concentration measurements, measure-

ment heights, instrumental techniques, and references. Chemical compound abbreviations are defined as follows: MVK

MACR (methyl

vinyl ketone and methacrolein), NO

(nitrogen oxides), O

(ozone), OH (hydroxyl radical), HO

(hydroperoxy radical), RO

(peroxy rad-

ical), CH

O (formaldehyde), ISOPOOH (isoprene hydroxy hydroperoxide), IEPOX (isoprene epoxydiol), IHN (isoprene hydroxy nitrate),

MTN (monoterpene hydroxy nitrate). IEPOX-SOA and 91Fac represent the isoprene-derived epoxydiol organic aerosols and monoterpene-

derived organic aerosols, respectively. The canopy height is 22.5 m.

Measurement

Height

Instrumental technique

Reference

(m)

Isoprene

34, 21,

PTR-QiTOF (Ionicon Analytik, GmbH)

Millet et al. (2018)

Monoterpenes

17, 13,

MVK

MACR

9, 5

A dual-channel custom-built chemiluminescence

Geddes and Murphy (2014),

instrument by Air Quality Design Inc.

Shi et al. (2021)

Model 205 (2B Technologies, Inc.) dual-beam

UV absorption instrument

OH, HO

Indiana University Laser-Induced Fluorescence–Fluorescence

Dusanter et al. (2009)

Assay by Gas Expansion (LIF-FAGE)

Ethane CHemical AMPlifier (ECHAMP) technique

Wood et al. (2016)

5, 17,

Harvard Fiber Laser-Induced

Hottle et al. (2009),

21, 30

Fluorescence (FILIF)

DiGangi et al. (2011),

instrument

Cazorla et al. (2015)

ISOPOOH, IEPOX, IHN

High-resolution GC chemical ionization mass

Vasquez et al. (2018b)

spectrometer (GC-HR-ToF-CIMS)

IHN, MTN

19.5

Iodide-adduct chemical ionization mass spectrometer

Ammonia (gas), PM

sulfate,

5.5

Ambient ion monitor–ion chromatography (AIM-IC,

Markovic et al. (2012)

nitrate, and ammonium

model 9000D, URG Corp., Chapel Hill, NC)

Organic aerosols, sulfate,

HR-ToF-AMS

Bui et al. (2021)

nitrate, ammonium

(Aerodyne Research Inc., USA)

IEPOX-SOA, 91Fac

HR-ToF-AMS (Aerodyne Research Inc., USA)

Bui et al. (2021)

and positive matrix factorization

mar and Sharan, 2012) and is sufficiently strong to produce

a realistic isoprene diurnal cycle under well-mixed condi-

tions (Fig. 3c). On the stagnant second day of the simula-

tion (23 July), the revised parameterization also reproduces

the end-of-day peak in isoprene but the modeled isoprene is

much lower than the observed isoprene (Fig. 3c).

FORCAsT predicts the air temperature based on mixing

of surface heat instead of prescribing the measured temper-

ature to nudge the model; therefore, the vertical mixing im-

pacts the temperature profiles (Fig. 3b). During the daytime

the canopy behaves as a heat source based on the leaf energy

balance, and the stronger mixing distributes the heat more

evenly throughout the model atmosphere. Therefore, a larger

reduces air temperature during the daytime in FORCAsT

(Fig. 3b). At night the canopy is cooler than the air aloft due

to the longwave radiation emission from the canopy, and the

low nocturnal mixing (

K <

5 m

−

, not shown) fails to

mix the warmer air down to the canopy, resulting in night-

time cooling of 8

◦

C during 0–6 model hours and a smaller

minimum temperature than the observation around 30 model

hours (Fig. 3b). The unrealistic nighttime cooling during 0–

6 model hours indicates that (i) heat capacity of leaves may

be important as a heat source at night, and/or (ii) nighttime

mixing is too low in FORCAsT. The idea of low nighttime

mixing is also supported by the overestimated nighttime iso-

prene concentrations (Fig. 3c), suggesting the need for a

better data-constrained nighttime mixing scheme. Note that

air temperature (Thomas, 2011) and evening isoprene decay

(Wei et al., 2020) also depend on horizontal advection, which

is not considered in this study. Attention should be on these

estimates when the homogeneity assumption is not met at the

site, particularly under stable conditions.

3.3 Isoprene updates to the gas-phase chemistry

mechanism

RCIM incorporates the current knowledge of isoprene chem-

istry under low-NO conditions, which is most notably man-

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6318

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

Figure 1.

Impacts of operator splitting on the vertical distribution of isoprene. The modeled and observed isoprene vertical profiles (normal-

ized by the concentration at the canopy) in the morning

(a)

, midday

(b)

, and evening

(c)

on 22 July 2016. Note that the observed isoprene

at 34, 21, 17, 13, 9, and 5 m is interpolated into model levels for comparison. The corresponding root mean square error (RMSE) of the

modeled isoprene profiles compared to the observations for the morning

(d)

, midday

(e)

, and evening

(f)

on 22 July 2016. The RMSEs are

also calculated separately in the canopy (In-cpy) and above the canopy (Abv. cpy).

Figure 2.

Impacts of operator splitting on species with chemi-

cal lifetimes ranging from seconds (OH; Di Carlo et al., 2004)

to months (CO; Holloway et al., 2000). The relative differences

between the two orders of operator (i.e.,

new

−

old

) for daytime

(10:00–16:00 local time) average hydroxyl radicals (OH), nitric ox-

ide (NO), formaldehyde (CH

O), isoprene, ozone (O

), and carbon

monoxide (CO) at two heights (0.8 and 1.6

z/h

ifested by the changes in the oxidation products of isoprene

(Fig. 4). Adding the RCIM isoprene chemistry does not drive

large changes in simulated isoprene concentrations (Fig. 4a),

which respond to OH, the isoprene emission rate, and mix-

Figure 3.

Impacts of the revised mixing parameterization on air

temperature and isoprene.

(a)

Modeled vertical profile of the eddy

diffusivity at 14:00 local time on 22 July 2016.

(b)

Modeled and

observed air temperature (

) at 46 m.

(c)

Modeled and observed

isoprene at 21 m.

ing, yet there is a substantial increase in MVK

MACR

(Fig. 4c). Generally the RCIM simulates MVK

MACR that

shows better agreement with the observations, although sim-

ulated concentrations are lower than observed on the second

day, which is likely due to underestimated isoprene (Fig. 4a).

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6319

The measured MVK

MACR peaks before isoprene on the

first day, suggesting horizontal advection and/or vertical

transport from the residual layer. Both RCIM and CACM1.0

underestimate formaldehyde (CH

O), a high-yield product

of isoprene oxidation, by over 50 % (Fig. 4j). Marvin et al.

(2017) show that existing isoprene mechanisms generally

underestimate CH

O by 17 %–33 %. This suggests missing

sources for CH

O, possibly heterogeneous conversion of

ISOPOOH on leaves, and/or unaccounted for VOC chemical

oxidation in the chemical mechanisms (Canaval et al., 2020;

DiGangi et al., 2011).

The modeled NO, NO

, and O

are very similar between

the two mechanisms (Fig. 4d–f). The timing of the early

morning NO maxima, caused by photolysis of NO

trans-

ported downward during the morning breakup of the noctur-

nal boundary layer (Seok et al., 2013), is captured by both

mechanisms. Both mechanisms overestimate NO

during 0–

6 model hours, suggesting missing sinks for NO

, possibly

the aqueous-phase reaction of NO

→

Adding RCIM also impacts HO

, with little change in

(Fig. 4h) but a substantial increase in OH (by 50 % on

the second day; Fig. 4g). Both mechanisms have HO

re-

generation from the H-shift isomerization of isoprene per-

oxy radicals (ISOPOO), which is important to sustain the

OH concentrations under low-NO conditions (Peeters et al.,

2009; Bates and Jacob, 2019). CACM1.0 recycles 1.0 HO

via the 1,6-H shift, while RCIM recycles 1.5 OH

0.7 HO

through the 1,6-H shift pathway and 1.0 OH through the 1,5-

H shift pathway (3.2 equivalents to HO

). Note that the 1,6-H

shift pathway dominates the 1,5-H shift pathway by a fac-

tor of 8 in RCIM. In addition, the ISOPOO concentrations

in CACM1.0 are lower than in RCIM by 50 % (Fig. 5a).

Therefore, the combination of higher ISOPOO and OH re-

generation efficiency leads to larger OH concentrations in

RCIM. Overall, the modeled OH is in good agreement with

the campaign-average measurements (Fig. 4g).

RCIM predicts higher RO

(HO

) than CACM1.0

(Fig. 4i). Because both mechanisms predict similar HO

which is discussed in the previous paragraph, here we fo-

cus on RO

that is predominantly composed of ISOPOO and

monoterpene-derived RO

(MT-RO

). Daytime (10–16 and

34–40 model hours) RCIM-RO

is higher than CACM1.0-

, mainly due to differences in ISOPOO (Fig. 5a), ac-

counting for 60 % and 40 % of daytime RO

in RCIM and

CACM1.0, respectively. Nighttime (0–6 and 24–30 model

hours) RCIM-RO

comprises 16 % ISOPOO and 84 % MT-

. The majority (87 %) of the nighttime MT-RO

is from

-initiated reactions and the remaining (13 %) is O

initiated. The observed NO

is below the limit of detec-

tion of the instrument (1.4 pptv) on the 2 simulated days

(not shown). While simulated nighttime NO

4 pptv

(roughly 0.4 pptv) and monoterpenes agree with observa-

tions, the nighttime RO

is still overestimated, suggesting

missing sinks for nighttime RO

. One potential sink is re-

action with NO. Using the observed NO to constrain the

model reduces the nighttime RO

by roughly 30 % (not

shown), but this is not sufficient to reproduce the observed

. Another possible sink is the accretion reactions of

monoterpene-derived RO

. The most recent accretion re-

action rates for OH-initiated and O

-initiated monoterpene

are included in FORCAsT (3

−

and 9

−

molecule

−

, respectively; Berndt et al., 2018).

However, laboratory updates on accretion rates for NO

initiated monoterpene RO

that dominate at night are not

available. Finally, dry deposition for RO

is not included in

the model due to lack of data, and this has the potential to be

an additional nighttime sink. Overall, the results suggest that

a better understanding of nighttime sinks of RO

is needed,

including chemical losses and dry deposition.

Because of the numerous mechanism changes for low-

isoprene chemistry, large deviations between RCIM and

CACM1.0 occur for the isoprene oxidation products, includ-

ing ISOPOOH (isoprene hydroxy hydroperoxide), IEPOX

(isoprene epoxydiol), and IHN (isoprene hydroxy nitrate)

(Fig. 4k–n). ISOPOOH, formed by the reaction of HO

with

ISOPOO, is increased by a factor of 2 in RCIM through-

out most of the diurnal cycle (Fig. 4k) due to a combination

of higher ISOPOO concentrations and a higher fraction of

ISOPOO following the HO

pathway (Fig. 5). IEPOX, pre-

dominantly produced by the reaction of ISOPOOH

OH, is

slightly higher in RCIM (Fig. 4l).

One important difference between the two mechanisms is

that the ISOPOO

−

NO reaction pathway increases in RCIM

(29 %) compared to CACM1.0 (11 %) (Fig. 5b, c). From a

low-NO study in the Amazon (NO

1 ppbv), Liu et al.

(2016) found that the ratio of the HO

pathway to NO path-

way is about unity. For these UMBS simulations with low-

NO conditions, RCIM shows a similar NO

ratio of

1.3, in contrast to the ratio of 3 in CACM1.0 (Fig. 5b, c).

The higher NO pathway percentage in RCIM consequently

increases IHN concentrations (Fig. 4m, n). The in-canopy

IHN (19.5 m) measured by an iodide-adduct chemical ioniza-

tion mass spectrometer (CIMS) is well reproduced by RCIM,

while the above-canopy IHN (32 m) by GC-HR-ToF-CIMS

is overestimated by a factor of 2. Because we do not simu-

late strong vertical gradients in isoprene (supported by ob-

servations; e.g., Wei et al., 2020), the IHN differences in the

observations are not likely due to vertical mixing or horizon-

tal advection. Vasquez et al. (2020) show that the 1,2-IHN

isomer undergoes a rapid hydrolysis loss not experienced

by the 4,3-IHN isomer. The in-canopy (19.5 m) measure-

ments by the iodide-adduct CIMS may be subject to uncer-

tainty in the hydrolysis-related isomer sensitivity (Fig. 4m,

n). Vasquez et al. (2020) estimate a condensed-phase hy-

drolysis coefficient of 4

M atm

−

for 1,2-IHN to

match their observed 1,2-IHN to 4,3-IHN isomer ratio. Us-

ing this hydrolysis loss coefficient, our model is able to repro-

duce the isomer ratio. However, the resulting total IHNs (1,2-

IHN+4,3-IHN) are still significantly higher than the mea-

surements (dashed line in Fig. 4m, n). A hydrolysis coef-

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6320

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

Figure 4.

Comparison of FORCAsT 2.0 with the RCIM (red) and CACM1.0 (blue) mechanisms against observations in and above the canopy

(heights noted). Gas-phase concentrations are evaluated versus observations (grey), including isoprene (34 m), monoterpenes (34 m), methyl

vinyl ketone

methacrolein (MVK

MACR, 34 m), nitric oxide (NO, 29 m), nitrogen dioxide (NO

, 29 m), ozone (O

, 6 m), hydroxyl radi-

cals (OH, 32 m), hydroperoxy radicals (HO

, 32 m), total peroxy radicals (RO

, 30 m), formaldehyde (CH

O, 21 m), hydroxy hydroperoxide

isomers (1,2-ISOPOOH

4,3-ISOPOOH, 32 m), epoxydiol isomers (

cis

-IEPOX

rans

-IEPOX, 32 m), isoprene hydroxy nitrates (1,2-IHN

4,3-IHN, 32 m), IHN (19.5 m), and monoterpene hydroxy nitrates (MTN, 19.5 m) in sequence from

(a)

(o)

. Grey shaded areas denote

1 standard deviation of the data when available. Instrumental information is in Table 4. The red dashed lines in

(m)

and

(n)

denote modeled

IHN with a hydrolysis loss rate of 4

M atm

−

for 1,2-IHN.

ficient of 1

M atm

−

is required for the model to

match the measured 1,2-IHN concentrations at 32 m; how-

ever, the ratio is then underestimated because the model sim-

ulates more 4,3-IHN than observed (not shown), which may

be due to the greater production of the 4,3-IHN driven by

higher than observed HO

. Overall, the differences in IHN

observations and the measured–modeled discrepancies in the

IHN isomers suggest that large uncertainties still exist in the

simulated IHN concentrations.

The monoterpene hydroxy nitrates (MTNs) in RCIM are

slightly lower than that in CACM1.0 at night (Fig. 4o),

likely driven by lower nighttime monoterpene concentra-

tions in RCIM (Fig. 4b). Overall, the MTN is underesti-

mated by both RCIM and CACM1.0, and this could be af-

fected by the rapid monoterpene OH RO

accretion reac-

tions, which decrease the available monoterpene RO

for

generating MTN during the daytime. Additionally, the dry

deposition of MTN is prescribed to be 50 % of IHNs, which

lack observational constraints and may result in an overcon-

sumption of MTN. To further evaluate the representation of

the low-NO chemistry of isoprene by the two mechanisms,

we compare the ratio of ISOPOOH to MVK

MACR. For

unpolluted regions, the reaction of ISOPOO with HO

is the

dominant pathway, with ISOPOOH isomers as the major ox-

idation products. Reaction of ISOPOO with NO dominates

in polluted regions, with major oxidation products MVK

and MACR. MVK

MACR can also be produced through

the isomerization and self-reaction of ISOPOO, with MVK

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

6321

Figure 5.

Fate of isoprene hydroxy peroxy radicals (ISOPOO) in

RCIM and in CACM1.0.

(a)

Diurnal profiles of ISOPOO concentra-

tions in RCIM and CACM1.0.

(b, c)

Fraction of daytime-averaged

(12:00–14:00 LT) production rates of ISOPOO through reactions

with peroxy radicals (RO

), hydroperoxy radicals (HO

), nitric ox-

ide (NO), and H-shift isomerization in RCIM and CACM1.0.

dominating over MACR by approximately a factor of 2.

Therefore, the ratio of ISOPOOH to MVK

MACR reflects

the contribution of the HO

oxidation pathway relative to

other oxidation pathways of ISOPOO. The consistency be-

tween the measured and modeled ISOPOOH fraction thus

reflects the skill of the isoprene mechanisms in represent-

ing the isoprene chemistry under low-NO conditions. The

most frequent (51 %) observed daytime NO levels when the

ISOPOOH data are available range from 20 to 40 pptv. The

measurement-based ratio of ISOPOOH to MVK

MACR de-

creases quickly and ranges from 0.15 to 0.05 with increas-

ing NO when NO

40 pptv (Fig. 6), indicating that the fate

of ISOPOO is highly sensitive to NO levels under low-NO

conditions. The simulated dependence of the ratio on NO

concentrations by RCIM is in good agreement with the ob-

served, although the simulated ratio decays faster than the

observations. In contrast, the ratios simulated by CACM1.0

are higher and decrease slower with increasing NO concen-

trations, partly due to lower OH (Fig. 6) that reacts faster

with ISOPOOH than with MVK.

In summary, the fate of ISOPOO in RCIM is different than

in CACM1.0 (Fig. 5), leading to differences in isoprene oxi-

dation products. The HO

pathway ratio is higher in RCIM,

leading to higher ISOPOOH concentrations, and the NO

pathway of ISOPOO is also enhanced in RCIM, resulting in

higher MVK

MACR and IHN concentrations. This comes

at the expense of reducing the H-shift isomerization pathway

in RCIM, yet OH concentrations are still enhanced due to

higher HO

recycling efficiency in RCIM. Overall, isoprene

Figure 6.

Daytime (08:00–17:00 local time) measurement-based

and modeled ratio of ISOPOOH (isoprene hydroxy hydroperox-

ides) to MVK

MACR (methyl vinyl ketone

methacrolein).

Measurement-based ratios of ISOPOOH to MVK

MACR are from

23–27 July 2016 (Grey dots). The fitted exponential curve is for

measurement-based ratios (black line), RCIM ratios (red line),

CACM1.0 ratios (blue line), RCIM OH concentrations (red dash

line), and CACM1.0 OH concentrations (blue dash line). The yel-

low patch denotes the most frequently observed nitric oxide (NO)

levels (20–40 pptv) during the period of 23–27 July 2016.

oxidation products in RCIM generally compare well with

observations, and RCIM captures the decreasing trend of

ISOPOOH

MVK

MACR

)

with increasing NO concentra-

tions, indicating a better representation of low-NO chemistry.

One aspect to note about the one-dimensional framework is

that it prioritizes local processes (emissions, chemistry, and

deposition) over that of advection. If we assume a horizon-

tally homogeneous canopy as utilized in flux analysis, then

advection should not be a factor. However, observations are

clear that the site can be influenced by long-range transport,

such as the advection of high-NO

conditions from southern

urban locations like Detroit, Milwaukee, and Chicago (e.g.,

Cooper et al., 2001; VanReken et al., 2015). One factor to

consider is if horizontal advection also influences the diurnal

cycle of isoprene and monoterpene oxidation products. For

example, above-canopy short-lived oxidation products (such

as ISOPOOH and IHN) are overestimated in the model com-

pared to the observations, whereas longer-lived species such

as IEPOX (Fig. 4k–m) show better agreement with the ob-

servations. Further analysis is needed that includes upwind

sources of these oxidation products, as is potential analysis

with wind direction and VOC sources.

3.4 Secondary organic aerosol formation

As discussed in Sect. 2.5, an isoprene-derived SOA (iSOA)

parameterization based on Griffin et al. (2005) is imple-

mented in FORCAsT 2.0, incorporating six new surrogate

groups for the isoprene-derived precursors such as IEPOX

https://doi.org/10.5194/gmd-14-6309-2021

Geosci. Model Dev., 14, 6309–6329, 2021

6322

D. Wei et al.: FORest Canopy Atmosphere Transfer (FORCAsT) 2.0

(isoprene epoxydiol) and tetrafunctionals. As described in

Ashworth et al. (2015), the partitioning of organic gases into

the aerosol phase is a function of the liquid water content

(Eq. 11).

Both monoterpene SOA (MNT-SOA) and iSOA mainly

form in the boundary layer (

∼

5 km; Fig. 7a, b) in response

to vertical mixing, vertical distributions of SOA precursors,

and meteorological conditions such as relative humidity. The

simulated boundary-layer-averaged MNT-SOA is roughly

1 μg m

−

. The MNT-SOA above the canopy (30 m)

is underestimated by a factor of 2–3 compared to the HR-

ToF-AMS positive matrix factorization (PMF) concentra-

tions (Fig. 3c), which is likely due in part to the underesti-

mated precursors such as MTN (Fig. 4o). Another possible

reason is that the partitioning efficiency for MNT-SOA is un-

derestimated due to phase separation in MPMPO. MPMPO

considers the two phases separately: the organic phase and

the aqueous phase. In reality, water uptake to organics in-

creases the partitioning medium, which in turn increases the

partitioning efficiency of condensable organic species (Pye

et al., 2017), and this could effectively reduce the simulated

MNT-SOA in MPMPO.

The majority of simulated iSOA (

99 %) is formed

through the aqueous phase is and therefore sensitive to the

liquid water content (Fig. 7b, e, f). The boundary-layer-

averaged iSOA is 1

0 μg m

−

(Fig. 7b). IEPOX-SOA

dominates the modeled iSOA formation (Fig. 8a, b), and

simulated concentrations agree well with the daytime HR-

ToF-AMS PMF-derived IEPOX-SOA at 30 m (Fig. 7d). The

modeled IEPOX-SOA shows a diurnal cycle that follows gas-

phase IEPOX concentrations (Fig. 4l), and the rapid increase

in IEPOX-SOA around 35 model hours results from the in-

crease in liquid water content (Fig. 7e). No such increase in

IEPOX-SOA is manifested around 6 model hours when liq-

uid water content is also high due to the absence of iSOA

precursors at the beginning of the simulation (Fig. 4l).

FORCAsT 2.0 simulates an in-canopy iSOA mass yield

of 7 %, which is slightly higher than values commonly used

in global models (0.9 %–6.8 %), to represent the ambient at-

mosphere (Carlton et al., 2009). In the model, the two dom-

inant iSOA surrogate groups are the tetrafunctionals (tetra-

SOA) and IEPOX (IEPOX-SOA), which combined account

for 86 % of the total iSOA in the canopy and 98 % in the

boundary layer (Fig. 8a, b). RCIM has a low yield of gly-

oxal from isoprene, and the GLYX-SOA component is 4 %

in the canopy and negligible in the boundary layer. The iSOA

from organic nitrates (IHN-SOA and C4-SOA) is more im-

portant (11 %) in the canopy likely due to higher NO lev-

els (Fig. 8a, b). Both IEPOX-SOA and tetra-SOA increase

with height (Fig. 8c, d) following the gas-phase precursors

(not shown). We note that dry deposition of tetrafunctionals

lacks data constraints, and in this study, the dry deposition

of these compounds is set to be the same as IHN (Table 2)

due to the common nitrate functional group. This leads to the

deposition of 60 % of the gas-phase tetrafunctionals in the

canopy and is a major source of uncertainties for estimation

of tetra-SOA. Using RCIM, Bates and Jacob (2019) report

an even higher iSOA mass yield of 25 %, with IEPOX-SOA,

organonitrates, and tetra-SOA each contributing about one-

third to the total mass. The differences in the iSOA yield

and composition between this study and Bates and Jacob

(2019) may derive from (i) the NO

levels, which are ex-

pected to significantly impact the iSOA precursors (i.e., gas-

phase IEPOX, tetrafunctionals, and organic nitrates), and/or

(ii) the different iSOA formation schemes.

There are a few points to note in our estimation for LWC.

The LWC

org

accounts for almost 40 % of total LWC on av-

erage (Fig. 7e), suggesting that water uptake to organics is

important for iSOA formation at the site. However, LWC

org

is highly sensitive to the hygroscopicity growth for organ-

ics (

org

), ranging from 0.01 to 0.5 for slightly to very hy-

groscopic organic species (Petters and Kreidenweis, 2007).

Our calculation gives a mean

org

of 0.08 for the simulation

period. Cerully et al. (2015) partitioned

org

for different or-

ganic aerosol compositions with 0

02 for monoterpene

SOA and 0

02 for isoprene SOA, respectively. To test

the impact of

org

on LWC estimation, we run a sensitivity

simulation with a higher

org

of 0.15, resulting in an increase

of 0.4 μg m

−

in LWC

org

on average and thus a relatively

small increase in total LWC (Fig. 7e). Therefore, the higher

org

of 0.15 leads to a negligible increase in iSOA formation.

Slade et al. (2019) derived different

org

values for daytime

(0.02) and nighttime (0.15) due to diurnal changes in organic

aerosol composition at the same study site, which would re-

quire a diurnal cycle of LWC

org

with lower values at night.

Because this result contrasts with current methods to derive

LWC, we were unable to test whether this would impact our

partitioning simulations, but we note that this large discrep-

ancy could have implications for the simulation of biogeni-

cally derived organic aerosol over the diurnal cycle. Finally,

we note that due to lack of measurements outside the canopy,

the vertical distribution of LWC is scaled by relative humid-

ity (RH) in the model. During the daytime (nighttime), LWC

is relatively constant within the boundary (surface) layer and

decreases rapidly with altitude above the boundary (surface)

layer (Fig. 7f). In general, this trend is consistent with the

meteorology soundings of RH from the nearest location (in

Gaylord, Michigan), but the soundings suggest that the day-

time boundary layer height and the nighttime surface layer

height are underestimated in the model despite the improved

vertical mixing formulation.

4 Conclusions

We update FORCAsT 1.0 to FORCAsT 2.0 by including

(i) splitting the integration of chemistry from emission and

dry deposition to provide more realistic representations of

vertical gradients in the forest canopy and to make the chem-

ical module more flexible for future chemical mechanism

Geosci. Model Dev., 14, 6309–6329, 2021

https://doi.org/10.5194/gmd-14-6309-2021