Decentralized charging coordination of electric vehicles using multi-population games

Decentralized Charging Coordination of Electric Vehicles using Multi-Population Games Juan Martinez-Piazuelo, Nicanor Quijano, Senior Member IEEE, Carlos Ocampo-Martinez, Senior Member IEEE Abstract— This paper addresses the decentralized charging coordination of a fleet of plug-in electric vehicles (PEVs). In particular, we cast the charging coordination task as a constrained multi-objective optimization problem, and we solve it using a novel receding horizon decentralized optimization method based on multi-population games. Our proposed method is able to coordinate the charging process of arbitrary fleets of PEVs, while satisfying hard operational constraints over the system’s variables. Our theoretical developments are illustrated through numerical simulations of various PEV-fleets of different sizes. I. INTRODUCTION The current world-wide environmental problems and the augmenting levels of pollution have motivated low-emission transportation systems like plug-in electric vehicles (PEVs). Such transportation systems provide an alternative to gasbased vehicles, and promise a significant reduction on the ecological footprint of the modern industry. However, the raising number of PEVs imposes additional challenges for the managers of the electric networks and for the future smart grids. For instance, it is known that the uncoordinated charge of large fleets of PEVs has non-negligible effects on the distribution grid. In fact, such uncoordinated charging can lead to power losses and overloads on the transformers or feeders of the grid. Therefore, it is paramount to develop coordination algorithms that can deal with the operational constraints of the electric networks and the PEVs. Several recent works have addressed the charging coordination problem both from centralized and decentralized perspectives. The centralized approaches typically assume that a central coordinator defines the charging schedule of each PEV, or that the information of each PEV is broadcast to all the other PEVs within the same distribution network. However, optimizing over large fleets of PEVs in a centralized way might lead to high computational costs, and, in practice, the PEVs might not be willing to share their information with other PEVs. For such reasons, the decentralized approaches have become more popular within the research community. For instance, the authors in [1] and [2] develop optimal decentralized valley-filling algorithms to fill the overnight demand valley. Alternatively, the authors in [3] and [4] propose decentralized coordination methods based on mean-field games and aggregative games, respectively, where Juan Martinez-Piazuelo and Nicanor Quijano are with the Department of Electrical and Electronics Engineering at Universidad de los Andes, Bogot´ a, Colombia. Carlos Ocampo-Martinez is with the Automatic Control Department at Universitat Polit` ecnica de Catalunya, Barcelona, Spain. Corresponding author: Juan Martinez-Piazuelo ([email protected]) the goal is to minimize the charging and battery degradation costs of each PEV within a large fleet. On the other hand, the authors in [5] develop a decentralized dynamic pricing method to coordinate the charging of the PEVs with some renewable energy generation. Finally, in [6], a decentralized gradient projection method is used to shave the peak of the demand of a fleet of PEVs, and hard constraints on the capacity of the feeder lines are considered. In contrast with the previous works, in this paper the coordination goal is to allocate the charging-load of each PEV to minimize the collective deviation from a given reference demand profile, as well as the deviation from the desired level of charge of each PEV, i.e., a multi-objective minimization problem. The reference demand profile might reflect for instance a peak shaving objective, or might be related to some estimated renewable energy generation or to an energy profile bought in a day-ahead market. Hence, most of the objectives of the aforementioned works can be considered under our coordination framework. Furthermore, similar to [6], in this work we consider hard constraints on the power capacity of the feeder lines, and we also consider some additional hard constraints related to the charging stations and the PEVs. Moreover, as in [3], we consider an aggregator-based communication scheme to reduce the communication costs of the decentralized approach. Consequently, the main contribution of this paper is the development of a novel receding horizon decentralized optimization method, based on the theory of multi-population games [7], which we use to solve the proposed coordination problem. Unlike other available works on multi-population games, in this paper we consider the existence of hard coupledconstraints over the multiple populations. Furthermore, we provide sufficient conditions to guarantee the asymptotic convergence of our proposed method under asynchronous updates of the PEVs, which depict a robust property of our decentralized algorithm. Lastly, all our theoretical results are illustrated through numerical experiments that emulate the charging coordination of PEV-fleets with different sizes. The remainder of this paper is organized as follows. First, Section II states the coordination goal as a decentralized constrained multi-objective optimization problem. Next, Section III presents our proposed method and provides our theoretical results1. Afterwards, Section IV depicts the application of the 1Due to space limitations, the complete proofs of the theoretical results can be found at: https://drive.google.com/drive/folders/ 1vNATfiJdr9ZmPRyA8qdHb1ThA6Yt4wMH?usp=sharing developed theory to some numerical experiments. Finally, Section V concludes the paper. II. PROBLEM STATEMENT AND SYSTEM DESCRIPTION In this paper, we consider a microgrid with multiple charging stations for PEVs, which seeks to coordinate the energy demands of its stations subject to some operational constraints. The goal of the coordination consists of two parts: i) to follow a reference energy profile provided by the manager of the microgrid; and ii) to charge the PEVs that are connected to the charging stations. Clearly, such coordination is subject to the energy and power constraints of the microgrid, its charging stations, and the connected PEVs. In consequence, the aforementioned coordination goal is cast as a constrained multi-objective optimization problem that can be solved in a receding horizon fashion. To define such optimization problem more formally, let T∞={t0, t1, t2,...,∞} be a discretization of the continuous time, where t0∈R≥0denotes the first operational time slot of the system, and ti∈R>0is given by ti=ti−1+∆t, for all i∈Z>0. Here, R≥0(R>0) is the set of non-negative (positive) real numbers; Z>0is the set of positive integers; and ∆t∈R>0is the (fixed) duration of each time slot of T∞, e.g., ∆t= 0.25 hours. Moreover, let t∈ T∞ denote the current operational time slot of the system, and let Tt={t, t + ∆t, t + 2∆t, . . . , t +Tt∆t}be the set of receding horizon time slots from time slot t, where Tt∈Z>0 is the optimization horizon for the optimization performed during the time slot t. At the beginning of every time slot t∈ T∞, the manager of the microgrid sets a reference profile dt= [dt, dt+∆t, . . . , dt+Tt∆t]>, where diis the desired energy (kWh) to be consumed by the entire microgrid during the i-th time slot (we assume that di∈R≥0, for all i∈ Tt). Since the microgrid is composed of a set of charging stations for PEVs, we let S={s0, s1· · · , sS}denote such set, where S∈R>0refers to the total number of stations. Throughout, we assume that each station s∈ S can charge at most one PEV at a time. Besides, when a PEV arrives at station s, it connects to the charging port and specifies its expected leaving time slot t(s) L∈ T∞, i.e., the PEV informs the station sof the future time slot, t(s) L> t, at which it expects to leave the station. To represent this information, we use the vector δ(s) t=hδ(s) t, δ(s) t+∆t, . . . , δ(s) t+Tt∆ti> , for all s∈ S, where δ(s) i= 1 if there is a PEV expected to be connected to station sat time slot i∈ Tt, and δ(s) i= 0 otherwise. Thus, if t(s) L> t +Tt∆t, the connected PEV is expected to be connected to station sduring all the time slots in Tt, and, hence δ(s) i= 1, for all i∈ Tt. In contrast, if no PEV is connected to station sat time slot t, then δ(s) i= 0, for all i∈ Tt. Furthermore, if a PEV remains connected beyond its expected departure time, i.e., t(s) L≤t, then it is assumed that δ(s) t= 1 and δ(s) i= 0 for all i>t. Therefore, for any current operational time slot t∈ T∞, we denote the set of active stations as St=ns∈ S δ(s) t= 1o. Additionally, we assume that each station is capable of measuring the state of charge (SOC) of the connected PEV. Thus, any active station is able to determine the maximum energy (kWh) that can be injected to the connected PEV to achieve a full SOC. Consequently, for a PEV connected to station s∈ Stat the current time slot t∈ T∞, we denote the remaining energy capacity as e(s) t∈R≥0, and we set e(s) t= 0 by default if no PEV is connected to station s. Under the aforementioned formulations, at any current time slot t∈ T∞, the coordination task is to schedule the power (kW) to be injected by each station s∈ Stat every time slot i∈ Tt. The goal of such scheduling is to minimize the collective energy deviation from the reference profile dt, as well as the energy deviation from the remaining energy capacity e(s) t, for all s∈ St. Clearly, such minimization is subject to the operational constraints of the entire system. Letting x(s) i∈R≥0denote the power (kW) that is scheduled to be injected by the station s∈ Stat the time slot i∈ Tt, and setting x(s)=hx(s) t, x(s) t+∆t,· · · , x(s) t+Tt∆ti> , for all s∈ St, and x=x(s0)>,x(s1)>,· · · ,x(s|St|)>> , the coordination problem at any current time slot t∈ T∞ can be cast as the receding horizon minimization given by min x g(x) := g1(x) + ωg2(x)(1a) with g1(x) = 1 2X i∈Tt X s∈St x(s) i∆t−di!2 , g2(x) = 1 2X s∈St X i∈Tt x(s) i∆t−e(s) t!2 , subject to x(s) i= 0,if δ(s) i= 0,∀i∈ Tt,∀s∈ St(1b) X i∈Tt x(s) i∆t≤e(s) t,∀s∈ St(1c) X s∈St x(s) i≤xi,max,∀i∈ Tt(1d) 0≤x(s) i≤x(s) lim ,∀i∈ Tt,∀s∈ St.(1e) Here, note that the cost function (1a) considers two quadratic deviations. Namely, g1(x)regards the collective deviation from the reference profile dt;g2(x)considers the total deviation from the remaining energy capacities of the PEVs; and ω∈R≥0is a weighting parameter to adjust the relative importance of each of these deviations. Throughout, we assume that ωis set by the manager of the microgrid and then we do not design it. Next, note that the constraint (1b) states that no power can be allocated to station s∈ St, at time slot i∈ Tt, if no PEV is expected to be connected to that station at that time slot. In contrast, the constraint (1c) establishes that no more than e(s) tenergy can be injected to the PEV connected to station s∈ St, and the constraint (1d) states that no more than xi,max ∈R>0power (kW) can be collectively injected by all stations at time slot i∈ Tt. Finally, the constraints in (1e) state that no power can be drawn from any PEV, i.e., we do not consider vehicle-togrid technology, and that the maximum power that can be ... Aggregator Charging station Charging station Charging station Charging station No PEV connected Fig. 1: Communication architecture of the microgrid. injected at any time by station s∈ Stis x(s) lim ∈R>0. In addition to the power and energy constraints defined in (1b)-(1e), we also consider some information-related constraints that depend on the system’s communication architecture. In this work, we assume that each station establishes its own power schedule in a decentralized fashion and under partial information. For such, we assume that the microgrid has the communication structure depicted in Fig. 1. In particular, we assume the existence of an aggregator node that can communicate with each charging station using bidirectional communication links, and we do not consider communication links between the stations. Due to privacy constraints, we assume that the information that the aggregator sends to the stations contains only aggregated variables of the entire system, e.g., aggregated power profiles. Hence, the individual information of each charging station is not explicitly disclosed over the microgrid. To design a decentralized algorithm that can deal with all the operational constraints of the system, in this paper we rely on the theory of multi-population games [7]. The motivation being that multi-population games offer a framework that allows the consideration of hard coupled constraints under decentralized decision-making architectures. In the following section we present our approach. III. DECENTRALIZED OPTIMIZATION USING MULTI-POPULATION GAMES A. A multi-population game approach For any current time slot t∈ T∞, consider a society composed of a set of populations, St, that represents the set of active stations of the microgrid, and let each population s∈ Stbe composed of a large, finite, and constant number of agents that interact strategically with their population peers. Thus, each population s∈ Stis characterized by a constant mass of agents that is denoted as m(s) t∈R>0. Note that the assumption of a constant mass m(s) timplies that there is no mortality of agents and that agents cannot leave their population. At the beginning of the current time slot t, the agents of each population engage in an iterative interaction process that lasts the whole duration of the time slot t. Such iterative interaction process is characterized by the strategic decisions that the agents make. More precisely, at any iteration k∈Z≥0, each agent within population s∈ St is allowed to select a strategy from a finite set of strategies Vt, which is assumed equal for all s∈ St. Moreover, the goal of each agent is to select the strategy that leads to the highest possible payoff within the strategic framework. In this work, at any time slot t∈ T∞, the set of strategies is given by Vt=Tt∪ {`}, which is the union of the set of receding horizon time slots, Tt, and an additional fictitious strategy, `, whose meaning will be explained shortly. Under such definition, at any iteration k, the mass of agents of population s∈ Stthat selects the strategy i∈ Vtis denoted as x(s) i[k]∈R≥0, and can be interpreted as the value of the optimization variable x(s) iat the optimization iteration k, i.e., such portion of agents is interpreted as the power (kW) that is scheduled to be injected by the station s∈ Stat the time slot i∈ Tt. The goal is thus to design the interaction game such that the population portions x(s) i[k], for all i∈ Vtand all s∈ St, converge to a minimizer of (1) as k→ ∞, i.e., convergence in the asymptotic sense. In order to do so, we establish that for all i∈ Vt, all s∈ St, and all k∈Z≥0, the dynamical evolution of x(s) i[k]is determined by x(s) i[k+1] = x(s) i[k]+(s) t[k]X j∈Vthρ(s) ji x(s) j[k]−ρ(s) ij x(s) i[k]i, (2) where (s) t[k]∈R≥0is the step size for the updates of population s∈ St, and ρ(s) ij ∈R≥0is a revision protocol [7] that determines the incentive that is given to an agent of population s∈ Stto switch from strategy i∈ Vtto strategy j∈ Vt. Notice that to solve (1) using iterative updates of the form (2), the revision protocols ρ(s) ij must reflect the cost function (1a), and must consider the optimization constraints (1b)-(1e) as well as the informationrelated constraints given by the communication structure of the microgrid. We discuss the design of such revision protocols in Section III-B. However, before doing so, let us explain the meaning of the fictitious strategy `. Note that, under the aforementioned interpretations, the population portions selecting the fictitious strategy `can be thought as slack variables of the optimization problem, i.e., x(s) `do not appear in (1) for any s∈ St. Moreover, notice that from the assumption that m(s) tis constant for all s∈ St, at any iteration k, it holds that Pi∈Vtx(s) i[k] = m(s) t, and, in consequence, it holds that Pi∈Ttx(s) i[k]≤m(s) t. Thus, if we let m(s) t=e(s) t∆t−1, for all s∈ St, the fictitious strategy `allows us to deal with the inequality constraints in (1c). Remark 1: Since (s) t[k]belongs to R≥0, we have that the step sizes can be zero. These properties allow us to consider heterogeneous step sizes and asynchronous updates over the different populations, which are desirable traits for distributed algorithms. To model such cases, note that (s) t[k] can always be expressed as (s) t[k] = t[k]φ(s)[k], for all s∈ St, where t[k]∈R≥0is a global upper bound on the possible step sizes, and φ(s)[k]∈R≥0is a scaling factor for the s-th population and satisfies 0≤φ(s)[k]≤1, for all k. Therefore, 0≤(s) t[k]≤t[k], for all s∈ St. B. Discrete-time coupled Smith dynamics Several revision protocols have been proposed in the literature and their application to (2) have led to different population dynamics [7], [8], [9], [10]. Inspired by [10], in this section we propose a novel revision protocol considers all the constraints of the optimization problem (1), and generates some novel multi-population dynamics, here referred to as the discrete-time coupled Smith dynamics (DT-CSDs). More formally, we propose the revision protocol ρ(s) ij =δ(s) jα(s) jβjhf(s) j−f(s) ii+,∀i, j ∈ Vt,∀s∈ St, (3) with [·]+ . = max{·,0}, and where α(s) j="1−x(s) j[k] x(s) lim #+ ,∀j∈ Tt,∀s∈ St,∀k(4a) βj="1−1 xj,max X z∈St x(z) j[k]#+ ,∀j∈ Tt,∀k(4b) δ(s) `=α(s) `=β`= 1,∀s∈ St(4c) f(s) i=−∂g1(x[k]) ∂x(s) i −ω∂g2(x[k]) ∂x(s) i ,∀i∈ Vt,∀s∈ St. (4d) First, note that the term δ(s) jis included in (3) to consider the constraints (1b). Note that if a PEV is not expected to be connected to station s∈ Stat the time slot j∈ Vt, then the agents in population s∈ Sthave no incentive to select such strategy. Second, the term (4a) is used to consider the upper bound in the constraints (1e). Notice that if a strategy j∈ Vthas reached its maximum capacity within its population, i.e., x(s) lim ∈R>0, then no agent in population s∈ Sthas incentives to select the strategy j∈ Vt. Similarly, the term (4b) is included to deal with the constraint (1d). In this case, if the strategy j∈ Vthas reached its maximum coupled capacity over all populations, which is xj,max ∈R>0, then no agent in the entire society has incentives to select the strategy j∈ Vt. Next, given that the fictitious strategy, `, is not related to an actual time slot, we do not constrain the upper magnitude of the mass of agents that select it, and we assume that it can always be selected regardless of the PEV connection status. Thus, we set the corresponding values as in (4c). Finally, the term f(s) iin (4d) is used to link the revision protocols to the cost function (1a). The term f(s) iis the fitness or payoff of the strategy i∈ Vtwithin population s∈ St, and, as shown in (4d), it depends on the overall society state x[k]. Here, resembling the definitions presented in Section II, x[k] = x(s)[k]∈Rnt ≥0, where nt=|Vt||St| and, for all s∈ St,x(s)[k] = hx(s) i[k]i∈R|Vt| ≥0. Hence, for all k∈Z≥0, the vector x[k]determines the distribution of the agents of the entire society over the set of strategies Vt. To derive the DT-CSDs, we replace the proposed revision protocol (3) into (2). By doing so, we obtain x(s) i[k+ 1] = x(s) i[k] + (s) t[k]X j∈Vt θ(s) ij f(s) i−f(s) j,(5) for all i∈ Vtand all s∈ St, and with θ(s) ij =         ϕ(s) ix(s) j[k],if f(s) i> f(s) j, ϕ(s) jx(s) i[k],if f(s) i< f(s) j, 1 2ϕ(s) ix(s) j+1 2ϕ(s) jx(s) i,if f(s) i=f(s) j, (6) with ϕ(s) i:= δ(s) iα(s) iβi, for all i, j ∈ Vtand all s∈ St. Moreover, notice that the DT-CSDs for the entire society can be written as x[k+ 1] = x[k] + t[k]Lf ,(7) where t[k]is the global upper bound for the step sizes defined (c.f., Remark 1); f=f(s)∈Rntis the society fitness vector, with f(s)=hf(s) t, f(s) t+∆t, . . . , f(s) t+Tt∆t, f(s) `i> , for all s∈ St; and L∈Rnt×ntis a block-diagonal matrix where the block (s, s)is given by φ(s)[k]L(s), for all s∈ St, and the matrix L(s)∈R|Vt|×|Vt|has elements l(s) ii =Pj∈Vt,j6=iθ(s) ij and l(s) ij =−θ(s) ij , for all i, j ∈ Vtwith i6=j, and all s∈ St. Now that we have presented the proposed dynamics, we proceed to show that the DT-CSDs can indeed be applied to solve the constrained optimization problem (1) in a decentralized receding horizon fashion. However, for the sake of clarity, let us first explain the steps of our proposed algorithm. The following steps are repeated for all t∈ T∞: 1) Before the beginning of the current time slot t:The aggregator broadcasts a signal to all s∈ S requesting the expected PEV connection status for time slot t. Then, each station s∈ S sends its corresponding status back to the aggregator, and the aggregator gathers such signals and updates the set of active stations St. 2) At the beginning of the current time slot t:Assuming that there is at least one active station, the aggregator computes the vector z∈R|Vt|, whose elements are z`= 0 and zi=xi,max/|St|, for all i∈ Tt. Afterwards, the aggregator broadcasts z,dt, and ω, to all s∈ St. Then, with the received vector, z, each station s∈ St computes its own initial population state as follows: x(s) i[0] = δ(s) imin (zi, γ(s)x(s) lim ,e(s) t ∆t|Vt|),∀i∈ Tt x(s) `[0] = e(s) t ∆t−X i∈Tt x(s) i[0], where γ(s)∈R>0satisfies 0< γ(s)<1, for all s∈ St. Note that such initial society state, x[0], satisfies (1b)-(1e) without binding on (1c)-(1e). Next, each station s∈ Stapplies the power x(s) t[0] to the connected PEV, and sends its corresponding x(s)[0] to the aggregator. Finally, the aggregator computes and broadcasts the aggregated power profile for all i∈ Tt, i.e., Ps∈Stx(s) i[0] for all i∈ Tt. 3) During the current time slot t:The following procedure is repeated until the termination of the time slot t. At each optimization iteration k, each station s∈ St computes ˆ x(s)[k] = L(s)f(s), and sends the vectors x(s)[k],ˆ x(s)[k], and f(s)to the aggregator. Then, using the received information, the aggregator computes and broadcasts the aggregated power profile for all i∈ Tt, i.e., Ps∈Stx(s) i[k]for all i∈ Tt, as well as a global bound for (s) t[k]for all s∈ St. Such global bound is computed according to Theorems 1 and 2, which are presented in Sections III-C and III-D, respectively. Finally, each station applies the charging power x(s) t[k] to the connected PEV, and computes x(s)[k+ 1] using (5) for all i∈ Vt. C. Information-related constraints and invariant set analysis In this section, we show that the proposed population dynamics satisfy the information-related constraints of the considered system, and we provide sufficient conditions to guarantee the satisfaction of all the constraints in (1b)-(1e). Let us begin with the information-related constraints. First, note that from (4d) and (1a) we have f(s) i=di+ωe(s) t−X z∈St x(z) i[k]∆t−ωX i∈Tt x(s) i[k]∆t, f(s) `= 0, for all i∈ Ttand all s∈ St. Clearly, with the aid of the aggregator node, f(s)can be computed at all s∈ St. Moreover, note that besides system-level terms as dtand ω, the only information that the aggregator broadcasts is the aggregated power profile Pz∈Stx(z) i[k], for all i∈ Tt. Hence, the information provided by the aggregator does not disclose individual information of any station. Second, note that the only other term in (3) that requires system-level information is βj, for all j∈ Tt. However, in this case the required information is again the aggregated power profile provided by the aggregator. Thus, the proposed DT-CSDs do preserve the information-related constraints of the system. To analyze the remaining constraints, consider the feasible set Ft=nx[k]∈Rnt ≥0x[k]satisfies (1b)-(1e)o, and its relative interior defined as int (Ft) = x[k]∈Rnt >0 x[k]satisfies (1b)-(1e) and (1c)-(1e) are not binding . (8) Using these definitions, the following results provide sufficient conditions to guarantee the forward invariance of int (Ft)under the DT-CSDs. That is, x[k]∈int (Ft)implies that x[k+ 1] ∈int (Ft), for all k. Assumption 1: All the population portions of the same population are updated at the same time. Lemma 1: Consider a population s∈ Stthat updates its population portions using (5), and let Assumption 1 hold. Then, m(s) t=Pi∈Vtx(s) i[0] = Pi∈Vtx(s) i[k], for all k≥0. Theorem 1: Let Assumption 1 hold, and let the step size of each population s∈ Stbe bounded by 0≤(s) t[k]<min nt,c[k], (s) t,d [k], (s) t,p [k]o,∀k, ∀t∈ T , (9) with t,c[k] = min i∈Tt     xi,max −Ps∈Stx(s) i[k] max Ps∈Sthˆx(s) i[k]i+, η    ,(10a) (s) t,d [k] = min i∈Tt  x(s) lim −x(s) i[k] max nˆx(s) i[k], ηo ,∀s∈ St,(10b) (s) t,p [k] = min i∈sppδ(s) t∪{`}  x(s) i[k] max n−ˆx(s) i[k], ηo ,∀s∈ St, (10c) where ˆx(s) i[k] = Pj∈Vtθ(s) ij f(s) i−f(s) j, for all i, j ∈ Vt, all s∈ St, and all k∈Z≥0;η∈R>0is a small positive constant that is included to avoid divisions over zero; and spp δ(s) tdenotes the support of δ(s) t, for all s∈ St. Then, int (Ft)is forward invariant under the DT-CSDs of (5). D. Convergence analysis In this section, we present the convergence analysis of the proposed dynamics and derive sufficient conditions to guarantee asymptotic convergence to a minimizer of (1). Lemma 2: Let x[k]∈ Ft. Then, L0at iteration k. Theorem 2: Consider the optimization problem (1) and the DT-CSDs (7). Let X∗⊂ Ftdenote the set of minimizers of (1), i.e., X∗={x∗∈ Ft|g(x∗)≤g(x),∀x∈ Ft}. Moreover, suppose that Assumption 1 holds, and let the following conditions be satisfied: 1) The initial society state belongs to the interior of the feasible set Ft, i.e., x[0] ∈int (Ft). 2) It holds that X∗∩int (Ft)6=∅, where ∅is the empty set, i.e., there is at least one interior minimizer of (1). 3) For every population s∈ St, the step size (s) t[k]is bounded by 0≤(s) t[k]≤min (s) thm1[k],f>Lf max {f>LHLf , η}, (11) for all kand all t∈ T . Here, (s) thm1[k]∈R>0is any positive step size that satisfies (9); H∈Rnt×ntis the Hessian of (1a); and η∈R>0is a small positive constant that prevents the division over zero. 4) As k→ ∞, it holds that (s) t[k]6= 0 infinitely often with probability one, for all s∈ Stand all x[k]/∈ X∗. Then, the DT-CSDs converge, in the asymptotic sense and with probability one, to an interior minimizer of (1). Remark 2: Due to the forward invariance of int (Ft), the convergence of (7) to an optimal boundary solution of (1) is harder to analyze, and so we leave it for a future work. Notice, however, that regardless of the location of the optimal solutions, the (constrained) asymptotic convergence to an interior minimizer is guaranteed for any number of charging stations. The convergence rate, on the other hand, might indeed depend on the number of active stations. Nonetheless, we leave such research for a future work as well. IV. NUMERICAL EXPERIMENTS To illustrate our developments, we consider four microgrids with different number of stations and power constraints. More precisely, we let S= 100,S= 50,S= 10, and S= 10, for experiments 1,2,3, and 4, respectively. For all experiments we let x(s) lim = 3.7kW, for all s∈ S, and we set xi,max = 2.96SkW, for all i∈ Tt, for experiments 1-3, and xi,max = 0.74SkW, for all i∈ Tt, for experiment 4. The constraints of experiments 1-2are selected so that the optimal solutions lie within int (Ft), while the constraints of the experiments 3-4lead to boundary solutions. Besides, the constraints of the experiment 4are set to make the solution of experiment 3unfeasible. Without loss of generality, the arrival and departure time slots of each PEV are sampled from the normal distributions N18:00,(4∆t)2 and N6:00,(4∆t)2, respectively. Here, ∆t= 0.25 hours. Moreover, the actual departure time slots are sampled from a normal distribution centered at the expected departure time slot and with a variance of ∆t2. Hence, there is uncertainty on the expected departure time slots. As illustration, Fig. 2.a depicts the number of active stations over time. For all experiments, the optimization horizon is Tt= 10 hours, for all t∈ T∞, and the reference energy profile is the one shown in Fig. 2.b. Furthermore, the initial remaining energy capacity of each PEV is sampled from N10kWh,(5kWh)2, and we set ω= 10. In addition, we assume that each charging station has a 97% charging efficiency, i.e., 3% of the injected power is lost. Finally, for all experiments we set η= 10−10, γ(s)= 0.9, for all s∈ S, and we assume that the charging stations update their optimization variables every three seconds with an independent probability of 0.7, i.e., each time slot has 300 optimization iterations and not all stations update at all iterations. As illustration, Fig. 2.b depicts the overall energy consumption of the microgrids over the considered time window. In particular, note that, in all cases, the peak demands are shifted towards the time slots with the higher reference profile. Moreover, for experiment 1the fleet of PEVs is charged up to 87.6% of its maximum collective energy capacity, while for experiments 2,3, and 4, such values are 96.3%,98.8%, and 94.4%, respectively. Furthermore, using a centralized optimization method, it is verified that in all the presented cases the power allocation at the end of each time slot lies within a 0.004% gap of an optimal solution. V. CONCLUDING REMARKS This work has proposed a novel method for decentralized constrained optimization and has illustrated its application to the charging coordination of multiple PEVs. The proposed method can deal with several hard constraints on the optimization variables, and its asymptotic convergence to optimality is guaranteed under certain assumptions. REFERENCES [1] H. Ito, “Disturbance and delay robustness guarantees of gradient systems based on static noncooperative games with an application to feedback control for pev charging load allocation,” IEEE Transactions on Control Systems Technology, vol. 21, pp. 1374–1385, August 2013. 17 18 19 20 21 22 23 0 1 2 3 4 5 6 7 Simulation time (from 17:00 to 7:00) 10 50 100 Number of active stations Exp. 1 Exp. 2 Exps. 3-4 (a) 17 18 19 20 21 22 23 0 1 2 3 4 5 6 7 Simulation time (from 17:00 to 7:00) 0 5 10 15 20 25 Collective energy consumption (kWh) Exp. 1 Exp. 2 Exp. 3 Exp. 4 Reference profile (b) Fig. 2: Numerical experiments over the window 17:00-7:00. (a) Number of active stations. (b) Collective energy profiles. [2] H. Xing, M. Fu, Z. Lin, and Y. Mou, “Decentralized optimal scheduling for charging and discharging of plug-in electric vehicles in smart grids,” IEEE Transactions on Power Systems, vol. 31, pp. 4118–4127, Sep. 2016. [3] M. A. Tajeddini and H. Kebriaei, “A mean-field game method for decentralized charging coordination of a large population of plug-in electric vehicles,” IEEE Systems Journal, vol. 13, pp. 854–863, March 2019. [4] S. Grammatico, “Exponentially convergent decentralized charging control for large populations of plug-in electric vehicles,” in Proceedings of the 55th IEEE Conference on Decision and Control (CDC), pp. 5775–5780, December 2016. [5] Y. Yang, Q. Jia, X. Guan, X. Zhang, Z. Qiu, and G. Deconinck, “Decentralized EV-based charging optimization with building integrated wind energy,” IEEE Transactions on Automation Science and Engineering, vol. 16, pp. 1002–1017, July 2019. [6] Z. Ma, N. Yang, S. Zou, and Y. Shao, “Charging coordination of plugin electric vehicles in distribution networks with capacity constrained feeder lines,” IEEE Transactions on Control Systems Technology, vol. 26, pp. 1917–1924, Sep. 2018. [7] W. H. Sandholm, Population games and evolutionary dynamics. MIT press, 2010. [8] N. Quijano, C. Ocampo-Martinez, J. Barreiro-Gomez, G. Obando, A. Pantoja, and E. Mojica-Nava, “The role of population games and evolutionary dynamics in distributed control systems: The advantages of evolutionary game theory,” IEEE Control Systems Magazine, vol. 37, pp. 70–97, Feb 2017. [9] J. Barreiro-Gomez, G. Obando, and N. Quijano, “Distributed population dynamics: Optimization and control applications,” IEEE Transactions on Systems, Man, and Cybernetics: Systems, vol. 47, pp. 304– 314, Feb 2017. [10] J. Barreiro-Gomez, F. D¨ orfler, and H. Tembine, “Distributed robust population games with applications to optimal frequency control in power systems,” in Proceedings of the 2018 American Control Conference (ACC), pp. 5762–5767, June 2018.