An Integrated Air Quality Model and Optimization Model for Regional Economic and Environmental Development : A Case Study of Tangshan , China

An interval programming optimization model was formulated to develop effective and feasible regional economic structure adjustment plan using Tangshan Municipality in China as a case study. The optimization model was coupled with a WRF-CAMx-PSAT air quality simulation system through the estimated industrial emission sensitivity coefficients and equivalent coefficients for PM2.5 concentrations. Seven categories of industries were examined, and the results indicated that industries with higher emission sensitivity coefficients should be given priority for control. The effectiveness of the obtained optimal schemes was further assessed by the air quality simulation system. It indicated that PM2.5 concentrations in Tangshan would decrease by [33.5%, 39.3%] than those in 2013. This study provided an effective method framework for industries to maximize profits while meeting certain air quality constraints under uncertainty through the coupling of air quality simulation and optimization models.


INTRODUCTION
Ambient air is crucial for people's daily life but is confronted with serious pollution with the development of the economy.As reported by the World Health Organization (WHO) in its global haze report of 2014, only 12% of the world's population live in air that satisfies WHO air quality standards (i.e., 10 µg m -3 of annual average PM 2.5 concentration), and about half of the cities investigated have experienced air pollution with pollutant concentrations of up to 2.5 times more than WHO standards (http://www.who.int/phe/health_topics/outdoorair/databases/cities-2014/en/).In particular, it is reported in the WHO report that the cities with the most seriously polluted air are all from developing countries.Taking for an example China, an emerging market, it was reported by the environmental state bulletins of the relevant governments that the number of days with air pollution in 2014 was 175, 197 and 264 days in Beijing, Tianjin, and Shijiazhuang, respectively.To improve air quality the central government of China released a plan in September 2013 titled "National Activity Plan for Air Pollution Prevention and Control" (NAPAPPC).Other provincial governments in China also released their activity plan for control of air pollution successively including air quality targets in 2017, detail implementation plans and safeguard measures.For example, it is reported that the governments of Beijing, Tianjin, Shijiazhuang and Tangshan set the targets of the PM 2.5 concentrations in 2017 to be no more than 60 µg m -3 , 60 µg m -3 , 76 µg m -3 and 79 µg m -3 , respectively (Wen et al., 2016a).The concentrations would decrease by 25%, 25%, 33%, and 33% than those in 2012, respectively.In general, air pollution, which possesses both local and regional characteristics due to large-scale atmospheric motion, is affected by various factors including meteorological conditions and emission sources (Chen et al., 2008).To control air pollution and maintain good air quality, it is of great significance to study these factors, especially emission sources.Thus, much research has focused on the study of emission sources and has attached importance to the identification of emission source contribution and control of emission sources (Wen et al., 2016b).
A variety of simulation models have been used to investigate both the characteristics and contributions of emission sources to air quality (Olson et al., 2009;Huang et al., 2010;Dumanoglu et al., 2014;Wang et al., 2015).
For example, Vedantham et al. (2014) developed a new source-type identification method (i.e., Reduction and Species Clustering Using Episodes) to examine the tenuous association between some species and emission source types.Li et al. (2013) proposed a MM5-CMAx-PSAT modeling approach to analyze the variation of emission contribution to the urban air quality of Tangshan in China.However, further insights into the effect of emission sources on air quality should be obtained through emission source sensitivity analysis.The emission source sensitivity analysis has also been examined by using the combination of different models (Koo et al., 2009;Tang et al., 2012).Burr and Zhang (2011) made a comparison of the source apportionment simulation of PM 2.5 over Eastern U.S. using CAMx/PSAT and CMAQ source sensitivity simulations.Baker and Kelly (2014) applied source apportionment and source sensitivity approaches in a photochemical grid model to isolate the impacts of a specific facility and compared these single source impacts with in-plume measurements.
Simulation models can identify emission source contributions and then screen the sensitive emission sources (Wen et al., 2016a), and often serve for air quality management planning.In addition, optimization models, which are associated with various uncertainties and complexities, have also been applied to air quality management planning.For example, Li et al. (2006) proposed a hybrid two-stage fuzzy-stochastic robust programming (TFSRP) model for the planning of an air-quality management system.Xie et al. (2014) combined interval parameter programming (IPP), stochastic robust optimization (SRO), and multistage stochastic programming (MSP) in their robust optimization model for regional energy system management and gained different electricity generation schemes with varied system cost and system-failure risk.
Some researchers have for decades combined simulation and optimization models in order to gain better planning results.Herrera et al. (2013) designed an optimization model of process plants based on mixed integer linear programming and multimodal transport simulation; Yun et al. (2013) developed energy load-adjusted emissions estimates and methodologies for their inclusion in air quality forecasting models to eventually facilitate dynamic intervention strategies; Biagi et al. (2016) employed the combined numerical solver to obtain the optimal CO 2 injection rate via a series of simulations.McDonald-Buller et al. (2016) integrated power system and air quality modeling to evaluate time-differentiated emissions price signals on high ozone days and gained that time-differentiated pricing aimed at reducing ozone concentrations had particulate matter reduction cobenefits.In general, the combination of simulation models and optimization models has always been a challenging task, and it is difficult to couple simulation models with optimization models due to the numerous and complicated simulation processes required.The objective of this study was thus to develop a combined simulation and optimization modeling framework for industries in order to maximize total net profits while meeting certain air quality constraints and then apply it to Tangshan City, China.It aimed at simplifying the air quality simulation process, addressing the optimization model's uncertainties and investigating the effect of air quality on regional planning.A number of tasks were included: (i) formulating an optimization model for regional economic structure adjustment; (ii) utilizing the interval programming method to address uncertainties; (iii) introducing the PM 2.5 equivalent coefficients to express the relationship between the primary particulate matter (PPM) and the precursors of secondary particulate matter (SPM); (iv) coupling an air quality model with an optimization model through the sensitivity coefficients obtained from the WRF-CAMx-PSAT model; and (v) assessing the effectiveness of optimization schemes to provide more valuable alternatives for decision-makers.

OVERVIEW OF CASE STUDY AREA
Tangshan Municipality, a metropolis at the geographical coordinates of 118°14'3"-118°49'45" East and 39°34'39"-39°58'25" North, is the economic center of Hebei province in northern China.It is situated in the central region of the Bohai Bay Economic Belt with the Bohai Sea to the south, Yanshan Mountain to the north, Qinhuangdao city to the east, and Beijing Municipality and Tianjin Municipality to the west.It has a total population of 7.5 million (2013) and covers a total area of 17,912 km 2 including 13,472 km 2 of land area with 12 districts and 4,440 km 2 of sea area, as shown in Fig. 1 (Tangshan statistical yearbook in 2014).It also has rich mineral and energy resources.For example, Tangshan possesses about 6.2 billion tonnes of iron ore reserves, one billion tonnes of oil reserves, and a 100 years' history of coal mining.
As one of the biggest industrial centers in northern China, Tangshan Municipality possesses many heavy industries, such as iron/steel industries, electric power industries, and coking industries (Tangshan statistical yearbook in 2014).Although these industries have made a significant economic contribution to regional development, they have also posed great challenges to the protection of the environment.For example, the number of days in Tangshan in 2013 with air quality satisfying the national ambient air quality standard was 6% less than the average of Hebei province in the same year.According to the government statistics annual report, the observed SO 2 , NO 2 and CO concentrations in 2013 were two times higher than the averages of Hebei province.Moreover, haze weather has been occurring frequently in recent years in Tangshan due to the high emission of pollutants from various industries.These have not only caused great deterioration in the air quality of Tangshan but have also affected the air quality of Beijing and Tianjin through regional atmospheric transport.In order to effectively improve air quality, many measures have been taken to ensure the coordinated and harmonious development of the Beijing-Tianjin-Hebei region, such as the elimination of pollution-causing industries, restrictions on traffic, and the improvement of emission standards.However, it is a challenge for decision-makers to understand how to develop a best regional economy development with good air quality.

Interval Linear Programming Model
An air quality management system is associated with various uncertainties.The Interval-parameter Linear Programming (ILP) model is one of the methods used to effectively address uncertainty in optimization models.It allows the interval information to be directly communicated into the optimization process and the resulting solution (Huang et al., 1995).The ILP model can be formulated as follows (Huang et al., 1993): where , and R ± denote a set of interval numbers.The solution for model (1) can be obtained through a two-step method (Huang et al., 1995).Many coefficients in the later optimization model (e.g., pollutants emission amount, product prices, and operation expense) are uncertain, which can't be expressed by a certain number.To express the coefficients better, interval numbers (lower bound and upper bound) are applied here.So the optimization model, which is Interval-parameter Linear Programming (ILP) model, would be solved by the method of two-step algorithm.

WRF-CAMx-PSAT Air Quality Modeling System
The Weather Research and Forecasting (WRF) model is a mesoscale numerical weather prediction system that is frequently used to provide meteorological inputs for many air quality modeling systems (Wang et al., 2010).The Comprehensive Air Quality model with extensions (CAMx) is a Eulerian photochemical dispersion model, which has been widely used for pollution simulation research, especially with Particulate Source Apportionment Technology (PSAT) (Wagstrom et al., 2008;Yarwood et al., 2008;Li et al., 2013).The WRF (version 3.3) model was applied and configured using three-level nested modeling domains (112-120°E, 37-43°N), as shown in Fig. 3, where domain 1 has a spatial resolution of 27 × 27 km and has been established with a dimension of 51 × 46 grid cells, and domain 2 has a spatial resolution of 9 × 9 km and has been established with a dimension of 108 × 99 grid cells.The 3-D first-guess meteorological fields for modeling, which are of great  importance for regional air quality simulation, were obtained from the Global Tropospheric Analyses datasets provided by the US National Center for Environmental Prediction (NCEP FNL data).The emission inventory of Beijing-Tianjin-Hebei in 2013 was collected and calculated based on our group's long-term research studies (Cheng et al., 2012).Six emission source types were classified including metallurgical industry emissions, cement industry emissions, electricity industry emissions, thermal industry emissions (including residential emissions), traffic emissions, and other industries emissions (including pollutant emissions from industries except the former five emission source types) (Huang et al., 2010;Zhou et al., 2012aZhou et al., , 2014Zhou et al., , 2015)).Thermal industry emissions in the inventory here included residential emissions by use of coal.However, the thermal industry emissions did not give complete residential emissions statistics, which caused some uncertainties in the inventory.The outputs of WRF were used for the input meteorological fields for CAMx, and PSAT was introduced into CAMx to provide source apportionment of PM 2.5 from different emission source categories and regions.Consequently, the hourly concentrations of primary and secondary PM 2.5 in each grid and their contributions would be calculated.The hourly concentrations and contributions of Tangshan Municipality, which occupies three grids in the 9 × 9 km grids, can thus be calculated.In addition, the hourly transport concentrations and contributions to Tangshan from boundary areas can be obtained by the meteorological mechanisms option and source regions partitions.ACM 2 (Asymmetric Convective Model) was chosen for Planetary Boundary Layer Scheme, the biases of which were studied by many researchers (Perez et al., 2006;Hu et al., 2010;Zhang and Yin, 2013;Dawn and Mandal, 2014).Five boundary regions, which were Beijing, Tianjin, North of Hebei, South of Hebei, and Shandong, are selected for transport concentrations and contributions.

PM 2.5 Equivalent Coefficient
PM 2.5 consists of both primary and secondary fine particulates.The primary PM 2.5 , such as Element Carbon (EC) and Primary Organic Aerosol (POA), is emitted directly from the emission source.The secondary PM 2.5 (e.g., SO 4 2-, NO 3 -, and Secondary Organic Aerosol (SOA)) are indirectly transferred from their precursors (e.g., SO 2 , NO x , and VOC s ) (Guerra et al., 2014).Thus, it is of great importance to investigate the relationship between PM 2.5 and their precursors (e.g., SO 2 , NO x , VOC s ) in order to further quantify the emission source contributions.PM 2.5 equivalent coefficients were introduced to establish such relationship, as follows.
where PPM 2.5 denotes emission amount of primary PM 2.5 ; α, β and η denote equivalent coefficients of SO 2 , NO x and VOC s , which are defined as equivalent rates of PPM per unit emission amount of SO 2 , NO x and VOC s , respectively.Meanwhile, the CAMx-PSAT model, which was applied to simulate PM 2.5 composition concentrations, was used for calculating the PM 2.5 equivalent coefficients α, β, η based on the assumptions: (i) SO 4 2-, NO 3 -, SOA are all converted from SO 2 , NO x , VOCs, respectively, (ii) the conversion proportion of SO 2 , NO x , and VOC s to SO 4 2-, NO 3 -, and SOA remains the same for different emission sources, and (iii) the PM 2.5 , PPM 2.5 and the precursors are from local emission.The equivalent coefficients are calculated as follows.
where C p (µg m -3 ) denotes the concentration of PM 2.5 components related to pollutant p (p = 1 for PPM 2.5 , p = 2 for SO 2 , p = 3 for NO x , p = 4 for VOC s ) based on modeling simulations; Q p (t a -1 ) denotes the emission amount of pollutant p; δ p (µg (m 3 t a -1 ) -1 ) is contribution concentration to PM 2.5 components per unit emission amount of pollutant p, and it can be used to calculate the emission source sensitivity coefficient as follows.
' , ' , ' , where δ' p,i (µg (m 3 t a -1 ) -1 ) is the emission source sensitivity coefficient, which can be defined as contribution concentration to PM 2.5 components per unit emission amount of pollutant p from emission source i; C' p,i (µg m -3 ) denotes the concentration of PM 2.5 components related to pollutant p from emission source i; Q' p,i (t a -1 ) denotes the emission amount of pollutant p from emission source i.

Optimization Model Construction
Consistent with the emission source categories in the inventory, seven typical industries in Tangshan's economic structures were considered in the optimization model including iron/steel industry, cement industry, electricity industry, coking industry, thermal industry, traffic industry, and other industries.The metallurgical industry in the emission inventory was classified into iron/steel industry and coking industry.Other industries referred to the industries excluding iron/steel industry, cement industry, electricity industry, coking industry, thermal industry, and traffic.In order to easily model them in the optimization model, other industries were classified into three typical industries, whose products were paper, ceramics, and gypsum, respectively.The optimization objective was to maximize total net profits from all the industries (e.g., profits, operating expense, pollutants treatment expenses) while meeting certain constraints.Measures referred in the air pollution control plan of Tangshan from 2013 to 2017 were for developing the constraints, such as coal consumption, production capacity, pollutant mitigation goal, and the demand and supply of energy and resources.The interval linear programming optimization model was then constructed to address the uncertainties in optimization parameters, such as the capacity utilization and emission coefficients of the industry.The optimization objective is as follows.
Eq. ( 3) indicates the maximum net profits from all the industries including two parts, the benefit and the cost, where f ± is the net profit (RMB); cost is the total cost of operation and pollution treatment for all the industries (RMB); benefit is the total profits from the products of the industries (RMB); i, j represent the types of industry and product, respectively (shown in Table 1); r represents pollutants, where r = 1 for SO 2 , r = 2 for NO x , r = 3 for PM, r = 4 for VOCs; OPU  denote coal, electric power, natural gas, coking and oil consumption amount per unit product j for industry i, respectively; PC ± , PE ± , PG ± , PJ ± , PO ± denote the annual mean price of coal, electric power, natural gas, coking and oil, respectively; , i j PA  represents amount of product j for industry i; , , i j r CIP  denotes treatment expense per unit reduction of pollutant r for product j of industry i; , i j CP  denotes the annual mean price of product j for industry i.
The optimization constraints are described below.(1) Total amount of coal consumption constraints , , 1 1 Eq. (3-c) is about the constraints on total coal consumption.It indicates that the total coal consumption of all the industries in 2017 would be less than the constrained coal consumption of 2017 in the air pollution control plan of Tangshan from 2013 to 2017, where OTC represents the amount of coal consumption amount in 2013 (10 4 t); and TC ± represents the decrease of coal consumption from 2013 to 2017 as reported in the air pollution control plan of Tangshan from 2013 to 2017 (10 4 t).
(2) Production capacity constraints , , , , , Eq. (3-d) is about constraints on the production capacity of every industry.It indicates that the production capacity of every industry would satisfy the social demands for their production, where OAS i,j represents the amount of product j for industry i in 2013; GAS i,j represents the production capacity reduction of product j for industry i from 2013 to 2017 as reported in the air pollution control plan of Tangshan from 2013 to 2017; and , i j   represents the capacity utilization (i.e., production amount per unit production capacity) of product j for industry i.
(3) Environmental constraints , ,, 1 1 1, 2 Eq. (3-e) is about constraints on SO 2 and NO x emission.It indicates that pollutants emission amount in 2017 would be less than the constrained pollutant emission amount of 2017 in the air pollution control plan of Tangshan from 2013 to 2017, where , , i j r EA  represents the emission amount of pollutant r per unit product j of industry i; r TS  represents the emission amount of pollutant r in 2013, (t); and r GS  is the mitigation amount of pollutant r from 2013 to 2017, (t)., ,, , 0 1 1 1 Table 1.Detail information of industry and product.
Thermal Industry j = 1 Heating Power i = 5 Coking Industry j = 1 Coke i = 6 Traffic Industry j = 1 Heavy Duty Vehicles j = 2 Medium Duty Vehicles j = 3 Light Duty Vehicles j = 4 Motorcycle i = 7 Other Industries j = 1 Paper j = 2 Ceramics j = 3 Gypsum Eq. (3-f) is about constraints on PM 2.5 concentration.It indicates that PM 2.5 concentration in 2017 caused by the industries would be lower than the PM 2.5 concentration target C 0 (79 µg m -3 ), where EC r represents equivalent coefficient related to pollutant r; SC i,r represents sensitivity coefficient of pollutant r for industry i, µg (m 3 t a -1 ) -1 , which can be calculated by δ' p,i ; and C 0 denotes the goal PM 2.5 concentration for 2017 (79 µg m -3 ).( 4) Demand constraint M 2 6, 1 (1 ) Eq. (3-g) is about constraints on vehicle amount and it indicates that the vehicle amount in 2017 would not only satisfy social demands but remain the lowest rate of increase, where TT ± denotes the total amount of vehicles in 2015 (vehicles); and τ denotes the maximum rate of increase for vehicles in Tangshan.
Eq. (3-h) is about constraints on electricity demands and it indicates that electricity demands for industries in 2017 would be more than the minimum electricity demands, where IE ± denotes electricity demands for satisfying basic industry development in 2014, (GWH); and σ denotes the minimum rate of increase for industrial electricity demands.

Data Collection
The economic data for the optimization model including product prices and operating expense for pollutants treatment per unit product, were obtained based on annual industry reports, annual government work reports, statistical bulletins and previous studies (Wen, 2015;Zhen et al., 2016).Operation expenses and product prices for different industries were shown in Table 2.In addition, some data about industry pollutant emission and pollutant treatment cost were obtained from Tangshan government through long-term co-operation on programs.The price per unit product was obtained by calculating the average price of the product in recent years.Pollutant emission coefficients of different industries were obtained based on the documents released by Ministry of Environmental Protection of the People's Republic of China and literature references (Zheng et al., 2009;Huo et al., 2010Huo et al., , 2012;;Huang et al., 2011;Sun et al., 2011;Zhao et al., 2012;Zhou, 2012a;Lang et al., 2012), as shown in Table 3.Additionally, data of pollutant emission amount per vehicle in the traffic industry was calculated by the annual mean mileage of various vehicles, fuel consumption per mileage, and pollutant emission coefficients per unit fuels.

Air Quality Model Verification
The normalized mean gross error (NME) and correlation coefficients (RC) were applied to assess the performance of the CAMx-PSAT simulation (USEPA 2007).The simulated concentrations of PM 2.5 in Tangshan in 2013 were compared with the observed mean concentrations of PM 2.5 which were obtained from six normal state-controlled monitoring sites in Tangshan.To make it clear, we only compared the data of four representative months (i.e., January, April, July and October) (shown in Fig. 4 and Table 4).Generally, the simulated results were lower than the observed value.This might be because the emission inventory of Tangshan were  (Wen, 2015;Zhen et al., 2016).underestimated (e.g., incomplete data collection of industries' information and underestimated pollutants emission coefficients).In addition, there might be some significant problems exiting in simulation models, which caused the underestimation.For example, the heterogeneous reactions, which are not well captured in the model, have been proposed to be important in the troposphere and contributed a great deal to model accuracy (Sarwar et al., 2008;Wang et al., 2012).Mueller et al. (2006) proposed that heterogeneous reactions were very important for the model and a model must correctly simulate cloud cover to avoid serious bias.Zhang et al. (2011) incorporated nitrous acid, direct emissions, two heterogeneous reactions, and two surface photolysis reactions into the model to improve the model performance and obtained reasonably good results.The NMB value for PM 2.5 in four months were 25%, 27%, 21%, and 23%, respectively, which is a bit larger (Wen et al., 2016a).This might be because there were uncertainties not only in emission inventories but also in the meteorology simulation, which would affect the accuracy of simulated PM 2.5 concentrations.However, the RC value were better than those in some of previous studies (e.g., Zhou et al., 2012c;  Lang et al., 2013; Wen et al., 2016a).Compared with other related studies, the simulation results in this study were acceptable (Liu et al., 2010; Li et al., 2013; Wu et al., 2013;  Zhang et al., 2013b).Sensitivity Coefficients for Different Industries PM 2.5 concentrations and contribution percentages of various industries to Tangshan Municipality were obtained from CAM x /PSAT simulation (shown in Fig. 5).It was found that the iron/steel industry and coking industry, which were the pillar industries of Tangshan Municipality's economy, were the primary emission sources (Li et al., 2013).This illustrated that pollutant mitigation and economic development should be better balanced in the future for the harmonious development of both the environment and the economy.Additionally, other industries contributed a big percentage (22.8%) to the PM 2.5 concentration of Tangshan Municipality, indicating that the control of other industries would not be neglected.Moreover, each of traffic, electricity industry, thermal industry, and cement industry accounted for over 10% of the PM 2.5 concentration.
The sensitivity coefficients for different industries (shown in Table 5) were obtained based on the CAM x /PSAT simulation and pollutant emission amount (shown in Table 6).It was found that the sensitivity coefficients from traffic, thermal industry, and other industries were significantly greater than those from the other four industries.This indicated that the PM 2.5 concentration contributions of the three industries per unit pollutants emission were larger than those of the four industries.The sensitivity coefficients of other industries were the biggest for each pollutant, and they were significantly greater than those of the iron/steel industry and the coking industry.As compared with the results shown in Fig. 5, it was obvious that the control of other industries had a more significant effect on air quality improvement.This might be because other industries were located in regions which were more sensitive to air quality.
The sensitivity coefficients, which were required information for optimization model in Eq. (3-f), were of great significance in the combination of simulation model and optimization model.The accuracy of the coefficients would directly affect the optimized results.However, there were some uncertainties in the coefficients.According to the equations in section "3.3 PM 2.5 equivalent coefficient", it was obvious that not only the PM 2.5 equivalent coefficients but also the sensitivity coefficients were closely related to pollutants emission inventory and simulated pollutant concentrations.Therefore, an accurate inventory is of great importance.However, there are several hundred industrial facilities in the study region, and the accuracy of the surveyed information could not be ensured.In addition, the pollutant emission coefficients of industries released by China National Environmental Monitoring Centre could not reflect all of the emission characteristics and therefore should be localized.Pollutants emission amount were various on the different operation conditions.Both of these could cause the uncertainties in the inventory.Another important factor might cause the uncertainties of the coefficients was the nonlinearity.The non-linear relationship between pollutant concentrations to the changes precursors' emissions was complex and a challenge to deal with (Xing et al., 2010;  Zhao et al., 2015).To simplify the optimization model, we assume a linear relationship between them.This would cause uncertainties and needed further studies to improve the sensitivity coefficients.Moreover, the uncertainties of various modeling parameters, such as background emission inventory input, boundary conditions, and initial conditions, were unavoidable due to the complex interactions of atmospheric chemistry and physics (Zhou et al., 2012b).All of the above would contribute greatly to the uncertainties of the coefficients, which could vary a great deal across regions.For example, the sensitivity coefficients in Tangshan could be 0.0007-0.8956µg m -3 (t d -1 ) -1 (Zhou et al., 2012b) and the equivalent coefficients in Beijing-Tianjin-Hebei could be 0.049-0.27(Wen, 2015).

Air Quality Targets
The WRF-CAMx-PSAT model was applied to simulate and evaluate the effectiveness of the obtained optimal schemes.The PM 2.5 concentration target in 2017 is 79 µg m -3 (Government of Tangshan 2013).When the value of optimization results were put into the left-hand side of Eq. (3-f), the PM 2.5 concentration for Tangshan calculated would be about [71.8, 78.6] µg m -3 , which was about [33.5%, 39.3%] less than that in 2013.It is obvious that the optimized schemes can satisfy the air quality target.However, this is based on the current situation of the lower capacity utilization of industries (less than 70%), which will not only cause more wasting of resources but also lead to more emission of pollutants.If more compulsory measures for higher levels of capacity utilizations were implemented by governments, then more industrial output would be produced with the same energy consumption.Less energy would be consumed for the same social demand, thereby leading to less amount of pollutants emission and more resources conservation.

Optimization Schemes
Fig. 6 presents the optimal schemes of economic structure adjustment for different industries in Tangshan Municipality.Two optimal schemes with relaxed constraints and with tightened constraints were obtained through the solution method of Huang et al. (1995).A large range of solutions were obtained when we used relaxed constraints while a small range of solutions were obtained when we used relaxed constraints.Taking Eq. (1-a) for an example, when the right-hand side value of the equation is B + , the constraint is relaxed.When the right-hand side value of the equation is B -, the constraint is tightened.We found that when the constraints were tightened, the iron/steel industry, the coking industry, the traffic, the thermal industry, and other industries obtained lower bound value.When the constraints were relaxed, the cement industry and the electricity industry obtained upper bound value.This might be because (i) the traffic, the thermal industry, and other industries obtained bigger sensitivity coefficients shown in Table 5, which contributed a great deal in PM 2.5 targets realization for Eq.(3-f); (ii) the cement industry and the electricity industry obtained relatively smaller sensitivity coefficients in the constraints and relatively bigger economic parameters in the objective function(i.e.parameters of net profits per unit product) than those of the traffic, the thermal industry, and other industries; and (iii) the iron/steel industry and the coking industry obtained relatively smaller economic parameters than those of the cement industry and the electricity industry.It was also found that a large interval existed in some optimized production, such as [5.7 × 10 6 , 28.7 × 10 6 ] GJ for heating power, [12.8 × 10 6 , 22.7 × 10 6 ] t for coke, and [7.4 × 10 6 , 23.1 × 10 6 ] t for gypsum.This might be because these products were associated with a large interval in the optimization modeling parameters owing to limitation in data collection of local real situation, such as pollutant emission coefficients, energy consumption per unit product, and pollution treatment expense per unit product.The solution for ILP was obtained through a two-step method (Huang et al., 1995), which has some limitations and needs further improvement (Zou et al., 2010; Zhang et  al., 2013a).Compared with production in 2013 for different Note: the units of the production for clinker, cement, coke, paper, ceramics, gypsum, heating power, electric power, steel, iron, heavy duty vehicles, medium duty vehicle, light duty vehicles and motorcycles are 10 6 t, 10 6 t, 10 6 t, 10 6 t, 10 5 t, 10 6 pieces, 10 6 GJ, 10 3 GWH, 10 6 t, 10 6 t, 10 3 vehicles, 10 3 vehicles, 10 4 vehicles and 10 4 vehicles, respectively.industries (shown in Table 6), all the optimal production, except vehicle population, was less than that of 2013.For example, the optimal amount of clinker, cement, and coke would be about [28.1%, 37.4%], [30.0%, 36.4%], and [19.8%, 54.8%] less than that of 2013, respectively.This indicates an economic structure adjustment should be implemented in Tangshan to reduce the scales of industries and to ensure good air quality and economic development.Nevertheless, Tangshan was faced with serious problems in economic structure adjustment.This might be because (1) industries in Tangshan possessed excessive production activity; and (2) some production processes had lower capacity utilizations and lower pollutant treatment efficiency.In fact, the capacity utilization of industries in Tangshan was lower than the international standard.According to the statistical data from the government, it is about 65%-80% in iron/steel industry and cement industry.Moreover, in some regions, some factories would prefer to be in debt for a long time rather than stop production, which puts a large pressure on local resources conservation and environment protection.Thus, mandatory measures for eliminating backward production processes and factories have been taken in Tangshan Municipality, especially between 2013 and 2017.Additionally, although a decrease of the vehicle population would contribute to better air quality, the optimal vehicle population increased by about [3.5%, 14.9%] than that in 2013.The annual mean increasing rate is [0.9%, 3.7%], which is lower than the lowest rate of increase from 2009 to 2014, as shown in Fig. 7.This might be because with the development of the economy demand for vehicles is becoming stronger.However, this trend would pose a pressure on the local government to implement some kind of vehicle control policy, such as eliminating yellow label cars (namely vehicles with heavy pollutant emissions) as well as improving vehicle oil quality.
The optimal and original (2013) pollutant emission amounts of the various industries are shown in Fig. 8.Although the emission sensitivity coefficient of iron/steel industry was the lowest, it was found that the change in percentage of pollutant emission amount in iron/steel industry was significantly greater than that in the other six industries.This is because in the crucial activity plan for the control of air pollution from 2013 to 2017 the iron/steel industry was compelled to eliminate a large amount of capacities owing to backward production process and energy waste.In addition, the cement industry also had a large percentage change in amounts of SO 2 , PM and VOC s emission.Other industries possessed the biggest change percentage of pollutant emission amount among all the industries excluding iron/steel and cement industry.This is because the largest sensitivity coefficient was obtained by other industries.In order to satisfy air quality and obtain better regional economic benefits, other industries should implement a relatively large change in pollutant emission reduction.Moreover, traffic and thermal industry also possessed relatively large change percentages of pollutant emission amount owing to their relatively large sensitivity coefficients.As a result, industries with large sensitivity coefficients should be given special attention by the government for regional economic and environmental development.

Assessment
The optimized schemes obtained from the proposed modeling framework can not only maximize the total net profits from all the industries for the planning region but also contribute to achieving air quality targets.Because the simulation model can be applied for any region and the optimization model is suited to study areas with the same types of industry, the proposed modeling framework is applicable to other regions.However, some shortcomings still exist in the research that need further study in the future.This study was based on a few assumptions including (i) equivalent coefficients for different emission sources between PM 2.5 and SO 2 , between PM 2.5 and NO x , and between PM 2.5 and VOC s remain invariable, respectively; (ii) the yearly variation of the equivalent coefficients was neglected; (iii) all the SO 4 2-, NO 3 -and SOA were transferred from the precursors SO 2 , NO x and VOC s , respectively, neglecting the complex and nonlinear processes; (iv) the nonlinear relationships between pollutant emissions and pollutant concentrations were neglected.As a result, more sampling and analysis are required in future studies for more reliable source apportionment, emission sensitivity coefficients estimation, and nonlinear relationships between pollutant emissions and pollutant concentrations.
Furthermore, two environmental constraints, which were emission constraints and concentration constraints, were considered in the framework.They were connected though they were different.When the pollutant emission constraints were considered without considering the PM 2.5 concentration   constraints, the value of PM 2.5 concentrations calculated by the optimization model were less than those of PM 2.5 concentration targets.It indicated that only the emission constraints might not be enough for achieving air quality targets.However, if the emission constraints for ever region or some typical regions of Thanshan are considered, the PM 2.5 concentration targets might be achieved.Future work for typical regions selection and optimal emission amounts for every region would be more meaningful, and more valuable for decision-makers.Additionally, inter-regional pollution transport was neglected in the optimization model, because the pollutant emissions of the boundary regions and the meteorological conditions in 2017 were uncertain.To avoid the uncertainties from the pollutant concentration prediction of inter-regional pollution transport, inter-regional pollution transport was not considered in the environmental constraints of the optimization model.It indicated that constraints in Eq. (3-f) were relaxed constraints.Therefore, the optimal production of industries might be bigger than those under the condition considering inter-regional pollution transport impact.
Secondary ions, which are important compositions of secondary particulate matter, are significant contributors to secondary particulate matter concentrations.Meanwhile, there exist nonlinearities and uncertainties in the relationship between secondary particulate matter and their precursor gases owing to the complicated atmospheric physical and chemical processes.Thus, it is necessary to investigate the relationship between secondary ions and precursor gases.Fu et al. (2014) investigated the relationship between secondary ions and their precursors in the Pearl River Delta region in China and found that each 1% SO 2 concentration decrease would lead to 0.59% SO 4 2-concentration decrease and each 1% NO x concentration decrease would lead to 0.97% NO 3 -concentration decrease.Wen et al. (2015) also investigated such relationship in the Beijing-Tianjin-Hebei region during the 2014 APEC summit and found that each 1% SO 2 concentration decrease would lead to [0.6%, 0.8%] SO 4 2concentration decrease and each 1% NO x concentration decrease would lead to [0.6%, 0.7%] NO 3 surface concentration decrease.In our study, the results demonstrated that each 1% SO 2 concentration decrease would lead to about [0.7%, 0.8%] SO 4 2-concentration decrease and each 1% NO x concentration decrease would lead to about [0.7%, 0.8%] NO 3 -concentration decrease in Tangshan.Compared with previous research, the mean simulated results were slightly higher than those reported in Wen et  al. (2015).In addition, Fu's results relevant to SO 2 were slightly lower than Wen's, and Fu's results relevant to NO x were slightly higher than Wen's.In general, the difference is not very large.Because the process of secondary ions transformation from their precursors is a complex chemical process and is closely relevant to meteorological conditions, the atmospheric oxidation and environmental power of hydrogen, the secondary ions transformation possessed to some degree seasonal characteristics.Although the relationship between secondary ions and their precursors is complex and uncertain, researches in this field are necessary and important not only for the PM 2.5 source apportionment but for the control of hybrid air pollution.

CONCLUSION
The emission sensitivity coefficients of typical industries in Tangshan Municipality were obtained based on WRF-CAM x -PSAT air quality simulation.The simulation results were then applied to the interval programming optimization model for regional economic structure adjustment planning under uncertainty.The optimization results indicated that the industries, especially other industries which had bigger sensitivity coefficients, should be given priority for emission control.The optimal annual increase rate of vehicles population in Tangshan should be controlled in range [0.9%, 3.7%].If the optimal economic structure adjustment schemes were implemented, the PM 2.5 concentration of Tangshan in 2017 would decrease by about 37.7% than that in 2013.It would satisfy the air quality goal of Tangshan, indicating that the obtained optimal schemes were effective and reliable.This study provided an effective method framework for regional economic and environmental planning under uncertainty by integrating advanced air quality simulation with system optimization model.

Fig. 1 .
Fig. 1.Schematic diagram of Tangshan Municipality and its main industrial distributions.

Fig. 2 .
Fig. 2. Schematic diagram of the combination of air quality simulation model and optimization model.

Fig. 5 .
Fig. 5. PM 2.5 concentration contribution and contribution percentage for different industries.

Table 2 .
Operation expense and product prices for different industries.
a based on the product prices during the past five years; b denotes no operation cost in this industry or without pollutants treatment; c based on annual government work reports, statistical bulletins and previous studies

Table 3 .
Pollutant emission coefficient of different industries.

Table 6 .
Production in 2013 for different industries.