Sistemi Multivariabile

Jan 31, 2024

Introduzione

La notazione di questo paragrafo differisce da quello precedenti. L’operatore di derivazione viene indicato con uno (o n) punti sopra la funzione che viene derivata:

\dot{x}(t) \equiv \frac{d}{dt} x(t) \quad ; \quad \ddot{x}(t) \equiv \frac{d^2}{dt^2} x(t)\quad ; \quad\dots

Le variabili misurabili sono le caratteristiche di un sistema che si possono esprimere in modo numerico. Il modello matematico di un sistema è l’insieme di relazioni che legano le variabili misurabili del suddetto sistema.

Classificazione dei sistemi

Sistema Orientato

Un sistema orientato distingue le proprie variabili misurabili in input (cause) ed output (effetti).

\begin{matrix}\xrightarrow{u_1} \\ \dashrightarrow \\ \xrightarrow{u_r}\end{matrix} \;\boxed{\eq{\\\quad&\text{sistema}\quad \\ &\text{orientato}\\&}}\; \begin{matrix}\xrightarrow{y_1} \\ \dashrightarrow \\ \xrightarrow{y_m}\end{matrix}

Prendendo in prestito un paragrafo da Segnali Fondamentali, si può ricordare che:

\begin{equation} \begin{split} x : & \quad A \longrightarrow B \\ & \quad t \longmapsto x(t) \end{split} \end{equation}

Il segnale $x$ è:

tempo-continuo se $t \in \R$ , ovvero se $A \equiv \R$ ,
tempo-discreto se $t \in \Z$ , ovvero se $A \equiv \Z$ .

Inoltre, in base al codominio B, il segnale $x$ è:

a valori reali (continui) se $B \sube \R$
a valori discreti se $B = \{b_1,\dots,b_M\}\sub \Z$

Un segnale analogico è tempo-continuo e a valori continui (reali). Un segnale digitale è tempo-discreto e a valori discreti (interi). Un segnale tempo-discreto, per convenzione, si scrive:

\begin{equation} \begin{split} x : & \quad A \longrightarrow B \\ & \quad n \longmapsto x[n] \quad \forall n \in \Z \end{split} \end{equation}

Sistema Discreto

In un sistema discreto, i segnali sono associati ad una successione di numeri interi $t_k$ detti passi (step). I passi rappresentano istanti temporali multipli di un periodo di campionamento fissato a priori: $t_k = k T_s$ con $k \in \Z$ .

\xrightarrow{\;u(t_k)\;} \boxed{\text{ sistema orientato }} \xrightarrow{\;y(t_k)\;}

Di seguito si analizzano diversi aspetti che differenziano i sistemi multi-variabile.

Sistema Statico

Un sistema statico è detto privo di memoria. Un resistore elettrico è un sistema statico. Dato il valore di resistenza R e la funzione $u(t)$ associata alla corrente che vi scorre attraverso, la funzione $y(t)$ associata alla tensione ai capi del resistore non tiene conto dei valori passati della corrente:

y(t) = R \cdot u(u)

Sistema Dinamico

Il suo opposto è il sistema dinamico, il quale è detto dotato di memoria: l’output ad un dato istante dipende anche dall’input negli istanti precedenti a quello considerato. Bipoli elettrici quali il condensatore e l’induttore sono sistemi dinamici. Preso un condensatore con valore di capacità elettrica C, siano $u(t)$ la corrente che scorre attraverso il bipolo e $y(t)$ la tensione ai capi dello stesso, allora:

u(t)=C\dot{y}(t) \implies y(t)=y(t_0)+\frac{1}{C}\int_{t_0}^t u(t)dt\quad ;\;t_0 \leq t

Preso invece un induttore con induttanza L, si può scrivere:

u(t)=u(t_0)+\frac{1}{L}\int_{t_0}^t y(t)dt \implies y(t) = L \dot{u}(t) \quad ;\;t_0 \leq t

Lo stato di un sistema dinamico è la sua memoria ed è costituito da un vettore di dimensione n: ordine del sistema. Al vettore degli stati è associato lo spazio chiamato spazio degli stati. Lo spazio degli stati è $x(t)$ :

\begin{bmatrix}x_1\\ x_2 \\ \vdots \\ x_n \end{bmatrix} = x(t) \in \R^n

Nei modelli matematici dei circuiti elettrici, le variabili di stato indicano sempre la presenza di elementi che immagazzinano energia (ovvero condensatori e induttori). I resistori sono esclusi poiché dissipano potenza.

Un sistema dinamico nello spazio degli stati è definito dai seguenti insiemi:

T: insieme dei tempi
U: insieme degli ingressi
$\bold{U_f}$ : insieme delle funzioni di ingresso
X: insieme degli stati
Y: insieme delle uscite

All’istante $t\in \bold{T}=\R$ , nel campo dei modelli differenziali vettoriali:

x(t) = \begin{bmatrix}x_1(t) \\ x_2(t) \\ \vdots \\ x_n(t)\end{bmatrix} ;\; u(t) = \begin{bmatrix}u_1(t) \\ u_2(t) \\ \vdots \\ u_r(t)\end{bmatrix} ;\; y(t) = \begin{bmatrix}y_1(t) \\ y_2(t) \\ \vdots \\ y_m(t)\end{bmatrix}

Si ottengono i vettori di:

Stato: dimensione n, vettore $x(t) \in \bold{X}$ con $\bold{X} = \R^n$
Ingresso: dimensione r, vettore $u(t) \in \bold{U}$ con $\bold{U} = \R^r$
Uscita: dimensione m, vettore $y(t) \in \bold{Y}$ con $\bold{Y} = \R^m$

dato l’istante $t \in \bold{T}=\R$ .

Un evento è una coppia tempo-stato:

\big(t,\,x(t)\big) \in \bold{T} \times \bold{X}

La funzione transizione è:

\phi\big(t,t_0,x(t_0),u(\cdot)\big)

Il moto è l’insieme degli eventi definiti dalla funzione transizione per $t \in [t_0,\,t_1]$ .

La traiettoria è l’immagine di X della funzione transizione per $t \in [t_0,\,t_1]$ .

Lo stato di equilibrio temporaneo è lo stato $x\in X$ per il quale:

\exist u(\cdot) \in U_f \text{ che soddisfa } x=\phi(t,t_0,x,u(\cdot)) \quad \forall t \in [t_0,\,t_1]

Lo stato di equilibrio è lo stato di equilibrio temporaneo in $[t_0,\,t_1] \; \forall t_0,\,t_1 \in T$ con $t_0 < t_1$

Classificazione dei modelli

Lineare e Non-Stazionario

\sis{ \dot{x}(t) = A(t)x(t)+B(t)u(t) \\ y(t)= C(t)x(t)+D(t)u(t) }

Con $A(t)$ , $B(t)$ , $C(t)$ e $D(t)$ continue a tratti.

Lineare e Stazionario (LTI)

\sis{ \dot{x}(t) = Ax(t)+Bu(t) \\ y(t)= Cx(t)+Du(t) }

Non-Lineare e Non-Stazionario

\sis{ \dot{x}(t) = f\big(x(t),u(t),t\big) \\ y(t)=g\big(x(t),u(t),t\big) }

Non-Lineare e Stazionario

\sis{ \dot{x}(t) = f\big(x(t),u(t)\big) \\ y(t)=g\big(x(t),u(t)\big) }

Matrici caratteristiche

Tutti i sistemi descritti di seguito sono dinamici.

Sistema Lineare Non-stazionario

Dato un sistema lineare tempo continuo, $t \in \R$ :

\sis{ \dot{x}(t) = A(t)x(t)+B(t)u(t) \\ y(t)= C(t)x(t)+D(t)u(t) }

Oppure dato un sistema lineare tempo discreto, $k \in \Z$ :

\sis{ x(k+1)=A(k)x(k)+B(k)u(k)\\ y(k) = C(k)x(k)+D(k)u(k) }

Allora, preso il caso tempo continuo, $x(t)$ è un vettore di dimensione n. $A(t)$ è una matrice quadrata $n\times n$ :

A(t) = \begin{bmatrix}a_{11}(t) & \dots & a_{1n}(t) \\ \vdots & \cdots & \vdots \\ a_{n1}(t) & \dots & a_{nn}(t)\end{bmatrix}

$u(t)$ è un vettore di dimensione r. $B(t)$ è una matrice di dimensione $n\times r$ :

B(t) = \begin{bmatrix}b_{11}(t) & \dots & b_{1r}(t) \\ \vdots & \cdots & \vdots \\ b_{n1}(t) & \dots & b_{nr}(t)\end{bmatrix}

$y(t)$ è un vettore di dimensione m. $C(t)$ una matrice di dimensione $m\times n$ :

C(t) = \begin{bmatrix}c_{11}(t) & \dots & c_{1n}(t) \\ \vdots & \cdots & \vdots \\ c_{m1}(t) & \dots & c_{mn}(t)\end{bmatrix}

$u(t)$ è un vettore di dimensione r. $D(t)$ è una matrice di dimensione $m \times r$ :

D(t) = \begin{bmatrix}d_{11}(t) & \dots & d_{1r}(t) \\ \vdots & \cdots & \vdots \\ d_{m1}(t) & \dots & d_{mr}(t)\end{bmatrix}

Sistema Lineare Stazionario

Quando il sistema è sia lineare che tempo-invariante, ovvero LTI, le matrici caratteristiche diventano costanti:

\eq{ & \sis{ \dot{x}(t) = Ax(t)+Bu(t) \\ y(t)= Cx(t)+Du(t) } \\ \\ &\sis{ x(k+1)=Ax(k)+Bu(k)\\ y(k) = Cx(k)+Du(k) } }

Le matrici hanno ancora la stessa dimensione:

A = \begin{bmatrix}a_{11} & \dots & a_{1n} \\ \vdots & \cdots & \vdots \\ a_{n1} & \dots & a_{nn}\end{bmatrix} \\ \\ B = \begin{bmatrix}b_{11} & \dots & b_{1r} \\ \vdots & \cdots & \vdots \\ b_{n1} & \dots & b_{nr}\end{bmatrix} \\ \\ C = \begin{bmatrix}c_{11} & \dots & c_{1n} \\ \vdots & \cdots & \vdots \\ c_{m1} & \dots & c_{mn}\end{bmatrix} \\ \\ D = \begin{bmatrix}d_{11} & \dots & d_{1r}\\ \vdots & \cdots & \vdots \\ d_{m1} & \dots & d_{mr}\end{bmatrix}

Sistema Multi-Input Multi-Output

Un sistema SISO, ovvero Single-Input Single-Output ha i coefficienti m ed r uguali ad uno:

m = 1 \quad ; \quad r = 1

Un sistema MIMO, ovvero Multi-Input Multi-Output (o multivariabile) ha i coefficienti m ed r maggiori di uno:

m > 1 \quad ; \quad r > 1

Caratterizzazione LTI MIMO

La caratterizzazione di un sistema Lineare Tempo-Invariato Multi-Input Multi-Output, ovvero LTI MIMO, si differenzia sulla base della natura del tempo. Le funzioni $f$ e $g$ sono lineari rispetto ad $x$ e $u$ per ogni $t$ , sono esprimibili come prodotto tra matrici e vettori:

Il sistema tempo continuo è:

\sis{ \dot{x}(t) = Ax(t)+Bu(t) \\ y(t)= Cx(t)+Du(t) }

Il sistema tempo discreto è:

\sis{ x(k+1) = Ax(k)+Bu(k) \\ y(k)= Cx(k)+Du(k) }

Risposta di un sistema dinamico

Sostituendo la funzione di transizione dello stato $\phi()$ nella funzione di uscita, si ottiene la funzione di risposta:

\eq{y(t)&=g\Big(\phi\big(t,t_0,x(t_0),u(\cdot)\big),u(t),t \Big)=\\ &=\gamma\big(t,t_0,x(t_0),u(\cdot)\big)}

La $\dot{x}(t)$ è lo stato iniziale ed è la derivata 1° della funzione di transizione dello stato:

x(t)=\phi\big(t,t_0,x(t_0),u(\cdot)\big)

Si passa dal modello differenziale ingresso-stato-uscita al modello in cui sono fissate le condizioni iniziali su $x(t_0)$ .

Matrice di transizione

Dato un sistema LTI con input costantemente nullo $u(t)=0$ , allora:

\dot{x}(t) = Ax(t); \quad x(0)=x_0,\; x(t)\in \R^n

La soluzione dell’equazione differenziale vettoriale è:

x(t)=e^{At}x_0

Si definisce la funzione esponenziale della matrice A, detta anche matrice di transizione:

\eq{ & e^{At} = \sum_{i=0}^\infty A^i \frac{t^i}{i!} \\ \implies & \frac{d}{dt} e^{At} = Ae^{At}\;;\quad e^{At}_{t=0} = I }

I è la matrice quadrata identità.

x(t)=\phi(t,t_0,x(t_0),0)=e^{A(t-t_0)}x(t_0)

Trasformata della Matrice di Transizione

Dato lo stesso sistema LTI con input costantemente nullo, in forma matriciale:

\dot{X}(t) = AX(t); \quad X(t_0)=I

L’esponenziale di matrice rappresenta la matrice di transizione di un sistema dinamico, nonché la soluzione dell’equazione differenziale matriciale sopra descritta.

Applicando la trasformazione funzione di Laplace alla matrice di tansizione, si ottiene la matrice complessa $X(s)$ . Grazie al teorema della trasformata di derivata, si ottiene:

\eq{ & \mathcal{L}\Big[e^{At}\Big] = X(s)\\ \implies & AX(s)= sX(s)-I\\ \implies & I = X(s)\Big(sI-A\Big)\\ \implies & \mathcal{L}\Big[e^{At}\Big]=(sI-A)^{-1} }

Ogni elemento di $(sI-A)^{-1}$ è un rapporto tra due polinomi: a denominatore si trova il polinomio caratteristico di A (grado n), a numeratore ha grado massimo $n-1$ .

(sI-A)^{-1} = \frac{\text{agg}(sI-A)}{\text{det}(sI-A)}

Matrice di Trasferimento

Dato il modello nello spazio degli stati precedentemente descritto, con condizioni iniziali fissate:

\sis{ \dot{x}(t) = Ax(t)+Bu(t) \\ y(t)= Cx(t)+Du(t) }\quad ; \quad x(0)=x_0

Applicando la trasformazione di Laplace, si ottengono le funzioni:

\eq{ & \mathcal{L}[x(t)]=X(s)\;;\\ & \mathcal{L}[y(t)]=Y(s)\;;\\ & \mathcal{L}[u(t)]=U(s)\;; }

Il sistema si può riscrivere come:

\eq{ & \sis{ sX(s)-x_0 = AX(s)+BU(s) \\ Y(s) = CX(s)+DU(s) }\\ \implies & X(s) = (sI-A)^{-1}x_0 + (sI-A)^{-1}BU(s)\\ \implies & Y(s) = C(sI-A)^{-1}x_0 + \Big[C(sI-A)^{-1}B+D\Big]U(s) }

Si definisce la matrice di trasferimento del sistema:

G(s) = C(sI-A)^{-1}B+D

Caratterizza il comportamento ingresso-uscita del sistema in condizione di quiete, ovvero quando lo stato iniziale è nullo: $x_0 \equiv x(0)=0$ .

Y(s)=G(s)U(s)

La matrice di trasferimento è inoltre univoca e indipendente dalla scelta di variabili di stato.

Le proprietà che seguono sono valide se il modello nello spazio degli stati è completamente raggiungibile-controllabile e completamente osservabile-ricostruibile.

Ogni elemento della matrice di trasferimento è il rapporto di due polinomi in s: a denominatore si trova il polinomio minimo di A. Le radici del denominatore sono dette poli del sistema e corrispondono agli autovalori della matrice A. I polinomi al numeratore, negli elementi di $G(s)$ , hanno grado inferiore a quello del polinomio minimo di A se il sistema è puramente dinamico ( $D = 0$ ), altrimenti hanno grado minore o uguale a quello del polinomio minimo di A se il sistema non è puramente dinamico.

In un sistema SISO, la matrice di trasferimento diventa una funzione di trasferimento: funzione scalare, rapporto di polinomi. Le radici del numeratore di tale funzione sono dette zeri del sistema.

\Large \xrightarrow{\;u\;}\boxed{G(s)}\xrightarrow{\;y\;}

In un sistema MIMO con r ingressi e m uscite, la matrice di trasferimento (dimensione $m \times r$ ) contiene $m \cdot r$ differenti funzioni di trasferimento SISO.

Ad esempio: $m=r=2$

\underbrace{\begin{bmatrix}y_{1}(s) \\ y_{2}(s)\end{bmatrix}}_{Y(s)} = \underbrace{\begin{bmatrix}G_{11}(s) & G_{12}(s) \\ G_{21}(s) & G_{22}(s)\end{bmatrix}}_{G(s)} \cdot \underbrace{\begin{bmatrix}u_{1}(s) \\ u_{2}(s)\end{bmatrix}}_{U(s)}

Da cui si ricavano $m \cdot r$ sistemi a blocchi:

\Large \eq{ & \xrightarrow{\;u_1(s)\;}\boxed{G_{11}(s)}\xrightarrow{\;y_{11}(s)\;} \\ & \xrightarrow{\;u_2(s)\;}\boxed{G_{12}(s)}\xrightarrow{\;y_{12}(s)\;} \\ & \xrightarrow{\;u_1(s)\;}\boxed{G_{21}(s)}\xrightarrow{\;y_{21}(s)\;} \\ & \xrightarrow{\;u_2(s)\;}\boxed{G_{22}(s)}\xrightarrow{\;y_{22}(s)\;} \\ \normalsize \implies & \normalsize y_1(s) = y_{11}(s) + y_{12}(s)\\ \normalsize \implies & \normalsize y_2(s) = y_{21}(s) + y_{22}(s)\\ }

Espressione diretta

I modelli descrittivi dei sistemi LTI possono anche essere espressi con equazioni differenziali di ordine n, del tipo:

\sum_{j=0}^n a_j \frac{d^j y}{dt^j} = \sum_{j=0}^m b_j \frac{d^j u}{dt^j} \\ \implies a_0 y(t) + \dots + a_n \frac{d^n y}{dt^n} = b_0 u(t) + \dots + b_m \frac{d^m u}{dt^m}

L’indice n è il massimo ordine della derivata della funzione di uscita $y(t)$ . L’indice m è il massimo ordine della derivata dell’ingresso $u(t)$ .

Se il modello differenziale è espresso rispetto ad $u$ e $y$ , tolte le variabili di stato e in condizioni iniziali nulle, la funzione di trasferimento si ottiene direttamente con la trasformazione di Laplace ed il corrispondente teorema della derivata. Ogni termine si trasforma, grazie a Laplace, come segue:

\eq{ & \mathcal{L}[a_0 y(t)] = a_0 Y(s) \quad \dots \quad \mathcal{L} \bigg[a_n \frac{d^n y}{dt^n}\bigg]=a_n s^n Y(s)\\ & \mathcal{L}[b_0 u(t)] = b_0 U(s) \quad \dots \quad \mathcal{L} \bigg[b_m \frac{b^m u}{dt^m}\bigg]=b_m s^m U(s)\\ \implies & G(s) = \frac{Y(s)}{U(s)}=\frac{\sum_{j=0}^m b_j s^j}{\sum_{j=0}^n a_j s^j} =\frac{b_0 s^0 + \dots + b_m s^m}{a_0 s^0 + \dots + a_n s^n} }

Risposta impulsiva

La risposta impulsiva di un sistema LTI tempo-continuo puramente dinamico è data dalla matrice:

W(t) = Ce^{At}B

L’output è dato dall’integrale di convoluzione tra input e risposta impulsiva:

y(t) = \int_0^t W(t-\tau)\delta(\tau)d\tau = W(t)

Applicando la trasformazione funzionale di Laplace:

\mathcal{L}[W(t)] = C\mathcal{L}\Big[e^{At}\Big]B = C(sI-A)^{-1}B

Nei sistemi completamente raggiungibili-controllabili e completamente osservabili-ricostruibili, la trasformata della risposta impulsiva coincide con la matrice di trasferimento:

G(s) =\mathcal{L}[W(t)]

La risposta impulsiva di un sistema LTI tempo-discreto è data dalla matrice $W(t)=CA^{k-1}B$ . Si può determinare applicando input del tipo:

u(k) = \sis{ 1 \quad \textrm{ se } k = 0 \\ 0 \quad \textrm{ se } k \neq 0 }

Altre proprietà

Stabilità

N.B. Il simbolo $\iff$ significa “se e solo se”, abbreviabile con sse o iff.

Un sistema LTI tempo-continuo è:

semplicemente stabile $\iff$
1. tutti gli autovalori di A hanno parte reale $\leq 0$ ,
2. gli autovalori a parte reale nulla hanno molteplicità unitaria nel polinomio minimo di A
asintoticamente stabile $\iff$ tutti gli autovalori di A hanno parte reale negativa.

Un sistema LTI tempo-discreto è:

semplicemente stabile $\iff$
1. tutti gli autovalori di A hanno modulo $\leq 1$ ,
2. gli autovalori a modulo unitario hanno molteplicità unitaria nel polinomio minimo di A
asintoticamente stabile $\iff$ tutti gli autovalori di A hanno modulo $< 1$ .

Il sistema LTI ottenuto come cascata o come parallelo di due sistemi LTI è asintoticamente stabile sse lo sono entrambi i sistemi interconnessi.

Raggiungibilità e Controllabilità

Lo stato $x_1$ è raggiungibile da $x_0$ , ovvero $x_0$ è controllabile a $x_1$ , nell’intervallo $[t_0,\,t_1]$ con $t_0<t_1$ se:

\exist u(\cdot) \in U_f \quad \textrm{ t.c. } \quad x_1=x(t_1)=\phi(t_1,t_0,x(t_0),u(\cdot))

Un sistema LTI MIMO tempo-continuo è completamente raggiungibile se qualunque stato può essere raggiunto da $x=0$ in un tempo finito. Un sistema LTI MIMO tempo-discreto è completamente controllabile se $x=0$ può essere raggiunto da qualunque stato in un tempo finito.

Si definisce matrice di raggiungibilità di un sistema LTI la matrice $P = [B \;\; AB \;\; A^2B \;\; \dots \;\; A^{n-1}B]$ di dimensione $n \times (n\cdot r)$ . Il sistema è dunque completamente raggiungibile e completamente controllabile $\iff \textrm{rango}(P)=n$ .

In MATLAB, la matrice P, dati A e B, si calcola con il comando ctrb:

A = [1  1; 4 -2]
B = [1 -1; 1 -1]
P = ctrb(A,B) % controllability matrix [n x (n r)]
control_rank = rank(P);
n = size(P,1);

fprintf("»»» rank of P = %d\n", control_rank);
fprintf("»»» n of P = %d\n", n);

if n == control_rank
    fprintf("»»» P is fully reachabile and fully controllable\n");
else
    fprintf("»»» P is NOT fully reachabile and fully controllable\n");
end

Osservabilità e Ricostruibilità

Proprietà che descrivono la possibilità di determinare lo stato iniziale $x(t_0)$ e finale $x(t_1)$ avendo noto il comportamento ingresso-uscita $u[t_0,\,t_1],\;y[t_0,\,t_1]$ .

Un sistema è completamente osservabile in $[t_0,\,t_1]$ se la conoscenza del comportamento ingresso-uscita $u[t_0,\,t_1],\;y[t_0,\,t_1]$ permette di determinare univocamente lo stato iniziale $x(t_0) \; \forall u[t_0,\,t_1]$ .

Un sistema è completamente ricostruibile in $[t_0,\,t_1]$ se la conoscenza del comportamento ingresso-uscita $u[t_0,\,t_1],\;y[t_0,\,t_1]$ permette di determinare univocamente lo stato iniziale $x(t_1) \; \forall u[t_0,\,t_1]$ .

Per sistemi LTI si può considerare anche solo $t=t_1 - t_0$ .

Completa osservabilità $\implies$ completa ricostruibilità.

Si definisce matrice di osservabilità di un sistema LTI la matrice:

Q = [C \;\; CA \;\; CA^2 \;\; \dots \;\; CA^{n-1}]^T \equiv Q^T = [C^T \;\; A^T C^T \;\; (A^T)^2 C^T \;\; \dots \;\; (A^T)^{n-1}C^T]

Tale matrice ha dimensione $n \times (n\cdot m)$ . Il sistema è dunque completamente osservabile e completamente ricostruibile $\iff \textrm{rango}(Q^T)=n$ .

Si osservi che per Q è stata introdotto la forma trasposta per poter permettere un’analogia con la matrice P.

In MATLAB, la matrice Q^T, dati A e C, si calcola con il comando obsv:

A = [1  1; 4 -2]
C = [-4 1; 1 -1]
Q = obsv(A,C) % observability matrix [n x (n m)]
obsv_rank = rank(Q);
n = size(Q,1);

fprintf("»»» rank of Q = %d\n", obsv_rank);
fprintf("»»» n of Q = %d\n", n);

if n == obsv_rank
    fprintf("»»» Q is fully observable and fully reconstructable\n");
else
    fprintf("»»» Q is NOT fully observable and fully reconstructable\n");
end