Infrastructure systems are the foundation of our modern society. The Internet, power grids, and transportation networks are just some examples of the several critical systems that our current lifestyle relies on. Due to their large scale, the operational and fault-associated costs that these systems incur are both in the range of hundreds of millions of dollars to several billion dollars [1]. Therefore, operators continuously face the conflicting tasks of operating these systems as efficiently as possible and guaranteeing certain levels of security or robustness are maintained.
Traditionally, this balancing between efficiency and security is achieved by separating tasks across different timescales. Efficiency goals are met using optimization algorithms running at a slow time-scale, and stability/robustness goals are attained using fast time-scale controllers. In power systems, for example, generators are optimally scheduled by solving an economic dispatch optimization problem at a slow time-scale (every 5/15 minutes, hour, or day) [2], whereas fast time-scale controllers based on frequency measurements *This work was supported by the Army Research Office contract W911NF-17-1-0092, NSF grants (CNS 1544771, EPCN 1711188), and Johns Hopkins WSE startup funds.
Zachary E. Nelson and Enrique Mallada are with the Department of Electrical and Computer Engineering, The Johns Hopkins University, 3400 N. Charles Street, Baltimore, MD 21218, emails: {znel-son2,mallada}@jhu.edu [3] focus on preserving the system stability [4], not efficiency.
Unfortunately, the state of flux that these infrastructure systems currently experience due to the growing population, deployment of sensing and communication technologies, and sustainability trends, is pushing system operation towards its limits, and therefore rendering this traditional approach obsolete. Operating at maximum capacity does not leave room for the inefficiencies incurred by the timescale separation. Moreover, the limited coordination capabilities that today's controllers provide -when compared with optimization algorithms-does not allow the system to react to unprescribed events quickly. Motivated by this problem, this paper aims to remove the time-scale separation by building controllers that can simultaneously achieve steadystate optimality while preserving the system stability.
More precisely, this paper proposes a systematic design framework for feedback controllers that, given an LTI system and an unconstrained optimization problem, generates a family of nonlinear controllers that seek to drive the system towards the optimal solution of the optimization problem. Our equilibrium analysis gives a criterion for steady-state optimality of the closed-loop equilibrium. Furthermore, we leverage recent studies of optimization algorithms using Integral Quadratic Constraints (IQCs) [5], [6], [7] to provide sufficient conditions, based on Linear Matrix Inequalities [8], that guarantee global exponential asymptotic stability. The derived LMIs provide an explicit bound on the rate of convergence and allow us to design an algorithm that computes the maximum rate of convergence.
Our controllers have two main distinctive features. Firstly, they can be functionally separated into three components/modules: (i) an Estimator, that aims to estimate the system's state from the output; (ii) an Optimizer, that uses the estimated state to compute the drift direction necessary to achieve optimality or outputs zero when the estimated state is optimal; and (iii) a Driver (PI controller) that generates the necessary input to drive the system toward the optimal solution. Secondly, the Optimizer module can be implemented using one of many optimization algorithms, leading to a family of optimization-based nonlinear controllers. The only required conditions are that (i) in steady-state the output of the optimizer is zero if and only if its input (the estimated state) is the optimal solution of the optimization problem, and (ii) there exists an IQC that captures the input-output relationship of the optimizer. Related Work: Optimization-based control design for achieving optimal steady-state performance has a been a popular subject of research for more than three decades. It has been used in communication networks to reverse engineer TCP/IP congestion control protocols [9], [10] and provide a design framework for novel congestion control algorithms [11], distributed multi-path routing [12], [13], and admission control [14], and access control in wireless networks [15]. In the context of power systems, optimization-based control design has been used for the design of distributed controllers that can achieve efficient supply-demand balance [16], frequency restoration [17], congestion management [18], and economic steady-state optimality [19], [20], [21]. Some of these approaches have been further extended to more general settings such as [22] and [23]. In general, these solutions either require that the dynamical system to be optimized has a specific structure, such as being passive [18], [19], having primal-dual dynamics [17], [21], [23], or having direct access to a subset of the system state [22].
More recently, real-time optimization algorithms have been proposed as a method to mitigate the substantial fluctuations that renewable energy introduces in power networks. The solutions fall within two categories depending on whether the system dynamics are considered as perturbations of the optimization algorithms [24], [25], or the system is modeled as a set of nonlinear algebraic constraints with slowly time-varying parameters [26], [27]. Our work distinguishes from these works by explicitly modeling the system dynamics and simultaneously guaranteeing the stability of the dynamical system and convergence to the optimal solution. Notably, while our framework today does not include optimization constraints or nonlinearities in the system dynamics, extending our framework to incorporate these features is the subject of our current research. Paper Organization: The organization of the paper is as follows. Section II gives the reader preliminary tools that are necessary for the later analysis. Section III sets up the problem and discusses some of the challenges. Section IV proposes a design framework of controllers that address the challenges. Section V shows the systematic procedure for analyzing steady-state optimality and stability. Section VI considers a numerical example to illustrate the practicality of this approach. Lastly, Section VII summarizes the major points of the paper and suggests future work.
The following notation will be used throughout the remainder of the paper. The n × n identity matrix is denoted as I n . The m × n zero matrix is denoted as 0 m×n . The zero vector with length n is denoted as 0 n . The subscripts are removed when the dimensions are implied by context. A positive (semi) definite matrix P ∈ R n×n is denoted as P 0 ( 0). All norms || • || : R n → R are the standard
Given a nonlinear mapping φ : p → q, with p, q ∈ R n , and an input-output reference (p * , φ(p * )) ∈ R n × R n , we consider the class of pointwise IQCs.
Definition 1: The mapping φ is said to satisfy the pointwise IQC defined by
Next, we discuss two particular choices of the nonlinear map φ that are commonly used in optimization algorithms.
One source of nonlinearity that commonly arises in optimization algorithms is the gradient ∇f (p) of a function f : R n → R. In particular, characterizing the input-output properties of the gradient of a strongly convex function with a Lipschitz continuous gradient is of interest.
Definition 2: The gradient mapping ∇f :
Definition 3: The function f : R n → R n is said to be strongly convex if
Using these two properties, it is possible to show that ∇f satisfies the pointwise IQC (Q f , p * , ∇f (p * )) defined by the matrix
We refer the reader to [5] or [28] for a proof of this statement.
The second type of nonlinearity that will be used in this paper arises from the proximal mapping of a function. Definition 4: The proximal mapping Π ρf : R n → R n of the function f : R n → R with real parameter ρ > 0 is defined as
The optimality condition of ( 1) is:
From ( 2), the proximal mapping can be viewed as the composition of the gradient mapping with an affine operator, followed by an inversion operation:
A known result is that Π ρf satisfies the pointwise IQC (Q Π ρf , p * , Π ρf (p * )) defined by the matrix
This result can be derived by using Lemma 1 followed by an IQC for inversion operations [6]. The proximal mapping is also know as the Moreau envelope. We refer the reader to [29] for additional background.
The following lemma, whose proof can be found in [6], shows how to derive IQCs when a nonlinearity φ is composed with an affine map.
Lemma 1: (IQC for Affine Operations) Consider the nonlinear mapping φ that satisfies the pointwise IQC defined by (Q φ , p * , φ(p * )). Define the affine mapping ψ : R n → R n to be ψ(p) := S 2 p + S 1 φ(S 0 p) where S 0 , S 1 , S 2 ∈ R n×n and S 1 is invertible. Then, ψ satisfies the pointwise IQC defined by (Q ψ , p * , ψ(p * )), where
The following lemma is useful when deriving stability conditions in terms of a LMI. See [8] and [30] for details.
Lemma 2:
We will now show how the input-output properties of an IQC can be used to generate a sufficient stability condition for the feedback interconnection of a LTI system and nonlinearity φ.
Proposition 1: Consider a LTI system defined by the matrices
Suppose the LTI system has the nonlinearity φ : R m → R m as feedback so that q = φ( Ĉξ + Dq). Assume φ satisfies the pointwise IQC (Q φ , p * , φ(p * )) and the feedback interconnection is well-posed. 1 Then, the closed-loop equilibrium point ξ * ∈ R n has global exponential asymptotic stability of at least rate α if the LMI ÂT P + P Â + αP P
is feasible for some σ ≥ 0, α > 0, and P 0.
Proof: Assume that (3) is feasible. Let δξ := ξ -ξ * and δq := q -q * , where q * is the input that achieves equilibrium. Consider the quadratic function V (δξ) = (δξ) T P δξ, where P ∈ R n×n , P 0. Lyapunov theory states that if V satisfies:
then the equilibrium point has global exponential asymptotic stability of at least rate α [31].
Clearly, V (0) = 0 and V (δξ) > 0, ∀δξ = 0 hold because P 0. The radial unboundedness property similarly follows from P 0. Using the fact that Âξ * + Bq * = 0, the third 1 The definition of well-posedness is given in Section III.
property can be expressed as
Since ( 3) is feasible and ( 5) holds, it directly follows from Lemma 2 that (4) holds. Hence, the equilibrium ξ * has global exponential asymptotic stability of at least rate α.
The problem setup is illustrated in Fig. 1, where we consider a LTI system represented by a state-space model where x ∈ R n is the state, u ∈ R m is the input, and y ∈ R p is the output:
The input u(t) is the sum of a control signal r(t) ∈ R m and an unknown constant disturbance
Ideally, a disturbance on the measurement y(t) should also be considered and is a current topic of our research. Finally, the feedback operator Ψ(•) denotes the (possibly nonlinear) feedback control to be designed. Our goal is to, given the measurement y, design a control input r = Ψ(y) that drives the system (6) to a steady-state x * that is an optimal solution of a predefined optimization problem
where f : R n → R is a given cost function. Therefore, given the measurement y(t), the feedback Ψ(•) must produce a control r = Ψ(y) such that x(t) → X * , where X * is the set of optimal solutions to (7), i.e.,
Throughout this paper we make the following assumption.
Assumption 1: The cost function f (x) of the optimization problem ( 7) is continuously differentiable, strongly convex, and has a Lipschitz continuous gradient. This implies that the set X * is a singleton. 2Finally, we provide a concrete model for Ψ(•). As the optimality conditions for optimization problem (7) are in general nonlinear, the feedback controllers to be designed will be necessarily of the same type. Thus, we consider the nonlinear feedback Ψ using the nonlinear dynamics
where η ∈ R d is the state of the feedback dynamics, r ∈ R m is the output of the feedback dynamics, and the mappings
Remark 1 (Well-Posedness): From the feedthrough terms present in ( 6) and ( 8), it is possible a priori that the feedback interconnection is not well-posed. 3 A sufficient condition that prevents this problem is by enforcing that whenever D = 0, the map H depends only on η, i.e., r(t) = H(η(t)). We will further discuss this condition in Section IV.
There are several challenges associated to designing (8) such that in steady-state, x * ∈ X * .
• Lack of direct access to x(t): The system output matrix C is not necessarily invertible. Thus, recovering x(t) from y(t) is not straightforward. • Finding the solution x * ∈ X * : Finding the optimal solution to ( 7) is usually challenging or the cost function may change, giving not enough time to recompute x * . • Driving x(t) to x * ∈ X * : Even if one has access to the optimal solution x * , one then needs to design the correct r(t) that ensures that x(t) converges to it. Interestingly, some of these challenges can be easily handled using tools from control theory, such as recovering x(t) from y(t) or driving x(t) to x * . On the other hand, finding an optimal solution x * is the major goal within optimization theory. Therefore, when the timescale of the control and optimization tasks do not intersect, our problem can be easily solved using standard tools from control and optimization. However, when the timescale separation is no longer present, the problem becomes more challenging as there are no standard tools to address it. This problem is systematically addressed in the next section.
In this section we describe the proposed optimizationbased controllers, that combine tools from control and optimization, and leverage the IQC framework described in the preliminaries. The crux of our solution is a modularized architecture that breaks down the feedback dynamics (8) into three serial components that systematically addresses the challenges described in the previous section and allow a straightforward application of Proposition 1 to certify global exponential asymptotic stability. The proposed architecture is described in Fig. 2. The first component E : y(t) → z(t) is a state estimator that takes the output of the LTI system and produces a state estimate z(t). The second component ϕ : z(t) → e(t), referred to as the optimizer, takes the state estimate and produces a measurement of the optimality error or direction of desired drift e(t), which is required to be zero if and only if the input is in the set X * . The optimizer can be thought of as the part of optimization algorithm that dictates the direction of the next step. The third component D : e(t) → r(t), the driver, takes the optimality error and produces the input to the LTI system that ensures that the equilibrium satisfies e * = 0.
Remark 2: One of the advantages of the proposed architecture is the role independence of each component. This allows for a subsystem to be skipped if the functionality is not required. For example, in cases where y(t) = x(t) or the optimization problem uniquely depends on y(t), then the estimator block can be avoided.
In the remainder of this section, we describe the design requirements of each proposed component/subsystem and give examples on how to implement them.
The estimator component E : y(t) → z(t) is perhaps the simplest to design. Its goal to build an estimate of the state, z(t), from y(t). The design requirement of E is:
• A.1: If the dynamics of (6) are in equilibrium, then the dynamics of E are in equilibrium if and only if z * = x * . Therefore, an obvious choice for E is an observer/state estimator. The dynamics of E are given by
where L ∈ R n×p is a matrix to be designed. A standard argument for observers shows that the evolution of the error δx(t) := x(t) -z(t) is given by
Moreover, if ( 6) is observable, L can be chosen to satisfy
The optimizer ϕ has two design requirements:
• B.1: The optimizer must take the estimated state z(t) to produce a measure of optimality error or direction of drift e(t) such that e * = 0 if and only if z * ∈ X * . • B.2: The input-output characteristics of ϕ must be captured by an IQC (Q ϕ , z * , ϕ(z * )). For the purpose of this paper, we consider two possible solutions.
ϕ 1 : Gradient Descent. The first solution considered is the standard gradient descent mapping, i.e.,
It is straightforward to verify that requirement B.1 holds. The following lemma explicitly computes the IQC for ϕ 1 . Lemma 3: Assume the pointwise IQC (Q f , z * , ∇f (z * )) is satisfied. Then, the pointwise IQC (Q ϕ1 , z * , ϕ 1 (z * )) defined by the matrix
The IQC immediately follows from Lemma 1 with φ = ∇f , S 0 = I n , S 1 = -I n , and S 2 = 0 n×n . ϕ 2 : Proximal Tracking. Our second option for the optimizer block is inspired by the proximal mapping (1). It essentially computes the error between the input z(t) and the solution given by the proximal operator Π ρf (z(t)), that is,
The following proposition shows that (11) satisfies the first design requirement. Proposition 2: The mapping ϕ 2 satisfies B.1.
Proof: Assume that e * = 0. Then from (11), we have that
By the definition of X * , z * ∈ X * . Conversely, assume that z * ∈ X * . Since The next lemma computes the IQC that characterizes ϕ 2 . Lemma 4: Assume the pointwise IQC (Q f , z * , ∇f (z * )) is satisfied. Then, the pointwise IQC (Q ϕ2 , z * , ϕ 2 (z * )) defined by the matrix
is satisfied.
Proof: The IQC immediately follows from Lemma 1 with φ = Π ρf , S 0 = I n , S 1 = I n , and S 2 = -I n .
The last component of the proposed solution is in charge of generating the control signal r(t) that drives the system towards the optimal solution of (7). The design requirement for D is: It is simple to show that ėI = 0 if and only if e * = 0. In fact, this also shows that we only need an integrator to satisfy the design requirement. However, a PI controller provides better dynamic properties than a pure integrator and therefore we choose to add the proportional term.
This resulting interconnected system using the example design choices is shown in Fig. 3. Remark 3: Whenever the estimator subsystem E is included in the interconnection, the feedback interconnection will be well-posed because r(t) only depends on η(t). However, well-posedness is not guaranteed when D = 0, K P = 0, and there is no estimator subsystem because r(t) = K P ϕ(Cx(t) + D(w + r(t))) + K I e I (t) depends on itself. One simple solution to overcome this issue is add a module E such that z(t) = y(t) -Du(t).
The next section shows that indeed the integrated system is able to guarantee steady-state optimality under mild conditions and illustrates how the IQC framework can be leveraged to guarantee global exponential asymptotic stability.
The optimality analysis requires the following assumption. Assumption 2: The system is steady-state controllable. That is, given any steady-state x * , there exists an input u * such that Ax * + Bu * = 0.
Assumption 2 is in some sense necessary to ensure that the system can achieve an arbitrary steady-state. Although this assumption is stronger than the standard controllability assumption, we point out that while controllability is sufficient to drive x(t) towards any state x * in finite time, it does not requires that x(t) remains equal to x * .
This section will derive a sufficient condition for the global exponential asymptotic stability of the equilibrium point considered in the optimality analysis. For this analysis, it is useful to group the linear dynamics of E and D into the LTI system to essentially create a larger dimension LTI system. The resulting system is
) is feasible for some σ ≥ 0, α > 0, and P 0.
Proof: The proof is a direct application of Proposition 1 to the larger dimension LTI system.
Finally, we show how the LMI condition derived in Theorem 2 can be leveraged to compute the maximum convergence rate that the system can achieve. Our goal here is to solve the optimization problem: maximize σ≥0,α>0,P 0 α subject to (12).
(
The main challenge is that because α multiplies P in ( 12), the optimization problem is non-convex. However, for a fixed α, finding whether ( 12) is feasible can be done efficiently. Therefore, it is possible to implement a line search in α that finds the maximum value α max that satisfies (12).
The purpose of this section is to jointly compare the performance of the gradient and proximal optimizers on a simple MIMO system. We leave the more extensive numerical study for the journal version of this paper. The example we are interested in is a second order LTI system defined by
with x, u ∈ R 2 and y ∈ R. Since the output of the LTI system only has information about the first state, an estimator module is obviously needed. We choose the controller parameters as
and L = 1 1 .
We first find the solution of ( 13) as a function of ρ for the specified LTI system and controller. Several curves corresponding to different convexity and Lipschitz parameters are plotted in Fig. 4. With a sufficiently large ρ, ϕ 2 was able to achieve a larger α max than ϕ 1 when m = 0.75, but was not able to when m = 1.25. For both optimizer types, L = 1.25 resulted in a larger α max than L = 1.5. After choosing a sufficiently large ρ, the solution of ( 13) was plotted as a function of m and is shown in Fig. 5. The Lipschitz constant was chosen as multiples of m. The ϕ 1 optimizer resulted in a larger α max when L was chosen closer to m. Conversely, the ϕ 2 optimizer resulted in a larger α max when the multiple was chosen farther from m.
Let the cost function of the optimization problem be f (x) = 1 2
x T 1 1/6 1/6 2/3
x -17/3 4/3 x.
The Lipschitz constant is the larger eigenvalue of the quadratic matrix L ≈ 1.0690 and the strong convexity constant is the smaller eigenvalue m ≈ 0.5976. The optimal solution of the problem is x * ≈ 5.5652 0.6087 T .
The system's state trajectory for various optimizer choices is given in Fig. 6. The constant disturbance signal was initially set to zero, but at t = 75s was changed to w = 1 1
T . All trajectories were able to recover from the change in disturbance and there were cases when ϕ 2 outperformed ϕ 1 . It is particularly interesting that x 2 has the ability to reach the optimal solution despite the fact that it is not being measured.
This paper introduces a framework of nonlinear controllers whose purpose is to drive a given LTI system to the optimal solution of some predefined optimization problem. The controllers are composed of an estimator, optimizer, and driver. We give specific design requirements and possible design choices for each of these modules. Our analysis shows that under mild assumptions, we can guarantee optimal steady-state performance. Moreover, we give a sufficient LMI condition so that global exponential asymptotic stability of the optimal steady-state is guaranteed. Lastly, we present numerical illustrations that demonstrate how the design choices relate to the rate of exponential convergence. The main focus of future work includes further generalizing the proposed framework. In particular, we will consider constrained optimization problems, multiple distributed LTI systems, and performance tuning (i.e., how to choose controller parameters).
-norm. The Kronecker product of two matrices is denoted by the symbol ⊗.
Relaxing this assumption is desired and is a subject of future work.
The feedback interconnection of (6) and (8) is well-posed if u(t) and y(t) are uniquely defined for every choice of states x(t) and η(t).