domain: spi math

text:
Sina Baghali

Zhaomiao Guo

Impacts of Privately Owned Electric Vehicles on Distribution System Resilience: A Multi-agent Optimization Approach

Keywords: Distribution system resilience, electric vehicles, incentive design, multi-agent optimization, convex reformulation

We investigate the effects of private electric vehicles (EVs) on the resilience of distribution systems after disruptions. We propose a framework of network-based multi-agent optimization problems with equilibrium constraints (N-MOPEC) to consider the decentralized decision making of stakeholders in transportation and energy systems. To solve the highdimensional non-convex problem, we develop an efficient computational algorithm based on exact convex reformulation. Numerical studies are conducted to illustrate the effectiveness of our modeling and computational approach and to draw policy insights. The proposed modeling and computational strategies could provide a solid foundation for the future study of power system resilience with private EVs in coupled transportation and power networks.

Introduction

Power system resilience reflects the versatility of a system to withstand and recover rapidly from unexpected disruptions [START_REF] Khazeiynasab | Resilience analysis and cascading failure modeling of power systems under extreme temperatures[END_REF][START_REF] Lei | Routing and scheduling of mobile power sources for distribution system resilience enhancement[END_REF]. An increasing market penetration of private electric vehicles (EVs) provides new opportunities for enhancing power system resilience due to their mobility and fast regulating characteristics [START_REF] Liu | Ev dispatch control for supplementary frequency regulation considering the expectation of ev owners[END_REF]. In this paper, we will investigate the possible impacts of private EVs on assisting the restoration process of distribution systems (DSs).

Current research on DS restoration with EVs have focused on a centralized perspective without considering the noncooperative travel and charging behavior of private EVs, as well as other decentralized power system stakeholders [START_REF] Haggi | A review of smart grid restoration to enhance cyber-physical system resilience[END_REF]. For example, researchers in [START_REF] Lei | Routing and scheduling of mobile power sources for distribution system resilience enhancement[END_REF][START_REF] Xu | Resilience-oriented distribution system restoration considering mobile emergency resource dispatch in transportation system[END_REF] consider EVs as a type of mobile energy sources (MESs) along with mobile energy storage systems and mobile generators in DS restoration. Authors in [START_REF] Jamborsalamati | Enhancing power grid resilience through an iec61850-based ev-assisted load restoration[END_REF] investigate the centralized real time management of DS restoration with EV aggregators and distributed generators (DGs). In all the mentioned studies, a centralized entity, i.e., distribution system operator (DSO) is in full control of DS and restoration strategies, making decisions for DGs and EVs. However, centralized controlling techniques are proven to be challenging in terms of costly communication infrastructure and single-point failures [START_REF] Nejad | Distributed load restoration in unbalanced active distribution systems[END_REF]. In addition, privately-owned EVs, as well as other energy providers in DSs, may have their own objectives and can not be considered as manageable elements.

Moreover, due to large-scale private EVs participating in Vehicle-to-Grid (V2G) services, power and transportation systems become more interdependent. However, most of the existing studies on power system restoration with EVs take a power-system-centric perspective and completely neglect transportation systems [START_REF] Jamborsalamati | Enhancing power grid resilience through an iec61850-based ev-assisted load restoration[END_REF][START_REF] Momen | Using evs as distributed energy resources for critical load restoration in resilient power distribution systems[END_REF][START_REF] Sun | Optimal distribution system restoration using phevs[END_REF] or simplify the modeling of transportation systems by considering constant predefined travel time on routes [START_REF] Lei | Routing and scheduling of mobile power sources for distribution system resilience enhancement[END_REF][START_REF] Yao | Transportable energy storage for more resilient distribution systems with multiple microgrids[END_REF][START_REF] Yao | Resilient load restoration in microgrids considering mobile energy storage fleets: A deep reinforcement learning approach[END_REF]. This stream of literature overlooks the inherent relationship between traffic distribution and travel time that are critical especially when largescale EVs are available. As a consequence, spatial and temporal interdependence between transportation and power systems can not be properly investigated.

We address the aforementioned gaps by proposing a a framework of network-based multi-agent optimization problems with equilibrium constraints (N-MOPEC) in a coupled distribution and transportation system. The main contribution is two-fold: [START_REF] Khazeiynasab | Resilience analysis and cascading failure modeling of power systems under extreme temperatures[END_REF] The modeling framework captures the decentralized behavior of stakeholders during distribution system restoration process and the spatial and temporal interdependence between transportation and power systems, which allows for rigorous system analyses and optimal V2G incentives design. (2) To facilitate large-scale computation, we reformulate the multi-agent optimization problems as an exact convex optimization problem, which can be efficiently solved by commercial nonlinear solvers.

Mathematical Modeling

We consider the following stakeholders: DSO, EV drivers, charging station aggregator (CSA), and DG owners for distribution system resilience. The decentralized decision making of each stakeholder is modeled in this section. Throughout this paper, set I represents the DS nodes, and I CS , I DG , and I L are subsets of this set representing charging station nodes, DG nodes, and load nodes respectively. The time period of 24 hours are included in set T , and L represents the set of distribution lines. Transportation networks are presented as a directed graph G(N, A), within which EV drivers starting from set R select from charging station set S to provide grid restoration services and/or to charge.

DG Owners Modeling: Each DG i (∈ I DG ) determines its generation quantity p DG i,t for each time step t (∈ T ) to optimize its profits. Because individual DG does not have market power to influence the locational electricity prices ρ i,t , the decision making of all DG owners can be aggregated into a single optimization problem, as formulated in model [START_REF] Khazeiynasab | Resilience analysis and cascading failure modeling of power systems under extreme temperatures[END_REF]. Objective (1a) maximizes the profits of DG owners calculated as the total revenue i∈I DG t∈T ρ i,t p DG i,t subtracting the total production costs i∈I DG t∈T C i (p DG i,t ). We assume C i (•) to be a convex function with respect to p DG i,t [START_REF] Yan | Robust two-stage regional-district scheduling of multi-carrier energy systems with a large penetration of wind power[END_REF]. Constraint (1b) determines the upper and lower bounds ( PDG i,t /P DG i,t ) for power generation at each DG node i for time t. When DGs can be disconnected from the systems, P DG i,t = 0. DSO Modeling: One of the key responsibilities of a DSO after disruptions is to restore services as soon as possible. Given different characteristics of the loads (e.g., hospital, emergency responses) and limited resources available, the DSO may need to prioritize the restoring procedures. We assume that the DSO intends to maximize the importance of loads being served within the system while minimizing the cost of energy purchased, which can be formulated in model [START_REF] Lei | Routing and scheduling of mobile power sources for distribution system resilience enhancement[END_REF].

max p d ,v≥0, p s ,pf ,q f i∈I L t∈T ω i,t p d i,t - i∈I DG ∪I CS t∈T ρ i,t p s i,t (2a)
s.t. l∈L p f l,t • LT l,i - l∈L p f l,t • LF l,i = p d i,t -p s i,t , ∀i ∈ I, t ∈ T (2b) l∈L q f l,t • LT l,i - l∈L q f l,t • LF l,i = q d i,t -q s i,t , ∀i ∈ I, t ∈ T (2c) 0 ≤ p d i,t ≤ Pd i,t , ∀i ∈ I, t ∈ T (2d) q d i,t = ( Qd i,t / Pd i,t ) • p d i,t , ∀i ∈ I, t ∈ T (2e) p f 2 l,t + q f 2 l,t ≤ λ l,t • (S max l ) 2 , ∀l ∈ L, t ∈ T (2f) v FN l ,t -v TN l ,t ≤ (1 -λ l,t ) • K + 2 • (r l • p f l,t + x l • q f l,t ), ∀l ∈ L, t ∈ T (2g) v FN l ,t -v TN l ,t ≥ (λ l,t -1) • K + 2 • (r l • p f l,t + x l • q f l,t ), ∀l ∈ L, t ∈ T (2h) (V min i ) 2 ≤ v i,t ≤ (V max i ) 2 , ∀i ∈ I, t ∈ T (2i)
Objective (2a) maximizes the weighted sum of load demand i∈I l t∈T ω i,t p d i,t , where the weights ω i,t determine the priority of the loads p d i,t ; i∈I s t∈T ρ i,t p s i,t is the cost of energy p s i,t purchased from either DG/EV owners, who want to participate in the restoration services. Operational constraints are formulated in (2b)-(2i), which are adapted based on the Dist-Flow equations proposed in [START_REF] Baran | Network reconfiguration in distribution systems for loss reduction and load balancing[END_REF][START_REF] Guo | Impacts of integrating topology reconfiguration and vehicle-to-grid technologies on distribution system operation[END_REF]. Active and reactive power flow balances are modeled in (2b) and (2c), respectively. p f l,i /q f l,t are active/reactive line flow of distribution line l (∈ L) at time t (∈ T ), and q i,t is the reactive load demand at node i at time step t. LF l,i /LT l,i are mapping matrices where LF l,i /LT l,i equals to 1 if line l starts/connects from/to node i. Constraint (2d) limits the served load demand p d i,t to be less than the expected load demand Pd i,t . Constraint (2e) maintains the same power factor for restored loads as the ratio between expected active ( Pd i,t ) and reactive ( Qd i,t ) load demands. Constraint (2f) models line capacity S max l considering binary line status λ l,t , with 0 indicating line outage. Constraints (2g) and (2h) calculate the squared value of voltage v i,t at each node where FN l and TN l are connection incidence matrices representing the start and end nodes of line l, respectively, and r l /x l are resistance/reactance of line l. K is a big number to ensure that when line l is disconnected at t (i.e., λ l,t = 0), these two constraints are always satisfied. Constraint (2i) defines the acceptable voltage range

[V min i , V max i ].
CSA Modeling: Communicating and controlling EVs/charging stations individually is challenging for DSO. Therefore, we assume a private CSA is responsible to coordinate the power transactions between all charging stations and DS. CSA aims to maximize its profits [START_REF] Duan | Bidding strategies in energy and reserve markets for an aggregator of multiple ev fast charging stations with battery storage[END_REF] while maintaining required charging demand of each EV. The decision making of the CSA is formulated in model [START_REF] Liu | Ev dispatch control for supplementary frequency regulation considering the expectation of ev owners[END_REF].

Objective (3a) maximizes CSA's profits, calculated as revenue of selling electricity to the DSO ( are based on EV type e, which depends on their arriving state of charge (SOC) and charging needs. For example, an EV with higher arriving SOC (soc r,e,T arr r,e ) and lower charging needs should receive higher incentives due to their higher value to the DS restoration. We adapt the Ah-throughput counting model [START_REF] Peterson | Lithium-ion battery cell degradation resulting from realistic vehicle and vehicle-to-grid utilization[END_REF] to model the degradation cost of EVs (C deg,e i,r,t (p e i,r,t )) as a set of linear constraints. An identical strategy is adopted in [START_REF] Guo | Impacts of integrating topology reconfiguration and vehicle-to-grid technologies on distribution system operation[END_REF], which one can refer to for details. Constraints (3b)-(3f) specify SOC requirements that the CSA needs to fulfill. We discretize EVs into different homogeneous groups e based on their travel and charging characteristics, including arriving/departing time (T arr r,e /T EV Drivers Modeling: Travel and charging behavior should be considered in modeling EV drivers. We have modeled the charging behavior of EVs based on their SOC requirements as a part of CSA responsibilities in (3b-3e). In this section, we will model the routing and charging location choices of EVs in transportation system.

The utility function U e rs of a driver in type e selecting station s is formulated in (4) [START_REF] Guo | Impacts of integrating topology reconfiguration and vehicle-to-grid technologies on distribution system operation[END_REF][START_REF] Guo | A stochastic multi-agent optimization framework for interdependent transportation and power system analyses[END_REF]. EV drivers make charging location choices based on four aspects: locational attractiveness β 0,s , travel time -β 1 tt rs , charging cost/revenue from charging/discharging their stored energy β 2 α e rs , and a random term . In this study, we adopt a multinomial logit model for charging location choices, in which follows an extreme value distribution.

U e rs = β 0,s -β 1 tt rs + β 2 α e rs + (4)
The destination choice of EVs (q e rs ) and path travel time (tt rs ) are coupled, see utility function (4). To capture these couplings, we adapt the classic combined distribution and assignment (CDA) model [START_REF] Sheffi | Urban transportation networks[END_REF] to model their destination choices and route choices, as shown in [START_REF] Xu | Resilience-oriented distribution system restoration considering mobile emergency resource dispatch in transportation system[END_REF]. Notice that ∀τ ∈ T arr , we will have an CDA model to solve for the traffic patterns at time τ.

The objective function in (5a) consists of two parts: the first part is the summation of the area under all the link travel cost functions tt a (•) (e.g., Bureau of Public Roads (BPR) function); the second part is the entropy of traffic distribution q e rs (ln q e rs -1) and utility terms (excluding time) in (4). Objective (5a) is constructed in this form to guarantee the optimal solutions of ( 5) are consistent with the first Wardrop principal [START_REF] Wardrop | Some theoretical aspects of road traffic research[END_REF] and the multinomial logit destination choice assumption. For technical details to prove this claim, one can refer to [START_REF] Sheffi | Urban transportation networks[END_REF]. Constraint (5b) calculates link flows v by summing link flows of EVs x τ rs and conventional vehicles xτ rs traveling at time τ. Constraints (5c-5d) are the vehicle flow conservation at each node for the travel demand of EV (q e rs ) and conventional vehicle ( qτ rs ), respectively. A is the nodelink incidence matrix of transportation network, and E rs is the origin-destination (OD) incidence vector. Constraint (5e) guarantees the summation of EV traffic flow distribution to each s equals to the total EV travel demand from r, Q e r . The equilibrium travel time for each OD pair rs can be calculated as tt rs η τ rs,r -η τ rs,s , where η τ rs,n is the dual variable for (5c). q e rs E rs , ∀r ∈ R, s ∈ S, τ ∈ T arr (5c)

A xτ rs = qτ rs E rs , ∀r ∈ R, s ∈ S, τ ∈ T arr (5d) s∈S q e rs = Q e r , ∀r ∈ R, e ∈ E τ (5e)
Market Clearing Conditions: The hourly market clearing conditions can be stated as [START_REF] Jamborsalamati | Enhancing power grid resilience through an iec61850-based ev-assisted load restoration[END_REF]. In a stable market, the power purchased and supplied by DSO needs to be balanced with locational power generation and power load, respectively. (6a) guarantees that the total energy purchased by DSO is equal to the total energy generated at each node. (6b) enforces the balance between EV flow of type e demanded and supplied at each location s from r. Locational prices of electricity ρ i,t and EV incentives α e rs can be calculated from the dual variables for the market clearing conditions, respectively.

(ρ i,t ) p s i,t = p DG i,t + P CS i,t , ∀i ∈ I DG ∪ I CS , ∀t ∈ T (6a)
(α e rs ) q e rs = q e rs , ∀r ∈ R, s ∈ S, e ∈ E (6b)

Convex Reformulation

The decision making of each stakeholder and market clearing conditions presented in Sections 2. are interdependent and need to be solved simultaneously to achieve the equilibrium solutions. However, due to the non-convex nature of the N-MOPEC, directly solving it is challenging. In this section, we propose an exact convex reformulation that can recover the optimal primal and dual variables efficiently. We observe that models (1)∼( 3) and ( 5) are convex optimization problems with constraints completely separable. The objective functions of these models are almost separable except the multiplication terms of primal and dual variables in market clearing conditions [START_REF] Jamborsalamati | Enhancing power grid resilience through an iec61850-based ev-assisted load restoration[END_REF]. This type of problem can be reformulated by linearly combining all the objective functions and constraints and applying the reverse procedures of Lagrangian relaxation on the market clearing conditions (6) [START_REF] Guo | A stochastic multi-agent optimization framework for interdependent transportation and power system analyses[END_REF]. Accordingly, the N-MOPEC can be reformulated as a single convex optimization problem [START_REF] Nejad | Distributed load restoration in unbalanced active distribution systems[END_REF]. max (p d ,p DG ,v,x, x,q,q )≥0, p s ,pf,q f,p,p

Numerical Simulation

We test the proposed models and reformulation techniques with a four-node test system shown in 3), medium (SOC arr r,2 = 0.6), and high (SOC arr r,3 = 0.8). We assume that two traffic flows of 30 EVs/hour, with 10 EVs/hour from each group e, depart from transportation nodes 1 and 4 and arrive at random hours τ (∈ T arr ) to charging stations. We consider losing line (1-2) unexpectedly at t = 10, disconnecting the important DG at node 1 from the system. The line comes back on at t = 20. During this period, generation from DG at node 4 is not enough to serve load demand of all the nodes, and we will considered two levels of SOC dep = 0.7 and 0.5 to investigate EVs participation in DS restoration.

The load pickup graphs presented in Figure 2a, 2d shows that the system is able to provide the energy for the higher priority load at node 3. With lower SOC dep , more lower priority loads are picked up. The total load loss decrease from 2.316 pu to 1.056 pu when SOC dep decrease from 0.7 to 0.5. The power injection of charging stations in Figure 2b, 2e shows more participation of EVs with SOC dep = 0.5 during the disruption compared to SOC dep = 0.7. Energy prices

Discussion

We investigate the impact of private EVs on distribution system restoration in a N-MOPEC framework, where each stakeholder, including DG owners, DSO, EV drivers, and CSA, maximizes their own objectives in a coupled transportation and distribution system. We reformulate the multi-agent problems as an equivalent convex optimization problem, which can be efficiently solved by commercial nonlinear solvers. Numerical results on the small-scale test system show that participation of EVs helps to reduce the load loss during restoration process, and the CSA could provide different incentives to EV drivers based on their value for the distribution system restoration.

This work can be extended in multiple directions. First, the proposed modeling framework can be used to design optimal incentives for EVs to participate in distribution system restoration services and guide the planning and expansion of transportation and power systems. Second, stochastic modeling on renewable energy sources, EV arrival and departure time, and arrival SOC can be further examined. Third, decomposition based algorithms can be developed based on the convex reformulation to facilitate computation for extreme large-scale systems.

Figure 1 .

 1 Nodes 2 and 3 are the connecting points of the distribution and transportation systems representing the location of charging stations. Nodes 2 ∼ 4 are load nodes where loads at node 3 have higher priority (ω 3 = 60) compared to the other nodes (ω 1 = ω 2 = 50). All the distribution lines have the same capacity of S max l = 2 (pu), and transportation links have the capacity of 20 vehicles/hour. EVs are in one of the three groups based on their arrival SOCs: low (SOC arr r,1 = 0.

r 1 = 5 .75 e - 3 x 1 = 2 .93 e - 3 r 2 = 3 .076 e - 2 x 2 = 1 .567 e - 2 r 3 = 2 .284 e - 2 x 3 = 1 .

 15312323221232231 163 e -2

Figure 1 :Figure 2 :

 12 Figure 1: Four node test system

Figure 3 :

 3 Figure 3: Transportation system results: Charging station incentives for EVs (a) SOC dep = 0.7 (b) SOC dep = 0.5. EV traffic flows (c) SOC dep = 0.7 (d) SOC dep = 0.5

  CS i,t can have negative values representing that EVs are charging instead of discharging. In this case, ρ i,t p CS i,t would be the cost for CSA to charge EVs. Power prices ρ i,t and incentives α e rs provided to EV drivers are endogenously determined by market in the proposed modeling framework. Incentives α e rs

i∈I CS t∈T ρ i,t p CS i,t ) subtracting incentives ( r∈R s∈S e∈E α e rs q e rs ) and battery degradation compensation ( i∈I CS e∈E t∈T C deg,e i,r,t (p e i,r,t )) paid to EV drivers. Charging station energy provision p

  ) of the EVs from r. Battery capacity Cap e is considered to normalize the charged/discharged energy of EVs i∈I CS p e i,r,t . Constraint (3c) limits the maximum and minimum SOC (SOC e /SOC e ) for each groups of EVs based on battery specifications, and q e ri(s) is the flow of group e EVs departing from node r that selected charging station i (connected to node s in transportation system). The arrival SOC and minimum departure SOC of EVs are constrained in (3d) and (3e), respectively. Constraint (3f) determines the total power supply/demand (p CS i,t ) of charging station i ∈ I CS at time t by summing the normalized charging/discharging of EVs at each station where S base is the nominal capacity of DS.

									dep r,e ), arriving
	SOCs (SOC arr r,e ), and minimum departing SOC (SOC dep e ). Constraint (3b) models the dynamics of the aggregated SOC
	(soc r,e,t max p CS ,q ,soc	i∈I CS t∈T	ρ i,t p CS i,t -	r∈R s∈S e∈E	α e rs q e ri(s) -	i∈I CS r∈R t∈T e∈E	C deg,e i,r,t (p e i,r,t )	(3a)
	s.t.	soc r,e,t = soc r,e,t-1 -	p e i,r,t /Cap e , ∀r ∈ R, e ∈ E, t ∈ (T arr r,e , T dep r,e ]	(3b)
						i∈I CS
		i∈I CS	q e ri(s) SOC i∈I CS	q e ri(s) SOC e , ∀r ∈ R, e ∈ E, t ∈ [T arr r,e , T dep r,e ]	(3c)
		soc r,e,T arr r,e =	q e ri(s) SOC arr r,e , ∀r ∈ R, e ∈ E	(3d)
				i∈I CS
		soc r,e,T dep r,e	≥	q e ri(s) SOC dep e , ∀r ∈ R, e ∈ E	(3e)
				i∈N CS
		p CS i,t =	p e i,r,t /S
			r∈R e∈E

e ≤ soc r,e,t ≤ base , ∀i ∈ I CS , t ∈ [T arr r,e , T dep r,e ] (3f)