Control and Optimization in Smart Power Grids INCITE Seminar @ - PowerPoint PPT Presentation

National & international redispatch • unforeseen congestion Cost of ancillary services of German TSOs in mio. Euros or voltage problems 111.8 82.3 primary frequency 85.2 • manually re-dispatched control reserves 103.4 110.9 on a 15-minute timescale 371.9 267.1 secondary frequency 352.9 control reserves 227.6 154.8 104.2 67.4 Redispatch actions in the German tertiary frequency 15 811 156.1 transmission grid control reserves 106.0 50.2 in hours 27.0 68.3 reactive power 33.0 26.7 8 453 7 965 32.6 7 160 41.6 164.8 5 030 national & internat. 113.3 redispatch 185.4 411.9 1 588 2011 2012 2013 2014 2015 2010 2011 2012 2013 2014 2015 [Bundesnetzagentur, Monitoringbericht 2016] [Bundesnetzagentur, Monitoringbericht 2016] 15

Proposal: online optimization in closed loop prediction (load, generation) δ short-term real-time low-level dynamic scheduling planning operations controllers model optimization stage steady-state model 16

Proposal: online optimization in closed loop prediction (load, generation) δ short-term real-time low-level dynamic scheduling planning operations controllers model optimization stage steady-state model combining optimization & feedback control for real-time operation robust (feedback strategy) steady-state optimality fast response satisfaction of operational constraints disclaimer: no predictive optimization (only for static systems) focus today on real-time (no distributed) aspects 16

Proposal: online optimization in closed loop prediction (load, generation) δ short-term real-time low-level dynamic scheduling planning operations controllers model optimization stage steady-state model combining optimization & feedback control for real-time operation robust (feedback strategy) steady-state optimality fast response satisfaction of operational constraints disclaimer: no predictive optimization (only for static systems) focus today on real-time (no distributed) aspects lots of related work: [Bolognani et. al, A Survey of Distributed Optimization and Control 2015], [Dall’Anese and Simmonetto, Algorithms for Electric Power Systems Daniel K. Molzahn, ∗ Member, IEEE , Florian D¨ orfler, † Member, IEEE , Henrik Sandberg, ‡ Member, IEEE , 2016], [Gan and Low, 2016], ... Steven H. Low, § Fellow, IEEE , Sambuddha Chakrabarti, ¶ Student Member, IEEE , Ross Baldick, ¶ Fellow, IEEE , and Javad Lavaei, ∗∗ Member, IEEE 16

OVERVIEW 1. The power flow manifold, representations, and approximations 2. Projected gradient flow on the power flow manifold 3. Tracking performance and robustness of closed-loop optimization 4. Output feedback and state uncertainty 17

THE POWER FLOW MANIFOLD, REPRESENTATIONS, AND APPROXIMATIONS 18

Steady-state AC power flow model quasi-stationary dynamics → complex impedances and voltages sources: locally controlled → buses are PQ or PV or slack V θ loads: constant impedance, current, or PQ power (today) Ohm’s Law Current Law 12 13 11 8 7 10 9 6 3 AC power 4 2 AC power flow equations 5 line impedance nodal voltage line current current injection power flow power injections (all variables and parameters are -valued) 19

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization 20

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization • rectangular form: replace V k = e k + jf k and split real & imaginary parts → real-valued quadratic and useful for homotopy methods & QCQPs 20

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization • rectangular form: replace V k = e k + jf k and split real & imaginary parts → real-valued quadratic and useful for homotopy methods & QCQPs • matrix form: replace W kl = V k · V ∗ l where W is unit-rank p.s.d. Hermitian matrix → linear and useful for relaxations in convex optimization problems 20

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization • rectangular form: replace V k = e k + jf k and split real & imaginary parts → real-valued quadratic and useful for homotopy methods & QCQPs • matrix form: replace W kl = V k · V ∗ l where W is unit-rank p.s.d. Hermitian matrix → linear and useful for relaxations in convex optimization problems • polar form: replace V k = | V k | e j θ k and split real / imaginary parts → this is how power system engineers think: all specs on | V k | and d dt θ k 20

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization • rectangular form: replace V k = e k + jf k and split real & imaginary parts → real-valued quadratic and useful for homotopy methods & QCQPs • matrix form: replace W kl = V k · V ∗ l where W is unit-rank p.s.d. Hermitian matrix → linear and useful for relaxations in convex optimization problems • polar form: replace V k = | V k | e j θ k and split real / imaginary parts → this is how power system engineers think: all specs on | V k | and d dt θ k • branch flow: parameterized in flows: I k → l = y kl ( V k − V l ) and S k → l = V k I ∗ k → l → useful in radial networks: equations can be expressed in magnitudes only 20

Power flow representations • complex form: S k = P k + jQ k = � l ∈ N ( k ) y ∗ kl V k · ( V ∗ k − V ∗ l ) where y kl = 1 / z kl → complex-valued quadratic and useful for calculations & optimization • rectangular form: replace V k = e k + jf k and split real & imaginary parts → real-valued quadratic and useful for homotopy methods & QCQPs • matrix form: replace W kl = V k · V ∗ l where W is unit-rank p.s.d. Hermitian matrix → linear and useful for relaxations in convex optimization problems • polar form: replace V k = | V k | e j θ k and split real / imaginary parts → this is how power system engineers think: all specs on | V k | and d dt θ k • branch flow: parameterized in flows: I k → l = y kl ( V k − V l ) and S k → l = V k I ∗ k → l → useful in radial networks: equations can be expressed in magnitudes only • many variations, coordinate changes, convexifications, etc. → some problems become easier in different coordinates 20

A brief history of power flow approximations for computational tractability, analytic studies, & control/optimization design • DC power flow : polar form → ℜ ( Z ) = 0 , | V | = 1 , and linearization B. Stott, J. Jardim, & O. Alsac, DC Power Flow Revisited. IEEE TPS , 2009. → standard (but often poor) approximation for transmission networks 21

A brief history of power flow approximations for computational tractability, analytic studies, & control/optimization design • DC power flow : polar form → ℜ ( Z ) = 0 , | V | = 1 , and linearization B. Stott, J. Jardim, & O. Alsac, DC Power Flow Revisited. IEEE TPS , 2009. → standard (but often poor) approximation for transmission networks • linear coupled power flow : polar form → linearization for small angles/voltages → preserves losses and angles/voltages cross-coupling: suited for distribution 21

A brief history of power flow approximations for computational tractability, analytic studies, & control/optimization design • DC power flow : polar form → ℜ ( Z ) = 0 , | V | = 1 , and linearization B. Stott, J. Jardim, & O. Alsac, DC Power Flow Revisited. IEEE TPS , 2009. → standard (but often poor) approximation for transmission networks • linear coupled power flow : polar form → linearization for small angles/voltages → preserves losses and angles/voltages cross-coupling: suited for distribution • LinDistFlow : branch flow → parameterization | V | 2 coordinates and linearization M.E. Baran & F.F. Wu, Optimal sizing of capacitors placed on a radial distribution system. PES , 1988. → very useful for voltages in (radial) distribution networks 21

A brief history of power flow approximations for computational tractability, analytic studies, & control/optimization design • DC power flow : polar form → ℜ ( Z ) = 0 , | V | = 1 , and linearization B. Stott, J. Jardim, & O. Alsac, DC Power Flow Revisited. IEEE TPS , 2009. → standard (but often poor) approximation for transmission networks • linear coupled power flow : polar form → linearization for small angles/voltages → preserves losses and angles/voltages cross-coupling: suited for distribution • LinDistFlow : branch flow → parameterization | V | 2 coordinates and linearization M.E. Baran & F.F. Wu, Optimal sizing of capacitors placed on a radial distribution system. PES , 1988. → very useful for voltages in (radial) distribution networks • rectangular DC power flow : fixed-point expansion for small S 2 / V 2 slack S. Bolognani & S. Zampieri, On the existence and linear approximation of the power flow solution in power distribution networks. IEEE TPS , 2015. → works amazingly well in distribution and transmission 21

A brief history of power flow approximations for computational tractability, analytic studies, & control/optimization design • DC power flow : polar form → ℜ ( Z ) = 0 , | V | = 1 , and linearization B. Stott, J. Jardim, & O. Alsac, DC Power Flow Revisited. IEEE TPS , 2009. → standard (but often poor) approximation for transmission networks • linear coupled power flow : polar form → linearization for small angles/voltages → preserves losses and angles/voltages cross-coupling: suited for distribution • LinDistFlow : branch flow → parameterization | V | 2 coordinates and linearization M.E. Baran & F.F. Wu, Optimal sizing of capacitors placed on a radial distribution system. PES , 1988. → very useful for voltages in (radial) distribution networks • rectangular DC power flow : fixed-point expansion for small S 2 / V 2 slack S. Bolognani & S. Zampieri, On the existence and linear approximation of the power flow solution in power distribution networks. IEEE TPS , 2015. → works amazingly well in distribution and transmission • many variations, extensions, sensitivity and Jacobian methods, etc. 21

A unifying geometric perspective: the power flow manifold 1.2 node 1 node 2 1 0.8 y = 0 . 4 − 0 . 8 j v 2 0.6 v 1 = 1 , θ 1 = 0 v 2 , θ 2 0.4 p 1 , q 1 p 2 , q 2 1 0.5 1 0 0.5 q 2 0 -0.5 p 2 22

A unifying geometric perspective: the power flow manifold 1.2 node 1 node 2 1 0.8 y = 0 . 4 − 0 . 8 j v 2 0.6 v 1 = 1 , θ 1 = 0 v 2 , θ 2 0.4 p 1 , q 1 p 2 , q 2 1 • variables: all of x = ( | V | , θ, P , Q ) 0.5 1 0 0.5 q 2 0 • power flow manifold: M = { x : h ( x ) = 0 } -0.5 p 2 → submanifold in R 2 n or R 6 n (3-phase) 22

A unifying geometric perspective: the power flow manifold 1.2 node 1 node 2 1 0.8 y = 0 . 4 − 0 . 8 j v 2 0.6 v 1 = 1 , θ 1 = 0 v 2 , θ 2 0.4 p 1 , q 1 p 2 , q 2 1 • variables: all of x = ( | V | , θ, P , Q ) 0.5 1 0 0.5 q 2 0 • power flow manifold: M = { x : h ( x ) = 0 } -0.5 p 2 → submanifold in R 2 n or R 6 n (3-phase) � ∂ h ( x ) x ∗ = A T • normal space spanned by � x ∗ ∂ x 1.4 1.2 • tangent space A x ∗ ( x − x ∗ ) = 0 1 v 2 1.5 → best linear approximant at x ∗ 0.8 1 0.6 0.5 -1 -0.5 0 0 -0.5 0.5 q 2 1 -1 p 2 1.5 22

A unifying geometric perspective: the power flow manifold 1.2 node 1 node 2 1 0.8 y = 0 . 4 − 0 . 8 j v 2 0.6 v 1 = 1 , θ 1 = 0 v 2 , θ 2 0.4 p 1 , q 1 p 2 , q 2 1 • variables: all of x = ( | V | , θ, P , Q ) 0.5 1 0 0.5 q 2 0 • power flow manifold: M = { x : h ( x ) = 0 } -0.5 p 2 → submanifold in R 2 n or R 6 n (3-phase) � ∂ h ( x ) x ∗ = A T • normal space spanned by � x ∗ ∂ x 1.4 1.2 • tangent space A x ∗ ( x − x ∗ ) = 0 1 v 2 1.5 → best linear approximant at x ∗ 0.8 1 0.6 0.5 -1 -0.5 0 ∂ 2 h ( x ) • accuracy depends on curvature 0 -0.5 0.5 ∂ x 2 q 2 1 -1 p 2 1.5 → constant in rectangular coordinates 22

Accuracy illustrated with unbalanced three-phase IEEE13 ◦ exact solution ⋆ linear approximant Matlab/Octave code @ https://github.com/saveriob/1ACPF 23

Special cases reveal some old friends • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) 24

Special cases reveal some old friends • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) � � � � � � ℜ ( Y ) −ℑ ( Y ) | V | P ⇒ tangent space parameterization: = −ℑ ( Y ) ℜ ( Y ) θ Q gives linear coupled power flow [D. Deka, S. Backhaus, and M. Chertkov, ’15] 24

Special cases reveal some old friends • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) � � � � � � ℜ ( Y ) −ℑ ( Y ) | V | P ⇒ tangent space parameterization: = −ℑ ( Y ) ℜ ( Y ) θ Q gives linear coupled power flow [D. Deka, S. Backhaus, and M. Chertkov, ’15] ⇒ ℜ ( Y ) = 0 gives DC power flow(s): −ℑ ( Y ) θ = P and −ℑ ( Y ) E = Q 24

Special cases reveal some old friends • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) � � � � � � ℜ ( Y ) −ℑ ( Y ) | V | P ⇒ tangent space parameterization: = −ℑ ( Y ) ℜ ( Y ) θ Q gives linear coupled power flow [D. Deka, S. Backhaus, and M. Chertkov, ’15] ⇒ ℜ ( Y ) = 0 gives DC power flow(s): −ℑ ( Y ) θ = P and −ℑ ( Y ) E = Q 1.5 1 power flow manifold 0.5 linear coupled power flow p 2 0 DC power flow approximation -0.5 (neglects PV coupling) -1 2 0.6 1 0.8 0 1 -1 1.2 v 2 1.4 -2 ! 2 24

Special cases reveal some old friends cont’d • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) 25

Special cases reveal some old friends cont’d • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) ⇒ rectangular coordinates ⇒ rectangular DC flow [S. Bolognani and S. Zampieri, ’15] 25

Special cases reveal some old friends cont’d • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) ⇒ rectangular coordinates ⇒ rectangular DC flow [S. Bolognani and S. Zampieri, ’15] • nonlinear change to quadratic coordinates from | V k | to | V k | 2 ⇒ linearization gives (non-radial) LinDistFlow [M.E. Baran and F.F. Wu, ’88] 25

Special cases reveal some old friends cont’d • flat-voltage/0-injection point: x ∗ = ( | V | ∗ , θ ∗ , P ∗ , Q ∗ ) = ( 1 , 0 , 0 , 0 ) ⇒ rectangular coordinates ⇒ rectangular DC flow [S. Bolognani and S. Zampieri, ’15] • nonlinear change to quadratic coordinates from | V k | to | V k | 2 ⇒ linearization gives (non-radial) LinDistFlow [M.E. Baran and F.F. Wu, ’88] 1.4 1.3 power flow manifold 1.2 1.1 linear approximation 1 v 2 0.9 linear approximation 0.8 in quadratic coordinates 0.7 0.6 -1 1.5 0 1 0.5 1 0 -0.5 2 -1 25 p 2 q 2

Properties of power flow manifold that we will exploit 1.2 nonlinear power flow is smooth manifold 1 → coordinate-independent – no singularities 0.8 v 2 0.6 → better local linear approximations 0.4 → methods for manifold optimization/control 1 0.5 natural concept for closed-loop dynamics 1 0 0.5 q 2 0 -0.5 → M is attractive for grid dynamics p 2 → closed-loop trajectories x ( t ) live on M → control task: steer ˙ x ( t ) in tangent space 1.4 1.2 v 2 1 const.-rank linearization A x ∗ ( x − x ∗ ) = 0 1.5 0.8 1 → implicit – no input/outputs (no disadvantage) 0.6 0.5 -1 -0.5 0 → sparse – A x ∗ has the sparsity of the grid 0 -0.5 0.5 q 2 1 p 2 -1 1.5 → structure – elements of A x ∗ are local → S. Bolognani & F. Dörfler (2015) “Fast power system analysis via implicit linearization of the power flow manifold” 26

PROJECTED GRADIENT FLOW ON THE POWER FLOW MANIFOLD 27

AC power flow model, constraints, and objectives model (physical constraint): x ∈ M Ohm’s Law Current Law 12 13 11 8 7 10 9 AC power 6 3 4 AC power flow equations 2 5 line impedance nodal voltage line current current injection power injections power flow (all variables and parameters are -valued) operational constraints: generation capacity, voltage bands, no congestion objective: economic dispatch, minimize losses, distance to collapse, etc. control: state measurements and actuation via generator set-points 28

Ancillary services as a real-time OPF Real-time optimal power flow (OPF) � cost k ( P G • minimize cost of generation minimize k ) k ∈N P G + jQ G = P L + jQ L + diag( V ) Y ∗ V ∗ • satisfy AC power flow laws subject to P k ≤ P G k ≤ P k , Q k ≤ Q G • respect generation capacity k ≤ Q k ∀ k ∈ N • no over-/under-voltage V k ≤ | V k | ≤ V k ∀ k ∈ N • no congestion | P kl + jQ kl | ≤ S kl ∀ ( k , l ) ∈ E Y admittance matrix, P G k , Q G k power generation, P L k , Q L k load, { V k , V k , . . . } nodal limits, S kl line flow limit 29

Ancillary services as a real-time OPF Real-time optimal power flow (OPF) � cost k ( P G • minimize cost of generation minimize k ) k ∈N P G + jQ G = P L + jQ L + diag( V ) Y ∗ V ∗ • satisfy AC power flow laws subject to P k ≤ P G k ≤ P k , Q k ≤ Q G • respect generation capacity k ≤ Q k ∀ k ∈ N • no over-/under-voltage V k ≤ | V k | ≤ V k ∀ k ∈ N • no congestion | P kl + jQ kl | ≤ S kl ∀ ( k , l ) ∈ E Y admittance matrix, P G k , Q G k power generation, P L k , Q L k load, { V k , V k , . . . } nodal limits, S kl line flow limit Loads P L , Q L generator setpoints physical, steady-state power system Real-time (AC power flow equations) operation P G = P L + ℜ{ diag( V ) Y ∗ V ∗ } Q G = Q L + ℑ{ diag( V ) Y ∗ V ∗ } state measurements 29

Ancillary services as a real-time OPF Real-time optimal power flow (OPF) � cost k ( P G • minimize cost of generation minimize k ) k ∈N P G + jQ G = P L + jQ L + diag( V ) Y ∗ V ∗ • satisfy AC power flow laws subject to P k ≤ P G k ≤ P k , Q k ≤ Q G • respect generation capacity k ≤ Q k ∀ k ∈ N • no over-/under-voltage V k ≤ | V k | ≤ V k ∀ k ∈ N • no congestion | P kl + jQ kl | ≤ S kl ∀ ( k , l ) ∈ E Y admittance matrix, P G k , Q G k power generation, P L k , Q L k load, { V k , V k , . . . } nodal limits, S kl line flow limit Loads A control problem with P L , Q L challenging specifications on the closed-loop system: generator setpoints physical, steady-state 1. its trajectory x ( t ) must satisfy power system Real-time (AC power flow equations) the constraints at all times operation P G = P L + ℜ{ diag( V ) Y ∗ V ∗ } Q G = Q L + ℑ{ diag( V ) Y ∗ V ∗ } 2. it must converge to x ⋆ , the state solution of the AC OPF measurements 29

Ancillary services as a real-time OPF Real-time optimal power flow (OPF) � cost k ( P G • minimize cost of generation minimize k ) k ∈N P G + jQ G = P L + jQ L + diag( V ) Y ∗ V ∗ • satisfy AC power flow laws subject to P k ≤ P G k ≤ P k , Q k ≤ Q G • respect generation capacity k ≤ Q k ∀ k ∈ N • no over-/under-voltage V k ≤ | V k | ≤ V k ∀ k ∈ N • no congestion | P kl + jQ kl | ≤ S kl ∀ ( k , l ) ∈ E Y admittance matrix, P G k , Q G k power generation, P L k , Q L k load, { V k , V k , . . . } nodal limits, S kl line flow limit x = � | V | θ P Q � Prototype of real-time OPF grid state φ : R n → R objective function minimize φ ( x ) M ⊂ R n AC power flow equations subject to x ∈ K = M ∩ X X ⊂ R n operational constraints 30

Unconstrained optimization on the power flow manifold geometric objects: manifold M = { x : h ( x ) = 0 } tangent space T x M = ker h ( x ) objective φ : M → R Riemann metric g : T x M × T x M → R (degree of freedom) 31

Unconstrained optimization on the power flow manifold geometric objects: manifold M = { x : h ( x ) = 0 } tangent space T x M = ker h ( x ) objective φ : M → R Riemann metric g : T x M × T x M → R (degree of freedom) target state: local minimizer on the power flow manifold x ⋆ ∈ arg min x ∈M φ ( x ) 31

Unconstrained optimization on the power flow manifold geometric objects: manifold M = { x : h ( x ) = 0 } tangent space T x M = ker h ( x ) objective φ : M → R Riemann metric g : T x M × T x M → R (degree of freedom) target state: local minimizer on the power flow manifold x ⋆ ∈ arg min x ∈M φ ( x ) always feasible due to physics: trajectory remains on power flow manifold M 31

Unconstrained optimization on the power flow manifold geometric objects: manifold M = { x : h ( x ) = 0 } tangent space T x M = ker h ( x ) objective φ : M → R Riemann metric g : T x M × T x M → R (degree of freedom) target state: local minimizer on the power flow manifold x ⋆ ∈ arg min x ∈M φ ( x ) always feasible due to physics: trajectory remains on power flow manifold M linear approximant continuous-time gradient descent on M : Gradient of cost function 1. grad φ ( x ): gradient of cost function Projected gradient (& soft constraints) in ambient space ˙ x x ( t ) 2. Π x grad φ ( x ): projection of gradient on the linear approximant T x M 3. flow on manifold: ˙ x = − γ Π x grad φ ( x ) power flow manifold 31

Constraints: projected dynamical systems for feasibility Operational constraints Per specification, the trajectories need to satisfy operational constraints at all times. x ( t ) ∈ K = M ∩ X where M power flow manifold X operational constraints → ˙ x ( t ) must belong to a feasible cone , subset of the tangent space of M precisely: ˙ x ( t ) ∈ T x K ⊂ T x M , the inward tangent cone at x 32

Constraints: projected dynamical systems for feasibility Operational constraints Per specification, the trajectories need to satisfy operational constraints at all times. x ( t ) ∈ K = M ∩ X where M power flow manifold F : R n → R n vector field, K ⊂ R n closed domain X operational constraints Projected dynamical systems: → ˙ x ( t ) must belong to a feasible cone , � x , F ( x ) � x = Π K ˙ subset of the tangent space of M where precisely: ˙ x ( t ) ∈ T x K ⊂ T x M , the inward tangent cone at x Π K ( x , F ( x )) ∈ arg min v ∈ T x K � F ( x ) − v � g 32

Projected gradient descent on the power flow manifold x = Π K ( x , − grad φ ( x )) , ˙ x (0) = x 0 33

Projected gradient descent on the power flow manifold x = Π K ( x , − grad φ ( x )) , ˙ x (0) = x 0 • Does a solution trajectory exist for a non-convex K ? Is it unique ? • Are solution trajectories (asymptotically) stable˙ ? • Do solution trajectories converge to a minimizer of φ ? 33

Projected gradient descent on the power flow manifold x = Π K ( x , − grad φ ( x )) , ˙ x (0) = x 0 • Does a solution trajectory exist for a non-convex K ? Is it unique ? • Are solution trajectories (asymptotically) stable˙ ? • Do solution trajectories converge to a minimizer of φ ? Corollary (simplified) Let x : [0 , ∞ ) → K be a (Carathéodory-)solution of the initial value problem x = Π K ( x , − grad φ ( x )) , ˙ x (0) = x 0 . If φ has compact level sets on K , x ( t ) will converge to a critical point x ⋆ of φ on K . Furthermore, if x ⋆ is asymptotically stable then it is a local minimizer of φ on K . → Hauswirth, Bolognani, Hug, & Dörfler (2016) “Projected gradient descent on Riemanniann manifolds with applications to online power system optimization” 33

How to induce the projected gradient flow actuate Controlled system u static minimize u , x φ ( x ) feedback system optimizer h ( x ) = 0 subject to x ∈ K g ( x ) = u g ( x ) = u x measure the state x is uniquely determined by – the algebraic model h ( x ) = 0 describing the power flow equations – an algebraic input constraint g ( x ) = u 34

How to induce the projected gradient flow actuate Controlled system u static minimize u , x φ ( x ) feedback system optimizer h ( x ) = 0 subject to x ∈ K g ( x ) = u g ( x ) = u x measure the state x is uniquely determined by – the algebraic model h ( x ) = 0 describing the power flow equations – an algebraic input constraint g ( x ) = u steady state: the closed-loop system converges to the solution of the OPF closed-loop trajectory remains in K at all times → no need to solve the optimization problem numerically → no need to solve any power flow equation 34

From projected gradient flow to discrete-time feedback control � � x exo partition: x = x endo exogenous variables: inputs/disturbances (e.g., reactive injection Q k ) endogenous variables: determined by the physics (e.g., voltage V k ) 35

From projected gradient flow to discrete-time feedback control linear approximant � � x exo Gradient of cost function partition: x = x endo Projected gradient exogenous variables: x ( t ) inputs/disturbances (e.g., reactive injection Q k ) Retraction x ( t + 1) endogenous variables: determined by the physics power flow manifold (e.g., voltage V k ) 1. compute continuous feasible descent direction : d t = Π K ( x , − grad φ ( x ( t ))) x ( t + 1) = x ( t ) + α · d t 2. Euler integration step to compute new set-points : ˜ 3. actuate exogeneous variables (inputs) based on ˜ x endo ( t + 1) (note: x exo will be updated accordingly since h ( x ) = 0 holds implicitly by physics) 4. retraction step x ( t + 1) = R x ( t ) (˜ x ( t + 1)) ⇒ x ( t + 1) ∈ M (note: carried out by physics since M is attractive / use AC PF solver in simulations) 35

Simple illustrative case study actuate u static feedback system optimizer h ( x ) = 0 g ( x ) = u x measure Objective Value [$] 10 real time cost 5 global minimum 0 0 50 100 150 200 250 300 Bus voltages [p.u.] 1.05 1 0.95 0 50 100 150 200 250 300 Active power generation [MW] 2 Slack bus Gen A Gen B 1 0 0 50 100 150 200 250 300 36 iteration

TRACKING PERFORMANCE AND ROBUSTNESS OF CLOSED-LOOP OPTIMIZATION 37

The tracking problem the power system state is also affected by exogeneous inputs w t → because of these inputs, the state could leave the feasible region K → outside of K , the projected gradient flow is not defined 38

The tracking problem the power system state is also affected by exogeneous inputs w t → because of these inputs, the state could leave the feasible region K → outside of K , the projected gradient flow is not defined 38

The tracking problem the power system state is also affected by exogeneous inputs w t → because of these inputs, the state could leave the feasible region K → outside of K , the projected gradient flow is not defined w t U u static system feedback h ( x , w t ) = 0 optimizer g ( x ) = u x constraints satisfaction for non-controllable variables: K accounts only for hard constraints on controllable variables u (e.g., generation limits) gradient projection becomes input saturation (saturated proportional feedback control) soft constraints included via penalty functions in φ (e.g., thermal and voltage limits) → alternative method (not discussed today) is dualization (i.e., integral control) 38

Tracking performance Bus voltages [p.u.] G Generator 1.1 Synchronous Condensor C 1.05 S Solar Wind W W 1 0.95 G 1 0.9 S C 3 Branch current magnitudes [p.u.] C 2 1 G 2 C 1 0.5 0 Aggregate Load & Available Renewable Power [MW] 400 Active power injection [MW] Load Solar Wind 200 300 150 200 100 100 50 0 0 0 2 4 6 8 10 12 14 16 18 20 22 24 0 2 4 6 8 10 12 14 16 18 20 22 24 Time [hrs] Time [hrs] Gen1 Gen2 Solar Wind controller: penalty + saturation → Hauswirth, Bolognani, Dörfler, & Hug (2017) “Online Optimization in Closed Loop on the Power Flow Manifold” 39

Tracking performance Generator G Comparison C Synchronous Condensor S Solar Wind W W closed-loop feedback trajectory benchmark: feedforward OPF G 1 S C 3 C 2 (solution of an ideal OPF without computation delay) G 2 C 1 Generation cost Feedback OPF 1 , 500 Optimal cost practically exact tracking 1 , 000 + trajectory feasibility 500 + robustness to model mismatch 0 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 Time [hrs] 40

Trajectory feasibility The feasible region K = M ∩ X often has disconnected components . x ∗ K x 0 M 41

Trajectory feasibility The feasible region K = M ∩ X often has disconnected components . x ∗ K x 0 M feedback (gradient descent) → the closed-loop trajectory x ( t ) is guaranteed to be feasible → convergence of x ( t ) to a local minimum is guaranteed feedforward (OPF) – optimizer x ⋆ = arg min x ∈K φ ( x ) can be in different disconnected component → no feasible trajectory exists: x 0 → x ⋆ must violate constraints 41

Illustration of trajectory feasibility 5-bus example known to have two Objective Value [$] 1200 disconnected feasible regions: Feedback Feed-forward 1000 2 0.55+j0.90 0.04+0.09 800 0 1000 2000 3000 4000 5000 4 Voltage Levels [p.u.] 0.07+j0.09 1.05 G 0.06+j0.1 1 1 0.95 0 5 1 0 1000 2000 3000 4000 5000 . G 0 j + 5 Active Power Generation P [MW] 0 . 0.55+j0.90 0 300 200 3 Gen1 100 Gen2 0 1000 2000 3000 4000 5000 [0s,2000s]: separate feasible regions Reactive Power Generation Q [MVAR] 200 [2000s,3000s]: loosen limits on Gen1 Gen2 100 reactive power Q 2 → regions merge Q5min 0 [4000s,5000s]: tighten limits on Q 2 0 1000 2000 3000 4000 5000 → vanishing feasible region 42

Robustness to model mismatch Intuition in 2D case: cost on x 1 , soft penalty for constraint x 2 ≤ ¯ x 2 , actuation on x 1 x 2 ¯ x 2 x 1 = u 43

Robustness to model mismatch Intuition in 2D case: cost on x 1 , soft penalty for constraint x 2 ≤ ¯ x 2 , actuation on x 1 x 2 x 2 ¯ ¯ x 2 x 2 u ∗ x 1 = u x 1 = u ↑ feedforward (OPF) model-based approach: model mismatch directly affects the decision u ⋆ 43

Robustness to model mismatch Intuition in 2D case: cost on x 1 , soft penalty for constraint x 2 ≤ ¯ x 2 , actuation on x 1 x 2 x 2 ¯ ¯ x 2 x 2 u ∗ x 1 = u x 1 = u x 2 ↑ feedforward (OPF) model-based approach: model ¯ mismatch directly affects the decision u ⋆ x 2 ← feedback (gradient descent) grad φ is orthogonal to the tangent plane x 1 = u 43

Illustration of robustness to model mismatch IEEE 30-bus test system w t controller: G Generator Synchronous Condensor C U saturation of Solar S u Wind W W generation static system constraints feedback h ( x , w t ) = 0 G 1 optimizer S g ( x ) = u C 3 penalty for C 2 operational x G 2 constraints C 1 no automatic re-dispatch feedback optimization f − f ∗ � v − v ∗ � f − f ∗ � v − v ∗ � model uncertainty feasible ? feasible ? loads ± 40% no 94.6 0.03 yes 0.0 0.0 line params ± 20% yes 0.19 0.01 yes 0.01 0.003 2 line failures no -0.12 0.06 yes 0.19 0.007 44

Illustration of robustness to model mismatch IEEE 30-bus test system w t controller: G Generator Synchronous Condensor C U saturation of Solar S u Wind W W generation static system constraints feedback h ( x , w t ) = 0 G 1 optimizer S g ( x ) = u C 3 penalty for C 2 operational x G 2 constraints C 1 no automatic re-dispatch feedback optimization f − f ∗ � v − v ∗ � f − f ∗ � v − v ∗ � model uncertainty feasible ? feasible ? loads ± 40% no 94.6 0.03 yes 0.0 0.0 line params ± 20% yes 0.19 0.01 yes 0.01 0.003 2 line failures no -0.12 0.06 yes 0.19 0.007 on-going work: observations can be made mathematically rigorous and quantified 44

OUTPUT FEEDBACK AND STATE UNCERTAINTY 45

Use real-time output measurements to reduce uncertainty w actuate u static feedback system optimizer h ( x , w ) = 0 g ( x ) = u y = y ( x ) output How to project the trajectory to K = M ∩ X when the state is partially known ? power flow manifold M : attractive manifold + robustness � operational constraints X : how to deal with state uncertainty ? 46

Use real-time output measurements to reduce uncertainty w actuate u static feedback system optimizer h ( x , w ) = 0 g ( x ) = u y = y ( x ) output How to project the trajectory to K = M ∩ X when the state is partially known ? power flow manifold M : attractive manifold + robustness � operational constraints X : how to deal with state uncertainty ? Chance constraints generally non-convex set of all u such that P [ x ∈ X w | y ( x ) = y ] ≥ 1 − ǫ where w is random and ǫ ∈ (0 , 1) is probability of constrained violation 46

Scenario approach to chance-constrained optimization chance constraint: P [ x ∈ X w ] ≥ 1 − ǫ where w is random and ǫ ∈ (0 , 1) → often intractable for complex (possibly unknown) distributions/constraints 47

Scenario approach to chance-constrained optimization chance constraint: P [ x ∈ X w ] ≥ 1 − ǫ where w is random and ǫ ∈ (0 , 1) → often intractable for complex (possibly unknown) distributions/constraints sample from distribution → deterministic constraints x ∈ X w ( i ) , i ∈ { 1 , . . . , N } convert stochastic constraint to large set of deterministic ones: X w ≈ � N i =1 X w ( i ) → # samples to approximate chance constraint depends on n , ε , and accuracy 47

Scenario approach to chance-constrained optimization chance constraint: P [ x ∈ X w ] ≥ 1 − ǫ where w is random and ǫ ∈ (0 , 1) → often intractable for complex (possibly unknown) distributions/constraints sample from distribution → deterministic constraints x ∈ X w ( i ) , i ∈ { 1 , . . . , N } convert stochastic constraint to large set of deterministic ones: X w ≈ � N i =1 X w ( i ) → # samples to approximate chance constraint depends on n , ε , and accuracy IEEE 13 grid with random demand and actuation (microgenerators & tap changers) feasible region with scenario approach 47

Scenario approach with real-time measurements scenario approach: stochastic constraint → large set of deterministic ones P w [ x ∈ X w ] ≥ 1 − ǫ → x ∈ X w ( i ) , i ∈ { 1 , . . . , N } 8 6 two sources of information on the unknown w 1. historical samples w ( i ) of prior distribution 4 → classic scenario-based approach 2 0 − 2 − 4 − 2 0 2 4 6 48

Scenario approach with real-time measurements scenario approach: stochastic constraint → large set of deterministic ones P w [ x ∈ X w ] ≥ 1 − ǫ → x ∈ X w ( i ) , i ∈ { 1 , . . . , N } 8 6 two sources of information on the unknown w 1. historical samples w ( i ) of prior distribution 4 → classic scenario-based approach 2 2. online measurements y from the system 0 − 2 − 4 − 2 0 2 4 6 48

Scenario approach with real-time measurements scenario approach: stochastic constraint → large set of deterministic ones P w [ x ∈ X w ] ≥ 1 − ǫ → x ∈ X w ( i ) , i ∈ { 1 , . . . , N } 8 6 two sources of information on the unknown w 1. historical samples w ( i ) of prior distribution 4 → classic scenario-based approach 2 2. online measurements y from the system 0 → use measurements to reduce uncertainty? − 2 − 4 − 2 0 2 4 6 48

Scenario approach with real-time measurements scenario approach: stochastic constraint → large set of deterministic ones P w [ x ∈ X w ] ≥ 1 − ǫ → x ∈ X w ( i ) , i ∈ { 1 , . . . , N } 8 6 two sources of information on the unknown w 1. historical samples w ( i ) of prior distribution 4 → classic scenario-based approach 2 2. online measurements y from the system 0 → use measurements to reduce uncertainty? − 2 − 4 − 2 0 2 4 6 re-sampling solution: scenario approach based on conditional distribution → high computational demand, large memory footprint, not suited for real time 48

Scenario approach with real-time measurements scenario approach: stochastic constraint → large set of deterministic ones P w [ x ∈ X w ] ≥ 1 − ǫ → x ∈ X w ( i ) , i ∈ { 1 , . . . , N } 8 6 two sources of information on the unknown w 1. historical samples w ( i ) of prior distribution 4 → classic scenario-based approach 2 2. online measurements y from the system 0 → use measurements to reduce uncertainty? − 2 − 4 − 2 0 2 4 6 re-sampling solution: scenario approach based on conditional distribution → high computational demand, large memory footprint, not suited for real time today: online computation of posterior distribution after measurement 48

Linear case linear grid model x = Au + Bw polytopic constraints Cx ≤ z linear measurement y = Hw 49

Linear case linear grid model x = Au + Bw polytopic constraints Cx ≤ z linear measurement y = Hw Approximate conditioning affine transformation: w y = w + K ( y − Hw ) ˆ where K = Σ H ⊤ � H Σ H ⊤ � − 1 → projection of uncertainty in the subspace { y = Hw } → Gaussian case: recovers the conditional distribution 49

Linear case Bimodal distribution Mean Variance Skewness Kurtosis linear grid model True posterior 3.35 4.23 -0.74 2.00 x = Au + Bw Gaussian approximation 3.20 3.57 0 3 Affine transformation 3.20 3.57 -0.54 2.35 polytopic constraints 0 . 4 Cx ≤ z 10 linear measurement 5 0 . 2 y = Hw 0 0 Approximate conditioning − 5 − 10 − 5 0 5 10 − 5 0 5 10 affine transformation: Annular distribution Mean Variance Skewness Kurtosis True posterior -0.6 32.9 0 1.08 w y = w + K ( y − Hw ) ˆ Gaussian approximation -0.6 17.8 0 3 Affine transformation -0.6 17.8 0 1.60 where K = Σ H ⊤ � H Σ H ⊤ � − 1 0 . 2 10 0 . 15 5 → projection of uncertainty 0 . 1 0 in the subspace { y = Hw } 0 . 05 − 5 → Gaussian case: recovers 0 − 10 the conditional distribution − 10 − 5 0 5 10 − 10 0 10 49