Chance-Constrained Shrunken-Primal-Dual Subgradient (CC-SPDS) Approach for Decentralized Electric Vehicle Charging Control

TL;DR

This paper introduces a decentralized EV charging control framework using a novel chance-constrained shrunken-primal-dual subgradient algorithm to optimize valley-filling while respecting network constraints and individual needs.

Contribution

It presents a new CC-SPDS algorithm for decentralized control under stochastic constraints, ensuring convergence and effectiveness in distribution network management.

Findings

The CC-SPDS algorithm effectively manages EV charging with chance constraints.

The method achieves valley-filling while satisfying network and individual constraints.

Convergence of the proposed algorithm is verified through simulations.

Abstract

In this paper, we develop a chance-constrained decentralized electric vehicle (EV) charging control framework to achieve "valley-filling" meanwhile meeting individual charging requirements and satisfying distribution network constraints. The control design is formulated as an optimization problem with a stochastic non-separable objective function and globally coupled chance constraints. We propose a novel chance-constrained shrunken-primal-dual subgradient (CC-SPDS) algorithm to support the chance-constrained decentralized control scheme and verify its efficacy and convergence with a representative distribution network model.

Equations

$$x_{ı,\hat{ı}}(T+1)=x_{ı,\hat{ı}}(T)+B_{ı,\hat{ı}}{u}_{ı,\hat{ı}}(T),$$

$$x(T+1)=x(T)+\sum _{ı=1}^h\sum _{\hat{ı}=1}^{n_{ı}}{B}_{ı,\hat{ı}}^c{u}_{ı,\hat{ı}}(T),$$

$$\mathcal{X}(k)$$

$${\mathcal{U}}_{ı,\hat{ı}}(k)$$

$$x(k+K|k)=x(k)+\sum _{ı=1}^h\sum _{\hat{ı}=1}^{n_{ı}}{\mathcal{B}}_{ı,\hat{ı}}^l{\mathcal{U}}_{ı,\hat{ı}}(k)=0,$$

$$\mathbf{V}(T)={\mathbf{V}}_0-2\mathbf{R}p(T)-2\mathbf{X}q(T),$$

$$\mathbf{V}(T)$$

$${\mathbf{V}}_0$$

$$p(T)$$

$$q(T)$$

$$\mathbf{R}\in {\mathbb{R}}^{h\times h},{\mathbf{R}}_{ıȷ}=\sum _{(\hat{ı},\hat{ȷ})\in {\mathbb{S}}_{ı}\cap {\mathbb{S}}_{ȷ}}{r}_{\hat{ı}\hat{ȷ}},$$

$$\mathbf{X}\in {\mathbb{R}}^{h\times h},{\mathbf{X}}_{ıȷ}=\sum _{(\hat{ı},\hat{ȷ})\in {\mathbb{S}}_{ı}\cap {\mathbb{S}}_{ȷ}}{x}_{\hat{ı}\hat{ȷ}},$$

$$p_{ı}(T)$$

$$q_{ı}(T)$$

$${\mathbf{V}}_b(T)=2\mathbf{R}G{p}^b(T)+2\mathbf{X}G{q}^b(T)=\mathbf{H}{p}^b(T),$$

$$\mathbf{V}(T)={\mathbf{V}}_0-{\mathbf{V}}_b(T)-2\mathbf{R}G\bar{P}u(T),$$

$$y(T)=y_d(T)+Du(T),$$

$${\mathcal{Y}}_k={\mathcal{Y}}_{dk}+\sum _{ı=1}^h\sum _{\hat{ı}=1}^{n_{ı}}{\mathcal{D}}_{ı,\hat{ı}}{\mathcal{U}}_{ı,\hat{ı}}(k),$$

$${\mathcal{Y}}_k$$

$${\mathcal{Y}}_{dk}$$

$${\mathcal{D}}_{ı,\hat{ı}}$$

$$\mathcal{Y}_{dk}{\sim}\mathcal{N}\left(\left[\begin{array}[]{c}\boldsymbol{V}_{0}-\tilde{\boldsymbol{\mu}}^{p}(k)\\ \vdots\\ \boldsymbol{V}_{0}-\tilde{\boldsymbol{\mu}}^{p}(k+K-1)\end{array}\right],\left[\begin{array}[]{ccc}\tilde{\boldsymbol{\Sigma}}^{p}&{}\hfil&{}\hfil\\ {}\hfil&\ddots&{}\hfil\\ {}\hfil&{}\hfil&\tilde{\boldsymbol{\Sigma}}^{p}\end{array}\right]\right).$$

$$\underset{\mathcal{U}}{\min }$$

$${\mathcal{U}}_{ı,\hat{ı}}\in {\mathbb{U}}_{ı,\hat{ı}},\forall ı=1,\dots ,h,\hat{ı}=1,\dots ,n_{ı},$$

$$d(\mathcal{U})\triangleq {\underline{\nu }}^2{\hat{\mathbf{V}}}_0-{\mathcal{Y}}_{dk}-\sum _{ı=1}^h\sum _{\hat{ı}=1}^{n_{ı}}{\mathcal{D}}_{ı,\hat{ı}}{\mathcal{U}}_{ı,\hat{ı}}\le 0,$$

$$\mathcal{U}$$

$${\mathbb{U}}_{ı,\hat{ı}}$$

$$\underset{\mathcal{U}}{\min }$$

$${\mathcal{U}}_{ı,\hat{ı}}\in {\mathbb{U}}_{ı,\hat{ı}},\forall ı=1,\dots ,h,\hat{ı}=1,\dots ,n_{ı},$$

$$\text{Pr}\left([z{]}_{\kappa }\ge [\hat{\mathcal{Y}}{]}_{\kappa }\right)\ge \delta ,\forall \kappa =1,\dots ,K,$$

$$d_{\kappa }(\mathcal{U})=-F_{[\hat{\mathcal{Y}}{]}_{\kappa }}([z{]}_{\kappa })+\delta \le 0,\forall \kappa =1,\dots ,K,$$

$$\mathcal{L}(\mathcal{U},\mathbf{\lambda })=\mathbb{E}\left[\mathcal{F}(\mathcal{U})\right]+\sum _{\kappa =1}^K{\lambda }_{\kappa }{d}_{\kappa }(\mathcal{U}).$$

$${\mathcal{U}}_{ı,\hat{ı}}^{(\ell +1)}={\Pi }_{{\mathbb{U}}_{ı,\hat{ı}}}\left(\frac{1}{{\tau }_{\mathcal{U}}}{\Pi }_{{\mathbb{U}}_{ı,\hat{ı}}}\left({\tau }_{\mathcal{U}}{\mathcal{U}}_{ı,\hat{ı}}^{(\ell )}-\alpha {\nabla }_{{\mathcal{U}}_{ı,\hat{ı}}}\mathcal{L}({\mathcal{U}}^{(\ell )},{\mathbf{\lambda }}^{(\ell )})\right)\right),$$

$${{\lambda }_{\kappa }}^{(\ell +1)}={\Pi }_{\mathbb{D}}\left(\frac{1}{{\tau }_{{\lambda }_{\kappa }}}{\Pi }_{\mathbb{D}}\left({\tau }_{{\lambda }_{\kappa }}{{\lambda }_{\kappa }}^{(\ell )}+\beta {\nabla }_{{\lambda }_{\kappa }}\mathcal{L}({\mathcal{U}}^{(\ell )},{\mathbf{\lambda }}^{(\ell )})\right)\right),$$

$$\mathbb{E}\left[\mathcal{F}(\mathcal{U})\right]$$

$$\leftrightarrow {\mathbb{E}}^{\mathsf{T}}\left[P_b\right](\tilde{P}\mathcal{U})+\frac{1}{2}{\Vert \tilde{P}\mathcal{U}\Vert }_2^2+\frac{\rho }{2}{\Vert \mathcal{U}\Vert }_2^2,$$

$${\nabla }_{{\mathcal{U}}_{ı,\hat{ı}}}\mathbb{E}[\mathcal{F}({\mathcal{U}}^{(\ell )})]=$$

$$+\rho {\mathcal{U}}_{ı,\hat{ı}}^{(\ell )},$$

Full text

I Introduction

Electric vehicles (EVs) have shown their effectiveness in increasing energy conversion efficiency and reducing GHG emissions [1]. Besides, with proper control, the charging process of a large fleet of EVs can be leveraged for provisioning grid-level services [2]. Extensive studies have been conducted on the modeling and control problems of EVs for “valley-filling” [3], load balancing [4], and frequency regulation [5]. The potential of EV charging control for facilitating renewable energy integration and providing additional services is discussed in [6].

Apart from the benefits, charging a large fleet of EVs without proper control would have negative impacts on distribution networks. Lopes et al. [7] surveyed the challenges of integrating EVs into the mid- and low-voltage distribution networks, including voltage sag, transformer overloading, network congestion, and increased power loss. The negative impacts caused by non-network-aware controlled charging were verified in [8]. Hence, an impact-free EV integration necessitates designing EV charging control in the context of distribution network constraints.

Recently, an increasing amount of attention is being paid to controlling EV charging for distribution network impact alleviation and grid service provisioning. Considering network impacts, Richardson et al. [9] optimized EV charging profiles to minimize the total power consumption. In [10], a real-time control for power loss minimization was studied. Authors in [11] and [12] separately developed decentralized EV charging controllers to achieve valley-filling under voltage constraints. Studies in a large amount of literature have shown that one of the most effective services EVs could provide is to fill the load valley at night [13]. In our paper, we target at valley-filling meanwhile meeting both local charging needs and distribution network constraints on voltage.

Most existing EV charging control algorithms assume a powerful central controller that can handle large computational loads, [14, 9, 10, 15]. As the EV deployment in the distribution network keeps increasing, centralized approaches however do not scale well and become unrealistic. In this paper, we aim to design an optimization-based decentralized control scheme that does not require a communication network among EV chargers. Achieving this is notably challenging, because firstly the valley-filling objective constructs a coupled and non-separable function, and secondly the network voltage constraints strongly couple the individual charging powers in form of linear inequalities [12]. For this type of problem, Koshal et al. [16] developed a regularized primal-dual subgradient approach, however convergence errors exist due to the regularization of the Lagrangian. In [11], the authors proposed an Alternating Direction Method of Multipliers (ADMM)-based decentralized algorithm. Though proved to be practically applicable, the required two-layer communication network complicates the communication and poses computing burdens to all buses. Additionally, ADMM-based algorithms often face the problem of a large number of iterations. In our previous work [12], we proposed a shrunken-primal-dual subgradient (SPDS) approach that can eliminate convergence errors, reduce the number of iterations, and alleviate communication loads.

Unfortunately, the aforementioned algorithms are not capable of handling uncertainties in EV charging (e.g., random baseline load and random customer’s behaviors). Among the limited literature, Liu et al. [17] proposed an SPDS-based event-triggering decentralized control algorithm to handle EV drivers’ random arrivals and departures, where customers’ baseline loads are assumed perfectly predicted. Hassan et al. [18] developed a chance-constrained ADMM for decentralized distributed energy resources control considering random solar generation. However, decentralized implementation of chance constraints on voltage was treated deterministically, leaving that problem open. Handling chance constraints on voltage in a decentralized fashion is notably challenging because: (i) Voltage constraints are strongly coupled; (ii) Gradient of the chance constraint function is ill-conditioned, leading to zero updates in iteration-based algorithms. In this paper, we aim to leverage the outstanding convergency, optimality, and speed of SPDS, to design a chance-constrained SPDS (CC-SPDS) to handle decentralized EV charging control problems that have chance constraints on nodal voltages. To the author’s best knowledge, this has never been attempted before.

The main contribution of this paper is two-fold. First, this is the first paper studying decentralized EV charging control for provisioning valley-filling under chance constraints due to uncertain baseline loads. Second, a novel CC-SPDS algorithm is developed to solve a class of optimization problems consisting of non-separable objective functions and strongly coupled network chance constraints.

II Preliminaries & Problem formulation

II-A EV charging model

For a radial distribution network with $h$ nodes, let $n_{\imath}$ denote the number of EVs connected at the $\imath$th node and $n=\sum_{\imath=1}^{h}n_{\imath}$ denote the total number of EVs. The individual EV charging dynamics can be represented as [17]

[TABLE]

where the subscript ${}_{\imath,\hat{\imath}}$ denotes the $\hat{\imath}$th EV connected at the $\imath$th node, $T$ is the discrete-time index, $x_{\imath,\hat{\imath}}(T)$ denotes the energy remaining to be charged, $B_{\imath,\hat{\imath}}{=}-\eta_{\imath,\hat{\imath}}\Delta t\bar{P}_{\imath,\hat{\imath}}$, $\eta_{\imath,\hat{\imath}}$ is the charging efficiency, $\Delta t$ is the sampling interval, $\bar{P}_{\imath,\hat{\imath}}$ is the maximum charging power, and $0{\leq}u_{\imath,\hat{\imath}}(T){\leq}1$ is the control signal.

Let $k_{\imath,\hat{\imath}}$ and $k_{\imath,\hat{\imath}}{+}K_{\imath,\hat{\imath}}$ denote individual plug-in time and designated charging deadlines; let $[k,k{+}K]$ denote the valley-filling service period. In this paper, we assume that individual charging periods $[k_{\imath,\hat{\imath}},k_{\imath,\hat{\imath}}{+}K_{\imath,\hat{\imath}}]$ are contracted to fully cover the valley-filling period. This assumption can be easily relaxed. The charging efficiencies, initial state of charge, maximum charging powers, and battery capacities are heterogeneous. Augmenting all connected EVs, we have the system dynamics represented as

[TABLE]

where $x(T)=[x_{1,1}(T)\cdots x_{1,n_{1}}(T)\cdots x_{h,1}(T)\cdots x_{h,n_{h}}(T)]^{\mathsf{T}}$, $u(T)=[u_{1,1}(T)\cdots u_{1,n_{1}}(T)\cdots u_{h,1}(T)\cdots u_{h,n_{h}}(T)]^{\mathsf{T}}$, and $B_{\imath,\hat{\imath}}^{c}$ is the $(\imath,\hat{\imath})$th column of the diagonal matrix $\text{diag}\left\{B_{\imath,\hat{\imath}}\right\}\in\mathbb{R}^{n\times n}$, $\hat{\imath}=1,\ldots,n_{\imath},~{}\imath=1,\ldots,h$.

Augmenting the system state $x(T)$ and individual control signal $u_{\imath,\hat{\imath}}(T)$ in (1) along the valley-filling period $[k,k+K]$, we have

[TABLE]

Herein, $x(k+\kappa|k)$, $\kappa=1,\cdots,K$, is the system state at time $k+\kappa$ predicted at time $k$; $u_{\imath,\hat{\imath}}(k+\kappa-1|k)$, $\kappa=1,\cdots,K$, is the control signal at time $k+\kappa-1$ predicted at time $k$.

To meet all drivers’ charging requirements, the energy remaining to be charged at the end of valley-filling, i.e., the $K$th vector block in $\mathcal{X}(k)$, must satisfy

[TABLE]

where $\mathcal{B}_{{\imath,\hat{\imath}}}^{l}=[B_{{\imath,\hat{\imath}}}^{c}~{}B_{{\imath,\hat{\imath}}}^{c}~{}\cdots~{}B_{{\imath,\hat{\imath}}}^{c}]\in\mathbb{R}^{n\times K}$.

II-B Distribution network model

In this paper, we consider radial distribution networks. Let $\mathbb{H}=\{\imath|\imath=1,\ldots,h\}$ denote the set of downstream nodes and let $\mathbb{S}$ denote the set of all downstream line segments. The feeder head decouples interactions in the downstream distribution system from the rest of the grid and maintains its own voltage magnitude $\left|V_{0}\right|$. To better illustrate the EV charging impacts and how controlled EV charging can help alleviate the impacts, we do not consider any distributed energy resources (DERs), voltage regulations, or reactive power supplies.

At time $T$, let $\left|V_{\imath}(T)\right|$ denote the voltage magnitude at Node $\imath$; let $p_{\imath}(T)$ and $q_{\imath}(T)$ denote the real and reactive power consumption at Node $\imath$; and with a slight abuse of notations, let $r_{\imath\jmath}+ix_{\imath\jmath}$ denote the impedance of the line segment $(\imath,\jmath)$. By omitting line losses, the LinDistFlow model of this distribution network is derived as [19, 20, 12]

[TABLE]

where

[TABLE]

and

[TABLE]

where $\mathbb{S}_{\imath}$ ($\mathbb{S}_{\jmath}$) is the set containing downstream line segments connecting the feeder head and Node $\imath$ (Node $\jmath$) [20].

At each node and each house, real and reactive loads can be separated into baseline loads and EV charging loads, indicating

[TABLE]

where $p_{\imath,\hat{\imath}}^{b}(T)$ and $q_{\imath,\hat{\imath}}^{b}(T)$ denote the baseline real and reactive power consumption, respectively, and $p_{\imath,\hat{\imath}}^{EV}(T)$ and $q_{\imath,\hat{\imath}}^{EV}(T)$ denote the real and reactive EV charging power, respectively. As it is unrealistic to assume perfect forecast of the baseline load, a realistic problem formulation must consider $p_{\imath,\hat{\imath}}^{b}(T){\sim}\mathcal{N}(\mu^{p}_{\imath,\hat{\imath}}(T),\sigma^{p}_{\imath,\hat{\imath}}(T))$ and $q_{\imath,\hat{\imath}}^{b}(T){\sim}\mathcal{N}(\mu^{q}_{\imath,\hat{\imath}}(T),\sigma^{q}_{\imath,\hat{\imath}}(T))$, where $\mu^{p}_{\imath,\hat{\imath}}(T)$ and $\mu^{q}_{\imath,\hat{\imath}}(T)$ denote the means of the normal distributions, and $\sigma^{p}_{\imath,\hat{\imath}}(T)$ and $\sigma^{q}_{\imath,\hat{\imath}}(T)$ denote the standard deviations of the normal distributions. In this paper, we make the following three assumptions:

• A1

At each time $T$, power consumption of all houses follows independent normal distributions. The standard derivation $\sigma^{p}_{\imath,\hat{\imath}}(T)$ is assumed to be a constant $\sigma^{p}$ at all time.

• A2

EVs do not consume reactive power, i.e., $q_{\imath,\hat{\imath}}^{EV}(T)=0$.

• A3

Power factor at each house is a constant, which implies $q_{\imath,\hat{\imath}}^{b}(T)=\gamma p_{\imath,\hat{\imath}}^{b}(T)$.

Remark 1: In A1, other distributions can also be used, but normal distribution can better reflect customers’ random behaviors. Constant standard deviation is assumed for the simplicity of presentation. A2 follows as no fourth-quadrant charger is considered. A3 follows by assuming power factor correction equipment is available at each node. $\blacksquare$

Let $\hat{p}^{b}(T)=[p_{1,1}^{b}(T)\cdots p_{h,n_{h}}^{b}(T)]^{\mathsf{T}}\in\mathbb{R}^{n}$ contain baseline power at all $n$ houses, then $\hat{p}^{b}(T){\sim}\mathcal{N}\left(\boldsymbol{\mu}^{p}(T),\boldsymbol{\Sigma}^{p}\right)$, where $\boldsymbol{\mu}^{p}(T)=[\mu^{p}_{1,1}(T)\cdots\mu^{p}_{h,n_{h}}(T)]^{\mathsf{T}}$ and $\boldsymbol{\Sigma}^{p}=\text{diag}\left\{(\sigma^{p})^{2}\right\}\in\mathbb{R}^{n\times n}$. Further let $p^{b}(T)=G\hat{p}^{b}(T)=[p_{1}^{b}(T)~{}\cdots~{}p_{h}^{b}(T)]^{\mathsf{T}}$ denote the nodal aggregated baseline power, where $G=\text{diag}\left\{G_{\imath}\right\}\in\mathbb{R}^{h\times n}$, $\imath=1,\ldots,h$, $G_{\imath}={\mathbb{1}}_{n_{\imath}}^{\mathsf{T}}$ is the nodal aggregation vector. Since (3) is linear, the voltage drop caused by the baseline load can be represented as

[TABLE]

where $\boldsymbol{H}=2(\boldsymbol{R}+\gamma\boldsymbol{X})G$, yielding

[TABLE]

where $\bar{P}=\text{diag}\left\{\bar{P}_{\imath,\hat{\imath}}\right\}$, $\hat{\imath}=1,\ldots,n_{\imath}$, $\imath=1,\ldots,h$ and $\boldsymbol{V}_{b}(T)\sim\mathcal{N}(\tilde{\boldsymbol{\mu}}^{p}(T),\tilde{\boldsymbol{\Sigma}}^{p})$ with $\tilde{\boldsymbol{\mu}}^{p}(T)=\boldsymbol{H}\boldsymbol{\mu}^{p}(T)$ and $\tilde{\boldsymbol{\Sigma}}^{p}=\boldsymbol{H}\boldsymbol{\Sigma}^{p}\boldsymbol{H}^{\mathsf{T}}$.

Let $D\in\mathbb{R}^{h\times n}$ denote $-2\boldsymbol{R}G\bar{P}$, $y_{d}(T)$ denote $\boldsymbol{V}_{0}-\boldsymbol{V}_{b}(T)$, and $y(T)$ denote $\boldsymbol{V}(T)$, we have

[TABLE]

where $y_{d}(T){\sim}\mathcal{N}(\boldsymbol{V}_{0}-\tilde{\boldsymbol{\mu}}^{p}(T),\tilde{\boldsymbol{\Sigma}}^{p})$. The system output augmented along the valley-filling period is written as

[TABLE]

where

[TABLE]

and

[TABLE]

Herein, $\oplus$ denotes the direct matrix sum.

II-C Deterministic valley-filling

Suppose $p^{b}_{\imath,\hat{\imath}}(T)$ is deterministic or can be accurately forecasted. According to [12, 17], optimal EV charging sequences for valley-filling under nodal voltage constraints are obtained by solving

[TABLE]

where the time index $k$ is dropped for simplicity, $\hat{\boldsymbol{V}}_{0}=\left[{\boldsymbol{V}}_{0}^{\mathsf{T}}\cdots{\boldsymbol{V}}_{0}^{\mathsf{T}}\right]^{\mathsf{T}}\in\mathbb{R}^{nK}$ and

[TABLE]

Herein, $P_{b}=\left[\left\|p^{b}(k)\right\|_{1}\cdots\left\|p^{b}(k+K-1)\right\|_{1}\right]^{\mathsf{T}}$ denotes the aggregated baseline load profile at the feeder head, $\tilde{P}$ is the EV charging power aggregation matrix, $\rho/2\left\|\mathcal{U}\right\|_{2}^{2}$ is for battery state of health protection, and $\underline{\nu}$ is the lower bound of the transformer service range.

II-D Chance-constrained valley-filling

In more realistic situations, where baseline load cannot be accurately forecasted, we have normally distributed $p_{\imath,\hat{\imath}}^{b}(T){\sim}\mathcal{N}(\mu^{p}_{\imath,\hat{\imath}}(T),\sigma^{p})$. Under this circumstance, the deterministic problem (4) becomes

[TABLE]

where $[\hat{\mathcal{Y}}]_{\kappa}{\sim}\mathcal{N}((\underline{\nu}^{2}-1)\boldsymbol{V}_{0}+\tilde{\boldsymbol{\mu}}^{p}(k+\kappa-1),\tilde{\boldsymbol{\Sigma}}^{p})\in\mathbb{R}^{h}$ is the $\kappa$th block of $\hat{\mathcal{Y}}=\underline{\nu}^{2}\hat{\boldsymbol{V}}_{0}-\mathcal{Y}_{dk}$, $[z]_{\kappa}=[\sum_{\imath=1}^{h}\sum_{\hat{\imath}=1}^{n_{\imath}}\mathcal{D}_{\imath,\hat{\imath}}\mathcal{U}_{\imath,\hat{\imath}}]_{\kappa}\in\mathbb{R}^{h}$ is the $\kappa$th block of $\sum_{\imath=1}^{h}\sum_{\hat{\imath}=1}^{n_{\imath}}\mathcal{D}_{\imath,\hat{\imath}}\mathcal{U}_{\imath,\hat{\imath}}$, and $1-\delta$ is the allowed chance of violating voltage limits. For simplicity of presentation, we use $\mathbb{E}$ to denote the expectation operator $\mathbb{E}_{p^{b}_{\imath,\hat{\imath}}{\sim}\mathcal{N}(\cdot,\cdot)}$ hereinafter. Designing decentralized algorithms for stochastic problem (5) is more complicated than for problem (4) as the coupled objective function becomes an expectation and the coupled inequality constraints become coupled chance constraints.

III Chance-constrained SPDS design

Let the multivariate normal cumulative distribution function (MVNCDF) of $[\hat{\mathcal{Y}}]_{\kappa}$ be denoted by $F_{[\hat{\mathcal{Y}}]_{\kappa}}(\cdot)$, then the chance constraints in problem (5) become

[TABLE]

yielding the unconstrained Lagrangian of problem (5) as

[TABLE]

It is proved in [12] that a convex, strongly coupled optimization problem can be efficiently solved by SPDS by following

[TABLE]

where $\Pi_{\cdot}(\cdot)$ is the Euclidean projection, $\ell$ denotes the iteration number, $\alpha$ and $\beta$ are the update step sizes for primal and dual variables, respectively, $0<\tau_{\mathcal{U}},\tau_{\lambda_{\kappa}}<1$ are the shrinking parameters, and $\mathbb{D}$ is the constraint set for dual variables. The definition of $\mathbb{D}$ is referred to [12, 17]. It has been shown in [12] that the convergence and optimality are guaranteed if $\alpha$ and $\beta$ are appropriately chosen for a designated tuple $(\tau_{\mathcal{U}},\tau_{\lambda_{\kappa}})$.

Nevertheless, SPDS in (7) cannot be directly applied to the chance-constrained problem (5) due to two challenges: (i) The feasible region is non-convex due to the non-convexity of $F_{[\hat{\mathcal{Y}}]_{\kappa}}(\cdot)$; (ii) The gradient $\nabla_{\mathcal{U}_{\imath,\hat{\imath}}}d_{\kappa}(\mathcal{U})$ is bounded and in most cases it realizes the value of [math], implying that updates of the primal variables can barely capture constraints. (This will be illustrated later). In the following, we will present the details of CC-SPDS and show via simulations that it can handle the non-convex chance-constrained problem (5). Rigorous proof of convergency and optimality gap will be provided in a journal extension of this paper.

The objective function of problem (5) can be rewritten as

[TABLE]

and its gradient at $\mathcal{U}^{(\ell)}$ can be calculated as

[TABLE]

where inside the bracket of the first term on the right hand side is the mean of the total load profile at the $\ell$th iteration. This is a universal information broadcasted to all EVs in each iteration. To obtain the complete form of $\nabla_{\mathcal{U}_{\imath,\hat{\imath}}}\mathcal{L}(\mathcal{\mathcal{U},\boldsymbol{\lambda}})$ we first introduce the following theorem.

Theorem 1[21] Let $\xi{\sim}\mathcal{N}(\boldsymbol{\mu},\boldsymbol{\Sigma})$ with some positive definite covariance matrix $\boldsymbol{\Sigma}=(\sigma_{ij})\in\mathbb{R}^{s\times s}$. Then, the distribution function $F_{\xi}$ is continuously differentiable at any $z\in\mathbb{R}^{s}$ and

[TABLE]

Here, $f_{\xi_{j}}$ denotes the 1-D normal density of the component $\xi_{j}$, $\xi_{\bar{j}}$ is an $(s-1)$-D normal random vector distributed according to $\xi_{\bar{j}}{\sim}\mathcal{N}(\hat{\boldsymbol{\mu}},\hat{\boldsymbol{\Sigma}})$, $\hat{\boldsymbol{\mu}}$ results from the vector $\boldsymbol{\mu}+\sigma_{jj}^{-1}(z_{j}-\mu_{j})\boldsymbol{\sigma}_{j}$ by deleting component $j$, and $\hat{\boldsymbol{\Sigma}}$ results from the matrix $\boldsymbol{\Sigma}-\sigma_{jj}^{-1}\boldsymbol{\sigma}_{j}\boldsymbol{\sigma}_{j}^{\mathsf{T}}$ by deleting row $j$ and column $j$, where $\boldsymbol{\sigma}_{j}$ refers to column $j$ of $\boldsymbol{\Sigma}$. Moreover, $\hat{\boldsymbol{\Sigma}}$ is positive definite. $\blacksquare$

By applying Theorem 1, we can write

[TABLE]

where $G_{\kappa}=[\boldsymbol{0}\cdots\boldsymbol{0}~{}\underbrace{\boldsymbol{I}_{h}}_{\kappa\text{th}}\ ~{}\boldsymbol{0}\cdots\boldsymbol{0}]$, $\boldsymbol{J}(\cdot)$ denotes the Jacobian matrix, $[z]_{\kappa,j}$ denotes the $j$th component of $[z]_{\kappa}$, and $[z]_{\kappa,\bar{j}}$ is obtained by removing the $j$th component of $[z]_{\kappa}$.

Though theoretically sound, directly applying $\nabla_{\mathcal{U}_{\imath,\hat{\imath}}}d_{\kappa}(\mathcal{U})$ in (8) to the primal update equation (7a) is ill-conditioned. For example in the independent bivariate normal distribution case, the gradient of the cumulative function is $\nabla F_{x,y}(X,Y)=[f_{x}(X)F_{y}(Y)~{}f_{y}(Y)F_{x}(X)]^{\mathsf{T}}$.

As illustrated in Figure 1, there are four cases where $f_{x}(X)F_{y}(Y)$ is extremely small or zero: (i & ii) $X$ is far away from $\mu_{x}$ to the right and $Y$ is far away from $\mu_{y}$ to the left (and right); (iii & iv) $X$ is far away from $\mu_{x}$ to the left and $Y$ is far away from $\mu_{y}$ to the left (and right); Translated into (6), these four cases mean: (i & ii) Constraint is satisfied with a high probability; No (or slight) updates w.r.t. the constraint are needed in the iteration; (iii & iv) Constraint is not satisfied with a high probability; Updates w.r.t. the constraint are needed in the iteration. However, in the last two cases, updates cannot be realized as the updating direction (gradient) is [math]. To overcome this ill-conditioned gradient, we define an adjusted gradient as

[TABLE]

where

[TABLE]

and $\hat{f}_{y}(Y)$, $\hat{F}_{x}(X)$ are defined in the same way. The adjusted density function and cumulative function are illustrated in Figure 1 as the dashed curves which only adjust the left half planes. These adjustments guarantee effective updates once chance constraints are violated. Following this idea, we have ${\nabla}_{\mathcal{U}_{\imath,\hat{\imath}}}d_{\kappa}(\mathcal{U})$ in (8) adjusted to

[TABLE]

and consequently the primal update equation (7a) becomes

[TABLE]

where

[TABLE]

The dual update equation (7b) can be readily calculated as

[TABLE]

Together (10) and (7b) construct the CC-SPDS.

Remark 2: In EV charging control, the adjusted cumulative gradient vector on the right hand side of (9) is broadcasted to all EVs as a public key. Each EV uses its private key $\mathcal{D}_{\imath,\hat{\imath}}$ to interpret the public key, so that the control sequence update can fully acknowledge the global chance constraints. $\blacksquare$

Remark 3: CC-SDPS, especially the approach of adjusting the gradient of MVNCDFs, can be applied to many other popular distributions, e.g., Poisson distribution and uniform distribution, with $\hat{f}$ and $\hat{F}$ properly defined. $\blacksquare$

IV Simulations

The simplified IEEE 13 Node Test Feeder used in [12] is adopted in this paper. In specific, Node 632 and Node 671 are not connected to any load, and each of the other nodes is connected with $70$ houses equipped with level-2 chargers. Battery capacities are uniformly distributed in $[18,20]$ kWh. Initial and designated SOCs are uniformly distributed in $[0.3,0.5]$ and $[0.7,0.9]$, respectively. Primal step size $\alpha=1\times 10^{-11}$; dual step size $\beta=2$; $d_{\lambda_{\kappa}}=5\times 10^{5}$ in $\mathbb{D}$; shrinking parameters are $\tau_{\mathcal{U}}=\tau_{\lambda_{\kappa}}=0.974$. The above parameters were empirically chosen to accelerate the convergence speed. Initial values of $\mathcal{U}_{i}^{(0)}$ and $\lambda_{\kappa}^{(0)}$ are set to zero vectors. The voltage lower bound is chosen to $\underline{\nu}=0.954$, which is slightly higher than the $0.95$ p.u. in ANSI C84.1 service standard, to compensate for the discarded line losses in LinDistFlow model and possible voltage violations. For the normal distribution of baseline load at each house, the mean is obtained by scaling the data from Southern California Edison [22]; the constant standard deviation $\sigma^{p}=400$W. The contracted valley-filling service period is from 19:00 to 8:00 next day. The sampling time is 15 minutes. The allowed voltage constraint violation chance $1{-}\delta$ is set to $10\%$.

Figure 2 shows the valley-filling performance under the proposed CC-SPDS decentralized framework, where convergence of the algorithm can be readily revealed.

In total 50 iterations were executed and each iteration took about 0.16 second. Comparing with the purely deterministic case discussed in [12, 23] where the converged total load profile is flat, the converged total load profile in the stochastic case presents slight variations. This is because the objective function only captures the mean of the baseline load profile, i.e., CC-SPDS fills the mean baseline load valley, and once random behaviors realized the valley-filling performance is slightly compromised. Figure 2 also shows the total load profile resulting from directly applying deterministic SPDS to the stochastic case, i.e., solving (4) by using the mean baseline loads in both the objective function and voltage constraints. The difference between CC-SPDS and SPDS w.r.t. valley-filling performance is neglectable, as both of them are using mean baseline load in the objective function.

Figure 3 shows the nodal voltage magnitudes of all nodes during the valley-filling period.

Under baseline loads, all nodal voltage magnitudes are above $0.954$ p.u. lower bound. Under the proposed CC-SPDS-based decentralized EV charging control framework, only few nodes occasionally exceed the bound during 1:00 to 5:00. Under the SPDS-based control, voltage magnitudes are continuously violating the limit during 23:00 to 5:30 and are much lower than those of the CC-SPDS control. This observation reveals the effectiveness of CC-SPDS.

Figure 4 presents the statistics of voltage violations.

The upper plot reveals that in the total 52 time slots during the valley-filling period, CC-SPDS can suppress 50 (96%) of them to have less than 10% (only one) of all nodes violating 0.954 p.u., and 2 (4%) of them to have 25% (three) of all nodes violating 0.954 p.u. At each time slot, percentage of the total nodes violating the limit can be treated as a measure of the chance constraint. The 25% violation is as expected since the chance constraints only reflect a probability which cannot be fully reflected by a single sample. Comparing with CC-SPDS, SPDS has the histogram more distributed to the right, leading to more violation events and greater violation magnitude. In the lower plot, boxes represent the 25th and 75th percentiles and legs with bars represent the extreme values. It can be concluded that under the control of CC-SPDS: Nodal voltage magnitudes have less variations; The 25th percentiles are guaranteed to stay above 0.954 p.u.; Number and magnitudes of violations are less than those of SPDS. These observations further reveal the effectiveness of CC-SPDS.

V Conclusion

This paper developed a novel decentralized EV charging control framework where distribution network constraints were formulated as chance constraints due to customers’ random behaviors. A new decentralized optimization algorithm – CC-SPDS was developed to solve the chance-constrained valley-filling problem as well as generic optimization problems that have coupled objective functions and coupled global chance constraints. Convergence and effectiveness of CC-SPDS were verified via simulations. CC-SPDS can be generalized to other distributions and extended to consider local chance constraints.