Linear–quadratic optimal steady state controllers for engineering students and practicing engineers

Abstract

This paper is an overview of fundamental linear–quadratic optimal control techniques used for linear dynamic systems. The presentation is suitable for undergraduate and graduate students and practicing engineers. The paper can be used by class instructors as supplemental material for undergraduate and graduate control system courses. The paper shows how to find the solution to a dynamic optimization problem: optimize an integral quadratic performance criterion along trajectories of a linear dynamic system over an infinite time period (steady-state linear–quadratic optimal control problem). The solution is obtained by solving a static optimization problem. All derivations done in the paper require only elementary knowledge of linear algebra and state space linear system analysis. The results are presented also for the observer-driven linear–quadratic steady-state optimal controller, output feedback-based linear–quadratic optimal controller, and the Kalman filter-driven linear–quadratic stochastic optimal controller. Having full understanding of derivations of the linear–quadratic optimal controller, observer-driven linear–quadratic optimal controller, optimal linear–quadratic output feedback controller, and optimal linear–quadratic stochastic controller, students and engineers will feel confident to use these controllers in numerous engineering and scientific applications. Several optimal linear–quadratic control case studies involving models of real physical systems, with the corresponding Simulink block diagrams and MATLAB codes, are included in the paper.

Keywords

Linear–quadratic optimal controller output feedback optimal controller linear observer Kalman filter linear–quadratic stochastic optimal controller

Introduction

Optimal control theory is considered as a graduate level topic and taught as a graduate level engineering course. Optimal control is hardly mentioned to undergraduate engineering students. The main reason is certainly the use of advanced mathematics needed to formulate the optimal control problem and derive corresponding optimal controllers. Either the calculus of variations or dynamic programing techniques can be used to solve the general problem of optimal control theory: minimize a performance criterion along dynamic system trajectories. The dynamic programing technique was solely developed by Bellman in the middle of the 1950s.¹ Bellman formulated his famous principle of optimality. The calculus of variations,² a mathematical discipline that was rapidly developing at the beginning of the past century, was used by Russian mathematicians headed by Pontryagin to derive the well-known minimum principle,³ also known as the Pontryagin minimum principle. Both approaches, dynamic programing and calculus of variations, were essential for the development of optimal control theory. Due to their mathematical complexity, these approaches are out of reach of undergraduate students. They are commonly taught in graduate courses on optimal control theory and its applications. There are several optimal control theory textbooks appropriate to teach optimal control to engineering graduate students.^4–15

A particular class of optimal controllers will be derived in this paper and presented in sufficient detail such that it can be fully understood by undergraduate students. This class of controllers, linear–quadratic (LQ) optimal controllers (linear system and quadratic performance criterion), has the greatest potential for actual applications among all optimal controllers. The problem formulation of the LQ optimal controller will be done over an infinite time interval (steady state), and all steps in the process of deriving the optimal LQ-controller will be analytically justified. To that end, only the knowledge of elementary linear algebra and basic state space results about linear dynamic systems will be needed—the knowledge that most of the engineering students receive during their undergraduate studies. In general, undergraduate students have knowledge and understanding of linear dynamic systems and corresponding linear system control techniques. Formulating an integral quadratic performance criterion that has to be optimized along trajectories of a linear dynamic system is not difficult, and it should be easily understood by undergraduate students as well as by practicing control engineers. In the last part of the paper, we will extend the deterministic formulation of the LQ optimal controller to linear stochastic systems mostly by paralleling results obtained for deterministic controllers, and without a need for rigorous analysis of stochastic processes and linear stochastic systems, but simply by indicating dual implementation strategies and dual results.

The material presented in this paper has been used, in various forms, by the author and her colleagues over several years at several academic institutions: Rutgers University, Villanova University, California State University, Los Angeles, American University of Sharjah, University of Belgrade, and Lafayette College. The original version of the first part of this paper that deals only with the zero- and non-zero set point LQ full-state optimal controllers was presented by the author at the American Control Conference.¹⁶ In the second part of the paper, the results from the first part are extended to the design of observer-driven LQ optimal controllers. Observer-driven controllers are used when not all state variables are available for feedback so that the state variables have to be estimated using observers. Observers are dynamic systems that have to be designed by control engineers. The third part of the paper considers the LQ optimal output feedback controller. In some applications, the use of linear observers can be avoided, and instead the linear system output that is instantly available for feedback can be used. Both the zero- and non-zero set point output feedback controllers are considered. The forth part of the paper considers the optimal LQ stochastic controller (known also as LQG controller, where G stands for Gaussian), where the Kalman filter plays a role of an observer. The presentation is done in the continuous-time domain. Dual results can be presented in the discrete-time domain by paralleling the continuous-time domain presentation of this paper.

LQ optimal control problem formulation

A time-invariant linear continuous-time system in state space form is defined by

\frac{dx (t)}{dt} = Ax (t) + Bx (t), x (0) = x_{0}

(1)

where

x (t)

is the state space vector of dimension

n

, and

u (t)

is the system input vector of dimension

m

. Matrices

A

and

B

are constant and of appropriate dimensions. An integral quadratic performance criterion to be optimized with respect to the system input is associated with equation (1)

J = \frac{1}{2} \int_{0}^{\infty} [x^{T} (t) R_{1} x (t) + u^{T} (t) R_{2} u (t)] d t, R_{1} = R_{1}^{T} \geq 0, R_{2} = R_{2}^{T} > 0

(2)

In the performance criterion, the matrices $R_{1}$ and $R_{2}$ are penalty matrices, respectively, positive semidefinite (all the eigenvalues of $R_{1}$ are in the closed right-half of the complex plane (including the imaginary axis)) and positive definite (all the eigenvalues of $R_{2}$ are strictly in the right-half complex plane (excluding the imaginary axis)). It follows from equation (2) that the optimization is done over an infinite time interval, and the corresponding controller is also known as the steady-state optimal LQ controller. The solution obtained at steady state leads to a constant feedback gain that preserves the time-invariance of the considered closed-loop (feedback) system. If the optimization is required over a finite time period (integral in equation (2) is defined on the interval $[0, t {}_{f}]$ with $t_{f} < \infty$ ), the obtained feedback gain will be time varying and the closed-loop system will become a time-varying linear system. Since time-varying linear systems are much more complex than time-invariant systems, this will greatly complicate analysis and implementation of such optimal solutions. For example, it is well known that stability of linear time-varying systems cannot be determined in terms of time-varying system eigenvalues. It is possible that all eigenvalues of such a system are in the left-half complex plane for all times, but the time-varying linear system is unstable.¹⁷ Moreover, implementation of feedback loops with time-varying gains is much more complex than corresponding implementation with constant gains. A particular complexity issue that one is faced with while solving the finite time optimal LQ control problem, and the corresponding Riccati differential equation, was considered in a recent paper.¹⁸ To avoid dealing with time-varying linear systems, in the case of finite-time LQ optimization, the use of a suboptimal solution obtained via the infinite-time (steady state) LQ optimization was discussed and analyzed in Radisavljevic and Koskie.¹⁹ Optimization problem (equations (1) and (2)) is formulated naturally in the time domain. The control system optimization problem can be also defined in the frequency domain as discussed in Iorga et al.²⁰

The choice of a quadratic performance criterion in equation (2) is pretty reasonable. Namely, the goal is that equation (2) be minimized at every time instant, meaning that the square of the state variables and the “square” of the control input signals should be minimized. Hence, the absolute goal is to regulate these variables to zero. This is very realistic when the considered state space model represents system deviations from the nominal trajectories, in which case regulating state space variables to zero means bringing the system back to the nominal trajectories. Such state space models are usually obtained by performing linearization of nonlinear systems with respect to their nominal trajectories and nominal control inputs. Minimizing the square of the system input, large signals that potentially enter the system will be avoided. In the case when the state space model represents dynamics of a pure linear system, the regulation of the system state space trajectories to zero might not be the desired goal. Fortunately, there is a technique that allows regulation of state space trajectories of a linear system to non-zero values in an optimal manner leading to so-called set-point optimal LQ controller. Such LQ optimal controllers basically have the same problem formulation as equations (1) and (2) with a minor modification in the implementation of the system closed-loop in order to achieve the desired non-zero set point. This modification will be presented in the follow-up of the paper.

The penalty matrices are chosen by control engineers using their experience dealing with particular control problems, and very often they are chosen as diagonal matrices. In the case when they are diagonal matrices, the engineering experience indicates that the more weight we put on a given state space variable or a given input variable, the more important that variable is. For example, for a system of order $n = 2$ , the choice $R_{1} = diag {10, 1}$ indicates that for this particular optimization problem, it is much more important that the state variable $x_{1} (t)$ is as close as possible to zero than the state variable $x_{2} (t)$ . If there is no need to regulate a particular state variable to near zero, the corresponding weight should be equal to zero. For example, the choice $R_{1} = diag {10, 0}$ indicates that the state variable $x_{2} (t)$ can take arbitrary (subject to system dynamics constraints) values in the optimization process. The choice of penalty matrices gives a lot of freedom to control engineers to achieve their goal of optimally controlling a system. By increasing diagonal elements in matrix $R_{2}$ , we will get the optimal solution and the corresponding controller that has smaller input signals. In order to have a unique solution to the optimization problem defined in equations (1) and (2), the matrix $R_{2}$ must be nonsingular. This will become obvious from the expression that gives the optimal solution for the considered optimization problem to be derived in the next section.

The solution to the optimization problem in which the performance criterion (equation (2)) is minimized subject to constraint (equation (1)) will be sought in terms of linear full-state feedback. Namely, it is assumed that that full-state feedback is available so that the feedback controller can be defined as

u (x (t)) = - Fx (t)

(3)

where F is the constant feedback gain to be determined in the optimization procedure. The choice given in equation (3) is very natural for control engineers. It should be noted that using feedback control in linear dynamic systems has several good features, for example, system stabilization, improvement of system output transient and steady state responses, reduction of system sensitivity to external disturbances, and reduction of system sensitivity to system parameter variations. The optimization problem, defined by equations (1) to (3), was solved in 1960 by Kalman in his paper.²¹ The corresponding LQ optimal controller is called in the control literature the Kalman regulator (controller).

In the next section, it will be shown how the solution to the optimization problem defined in equations (1) to (3) can be obtained using the elementary knowledge of undergraduate linear algebra and state space analysis. As a matter of fact, the state feedback controller given in equation (3) is a simple constant feedback controller (assuming that all state variables are available for feedback, full-state feedback). It is much simpler for implementation than dynamic controllers that require feedback signal differentiation or signal integration, which requires processing of feedback signals though another linear dynamic system.

Having full understanding of the problem formulation and derivations of the solution for the LQ optimal controller, undergraduate students and practicing engineers will feel confident to use this optimal controller and its variants in numerous engineering and scientific applications. Moreover, similar derivations will be used for the observer-driven LQ optimal controller, and the output feedback LQ optimal controller to be presented in the second and third parts of this paper.

Solution of the LQ optimal control problem

The solution to the optimization problem defined in equations (1) to (3) will be obtained using a constrained static optimization technique. To that end, four preliminary and simple mathematical results will be derived and explained in detail using basic knowledge of linear algebra and common results of state space analysis. Under state feedback control (equation (3)), equations (1) and (2) become

\frac{dx (t)}{dt} = (A - BF) x (t), x (0) = x_{0}

(4)

J = \frac{1}{2} \int_{0}^{\infty} x {}^{T}{(t)} (R_{1} + F^{T} R_{2} F) x (t) d t

(5)

The closed-loop state space trajectory can be found from equation (4) in terms of the matrix exponential⁸ using the well-known state space result^8,15

x (t) = e^{(A - BF) t} x (0)

(6)

Substituting equation (6) into equation (5) produces

J = \frac{1}{2} x^{T} (0) {\int_{0}^{\infty} e^{{(A - BF)}^{T} t} (R_{1} + F^{T} R_{2} F) e^{(A - BF) t} d t} x (0)

(7)

The integral in equation (7) can be evaluated in terms of solution of the following algebraic Lyapunov equation^17,22

{(A - BF)}^{T} P + P (A - BF) + R_{1} + F^{T} R_{2} F = 0

(8)

using only matrix differentiation and integration. This fact, called Result 1, is the first basic result, out of four mathematical results that will be needed to complete the optimization process.

Result 1. It can be shown using simple algebra that the solution of the algebraic Lyapunov equation

M^{T} P + PM + N = 0, N = N^{T}

Assuming that matrix M is an asymptotically stable is given by

P = \int_{0}^{\infty} e^{M^{T} t} N e^{Mt} d t

This can be proved by observing that

\frac{d}{dt} {e^{M^{T} t} N e^{M t}} = M^{T} e^{M^{T} t} N e^{M t} + e^{M^{T} t} N e^{M t} M

and integrating both sides of this expression from 0 to

\infty

. The left-hand side produces

e^{M^{T} t} N e^{M t} |_{t = 0}^{t = \infty} = e^{M^{T} \infty} N e^{M \infty} - e^{M^{T} 0} N e^{M 0} = - N

The right-hand side yields

M^{T} \int_{0}^{\infty} e^{M^{T} t} N e^{M t} d t + \int_{0}^{\infty} e^{M^{T} t} N e^{M t} d t M = M^{T} P + PM

- N = M^{T} P + PM

which verifies the stated result. Note that in the above derivations, one of the main properties of the state transition matrix (matrix exponential) was used, that is

\frac{d}{dt} {e^{M t}} = M e^{M t} = e^{M t} M

e^{M 0} = I, e^{M \infty} = 0 for asymptotically stable

Note that in the passage from equations (7) to (8), the matrix $A - BF$ plays the role of the matrix M, and it is an asymptotically stable matrix as a closed-loop system matrix. The performance criterion defined in equation (7) can be now evaluated as

J = \frac{1}{2} x^{T} (0) Px (0) = \frac{1}{2} trace {Px (0) x^{T} (0)}

(9)

where the trace operator represents the sum of all diagonals elements of a square matrix. The right-hand side of equation (9) can be easily justified as demonstrated in Result 2.

Result 2. The trace operator property states that $trace (MN) = trace (NM)$ , assuming that matrices $M$ and $N$ are of appropriate dimensions so that both products $MN$ and $NM$ are defined. Using this property, the following equality $trace {x^{T} (0) Px (0)} = trace {Px (0) x^{T} (0)}$ holds. The fact that $x^{T} (0) Px (0)$ is a scalar, implies the next result $x^{T} (0) Px (0) = trace {x^{T} (0) Px (0)}$ . Students can demonstrate this property easily starting with a simple two-by-two example and generalizing the result obtained to a $n \times n$ matrix.

To avoid the dependence of the performance criterion on the system initial conditions, it is common to assume that, in general, the system initial conditions are uniformly distributed on the unit sphere, so that the performance criterion is defined by

J = \frac{1}{2} trace {P}

(10)

With the above derivations, the dynamic optimization problem: minimize equation (2) subject to equations (1) and (3) is converted into the static optimization problem: minimize equation (10) subject to constraint given in equation (8). In this optimization problem, the unknown constant matrix gain $F$ has to be determined.

Derivations of the optimal solution via static optimization

The solution to the defined static optimization problem: minimize equation (10) subject to equation (8) can be obtained using the classic static optimization method to minimize a function subject to an algebraic constraint, which requires introduction of a Lagrange multiplier and the use of partial derivatives to obtain the necessary conditions for minimum. For the defined optimization problem, the Lagrangian is formed as

L = trace {P + [{(A - BF)}^{T} P + P (A - BF) + R_{1} + F {}^{T}R {}_{2}F] L} = L (P, L, F)

(11)

where

L

plays the role of a Lagrange multiplier. The necessary conditions for a minimum are given by

\begin{array}{l} \frac{\partial L (P, L, F)}{\partial F} = 0 \\ \frac{\partial L (P, L, F)}{\partial L} = 0 \\ \frac{\partial L (P, L, F)}{\partial P} = 0 \end{array}

(12)

At this point, an explanation of the matrix operation stated in equation (12) is needed. Note that $L$ is a scalar function of matrix variables and that $F, L, P$ are matrices. The above partial derivatives are known as gradient matrices²³ and defined by a matrix whose $i, j$ elements are given by (for example, in the case of $\partial L / \partial P$ )

{[\frac{\partial L}{\partial P}]}_{ij} = {[\frac{\partial L}{\partial p_{ij}}]}_{ij}, i, j = 1, 2, \dots, n, p_{ij} = {[P]}_{ij}

(13)

The known formulas for the gradient matrix derivatives can be found in Athans²³ and Bernstein.²⁴

\begin{array}{l} \frac{\partial}{\partial X} trace {M ​ X} = M^{T} \\ \frac{\partial}{\partial X} trace {M ​ X^{T}} = M \\ \frac{\partial}{\partial X} trace {M ​ X ​ N ​ X} = M^{T} X^{T} N^{T} + N^{T} X^{T} M^{T} \\ \frac{\partial}{\partial X} trace {M ​ X ​ N ​ X^{T}} = M^{T} X ​ N^{T} + M ​ X ​ N \end{array}

(14.a)

Using the results $trace (MN) = trace (NM)$ and $trace (M) = trace (M^{T})$ , we have the following variants of equation (14.a)

\begin{array}{l} \frac{\partial}{\partial X} trace {X ​ M^{T}} = M \\ \frac{\partial}{\partial X} trace {M ​ X ​ N} = M^{T} N^{T} \\ \frac{\partial}{\partial X} trace {M ​ X^{T} N} = N ​ M \end{array}

(14.b)

These formulas can be verified using a simple $2 \times 2$ example and generalizing it to the $n \times n$ case. The formulas were established and used in dynamic optimization for the first time in Athans.²³ One of them is verified in Result 3.

Result 3. Consider the scalar function of a matrix variable defined in the first formula in equation (14.a), that is

trace {M ​ X} = trace {\begin{matrix} [\begin{array}{l} m_{11} m_{12} \dots m_{1 n} \\ m_{21} m_{22} \dots m_{2 n} \\ ⋮ ⋮ ⋱ ⋮ \\ m_{n 1} m_{n 2} \dots m_{nn} \end{array}] & [\begin{matrix} x_{11} & x_{12} & \dots & x_{1 n} \\ x_{21} & x_{22} & \dots & x_{2 n} \\ ⋮ & ⋮ & ⋱ & ⋮ \\ x_{n 1} & x_{n 2} & \dots & x_{nn} \end{matrix}] \end{matrix}}

= trace {M ​ X} = {[\begin{matrix} \sum_{j = 1}^{n} m_{1 j} x_{j 1} & * & \dots & * \\ * & \sum_{j = 1}^{n} m_{2 j} x_{j 2} & \dots & * \\ ⋮ & ⋮ & ⋱ & ⋮ \\ * & * & \dots & \sum_{j = 1}^{n} m_{n ​ j} x_{jn} \end{matrix}]} = \sum_{i = 1}^{n} \sum_{j = 1}^{n} m_{ij} x_{ji}

where * denotes elements irrelevant for the matrix trace. The first formula in equation (14.a) can be verified as follows

{[\frac{\partial}{\partial X} trace {M ​ X}]}_{ij} = \frac{\partial}{\partial x_{ij}} (\sum_{i = 1}^{n} \sum_{j = 1}^{n} m_{ij} x_{ji}) = m_{ji} = {[M^{T}]}_{ij}

Using the above formulas, it can be easily shown that equations (11) and (12) imply

\frac{\partial L}{\partial F} = 0 = - B^{T} P L - {(LPB)}^{T} + R_{2} F L + R_{2} F L \Rightarrow F_{opt} = R_{2}^{- 1} B^{T} P

(15)

which gives the expression for the optimal feedback gain that minimizes the performance criterion. Formulas (11) and (12) also produce

\frac{\partial L}{\partial L} = 0 = {(A - BF)}^{T} P + P (A - BF) + R_{1} + F^{T} R_{2} F \Rightarrow A^{T} P + PA + R_{1} - P B^{T} R_{2}^{- 1} B P = 0

(16)

Leading to the algebraic Riccati equation obtained after the substitution of the optimal feedback gain defined in equation 15). The third necessary condition for optimality in formulas (12) and (11) give

\frac{\partial L}{\partial P} = 0 = I + {[L {(A - BF)}^{T}]}^{T} + {[L (A - BF)]}^{T} \Rightarrow {(A - BF)}^{T} L + L ​ (A - BF) + I = 0

(17)

The last equation in (17) is another Lyapunov algebraic equation. It can be used to determine the system asymptotic stability via the well-known direct Lyapunov stability method.²²

Comment: Note that if the initial condition is kept in the expression for the performance criterion as given in equation (9), equation 17) becomes

{(A - BF)}^{T} L + L (A - BF) + x_{0} x_{0}^{T} = 0

Since equations (15) and (16) are not affected by the system initial conditions, the expression for the optimal feedback gain is obviously not affected by the choice of the system initial conditions.

All important formulas established for the LQ optimal controller and its optimal feedback gain, optimal feedback system, and optimal performance criterion are now summarized as

\begin{array}{l} u_{opt} (x_{opt} (t)) = - R_{2}^{- 1} B^{T} P ​ x_{opt} (t) = - F_{opt} x_{opt} (t) A^{T} P + PA + R_{1} - PB R_{2}^{- 1} B^{T} P = 0 \\ \frac{d x_{opt} (t)}{dt} = (A - B F_{opt}) ​ x_{opt} (t), x_{opt} (0) = x_{0} \\ J_{opt} = \frac{1}{2} x_{0}^{T} P x_{0} = \frac{1}{2} trace {P x_{0} x_{0}^{T}} \end{array}

(18)

Note that the optimal feedback gain must be such that the closed-loop system is asymptotically stable, otherwise the performance criterion will become infinity. Since the algebraic Riccati equation is nonlinear and in the matrix form, it has many solutions. It was shown in the original work of Kalman²¹ that for a controllable and observable system, the algebraic matrix Riccati equation has a unique solution that stabilizes the closed-loop systems, that is, the matrix $A - B F_{opt}$ is asymptotically stable. Even more, the unique stabilizing solution is positive definite, $P > 0$ . These facts are formally provided as the next result.

Result 4. Under the assumption that the system is controllable and observable, the unique positive definite solution $P > 0$ of the algebraic Riccati equation exists such that the closed-loop system matrix $(A - B F_{opt})$ is asymptotically stable.²¹ The controllability is apparently with respect to the pair $(A, B)$ . The observability is required with respect to the pair $(A, C)$ where $R_{1} = C^{T} C$ ( $C$ is the Cholesky factor that can be obtained via the Cholesky decomposition, which states that any positive semidefinite matrix $M \geq 0$ can be written as $M = N^{T} N$ ). One way of interpreting and understanding the observability condition imposed on the pair $(A, C)$ is to assume that in the quadratic performance criterion (equation (2)), instead of the term $x^{T} (t) R_{1} x (t)$ , the “square” of the system output is minimized, that is, $y^{T} (t) y (t)$ is included in equation (2) with $y (t) = Cx (t)$ , which leads to $y^{T} (t) y (t) = x^{T} (t) C^{T} C x (t)$ . In that case, the requirement $(A, C)$ is observable seems to be very natural. Hence, before solving the algebraic Riccati equation (16), first a test for controllability and observability has to be performed, which can be done using MATLAB to find the controllability and observability matrices via the MATLAB functions ctrb and obsv, respectively, and examining their rank via the MATLAB function rank. After the initial work of Kalman, it was shown by Wonham²⁵ that controllability and observability conditions imposed, respectively, on $(A, B)$ and $(A, C)$ can be relaxed to the stabilizability condition of $(A, B)$ and detectability condition of $(A, C)$ . Stabilizability means controllability of only unstable system state space variables, and detectability means observability of only unstable system state space variables.²⁴ Hence, applicability of the LQ optimal controller requires the following assumption.

Assumption 1: The pair $(A, B)$ is controllable and the pair $(A, C)$ is observable.

There are several methods for testing controllability of linear time invariant systems.¹⁷ For example, controllability can be tested either by examining whether the rank of the controllability matrix¹⁷ is equal to the system order, equal to n, that is

rank {[\begin{matrix} B & AB & A^{2} B & \dots & A^{n - 1} \end{matrix} B]} = n

(19)

or testing n times the rank of n simpler matrices

rank [λ_{i} I - A] = n, \forall λ_{i}, i = 1, 2, \dots, n, \det (λ_{i} I - A) = 0

(20)

and requiring that all of them have ranks equal to n for all eigenvalues of matrix

A

. Similarly, observability can be tested either by examining the rank of the observability matrix¹⁷

rank {[\begin{matrix} C \\ CA \\ C A^{2} \\ ⋮ \\ C A^{n - 1} \end{matrix}]} = n

(21)

or the ranks of the following n matrices¹⁷

rank {[\begin{matrix} C \\ λ_{i} I - A \end{matrix}]} = n, \forall λ_{i}, i = 1, 2, \dots, n, \det (λ_{i} I - A) = 0

(22)

A weaker form of Assumption 1 can be given in terms of stabilizability and detectability conditions.^17,25

Assumption 1a: The pair $(A, B)$ is stabilizable, and the pair $(A, C)$ is detectable.²⁵

A linear system is both stabilizable and detectable if the following rank conditions are satisfied for all unstable eigenvalues of matrix $A$

rank [λ_{i} I - A] = n, \forall Re {λ_{i} (A) > 0}, i = 1, 2, \dots, n, \det (λ_{i} I - A) = 0

(23)

rank {[\begin{matrix} C \\ λ_{i} I - A \end{matrix}]} = n, \forall Re {λ_{i} (A) > 0}, i = 1, 2, \dots, n, \det (λ_{i} I - A) = 0

(24)

Note that when the system matrix $A$ has all eigenvalues strictly in the left half complex plane (asymptotically stable linear system), then the system is both stabilizable and detectable so that another simple assumption for the existence of the positive semi-definite stabilizing solution of the algebraic Riccati equation (16), and hence of the LQ optimal control problem is given as follows.

Assumption 1 b: The system matrix $A$ is asymptotically stable.

In summary, the optimal feedback gain that minimizes equation (2) subject to equation (1) is given by equation (15), where the matrix $P$ satisfies the algebraic Riccati equation (16). The optimal (minimal) value for the performance criterion is given by equation (10) with $P$ obtained from equation (16). The solution of the algebraic Riccati equation and optimal feedback gain can be obtained using the MATLAB function lqr as [Fopt,P]=lqr(A,B,R1,R2), where lqr stands for the LQ regulator. Note that the optimal controller can be easily implemented since it is a constant state feedback proportional controller. This block diagram of this controller implementation is presented in Figure 1.

Note that when the LQ optimal feedback controller is applied to the linear system, its dynamics is described by $d x_{opt} (t) / dt = (A - B F_{opt}) x_{opt} (t)$ . Due to the fact that the optimal feedback matrix is asymptotically stable, all state variables tend to zero, so that the considered LQ optimal controller is also known as the zero-set point LQ optimal controller. It importance is in cases when the system dynamics represents the linearized dynamics of a nonlinear system, in which case $x_{opt} (t) \to 0$ and $u_{opt} (t) \to 0$ means that the system is going back to its nominal trajectory and its nominal control input. Applications of such LQ optimal controllers are several in engineering areas, especially in avionics and robotics. In Subsection “Optimal deterministic LQ non-zero set-point controller”, we will define and present a solution to the non-zero set optimal LQ controller.

Case study: LQ controller design for an F-15 aircraft

The linearized mathematical model of an F-15 aircraft longitudinal dynamics for supersonic flight conditions is given by^26,27

A = [\begin{matrix} - 0.00819 & - 25.70839 & 0 & - 32.17095 \\ - 0.00019 & - 1.27626 & 1 & 0 \\ 0.00069 & 1.027626 & - 2.40523 & 0 \\ 0 & 0 & 1 & 0 \end{matrix}], B = [\begin{matrix} - 6.80939 \\ - 0.14968 \\ - 14.06111 \\ 0 \end{matrix}]

where the components of the fourth dimensional state vector represent the following physical quantities:

x_{1} (t)

: velocity in feet/sec,

x_{2} (t)

: angle of attack in radians,

x_{3} (t)

: pitch rate in radians per second, and

x_{4} (t)

: pitch angle in radians. The eigenvalues of

A

are

- 0.6835, - 3.0036, - 0.0013 \pm j 0.1037

, which indicates that the matrix is

A

asymptotically stable so that Assumption 1 b is satisfied and the corresponding LQ optimal controller can be designed.

The MATLAB program needed for the design of a LQ optimal steady state controller is given by

≫ % MATLAB program for the LQ optimal controller design

≫ % input matrices A, B, C, R1 = C′*C, R2, and vector x0

≫ eig(A); % finds the eigenvalues of A

≫ CM=ctrb(A,B); rank(CM); % controllability test

≫ OM=obsv(A,C); rank(OB); % observability test

≫ [Fopt, P]=lqr(A,B,R1,R2); % solves the algebraic Riccati equation and finds the corresponding LQ optimal feedback gain

This program is demonstrated on the F-15 aircraft example using the quadratic performance criterion penalty matrices as $R_{1} = C^{T} C = I_{4}, R_{2} = 1$ , indicating that all state variables and the control input have the same penalties, and with all initial conditions set to 1. The following results were obtained for the optimal feedback gain and for the optimal performance

F_{opt} = [\begin{matrix} 0.9980 & - 5.5322 & - 1.5484 & - 8.0155 \end{matrix}], J = 9.0046

MATLAB has no function for testing stabilizability and detectability, but that can be achieved via the rank test presented in formulas (23) and (24) using the MATLAB functions eig and rank.

Optimal deterministic LQ non-zero set-point controller

In many applications, the goal is to regulate either all or some state variables to constant values. These applications include, among others, cruise control that regulates velocity of a moving vehicle at a constant value, temperature control, and voltage control in electrical circuits. Denote the $n_{1}$ -dimensional vector ( $n_{1} < n$ ) of controlled state variables $z (t)$ and relate it with the state variables by

z (t) = C_{z} x (t)

(25)

The goal is that $z (t)$ is as close as possible to a certain constant vector, say $z (t) \approx r$ , at least at steady state, that is

z_{ss} = C_{z} x_{ss} = r = const

(26)

From the original linear system equation (1), the steady state values of the state space variables satisfy

0 = A x_{ss} + B u_{ss}

(27)

From the last two equations, the steady state values can be found for the system state variables and the system input by solving linear algebraic equations

{[\begin{matrix} A & B \\ C_{z} & 0 \end{matrix}]}^{(n + n_{1}) \times (n + m)} {[\begin{matrix} x_{ss} \\ u_{ss} \end{matrix}]}^{(n + m) \times 1} = {[\begin{matrix} 0 \\ r \end{matrix}]}^{(n + n_{1}) \times 1}

(28)

This system produces solutions for $x_{ss}$ and $u_{ss}$ if the following holds

rank {[\begin{matrix} A & B \\ C_{z} & 0 \end{matrix}]}^{(n + n_{1}) \times (n + m)} = n + n_{1}

(29)

The unique solutions for $x_{ss}$ and $u_{ss}$ will be obtained if $m = n_{1}$ , which is the case when the dimension of the system input is equal to the dimension of the controlled variables.

Having obtained the values for $x_{ss}$ and $u_{ss}$ , the change of variables can be introduced as

\bar{x} (t) = x (t) - x_{ss}, \bar{u} (t) = u (t) - u_{ss}

(30)

Using equations (1) and (30), a new linear system can be derived, which in the translated coordinates has the form

\frac{d \bar{x} (t)}{dt} = A \bar{x} (t) + B \bar{u} (t) + ︸ \underset{= 0}{A x_{ss} + B u_{ss}}, \bar{x} (t_{0}) = x_{0} - x_{ss}

(31)

For the linear system (equation (31)), the performance criterion (equation (2)) can be used with the penalty matrices chosen as needed. The optimal feedback controller (for the zero set-point controllers) is

{\bar{u}}_{opt} (x (t)) = - F_{opt} {\bar{x}}_{opt} (t)

(32)

Going back to the original coordinates, it follows that

\begin{array}{l} u_{opt} (t) = {\bar{u}}_{opt} (t) + u_{ss} = - F_{opt} {\bar{x}}_{opt} (t) + u_{ss} = - F_{opt} (x_{opt} (t) - x_{ss}) + u_{ss} \\ = - F_{opt} x_{opt} (t) + F x_{ss} + u_{ss} \end{array}

(33)

The corresponding block diagram is presented in Figure 2.

Case study: Non-zero set point LQ controller design for a voltage regulator

The corresponding state space voltage regulator mathematical model can be found in Kokotovic et al.²⁸ The design goal is to regulate in an optimal manner the steady state controlled output voltage to a constant value, say $z_{ss} = C_{z} x_{ss} = 6 V$ . The problem state space matrices are given by²⁸

A = [\begin{matrix} - 0.2 & 0.5 & 0 & 0 & 0 \\ 0 & - 0.5 & 1.6 & 0 & 0 \\ 0 & 0 & - 10 / 7 & 60 / 7 & 0 \\ 0 & 0 & 0 & - 2.5 & 7.5 \\ 0 & 0 & 0 & 0 & - 1 \end{matrix}], B = [\begin{matrix} 0 \\ 0 \\ 0 \\ 0 \\ 3 \end{matrix}], C_{z} = [\begin{matrix} 1 & 0 & 0 & 0 & 0 \end{matrix}]

The performance criterion penalty matrices and system initial conditions are given within the MATLAB code and presented in the follow-up of this section.

% MATLAB evaluation of non-zero steady state values

r = 6; Cz=[1 0 0 0 0]; M=[A B;Cz 0];

sol=inv(M)*[zeros(5,1);r];

xss=sol(1:5,1:1); uss=sol(6:6,1:1)

% Simulink state space block diagram data needed for full state feedback

C=eye(5); D=zeros(5,1); % for the purpose of obtaining full state feedback

R1 = Cz′*Cz; R2 = 1;

CM=ctrb(A,B); rank(CM); OM=obsv(A,C); rank(OM);

[Fopt,P]=lqr(A,B,R1,R2)

ev_cl=eig(A-B*Fopt); x0=[2; 1; 0; 0; 0]; Jopt = 0.5*x0′*P*x0

The corresponding Simulink block diagram is presented in Figure 3. For the optimal value of the performance criterion, we have obtained $J_{opt} = 14.8161$ . The optimal feedback gain is found to be

F = [\begin{matrix} 0.7805 & 0.4013 & 0.2718 & 0.6081 & 1.4419 \end{matrix}]

The steady state values of the system state space variables are obtained from equation (28) as

x_{ss} = [\begin{matrix} 6.0000 & 2.4000 & 0.7500 & 0.1250 & 0.0417 \end{matrix}], u_{ss} = 0.0139

The linear system controlled output is presented in Figure 4. It can be seen from Figure 4 that after six seconds, the output voltage is regulated to the desired constant value of 6 V.

Observer-driven LQ optimal controller

There is an implementational difficulty with the LQ optimal controller considered in the previous sections. This difficulty is in the fact that all state space variables must be available for feedback control, see equation (3), which in the case of higher-order dimensional systems creates a lot of feedback loops and might be very impractical. Sometimes, it is not feasible to feedback some state variables, and the output signal that represents a linear combination of the state variables, $y (t) = Cx (t)$ , has to be used. In practice, the dimension of the output signal is much smaller than the dimension of the state space variables, that is, $\dim {y (t)} = c < n = \dim {x (t)}$ . In such a case, under certain conditions, an observer, a dynamic system driven by the input and output signals of the system, can be designed with the goal to estimate at all times all the system state variables. These days, the basic steps of the observer design can be found in several standard undergraduate control engineering textbooks.^29–36

An observer is an artificial system designed by a control engineer whose purpose is to provide information about the estimates of the state variables at all times. Being a dynamic system, an observer is either built using capacitors and resistors what electrical engineers would do, or built using masses, springs, and frictional elements what mechanical engineers would prefer, or constructed using simply a personal computer that simulates and solves the corresponding control system differential equation what anybody with basic knowledge of differential equations can do. It takes some time for an observer to provide estimates of the state space variables, so that the obtained state estimate, say $\hat{x} (t)$ will only approximate $x (t)$ , that is, $\hat{x} (t) \approx x (t)$ and hopefully in the long run $\hat{x} (t)$ will tend asymptotically to $x (t)$ . Using $\hat{x} (t)$ instead of $x (t)$ for the feedback signal of the LQ controller will generate an error, comparing to the solution obtained using the full state feedback (instantaneous) LQ optimal controller. The corresponding optimal performance loss in the case of an observer-driven optimal LQ controller will be considered at the end of this section, after the main steps used in the observer design are introduced.

Observer design

The theory of observers originated in the work of Luenberger in the middle of the 1960s.^37–39 According to Luenberger, any system driven by the output of the given system can serve as an observer for that system. Consider a linear dynamic system represented in state space form by

\begin{array}{l} \dot{x} (t) = A x (t) + B u (t), x (t_{0}) = x_{0} = unknown \\ y (t) = C x (t) \end{array}

(34)

where

x (t)

is the state space vector of dimension

n

u (t)

is the system input vector of dimension

m

, and

y (t)

is the system measured output of dimension

c

. Matrices

A

B

, and

C

are constant and of appropriate dimensions. It has been assumed that matrix

C

has full rank equal to

c

(so that there are no redundant measurements). The system output variables are available at all times, and that information can be used to construct an artificial dynamic system of the same order as the system under consideration, the system of order n that will estimate the system state space variables at all times. Since matrices

A, B, C

are known, it is rational to postulate that an observer of equation (34) has exactly the same state space form as the original system, that is

\begin{array}{l} \dot{\hat{x}} (t) = A \hat{x} (t) + B u (t), \hat{x} (t_{0}) = {\hat{x}}_{0} \\ \hat{y} (t) = C \hat{x} (t) \end{array}

(35)

Note that equations (34) and (35) differ only in the initial conditions, and that the observer initial condition can be arbitrarily chosen. If we compare the outputs $y (t)$ and $\hat{y} (t)$ , they will be different, $y (t) \neq \hat{y} (t)$ , since in the first case the initial condition of equation (34) is unknown, and the initial condition for the proposed observer (equation (35)) is arbitrarily chosen by the designer and hence, in general, different from the system initial condition, $x (t_{0}) \neq \hat{x} (t_{0})$ . The difference between these two outputs generates an error signal

y (t) - \hat{y} (t) = Cx (t) - C \hat{x} (t) = Ce (t)

(36)

which can be used as the feedback signal to the observer such that the observation (estimation) error

e (t) = x (t) - \hat{x} (t)

is reduced. The observer that takes into the account this information about the observation error is given by^37–39

\begin{array}{l} \dot{\hat{x}} (t) = A \hat{x} (t) + Bu (t) + K (y (t) - \hat{y} (t)) = A \hat{x} (t) + Bu (t) + KCe (t) \\ \hat{y} (t) = C \hat{x} (t) \end{array}

(37)

where the matrix

K

represents the observer gain that has to be chosen such that the observation error tends to zero as time increases. Differential equation (37) defines the Luenberger observer that represents a copy of the original system with the output injection that carries information about the observation error.

From equations (34) and (37), the observation error dynamics is given by

\dot{e} (t) = (A - KC) e (t)

(38)

If the observer gain $K$ is chosen such that the feedback matrix is asymptotically stable, then the observation error $e (t)$ will decay to zero for any initial condition $e (t_{0})$ . This can be achieved if the pair $(A, C)$ is observable. By taking the transpose of the observation error feedback matrix, that is, $A^{T} - C^{T} K^{T}$ , it can be seen that the pair $(A^{T}, C^{T})$ must be controllable in order to place the observer eigenvalues in the left half complex plane and make the observation error dynamics equilibrium point at the origin asymptotically stable. It is well known that controllability of the pair $(A^{T}, C^{T})$ is equal to observability of the pair $(A, C)$ ,¹⁷ which is the main condition needed for the observer design.

Assumption 2: The pair $(A, C)$ is observable.

For a small observation error, the observer should be chosen to be about 10 times faster than the system. This can be achieved by setting the minimal real part of the observer eigenvalues to be 10 times bigger than the maximal real part of the closed-loop system eigenvalues

{| Re {λ_{\min} (A - KC)} |}_{observe r} > 10 {| Re {λ_{\max} (A - BF)} |}_{system}

(39)

Theoretically, an observer can be made arbitrarily fast by pushing its eigenvalues far to the left in the complex plane, but very fast observers generate noise, and this is not desirable. Moreover, placing observer eigenvalues far to the left in the complex plane requires large values for the observer feedback matrix gain $K$ , which has a consequence of generating large signals that might be harmful to the system (and to the observer), so that in practice the ratio of 10 given in formula (39) is replaced by 5 or 6.

Separation principle

It is important to point out that the system-observer structure preserves the closed-loop system eigenvalues that would have been obtained if the linear perfect state feedback control had been used. The system (equation (34)) under the perfect state feedback control, that is, $u (t) = - Fx (t)$ has the closed-loop form

\dot{x} (t) = (A - BF) x (t)

(40)

The eigenvalues of matrix $A - BF$ are the closed-loop system eigenvalues under perfect state feedback. In the case of the system-observer structure, the actual control input applied to both the system and the observer is

u (t) = - F \hat{x} (t) = - Fx (t) + Fe (t)

(41)

From equations (34), (38), and (41), the following augmented system can be formed

[\begin{matrix} \dot{x} (t) \\ \dot{e} (t) \end{matrix}] = [\begin{matrix} A - BF & BF \\ 0 & A - KC \end{matrix}] [\begin{matrix} x (t) \\ e (t) \end{matrix}] = A_{xe} [\begin{matrix} x (t) \\ e (t) \end{matrix}]

(42)

with the obvious definition of matrix

A_{xe}

. Since the state matrix of this augmented system is upper block triangular, its eigenvalues are the union of the eigenvalues of matrices

A - BF

and

A - KC

. On the other hand, the original system-observer configuration defined by equations (34), (37), and (41) has the state space form

[\begin{matrix} \dot{x} (t) \\ \dot{\hat{x}} (t) \end{matrix}] = [\begin{matrix} A & - BF \\ KC & A - KC - BF \end{matrix}] [\begin{matrix} x (t) \\ \hat{x} (t) \end{matrix}] = A_{x \hat{x}} [\begin{matrix} x (t) \\ \hat{x} (t) \end{matrix}]

(43)

The augmented matrix $A_{x \hat{x}}$ does not clearly indicate what the closed-loop eigenvalues are. However, a very simple relationship among $x (t), e (t),$ and $\hat{x} (t)$ can be obtained from the definition formula of the estimation error as

[\begin{matrix} x (t) \\ e (t) \end{matrix}] = [\begin{matrix} I & 0 \\ I & - I \end{matrix}] [\begin{matrix} x (t) \\ \hat{x} (t) \end{matrix}] = T [\begin{matrix} x (t) \\ \hat{x} (t) \end{matrix}]

(44)

with the matrix

T

being nonsingular. The augmented system (equation (42)) represented in

xe

-coordinates and the augmented system (equation (43)) represented in

x \hat{x}

-coordinates are related via the similarity transformation

T

defined (equation (44)). One of the main properties of the similarity transformation is that it preserves the system eigenvalues,¹⁷ hence matrices

A_{xe}

and

A_{x \hat{x}}

have the same eigenvalues, that is,

λ (A - BF)

and

λ (A - KC)

, are also the eigenvalues in the

x \hat{x}

-coordinates of the augmented observer-system configuration. This fact that the system-observer configuration has the closed-loop eigenvalues separated into the original system closed-loop eigenvalues under perfect state feedback and the actual observer closed-loop eigenvalues is known as the separation principle.

Observer implementation

In practice, the observer is implemented as a linear dynamic system driven by the original system input and output signals, that is, $u (t)$ and $y (t)$ . Eliminating $\hat{y} (t)$ from equation (37), the following observer form, convenient for practical implementation, is obtained

\dot{\hat{x}} (t) = (A - KC) \hat{x} (t) + Bu (t) + Ky (t)

(45)

The corresponding block diagram of the system-observer configuration is presented in Figure 5. Details about MATLAB/Simulink block diagram observer implementation can be found in Radisavljevic-Gajic.⁴⁰

In this subsection, some basic facts are indicated about the observer implementation using Simulink. Since the state space form in Simulink allows only for one vector input and one vector output, the observer with one augmented input can be represented as

\dot{\hat{x}} (t) = (A - KC) \hat{x} (t) + [\begin{matrix} B & K \end{matrix}] [\begin{matrix} u (t) \\ y (t) \end{matrix}]

(46)

Using the given dimensions of the state, input, and output variables, respectively, given by $n$ , $m$ , and $c$ , the system matrices for the state space Simulink blocks should be defined as follows.⁴⁰ For the system, the matrices are respectively: A, B, C, D=zeros(c,m). Note that in the original state space form, the input signal is not present in the system measurements so that $D = 0$ . The system initial condition is unknown, but in order to be able to run simulation, it has to be defined as any vector of dimension $n$ . For the observer, the Simulink state space block matrices are respectively set as

A_{obs} = A - K * C, B_{obs} = [B K], C_{obs} = eye (n), D_{obs} = zeros (n, m + c)

(47)

The observer output matrix is set to identity since all state variables are supposed to be seen on the observer output. Since an observer is an artificial system built by a control engineer, all state variables should be available on its output, and eventually used for feedback. The corresponding MATLAB/Simulink observer implementation block diagram is presented in Figure 6. Figure 7 provides MATLAB/Simulink state space blocks input data for both the system and the observer.

Note that the observer initial condition has to be specified. In general, the observer initial condition can be any vector whose dimension is equal to the system dimension. Engineering experience can help to choose observer’s initial condition, but no general guidelines exist. Possible choices for the observer initial condition can be determined as follows. In practice, engineers usually set $\hat{x} (0) = 0$ ,³⁸ so that $e (0) = x (0) - \hat{x} (0) = x (0)$ .

It was recommended in Johnson,⁴¹ that partial information about the system initial condition available in the system measurements at the initial time, that is, $y (0) = Cx (0)$ be used to obtain a value for $\hat{x} (0)$ . In such a case, $\hat{x} (0)$ can be obtained via the least square method as $\hat{x} (0) = {(C^{T} C)}^{- 1} C^{T} y (0)$ . Note $y (0)$ is known measured signal, and that $rank {C^{T} C} = c < n$ so that the generalized inverse has to be used in the given least square formula.

Another technique for the observer design can be found in Tsui⁴² and Ali et al.⁴³ In this paper, we have followed the classic approach of Luenberger, and the observer presented is often known in the literature as the Luenberger observer.

Observer-driven LQ controller optimal performance

In this section, the quadratic performance criterion defined in equation (2) will be evaluated for the observer-driven controller defined in equation (41), that is, $u (t) = - F \hat{x} (t) = - Fx (t) + Fe (t)$ , whose optimal feedback gain $F$ is determined in equations (15) and (16). Using the observer-driven controller (equation (41)), the quadratic performance criterion (equation (2)) becomes

J = \frac{1}{2} \int_{0}^{\infty} [x^{T} (t) R_{1} x (t) + {\hat{x}}^{T} (t) F^{T} R_{2} F \hat{x} (t)] d t

(48)

This can be also expressed in terms of observer estimation error, $\hat{x} (t) = x (t) - e (t)$ , as

J = \frac{1}{2} \int_{0}^{\infty} [x^{T} (t) (R_{1} + F^{T} R_{2} F) x (t) - 2 x^{T} (t) F^{T} R_{2} F e (t) + e^{T} (t) F^{T} R_{2} F e (t)] d t

(49)

or in a more compact form

\begin{array}{l} J = \frac{1}{2} \int_{0}^{\infty} {{[\begin{matrix} x (t) \\ e (t) \end{matrix}]}^{T} [\begin{matrix} R_{1} + F^{T} R_{2} F & - F^{T} R_{2} F \\ - F^{T} R_{2} F & F^{T} R_{2} F \end{matrix}] [\begin{matrix} x (t) \\ e (t) \end{matrix}]} d t \\ = \frac{1}{2} \int_{0}^{\infty} {{[\begin{matrix} x (t) \\ e (t) \end{matrix}]}^{T} R [\begin{matrix} x (t) \\ e (t) \end{matrix}]} d t = \frac{1}{2} {[\begin{matrix} x (0) \\ e (0) \end{matrix}]}^{T} \int_{0}^{\infty} e^{A_{xe}^{T} t} R e^{A_{xe} t} d t [\begin{matrix} x (0) \\ e (0) \end{matrix}] \end{array}

(50)

with the obvious definition of the matrix

R

, and the solution for the augmented vector

[\begin{matrix} x (t) \\ e (t) \end{matrix}]

obtained from equation (42). According to Result 1, this integral can be evaluated along trajectories of the augmented linear system (equation (42)) via the solution of the following algebraic Lyapunov equation

{PA}_{xe} + A_{xe}^{T} P + R = 0

(51)

Partitioning the symmetric matrix $P$ as

P = [\begin{matrix} P_{11} & P_{12} \\ P_{12}^{T} & P_{22} \end{matrix}]

(52)

and performing block matrix multiplications in equation (51) using expressions for

A_{xe}

and

P

, a system of three matrix algebraic equations is obtained as

P_{11} (A - BF) + {(A - BF)}^{T} P_{11} + R_{1} + F^{T} R_{2} F = 0

(53.a)

P_{12} (A - KC) + {(A - BF)}^{T} P_{12} + P_{11} B F - F^{T} R_{2} F = 0

(53.b)

P_{22} (A - KC) + {(A - KC)}^{T} P_{22} + P_{12}^{T} B F + F^{T} B^{T} P_{12} + F^{T} R_{2} F = 0

(53.c)

We see that equation (53.a) is identical to equation (16), hence $P_{11} = P$ . If the solution of equation (53.a) is used for the optimal feedback gain derived in equation (15), that is, $F = R_{2}^{- 1} B^{T} P_{11}$ , equation (53.a) will become the algebraic Riccati equation (already derived in equation (16)), in this case given by

P_{11} A + A^{T} P_{11} + R_{1} - P_{11} B R_{2}^{- 1} B^{T} P_{11} = 0

(54)

with

F = R_{2}^{- 1} B^{T} P_{11}

, equation (53.b) becomes

P_{12} (A - KC) + {(A - BF)}^{T} P_{12} = 0

(55)

Since $A - KC$ and $A - BF$ are asymptotically stable, the homogeneous algebraic equation (55) has the unique solution $P_{12} = 0$ . With this result, $P_{12} = 0$ , and the optimal feedback gain $F = R_{2}^{- 1} B^{T} P$ , equation (53.c) produces

P_{22} (A - KC) + {(A - KC)}^{T} P_{22} + P_{11} B R_{2}^{- 1} B^{T} P_{11} = 0

(56)

The performance criterion under the observer-driven optimal controller is now given by

\begin{array}{l} {\hat{J}}_{opt} = \frac{1}{2} {[\begin{matrix} x (0) \\ e (0) \end{matrix}]}^{T} P [\begin{matrix} x (0) \\ e (0) \end{matrix}] = \frac{1}{2} {[\begin{matrix} x (0) \\ e (0) \end{matrix}]}^{T} [\begin{matrix} P_{11} & 0 \\ 0 & P_{22} \end{matrix}] [\begin{matrix} x (0) \\ e (0) \end{matrix}] \\ = \frac{1}{2} x^{T} (0) P_{11} x (0) + \frac{1}{2} e^{T} (0) P_{22} e (0) \\ = J_{opt} + J^{loss} = J_{opt} + \frac{1}{2} e^{T} (0) P_{22} e (0) \end{array}

(57)

It can be seen from equation (57) that when $e (0) \to 0$ , then ${\hat{J}}_{opt} \to J_{opt}$ .

At the end of this section, it should be emphasized that there are numerous applications of observers and observer-driven controllers in almost all areas of engineering and sciences.^39,40 It is also interesting to mention that a recent paper⁴⁴ gives a broader recognition to observers by showing that almost all full-order dynamic controllers (compensators) are observer-based controllers. To conclude this section, it is important to indicate that instead of the full-order observer of dimension n, an observer of the reduced-order of dimension $n - c$ ,^17,40 could have been used since the measurement formula $y (t) = C x (t)$ under assumption that matrix $C $ has full rank equal to $c$ , provides at all times $c$ independent algebraic equations for $n$ unknowns of the vector $x (t)$ . Hence, a reduced-order observer-based LQ controller can be designed for engineering applications. The LQ optimal performance loss for a reduced-order observer-based controller has been recently derived in Radisavljevic-Gajic and Milanovic.⁴⁵

Observer-driven optimal LQ controller design case study: F-15 aircraft

The mathematical model of an F-15 aircraft is presented in Section “Case study: LQ controller design for an F-15 aircraft”, where it’s state space matrices $A$ and $B$ are specified. The matrix $C$ is given by⁴⁶

C = [\begin{matrix} 1 & 0 & 1 & 0 \\ 0 & 1 & 0 & 1 \end{matrix}]

The system feedback gain is obtained via the LQ optimal control theory with the penalty matrices chosen as $R_{1} = I_{4}$ and $R_{2} = 1$ . The system closed-loop eigenvalues under the optimal feedback gain

F_{opt} = [\begin{matrix} 0.9975 & - 6.4115 & - 2.1094 & - 10.0049 \end{matrix}]

are given by

λ (A - BF) = {\begin{matrix} - 14.3778 & - 6.2164 \pm j 4.1813 & - 0.7094 \end{matrix}}

The optimal full-state LQ controller performance is given by $J = 92.3487$ . Since, in practice, the observer has to be 5–6 times faster than the system, the observer closed-loop eigenvalues are chosen as

λ (A - KC) = {\begin{matrix} - 100 & - 110 & - 120 & - 130 \end{matrix}}

The MATLAB function “place” was used to find the observer feedback gain that places the observer closed-loop eigenvalues in the desired locations. The following result was obtained for the observer feedback gain

K = 10^{3} [\begin{matrix} - 2.5263 & - 8.9395 \\ 3.6578 & 4.4145 \\ 2.7563 & 8.9239 \\ - 3.6490 & - 4.1883 \end{matrix}]

The system initial conditions (that are unknown, but have to be specified in order to run the Simulink block diagram) are chosen as $x (0) = {[\begin{matrix} 2 & - 3 & 1 & - 2 \end{matrix}]}^{T}$ , and the observer initial conditions are chosen as $\hat{x} (0) = {[\begin{matrix} 0 & 0 & 0 & 0 \end{matrix}]}^{T}$ .

The observer MATLAB/Simulink block diagram is presented in Figure 6 with the input data for the observer state space blocks generated using equation (47) and Figure 7. Performing simulation using the block diagram in Figure 6, the observer output error defined by $e_{y} (t) = y (t) - \hat{y} (t)$ is presented in Figure 8. It can be seen from this figure that the observation error tends to zero very rapidly owning to the fact that the observer closed-loop eigenvalues, $λ (A - KC)$ , are placed far to the left in the complex plane.

The corresponding MATLAB code used is given by

% input A, B, C, R1, R2;

[Fopt,P]=lqr(A,B,R1,R2);

ev_cl=eig(A-B*Fopt) % closed-loop eigenvalues under the optimal LQ feedback control

x0=[2; −3; 1; −2]; Jopt = 0.5*x0′*P*x0; % optimal full-state LQ control performance

lambda_obs_desired=[−100 −110 −120 −130] % observer closed loop eigenvalues

KT=place(A′,C′,lambda_obs_desired); K=KT′

Cobs=eye(4); Dobs=zeros(4,3);

% Optimal performance loss due to use of the observer

P22 = lyap((A-K*C)′,P*B*inv(R2)*B′*P)

x0hat=[0;0;0;0]; e0 = x0-xohat; J_loss = 0.5*e0′*P22*e0

The optimal LQ controller performance loss due to the use of the observer, (equation (57)) is given by $J^{loss} = 530.8158$ . This result indicates that the optimal performance loss can be significant when an observer is used to provide estimates of state space variables to be used for the LQ optimal controller.

LQ output feedback optimal controller

Another approach in LQ optimization is that instead of using the system state space variables, the system output variables are used for feedback,^11,47 that is $y (t) = Cx (t)$ . In such a case, the optimal static output feedback controller is sought in the form

u (y (t)) = - Fy (t) = - FCx (t)

(58)

In contrast to the static output feedback controller that uses for feedback the output of the original system, the dynamic output feedback controller uses the output of another dynamic system (driven by the original system input and the original system output),^47,48 to generate the feedback signal. The block diagram for the output feedback controller is presented in Figure 9.

In practice, the dimension of $y (t)$ is usually much smaller than the dimension of the system state space variables. Derivations of the corresponding optimal output feedback controller are presented below following the same kind of derivations as was done in the case of the LQ optimal controller by converting the dynamic optimization problem into a static optimization problem and taking into account that the optimal control signal is constrained to $u (y (t)) = - Fy (t)$ .

Derivations of the optimal output feedback controller

The linear system (equation (1)) under the output feedback control (equation (58)) is given by

\frac{dx (t)}{dt} = (A - BFC) x (t), x (t_{0}) = x_{0}

(59)

Elimination of the control variable from the quadratic performance criterion (equation (2)) leads to

J = \min_{F} {\frac{1}{2} \int_{0}^{\infty} x^{T} (t) [R_{1} x (t) + C^{T} F^{T} R_{2} F C] x (t) d t}

(60)

It is known from Result 1 that the value of the above integral is given by

J_{ss}^{opt} = \frac{1}{2} x^{T} (0) Kx (0)

(61)

where the matrix

K

satisfies the algebraic Lyapunov equation

{(A - BFC)}^{T} K + K (A - BFC) + R_{1} + C^{T} F^{T} R_{2}^{- 1} F C = 0

(62)

The original dynamic optimization problem is now converted into a static optimization problem, in which the problem is to find the matrix $F$ that minimizes (equation (60)) along the static constraint (equation (62)). As done in the previous parts of this paper, a Lagrange (matrix) multiplier is introduced and the corresponding matrix derivatives are taken to provide the necessary conditions for a minimum. The Lagrangian for equations (60) and (62) is given by

Λ = \frac{1}{2} trace {P x_{0} x_{0}^{T}} + \frac{1}{2} trace {[{(A - BFC)}^{T} P + P (A - BFC) + (R_{1} + C^{T} F^{T} R_{2}^{- 1} F C)] L}

(63)

The necessary conditions for minimum are given in equation (12). Using equations (12) and (14), the expression for the optimal gain and for the $K - L$ equations (similar to the methodology previously presented) are derived as

F = R_{2}^{- 1} B^{T} K L C^{T} {(CL C^{T})}^{- 1}

(64)

and

{(A - BFC)}^{T} K + K (A - BFC) + R_{1} + C^{T} F^{T} R_{2} F C = 0

(65)

(A - BFC) L + L {(A - BFC)}^{T} + x_{0} x_{0}^{T} = 0

(66)

An algorithm can be found in Moerder and Calise⁴⁹ for solving equations (64) to (66) in terms of decoupled Lyapunov equations as follows

Step 1:

Choose F^{(0)} such that A - B F^{(0)} C is stable

(67)

Step 2: Solve the algebraic Lyapunov equations

{(A - B F^{(i)} C)}^{T} K^{(i + 1)} + K^{(i + 1)} (A - B F^{(i)} C) + R_{1} + C^{T} F^{(i)}^{T} R_{2} F^{(i)} C = 0

(68)

(A - B F^{(i)} C) L^{(i + 1)} + L^{(i + 1)} {(A - B F^{(i)} C)}^{T} + x_{0} x_{0}^{T} = 0

(69)

The solution of the algebraic Lyapunov equations (68) and (69) can be found by using the MATLAB function lyap. For example, to solve equation (69) in MATLAB, Li=lyap(A-B*F0*C,x0*x0’) can be used.

Step 3: Update the value for the feedback gain

F^{(i + 1)} = R_{2}^{- 1} B^{T} K^{(i + 1)} L^{(i + 1)} C^{T} {(C L^{(i + 1)} C^{T})}^{- 1}

(70)

This algorithm converges to a local minimum under a very nonrestrictive assumption. In the case when convergence does not happen due to the fact that the minimum is overshot, the updated value for $F$ should be calculated using the formula⁴⁷

F_{\mod}^{(i + 1)} = F^{(i)} + α (F^{(i + 1)} - F^{(i)})

(71)

with

α \in (0, 1]

. Since the optimal value for the performance criterion is given by

J_{opt} = \frac{1}{2} trace {K x_{0} x_{0}^{T}}

(72)

the algorithm convergence can be monitored by requiring that at each iteration, the following holds

J^{(i + 1)} = \frac{1}{2} trace {K^{(i + 1)} x_{0} x_{0}^{T}} < J^{(i)} = \frac{1}{2} trace {K^{(i)} x_{0} x_{0}^{T}}

(73)

If condition (73) is not satisfied, either the local minimum is reached or a lower value for the parameter α should be chosen.

Non-zero set-point output feedback controller

Similarly, as was done for non-zero set point full-state feedback, a non-zero set-point output feedback can be constructed as follows. Assume that the system output $y (t)$ has to be as close as possible to a certain constant vector, say $y (t) \approx r$ , and at steady state to satisfy

y_{ss} = C x_{ss} = r

(74)

From the original system, equation (1), the steady state values of the state space and control variables satisfy

0 = A x_{ss} + B u_{ss}

(75)

From the last two equations, the steady state values for the system state variables and the system input can be found by solving linear algebraic equations

{[\begin{matrix} A & B \\ C & 0 \end{matrix}]}^{(n + p) \times (n + m)} {[\begin{matrix} x_{ss} \\ u_{ss} \end{matrix}]}^{(n + m) \times 1} = {[\begin{matrix} 0 \\ r \end{matrix}]}^{(n + p) \times 1}

(76)

This system produces solutions for $x_{ss}$ and $u_{ss}$ if the following holds

rank {[\begin{matrix} A & B \\ C & 0 \end{matrix}]}^{(n + p) \times (n + m)} = n + p

(77)

The unique solutions for $x_{ss}$ and $u_{ss}$ will be obtained if $m = p$ , which is the case when the number of the system inputs is equal to the number of system outputs. Having obtained the values for $x_{ss}$ and $u_{ss}$ , the value for $y_{ss} = C x_{ss}$ is known. Introducing the change of variables as

\bar{x} (t) = x (t) - x_{ss}, \bar{u} (t) = u (t) - u_{ss}, \bar{y} (t) = y (t) - y_{ss}

(78)

It follows

\begin{array}{l} \frac{d \bar{x} (t)}{dt} = A \bar{x} (t) + B \bar{u} (t), \bar{x} (t_{0}) = x_{0} - x_{ss} \\ \bar{y} (t) = C \bar{x} (t) \\ \bar{u} (t) = F \bar{y} (t) = FC \bar{x} (t) \end{array}

(79)

The optimal solution to the output feedback control problem defined by equations (60) and (79) leads to the same set of algebraic equations (64) to (66), which can be solved using the algorithm given in equations (67) to (71). This leads to ${\bar{u}}_{opt} (t) = F_{opt}^{out} {\bar{y}}_{opt} (t)$ . Going back to the original coordinates, the optimal output feedback control signal is obtained as

\begin{array}{l} u_{opt} (t) = {\bar{u}}_{opt} (t) + u_{ss} = - F_{opt}^{out} {\bar{y}}_{opt} (t) + u_{ss} = - F_{opt}^{out} (y_{opt} (t) - y_{ss}) + u_{ss} \\ = - F_{opt}^{out} y_{opt} (t) + F_{opt}^{out} y_{ss} + u_{ss} \end{array}

(80)

The corresponding block diagram is represented in Figure 10.

Extension to the stochastic LQG optimal control problem

Three fundamental results of the LQ optimal controller and observer-driven LQ optimal controller, derived in the previous sections for deterministic linear systems, can be extended to stochastic linear systems. Namely, the formula for the optimal feedback gain (equation (15)), the separation principle formula (equation (30)), and the observer structure (equation (31)) will hold also in the stochastic case taking into account specifics of stochastic processes and linear stochastic systems.^8,15 Linear stochastic systems driven by general stochastic processes may have very complex behavior so that they hardly can be studied analytically. For that reason, the study in this paper is limited to the Gaussian, zero-mean, stationary, and white noise stochastic processes, which is almost always the case with all linear stochastic systems considered in control engineering theory and applications. The corresponding LQ optimal control problem for stochastic systems is known in the literature as the linear–quadratic–Gaussian (LQG) optimal control problem.

Stochastic LQG optimal control problem formulation

A linear dynamic stochastic system is assumed to be driven (disturbed) by a Gaussian zero-mean white noise stochastic process $w (t)$ with the system measurements corrupted by another Gaussian zero-mean white noise stochastic process $v (t)$ . Since, in general, $w (t)$ and $v (t)$ come from different sources, it is assumed that they are uncorrelated stochastic processes,⁵⁰ which in this case means $E {w (t) v^{T} (t)} = 0$ for all times. The assumptions that the white noise stochastic processes are zero-mean and that they are mutually uncorrelated are only the simplifying assumptions, and they are not crucial for the presented method.

The linear stochastic system under consideration, and all needed statistical information about introduced stochastic processes, are defined by

\begin{array}{l} \dot{x} (t) = A ​ x (t) + B ​ u (t) + Gw (t), E {x (0)} = {\bar{x}}_{0}, Var {x (0)} = Q_{0} \\ E {w (t)} = 0, E {w (t) w^{T} (t + τ)} = W δ (τ), W \geq 0 \\ y (t) = C ​ x (t) + v (t), E {v (t)} = 0, E {v (t) v^{T} (t + τ)} = V δ (τ), \\ V > 0, E {v (t) w^{T} (t + τ)} = 0 \end{array}

(81)

where

x (t) \in R^{n}

are state space variables,

u (t) \in R^{m}

are system inputs,

y (t) \in R^{p}

are system outputs, and

w (t) \in R^{q}

and

v (t) \in R^{p}

are mutually uncorrelated, Gaussian, zero-mean, stationary, white noise stochastic processes.

W

and

V

represent constant spectral densities of white noise stochastic processes

w (t)

and

v (t)

. Namely, the magnitude spectra of these white noise stochastic processes are constant and equal to

W

and

V

In the linear stochastic control problem, the optimal performance has to be defined appropriately taking into account stochastic nature of variables involved. Due to the fact that the system and measurements are disturbed (corrupted) by white noise stochastic processes at all times, the state space variables always have some non-zero values. If the quadratic performance criterion is defined as in equation (2) with the integral’s upper limit set to infinity (which corresponds to steady state optimization), the performance criterion will be equal to infinity. Such a performance criterion will produce a finite value only if some kind of averaging is introduced. In addition, due to the stochastic nature of the problem, the expected value operator has to be used. For these reasons, the quadratic performance criterion for linear stochastic systems is defined by

J_{stoch} = \frac{1}{2} E {\lim_{T \to \infty} [\frac{1}{T} \int_{0}^{T} (x^{T} (t) R_{1} x (t) + u^{T} (t) R_{2} u (t)) d t]} R_{1} = R_{1}^{T} \geq 0, R_{2} = R_{2}^{T} > 0

(82)

Clarification of assumptions and justification of the use of white noise

The Gaussian assumption means that the complete system statistics is described by the mean and variance (all higher order moments, third, fourth, are equal to zero).^8,50 Moreover, having a Gaussian input into a linear stochastic system produces the Gaussian system state space variables and the Gaussian system output.⁵⁰ The fact that the white noise processes is stationary (the concept that corresponds to time invariance for deterministic systems) means that the white noise spectral densities W and V are constant (not time functions). Note that for a stationary white noise processes, the spectrum (Fourier transform $F {W δ (τ)} = W$ ) is constant at all frequencies. The assumption that $w (t)$ and $v (t)$ are white noise stochastic processes have the following justification

The white noise assumption facilitates analytical derivations;

White noise either represents (or at least very well approximates) some real physical processes such as white light, electron thermal noise, wind, unevenness of the road, inaccuracies of instruments and devices, …;

White noise can be used as the worst case scenario for other stochastic processes;

White noise approximation is a realistic assumption for wideband stochastic processes (they have constant spectral density in a very broad range of frequencies). When they are approximated by white noise processes, a real physical system (that in general has a narrow magnitude spectrum) cannot distinguish between the wideband noise and the white noise (the system does not “see”, does not detect, very high frequencies).

Kalman filter

The first step in control of the linear stochastic system defined in equation (81) is to estimate its state space variables. This can be done in an optimal manner by minimizing the mean and the variance of the estimation error $e (t) = x (t) - {\hat{x}}_{KF} (t)$ using the famous Kalman filter^51,52 (see also Kwakernaak and Sivan,⁸ Sinha,¹⁵ Gelb,⁵³ Simon,⁵⁴ and Grewal and Andrews⁵⁵). It is interesting to observe that the Kalman filter has the same structure as the deterministic observer defined in equation (45), that is

{\dot{\hat{x}}}_{KF} (t) = (A - K_{KF} C) {\hat{x}}_{KF} (t) + Bu (t) + K_{KF} y (t), E {{\hat{x}}_{KF} (0)} = {\bar{x}}_{0} = E {x (0)}

(83)

with the Kalman filter gain obtained via optimization that minimizes the variance of the estimation error that is

\min_{K_{KF}} {E [(e (t) - \bar{e} (t)) {(e (t) - \bar{e} (t))}^{T}]}, e (t) = x (t) - {\hat{x}}_{KF} (t)

(84)

Note that the choice of the Kalman filter initial condition specified in equation (83), and the fact that white noise processes are assumed to be zero-mean, make the mean value of the estimation error equal to zero at all times, that is $E {e (t)} = \bar{e} (t) = 0, t \geq 0$ . This result can be shown as follows. From equations (81) and (83), the estimation error dynamics is given by

\begin{array}{l} \dot{e} (t) = (A - K_{KF} C) e (t) + Gw (t) - K_{KF} v (t), e (0) = x (0) - {\hat{x}}_{KF} (0) \\ � e (t) = e^{(A - K_{KF} C) t} e (0) + \int_{0}^{t} e^{(A - K_{FC} C) (t - τ)} [Gw (τ) - K_{KF} v (τ)] d τ \end{array}

(85)

Taking the expected value of both sides with the initial conditions specified in equations (81) and (83) produces

\begin{array}{l} E {\dot{e} (t)} = (A - K_{KF} C) E {e (t)} + G E {w (t)} - K_{KF} E {v (t)}, \\ E {e (0)} = E {x (0)} - E {{\hat{x}}_{KF} (0)} = {\bar{x}}_{0} - {\bar{x}}_{0} = 0 \\ E {\dot{e} (t)} = (A - K_{KF} C) E {e (t)}, E {e (0)} = 0 \\ \Rightarrow E {e (t)} = \bar{e} (t) = e^{(A - K_{KF} C) t} \bar{e} (0) = 0, \forall t � 0 \end{array}

(86)

The result $\bar{e} (t) = 0, \forall t$ justifies why the Kalman filter initial condition in equation (83) is set to be equal to the expected value of the system initial condition.

By forming the product $e (t) e^{T} (t)$ , and finding its expected value $E {e (t) e^{T} (t)}$ , the variance for the estimation error at steady state denoted by Q can be found using equation (85) and taking into account the result established in equation (86)

\begin{array}{l} Q = \lim_{t \to \infty} E {[e (t) - \bar{e} (t)] [e (t) - \bar{e} (t)]} = \lim_{t \to \infty} E {e (t) e^{T} (t)} \\ = E {[\int_{0}^{\infty} e^{(A - K_{KF} C) (t - τ)} [Gw (τ) - K_{KF} v (τ)] d τ] \times {[\int_{0}^{\infty} e^{(A - K_{FC} C) (t - σ)} [Gw (σ) - K_{KF} v (σ)] d σ]}^{T}} \\ = \int_{0}^{\infty} e^{(A - K_{KF} C) (t - τ)} [GW G^{T} + K_{KF} V K {}_{KF}] e^{{(A - K_{FC} C)}^{T} (t - τ)} d τ \end{array}

(87)

The assumption that the white noise processes $w (t)$ and $v (t)$ are uncorrelated has been used in simplification of the above integral. Also, one of the integrals is eliminated due to integration of $GW G^{T} δ (τ - σ)$ and $K_{KF} V K_{KF}^{T} δ (τ - σ)$ terms. Using Result 1, derived in the first section of the paper, it follows from equation (87) that the estimation error steady state variance satisfies the following algebraic Lyapunov equation

(A - K_{KF} C) Q + Q (A - K_{KF} C) + GW G^{T} + K_{KF} V K_{KF}^{T} = 0

(88)

This equation is dual to the optimal regulator, equation (8).

To find the optimal Kalman filter gain, we go through the optimization process presented in the first section of this paper. The optimization goal in the Kalman filtering is to minimize the variance of the estimation error subject to constraint (equation (88)). This minimization problem at steady state can be formulated as a static optimization problem, similar to what was done in the optimal LQ controller problem (see dual formula, equation (11)) by forming the corresponding Lagrangian as

L = trace {Q + [(A - K_{KF} C) Q + Q (A - K_{KF} C) + GW G^{T} + K_{KF} V K_{KF}^{T}] L} = L (Q, L_{KF}, K_{KF})

(89)

where

L_{KF}

is a Lagrange multiplier for the Kalman filter optimization problem. The necessary conditions for a minimum, defined in equation (12) for the controller, now in the case of the Kalman filter, become

\begin{matrix} \frac{\partial L (Q, L_{KF}, K_{KF})}{\partial K_{KF}} = 0 \\ \frac{\partial L (Q, L_{KF}, K_{KF})}{\partial L_{KF}} = 0 \\ \frac{\partial L (Q, L_{KF}, K_{KF})}{\partial Q} = 0 \end{matrix}

(90)

Using the formulas from Appendix 1, conditions (90) lead to the following algebraic equations

\frac{\partial L (Q, L_{KF}, K_{KF})}{\partial K_{KF}} = - L_{KF} Q C^{T} - L_{KF} Q C^{T} + 2 L_{KF} K_{KF} V = 0 \Rightarrow K_{KF}^{opt} = Q C^{T} V^{- 1}

(91)

\frac{\partial L (Q, L_{KF}, K_{KF})}{\partial L_{KF}} = (A - K_{KF} C) Q + Q (A - K_{KF} C) + GW G^{T} + K_{KF} V K_{KF}^{T} = 0

(92)

\frac{\partial L (Q, L_{KF}, K_{KF})}{\partial Q} = I + {(A - K_{KF} C)}^{T} L_{KF} + L_{KF} (A - K_{KF} C) = 0

(93)

Note that since $A - K_{KF} C$ is the Kalman filter feedback matrix, it is asymptotically stable so that Lyapunov algebraic equation (93) implies a positive definite (nonsingular) solution for $L_{KF}$ This facilitates the solution for the optimal Kalman filter gain $K_{KF}^{opt}$ as given in equation (91). Substituting the optimal Kalman filter gain $K_{KF}^{opt}$ in equation (92), it produces

AQ + Q A^{T} + GW G^{T} - Q C^{T} V^{- 1} C Q = 0

(94)

which is the filter algebraic Riccati equation. This equation is dual to the controller algebraic Riccati equation (16). The matrix

Q

representing the optimal steady state variance of the estimation error is a positive semidefinite stabilizing solution of equation (94). Such a solution exists under a dual assumption stated in Result 4 for the controller algebraic Riccati equation. Comparing equations (91), (94), (15), and (16), the duality between the optimal Kalman filter and the optimal LQ controller (also known as the Kalman controller or Kalman regulator) can be established as

A \to A^{T}, B \to C^{T}, R_{1} \to GW G^{T}, R_{2} \to V

(95)

For the existence of the positive definite (semi-definite) solution of the controller algebraic Riccati equation, in Result 4, the controllability (stabilizability) of the pair $(A, B)$ and observability (detectability) of the pair $(A, C)$ with $R_{1} = C^{T} C$ were required. For the filter algebraic Riccati equation (94), these conditions translate into controllability (stabilizability of the pair $(A^{T}, C^{T})$ , which is equivalent to observability (detectability) of the pair $(A, C)$ , and observability (detectability) of the pair $(A^{T}, H^{T})$ where $H^{T} H = G^{T} W G$ ( $H$ is Cholesky factor of $G^{T} W G$ ), which is equivalent to controllability (stabilizability) of the pair $(A, H)$ .

Assumption 3: The pair $(A, H)$ is controllable (stabilizable) and the pair $(A, C)$ is observable (detectable).

Optimal LQG controller

Having obtained the optimal estimates of the state space variables ${\hat{x}}_{KF} (t)$ using the Kalman filter, they can be used for feedback to control linear stochastic system (equation (81)). It is important that the separation principle holds in the linear stochastic system case as well,^8,15 so that the optimal feedback controller has form

u_{opt} (t) = - F_{opt} {\hat{x}}_{KF} (t)

(96)

That is, the optimal state estimation and the optimal LQ controller designs can be separated. The optimal LQ controller stage was solved in the first part of this paper, and the corresponding results obtained in equations (15) and (16) can be used here as well, that is

F_{opt} = R_{2}^{- 1} B^{T} P

(97)

Justification of the separation principle and optimal performance evaluation can be done as following. Eliminate $x (t)$ from equation (82) using the expression for the estimation error defined in equation (84), which produces

\begin{array}{l} J_{s t o c h} = E {\lim_{T \to \infty} \frac{1}{T} \int_{0}^{T} ({\hat{x}}_{K F}^{T} (t) R_{1} {\hat{x}}_{K F} (t) + u^{T} (t) R_{2} u (t) + e^{T} (t) R_{1} e (t) + 2 {\hat{x}}_{K F}^{T} (t) R_{1} e (t)) d t} \\ = E {\lim_{T \to \infty} \frac{1}{T} \int_{0}^{T} ({\hat{x}}_{K F}^{T} (t) R_{1} {\hat{x}}_{K F} (t) + u^{T} (t) R_{2} u (t)) d t} \\ + E {\lim_{T \to \infty} \frac{1}{T} \int_{0}^{T} e^{T} (t) R_{1} e (t) d t} = J_{u} + J_{e} \end{array}

(98)

In the last formula, the term $E {{\hat{x}}^{T} (t) R_{1} e (t)} = R_{1} E {{\hat{x}}^{T} (t) e (t)} = 0$ since the optimal estimate $\hat{x} (t)$ is orthogonal to the estimation error $e (t)$ , which is one of the basic results of the optimal mean-square estimation theory known as the orthogonality principle.^8,50 Hence, the stochastic quadratic performance criterion is represented as $J_{stoch} = J_{u} + J_{e}$ , where the first term depend only on the system control and the second term depends only on the estimation error, which in fact represents the separation principle. At steady state, eliminating from equation (98) the control input given by equation (96) and using Result 2, established in the first section of this paper (see formula (9)), the above performance criterion becomes

J_{stoch}^{opt} = tr {(R_{1} + F_{opt}^{T} R_{2} F_{opt}) E [{\hat{x}}_{KF} (t) {\hat{x}}_{KF} (t)]} + tr {R_{1} E [e^{T} (t) e (t)]}

(99)

It can be shown after simple but lengthy calculations⁸ that the optimal performance value (equation (99)) is given in terms of the quantities that have already been defined and determined as

J_{stoch}^{opt} = tr {Q R_{1} + P K_{opt} V K_{opt}^{T}} = tr {PGW G^{T} + Q F_{opt}^{T} R_{2} F_{opt}}

(100)

Case study: LQG control of an F-15 aircraft

The aircraft matrices $A, B, C$ and penalty matrices $R_{1}, R_{2}$ are given in Sections “Case study: LQ controller design for an F-15 aircraft” and “Observer-driven optimal LQ controller design case study: F-15 aircraft”. The matrix $G$ , according to Gajic and Lim⁴⁶ is $G = B$ . The spectral densities of the white noise stochastic processes are given by⁴⁶

W = 0.000315, V = [\begin{matrix} 0.00686 & 0 \\ 0 & 40 \end{matrix}]

The following MATLAB code can be used to find the optimal system feedback gain, the optimal Kalman filter gains, as well as the optimal performance value defined by formula (100).

[Fopt,P]=lqr(A,B,R1,R2); % lqr stands for the LQ regulator

[Kopt,Q]=lqe(A,G,C,W,V); lqe stands for the LQ estimator

J1opt=trace(Q*R1+P*Kopt*V*Kopt')

J2opt=trace(P*G*W*G'+Q*Fopt'*R2*Fopt); % J2opt = J1opt

The obtained optimal performance value for the considered LQG controller is $J_{stoch}^{opt} = 0.5014$ .

Educational effectiveness of this paper

There is no single publication that studies LQ optimal full-state or partial-state feedback controllers, and/or observer and/or Kalman filter based LQ optimal controllers in a unique way presented in this educational paper. Neither the material presented in this paper can be obtained by taking particular sections from different textbooks and/or different journal and conference papers. All derivations done in the paper follow from each other and they represent a unit. The presentation of the material in this paper comes from the author’s many years of teaching experience at various schools, where this material was taught to both undergraduate and graduate students and tested via exams, homework, term papers, projects, and laboratories. In this paper, difficult concepts, such as Kalman filtering and LQG control, and dynamic optimization were fully derived and explained using elementary knowledge of linear algebra and differential equations that all junior engineering students and many students in sciences should possess. Equipped with such knowledge, the students and practicing engineers via self-study can get full knowledge and understanding of several types of LQ optimal controllers. The most essential ones are considered in this paper. Other variants, such as constant disturbance LQ optimal controllers, could have been easily derived from the presented controllers, or fully understood by reading a book or journal section on that type of controllers. The results derived in this paper are in continuous-time. Their dual discrete-time counterparts can be easily derived by paralleling the techniques presented in this paper.

By reading this paper, students and engineers with solid mathematical background, even those who are not familiar with control systems may grasp the essence of the design of LQ optimal controllers. This is particularly important since these days control systems find applications in numerous areas such as energy systems (solar and fuel cells, wind turbines, power systems), wireless communications, computer networks and Internet, automated highways, autonomous aerial vehicles, autonomous robots, biomedical systems, machine learning algorithms (where the enforcement learning technique is in fact the approximate dynamic programing method of optimal control), and artificial intelligence (where the Kalman filter in the discrete-time domain represents one important iterative technique).

Conclusion

It has been shown how elementary knowledge of linear algebra and state space linear system analysis can be used to completely derive and understand steady state solutions and implementation techniques for some standard optimal control problems: LQ controller, observer-driven LQ controller, Kalman filter driven LQ stochastic controller, and output feedback LQ controller (including non-zero set point optimal LQ controllers) so that these optimal control techniques can be taught to undergraduate students and practicing control engineers.

Footnotes

Declaration of conflicting interests

The author declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author received no financial support for the research, authorship, and/or publication of this article.

Appendix 1: Some linear algebra results used in the paper

References

Bellman

Dynamic programming. Princeton, NJ: Princeton University Press, 1957.

Cossu

An introduction to optimal control. Appl Mech Rev 2014; 66: 024801–1-024801-15.

Pontryagin

Boltyanskii

Gamikrelidze

, et al. The mathematical theory of optimal processes. Geneva, Switzerland: Inderscience Publishers, 1962.

Athans

Falb

Optimal control. New York, NY: McGraw-Hill, 1966, republished by Dover Publications in 2006.

Anderson

Moore

Linear optimal control. Englewood Cliffs, NJ: Prentice Hall, 1970.

Anderson

Moore

Optimal control: Linear quadratic methods. Englewood Cliffs, NJ: Prentice Hall, 1990, republished by Dover Publications in 2007.

Kirk

Optimal control theory: An introduction. Englewood Cliffs, NJ: Prentice-Hall, 1970, republished by Dover Publications in 2004.

Kwakernaak

Sivan

Linear optimal control systems. New York, NY: Wiley Interscience, 1972.

Bryson

Applied optimal control. New York, NY: Hemisphere, 1975.

10.

Sage

White

Optimum systems control. Englewood Cliffs, NJ: Prentice Hall, 1977.

11.

Lewis

Optimal control. New York, NY: John Wiley and Sons, 1986.

12.

Stengel

Optimal control and estimation. Mineola, NY: Dover Publications, 1994.

13.

Naidu

Optimal control systems. Boca Raton, FL: CRC Press, 2001.

14.

Hull

Optimal control theory and applications. Berlin, Germany: Springer, 2003.

15.

Sinha

Linear systems: Optimal and robust control. Boca Raton, FL: CRC Press, 2007.

16.

Radisavljevic-Gajic

Teaching linear-quadratic optimal control to undergraduate students. In: Proceedings of the American control conference, June 30-July 2, 2010, pp. 378–383, Baltimore, USA

17.

Chen

C-T.

Linear system theory and design. Oxford, UK: Oxford University Press, 2013.

18.

Prach

Tekinalp

Bernstein

Infinite-horizon linear-quadratic control by forward propagation of the differential Riccati equation. IEEE Control Systems Magazine, April 2015, pp. 78–93.

19.

Radisavljevic

Koskie

Suboptimal strategy for finite-time linear-quadratic finite optimal control problem. IET Control Theory Appl 2012; 6: 1516–1521.

20.

Iorga

Baruh

Ursu

A review of H_∞ robust control of piezoelectric smart structures. Appl Mech Rev 2008; 61: 040802.

21.

Kalman

Contribution to the theory of optimal control. Boletin Sociedad Matematica Mexicana 1960; 5: 102–119.

22.

Gajic

Qureshi

The lyapunov matrix equation in systems stability and control. San Diego, CA: Academic Press, 1995, republished by Dover Publications, 2008.

23.

Athans

The matrix minimum principle. Inf Control 1968; 11: 592–606.

24.

Bernstein

DS.

Matrix mathematics. Princeton, NJ: Princeton University Press, 2009.

25.

Wonham

On the matrix Riccati equation of stochastic control. SIAM J Control 1968; 6: 681–697.

26.

Brumbaugh

An aircraft model for the AIAA controls design challenge. AIAA J Guid Control Dyn 1994; 17: 747–752.

27.

Shoming

Sznaier

Mixed H2/H∞ control of multimodel plants. AIAA J Guid Control Dyn 1995; 18: 525–531.

28.

Kokotovic

Khalil

O’Reilly

Singular perturbation methods in control: Analysis and design. Philadelphia, PA: SIAM Publishers, 1999.

29.

Nise

Control systems engineering. Hoboken, NJ: John Wiley and Sons, 2008.

30.

Stefani

Savant

Shahian

, et al. Design of feedback control systems. Orlando, FL: Saunders College Publishing, 1994.

31.

Franklin

Powel

Emami-Naeini

Feedback control of dynamic systems. Upper Saddle River, NJ: Prentice Hall, 2002.

32.

Dorf

Bishop

Modern control systems. Upper Saddle River, NJ: Pearson Education, 2005.

33.

Ogata

Modern control engineering. Upper Saddle River, NJ: Prentice Hall, 2002.

34.

Shinners

Modern control systems theory and design. New York, NY: John Wiley and Sons, 1992.

35.

Gajic

Lelic

Modern control systems engineering. London, UK: Prentice Hall International, 1996.

36.

Friedland

Control systems design: An introduction to state space methods. Mineola, NY: Dover Publications, 2005.

37.

Luenberger

Observing the state of a linear system. IEEE Trans Mil Electron 1964; 8: 74–80.

38.

Luenberger

Observers for multivariable systems. IEEE Trans Automat Contr 1966; 11: 190–197.

39.

Luenberger

An introduction to observers. IEEE Trans Automat Contr 1971; 16: 596–602.

40.

Radisavljevic-Gajic

Full-order and reduced-order observer implementations in MATLAB/SIMULINK. IEEE Control Systems Magazine: Lecture Notes, October 2015, pp. 91–101.

41.

Johnson

CD.

Optimal initial condition for full-order observers. Int J Control 1988; 48: 857–864.

42.

Tsui

CC.

Observer design: A survey. Int J Autom Comput 2015; 12: 50–61.

43.

Ali

Hoang

Hussain

, et al. Review and classification of recent observers applied in chemical process systems. Comput Chem Eng 2015; 76: 27–41.

44.

Rahman

Bernstein

Are all full-order dynamic compensators observer based? IEEE Control Systems Magazine, February 2017, pp. 116–122.

45.

Radisavljevic-Gajic

Milanovic

Loss of optimal performance of the finite-horizon continuous-time linear-quadratic controller driven by a reduced-order observer. J Dyn Sys Meas Control 2018; 140: 061014.

46.

Gajic

Lim

M-T.

Optimal control of singularly perturbed linear systems and applications. Didcot, UK: Taylor and Francis, 2001.

47.

Levine

Athans

On the determination of the optimal constant output feedback gains for linear multivariable systems. IEEE Trans Automat Contr 1970; 15: 44–48.

48.

How

JH.

Dynamic output feedback control architectures matter.

IEEE Control Systems Magazine, December 2017, pp. 88–92.

49.

Moerder

Calise

Convergence of a numerical algorithm for calculating optimal output feedback gains. IEEE Trans Automat Contr 1985; 30: 900–903.

50.

Papoulis

Pillai

Random variables, probability and stochastic processes. New York, NY: McGraw-Hill, 2012.

51.

Kalman

A new approach to linear filtering and prediction problems. J Basic Eng 1960; 82: 35–45.

52.

Kalman

Bucy

New results in linear filtering and prediction theory. J Basic Eng 1961; 83: 95–108.

53.

Gelb

Applied optimal estimation. Cambridge, MA: MIT Press, 1974.

54.

Simon

Optimal state estimation, kalman, H infinity, and nonlinear approaches. New York, NY: Wiley-Interscience, 2006.

55.

Grewal

Andrews

Kalman filtering: Theory and practice with MATLAB. New York, NY: Wiley-IEEE Press, 2014.