Recover Continuous time Linear System From Discrete time

CONTROLLABILITY, OBSERVABILITY, AND DISTANCE TO UNCONTROLLABILITY

BISWA NATH DATTA , in Numerical Methods for Linear Control Systems, 2004

6.3 OBSERVABILITY: DEFINITIONS AND BASIC RESULTS

In this section we state definitions and basic results of observability. The results will not be proved here because they can be easily proved by duality of the results on controllability proved in the previous section.

6.3.1 Observability of a Continuous-Time System

The concept of observability is dual to the concept of controllability.

Definition 6.3.1.

The continuous-time system (6.2.1) is said to be observable if there exists t ₁ > 0 such that the initial state x(0) can be uniquely determined from the knowledge of u(t) and y(t) for all t, 0 ⩽ t ⩽ t ₁.

Remark

•

The observability of the system (6.2.1) is often referred to as the observability of the pair (A, C).

Analogous to the case of controllability, we state the following criteria of observability.

Theorem 6.3.1. Criteria for Continuous-Time Observability.

The following are equivalent:

(i)

The system (6.2.1) is observable.

(ii)

The observability matrix

$O_{\text{M}}=\left(\begin{array}{l}C\hfill \\ CA\hfill \\ C{A}^2\hfill \\ ⋮\hfill \\ C{A}^{n-1}\hfill \end{array}\right)$

has full rank n.

(iii)

The matrix

$W_{\text{O}}={\int }_0^{t_1}{e}^{A^T\tau }{C}^{\text{T}}C{e}^{{\text{A}}^{\tau }}d\tau$

is nonsingular for any t ₁ > 0.

(iv)

The matrix

$\left[\begin{array}{c}\lambda I-A\\ C\end{array}\right]$

has rank n for every eigenvalue λ of A.

(v)

None of the eigenvectors of A is orthogonal to the rows of C, that is, if (λ, y) is an eigenpair of A, then Cy ≠ 0.

(vi)

There exists a matrix L such that the eigenvalues of A + LC can be assigned arbitrarily, provided that the complex eigenvalues occur in conjugate pairs.

We only prove (iii) ⇔ (i) and leave the others as an exercise ( Exercise 6.6 ).

Theorem 6.3.2.

The pair (A, C) is observable if and only if the matrix W _O is nonsingular for any t ₁ > 0.

Proof. First suppose that the matrix W _O is nonsingular. Since y(t) and u(t) are known, we can assume, without any loss of generality, that u(t) = 0 for every t. Thus,

This gives

Thus, x(0) is uniquely determined and is given by $x(0)=W_{\text{O}}^{-1}{\int }_0^{t_1}{e}^{A^{\text{T}}\tau }{C}^{\text{T}}y(\tau )d\tau$ .

Conversely, if W _O is singular, then there exists a nonzero vector z such that W _O z = 0, which in turn implies that Ce ^Aτ z = 0. So, y(τ) = Ce ^Aτ(x(0) + z) = Ce ^Aτ x(0).

Thus, x(0) cannot be determined uniquely, implying that (A, C) is not observable.

Component observability. As in the case of controllability, we can also speak of component observability when all the states are not observable by certain output. The rank of the observability matrix

where c _j , the jth row of the output matrix C, determines the number of states that are observable by the output y _j (t).

6.3.2 Observability of a Discrete-Time System

Definition 6.3.2.

The discrete-time system (6.2.6) is said to be observable if there exists an index N such that the initial state x ₀ can be completely determined from the knowledge of inputs u ₀, u ₁,…, u _N-1, and the outputs y ₀, y ₁,…, y _N .

The criteria of observability in the discrete-time case are the same as in the continuous-time case, and therefore, will not be repeated here.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780122035906500100

REVIEWS

Wing-Kuen Ling , in Nonlinear Digital Filters, 2007

Controllability and observability

For a continuous time system, assume that x(0) = 0. ∀x ₁, if ∃t ₁ > 0 and u(t) such that x(t ₁) = x ₁, then the continuous time system is said to be reachable. Similarly, for a discrete time system, assume that x(0) = 0. ∀x ₁, if ∃n ₁ > 0 and u(n) such that x(n ₁) = x ₁, then the discrete time system is said to be reachable. For a continuous time system, if ∀x ₀, x ₁, ∃t ₁ > 0 and u(t) such that x(0) = x ₀ and x(t ₁) = x ₁, then the continuous time system is said to be controllable. Similarly, for a discrete time system, if ∀x ₀, x ₁, ∃n ₁ > 0 and u(n) such that x(0) = x ₀ and x(n ₁ = x ₁, then the discrete time system is said to be controllable. For LTI systems, the set of reachable state is R(|BAB… AⁿB|), where R(A) is defined as the range of A, that is R(A) = {y : y = A x}. Also, the LTI systems are controllable if and only if R(A) = R ⁿ Or in other words, rank(|BAB… Aⁿ _B |) = n.

For a continuous time system, ∀x ₁, if ∃t ₁ > 0 and u(t) such that x ₁ can be determined from y(t) for t > t ₁, then the continuous time system is said to be observable. Similarly, for a discrete time system, ∀x ₁, if ∃n ₁ > 0 and u(n) such that x ₁ can be determined from y(n) for n > n ₁, then the discrete time system is said to be observable. For LTI systems, the set of unobservable state is

$\text{N}\left(\left[\begin{array}{l}\text{C}\hfill \\ \text{CA}\hfill \\ ⋮\hfill \\ {\text{CA}}^{n-1}\hfill \end{array}\right]\right),$

where N(A) is defined as the null space of the kernel of A, that is N(A) ≡ {x: Ax = 0}. Also, the LTI systems are observable if and only if

$\text{N}\left(\left[\begin{array}{l}\text{C}\hfill \\ \text{CA}\hfill \\ ⋮\hfill \\ {\text{CA}}^{n-1}\hfill \end{array}\right]\right)=\left\{0\right\}.\text{Or\hspace{0.5em}in\hspace{0.5em}other\hspace{0.5em}words},rank\left(\left[\begin{array}{l}\text{C}\hfill \\ \text{CA}\hfill \\ ⋮\hfill \\ {\text{CA}}^{n-1}\hfill \end{array}\right]\right)=n.$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780123725363500028

STABILITY, INERTIA, AND ROBUST STABILITY

BISWA NATH DATTA , in Numerical Methods for Linear Control Systems, 2004

7.10 SUMMARY AND REVIEW

The stability of the system:

or that of

is essentially governed by the eigenvalues of the matrix A.

Mathematical Criteria of Stability

The continuous-time system $\dot{x}(t)=Ax(t)$ is asymptotically stable if and only if the eigenvalues of A are all in the left half plane ( Theorem 7.2.1 ). Similarly, the discrete-time system x(k + 1) = Ax(k) is asymptotically stable if and only if all the eigenvalues of A are inside the unit circle. ( Theorem 7.3.1 ). Various Lyapunov stability theorems (Theorems 7.2.3–7.2.9 , and Theorem 7.3.2 ) have been stated and proved.

The Inertia of a Matrix

Two important inertia theorems ( Theorems 7.4.1 and 7.4.2 ) and the classical Sylvester Law of Inertia have been stated. These inertia theorems generalize the Lyapunov stability results.

Methods for Determining Stability and Inertia

The Characteristic Polynomial Approach and the Matrix Equation Approach are two classical approaches for determining the stability of a system and the inertia of a matrix. Both these approaches have some computational drawbacks.

The zeros of a polynomial may be extremely sensitive to small perturbations. Furthermore, the numerical methods to compute the characteristic polynomial of a matrix are usually unstable.

The most numerically effective method (The Schur method, described in Chapter 8), for solving a Lyapunov matrix equation is based on reduction of the matrix A to RSF, and the RSF displays the eigenvalues of A or the eigenvalues can be trivially computed out of this form.

Thus, the characteristic equation approach is not numerically viable and the matrix equation approach for stability and inertia is counterproductive.

Hence, the most numerically effective approach for stability and inertia is the eigenvalue approach: compute all the eigenvalues of A.

By explicitly computing the eigenvalues, one, however, gets much more than what is needed for stability and inertia. Furthermore, since the eigenvalues of a matrix can be very sensitive to small perturbations, determining the inertia and stability by computing explicitly the eigenvalues can be misleading.

An implicit matrix equation approach ( Algorithm 7.5.1 ), which does not require computation of eigenvalues nor explicit solution of any matrix equation has been described. Algorithm 7.5.1 is about three times faster than the eigenvalue method (According to the flop-count).

Distance to an Unstable System

Given a stable matrix A, the quantity β(A) defined by

$\beta (A)=min\left\{\right|\left|E\right|{|}_F\text{such}\text{that}A+E\in U\},$

where U is the set of n × n matrices with at least one eigenvalue on the imaginary axis, is the distance to the set of unstable matrices.

A bisection algorithm ( Algorithm 7.6.1 ) based on knowing if a certain Hamiltonian matrix (the matrix (7.6.3)) has a purely imaginary eigenvalue, is described. The algorithm is based on Theorem 7.6.1 , which displays a relationship between a spectral property of the Hamiltonian matrix and the quantity β(A).

The discrete-analog of β(A) is defined to be

$\gamma (A)=min\left\{\right|\left|E\right||\text{for}\text{some}\theta \in ℝ;e^{i\theta }\in \Omega (A+E)\}.$

An analog of Theorem 7.6.1 ( Theorem 7.6.3 ) is stated and a bisection algorithm ( Algorithm 7.6.2 ) based on this theorem is described.

Robust Stability

Given a stable matrix A, one naturally wonders if the matrix A + E remains stable, where E is a certain perturbed matrix. Two bounds for E guaranteeing the stability of the perturbed matrix (A + E) are given, in terms of solutions of certain Lyapunov equations ( Theorems 7.7.1 and 7.7.2 ).

Stability Radius

Section 7.8 deals with the structured stability radius. If the perturbations are of the form BΔC, where Δ is an unknown perturbation matrix, then it is of interest to know the size of smallest Δ (measured using 2-norm) that will destabilize the perturbed matrix A + BΔC. In this context, the concept of stability radius is introduced, and formulas both for the complex stability radius ( Theorem 7.8.1 ) and the real stability radius are stated.

H₂-Norm

The H ₂-norm of a stable transfer, transfer function measures the steady-state covariance of the output response y = Gv to the white noise inputs v. An algorithm ( Algorithm 7.2.1 ) for computing the H ₂-norm, based on computing the controllability or observability Grammian via Lyapunov equations is given.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780122035906500112

Continuous-Time Systems

Luis F. Chaparro , Aydin Akan , in Signals and Systems Using MATLAB (Third Edition), 2019

2.3.2 Time-Invariance

A continuous-time system $\mathcal{S}$ is time-invariant if whenever for an input $x(t)$ , with a corresponding output $y(t)=\mathcal{S}[x(t)]$ , the output corresponding to a shifted input $x(t\mp \tau )$ (delayed or advanced) is the original output equally shifted in time, $y(t\mp \tau )=\mathcal{S}[x(t\mp \tau )]$ (delayed or advanced). Thus

(2.9) $\begin{array}{ccc}\hfill x(t)& \hfill \Rightarrow \hfill & y(t)=\mathcal{S}[x(t)],\hfill \\ \hfill x(t\mp \tau )& \hfill \Rightarrow \hfill & y(t\mp \tau )=\mathcal{S}[x(t\pm \tau )].\hfill \end{array}$

That is, the system does not age—its parameters are constant.

A system that satisfies both the linearity and the time-invariance is called Linear Time-Invariant or LTI.

Remarks

1.

It should be clear that linearity and time-invariance are independent of each other. Thus, one can have linear time-varying, or non-linear time-invariant systems.

2.

Although most actual systems are, according to the above definitions, non-linear and time-varying, linear models are used to approximate (around an operating point) the non-linear behavior, and time-invariant models are used to approximate (in short segments) the system's time-varying behavior. For instance, in speech synthesis the vocal system is typically modeled as a linear time-invariant system for intervals of about 20 ms attempting to approximate the continuous variation in shape in the different parts of the vocal system (mouth, cheeks, nose, etc.). A better model for such a system is clearly a linear time-varying model.

3.

In many cases time invariance can be determined by identifying—if possible—the input, the output and letting the rest represent the parameters of the system. If these parameters change with time, the system is time-varying. For instance, if the input $x(t)$ and the output $y(t)$ of a system are related by the equation $y(t)=f(t)x(t)$ , the parameter of the system is the function $f(t)$ , and if it is not constant, the system is time-varying. Thus, the system $y(t)=Ax(t)$ , where A is a constant, is time invariant as can be easily verified. But the AM modulation system given by $y(t)=\cos ({\mathrm{\Omega }}_{0}t)x(t)$ is time-varying as the function $f(t)=\cos ({\mathrm{\Omega }}_{o}t)$ is not constant.

The Beginnings of Radio

The names of Nikola Tesla (1856–1943) and Reginald Fessenden (1866–1932) are linked to the invention of amplitude modulation and radio [65,78,7]. Radio was initially called "wireless telegraphy" and then "wireless" [26].

Tesla was a mechanical as well as an electrical engineer, but mostly an inventor. He has been credited with significant contributions to electricity and magnetism in the late 19th and early 20th century. His work is the basis of the alternating current (AC) power system, and the induction motor. His work on wireless communications using the "Tesla coils" was capable of transmitting and receiving radio signals. Although Tesla submitted a patent application for his basic radio before Guglielmo Marconi, it was Marconi who was initially given the patent for the invention of radio (1904). The Supreme Court in 1943 reversed the decision in favor of Tesla [6].

Fessenden has been called the "father of radio broadcasting." His early work on radio led to demonstrations on December 1906, of the capability of point-to-point wireless telephony, and what appears to be the first audio radio broadcasts of entertainment and music ever made to an audience (in this case shipboard radio operators in the Atlantic). Fessenden was a professor of electrical engineering at Purdue University and the first chairman of the electrical engineering department of the University of Pittsburgh in 1893.

AM Communication System

Amplitude modulation (AM) communication systems arose from the need to send an acoustic signal, a "message", over the airwaves using a reasonably sized antenna to radiate it. The size of the antenna depends inversely on the highest frequency present in the message, and voice and music have relatively low frequencies. A voice signal typically has frequencies in the range of 100 Hz to about 5 kHz (the frequencies needed to make a telephone conversation intelligible) while music typically displays frequencies up to about 22 kHz. The transmission of such signals with a practical antenna is impossible. To make the transmission possible, modulation was introduced, i.e., multiplying the message $m(t)$ by a periodic signal such as a cosine $\cos ({\mathrm{\Omega }}_{0}t)$ , the carrier, with a frequency ${\mathrm{\Omega }}_0$ much larger than those in the acoustic signal. Amplitude modulation provided the larger frequencies needed to reduce the size of the antenna. Thus $y(t)=m(t)\cos ({\mathrm{\Omega }}_{0}t)$ is the signal to be transmitted, and the effect of this multiplication is to change the frequency content of the input.

The AM system is clearly linear, but time-varying. Indeed, if the input is $m(t-\tau )$ , the message delayed τ seconds, the output would be $m(t-\tau )\cos ({\mathrm{\Omega }}_{0}t)$ which is not $y(t-\tau )=\cos ({\mathrm{\Omega }}_0(t-\tau ))m(t-\tau )$ as a time-invariant system would give. Fig. 2.4 illustrates the AM transmitter and receiver. The carrier continuously changes independent of the input and as such the system is time-varying.

Figure 2.4. AM modulation: transmitter and receiver.

FM Communication System

In comparison with an AM system, a frequency modulation (FM) system is non-linear and time-varying. An FM modulated signal $z(t)$ is given by

$z(t)=\cos ({\mathrm{\Omega }}_{c}t+{\int }_{-\infty }^{t}m(\tau )d\tau ),$

where $m(t)$ is the input message.

To show the FM system is non-linear, assume we scale the message to $\gamma m(t)$ , for some constant γ, the corresponding output is given by

$\cos ({\mathrm{\Omega }}_{c}t+\gamma {\int }_{-\infty }^{t}m(\tau )d\tau )$

which is not the previous output scaled, i.e., $\gamma z(t)$ , thus FM is a non-linear system. Likewise, if the message is delayed or advanced, the output will not be equally delayed or advanced—thus the FM system is not time-invariant.

Vocal System

A remarkable system that we all have is the vocal system (see Fig. 2.5). The air pushed out from the lungs in this system is directed by the trachea through the vocal cords making them vibrate and create resonances similar to those from a musical wind instrument. The generated sounds are then muffled by the mouth and the nasal cavities resulting in an acoustic signal carrying a message. Given the length of the typical vocal system, on average for adult males about 17 cm and 14 cm for adult females, it is modeled as a distributed system and represented by partial differential equations. Due to the complexity of this model, it is the speech signal along with the understanding of the speech production that is used to obtain models of the vocal system. Speech processing is one of the most fascinating areas of electrical engineering.

Figure 2.5. Vocal system: principal organs of articulation. Model for speech production.

A typical linear time-invariant model for speech production considers segments of speech of about 20 ms, and for each a low-order LTI system is developed. The input is either a periodic pulse for the generation of voiced sounds (e.g., vowels) or a noise-like signal for unvoiced sounds (e.g., the /sh/ sound). Processing these inputs gives speech-like signals. A linear time-varying model would take into consideration the variations of the vocal system with time and it would be thus more appropriate.

Example 2.4

Consider constant linear capacitors and inductors, represented by ordinary differential equations

$\frac{d{v}_c(t)}{dt}=\frac{1}{C}i(t),\frac{d{i}_L(t)}{dt}=\frac{1}{L}v(t),$

with initial conditions $v_c(0)=0$ and $i_L(0)=0$ , under what conditions are these time-invariant systems?

Solution: Given the duality of the capacitor and the inductor, we only need to consider one of these. Solving the ordinary differential equation for the capacitor we obtain (the initial voltage $v_c(0)=0$ ):

$v_c(t)=\frac{1}{C}{\int }_0^{t}i(\tau )d\tau .$

Suppose then we delay the input current $i(t)$ by λ s. The corresponding output is given by

(2.10) $\frac{1}{C}{\int }_0^{t}i(\tau -\lambda )d\tau =\frac{1}{C}{\int }_{-\lambda }^{0}i(\rho )d\rho +\frac{1}{C}{\int }_0^{t-\lambda }i(\rho )d\rho$

by changing the integration variable to $\rho =\tau -\lambda$ . For the above equation to equal the voltage across the capacitor delayed λ s, given by

$v_c(t-\lambda )=\frac{1}{C}{\int }_0^{t-\lambda }i(\rho )d\rho ,$

we need that $i(t)=0$ for $t<0$ , so that the first integral in the right expression in Equation (2.10) is zero. Thus, the system is time-invariant if the input current $i(t)=0$ for $t<0$ . Putting this together with the condition on the initial voltage $v_c(0)=0$ , we can say that for the capacitor to be a linear and time-invariant system it should not be initially energized—i.e., the input $i(t)=0$ for $t<0$ and no initial voltage across the capacitor, $v(0)=0$ . Similarly, using duality, for the inductor to be linear and time-invariant, the input voltage across the inductor $v(t)=0$ for $t<0$ (to guarantee time-invariance) and the initial current in the inductor $i_L(0)=0$ (to guarantee linearity). □

Example 2.5

Consider the RLC circuit in Fig. 2.6 consisting of a series connection of a resistor R, an inductor L and a capacitor C. The switch has been open for a very long time and it is closed at $t=0$ , so that there is no initial energy stored in either the inductor or the capacitor and the voltage applied to the elements is zero for $t<0$ . Obtain an equation connecting the input $v(t)$ and output the current $i(t)$ .

Figure 2.6. RLC circuit with a switch for setting the initial conditions to zero.

Solution: Because of the presence of the capacitor and the inductor, both capable of storing energy, this circuit is represented by a second order ordinary differential equation with constant coefficients. According to Kirchhoff's voltage law:

$v(t)=Ri(t)+\frac{1}{C}{\int }_0^{t}i(\tau )d\tau +L\frac{di(t)}{dt}.$

To get rid of the integral we find the derivative of $v(t)$ with respect to t:

$\frac{dv(t)}{dt}=R\frac{di(t)}{dt}+\frac{1}{C}i(t)+L\frac{d^{2}i(t)}{d{t}^2},$

a second order ordinary differential equation, with input the voltage source $v(t)$ , and output the current $i(t)$ . Because the circuit is not energized for $t<0$ , the circuit represented by the above ordinary differential equation is linear and time-invariant. □

Remark

An RLC circuit is represented by an ordinary differential equation of order equal to the number of independent inductors and capacitors in the circuit. If two or more capacitors (two or more inductors) are connected in parallel (in series), sharing the same initial voltage across them (same initial current) we can convert them into an equivalent capacitor with the same initial voltage (equivalent inductor with the same current).

A system represented by a linear ordinary differential equation, of any order N, having constant coefficients, and with input $x(t)$ and output $y(t)$ :

(2.11) $a_{0}y(t)+a_1\frac{dy(t)}{dt}+\cdots +\frac{d^{N}y(t)}{d{t}^N}=b_{0}x(t)+b_1\frac{dx(t)}{dt}+\cdots +b_M\frac{d^{M}x(t)}{d{t}^M}t>0$

is linear time-invariant if the system is not initially energized (i.e., the initial conditions are zero, and the input $x(t)$ is zero for $t<0$ ). If one or more of the coefficients $\{a_i,b_i\}$ are functions of time, the system is time-varying.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128142042000120

Techniques in Discrete and Continuous Robust Systems

Ahmad A. Mohammad , J.A. De Abreu-Garcia , in Control and Dynamic Systems, 1996

3 Theoretical Aspects of the LBCT

This is best summarized in the following theorem:

Theorem: For any given stable, minimal balanced DTS (A_db, B_db, C_db, D_d, Σ, T), with distinct second order modes, there exists a unique, minimal, stable, and balanced CTS model (A_cb, B_cb, C_cb, D_c  =   D_d, Σ), where

${\mathrm{B}}_{\mathrm{c}\mathrm{b}}={\mathrm{B}}_{\mathrm{d}\mathrm{b}}/\sqrt{\mathrm{T}},{\mathrm{C}}_{\mathrm{c}\mathrm{b}}={\mathrm{C}}_{\mathrm{d}\mathrm{b}}/\sqrt{\mathrm{T}}\text{,}$

and A_cb is the unique solution of the CT Lyapunov equations

$\begin{array}{l}{\mathrm{A}}_{\mathrm{c}\mathrm{b}}\sum +\sum {\mathrm{A}}_{\mathrm{c}\mathrm{b}}^{\mathrm{T}}=-{\mathrm{B}}_{\mathrm{c}\mathrm{b}}{\mathrm{B}}_{\mathrm{c}\mathrm{b}}^{\mathrm{T}}=-{\mathrm{T}\mathrm{B}}_{\mathrm{d}\mathrm{b}}{\mathrm{B}}_{\mathrm{d}\mathrm{b}}^{\mathrm{T}}\\ {\mathrm{A}}_{{}_{\mathrm{c}\mathrm{b}}}^{\mathrm{T}}\sum +\sum {\mathrm{A}}_{\mathrm{c}\mathrm{b}}=-{\mathrm{C}}_{{}_{\mathrm{c}\mathrm{b}}}^{\mathrm{T}}{\mathrm{C}}_{\mathrm{c}\mathrm{b}}=-{\mathrm{T}\mathrm{C}}_{{}_{\mathrm{d}\mathrm{b}}}^{\mathrm{T}}{\mathrm{C}}_{\mathrm{d}\mathrm{b}}\text{.}\end{array}$

Moreover, the CTS approximation of the DTS model has the following properties:

1)

The CTS model is minimal and stable if the DTS model is minimal and stable.

2)

There exists a solution for Eqs.(84) – (85), for a fixed choice of T, B, and C.

3)

The solution of 2) is unique.

4)

The |H|₂ norm of the DTS is equal to the |H|₂ norm of the CTS (notice that we only deal with the strictly proper part).

5)

The initial error in the step response is identically zero.

Proof:

(1)

Guaranteed by Lyapunov equations.

(2, 3)

Proved in section 2.

(4)

A direct consequence of matching the Hankel Singular values of the DT and CT systems.

(5)

First, it is noted that both systems share the same feedthrough term D. Moreover, from the initial value theorem, the initial response of the strictly proper part is identically zero.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S0090526796800673

State–space representation

M. Sami Fadali , Antonio Visioli , in Digital Control Engineering (Third Edition), 2020

7.9 Similarity transformation

Any given linear system has an infinite number of valid state–space representations. Each representation corresponds to a different set of basis vectors in state–space. Some representations are preferred because they reveal certain system properties, whereas others may be convenient for specific design tasks. This section considers transformation from one representation to another.

Given a state vector x(k) with state–space representation of the form of Eq. (7.59), we define a new state vector z(k)

(7.99) $\mathbf{\text{x}}\left(k\right)=T_r\mathbf{\text{z}}\left(k\right)\iff \mathbf{\text{z}}\left(k\right)=T_r^{-1}\mathbf{\text{x}}\left(k\right)$

where the transformation matrix T _r is assumed invertible. Substituting for x(k) from Eq. (7.99) in the state Eq. (7.59) and the output Eq. (7.67) gives

(7.100) $\begin{array}{c}{T}_r\mathbf{\text{z}}\left(k+1\right)=A_d{T}_r\mathbf{\text{z}}\left(k\right)+B_d\mathbf{\text{u}}\left(k\right)\\ \mathbf{\text{y}}\left(k\right)=C{T}_r\mathbf{\text{z}}\left(k\right)+D\mathbf{\text{u}}\left(k\right)\end{array}$

Premultiplying the state equation by $T_r^{-1}$ gives

(7.101) $\mathbf{\text{z}}\left(k+1\right)=T_r^{-1}{A}_d{T}_r\text{z}\left(k\right)+T_r^{-1}{B}_d\mathbf{\text{u}}\left(k\right)$

Hence, we have the state–space quadruple for the state vector z(k)

(7.102) $\left(\text{A},\text{B},\text{C},\text{D}\right)=\left(T_r^{-1}{A}_d{T}_r,T_r^{-1}{B}_d,C{T}_r,D\right)$

Clearly, the quadruple for the state vector x(k) can be obtained from the quadruple of z(k) using the inverse of the transformation T_r (i.e., the matrix $T_r^{-1}$ ).

For the continuous-time system of Eq. (7.3), a state vector z(t) can be defined as

(7.103) $\mathbf{\text{x}}\left(t\right)=T_r\mathbf{\text{z}}\left(t\right)\iff \mathbf{\text{z}}\left(t\right)=T_r^{-1}\mathbf{\text{x}}\left(t\right)$

and substitution in Eq. (7.3) yields Eq. (7.102). Thus, a discussion of similarity transformation is identical for continuous-time and discrete-time systems. We therefore drop the subscript d in the sequel.

Example 7.20

Consider the point mass m driven by a force f of Example 7.1, and determine the two state–space equations when m  =   1 and

1.

the displacement and the velocity are the state variables and

2.

the displacement and the sum of displacement and velocity are the state variables.

Solution

1. From the first principle, as shown in Example 7.1, we have that the equation of motion is

$m\ddot{y}=f$

and in case 1, by considering x ₁ as the displacement and x ₂ as the velocity, we can write

$\begin{array}{c}{\dot{x}}_1\left(t\right)=x_2\left(t\right)\\ {\dot{x}}_2\left(t\right)=u\left(t\right)\\ y\left(t\right)=x_1\left(t\right)\end{array}$

Thus, we obtain the realization

$A\left[\begin{array}{cc}1& 0\\ 0& 0\end{array}\right]B=\left[\begin{array}{c}0\\ 1\end{array}\right]C=\left[\begin{array}{cc}1& 0\end{array}\right]D=0$

2.

We express the new state variables, z ₁ and z ₂, in terms of the state variable of part 1 as

$\begin{array}{l}{z}_1\left(t\right)=x_1\left(t\right)\\ z_2\left(t\right)=x_1\left(t\right)+x_2\left(t\right)\end{array}$

This yields the inverse of the transformation matrix

$T_r^{-1}=\left[\begin{array}{cc}1& 0\\ 1& 1\end{array}\right]$

Using Eq. (7.102) gives the realization

$\begin{array}{l}\text{A}=T_r^{-1}{A}_d{T}_r=\left[\begin{array}{cc}-1& 1\\ -1& 1\end{array}\right]\text{B}=T_r^{-1}{B}_d=\left[\begin{array}{c}0\\ 1\end{array}\right]\\ \text{C}=\left[\begin{array}{cc}1& 0\end{array}\right]\text{D}=0\end{array}$

To obtain the transformation to diagonal form, we recall the expression

(7.104) $A=V\mathit{\Lambda }{V}^{-1}\iff \mathit{\Lambda }=V^{-1}AV$

where V is the modal matrix of eigenvectors of A and Λ  =   diag{λ ₁, λ ₂, …, λ _n } is the matrix of eigenvalues of A. Thus, for A   = Λ in Eq. (7.102), we use the modal matrix of A as the transformation matrix. The form thus obtained is not necessarily the same as the diagonal form obtained in Section 8.5.3 from the transfer function using partial fraction expansion. Even though all diagonal forms share the same state matrix, their input and output matrices may be different.

Example 7.21

Obtain the diagonal form for the state–space equations

$\left[\begin{array}{c}{x}_1\left(k+1\right)\\ x_2\left(k+1\right)\\ x_3\left(k+1\right)\end{array}\right]=\left[\begin{array}{ccc}0& 1& 0\\ 0& 0& 1\\ 0& -0.04& -0.5\end{array}\right]\left[\begin{array}{c}{x}_1\left(k\right)\\ x_2\left(k\right)\\ x_3\left(k\right)\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]u\left(k\right)$

$y\left(k\right)=[\begin{array}{ccc}1& 0& 0\end{array}]\left[\begin{array}{c}{x}_1\left(k\right)\\ x_2\left(k\right)\\ x_3\left(k\right)\end{array}\right]$

Solution

The eig command of MATLAB yields the eigenvalues and the modal matrix whose columns have unity norm

$\mathit{\Lambda }=\text{diag}\left\{\begin{array}{ccc}0,& -0.1,& -0.4\end{array}\right\}V=\left[\begin{array}{ccc}1& -0.995& 0.9184\\ 0& 0.0995& -0.36741\\ 0& -0.00995& 0.1469\end{array}\right]$

The state matrix is in companion form and the modal matrix is also known to be the Van der Monde matrix (column norms need not be unity):

$V=\left[\begin{array}{ccc}1& 1& 1\\ {\lambda }_1& {\lambda }_2& {\lambda }_3\\ {\lambda }_1^2& {\lambda }_2^2& {\lambda }_3^2\end{array}\right]=\left[\begin{array}{ccc}1& 1& 1\\ 0& -0.1& -0.4\\ 0& 0.01& 0.16\end{array}\right]$

The MATLAB command for similarity transformation is ss2ss. It requires the inverse $T_r^{-1}$ of the similarity transformation matrix T_r . Two MATLAB commands transform to diagonal form. The first requires the modal matrix of eigenvectors V and the system to be transformed

$\gg \mathbf{s}\_\mathbf{d}\mathbf{i}\mathbf{a}\mathbf{g}=\mathbf{s}\mathbf{s}\mathbf{2}\mathbf{s}\mathbf{s}\left(\mathbf{s}\mathbf{y}\mathbf{s}\mathbf{t}\mathbf{e}\mathbf{m},\mathbf{i}\mathbf{n}\mathbf{v}\left(\mathbf{v}\right)\right)$

The second uses similarity transformation but does not require the transformation matrix

$\gg \mathbf{s}\_\mathbf{\text{d}}\mathbf{\text{i}}\mathbf{\text{a}}\mathbf{\text{g}}=\mathbf{c}\mathbf{a}\mathbf{n}\mathbf{n}\mathbf{o}\mathbf{n}\left(\mathbf{s}\mathbf{y}\mathbf{s}\mathbf{t}\mathbf{e}\mathbf{m},{}^{\text{'}}\mathbf{\text{modal}}^{\text{'}}\right)$

The two commands yield

$\begin{array}{l}{\text{A}}_{\text{t}}=diag\left\{\begin{array}{ccc}0,& -0.1,& -0.4\end{array}\right\},{\text{B}}_{\text{t}}={\left[\begin{array}{ccc}25& 33.5012& 9.0738\end{array}\right]}^T\\ {\text{C}}_{\text{t}}=\left[\begin{array}{ccc}1.0000& -0.9950& 0.9184\end{array}\right]{\text{D}}_{\text{t}}=0\end{array}$

For complex conjugate eigenvalues, the command canon yields a real realization but its state matrix is not in diagonal form, whereas ss2ss will yield a diagonal but complex matrix.

7.9.1 Invariance of transfer functions and characteristic equations

Similar systems can be viewed as different representations of the same systems. This is justified by the following theorem.

Theorem 7.1

Similar systems have identical transfer functions and characteristic polynomials.

Proof

Consider the characteristic polynomials of similar realizations (A, B, C, D) and (A ₁, B ₁, C ₁, D):

$\begin{array}{c}\det \left(z{I}_n-A_1\right)=\det \left(z{I}_n-T_r^{-1}A{T}_r\right)\\ =\det \left[T_r^{-1}\left(z{I}_n-A\right)T_r\right]\\ =\det \left[T_r^{-1}\right]\det \left(z{I}_n-A\right)\det \left[T_r\right]=\det \left(z{I}_n-A\right)\end{array}$

where we used the identity $\det \left[T_r^{-1}\right]\times \det \left[T_r\right]=1$ .

The transfer function matrix is

$\begin{array}{c}{G}_1\left(s\right)=C_1\left[\mathit{z}{\mathit{I}}_n-A_1\right]B_1+D_1\\ =C{T}_r\left[z{I}_n-T_r^{-1}A{T}_r\right]T_r^{-1}B+D\\ =C\left[T_r\left(z{I}_n-T_r^{-1}A{T}_r\right)T_r^{-1}\right]B+D\\ =C\left[z{I}_n-A\right]B+D=G\left(s\right)\end{array}$

where we used the identity (A B C)⁻¹  = C ^{– 1} B ⁻¹ A ⁻¹.

Clearly, not all systems with the same transfer function are similar, due to the possibility of pole-zero cancellation. Systems that give rise to the same transfer function are said to be equivalent.

Example 7.22

Show that the following system is equivalent to the system shown in Example 7.19 (2).

$\begin{array}{c}x\left(k+1\right)=0.8187x\left(k\right)+9.0635\times {10}^{-2}u\left(k\right)\\ y\left(k\right)=x\left(k\right)\end{array}$

Solution

The transfer function of the system is

$G\left(z\right)=\frac{9.0635\times {10}^{-2}}{z-0.8187}$

which is identical to the reduced transfer function of Example 7.19 (2).

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128144336000077

State estimation strategy for continuous-time systems with time-varying delay via a novel L-K functional

Feisheng Yang , ... Peipei Kang , in Control Strategy for Time-Delay Systems, 2021

5.5.2 3rd order example

Consider the 3D continuous-time system in (5.1) with time-varying delay satisfying (5.2) and

$A=\left[\begin{array}{ccc}\hfill 1.06\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 1.42\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 0.88\hfill \end{array}\right],A_d=\left[\begin{array}{ccc}\hfill -0.32\hfill & \hfill 0.85\hfill & \hfill -1.36\hfill \\ \hfill 1.10\hfill & \hfill 0.41\hfill & \hfill -0.50\hfill \\ \hfill 0.42\hfill & \hfill 0.82\hfill & \hfill -0.95\hfill \end{array}\right],C=\left[\begin{array}{ccc}\hfill 1\hfill & \hfill 0.5\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill -0.5\hfill & \hfill 0.6\hfill \end{array}\right],D=\left[\begin{array}{ccc}\hfill 0\hfill & \hfill 1\hfill & \hfill 0.2\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 0.5\hfill \end{array}\right],L_0=\left[\begin{array}{ccc}\hfill 0\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 0\hfill \end{array}\right],L_1=\left[\begin{array}{ccc}\hfill 1\hfill & \hfill 0\hfill & \hfill 0.5\hfill \\ \hfill 1\hfill & \hfill 0\hfill & \hfill 1\hfill \\ \hfill 0\hfill & \hfill -1\hfill & \hfill 1\hfill \end{array}\right],B_1=\left[\begin{array}{ccc}\hfill 0.2\hfill & \hfill 0.2\hfill & \hfill 0.2\hfill \end{array}\right],B_2=\left[\begin{array}{cc}\hfill 0.4\hfill & \hfill -0.3\hfill \end{array}\right].$

Given different $H_{\infty }$ performance indices γ, the maximum allowable delay upper bounds are given in Table 5.3, and the $H_{\infty }$ performance indices ${\gamma }_{\min }$ for given different time delay upper bounds τ are calculated in Table 5.4.

Table 5.3. Maximum allowable delay upper bounds τ for different γ.

τ	γ
Method	10	20	30	40	50
Theorem 5.2	0.31	0.33	0.34	0.37	0.39

Table 5.4. Minimum H _∞ performance index γ _min for different τ.

γ min	τ
Method	0.3	0.5	0.7	0.8	1.0
Theorem 5.2	7	48	58	62	85

For the simulation verification, we design a suitable $H_{\infty }$ state estimator for the system under study. Choosing the time-varying delay $h(t)=0.5+\sin (t)$ as shown in Fig. 5.6 and applying Theorem 5.2, we obtain the following $H_{\infty }$ estimator gain matrix:

Figure 5.6. The curve for time-varying delay $0.5+\sin (t)$ .

$K=\left[\begin{array}{cc}\hfill 0.1059\hfill & \hfill 14.9374\hfill \\ \hfill 2.8186\hfill & \hfill -21.0109\hfill \\ \hfill 0.9271\hfill & \hfill -7.7306\hfill \end{array}\right].$

By choosing the initial condition as $e(0)={[31.5-1.5]}^T$ Fig. 5.7 depicts the responses of the estimation errors associated with the above estimator. Clearly, the designed $H_{\infty }$ state estimator produces pretty good estimation on the system states. The simulation results further confirm the effectiveness of Theorem 5.2 to guaranteed $H_{\infty }$ performance state estimation of the continuous-time system with time-varying delay.

Figure 5.7. The state estimation errors in the 3rd order example.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780323853477000109

12th International Symposium on Process Systems Engineering and 25th European Symposium on Computer Aided Process Engineering

Muxin Sun , ... Efstratios N. Pistikopoulos , in Computer Aided Chemical Engineering, 2015

5.2 Comparison with multiparametric programming

By discretising the continuous time system, the optimization problem can be formulated as conventional multiparametric programming problem (Pistikopoulos, 2012). Critical regions can be obtained by solving mpQP, where input is assumed constant for each time interval, are shown in Figure 4. Critical regions for mpQP are approximations of the regions for mpDO, the smaller the sample time, the more closely the critical regions approach to the exact solutions in Figure 2. A special region in both cases is the middle one, namely CR08 (a) or CR15 (b), which represents the case that no constraint is active. The others are all with the first couple of time intervals containing active constrains, with same solution structure as CR02 or CR03 in Figure 2. This unifies the definition of critical regions for different branches of multiparametric programing and also provides a way to explore the solution structures of mpDO by solving approximate mpQP.

Figure 4. Critical Regions with sampling time (a) Ts   =   0.1 and (b) Ts   =   0.05.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444635778500978

The Z-transform

Luis F. Chaparro , Aydin Akan , in Signals and Systems Using MATLAB (Third Edition), 2019

10.6 State Variable Representation

Modern control theory uses the state-variable representation of systems, whether continuous- or discrete-time. In this section we introduce the discrete-time state-variable representation of systems. In many respects, it is very similar to the state-variable representation of continuous-time systems.

State variables are the memory of a system. In the discrete-time, just like in the continuous-time, knowing the state of a system at a present index n provides the necessary information from the past that together with present and future inputs allows us to calculate the present and future outputs of the system. The advantage of a state-variable representation over a transfer function is the inclusion of initial conditions in the analysis and the ability to use it in multiple-input and multiple-output systems.

Assume a discrete-time system is represented by a difference equation (which could have been obtained from an ordinary differential equation representing a continuous-time system) where $x[n]$ is the input and $y[n]$ the output:

(10.42) $y[n]+a_{1}y[n-1]+\cdots +a_{N}y[n-N]=b_{0}x[n]+b_{1}x[n-1]+\cdots +b_{M}x[n-M],n\ge 0,$

where $M\le N$ and the initial conditions of the system are $\{y[i],-N\le i\le -1\}$ . As in the continuous-time, the state-variable representation of a discrete-time system is not unique. To obtain a state-variable representation from the difference equation in (10.42) let

(10.43) $\begin{array}{ccc}\hfill & \hfill \hfill & v_1[n]=y[n-1],\hfill \\ \hfill & \hfill \hfill & v_2[n]=y[n-2],\hfill \\ \hfill & \hfill \hfill & ⋮\hfill \\ \hfill & \hfill \hfill & v_N[n]=y[n-N].\hfill \end{array}$

We then obtain from the difference equation the state-variable equations:

(10.44) $\begin{array}{ccc}\hfill & \hfill \hfill & v_1[n+1]=y[n]=-a_1{v}_1[n]-\cdots -a_N{v}_N[n]+b_{0}x[n]+\cdots +b_{M}x[n-M],\hfill \\ \hfill & \hfill \hfill & v_2[n+1]=v_1[n],\hfill \\ \hfill & \hfill \hfill & ⋮\hfill \\ \hfill & \hfill \hfill & v_N[n+1]=v_{N-1}[n],\hfill \end{array}$

and the output equation

$y[n]=-a_1{v}_1[n]-\cdots -a_N{v}_N[n]+b_{0}x[n]+b_{1}x[n-1]+\cdots +b_{M}x[n-M].$

These state and output equations can then be written in a matrix form

(10.45) $\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{v}[n+1]=\mathbf{A}\mathbf{v}[n]+\mathbf{B}\mathbf{x}[n],\hfill \\ \hfill & \hfill \hfill & y[n]={\mathbf{c}}^T\mathbf{v}[n]+{\mathbf{d}}^T\mathbf{x}[n],n\ge 0,\hfill \end{array}$

by appropriately defining the matrices A and B, the vectors c and d as well as the state vector $\mathbf{v}[n]$ and the input vector $\mathbf{x}[n]$ .

Fig. 10.13 displays the block diagrams for a delay, a constant multiplier and an adder that are used in obtaining block diagrams for discrete-time system. Different from the continuous-time system representation, the representation of discrete-time systems does not require the equivalent of integrators, instead it uses delays.

Figure 10.13. Block diagrams of different components used to represent discrete-time systems (top to bottom): delay, constant multiplier and adder.

Example 10.27

A continuous-time system is represented by the ordinary differential equation

$\frac{d^{2}y(t)}{d{t}^2}+\frac{dy(t)}{dt}+y(t)=x(t)t\ge 0.$

To convert this ordinary differential equation into a difference equation we approximate the derivatives by

$\begin{array}{ccc}\hfill & \hfill \hfill & \frac{dy(t)}{dt}\approx \frac{y(t)-y(t-T_s)}{T_s},\hfill \\ \hfill & \hfill \hfill & \frac{d^{2}y(t)}{d{t}^2}\approx \frac{y(t)-2y(t-T_s)+y(t-2T_s)}{T_s^2}.\hfill \end{array}$

Find the difference equation when $T_s=1$ and $t=n{T}_s$ , and find the corresponding state and output equations.

Solution: Replacing the approximations of the derivatives for $T_s=1$ and letting $t=n{T}_s=n$ we obtain the difference equation

$\begin{array}{ccc}\hfill & \hfill \hfill & (y[n]-2y[n-1]+y[n-2])+(y[n]-y[n-1])+y[n]=x[n]\text{or}\hfill \\ \hfill & \hfill \hfill & y[n]-y[n-1]+\frac{1}{3}y[n-2]=\frac{1}{3}x[n],\hfill \end{array}$

which can be realized as indicated by the block diagram in Fig. 10.14. Because the second-order difference equation is realized using two unit delays, this realization is called minimal.

Figure 10.14. Block diagram of y[n]−y[n − 1]+(1/3)y[n − 2]=(1/3)x[n] displaying the state variables v ₁[n]=y[n − 2] and v ₂[n]=y[n − 1].

Letting the outputs of the delays be the state variables

$v_1[n]=y[n-2],v_2[n]=y[n-1],$

we then have the following matrix equation for the state:

$\left[\begin{array}{c}\hfill v_1[n+1]\hfill \\ \hfill v_2[n+1]\hfill \end{array}\right]=\underset{\mathbf{A}}{\underbrace{\left[\begin{array}{cc}\hfill 0& \hfill 1\\ \hfill -1/3& \hfill 1\end{array}\right]}}\left[\begin{array}{c}\hfill v_1[n]\hfill \\ \hfill v_2[n]\hfill \end{array}\right]+\underset{\mathbf{b}}{\underbrace{\left[\begin{array}{c}\hfill 0\hfill \\ \hfill 1/3\hfill \end{array}\right]}}x[n].$

The output equation is then

$\begin{array}{ccc}\hfill & \hfill \hfill & y[n]=-\frac{1}{3}{v}_1[n]+v_2[n]+\frac{1}{3}x[n]\text{or in matrix form}\hfill \\ \hfill & \hfill \hfill & y[n]=\underset{{\mathbf{c}}^T}{\underbrace{\left[\begin{array}{cc}\hfill -\frac{1}{3}\hfill & \hfill 1\hfill \end{array}\right]}}\left[\begin{array}{c}\hfill v_1[n]\hfill \\ \hfill v_2[n]\hfill \end{array}\right]+\underset{d}{\underbrace{\left[\frac{1}{3}\right]}}x[n].\hfill \end{array}$

The state variables are not unique. Indeed we can use an invertible transformation matrix F to define a new set of state variables

$\mathbf{w}[n]=\mathbf{F}\mathbf{v}[n]$

with a matrix representation for the new state variables and the output given by

$\begin{array}{ccc}\hfill \mathbf{w}[n+1]& \hfill =\hfill & \mathbf{F}\mathbf{v}[n+1]=\mathbf{F}\mathbf{A}\mathbf{v}[n]+\mathbf{F}\mathbf{b}x[n]\hfill \\ \hfill & \hfill =\hfill & \mathbf{F}\mathbf{A}{\mathbf{F}}^{-1}\mathbf{w}[n]+\mathbf{F}\mathbf{b}x[n],\hfill \\ \hfill y[n]& \hfill =\hfill & {\mathbf{c}}^T\mathbf{v}[n]+\mathbf{d}x[n]={\mathbf{c}}^T{\mathbf{F}}^{-1}\mathbf{w}[n]+\mathbf{d}x[n].\text{□}\hfill \end{array}$

Example 10.28

A discrete-time system is represented by the difference equation

$y[n]-y[n-1]-y[n-2]=x[n]+x[n-1]$

where $x[n]$ is the input and $y[n]$ the output. Obtain a state-variable representation for it.

Solution: Notice that this difference equation has $x[n-1]$ as well as $x[n]$ in the input. The transfer function corresponding to the difference equation is

$H(z)=\frac{Y(z)}{X(z)}=\frac{1+z^{-1}}{1-z^{-1}-z^{-2}},$

i.e., it is not a "constant-numerator" transfer function. A direct realization of the system in this case will not be minimal. Indeed, a block diagram obtained from the difference equation (see Fig. 10.15) shows that three delays are needed to represent this second-order system. It is important to realize that this representation despite being non-minimal is a valid representation. This is different from an analogous situation in the continuous-time representation where the input and its derivatives are present. Differentiators in the continuous-time representation are deemed not acceptable, while delays in the discrete-time representations are.   □

Figure 10.15. Block diagram of y[n]−y[n − 1]−y[n − 2]=x[n]+x[n − 1] displaying the state variables v ₁[n] and v ₂[n], and corresponding to a non-minimal realization.

Solution of the State and Output Equations

The solution of the state-variable equations

$\mathbf{v}[n+1]=\mathbf{A}\mathbf{v}[n]+\mathbf{B}\mathbf{x}[n],n\ge 0,$

can be obtained recursively:

$\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{v}[1]=\mathbf{A}\mathbf{v}[0]+\mathbf{B}\mathbf{x}[0],\hfill \\ \hfill & \hfill \hfill & \mathbf{v}[2]=\mathbf{A}\mathbf{v}[1]+\mathbf{B}\mathbf{x}[1]={\mathbf{A}}^2\mathbf{v}[0]+\mathbf{AB}\mathbf{x}[0]+\mathbf{B}\mathbf{x}[1],\hfill \\ \hfill & \hfill \hfill & ⋮\hfill \\ \hfill & \hfill \hfill & \mathbf{v}[n]={\mathbf{A}}^n\mathbf{v}[0]+\sum _{k=0}^{n-1}{\mathbf{A}}^{n-1-k}\mathbf{B}\mathbf{x}[k],\hfill \end{array}$

where ${\mathbf{A}}^0=\mathbf{I}$ , the identity matrix. The complete solution is then obtained,

(10.46) $y[n]=\underset{y_{zi}[n]\text{zero-input response}}{\underbrace{{\mathbf{c}}^T{\mathbf{A}}^n\mathbf{v}[0]}}+\underset{y_{zs}[n]\text{zero-state response}}{\underbrace{\sum _{k=0}^{n-1}{\mathbf{c}}^T{\mathbf{A}}^{n-1-k}\mathbf{B}\mathbf{x}[k]+\mathbf{d}\mathbf{x}[n]}}.$

By definition of the state variables in Equation (10.43) the initial conditions of the state variables coincide with the initial conditions of the system,

(10.47) $v_1[0]=y[-1],v_2[0]=y[-2],\cdots v_N[0]=y[-N].$

Using the Z-transform we can obtain a close solution to the state and output equations. Indeed, if the state and output equations are in their matrix form

(10.48) $\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{v}[n+1]=\mathbf{A}\mathbf{v}[n]+\mathbf{B}\mathbf{x}[n],\hfill \\ \hfill & \hfill \hfill & y[n]={\mathbf{c}}^T\mathbf{v}[n]+{\mathbf{d}}^T\mathbf{x}[n],n\ge 0,\hfill \end{array}$

calling the Z-transforms of the state variables $V_i(z)=\mathcal{Z}(v_i[n])$ , for $i=1,\cdots ,N$ ; of the input $X_m(z)=\mathcal{Z}(x[n-m])$ , for $m=0,\cdots ,M$ , and of the output $Y(z)=\mathcal{Z}(y[n])$ we have the following matrix expression for the Z-transforms of the state equations in (10.48):

$z\mathbf{V}(z)-z\mathbf{v}[0]=\mathbf{A}\mathbf{V}(z)+\mathbf{B}\mathbf{X}(z)\text{or}(z\mathbf{I}-\mathbf{A})\mathbf{V}(z)=z\mathbf{v}[0]+\mathbf{B}\mathbf{X}(z)$

using the Z-transform of $v_i[n+1]$ , and $\mathbf{v}[0]$ is the vector of initial conditions of the state variables and I the unit matrix. Assuming that the inverse of $(z\mathbf{I}-\mathbf{A})$ exists, i.e., $\text{det}(z\mathbf{I}-\mathbf{A})\ne 0$ , we can solve for $\mathbf{V}(z)$ :

(10.49) $\begin{array}{ccc}\hfill \mathbf{V}(z)& \hfill =\hfill & {(z\mathbf{I}-\mathbf{A})}^{-1}z\mathbf{v}[0]+{(z\mathbf{I}-\mathbf{A})}^{-1}\mathbf{B}\mathbf{X}(z)\hfill \\ \hfill & \hfill =\hfill & \frac{\text{Adj}(z\mathbf{I}-\mathbf{A})}{\text{det}(z\mathbf{I}-\mathbf{A})}z\mathbf{v}[0]+\frac{\text{Adj}(z\mathbf{I}-\mathbf{A})}{\text{det}(z\mathbf{I}-\mathbf{A})}\mathbf{B}\mathbf{X}(z)\hfill \end{array}$

by expressing the inverse of a matrix in terms of its adjoint and determinant. We can then obtain the Z-transform of the output as

(10.50) $Y(z)=\frac{{\mathbf{c}}^T\text{Adj}(z\mathbf{I}-\mathbf{A})}{\text{det}(z\mathbf{I}-\mathbf{A})}z\mathbf{v}[0]+[\frac{{\mathbf{c}}^T\text{Adj}(z\mathbf{I}-\mathbf{A})}{\text{det}(z\mathbf{I}-\mathbf{A})}\mathbf{B}+\mathbf{d}]\mathbf{X}(z).$

If the initial conditions $\mathbf{v}[0]$ are zero, and the input is $x[n]$ , then we find that the transfer function is given by

(10.51) $H(z)=\frac{Y(z)}{X(z)}=\frac{{\mathbf{c}}^T\text{Adj}(z\mathbf{I}-\mathbf{A})}{\text{det}(z\mathbf{I}-\mathbf{A})}\mathbf{b}+d.$

Example 10.29

Consider the state-variable representation of the system in Example 10.27 with matrices:

$\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{A}=\left[\begin{array}{cc}\hfill 0\hfill & \hfill 1\\ \hfill -1/3\hfill & \hfill 1\end{array}\right],\mathbf{b}=\left[\begin{array}{c}\hfill 0\hfill \\ \hfill 1/3\hfill \end{array}\right],\hfill \\ \hfill & \hfill \hfill & {\mathbf{c}}^T=\left[\begin{array}{cc}\hfill -\frac{1}{3}\hfill & \hfill 1\hfill \end{array}\right],d=\left[\frac{1}{3}\right].\hfill \end{array}$

Determine the transfer function of the system.

Solution: Instead of finding the inverse ${(z\mathbf{I}-\mathbf{A})}^{-1}$ we can use Cramer's rule to find $Y(z)$ . Indeed, writing the state equations with zero initial conditions

$\underset{(z\mathbf{I}-\mathbf{A})}{\underbrace{\left[\begin{array}{cc}\hfill z\hfill & \hfill -1\hfill \\ \hfill 1/3\hfill & \hfill z-1\hfill \end{array}\right]}}\underset{\mathbf{V}(z)}{\underbrace{\left[\begin{array}{c}\hfill V_1(z)\hfill \\ \hfill V_2(z)\hfill \end{array}\right]}}=\underset{\mathbf{b}X(z)}{\underbrace{\left[\begin{array}{c}\hfill 0\hfill \\ \hfill X(z)/3\hfill \end{array}\right]}}$

according to Cramer's rule we have

$\begin{array}{ccc}\hfill V_1(z)& \hfill =\hfill & \frac{\text{det}\left[\begin{array}{cc}\hfill 0\hfill & \hfill -1\hfill \\ \hfill X(z)/3\hfill & \hfill z-1\hfill \end{array}\right]}{\mathrm{\Delta }(z)}=\frac{X(z)/3}{\mathrm{\Delta }(z)},\hfill \\ \hfill V_2(z)& \hfill =\hfill & \frac{\text{det}\left[\begin{array}{cc}\hfill z\hfill & \hfill 0\hfill \\ \hfill 1/3\hfill & \hfill X(z)/3\hfill \end{array}\right]}{\mathrm{\Delta }(z)}=\frac{zX(z)/3}{\mathrm{\Delta }(z)},\hfill \\ \hfill \mathrm{\Delta }(z)& \hfill =\hfill & z^2-z+1/3,\hfill \end{array}$

and the Z-transform of the output is obtained after replacing the Z-transforms of the state variables:

$\begin{array}{ccc}\hfill Y(z)& \hfill =\hfill & \left[\begin{array}{cc}\hfill -1/3\hfill & \hfill 1\hfill \end{array}\right]\left[\begin{array}{c}\hfill V_1(z)\hfill \\ \hfill V_2(z)\hfill \end{array}\right]+\frac{X(z)}{3}\hfill \\ \hfill & \hfill =\hfill & \frac{-1/9+z/3+(z^2/3-z/3+1/9)}{z^2-z+1/3}X(z)=\frac{z^2/3}{z^2-z+1/3}X(z),\text{so that}\hfill \\ \hfill H(z)& \hfill =\hfill & \frac{Y(z)}{X(z)}=\frac{z^2/3}{z^2-z+1/3}=\frac{1/3}{1-z^{-1}+z^{-2}/3}.\text{□}\hfill \end{array}$

Example 10.30

The transfer function of a LTI discrete-time system is

$H(z)=\frac{z^{-2}}{1-0.5z^{-1}+0.5z^{-2}}.$

Obtain a minimal state-variable realization (i.e., a realization that only uses two delays, corresponding to the second-order system). Determine the initial conditions of the state variables in terms of initial values of the output.

Solution: The given transfer function is not "constant-numerator" transfer function and the numerator indicates the system has delayed inputs. Letting $x[n]$ be the input and $y[n]$ the output with Z-transforms $X(z)$ and $Y(z)$ , respectively, we factor $H(z)$ as follows:

$H(z)=\frac{Y(z)}{X(z)}=\underset{Y(z)/W(z)}{\underbrace{z^{-2}}}\times \underset{W(z)/X(z)}{\underbrace{\frac{1}{1-0.5z^{-1}+0.5z^{-2}}}}$

so that we have the following equations:

$\begin{array}{ccc}\hfill & \hfill \hfill & \text{w}[n]=0.5\text{w}[n-1]-0.5\text{w}[n-2]+x[n],\hfill \\ \hfill & \hfill \hfill & y[n]=\text{w}[n-2].\hfill \end{array}$

The corresponding block diagram in Fig. 10.16 is a minimal realization as it only uses two delays. As shown in the figure, the state variables are defined as

Figure 10.16. Minimal realization of a system with transfer function H(z)=z ⁻²/(1 − 0.5z ⁻¹ + 0.5z ⁻²).

$v_1[n]=\text{w}[n-1],v_2[n]=\text{w}[n-2].$

The state and output equations are in matrix form

$\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{v}[n+1]=\left[\begin{array}{cc}\hfill 1/2\hfill & \hfill -1/2\\ \hfill 1\hfill & \hfill 0\end{array}\right]\mathbf{v}[n]+\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 0\hfill \end{array}\right]x[n],\hfill \\ \hfill & \hfill \hfill & y[n]=\left[\begin{array}{cc}\hfill 0\hfill & \hfill 1\hfill \end{array}\right]\mathbf{v}[n].\hfill \end{array}$

Expressing the state variables in terms of the output:

$v_1[n]=\text{w}[n-1]=y[n+1],v_2[n]=y[n],$

we then have $v_2[0]=y[0]$ and $v_1[0]=y[1]$ . The reason for these initial conditions is that in the difference equation

$y[n]-0.5y[n-1]+0.5y[n-2]=x[n-2]$

the input is delayed two samples and so the initial conditions should be $y[0]$ and $y[1]$ .   □

Example 10.31

Consider the state-variable representation of a system having the matrices

$\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{A}=\left[\begin{array}{cc}\hfill 0\hfill & \hfill 1\hfill \\ \hfill 1/8\hfill & \hfill 1/4\hfill \end{array}\right],\mathbf{b}=\left[\begin{array}{c}\hfill 0\hfill \\ \hfill 1\hfill \end{array}\right],\hfill \\ \hfill & \hfill \hfill & {\mathbf{c}}^T=\left[\begin{array}{cc}\hfill 1\hfill & \hfill 0\hfill \end{array}\right],d=\left[0\right].\hfill \end{array}$

Determine the zero-input response of the system for arbitrary initial conditions $\mathbf{v}[0]$ .

Solution: The zero-input response in the time-domain is

$y[n]={\mathbf{c}}^T{\mathbf{A}}^n\mathbf{v}[0]$

requiring one to compute the matrix ${\mathbf{A}}^n$ . To do so, consider

$\mathcal{Z}[{\mathbf{A}}^n]=\sum _{n=0}^{\infty }{\mathbf{A}}^n{z}^{-n}={(\mathbf{I}-\mathbf{A}{z}^{-1})}^{-1},\text{det}(\mathbf{I}-\mathbf{A}{z}^{-1})\ne 0,$

which can be verified by pre-multiplying by $\mathbf{I}-\mathbf{A}{z}^{-1}$ the above equation, giving

$[\mathbf{I}-\mathbf{A}{z}^{-1}]\sum _{n=0}^{\infty }{\mathbf{A}}^n{z}^{-n}=\sum _{n=0}^{\infty }{\mathbf{A}}^n{z}^{-n}-\sum _{n=0}^{\infty }{\mathbf{A}}^{n+1}{z}^{-(n+1)}=\mathbf{I}$

or that the infinite summation is the inverse of $(\mathbf{I}-\mathbf{A}{z}^{-1})$ . In positive powers of z the Z-transform of ${\mathbf{A}}^n$ is

$\mathcal{Z}[{\mathbf{A}}^n]={(z^{-1}[\mathbf{I}z-\mathbf{A}])}^{-1}=z{(\mathbf{I}z-\mathbf{A})}^{-1}.$

Now, using the fact that for a $2\times 2$ matrix

$\mathbf{F}=\left[\begin{array}{cc}\hfill a\hfill & \hfill b\hfill \\ \hfill c\hfill & \hfill d\hfill \end{array}\right]\Rightarrow {\mathbf{F}}^{-1}=\frac{1}{ad-bc}\left[\begin{array}{cc}\hfill d\hfill & \hfill -b\hfill \\ \hfill -c\hfill & \hfill a\hfill \end{array}\right]$

provided that the determinant of the matrix, $ad-bc\ne 0$ , we then have

$\begin{array}{cc}\hfill & \mathbf{I}-\mathbf{A}{z}^{-1}=\left[\begin{array}{cc}\hfill 1\hfill & \hfill -z^{-1}\hfill \\ \hfill -z^{-1}/8\hfill & \hfill 1-z^{-1}/4\hfill \end{array}\right]\hfill \\ \hfill & \Rightarrow {(\mathbf{I}-\mathbf{A}{z}^{-1})}^{-1}=\frac{1}{1-z^{-1}/4-z^{-2}/8}\left[\begin{array}{cc}\hfill 1-z^{-1}/4\hfill & \hfill z^{-1}\hfill \\ \hfill z^{-1}/8\hfill & \hfill 1\hfill \end{array}\right].\hfill \end{array}$

We then need to determine the inverse Z-transform of the four entries to find ${\mathbf{A}}^n$ . If we let

$\begin{array}{ccc}\hfill & \hfill \hfill & P(z)=\frac{1}{1-z^{-1}/4-z^{-2}/8}=\frac{2/3}{1-0.5z^{-1}}+\frac{1/3}{1+0.25z^{-1}}\text{with inverse Z-transform}\hfill \\ \hfill & \hfill \hfill & p[n]=(\frac{2}{3}{0.5}^n+\frac{1}{3}{(-0.25)}^n)u[n]\hfill \end{array}$

then

${\mathbf{A}}^n=\left[\begin{array}{cc}\hfill p[n]-0.25p[n-1]\hfill & \hfill p[n-1]\hfill \\ \hfill 0.125p[n-1]\hfill & \hfill p[n]\hfill \end{array}\right],$

which is used in finding the zero-input response

$\begin{array}{ccc}\hfill y[n]& \hfill =\hfill & \left[\begin{array}{cc}\hfill 1\hfill & \hfill 0\hfill \end{array}\right]\left[\begin{array}{cc}\hfill p[n]-0.25p[n-1]\hfill & \hfill p[n-1]\hfill \\ \hfill 0.125p[n-1]\hfill & \hfill p[n]\hfill \end{array}\right]\mathbf{v}[0]\hfill \\ \hfill & \hfill =\hfill & (p[n]-0.25p[n-1])v_1[0]+p[n-1]v_2[0].\text{□}\hfill \end{array}$

Canonical Realizations

Just like in the continuous case, discrete-time state-variable realizations have different canonical forms. In this section we will illustrate the process for obtaining a parallel realization to demonstrate the similarity.

The state-variable parallel realization is obtained by realizing each of the partial fraction expansion components of a transfer function. Consider the case of simple poles of the system, so that the transfer function is

$H(z)=\frac{Y(z)}{X(z)}=\sum _{i=1}^N\underset{H_i(z)}{\underbrace{\frac{A_i}{1-{\alpha }_i{z}^{-1}}}}$

so that

$Y(z)=\sum _{i=1}^N\underset{Y_i(z)}{\underbrace{H_i(z)X(z)}}.$

The state variable for $Y_i(z)(1-{\alpha }_i{z}^{-1})=X(z)$ or $y_i[n]-{\alpha }_i{y}_i[n-1]=x[n]$ would be $v_i[n]=y_i[n-1]$ so that

$v_i[n+1]=y_i[n]={\alpha }_i{v}_i[n]+x[n],i=1,\cdots ,N.$

For the whole system we then have

$\begin{array}{ccc}\hfill & \hfill \hfill & v_1[n+1]={\alpha }_1{v}_1[n]+x[n],\hfill \\ \hfill & \hfill \hfill & v_2[n+1]={\alpha }_2{v}_2[n]+x[n],\hfill \\ \hfill & \hfill \hfill & \cdots \hfill \\ \hfill & \hfill \hfill & v_N[n+1]={\alpha }_N{v}_N[n]+x[n],\hfill \end{array}$

and the output is

$y[n]=\sum _{i=1}^N{y}_i[n]=\sum _{i=1}^N{\alpha }_i{v}_i[n]+Nx[n]$

or in matrix form

$\begin{array}{ccc}\hfill & \hfill \hfill & \mathbf{v}[n+1]=\mathbf{A}\mathbf{v}[n]+\mathbf{b}x[n],\hfill \\ \hfill & \hfill \hfill & y[n]={\mathbf{c}}^T\mathbf{v}[n]+dx[n],\hfill \\ \hfill & \hfill \hfill & \mathbf{A}=\left[\begin{array}{cccc}\hfill {\alpha }_1\hfill & \hfill 0\hfill & \hfill \cdots \hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill {\alpha }_2\hfill & \hfill \cdots \hfill & \hfill 0\hfill \\ \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill & \hfill ⋮\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill \cdots \hfill & \hfill {\alpha }_N\hfill \end{array}\right],\mathbf{b}=\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 1\hfill \\ \hfill ⋮\hfill \\ \hfill 1\hfill \end{array}\right],\hfill \\ \hfill & \hfill \hfill & {\mathbf{c}}^T=\left[\begin{array}{cccc}\hfill 1\hfill & \hfill 1\hfill & \hfill \cdots \hfill & \hfill 1\hfill \end{array}\right],d=N.\hfill \end{array}$

In general, whenever the transfer function is factored in terms of first- and second-order systems, with real coefficients, a parallel realization is obtained by doing a partial fraction expansion in terms of first- and second-order components and realizing each of these components. The following example illustrates the procedure.

Example 10.32

Obtain a parallel state-variable realization from the transfer function

$H(z)=\frac{z^3}{(z+0.5)[{(z-0.5)}^2+0.25]}.$

Show a block realization of the corresponding state and output equations.

Solution: If we obtain a partial fraction expansion of the form

$H(z)=\frac{1}{(1+0.5z^{-1})[{(1-0.5z^{-1})}^2+0.25z^{-2}]}=\frac{A_1}{1+0.5z^{-1}}+\frac{A_2+A_3{z}^{-1}}{1-z^{-1}+0.5z^{-2}}$

where

$A_1=H(z)(1+0.5z^{-1}){|}_{z^{-1}=-2}=\frac{1}{1-z^{-1}+0.5z^{-2}}{|}_{z^{-1}=-2}=\frac{1}{5}$

and then

$\begin{array}{ccc}\hfill \frac{A_2+A_3{z}^{-1}}{1-z^{-1}+0.5z^{-2}}& \hfill =\hfill & \frac{1}{(1+0.5z^{-1})(1-z^{-1}+0.5z^{-2})}-\frac{1/5}{1+0.5z^{-1}}\hfill \\ \hfill & \hfill =\hfill & \frac{1-0.2(1-z^{-1}+0.5z^{-2})}{(1+0.5z^{-1})(1-z^{-1}+0.5z^{-2})}=\frac{0.8-0.2z^{-1}}{1-z^{-1}+0.5z^{-2}}\hfill \end{array}$

so that $A_2=0.8$ and $A_3=-0.2$ . This allows us to write the output as

$Y(z)=\underset{Y_1(z)}{\underbrace{\frac{0.2X(z)}{1+0.5z^{-1}}}}+\underset{Y_2(z)}{\underbrace{\frac{(0.8-0.2z^{-1})X(z)}{1-z^{-1}+0.5z^{-2}}}}$

or for two difference equations

$\begin{array}{ccc}\hfill & \hfill \hfill & y_1[n]+0.5y_1[n-1]=0.2x[n],\hfill \\ \hfill & \hfill \hfill & y_2[n]-y_2[n-1]+0.5y[n-2]=0.8x[n]-0.2x[n-1],\hfill \end{array}$

for which we need to obtain state variables. Minimal realizations can be obtained as follows. The first system is a "constant-numerator" system, in terms of negative powers of z, and a minimal realization is obtained directly. The transfer function of the second difference equation can be written

$Y_2(z)=\underset{W(z)}{\underbrace{\frac{X(z)}{1-z^{-1}+0.5z^{-2}}}}\times \underset{Y_2(z)/W(z)}{\underbrace{(0.8-0.2z^{-1})}},$

from which we obtain the equations

$\begin{array}{ccc}\hfill & \hfill \hfill & \text{w}[n]-\text{w}[n-1]+0.5\text{w}[n-2]=x[n],\hfill \\ \hfill & \hfill \hfill & y_2[n]=0.8\text{w}[n]-0.2\text{w}[n-1].\hfill \end{array}$

Notice that the realization of these equations only require two delays. If we let

$\begin{array}{c}\hfill v_1[n]=\text{w}[n-1],\\ \hfill v_2[n]=\text{w}[n-2],\end{array}$

we get

$\begin{array}{ccc}\hfill & \hfill \hfill & v_1[n+1]=\text{w}[n]=v_1[n]-0.5v_2[n]+x[n],\hfill \\ \hfill & \hfill \hfill & v_2[n+1]=v_1[n],\hfill \\ \hfill & \hfill \hfill & y_2[n]=0.8(v_1[n]-0.5v_2[n]+x[n])-0.2v_1[n]=0.6v_1[n]-0.4v_2[n]+0.8x[n],\hfill \end{array}$

so that

$\begin{array}{ccc}\hfill & \hfill \hfill & \left[\begin{array}{c}\hfill v_1[n+1]\hfill \\ \hfill v_2[n+1]\hfill \end{array}\right]=\left[\begin{array}{cc}\hfill 1\hfill & \hfill -0.5\hfill \\ \hfill 1\hfill & \hfill 0\hfill \end{array}\right]\left[\begin{array}{c}\hfill v_1[n]\hfill \\ \hfill v_2[n]\hfill \end{array}\right]+\left[\begin{array}{c}\hfill 1\hfill \\ \hfill 0\hfill \end{array}\right]x[n],\hfill \\ \hfill & \hfill \hfill & y_2[n]=\left[\begin{array}{cc}\hfill 0.6\hfill & \hfill -0.4\hfill \end{array}\right]\left[\begin{array}{c}\hfill v_1[n]\hfill \\ \hfill v_2[n]\hfill \end{array}\right]+0.8x[n].\hfill \end{array}$

The following script shows how to obtain the state-variable representation for the two components of $H(z)$ using the function tf2ss. The second part of the script shows that by using the function ss2tf we can get back the two components of $H(z)$ . Notice that these are the same functions used in the continuous case. Fig. 10.17 shows the block diagram of the parallel realization.   □

Figure 10.17. Block diagram of a parallel state-variable realization of the system with H(z)=z ³/(z + 0.5)[(z−0.5)² + 0.25]. Notice this is a minimum realization with three delays corresponding to the third-order system.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128142042000211

Impulsive MPC schemes for biomedical processes: Application to type 1 diabetes

Pablo S. Rivadeneira , ... A.H. González , in Control Applications for Biomedical Engineering Systems, 2020

2 Dynamic systems with short-duration inputs

In many biomedical applications, the medication is administered to the patient in a discrete way, in the form of small pulses or impulses. For instance, in the AP problem, the insulin can be injected to the patient in both ways.

If a fixed time period is assumed, then it can be argued that an insulin dose is administered by a zero-order hold which maintains their rate all the period long, and changes its value only at the beginning of a new periods. This is in fact what is made in most of the approaches, although this modeling of the input is not realistic since insulin pumps are not able to do that. Actually, insulin doses are administered by pulses that covers only a short interval of these periods. If these pulses are extremely short in comparison with the period, they approach impulses (i.e., pulses of infinitesimal duration). These latter schemes work, for instance, to represent diseases that are controlled by taking pills one or two times per day. However, in the AP, it is not so clear if the (probably fast) insulin injection is short enough to have a good representation by the impulsive approaches. In any case, having a general pulse representation that account for both the zero-order hold input and impulsive input as particular cases is desired.

Consider the following affine continuous-time system:

(1) $\begin{array}{ll}\hfill \dot{x}(t)& =Ax(t)+B_{u}u(t)+B_{r}r(t)+E,x(0)=x_0,\hfill \\ \hfill y(t)& =Cx(t),\hfill \end{array}$

where x represents the state, u the control action, y the output of the system, and r the disturbances, with A, B _u , C, and B _r being the corresponding matrices, respectively, and E being a constant term. If the drug infusion u is injected to the system only at certain time instants given by kT (and quickly enough), where T is the fixed period, and $k\in \mathbb{N}$ , then it is possible to work under the impulsive system framework recently developed by Rivadeneira et al. (2018) (see Fig. 2 for a schematic plot). Formally, assume that the input is given by

Fig. 2. Impulsive scheme.

$\begin{array}{l}\hfill u(t)=u(kT)\delta (t-kT),t\in [kT,(k+1)T),k\in \mathbb{N},\end{array}$

where δ(t) is a generalized function (or distribution), Dirac delta, ^a that fulfills δ(0) = ∞, δ(t) = 0 for any t  ≠   0, and

$\begin{array}{l}\hfill {\int }_{-\infty }^{\infty }g(\zeta )\delta (\zeta )d\zeta =g(0),\end{array}$

for all continuous, compactly supported functions g.

The solution of Eq. (1) at each periods T can be divided into two parts, the first one, describing the system in the period [kT, kT + ΔT], and the second one describing the system in the period (kT + ΔT, (k + 1)T), for a positive and arbitrary smallΔ, 0 < Δ < 1:

$\begin{array}{l}\phi (t;x(kT),u,r)=e^{A(t-kT)}x(kT)+{\int }_{kT}^t{e}^{A(t-\zeta )}{B}_{u}u(\zeta )\delta (\zeta -kT)d\zeta \hfill \\ \hfill & +{\int }_{kT}^t{e}^{A(t-\zeta )}{B}_{r}r(\zeta )d\zeta +{\int }_{kT}^t{e}^{A(t-\zeta )}d\zeta E,\hfill \end{array}$

for t ∈ [kT, kT + Δ], and

$\begin{array}{ll}\phi (t;x(kT+\mathrm{\Delta }T),u,r)\hfill & =e^{A(t-kT)}x(kT+\mathrm{\Delta }T)+{\int }_{kT+\mathrm{\Delta }T}^t{e}^{A(t-\zeta )}{B}_{r}r(\zeta )d\zeta +{\int }_{kT+\mathrm{\Delta }T}^t{e}^{A(t-\zeta )}d\zeta E,\hfill \end{array}$

for $t\in \left(kT+\mathrm{\Delta },(k+1)T\right)$ .

Now, if we consider the limits of this solution for Δ → 0, it follows that

(2) $\begin{array}{l}\hfill x(k{T}^{+})\triangleq \underset{\mathrm{\Delta }\to 0}{lim}\phi (t;x(kT),u,r)=x(kT)+B_{u}u(kT),\end{array}$

where $x(k{T}^{+})\triangleq \underset{\mathrm{\Delta }\to 0}{lim}x(kT+\mathrm{\Delta }T)$ , and

$\begin{array}{l}x(t)\triangleq \phi (t;x(k{T}^{+}),u,r)=e^{A(t-kT)}x(k{T}^{+})+{\int }_{kT}^t{e}^{A(t-\zeta )}{B}_{r}r(\zeta )d\zeta +{\int }_{kT}^t{e}^{A(t-\zeta )}d\zeta E,\hfill \end{array}$

for t ∈ (kT, (k + 1)T). This latter solution is the one corresponding to the continuous-time system

(3) $\begin{array}{l}\hfill \dot{x}(t)=Ax(t)+B_{r}r(t)+E,x(0)=x_0,\end{array}$

for the aforementioned period of time.

This way, merging Eqs. (2), (3) into a single model, the following impulsive system can be obtained:

(4) $\begin{array}{l}\left\{\begin{array}{ll}\dot{x}(t)=Ax(t)+B_{r}r(t)+E,x(0)=x_0,& t\ne kT,\\ x(k{T}^{+})=x(kT)+B_{u}u(kT),k\in \mathbb{N},\end{array}\right.\hfill \end{array}$

where x(kT ⁺) denotes the limits of x(t) when t approaches kT from the right. In Eq. (4), the latter equation describes the discontinuities that the impulsive input causes into the system, while the former describes the free response (only affected by the disturbance r) between the discontinuities.

Remark 1

Note that the affine impulsive system (4) has no formal equilibrium points. In fact, there is no triplet $(x_s^i,u_s^i,r_s^i)$ fulfilling the condition:

(5) $\begin{array}{ll}\hfill & A{x}_s^i+B_r{r}_s^i+E=0,\hfill \end{array}$

(6) $\begin{array}{ll}\hfill & x_s^i=x_s^i+B_u{u}_s^i,\text{(no discontinuity)},\hfill \end{array}$

because the first equation has no solution, for $r_s^i=0$ (vanishing disturbance).

An interpretation of the latter situation is that it is not possible to maintain the system in a given point by only applying impulsive inputs, even when this input is null. However, as it will be shown in the following section, it is possible to maintain the system (switching) inside a given region. The conditions for that are: (i) to find states before and after the discontinuity that, although different, remain constant and (ii) to ensure that the transient state trajectories between these states remain inside a given set.

2.1 Underlying discrete-time subsystem

A natural way to obtain a discretization (discrete-time system) of the continuous-time impulsive system (4) is by sampling it with a sampling time given by the period of the impulses, T. According to the results presented by Rivadeneira et al. (2015, 2018), the idea is to characterize two discrete-time subsystems describing the evolution of the states x(kT) and x(kT ⁺); that is, the evolution of the states just before and after the discontinuities, respectively. If the disturbance r is assumed to be piece-wise constant, that is, r(t) = r(kT), for t ∈ [kT, (k + 1)T), $k\in \mathbb{N}$ , the subsystems are as follows:

(7) $\begin{array}{ll}\hfill x^{•}(k+1)& =A^{•}{x}^{•}(k)+B_u^{•}{u}^{•}(k)+B_r^{•}{r}^{•}(k)+E^{•},x^{•}(0)=x(0)=x_0,\hfill \end{array}$

(8) $\begin{array}{ll}\hfill x^{°}(k+1)& =A^{°}{x}^{°}(k)+B_u^{°}{u}^{°}(k)+B_r^{°}{r}^{°}(k)+E^{°},x^{°}(0)=x(0^{+})=x_0+B_{u}u(0),\hfill \end{array}$

where $x^{•}(k)\triangleq x(kT)$ , $x^{°}(k)\triangleq x(k{T}^{+})$ , A ^• = A° = e ^AT , $B_u^{•}=e^{AT}{B}_u$ , and $B_u^{°}=B_u$ . Furthermore, the input and disturbance of each subsystem are such that u°(k + 1) = u ^•(k) = u(kT), r°(k) = r ^•(k) = r(kT), and $E^{•}=E^{°}\triangleq {\int }_0^T{e}^{A\zeta }d\zeta E$ , $B_r^{•}=B_r^{°}\triangleq {\int }_0^T{e}^{A\zeta }d\zeta B_r$ .

These two subsystems are not only useful to describe the evolution of the impulsive system (4) at the sampling times kT and kT ⁺—which is necessary to implement an MPC—but also to characterize the equilibrium regions of Eq. (4) accounting for both the discontinuities and the free responses between them.

2.2 Extended equilibrium of the impulsive system

The idea now is to find an equilibrium set for Eq. (4), based on the underlying discrete-time subsystems (7), (8). To do that, we define first the equilibrium pairs of Eqs. (7), (8)— $(x_s^{•},u_s^{•})$ and $(x_s^{°},u_s^{°})$ , respectively—as the ones fulfilling the conditions ^b :

(9) $\begin{array}{l}\hfill x_s^{•}=A^{•}{x}_s^{•}+B_u^{•}{u}_s^{•}+E^{•},\\ \hfill x_s^{°}=A^{°}{x}_s^{°}+B_u^{°}{u}_s^{°}+E^{°}.\end{array}$

Given that u°(k + 1) = u ^•(k) by definition (and E° = E ^•), we have that at steady state it is $u_s^{•}=u_s^{°}$ . This means that one equilibrium input (defined as $u_s^{•}(k)$ for simplicity) has to fulfill both equilibrium equations, that is, the last equation can be written as

(10) $\begin{array}{l}\hfill x_s^{°}=A^{°}{x}_s^{°}+B_u^{°}{u}_s^{•}+E^{•}.\end{array}$

This way, given that Eqs. (7), (8) are subsystems describing one single impulsive system, we need to find triplets $(x_s^{•},x_s^{°},u_s^{•})\in \mathcal{X}\times \mathcal{X}\times \mathcal{U}$ fulfilling both Eqs. (9), (10). Furthermore, a final condition that the equilibrium triplet $(x_s^{•},x_s^{°},u_s^{•})$ must fulfill to be a feasible extended equilibrium is that the free responses corresponding to them (orbits) must be feasible:

$\begin{array}{l}\hfill o_s(x_s^{•},u_s^{•})\in \mathcal{X},\end{array}$

where

$\begin{array}{l}\hfill o_s(x_s^{•},u_s^{•})\triangleq \left\{\phi (t;x_s^{•},u_s^{•},0),t\in (kT,(k+1)T],k\in \mathbb{N}\right\},\end{array}$

and $\phi (t;x_s^{•},u_s^{•},0)=e^{At}(x_s+B_u{u}_s^{•})+{\int }_{kT}^t{e}^{A(t-\zeta )}d\zeta E$ (i.e., φ(⋅) it is the solution of the system (3), for an impulsive input of value ( $u_s^{•}$ ) and considering no disturbances).

Summarizing, we can characterize the extended equilibrium set, ${\mathcal{X}}_s^{•}$ , as ^c :

(11) $\begin{array}{l}\hfill {\mathcal{X}}_s^{•}\triangleq \left\{x_s^{•}\in \mathcal{X}:\exists u_s^{•}\in \mathcal{U}\text{such that}{A}^{•}{x}_s^{•}+B_u^{•}{u}_s^{•}+E^{•}=x_s^{•},o_s(x_s^{•},u_s^{•})\in \mathcal{X}\right\}\end{array}$

Furthermore, the set of equilibrium inputs—that is, the set of inputs for which an equilibrium state exists—is denoted as ${\mathcal{U}}_s^{•}$ . Note that this latter impulsive equilibrium is different from the input equilibrium defined in Section 4.2. For details in the procedure to compute such an equilibrium sets ${\mathcal{X}}_s^{•}$ and ${\mathcal{U}}_s^{•}$ , in a general context, see Rivadeneira et al. (2018).

2.3 Pulse input scheme

Here, a pulse input scheme (i.e., when the inputs are pulses of a given duration, instead of impulses) is introduced to be consistent with the real application of the control action that an actuator can supply in some biomedical cases, in which impulsive inputs are not allowed. See Fig. 3 for a schematic plot of the input pulses and state responses. To properly represents model (1) under the pulse scheme, some modifications need to be made on the underlying subsystems of Section 2.1. Specifically, matrix $B_u^{•}$ must be redefined as

Fig. 3. Pulse scheme.

$\begin{array}{l}\hfill B_u^{•}\triangleq e^{A(T-\mathrm{\Delta }T)}{\int }_0^{\mathrm{\Delta }T}{e}^{A(\mathrm{\Delta }T-\zeta )}d\zeta B_u,\end{array}$

where ΔT represents now the input pulse duration. Furthermore, to characterize the extended equilibrium, matrix $B_u^{°}$ takes the form

$\begin{array}{l}\hfill B_u^{°}\triangleq {\int }_0^{\mathrm{\Delta }T}{e}^{A(\mathrm{\Delta }T-\zeta )}d\zeta B_u.\end{array}$

This way, the extended equilibrium sets ${\mathcal{U}}_s^{•}$ and ${\mathcal{X}}_s^{•}$ are different from the ones obtained for the impulsive case. Clearly, the latter representation is a generalization of the former, in such a way that the impulsive representation is a particular case of the pulsive one, when the pulse duration tends to zero.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128174616000032

greengrecelf.blogspot.com

Source: https://www.sciencedirect.com/topics/engineering/continuous-time-system

Share this post