Vertical Distribution of PM 2.5 Transport Flux in Summer and Autumn in Beijing

Temporal and spatial distribution characteristics of PM 2.5 transport fluxes between Beijing and its adjacent cities in July and October 2016 were estimated quantitatively using the cross-boundary transport flux method based on the WRF-CMAQ model in this study. Vertical distribution of PM 2.5 net fluxes was revealed, and hourly evolution characteristics of PM 2.5 transport fluxes in typical heavy pollution processes were illustrated. We verified significant seasonal differences in PM 2.5 transport fluxes between Beijing and its adjacent cities during the study period. PM 2.5 total inflow and outflow fluxes in autumn were 3187 t d – 1 and – 2721 t d – 1 , which were significantly higher than that of 2134 t d – 1 and – 2172 t d – 1 in summer. In autumn, Beijing received more PM 2.5 from adjacent cities, while the results in summer were opposite. Maximum net fluxes appeared at 600 – 800 m and 1000 – 1260 m above the ground in summer and autumn, respectively. The vertical distribution characteristics of pollution days and clean days were consistent, both of which show a net inflow from Baoding and Langfang to Beijing, while Beijing was in a state of net outflow to Chengde. During the whole heavy pollution process, the evolution characteristics of PM 2.5 flux at low and high altitudes were consistent, and the intensity of the latter was 2.15 – 5.30 times of the former. Meanwhile, local discharging was more likely to cause extreme PM 2.5 heavy pollution, and a better peak clipping effect can be achieved when emission reduction measures are initiated 1 to 2 days before heavy pollution. The results can provide scientific support to put forward effective joint control measures and obtain insights into the evolutionary mechanism of haze episodes in Beijing.


INTRODUCTION
Currently, China's air pollution presents regional compound pollution characterized by PM2.5 and O3 pollution. Its essence is the coupling of multiple atmospheric environmental problems caused by the emission of a variety of air pollutants (Hua et al., 2021;Yan et al., 2020). Affected by the atmospheric circulation, multiple pollutants were transported across the boundary, resulting in a series of atmospheric pollution , such as the rise of PM2.5 and O3 concentration levels in contiguous areas, haze, photochemical smog, etc. PM2.5 has small particle sizes, which are rich in a large number of toxic and harmful chemical components (Zhang et al., 2022). In addition, PM2.5 can stay in the atmosphere for a long time  and be transported over a long distance . Therefore, it has a significant impact on air quality (Vega et al., 2021), atmospheric visibility , human health (Qi et al., 2016), climate changes (Zhong et al., 2019), and other aspects. China's air quality has improved significantly since 2013 when the Chinese government implemented the "Air Pollution Prevention and Control Action Plan". However, regional heavy pollution incidents dominated by PM2.5 still occurred frequently, especially in the Beijing-Tianjin-Hebei region (Li and Han, 2016), Yangtze River Delta , Pearl River Delta , and Sichuan-Chongqing (Li et al., 2017a), which have attracted widespread attention from domestic and foreign scholars. Therefore, it is imperative to quantitatively identify the quantitative characteristics of transboundary transport of PM2.5 pollution and to elucidate the contribution of local emissions and regional transport.
Technical means of studying PM2.5 cross-boundary transport at home and abroad mainly included mathematical statistics, stereo observation, and modeling simulation. Ge et al. (2018) developed a mathematical-statistical approach to calculate the relative contributions of minimal local emission (MinLEC) and the maximum regional transport (MaxRTC), which is based on the measurement results of several air pollutants concentrations near the surface of the North China Plain. They estimated impact of local emissions and regional transport on air quality, and illustrated that the average annual MinLEC (MaxRTC) contributions of NO2, SO2, PM2.5, and CO in the North China Plain were 61.7% (30.7%), 46.6% (48%), 52.1% (40.2%) and 35.8% (45.5%), respectively. Lv et al. (2017) have studied heavy haze on the sixth ring road in Beijing based on a vehicle-based mobile lidar and found that the average inflow flux from the southwest on the day of heavy pollution was about 308 µg m -2 s -1 . This showed a significant transport path from the southwest to the northeast, demonstrating a greater impact on the formation of this heavy pollution along the southwest transport path. In addition, the air quality model has become one of the most widely used and most important methods for studying the laws of cross-boundary transport of pollution, which can reproduce the physical and chemical processes of atmospheric pollutants in the atmosphere, such as the transmission, transformation, diffusion, and sedimentation. At present, scholars at home and abroad have quantitatively identified the contribution of atmospheric PM2.5 pollution transport near the surface in multiple regions, indicated the importance of regional transport to the formation of heavy PM2.5 pollution in multiple cities/regions, based on the source apportionment techniques in numerical models, such as the Particulate Matter Source Apportionment Technology (PAST) in CAMx (Comprehensive Air Quality Model with Extension) model , the Integrated Source Apportionment Model (ISAM) in CMAQ (Community Multiscale Air Quality) Model (Dong et al., 2020), the on-line tracer-tagged technology in NAQPMS (a threedimensional regional Eulerian chemical transport model developed by the Institute of Atmospheric Physics) model , the source sensitivity (brute force zero-out) approach (Baker et al., 2016), and the process analysis (Zhang et al., 2019a). However, the above-mentioned studies mostly focused on the contribution of near-surface PM2.5 local emissions and regional transport and lacked the quantitative characteristics of the transport of PM2.5 pollution at different ground clearances in the vertical direction.
Transport flux method has been widely used to study the characteristics of air pollutants crossboundary transport and can quantitatively evaluate the transportation of pollutants between neighboring cities and regions by calculating the transport flux of pollutants at different vertical heights, based on wind field and pollutant concentration outputted by meteorological and air quality models (Yao et al., 2021). At present, most studies focused on the vertical distribution characteristics of PM2.5 transport flux in typical seasons (Chang et al., 2018). However, there were relatively few studies focusing on typical heavy pollution processes, elucidating the hourly evolution laws and vertical distribution characteristics of PM2.5 transport flux, and the difference of PM2.5 transport flux between clean days and heavy pollution days need to be further investigated.
Based on the implementation of regional air pollution prevention and control policies by the Chinese government, the PM2.5 pollutions had been greatly improved during 2016-2022. However, heavy pollution has not been completely eliminated, especially in autumn when the wind is high and the air is crisp, heavy pollution events still occur occasionally. Therefore, we quantitatively identified the vertical distribution characteristics of PM2.5 transport flux between Beijing and neighboring cities during July and October 2016, based on the WRF-CMAQ model, analyzed the difference in flux intensity under different seasons and pollution levels, and clarified the hourly evolution characteristics of PM2.5 transport flux in typical heavy pollution processes. Flux intensity calculations can provide scientific support to formulate scientific and effective air pollution prevention and control strategies and PM2.5 heavy pollution response plans in China.

Model Configuration
In this study, the Weather Research and Forecasting model (WRF, version 3.7.1) (Yang et al., 2018) and the Comprehensive Air Quality Model with Extension model (CMAQ, version 5.0.1) (Zhang et al., 2019b) were employed. The study period was from July 5 to 30 and October 5 to 30, 2016, with the simulation time interval of 1 h. Among them, the parameter setting of WRF model can be found in Table 1. It is noteworthy that the cumulus cloud parameterization scheme is adopted for the Domain 1 and 2, but not for the Domain 3. δ coordinate represents the height parameter set in the vertical direction by the WRF and CMAQ model. Both models set 28 air pressure layers in the vertical direction of the simulation area, with the corresponding δ coordinates of 1. 000, 0.994, 0.988, 0.981, 0.969, 0.956, 0.944, 0.926, 0.902, 0.881, 0.852, 0.828, 0.796, 0.754, 0.704, 0.648, 0.589, 0.546, 0.495, 0.445, 0.387, 0.287, 0.187, 0.136, 0.091, 0.061, 0.020 and 0.000. The corresponding height of each δ coordinate have been summarized in Table 2.
The horizontal simulation grids of the two models used the Lambert projection coordinate system, with three nested domains set up, and the grid resolution remained the same. The boundary layer of the outer grid adopted the default setting of the CMAQ model, and the boundary layer of the inner grid was determined by outer layer simulation. The central point for the CMAQ model was set at the coordinate (39.02°N, 116.47°E). Domain 1 covered most of eastern and central China, with a horizontal resolution of 36 km; Domain 2 covered the Beijing-Tianjin-Hebei (BTH) region, Henan, Shanxi, Shandong, parts of the Inner Mongolia Autonomous Region, Jiangsu, Hubei, Shaanxi, Anhui, Liaoning, Zhejiang, and Shanghai, with a horizontal resolution of 12 km; Domain 3 covered Beijing and neighboring cities, with a horizontal resolution of 4 km, as shown in Fig. 1. Since the initial conditions used in the inner and outer layer simulations were all default values, the actual simulation period starts on the 1 st of the month. The first four days were to ensure that the model simulation was initialized and was not used for subsequent analysis.

Data Collection
In this study, the hourly meteorological parameters at the surface, including temperature (T), wind speed (WS), wind direction, and relative humidity (RH), were obtained from the Chinese MICAPS (Meteorological Information Combine Analysis and Process System) (Wang et al., 2017b). The hourly data of surface PM2.5 mass concentration were acquired from the website (https://www.aqistudy.cn/historydata) provided by the China Air Quality Online Monitoring and Analysis Platform, which were used to verify the accuracy of WRF and CMAQ model simulations. The initial weather field and boundary conditions adapted global troposphere analysis data based on the NCEP (National Centers for Environmental Prediction), with a resolution of 1° × 1° and a temporal resolution of 6 h. We employed a high-resolution county-level pollution emission inventory for the BTH region, which was developed by the Beijing University of Technology (Zhou et al., 2017). Emissions data outside the BTH region were obtained from the multi-resolution emission inventory for China (http://www.meicmodel.org/). The pollutants in the emission inventory include SO2, CO, NOx, NH3, VOCs, PM2.5, PM10, and others.

PM2.5 Flux Calculation
The PM2.5 flux stands for the mass of PM2.5 that flows through a special vertical surface during a certain period (Zhang et al., 2019c), which was used to quantitatively analyze the characteristics of PM2.5 cross-boundary transport at different heights between the target region and the surroundings. We took Beijing as the target city, determined the boundary line and its corresponding vertical section according to China's administrative divisions, whose adjacent cities consisted Zhangjiakou, Baoding, Langfang (S), Langfang (N), Tianjin, and Chengde. The WRF and CMAQ model belongs to the Euler three-dimensional grid mode, the simulation area can be regarded as a three-dimensional box composed of multiple grids (Thongthammachart et al., 2021). Each grid has a fixed threedimensional position and generates the wind field and PM2.5 concentration in the corresponding grid. Therefore, the outflow data of the model can be used to calculate the transport flux of PM2.5 on a specific section, after the section to be calculated was discretized according to the meshing method in the CMAQ model. The PM2.5 inflow flux stands for the mass of PM2.5 that flows into the target city from surrounding cities during a certain period, with represents of "+". The PM2.5 outflow flux stands for the mass of PM2.5 that flows from the target city to surrounding cities during a certain period, with represents of "-". The PM2.5 net flux stands for the vector sum of the two fluxes. The calculation formula of PM2.5 flux is indicated as below (Eq. (1)): where Flux is the PM2.5 transport flux (t d -1 ); h is the top height, with set to 12 layers in this study; l is the boundary line between two adjacent cities; L is the boundary length (m); Hk is the altitude difference between layer k and k + 1 (m); c is the PM2.5 concentration (µg m -3 ); v is the wind vector (° and m s -1 ); and n is the normal vector corresponding to the vertical grid cell; 8.64 × 10 -8 is the unit conversion coefficient that converts µg s -1 into t d -1 .

Model Evaluation
To evaluate the performance of the WRF and CMAQ model, we selected correlation coefficient (COR), normalized mean bias (NMB), and normalized mean error (NME) as statistical parameters, which refer to U.S. Environmental Protection Agency (U.S. EPA)'s evaluation standards (U.S. EPA, 2007), calculated as follows (Eq. (2); Eq. (3); Eq. (4)): where Sim(i) stands for the simulations of meteorological factors or PM2.5 mass concentration, Obs(i) stands for the observations of meteorological factors or PM2.5 mass concentration, n stands for the number of samples.
The meteorological observation data of two first-level stations in Beijing and Shijiazhuang were selected to verify the accuracy of the WRF model, including temperature at 2 m (T2), relative humidity (RH), and wind speed at 10 m (WS10). The verification results are shown in Fig. 2 and Fig. 3. On the whole, the simulations of T2, RH, and WS10 were relatively consistent with the observations, with a significantly positive correlation. For T2, the simulations of the two cities in July and October were generally higher than observations, with CORs of 0.74-0.93, NMBs ranging from 0.05% to 1.35% and NMEs ranging from 0.35% to 1.35%, so as to explain the WRF model can better reproduce the temporal and spatial distribution characteristics of air temperature. In terms of RH, the simulations were underestimated, with NMBs of -43.09% to -18.99%, NMEs of 21.27%-43.09%, but the CORs were both above 0.54. For WS10, the CORs of Beijing (Shijiazhuang) in July and October were 0.41 (0.63) and 0.67 (0.58), respectively, which were lower than the CORs of T2 and the CORs of RH. Simulations were generally higher than observations, the simulation deviations of Shijiazhuang were relatively greater than that of Beijing. NMBs of Beijing (Shijiazhuang) during July and October were 14.33% (47.46%) and 26.91% (74.87%), respectively, and NMEs of Beijing (Shijiazhuang) during July and October were 30.86% (49.63%) and 42.63% (75.73%). Notably, the simulations of near-surface wind speed were generally greater than previous model studies , with reasons that may include differences in the local terrain, insufficient data resolution of the urban underlying surface, etc. (Han et al., 2018).
For PM2.5 mass concentration, the simulations and the observations of the two cities had better consistency, with CORs both above 0.67. And the simulations were generally lower than 6 of 15 Volume 22 | Issue 8 | 220143 the observations, with NMBs and NMEs varied between -35.74% to -5.44% and 31.98%-40.36%, respectively. There were varying degrees of underestimation of PM2.5 mass concentration during heavy pollution days, which was similar to the other studies and may be the main reason for the larger overall simulation error. This may be related to many factors, such as the uncertainty of the emission inventory, the imperfect secondary pollutants conversion mechanism in the model, the horizontal and vertical resolution, and the WRF simulation error (Gao et al., 2016;Liu et al., 2017).
In summary, WRF and CMAQ models could accurately reproduce the evolution characteristics of various meteorological elements and PM2.5 concentration in the simulation region, and the simulation error was within an acceptable range, which could be further used for the simulation study of PM2.5 transport flux.

Seasonal difference
There were significant differences in the total transport flux of PM2.5 between Beijing and adjacent cities in July and October 2016, as shown in Fig. 4. Overall, based on the analysis of PM2.5 and wind speed simulation error, the PM2.5 concentration in July was underestimated more than that in October, meanwhile the wind speed simulations in July were less overestimated than those in October. These simulation errors may have increased the difference in the intensity of the inflow and outflow flux between July and October. The PM2.5 average daily total inflow flux (2134 t d -1 ) was equivalent to the total outflow flux (-2172 t d -1 ), resulting in a small total net flux (-38 t d -1 ), which showed that the transport between Beijing and surrounding cities was basically in balance. The total inflow flux from neighboring cities to Beijing was in the order of   Baoding > Langfang (N) > Langfang (S) > Zhangjiakou > Tianjin > Chengde, the total outflow flux from Beijing to surrounding cities was in the order of Chengde > Zhangjiakou > Langfang (N) > Baoding > Langfang (S) > Tianjin. The PM2.5 total net flux in October was relatively greater (467 t d -1 ), which was about 12 times that in July. As for July, Beijing showed a net output to Chengde, while Baoding and Langfang showed a net input to Beijing. In October, both the total inflow flux   (3187 t d -1 ) and outflow flux (-2721 t d -1 ) were greater than the flux in July. As in July, the total inflow flux from Baoding to Beijing was the strongest in October, and the total outflow flux from Beijing to Chengde was the strongest in October.

Analysis of the difference between clean days and polluted days
In this study, the periods when the average daily concentration of PM2.5 was not higher than 75 µg m -3 were defined as clean days, the periods when the average daily concentration of PM2.5 was higher than 75 µg m -3 were defined as pollution days. On the whole, the simulation errors of PM2.5 concentration in clean days and wind speed in polluted days were relatively smaller. While PM2.5 concentration in polluted days were underestimated and wind speed in clean days were overstated. The above results may lead to the increase of inflow and outflow flux intensity difference between clean days and polluted days. As shown in Fig. 5, the total inflow flux from Langfang (N) and Tianjin to Beijing were relatively larger, with the flux of 994 and 771 t d -1 , respectively. The total outflow flux from Beijing to Chengde and Zhangjiakou was relatively larger, with the flux of -914 and -1080 t d -1 , respectively. This showed that the transport contribution  of high concentrations-PM2.5 from external sources in Beijing in summer was mainly affected by the southeast monsoon climate, and the transport pathway along the southeast-northwest direction was relatively significant. In addition, the contributors of the total inflow and outflow flux of clean days in July were different from the pollution days. The major contributors were Baoding and Chengde, respectively, with the total net flux of 414 and -737 t d -1 , respectively, which indicated that PM2.5 in Beijing was mainly transported along the southwest-northeast direction during clean days in July. Contrary to the above, Baoding and Langfang (S) became the main contributors of the total inflow flux during polluted days in October, with the flux of 1645 and 1310 t d -1 , respectively. Chengde became the main receiver of the total outflow flux, with the flux of -2478 t d -1 , while the main contributors of the total net flux were Baoding, Langfang (S), and Chengde. Thus, higher concentration may be affected by the southwest-northeast transport pathway in October in Beijing. Compared to the polluted days, the total inflow flux, outflow, and net flux during clean days were significantly lower, were 0.20-0.46 times the total flux of polluted days, and fluxes of the main contributing cities were relatively close.

Seasonal difference
The vertical distribution characteristics of PM2.5 net flux in Beijing in July and October 2016 are shown in Fig. 6. On the one hand, it showed that the PM2.5 net flux varied with the height above the ground. On the other hand, it also revealed that there were significant differences in the direction of net transport flux in different seasons. The total net flux of PM2.5 at each height above the ground in July was negative, except for those below 153 m, revealing Beijing became a "source" that released more outflow than the inflow from the surroundings. Contrary to July, the total net flux of PM2.5 at each height above the ground in October was positive, indicating Beijing acts as a role of "convergence" city that received the inflow from surroundings exceeded the outflow released by itself. In July, the PM2.5 total net flux below 153m was expressed as a net inflow and mainly from Langfang (N) and Tianjin, with flux between 16 t d -1 -20 t d -1 and 16 t d -1 -18 t d -1 , respectively. The PM2.5 total net flux that ranged from 153 m to 1261 m expressed as a net outflow appeared as a bimodal distribution and reached the maximum (-21 t d -1 ) at approximately 600-800 m. Langfang and Baoding were the main contributors of the net inflow flux. The net flux transported from Baoding to Beijing increased significantly with height (from -1 t d -1 to 68 t d -1 ), so as revealed the significant influence of pollutants transportation from Baoding to Beijing in the southwest transport path. In addition, Beijing had a continuous net outflow flux to Chengde and Zhangjiakou, especially Chengde, where received the net outflow flux reached -120 t d -1 .
In October, the PM2.5 total net flux below 459 m changed relatively stable (20-36 t d -1 ), then the flux increased significantly with height and reached the maximum (70 t d -1 ) at approximately  1000-1260 m. Similar to July, Langfang and Baoding were the main contributors of the net inflow flux, and Chengde was the main receiver of the net outflow flux. Generally, Beijing was greatly affected by regional transport from Baoding and Langfang, which had a greater influence on the outflow to Chengde.

Difference between clean and polluted days
The differences in vertical distribution characteristics of PM2.5 net flux during the clean days and polluted days are shown in Fig. 7. As shown in Figs. 7(c) and 7(f), the PM2.5 total net flux in July and October was shown in the vertical direction as follows: near-surface (+), mid-altitude (-), high altitude (+), and the sum of flux in polluted days was 1.71 times, 6.80 times, and 2.04 times that of the clean days, respectively. On the whole, the main contributor cities to the PM2.5 net inflow flux and net outflow flux in polluted days and clean days in the vertical direction were similar. The net inflow flux mainly came from Baoding and Langfang, with the peak in polluted days (clean days) of 90 t d -1 (58 t d -1 ) and 108 t d -1 (46 t d -1 ), respectively. Chengde was the main receiver of the net outflow flux, with the peak in polluted days (clean days) of -224 t d -1 (-102 t d -1 ).
For the polluted days in July, Langfang (N) and Tianjin were the main contributors of the PM2.5 net inflow flux, and the net outflow flux was mainly transported to Zhangjiakou, Chengde, Langfang (S), and Baoding. The total net flux appeared at a relatively high level between 252 m-817 m, ranging from -107 t d -1 to 125 t d -1 . But Beijing mainly received the net inflow from Baoding and Langfang and transported the net outflow to Chengde in clean days, with a relatively lower level of the total net flux (-3 to 25 t d -1 ). In October, there was no significant difference in the net flux between polluted days and clean days. Chengde was the main receiving city of the net outflow flux, Baoding and Langfang were the main contributing cities of the net inflow flux. It is noteworthy that the net flux in polluted days in October was significantly higher than that in clean days. Taking Baoding and Chengde, the cities with the largest net inflow and net outflow respectively, as examples: the maximum net flux in the vertical direction of Baoding and Chengde in polluted days was 225 t d -1 and -331 t d -1 respectively, the corresponding maximum net flux in clean days was 47 t d -1 and -85 t d -1 respectively, and the former was 4.82 times and 4.63 times the latter, respectively, so as to explain that the mutual transportation between Beijing and surrounding cities during the PM2.5 pollution period was relatively significant.

Hourly Evolution of PM2.5 Transport Fluxes during a Heavy Pollution Process
To further clarify the evolution characteristics of the PM2.5 transport flux at different stages of the heavy pollution process, we took a period in Beijing from 0:00 on the 17 th to 23:00 on the 21 st in October 2016 as an example and divided the heavy pollution into three stages: Stage I (0:00 on the 17 th -16:00 on the 18th), Stage II (17:00 on the 18 th -7:00 on the 20 th ), Stage III (8:00 on the 20 th -23:00 on the 21 st ). In the vertical direction, it is divided into low-altitude (below 400 m) and high-altitude (above 400 m). The time-wise changes of PM2.5 inflow, outflow, and net flux in the three stages of heavy pollution were analyzed, as shown in Fig. 8. On the whole, the temporalspatial evolution characteristics of PM2.5 transport flux in low altitude and high altitude were relatively consistent. The total inflow, outflow, and net flux of high altitude PM2.5 were 2.64, 2.15, and 5.13 times that of low altitude, respectively. Results of the low altitude were analyzed as follows: the average flux of the total inflow, outflow, and net of PM2.5 in the entire heavy pollution process was 1380 t d -1 , -1470 t d -1 , and -90 t d -1 , respectively. The total inflow and outflow fluxes were 1.32 and 1.61 times the monthly average flux of the low altitude in October, respectively, so as explained that the influence of the transportation during heavy pollution process has increased compared to the entire October. In the early stage of heavy pollution (Stage I), PM2.5 total inflow and outflow flux increased slowly at a low level, with the average flux of 790 t d -1 and -610 kt d -1 , respectively. While Beijing accepted the inflow flux mainly from Baoding, it also continues to release outflow flux to Chengde. The total net flux was basically positive, indicating that adjacent cities had a continuous net input to Beijing, which led to a slow increase in Beijing's PM2.5 concentration, from 8 µg m -3 at 0:00 on the 17 th to 109 µg m -3 at 12:00 on the 18 th , which reached the first peak of the heavy pollution process. Fig. 8. Hourly variations of PM2.5 inflow, outflow, and net flux at the lower-and higher-altitude between Beijing and its adjacent cities during the whole haze process.
Subsequently, PM2.5 concentration remained at a high level into the Stage II of heavy pollution, and the PM2.5 total inflow and outflow flux were significantly enhanced compared with Stage I, with the flux of 2720 t d -1 and -3050 t d -1 , respectively. However, the net flux was -330 t d -1 , so as indicated that the amount of PM2.5 outputted from Beijing in stage II was higher than the amount inputted to itself by adjacent cities. From 20:00 on the 18th to 15:00 on the 19 th , the PM2.5 total net flux basically remained positive, with an average flux of 840 t d -1 , of which the net inflow flux mainly came from Baoding and Langfang (S), with an average flux of 701 and 752 t d -1 , respectively. The net outflow flux was mainly transported to Chengde with a mean flux of -1221 t d -1 . As a result, PM2.5 concentration maintained a continuous increase from 99 µg m -3 to 183 µg m -3 , indicating that regional transportation played an important role in heavy pollution before it reached the extreme concentration. Subsequently, the PM2.5 total net flux showed a high negative value, especially the net flux reached -8731 t d -1 at 3:00 on the 20 th , where Baoding was the largest recipient of outflow flux with the highest net flux of -5700 t d -1 , and Tianjin and Langfang (N) were the largest contributors of inflow flux with the highest net flux of 2842 t d -1 . It is noteworthy that PM2.5 concentration didn't decrease due to the transportation (local to adjacent cities), but continued to rise, and reached the second (218 µg m -3 ) and third concentration peaks (225 µg m -3 ) at 19:00 on the 19 th and 2:00 on the 20 th respectively, thus indicating that this peak was mainly due to the contribution of local emissions in Beijing. In other words, local emissions have a greater potential to cause extreme PM2.5 heavy pollution to occur, while external transport tends to have a significant impact before heavy pollution peaks, which was relatively consistent with the findings of Li et al. (2017b) and Tang et al. (2015). Therefore, in the 1-2 days before the PM2.5 concentration reached the peak of heavy pollution, strengthening emergency emission reduction measures in neighboring cities will achieve a better peak shaving effect. After 2:00 on the 20 th , influenced by northerly winds, clean air from Chengde transported into Beijing, and Beijing continued to export to Baoding and Langfang (S), and the average net flux of Baoding was as high as -4000 t d -1 , resulting in a sharp drop in Beijing PM2.5 concentration to 42 µg m -3 within 5 hours. In Stage III, PM2.5 total inflow, outflow, and net flux were relatively lower, with mean fluxes of 682, -831, and 149 t d -1 , respectively. PM2.5 concentrations continued to remain a low-level, basically below 40 µg m -3 , marking the basic end of this heavy pollution.
Note that the characteristics of PM2.5 flux at the higher-altitude layer were almost similar to that of the lower-altitude layer ( Fig. 8(b)). This demonstrated that wind direction was more uniform at various altitudes, which were not conductive to the vertical diffusion of haze pollution. Additionally, there were large positive inflows from Zhangjiakou at the higher-altitude layer during stage II, which was mainly due to the strong winds rather than high PM2.5 concentration. On the whole, the total inflow, outflow and net fluxes at the higher-altitude layer were probably two times the average at the lower-altitude layer, indicating that transport contribution at the higher-altitude layer may have a greater influence on regional air quality compared to that of the lower-altitude layer.

CONCLUSION
Significant transport effects existed between Beijing and its adjacent cities in July and October, and seasonal differences were evident, with both PM2.5 total inflow and outflow flux intensities higher in autumn than in summer. Beijing exported more impact outward than received external inflow in summer and vice versa in autumn. The vertical distribution of PM2.5 net fluxes was significant, and maximum net fluxes in July and October appeared at 817 m and 1261 m above the ground, respectively. Both Baoding and Langfang showed a net inflow to Beijing in the vertical direction, while Beijing showed a net outflow to Chengde.
Total average PM2.5 net fluxes on polluted days and clean days in July and October were consistently distributed vertically, and the former was 1.71 to 6.80 times higher than the latter. Spatial and temporal distributions of PM2.5 transport fluxes at low and high altitudes during heavy pollution were consistent, where the average PM2.5 total inflow, outflow, and net fluxes at low altitude were 1380 t d -1 , -1470 t d -1 , and -90 t d -1 , respectively, and the flux at high altitude was 2.15 to 5.30 times higher than at low altitude. Meanwhile, we found that local discharging was more likely to cause extreme PM2.5 heavy pollution, and a better peak clipping effect can be achieved when emission reduction measures are initiated 1 to 2 days before heavy pollution.