Introduction

Aerosols play an important role in the radiative balance of the Earth’s atmosphere system, making it crucial to quantify the changes in radiation and heat budget caused by aerosols for a better understanding of synoptic or climate variations1. Aerosol feedback can reduce solar radiation and planetary boundary layer (PBL) wind speeds2,3,4, alter horizontal and vertical wind speeds5,6, and cause abnormal cyclonic circulations7,8,9, and pose a threat to human health10. The climatic effects of aerosols have raised widespread concerns in the scientific communities11,12,13,14. However, there are few studies based on dynamics, especially lacking in research on whether it will lead to the appearance of secondary circulation (SC). SC plays a vital role in atmospheric motion15,16, and it is always associated with synoptic processes and microclimate to some extent17,18,19. In addition, due to the complexity of the synoptic or climatic system and the interactions between atmospheric pollutants and synoptic/climatic systems, assessing the aerosol impact on climate is still a challenge under the current polluted air conditions20,21,22,23.

Here, we propose a two-way interaction dynamical mechanism between the SCs and PM concentrations. Figure 1 presents a schematic view of the SCs evolution (including the process of reinforced SC forming and SC occurrence) influenced by aerosol feedback. Suppose that the region is under the westerly wind regime, since heavy aerosol pollution tends to reduce wind speed, the closer to the center of heavy pollution, the weaker the wind speed is. Consequently, the upstream area of the heavy pollution region forms a convergence zone, and the downstream area becomes a divergence zone, which would generate vertical motions to maintain mass conservation. Under favorable synoptic-scale background meteorological conditions, the downward and upward airflow strengthen further, leading to the formation of a convergence zone above the divergence zone. Eventually, a strong and well-structured clockwise rotating SC and a weaker counterclockwise rotating SC are formed. But when the aerosol feedback impact is not very strong, often only the well-structured clockwise rotating SC appears.

Fig. 1: Schematic view of the mesoscale SCs evolution influenced by aerosol feedback effect in the atmospheric boundary layer.
figure 1

The green dashed circles are reinforced SCs, the white dashed circle is SC, the blue arrow lines denote wind speeds, and the dashed read arrow lines indicate the wind speed difference between feedback and non-feedback effect. The solid green thick arrow lines represent the downward and upward motion induced by low atmospheric boundary layer (ABL) divergence or convergence. The shaded color region stands for PM concentrations and ∇Vh is horizontal divergence \(({{\rm{s}}}^{-1},{\rm{div}}=-{\nabla }_{h}\vec{V})\). Convergent and divergent areas caused by the feedback effects are marked as shallow green oval areas.

Full size image

Beijing-Tianjin-Hebei is one of the most severely polluted areas in China24,25. The emission actions in Beijing–Tianjin–Hebei (BTH) and its surrounding regions during the 2014 Asia-Pacific Economic Cooperation (APEC) Summit in Beijing, resulted in a significant improvement in air quality, and the air tended to be cleaner26,27. This provided unique data and scenarios for better understanding of the impact of particulate matter28,29 on SC. To further verify the reliability of the results, two emission reduction events in 2015 and 2017 in the BTH30 were also selected. Due to the implementation of emission policies in recent years, PM2.5 concentrations have been declining31,32. Therefore, January and July of 2014, 2017, and 2020 were selected as the validation of long-term emission reduction policies, rather than emission reduction events. WRF-Chem model was used to simulate the influence of aerosol feedback on SC, further details of model setup in methods.

Results and discussion

Association between PM2.5 and SCs

On the heavily polluted day closest before APEC without emission reduction (BAPEC, Supplementary Table 2), in the morning of October 30th (Fig. 2a), the planetary boundary layer height (PBLH) was low and the ground-level PM2.5 concentrations were high. A clockwise rotating SC can be identified in the vicinity of Beijing with a horizontal scale of about 100 km and a vertical scale of about 1 km. Due to the upward transport by the SC, a high PM2.5 plum occurred in the western side of SC. In the afternoon (Fig. 2b), SC still maintained, on a larger spatial scale compared to the morning, PM2.5 concentrations continued to rise as they were not easily dispersed. By evening (Fig. 2c), the SC weakened and shifted eastward, resulting in a decrease in ground-level concentrations in the Beijing area. At night (Fig. 2d), SC further weakened and PM2.5 concentration continued to decrease. On the most severe heavily polluted day during APEC with emission reduction (DAPEC), in the morning of November 4th (Fig. 2e), there was a SC, but its structure was not as clear as that of the morning in BAPEC, and PM2.5 concentrations were low. In the afternoon (Fig. 2f), only one smaller and weaker scale SC occurred compared to the morning. PM2.5 was notably smaller than those in the BAPEC. In the evening (Fig. 2g) and at night (Fig. 2h), SC was not detectable, and PM2.5 concentrations were significantly lower than that of BAPEC. On the heavily polluted day closest after APEC without emission reduction (AAPEC), in the morning of November 15th (Fig. 2i), a smaller SC can be discerned. In the afternoon (Fig. 2j), a well-defined and larger-scale SC occurred. In the evening (Fig. 2k), the SC further intensified. But at night (Fig. 2l), SC weakened and PM2.5 concentration decreased.

Fig. 2: Feedback simulation of vertical cross-section of modeled hourly resultant vector winds (U-W, m s−1) and PM2.5 concentrations in three representative heavily polluted days.
figure 2

a 08:00 LST, October 30th in BAPEC, b same as a but for 14:00 LST, c same as a but for 20:00 LST, d same as a but for 02:00 LST, October 31st, e 08:00 LST November 4th in DAPEC, f same as e but for 14:00 LST, g same as e but for 20:00 LST, h same as e but for 02:00 LST, November 5th, i 08:00 LST November 15th in AAPEC, j same as i but for 14:00 LST, k same as i but for 20:00 LST, l same as i but for 02:00 LST, November 16th. The red solid line denotes the PBLH. The white dashed circles highlight SCs, and the green dashed line is the region to the flat area (east of the 47th grid point) and below 2 km altitude in the vertical direction. The three black dots from left to right on the bottom of each figure stand for the positions of BJP, BTP, and TSP, respectively (Supplementary Table 3).

Full size image

The mesoscale SC mainly occurs in the region to the flat area and below 2 km altitude in the vertical direction. The SC during DAPEC was significantly weakened with the evolution from its occurrence in the morning, weakening at afternoon, and vanishing in the evening and night. In contrast, during BAPEC and AAPEC, SC also strengthened throughout the daytime. PM2.5 pollution during DAPEC was significantly lower than that during BAPEC and AAPEC, indicating a strong relationship between the PM2.5 levels and mesoscale SCs.

We use the time levels to represent how many times the SC could be identified in specific hour in specific period. The time levels of SCs (Fig. 3a) during DAPEC significantly decrease to ~37.3% of the periods during BAPEC and AAPEC (Supplementary Table 9). The simulated PM2.5 (Fig. 3b) during DAPEC are ~47.6% of BAPEC and AAPEC, and the measured concentration ratio at 47.7%. We further used the mean value of U-W vorticity (Vorj, 10−3 s−1) on vertical profile to quantify the occurrence and strength of SC and set a threshold value of 3.0 × 10−3 s−1 to measure changes in SC. We found that the absolute value of Vorj (Vorj_Abs) averaged over DAPEC reduced by 52.8% from BAPEC and AAPEC (Fig. 3c). And the result reveals that the accuracy of SC occurrence determined by vorticity (Fig. 3d) compared to manual identification was 79.2%. During DAPEC, PM2.5 concentrations, its precursors (e.g., NO, SO2, and TOL) and Vorj_Abs all decreased compared to the BAPEC and AAPEC, and their trends over time were quite similar (Supplementary Fig. 2).

Fig. 3: Statistical data on variables related to secondary circulation in the Feedback simulation.
figure 3

a Time levels with the occurrence of mesoscale SCs, during BAPEC, DAPEC, and AAPEC at 08:00 LST, 14:00 LST, 20:00 LST, and 02:00 LST, b same as a but for simulated PM2.5 concentrations, c same as a but for simulated Vorj_Abs, d same as a but for the number of hours with the occurrence of mesoscale SCs judged by Vorj, e time levels with the occurrence of mesoscale SCs judged by Vorj, before, during and after the emission reduction period in 2014 APEC, 2015 and 2017 events, f same as e but for simulated PM2.5 concentrations, g same as e but for simulated Vorj_Abs, h same as e but for the numbers of easterly U-winds, i same as e but for the frequency of observed easterly U-winds, j time levels with the occurrence of mesoscale SCs judged by Vorj, in January and July of 2014, 2017 and 2020, k same as j but for simulated PM2.5 concentrations, l same as j but for simulated Vorj_Abs, m same as j but for the numbers of easterly U-winds, n same as j but for the frequency of observed easterly U-winds.

Full size image

Compared to the feedback simulation, the time levels of SCs are slightly lower in the non-feedback simulation (Supplementary Fig. 3, Table 10), and the largest variation in SC numbers occurred during the afternoon, which should be related to the stronger solar radiation in that afternoon (Supplementary Table 11). Since aerosol feedback acts directly to alter radiation33,34, the stronger incoming radiation in afternoon leads to a greater forcing aerosol effect on SC, and it should even induce the occurrence of SC. The results indicate that SC is likely to be reinforced during periods with high PM2.5 concentrations and stronger radiation.

During the other two emission reduction events in 2015 and 2017 (Supplementary Fig. 4), PM2.5 concentrations and the number of identified SCs significantly decreased during the emission reduction period, compared to the before and after periods without emission reduction. In addition, PM2.5 precursors concentrations (e.g., NO, SO2, and TOL) and Vorj_Abs are all decreased during the emission reduction period, and their trends were also quite similar (Supplementary Figs. 5 and 6). For all the times of three emission reduction events (Fig. 3e–i), the time levels of SCs occurrence measured by Vorj, PM2.5 concentrations, and Vorj_Abs were significantly decreased during emission reduction period. Moreover, the surface U-wind direction also shows changes corresponding to the number of SCs or U-W vorticity. The main manifestation is that after the strengthening of SCs by aerosol feedback, the frequency of easterly U-winds increases, except for the D_15ER period of observed winds, which should be related to the background wind field not being dominated by westerly winds. Since in recent years, most of PM2.5 are secondary in Beijing. Moreover, we also noticed that the average proportion of secondary aerosols for the all times during the 2014 APEC events was 59.2%, and the proportion above 70% often occurred (Supplementary Fig. 7), the correlation coefficient between Vorj and secondary PM2.5 (R = 0.239, which passed the two-tailed significance T-tests, P < 1.0 × 10−10) exceeds the correlation coefficient with total PM2.5 (0.227). However, the proportion of secondary PM2.5 in the 2015 and 2017 reduction events was not high, which should be due to the higher proportion as winter approached23,35. In winter, secondary aerosols may have a better correlation with SC.

From 2014 to 2020 when emission reduction actions were taken36, the emissions of PM2.5 and its precursors (SO2, NOx, VOC, and NH3) have also shown a downward trend, and the observed concentration of PM2.5 has also decreased (Supplementary Fig. 8). The time levels of SCs occurrence by Vorj, PM2.5 concentrations, and Vorj_Abs were decreased in January and July from 2014 to 2017 and then to 2020 (Fig. 3j–n), and PM2.5 precursors concentrations (e.g., NO, SO2, and TOL) also decreased (Supplementary Figs. 9 and 10). The surface U-wind direction also shows changes corresponding to the number of SCs or U-W vorticity (Fig. 3m, n). But the relationship between simulated wind and SC is not better than observed, that should be because the input data for WRF are once every 6 h, which cannot more accurately reflect the hourly variation of surface wind like the observed hourly data. Additionally, some inconsistencies should originate from the background wind field.

Potential mechanism of aerosol feedback influence on SCs

During the heavily polluted days (Fig. 4a–k), reinforced SCs caused by aerosol feedback can be clearly discerned at each time. Comparing DAPEC with the other two periods, the vertical wind speed is significantly lower, resulting in weaker reinforced SCs in wind speed and smaller scales in spatial structure. When the westerly wind (Fig. 2f–i) prevails, the reinforced SCs result in convergence and divergence in the downwind areas of high PM2.5 concentrations (upper right in Fig. 4f–i). When the predominant wind direction is easterly, due to the blocking effect of Beijing’s Western Hills on aerosols, high concentrations primarily accumulate at the foothills of the Western Hills and do not significantly transport downstream37,38 (Fig. 4a–c, g, k), resulting in a strong convergence zone in the upstream region due to the reduction in wind speed, along with the forced upward flow due to the blocking effect of the mountains. This leads to the occurrence of a strong reinforced SC in the upstream region (upper right in Fig. 4a–c, g, k). However, in the downstream region where the high concentration decreases, which is difficult for a mesoscale SC with a certain spatial scale to be reinforced. It means that the occurrence of reinforced mesoscale SC corresponds well with PM2.5 concentrations. The reinforced SC generally appears in the upper right of areas with high PM2.5 concentrations, which is in line with the mechanism for the impact of aerosol feedback on SCs (Fig. 1).

Fig. 4: Vertical cross-sections of differences of hourly U-W wind vectors between feedback and non-feedback simulations and PM2.5 concentrations from feedback simulations in three representative heavily polluted days.
figure 4

a 08:00 LST, October 30th in BAPEC, b same as a but for 14:00 LST, c same as a but for 20:00 LST, d same as a but for 02:00 LST, October 31st, e 08:00 LST November 4th in DAPEC, f same as e but for 14:00 LST, g same as e but for 20:00 LST, h same as e but for 02:00 LST November 5th, i 08:00 LST November 15th in AAPEC, j same as i but for 14:00 LST, k same as i but for 20:00 LST, l same as i but for 02:00 LST November 16th. The blue dashed circles highlight reinforced SCs. The three black dots on the abscissa denote, respectively, BJP, BTP, and TSP from left to right (Supplementary Table 3).

Full size image

For the average over BAPEC and AAPEC, high PM2.5 concentrations (Supplementary Fig. 11) correspond well with reinforced mesoscale SCs (Supplementary Fig. 12), but it is not straightforward to identify well-structured reinforced mesoscale SCs during DAPEC. Moreover, compared to the hourly U-W differences (Fig. 4), the SC structure averaged over the three periods (Supplementary Fig. 12) is more complete, with larger spatial scales and better correspondence with PM2.5 concentrations. This is likely because the mean U and W winds averaged over the three periods filter out some incidental factors that may disturb SC, focusing mainly on the relationship between aerosol feedback and SCs.

Relationship between SC and synoptic background

Different synoptic background wind fields could result in different characteristics of reinforced SCs. When westerly flow prevailed (Fig. 4g, Supplementary Fig. 12a), there was a better match between the reinforced SCs and the conceptual SC (Fig. 1). When easterly background flow prevailed (Fig. 2b, c), the updrafts in the reinforced SCs were often stronger than the downdrafts. The reinforced SCs were primarily located near the top of PBL and were associated with areas with high PM2.5 concentrations. The SC occurred is the result of the combined effect of background and reinforced SC. The synoptic charts (500-hPa, 850-hPa, surface, and surface wind field map, Supplementary Figs. 13–18) indicate that the synoptic situation were not the primary cause of SCs and its vertical motion at the low atmosphere (Supplementary Information provide details). Although trough and surface convergence should contribute to high concentrations and evolution of SCs37,39,40. Comparing the feedback and non-feedback simulations, we can identify similarities of hourly wind fields as well as geopotential heights and temperatures (Supplementary Figs. 19–22). As some studies have previously pointed out, the impact of aerosol feedback on wind field is about an order of magnitude smaller than that of the background wind field2,41,42. This is consistent with the average situation (Supplementary Fig. 12) of our research. However, on hourly basis, the impact of aerosol feedback on wind field is of the same order of magnitude as the background wind field (Figs. 2 and 4). Furthermore, the R between Vorj and meteorological scalar fields (e.g., Div, wind speed and its gradient and temperature and its gradient, etc. Supplementary Tables 12–16) are notably lower than those for PM2.5 and Vorj. Notably, the R between Vorj below 1 km and PM2.5 reaches up to 0.717 (P < 0.001) in afternoon. This suggests that while synoptic background does influence the ___location and characteristics of SC, aerosol impact should be the primary factor associated with SC occurrence.

The statistical results (Fig. 3h, i, m, n) show that the influence of background meteorological fields on SCs is also very important. In the BTH area, the upper atmosphere (e.g., 500-hpa, Supplementary Figs. 13, 14, 23, 24) is mainly dominated by westerly winds43,44, but the prevailing wind direction near the top of the boundary layer (e.g., 850-hpa, Supplementary Figs. 15, 25, 26) and on the surface may not necessarily be mainly westerly (Supplementary Figs. 16–18, 27–31). When the dominant wind direction of the background is westerly (Supplementary Figs. 16, 27d–f, 30a–c), the SC generated by aerosol feedback is clockwise, resulting in an increase in the frequency of surface easterly U-winds, and when easterly winds prevail on the surface (Supplementary Figs. 27b, d, 30e, f), it will result in a decrease in the frequency of surface easterly U-winds. For example, westerly winds are prevalent on the surface and 500-hpa in 2014 APEC and 2017 emission reduction events, during the emission reduction period, SC decreases and the frequency of easterly U-winds decreases. For January 2014, 2017 and 2020, modeled wind direction is easterly but observed as westerly, so the frequency of easterly U-winds increases in modeled but decreases in observed. Although there is a change in wind direction on the surface opposite to the background, which is not as good as the correlation between PM2.5 and SC. However, for most periods, the observed wind corresponds better to the simulated SC than the simulated wind, which should be attributed to the fact that the observed wind speed contains more information about hourly variations.

In summary, heavy PM2.5 pollution-induced aerosol feedback significantly enhances the evolution of mesoscale SCs, especially when solar radiation is strongest. It is primarily the aerosol feedback that serves as the main factor driving SC occurrences, although the synoptic background influences the ___location and characteristics of SC, and the frequency of surface U-winds opposite to the prevailing background wind associated with the change of SCs by aerosol feedback. SCs are constrained by surface, meteorological, dynamic, and thermodynamic conditions within the PBL45. The locations and sizes of SCs associated with aerosol concentration distribution under aerosol feedback. Since our proposed two-way interaction dynamical mechanism between the SCs and PM concentrations, is based on the universal phenomenon of aerosol reducing wind speed and constrained by the law of conservation of air mass, the relevant conclusions should also have some reference significance in other places.

However, there are still uncertainties in the impact of aerosol feedback on weather or climate, such as aerosol-cloud interactions46, the uncertainties in pollutant emission inventories47 and WRF input data which contain limited hourly meteorological variations. Additionally, the mechanism we proposed is the main process by which aerosols affect SCs, rather than the full chain of events associated with the feedback. The non-feedback simulation suggests that at least half of SCs should be attributed to the aerosol feedback effects during APEC non-emission reduction period. The observed PM2.5 and wind direction have a better correspondence with the simulated SC or U-W vorticity than the simulated PM2.5 and wind direction. These imply that the WRF-Chem model has shortcomings and might underestimate the impact of aerosol feedback48 on SCs. Additionally, since the observed atmosphere is already polluted rather than clean air, the measured meteorological conditions might be subject to polluted air. As a result, the impact of anthropogenic aerosols on mesoscale SCs are likely to be more severe than the present research results. Despite significant efforts to improve air quality, air pollution in China continues to be a major environmental concern49. The emergence of SCs caused by aerosol feedback may disturb the weather system, especially under heavy pollution conditions, which may increase the uncertainty of weather simulation. There should be many aspects that need to be studied in the two-way interactions between aerosols and meteorology and associated feedback6,11,48,50. The better we understand aerosol feedback, the more effectively we can mitigate air pollution as well as mitigate climate forces and adapt to climate change.

Methods

Data introduction

The final operational global analysis data from national centers for environmental prediction (NCEP-FNL) with a resolution of 1° × 1° and a 6-h interval, were used to generate boundaries and initial conditions for WRF-Chem. The MEIC emission inventory for China was utilized, which has a resolution of 0.25° × 0.25° lat/lon51. The present study also incorporated MEGAN biogenic emissions data52 and biomass burning emissions data14. The initial conditions and boundary conditions for pollutant concentrations in each simulation were provided by the global model MOZART53. The meteorological observation data are sourced from National Oceanic and Atmospheric Administration Integrated Surface Database (ISD-Lite). Pollutant monitoring data are obtained from the National Urban Air Quality Real-time Release Platform (https://air.cnemc.cn:18007/), which collects hourly PM2.5 concentration data. A total of 48 pollution monitoring sites within the model Domain 3 were selected and the site locations and names are provided in Supplementary Table 1.

During the APEC, extraordinary emission reduction measures were implemented in key regions (including BTH region and the western part of Shandong Province, including Liao Cheng, Dezhou, Jinan, Binzhou, Dongying, and Zibo, as delineated by the blue boundary line in Supplementary Fig. 1a). Significant changes in pollutant emissions occurred during the DAPEC period. The emission reduction ratios for various pollutants during the DAPEC period were referred to relevant literatures27,54. The emission for PM2.5, PM10, SO2, NOx, NO2, CO, and VOCs in Beijing were decreased by 39.2%, 66.6%, 61.6%, 49.6%, 39.9%, 28.8%, and 33.6%, respectively. Due to the lack of comprehensive emission reduction data for other key regions apart from Beijing, a rough estimate was made by assuming that the emission reduction ratios in Tianjin were 90% of those in Beijing, and for other regions, they were 80% of those in Beijing. In the emission reduction action from August 20 to September 3, 2015, the emission reduction ratios for PM2.5, PM10, SO2, NOx, NO2, CO, and VOCs in Beijing were 73.0%, 69.0%, 46.0%, 52.0%, 52.0%, 38.0%, and 53.0%, respectively. In the emission reduction action from May 5 to 15, 2017, the emission reduction ratios for PM2.5, PM10, SO2, NOx, NO2, CO, and VOCs in Beijing were 40.0%, 25.0%, 25.0%, 25.0%, 25.0%, 20.0%, and 25.0%, respectively. The emission reduction ratios of Tianjin were 90% of those in Beijing, and for other regions, they were 80% of those in Beijing for both emission reduction events.

Model setup

The WRF-Chem 4.1.3 version, a regional air quality model that fully couples meteorology and chemistry online55, was used in this study. The Lambert projection was employed, and the model setup consisted of three nested domains (Supplementary Fig. 1) with grid nesting ratios of 1:3:3. The horizontal resolutions and grid sizes for each ___domain were as follows: Domain 1: 27 km, 164 × 140 grids; Domain 2: 9 km, 156 × 153 grids; Domain 3: 3 km, 123 × 123 grids. Domain 1 covers the mainland of China, Domain 2 encompasses the North China Plain, Shanxi Province, parts of Shaanxi Province, and Inner Mongolia, and Domain 3 encircles Beijing, Tianjin, and surrounding areas. The model atmosphere has 28 vertical layers.

The physics parameterizations used in the model include the Morrison (2 moments) microphysics scheme56, the Improved Rapid Radiative Transfer Model for General Circulation Models for longwave and shortwave radiation57, the Monin-Obukhov surface layer scheme from MM558, the Unified Noah land-surface scheme59, the MYNN 2.5 level Turbulence Kinetic Energy boundary layer scheme60,61, and the Grell 3D ensemble cumulus parameterization with the Xu-Randall method for calculating cloud optical thickness62. The chemical module employed the MOZART mechanism63, and the MOSAIC module with the KPP library is used for aerosol chemistry processes64.

The simulation period extends from October 18th, 2014, 00:00 UTC to November 25th, 2014, 00:00 UTC, divided into three simulation periods, each lasting 14 days, in 2014 APEC. In the 2015 and 2017 emission reduction events, the simulation also divided into three simulation periods, each lasting 17 or 9 days, the period extends from August 3th, 2015, 00:00 UTC to September 19th, 2015, 00:00 UTC, and from April 30th, 2017, 00:00 UTC to May 23th, 2017, 00:00 UTC, respectively. The starting time of each subsequent period is 2 days prior to the end time of the previous period. The first two days of each simulation period are designed as the spin-up period to allow the model’s meteorological processes to reach a quasi-steady state. In the simulations conducted in January and July 2014, 2017, and 2020, each simulation lasted for 34 days, with the first three days used as the spin-up period and the remaining 31 days of each month used as the analysis data. The direct feedback of aerosols refers to their scattering and absorption of radiation, while the indirect feedback refers to the influence of aerosols on radiation after participating in cloud processes. In the feedback configuration, the bi-directional feedback mechanism between meteorology and pollutants was turned on, considering both direct and indirect effects of aerosols. While in the non-feedback configuration, both the direct and indirect feedback effects were turned off, and the cloud droplet concentration was set to a constant value of 250 cm−3 2,41,65.

Model simulations

This study investigates the impact of aerosol feedback on SC by comparing the simulation results of feedback and non-feedback configuration. The period from 08:00 LST on November 1st to 20:00 LST on November 12th, 2014, is considered as the APEC emission reduction period (referred to as DAPEC), while the period before APEC with the same duration (from 08:00 LST on October 20th to 20:00 LST on October 31st) is referred to as BAPEC, and the period after APEC with the same duration (from 08:00 LST on November 13th to 20:00 LST on November 24th) is referred to as AAPEC, respectively (Supplementary Table 2). Periods of other two emission reduction events are also listed in Supplementary Table 2. Considering that heavy air pollution occurs typically at the city scale, and the SC subject to aerosol feedback also should mainly occur at the city scale, the innermost nested ___domain (Domain 3) is selected for the study.

To visualize the SC, this study uses vertical cross-section diagrams that span Beijing and Tangshan cities along the Y-axis. The U-W composite wind vector profile is visualized by multiplying the vertical velocity W by 20 and then combining it with the U wind component (zonal wind) to highlight the influence of vertical velocity and better identify the boundary layer SC structure. We estimated the differences of meteorological variables and PM2.5 concentrations between feedback and non-feedback simulations at each time step. Three model grids were selected to plot wind profile: one located in the center of Beijing (referred to as BJP or the urban area), one in the center of Tangshan (referred to as TSP or the urban area), and one located between Beijing and Tangshan (referred to as BTP or the suburban area). Of these grids, BTP is located 75 km east of BJP, and TSP is located 75 km east of BTP. Supplementary Fig. 1b shows the positions of the three grid points, and detailed information about their locations can be found in Supplementary Table 3. Detailed model verifications are presented in the Supporting Information (Supplementary Tables 4–7). Model validation shows that WRF-Chem demonstrates strong simulation capabilities in both feedback and non-feedback simulations66,67, the simulation of PM2.5 and meteorological factors under feedback configuration is often better than that under non-feedback configuration, which also indicates the importance of aerosol feedback.

U-W vorticity, divergence, and gradient

Vorticity is defined as the curl of the velocity field and serves as a microscopic measure of fluid rotation characteristics. The spatial distribution of vorticity is referred to as the vorticity field, which is a vector field. Since SC mainly occurs in the flat area east of the Western Hills of Beijing in our study area, to eliminate the influence of the Western Hills on Vorj (Formulas in Supplementary Table 8), the U-W plane’s average Vorj value was selected in the region to the east of the 47th grid point and below 2 km altitude in the vertical direction. To determine whether SC occurs by Vorj, we selected regions on the U-W plane with relatively flat topography. We consider the absolute mean value of Vorj greater than or equal to at threshold value of 3.0 × 10−3/s as the criterion for judging its presence. If the absolute value exceeds the threshold, secondary circulation is considered to be present. Otherwise, no SC occurs. The U–V divergence and vorticity and horizontal gradient of a scalar field were also calculated for the statistical testing of variables related to SC (Formulas in Supplementary Table 8).

The U-W vorticity (Vorj) was estimated using a central differencing scheme with second-order accuracy, where i and j represent ith grid along the east-west and jth level along vertical directions, respectively. mapfac_u and mapfac_w are the map factors for U and W. The vertical resolution ∆z is 55 m, and the horizontal resolution ∆x is 3000 m. Positive value of Vorj (1/s) indicates clockwise circulation and negative value of Vorj represents counterclockwise circulation. For the convenience of writing, the unit is set as 10−3/s in this paper. The U-V divergence (Div, 10−5/s) and vorticity (Vor 10−5/s) were also estimated using a central differencing scheme with second-order accuracy, where mapfac_v is the map factors for V, and the horizontal resolution ∆y is 3000 m. The horizontal gradient of a scalar field was calculated using a central differencing scheme with second-order accuracy, where mapfac_f is the map factors for f, these equations are used to calculate the gradients of meteorological scalar fields such as wind speed, temperature, and pressure (Formulas in Supplementary Table 8).