Exploring Regional Fine Particulate Matter (PM2.5) Exposure Reduction Pathways Using an Optimal Power Flow Model: The Case of the Illinois Power Grid

This work develops an exposure-based optimal power flow model (OPF) that accounts for fine particulate matter (PM2.5) exposure from electricity generation unit (EGU) emissions. Advancing health-based dispatch models to an OPF with transmission constraints and reactive power flow is an essential development given its utility for short- and long-term planning by system operators. The model enables the assessment of the exposure mitigation potential and the feasibility of intervention strategies while still prioritizing system costs and network stability. A representation of the Illinois power grid is developed to demonstrate how the model can inform decision making. Three scenarios minimizing dispatch costs and/or exposure damages are simulated. Other interventions assessed include adopting best-available EGU emission control technologies, having higher renewable generation, and relocating high-polluting EGUs. Neglecting transmission constraints fails to account for 4% of exposure damages ($60 M/y) and dispatch costs ($240 M/y). Accounting for exposure in the OPF reduces damages by 70%, a reduction on the order of that achieved by high renewable integration. About 80% of all exposure is attributed to EGUs fulfilling only 25% of electricity demand. Siting these EGUs in low-exposure zones avoids 43% of all exposure. Operation and cost advantages inherent to each strategy beyond exposure reduction suggest their collective adoption for maximum benefits.


Network Equations
A brief background for the optimal power flow (OPF) formulation will be provided here. The purpose behind this background is to provide those who may not be entirely familiar with OPF models with an accessible introduction to their fundamental theoretical components. For a more detailed overview, interested readers are referred to Frank & Rebbenack, 1 which this formulation adopted from. Table S1 summarizes the model's variables and parameters.
Mathematical models are used for the formal and systematic representation of power systems as well as for solving the OPF problem. A general setup of the OPF problem and necessary background will be presented here. The most common representation of power systems is through an undirected graph. Within this study, a general power systems network will be represented by Η = ( , ), where represents the set of nodes (i.e., buses) and represents the set of edges within the network. Each edge represents a ℎ ( , ) that connects an upstream ( ) to a downstream ( ). Branches can be physical representations of different power systems equipment such as lines or transformers. Relevant ℎ ( , ) parameters that characterize power systems equipment include the impedance (ℤ ), which is broken up into ( ) and ( ) as per (1). Within this notation, roman letters (e.g., ℤ, , , ) will be used to signify phasor quantities.

ℤ = + ( )
The current on a line is generally represented by ( ). This can be used to represent the general form of Ohm's law to calculate the voltage drop ( ) across a circuit element, which is shown in (2). In this general formulation S3 = ℤ ( ) Therefore, the current leaving (or injected into) ( ) can be represented using (3).

= ℤ ( )
In power flow problems, it is often more convenient to replace the impedance (ℤ ) of a line by its admittance ( ) as per (4).
where ( ) represents the ℎ and ( ) represents the ℎ . Therefore, Ohm's law can be rewritten as shown in Equation (5).

= ( )
This allows one to write Ohm's law for an entire network, and it can be represented compactly in matrix notation as per (6), where ̅ = ( 1 , … , ) is the vector of phasor source currents that are injected into each bus, ̅ = ( ̅ 1 , … , ̅ ) is the phasor bus voltage vector, and ( ̿ ) is the complex bus admittance matrix.

Power Flow Equations
The power flow problem is a feasibility problem that seeks to determine a solution to the network equations, without a specific objective function. The power flow equations are thus used as constraints in the OPF formulation to ensure that the optimal solution to the OPF problem is a S4 feasible one. To derive the power flow equations, we begin by assuming that each node within the network will represent a bus where power can be generated or consumed. Each ( ∈ ) will be characterized by the following four variables: 1. Voltage magnitude: | | 2. Voltage phase angle: 3. Net real power injection: 4. Net reactive power injection: Solving the power flow problem requires developing a relationship between all four variables at all points within a network. Doing so requires developing a system of power flow equations for the network which will have the form shown in (7) and (8). The final form of the equations will be shown in (14) and (15) after they are derived in this section. Equations (7) and (8) represent power injections into each , which also equate to the difference between power generation ( , ) and power load ( , ) at each bus. Next, we move towards specifying the exact form of equations (7) and (8). This is achieved by transforming the current flows within the network into power flows. Using the network equations The equations can also be written in rectangular form by taking the real (14) and imaginary parts (15), separately. These represent the full form of the power flow equations that will be used as constraints in the OPF to ensure the feasibility of the OPF solution.

AC Optimal Power Flow
The OPF can be formulated as an optimization problem when the power flow equations derived in (14) and (15) are combined with an objective function. Objective functions of OPF problem typically minimize total costs of electricity generation while also ensuring that the system is operating within safety limits. The optimization variables within the OPF formulation are the real power injections at generator buses (Γ ⊆ ), while voltage magnitudes and angles are dependent state variables used to formulate the constraints. The AC OPF problem for a specific time interval (e.g., 1-hour) can be formulated as shown in (16).
S7 (16) min ∑ ( ) ∈Γ s. t. Where: , ∈ is the power produced by the ith generator and represents the optimization variable of the model.
( ) is the cost to operate unit i at the output level .
The OPF will minimize system costs by controlling the power produced at all the generators subject to the network constraints. The constraints can be interpreted as follows: • The first two equality constraints denote the power flow equations derived in the previous subsection that must be satisfied. AC power flow equations must be satisfied by the converged to solution in order for it to be a feasible operating point of the power S8 network. The load real ( ) and reactive power ( ) are given and fixed for a specific time-period (i.e., 1-hr time interval in this instance).
• The third constraint represents the real power limits at the buses. ( , ) is usually set as 0 unless a generator has a fixed minimum generation requirement. Renewable sources are set to their historic hourly output by fixing their real power generation limits.
• The fourth constraint represents the reactive power limits at the buses. Reactive power loads and minimum and maximum limits are based on power factor data approximated from other synthetic networks. [2][3][4] • The fifth constraint represents the branch limits, which sets an upper limit on the magnitude of the current that can pass by any branch within the network. This constraint is derived using Ohm's law, where branch current is dependent on the voltages between the branch and the branch admittance. Max loading is also based on parametrized line data. [2][3][4][5] The voltages in the constraint are expanded and written in rectangular form to match the rest of the formulation, but it can be compacted and written in phasor form as follows: • The sixth constraint represents the maximum and minimum allowable bus voltage magnitude limits. For this formulation, | | is set to 0.98 p.u. while | | is set to 1.04 p.u.
• The seventh constraint represents the maximum and minimum allowable bus voltage angle limits. For this formulation, is set to −45° while is set to +45°.

DC Optimal Power Flow
To quantify the additional costs and damages that arise from accounting for network constraints and reactive power flow, we employ a DC power flow approximation to the network. A DC power flow approximation is a linearized version of the AC power flow equations. The following assumptions hold for the DC power flow approximation: • There is no reactive power flow on the network.
• The transmission system is assumed to be lossless with branch resistances ( = 0).
• All bus voltages in the network are at unity (| | ≈ 1.0).
Applying these assumptions to equation (15) leads to the following approximation for real power transfer: Thus, the DC approximate optimization program is as follows:

Network Development
Details regarding network development, assumptions, and data sources are provided below.

Transmission Line and Substation Data
Transmission line and substation data are obtained from the Homeland Infrastructure

Load Demand Data
Running an OPF requires knowing the spatial distribution of load within the network. To our knowledge, the exact distribution of network load is not publicly available. Therefore, several steps are taken to obtain a representative, spatially and temporally resolved load profile for this network. It is assumed that power load demand across the network is consumed at the substation level, specifically at substations where voltage levels are stepped down to distribution level S13 voltages. Substations which meet that voltage criteria are identified and are designated as load buses. A load weight is then assigned to each of the load buses. The load weight represents the percent of total load in the network consumed at that bus. It is quantified based on the location of the load buses within the network with respect to Illinois' Electric Retail Service Territories. 9 Each load bus is assigned to the Retail Service Territory(ies) it lies within. A load weight for each territory is developed based on the ratio of its total electric sales relative to the total electricity sales of all territories within Illinois. It is assumed that the utilities that lie within Illinois' Retail Service Territories, in addition to the external balancing grid, meet all the electric load demand in Illinois. The load weight is then equally distributed amongst all load buses within each territory. Since certain load buses lie within two or more territories, their respective load weights from each territory are summed to obtain a total load weight for each bus. Figure 1(b) in the main text shows the distribution of load weights to the network substations.
The next step involves obtaining the total annual electric energy demand for the state of Illinois.
Specifically, 2019 annual electric energy consumption data is used. 10,11 The temporal variation in total system load demand for each 24-hr period of our simulation was generated based on hourly load data obtained from Midcontinent Independent System Operator (MISO), which most of the state's power utilities are operated by. 12 Once the total hourly load for the entire network is known, it is then spatially distributed to the load buses based on their corresponding load weights. The load buses are also augmented with reactive power data based on power factor data approximated from other synthetic networks. 2,4,8

EGU Data
EGU data are obtained from publicly available eGRID2019. 10 Each EGU is designated as a generator bus and is then connected to the network based on its location. Since voltage ratings of S14

Exposure Reduction Strategies Analyzed
Details regarding the development and assumptions of the different strategies are provided below.

Post Combustion NOx and SOx Emissions Control
In  to retrofit each EGU with emission control technologies was assumed to be possible without further investigation into any necessary installation prerequisites.

Higher Renewable Generation
The next strategy models the impacts of an increase in the renewable generation within the network, specifically from wind and solar sources. The 2019 baseline renewable generation is around 7.5% (eGRID). Three different scenarios are adopted to analyze the effects of integrating more renewables within the network. The first two are based on Illinois' renewable portfolio standards which assume 25% and 50% renewable generation by 2026 and 2040, respectively. 19 A third hypothetical scenario assumes 75% renewable generation by 2050 to model the effects of an aggressive renewable adoption. The expected increase in electricity demand (around a 1% increase annually) by the given years is factored in to reflect the amount of fossil-fuel based generation that will still be taking place during those years. 20 This allows the results to reflect the actual exposure damages that would still take place within those years while still running the OPF based on current load demand, allowing this scenario to be comparable to the others. To make the renewable scenarios comparable to the others, it is assumed that the expected increase in renewable generation capacity is adopted by the current grid conditions, so that all other factors (e.g., population numbers, underlying health risks, electricity demand, grid infrastructure) remains the same. Thus, we are not modeling future scenarios, rather what exposure would take place if future renewable generation capacity was present today. The expected renewable share S17 split between wind and solar is also accounted for based on EIA projections. 20 The increase in renewable generation within the OPF is modeled by increasing generation and capacity from existing solar and wind plants. Transmission equipment capacity is adjusted to handle the expected increased load from these sources. The capital costs of installation of the plants are quantified following a methodology outlined in Sergi et al. 21,22 The total required annual generation from these plants is quantified followed by the required installed capacity to meet generation demand using the plants' capacity factors. Capacity factors for wind and solar sources are derived the National Renewable Laboratory's (NREL) Wind Integration National Dataset (WIND) 23 and Global Horizontal Irradiance (GHI) 24 databases, respectively. To control for the rising intermittency issues from increased renewable generation capacity, installed utility-scaled battery storage capacity is also modeled following Sergi et al. 21,22 Additional battery storage leads to an increase in both operation and capital costs.

Relocation of High Polluting Plants to Low iF Zones
The final strategy models relocating high polluting EGUs to lower intake fraction ( actually host these EGUs from a regulatory perspective is beyond the scope of this study.

Linear Exposure Damages Factor
Given that nonlinear concentration-response functions are computationally difficult to integrate into the OPF model, linearized PM2.5 intake effect and severity factor parameters 26 Figure S1 shows a map of the exposure domain and network EGUs. 10 Figure S1: Exposure domain considered in the exposure-base optimal power flow (OPF) model in this study. Exposure is measured at the census tract level. The domain includes the state of S21 Illinois and its 10 neighboring states to capture the impacts of secondary PM2.5 exposure:

Exposure Domain
Arkansas, Indiana, Iowa, Kentucky, Michigan, Minnesota, Missouri, Ohio, Tennessee, and Wisconsin. Illinois based electricity generation units (EGUs) disaggregated by plant type (which represent the emissions sources in this study) are also highlighted. Figure S2: Annual electric generation [TWh/y] disaggregated by EGU plant type.