3.5: Review of Model-Based Control Design and Gain Tuning Methods

Learning Objectives

After completing this section, you will be able to:

1. Design model-based controllers utilizing pole assignment, LQR (linear quadratic regulation), and LQG (linear quadratic Gaussian) methods.
2. Design OCC (output covariance constrained) controllers.
3. Describe model-based gain tuning method based on the OCC control design process.

This section reviews various model-based control design methods for mechatronic systems. These include pole assignment [18], linear quadratic regulation (LQR) [5], and linear quadratic Gaussian (LQG) methods [4], as well as output covariance constrained (OCC) control [6] design methods, where the OCC control methodology can be used for gain tuning for model-based control design method. Note that gain tuning is very important for model-based control because no other methodology exists that can be used directly for model-based control gain tuning (see Chapter 6 of reference [19] for details).

For model-based control, we use the following state-space system model

\[G\left(s\right):\ \left\{\ \begin{array}{c} {\dot{x}}_p=A_px_p+B_pu \\ y_p=C_px_p+H_pu \end{array} \right.\label{3.74}\]

where \(x_p\in {\mathbb{R}}^n\) is the system (plant) state vector with dimension \(n\) (also called system order), \(u\in {\mathbb{R}}^{n_u}\) is the system control vector with dimension \(n_u\), \(y_p\in {\mathbb{R}}^{n_y}\) is the system output vector with dimension \(n_y\), and \(A_p\), \(B_p\), \(C_p\), and \(H_p\) are state-space system model matrices with proper dimensions. This state-space model can be also converted into a transfer function matrix as below

\[G\left(s\right)=C_p{\left(sI-A_p\right)}^{-1}B_p+H_p \label{3.75}\]

with \(Y_p\left(s\right)=G\left(s\right)U(s)\), where \(Y_p(s)\) and \(U(s)\) are Laplace transforms of \(y_p(t)\) and \(u(s)\), respectively. A set of system model differential equations can be easily converted into a state-space model; Example 3.15 provides a demonstration.

Example 3.15

Consider Example 2.4 in Chapter 2 with the following differential equation.

\[\ddot{z}\left(t\right)+3\dot{z}\left(t\right)+2z(t)=f(t) \label{3.76}\]

Letting \(x={\left[x_1,\ x_2\right]}^T={\left[z,\ \dot{z}\right]}^T\) and \(u=f\), we have the following.

\[\begin{array}{c} {\dot{x}}_1=x_2 \\ {\dot{x}}_2=-2x_1-3x_2+u \end{array}\label{3.77}\]

Based on these differential equations, the system described by equation \ref{3.76} can be put into the form of equation \ref{3.74} with the following system matrices

\[A_p=\left[ \begin{array}{cc} 0 & 1 \\ -2 & -3 \end{array} \right], B_p=\left[ \begin{array}{c} 0 \\ 1 \end{array} \right], C_p=\left[ \begin{array}{cc} 1 & 0 \end{array} \right], H_p=0\label{3.78}\]

where the transfer function from \(F\left(s\right)=L[f\left(t\right)]\) to \(Z\left(s\right)=L[z\left(t\right)]\) is shown below.

\[G\left(s\right)=\dfrac{Z(s)}{F(s)}=C_p{\left(sI-A_p\right)}^{-1}B_p+H_p=\dfrac{1}{s^2+3s+2} \label{3.79}\]

Before proceeding with model-based control, we need the following definitions for system controllability (stabilizability) and observability (detectability).

Definition 3.4

Consider a linear time-invariant (LTI) dynamic system, as defined in equation \ref{3.74}, or the pair \((A_p,\ B_p)\). The LTI system is said to be controllable if for any initial state \(x_p\left(0\right)=x_0\), given the final state \(x_1=x_p(t_1)\), and for \(t_1>t_0\), there exists a piecewise continuous input \(u(\cdot )\) such that the solution of the LTI system satisfies \(x_p\left(t_1\right)=x_1\). Otherwise, the system is said to be uncontrollable.

An LTI dynamic system or the pair (\(A_p,\ B_p\)) is said to be stabilizable if there exists a state feedback control, \(u=Fx_p\), such that the feedback system is stable, That is, \(A_p+B_pF\) is stable or \(Re\left[eig(A_p+B_pF)\right]<0\), where \(eig(\cdot )\) denotes eigenvalues of a square matrix and \(Re\left[\cdot \right]\) represents the real part of the complex number.

Theorem 3.5 [20]

For the system defined in equation \ref{3.74}, the following statements are equivalent.

1. \((A_p,\ B_p)\) is controllable.
2. The controllability matrix \(\mathcal{C}=\left[B_p,\ A_pB_p,\ A^2_pB_p,\ \cdots ,\ A^{n-1}_pB_p\right]\) has full row rank.
3. The matrix \(\left[A_p-λI,\ B_p\right]\) has full row rank for all \(λ\) (eigenvalues of \(A_p\)) in \(\mathbb{C}\) (complex plane).

In addition, for the system defined in equation \ref{3.74}, the following statements are equivalent.

1. \((A_p,\ B_p)\) is stabilizable.
2. The matrix \(\left[A_p-λI,\ B_p\right]\) has full row rank for all \(Re(λ)\ge 0\) (eigenvalues of \(A_p\)).
3. There exists a matrix \(F\) such that \(A_p+B_p F\) is stable, i.e., \(\operatorname{Re}\left(\lambda\left(A_p+B_p F\right)\right)<0\).

Definition 3.5

An LTI dynamic system denoted by the pair \((C_p,{\ A}_p)\), defined in equation \ref{3.74}, is said to be observable if for any \(t_1>0\), the initial state \(x_p\left(0\right)=x_0\) can be determined from the history of the input \(u(t)\) and the output \(y(t)\) over the interval of \([0,\ t_1]\). Otherwise, the system is said to be unobservable.

An LTI dynamic system denoted by the pair \((C_p,{\ A}_p)\) is said to be detectable if there exists a matrix \(L\) such that \(A_p-LC_p\) is stable. That is, \(Re[eig\left(A_p-LC_p\right)]\ <\ 0\), where \(eig(\cdot )\) denotes eigenvalues of a square matrix.

Theorem 3.6 [20]

For the system defined in equation \ref{3.74}, the following statements are equivalent.

1. \((C_p,{\ A}_p)\) is observable.
2. The observability matrix \(\mathcal{O}={\left[C^T,\ A^T_pC^T_p,\ {\left(A^T_p\right)}^2C^T_p,\ \cdots ,\ {\left(A^T_p\right)}^{n-1}C^T_p\right]}^T\) has full row rank.
3. The matrix \(\left[A^T_p-λI,\ C^T_p\right]\) has full row rank for all \(λ\) (eigenvalues of \({\ A}_p\)) in \(\mathbb{C}\) (complex plane).

In addition, for the system defined in \ref{3.74}, the following statements are equivalent.

1. \((C_p,{\ A}_p)\) is detectable.
2. The matrix \(\left[A^T_p-λI,\ C^T_p\right]\) has full row rank for all \(Re(λ)\ge 0\) (eigenvalues of \({\ A}_p\)).
3. There exists a matrix \(L\) such that \(A_p-L C_p\) is stable, i.e., \(\operatorname{Re}\left(\lambda\left(A_p-L C_p\right)\right)<0\).

Now, we are ready to discuss the pole assignment with state feedback. In Section 3.3, we discussed the relationship between system performance indices (e.g., overshoot and settling time) and the closed-loop system pole locations. If the system is stable, a small \(\left|Re[λ_i]\right|\) produces slow response convergence (\(λ_i\) is the \(i^{th}\) eigenvalue of \({\ A}_p\)). A large \(\left|Re[λ_i]\right|\) leads to fast convergence. A real \(λ_i\) produces a monotonic response, whereas a pair of complex \(λ_i\) and \(λ_{i+1}\) result in an oscillatory response. In addition, for a pair of complex \(λ_i\) and \(λ_{i+1}\), the response becomes more oscillatory as the ratio \(\left|Re\left[λ_i\right]/Im[λ_i]\right|\) decreases. This last can also be expressed as \({{\mathrm{tan}}^{-1} \left(\left|Re\left[λ_i\left({\ A}_p\right)\right]/Im[λ_i({\ A}_p)]\right|\right)\ }\). Therefore, in order to achieve required system performance, there is a motivation to assign the closed-loop poles to desired locations using state feedback control.

Consider the state feedback control \(u=Fx+v\), where \(F\) is \(n_u\times n\) matrix called state feedback gain and \(v(t)\) is an \(n_u\) dimensional vector. The resulting closed-loop system of \ref{3.74} is shown below.

\[G_{CL}\left(s\right):\ \left\{\ \begin{array}{c} {\dot{x}}_p=(A_p+B_pF)x_p+B_pv \ \ y_p=(C_p+H_pF)x_p+H_pv\ \end{array} \right.\label{3.80}\]

The stability and transient response of the closed-loop system are determined by the eigenvalues of \((A_p+B_pF).\ \)Theorem 3.7 states the condition for relocating these open-loop poles to the desired locations. The uncontrollable eigenvalues \(λ_i\) of \((A_p,\ B_p)\), which makes matrix \(\left[A_p-λ_iI,\ B_p\right]\) no longer full row rank, cannot be relocated by state feedback control. That is, if you have any uncontrollable unstable eigenvalue \(λ_i\) of \((A_p,\ B_p)\), the closed-loop system cannot be stabilized by a state feedback control.

Theorem 3.7 [18]

For the dynamic system defined in equation \ref{3.74}, let \((A_p,\ B_p)\) be controllable and \(λ_1\), \(λ_2\), …, \(λ_n\) be a given set of complex numbers with complex numbers being conjugate pairs. Then, there exists a real \(n_u\times n\) matrix \(F\) such that eigenvalues of \((A_p+B_pF)\) are \(λ_1\), \(λ_2\), …, \(λ_n\).

For single-input and single-output (SISO) systems, the system transfer function can be expressed by the following transfer function.

\[G\left(s\right)=\dfrac{b_{n-1}s^{n-1}+\cdots +b_1s+b_0}{s^n+a_{n-1}s^{n-1}+\cdots +a_1s+a_0}+d \label{3.81}\]

The associated state-space realization results in the following system matrices.

\[A_p=\left[ \begin{array}{cc} \begin{array}{ccc} 0 \ & 1 \ & 0 \\ 0 \ & 0 \ & 1 \\ \vdots \ & \vdots \ & \vdots \end{array} & \begin{array}{cc} \cdots \ & 0 \ \\ \cdots \ & 0 \ \\ \ddots \ & 0 \ \end{array} \\ \begin{array}{ccc} 0 & 0 & 0 \\ -a_0 & -a_1 & {-a}_2 \end{array} & \begin{array}{cc} \cdots & 1 \\ \cdots & {-a}_{n-1} \end{array} \end{array} \right], B_p=\left[ \begin{array}{c} \begin{array}{c} 0 \\ 0 \\ \vdots \end{array} \\ \begin{array}{c} 0 \\ 1 \end{array} \end{array} \right], C^T_p=\left[ \begin{array}{c} b_0 \\ b_1 \\ \begin{array}{c} \vdots \\ b_{n-2} \\ b_{n-1} \end{array} \end{array} \right], H_p=d\label{3.82}\]

For a given set of desired poles \(λ_1\), \(λ_2\), …, \(λ_n\), we can define \(q\left(s\right)\) as below.

\[q\left(s\right)=\prod^n_{i=1}{\left(s-λ_i\right)=}s^n+{\overline{a}}_{n-1}s^{n-1}+\cdots +{\overline{a}}_1s+{\overline{a}}_0 \label{3.83}\]

Then, the following state feedback gain \(F\) results in the closed-loop system poles of \(λ_1\), \(λ_2\), …, \(λ_n\).

\[F=\left[ \begin{array}{cc} \begin{array}{ccc} a_0-{\overline{a}}_0 & a_1-{\overline{a}}_1 & a_2-{\overline{a}}_2 \end{array} & \begin{array}{cc} \cdots & a_{n-1}- \end{array} \end{array} {\overline{a}}_{n-1}\right]\label{3.84}\]

Note that

\[A=A_p+B_pF=\left[ \begin{array}{cc} \begin{array}{ccc} 0 \ & 1 \ & 0 \\ 0 \ & 0 \ & 1 \\ \vdots \ & \vdots \ & \vdots \end{array} & \begin{array}{cc} \cdots \ & 0 \ \\ \cdots \ & 0 \ \\ \ddots \ & 0 \ \end{array} \\ \begin{array}{ccc} 0 & 0 & 0 \\ -{\overline{a}}_0 & -{\overline{a}}_1 & {-\overline{a}}_2 \end{array} & \begin{array}{cc} \cdots & 1 \\ \cdots & {-\overline{a}}_{n-1} \end{array} \end{array} \right]\label{3.85}\]

Example 3.16

Consider the system defined in equation \ref{3.74} with system matrices \(A_p\) and \(B_p\) shown below.

\[A_p=\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 1 \end{array} \right], B_p=\left[ \begin{array}{c} 0 \\ 0 \\ 1 \end{array} \right]\label{3.86}\]

The system characteristic polynomial is \({\mathrm{det} \left(sI-A_p\right)\ }=s^3-s^2+0s-1\) with eigenvalues of \(A_p\) as \(s_1=1.4656\) and \(s_{2,3}=-0.2328\pm 0.7926j\). The open-loop system is not stable because \(s_1>0\). We would like to design a state feedback gain \(F\) to assign the eigenvalues of \((A_p+B_pF)\) to \(λ_1=-1\) and \(λ_{2,3}=-1\pm j\) so that the closed-loop system is stable. Based on the derivation above, the desired characteristic polynomial is shown below.

\[q\left(s\right)=(s+1)(s^2+2s+2)=s^3+3s^2+4s+2 \label{3.87}\]

Therefore, the feedback gain assigning eigenvalues of \((A_p+B_pF)\) to \(λ_1=-1\) and \(λ_{2,3}=-1\pm j\) can be calculated as below.

\[F=\left[ \begin{array}{ccc} -1-2 & 0-4 & -1-3 \end{array} \right]=\left[ \begin{array}{ccc} -3 & -4 & -4 \end{array} \right]\label{3.88}\]

For multi-input and multi-output (MIMO) systems, the formulae for calculating the feedback gain \(F\) can be found in [18]; the solution may not be unique. Consider MATLAB [21] command F=place(Ap,Bp,p), where \(Ap=A_p\) and \(Bp=B_p\). It can be used to calculate the state feedback gain, where matrices \(A_p\) and \(B_p\) are system matrices, and \(p={\left[ \begin{array}{cc} \begin{array}{ccc} λ_1, & λ_2 & \cdots \end{array} & λ_n \end{array} \right]}^T\) is a vector of poles to be assigned.

In cases where not all system states are available for feedback control, a state observer can be designed if the pair \((C_p,{\ A}_p)\) is detectable, and the estimated states can be used in place of measured states. The state observer [22] is in the following form

\[\begin{array}{c} \dot{\hat{x}}=A_p\hat{x}+B_pu+L(y_p-\hat{y}) \\ \hat{y}=C_p\hat{x}+H_pu \end{array}\label{3.89}\]

where \(\hat{x}\) is the estimated state vector, \(\hat{y}\) is estimated system output, and \(L\) is the estimation gain matrix. To make sure the estimated state vector \(\hat{x}\) converges to the original state vector \(x\), the matrix \((A_p-LC_p)\) will be stable, i.e., all eigenvalues of \((A_p-LC_p)\) have negative real parts or \(Re\left[A_p-LC_p\right]<0\). Note that the eigenvalues of \((A_p-LC_p)\) can be freely assigned if \((C_p,A_p)\) is observable. The problem of assigning observer eigenvalues of \((A_p-LC_p)\) can be converted to a pole assignment problem. Subtracting the estimated state equation in \ref{3.89} by the state equation in \ref{3.74} leads to the following.

\[\dot{e}={\dot{x}}_p-\dot{\hat{x}}=(A_p-LC_p)\left(x-\hat{x}\right)=(A_p-LC_p)e \label{3.90}\]

Because \((A_p-LC_p)\) is stable, as time goes to infinity, the error vector \(e\) converges to zero. Now, let \(u=F\hat{x}+v\) and \(x_{cl}={\left[ \begin{array}{cc} x^T, & e^T \end{array} \right]}^T\). The resulting closed-loop system of \ref{3.74} with state estimator defined in equation \ref{3.89} can be reorganized into the following form.

\[\dot{x}=\left[ \begin{array}{c} {\dot{x}}_p \\ \dot{e} \end{array} \right]=\left[ \begin{array}{cc} A_p+B_pF & -B_pF \\ 0 & A_p-LC_p \end{array} \right]\left[ \begin{array}{c} x_p \\ e \end{array} \right]+\left[ \begin{array}{c} B \\ 0 \end{array} \right]v=Ax+Bv\label{3.91}\]

The closed-loop system is stable because eigenvalues of the closed-loop system matrix \(A_{cl}\) are those from matrices \(A_p+B_pF\) and \(A_p-LC_p\) that are stable. Therefore, as long as the eigenvalues of \(A_p+B_pF\) and \(A_p-LC_p\) are assigned to be stable, the closed-loop system based on the estimated state vector is asymptotically stable.

During the state estimator design, we normally want the convergence rate of estimated state vector (associated with poles of \(A_p-LC_p\)) to be much faster than the convergence rate of closed-loop system dynamics (associated with poles of \(A_p+B_pF\)). In other words, it is recommended that the minimal absolute value of real parts of the assigned estimator poles shall be at least a few times larger than the maximum absolute value of real parts of the assigned closed-loop system poles. In this case, the state estimation error will have minimal effect on the closed-loop system response, compared with the case using state feedback.

In the LQG (linear quadratic Gaussian) control to be discussed later, the estimation gain \(L\) will be optimized based on intensity matrices of system and measurement noises. For a given level of system noise, the lower the level of measurement noise (in terms of maximum singular value of the measurement noise intensity matrix), the larger the optimized estimation gain \(L\) in term of its maximum singular value, and vice versa. Note that large estimation gain \(L\) is able to move the estimator poles to further left on the complex plane, resulting in fast estimation convergence.

In addition to the pole assignment method based on the system state-space model, linear optimal control techniques are widely used for design optimal controllers. These include linear quadratic regulation (LQR), linear quadratic Gaussian (LQG), output covariance constraint (OCC), and linear parameter-varying (LPV) control. Note that LPV control is not addressed in this book; see references [9], [10], [11], [12], [13] for more details. Let’s start with LQR control.

LQR Control Problem

Consider the following system

\[{\dot{x}}_p=A_px_p+B_pu, x_p\left(0\right)=x_0 \label{3.92}\]

where system state vector \(x_p\), control input vector \(u\), and system matrices \(A_p\) and \(B_p\) are defined earlier. In addition, \(x_0\) is an arbitrary initial system state vector. Note that \(x_0\) is deterministic but can be any given vector with compatible dimension.

Assuming that the pair \((A,\ B)\) is stabilizable, find an optimal controller, \(u\left(t\right)=Fx_p(t)\), that stabilizes the system defined in equation \ref{3.92} and minimizes the following cost function

\[J_{LQR}=\int^{\infty }_0{\left\{x^T_p\left(t\right)Qx_p\left(t\right)+u^T\left(t\right)Ru(t)\right\}dt} \label{3.93}\]

for any given initial system state vector \(x_p\left(0\right)=x_0\) with \(R>0\) (positive definite) and \(Q=C^TC\ge 0\) (positive semidefinite).

From references [5], [6], and [20], the global optimal solution exists if \((A_p,{\ B}_p)\) is stabilizable and \((C,\ A_p)\) is detectable, and the resulting optimal state feedback control is in the following form

\[u\left(t\right)=Fx_p\left(t\right), F=-R^{-1}B^T_pK \label{3.94}\]

where \(K\ge 0\) is the positive semidefinite solution of the following Riccati equation.

\[0=A^T_pK+KA_p-KB_pR^{-1}B^T_pK+Q \label{3.95}\]

Note that the above Riccati equation has a positive semidefinite solution if \((A_p,{\ B}_p)\) is stabilizable and \((C,\ A_p)\) detectable; see [6] and [20] for details.

Remark 3.5

The optimal control law defined in equation \ref{3.94} is global optimal for any given initial state vector. In MATLAB, the optimal LQR control can be designed using the command F=-lqr(Ap,Bp,Q,R,0), where \(Ap=A_p\) and \(Bp=B_p\), by providing required system and weighting matrices. The requirement of \((C,\ A_p)\) detectable is a sufficient condition; it can be reduced down to the condition that \((C,\ A_p)\) has no unobservable mode on the imaginary axis [20]. Although the LQR controller is global optimal, it is optimal with respect to the cost function (equation 3.93) with given weighting matrices \(R>0\) and \(Q\ge 0\). How to select the weighting matrices \(Q\) and \(R\) to achieve desired closed-loop system performance is an open question. The output covariance constraint (OCC) control design generates an optimal \(Q\ge 0\) by converting the problem into a constrained optimization problem [6].

Example 3.17

Consider the system defined in equation \ref{3.92} with system matrices \(A_p\) and \(B_p\) shown below, and with output weighting matrix \(Q=10I_3\) and control weighting \(R=1\), where \(I_3\) is an identity matrix with dimension of three.

\[A_p=\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 1 \end{array} \right], B_p=\left[ \begin{array}{c} 0 \\ 0 \\ 1 \end{array} \right] \label{3.96}\]

An LQR controller can be designed using the MATLAB command F=-lqr(Ap,Bp,Q,R,0), where (\(Ap=A_p\) and \(Bp=B_p\)) with state feedback control gain listed below.

\[F=\left[ \begin{array}{ccc} -4.3166 & -6.3473 & -5.8677 \end{array} \right] \label{3.97}\]

The resulting closed-loop system poles of \((A_p+B_pF)\) are \(λ_1=-3.2138\) and \(λ_{2,3}=-0.827\pm 0.59j\).

Next, we will discuss the LQG control problem, where the following system model is considered with notations used in [6].

\[\left\{ \begin{array}{c} {\dot{x}}_p=A_px_p+B_pu+D_pw_p \\ y_p=C_px_p \ \\ z=M_px_p+v. \ \end{array} \right.\label{3.98}\]

where system state vector \(x_p\), control input vector \(u\), and system matrices \(A_p\) and \(B_p\) are defined in equation \ref{3.74}. In addition, \(w_p\) and \(v\) are system and measurement noises, respectively; \(D_p\), \(C_p\), and \(M_p\) are system noise input, performance output, and measurement output matrices, respectively; and \(z\in {\mathbb{R}}^{n_z}\) is the system measurement vector.

Note that both \(w_p\) and \(v\) are uncorrelated zero-mean white noises with intensity matrices \(W_p>0\) and \(V>0\). That is, letting \(\mathcal{E}\) be an expectation operator, we have the following.

\[\begin{array}{cc} \mathcal{E}\left[w_p(t)\right]=0, & \mathcal{E}\left[w_p(t)w^T_p(t-τ)\right]=W_pδ(τ) \ \mathcal{E}\left[v(t)\right]=0, & \mathcal{E}\left[v(t)v^T(t-τ)\right]=Vδ\left(τ\right) \ \end{array} \label{3.99}\]

For a strictly proper output feedback dynamic stabilizing control law given below,

\[\left\{ \begin{array}{c} {\dot{x}}_c=A_cx_c+B_cz \\ u=C_cx_c \ \end{array} \right.\label{3.100}\]

The resulting closed-loop system is

\[\left\{\ \begin{array}{c} \dot{x}=Ax+Dw \ \\ y=\left[ \begin{array}{c} y_p \\ u \end{array} \right]=\left[ \begin{array}{c} C_y \\ C_u \end{array} \right]x=Cx \end{array} \right.\label{3.101}\]

where \(x={\left[ \begin{array}{cc} x^T_p & x^T_c \end{array} \right]}^T\) and \(w={\left[ \begin{array}{cc} w^T_p & v^T \end{array} \right]}^T\) and system matrices \(A\), \(D\), and \(C\) are defined below.

\begin{equation}
A=
\begin{bmatrix}
A_p & B_pC_c\\
B_cM_p & A_c
\end{bmatrix},
\quad
D=
\begin{bmatrix}
D_p & 0\\
0 & B_c
\end{bmatrix},
\quad
C=
\begin{bmatrix}
C_y\\
C_u
\end{bmatrix}
=
\begin{bmatrix}
C_p & 0\\
0 & C_c
\end{bmatrix}
\label{3.102}
\end{equation}

For the closed-loop system defined in equation \ref{3.101}, let \(W=block\ diag\left[W_p,\ V\right]\) be the closed-loop covariance matrix. Then the closed-loop system state covariance matrix \(X\) can be defined by the positive semidefinite solution of the following Lyapunov equation.

\[0=AX+XA^T+DWD^T \label{3.103}\]

Note that \(X\) is also the controllability Gramian of the pair \(\left(A,\ DW^{0.5}\right).\ \)Let \({\mathcal{E}}_{\infty }[\cdot ]\triangleq {\mathop{\mathrm{lim}}_{t\to \infty } \mathcal{E}[\cdot ]\ }\) and we are ready to define the LQG control problem below.

LQG Control Problem

Assuming that the pair \((A_p,{\ B}_p)\) is stabilizable and \((M_p, A_p)\) is detectable, find an optimal full-order dynamic controller in the form defined in equation \ref{3.100} that stabilizes the closed-loop system defined in equation \ref{3.101} and minimizes the following cost function

\[J_{LQG}={\mathcal{E}}_{\infty }\left[{y^T_p(t)Qy}_p(t-τ)\right]+{\mathcal{E}}_{\infty }\left[u^T\left(t\right)Ru(t-τ)\right]=trace\left(QC_yXC^T_y+RC_uXC^T_u\right) \label{3.104}\]

for zero mean white noise input with covariance matrix \(W=block\ diag\left[W_p,\ V\right]\), where \(X\ge 0\) is defined in equation \ref{3.103}, \(R>0\) (positive definite), and \(C^T_pQC_p\ge 0\) (positive semidefinite).

From references [5], [6], and [20], the global optimal solution exists if \((A_p,{\ B}_p)\) is stabilizable and \((M_p, A_p)\) is detectable, and the optimal full order dynamic controller is in the form of \ref{3.100}, where

\[A_c=A_p+B_pF-LM_p, B_c=L=\tilde{X}M^T_pV^{-1}, C_c=F=-R^{-1}B^TK \label{3.105}\]

where \(\tilde{X}\ge 0\) and \(K\ge 0\) are the positive semidefinite solutions of the following Riccati equations.

\[\begin{array}{c} 0=A_p\tilde{X}+\tilde{X}A^T_p-\tilde{X}M^T_pV^{-1}M_p\tilde{X}+D_pW_pD^T_p \\ 0=A^T_pK+KA_p-KB_pR^{-1}B^T_pK+C^T_pQC_p \end{array}\label{3.106}\]

Note that these Riccati equations have a positive semidefinite solution \(\tilde{X}\ge 0\) if \(({M_p,\ A}_p)\) is detectable and \((A_p,\ D_p)\) is stabilizable. They have a positive semidefinite solution \(K\ge 0\) if \((A_p,\ B_p)\) is stabilizable and \((Q^{0.5}C_p,\ A_p)\) is detectable. The derivation can be found in references [6] and [20] for a generalized case. Note that \(\tilde{X}\ge 0\) is also called the state estimation error covariance matrix.

Remark 3.6

Note that the optimal control cost defined in equation \ref{3.104} is the summation of weighted closed-loop system output variance and weighted control variance. The optimal control law defined in equation \ref{3.105} is global optimal.

Remark 3.7te

In MATLAB, the optimal LQG control gains \(F\) and \(L\) in equation \ref{3.105} can be designed using MATLAB command \(lqr\) by providing required system and weighting matrices. Examples include F=-lqr(Ap,Bp,Qp,R,0)(\(Ap=A_p,\) \(Bp=B_p\), and \(Qp=C^T_pQC_p\)) and Lt=-lqr(Apt,Bpt,Wp,V,0 (Ap=A_p,\) \(Bp=B_p\), \(Wp={D_pW_pD}^T_p\), \(Lt=L^T\)). Similar to the LQR case, the requirement of \((Q^{0.5}C_p,\ A_p)\) detectable is a sufficient condition. It can be reduced down to the condition that \((Q^{0.5}C,\ A_p)\) has no unobservable mode on the imaginary axis (see [20] for details). Similarly, the requirement of \((A_p,\ D_p)\) stabilizable is a sufficient condition; it can be relaxed to the condition that \((A_p,\ D_p)\) has no unstabilizable mode on the imaginary axis.

Example 3.18

Consider the system defined in equation \ref{3.98} with system matrices \(A_p\), \(B_p\), \(D_p\), \(C_p\), and \(M_p\). These are shown below with output weighting matrix \(Q=10I_2\), control weighting \(R=1\), system input noise intensity matrix \(W_p=1\), and measurement noise intensity matrix \(V=0.001I_2\), where \(I_2\) is an identity matrix with dimension of two.

\[A_p=\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 1 \end{array} \right], B_p=D_p=\left[ \begin{array}{c} 0 \\ 0 \\ 1 \end{array} \right], \ C^T_p=M^T_p=\left[ \begin{array}{cc} 1 & 0 \\ 0 & 1 \\ 0 & 0 \end{array} \right]\label{3.107}\]

An LQG controller can be designed using the MATLAB commands F=-lqr(Ap,Bp,Qp,R,0)(\(Ap=A_p,\) \(Bp=B_p\), and \(Qp=C^T_pQC_p\)) and Lt=-lqr(Apt,Bpt,Wp,V,0 (Ap=A_p,\) \(Bp=B_p\), \(Wp={D_pW_pD}^T_pQ\), \(Lt=L^T\)) with dynamic feedback control gain \(F\) and the estimation gain \(L\) listed below.

\[\begin{array}{c} F=-\left[ \begin{array}{ccc} 4.3166 & 5.5794 & 4.4869 \end{array} \right] \\ L^T=\left[ \begin{array}{ccc} 0.9997 & 1.0246 & 1.2449 \\ 1.0246 & 8.9503 & 40.5793 \end{array} \right] \end{array}\label{3.108}\]

The resulting closed-loop system poles are the combination of these from \((A_p+B_pF)\ (λ_1=-1.1441\) and \(λ_{2,3}=-1.1714\pm 1.2356j\)) and these from \((A_p-{LM}_p)\ (λ_4=-1.0015\) and \(λ_{5,6}=-3.9743\pm 3.9745j\)).

One of the key barriers of applying advanced control theories (such as LQR and LQG) to solve practical control problems is that no systematic approach exists for tuning model-based control gains to optimize the system performance. Difference (or modeling error) always exists between the control design model and the physical system. Due to the modeling error, extremely high gain controllers may destabilize the closed-loop physical system [7]. Therefore, it is important to have a method to design a series of model-based controllers with their control gains varying from low to high. In that way, we can start testing with the low gain controller and gradually increase the control gain experimentally until the closed-loop system performance is no long improving. This result would indicate that the best closed-loop system has been achieved; see reference [7] for an application.

To further improving the closed-loop system performance, the control design model needs to be further improved to reduce modeling error. This process of model-based control tuning is motivated to form the following control problem called an output covariance constraint (OCC) control problem (introduced in reference [6] along with a guaranteed convergence control design algorithm). Before presenting the OCC control problem, partition the performance output \(y_p\) from equations \ref{3.98} and \ref{3.101} into \(y_p\cong {\left[y^T_1,\ y^T_2,\ \cdots ,y^T_m\right]}^T\), where \(y_i=C_ix=C_{p,i}x_p\in {\mathbb{R}}^{m_i}\ (C^T=\left[C^T_1,\ C^T_2,\dots ,C^T_m\right]\) and \(C^T_p=\left[C^T_{p,1},\ C^T_{p,2},\dots ,C^T_{p,m}\right]\)) for \(i=1,2,\ \dots ,\ m\) \((\sum^m_{i=1}{m_i}=n_y)\).

OCC Control Problem [6]

Assuming that the pair \((A_p,{\ B}_p)\) is stabilizable and \((M_p, A_p)\) is detectable. Find an optimal full-order dynamic controller in the form defined in equation \ref{3.100} that stabilizes the closed-loop system defined in equation \ref{3.101} and minimizes the following cost function

\[J_{OCC}=\left[u^T\left(t\right)Ru(t-τ)\right] = \text{trace}\left(RC_uXC^T_u\right) \label{3.109}\]

subject to output performance constraints

\[{\mathcal{E}}_{\infty }\left[y_i(t)y^T_i(t-τ)\right]=\ C_iXC^T_i\le {\overline{Y}}_i, i=1,2,\dots ,m \label{3.110}\]

for zero mean white noise input with intensity matrix \(W=block\ diag\left[W_p,\ V\right]\), where \(X\ge 0\) is defined in equation \ref{3.103}, \(R>0\) (positive definite), and \({\overline{Y}}_i>0\ \left(i=1,2,\dots ,m\right)\) (positive definite) are constraint matrices on the output covariance matrices.

Remark 3.8

Reference [6] indicates that the solution of the OCC control problem stated above is an LQG controller with a special choice of output weighting matrix \(Q=block\ diag[Q_1,\ Q_{2,\ \dots ,\ }Q_m]\ge 0\), \(Q_i\in {\mathbb{R}}^{m_i\times m_i} (i=1,\ 2,\ \dots ,\ m\)) such that the following conditions are satisfied.

\[0=\left(Y_i-{\overline{Y}}_i\right)Q_i=\left(C_iXC^T_i-{\overline{Y}}_i\right)Q_i, Y_i\le {\overline{Y}}_i, \ (i=1,\ 2,\ \dots ,\ m) \label{3.111}\]

The proof of this remark can be found in reference [6].

An Iterative OCC Control Design Algorithm [6]

Now, the question is how to find the weighting matrix \(Q\) satisfying conditions stated in equation \ref{3.111} in the form of \(=block\ diag[Q_1,\ Q_{2,\ \dots ,\ }Q_m]\ge 0\), \(Q_i\in {\mathbb{R}}^{m_i\times m_i}\ (i=1,\ 2,\ \dots ,\ m\)). Fortunately reference [6] provides an iterative algorithm to find such a weighting matrix \(Q\) with guaranteed convergence if the optimal OCC solution exists; see reference [6] for details.

1. Assume this system has matrices \(A_p,\ B_p,\ D_p,C_p\) and \(M_p\). It has system and measurement noise intensity matrices \(W_p\) and \(V\) and control weighting matrix \(R>0.\ \) In addition, output covariance bound matrix \(\overline{Y}\ =block\ diag\ \left[{\overline{Y}}_1,{\overline{Y}}_2,\dots ,{\overline{Y}}_m\right]>0,\)initial weighting matrix \(Q\left(0\right)=block\ diag\ [Q_1(0),\ Q_2(0),\ .\ .\ .\ ,\ Q_m(0)]>0\), convergence criterion \(ε>0\), and constants \(α>0,\ 0\ <β<\ 1\). Let \(j\ =\ 0\), and proceed to Step 2.
2. Solve the following Riccati equation for \(\tilde{X}\ge 0\) and let \(L=\tilde{X}M^T_pV^{-1}\) \[0=A_p\tilde{X}+\tilde{X}A^T_p-\tilde{X}M^T_pV^{-1}M_p\tilde{X}+D_pW_pD^T_p \label{3.112}\]
3. Compute \(K(j)\ge 0\) by solving the following Riccati equation and letting \(F\left(j\right)=-R^{-1}B^T_pK(j)\). \[0=A^T_pK(j)+K{\left(j\right)A}_p-K(j)B_pR^{-1}B^T_pK(j)+C^T_pQ(j)C_p \label{3.113}\]
4. Compute \(X_p(j)\) by solving the following Lyapunov equation. \[0=\left(A_p+B_pF(j)\right)X_p(j)+X_p(j){\left(A_p+B_pF(j)\right)}^T+LVL^T \label{3.114}\]
5. Compute the output covariance matrices. \[Y_i\left(j\right)=C_i\left[ \begin{array}{c} I \\ 0 \end{array} \right]\left\{\tilde{X}+X_p\left(j\right)\right\}\left[ \begin{array}{cc} I & 0 \end{array} \right]C^T_i=C_{p,i}\left\{\tilde{X}+X_p\left(j\right)\right\}C^T_{p,i}, \ i=1,\ 2,\ \dots ,\ m \label{3.115}\]
6. If the following conditions are satisfied, \[\sum^m_{i=1}{\left\|\left[Y_i\left(j\right)-{\overline{Y}}_i\right]Q_i(j)\right\|\le ε, } \ Y_i\left(j\right)\le {\overline{Y}}_i, \ i=1,\ 2,\ \dots ,\ m \label{3.116}\] the OCC algorithm is converged, and the resulting dynamics controller matrices in the form of \ref{3.100} are listed in equation \ref{3.117}. Otherwise, proceed to Step 7. \[{ A}_c=A_p+B_pF\left(j\right)-LM_p,\ {\ B}_c=L=\tilde{X}M^T_pV^{-1},\ { C}_c=F(j)=-R^{-1}B^TK(j) \label{3.117}\]
7. Update output weighting matrix \(Q(j)\) as follows. \[Q(j+1)\ =\ βQ(j)+(1-β)\mathcal{P}[Q\left(j\right)+α\{Y(j)\ -\ \overline{Y}\}]\mathrm{\ }) \label{3.118}\]

where \(Y\left(j\right)=block\ diag\ [Y_1(j),\ Y_2(j),\ .\ .\ .\ ,\ Y_m(j)]\) and \(\mathcal{P}[\cdot ]\) is an operator that projects any symmetric matrix to the closest positive semidefinite matrix. Let \(j=j+1\), and proceed to Step 3.

Remark 3.9

The norm used in equation \ref{3.116} can be either 2-norm (maximum singular value of the matrix) or Frobenius norm (square root of summation of all squared entries). Consider the projection operator \(\mathcal{P}[\cdot ]\).For a real symmetric matrix \(M\), real Schur decomposition leads to \(M=\left[U_1,\ U_2\right]\left\{block\ diag[E_p,\ E_n]\right\}{\left[U_1,\ U_2\right]}^T\), where \(E_p\) and \(E_n\) are diagonal matrices containing strictly positive and nonpositive eigenvalues of \(M\) in decreasing order, respectively, and \(\left[U_1,\ U_2\right]\) is a unitary orthogonal matrix. Then, \(\mathcal{P}[\cdot ]\) can be defined below.

\[\mathcal{P}[M]=\left\{ \begin{array}{cc} 0, & M\le 0 \ \\ U_1E_pU^T_1, & otherwise \end{array} \right.\label{3.119}\]

Remark 3.10

Theorems 2.2 and 2.5 in reference [6] state if the OCC algorithm is converged, it provides a global optimal solution to the OCC problem. In addition, if the optimal solution exists, for a small enough parameter \(α>0\), the algorithm convergence is guaranteed.

With the OCC control problem introduced along with its iterative control design algorithm, we are now ready to discuss how to use the OCC control design method for control gain tuning to achieve the best performance for the given control design model. Note that \({\tilde{Y}}_i{=C}_{p,i}\tilde{X}C^T_{p,i}\ (i=1,\ 2,\ \dots ,\ m\)) provides the minimal output covariance matrices that can be achieved.

If the open-loop system is stable, calculate the open-loop system state covariance matrix \(X_o\) by solving the following Lyapunov equation. If it is not stable, find a low gain stabilizing controller to stabilize it. For instance, design an LQG controller with extremely small output weighting matrix \(Q\) to obtain a stabilizing controller.

\[0=A_pX_o+X_oA^T_p+D_pW_pD^T_p \label{3.120}\]

Then, the open-loop output covariance matrices can be calculated using the following formula.

\[{\hat{Y}}_i{=C}_{p,i}X_oC^T_{p,i}, i=1,\ 2,\ \dots ,\ m \label{3.121}\]

We can generate a sequence of output covariance constraint by varying parameter \(0<γ<1\) as follows.

\[{\overline{Y}}_i\left(l\right)={\tilde{Y}}_i+γ(l)\left({\hat{Y}}_i-{\tilde{Y}}_i\right), i=1,\ 2,\ \dots ,\ m \label{3.122}\]

For example, letting \(γ\left(l\right)=0.1\times l-0.05\ (l=1,\ 2,\ \dots ,\ 10\)) generates 10 upper bound covariance matrices that can be used to design 10 OCC controllers with gains from low to high with \(l=10,\ 9,\ \dots ,\ 1\). Now Example 3.19 demonstrates this control design process.

Example 3.19

Consider the system defined in equation \ref{3.98} with system matrices \(A_p\), \(B_p\), \(D_p\), \(C_p\), and \(M_p\) shown below, where control weighting is \(R=1\), system input noise intensity is \(W_p=1\), and the measurement noise intensity matrix is \(V=0.001I_2\).

\[A_p=\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & 0 & 1 \\ 1 & 0 & 1 \end{array} \right], B_p=D_p=\left[ \begin{array}{c} 0 \\ 0 \\ 1 \end{array} \right], \ C^T_p=M^T_p=\left[ \begin{array}{cc} 1 & 0 \\ 0 & 1 \\ 0 & 0 \end{array} \right]\label{3.123}\]

We will first design an OCC controller using the algorithm proposed before Remark 3.9. The constraints on output variances for outputs \(y_1\) and \(y_2\ (m_1=m_2=1\) and \(m=2\)) are listed below.

\[{\overline{Y}}_1=0.2073, \ {\overline{Y}}_2=0.1128 \label{3.124}\]

For \(ε=0.0001\) with \(β=α=0.5\), the OCC algorithm converges in 51 iterations with the following controller matrices defined in equation \ref{3.100}

\[A_c=\left[ \begin{array}{ccc} -0.9997 & -0.0246 & 0 \\ -1.0246 & -8.9503 & 1.0000 \\ -3.6393 & -49.0405 & -4.2335 \end{array} \right],\ B_c=\left[ \begin{array}{cc} 0.9997 & 1.0246 \\ 1.0246 & 8.9503 \\ 1.2449 & 40.4793 \end{array} \right],\ C^T_c=\left[ \begin{array}{c} -3.3944 \\ -8.4612 \\ -5.2335 \end{array} \right]\label{3.125}\]

and the output variances \(Y_1={\overline{Y}}_1=0.2073\) and \(Y_2={\overline{Y}}_2=0.1128\) with converged output weight matrix \(Q\) as below.

\[Q=\left[ \begin{array}{cc} 4.7332 & 0 \\ 0 & 53.3180 \end{array} \right]>0\label{3.126}\]

Note that \(Q\) is positive definite because achieved variances equal the provided upper bounds, respectively, and the resulting \(Q\) is different from the given \(Q=10I_2\) of LQG control. The resulting closed-loop system poles are the combination of these from \((A_p+B_pC_c)\ (λ_1=-0.3346\), \(λ_{2,3}=-1.9495\pm 1.8320j\)) and these from \((A_p-{LM}_p)\ (λ_4=-1.0015\) and \(λ_{5,6}=-3.9743\pm 3.9745j\)). The system estimation and control gains (\(L\) and \(F\)) are listed in equation \ref{3.127}. We can see that the estimation gain \(L\) is the same as that from the LQG control example (Example 3.18) because the same input and measurement noise intensity matrices are used. In addition, the control gain \(F\) is quite different from the LQG one because the resulting output weighting matrix \(Q\) is quite different from the LQG one (\(Q=10I_2\)).

\[\begin{array}{c} F=-\left[ \begin{array}{ccc} 3.3944 & 8.4612 & 5.2335 \end{array} \right] \\ L^T=\left[ \begin{array}{ccc} 0.9997 & 1.0246 & 1.2449 \\ 1.0246 & 8.9503 & 40.5793 \end{array} \right] \end{array} \label{3.127}\]

Now, we will discuss the tradeoff relationship between control and output variances. The open-loop system is not stable. Therefore, a stabilizing LQR controller is needed with very small weighting matrix \(Q={10}^{-6}I_2\) and \(R=1\), and using MATLAB command Fs=-lqr(Ap,Bp,Qp,R,0)(\(Ap=A_p,\) \(Bp=B_p\), and \(Qp=C^T_pQC_p\), \(Fs=F_s\),) to obtain the stabilizing controller below.

\[F_s=-\left[ \begin{array}{ccc} 2.0000 & 1.3647 & 2.9311 \end{array} \right]\label{3.128}\]

The corresponding stabilizing output covariance matrix \(\hat{Y}\) and lower bound output covariance matrix \(\tilde{Y}\), defined in equation \ref{3.122}, are shown below. In this case, the stabilizing controller is used only to calculate an upper bound for output variance bound matrix \(\hat{Y}\) below.

\[\hat{Y}=\left[ \begin{array}{cc} 0.5904 & 0 \\ 0 & 0.3057 \end{array} \right]>0,\  \tilde{Y}=\left[ \begin{array}{cc} 0.0009997 & 0.0010246 \\ 0.0010246 & 0.0089503 \end{array} \right]\label{3.129}\]

As discussed earlier, assume fifteen-dimensional vector \(γ\) as below. Note that as \(γ\) reduces, the resolution for \(γ\) should become fine because the sensitivity of designed controller gain goes high as \(γ\) reduces.

\[γ^T=\left[0.95,\ 0.85,\ 0.75,\ 0.65,\ 0.55,\ 0.45,\ 0.35,\ 0.3,\ 0.25,\ 0.20,\ 0.15,\ 0.125,\ 0.1,\ 0.075,\ 0.05\right] \label{3.130}\]

Using the formulae from equation \ref{3.122}, fifteen constraints of \({\overline{Y}}_i\left(l\right)\ (i=1,2\), \(l=1,2,\ \dots ,\ 15\)) can be generated based upon \(\hat{Y}\) and \(\tilde{Y}\) defined in equation \ref{3.126}. In addition, fifteen OCC controllers with gains low to high can be designed. In fact, the OCC controller shown earlier is the seventh OCC controller with \(γ\left(7\right)=0.35\). With the fifteen designed OCC controllers, we can calculate the unweighted OCC control and output performance costs (\({\mathcal{E}}_{\infty }\left[u(t)u^T(t-τ)\right]\) and \({\mathcal{E}}_{\infty }\left[y_1(t)y^T_1(t-τ)\right]+{\mathcal{E}}_{\infty }\left[y_2(t)y^T_2(t-τ)\right]\)). They are denoted as \({\mathcal{E}}_{\infty }\left[u^Tu\right]\) and \({\mathcal{E}}_{\infty }\left[y^Ty\right]\), respectively, in Figure 3.34 as horizontal axis and vertical axis.

Note that the blue line with circles in Figure 3.34 is the tradeoff relationship between \({\mathcal{E}}_{\infty }\left[u^Tu\right]\) and \({\mathcal{E}}_{\infty }\left[y^Ty\right]\) without any system modeling error. That is, the designed controller is applied to the control design model without any modeling error. Under this ideal situation, the closed-loop system performance \({\mathcal{E}}_{\infty }\left[y^Ty\right]\) continuous improving as the control energy \({\mathcal{E}}_{\infty }\left[u^Tu\right]\) increases. In practice, however, there are always modeling errors, which could lead to a quite different system performance. Consider the actual physical system model described by modified \(A_p\) and \(B_p\) as below (\(D_p,\ C_p\), and \(M_p\) remain unchanged). Uncertainties are marked in red.

\[A_p=\left[ \begin{array}{ccc} 0 & 1 & 0 \\ 0 & \textcolor{red}{-1.0} & 1 \\ 1 & \textcolor{red}{0.2} & 1 \end{array} \right], B_p=\left[ \begin{array}{c} 0 \\ \textcolor{red}{-0.15} \\ 1 \end{array} \right]\ \label{3.131}\]

Fifteen closed-loop systems were formed using the actual physical model in equation \ref{3.131} based upon the fifteen designed OCC controllers with gains from low to high. Both control and output variances, \({\mathcal{E}}_{\infty }\left[u^Tu\right]\) and \({\mathcal{E}}_{\infty }\left[y^Ty\right]\), are calculated for these fifteen closed-loop systems. They are plotted in Figure 3.34 in red line with stars. We can see that the output variance is higher than those from the closed-loop system evaluated with the design model with increased control variance.

In addition, beyond certain control variance \({\mathcal{E}}_{\infty }\left[u^Tu\right]\) level (around 33), the output performance in terms of output variance \({\mathcal{E}}_{\infty }\left[y^Ty\right]\) gets worse with increased \({\mathcal{E}}_{\infty }\left[y^Ty\right]\). This indicates that due to the modeling error, the best closed-loop system performance in terms of output variance \({\mathcal{E}}_{\infty }\left[y^Ty\right]\) is achieved at a certain level of control effort in terms of \({\mathcal{E}}_{\infty }\left[u^Tu\right]\). Beyond that level, increasing control effort will not improve the closed-loop system performance. Instead, it could make the resulting closed-loop system performance worse or unstable. That is, practically, control effort cannot be increased to a very high level.

Remark 3.11

From Example 3.19, it is important to use the OCC control design algorithm to find a family of controllers with control gains from low to high. As a model-based control gain tuning method, we can evaluate each controller by starting from the lowest gain and gradually increasing the gain until the system performance gets worse. Using this process, we can find the best performance controller from the family of designed OCC controllers. References [7] and [39] address this control tuning method with application to the Hubble Telescope [7] and JPL LSCL space antenna [39].

---