Amending the algorithm of aerosol–radiation interactions in WRF-Chem (v4.4)

WRF-Chem (Weather Research and Forecasting model coupled with Chemistry) is widely used to assess regional aerosol radiative feedback. However, in the current version, aerosol optical properties are only calculated in four shortwave bands, while only two of them are used to interpolate optical properties towards the 14 shortwave bands used in the Rapid Radiative Transfer Model for Global Climate Models (RRTMG) scheme. In this study, we use a “Resolved” algorithm to estimate aerosol radiative feedback in WRF-Chem, in which aerosol optical properties are calculated in all 14 shortwave bands. The impacts of changing this calculation algorithm are then evaluated. The simulation results of aerosol optical properties are quite different using the new Resolved algorithm, especially for dust aerosols. The alteration of aerosol optical properties results in considerably different aerosol radiative effects: the dust radiative forcing in the atmosphere simulated by the Resolved algorithm is about 2 times larger than the original “Interpolated” algorithm. The dust radiative forcing at the top of the atmosphere (TOA) simulated by the Interpolated algorithm is negative in the Sahara region, while the Resolved algorithm simulates positive forcing at TOA and can exceed 10 W m⁻² in the Sahara, which is more consistent with previous studies. The modification also leads to changes in meteorological fields due to alterations in radiative feedback effects of aerosols. The near-surface temperature is changed due to the difference in the radiation budget at the bottom of the atmosphere (BOT) and the heating effects of aerosols at the surface. Furthermore, the amendment of the algorithm partially corrects the wind field and temperature simulation bias compared to the reanalysis data. The difference in planet boundary layer height can reach up to ∼100 m in China and ∼200 m in the Sahara, which also results in a greater surface haze. The results show that correcting the estimation algorithm of aerosol radiative effects is necessary in WRF-Chem.

1 Introduction

Aerosol–radiation interactions and their impacts on meteorological processes and the aerosol cycle have been proven to be important (e.g., Ackerman, 1977; Dickerson et al., 1997; Jacobson, 1998; Zhao et al., 2010, 2011, 2012, 2013, 2014; Myhre et al., 2013; Bender, 2020; Bellouin et al., 2020; Huang and Ding, 2021). As the two-way interactions between aerosol and meteorological fields are complex, a fully coupled “online” meteorology–chemistry model is a necessary tool to account for these feedbacks in simulating aerosol concentrations and meteorological fields. WRF-Chem (the Weather Research and Forecasting model coupled with Chemistry) is one of the most widely used atmospheric models that considers aerosol–radiation interactions for investigating regional aerosol life cycles and radiative impacts (e.g., Zhao et al., 2010, 2011, 2013, 2014; Jiang et al., 2012; Ding et al., 2013; Wu et al., 2013; Gao et al., 2014; Chen et al., 2014; Zhong et al., 2016; Liu et al., 2016; Huang et al., 2016; Petäjä et al., 2016; Du et al., 2020, 2023; Zhang et al., 2020; Wang et al., 2022; Chen et al., 2022; Sharma et al., 2023; Wei et al., 2023). WRF-Chem is capable of performing regional-scale simulations at high spatial resolution. This allows for a detailed representation of aerosol and radiation processes at regional scales. By incorporating the interactions between aerosols and radiation, WRF-Chem can provide insights into the impacts of aerosols on regional weather patterns, climate, air quality, and the energy balance. Therefore, it is crucial for the modeling community to appropriately simulate aerosol optical properties and aerosol–radiation interactions in WRF-Chem.

In WRF-Chem, aerosol optical properties (i.e., aerosol optical depth (AOD), single-scattering albedo (SSA), and the asymmetry factor) for shortwave bands are first computed at four wavelengths of 300, 400, 600, and 999 nm following the method described in previous studies (Fast et al., 2006; Barnard et al., 2010; Zhao et al., 2013; see Sect. 2.2 for more details). Afterwards, these aerosol optical properties are used in radiative transfer schemes such as the Rapid Radiative Transfer Model for Global Climate Models (RRTMG; Mlawer et al., 1997; Iacono et al., 2000; Zhao et al., 2011). RRTMG calculates radiative fluxes and heating rates in 14 shortwave bands. However, due to the aerosol optical properties only being shortwave calculated for four spectral bands (as mentioned above), the RRTMG scheme interpolates the values at these 4 wavelengths to create the 14 wavelengths used here. For AOD, the scheme obtains the values for all 14 shortwave bands using the Ångström exponent (AE, Ångström, 1929) based on AOD at 400 and 600 nm. For SSA and the asymmetry factor, simple linear interpolation is applied (Barnard et al., 2010; Zhao et al., 2011). These interpolation methods for taking aerosol optical properties from 2 bands and creating 14 bands may lead to significant errors in estimating aerosol radiative forcing and subsequently simulating aerosol radiative feedback on meteorological fields.

Therefore, in this study we amend the aforementioned approach by directly calculating aerosol optical properties for all 14 RRTMG shortwave bands and examine the difference between the new algorithm and the original algorithm when simulating multi-band aerosol optical properties, radiative forcing, and aerosol impacts on meteorological fields. The paper is organized as follows. Section 2.1 briefly introduces WRF-Chem, while Sect. 2.2 describes the algorithms used for simulating aerosol–radiation interactions in WRF-Chem. The differences seen when simulating aerosol optical properties and radiative impacts with the two algorithms are investigated in Sect. 3. The conclusion and summary are in Sect. 4.

2 Methodology

2.1 WRF-Chem

WRF-Chem is a version of the WRF (Skamarock et al., 2021) that simulates trace gases and particulates simultaneously with the meteorological fields (Grell et al., 2005). The Model for Simulating Aerosol Interactions and Chemistry (MOSAIC) aerosol model (Zaveri et al., 2008) and photochemical Carbon Bond Mechanism version Z (CBM-Z) (Zaveri and Peters, 1999) implemented by Fast et al. (2006) into WRF-Chem, which include complex treatments of aerosol radiative properties and photolysis rates, are used in this study. Since this study only focuses on the amendment of the current algorithm of aerosol optical properties and its radiative feedback in WRF-Chem, more details about the physics and chemistry schemes in WRF-Chem are not given here and can instead be found in previous studies (Zhao et al., 2011, 2013).

2.2 Amendment of the aerosol–radiation interaction algorithm in WRF-Chem

In current (v4.4) and previous versions of WRF-Chem, aerosol optical properties such as extinction, SSA, and the asymmetry factor for scattering are computed as a function of wavelength and three-dimensional position. Currently, the methodology described by Ghan et al. (2001) is applied to compute the extinction efficiency Q_e and the scattering efficiency Q_s in WRF-Chem. In the model, the full Mie calculation is performed only once to obtain a table of seven sets of Chebyshev expansion coefficients, and later the full Mie calculations are skipped, while Q_e and Q_s are calculated using bilinear interpolation over the Chebyshev coefficients stored in the table. The detailed method of the computation of aerosol optical properties in the model is similar to the description given in Fast et al. (2006), Barnard et al. (2010), and Zhao et al. (2013). The Optical Properties of Aerosols and Clouds (OPAC) dataset (Hess et al., 1998) is used for the shortwave (SW) and longwave (LW) refractive indices of dust aerosols. Radiative feedback of aerosols is coupled with RRTMG (Mlawer et al., 1997; Iacono et al., 2000) for both longwave (LW) and shortwave (SW) radiation as in Zhao et al. (2011). AOD and direct radiative forcing of aerosols are diagnosed following the methodology of Zhao et al. (2013). In this methodology, the calculation of aerosol optical properties and the radiative transfer scheme is performed multiple times with the mass of one or more aerosol species (i.e., the mass of an individual or a group of aerosol species), and its associated water aerosol mass is also removed from the calculation each time. After this diagnostic procedure, the optical properties (e.g., AOD) and direct radiative forcing for aerosols can be estimated by subtracting the optical properties and direct radiative forcing from the diagnostic iterations from those estimated following the standard procedure for all the aerosol species. It can be described as follows:

$$\begin{array}{}\text{(1)}& {\displaystyle }{\text{AOD}}_{\left[\mathrm{species}i\right]}={\text{AOD}}_{\left[\mathrm{all}\mathrm{species}\right]}-{\text{AOD}}_{\left[\mathrm{without}\mathrm{species}i\right]}\text{(2)}& {\displaystyle }\begin{array}{rl}{\text{AForcing}}_{\left[\mathrm{species}i\right]}& ={\text{Forcing}}_{\left[\mathrm{all}\mathrm{species}\right]}\\ & -{\text{Forcing}}_{\left[\mathrm{without}\mathrm{species}i\right]}\end{array}.\end{array}$$

Currently, the aerosol optical properties for shortwave spectral bands (roughly 200–4000 nm) are only calculated for four shortwave bands centered at 300, 400, 600, and 999 nm following the method above. When coupled with the RRTMG shortwave radiative transfer scheme (RRTMG-SW) in WRF-Chem, the aerosol optical properties (i.e., AOD, SSA, and the asymmetry factor) need to be interpolated from the values at the 4 wave bands to those values at the 14 wave bands (from 232 to 3462 nm) used in the RRTMG shortwave scheme. The interpolation of AOD is based on the AE that is derived from the values at 400 and 600 nm:

$$\begin{array}{}\text{(3)}& \mathit{\alpha }=-\ln \left({\displaystyle \frac{{\mathit{\tau }}_{\mathrm{600}}}{{\mathit{\tau }}_{\mathrm{400}}}}\right)/\ln \left({\displaystyle \frac{\mathrm{600}}{\mathrm{400}}}\right),\end{array}$$

where α is the Ångström exponent (AE), while τ₄₀₀ and τ₆₀₀ are the AOD over 400 and 600 nm bands, respectively.

Given the AE, AOD over the wave band centered at λ_i is calculated as follows:

$$\begin{array}{}\text{(4)}& {\mathit{\tau }}_i={\mathit{\tau }}_{\mathrm{600}}{\left({\displaystyle \frac{{\mathit{\lambda }}_i}{\mathrm{600}}}\right)}^{-\mathit{\alpha }},(i=\mathrm{1},\mathrm{2},\mathrm{\dots },\mathrm{14}).\end{array}$$

The SSA and asymmetry factor values for other wave bands are linearly interpolated from the values at the four wave bands. This current algorithm of calculating optical properties and coupling with the shortwave radiation transfer scheme is referred to as the “Interpolated” algorithm for the rest of this paper.

In this study, we implement the “Resolved” algorithm, wherein the aerosol optical properties over the 14 shortwave bands are calculated and coupled with RRTMG-SW directly. While our study focused on the sectional representation of aerosols, our amendment of the algorithm is theoretically valid for bulk, sectional, and modal representations in WRF-Chem. However, we have not tested the differences between the Interpolated and Resolved methods for bulk and modal aerosol representations in this study. The difference between Interpolated algorithm and Resolved algorithm is defined as the bias due to the interpolation of aerosol optical properties for the radiation transfer scheme. The biases of the Interpolated algorithm and its impacts are investigated.

Figure 1 Simulation domains in (a) China and (b) the Sahara. The dashed boxes in panels (a) and (b) represent regions dominated by anthropogenic aerosols and dust aerosols, respectively. Spatial distributions of topography height are also shown. Locations of selected AERONET stations over dust-dominant areas in the Sahara are also denoted in panel (b).

2.3 Numerical experiments

In this study, four sets of experiments are conducted over two domains, as shown in Fig. 1. One covers China (8–55° N, 59–146° E) and represents a region with complex aerosol sources including large anthropogenic aerosol mass loading and natural dust (in the northwestern part of the domain). The other covers the Sahara (3.5° S–42° N, 24° W–44° E) and represents the region with the largest natural dust aerosol mass loading in the world. In addition, Fig. 1a and b also delineate regions dominated by anthropogenic aerosols and dust aerosols, respectively, using dashed boxes (referred to as anthro-dominant and dust-dominant regions, respectively). Over both domains, two sets of experiments, one with the Interpolated algorithm and the other with the Resolved algorithm, are conducted. The simulations are performed at 50 km × 50 km horizontal resolution with 120×100 grid cells and 40 vertical layers up to 100 hPa. The model time step is set at 150 s, and the aerosol optical properties are updated every 30 min in the model. The experiments are conducted for January and July 2015, representing boreal winter and summer, which ran from 25 December 2014 to 31 January 2015 and from 25 June 2015 to 31 July 2015, respectively. Only the results during January and July are used in the analysis to minimize the impact from the initial chemical conditions. In this study, all analysis results represent monthly averages for January and July 2015 unless otherwise stated. The meteorological initial conditions are derived from the European Centre for Medium-Range Weather Forecasts (ERA5) reanalysis dataset at approximately 25 km horizontal resolution and a 6 h temporal interval (Dee et al., 2011). The chemical lateral boundary conditions are from a quasi-global simulation with 360×145 grid cells (67.5° S–77.5° N, 180° W–180° E) at a 1° × 1° horizontal resolution (Zhao et al., 2013; Hu et al., 2016). Aside from the aforementioned sets of experiments, four sets of sensitive experiments with aerosol radiative feedback disabled are also conducted to examine the aerosol radiative feedback effects on meteorological fields. We also conducted a detailed analysis of the computational costs associated with both the Resolved method and the Interpolated method. For the Resolved method, the calculation of aerosol optical properties takes 11 458.6 s, which accounts for 9.3 % of the total simulation runtime of 122 717 s. In contrast, the Interpolated method requires 5223.05 s for the same calculations, which accounts for 4.9 % of its total runtime of 107 615 s. Therefore, the Resolved method has an approximately 2.19 times higher computational cost for aerosol optical property calculations compared to the Interpolated method. This difference in calculation time translates to an additional 14 % in total simulation runtime when using the Resolved method. It is worth noting that this increase in computational cost is less than the 3.5 times one might expect from increasing the number of shortwave bands from 4 to 14. This is because the aerosol optical property calculation process includes both shortwave and longwave calculations. The original WRF-Chem has already used the Resolved method to calculate the longwave component; therefore, we only modified the shortwave calculations from 4 bands to 14 bands.

Table 1 The lower and upper bound and volume fraction for the eight aerosol size bins used in the model.

Anthropogenic emission for the domain covering China is from the Multi-resolution Emission Inventory for China (MEIC) at 0.1° × 0.1° horizontal resolution for 2015 (Li et al., 2017a, b), while the one for the domain covering Sahara is obtained from the Hemispheric Transport of Air Pollution version 2 (HTAPv2) at 0.1° × 0.1° horizontal resolution for 2010 (Janssens-Maenhout et al., 2015). The dust emission flux is calculated with the GOCART dust emission scheme (Ginoux et al., 2001), and the size distribution of emitted dust particles follows a theoretical expression based on the physics of scale-invariant fragmentation of brittle materials derived by Kok (2011). The details of the aerosol size distributions from approximately 0.04 to 10 µm are given in Table 1. More details about the dust emission scheme in WRF-Chem can be found in Zhao et al. (2010, 2013). Biomass burning emissions are obtained from the Fire Inventory from NCAR (FINN) with hourly temporal resolution and 1 km horizontal resolution (Wiedinmyer et al., 2011). Sea salt emissions follow Zhao et al. (2013) and include correction of particles with a radius less than 0.2 µm (Gong, 2003) and dependence on sea surface temperature (Jaeglé et al., 2011). The detailed parameterization schemes for physical and chemical processes in WRF-Chem that are used in the study are summarized in Table 2.

Table 2 Physical and chemical options of WRF-Chem used in this study.

2.4 Dataset

To evaluate the modeling results, multi-spectral AOD measurements are required, which are critical for this study comparing the Interpolated and Resolved methods across wavelengths. Therefore, we retrieved total AOD from the AERONET network (Holben et al., 1998). Although satellite AOD products such as the Moderate Resolution Imaging Spectroradiometer (MODIS) have greater spatial coverage, the number of wavelengths is limited. Comparison at only one or very few wavelengths would not effectively demonstrate the distinctions between the two algorithms examined in this work. In this study, the monthly mean AOD from the AERONET Version 3 direct-sun algorithm Level 2.0 dataset is used. A subset of stations dominated by dust is used here. The selected stations need to meet the following conditions: (1) to reduce the impact of oceanic aerosols, only the sites that are located on land are used; (2) to allow for comparisons with the simulation results, the sites must contain data for January and July in 2015; and (3) the monthly average AE (500–870 nm) at each site should be less than 0.8 because the lower the AE is, the larger the dust fraction will be (Dubovik et al., 2002).

Figure 2 Simulated AOD as a function of wavelength averaged for January and July 2015. (a) AOD averaged over the anthro-dominant region in eastern China (as shown in Fig. 1). (b) AOD averaged over the dust-dominant region in the Sahara (as shown in Fig. 1). The blue line and red line represent the Interpolated method and the Resolved method, respectively.

Figure 3 Comparison of total AOD from simulations and AERONET observations in dust-dominant areas. The results are averaged for both January and July 2015. The blue line and red line represent the Interpolated method and the Resolved method, respectively. The AERONET AOD values are indicated by black dots in each panel. The simulation results are obtained from the grid box closest to the AERONET stations.

3 Results

3.1 Impacts on aerosol optical properties

AOD is one of the key optical properties of aerosol. In Fig. 2 we show the simulated AOD values over the anthro-dominant region in China and over the dust-dominant region in the Sahara that are measured at different wavelengths using the two algorithms. The results presented are spatial averages over the regions delineated by the dashed boxes in Fig. 1, as well as temporal averages for the months of January and July 2015. The blue line represents the simulated AOD for 14 shortwave bands used by the RRTMG radiation scheme with the Interpolated algorithm based on the calculated values for 400 and 600 nm wavelengths. The red line represents the simulated AOD for 14 shortwave bands directly calculated with the Resolved algorithm. As seen in Fig. 2a, both the Resolved and Interpolated algorithms produce similar exponential decaying trends in AOD for regions dominated by anthropogenic aerosols. The AOD values calculated with the Resolved algorithm are slightly higher than those obtained with the Interpolated algorithm. These trends indicate that Ångström's theory is applicable to a certain extent in these areas. However, results in Fig. 2b illustrate a significant impact of algorithm modification on the simulation of AOD in regions dominated by dust aerosols. Both the Resolved and Interpolated algorithms calculated an upward trend at 400 and 600 nm wavelengths. However, since the Interpolated algorithm only includes information from these two wavelengths, it results in an exponential increase across all bands. Conversely, the Resolved algorithm reveals a fluctuating downward trend for dust aerosols across all wavelengths. This discrepancy leads to a substantial difference at longer wavelengths, with the maximum divergence exceeding 50 %. To assess whether the Resolved algorithm is more accurate, we compared our simulated AOD results with the AERONET observations in dust-dominant regions of the Sahara, as shown in Fig. 3. The selected AERONET stations follow the algorithm described in Sect. 2.4. The geographical locations of these stations are marked in Fig. 1b. AERONET observations demonstrate a gradual decline in dust AOD with increasing wavelength, with minimal variation between 400 and 600 nm. Although both the Resolved and Interpolated algorithms incorrectly predict an increasing trend in AOD between 400 and 600 nm, the Resolved method's more detailed spectral optical calculations offer the opportunity to correct this discrepancy at other wavelengths. Conversely, the Interpolated algorithm's erroneous trend at these wavelengths, when applied in the AE theory, propagates and amplifies discrepancies at longer wavelengths. These findings underscore the need for cautious application of the AE theory, particularly when dealing with aerosols such as dust, whose optical properties exhibit minimal variation at the selected wavelengths (400 and 600 nm in WRF-Chem). In such cases, small uncertainties can lead to a reversal in the AOD–wavelength relationship, potentially resulting in significant simulated errors in spectral extrapolation. Additionally, the simulated total AOD was also compared with the AERONET results in the anthro-dominant region in China (see Fig. S1 in the Supplement). The results demonstrate that both the Interpolated and Resolved algorithm simulations exhibit significant decreasing trends in total AOD compared to those seen in the AERONET data. Therefore, the discrepancies introduced by the application of the AE theory have a much smaller effect on simulated results, which further indicates that the amendments to the algorithm have less impact on simulated AOD resulting from anthropogenic aerosols compared to dust aerosols.

Figure 4 Simulated SSA as a function of wavelength averaged for January and July 2015. (a) SSA averaged over the anthro-dominant region in eastern China (as shown in Fig. 1). (b) SSA averaged over the dust-dominant region in the Sahara (as shown in Fig. 1). The blue line and red line represent the Interpolated method and the Resolved method, respectively.

Figure 5 Simulated AAOD as a function of wavelength averaged for January and July 2015. (a) AAOD averaged over the anthro-dominant region in eastern China (as shown in Fig. 1). (b) AAOD averaged over the dust-dominant region in the Sahara (as shown in Fig. 1). The blue line and red line represent the Interpolated method and the Resolved method, respectively.

Another crucial aerosol optical property is SSA, which is illustrated in Fig. 4. The amendment of the algorithm results in smaller anthropogenic and dust SSA over a considerable range of wavelengths, which could result in an overall larger absorption effect. At wavelengths smaller than the 600 nm waveband, the SSA simulated by the Resolved algorithm generally follows a linear function, as suggested by Interpolated algorithm in both regions. However, the value of SSA no longer increases with increasing wavelength when it reaches around 600 nm. Moreover, beyond 2000 nm in wavelength, the Resolved algorithm starts to decrease, whereas the Interpolated algorithm does not simulate this pattern. Figure 5 shows the aerosol absorption optical depth (AAOD) as a function of wavelength. The AAOD is calculated by subtracting the scattering radiation from the extinction radiation (AOD) by aerosols, which can be described as $\mathrm{AAOD}=\text{AOD}\cdot (\mathrm{1}-\text{SSA})$. AAOD can represent the absorption effects caused by aerosols in the atmosphere. Although the Resolved algorithm simulates smaller anthropogenic AOD in China and dust AOD in the Sahara, the Resolved algorithm simulates larger AAOD at all 14 bands due to its simulated smaller SSA. In WRF-Chem, the water content of aerosols is also considered while calculating optical properties. Water exhibits a peak in absorption at approximately 2900 nm, which significantly influences the overall aerosol optical properties at this wavelength. As illustrated in Fig. S2, there is a pronounced peak in the imaginary part of the refractive index of water at around 2900 nm in the Resolved method. The imaginary part of the refractive index is directly related to absorption, with higher values indicating stronger absorption. This peak in water absorption leads to the observed increase in AAOD (decrease in SSA) at this wavelength. In contrast, the Interpolated method, which relies solely on optical property information at 400 and 600 nm, fails to capture this crucial spectral feature. Consequently, the Interpolated method is unable to accurately represent the complex wavelength-dependent optical properties of aerosols. In addition, the lower SSA values and higher AAOD values obtained from the Resolved method compared to the Interpolated method in the Sahara region are primarily due to the spectral variation in dust optical properties, particularly at shorter wavelengths. In the OPAC dataset, dust has a larger imaginary part of the refractive index at wavelengths shorter than ∼600 nm, indicating stronger absorption in this spectral range (also as shown in Fig. S2). However, in the latest version of WRF-Chem (Interpolated method), the imaginary part of the refractive index for dust is set to a constant value of 0.003 at all four shortwave bands. The Resolved method simulates higher AAOD using the wavelength-dependent refractive indices, especially at shorter wavelengths. This results in significant differences in SSA and AAOD between the two methods, particularly in dust-dominated regions like the Sahara.

The differences in AAOD indicate that the algorithm modifications lead to a larger absorption effect from aerosols, which will be illustrated in more detail in the discussions below (see Sect. 3.2 and 3.3).

Figure 6 Radiative forcing of all aerosols in China at TOA, BOT, and in the ATM averaged for January and July 2015. Panels (a)–(c) and (d)–(f) show the results using Interpolated method and Resolved method, respectively. Panels (g)–(i) show the differences between the Interpolated method and the Resolved method. The average results of anthro-dominant areas (dashed boxes) are shown in the top-right corner.

Figure 7 Radiative forcing of aerosols in the Sahara at TOA, BOT, and in the ATM averaged for January and July 2015. Panels (a)–(c) and (d)–(f) show the results using Interpolated method and Resolved method, respectively. Panels (g)–(i) show the differences between the Interpolated method and the Resolved method. The average results of dust-dominant areas (dashed boxes) are shown in the top-right corner.

3.2 Impacts on radiative forcing of aerosols

As discussed above, there are significant differences in the aerosol optical properties computed based on these two algorithms, which may also lead to biases in simulating aerosol radiative forcings. Radiative forcing is defined as the perturbation of radiative fluxes at the top of the atmosphere (TOA) and the bottom of the atmosphere (BOT), as well as the perturbation of radiative heating and cooling in the atmosphere (ATM) if a specific aerosol species is removed. It should be noted that the aerosol radiative forcing calculated in this section refers to the change in radiative fluxes resulting from the removal of aerosols in a single experiment (see Sect. 2.2), excluding the perturbations to other meteorological variables caused by the removal of aerosol radiative effects as introduced by Zhao et al. (2013). In this study, the net downward radiative flux at TOA and BOT is considered positive, while upward flux is considered negative. In contrast, the heating effect of radiative flux within the ATM is considered positive, while the cooling effect is considered negative. Figure 6 illustrates the spatial distribution of aerosol radiative forcing computed using the Interpolated and Resolved algorithms and the differences between the two algorithms for TOA, BOT, and the ATM in China. Figure 7 shows the aerosol radiative forcing in the Sahara. In both China and the Sahara, the Resolved algorithm simulates more aerosol “warming” effects in the ATM, more negative radiative forcing at BOT, and smaller negative forcing at TOA compared to the Interpolated algorithm. These discrepancies can be explained by the stronger aerosol absorption effects computed from the Resolved algorithm, as discussed in Sect. 3.1.

As shown in Fig. 6, in regions predominantly influenced by anthropogenic aerosols, the algorithm amendment primarily affects the aerosol radiative forcing in the ATM. As discussed in Sect. 3.1, the Resolved algorithm is capable of simulating stronger aerosol absorption effects compared to the Interpolated algorithm, resulting in a stronger warming effect in the ATM (approximately 30 % enhancement). Additionally, the algorithm amendment introduces more pronounced “cooling” effects at BOT due to the aerosol radiative perturbations. The combined effects in the ATM and at BOT contribute to a relatively small impact of algorithm amendment on aerosol radiative forcing at TOA. In contrast, in the northwestern part of China, which includes the Gobi Desert and the Taklimakan Desert (two major dust source regions), aerosols are predominantly composed of dust. In these regions, the impacts of algorithm amendment on radiative forcing are more prominent compared to the impacts in anthro-dominant regions. In the Sahara region, the Resolved algorithm predicts a much greater dust warming effect in ATM compared to the Interpolated algorithm (approximately 140 % higher). Furthermore, the Interpolated algorithm only simulates negative dust radiative forcing at TOA in dust-dominant areas. In contrast, the Resolved algorithm can simulate positive forcing at TOA and can exceed 10 W m⁻² in the Sahara due to the radiative absorption of dust when located over highly reflective surfaces. The positive aerosol radiative forcing at TOA in the Sahara region simulated by the Resolved algorithm is notably more consistent with previous studies (e.g., Albani et al., 2014; Feng et al., 2022, 2023) compared to the Interpolated algorithm. Furthermore, our previous study (Feng et al., 2023) initially employed the original optical calculation algorithm in WRF-Chem (referred to as the Interpolated algorithm), which yielded a dust direct radiative forcing at TOA of −0.75 W m⁻². This value significantly deviated from the observationally constrained estimate of −0.20 W m⁻² proposed by Kok et al. (2017). Subsequently, in final version of Feng et al. (2023) implemented the Resolved method, which resulted in a substantially improved estimate of −0.27 W m⁻². This marked improvement in alignment with observational constraints strongly suggests that the Resolved method demonstrates better performance in simulating radiative forcing, particularly for dust aerosols. Apart from TOA and ATM, the algorithm amendment also results in a stronger negative radiative forcing (approximately 50 % higher) of dust aerosols at BOT, which may lead to cooling effects at the Earth's surface. In summary, the modified algorithm has an effect on anthropogenic aerosol radiative forcing, although this effect is relatively small. However, its impact on dust aerosols is significantly pronounced, to the extent of yielding divergent outcomes at TOA. Hence, to enhance the accuracy of aerosol radiative feedback simulation in the model, the amendment of this algorithm is imperative.

Figure 8 Shortwave heating profile of aerosols averaged for January and July 2015 (a) over the anthro-dominant area in China and (b) over dust-dominant area in the Sahara.

Figure 8 displays the vertical profiles of the shortwave aerosol heating rates calculated by the two algorithms over anthro-dominant areas in China and dust-dominant areas in the Sahara, respectively. In both regions, the shortwave heating effect induced by aerosols is strongest at the surface and decreases rapidly as the altitude approaches approximately 500 hPa. Beyond this altitude, the aerosol-induced heating effect stabilizes at a relatively constant value, maintaining a positive heating effect throughout the atmospheric column. In the areas predominantly affected by anthropogenic aerosols, the Resolved algorithm simulates a more pronounced aerosol heating effect, displaying a vertical trend similar to that produced by the Interpolated algorithm. In the dust-dominant area, however, the impact of algorithmic amendment on the simulation of aerosol heating effects is significantly greater than that for anthropogenic aerosols, with the heating rate near the surface exceeding twice the result simulated by the Interpolated algorithm. This discrepancy diminishes near the altitude of 500 hPa but remains approximately 80 % stronger than that predicted by the Interpolated algorithm. This also explains the impact of algorithm amendment on the simulation of radiative forcing in ATM, as discussed previously. The results of the heating profile also illustrate that the effects of algorithm amendment vary significantly among different types of aerosols.

It should be noted that the refractive indices of dust can vary significantly depending on mineral composition and source region. Our study uses a set of refractive indices that represent a global average, which may not capture the full range of variability observed in different regions. Changing these values could affect the magnitude of the radiative forcing differences between the Interpolated and Resolved simulations. We have conducted additional sets of experiments with smaller radiative indices for dust aerosols in Resolved simulations. The results show that with less absorbing dust (smaller imaginary part of refractive indices), the differences in the radiative effects could become smaller (not shown here). Theoretically, one could identify a set of aerosol optical properties that would yield similar aerosol radiative forcing results for the Resolved and Interpolated algorithms. However, such a set might not accurately reflect the true physical properties of the aerosols. Despite these considerations, it is important to note that the spectral variation in optical properties would be better captured by the Resolved method compared to the Interpolated method, as discussed in Sect. 3.2. This provides a more robust foundation for accurately representing the complex interactions between aerosols and radiation across various wavelengths. The advantages of the Resolved method in capturing these spectral variations persist regardless of the specific refractive indices employed in the simulations.

3.3 Impacts on radiative effects of aerosols

The impacts of the algorithm amendment on radiative forcing of aerosols can further influence the radiative feedback of aerosols on meteorological fields, such as temperature, wind field, and planetary boundary layer (PBL) height, as discussed below. An additional set of experiments are conducted with the radiative feedback of aerosols disabled. The differences in the simulation results between the two sets of experiments (one with aerosol radiative feedbacks enabled and the other with them disabled) are used assess the radiative impacts of aerosols on meteorological fields.

Figure 9 (a, d) Spatial distribution of aerosol radiative effects on 2 m temperature from the Interpolated experiments in China and the Sahara, respectively. (b, e) Spatial distribution of aerosol radiative effects on 2 m temperature from the Resolved experiments in China and the Sahara, respectively. (c, f) The impacts of algorithm amendment on the simulated aerosol radiative effects (difference in aerosol radiative effects between Resolved and Interpolated experiments) on 2 m temperature in China and the Sahara, respectively. The average results of anthro-dominant (dashed boxes in a, b, and c) and dust-dominant areas (dashed boxes in d, e, and f) are shown in the top-right corner. The results are averaged for January and July 2015.

Figure 10 The same as Fig. 9 but for the net shortwave radiative fluxes at the bottom of the atmosphere (positive denotes downward).

Figure 9 shows the aerosol radiative effects on 2 m temperature from the Interpolated and the Resolved experiments, as well as the differences between the two algorithms. From Fig. 9a and b, it can be observed that in regions dominated by anthropogenic aerosols in China, aerosol radiative effects lead to surface cooling in both algorithms. In the Tibetan Plateau region, it results in surface warming. According to the results in Fig. 6c and f, aerosols exhibit a cooling effect on the surface due to their direct radiative effect in all regions. This difference indicates that the radiative effects of aerosols on near-surface temperature are influenced not only by radiative forcing (direct radiative effect) but also by other factors such as aerosol–cloud interactions, aerosol heating effects, and radiative effects on the surface radiative fluxes, as discussed below. Figure 10 illustrates the impacts of aerosols on the shortwave radiative fluxes at BOT. Compared with Fig. 6, these effects include not only the direct radiative effect (DRE) of aerosols but also their indirect consequences. The alterations in radiative forcing induced by aerosols can lead to effects on the atmospheric energy balance, which in turn influence other meteorological processes such as cloud formation and thus further affect radiation. Figure 10a and b exhibit similar spatial distribution characteristics compared to Fig. 9a and b, indicating that the radiative effects of aerosols on 2 m temperature are primarily influenced by their effects on surface radiative fluxes. In addition, 2 m temperature is also affected by the aerosol heating near the surface. As shown in Fig. 8a, the new Resolved algorithm simulates a stronger heating near the surface compared to the Interpolated algorithm. This proves that the average effect of aerosol on 2 m temperature is slightly larger with the Resolved algorithm than with the Interpolated algorithm, despite the average surface radiative cooling being stronger in the Resolved algorithm.

In the Sahara region, the effects of aerosol on radiative fluxes at the surface are similar to the aerosol radiative forcing (Fig. 7c, f), showing a cooling effect at the surface throughout the domain (Fig. 9d, e). It is noteworthy that although the effects of aerosol on radiative fluxes at the surface are negative in the Sahara region, the effects on 2 m temperature with both algorithms could exhibit a warming effect. This is due to the aerosol heating effect in the atmosphere near the surface, which also leads to stronger warming effect of aerosols at the surface from the Resolved algorithm in some areas. In both regions, the aerosol effects simulated with the new Resolved algorithm reduce the simulation biases in 2 m temperature with the Interpolated algorithm to some extent compared to the ERA5 data, particularly over the Sahara region (Fig. S3).

Figure 11 The same as Fig. 9 but for wind fields (vectors) and geopotential height (shaded contour) at 850 hPa.

The effects of aerosol on wind and geopotential height at 850 hPa simulated by both algorithms are illustrated in Fig. 11. The results indicate that the aerosol radiative effects simulated with the Interpolated algorithm lead to a significant decrease in geopotential height in eastern China, whereas the effects are relatively small with the Resolved algorithm. Interestingly, region the effects of aerosol on geopotential heights and wind fields at 850 hPa simulated by the Interpolated algorithm are small over the Sahara, while the effects simulated by the Resolved algorithm are much larger. This may be due to the significantly larger reduction and warming effects from the Resolved algorithm over the Sahara region. Distinct mechanisms of the effects of aerosol on wind fields and geopotential heights from different types of aerosols over different climate regimes deserve further investigation in the future. In both regions, the algorithm amendment results in significant differences in the effects of aerosol on wind fields, which partially reduces the biases in simulated cyclonic wind circulation compared to the ERA5 reanalysis along the southeastern coastal region of China (Fig. S4).

Figure 12 The same as Fig. 9 but for PBL height.

Previous studies have also highlighted the important role of the effects of aerosol on the development of the PBL and hence on the air quality near the surface (e.g., Liu et al., 2016; Wilcox et al., 2016; Yang et al., 2017). Aerosols could reduce the near-surface temperature, heat the upper atmosphere (as discussed above), and therefore suppress PBL development (Huang et al., 2018). Figure 12 illustrates the aerosol radiative effects on PBL heights and their differences between the simulations with the Resolved and Interpolated algorithms. In the regions dominated by anthropogenic aerosols in China, the heating rates in the upper atmosphere or around the top of the PBL due to absorbing aerosols are relatively small (Fig. 8a). Therefore, the differences in the effects of aerosol on PBL height between the two algorithms are primarily associated with their differences in 2 m temperature (Fig. 9). However, in the Sahara region, the difference in the effects of aerosol on 2 m temperature and the heating rates at the upper atmosphere or around the top of the PBL are significant between the two algorithms (Figs. 8 and 9). Overall, the effects of aerosol suppress PBL development in most areas of Sahara, leading to a decrease in PBL height. Consequently, the algorithm amendment significantly affects the effects of aerosol on PBL development in both regions, with an average reduction of ∼20 m and up to ∼100 m in China, with an average reduction of ∼40 m and up to ∼200 m in the Sahara.

Figure 13 The same as Fig. 9 but for particulate matter with diameters of 10 µm and smaller (PM₁₀) at surface.

As previously mentioned, aerosol radiative effects on the height of the PBL could concurrently affect air pollutant concentrations within the PBL (Ding et al., 2016). Figure 13 illustrates the aerosol radiative effects on PM₁₀ (particulate matter with diameters 10 µm and smaller) at the surface and the differences between the simulations with the two algorithms. In the China region, Fig. 13 shows similar spatial patterns to those in Fig. 12, confirming that lower PBL can raise the surface PM concentration. In the Sahara region, while aerosol radiative effects generally lead to an increase in surface PM₁₀ concentrations, there are still areas with a reduction, particularly for the simulations with the Interpolated algorithm. The Resolved algorithm results in significant differences in the effects on surface PM₁₀ in the Sahara. Please note that the impacts over the Sahara are more complex because aerosol radiative effects could affect both PBL heights and emissions (through near-surface wind) and hence the near-surface mass concentrations.

4 Summary and discussion

Aerosol–radiation interaction can have important impacts on meteorological processes and the aerosol cycle. WRF-Chem, a fully coupled online meteorology–chemistry model, has been widely used to investigate the impacts of aerosol–radiation interaction at regional scales. In this study, the original Interpolated algorithm for calculating aerosol optical properties and radiative effects in WRF-Chem is re-examined against the Resolved algorithm implemented in this study. Two domains are selected for investigating the difference between the two algorithms, with one covering China, which represents a region with complex aerosol sources including large anthropogenic aerosol mass loading and also natural dust over the northwest, and the other covering the Sahara, which represents the region with the largest natural dust aerosol mass loading in the world.

The discrepancies between the two algorithms show distinct regional characteristics. In China, where anthropogenic sources dominate the aerosol composition, the differences between the Resolved and Interpolated algorithms are relatively small, which could potentially be indicative of similar patterns in other regions dominated by anthropogenic aerosols. In contrast, the Sahara, which is dominated by dust aerosols, exhibits significant differences between the two algorithms. The Resolved algorithm results of dust AOD show a general downward trend with increasing wavelength rather than an upward trend calculated by the AE. The maximum difference between the two algorithms can reach about 50 %. Further comparison with AERONET observations reveals that the Resolved algorithm's AOD simulations are in better agreement with the measured values at dust-dominant stations. This suggests that the Resolved approach can more accurately capture the optical properties of dust aerosols. The Resolved algorithm also simulates smaller SSA than the Interpolated algorithm. As it is affected by these two factors (AOD and SSA), the Resolved algorithm simulates larger AAOD than the Interpolated algorithm, resulting in larger aerosol heating effects.

The impacts of algorithm amendment on aerosol radiative forcing are different depending on the aerosol type and region. For areas in China with high concentrations of anthropogenic aerosols, the Resolved algorithm enhances the aerosol radiative absorption in the atmosphere by about 30 % compared to the Interpolated algorithm. It also introduces a larger cooling effect at the surface. The impact on the radiative forcing at the top of the atmosphere is small. In the areas dominated by dust aerosol, the impacts of algorithm amendment are substantially larger. The Resolved algorithm predicts a ∼140 % higher warming effect in the atmosphere from dust aerosol compared to the Interpolated algorithm. Moreover, the Resolved algorithm can simulate positive radiative forcing at the top of atmosphere (exceed 10 W m⁻²) in dust-dominant areas, aligning better with previous studies constrained by observations, against the negative values with the Interpolated algorithm. The algorithm amendment also causes a roughly 50 % larger negative radiative forcing at the surface from dust, leading to stronger surface cooling.

In addition to the impacts on aerosol optical properties and radiative forcing, the impacts of aerosol radiative effects on meteorological fields are also investigated. Both algorithms simulate that aerosols reduce (increase) near-surface temperature in the anthro-dominant areas (Tibetan Plateau) of China. The Resolved algorithm leads to a slightly larger increase in near-surface temperature in China than the Interpolated algorithm. Over the Sahara region, both algorithms simulate a dominant cooling effect on near-surface temperature over most of the region but with a warming effect over some areas. The Resolved algorithm leads to stronger effects in both cooling or warming areas. The difference in the effects of aerosol on near-surface temperature between the two algorithms can be explained from their difference in simulating net radiative fluxes at the surface that can result from both direct (radiation) and indirect (cloud) aerosol effects. The experiment with the Resolved algorithm better simulates the 2 m temperature as compared to the ERA5 reanalysis data than the one using the Interpolated algorithm. In addition, the algorithm amendment also leads to different the effects of aerosol on the wind fields and geopotential height over both regions, with a larger impact present over the Sahara compared to over China. Both algorithms simulate the aerosol radiative effects to suppress PBL development and thus reduce PBL height. The algorithm amendment leads to a further reduction in PBL height of ∼20 m in China and ∼40 m in the Sahara on domain average. The enhancement of aerosol radiative effects in terms of reducing the PBL height in the Resolved algorithm leads to more accumulation of surface concentrations of PM₁₀.

Please note that the impacts of aerosol–radiation interaction on meteorological fields and chemical species are not only through direct effects on radiation but also through indirect effects on cloud and then precipitation and winds. For example, some differences between the two algorithms in simulating surface PM₁₀ concentration may not be fully explained by their difference in radiative fluxes but instead are contributed to by other factors such as their induced changes in surface wind-driven emissions and precipitation-driven wet removals (Feng et al., 2023). The mechanisms driving the differences between the impacts of the two algorithms deserves further investigation in the future. This study underscores the importance of refining the algorithm of aerosol–radiation interaction for simulating the effects of aerosol on weather and climate more accurately. It cautions against the usage of the original Interpolated algorithm in WRF-Chem for simulating aerosol optical properties and their impacts on meteorological fields, which has some biases, particularly for regions with a large contribution from dust. It is necessary to update the model to use the new Resolved algorithm proposed in this study in future.

It is important to note that the primary aim of our study is to demonstrate the importance of using a spectrally resolved method for calculating aerosol optical properties rather than to improve the overall model performance. Many factors can affect the simulation of meteorological fields and radiative processes beyond the optical properties methods we are investigating in this study. For example, while our study employs a 50 km grid resolution, which is suitable for investigating aerosol–radiation interactions, higher resolution could enhance the simulation of aerosol emission, deposition, and transport processes, potentially leading to a more accurate representation of aerosol distributions and their radiative effects (Feng et al., 2023; Tan et al., 2015; Tao et al., 2020). In addition to model resolution, the calibration factors for aerosol emission rates, the representation of aerosol size distributions, the quality and accuracy of input data, and the selection and implementation of parameterization schemes for various physical processes could all introduce uncertainties and potential impacts on the simulated results. Given this complexity, direct comparisons of the Resolved and Interpolated methods with observations may not provide a conclusive assessment of whether our modifications improve the model's overall simulation abilities. Therefore, we did not evaluate the model's simulation results of meteorological fields from two methods by comparing them with more observation results other than the ERA5 reanalysis dataset. Our results show that the Resolved method can capture complex relationships between aerosol optical properties and wavelengths that the Interpolated method may miss, particularly for dust-dominated regions and at specific wavelengths where water contents have a significantly larger absorption than other wavelengths. In addition, our results demonstrate that the amendment of algorithms can significantly affect the simulation results of meteorological fields. While these changes may not necessarily lead to better agreement with observations in all cases, they give the model more potential to improve simulation abilities by more accurately representing the underlying physical processes.