Distributed State Estimation for AC Power Systems using Gauss-Newton ALADIN

TL;DR

This paper introduces a distributed state estimation algorithm for AC power systems using a tailored ALADIN method with Gauss-Newton approximations, improving efficiency in handling non-convex power flow equations.

Contribution

It develops a novel distributed state estimator based on ALADIN with Gauss-Newton Hessian approximations for efficient nonlinear power system estimation.

Findings

Successfully applied to IEEE 30-Bus system

Demonstrates improved computational efficiency

Outperforms existing distributed estimation methods

Abstract

This paper proposes a structure exploiting algorithm for solving non-convex power system state estimation problems in distributed fashion. Because the power flow equations in large electrical grid networks are non-convex equality constraints, we develop a tailored state estimator based on Augmented Lagrangian Alternating Direction Inexact Newton (ALADIN) method, which can handle the nonlinearities efficiently. Here, our focus is on using Gauss-Newton Hessian approximations within ALADIN in order to arrive at at an efficient (computationally and communicationally) variant of ALADIN for network maximum likelihood estimation problems. Analyzing the IEEE 30-Bus system we illustrate how the proposed algorithm can be used to solve highly non-trivial network state estimation problems. We also compare the method with existing distributed parameter estimation codes in order to illustrate its performance.

Figures

Tables

Table 1. TABLE I: Relative a-posteriori standard deviations at selected nodes.

Relative Standard Deviation of the A-Posteriori Distribution
Bus$\#$	$\theta$	$v$	$p$	$q$
1	$\ast$	$0.12\%$	$0.16\%$	$2.65\%$
8	$0.32\%$	$0.12\%$	$0.16\%$	$0.23\%$
13	$0.65\%$	$0.17\%$	$0.05\%$	$0.79\%$
20	$0.25\%$	$0.18\%$	$0.46\%$	$14.16\%$
30	$0.33\%$	$0.23\%$	$0.28\%$	$2.95\%$
AVG	$0.49\%$	$0.17\%$	$0.98\%$	$2.75\%$

Equations

$$Y_{k,l}=\left\{\begin{array}[]{ll}\sum\limits_{l\neq k}\left(g_{k,l}+jb_{k,l}\right)&\text{if}\;k=l,\\[7.11317pt] -\left(g_{k,l}+jb_{k,l}\right)&\text{if}\;k\neq l,\end{array}\right.$$

$$x_k^{\top }=({\theta }_k{v}_k{p}_k{q}_k{)}^{\top }\in {\mathbb{R}}^4,$$

$$0=p_k-v_k\sum _{l\in \mathcal{N}}{v}_l(G_{k,l}\cos ({\theta }_{k,l})+B_{k,l}\sin ({\theta }_{k,l})),$$

$$0=q_k-v_k\sum _{l\in \mathcal{N}}{v}_l(G_{k,l}\sin ({\theta }_{k,l})-B_{k,l}\cos ({\theta }_{k,l})),$$

$$\begin{array}{rl}{f}_p(x_k,x_l)=& v_k[v_k{g}_{k,l}-v_l{g}_{k,l}\cos ({\theta }_{k,l})]\\ & -v_k[v_l{b}_{k,l}\sin ({\theta }_{k,l})],\\ f_q(x_k,x_l)=& -v_k[v_k{b}_{k,l}-v_l{b}_{k,l}\cos ({\theta }_{k,l})]\\ & +v_k[v_l{g}_{k,l}\sin ({\theta }_{k,l})],\\ f_i(x_k,x_l)=& \frac{f_{pt}(x_k,x_l{)}^2+f_{qt}(x_k,x_l{)}^2}{v_k^2},\end{array}$$

$${\mathcal{F}}_{\mathcal{N},\mathcal{L}}(x):=\left(\begin{array}{c}{\Sigma }_k^{\frac{1}{2}}(x_k-{\hat{x}}_k{)}_{k\in \mathcal{N}}\\ W_{k,l}^{\frac{1}{2}}{\left(f(x_k,x_l)-{\hat{w}}_{k,l}\right)}_{(k,l)\in \mathcal{L}}\end{array}\right)$$

$$\begin{array}{rl}\underset{x}{\min }& \Vert {\mathcal{F}}_{\mathcal{N},\mathcal{L}}(x){\Vert }_2^2\\ & \text{s.t. }\text{eqref}eq:PFeq\text{for all }k\in \mathcal{N}.\end{array}$$

$${\theta }_k={\theta }_l,v_k=v_l,\text{ for all }(k,l)\in \mathcal{A},$$

$$\underset{z}{\min }\sum _{i\in \mathcal{R}}\Vert \mathcal{F}$$

$$\text{subject to}\sum _{i\in \mathcal{R}}{A}_i{z}_i$$

$${\mathcal{H}}_i(z_i)$$

$$\begin{array}{rl}\underset{y_i}{\min }& {\Vert {\mathcal{F}}_i(y_i)\Vert }_2^2+{\lambda }^{\top }{A}_i{y}_i+\frac{\rho }{2}{\Vert y_i-z_i\Vert }_2^2\\ \text{s.t.}& {\mathcal{H}}_i(y_i)=0\mid {\kappa }_i^d\end{array}$$

$$b_i={\mathcal{F}}_i(y_i),B_i=\nabla {\mathcal{F}}_i(y_i{)}^{\top },C_i=\nabla {\mathcal{H}}_i(y_i{)}^{\top }$$

$$\Vert \sum _{i\in \mathcal{R}}{A}_i{y}_i\Vert \le ϵ\text{ and }{\Vert y_i-z_i\Vert }_{\infty }\le ϵ.$$

$$\begin{array}{rl}\underset{\Delta y}{\min }& \sum _{i\in \mathcal{R}}{\Vert B_i\Delta y_i\Vert }_2^2+2\Delta y_i^{\top }{B}_i{b}_i\\ \text{s.t.}& \sum _{i\in \mathcal{R}}{A}_i(y_i+\Delta y_i)=0\mid {\lambda }^{\mathrm{QP}},\\ & C_i\Delta y_i=0\mid {\kappa }_i^{\mathrm{QP}}i\in \mathcal{R}.\end{array}$$

$$\left\|\left[\begin{array}[]{c}y-z^{*}\\ \kappa^{d}-\kappa^{*}\end{array}\right]\right\|_{2}\leq\alpha\left\|\left[\begin{array}[]{c}z-z^{*}\\ \lambda-\lambda^{*}\end{array}\right]\right\|_{2}$$

$$\left\|\left[\begin{array}[]{c}z^{+}-z^{*}\\ \lambda^{+}-\lambda^{*}\end{array}\right]\right\|_{2}\leq\gamma\left\|\left[\begin{array}[]{c}y-z^{*}\\ \lambda-\lambda^{*}\\ \kappa^{d}-\kappa^{*}\end{array}\right]\right\|_{2}\;.$$

$$\left\|\left[\begin{array}[]{c}z^{+}-z^{*}\\ \lambda^{+}-\lambda^{*}\end{array}\right]\right\|_{2}\leq\gamma(\alpha+1)\left\|\left[\begin{array}[]{c}z-z^{*}\\ \lambda-\lambda^{*}\end{array}\right]\right\|_{2}$$

$$\sum _{i\in \mathcal{R}}6|{\mathcal{N}}_i|+16|{\mathcal{N}}_i{|}^2+2|\mathcal{A}|$$

$${\Sigma }_k=\mathrm{cov}({\hat{x}}_k{)}^{-1}=\left[\begin{array}{cccc}10^4& 0& 0& 0\\ 0& 10^5& 0& 0\\ 0& 0& 10^4& 0\\ 0& 0& 0& 10^4\end{array}\right],$$

$$W_{k,l}=\mathrm{cov}({\hat{w}}_k{)}^{-1}=\left[\begin{array}{ccc}10^4& 0& 0\\ 0& 10^4& 0\\ 0& 0& 10^4\end{array}\right]$$

$${\left[\begin{array}{c}I\\ 0\\ 0\end{array}\right]}^{\top }{\left[\begin{array}{ccc}{B}_{\mathcal{R}}^{\top }{B}_{\mathcal{R}}& C_{\mathcal{R}}^{\top }& A_{\mathcal{R}}^{\top }\\ C_{\mathcal{R}}^{\top }& 0& 0\\ A_{\mathcal{R}}^{\top }& 0& 0\end{array}\right]}^{-1}\left[\begin{array}{c}I\\ 0\\ 0\end{array}\right],$$

Full text

I Introduction

State estimation is of increasing importance in modern electricity transmission and distribution systems. Due to the integration of renewable energy systems, effective grid operation often requires detailed knowledge of the system state. High-accuracy measurement devices are usually costly. Hence it is relevant to consider all available information and also cost-effective (including possibly inaccurate) measurement devices for determining the power system’s state.111We remark that in power systems the notion of state variables differs slightly from control. Hence here we refer to a set of variables defined as the solution of a stationary nonlinear system of equations. A standard method to solve the arising Power System State Estimation (PSSE) problem is via weighted nonlinear least squares [27, 1, 22].

Centralized formulations of AC PSSE—i.e. considering the full AC power flow equations—have a long history and can be traced back to [25]. AC PSSE is in general hard to solve as it is usually formulated as nonlinear least squares problem yielding a large-scale non-convex optimization problem. Different formulations including polar vs. rectangular coordinates and algorithms with different Jacobian approximations, including exact Jacobian [27], $p$-$q$ decoupled Jacobian [13], and Gauss-Newton approximation [22], have been considered for the AC PSSE problem.

The non-convexity of the power flow equations makes large-scale PSSE problems difficult to solve. Hence several authors have considered convex formulations and relaxations of this problem. This includes DC approximations [26, 21] and more recently SDP relaxations [32, 19]. However, as voltage and reactive power are often needed (especially in context of renewables) the practical usefulness of DC approximations is in general limited.

As power systems are large scale and as limited information exchange is desirable, distributed approaches have also been considered for AC and DC state estimation [25, 2, 16].222Sometimes distributed approaches are also called hierarchical depending on the amount of central coordination [11]. The DC case is considered in [26, 31, 28]. In a distributed setting, AC PSSE is even harder to solve than in the centralized setting. The reason being that so far there are only a few algorithms for general distributed non-convex optimization [5, 14]. Classical methods for distributed AC state estimation often exploit the sparsity pattern of the Jacobian of the measurement equations. These works started already in the 1970s with the seminal paper [25], continued with [2, 15]; and can also be found today [17]. Recently [20] proposed a distributed Gauss-Newton approach using matrix splitting techniques with promising results. Exploiting sparsity one obtains distributed methods where a coordinator typically solves preferably simple coordination problems.

A second and more recent line of research applies distributed optimization techniques coming from convex optimization to AC PSSE. These approaches include the auxiliary problem principle [7] and the popular Alternating Direction of Multipliers Method (ADMM) [16]. An algorithm based on gossiping techniques can be found in [31]. These methods usually have an advantage over splitting techniques—they are typically decentralized, i.e. they avoid central coordination and communicate based on neighborhood information only. However, despite working well for many cases, these methods usually have limited convergence guarantees for AC PSSE. A recent overview on distributed AC state estimation can be found in [11].

In this paper, we follow a different route tailoring the Augmented Lagrangian Alternating Direction Inexact Newton (ALADIN) method [14] to AC PSSE problems. After introducing the problem at hand in Section II, we explain how to exploit its distributed structure in Section III. Section IV introduces the main algorithmic contribution of this paper; i.e. we construct a variant of ALADIN based on generalized Gauss-Newton Hessian approximation. Such Gauss-Newton approximation based approaches have been analyzed exhaustively for unconstrained nonlinear least-squares problems [23, Chapter 10.3] as well as in the context of centralized parameter estimation of constrained problem as analyzed in [3, 4] and [24]. However, one interesting contribution of Section IV is analyzing such Hessian approximations in the context of distributed least squares estimation with the ALADIN framework proving a local convergence proof. The main contribution of this paper is presented in Section V, where we not only illustrate the performance of the proposed algorithm on the IEEE 30-bus system in comparison to the widely used ADMM algorithm [5], but also elaborate on the communication effectiveness of the proposed scheme.

II Power System State Estimation

This section briefly reviews the main physical relations in electrical grids, including power flow equations, and introduces the AC PSSE problem.

II-A Grid Model

This paper considers a power system $(\mathcal{N},\mathcal{L},Y)$ consisting of a set of buses $\mathcal{N}=\{1,\dots,N\}$, a set of transmission lines $\mathcal{L}\subseteq\mathcal{N}\times\mathcal{N}$, and a sparse, complex-valued admittance matrix $Y=G+jB\in\mathbb{C}^{N\times N}$ with $j=\sqrt{-1}$ . The admittance matrix is defined by

[TABLE]

where $g_{k,l}$ denotes the line conductance and $b_{k,l}$ denotes the line susceptance for all transmission lines $(k,l)\in\mathcal{L}$. For $(k,l)\notin\mathcal{L}$, we have $g_{k,l}=b_{k,l}=0$.

To each node in the grid, we assign a state as

[TABLE]

where $\theta_{i}$ is the voltage angle, $v_{i}$ is the voltage magnitude and $p_{i},q_{i}$ are the net active and reactive power at node $i\in\mathcal{N}$. The state of the grid is then defined as $x^{\top}=(x_{1},\dots,x_{N})^{\top}\in\mathbb{R}^{4N}$. The grid physics are described by the the power flow equations in polar form as

[TABLE]

for all nodes $k\in\mathcal{N}$ with $\theta_{k,l}=\theta_{k}-\theta_{l}$, cf. [1]. Note that $G_{k,l}$ and $B_{k,l}$ refer to the real and imaginary pats of the entries of the admittance matrix $Y$.

II-B Measurement Functions

PSSE aims at determining the steady state, $x$, of the grid using a given set of measurements. In general one considers two types of measurements: firstly one can directly measure the system states $x_{k}$ (or parts thereof) at the nodes. And secondly, one can attempt to measure the power flows at the transmission lines, which depend on the state of the grid at neighboring nodes. In order to arrive at a model that allows us to take the second type of measurements into account, we introduce so-called measurement functions, which relate the nodal states to the power flows over the transmission lines. These functions are given by

[TABLE]

where $f_{p},f_{q},f_{i}:\mathbb{R}^{4}\times\mathbb{R}^{4}\rightarrow\mathbb{R}$ denote the active power, respectively, reactive power and the current in the transmission line $(k,l)\in\mathcal{L}$. The complete vector-valued measurement function $\mathcal{F}_{\mathcal{N},\mathcal{L}}:\mathbb{R}^{4|\mathcal{N}|}\rightarrow\mathbb{R}^{4|\mathcal{N}|+3|\mathcal{L}|}$ is then given by

[TABLE]

where $f=(f_{p},\;f_{q},\;f_{i})^{\top}$ is used. Moreover, we use the shorthand $\hat{w}_{k,l}=(\hat{p}_{k,l},\;\hat{q}_{k,l},\;\hat{i}^{2}_{k,l})^{\top}$ to collect the actual measurements of the active and reactive power, $\hat{p}_{k,l}$ and $\hat{q}_{k,l}$, as well as the measurement $\hat{i}_{k,l}$ of the current at the transmission line $(k,l)\in\mathcal{L}$. In this context, it is assumed that approximations of the inverse of variance matrices of the associated measurement errors, $\Sigma_{k}\in\mathbb{S}^{4}$ and $W_{k,l}\in\mathbb{S}^{3}$, are given by positive semi-definite matrices [3].

II-C * Maximum Likelihood State Estimation*

The above model is used to formulate the AC PSSE problem of interest as the following nonlinear least-squares optimization problem

[TABLE]

Here, the underlying assumption is that the measurement errors have Gaussian probability distributions. This way (4) can be interpreted as a maximum likelihood parameter estimation problem [4, 29], recalling that the inverse variance matrices of the measurement errors, $\Sigma_{k}\in\mathbb{S}^{4}$ and $W_{k,l}\in\mathbb{S}^{3}$, are assumed to be given.

Remark 1

Although the theoretical properties of nonlinear least-squares optimization problems are rather well-understood [18], Problem (4) is in general a large-scale non-convex optimization problem with non-convex objective and non-convex constraint set over the complete electrical grid. As it turns out, AC PSSE problems are rather callenging to solve in practice. In particular, there might be multiple local minima and numerical algorithms might converge to one or the other minimum depending on the initialization [6].

III Distributed Least Squares Estimation

This section outlines the main strategy for breaking large-scale electrical grid networks into smaller sub-regions, thereby revealing the distributed structure of AC PSSE problems.

III-A Problem Decomposition

In order to solve (4) in distributed fashion, we reformulate (4) in affinely-coupled separable from [14]. To this end, we recall the partitioning method from [9], for alternative partitioning schemes see [10, 11, 16]. Figure 1 depicts the whole IEEE 30-Bus network as well as the partitioning strategy used throughout this paper.

We first divide the bus set $\mathcal{N}$ into several node sets $\mathcal{N}_{i}^{0}$, one for each region $\mathcal{R}=\{1,\dots R\}$ such that $\underset{i\in\mathcal{R}}{\cup}\mathcal{N}^{0}_{i}=\mathcal{N}$ and $\mathcal{N}^{0}_{i}\cap\mathcal{N}^{0}_{j}=\emptyset$ for all $i,j\in\mathcal{R}$ with $i\neq j$. At each transmission line connecting two adjacent regions, i.e. all $(m,n)\in\mathcal{L}$ with $m\in\mathcal{N}_{i}$ and $n\in\mathcal{N}_{j}$, $i\neq j$, we introduce an auxiliary bus pair $(k,l)$ and we collect all auxiliary bus pairs in set $\mathcal{A}$. The set of auxiliary buses of region $i\in\mathcal{R}$ are denoted as $\mathcal{A}_{i}$. Finally, we combine all auxiliary nodes and original nodes belonging to one region in combined node sets $\mathcal{N}_{i}=\mathcal{N}_{i}^{0}\cup\mathcal{A}_{i}$. The line set connecting original nodes with each other and all auxiliary nodes for region $i\in\mathcal{R}$ is denoted by $\mathcal{L}_{i}$. We assume a decomposition in the middle of each transmission line connecting two regions. This leads to new line admittances $y_{m,k}=2y_{m,n}$ and $y_{n,l}=2y_{m,n}$, respectively. The partitioning strategy is graphically illustrated in Figure 2 and Figure 1.

In order to resemble the original physical properties of the grid model, we introduce the consensus constraints333Note that there exist different variants of coupling. Here we only couple voltages as we observed best performance in this case. However additional active/reactive power coupling is also possible [8], cf. [11] for an overview of different coupling methods.

[TABLE]

a particular class of affine equality constraints.

III-B Distributed Formulation of the AC PSSE Problem

This section reformulates (4) in affinely coupled separable form, as required in the context of distributed optimization [14]. To this end, we introduce state vectors $z_{i}$ for all regions $i\in\mathcal{R}$ such that $z_{i}^{\top}=(\theta_{i}\quad v_{i}\quad p_{i}\quad q_{i})_{i\in\mathcal{N}_{i}}\in\mathbb{R}^{4|\mathcal{N}_{i}|}$. Summarizing the measurement equations for all nodes $\mathcal{N}_{i}$ and transmission lines $(k,l)\in\mathcal{L}_{i}$ in each region $i\in\mathcal{R}$, i.e. $\mathcal{F}_{\mathcal{N}_{i},\mathcal{L}_{i}}$ yields

[TABLE]

with $z^{\top}=(z_{1}^{\top},\dots,z_{R}^{\top})\in\mathbb{R}^{4|\mathcal{N}|}$, (6b) contains equations (5) by appropriate choice of coupling matrices $A_{i}\in\mathbb{R}^{2|\mathcal{A}|\times 4|\mathcal{N}|}$ and $\lambda$ denotes Lagrange multipliers assigned to (6b). Subsequently we denote $\mathcal{F}_{\mathcal{N}_{i},\mathcal{L}_{i}}$ as $\mathcal{F}_{i}$ for simplicity. The equality constraint $\mathcal{H}_{i}$ collects the power flow equations $\eqref{eq:PFeq}$ for all $i\in\mathcal{R}$.

IV Distributed Optimization Algorithm

This section introduces a Gauss-Newton Hessian approximation based variant of the distributed optimization algorithm ALADIN [14], which is tailored to AC PSSE problems in nonlinear least-squares form.

IV-A Main Algorithmic Steps

Algorithm 1 outlines a variant of ALADIN for solving (6). Similar to the traditional ALADIN procedure, there are two main steps: A parallelizable Step 1) and a consensus Step 3). In Step 1), decoupled NLPs (7a) are solved followed by a sensitivity evaluation (7b)—both in parallel. Note that due to the Gauss-Newton Hessian approximation, we only need to compute first-order derivatives. This way the computational burden and communication overhead is reduced significantly compared with standard ALADIN.

IV-B Local Convergence Analysis

Let $(z^{*},\lambda^{*},\kappa^{*})$ denote a primal-dual locally optimal solution of (6), where $\lambda^{\star}$ denotes the multiplier of the linear coupling constraints and $\kappa^{\star}$ the multiplier that is associated with the power-flow equations. In the following, we introduce the following regularity assumption on the physical power flow equations.

Assumption 1

The Jacobian of the power flow equations (2) with respect to all states $x$ of the network at $z^{\star}$ has full row-rank.

Notice that a detailed discussion of mathematical conditions under which the linear inpendendence constraint qualification (LICQ) condition in Assumption 1 is satisfied for power flow networks can be found in [12], where it is also discussed why this assumption is essentially satisfied for all power-flow networks of practical interest.

Proposition 1

If Assumption 1 holds, then the LICQ condition for the decoupled NLPs (7a) as well as the coupled QP (7d) is satisfied, that is, the matrix $[A^{\intercal}\;\;C^{\intercal}]^{\intercal}$ has full row rank.

*Proof.  *The proof of the proposition follows from the fact that the Jacobian of the consensus constraint has—by definition—full rank, as this constraint enforces linear coupling between neighboring regions. Because the power-flow equations are local (decoupled) in the reformulated problem (6), they satisfy the decoupled LICQ conditions (since Assumption 1 holds), and, additionally, cannot possibly be redundant to the coupling constraints. $\blacksquare$

In order to further ensure that any local solution $(z^{*},\lambda^{*},\kappa^{*})$ is a regular KKT point of (6), the following proposition is introduced.

Proposition 2

If the residuum $\sum_{i}\|\mathcal{F}_{i}(z_{i}^{\star})\|_{2}^{2}$ in the optimal solution is sufficiently small and if the matrices $\Sigma_{k}$ are positive definite, then the second order sufficient optimality condition (SOSC) is satisfied for (6) at $z^{\star}$ and the Gauss-Newton Hessian approximation, $\nabla\mathcal{F}_{i}(z_{i})\nabla\mathcal{F}_{i}(z_{i})^{\top}\succ 0$ is positive definite in a local neighborhood of $z^{\star}$.

*Proof.  *The statement of the above proposition is well-known in the context of Gauss-Newton SQP methods and a formal proof can be found in [3]. We remark that the conditions therein are indeed satisfied if the matrices $\Sigma_{k}$ are positive definite, as this condition trivially ensures identifiability of all measured states. $\blacksquare$

Note that the conditions of the above proposition are satisfied in practice if the model-data mismatch is small—but it can be violated otherwise.

Theorem 1

Let Assumption 1 be satisfied and let the residuum $\sum_{i}\|\mathcal{F}_{i}(z_{i}^{\star})\|_{2}^{2}$ at the local minimizer be sufficiently small such that Propsition 2 is applicable. Then the iterates $(z,\lambda)$ locally converge to $(z^{*},\lambda^{*})$ achieving a locally linear convergence rate.

*Proof.  *Propositions 1 and 2 ensure that minimizers of the decoupled NLPs (7a) are regular KKT points in a neighborhood of the optimal solution. Hence we can apply Lemma 3 in [14] to show that the solution $(y,\kappa^{d})$ of the decoupled NLP satisfies

[TABLE]

for a constant $\alpha<\infty$. Furthermore, in [14] it has been shown that the consensus step of ALADIN is locally equivalent to one SQP iteration. Thus, as we employ a Gauss-Newton Hessian approximation, we have

[TABLE]

Next, recall that the linear convergence rate of Gauss-Newton methods is locally proportional to the least-squares residuum at the optimal solution. In other words, we have $\gamma=\mathbf{O}(\sum_{i}\|\mathcal{F}_{i}(z_{i}^{\star})\|_{2}^{2})$, as proven in [3]. Thus, as long as the least-squares residuum is sufficiently small, it holds that

[TABLE]

with $\gamma(\alpha+1)<1$. This finishes the proof. $\blacksquare$

Note that Theorem 1 establishes local convergence of Algorithm 1 only. Thus, if one has poor initial guesses for the state, the proposed method needs to be augmented by a globalization routine, as discussed in [14]. This is subject to future work.

IV-C Communication Overhead

Step 3 of Algorithm 1 communicates between different regions. The forward communication collects matrices $B_{i}^{\top}B_{i}$, $C_{i}$ and vectors $B_{i}b_{i}$, $A_{i}y_{i}$ such that there are in total

[TABLE]

floats that need to be uploaded. The download phase after (7d) is solved sends the dual update $\lambda^{+}$ and local direction $\Delta y_{i}$ to each region, which requires $2|\mathcal{A}|+4|\mathcal{N}_{i}|$ floats in total for each region $i\in\mathcal{R}$.

V Numerical Example

In this section, we illustrate the performance of Algorithm 1 drawing upon the 30-bus system shown in Figure 1.

V-A Implementation and Data

The problem data is obtained from the MATPOWER dataset [33], although in our case study shunt elements are neglected. The system is partitioned into four regions $\mathcal{R}=\{1,2,3,4\}$ which are linked by $|\mathcal{A}|=8$ pairs of auxiliary nodes. We use nodal measurements and line measurements for all original nodes $k\in\mathcal{N}^{0}$ and all lines connecting original nodes $(k,l)\in\mathcal{N}_{i}\setminus\mathcal{A}_{i}\times\mathcal{N}_{i}\setminus\mathcal{A}_{i}$ and all $i\in\mathcal{R}$.

The measurements in our case study have been obtained by running a realistic scenario simulation in MATPOWER. During this simulation, we have introduced an additional Gaussian white noise with zero mean and a relative error variance of $10^{-4}$ the states $\theta_{k},p_{k},q_{k}$. In addition, the relative variance of the noise added to the voltage magnitude $v_{k}$ has been set to $10^{-5}$. Notice that such a white noise as been added for all $k\in\mathcal{N}^{0}$. The flows over transmission lines $p_{k,l},q_{k,l}$ and $i_{k,l}$ are subject to a relative variance of $10^{-5}$, which is a standard value that is often used in the context of PSSE [30]. The considered associated weighting matrices of the least-squares objective are

[TABLE]

for all nodes $k\in\mathcal{N}^{0}$ and

[TABLE]

for all $(k,l)\in\mathcal{N}_{i}\setminus\mathcal{A}_{i}\times\mathcal{N}_{i}\setminus\mathcal{A}_{i}$. It can be checked numerically that this choice ensures that the local convergence conditions of Theorem 1 are indeed satisfied.

V-B Numerical Comparison of ADMM and ALADIN

The implementation of Algorithm 1 relies on Casadi-v3.4.5 with IPOPT and MATLAB 2018a. The tuning parameters in Algorithm 1 are set to $\rho=10^{4}$ and $\epsilon=10^{-4}$. Moroever, in order to assess the numerical performance of the proposed Gauss-Newton ALADIN algorithm, we compare our implementation with a standard implementation of ADMM, where the augmented Lagrangian parameter is set to $\rho^{\text{ADM}}=10^{4}$, too. Note that ADMM does provide convergence guarantees for general non-convex problems. Indeed counter-examples where ADMM is divergent are given in [14]. However, it turns out that both ADMM and Gauss-Newton ALADIN converge for this particular PSSE case study.

Figure 3 shows the convergence of voltage angles, voltage magnitudes and active/reactive powers over the iteration index $k$. One can observe that ALADIN converges at a fast linear rate while ADMM converges quite slower to lower accuracy. Figure 4 shows the convergence of the corresponding consensus constraint violation $\|Ax^{k}-b\|_{\infty}$, which can be interpreted as the degree of matching of the voltage angles and magnitudes at auxiliary nodes according to (5). The fast linear convergence of ALADIN can also here be witnessed; that is, at least for this PSSE problem, our numerical results indicate a much better convergence of ALADIN when compared to ADMM. However, as shown in [8], one should keep in mind that ALADIN has a higher per-step communication overhead and complexity compared to ADMM.

Remark 2** **(Effects of scaling on convergence)

It is well-known that ADMM is rather scaling dependent [5]. In our numerical implementation all variables have already been scaled before running the ADMM routine, but, of course, we cannot exclude that it is possible to further improve the performance of ADMM by developing more sophisticated scaling or pre-conditioning strategies. However, one of the key advantages of Algorithm 1 versus ADMM is that it works “out of the box”; that is, there is no pre-conditioning or scaling needed, as Gauss-Newton methods are naturally invariant with respect to scaling [3, 14].

V-C A Posteriori Error Analysis

As for Bayesian inference or maximum likelihood estimation problem, there arises also in PSSE the question what can be said about the quality of the a posteriori distribution of the parameter estimate. At this point, we rely on a mature body of literature on nonlinear parameter estimation theory as reviewed in [18]. Therein, it has been proven that the inverse of the Fisher information matrix of a nonlinear least-squares parameter estimation problem is as a lower bound of the a-posteriori parameter estimation variance matrix. This relation is also known as Cramér-Rao bound. We refer to [18, 29] for an in-depth discussion and further references. Note that the inverse Fisher information matrix of the state estimate of the $i$-th subregion is in our context given by

[TABLE]

where $B_{\mathcal{R}}=\mathrm{diag}(B_{i})_{i\in\mathcal{R}}$ denotes the derivative of $\mathcal{F}$ while $A_{\mathcal{R}}$ and $C_{\mathcal{R}}$ denote the associated constraint Jacobian matrices, all evaluated a-posteriori at the optimal solution. A detailed derivation of this expression for the inverse Fisher information matrix in the context of constrained Gauss-Newton methods can be found in [4], see also [29]. The square-roots of selected diagonal elements of the above matrix relative to the nominal value of the associated parameter estimate can be found in Table I.

As there are $30$ busses in total not all values are shown. However, the last line lists the average values for the standard deviation over the whole network. Since most of the relative errors are below $1\%$, it can certainly be claimed that no over-fitting effects are visible in our PSSE case study. Nevertheless, some of the a-posteriori standard deviations of the reactive power estimates are around $14\%$, which indicates that some of the states of this 30-bus power system are rather difficult to estimate from the measurement data. In fact, these results suggest that a more detailed analysis of the parameter estimation accuracy in PSSE, as well as the optimization of sensor locations in power grids are an interesting future research direction.

VI Conclusion & Outlook

This work has introduced a distributed state estimation algorithm for non-convex AC PSSE problems based on ALADIN and a Gauss-Newton Hessian approximation. A local convergence condition for this algorithm has been given in Theorem 1. Moreover, we have illustrated the promising convergence behavior of Gauss-Newton compared to state-of-the-art ADMM methods by analyzing a highly non-trivial IEEE $30$-bus power grid.