Interaction between aerosol and thermodynamic stability within the PBL during the wintertime over the North China Plain: Aircraft observation and WRF-Chem simulation

. Aerosol-planetary boundary layer (PBL) interaction has been proposed as a key mechanism for stabilizing the atmosphere and exacerbating surface air pollution. Although the understanding of this process has progressed enormously, its magnitude and impact remain uncertain and vary widely concerning aerosol types, vertical distributions, synoptic conditions, etc. In this study, our primary interest is to distinguish the aerosol-PBL interaction of absorbing and scattering 20 aerosols under contrasting synoptic patterns and aerosol vertical distributions. Detailed in-situ aircraft (KingAir-350) measurements and online coupled model Weather Research and Forecasting with Chemistry (WRF-Chem) simulations are explored over the North China Plain (NCP). Furthermore, a long-term PBL stability trend from 1980 to 2020 over the NCP is also investigated. The aircraft measurements and surface observations show that the surface air pollution over the Baoding City on 3 January is heavier than that on 4 January, 2020. In addition, the aerosols are restricted to the low layer on 3 25 January, whereas the aerosols mix more homogeneous upwards on 4 January. Thereupon, we focus on the two days with distinct synoptic circumstances, PBL stability, and aerosol vertical distributions over the NCP. According to the WRF-Chem modelling, the synoptic pattern over the Baoding City differs between the two days. The prevailing wind direction is opposite with a southwest wind on 3 January and a northeast wind on 4 January. The results indicate that the synoptic condition may affect the PBL thermal structure, thus affecting the aerosol vertical distribution. Additionally, the sensitive numerical experiments reveal that the light-absorbing and light-scattering aerosols have different effects on altering the PBL thermal structure. The inhibition effect of scattering aerosols on the PBL appears to be independent of the aerosol height distribution and solely depends on its concentration. However, aerosol-PBL feedback of absorbing aerosols is highly dependent on its vertical distribution. Our analysis highlights that we should principally concentrate on controlling the emissions of scattering aerosols under the stable stratification while cooperating to control the emissions of scattering and 35 absorbing aerosols in an unstable stratification. Moreover, the long-term inter-annual variation in PBL stability shows a strong correlation with the East Asian Winter Monsoon, which seems to be valuable in determining which pollutants to target in different monsoon years and attaining more precise air pollution control. Based on the numerical simulations and observational constraints, a concept scheme description has been concluded to deepen our recognition of the interactions between thermodynamic stability and aerosols within the PBL over the NCP region. The sensitive numerical experiments reveal that the absorbing and scattering aerosols have different effects when modifying the PBL thermal structure, and the AREs are also influenced by the aerosol vertical distribution. Aerosol-PBL feedback of absorbing aerosols is highly dependent on its vertical distribution. In the case of stable weather patterns under the influence of the south wind, the absorbing aerosols concentrated in the low layer have a strong radiative heating effect on the atmosphere (stove effect). This stove effect disturbs the stratification stability of the low layer, promotes the development 415 of the PBL in coordination with the surface sensible heat flux, and transports heat and pollutants to the upper layer through turbulence. In the case of unstable weather patterns under the influence of the north wind, the absorbing aerosols in the upper atmosphere heat the atmosphere (dome effect) and inhibit the development of the PBL, which further leads to the accumulation of surface pollutants. The inhibition effect of scattering aerosols on the PBL is independent of the aerosol height distribution but only depends on its concentration. Scattering aerosols reflect solar radiation back to the space 420 (umbrella effect), thus weakening the surface sensible heat flux and upward heat transfer, inhibiting the development of the PBL and further accumulating surface pollutants. In summary, the umbrella, stove and dome effects of the scattering and absorbing aerosols under different synoptic patterns are illustrated in Fig. 15. The conclusions enlighten that we should mainly control the emissions of scattering aerosols under the stable stratification, and cooperate to control the emissions of scattering and absorbing aerosols in an unstable stratification, especially to prevent the accumulation of surface pollutants 425 caused by the dome effect of absorbing aerosols.

impeding the PBL development, which is known as the umbrella effect (Ma et al., 2020;Shen and Zhao, 2020). On the other hand, light-absorbing aerosols can trap solar energy and heat the atmosphere, which has different effects depending on the aerosol vertical distribution. For one thing, the absorbing aerosols in the low layer of the atmosphere strongly heat the nearsurface layer to form a well-mixed PBL, which is recognized as the stove effect (Ma et al., 2020). For another, the absorbing aerosols in the upper layer of the atmosphere strengthen the PBL inversion intensity and weaken its development, which is 70 noticed as the dome effect .
In addition, the large-scale synoptic patterns disturb the sensitivity of the aerosol-PBL feedback mechanism. PBLH is not only influenced by the local-scale surface properties (e.g., surface cover types, sensible heat fluxes, latent heat fluxes, etc.), but also regulated by large-scale synoptic forcing (Hu et al., 2016;Miao et al., 2020). Synoptic conditions can alter the PBL thermodynamic stability through cold or warm advection, thus affecting the aerosol vertical dispersion (Zhang, J. et al., 75 2012;Ye et al., 2016;Miao et al., 2019). Furthermore, due to the prevailing south winds, aerosols emitted in the southern regions can be carried to NCP via cross-region transport, and the increased aerosols eventually modulate the PBL thermodynamic stability (Du et al., 2020;Ma et al., 2020). On the inter-annual scale, previous studies reported that the wintertime air quality is closely connected with the East Asian Winter Monsoon (EAWM), and the strengthening (weakening) of EAWM is typically associated with changes in representative circulation patterns, which can improve 80 (worsen) regional air quality (Niu et al., 2010;Zhang, Y. et al., 2016).
According to the aforementioned efforts, most of them primarily concentrated on the aerosol-PBL interactions neglecting the simultaneous influences of the synoptic condition, the aerosol type and vertical distribution. With numerical simulations and observational constraints, the roles of synoptic pattern, PBLH, aerosol type and vertical distribution in aerosol-PBL interactions warrant further investigation. To address this issue, this study attempts to analyze the effects of 85 synoptic forcing and aerosol types on the PBL thermodynamic stability, aerosol vertical distribution, and aerosol-PBL interaction by using in-situ aircraft measurements, surface observations, and on-line coupled model Weather Research and Forecasting with Chemistry (WRF-Chem) simulations. On the 3 rd and 4 th of January, 2020, two days with distinct synoptic conditions, PBL stability, and aerosol vertical distributions over NCP are analyzed. Five parallel numerical experiments are conducted to investigate the AREs of scattering and absorbing aerosols under different aerosol vertical distributions. 90 Moreover, the long-term inter-annual variation in PBL stability from 1980 to 2020 over the NCP region is examined, and its driving factors are figured out. The conclusions appear to be beneficial in determining which pollutants to target in different weather conditions and achieving more precise controls of air pollution.
The remainder of this paper is organized as follows. Section 2 describes the data and methods used in this study.
Section 3 presents and discusses the aerosol-PBL interactions, as well as the factors of synoptic forcing and aerosol types. 95 The conclusions and summary are provided in section 4.  Fig. 1e) were carried out to measure the vertical profiles on January 3 and 4, 2020. The terrain height and three-dimensional flight routes are shown in Fig. 1 during the study periods. The flight routes on the two days were almost identical, which is conducive to the comparative analysis. The 105 KingAir-350 took off and climbed to about 3000 m over the Shahe airport (40.1° N,116.3° E), which is located in Beijing City. Afterward, the aircraft proceeded towards the southwest for about 150 km, arriving in Baoding City and performing vertical measurements ranging from around 650 m to 3000 m.  In this study, the aerosol number concentrations and size distributions were measured by a passive cavity aerosol 115 spectrometer probe (PCASP) mounted on the wingtip of the KingAir-350 aircraft. PCASP detected aerosol size distributions ranging from 0.1 to 3.0 μm at 1 Hz. Given the substantial noise of the first channel, we considered all the spectral bands except the first one, i.e. 0.11 to 3.0 μm, in this study (Twohy et al., 2005). In addition, an integrated meteorological measurement system (AIMMS) was employed to measure the aircraft location and ambient meteorological conditions (e.g. temperature, wind, etc.). The detailed information of aircraft observation is summarized in Table 1. 120

Advantech Research
Inc.

Ground-based observations
Surface meteorological observations at three ground-based automatic weather stations (AWS) in Beijing City were collected for model meteorological validation, which includes the hourly datasets of surface pressure, 2 m temperature, 10 m wind speed, etc. In addition, air quality observations, including hourly ground-based PM2.5 concentration measurements at 125 six air quality monitoring sites (AQMS) in Baoding City, were used in this study to validate the model performance. The observational air pollutant data was collected from the China National Environmental Monitoring Center (http://www.cnemc.cn/). Ground-based observations on January 3 and 4, 2020 were applied in this study, and the locations of AWS and AQMS are shown in Fig. 1b.

Reanalysis data 130
The data of isobaric temperature taken from the European Centre for Medium-Range Weather Forecasts (ECMWF) ERA5 reanalysis dataset was used to summarize the long-term thermal structure of PBL in North China Plain. The temperature data of seven isobaric surfaces between 1000 hPa and 850 hPa were considered to be within the PBL. The sea level pressure from the ERA5 reanalysis product was used to calculate and describe the intensity of the EAWM and Siberian High (SH), which are demonstrated to be the aspects that vigorously govern the PBL thermal structure over the NCP in this 135 https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License. study. Here, the monthly average data during the wintertime (December, January, February) from 1980 to 2020 with a grid of 0.25°×0.25° was adopted.

WRF-Chem model
In this study, we simulated the aerosol concentrations and meteorological conditions by using the WRF-Chem model version 3.9, which is an online-coupled three-dimensional Eulerian chemical transport model considering complex physical 140 and chemical processes (Grell et al., 2005). The three nesting domains with horizontal resolutions of 27, 9 and 3 km, respectively, are shown in Fig. 1a. There were 29 vertical layers extended from the ground to the top pressure of 50 hPa in the model, with more than 15 layers located below 3 km to fully describe the vertical structure of the PBL. The simulations were carried out from 0000 UTC on January 2 to 1800 UTC on January 4, 2020, and the first 16 hours were set as the model spin-up time. 145 The 6-hour National Centers for Environmental Prediction (NCEP) global final (FNL) reanalysis fields with a spatial resolution of 1°×1° were input for the model initial and lateral boundary meteorological conditions. Anthropogenic emissions were adopted from the Multi-resolution Emission Inventory for China (MEIC) developed by Tsinghua University (http://meicmodel.org), whose emission sources were classified into five sectors: industrial process, power plants, residential combustion, on-road mobile sources, and agricultural activities (Li, M., Liu, H. et al., 2017b;Zheng et al., 2018). The Model 150 of Emissions of Gases and Aerosols from Nature (MEGAN) was used online in the simulation to calculate the biogenic emissions (Guenther et al., 2006).
The model physical configurations we used in the simulation included the Morrison double-moment microphysics scheme (Morrison et al., 2005), the RRTMG radiation scheme (Iacono et al., 2008), the Yonsei University (YSU) boundary layer scheme (Noh et al., 2003) and the Noah land surface scheme (Ek et al., 2003). Model for Simulating Aerosol 155 Interactions and Chemistry (MOSAIC) was applied as the sectional aerosol scheme, with aerosols specified as 8 size sections (bins) ranging from 39 nm to 10 μm (Zaveri et al., 2008).
To investigate the interactions between aerosol and PBL thermodynamic stability over the NCP, five parallel numerical experiments were conducted as presented in Table 2. In the control experiment (denoted as the EXP_Ctrl), the anthropogenic emissions were performed at the normal level, and the aerosol optical properties were calculated at each time step and then 160 coupled with the short-and long-wave radiative transfer model. In the EXP_WoF, the emissions were the same as the EXP_Ctrl, except that ARE was not included. In the other three experiments, ARE was taken into consideration, but the emissions varied. Specifically, the BC emission was subtracted in the EXP_WFexBC, the BC emission was 20 times higher in the EXP_WF20BC, and the total aerosol emissions were strengthened to 20 times in the EXP_WF20Aer.

Experiment Description
EXP_Ctrl Control experiment, with the aerosol radiative feedback.

EXP_WoF
Without any aerosol radiative feedback.

EXP_WFexBC
With the aerosol feedback excluding the contribution of BC.

EXP_WF20BC
With the aerosol feedback including 20 times the BC emission.

EXP_WF20Aer
With the aerosol feedback including 20 times the total aerosol emission.

Statistical methods for model validation 170
Correlation coefficient (R), mean bias (MB), normalized mean bias (NMB), mean error (ME) and root mean square error (RMSE) were applied to assess the WRF-Chem model veracity in simulating meteorological parameters and air pollutants against the ground-based observations with the following equations (Itahashi et al., 2015;Granella et al., 2021): where i S and i O are the simulated and observed parameters, respectively, n is the total number of the values used for validation, and S and O are the averages of the simulation and observation, respectively. 180

Validations of meteorological parameters and aerosol concentrations
Considering that the weather conditions play important roles in pollutant transport, aerosol dispersion, and chemical reactions in the atmosphere, it is crucial to validate the simulated meteorological parameters with the observations. Fig. 2 shows the validation of meteorological parameters between the modeled and observed results, and model evaluations suggest 185 https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License.
that the WRF-Chem model is able to simulate the weather characteristics in the NCP region. The simulated temporal series and magnitudes of temperature, pressure and wind speed generally agree well with the ground-based and vertical aircraft observations. The accuracy of simulation during the daytime is superior to that at night, which is instrumental in the intention of discussing the PBL thermodynamic stability and ARE in the daytime.
The statistical evaluations of model performance during the daytime from 1000LST to 1800LST on January 3 and 4, 190 2020 are concluded in Table 3. Underestimation of temperature exists in both the surface level and vertical profile, though the R between the ground-based observation and EXP_Ctrl is 0.87 that indicates the model well describes the temperature variation. The surface pressure is well simulated with an MB of 0.81 and a small NMB of 0.08%. Slight overestimation of wind speed (0.40 m s -1 ) in the simulation might be due to the unresolved topography in the WRF model, which is also demonstrated in the previous studies (Jiménez et al., 2013;Li, M., Wang, T. et al., 2021a). 195    In addition, the EXP_Ctrl modeled surface PM2.5 mass concentrations in Baoding City also compare well with the 205 ground-based measurements, especially during the daytime (Fig. 2). The statistical validations indicate an R of 0.79, a total MB of -4.91 μg m -3 , and an NMB of 4.69% during the daytime (Table 3). Both simulation and observation display a high level of air pollution on 3 January, and good air quality on 4 January 4. Therefore, in this study, we consider that the WRF-Chem simulation is in line with the observation and can capture the general distribution and variation in air pollutants.

Aerosol vertical characteristics and PBL thermal structures
Despite the large-scale weather condition, the local-scale PBL thermal structure is also reported as a critical aspect that substantially affects the surface aerosol loading (Zhang, Q., Ma, X. et al., 2009). The PBL thermodynamic stability dictates 215 the PBLH, thereby dominating the vertical dissipation of surface pollutants to some extent. Therefore, the detailed synoptic condition, PBL thermal structure, and their influences on the aerosol vertical characteristics are discussed in the following subsections.

Vertical thermal structure and synoptic condition
The temporal evolution of the EXP_Ctrl simulated temperature profile is presented in Fig. 4, which reflects the 220 variation in the thermal structure and stability during the study period. The results indicate that the thermal structures differ significantly between the two days. It is noticed that the temperature on 3 January is warmer than that on 4 January, suggesting different wintertime synoptic characteristics. Furthermore, the average lapse rate below 1.5 km is lower on 3 January, with lapse rate values of 3.62 ℃ km -1 and 5.38 ℃ km -1 during the aircraft observation periods on 3 and 4 January,  The thermal structure of the atmosphere is closely related to and highly sensitive to the synoptic condition over the NCP, especially during the wintertime Miao et al., 2019). Figs. 5a and 5b display the surface synoptic patterns 230 during the aircraft observation periods on 3 and 4 January, respectively. In the adjacent regions of the Baoding City, the wind fields are practically opposite in these two periods, and the contrast is also intuitively shown in Fig. 5c with the vertical profile difference in the meridional wind speed. The prevailing surface wind direction on 3 January is southwest, while on 4 January is northeast. The south wind at the low layer warms the atmosphere below 1.5 km on 3 January, whereas the north wind leads to a cooling effect on 4 January, thus there appears a large temperature difference manifested in Fig. 4. 235 Furthermore, the south wind warms the atmosphere, particularly the aloft layer at about 1 km altitude, leading to a stable thermal stratification below 1 km and a low PBLH on 3 January.

Observed and simulated aerosol vertical distributions
As demonstrated in Fig. 6, contrasting aerosol vertical distributions are observed with aircraft during the daytimes of January 3 and 4, 2020. The accumulation mode particles (0.1-1 μm) dominate the aerosol number concentrations during the 245 study period, which is primarily emitted by anthropogenic pollutants, such as BC aerosol loadings . On 3 January, larger amounts of aerosol particles are constrained to the low layer of the atmosphere (below 650 m), while the pollutants decrease sharply above 650 m and gradually tend to a stable low value. On 4 January, the aerosol number concentrations below 650 m are lower than that on 3 January, but higher when the height is above 1000 m. This contrast is attributed to the differences in atmospheric stability (Su et al., 2020). The great difference in lapse rate shown in Fig. 4 leads  250 to the disparity in atmospheric stability and aerosol dispersion ability. On 4 January, the more unstable thermal stratification promotes the dispersion of the near-surface pollutants, and meanwhile, exacerbating the pollution on the elevated layer. The temporal evolution of the EXP_Ctrl simulated vertical profiles of total aerosol number concentrations and BC mass concentrations are shown in Fig. 7. Both profiles have similar features, indicating that anthropogenic fine particles like BC, are the dominant pollutants, which is consistent with the particle size observations exhibited in Figs. 6b and 6c. The diurnal variation in the aerosol profile is strongly driven by the evolution of PBLH, which is well examined in the previous studies 260 (Qu et al., 2017;Huang et al., 2018). The pollutants are well constrained to the near-surface by the stable PBL at nighttime, while the vertical distribution of aerosol particles is more homogenous by turbulent mixing in the afternoon. In addition, on 4 January, the vertical mixing of aerosols is more homogenous, transporting aerosols to a higher layer, which implies a more unstable thermal stratification and a higher PBLH.  The contrasting aerosol vertical distributions between the two days are interpreted from two aspects. From the perspective of aerosol regional transport, the air quality is remarkably altered by the synoptic condition (Zhang, Y. et al., 270 2016;An et al., 2019). The accumulation of pollutants over NCP is contributed by regional transport from the polluted southwest region via the northward movement of the air mass, whereas the cessation of pollution over NCP is caused by prevailing northern air mass and then the pollution hotspot moves to the south (Luo et al., 2021). Therefore, the observed and simulated high surface pollution loading on 3 January is externally affected by the unfavorable synoptic condition, which is consistent with the south wind demonstrated in Fig. 5. From the perspective of aerosol vertical distribution, the surface air 275 pollution loading is efficiently driven by the thermal structure and stability of the PBL (Miao et al., 2020;Li, Q. et al., 2021b). On 3 January, the prevailing surface south wind leads to a stable stratification, restraining aerosol vertical dispersion https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License. and finally exacerbating the surface air pollution. Besides, the north wind condition on 4 January has the opposite effect, which reduces the surface aerosol loading.
In conclusion, the stable stratification related by the prevailing surface south wind limits the PBL development and 280 exacerbates the surface air pollution by restraining the aerosols to the low layer. Conversely, the unstable stratification generated by surface north wind promotes the PBL turbulent mixing and lightens the surface aerosol concentration via vertical diffusion.

ARE on PBL thermal structure 285
PBL structure modulates the aerosol vertical distribution by turbulent mixing, and inversely, the suspended aerosols may also modify the PBL thermal structure by ARE to some extent. Fig. 8 presents the temperature variation affected by ARE, and it is noticeable that light-absorbing aerosols heat the atmosphere while light-scattering aerosols contribute to a cooling effect. On 3 January, aerosols are constrained to the low layer during the daytime, and the overall ARE results in a warming effect below 1 km. The conspicuous heating suggests that the leading role of absorptive BC aerosols on the 290 variation in PBL thermal structure when coexisting with other types of scattering aerosols on 3 January. On 4 January, due to the strong turbulence mixing, aerosols are carried to the aloft layer. The overall ARE results in a cooling effect below 0.6 km and a warming effect between 0.6 and 1.5 km. When the aerosols are transported to a high altitude, the absorbing and scattering particles inhibit the incident solar radiation from reaching the low layer and exhibit a cooling effect, while the absorbing aerosols heat the upper layer. The contrasting aerosol vertical distribution caused by the varied PBL thermal 295 structure between the two daytimes results in different AREs, therefore, it is quite crucial to examine the aerosol vertical structure when evaluating the ARE.  The simulated ARE modify the thermal structure with a temperature variation of less than 0.5 ℃ shown in Fig. 8, to amplify the signal of ARE when discussing its interactions with PBL, we further performed two experiments, one with strong overall aerosol loading (EXP_WF20Aer) and the other only with strong BC aerosol loading (EXP_WF20BC). The simulated thermal structures are presented in Fig. 9. The amplified simulations reveal a temperature variation of less than 305 2.5 ℃, which is consistent with the altered magnitude of aerosol emissions. ARE by only absorptive BC aerosols heat the low layer more remarkable when compared with the role of overall aerosols, suggesting that the coexists of absorbing and scattering aerosols partially conceal the signal of ARE. This result is owing to the opposite optical roles of the two aerosol types in the low layer, as well as the strong inhibition of the incident solar radiation by the aerosols in the upper layer.

ARE on the PBLH, lapse rate, and surface pollution 315
Since ARE can modulate the PBL thermal structure, the related impacts caused by the thermal variation are discussed here. Fig. 10 shows the time series of the simulated lapse rate, PBLH, and surface carbon monoxide (CO) concentration under five parallel experiments.
Generally, on 3 January, the PBLH is lower and the surface air pollution is heavier than those on 4 January. Influenced by the southwest winds, the warmer air mass can be guided to NCP (Fig. 5a), enhancing the thermal stability and restraining 320 https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License. the growth of the PBL (Fig. 10b). On the one hand, from the local-scale view, the lower PBLH constrains the air pollution to the near-surface; on the other hand, from the regional transport view, the southwest prevailing winds can transport the pollutants emitted from upstream plain regions to NCP, which further worsen the surface air quality (Fig. 10c). However, the relatively higher PBLH and the northeast winds on 4 January are able to reduce the surface air pollution and inversely carry the pollutants to the higher altitude with the strong turbulent mixing. The disparate synoptic conditions regulate the vertical 325 aerosol structures and further affect the ARE and other related variables.
The distinct aerosol types contribute to varied AREs when different aerosol vertical distributions are dominant. The lapse rate and PBLH are closely associated with the vertical thermal structure. A warmer low layer and a colder upper layer induce a higher lapse rate and PBLH, while a colder low layer and a warmer upper layer decrease the lapse rate and PBLH.
On 3 January, the pollutants are restricted to the low layer, the absorbing and scattering aerosols play the opposite roles. The 330 absorbing aerosols heat the low layer and invigorate the turbulent mixing (stove effect), leading to a higher lapse rate and PBLH. On the contrary, the scattering aerosols cool the low layer and inhibit the growth of the PBL (umbrella effect), causing a lower lapse rate and PBLH (Figs. 10a and 10b). On 4 January, the aerosol vertical profile differs from that on 3 January, with a huge number of pollutants being carried to the high altitude by the unstable stratification. The absorbing and scattering aerosols both prevent the solar radiation to the low layer (umbrella effect), while the absorbing aerosols also heat 335 the upper layer (dome effect), generating a lower lapse rate and PBLH (Figs. 10a and 10b).
The variations in lapse rate and PBLH further alter the surface air pollution by the thermodynamic and dynamic effect, and high PBLH ordinarily lightens surface air pollution. Due to the CO is a relatively stable pollutant that is hardly interactions with other pollutants and the main sources of CO are primary emissions from anthropogenic pollution and open biomass burning (Chi et al., 2013), we use CO concentration to evaluate the variation in surface air pollution caused by the 340 ARE in this study. The surface CO concentration in the EXP_WF20Aer is divided by a factor of 20 to match other experiments since the emissions are magnified 20 times. The results presented in Fig. 10c indicate that the surface air pollution variation is closely related to the PBLH variation and driven by the aerosol vertical dispersion. When the stable stratification confines the absorbing aerosols to the low layer, the stove effect destroys the stability, enhancing the thermal convection and decreasing the surface air pollution. The absorbing aerosols in the upper layer or the scattering aerosols in 345 both the low and upper layers stabilize the atmosphere, deteriorating the surface air quality.

ARE under different PBL thermodynamic stabilities
It is noticeable that different aerosol vertical distributions between the two days contribute to distinct AREs due to the PBL thermal structure and stability differences. In particular, the absorptive BC aerosols have both stove and dome effects, which are sensitive to the PBL thermal structure. Here, we further analyze the spatial distribution of atmospheric stability during the daytime of the aircraft observation periods on 3 January. The atmospheric stability is divided into two categories 355 following the method developed by Wolyn and McKee (1989): a stable stratification known as the deep stable layer (DSL) and an unstable stratification known as the non-deep stable layer (non-DSL). The DSL is identified supposing that 65% of the lowest 1.5 km has a lapse rate of 2.5 ℃ km -1 or less, and the DSL definition examines a layer much deeper than a typical nocturnal inversion (Wolyn and McKee, 1989). The spatial distributions of the DSL and non-DSL are displayed in Fig. 11. The impacts of PBL thermodynamic stability on the BC ARE are compared between DSL and non-DSL conditions demonstrated in Fig. 12. Under the non-DSL condition, the BC dome effect is dominant, the BC aerosols in the upper layer absorb the solar radiation, stabilize the PBL and decrease the PBLH (Fig. 12a), and further increase the surface air pollution 365 (Fig. 12c). Under the DSL condition, the BC stove effect plays the dominant role, the strong BC concentrations in the low layer contribute to the development of the PBL and increase the PBLH (Fig. 12b), and further decrease the surface air pollution (Fig. 12d).

Figure 12: Scatter plots of the correlations between EXP_Ctrl and EXP_WF20BC for (a-b) PBLH and (c-d) CO concentration
370 under non-DSL condition and DSL condition, respectively. The PBLHs lower than 500 m are discarded to avoid the incorrect simulation.

Long-term variation in PBL thermodynamic stability over NCP
The linkage between synoptic pattern and PBL thermal structure is analyzed from a climatological view, which is critical for the examination of aerosol vertical distribution and ARE. The potential impact of inter-annual climate variability 375 on the PBL thermodynamic stability is explored in this study over the NCP during the wintertime of 1980-2013. Following the method developed by Wu and Wang (2002), the index of SH and the index of EAWM are used to analyze the possible factors that influence the PBL thermal structure over the NCP. The index of SH and EAWM are the two variables describing the meridional wind intensity over the coast of East Asia, which are crucial to the PBL thermal structure. The PBL lapse rate is calculated between 1000 hPa and 850 hPa. 380 Fig. 13 suggests that the relationships between SH index, EAWM index, and PBL lapse rate are highly correlated.
Weak EAWM and SH with the prevailing south wind stabilize the PBL and constrain the aerosols to the near-surface, while https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License. the strong EAWM and SH play the opposite role. The result indicates that the inter-annual variability of the EAWM and SH can influence surface air quality, aerosol vertical distribution, and ARE through its reflection on the intensity of PBL thermodynamic stability over the NCP region. 385

Conclusions and Summary 390
In this study, the complex relationships between large-scale synoptic patterns, local PBL thermal structures, aerosol vertical distributions, and AREs of different aerosol types are systematically investigated by the combined aircraft observations, surface measurements, reanalysis data, and WRF-Chem simulations. The in-situ aircraft observations were carried out over the Beijing and Baoding Cities during the daytimes on January 3 and 4, 2020, and the numerical simulations were performed over NCP from 3 to 4 January. By and large, the meteorological variables and air pollutants are well 395 simulated when validated by the observations.
Observations show that the surface air pollution over the Baoding City on 3 January is heavier than that on 4 January, and the aerosols are constrained to the low layer on 3 January while the aerosols mix more homogeneous vertically on 4 January. The numerical simulations suggest that the synoptic pattern over the Baoding City differs between the two days, and the prevailing wind direction is opposite with a southwest wind on 3 January and a northeast wind on 4 January. The 400 synoptic condition may affect the PBL thermal structure, thus affecting the vertical dispersion of aerosols. The south winds flow the warm and polluted air to the NCP and stabilize the PBL, hence exacerbating the surface air pollution. The north winds flow the cold and clean air to the NCP, establishing an unstable stratification and transporting the pollutants to the upper layer.
The suspended aerosols may also modify the PBL thermodynamic stability through ARE to some extent, and the 405 synoptic condition can modulate the sensitivity of PBL-aerosol feedback by influencing the PBL thermal structure and aerosol vertical distribution. The interactions among the synoptic pattern, PBL thermal structure, and aerosol vertical distribution are demonstrated in Fig. 14.

430
In addition, the factors affecting the wintertime PBL stability over the NCP area are further analyzed from the climatological perspective, which is mainly controlled by the intensity of EAWM and SH. A strong winter monsoon creates an unstable PBL stratification, whereas a weak winter monsoon stabilizes the PBL. Combined with the umbrella, stove, and dome effects of scattering and absorbing aerosols under different synoptic conditions, this finding aids us in determining which pollutants to target in different monsoon years and achieving more precise air pollution control. 435 https://doi.org/10.5194/acp-2021-769 Preprint. Discussion started: 12 October 2021 c Author(s) 2021. CC BY 4.0 License.
It should be noted that only two days of pollution episodes are examined in this study. While the model simulations produce evidence that is compatible with the observations and substantiate the observational analysis, additional measurements should be carried out in the near future.

Code/Data availability
The surface observational air pollutant data are collected from http://www.cnemc.cn. The meteorological data taken from the 440 ERA5 reanalysis dataset are obtained from https://www.ecmwf.int/en/forecasts/datasets/reanalysis-datasets/era5. The NCEP-FNL global reanalysis dataset used in the WRF-Chem model is available from https://rda.ucar.edu/datasets/ds083.2. The MEIC anthropogenic emission inventories are available from www.meicmodel.org, and for more information, please contact Q. Zhang (qiangzhang@tsinghua.edu.cn). The WRF-Chem model version 3.9 is available from https://www2.mmm.ucar.edu/wrf/users/download/get_source.html. The aircraft data and surface meteorological data are 445 available from the corresponding author. The codes that support the findings of this study are available from the corresponding author upon reasonable request.