COMPTE RENDU DU TP D’AUTOMATIQUE DES SYSTEMES LINEAIRES MULTIVARIABLES

PARTIE I : D’UN PORTIQUE DE STRANSPORT

Dans cette première partie, nous allons réaliser la commande d’un système de transport utilisé pour
le déchargement des marchandises dans un port. Il est représenté a la figure 1.

Son principe de fonctionnement est le suivant :

 Un chariot mobile de masse M repéré par son abscisse xl et se déplaçant sur le portique sous
l'action d'une force f(t) développée par un moteur électrique. Le chariot est soumis à des
forces de frottement visqueux ;
 Une charge de masse m, suspendue à l'extrémité d'un câble de longueur l. La position de la
charge est repérée par l'angle  d'inclinaison du câble par rapport à la verticale.

Notre objectif tout au long de cette partie sera de réguler la distance d fait par le chariot ainsi que
l’angle θ que fait le câble liant le chariot et la charge par rapport à la verticale, le tout en prenant
comme consigne, la force f(t) appliqué au chariot.

Notre système possède alors une entrée (f(t)) et deux sortie (d et θ)

1) Modélisation et étude du principe de transport

( M +m ) . d̈ ( t )+ m. l. θ̈ ( t )+ α . ḋ ( t ) =f (t )
l . θ̈ ( t ) + d̈ ( t ) =−g .θ (t)
On fait les choix suivants :
 []
d (t)
ḋ (t)
 Le vecteur d’état : x (t )=
θ(t)
θ̇(t)

 Le vecteur de sortie est le suivant : y(t)=

[ ]
d (t)
θ(t)
1-1)Expliquons le choix du vecteur d’état :

Nous savons que le choix du vecteur d’état n’est pas unique et dépend généralement de plusieurs
paramètres tels que la disponibilité de leurs capteurs, la structure du système (observabilité et
commandabilité) et même leur pertinence. Dans notre cas, le vecteur d’état contient les deux sortie
(choix logique) mais aussi leurs dérivées premières car il sera judicieux de connaitre la valeur de la
vitesse linéaire de déplacement du chariot et la vitesse de rotation de la charge afin de synthétiser
une commande qui prend en compte non seulement la position du chariot et de la charge, mais aussi
leurs vitesses de déplacement.

1-2) Ecrivons les équations d’état lié à ce système :

Nous savons que la représentation d’état système est généralement de la dorme suivante :

{ ẋ ( t )= A . x ( t ) + B .u (t)
y ( t )=C . x ( t )+ D . u(t )

D’après les équations 1 et 2, on peut déduire la représentation d’état suivante pour notre
système :

[ ][ ] [ ][ ]
0 1 0 0 0
ḋ (t) −a m d (t ) 1
0 ×g O
d̈ (t) M M ḋ (t )
= × + M × f (t)
θ̇(t) 0 0 0 1 θ (t) 0
θ̈(t) 0 a −( M +m) θ̇ (t) −1
×g 0
M ×l M ×l M ×l

[ ][]
d (t )

[ ][
d (t) 1 0 0 0
=
θ(t) 0 0 1 0
×
ḋ (t ) 0
θ(t) 0 ]
+ × f ( t)

θ̇(t)

Apres application numérique, on a :

 {
[ ][ ][ ][ ]
ḋ (t) 0 1 0 0 d (t) 0
d̈ (t) 0 −0.0400 1 0 ḋ (t) 10−4
= × + × f (t)
θ̇ (t) 0 0 0 1 θ(t ) 0
−4
θ̈ (t) 0 0.0011 −0.3143 0 θ̇(t ) −0.0286 ×10

[ ][]
d (t )

[ ][
d (t) 1 0 0 0
θ(t )
=
0 0 1 0
× ]
ḋ (t ) 0
θ (t) 0
+ × f (t )

θ̇ (t)

On constate par identification que :

[ ] [ ]
0 1 0 0 0
A=
0 −0.0400
0 0
1
0
0
1
; B=
10−4
0
;C= [
1 0 0 0
0 0 1 0
et D =
0
0 ] []
0 0.0011 −0.3143 0 −0.0286 ×10
−4

NB : notre système est dit propre car D est une matrice nulle.

2) Analyse du système en boucle ouverte

Avant toute étude, nous allons tout d’abord définir notre système (syst_detat) sous Matlab grâce à la
fonction « ss » qui prend en argument les matrice A,B,C et D.

Le résultat est le suivant :

2-1) Etudions :

a) La stabilité

Un système continu est par définition stable si les valeurs propre de sa matrice d’état (A) sont a
partie réelles positives.

Nous allons calculer ces modes propres en utilisant la commande «pole(Syst_detat) » sous Matlab.

Le résultat est le suivant :

Notre système possède 4 pôles (Logique vue qu’on a 4 variables d’état). Si les trois dernies pôles sont
à partie réels négatifs, le premier pôle est nul (à la limite de la stabilité). On peut donc conclure que
notre système est à la limite de la stabilité.

b) La commandabilité
Selon Kalman, un système est commandable si la matrice de commandabilité (Com=
( B AB A 2 B … .. An− A B )) est de rang plein et égal a n (l’ordre du système).
Nous allons calculer sous Matlab la matrice de commandabilité grâce à la fonction « ctrb » et on
interrogera son rang grâce à la commande « rank ».

Le résultat est le suivant :

Notre système bel et bien commandable.

c) L’observabilité

Selon Kalman, un système est observable si la matrice d’observabilité (Cob=

T
( C CA CA 2 ….. CA n− A ) ) est de rang plein et égal a n (l’ordre du système).

Nous allons calculer sous Matlab la matrice d’observabilité grâce à la fonction « obsv» et on
interrogera son rang grâce à la commande « rank ».

Le résultat est le suivant :

Notre système est observable.

[]
0
0
2-2) Traçons la réponse la réponse de notre système si : f(t)=0 ; et que x ( 0 )=
0.1
0

Il s’agit ici de déterminer la réponse en régime libre.

Pour cela, on utilisera la fonction « initial » qui prend en argument le système et le vecteur d’état a
l’instant t= t0.

Le résultat est le suivant :

Si on réalise un zoom vers les premières secondes, on constate qu’a t=0s, d=0 m et θ=0.1 rad.

2-3) représentons notre système sous Simulink.

a) En utilisant le block State space

Le schéma est le suivant :
 La réponse du système ainsi que la configuration du bloc state space sont donné dans la
figure suivante :

L’on constate que le système converge mais avec des oscillations, il est a la limite de la stabilité.

b) En utilisant le bloc à base d’intégrale

Le schéma issu des équations d’état est le suivant :

La réponse du système est la suivante :

Cette courbe est identique a celle obtenu par le bloc state space et démontre toujours que notre
système est à la limite de la stabilité.

3) Commande par retour d’état.

L’on souhaite améliorer les performances de ce système (statique et dynamique). Pour cela, on
applique un rebouclage afin de lui appliquer une commande du style f(t)=v(t)-Lx(t) avec L la matrice
de contre réaction d’état. On aimera que le système bouclé ait les valeurs propres suivantes :

{
−2+2 j
−2−2 j
−5
−10

3-1) En utilisant le passage par une forme commandable, calculons L

Cette méthode de calcul de L passe par le calcul d’une matrice M tel que x (t )=M . ^x ( t )

Avec : ^x (t ) le vecteur d’état dans la forme commandable et M la matrice de valeur propre

permettant le passage de la forme quelconque vers la zone commandable.

L’équation d’état sous la forme commandable s’écrit alors :

x^˙ ( t )=~
A . ^x ( t ) + ~
{ ~ ~
B .u (t)
y ( t )=C . x^ ( t )+ D . u(t )

Avec ~ −1 ~ −1 ~ ~
A=M . A . M ; B=M . B ; C=C . M et D=D
L’on calculera alors une peudo-matrice de contre réaction dans l’espace commandable appelé L com et
grâce a la matrice M, l’on retrouvera L la matrice de contre réaction d’état dans la forme quelconque
du système.

Avec Lcom=[l com 3 l com 2 l com 1 l com 1 ]

L’on comprend donc que le préalable a la détermination de L est le calcul de la matrice M.

M est généralement une matrice carrée de nxn, dans notre cas notre matrice sera donc 4x4.

M =[ m 4 m 3 m2 m1 ]
m4 =B

Et m n− j =A . mn− j+1 +a n− j . B j Avec j∈[1; n−1]

i '
Avec a i≤coéfficient lié à P dans l équation caractéristique de système en boucle ouverte

Les calculs de M et M-1 ce sont fait sous Matlab par le code suivant :

Le résultat est le suivant :

Les coefficients de la matrices Lcom se calculent de la manière suivante :

l com i=β i−ai

i '
Avec β i≤coéfficient lié à P dans l équation caractéristique de système en boucle fermé
Le code permettant le calcul de la matrice Lcom est le suivant :
Le résultat est le suivant :

Lcom= [ 18.96 117.6857 319.9886 4000 ]

Le calcul de L se fait de la manière suivante : L = M-1.Lcom

Tout calcul fait sous matlab, on a :

La matrice de contre réaction L nécessaire pour atteindre nos objectifs est donc :

L=10 . [ 0.14 0.112 4.4881 3.8535 ]

8

3-2) Vérifions la véracité de notre calcul de la question 3-1) en utilisant la fonction « Place » sous
Matlab.

Place est une fonction de Matlab qui calcul la matrice de contre réaction d’état et prend en argument
le système d’état en boucle ouverte et le vecteur contenant les pôles souhaité en boucle fermé.

Le résultat est le suivant :

Lplace=10 . [ 0.14 0.112 4.4881 3.8535 ]

8

On constate donc que les résultats sont identiques.

3-3) Traçons la réponse de notre système bouclé avec des conditions initiales nulles.

Pour cela, nous allons définir un nouveau système dans l’espace d’état appelé Syst_detat_bf qui
représente notre système bouclé.
Pour rappel, l’équation d’état du système bouclé est de la forme suivante :

{ ẋ ( t ) =( A−B . L ) . x ( t ) +B . K . v (t)
y ( t )= (C−D . L ) . x ( t ) + D . K . v (t)

Notre système étant propre, on a :

{ ẋ ( t )=( A−B . L ) . x ( t ) + B . K . v (t)

y ( t )=C . x ( t )

Nous supposerons que K=1 pour notre simulation, le code Matlab est le suivant :

Le résultat est le suivant :

L’on peut constater que :

 La valeur finale de d est 7.12x10-8 m

3-4) Déterminons K pour avoir une érreur statique nulle sur la position du chariot.

f (+∞) 1
K= =
d (+ ∞) 7.12 x 10−8
3-5) Pour cette nouvelle valeur de K, traçons la réponse du système avec Simulink.
L’on peut constater sur ces courbes que :

 La dynamique du système n’a pas changé malgré le changement des valeurs de k.

 La distance d est à la valeur 1m en régime permanent ;
 La position θ est nulle à t=0 mais débute par un angle négatif car dès le démarrage du
chariot, la charge aura un léger retard, ce qui causera cet angle négatif.

4) Observateur d’état d’ordre plein.

Un observateur d’état en automatique a pour but de synthétiser un vecteur d’état identique au

vecteur d'état réel. Il s’agit d’un « capteur logiciel ». Plusieurs causes peuvent nous contraindre a
faire recours a un observateur pour la commande d’un système tels que (l’inexistence du capteur
nécessaire pour la mesure de la variable d’état, du bruitage des informations issus des capteurs ou
même l’inaccessibilité de la variable d’état à mesurer). Il existe plusieurs observateurs d’état. Dans ce
cas, nous allons réaliser la synthèse d’un observateur de type Luenberger d’ordre Plein (on va
procéder à la reconstruction totale du vecteur d’état du système). L’observateur d’ordre plein
s’oppose à l’observateur d’ordre réduit ou l’on ne reconstitue pas la totalité du vecteur d’état.

4-1) Choisissons les pôles de nôtres observateurs.

Nous avons choisi pour nôtres observateur les pôles suivants :

P1=−10 P2=−20 P3=−30 Et P4 =−40

Le choix de ses pôles a été fait avec les contraintes suivantes :

 L’observateur doit être stable ;

 Sa dynamique doit être meilleur que celle du système a commander en boucle ouverte. En
effet, la commande se faisant par retour de l’état reconstitué, il est évident de constater que
la commande sera erronée si la valeur du vecteur d’état l’est également. Pendant son régime
transitoire, le vecteur d’état reconstitué est différent du vecteur d’état réel, d’où la nécessité
de lui donner un temps de réponse court.

4-2) En utilisant la méthode d’Ackerman, déterminons la valeur de la matrice de gain de

l’observateur d’état.
 La méthode d’Ackerman est une méthode de synthèse des observateurs pour les systèmes
SISO qui se base sur :

 Le principe de système dual ;

 Le thèorème de Caley-Hamilton

Notre observateur possèdera une matrice de contreréaction K obs de dimension 4x2 (on a 4 variables
d’état et 2 sorties). Afin d’utiliser la méthode d’Ackerman, nous allons nous ramener à un cas à une
sortie s(t)=r(t).y(t) Avec r(t)=[r1 r2 r3….rp]

Le vecteur r est fixé par nos soins afin d’attribuer les mêmes importances aux sorties du système.

Dans notre cas, nous avons pris r =[ 1 10 ]

Ainsi, K obs =r.g avec g(t)=[g1 g2 g3….gn] à calculer en se ramenant a un cas monovariable du genre :

x^˙ (t )=( A−g . c 1 ) . ^x ( t ) + B .u ( t ) + g . s (t)

La formule d’Ackerman pour le calcul de la matrice de gains de l’observateur est la suivante :

K gain obs =[ 0 … .. 1 ] ¿ ¿ Avec :

T

 ¿ ¿ : inverse de la matrice de commandabilité du système dual.

 ϕ bf ( p ) : Équation caractéristique de l’observateur.

Le code de calcul de la matrice de gain de l’observateur d’état d’ordre pleint de notre système est le
suivant :

Le résultat est le suivant :

 [ ]
0.0379 0.3791
7 0.1098 1.0982
K obs =10 .
−0.0038 −0.0379
−0.0109 −0.1095

4-3) Simulons l’ensemble sous Matlab.

Il est a noté que :

 La synthèse de la commande et de l’observateur sont indépendant, par conséquent, notre

matrice de contre réaction d’état calculé en 3-1 est toujours valable.
 L’observateur est sous forme d’un système d’état bouclé avec la matrice Kobs

Le schéma de simulation est le suivant :

Le bloc Observateur de rang plein est rempli de la manière suivante :

Les courbes obtenues sont les suivantes :

Cette courbe démontre que notre observateur converge rapidement et est précis car les courbes
représentantes de chaque variable d’état réel du système se confond parfaitement avec celle
calculées par l’observateur.
Cette courbe démontre également l’exactitude de notre observateur car le vecteur de sortie du
système est identique au vecteur de sortie estimé.

 Problématique de l’observateur de Luenberger :

Cet observateur est certes bon, mais il est très vulnérable au bruit (réalité physique de tout système
physique). Nous allons dans l’ la simulation suivante rajouté du bruit de mesure (usure des capteurs
de sortie ou détérioration des filtres faisant partie du système de conditionnement des dits capteurs)
et du bruit de commande dut par exemple a la faible précision numérique de notre modèle.

Dans la simulation suivante, l’on a rajouté au système de la question 4-3 un bruit de mesure de
moyenne nulle et de variance 0.0001. Les courbes du vecteur d’état sont les suivantes :
L’on constate qu’il y’a une différence notoire entre le vecteur d’état réel et celui estimé. Ceci est la
limite de l’observateur de Luenberger, il est très peux robuste dès lors qu’il y a apparition du bruit.

5) Observateur de Kalman discret

L’observateur de Kalman permet de résoudre le problème de la robustesse vis-à-vis du bruit de

l’observateur de Luenberger. En effet, ce filtre permet de synthétiser un vecteur d’état fidèle a celui
du système réel même en présence de bruit de mesure et de commande. Ce choix se fait au travers
des matrice Q et R en fonction de la confiance qu’on accorde au model (R) et à la mesure (Q). Dans
notre cas, on réalisera un filtre de Kalman discret.

5-1)Discrétisons le système avec Te=0.01s

Nous réaliserons cette opération sous Matlab en utilisant la fonction «c2d » qui prend en argument
le système d’état continu et la période d’échantillonnage.

Le résultat est à la figure suivante :

 { [ ] [ ]
−3 −5 −7 −9
1 9.998 ×10 4.999× 10 1.666 × 10 4.999 ×10
−3
0 0.9996 9.998 × 10 4.999 ×10−5 9.998 ×10−7
x ( k + 1 )= −8 −2 × x (k )+ −10 × f (k )
0 5.714 ×10 1 10 −1.428 ×10
0 1.143 ×10−5 −3.143 × 10−3 1 −2.857 ×10−8

[ ][
d (t ) 1 0 0 0
θ(t)
=
0 0 1 0 ]
0
[]
× x (k)+ × f (k )
0

5-2)Calculons Kf la matrice de gain de l’observateur de Kalman.

 On a plus confiance en la mesure qu’au modèle

Le code permettant le calcul de Kf dans ce cas est :

Le résultat est le suivant :

[ ]
1
2.35×−10 0
0.9×−10 3 3.2×−103
k f 1=
0 2.35 ×−101
4
0 −10

Dans la simulation suivante, l’on va rajouter au système de la question 4-3 un bruit de mesure de
moyenne nulle et de variance 0.0001. L’observateur de Kalman et celui de Luenberger ont des
structures identiques, ainsi l’on modifiera uniquement les gains de l’observateur dans les schémas de
simulation.

Les courbes de simulation pour le cas 1 sont les suivantes :

 On a plus confiance au modèle

Le code permettant le calcul de Kf dans ce cas est :

PARTIE II : REGULATION DU NIVEAU D’EAU ET DE TEMPERATURE D’UNE CUVE

Le système est constitué d’un bac constitué d’un réservoir de section S et équipé d’une résistance
chauffante.

Les équations de fonctionnements du système sont les suivants :

dh −Q e 0 1
= h+ qe
dt 2 S H 0 S

d Θ −Θ s 0−Θe Q e 0 ( Θ s 0−Θe ) Qe0

= qe + pu− Θ
dt S H0 S H 0 Pu0 S H0
 1) Modélisation

Modélisons le système :

{
[ ][ ][
−Q e0
ḣ(t)
θ̇(t)
= 2 S H0
0
0

−0.033
×
][
h(t)
+
1 0
θ(t) 0.05 0.05
×][ ]
qe
pu

[ ][ ][ ][ ][ ]
d (t) 1 0
=
θ (t) 0 1
×
h(t ) 0 0
+
θ(t ) 0 0
×
qe
pu

Le point de fonctionnement est fixé au couple suivante : (H0, Pu0) = (0.6m,20Kw)

Notre système peut être mis sous la forme d’état suivante :

{
[ ][ ][ ][ ][ ]
ḣ(t) −0.0167 0 h(t)
= × + 1 0 × qe
θ̇(t) 0 −0.033 θ(t) 0.05 0.05 pu

[ θ(t)] 0 1 [ θ(t)]+[ 0 0] ×[ qepu]

=[
d (t) 1 0
] ×
h(t) 0 0

2) Analyse du système
2-1) Système échantillonné dans l’espace d’état

Donnons la représentation échantillonnée du système (Te=6min)

{[ ][ ][][ ][ ]
ḣ ( kTe+1 ) 0.002479 0 h + 59.85 0 × qe
= −6 ×
θ̇ ( kTe+1 ) 0 6.144 × 10 θ 1.5 1.5 pu

[
d ( kTe+1 )
θ ( kTe+1 ) ][ ][][ ][ ]
=
1 0
0 1
h 0 0
× +
θ 0 0
×
qe
pu

2-1)Système échantillonné dans l’espace d’état

Afin de convertir le système continu en analogique en numérique avec une période

d’échantillonnage de Te = 6 min.

Le code Matlab est le suivant :

…………..

{ [ ][ ] [][ ][ ]
ḣ ( kTe+1 ) −0.01667 0
= ×h + 1 0 × qe
θ̇ ( kTe+1 ) 0 −0.0333 θ 0.05 0.05 pu

[ ][ ][][ ][ ]
d ( kTe+1 ) 1 0
θ ( kTe+1 )
=
0 1
h 0 0
× +
θ 0 0
×
qe
pu

2-2) Étudions la stabilité du système d’après l’Lyapunov

Nous avons fixé la matrice Q= [ 60 06] (matrice diagonale et positive).

Nous avons déterminé la matrice P grâce à la commande Matlab suivante :

Le système est stable si et seulement si la matrice P est positive (Ses valeurs propres sont positives).
Le code pour obtenir P et déterminer ses valeurs propres est le suivant :

2-3)Étudions la commandabilité et l’observabilité du système grâce au gramien

2-3-1) La commandabilité

Le gramien de commandabilité nous permet de nous prononcer sur la commandabilité du système

en fonction de l’énergie a consommé pour réaliser une telle opération.

Le code nécessaire pour réaliser cela est le suivant :

…………..

Le résultat est le suivant :

……

2-3-2) L’observabilité

Le gramien d’observabilité nous permet de nous prononcer sur l’observabilité du système en

fonction de l’énergie a consommé pour réaliser une telle opération.

Le code nécessaire pour réaliser pour cela est le suivant :

…………..

Le résultat est le suivant :

……

2-3-3) Déterminons la matrice de transfère du système.

Le code Matlab permettant de l’obtenir est le suivant :

….

Le résultat est le suivant :

………

Il faut noter que le code Matlab donne les fonctions de transfert par entré.

2-3-3) Déterminons la matrice de transfère du système.