Beijing Climate Center Earth System Model version 1 (BCC-ESM1): model description and evaluation of aerosol simulations

The Beijing Climate Center Earth System Model version 1 (BCC-ESM1) is the first version of a fully coupled Earth system model with interactive atmospheric chemistry and aerosols developed by the Beijing Climate Center, China Meteorological Administration. Major aerosol species (including sulfate, organic carbon, black carbon, dust, and sea salt) and greenhouse gases are interactively simulated with a whole panoply of processes controlling emission, transport, gas-phase chemical reactions, secondary aerosol formation, gravitational settling, dry deposition, and wet scavenging by clouds and precipitation. Effects of aerosols on radiation, cloud, and precipitation are fully treated. The performance of BCC-ESM1 in simulating aerosols and their optical properties is comprehensively evaluated as required by the Aerosol Chemistry Model Intercomparison Project (AerChemMIP), covering the preindustrial mean state and time evolution from 1850 to 2014. The simulated aerosols from BCC-ESM1 are quite coherent with Coupled Model Intercomparison Project Phase 5 (CMIP5)-recommended data, in situ measurements from surface networks (such as IMPROVE in the US and EMEP in Europe), and aircraft observations. A comparison of modeled aerosol optical depth (AOD) at 550 nm with satellite observations retrieved from the Moderate Resolution Imaging Spectroradiometer (MODIS) and the Multi-angle Imaging SpectroRadiometer (MISR) and surface AOD observations from the AErosol RObotic NETwork (AERONET) shows reasonable agreement between simulated and observed AOD. However, BCC-ESM1 shows weaker upward transport of aerosols from the surface to the middle and upper troposphere, likely reflecting the deficiency of representing deep convective transport of chemical species in BCC-ESM1. With an overall good agreement between BCC-ESM1 simulated and observed aerosol properties, it demonstrates a success of the implementation of interactive aerosol and atmospheric chemistry in BCC-ESM1.

1 Introduction

Atmosphere is a thin gaseous layer around the Earth, consisting of nitrogen, oxygen, and a large number of trace gases including important greenhouse gases (GHGs) such as water vapor, tropospheric ozone (O₃), carbon dioxide (CO₂), methane (CH₄), nitrous oxide (N₂O), and chlorofluorocarbons (CFCs). Besides gaseous components, atmosphere also contains various aerosols, which are important for cloud formation and radiative transfer. Atmospheric trace gases and aerosols are actually interactive components of the climate system. Their inclusion in global climate models (GCMs) is a significant enhancement for most state-of-the-art climate models (Lamarque et al., 2013; Collins et al., 2017). Early attempts at coupling global climate dynamics with atmospheric chemistry can be traced back to the late 1970s, when 3D transport of ozone and simple stratospheric chemistry were first incorporated into a GCM to simulate global O₃ production and transport (e.g., Cunnold et al., 1975; Schlesinger and Mintz, 1979). Since the mid-1980s, a large number of online global climate–chemistry models have been developed to address issues of the Antarctic stratospheric O₃ depletion (e.g., Cariolle et al., 1990; Austin et al., 1992; Solomon, 1999), tropospheric O₃ and sulfur cycle (e.g., Feichter et al., 1996; Barth et al., 2000), tropospheric aerosol, and its interactions with cloud (e.g., Chuang et al., 1997; Lohmann et al., 2000; Ghan and Easter, 2006; Jacobson, 2012). Aerosols and chemically reactive gases in the atmosphere exert important influences on global and regional air quality and climate (Collins et al., 2017).

Since 2013, the Beijing Climate Center (BCC), China Meteorological Administration, has continuously developed and updated its fully coupled GCM, the Beijing Climate Center Climate System Model (BCC-CSM) (Wu et al., 2013, 2014, 2019). BCC-CSM version 1.1 was one of the comprehensive carbon–climate models participating in phase five of the Coupled Model Intercomparison Project (CMIP5; Taylor et al., 2012). When forced by prescribed historical emissions of CO₂ from the combustion of fossil fuels and land-use change, BCC-CSM1.1 successfully reproduced the trends in observed atmospheric CO₂ concentration and global surface air temperature from 1850 to 2005 (Wu et al., 2013). During recent years, BCC-CSM1.1 has been used in numerous investigations on soil organic carbon changes (e.g., Todd-Brown et al., 2014), ocean biogeochemistry changes (e.g., Mora et al., 2013), and carbon–climate feedbacks (e.g., Arora et al., 2013; Hoffman et al., 2014). BCC-CSM includes the main climate–carbon cycle processes (Wu et al., 2013), and the global mean atmospheric CO₂ concentration is calculated from a prognostic equation of CO₂ budget taking into account global anthropogenic CO₂ emissions and interactive land–atmosphere and ocean–atmosphere CO₂ exchanges.

In recent years, BCC has put large efforts into developing a global climate–chemistry–aerosol fully coupled Earth system model (the Beijing Climate Center Earth System Model version 1, BCC-ESM1) on the basis of BCC-CSM2 (Wu et al., 2019). The objective is to interactively simulate global aerosols (e.g., sulfate, black carbon) and the main greenhouse gases (e.g., O₃, CH₄, N₂O, and CO₂) in the atmosphere and to investigate feedbacks between climate and atmospheric chemistry. BCC-ESM1 is at the point of being publicly released, and it is actively used by BCC for several Coupled-Model Intercomparison Project 6 (CMIP6)-endorsed research initiatives (Eyring et al., 2016), including the Aerosol Chemistry Model Intercomparison Project (AerChemMIP; Collins et al., 2017) and the Coupled Climate–Carbon Cycle Model Intercomparison Project (C4MIP; Jones et al., 2016).

The purpose of this paper is to evaluate the performance of BCC-ESM1 in simulating aerosols and their optical properties in the 20th century. The description of BCC-ESM1 is presented in Sect. 2. The experimental protocol is given in Sect. 3. Section 4 presents the evaluations of aerosol simulations with comparisons to CMIP5-recommended data (Lamarque et al., 2010) and data obtained from both global surface networks and satellite observations. The regional and global characteristics compared to observations and estimates from other studies are analyzed. Simulations of aerosol optical properties in the 20th century are also analyzed in Sect. 4. Conclusions and discussions are summarized in Sect. 5. Information about code and data availability are found at the end of the paper.

Table 1 Chemical species considered in BCC-AGCM3-Chem. Species marked with star (^*) denote those added in BCC-ESM1 apart from the 63 species used in MOZART2. In the column of surface emission, interactive surface emissions are considered for sea salt and dust.

2 Model description

BCC-ESM1 is an Earth system model with interactive chemistry and aerosol components, in which the atmospheric component is the BCC Atmospheric General Model version 3 (Wu et al., 2019) with interactive atmospheric chemistry (hereafter BCC-AGCM3-Chem), the land component BCC Atmosphere and Vegetation Interaction Model version 2.0 (hereafter BCC-AVIM2.0), the ocean component Modular Ocean Model version 4 (MOM4)-L40, and sea ice component (sea ice simulator, SIS). Different components of BCC-ESM1 are fully coupled and interact with each other through fluxes of momentum, energy, water, carbon, and other tracers at their interfaces. The coupling between the atmosphere and the ocean is done every hour.

The atmospheric component BCC-AGCM3-Chem is able to simulate global atmospheric composition and aerosols from anthropogenic emissions as forcing agents. Its resolution is T42 (approximately $\mathrm{2.8125}{}^{\circ }\times \mathrm{2.8125}{}^{\circ }$ transformed spectral grid). The model has 26 levels in a hybrid sigma–pressure vertical coordinate system with the top level at 2.914 hPa. Details of the model physics are described in Wu et al. (2019). The BCC-AGCM3-Chem combines 66 gas-phase chemical species and 13 bulk aerosol compounds as listed in Table 1. Apart from 3 gas-phase species of dimethyl sulfide (DMS), sulfur dioxide (SO₂), and ammonia (NH₃), the other 63 gas-phase species are the same as those in the “standard version” of MOZART2 (Model for Ozone and Related chemical Tracers, version 2), a global chemical transport model for the troposphere developed by the National Center for Atmospheric Research (NCAR) driven by meteorological fields from either climate models or assimilations of meteorological observations (Horowitz et al., 2003). Advection of all tracers in BCC-AGCM3-Chem is performed through a semi-Lagrangian scheme (Williamson and Rasch, 1989), and vertical diffusion within the boundary layer follows the parameterization of Holtslag and Boville (1993). The gas-phase chemistry of the 63 MOZART2 gas-phase species as listed in Table 1 is treated in the same way as that in the standard version of MOZART2 (Horowitz et al., 2003), and there are 33 photolytic reactions and 135 chemical reactions involving 30 dry deposited chemical species and 25 soluble gas-phase species. Dry-deposition velocities for the 15 trace gases including O₃, carbon monoxide (CO), CH₄, formaldehyde (CH₂O), acetic acid (CH₃OOH), hydrogen peroxide (H₂O₂), nitrogen dioxide (NO₂), nitric acid (HNO₃), polyacrylonitrile (PAN), acetone (CH₃COCH₃), peroxyacetic acid (CH₃COOOH), acetaldehyde (CH₃CHO), methylglyoxal (CH₃COCHO), nitric oxide (NO), and pernitric acid (HNO₄) are not computed interactively and are directly interpolated from MOZART2 climatological monthly mean deposition velocities (https://en.wikipedia.org/wiki/MOZART_(model), last access: 1 May 2014), which are calculated offline (Bey et al., 2001; Shindell et al., 2008) using a resistance-in-series scheme originally described in Wesely (1989). The dry-deposition velocities for the other 15 species including peroxy acetyl nitrate (PAN), methyl nitroacetate (ONIT), organic nitrates (ONITR), ethyl alcohol (C₂H₅OH), organic hydroxiperoxide (POOH), ethyl hydroperoxide (C₂H₅OOH), propylhydroperoxide (C₃H₇OOH), methylene glycol mono acetate (ROOH), glycolaldehyde (GLYALD), acetol (HYAC), methanol (CH₃OH), propanoic acid (MACROOH), isoprene hydroxy hydroperoxide (ISOPOOH), carboxylic acid (XOOH), formaldehyde (HYDRALD), and hydrogen (H₂) are calculated using prescribed deposition velocities of O₃, CO, CH₃CHO, or land surface type and surface temperature following the MOZART2 (Horowitz et al., 2003). Wet removal by in-cloud scavenging for 25 soluble gas-phase species in the standard version of MOZART2 uses the parameterization of Giorgi and Chameides (1985) based on their temperature-dependent effective Henry's law constants. In-cloud scavenging is proportional to the amount of cloud condensate converted to precipitation, and the loss rate depends on the amount of cloud water, the rate of precipitation formation, and the rate of tracer uptake by the liquid phase water. Other highly soluble species such as HNO₃, H₂O₂, ONIT, ISOPOOH, MACROOH, XOOH, and lead (Pb-210) are also removed by below-cloud washout as calculated using the formulation of Brasseur et al. (1998). Below-cloud scavenging is proportional to the precipitation flux in each model layer and the loss rate depends on the precipitation rate. Vertical transport of gas tracers and aerosols due to deep convection is not yet included in the present version of BCC-AGCM3-Chem, a process which is regarded as a part of the deep convection and occurs generally in a small spatial region on a GCM box with a low resolution (2.8^∘ lat. × 2.8^∘ long.). Another consideration is that a large uncertainty exists in treating transport of those water-soluble tracers by deep convection. But this effect will be involved in the next version of BCC model.

The BCC-AVIM2.0 is the land model with a terrestrial carbon cycle. It is described in detail in Li et al. (2019) and includes biophysical, physiological, and soil carbon–nitrogen dynamical processes. The terrestrial carbon cycle operates through a series of biochemical and physiological processes on photosynthesis and respiration of vegetation. Biogenic emissions from vegetation are computed online in BCC-AVIM2.0 following the algorithm of the Model of Emissions of Gases and Aerosols from Nature version 2.1 (MEGAN2.1; Guenther et al., 2012).

The oceanic component of BCC-ESM1 is the Modular Ocean Model version 4 with 40 levels (hereafter MOM4-L40) and the sea ice component SIS. MOM4-L40 uses a tripolar grid of horizontal resolution with 1^∘ longitude by 1∕3^∘ latitude between 30^∘ S and 30^∘ N, 1^∘ longitude by 1^∘ latitude from 60^∘ S and 60^∘ N poleward, and 40 z levels in the vertical. Carbon exchange between the atmosphere and the ocean are calculated online in MOM4-L40 using a biogeochemistry module that is based on the protocols from the Ocean Carbon Cycle Model Intercomparison Project–Phase 2 (OCMIP2, http://www.ipsl.jussieu.fr/OCMIP/phase2/, last access: 1 August 2011). SIS has the same horizontal resolution as MOM4-L40 and three layers in the vertical, including one layer of snow cover and two layers of equally sized sea ice. Details of oceanic component MOM4-L40 and sea ice component SIS that are used in BCC-ESM1 may be found in Wu et al. (2013, 2019).

In the following subsections, we will describe the treatments in BCC-ESM1 for 3 gas-phase species of DMS, SO₂ and NH₃, 13 prognostic aerosol species including sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$), 2 types of organic carbon (hydrophobic OC1, hydrophilic OC2), 2 types of black carbon (hydrophobic BC1, hydrophilic BC2), 4 categories of soil dust (DST01, DST02, DST03, DST04), and 4 categories of sea salt (SSLT01, SSLT02, SSLT03, SSLT04). Concentrations of all aerosols in BCC-ESM1 are mainly determined by advective transport, emission, dry deposition, gravitational settling, and wet scavenging by clouds and precipitation, except for ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, the gas-phase and aqueous-phase conversion of which from SO₂ are also considered. The present version of aerosol scheme belongs to a bulk aerosol model and mainly refers to the scheme of CAM-Chem (Lamarque et al., 2012), but the nucleation and coagulation of aerosols are still ignored.

Table 2 Gas-phase chemical reactions for NH₃ and bulk aerosols precursors following CAM-Chem (Lamarque et al., 2012). The reaction rates (s⁻¹) refer to Tie et al. (2001) and Sander et al. (2003) and Cooke and Wilson (1996). Temperature (T) is expressed in kelvin, air density (M) in molecule per cubic centimeter, and k_i and k_o in cubic centimeters per molecule per second.

2.1 SO₂, DMS, NH₃, and sulfate

SO₂ is a main sulfuric acid precursor to form aerosol sulfate ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$. Conversions of SO₂ to ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ occur by gas-phase reactions (Table 2) and by aqueous-phase reactions in cloud droplets. The dry-deposition velocity of SO₂ follows the resistance-in-series approach of Wesely (1989) using the formula $W_{{\mathrm{SO}}_{\mathrm{2}}}=\mathrm{1}/(r_a+r_b+r_c)$, in which r_a, r_b, and r_c are the aerodynamic resistance, the quasi-laminar boundary layer resistance, and the surface resistance, respectively and they are interactively computed in each model time step. The loss rate of SO₂ due to wet deposition is computed following the scheme in the global Community Atmosphere Model (CAM) version 4, the atmospheric component of the Community Earth System Model (Lamarque et al., 2012).

The sources of SO₂ mainly come from fuel combustion, industrial activities, and volcanoes. SO₂ can also be formed from the oxidation of DMS as listed in Table 2, in which their reaction rates follow CAM-Chem (Lamarque et al., 2012). The main source of DMS is from oceanic emissions via biogenic processes. It is prescribed with the climatological monthly data that are extracted from the MOZART2 package (https://www2.acom.ucar.edu/gcm/mozart-4, last access: 10 February 2018). ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ is one of the prognostic aerosols in BCC-AGCM3-Chem. Its treatment follows CAM4-Chem (Lamarque et al., 2012). It is produced primarily by the gas-phase oxidation of SO₂ (in Table 2) and by aqueous-phase oxidation of SO₂ in cloud droplets. The gas-phase reactions, rate constants, and gas–aqueous equilibrium constants are given by Tie et al. (2001). The heterogeneous reactions of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ occur on all aerosol surfaces. Their treatment follows a bulk aerosol model (BAM) used in CAM4 (Neale et al., 2010). The heterogeneous reactions depend strongly on pH values in clouds, which are calculated from the concentrations of SO₂, HNO₃, H₂O₂, NH₃, O₃, HO₂, and ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$. NH₃ is a gas tracer apart from MOZART2 (Table 1). Its sources include aircraft and surface emissions due to anthropogenic activity, biomass burning, and biogenic emissions from land soil and ocean surfaces (Table 4). ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ is assumed to be all in aqueous phase due to water uptake, although Wang et al. (2008a) showed that ∼ 34 % of sulfate particles are in solid phase globally due to the hysteresis effect of ammonium sulfate phase transition. However, in terms of radiative forcing, consideration of the solid sulfate formation process lowers the sulfate forcing by ∼ 8 % as compared to a consideration of all sulfate particles in the aqueous phase (Wang et al., 2008b). Future model development may consider the life cycle of NH₃. The sulfate in- and below-cloud scavenging follows Neu and Prather (2011). Washout of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ is set to 20 % of the washout rate of HNO₃ following Tie et al. (2005) and Horowitz (2006). The dry-deposition velocity of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ is also calculated by the resistance-in-series approach.

2.2 Aerosols of organic carbon and black carbon

BCC-AGCM3-Chem treats two types of organic carbon (OC), i.e., water-insoluble tracer OC1 and water-soluble tracer OC2, and two types of black carbon (BC), i.e., water-insoluble tracer BC1 and water-soluble tracer BC2. As shown in Table 2, hydrophobic BC1 and OC1 can be converted to hydrophilic BC2 and OC2 with a constant rate of $\mathrm{7.1}\times {\mathrm{10}}^{-\mathrm{6}}$ s⁻¹ (Cooke and Wilson, 1996). The four tracers of organic carbon and black carbon are mainly from emissions including both fossil fuel and biomass burning and are from the CMIP6 data package (https://esgf-node.llnl.gov/search/input4mips/, last access: 10 March 2019; Hoesly et al., 2018). Beside anthropogenic and biomass burning emissions, hydrophilic organic carbon OC2 can also come from natural biogenic volatile organic compound (VOC) emissions. Dry-deposition velocities for all the four OC and BC tracers are set to 0.001 m s⁻¹. OC2 and BC2 are soluble aerosols, and their sinks are primarily governed by wet deposition. Their in- and below-cloud scavenging follows the scheme of Neu and Prather (2011).

Table 3 Size and density parameters of bulk aerosols.

2.3 Sea salt aerosols

As shown in Table 3, sea salt aerosols in the model are classified into four size bins (0.2–1.0, 1.0–3.0, 3.0–10, and 10–20 µm) in diameter. They originate from oceans and are calculated online by BCC-ESM1. The upward flux F_sea-salt of sea salt productions for four bins is proportional to the 3.41 power of the wind speed u_10 m at 10 m height near the sea surface (Mahowald et al., 2006) and is expressed as

$$\begin{array}{}\text{(1)}& F_{\text{sea-salt}}=S\cdot (u_{\text{10\hspace{0.17em}m}}{)}^{\mathrm{3.41}},\end{array}$$

where S is a scaling factor and set to $\mathrm{4.05}\times {\mathrm{10}}^{-\mathrm{15}}$, $\mathrm{4.52}\times {\mathrm{10}}^{-\mathrm{14}}$, $\mathrm{1.15}\times {\mathrm{10}}^{-\mathrm{13}}$, and $\mathrm{1.20}\times {\mathrm{10}}^{-\mathrm{13}}$ for four size bins of sea salt aerosols in BCC-ESM1, respectively.

Dry deposition of sea salt depends on the turbulent deposition velocity in the lowest atmospheric layer using aerodynamic resistance and the friction velocity and the settling velocity through the whole atmospheric column for each bin of sea salt. The turbulent deposition velocity and settling velocity depend on particle diameter and density (listed in Table 3). In addition, the fact that the size of sea salt changes with humidity is also considered. The wet deposition of sea salt follows the scheme for soluble aerosols used in CAM4 and depends on prescribed solubility and size-independent scavenging coefficients.

2.4 Dust aerosols

Dust aerosols behave in a similar way to sea salt. Their variations involve three major processes: emission, advective transport, and wet or dry depositions. The dust emission is based on a saltation–sandblasting process and depends on wind friction velocity, soil moisture, and vegetation or snow cover (Zender et al., 2003). The vertical flux of dust emission is corrected by a surface erodible factor at each model grid cell which has been downloaded from the NCAR website (https://svn-ccsm-inputdata.cgd.ucar.edu/trunk/inputdata/atm/cam/dst/, last access: 1 May 2014). Soil erodibility is prescribed by a physically based geomorphic index that is proportional to the runoff area upstream of each source region (Albani et al., 2014). Like sea salt, dry deposition of dust aerosols includes gravitational and turbulent deposition processes, while wet deposition results from both convective and large-scale precipitation and is dependent on prescribed size-independent scavenging coefficients.

2.5 Effects of aerosols on radiation, clouds, and precipitation

The mass mixing ratios of bulk aerosols are prognostic variables in BCC-ESM1 and directly affect the radiative transfer in the atmosphere with their treatments following the NCAR Community Atmosphere Model (CAM3; Collins et al., 2004). Indirect effects of aerosols are taken into account in the present version of BCC-AGCM3-Chem (Wu et al., 2019). Aerosol particles act as cloud condensation nuclei and exert influence on cloud properties and precipitation and ultimately impact the hydrological cycle. Prognostic aerosol masses are used to estimate the liquid cloud droplet number concentration N_cdnc (cm⁻³) in BCC-AGCM3-Chem. N_cdnc is explicitly calculated using the empirical function suggested by Boucher and Lohmann (1995) and Quaas et al. (2006):

$$\begin{array}{}\text{(2)}& N_{\text{cdnc}}=\exp \left[\mathrm{5.1}+\mathrm{0.41}\ln \left(m_{\text{aero}}\right)\right],\end{array}$$

where m_aero (µg m⁻³) is the total mass of all hydrophilic aerosols,

$$\begin{array}{}\text{(3)}& m_{\text{aero}}=m_{\text{SS}}+m_{\text{OC}}+m_{{\mathrm{SO}}_{\mathrm{4}}}+m_{{\mathrm{NH}}_{\mathrm{4}}{\mathrm{NO}}_{\mathrm{2}}},\end{array}$$

i.e., the first bin of sea salt (m_SS), hydrophilic organic carbon (m_OC), sulfate ($m_{{\mathrm{SO}}_{\mathrm{4}}}$), and ammonium nitrite (NH₄NO₂). A dataset of NH₄NO₂ from NCAR CAM-Chem (Lamarque et al., 2012) is used in our model.

N_cdnc is an important factor in determining the effective radius of cloud droplets for radiative calculation. The effective radius of cloud droplets r_el is estimated as

$$\begin{array}{}\text{(4)}& r_{\text{el}}=\mathit{\beta }\cdot r_{\text{l,vol}},\end{array}$$

where β is a parameter dependent on the droplets' spectral shape and follows the calculation proposed by Peng and Lohmann (2003):

$$\begin{array}{}\text{(5)}& \mathit{\beta }=\mathrm{0.00084}{N}_{\text{cdnc}}+\mathrm{1.22}.\end{array}$$

r_l,vol is the volume-weighted mean cloud droplet radius:

$$\begin{array}{}\text{(6)}& r_{\text{l,vol}}={\left[\left(\text{3LWC}\right)/\left(\mathrm{4}\mathit{\pi }{\mathit{\rho }}_{\text{w}}{N}_{\text{cdnc}}\right)\right]}^{\mathrm{1}/\mathrm{3}},\end{array}$$

where ρ_w is the liquid water density and LWC the cloud liquid water content (g cm⁻³).

Aerosols also exert impacts on precipitation efficiency (Albrecht, 1989), which is taken into account in the parameterization of non-convective cloud processes. There are five processes that convert condensate to precipitate: auto-conversion of liquid water to rain, collection of cloud water by rain, auto-conversion of ice to snow, collection of ice by snow, and collection of liquid by snow. The auto-conversion of cloud liquid water to rain (PWAUT) is dependent on the cloud droplet number concentration and follows a formula that was originally suggested by Chen and Cotton (1987):

$$\begin{array}{}\text{(7)}& \text{PWAUT}=C_{\text{l,aut}}{\widehat{q}}_l^{\mathrm{2}}{\mathit{\rho }}_{\text{a}}/{\mathit{\rho }}_{\text{w}}{\left({\displaystyle \frac{{\widehat{q}}_{\text{l}}{\mathit{\rho }}_{\text{a}}}{{\mathit{\rho }}_w{N}_{\text{ncdc}}}}\right)}^{\mathrm{1}/\mathrm{3}}H\left(r_{\text{l,vol}}-r_{\text{lc,vol}}\right),\end{array}$$

where ${\widehat{q}}_{\text{l}}$ is the in-cloud liquid water mixing ratio, ρ_a and ρ_w are the local densities of air and water, respectively, and C_l,aut is a constant. H(x) is the Heaviside step function with the definition

$$\begin{array}{}\text{(8)}& H\left(x\right)=\left\{\begin{array}{ll}\mathrm{0},& x<\mathrm{0}\\ \mathrm{1},& x\ge \mathrm{0}\end{array}\right..\end{array}$$

r_lc,vol is the critical value of mean volume radius of liquid cloud droplets r_l,vol and is set to 15 µm.

The treatment of aerosol single-scattering (optical) properties (such as mass extinction efficiency, single-scattering albedo, and asymmetric factor) follows the lookup table approach in CAM (Collins et al., 2004). The optics for black and organic carbon, sea salt, and sea salt particles is assumed to be the same as the optics for soot and water-soluble aerosols in the Optical Properties of Aerosols and Clouds (OPAC) dataset (Hess et al., 1998). The optics for dust are derived by Mie calculations for the size distribution represented by each size bin (Zender et al., 2003). Similarly, for sulfate and nitrate particles, the same set of aerosol optical properties for ammonium sulfate is used and is taken from Wang et al. (2008b) with a treatment of aerosol hygroscopicity. The volcanic stratospheric aerosols are assumed to be comprised of 75 % sulfuric acid and 25 % water, as in Hess et al. (1998). For each model year, different aerosol types are assumed to be externally mixed in the calculation of bulk aerosol single-scattering properties that are in turn used in the radiative transfer calculations.

Table 4 Sources of emission data. MOZART2 data denote the standard tropospheric chemistry package for MOZART and contain surface emissions from the EDGAR 2.0 database (Olivier et al., 1996). ACCMIP data are downloaded from the IPCC ACCMIP emission inventory (http://accent.aero.jussieu.fr/ACCMIP.php, last access: 1 May 2014) and they vary from 1850 to 2000, in 10-year steps (Lamarque et al., 2010). CMIP6 data are from https://esgf-node.llnl.gov/search/input4mips/ (last access: 10 March 2019). Anthropogenic emission includes industrial and fossil fuel use, agriculture, ships, etc. Biomass burning includes vegetation fires, such as fuel wood and agricultural burning.

3 Experiment design for the 20th-century climate simulation

AerChemMIP (Collins et al., 2017) is endorsed by CMIP6 for documenting and understanding past and future changes in the chemical composition of the atmosphere and estimating the global-to-regional climate response from these changes. Modeling groups with full chemistry and aerosol models are encouraged to perform all AerChemMIP simulations (Collins et al., 2017). To assess the ability of our model to simulate aerosols (mean and variability), we have followed the historical simulation designed by CMIP6 (Eyring et al., 2016) which is called the “historical” experiment in the Earth System Grid Federation (ESGF). The historical experiment is forced with emissions evolving from 1850 to 2014 that include biomass burning emissions (Van Marle et al., 2017) and anthropogenic and open burning emissions (Hoesly et al., 2018; Feng et al., 2019). O₃ in the historical simulation is an interactive prognostic variable and feeds back on radiation, and the concentrations of other well-mixed greenhouse gases (WMGHGs), e.g., CH₄, N₂O, CO₂, CFC11, and CFC12, are prescribed using CMIP6 historical forcing data (Meinshausen et al., 2017). Although CH₄ and N₂O are prognostic variables in the chemistry scheme (Table 1), their prognostic values at each model step in the historical experiment are replaced by CMIP6 data (Meinshausen et al., 2017) throughout the model domain. The rest of the historical forcing data include (1) yearly global gridded land-use forcing datasets (Hurtt et al., 2011, 2017) and (2) solar forcing (Matthes et al., 2017). All these datasets were downloaded from https://esgf-node.llnl.gov/search/input4mips/ (last access: 10 March 2019). Climate feedback processes that involve changes to the atmospheric composition of reactive gases and aerosols may affect the temperature response to a given WMGHG concentration level.

3.1 Surface emissions

Surface emissions of chemical species from different sources are summarized in Table 4. They include anthropogenic emissions from fossil fuel burning and other industrial activities, biomass burning (including vegetation fires, fuel wood, and agricultural burning), biogenic emissions from vegetation and soils, and oceanic emissions. Most historical emissions from anthropogenic sources (surface, aircraft plus ship) and biomass burning from 1850 to 2014 are CMIP6-recommended data (available at https://esgf-node.llnl.gov/search/input4mips (last access: 10 March 2019). Anthropogenic or biomass burning sources of some tracers which are not included in the CMIP6 dataset (see Table 4), and anthropogenic emission of H₂ and N₂O are from the monthly climatological dataset provided by the MOZART-2 standard package. N₂O is a prognostic variable in BCC-ESM1, but it is replaced by CMIP6 prescribed concentration in the historical run. Other emissions including biomass burning (CH₃COCH₃) and anthropogenic emission (CH₃CHO, CH₃OH, and CH₃COCH₃) are from the Atmospheric Chemistry and Climate Model Intercomparison Project (ACCMIP) emission inventory (http://accent.aero.jussieu.fr/ACCMIP.php, last access: 1 May 2014) covering the period from 1850 to 2010 with 10-year intervals (see Table 4). Monthly lumped emissions of black carbon and organic carbon aerosols from 1850 to 2014 are downloaded from CMIP6-recommended data, but we used 80 % (for BC) and 50 % (for OC) of them in their hydrophobic forms (BC1 and OC1) and the rest in their hydrophilic forms (BC2 and OC2), following the work of Chin et al. (2002).

Five tracers of isoprene (C₃H₈O; ISOP), acetone (CH₃COCH₃), C₂H₄, C₃H₈, and monoterpenes (C₁₀H₁₆) in Table 1 belong to biogenic VOCs. As shown in Table 4, those VOC emissions are calculated online in BCC-ESM1 following the modeling framework of MEGAN2.1 (Guenther et al., 2012) using simple mechanistic algorithms to account for major known processes controlling biogenic emissions. MEGAN2.1 can provide a flexible scheme for estimating 16 tracers of biogenic emissions from terrestrial ecosystems including five VOC emissions used in BCC-ESM1 (Table 4). All the VOC emissions depend on current and past surface air temperature, solar flux, and the landscape types. Their calculation requires global maps of plant functional type (PFT) and leaf area index (LAI), which is a prognostic variable from the land model BCC-AVIM2. The effect of atmospheric CO₂ concentration on isoprene emissions is included. Overall, 10 % of the biogenic monoterpene emissions as calculated online with the MEGAN2.1 algorithm in BCC-AVIM2 is converted to hydrophilic organic carbon (OC2) to account for the formation of secondary organic aerosols following Chin et al. (2002) in this version of BCC-ESM1.

3.2 Volcanic eruptions, lightning, and aircraft emissions

As there is no stratospheric aerosol scheme in BCC-ESM1, concentrations of sulfate aerosol at heights from 5 to 39.5 km, which are of volcanic origin, are directly prescribed using the CMIP6-recommended data (Thomason et al., 2018) from 1850 to 2014. The effects of surface SO₂ emissions from volcanic eruption on the variation in SO₂ in the atmosphere and then on the variation in tropospheric ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ concentration are considered, and the SO₂ emissions from 1850 to 2014 are downloaded from the IPCC ACCMIP emission inventory (http://accent.aero.jussieu.fr/ACCMIP.php, last access: 1 May 2014). Aircraft emissions are provided for NO₂, CO, CH₄, NH₃, NO, SO₂, and aerosols of OC and BC (Table 1). The emissions of NO from lightning are calculated online in BCC-AGCM3-Chem following the parameterization in MOZART2, and the globally averaged mean during the period of 1850 to 2014 is 5.19 Tg (N) yr⁻¹, which is in agreement with observations within the range of 3 to 6 Tg (N) yr⁻¹ (Martin et al., 2002). The lightning frequency depends strongly on the convective cloud top height, and the ratio of cloud-to-cloud versus cloud-to-ground lightning depends on the cold cloud thickness from the level of 0^∘ to the cloud top (Price and Rind, 1992).

Figure 1 The time series of global and annual mean of (a) net energy budget at the top of the atmosphere (W m⁻²), (b) near-surface air temperature (K), and (c) sea surface temperature (K) in the last 450 years of the piControl simulation.

Figure 2 Same as in Fig. 1, but for the global burdens of (a) SO₂, (b) DMS, (c) OH, and (d–h) different aerosols in the troposphere (below 100 hPa). Units are teragrams.

3.3 Upper boundary of the atmosphere

As no stratospheric chemistry is included in the present version of BCC-AGCM3-Chem, it is necessary to ensure a proper distribution of chemically active stratospheric species. Concentrations of different tracers (O₃, CH₄, N₂O, NO, NO₂, HNO₃, CO, and N₂O₅) at the top two layers of the model are set to prescribed monthly climatological values, and concentrations from below the top two layers to the tropopause are relaxed at a relaxation time of 10 d towards the climatology. Climatological values of NO, NO₂, HNO₃, CO, and N₂O₅ at the top two layers are extracted from the MOZART2 data package available at the website (https://www2.acom.ucar.edu/gcm/mozart-4, last access: 1 February 2018), originating from the Study of Transport and Chemical Reactions in the Stratosphere (STARS; Brasseur et al., 1997). Concentrations for the other tracers (O₃, CH₄, and N₂O) at the top two model layers are the zonally averaged and monthly values from 1850 to 2014 derived from the CMIP6 data package.

Figure 3 Global distributions of annual mean mass burdens of sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$; a, f), organic carbon (OC; b, g), black carbon (BC; c, h), dust (d, i), and sea salt (e, j) aerosols in the whole atmospheric column. Panels (a–e) show the mean averaged for the last 100 years of BCC-ESM preindustrial piControl simulations, and (f–j) show the CMIP5 recommended aerosol concentrations in the year 1850 (the website at IIASA http://tntcat.iiasa.ac.at/RcpDb/, last access: 10 January 2012). Units: milligram per square meter.

3.4 The preindustrial model states

The preindustrial state of BCC-ESM1 is obtained from a piControl simulation of over 600 years in which all forcings including emissions data are fixed at 1850 conditions. The initial state of the piControl simulation itself is obtained through individual spin-up runs of each component of BCC-ESM1 in order for the piControl simulation to run stably and fast to reach its equilibrium. Figure 1a–c show the time series of global yearly means of the net energy budget at the top of the atmosphere (TOA), near-surface air temperature (TAS), and sea surface temperature (SST) from the piControl simulation for the last 450 years. It shows that the surface climate in BCC-ESM1 nearly reaches its equilibrium after 600 years of piControl simulation. The whole system in BCC-ESM1 fluctuates around +0.7 W m⁻² net energy flux at TOA without an obvious trend in 450 years (Fig. 1a). This level of TOA energy imbalance is close to the average imbalance (1.0 W m⁻²) among CMIP5 models (Wild et al., 2013). It means that there exists a surplus energy of +0.7 W m⁻² obtained by the whole system in BCC-ESM1, but this does not cause remarkable climate drift. The global mean TAS and SST stay around 288.1 K (Fig. 1b) and 295.05 K (Fig. 1c), respectively. During the last 450 years, there are (±0.2 K amplitude of TAS and SST) oscillations of centennial scale for the whole globe (Fig. 1b and c), which are certainly caused by internal variation in the system.

Figure 2a–c show the time series of global annual total burdens of SO₂, DMS, and OH in the troposphere (integrated from the surface to 100 hPa) in the last 450 years of the piControl simulation. Without any anthropogenic source, the SO₂ amount in the troposphere stays almost at the level of 0.0868 Tg in the 450 years of the piControl simulation. Tropospheric DMS varies around the value of 0.116 Tg. Tropospheric OH, as an important gas species oxidizing SO₂ to form ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ (Table 2), keeps at a stable level in the atmosphere. ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ also remains at a stable level of 0.556 Tg in the atmosphere in the whole period of the piControl simulation (Fig. 2d). The amounts of BC and OC in the troposphere vary around 0.0395 Tg and 0.275 Tg (Fig. 2e–f), respectively. Dust and sea salt aerosols are at the level of 22 Tg and 11.7 Tg (Fig. 2g–h), respectively. All those data are close to the global mean concentrations of 0.604 Tg ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, 0.046 Tg BC, 0.30 Tg OC, 22.18 Tg dust, and 11.73 Tg sea salt in 1850, which are estimated based on the CMIP5 prescribed data in 1850 (Lamarque et al., 2010).

Figure 3 shows the global spatial distributions of annual mean sulfate, organic carbon, black carbon, dust, and sea salt aerosols in the whole atmospheric column averaged for the last 100 years of the piControl simulation of BCC-ESM. We can compare them with CMIP5 recommended concentrations in the year 1850, regarded as the reference state at the preindustrial stage. At that time, there are fewer anthropogenic or biomass SO₂ emissions; the ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ over land is evidently smaller than that over oceans, especially over the tropical Pacific and Atlantic oceans, where DMS can be oxidized to SO₂ and then form ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$. There are several centers of high values of black carbon and organic carbon in East and South Asia, Europe, southeast America, and in the tropical rain forests in Africa and South America. They mainly result from biomass burning including vegetation fires, fuel wood, and agricultural burning. Dust aerosols are mainly distributed in North Africa, Central Asia, North China, and Australia, where arid and semiarid areas are located. Dust emitted from the Sahara can be transported to the tropical Atlantic by easterly wind. The sea salt aerosols are mainly distributed over the midlatitude southern oceans, the tropical southern Indian Ocean, and the tropical northern Pacific Ocean, where wind speeds near the sea surface are strong. As shown in Fig. 3, all the spatial distribution patterns of CMIP5-derived sulfate, black carbon, organic carbon, dust, and sea salt aerosols (Lamarque et al., 2010) are well simulated in BCC-ESM1. There are high spatial correlation coefficients (0.76 for sulfate, 0.77 for black carbon, 0.77 for organic carbon, 0.94 for dust, and 0.94 for sea salt) between CMIP5 data and BCC-ESM1 simulations. Relatively lower relations for sulfate, black carbon, and organic carbon are possibly due to different anthropogenic emission sources being used in BCC-ESM1 and to create CMIP5 data. Dust and sea salt belong to natural aerosols and depend on the land and sea surface conditions, so their spatial distributions are easy to capture and have relatively higher correlations between CMIP5 data and BCC-ESM1 simulations.

Figure 4 Zonal mean of the yearly mean concentration of the ozone column in the troposphere below 300 hPa to the ground from 1871 to 1999 for (a) BCC-ESM1 and (b) CMIP6 data. Unit: DU.

Figure 5 Vertical profiles of annual mean ozone concentrations from observations averaged for 2010–2014 in nine regions (black) and from the BCC-ESM1 simulations (red). The observations are derived from 41 global WOUDC sites.

4 Evaluation of O₃ and aerosol simulations in the 20th century

The rate of sulfate formation is dependent on the levels of oxidants in the troposphere. O₃ is an important oxidant. So, the evaluation of simulated tropospheric O₃ is helpful to understand the aerosol simulations. BCC-ESM1 is driven by most of the CMIP6-recommended emission data. As shown in Fig. 4, the zonal distributions of the total amounts of tropospheric O₃ below 300 hPa to the ground and their changes with time from 1850 to 2014 from the CMIP6-recommend dataset (Table 4) are well simulated by BCC-ESM1. Evident increasing trends since 1850 exist at almost every latitude, especially in the Northern Hemisphere where the contents of tropospheric O₃ are higher than in the Southern Hemisphere.

Figure 5 shows the vertical profiles of O₃ simulations with a comparison to global ozonesonde observations averaged for the monthly data over 2010–2014 from the World Ozone and Ultraviolet Radiation Data Centre (WOUDC; http://woudc.org/data.php, last access: 24 September 2019) in nine regions, which are averaged from 41 global WOUDC sites. The details of WOUDC data may be found in Lu et al. (2019). As shown in Fig. 5, BCC-ESM1 captures the observed ozone vertical structure well in all regions. In the lower and middle troposphere (i.e., below 6 km), the model typically shows positive bias within 5 ppbv for the Southern Hemisphere and 10 ppbv for the northern midlatitudes, similar to those simulated by many other global atmospheric chemical models (Young et al., 2013, 2018). The model has larger ozone overestimation in the upper troposphere and stratosphere in most regions, at least partly due to the use of prescribed stratospheric ozone as upper boundary conditions and/or errors in modeling ozone exchange between the stratosphere and the troposphere. The global tropospheric ozone burden derived from our simulation is 335 Tg averaged over 2010–2014, consistent with a recent assessment from multi chemistry models (Young et al., 2018).

Figure 6 Global annual anthropogenic, natural, and total emissions of SO₂, organic carbon (OC), and black carbon (BC) in the BCC-ESM1 historical simulation. All the biomass burning emissions are included in natural emissions in (a–c). Units: teragrams per year.

4.1 Global aerosol trends

Figure 6a–c show the time series of global total emissions of SO₂, OC, and BC to the atmosphere from natural and anthropogenic sources. Emissions of SO₂ are largely due to industrial production. From 1850 to 1915, SO₂ emissions increased year by year as the industrial revolution intensified and expanded. But from 1915 to 1945, the increased trend in SO₂ emissions became slower as the First and the Second World Wars broke out. After that period, with growing industrial production, SO₂ emissions increased again and reached a maximum around the end of 1970s. With a substantial decrease in SO₂ emissions in Europe and the United States, the global SO₂ emissions has been decreasing since the 1980s despite the rapid increase in SO₂ emissions in South and East Asia as well as in developing countries in the Southern Hemisphere in recent years (Liu et al., 2009). The OC and BC emissions have substantially increased since the 1950s just after the Second World War. The global total OC emission in 2010 was nearly twice as high as that in the preindustrial period (before the year 1850) and increased by 18 Tg yr⁻¹. Anthropogenic black carbon emissions increased from 1 Tg yr⁻¹ in 1850 to nearly 8 Tg yr⁻¹ in 2010.

Figure 7 The time series of global yearly amounts of (a) SO₂ and DMS and (b–f) aerosols in the whole atmosphere column from the CMIP6 historical simulations of BCC-ESM1 (black lines) and the CMIP5-recommended aerosol masses (red lines). The yearly CMIP5 data are interpolated from the time series at a 10-year interval. Units: teragrams.

Anthropogenic SO₂, OC, and BC emissions strongly affect the variations in atmospheric concentrations of sulfate, OC, and BC. The global $\mathrm{0.5}{}^{\circ }\times \mathrm{0.5}{}^{\circ }$ gridded data of CMIP5-recommended aerosol masses with a 10-year interval from 1850 to 2000 (Lamarque et al., 2010) provides an important reference to evaluate the aerosol simulations in BCC-ESM1. As shown in Fig. 7b–f, the annual total aerosol burdens of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, OC, and BC in the whole atmosphere column as simulated by the BCC-ESM1 20th-century historical simulation are generally consistent with the values derived from CMIP5-recommended aerosols concentrations. Due to increasing SO₂ emissions from 1850 to the present day (Fig. 6), the global SO₂ burden in the atmosphere increased from 100 Tg in the 1850s to 200 Tg in the 1980s (Fig. 7a) and has a high correlation coefficient of 0.996 with the anthropogenic emissions (Fig. 6a), as the lifetime of SO₂ is short. The burden directly followed the emission. DMS in the atmosphere is oxidized by OH and NO₃ to form SO₂ (Table 2). Its natural emissions from oceans from 1850 to 2010 in the model are the climatological monthly means (Dentener et al., 2006) from MOZART2 data package. As shown in Fig. 7a, the global amount of DMS in the whole atmosphere was about 0.12 Tg during 1850–1900 and decreased to 0.055 Tg in 2010. This decrease trend maybe partly results from the speeded rate of DMS oxidation with global warming, and the loss of DMS gradually exceeds the source of ocean DMS emission to cause a net loss of DMS in the atmosphere since the 1910s. Largely driven by SO₂ anthropogenic emissions, the sulfate burden shows three different stages from 1850 to the present. In the first period from the 1850s to the 1900s, the sulfate burden had a weak linear increase. It increased significantly in the second stage from the 1910s to the 1940s and then exploded from the 1950s until the middle 1970s and early 1980s. The sulfate burden then remained nearly stable and even showed slight decreases as seen from the CMIP5 data. As for global BC and OC burdens, BCC-ESM1 results have shown continuous increases since the 1850s, especially from 1950 to the present. From the 1910s to the 1940s, the CMIP5 data showed a slight decrease in BC and OC burdens in the atmosphere.

The dust and sea salt aerosols in the atmosphere are largely determined by the atmospheric circulations and states of the land and ocean surface. We can see that the global dust burden in the atmosphere showed an evident increase from 1980 to 2000, which could be partly caused by evident global warming since 1980 and increasing soil dryness resulting in more surface dust to be released into the atmosphere. Their details will be explored in second part.

4.2 Global aerosol budgets

We further evaluate global aerosol budgets by comparing a 10-year average of BCC-ESM results from 1990 to 2000 with various studies for sulfate, BC, OC, sea salt, and dust. Their annual total emissions, average atmospheric mass loading, and mean lifetimes are listed in Tables 5 and 6. It is worth emphasizing that the global mean total source and sink for each type of aerosols in BCC-ESM1 are almost balanced.

Table 5 Global budgets for DMS, SO₂, and sulfate in the period of 1991 to 2000. Units for sources and sinks: teragrams (S) per year; for burden: teragrams; for lifetime: days.

^a Liu et al. (2005). ^b Textor et al. (2006). ^c Values derived from CMIP5 prescribed aerosol masses averaged from 1991 to 2000. ^d Liu et al. (2012). ^e Matsui and Mahowald (2017). ^f Tegen et al. (2019). Values of DMS, SO₂, and sulfate burdens from ^dLiu et al. (2012) are transferred from teragrams S to teragrams (species) for unit consistency.

The global DMS emission from the ocean is 27.4 Tg (S) yr⁻¹ in BCC-ESM. This emission in BCC-ESM is nearly balanced by the gas-phase oxidation of DMS to form SO₂. The DMS burden is 0.12 Tg with a lifetime of 0.78 d, which is within the range of other models reported in the literature. As shown in Table 5, the total SO₂ production averaged for the period of 1991 to 2000 is 76.93 Tg (S) yr⁻¹. A rate of 13.2 Tg (S) yr⁻¹ (about 17 %) SO₂ is produced from the DMS oxidation; only 0.1 Tg (S) yr⁻¹ SO₂ is from airplane emissions to the atmosphere, and the rest (63.63 Tg (S) yr⁻¹, almost 82.7 %) is from anthropogenic activities and volcanic eruption at the surface. The amount of SO₂ produced from the DMS oxidation is in the range of other works (10.0 to 24.7 Tg (S) yr⁻¹) reported in Liu et al. (2005). All the SO₂ production is balanced by SO₂ losses by dry and wet deposition and by gas- and aqueous-phase oxidation. Half of its loss (38.74 Tg (S) yr⁻¹) occurs via its aqueous-phase oxidation to form sulfate. Other losses through dry and wet depositions and gas-phase oxidation to form ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ are also important (Table 2). All the sinks are in the range from the literature (Liu et al., 2005). The global burden of SO₂ in the atmosphere is 0.48 Tg with a lifetime of 1.12 d, consistent with values in literature (Liu et al., 2005).

Sulfate aerosol is mainly produced from aqueous-phase SO₂ oxidation (38.73 Tg (S) yr⁻¹) and partly from gaseous-phase oxidation of SO₂ (10.32 Tg (S) yr⁻¹) and is largely lost by wet scavenging (49.06 Tg (S) yr⁻¹). The total ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$ production in BCC-ESM is at the lower range of values in other models reported in Textor et al. (2006). Its global burden is 1.89 Tg and the lifetime is 4.69 d, which is within the range of 1.71 to 2.43 Tg and 3.3 to 5.4 d in the literature (Textor et al., 2006; Liu et al., 2012, 2016; Matsui and Mahowald, 2017; Tegen et al., 2019; the value is derived from CMIP5 data).

Sources of BC and OC are mainly from anthropogenic emissions. Based on the CMIP6 data, there are, on average, 7.22 Tg yr⁻¹ BC and 13.91 Tg yr⁻¹ OC from fossil and bio-fuel emissions and 18.38 Tg yr⁻¹ OC from natural emission during the period of 1991 to 2000. Most of these are scavenged through convective and large-scale rainfall processes. The rest returns to the surface by dry deposition. The simulated global BC and OC burdens are 0.13 and 0.62 Tg, respectively (Table 6), all close to the values of 0.114 Tg BC and 0.69 Tg OC derived from the CMIP5 data and within the range of 0.11–0.26 Tg BC (Textor et al., 2006; Matsui and Mahowald, 2017; Tegen et al., 2019) and less than the values of 1.25–2.2 Tg OC in other literature (Textor et al., 2006; Tegen et al., 2019). The simulated BC and OC lifetimes are 6.6 and 5.0 d, respectively, and are close to the recent values of 5.0–7.5 d BC and 5.4–6.6 d OC in the literature (Matsui and Mahowald, 2017; Tegen et al., 2019).

Table 6 Same as Table 5, but for global budgets for black carbon, organic carbon, dust, and sea salt. Units for sources and sinks: teragrams per year; for burden: teragrams; for lifetime: days.

Notes: References denote ^a for Liu et al. (2005), ^b for Textor et al. (2006), ^c derived from CMIP5 prescribed aerosol masses averaged from 1991 to 2000, ^d for Liu et al. (2012), ^e for Liu et al. (2016), ^f for Ginoux et al. (2001), ^g for Matsui and Mahowald (2017), and ^h for Tegen et al. (2019).

The emissions of dust and sea salt are mainly determined by winds near the surface. The annual total dust emission in BCC-ESM1 is 2592 Tg yr⁻¹, higher than the AeroCom multi-model mean (1840 Tg yr⁻¹; Textor et al., 2006) but comparable to other studies (Chin et al., 2002; Liu et al., 2012; Matsui and Mahowald, 2017). The average dust loading is 22.93 Tg, lower than the value of 35.9 Tg in Ginoux et al. (2001) but slightly higher than the value of 20.41 Tg derived from CMIP5 data. The average lifetime for dust particles is 3.23 d, which is shorter than the AeroCom mean (4.14 d) and the value of 3.9 d in a recent study (Matsui and Mahowald, 2017). The simulated sea salt emission is 4667.2 Tg yr⁻¹, slightly lower than the simulated value in Liu et al. (2012) and substantially lower than the AeroCom mean (16 600 Tg yr⁻¹; Textor et al., 2006). The simulated sea salt burdens are 11.89 Tg and close to the CMIP5 data. Their averaged lifetimes are 0.93 d and close to the value in the recent of Matsui and Mahowald (2017) but longer than the AeroCom mean (0.41 d; Textor et al., 2006).

4.3 Global aerosol distributions in the present day

Figures 8–12 show December–January–February (DJF) and June–July–August (JJA) mean column mass concentrations of sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$), OC, BC, dust, and sea salt aerosols averaged for the period of 1991–2000. Here, BCC-ESM1 simulated results are compared with the CMIP5-recommended data for the same period. Unlike the preindustrial level of sulfate shown in Fig. 2, sulfate concentrations in the present day (Fig. 8) are strongly influenced by anthropogenic emissions and have maximum concentrations in the industrial regions (e.g., East Asia, Europe, and North America). Their seasonal variations are distinct and are characterized by high concentrations in boreal summer and low concentrations in boreal winter. These spatial distributions simulated by BCC-ESM1 show good consistency with the CMIP5 data, with spatial correlation coefficients in DJF and JJA reaching 0.92 and 0.83 (Fig. 13), respectively. In DJF, the deviation in the spatial pattern in BCC-ESM1 from the CMIP5 data is less, but it is more in JJA (Fig. 13).

Figure 8 December–January–February (DJF; a, c) and June–July–August (JJA; b, d) mean sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$) aerosol column mass concentrations averaged for the period of 1971–2000. Panels (a, b) show the historical simulations of BCC-ESM1 and (c, d) the CMIP5-recommended data. Units: milligram per square meter.

Figure 9 The same as in Fig. 8, but for organic carbon (OC) aerosol column mass concentrations. Units: milligram per square meter.

Figure 10 The same as in Fig. 8, but for black carbon (BC) aerosol. Units: milligram per square meter.

Figure 11 The same as in Fig. 8, but for dust aerosol. Units: milligram per square meter.

Figure 12 The same as in Fig. 8, but for sea salt (SSLT) aerosol. Units: milligram per square meter.

Unlike sulfate whose maximum concentrations are mainly distributed between 60^∘ N and the Equator, peaking concentrations of BC and OC as shown in Figs. 9 and 10 are located near the tropics in the biomass burning regions (e.g., the Maritime Continent, Central Africa, South America), and their seasonal variations from DJF to JJA are evidently weaker than those of sulfate except in South America. In boreal summer, there are centers of high values in the industrial regions in the Northern Hemisphere midlatitudes (i.e., East Asia, South Asia, Europe, and North America). These main features of spatial and seasonal variations in CMIP5 data are well captured by BCC-ESM1, and the BCC-ESM1 vs. CMIP5 spatial correlation coefficients (Fig. 13) are 0.90 (OC in DJF), 0.91 (BC in DJF), 0.91 (OC in JJA), and 0.92 (BC in JJA). There are fewer deviations in spatial pattern for OC in DJF and JJA but larger deviation for BC from CMIP5 data (Fig. 13).

As shown in Fig. 11, dust concentrations in the atmosphere show the largest values over strong source regions such as northern Africa, southwest and Central Asia, and Australia and over their outflow regions such as the Atlantic and the western Pacific. In DJF, the CMIP5 data show centers of high concentrations over East Asia and central North America, but both centers are missing in BCC-ESM1. However, these two high-value centers in the CMIP5 data may not be true, since frozen soils in these areas in winter lead to unfavorable conditions for soil erosion by winds. The spatial correlation coefficients between CMIP5 and BCC-ESM1 remain high: 0.95 in JJA and 0.88 in DJF (Fig. 13). Small deviations in spatial pattern for dust simulations in BCC-ESM1 show a lower magnitude of dust maximums against CMIP5 data (Fig. 13).

Figure 13 Taylor diagram for the global aerosol climatology (1971–2000) of sulfate, organic carbon, black carbon, dust, and sea salt averaged for December–January–February (DJF) and June–July–August (JJA) and annually. The radial coordinate shows the standard deviation in the spatial pattern, normalized by the observed standard deviation. The azimuthal variable shows the correlation of the modeled spatial pattern with the observed spatial pattern. Analysis is for the whole globe. The reference dataset is the CMIP5 prescribed dataset.

As shown in Fig. 12, high sea salt concentrations are generally found over the storm track regions over the oceans, e.g., the midlatitudes in the northern oceans in DJF and the Southern Ocean in JJA, where wind speeds and thus sea salt emissions are higher. In addition, there is a belt of high sea salt concentrations in the subtropics of both hemispheres where precipitation scavenging is weak. Their spatial distributions in BCC-ESM1 are consistent with the CMIP5 data with correlation coefficients of 0.92 in JJA and 0.90 in DJF (Fig. 13). The spatial deviations in sea salt are much closer to CMIP5 data than those of sulfate, OC, BC, and dust distributions (Fig. 13).

Figure 14 Latitude–pressure distributions of zonally averaged annual mean sulfate, organic carbon, black carbon, dust, and sea salt aerosol concentrations for the period of 1971–2000. Left panels show the CMIP6 historical simulation of BCC-ESM1 and right panels the CMIP5 recommendation data. Units: micrograms per cubic meter.

Figure 15 Scatterplots showing observed versus simulated multi-year averaged annual mean sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$), organic carbon (OC), and black carbon (BC) mixing ratios at the IMPROVE and EMEP network sites. Observations are averages over the available years 1990–2005 for IMPROVE sites and 1995–2005 for EMEP sites. Simulated values are those at the lowest layer of BCC-ESM1.

Figure 16 Comparison of modeled black carbon (BC) aerosol (red lines) with observations from HIPPO aircraft campaigns over the Pacific Ocean (black symbols; bars represent the full data range). Observations from different HIPPO campaigns were averaged over 5^∘ latitude bins and three different altitude bands (left column: 1–2 km; middle column: 4–5 km; and right column: 7–8 km) along the flight track over the Pacific Ocean. Model results were sampled along the flight track and then averaged over the abovementioned regions for comparison.

Figure 14 shows vertical distributions of zonally averaged annual mean concentrations of sulfate, organic carbon, black carbon, dust, and sea salt aerosols in the period of 1991–2000. Both BCC-ESM1 and CMIP5 results show that strong sulfur, OC, and BC emissions in the industrial regions of the Northern Hemisphere midlatitudes can rise upward and be transported towards the North Pole in the mid-troposphere to upper troposphere. Most OC, BC, and dust aerosols are confined to below 500 hPa, while sulfate can be transported to higher altitudes. Sea salt aerosols are mostly confined to below 700 hPa, as the particles are large in size and favorable for wet removal and gravitational settling towards the surface. It can be seen that BCC-ESM1 tends to simulate less upward transport of aerosols than the CMIP5 data, likely reflecting the omission of deep convection transport of tracers in BCC-ESM1.

The CMIP5 data used here are mainly from model simulations. We will further evaluate the BCC-ESM1 model results with ground observations. Annual mean ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, BC, and OC aerosol observations from the Interagency Monitoring of Protected Visual Environments (IMPROVE) sites over 1990–2005 in the United States (http://vista.cira.colostate.edu/IMPROVE/, last access: 5 October 2019) and from the European Monitoring and Evaluation Programme (EMEP) (http://www.emep.int, last access: 5 October 2019) sites over 1995–2005 are used. As shown in Fig. 15a and b, the BCC-ESM simulated sulfate concentrations are in general comparable to the EMEP observations in Europe but are systematically about 1 µg m⁻³ higher than the US IMPROVE observations. As for BC, there are large model biases at both European and US sites (Fig. 15c and d); in particular, BCC-ESM overestimates BC concentrations at the IMPROVE sites. The observed OC concentrations are slightly overestimated for IMPROVE sites but systematically underestimated for EMEP sites. Some statistical features for simulated concentrations versus EMEP and IMPROVE observations are listed in Table 7. These comparisons are overall fairly reasonable considering the uncertainties in emissions and the coarse model resolution.

Table 7 Observed versus simulated concentrations of sulfate (${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$), organic carbon (OC), and black carbon (BC) for the regional mean and spatial standard deviation and minimum and maximum values at HIPPO aircraft observations (BC only), and IMPROVE and EMEP network sites and the spatial correlation between observed and simulated multi-year averaged annual means. Simulated values are selected for the same locations and the same valid observation time. The data used are the same as those in Fig. 12.

We then evaluate the simulated BC concentrations from BCC-ESM1 with the HIAPER (High-Performance Instrumented Airborne Platform for Environmental Research) Pole-to-Pole Observations (HIPPO) (Wofsy et al., 2011). The HIPPO campaign provided observations of black carbon concentration profiles over the Pacific Ocean and North America between 2009 and 2011. Following Tilmes et al. (2016), model results here are sampled along the HIPPO flight tracks and then averaged to different latitude and altitude bands for comparison. As shown in Fig. 16, BCC-ESM1 and HIPPO aircraft observations shows reasonable agreement in terms of the spatial distributions and seasonal variations in BC levels. BCC-ESM1 generally reproduces the observed hemispheric gradients of BC, i.e., the larger burden in the NH compared to the SH, consistent with Figs. 10 and 14. The mean value of modeled results along the flight track is 11.1 ng kg⁻¹, comparable to the 8.2 ng kg⁻¹ of the HIPPO observations. The model shows large overestimations of BC observations over the tropics, which is also found in the CAM4-chem global chemical model (Tilmes et al., 2016).

Figure 17 Global distribution of annual mean AOD simulated in BCC-ESM1 compared with the MISR and MODIS data for the year 2008.

Figure 18 Observed versus simulated annual means of AOD at AERONET sites. Each data point represents the mean averaged for available monthly values of AOD. The dot sizes denote the magnitudes of AOD at sites. The spatial correlation is 0.56.

4.4 Aerosol optical properties

Aerosol optical depth (AOD) is an indicator of the reduction in incoming solar radiation (at a particular wavelength) due to scattering and absorption of sunlight by aerosols. In this study, we calculate the AOD at 550 nm for all aerosols including sulfate, BC, organic carbon, sea salt, and dust as the product of aerosol dry mass concentrations, aerosol water content, and their specific extinction coefficients. The total AOD is calculated by summing the AOD in each model layer for each aerosol species using the assumption that they are externally mixed. The AOD observations retrieved from the Moderate Resolution Imaging Spectroradiometer (MODIS) and the Multi-angle Imaging SpectroRadiometer (MISR) over the period of 1997–2003 and from the AErosol RObotic NETwork (AERONET) over the period of 1998–2005 (http://aeronet.gsfc.nasa.gov, last access: 2 August 2019) are used to evaluate the averaged AOD at 550 nm in BCC-ESM. Figure 17 shows averages of MISR and MODIS AOD with corresponding averages from BCC-ESM. The BCC-ESM1 simulated AOD generally captures the spatial distribution of MISR and MODIS retrievals. The model overestimates AOD over east China. It also systematically underestimates the MODIS observations in the Southern Hemisphere but is closer to MISR observations. Figure 18 shows multi-year annual means of BCC-ESM1-simulated AOD values versus observations from AERONET over the period of 1998–2005. The basic pattern of modeled global AOD is similar to that of observations and their spatial correlation reaches 0.56. Large values of AOD are mainly distributed in land continents such as North African, South Asia, East Asia, Europe, and the eastern part of North America. Figure 19a–d present scatterplots of observed versus simulated multi-year monthly mean AOD at those sites of AERONET in Europe, North America, East Asia, and South Asia over the period of 1998–2005. Model-simulated monthly AOD generally agrees with observations within a factor of 2 for most sites. BCC-ESM slightly overestimates the AOD at European and North American sites. In those regions, BCC-ESM also slightly overestimates MODIS and MISR AOD observations (Fig. 17).

Figure 19 Scatterplots of observed versus simulated monthly mean AOD at AERONET sites in Europe, North America, East Asia, and South Asia over the period of 1998–2005.

5 Summary and discussions

This paper presents a primary evaluation of aerosols simulated in version 1 of the BCC-ESM1 with the implementation of the interactive atmospheric chemistry and aerosol based on the newly developed BCC-CSM2. Global aerosols (including sulfate, organic carbon, black carbon, dust, and sea salt) and major greenhouse gases (e.g., O₃, CH₄, N₂O) in the atmosphere can be interactively simulated when anthropogenic emissions are provided to the model. Concentrations of all aerosols in BCC-ESM1 are determined by the processes of advective transport, emission, gas-phase chemical reactions, dry deposition, gravitational settling, and wet scavenging by clouds and precipitation. The nucleation and coagulation of aerosols are ignored in the present version of BCC-ESM1. Effects of aerosols on radiation, cloud, and precipitation are fully included.

We evaluate the performance of BCC-ESM1 in simulating aerosols and their optical properties in the 20th century following CMIP6 historical simulation according to the requirement of the AerChemMIP. It is forced with anthropogenic emissions evolving from 1850 to 2014, but some WMGHGs such as CH₄, N₂O, CO₂, CFC11, and CFC12 are prescribed using CMIP6 prescribed concentrations (to replace prognostic values of CH₄ and N₂O from the chemistry scheme). Both direct and indirect effects of aerosols are considered in BCC-ESM1. Initial conditions of the CMIP6 historical simulation are obtained from a 600-year piControl simulation in the absence of anthropogenic emissions, which captures the preindustrial concentrations of ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, organic carbon (OC), black carbon (BC), dust, and sea salt aerosols well and is consistent with the CMIP5 recommended concentrations for the year 1850. With the CMIP6 anthropogenic emissions of SO₂, OC, and BC from 1850 to 2014 and their natural emissions implemented in BCC-ESM1, the model-simulated ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, BC, and OC aerosols in the atmosphere are highly correlated with the CMIP5-recommended data. The long-term trends in CMIP5 aerosols from 1850 to 2000 are also well simulated by BCC-ESM1. Global budgets of aerosols were evaluated through comparisons of BCC-ESM1 results for 1990–2000 with reports in various contributions to the literature for sulfate, BC, OC, sea salt, and dust. Their annual total emissions, atmospheric mass loading, and mean lifetimes are all within the range of values reported in the relevant literature. Evaluations of the spatial and vertical distributions of BCC-ESM1 simulated present-day ${\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}$, OC, BC, dust, and sea salt aerosol concentrations against the CMIP5 datasets and in situ measurements of surface networks (IMPROVE in the US and EMEP in Europe), and HIPPO aircraft observations indicate good agreement among them. The BCC-ESM1 simulates weaker upward transport of aerosols from the surface to the middle and upper troposphere (with reference to CMIP5-recommended data), likely reflecting a lack of deep convection transport of chemical species in the present version of BCC-ESM1. The AOD at 550 nm for all aerosols including sulfate, BC, OC, sea salt, and dust aerosols was further compared with the satellite AOD observations retrieved from MODIS and MISR and surface AOD observations from AERONET. The BCC-ESM1 model results are overall in good agreement with these observations within a factor of 2. All these comparisons demonstrate the success of the implementation of interactive aerosol and atmospheric chemistry in BCC-ESM1.

This work has only evaluated the ability of BCC-ESM1 to simulate aerosols. The variations in aerosols especially for sulfate are related to other gaseous tracers such as OH and NO₃ (Table 2), which are determined by the MOZART2 gaseous chemical scheme as implemented in BCC-ESM1, and require further evaluation. As the length of the text is limited, the other optical features of aerosols such as extinction coefficients, single-scattering albedo and asymmetry parameters, and even their feedbacks on radiation and global temperature change will be explored in a second part. O₃ is evaluated in this work. Other GHGs such as CH₄ and N₂O concentrations can be simulated when forced with emissions, and their simulations also needs to be evaluated in the future.