Singular Perturbation Methods

A strong control action can force a system to have fast and slow transients, i.e. to behave like a singularly perturbated system. In feedback systems, the strong control action is achieved by high feedback gain, e.g. as a result of a cheap optimal control problem. We give a brief account for the methods to be used in our analysis of the cheap control regulator. Basic results and singular pertubation techniques applied to different problems can be found in e.g. Kokotovič et al. [6].

Consider the singular perturbation model of finite-dimensional dynamic systems

(1) $\begin{equation*}\begin{array}{ccl}\dot x&=&f(x,z,\epsilon ,t),\quad x(t_0)=x^0,\quad x\in R^n,\\\epsilon \dot z&=&g(x,z,\epsilon ,t),\quad z(t_0)=z^0,\quad z\in R^m,\end{array}\end{equation*}$

where $f$ and $g$ are assumed to be sufficiently many times continuously differentiable. When we set $\epsilon =0$ , the dimension of the state-space of (1) is reduced from $n+m$ to $n$ as the second equation degenerates into the algebraic (or transcendental) equation

(2) $\begin{equation*}0=g(\bar{x},\bar{z},0,t),\end{equation*}$

where $\bar{x}$ and $\bar{z}$ belong to a system with $\epsilon =0$ . We say that the model (1) is in standard form if and only if the following assumption is satisfied.

Assumption 2.1: In a domain of interest, the equation (2) has $k\ge 0$ distinct real roots

(3) $\begin{equation*}z=\bar \phi_i(\bar x,t),\quad i=1,2,…,k.\end{equation*}$

This assumption ensures that a well-defined $n$ -dimensional reduced model will correspond to each root (3). Substituting (3) into the first equation of (1) gives the $i$ th reduced (quasi-steady-state) model

(4) $\begin{equation*}\dot {\bar x}=f(\bar x,\bar \phi_i(\bar x,t),0,t)=\bar f(\bar x,t),\quad\bar x(t_0)=x^0,\end{equation*}$

where we use the same initial condition for $\bar x(t)$ as for $x(t)$ .

Singular perturbation cause a multi-time-scale behavior of dynamic systems characterized by slow and fast transients in the system response. The slow response is approximated by (4) while the discrepancy between the response of the reduced model (4) and the full model (1) is the fast transient. In removing the variable $z$ from the original system (1) and replacing it by its quasi-steady-state $\bar z$ , we have to consider the fact that $\bar z$ is not free to start at $t_0$ from $z^0$ , the initial condition for $z$ . There may be a large discrepancy between the initial condition $z^0$ and the initial value of $\bar z$ ,

$\bar z(t_0)=\bar {\phi}(\bar x(t_0),t_0).$

Thus, $\bar z$ cannot be a uniform approximation of $z$ . The best approximation we can obtain is that

(5) $\begin{equation*}z=\bar z(t)+O(\epsilon )\end{equation*}$

will hold on an interval excluding $t_0$ , i.e. for all $t\in [t_1,T]$ , where $t_1>t_0$ . However, we can constrain the quasi-steady-state $\bar x$ to start from the prescribed initial condition $x^0$ , and therefore the approximation of $x$ by $\bar x$ may be uniform. In other words,

(6) $\begin{equation*}x=\bar x(t)+O(\epsilon )\end{equation*}$

will hold on an interval including $t_0$ , i.e. for all $t\in [t_0,T]$ on which $\bar x(t)$ exists. The approximation (5) establishes that during an initial “boundary layer” interval $[t_0,t_1]$ , the original variable $z$ approaches $\bar z$ and then, during $[t_1,T]$ , remains close to $\bar z$ . The speed of $z$ can be large, $\dot z=g/\epsilon$ . In fact, having set $\epsilon$ equal to zero in (1), we have made the transient instantaneous whenever
$g\not= 0$ . To see whether $z$ will escape to infinity during this transient or converge to its quasi-steady-state, we analyze $\epsilon \dot z$ , which may be finite, even when $\epsilon$ tends to zero and $\dot z$ tends to infinity. We change coordinates according to

$\epsilon \frac {dz}{dt}=\frac {dz}{d\tau }$

so that

$\frac {d\tau }{dt}=\frac 1{\epsilon },$

and use $\tau =0$ as the initial value at $t=t_0$ . We then have the new time variable

$\tau=\frac {t-t_0}{\epsilon },~\tau=0 {\rm~at}~t=t_0,$

so that, if $\epsilon$ tends to zero, $\tau$ tends to infinity even for fixed $t$ only slightly larger than $t_0$ . On the other hand, while $z$ and $\tau$ almost instantaneously change, $x$ remains very near its initial value $x^0$ . Consider the so-called boundary layer correction $\hat z=z-\bar z$ satisfying the boundary layer system

(7) $\begin{equation*}\frac {d\hat z}{d\tau}=g(x^0,\hat z(\tau )+\bar z(t_0),0,t_0)\end{equation*}$

with the initial condition $z^0-\bar z(t_0)$ , where $x^0$ and $t_0$ are fixed parameters. The solution $\hat z(t)$ of this initial value problem is used as a boundary layer correction of (5) for a possibly uniform approximation of $z$ :

(8) $\begin{equation*}z=\bar z(t)+\hat z(\tau )+O(\epsilon ).\end{equation*}$

Clearly, $\bar z$ is the slow transient of $z$ and $\hat z(\tau )$ is the fast transient of $z$ . For the corrected approximation (8) to converge, after a short period, to the slow approximation (5), the
correction term $\hat z(\tau )$ must decay as $\tau \rightarrow \infty$ to an $O(\epsilon )$ quantity. Note that in the slow time scale $t$ this decay is rapid since

$\frac {d\hat z(\tau )}{dt}=\frac {d\hat z(\tau )}{d\tau }\frac {d\tau }{dt}=\frac 1{\epsilon }\frac {d\hat z(\tau )}{d\tau }.$

The stability properties of the boundary layer system (7), which are crucial for the approximations (5), (6) and (8) to hold can be stated as two assumptions.

Assumption 2.2: The equilibrium $\hat z(\tau )=0$ of (7) is assyptotically stable uniformly in $x^0$ and $t_0$ , and $z^0-\bar z(t_0)$ belongs to its domain of attraction, so $\hat z(\tau )$ exists for $\tau \ge 0$ .

If this assumption is satisfied, then

$\lim_{\tau \to \infty}\hat z(\tau )=0$

uniformly in $x^0$ , $t_0$ . Thus, $z$ will come close to its quasi-steady-state $\bar z$ at some time $t_1>t_0$ . To ensure that $z$ stays close to $\bar z$ , we think as if any instant $t\in [t_1,T]$ can be the initial instant, and make the following assumption about the linearization of (7).

Assumption 2.3: The eigenvalues of $\partial g/\partial z$ evaluated, for $\epsilon =0$ , along $\bar x(t)$ , $\bar z(t)$ , have real parts smaller than a fixed negative number, i.e.

(9) $\begin{equation*}{\rm Re}(\sigma [\frac {\partial g}{\partial z}])\le -c<0.\end{equation*}$

Both assumptions desribe a strong stability property of the boundary layer system (7). If $z^0$ is assumed to be sufficiently close to $\bar z(t_0)$ , then assumption 2.3 encompasses assumption 2.2. We also note that the nonsingularity of $\partial g/\partial z$ implies that the root $\bar z(t)$ is distinct as required by assumption 2.1. We have the following fundamental result, the so called “Tikhonov’s theorem”.

Theorem 2.1: If assumptions 2.2 and 2.3 are satisfied, then the approximation (6), (8) is valid for all $t\in [t_0,T]$ , and there exists $t_1\ge t_0$ such that (5) is valid for all $t\in [t_1,T]$ .

Under certain circumstances, the assumption 2.3 can be relaxed to

${\rm Re}(\sigma [\frac {\partial g}{\partial z}])\not= 0.$

This relaxed form of assumption 3.2 can be used when dealing with optimal trajectories having boundary layers at both ends.

Now, consider the problem of finding a control $u(t)$ that steers the state $x,z$ of the singularly perturbated system

(10) $\begin{equation*}\begin{array}{ccl}\dot x&=&f(x,z,t,\epsilon ,u),\quad x\in R^n,\quad u\in R^r\\\epsilon \dot z&=&g(x,z,t,\epsilon ,u),\quad z\in R^m,\end{array}\end{equation*}$

from $x(0)=x^0$ , $z(0)=z^0$ to $x(1)=x^1$ , $(1)=z^1$ while minimizing the cost functional

(11) $\begin{equation*}J=\int_0^1 V(x,z,t,\epsilon ,u) dt\end{equation*}$

This problem can be simplified by neglecting $\epsilon$ in two ways. First, an optimality condition can be formulated for the exact problem and simplified setting $\epsilon =0$ . The result will be a reduced optimality condition. Second, by neglecting $\epsilon$ in the system (10) the same type of optimality condition can be formulated for the reduced system. When the obtained optimality conditions are identical, it is said that the reduced system is formally correct. First, we formulate a necessary optimality condition for the exact problem, using the Lagrangian

$L=V+\lambda_1^T (\dot x-f)+\lambda_2^T(\epsilon \dot z-g),$

where $\lambda_1$ and $\lambda_2$ are the adjoint variables associated with $x$ and $z$ respectively. The Lagrangian necessary condition

(12) $\begin{equation*}\begin{array}{ccl}\frac {\partial L}{\partial \lambda_1}&=&0\\\frac {\partial L}{\partial \lambda_2}&=&0\\\frac {\partial L}{\partial u}&=&0\\\frac {\partial L}{\partial x}-\frac d{dt}\frac {\partial L}{\partial \dot x}&=&0\\\frac {\partial L}{\partial z}-\frac d{dt}\frac {\partial L}{\partial \dot z}&=&0\end{array}\end{equation*}$

give

(13) $\begin{equation*}\begin{array}{ccl}\dot x&=&f,\quad x(0)=x^0\\\epsilon \dot z&=&g,\quad z(0)=z^0\\\nabla_u V&=&f_u^T\lambda_1+g_u^T\lambda_2\\\dot \lambda_1&=&\nabla_x L\\\epsilon \dot \lambda_2&=&\nabla_z L\end{array}\end{equation*}$

The reduced necessary condition is now obtained by setting $\epsilon =0$ in (13) and disregarding the requirement that $z(0)=z^0$ , which must be dropped because the second and fifth equation in (13)
degenerate to

(14) $\begin{equation*}\begin{array}{ccl}0&=&g\\0&=&\nabla_z L\end{array}\end{equation*}$

and are not free to satisfy arbitrary boundary conditions.

On the other hand, the reduced problem is obtained by setting $\epsilon =0$ in (10) and dropping the requirement that $z(0)=z^0$ and $z(1)=z^1$ , i.e. the reduced problem is defined as

(15) $\begin{equation*}\begin{array}{ccl}J&=&\int_0^1V(x,z,t,0,u) dt\\\dot x&=&f(x,z,t,0,u),\quad x(0)=x^0,\\0&=&g(x,z,t,0,u),\quad x(1)=x^1,\end{array}\end{equation*}$

Applying the Lagrangian necessary condition to this problem, we arrive at

(16) $\begin{equation*}\begin{array}{ccl}\dot x&=&f(x,z,t,0,u)\\0&=&g(x,z,t,0,u)\\\nabla_u V(x,z,t,0,u)&=&f_u^T(x,z,t,0,u)\lambda_1+g_u^T(x,z,t,0,u)\lambda_2\\\dot \lambda_1&=&\nabla_x L\\0&=&\nabla_z L\end{array}\end{equation*}$

subject to $x(0)=x^0$ , $x(1)=x^1$ . By setting $\epsilon =0$ in (13), and comparing the result with (16), we come to the following conclusion.

Lemma 2.1: The reduced problem (15) is formally correct.

Using the Hamiltonian form of (15), we define

$H=V+\lambda_1^Tf+\lambda_2^Tg$

and we get the equations

(17) $\begin{equation*}\begin{array}{ccl}\dot x&=&\nabla_{\lambda_1}H=f,\quad x(0)=x^0,\quad x(1)=x^1,\\\epsilon \dot z&=&\nabla_{\lambda_2}H=g,\quad z(0)=z^0,\quad z(1)=z^1,\\0&=&\nabla_uH\\\dot \lambda_1&=&-\nabla_xH\\\epsilon \dot \lambda_2&=&-\nabla_zH\end{array}\end{equation*}$

Note that we make use of $g$ rather than $g/\epsilon$ in the definition of $H$ . This means that the Hamiltonian adjoint variable associated with $z$ is $\epsilon \lambda_2$ . Applying the Hamiltonian condition to the problem

(18) $\begin{equation*}\begin{array}{ccl}J&=&\frac 1{2}\int_0^1(y^Ty+u^TRu)dt,\quad R>0\\\dot x&=&A_{11}x+A_{12}z+B_1u,\quad x(0)=x^0,\quad x(1)=x^1,\\\epsilon \dot z&=&A_{21}x+A_{22}z+B_2u,\quad z(0)=z^0,\quad z(1)=z^1,\\y&=&C_1x+C_2z.\end{array}\end{equation*}$

Using the Hamiltonian

$\begin{equation*}\begin{array}{ccl}H&=&y^Ty/2+u^TRu/2+\lambda_1^T(A_{11}x+A_{12}z+B_1u)+\\ &&\lambda_2^T(A_{21}x+A_{22}z+B_2u)\end{array}\end{equation*}$

we obtain

(19) $\begin{equation*}\begin{array}{ccl}u&=&-R^{-1}(B_1^T\lambda_1+B_2^T\lambda_2)\\\dot \lambda_1&=&-A_{11}^T\lambda_1-A_{21}^T\lambda_2-C_1^TC_1x-C_1^TC_2z\\\epsilon \dot \lambda_2&=&-A_{12}^T\lambda_1-A_{22}^T\lambda_2-C_2^TC_1x-C_2^TC_2z\end{array}\end{equation*}$

Substitution into (18) results in the boundary value problem

(20) $\begin{equation*}\begin{array}{ccl}\left(\begin{array}{c}\dot x\\c\\\dot \lambda_1\end{array}\right)&=&\left(\begin{array}{cc}A_{11}& -S_{11}\\ -Q_{11}& -A_{11}^T\end{array}\right)\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)+\\&+&\left(\begin{array}{cc}A_{12}& -S_{12}\\ -Q_{12}& -A_{21}^T\end{array}\right)\left(\begin{array}{c}z\\ \lambda_2\end{array}\right)\\\epsilon \left(\begin{array}{c}\dot z\\ \dot \lambda_2\end{array}\right)&=&\left(\begin{array}{cc}A_{21}& -S_{12}^T\\ -Q_{12}^T& -A_{12}^T\end{array}\right)\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)+\\&+&\left(\begin{array}{cc}A_{22}& -S_{22}\\ -Q_{22}& -A_{22}^T\end{array}\right)\left(\begin{array}{c}z\\ \lambda_2\end{array}\right)\end{array}\end{equation*}$

and the boundary conditions

$x(0)=x^0,\quad z(0)=z^0,\quad x(1)=x^1,\quad z(1)=z^1,$

where

$Q_{ij}=C_i^TC_j$

and

$S_{ij}=B_iR^{-1}B_j^T,~~i,j=1,2.$

Writing the system in a more compact form

(21) $\begin{equation*}\begin{array}{ccl}\left(\begin{array}{c}\dot x\\ \dot \lambda_1\end{array}\right)&=&F_{11}\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)+F_{12}\left(\begin{array}{c}z\\ \lambda_2\end{array}\right)\\\epsilon \left(\begin{array}{c}\dot z\\ \dot \lambda_2\end{array}\right)&=&F_{21}\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)+F_{22}\left(\begin{array}{c}z\\ \lambda_2\end{array}\right)\end{array}\end{equation*}$

we see that this system is in the familiar singularly perturbated form. To define its reduced system, we need $F_{22}^{-1}$ . However, anticipating that $F_{22}$ will be the Hamiltonian matrix for a fast subsystem, we
make the following assumption.

Assumption 2.4: The eigenvalues of $F_{22}(t)$ lie off the imaginary axis for all $t\in [0,1]$ , i.e. there exists a
constant $c$ such that

$\vert {\rm Re~}\sigma (F_{22})\vert \ge c>0,\quad t\in [0,1].$

Under this assumption $F_{22}^{-1}$ exists and the reduced necessary condition is

(22) $\begin{equation*}\begin{array}{ccl}\left(\begin{array}{c}\dot x\\ \dot \lambda_1\end{array}\right)&=&(F_{11}-F_{12}F_{22}^{-1}F_{21})\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)\\\left(\begin{array}{c}\dot z\\ \dot \lambda_2\end{array}\right)&=&-F_{22}^{-1}F_{21}\left(\begin{array}{c}x\\ \lambda_1\end{array}\right)\end{array}\end{equation*}$

subject to $x(0)=x^0$ . From lemma 2.1, we know that the same necessary condition is obtained by setting $\epsilon =0$ in (18), even if $A_{22}^{-1}$ does not exist.

Assumption 2.5: The solution of the reduced problem exists and is unique.

Theorem 2.2: Consider the time-invariant optimal control problem (18) where $R>0$ . If assumptions 2.4 and 2.5 hold and if

${\rm rank}[B_2,A_{22}B_2,…,A_{22}^{m-1}B_2]=m,$

then there exists an $\epsilon_0>0$ such that for all $\epsilon \in (0,\epsilon_0]$ and all $t\in [0,1]$ the optimal trajectory $x(t,\epsilon ),~z(t,\epsilon )$ and the corresponding optimal control $u(t,\epsilon )$ satisfy

(23) $\begin{equation*}\begin{array}{ccl}x(t,\epsilon )&=&x_s(t)+O(\epsilon )\\z(t,\epsilon )&=&z_s(t)+\eta_l(\frac t{\epsilon })+\eta_r(\frac {1-t}{\epsilon })+O(\epsilon )\\u(t,\epsilon )&=&u_s(t)+\nu_l(\frac t{\epsilon })+\nu_r(\frac {1-t}{\epsilon })+O(\epsilon ),\end{array}\end{equation*}$

where $x_s(t)$ , $z_s(t)$ , $u_s(t)$ is the optimal solution of the reduced problem (22), while $\eta_l$ , $\nu_l$ , $\eta_r$ , $\nu_r$ , are solutions to boundary layer systems. This theorem is simple to apply. First, solve the reduced problem (22). Next, we determine the initial and terminal conditions $\eta_l(0)$ and $\eta_r(0)$ from (23) at $t=0$ and $t=1$ as

$\eta_l(0)=z_0-z_s(0),~~\eta_r(0)=z_1-z_s(1).$

These values are the correctors of the discrepances in the boundary conditions for $z$ caused by the reduced solution $z_s(t)$ .

Consider the cheap control problem

(24) $\begin{equation*}\begin{array}{ccl}&&\min \frac 1{2}\int_0^T(\left(\begin{array}{cc}x^T& z^T\end{array}\right)Q\left(\begin{array}{c}x\\ z\end{array}\right)+\epsilon^2u^TRu)dt,\\\dot x&=&A_{11}x+A_{12}z,\quod x(0)=x^0,~x(T)=x^1,\\\dot z&=&A_{21}x+A_{22}z+u,\quod z(0)=z^0,~z(T)=z^1,\end{array}\end{equation*}$

where $\epsilon$ is small. There is no loss of generality in considering the system in (24), because every linear system with a full-rank control input matrix can be represented in this form after a similarity transformation. Although the system in (24) is not singularly perturbated, a large $u$ will force the $z$ -variable to act as a fast variable. The Hamiltonian function for this problem is

(25) $\begin{equation*}\begin{array}{ccl}H&=&\frac 1{2}\left(\begin{array}{cc}x^T& z^T\end{array}\right)Q\left(\begin{array}{c}x\\ z\end{array}\right)+\frac 1{2}\epsilon^2u^TRu+\\&&\left(\begin{array}{cc}p^T& q^T\end{array}\right)\left(A\left(\begin{array}{c}x\\ z\end{array}\right)+\left(\begin{array}{c}0\\ I\end{array}\right)u\right).\end{array}\end{equation*}$

For $\epsilon>0$ and $R>0$ , the optimal control is

(26) $\begin{equation*}u=-\frac 1{\epsilon^2}R^{-1}q.\end{equation*}$

Now we rescale the adjoint variable $q$ as $\bar q=q/\epsilon$ , transferring the singularity at $\epsilon =0$ from the functional to the dynamics, and write the adjoint equations in the form

(27) $\begin{equation*}\begin{array}{ccl}\dot p&=&-\nabla_x H,\\\epsilon \dot {\bar q}&=&-\nabla_z H.\end{array}\end{equation*}$

Combining (24), (26) and (27), we obtain the singularly perturbated boundary value problem

(28) $\begin{equation*}\begin{array}{ccl}\left(\begin{array}{c}\dot x\\ \dot p\end{array}\right)&=&\left(\begin{array}{cc}A_{11}& 0\\ -Q_{11}& -A_{11}^T\end{array}\right)\left(\begin{array}{c}x\\ p\end{array}\right)+\\&&+\left(\begin{array}{cc}A_{12}& 0\\ -Q_{12}& -\epsilon A_{21}^T\end{array}\right)\left(\begin{array}{c}z\\ \bar q\end{array}\right),\\x(0)&=&x^0,~x(T)=x^1,\\\epsilon \left(\begin{array}{c}\dot z\\ \dot {\bar q}\end{array}\right)&=&\left(\begin{array}\epsilon A_{21}& 0\\ -Q_{12}^T& -A_{12}^T\end{array}\right)\left(\begin{array}{c}x\\ p\end{array}\right)+\left(\begin{array}\epsilon A_{22}& -R^{-1}\\ -Q_{22}&-\epsilon A_{22}^T\end{array}\right)\left(\begin{array}{c}z\\ \bar q\end{array}\right),\\z(0)&=&z^0,~z(T)=z^1,\end{array}\end{equation*}$

which is a special case of the problem (20). If $Q_{22}>0$ , assumption 2.4 and 2.5 are satisfied, and we can set $\epsilon =0$ in (28) to obtain

(29) $\begin{equation*}\begin{array}{ccl}z&=&-Q_{22}^{-1}(Q_{12}^Tx+A_{12}^Tp)\\q&=&0.\end{array}\end{equation*}$

Substituting (29) into (28) gives the reduced slow system

$\begin{equation*}\begin{array}{ccl}\left(\begin{array}{c}\dot x\\ \dot p\end{array}\right)&=&\left(\begin{array}{cc}A_{11}-A_{12}Q_{22}^{-1}Q_{12}^T&-A_{12}Q_{22}^{-1}A_{12}^T\\ -Q_{11}+Q_{12}Q_{22}^{-1}Q_{12}^T&-A_{11}^T+Q_{12}Q_{22}^{-1}A_{12}^T\end{array}\right)\left(\begin{array}{c}x\\ p\end{array}\right),\\x(0)&=&x^0,~x(T)=x^1.\end{array}\end{equation*}$

A cheap control problem where $T=\infty$ , results in a problem which is easier to solve because the corresponding Riccati equation is algebraic, but we must require the existence of a unique positive semidefinite solution to that Riccati equation, which stabilizes the system.

Thus, consider the functional to be minimized

(30) $\begin{equation*}J=\frac 1{2}\int_0^{\infty}(\left(\begin{array}{cc}x^T& z^T\end{array}\right)\left(\begin{array}{cc}Q_{11}& Q_{12}\\ Q_{12}^T& Q_{22}\end{array}\right)\left(\begin{array}{c}x\\ z\end{array}\right)+\epsilon^2u^TRu)dt,\end{equation*}$

subject to the dynamics

(31) $\begin{equation*}\begin{array}{ccl}\dot x&=&A_{11}x+A_{12}z,\quad x(0)=x^0,~x(T)=x^1,\\\dot z&=&A_{21}x+A_{22}z+B_2u,\quad z(0)=z^0,~z(T)=z^1,\end{array}\end{equation*}$

where $B_2$ is nonsingular, $Q\ge 0$ , $R>0$ and $\epsilon$ is small. To regulate the output

$y=C_1x+C_2z$

we set $Q=C^TC$ . We will make the following assumptions.

Assumption 2.6: $Q_{22}$ is a positive-definite matrix.

Under this assumption, a two-time-scale decomposition of the problem
is possible.

Assumption 2.7: The algebraic Riccati equation corresponding to the problem (30)-(31) has a unique positive semidefinite stabilizing solution.

Under assumptions 2.6 and 2.7 we can apply theorem 2.2 to the problem (30)-(31) and approximate $x$ , $z$ and $u$ by the expressions in (23), respectively.