Lecture/week14

src/main.tex

@@ -24,7 +24,8 @@
     \input{sections/stab_and_perf_robustness.tex}
     \input{sections/modern_controller_synthesis.tex}
     \input{sections/h2_hinf_synthesis.tex}
-    \input{sections/nonlinear.tex}
+    \input{sections/nonlinear_systems.tex}
+    \input{sections/mpc.tex}
     \newpage
     \input{sections/Appendix.tex}
 

src/sections/h2_hinf_synthesis.tex

@@ -1,9 +1,21 @@
 \subsection{\texorpdfstring{$\mathcal{H}_2$}{H2} Synthesis}
+The dynamics of the standard setup for modern control synthesis are given by
+\begin{align*}
+    \quad \mathbf{x}(0) & =\mathbf{x}_{0}                                                                              \\
+    \dot{\mathbf{x}}(t) & =\quad \mathbf{A}\mathbf{x}(t)+\mathbf{B}_{w}\mathbf{w}(t)+\mathbf{B}_{u}\mathbf{u}(t)       \\
+    \mathbf{z}(t)       & =\quad \mathbf{C}_{z}\mathbf{x}(t)+\mathbf{D}_{zw}\mathbf{w}(t)+\mathbf{D}_{zu}\mathbf{u}(t) \\
+    \mathbf{y}(t)       & =\quad \mathbf{C}_{y}\mathbf{x}(t)+\mathbf{D}_{yw}\mathbf{w}(t)+\mathbf{D}_{yu}\mathbf{u}(t)
+\end{align*}
+LQR, LQE, LQG are applied to special cases of the standard setup.
 \subsubsection{LQR}
+
+\newpar{}
+\ptitle{Problem Statement}
+
 In the LQR problem we assume
 \begin{itemize}
     \item Full state feedback: $\mathbf{C_y}=\mathbf{\mathbb{I}}$
-    \item No disturbance input: $w=0$
+    \item No disturbance input: $w=0$ ($\mathbf{B}_w=\mathbf{D}_{zw}=\mathbf{D}_{yw}=\mathbf{0}$)
     \item $\mathbf{C}_z^{\mathsf{T}}\mathbf{D}_{zu}=0$
 \end{itemize}
 and try to find a control signal $u(t,x)$ that minimizes
@@ -26,13 +38,20 @@ \subsubsection{LQR}
     \|\mathbf{z}\|_2^2=\int_0^{+\infty}(\mathbf{x}^{\mathsf{T}}\mathbf{Q}\mathbf{x}+\mathbf{u}^{\mathsf{T}}\mathbf{R}\mathbf{u})\;dt
 \end{equation*}
 
-\ptitle{ARE}
-\noindent\begin{gather*}
-    \mathbf{A}^{\mathsf{T}}\mathbf{X}_F + \mathbf{X}_F \mathbf{A} - \mathbf{X}_F \mathbf{B}_u{\left(\mathbf{D}_{zu}^{\mathsf{T}}\mathbf{B}_{zu}\right)}^{-1} \mathbf{B}_u^{\mathsf{T}} \mathbf{X}_F + \mathbf{C}_z^{\mathsf{T}}\mathbf{C}_z = \mathbf{0}\\
+\newpar{}
+\ptitle{Solving for the Controller Matrix}
+
+One can show that minimizing $\|\mathbf{z}\|_2^2$ imposes the following ARE on $\mathbf{X}_F$:
+\begin{equation*}
+    \mathbf{A}^{\mathsf{T}}\mathbf{X}_F + \mathbf{X}_F \mathbf{A} - \mathbf{X}_F \mathbf{B}_u{\left(\mathbf{D}_{zu}^{\mathsf{T}}\mathbf{B}_{zu}\right)}^{-1} \mathbf{B}_u^{\mathsf{T}} \mathbf{X}_F + \mathbf{C}_z^{\mathsf{T}}\mathbf{C}_z = \mathbf{0}
+\end{equation*}
+Solving this ARE yields the controller
+\begin{equation*}
     \mathbf{F} = -{\left(\mathbf{D}_{zu}^{\mathsf{T}}\mathbf{D}_{zu}\right)}^{-1} \mathbf{B}_u^{\mathsf{T}}\mathbf{X}_F
-\end{gather*}
+\end{equation*}
 where $\mathbf{X}_F$ is symmetric and positive semidefinite.
 
+\newpar{}
 \ptitle{Technical Considerations}
 
 \begin{enumerate}
@@ -45,6 +64,7 @@ \subsubsection{LQR}
           Is for convenience.
 \end{enumerate}
 
+\newpar{}
 \ptitle{Remarks}
 
 \begin{itemize}
@@ -54,6 +74,10 @@ \subsubsection{LQR}
 
 
 \subsubsection{LQE}
+
+\newpar{}
+\ptitle{Problem Statement}
+
 In the LQE problem we assume
 \begin{itemize}
     \item That $\mathbf{u}$ takes the role uf the observer update
@@ -82,11 +106,16 @@ \subsubsection{LQE}
 \end{equation*}
 
 \newpar{}
-\ptitle{ARE}
-\noindent\begin{gather*}
-    \mathbf{A}\mathbf{Y}_L + \mathbf{Y}_L \mathbf{A}^{\mathsf{T}} - \mathbf{Y}_L \mathbf{C}_y^{\mathsf{T}}{\left(\mathbf{D}_{yw}\mathbf{D}_{yw}^{\mathsf{T}}\right)}^{-1} \mathbf{C}_y \mathbf{Y}_L + \mathbf{B}_w\mathbf{B}_w^{\mathsf{T}} = \mathbf{0}\\
-    \mathbf{F} = -\mathbf{Y}_L\mathbf{C}_y^{\mathsf{T}}{\left(\mathbf{D}_{yw}\mathbf{D}_{yw}^{\mathsf{T}}\right)}^{-1}
-\end{gather*}
+\ptitle{Solving for the Observer Matrix}
+
+One can show that the given problem statement on $\mathbf{u}$ imposes the following ARE on $\mathbf{Y}_L$:
+\begin{equation*}
+    \mathbf{A}\mathbf{Y}_L + \mathbf{Y}_L \mathbf{A}^{\mathsf{T}} - \mathbf{Y}_L \mathbf{C}_y^{\mathsf{T}}{\left(\mathbf{D}_{yw}\mathbf{D}_{yw}^{\mathsf{T}}\right)}^{-1} \mathbf{C}_y \mathbf{Y}_L + \mathbf{B}_w\mathbf{B}_w^{\mathsf{T}} = \mathbf{0}
+\end{equation*}
+The optimal observer matrix is then given by
+\begin{equation*}
+    \mathbf{L} = -\mathbf{Y}_L\mathbf{C}_y^{\mathsf{T}}{(\underbrace{\mathbf{D}_{yw}\mathbf{D}_{yw}^{\mathsf{T}}}_{=\mathbf{R}})}^{-1}    
+\end{equation*}
 where $\mathbf{Y}_L$ is symmetric and positive semidefinite.
 
 \newpar{}
@@ -130,6 +159,11 @@ \subsubsection{LQG}
 \end{equation*}
 where $\mathbf{Y}$ is the stabilizing solution to the corresponding ARE.
 
+\newpar{}
+\ptitle{Technical Considerations}
+
+Remember that LGQ combines LQR and LQE.
+
 \subsection{\texorpdfstring{$\mathcal{H}_\infty$}{H-infinity} Synthesis}
 \ptitle{Differentiation from $\mathcal{H}_2$ Synthesis}
 \begin{itemize}
@@ -151,27 +185,27 @@ \subsection{\texorpdfstring{$\mathcal{H}_\infty$}{H-infinity} Synthesis}
 \end{itemize}
 
 \subsubsection{Optimal \texorpdfstring{$\mathcal{H}_\infty$}{H-infinity} Synthesis}
-\begin{itemize}
-    \item In principle, we would like to find a controller K that minimizes the energy ($\mathcal{L}_2$) gain of the CL system i.e.
-          \begin{equation*}
-              \left\|\mathbf{T}_{zw}\right\|_{\mathcal{H}_\infty}=\sup_{w\neq0}\frac{\left\|\mathbf{z}\right\|_{\mathcal{L}_2}}{\left\|\mathbf{w}\right\|_{\mathcal{L}_2}}
-          \end{equation*}
-    \item However, the optimal controller(s) are such that
-          \begin{enumerate}
-              \item $\sigma_{\max}(\mathbf{T}_{zw}(j\omega))$ is a constant over all frequencies, the response does not roll off at high frequencies (\textit{Bode's integral law}).
-              \item The controller is not strictly proper. (The optimal controller is not unique)
-              \item Computing an optimal controller is numerically challenging.
-          \end{enumerate}
-\end{itemize}
+In principle, we would like to find a controller K that minimizes the energy ($\mathcal{L}_2$) gain of the CL system i.e.
+\begin{equation*}
+    \left\|\mathbf{T}_{zw}\right\|_{\mathcal{H}_\infty}=\sup_{w\neq0}\frac{\left\|\mathbf{z}\right\|_{\mathcal{L}_2}}{\left\|\mathbf{w}\right\|_{\mathcal{L}_2}}
+\end{equation*}
+However, the optimal controller(s) are such that
+\begin{enumerate}
+    \item $\sigma_{\max}(\mathbf{T}_{zw}(j\omega))$ is a constant over all frequencies, the response does not roll off at high frequencies (\textit{Bode's integral law}).
+    \item The controller is not strictly proper. (The optimal controller is not unique)
+    \item Computing an optimal controller is numerically challenging.
+\end{enumerate}
 
 \subsubsection[Suboptimal H-infinity Synthesis]{Suboptimal $\mathcal{H}_\infty$ Synthesis}
+\newpar{}
+\ptitle{Problem Statement}
 In practice one pursues a sub-optimal design, i.e., given $\gamma>0$, find a controller $K$ such that
 \begin{equation*}
     \|\mathbf{T}_{zw}\|_{\mathcal{H}_\infty}<\gamma
 \end{equation*}
 if one exists. This can be formulated mathematically as a cost function:
 \begin{equation*}
-    \|\mathbf{z}\|_{\mathcal{L}_2}^2-\gamma^2\|\mathbf{w}\|_{\mathcal{L}_2}^2
+    \|\mathbf{z}\|_{\mathcal{L}_2}^2-\gamma^2\|\mathbf{w}\|_{\mathcal{L}_2}^2 < 0
 \end{equation*}
 where we search for the smallest $\gamma$ such that the controller can achieve negative cost. We \textbf{need} $\boldsymbol{\gamma}\mathbf{\le1}$ for a stabilizing controller as otherways for larger gamma one could easily achieve negative cost even though the energy of the disturbance gets not damped at all! The formula can be understood as
 \begin{itemize}
@@ -180,6 +214,29 @@ \subsubsection{Optimal \texorpdfstring{$\mathcal{H}_\infty$}{H-infinity} Synthes
     \item The aforementioned formula becomes negative if the disturbance energy term is larger than the (hopefully damped) performance output $\mathbf{z}$. A small $\gamma$ gives an even better controller as it is more difficult to achieve negative cost then.
 \end{itemize}
 
+\newpar{}
+\ptitle{Conditions}
+
+For feasibility one must have
+\begin{enumerate}
+    \item $(\mathbf{A}_{ext}, \mathbf{B}_2)$ must be stabilizable
+    \item $(\mathbf{C}_{ext, y},\mathbf{A}_{ext})$ must be detectable
+    \item $\begin{bmatrix}\mathbf{A}-j\omega \mathbf{I}&\mathbf{B}_w\end{bmatrix},\begin{bmatrix}\mathbf{A}^{\prime}-j\omega \mathbf{I}&\mathbf{C}_z^{\prime}\end{bmatrix}$ must have full row rank
+\end{enumerate}
+\newpar{}
+A controller $\mathbf{K}$ fulfilling the cost inequality exists, only if 
+\begin{enumerate}
+    \item The following ARE has a solution for $\mathbf{X}_{\infty}$
+    \begin{equation*}
+        \mathbf{A}^{\prime}\mathbf{X}_{\infty}+\mathbf{X}_{\infty}\mathbf{A}+\mathbf{C}_{z}^{\prime}\mathbf{C}_{z}=\mathbf{X}_{\infty}(\mathbf{B}_{u}\mathbf{B}_{u}^{\prime}-\gamma^{-2}\mathbf{B}_{w}\mathbf{B}_{w}^{\prime})\mathbf{X}_{\infty}
+    \end{equation*}
+    \item The following ARE has a solution for $\mathbf{Y}_{\infty}$
+    \begin{equation*}
+        \mathbf{A}\mathbf{Y}_\infty + \mathbf{Y}_\infty \mathbf{A}^{\prime} + \mathbf{B}_w^{\prime} \mathbf{B}_w = \mathbf{Y}_\infty (\mathbf{C}_y \mathbf{C}_y^{\prime} - \gamma^{-2} \mathbf{C}_z \mathbf{C}_z^{\prime}) \mathbf{Y}_\infty
+    \end{equation*}
+    \item The matrix $\gamma^2 \mathbf{I} - \mathbf{Y}_{\infty} \mathbf{X}_{\infty}$ is positive definite
+\end{enumerate}
+
 \newpar{}
 \ptitle{Bisection Algorithm}
 
@@ -190,6 +247,20 @@ \subsubsection{Optimal \texorpdfstring{$\mathcal{H}_\infty$}{H-infinity} Synthes
     \item Repeat from step 2 until $\gamma_+-\gamma_-<\epsilon$.
     \item Return $\mathbf{K}_+$.
 \end{enumerate}
+\newpar{}
+Remarks
+\begin{itemize}
+    \item The conditions are only fulfilled if $\gamma\le 1$
+    \item For $\gamma>1$ use relaxed weights
+    \item The final controller gains are 
+    \begin{equation*}
+        \mathbf{F}_u=-\mathbf{B}_u^{\prime}\mathbf{X}_\infty,\quad \mathbf{F}_w=\frac1{\gamma^2}\mathbf{B}_w^{\prime}\mathbf{X}_\infty 
+    \end{equation*}
+    \item The final observer gain is 
+    \begin{equation*}
+        \mathbf{L}=-{(\mathbf{I}-\gamma^{-2}\mathbf{Y}_\infty \mathbf{X}_\infty)}^{-1}\mathbf{Y}_\infty \mathbf{C}_y^{\prime}
+    \end{equation*}
+\end{itemize}
 
 \newpar{}
 \ptitle{Simplified Setup}
@@ -206,13 +277,13 @@ \subsubsection{Applying Frequency Weights}
 \begin{itemize}
     \item In order for $\mathcal{H}_\infty$ synthesis to be possible, frequency weights must be \textbf{stable and proper}.
     \item Integrators (not asymptotically stable) can be approximated by very slow poles.
-    \item MATLAB: 
-    \begin{itemize}
-        \item \texttt{Pinf=augw(G,W1,0,W3)}\\
-        creates the augmented (weighted) system
-        \item \texttt{[Kinf, CL, GAM, INFO]=hinfsyn(Pinf,1,1)}\\
-        applies optimal $\mathcal{H}_\infty$ synthesis
-        \item \texttt{[K09, CL, GAM, I]=hinfsyn(Pinf,1,1,0.9)}\\
-        applies suboptimal $\mathcal{H}_\infty$ synthesis targeting $\gamma=0.9$. This can yield a more feasible controller (less control effort)
-    \end{itemize}
+    \item MATLAB:
+          \begin{itemize}
+              \item \texttt{Pinf=augw(G,W1,0,W3)}\\
+                    creates the augmented (weighted) system
+              \item \texttt{[Kinf, CL, GAM, INFO]=hinfsyn(Pinf,1,1)}\\
+                    applies optimal $\mathcal{H}_\infty$ synthesis
+              \item \texttt{[K09, CL, GAM, I]=hinfsyn(Pinf,1,1,0.9)}\\
+                    applies suboptimal $\mathcal{H}_\infty$ synthesis targeting $\gamma=0.9$. This can yield a more feasible controller (less control effort)
+          \end{itemize}
 \end{itemize}

src/sections/modern_controller_synthesis.tex

@@ -92,7 +92,7 @@ \subsection{Frequency Weights}
     \|\Delta(s)\|_{\mathcal{H}_\infty} < 1 \\
     \|M(s)\|_{\mathcal{H}_\infty} < 1
 \end{gather*}
-where $M(s)$ is the TF from the output of $\delta$ ti its input.
+where $M(s)$ is the TF from the output of $\Delta$ to its input.
 
 \subsubsection{First Order Weights}
 \ptitle{Lowpass}
@@ -297,7 +297,7 @@ \subsection{Youla's Q Parameterization}
 \begin{equation*}
     C_0(s)=-{X_0(s)}^{-1}Y_0(s)
 \end{equation*}
-is the TF when $Q(s)=0$ with controller gain $K_0$ and observer gail $L_0$
+is the TF when $Q(s)=0$ with controller gain $K_0$ and observer gain $L_0$
 \noindent\begin{align*}
     X_0(s) & = -K_0{(sI-A-L_0C)}^{-1}B   \\
     Y_0(s) & = -K_0{(sI-A-L_0C)}^{-1}L_0
@@ -309,7 +309,7 @@ \subsection{Youla's Q Parameterization}
 \end{align*}
 where $L_0$ is an observer gain that would stabilize the error dynamics.
 \newpar{}
-It can be show that the interconnection of the system is stable if the \textbf{Bezout identity} is fulfilled
+It can be shown that the interconnection of the system is stable if the \textbf{Bezout identity} is fulfilled
 \begin{equation*}
     D(s)X(s)-N(s)Y(s)=1
 \end{equation*}

src/sections/mpc.tex

@@ -0,0 +1,115 @@
+\section{Model Predictive Control (MPC)}
+
+For a discrete-time nonlinear control system of the form
+\begin{align*}
+    \mathbf{x}[k+1] & = f(\mathbf{x}[k], \mathbf{u}[k]) \\
+    \mathbf{y}[k]   & = h(\mathbf{x}[k],\mathbf{u}[k])
+\end{align*}
+a general description of the cost function and the constraints with the horizon length $H$ can be formulated as following
+\begin{align*}
+    \min_{u[0], \ldots, u[H-1]} J_{H}(\mathbf{x,u}) & := \sum_{k=0}^{H-1}g(\mathbf{x}[k],\mathbf{u}[k]) \\
+    \mathbf{x}[k+1]                                 & = f(\mathbf{x}[k]. \mathbf{u}[k])                 \\
+    \mathbf{x}[0]                                   & = \mathbf{x}_0                                    \\
+    \mathbf{x}[k]                                   & \in X                                             \\
+    \mathbf{u}[k]                                   & \in U                                             \\
+    k                                               & = 0,\ldots,H-1
+\end{align*}
+By solving these equations one can obtain the best control input $u$ for the next step. This is then done over and over again.
+
+\newpar{}
+
+The remaining costs from $k=H, \ldots, \infty$ are called the tail. If the prediction is too short-sighted (hence the tail is to large) the system can get unstable. To fix this, one can
+\begin{itemize}
+    \item shorten the tail (enlarge the horizon). This introduces more computational effort each step.
+    \item make the tail zero (or sufficiently small). See Section~\ref{mpc_terminal_constraint}.
+    \item replace the tail with an estimate. See Section~\ref{mpc_terminal_cost}.
+\end{itemize}
+
+The formulation becomes
+\begin{align*}
+    \min_{u[0], \ldots, u[H-1]} J_{H}(\mathbf{x,u}) & := \sum_{k=0}^{H-1}g(\mathbf{x}[k],\mathbf{u}[k]) \textcolor{blue}{+V(\textbf{x}[H])} \\
+    \mathbf{x}[k+1]                                 & = f(\mathbf{x}[k]. \mathbf{u}[k])                                                     \\
+    \mathbf{x}[0]                                   & = \mathbf{x}_0                                                                        \\
+    \mathbf{x}[k]                                   & \in X                                                                                 \\
+    \mathbf{u}[k]                                   & \in U                                                                                 \\
+    \textcolor{blue}{\mathbf{x}[H]}                 & \textcolor{blue}{\;\in X_H}                                                           \\
+    k                                               & = 0,\ldots,H-1
+\end{align*}
+
+\subsection{Terminal Constraints}\label{mpc_terminal_constraint}
+
+If the state is forced to an equilibrium at the end of the horizon it wil remain there and the tail will have no further effect. On the down side this may put excessive pressure on the actuators/control effort.
+
+\subsection{Terminal Cost}\label{mpc_terminal_cost}
+
+The best terminal cost estimate would be $V = J^*_\infty$ the cost from the end of the horizon up to infinity.
+
+Because this is not determinable we are looking for a \textbf{global Control Lyapunov Function (CLF)} such that
+\begin{equation*}
+    V(\mathbf{x}) \geq \mathbf{c}_v {\lVert \mathbf{x} \rVert}^2
+\end{equation*}
+and
+\begin{equation*}
+    \min_{u}(V(f(\mathbf{x,u})) + \mathbf{x}^T \mathbf{Qx + u}^T \mathbf{R u}) \leq 0 \quad \forall \; \mathbf{x}
+\end{equation*}
+which is an upper bound for the optimal cost
+\begin{equation*}
+    V(\mathbf{x}) \geq J^*_\infty (\mathbf{x})
+\end{equation*}
+If we ensure that the terminal cost of the next step $V(f(\mathbf{x,u}))$ is less or equal to the current terminal cost $V(\mathbf{x})$ we get a stabilizing control law.
+
+However, the MPC control law will have better performance -- if not optimal, i.e. % TODO: Maybe explain this in more detail
+\begin{equation*}
+    J^*_\infty(\mathbf{x}) \leq J^*_H(\mathbf{x}) \leq J^*_{CLF}(\mathbf{x})
+\end{equation*}
+
+An example of a CLF candidate may be given by the solution $\mathbf{P}$ of the Algebraic Ricatti Equation of the LQR problem, i.e., given $V(x) = \mathbf{x}^T \mathbf{Px}$.
+
+\subsection{Solving the Finite Horizon Optimization}
+
+For an unconstrained, differentiable and convex optimization problem one can
+\begin{itemize}
+    \item set the gradient to zero to find the stationary points.
+    \item use a gradient descent method.
+\end{itemize}
+
+In most cases one deals with limiting constraints. To incorporate these into the cost function a so called \textbf{barrier function} is added.
+\begin{equation*}
+    V(\mathbf{x,u}) - \log(-c(\mathbf{x,u}))
+\end{equation*}
+Note that the function $c$ is defined in a way such that
+\begin{equation*}
+    c(\mathbf{x,u}) \leq 0
+\end{equation*}
+Hence the barrier function is not defined for $(x,u)$ not satisfying the constraints, and it becomes very large when approaching the boundaries of the constraint set.
+
+To solve this kind of optimization problem the \textbf{barrier interior-point method} can be used.
+
+\subsubsection{Barrier Interior-Point Metod}
+Interior point methods solve \textbf{convex problems} of the form
+\begin{align*}
+    \min \quad         & f_0(x)                                 \\
+    \text{s.t.:} \quad & f_i(x) \leq 0 , \quad 0 = 1, \ldots, m \\
+                       & Ax = b
+\end{align*}
+where the scalar functions $f_0,f_1,\ldots,f_m$ are convex and twice differentiable.
+\newpar{}
+The problem is converted into a minimization with affine equality constraints by adding \textbf{barrier functions} in the following way
+\begin{align*}
+    \min \quad        & f_0(x) - \frac{1}{t}\sum_{i=1}^{m}\log(-f_i(x)) \\
+    \text{s.t.:}\quad & Ax=b
+\end{align*}
+
+As $t \to +\infty$, the objective function gets closer to $f_0$, but diverges very steeply at the boundary of the feasible set.
+\newpar{}
+The problem is then solved iteratively using gradient descent or Newton method for $t=\mu^k t_0$ with $t_0>0, \mu > 1$ where $k$ is the iteration variable. The solution of the pervious step is used as a starting point for the next step.
+
+% TODO: Add graphical example like the one drawn on the blackboard.
+
+\newpar{}
+\ptitle{Remarks:}
+\begin{itemize}
+    \item The sequence of intermediate optimal solutions forms what is known as the \textbf{cental path}, which is always contained in the interior of the feasible set.
+    \item This method terminates with $\mathcal{O}(\sqrt{m})$ time.
+\end{itemize}
+