Introduction

Dust is the main component of atmospheric aerosols, accounting for 75% of the global aerosol mass loading and 25% of the global aerosol optical depth (AOD)1. It is mainly caused by aeolian erosion over drylands and affects the energy budget of the Earth–atmosphere system and hydrological cycle through direct and indirect effects2,3,4,5. Moreover, dust particles can be enriched with germs6, organics7, and heavy metals8, posing a serious threat to human health and socioeconomic activities by affecting air quality9,10,11. After a long-distance transport, dust particles are deposited on the ocean surface, affecting the generation of biogenic aerosols over the ocean and changing the marine biogeochemical cycle and biological productivity12,13.

Aerosol–radiation interaction considerably impacts climate changes, weather processes, and air quality14. On the one hand, changes in meteorological processes can alter the transport and spatial distribution of dust15,16. On the other hand, dust affects the atmospheric system and air quality through radiative feedback17. Dust direct radiative forcing values for the Taklimakan Desert (TD) and Gobi Desert (GD) are reported as −3 and −7 W m−2 at the top of the atmosphere (TOA), −8 and −10 W m−2 at the surface, and +5 and +3 W m−2 in the atmosphere18. In addition, dust radiative forcing can considerably change atmospheric stratification and baroclinicity19. GD dust enhances the anthropogenic aerosol pollution in eastern China by changing the meteorological field20. The dust-induced pressure disturbance prompts secondary circulation, which reduces the wind speed at the low atmospheric level over the dust source region and enhances the wind speed at the downstream regions, resulting in dust reduction and enhancement in the upper source region and downstream regions, respectively21.

Located at the border between China and Mongolia, the GD is an important source of dust in East Asia22,23,24. As a typical plateau desert, the GD is characterized by a low average annual temperature and drastic temperature changes25. Dust storms occur abruptly in this region over wide areas but are short in duration. The dust-event frequency at the China–Mongolia border has shown a gradually increasing trend26,27, and GD dust is the main source of dust-type air pollution in the inland areas of China, particularly, Beijing–Tianjin–Hebei28. GD is characterized by flat terrain, and the entire troposphere is dominated by westerly winds. The combination of special terrain and background winds provides conditions for trans-Pacific eastward transmission of dust to North America29,30.

Although researchers have consistently discussed the importance of GD’s contributions to dust sources in East Asia28, the effect of GD dust on the weather and climate in this region has been poorly studied. As an important weather system, the Mongolian cyclone can generate strong surface winds that lift surface particles into the atmosphere and then produce severe dust storms31. The outbreak of dust storms in East Asia is often closely related to the activity of cold air and the Mongolian cyclone23,32. Cyclones with cold air are common in Mongolia, and the GD and TD generate large amounts of dust that is transported to the downstream regions33. Moreover, the westerly ridge and East Asian vortex at 500 hPa are key drivers of long-distance dust transport34,35. A severe dust event occurred in May 2019 directly affected 168 people in a county in Inner Mongolia, resulting in direct economic losses of 1.004 million RMB36. Herein, the May 2019 event is used as a typical case to further explore the relation between Mongolian cyclone and the dust radiative feedback in the GD. We use dynamic dust sources, the Weather Research and Forecasting (WRF) model coupled with Chemistry (WRF-Chem), and satellite retrievals to develop a comprehensive understanding of this event. Further, we explore the mechanism of dust radiative feedback in the GD on enhancing Mongolian cyclone and its continued effect on dust concentrations over downstream regions.

Results

Spatial–temporal dust distribution

MODIS composite image shows that this dust event was closely related to Mongolian cyclone (Fig. 1a). The obvious comma cloud system over northern China in the figure indicates that the Mongolian cyclone enhanced the easterly transport of mid-latitude dust aerosols via atmospheric circulation. At the same time, the simulated weather circulation also captures a deep vortex structure over northern China on May 13, 2019 (Fig. 1b). WRF–Chem effectively captured the maximum MODIS aerosol distribution in central and eastern Inner Mongolia and in northern Heilongjiang, with AOD values up to 1.8 (Fig. 1c, d). There are well-documented evidence for the connection of aerosol index (AI) and aerosol concentrations and optical properties37. It is a qualitative factor represent aerosol particles suspended in the atmosphere, mainly from desert dust, biomass burning and volcanic ash plumes, and it has already been used to identify dust aerosols38. Generally speaking, both OMI AI and WRF–Chem AOD show consistent zonal distribution in Inner Mongolia (Supplementary Fig. 3). The air in northern China was also severely polluted, with the AI peaking at 2.9. The AOD values shifted eastward owing to the lack of MODIS measurement in the main areas of the dust source in the event.

Fig. 1: Model evaluation of AOD spatiotemporal distribution.
figure 1

a MODIS Aqua composite image overlaid by OMI Aerosol Index (AI) on May 13 2019. The red dots indicate the AERONET observation stations including Irkutsk, AOE_Baotou, Beijing_PKU, and Anmyon. b Spatial patterns of the geopotential height (blue contour lines; unit: dagpm), temperature (red contour lines; unit: °C), and wind field (vectors; unit: m s−1) at 500 hPa on May 13, 2019, retrieved using WRF-Chem. Spatial distributions of the average AOD based on (c) MODIS (MYD08_D3) and (d) WRF-Chem simulations are shown for May 11–16, 2019. The wind field at 10 m (vectors; unit: m s−1) was obtained from (c) Final Operational Global Analysis (FNL) reanalysis and (d) WRF-Chem simulation. Daily variation of AOD from AERONET observations (Aeronet) at four sites (Anmyon (36.539°N, 126.330°E), Irkutsk (51.800°N, 103.087°E), AOE_Baotou (40.852°N, 109.629°E), Beijing_PKU (39.992°N, 116.310°E)) (e, f, g, h), MODIS on board Terra (MOD) and Aqua (MYD), and the WRF-Chem model (WRF-Chem) on May 11–30, 2019. The results associated with WRF-Chem are based on the EXP_CTRL experiment.

Full size image

To verify the numerical simulation of dust aerosols based on the WRF–Chem model, the simulated AOD on May 11–16 was evaluated using ground observations. The AOD recorded at Anmyon, Irkutsk, Academy of Optoelectronics (AOE)_Baotou, and Beijing_Peking University (PKU) sites derived from Terra, Aqua, and AERONET data were compared with the WRF–Chem results (Fig. 1). Simulations well reproduced the observed AOD in China. Specifically, Beijing suffered from serious air pollution during this dust event, with AOD values generally >0.5, whereas the WRF-Chem AOD was close to that of MODIS and AERONET (Fig. 1h). Although the simulated AOD in Baotou was overestimated compared with the observations, basic aerosol changes were noted (Fig. 1g). Moreover, the AOD simulations at the Anmyon station and the Irkutsk station were in good agreement with the observations (Fig. 1e, f).

Dust emission in the dust event was concentrated mainly over the GD (Supplementary Fig. 5). The maximum multiday average dust emission flux was as high as 182.1 µg m−2 s−1, and a small amount of dust particles was emitted from the TD. In the case of daily dust emission, peaks at 554.2 µg m−2 s−1 and 264.7 µg m−2 s−1 were simulated at the same grid (44.58408°N, 99.15584°E) on May 11 and 14, 2019, respectively, and the dust emission area extended eastward to 110°E. The occurrence of dust storms is usually closely associated to the invasion of cold air29. Cold air from Siberia is continuously transported to the GD, and the GD dust then affects the area over northeastern Mongolia via the southerly wind (Supplementary Fig. 6). In addition, Supplementary Fig. 6b shows that a weak cyclone was generated over northern Mongolia and the eastern region of Inner Mongolia on May 12, 2019. Mongolian cyclone gradually intensified and peaked on May 13, 2019, and a common cloud system is also found based on the satellite cloud images (Supplementary Fig. 7). The sea level pressure in the center of Mongolian cyclone on May 13 was lower than that on May 12. Moreover, the cyclonic wind field in northern Mongolia gradually dissipated on May 14 (Supplementary Fig. 6d). Therefore, it gradually intensified and peaked on May 13, 2019, and a comma cloud system is also found based on the composite image (Fig. 1a). Mongolian cyclone combined with southward cold air strengthened dust intensity in the dust event. The persistent maintenance of cyclones in central Siberia is conducive to continuous eastward transport of dust aerosols (Supplementary Fig. 6). Therefore, the combination of cold air and Mongolian cyclone promoted larger dust emission and dust loading on May 14 compared that on May 13. Moreover, the temperature field always lags behind the height field at different altitude levels, which is conducive to the development of trough–ridge systems. The large angle between the temperature and height fields in the middle and low atmospheric layers is beneficial for the southward movement of cold air from northern areas, and the transport of cold advection and the powerful cyclone in north of Mongolia enhances the development of dust storms (Supplementary Fig. 8). These phenomena enable the GD dust to continuously affect the downstream regions.

The dust transport flux at 850 hPa was concentrated over Inner Mongolia. With increasing altitude, the dust transport flux moved eastward to North China and Northeast China under westerly jets at 500 hPa (Fig. 2a, b). The spatial distribution of dust transport flux clearly shows that GD dust was transported to northeast of Mongolia, North China, and other downstream areas. Atmospheric uplift between 105°E and 110°E is conducive to the lifting of the dust, and then dust concentrations gradually spreads eastly through westerly winds (Fig. 2c). The GD and TD were the main dust emission areas in the dust event (Supplementary Fig. 5). To explore the relative contributions of dust and their uplifting capacity from these two sources28, we compared the dust concentrations and dust ratio in the middle and lower troposphere with those in the upper atmosphere (Fig. 2d, e). Specifically, the maximum GD dust concentrations at the middle and low troposphere (3–10 km) was 31.4 µg m−3, whereas that in high atmosphere (8–10 km) was 2.3 µg m−3 (Supplementary Fig. 11), with average values of 14.2 µg m−3 and of 0.5 µg m−3, respectively. Regarding TD dust, the maximum concentrations was 18.5 µg m−3 at the middle and low troposphere while that in the high layer was ~1.5 µg m−3, with average values of 6.4 and 0.2 µg m−3, respectively. The concentrations of GD dust in different layers of the atmosphere and the ratio of the dust concentrations in the upper layer to that in the lower layer (below 3 km) were both greater than those for TD dust. Specifically, the average dust ratios for the GD and TD were 3.88% and 2.97% at 3–10 km and 0.17% and 0.084% at 8–10 km, respectively. The larger dust concentrations and dust ratio for the GD indicate that GD has higher lofting of dust, and the contribution of GD to atmospheric dust concentrations was also greater than that of the TD.

Fig. 2: The basic physical processes in this dust event.
figure 2

Spatial distribution of dust transport flux (color contour, unit: µg m−2 s−1) at (a) 850 and (b) 500 hPa. c Average cross sections of dust concentrations (color contour, unit: µg m−3) and vertical speed (contour lines, unit: m s−1) during May 11–16, 2019, over domain 2 (solid and dashed lines represent ascent and descent motions, respectively). d Hourly average dust concentrations (unit: µg m−3) over the TD (35–45°N, 85–95°E; red solid line) and GD (35–45°N, 95–110°E; blue solid line) at 3–10 km altitude. e Ratio of dust concentrations (%) at 3–10 km against that measured below 3 km during May 11–16, 2019. The results are based on the EXP_CTRL experiment.

Full size image

Influence of dust radiative feedback on Mongolian cyclone

Dust affects the regional atmospheric thermal structure through direct radiative forcing, which in turn affects the regional climate39,40 and dust cycles41. The net radiative forcing values of dust at the TOA in North and Northeast China were mainly negative (~−20 W m−2), indicating that dust has a cooling effect in the Earth–atmosphere system in this region. As a typical absorbing aerosol, dust can warm the atmosphere40, and the net radiative forcing of dust in the atmosphere is mainly positive. The dust layer mainly appeared near 850 hPa, and the regional average dust concentrations reached the maximum at 875 hPa (729.12 μg m−3) (Fig. 3d). The influence of dust radiative feedback on the atmospheric heating rate varied with dust concentrations, and the atmospheric heating rate also reached 0.33 K day−1 at 700 hPa.

Fig. 3: Dust radiative forcing and dust-induced heating rate.
figure 3

Spatial distribution of dust net radiative forcing under all-sky (unit: W m−2) (a) at the top of atmosphere, (b) in the atmosphere, and (c) at the bottom of atmosphere (color contour; unit: W m−2). d Average vertical profiles for dust concentrations (blue solid line, unit: µg m−3) and dust-induced changes in the heating rate (red solid line, unit: K day−1). e Hourly heating rate (units: K day−1) induced by dust radiative feedback in the study area (domain 2) during May 11–16, 2019. The results are based on the difference between the two parallel experiments.

Full size image

Owing to the strong scattering and absorption of dust in the atmosphere as well as the influence of clouds, the amount of solar shortwave radiation reaching the surface was reduced. The net radiative forcing and shortwave radiative forcing of dust near the surface generally showed negative values (Fig. 3a, c). Atmospheric heating caused by dust exhibited considerable diurnal variation (Fig. 3e) and peaked at 1.22 K day−1 at 14:00 LT. Generally, dust aerosols change the dynamic structure by heating the atmosphere, and the enhanced wind speed enables more dust particles to be uplifted and transported, affording a positive feedback for dust transport.

Dust radiative feedback on the meteorological field (Fig. 4) facilitated the intensification of the Mongolian cyclone during the dust event and caused continuous accumulation of cold air in the middle and lower troposphere over northeast of Mongolia (Fig. 4a, b, c). The increase of low pressure always corresponds to the cold temperature anomaly, while the decrease of low pressure often corresponds to the warm temperature anomaly. Therefore, with increasing altitude, the low pressure gradually increased from 850 to 300 hPa, reaching maximum intensity near 300 hPa but gradually weakening above 300 hPa (Fig. 4g). The cold air accumulated at 700 hPa and 500 hPa (Fig. 4a, b) in this region also turned to warm ones at 200 hPa (Fig. 4c). Dust radiative feedback in the middle and low troposphere strengthens the zonal pressure and temperature gradient by enhancing the ridge of high pressure over northeast of China, intensifying the cold advection and Mongolian cyclone. In other words, dust radiative effect enhanced the temperature difference in the northeast of Mongolia, enhancing the geopotential height difference and cold advection, making the cyclone develop rapidly. Moreover, the wind field is also consistent with the variation in geopotential height field, which always presents cyclonic rotation in this region (Fig. 4d, e, f). The maximum potential vorticity induced by dust at 200 hPa exceeded 0.08 PVU. The potential vorticity decreases at 500 hPa and increases again at 700 hPa.

Fig. 4: Dust radiative feedback on enhancing Mongolian Cyclone.
figure 4

Spatial distribution of dust-induced temperature change (color contour; unit: °C) and geopotential height change (contour lines; unit: gpm) at (a) 700, (b) 500, and (c) 200 hPa and that of dust-induced potential vorticity (color contour; unit: PVU) and wind field (vectors; unit: m s−1) at (d) 700, (e) 500, and (f) 200 hPa. Vertical cross sections of dust-induced latitudinal anomalies of (g) temperature (color contour; unit: °C) and geopotential height (contour lines; unit: gpm) as well as (h) potential vorticity (color contour; unit: PVU) and wind field (vectors; unit: m s−1). The black dots and shadows represent statistical significance above the 90% level. The results are based on the difference between the two parallel experiments.

Full size image

Mongolian cyclone can be continuously maintained and strengthened via dust radiative feedback. Dust vertical distribution plays a key role in weather and climate effects42. To clarify the influence of dust radiative feedback on the temperature and pressure vertical structure (Fig. 4g), we analyzed the latitudinal anomalies section in the main cold air accumulation area (domain 1 in Supplementary Fig. 1). As shown in Fig. 4g, the colder air forms a wedge shape and descends deeply into the bottom of the warmer air, causing the warmer air to ascend along the colder ones. In the middle and back parts of the Mongolian cyclone, stable atmospheric structures were present with abundant colder (warmer) air accumulation in the lower (higher) layers. However, the front part of cyclone presented an opposite structure, with warmer (colder) air at lower (higher) altitudes. The development of the Mongolian cyclone was strengthened by the unstable and stable atmospheres present at the front and rear of the cyclone, respectively. Other’s research43 also support our results. It is reported that downward transfer of momentum from high altitudes is beneficial for intensifying Mongolian cyclone44. In the present study, the potential vorticity was maximum near 200 hPa and decreased at lower pressures (Fig. 4d, e, f). The cross section of the potential vortex with regional latitudinal anomalies shows that the maximum potential vortex (0.2 PVU) (where 1PVU = 1.0 × 10−6 m2 s−1 K kg−1) near 200 hPa experiences downward transfer of momentum under the action of the wind field (Fig. 4h). Therefore, the cyclone over northeast of Mongolia was strengthened.

We then used the thermodynamic equation to explore the reason for the continuous cold air accumulation in the middle and lower troposphere over northeast of Mongolia. The influence of dust radiative feedback on the cold air over this region originated mainly from the zonal wind (Fig. 5a) and zonal temperature gradient (Fig. 5e). When the zonal wind field was enhanced, dust radiative feedback increased the transport of atmospheric cold advection over Siberia, thus cooling the atmosphere over northeast of Mongolia. However, the disturbance of vertical atmospheric movement caused by the dust radiative feedback inhibited atmospheric cooling (Fig. 5h), which generally balances the thermodynamic equation. Dust radiative feedback yielded cyclone intensification over northeast of Mongolia, where atmospheric cooling and geopotential height reduction strengthened all layers of the cyclone. The strengthening of Mongolian cyclone further enhanced the upper westerly jets, forming a positive feedback process. This caused the GD dust aerosols to continue eastward transport, thus affecting northern China.

Fig. 5: Dust radiative feedback on the nine terms in the thermodynamic equation.
figure 5

Spatial distributions of nine terms (units: 10−2) of the thermodynamic equation (Eq. 10) in the enhanced Mongolian cyclone area: (a) \(-{U\text{'}}(\frac{\partial \bar{T}}{\partial x})\), (b) \(-\bar{U}(\frac{\partial {T\text{'}}}{\partial x})\), (c) \(-{U\text{'}}(\frac{\partial {T\text{'}}}{\partial x})\), (d) \(-{V\text{'}}(\frac{\partial \bar{T}}{\partial y})\), (e) \(-\bar{V}(\frac{\partial {T\text{'}}}{\partial y})\), (f) \(-{V\text{'}}(\frac{\partial {T\text{'}}}{\partial y})\), (g) \(\sigma ^{\prime} \bar{\omega }\), (h) σ¯ω′ and (i) \(\frac{\Delta Q}{{Cp}}\).

Full size image

Discussion

Based on two parallel experiments, the synoptic-scale effects of dust radiative forcing could be further explored. The effects of dust radiative feedback on the development and maintenance of Mongolian cyclone during a major dust event that occurred in May 2019 were investigated based on the dynamic dust sources coupled with the WRF-Chem model. Results showed that dust originating from the GD played a dominant role in the dust contribution over northern China. The intense Mongolian cyclone is responsible for the dust emission over the GD. With the influence of Mongolian cyclone, the GD dust continues to affect the downstream regions through westerly jets and is transported to northeast of Mongolia via southerlies. Dust radiative feedback shows a strong positive feedback on Mongolian Cyclone. By affecting zonal wind and temperature advection, the atmospheric cold anomaly in this region is enhanced. The intensified Mongolian cyclone is also sustained by enhancing the downward transmission of the high-level momentum of ~0.2 PVU and changing the stable atmospheric structure. After the Mongolian cyclone strengthened, westerly jets on its south side prevail and the GD dust is continuously transported to the downstream areas, causing higher dust concentrations in northern China (Fig. 6).

Fig. 6: Dust radiative feedback influence mechanism.
figure 6

Influence mechanism of dust radiative feedback over the GD on Mongolian cyclone intensification.

Full size image

Driven by the increase dust in Asia and North Africa, the global dust mass loading has increased 55 ± 30% since pre-industrial times45. Dust radiative feedback has recently attracted considerable scientific attention as it could change the atmospheric thermal structure, which in turn can affect the weather and climate change41,46. Dust aerosols influence the radiation budget of the Earth by absorption, scatting of the longwave and shortwave radiation47. The dust radiation effect greatly reduces the temperature deviation between the surface and the atmosphere48. The atmospheric heating rate increased with dust radiative effect, resulting an increase of about 0.2 K in the middle troposphere (about 600 hPa)49. Dust radiative effect can further affect circulation and weather and climate change by influencing the atmospheric thermal structure. Previous studies on dust radiative feedback mainly focused on changes in meteorological elements on the surface over dust source areas, such as the effect of dust radiative feedback on reducing near-surface wind speed and surface temperature50,51. Such studies form a foundation for the basic research on dust radiative feedback. Herein, we investigate the dust radiative feedback over the GD on the weather and climate effects of an intensified cyclone in northeast of Mongolia. In addition, strengthening of the cyclone facilitates further eastward transmission of dust aerosols over the GD.

As an important dust source in East Asia, the GD is the main contributor to air pollution in northern China28. The results of the present study indicate that dust radiative feedback over the GD can intensify Mongolian cyclone. The participation of dynamic sources in actual analysis resolves the poor dust simulation of numerical models in the Mongolian Plateau. However, owing to the uncertainty of emission inventory and the deficiency of numerical simulation in this area, the linkage between dust radiative feedback in the GD and Mongolian cyclone is worthy of further exploration in the future.

Methods

Model configuration and sensitivity simulations

To study the impact of dust aerosols on meteorology in China based on aerosol–meteorological interaction, a new-generation regional air quality model WRF-Chem52 coupled with an online meteorological and chemical model was used to simulate a severe dust event occurred in northern China in May 2019. WRF-Chem is a public model developed by National Center for Atmospheric Research and National Oceanic and Atmospheric Administration. Various physical and chemical processes including transport, wet and dry deposition, gas chemistry, radiation, and photolysis are considered in WRF-Chem, which has advantages in the numerical pollution simulation.

The meteorological initial fields of the WRF-Chem model were provided by the National Centers for Environmental Prediction final reanalysis data (NCEP/FNL) with a time interval of 6 h and a resolution of 1° × 1°. Related data of boundary conditions in the numerical simulation were output from the Community Atmosphere Model with Chemistry (CAM–Chem) model. The anthropogenic emission inventory was obtained from Emission Database for Global Atmospheric Research–Hemisphere Transport of Air Pollution (EDGAR–HTAP) global inventory at a horizontal resolution of 0.1° × 0.1° for 2010. EDGAR–HTAP provides specific inventory information for CH4, CO, SO2, NOX, nonmethane volatile organic compounds (NMVOCs), NH3, PM10, PM2.5, black carbon (BC), and organic carbon (OC). Biomass emission data were obtained from Model of Emission of Gases and Aerosol from Nature (MEGAN). Other detailed physical parameterization schemes used in this study are shown in Table 1.

Table 1 WRF-Chem configuration options for physical and chemical parameterizations.
Full size table

WRF-Chem model version 3.9.1 was adopted with Lambert projection and a one-way nested grid. This study used 42°N and 108°E as the center grid, with a spatial resolution of 151 × 191 grid points and a grid resolution of 30 km. Study areas included the TD and GD, which are the main dust source areas in China. Simulation period is 00:00 UTC May 10 to 00:00 UTC May 31, 2019, with the first day used as the spin-up period. The severe dust event occurred from May 11 to May 16, 2019, was selected as the main study period. The main simulation domain is displayed in Supplementary Fig. 1. Dust aerosol is mainly activated by GOCART dust emission scheme in the WRF-Chem. When the dust emission scheme is turned off, dust aerosol is also removed in the modeling. To investigate the feedback between dust-radiation interaction and meteorology, we conducted two parallel experiments with and without dust emission in the study. Two simulations that are respectively, designated as “EXP_CTRL” and “EXP_NODE” are carried out. The EXP_CTRL experiment is defined as a control simulation in which dust emission scheme is turned on and dust aerosols are allowed to provide it’s feedback to the radiation. The sensitivity experiment (“EXP_NODE”) is conducted by closing dust emission scheme and the direct and indirect feedback between dust aerosol and the radiation scheme. The difference between control and sensitivity result are considered as the synoptic-scale effects of dust radiative forcing53,54. This method is also widely used to explore the radiative forcing of different kinds of aerosol and its effects on meteorological fields55,56,57.

Dynamic dust source

To effectively describe dust emission, this study adopted the Georgia Tech/Goddard Global Ozone Chemistry Aerosol Radiation and Transport (GOCART) scheme, in which the dust emission flux (F) is calculated using the following formula:

$$F={CS}{s}_{p}{u}_{10m}^{2}\left({u}_{10m}-{u}_{t}\right)\left({u}_{10m}\ge {u}_{t}\right),$$
(1)

where C denotes the constant of the empirical derivation function, S is the source function of soil wind erosion, sp represents the fraction of each size bin in erodible dust, u10m denotes the horizontal wind speed at 10 m, and ut represents the threshold velocity of particle size, air density and soil moisture. More details about the coupling of the GOCART scheme with the WRF-Chem model can be found in Zhao et al. 58. When GOCART scheme is turned off, dust aerosol is removed in the modeling.

Traditional numerical models treat the potential dust source as static climatic distribution and ignore its dynamic activity, causing considerable uncertainty in the difference in seasonal dust emissions59. Dust sources used in previous GOCART simulations were based on average land covers from the Advanced Very High Resolution Radiometer (AVHRR) satellite, which has no temporal variance60. Therefore, a global-scale dynamic dust source has been developed using time-varying NDVI data61,62. The dust source function (S) is determined by surface bareness (B) and topographic features (H)61,62. Topographic differences are closely related to dust accumulation63. H represents the topographic differences in the grid relative to the surroundings, as shown in Eq. (2), where the lower H is, the easier the dust is to accumulate.

$$H={\left(\frac{{Z}_{\max }-{Z}_{i}}{{Z}_{\max }-{Z}_{\min }}\right)}^{5}.$$
(2)

In this equation, Z represents the terrain elevation in the grid, Zmax and Zmin correspond to the highest and lowest points in the surrounding areas (calculated grids) of the 10° × 10°, respectively, and Zi denotes the height of grid i. Power of 5 is used to improve the terrain contrast. Bare soil is indicated when the normalized difference vegetation index (NDVI) is <0.17 in each grid. Therefore, the surface bareness is calculated using the ratio of grid numbers below 0.17 (\({{\rm{N}}}_{{{ < }}{{0}}{{.}}{{17}}}\)) to the total grid number \({{\rm{N}}}_{{\rm{total}}}\):

$$B=\frac{{N}_{ < 0.17}}{{N}_{{total}}}.$$
(3)

However, the dynamic dust source function based on this calculation method is not accurate and the regions with perennial ice and snow cover at high latitudes also showing dust source function maximum. We also employed land cover from MODIS to constrain the calculated S to reproduce the spatial distribution of dust sources. Therefore, a global-scale dynamic dust source has been developed61,62.

Aerosol optical depth

The WRF-Chem can effectively simulate aerosol optical depth, which is vertically derived from the accumulation of the extinction coefficient. The aerosol extinction coefficient is determined by the following formula64,65:

$${b}_{{ext}}(\lambda ,X)=\mathop{\sum }\limits_{i=1}^{i=m}{Q}_{e}({x}_{i},X)\pi {r}_{i}^{2}n({r}_{i},X)$$
(4)

where \({{\rm{Q}}}_{{\rm{e}}}\) shows the exctinction efficiency, \({{\rm{r}}}_{{\rm{i}}}\) is the wet radius, \({{\rm{x}}}_{{\rm{i}}}\) is the size parameter, and \({\rm{n}}({{\rm{r}}}_{{\rm{i}}},{\rm{X}})\) is the number concentration (number per unit volume) associated with the section “i”.

AOD mainly depends on particle size distribution, composition, mixing state and hygroscopic property, and is calculated according to aerosol concentrations and aerosol optical properties66,67. It was obtained from WRF-Chem by vertical integration (from the ground to the top of the domain) of the aerosol extinction at 550 nm, which was obtained as the output of the WRF-Chem68.

$${AOD}={\int }_{0}^{\infty }{b}_{{ext}}(z){dz}$$
(5)

Aerosol radiative forcing

Direct radiative forcing of aerosols refers to radiation change caused by scattering and absorption. That is, direct radiative forcing of dust can be defined as the difference between the net radiative flux with and without dust, whereas atmospheric radiative flux refers to the radiative energy passing through any surface per unit time. Downward and upward radiation fluxes are characterized by positive and negative values, respectively. As calculated using the following equations, dust radiative forcing is always estimated at the TOA, in the atmosphere (ATM), and at the bottom of the atmosphere (BOA):

$${{DRF}}_{{ATM}}={{NF}}_{{ATM}}^{{with}}-{{NF}}_{{ATM}}^{{no}},$$
(6)
$${{DRF}}_{{TOA}}={{NF}}_{{TOA}}^{{with}}-{{NF}}_{{TOA}}^{{no}},$$
(7)
$${{DRF}}_{{BOA}}={{NF}}_{{BOA}}^{{with}}-{{NF}}_{{BOA}}^{{no}},$$
(8)
$${NF}={{NF}}_{\downarrow }-{{NF}}_{\uparrow },$$
(9)

where DRF represents direct radiative forcing and NF indicates the net radiative flux69,70,71.

Thermodynamic equations

The following thermodynamic equation were used to accurately evaluate the effect of dust radiative feedback on thermal structure in the meteorological field. We defined the experiment in EXP_NODE as the mean quantity in thermodynamic equation, and the difference between the two parallel experiments was considered as disturbance. In the actual analysis, each terms in the following equation were integrated from the surface to 300 hPa.

$$\begin{array}{l}{\int }_{Ps}^{p0}\left\{{\underbrace{-U^{\prime}\left(\frac{\partial \bar T}{\partial x}\right) -\bar {U}\left(\frac{\partial T^\prime} {\partial x}\right) -U^{\prime} \left(\frac{\partial T^\prime }{\partial x}\right)-V^{\prime} \left(\frac{\partial \bar T}{\partial y}\right)-\bar{V}\left(\frac{\partial T^\prime}{\partial y}\right)-V^{\prime} \left(\frac{\partial T^\prime }{\partial y}\right)+\sigma ^{\prime} \bar{\omega }+\bar{\sigma }\omega ^{\prime}}_{(adiabatic\,change)}}\right.\\ \left.+{\underbrace{\frac{\varDelta Q}{Cp}}_{(nonadiabatic\,change)}}\right\}\frac{dp}{g}=0\end{array}$$
(10)

The atmospheric apparent heat source Q1 is calculated using an inverse algorithm based on the thermodynamic equation72:

$$Q1=\,{C}_{p}\frac{\partial T}{\partial t}\,-{C}_{p}\left(\omega \sigma \,-\,V\,{{\cdot }}\,\nabla T\right).$$
(11)

Data

The MODIS cloud image provides effective records about the frequency and distribution of cloud. It can be obtained from the online website in png format. As a vital sensor onboard Terra and Aqua, the Moderate Resolution Imaging Spectroradiometer (MODIS) provides reliable global information on clouds, aerosols, land cover, and other parameters73 and always plays an important role in module evaluation. Its AOD products are at a wavelength of 550 nm with dark target and deep blue algorithms (MOD08 and MYD08) and contain daily degree grid average values of atmospheric parameters at a resolution of 1° × 1°. The aerosol index (AI) derived from Aura Ozone Monitoring Instrument (OMI) covers the period 2004–2023 and exhibits a horizontal resolution of 1° × 1°. The AI is extremely sensitive to ultraviolet (UV)-absorbing aerosols, such as smoke, mineral dust, and volcanic ash74. OMI is the successor to the total ozone mapping spectrometer (TOMS) instrument for monitoring ozone levels, air quality, and climate of the earth75. Moreover, the daily aerosol ground observations at sites are derived using the Aerosol Robotic Network (AERONET), which is a wide network of sun photometers located worldwide76.