Università degli Studi di Roma Tor Vergata
FACOLTÀ DI INGNEGNERIA
Corso di Laurea Magistrale in Ingegneria dell’automazione
Progetto per il corso di controllo dei processi
Analisi e controllo di uno scambiatore di calore
Candidato:
Relatore:
Francesco Ferrante
Prof Riccardo Marino
Matricola 0157498
Anno Accademico 2011–2012
francesco ferrante
A N A L I S I E C O N T R O L LO D I U N O S C A M B I ATO R E D I C A LO R E
Il matematico gioca un gioco in cui egli stesso inventa le regole. Il fisico gioca un gioco
in cui le regole sono fornite dalla Natura. Ma, con il passare del tempo, diventa sempre
più evidente che le regole che il matematico trova interessanti sono quelle che la
Natura ha scelto.
—Paul Adrien Maurice Dirac
INDICE
Notazione
I
1
2
3
II
4
5
6
12
Analisi dello scambiatore di calore
descrizione dell’impianto
14
15
modello matematico dello scambiatore di calore
2.1 Modellazione degli organi di misura e di attuazione
2.2 Analisi F.E.M
24
17
23
analisi numerica dello scambiatore di calore
29
3.1 Caratterizzazione input-output dello scambiatore di calore
29
m
3.1.1 Analisi del legame v − ϑlu
29
m
3.1.2 Analisi del legame ϑli − ϑlu
36
m
3.1.3 Analisi del legame ϑv − ϑlu
38
3.2 Analisi dell’impianto basata sul modello
42
3.2.1 Analisi del comportamento dinamico del processo attorno
al punto di lavoro nominale
44
3.2.2 Analisi del comportamento dinamico del processo attorno
ad alcuni punti di lavoro ammissibili
49
Controllo dello scambiatore di calore
controllo dello scambiatore di calore
60
4.1 Descrizione del problema di controllo
60
4.2 Problema di regolazione asintotica
61
4.2.1 Primi passi verso la soluzione del problema di regolazione e assunzioni sul problema
62
4.2.2 Regolazione con azione feed-back+feed-forward
68
4.2.3 Controllo Observer-Based
70
4.3 Metodi basati sulla risposta armonica per la regolazione
70
4.3.1 Controllo PID based
71
4.4 Valutazione delle prestazioni del sistema di controllo in piccolo
73
sintesi dei regolatori per lo scambiatore di calore
5.1 Verifica delle condizioni per la sintesi del controllore
5.2 Controllo PID dello scambiatore di calore
78
5.3 Controllo PID con rete compensatrice
84
5.4 Controllo LQG
89
5.5 Scelta del controllore più soddisfacente
97
alcuni risultati numerici
101
6.1 Controllo PID con azione in avanti
4
102
76
77
59
6.2
6.3
7
Controllo PID con rete compensatrice con azione in avanti
Controllo LQG con azione in avanti 115
conclusioni
121
a strumenti necessari per la sintesi dei regolatori
bibliografia
132
122
indice
109
5
ELENCO DELLE FIGURE
Figura 1
Figura 2
Figura 3
Figura 4
Figura 5
Figura 6
Figura 7
Figura 8
Figura 9
Figura 10
Figura 11
Figura 12
Figura 13
Figura 14
Figura 15
Figura 16
Scambiatore di calore
15
Tratto di condotta di lunghezza infinitesima considerato
per il bilancio energetico
18
Profili di temperatura in condizioni nominali
21
Grafico della mappa della velocità all’equilibrio
22
Piano (θli θv ) dei punti di equilibrio. In blu vengono indicati i punti al di fuori della regione ammissibile mentre in
rosso i punti appartenenti alla regione ammissibile
23
Caratteristica statica delle servovalvola
23
Andamento dell’indice J ( N )
27
Profili di temperatura del liquido e del rame in condizioni nominali confrontati con le relative approssimazione
ottenute dal modello ad elementi finiti
28
Risposta ad un gradino di velocità pari a 0,31 m s−1 a partire dalla condizione indicata in 3.1.1.1
30
Sollecitazione dell’impianto con un ingresso di velocità
v(t) = v̄ + 0.1 sin(0.6 t) [m s−1 ]
32
Sollecitazione dell’impianto con un ingresso di velocità
v(t) = v̄ + 0.2 sin(0.6 t) [m s−1 ]
33
m ( t ) a fronAndamento della variazione di temperatura ∆ ϑlu
te della variazione della velocità del fluido rispetto al valore nominale: ∆ v(t) = −0.0203 δ−1 (t − 20) [m s−1 ]
33
Andamento reale ed approssimato con il modello (F.O.P.D.T.)
m ( t ) rispetto al valore
della variazione di temperatura ∆ ϑlu
nominale a fronte della variazione della velocità del fluido rispetto al valore nominale ∆ v(t) = −0.0203 δ−1 (t −
20) [m s−1 ]
35
Andamento reale ed approssimato con il modello (S.O.R.P.P.D.T.)
m ( t ) rispetto al valore
della variazione di temperatura ∆ ϑlu
nominale a fronte della variazione della velocità del fluido rispetto al valore nominale ∆ v(t) = −0.0203 δ−1 (t −
20) [m s−1 ]
36
Andamento reale ed approssimato con i modelli (S.O.R.P.P.D.T.)
m ( t ) rie (F.O.P.D.T. della variazione di temperatura ∆ ϑlu
spetto al valore nominale a fronte della variazione della
velocità del fluido rispetto al valore nominale ∆ v(t) =
−0.0203 δ−1 (t − 20) [m s−1 ]
36
m ( t ) a fronAndamento della variazione di temperatura ∆ ϑlu
te della variazione della temperatura del fluido in ingresso rispetto al valore nominale ∆ θli = 5 δ−1 (t − 20) [◦ C ]
37
6
Figura 17
Figura 18
Figura 19
Figura 20
Figura 21
Figura 22
Figura 23
Figura 24
Figura 25
Figura 26
Figura 27
Figura 28
Figura 29
Figura 30
Figura 31
Figura 32
Figura 33
Elenco delle figure
Andamento reale ed approssimato con il modello (F.O.P.D.T.)
m ( t ) rispetto al valore
della variazione di temperatura ∆ ϑlu
nominale a fronte della variazione della temperatura del
fluido in ingresso rispetto al valore nominale ∆ θli (t) =
5 δ−1 (t − 20) [◦ C ]
38
m ( t ) a fronAndamento della variazione di temperatura ∆ ϑlu
te della variazione della temperatura del fluido in ingresso rispetto al valore nominale ∆ θv = 5 δ−1 (t − 20) [◦ C ]
39
Andamento reale ed approssimato con il modello (F.O.P.D.T.)
m ( t ) rispetto al valore
della variazione di temperatura ∆ ϑlu
nominale a fronte della variazione della temperatura del
fluido in ingresso rispetto al valore nominale ∆ θv (t) =
5 δ−1 (t − 20) [◦ C ]
39
Andamento reale ed approssimato con il modello (S.O.R.P.D.T.)
m ( t ) rispetto al valodella variazione di temperatura ∆ ϑlu
re nominale a fronte della variazione della temperatura
del vapore rispetto al valore nominale ∆ θv (t) = 5 δ−1 (t −
20) [◦ C ]
40
Andamento reale ed approssimato con il modello (F.O.P.D.T.)
e con il modello (S.O.R.P.P.D.T.) della variazione di temm ( t ) rispetto al valore nominale a fronte della
peratura ∆ ϑlu
variazione della temperatura del vapore rispetto al valore
nominale ∆ θv (t) = 5 δ−1 (t − 20) [◦ C ]
41
Diagramma di Bode della funzione di trasferimento Fy v (s)
46
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
46
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
47
Diagramma di Bode della funzione di trasferimento Fy v (s)
e della sua approssimante (3.3)
47
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
e della sua approssimante (3.5)
48
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
e della sua approssimante (3.7)
48
Diagramma di Bode della funzione di trasferimento Fy v (s)
per il punto di lavoro θli = 75 ◦C θv = 140 ◦C
49
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
per il punto di lavoro θli = 75 ◦C θv = 140 ◦C
50
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
per il punto di lavoro θli = 75 ◦C θv = 140 ◦C
50
Diagramma di Bode della funzione di trasferimento Fy v (s)
per il punto di lavoro θli = 50 ◦C θv = 140 ◦C
51
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
per il punto di lavoro θli = 50 ◦C θv = 140 ◦C
52
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
per il punto di lavoro θli = 50 ◦C θv = 140 ◦C
52
7
Figura 34
Figura 35
Figura 36
Figura 37
Figura 38
Figura 39
Figura 40
Figura 41
Figura 42
Figura 43
Figura 44
Figura 45
Figura 46
Figura 47
Figura 48
Figura 49
Figura 50
Figura 51
Figura 52
Figura 53
Elenco delle figure
Diagramma di Bode della funzione di trasferimento Fy v (s)
per il punto di lavoro θli = 75 ◦C θv = 105 ◦C
53
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
per il punto di lavoro θli = 75 ◦C θv = 105 ◦C
54
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
per il punto di lavoro θli = 75 ◦C θv = 105 ◦C
54
Diagramma di Bode della funzione di trasferimento Fy v (s)
per il punto di lavoro θli = 50 ◦C θv = 105 ◦C
55
Diagramma di Bode della funzione di trasferimento Fy w1 (s)
per il punto di lavoro θli = 50 ◦C θv = 105 ◦C
56
Diagramma di Bode della funzione di trasferimento Fy w2 (s)
per il punto di lavoro θli = 50 ◦C θv = 105 ◦C
56
Diagrammi di Bode della funzione di trasferimento Fy v (s)
per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma sono indicate con Fi (s) con
i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata
in (Tabella 1)
57
Diagrammi di Bode della funzione di trasferimento Fy w1 (s)
per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma sono indicate con Fi (s) con
i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata
in (Tabella 1)
57
Diagrammi di Bode della funzione di trasferimento Fy w2 (s)
per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma sono indicate con Fi (s) con
i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata
in (Tabella 1)
58
Schema di controllo adottato
65
Schema di controllo con azione integrale
68
Schema di controllo con azione avanti
68
Diagramma di Bode della risposta armonica del controllore PID
80
Diagramma di Bode della funzione d’anello L(s)
80
Diagramma di Nyquist della funzione d’anello L(s)
81
Diagrammi di Bode della funzione di sensitività complementare T (s) e della funzione H(s) T (s)
82
Diagramma di Bode del modulo della funzione di sensitività S(s)
82
Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s) con i ∈ {0, 1, 2, 3, 4} a
seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata in (Tabella 1)
83
Diagrammi di Bode della funzione di compensazione G (s)
Diagrammi di Bode del controllore C (s) G (s)
85
8
85
Figura 54
Figura 55
Figura 56
Figura 57
Figura 58
Figura 59
Figura 60
Figura 61
Figura 62
Figura 63
Figura 64
Figura 65
Figura 66
Figura 67
Figura 68
Figura 69
Figura 70
Figura 71
Figura 72
Figura 73
Figura 74
Figura 75
Figura 76
Figura 77
Elenco delle figure
9
Diagrammi di Bode della funzione di trasferimento P(s) a
confronto con la funzione d’anello L(s) = P(s) C (s) G (s)
86
Diagramma di Nyquist della funzione di anello L(s) nel
caso nominale
86
Diagramma di Bode della funzione sensitività nel caso
nominale
87
Diagramma di Bode della funzione T (s) H(s) nel caso
nominale
87
Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s) con i ∈ {0, 1, 2, 3, 4} a
seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata in (Tabella 1)
88
Valori singolari della funzione di trasferimento del processo linearizzato attorno al punto di lavoro nominale
90
Diagrammi di Bode della funzione di trasferimento del
processo nominale e della sua approssimazione di ordine
6
91
Diagrammi di Bode della funzione di trasferimento del
compensatore LQG con azione integrale
92
Valori singolari associati ai modi naturali presenti nel legame ingresso-uscita del controllore LQG
93
Diagramma di Bode di Cr (s) e di C (s) a confronto
93
Diagramma di Bode della funzione di trasferimento Cd (s)
94
Diagramma di Bode della funzione di trasferimento L(s)
94
Diagramma di Bode del modulo della funzione di sensitività S(s)
95
Diagramma di Bode della funzione di trasferimento H(s) T (s)
95
Diagramma di Nyquist della funzione d’anello L(s)
96
Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s) con i ∈ {0, 1, 2, 3, 4} a
seconda del punto di lavoro considerato, in maniera consistente alla numerazione riportata in (Tabella 1)
96
Diagrammi di Bode del modulo della funzione di sensitività per i tre controllori nel caso nominale
98
Diagrammi di Bode della funzione di trasferimento tra
la variazione del riferimento di temperatura e l’errore di
regolazione H(s) T (s) per i tre controllori nel caso nominale
99
Diagrammi di Bode delle funzioni di trasferimento C (s) e
Cd (z) a confronto
100
Andamento del riferimento y∗ (t)
101
Andamento delle variabili ϑli (t) e ϑv (t)
101
Andamento dell’uscita controllata y(t)
102
Andamento del controllo u(t)
102
Variazione delle temperatura di riferimento da 98 ◦ a 94 ◦
103
Figura 78
Figura 79
Figura 80
Figura 81
Figura 82
Figura 83
Figura 84
Figura 85
Figura 86
Figura 87
Figura 88
Figura 89
Figura 90
Figura 91
Figura 92
Figura 93
Figura 94
Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
104
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
105
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
106
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
107
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
108
◦
Variazione delle temperatura di riferimento da 98 a 94 ◦
109
Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
110
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
111
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
112
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
113
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
114
◦
Variazione delle temperatura di riferimento da 98 a 94 ◦
115
Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
116
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
117
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
118
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
119
Variazione della temperatura del fluido in ingresso e della
temperatura del vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
120
E L E N C O D E L L E TA B E L L E
Tabella 1
Tabella 2
Valore della temperatura del liquido in ingresso e del
vapore per i punti di lavoro scelti. Il punto di lavoro
nominale viene indicato con il numero 0.
44
Valore dei margini di guadagno e di fase ottenuti per
i punti di lavoro scelti con il controllo PID. Il punto di
lavoro nominale viene indicato con il numero 0.
83
10
Tabella 3
Tabella 4
Tabella 5
Elenco delle tabelle
Valore dei margini di guadagno e di fase ottenuti per i
punti di lavoro scelti con il controllo PID+ Rete compensatrice. Il punto di lavoro nominale viene indicato con il
numero 0.
88
Valore dei margini di guadagno e di fase ottenuti per i
punti di lavoro scelti con il controllo LQG. Il punto di
lavoro nominale viene indicato con il numero 0.
92
Valore dei margini di guadagno e di fase ottenuti per i
punti di lavoro scelti con il controllo LQG. Il punto di
lavoro nominale viene indicato con il numero 0.
98
11
N OTA Z I O N E
(
δ−1 (t − t0 ) :=
1
0
t ≥ t0
t < t0
Funzione gradino
C
Insieme dei numeri complessi
C+
Insieme dei numeri complessi con parte reale non negativa
R
Insieme dei numeri reali
R+
Insieme dei numeri reali non negativi
Rn
Vettori reali di dimensione n
Rn × m Matrici reali di dimensione n × m
In
|| u||
s
Matrice identità di dimensione n
√
u T u Norma euclidea
Variabile complessa del dominio di Laplace
T.H.D.% Distorsione armonica percentuale (Total Harmonic Distortion)
x T , ( A T ) Vettore x trasposto (matrice A trasposta)
C k ( D ) Insieme delle funzioni continue e con le prime k derivate continue nell’insieme D
L∞
Insieme delle funzioni reali limitate
12
INTRODUZIONE
I
n questa trattazione si vuole affrontare il problema della analisi e del
controllo di uno scambiatore di calore. La trattazione verrà affrontata dividendo il testo in due parti: La prima relativa all’analisi del processo
considerato mentre la seconda dedicata alla sintesi del sistema di controllo. In
particolare nella prima parte verrà derivato un modello matematico del processo basato sulle equazioni che governano il fenomeno fisico portante, successivamente su tale modello verrano condotte delle analisi atte a caratterizzare la
dinamica del processo fisico e il comportamento ingresso uscita del medesimo.
Nella seconda parte invece verrà formalizzato il problema di controllo dello
scambiatore di calore e verranno presentati dei risultati che condurranno alla
soluzione del problema. Infine si mostreranno i risultati ottenuti dal controllo
dello scambiatore di calore.
13
Parte I.
Analisi dello scambiatore di calore
14
1
I
D E S C R I Z I O N E D E L L’ I M P I A N T O
n termotecnica con il termine scambiatore di calore si fa riferimento
ad un componente atto a realizzare uno scambio di energia termica
tra due fluidi a temperature diverse. In generale gli scambiatori sono sistemi
aperti che operano senza scambio di lavoro, ovvero presentano un flusso costante di fluido e una distribuzione di temperatura a regime costante. Nella realtà
industriale esistono molteplici metodologie per effettuare lo scambio termico
tra due fluidi ossia esistono diverse tecnologie per realizzare uno scambiatore
di calore.
Nella trattazione che ci appresta a presentare verrà presa in considerazione una
particolare tipologia di scambiatore di calore detto Scambiatore a superficie. In
questo tipo di scambiatore le correnti assorbono il calore dalle superfici con
le quali sono separati, inoltre in tali dispositivi si riconoscono sempre due lati
contenenti i fluidi. In particolare la struttura che verrà presa in considerazione in questo ambito è mostrata in (Figura 1). In questa struttura il fluido ad
Figura 1.: Scambiatore di calore
un temperatura ϑlim , acqua nel caso in esame, scorre con velocità v all’interno
di una serpentina in rame la quale è immersa in una camera contenente del
vapore saturo a temperatura ϑv ed infine il fluido sbocca dalla serpentina con
una temperatura ϑlu . In particolare nell’impianto in esame la velocità del fluido
v può essere variata, agendo sulla apertura di una servovalvola calettata sulla
serpentina stessa, mentre la temperatura del fluido in uscita ϑlu può essere misurata attraverso una termocoppia installata allo sbocco della condotta dello
scambiatore di calore. E’ evidente che a seconda del valore delle variabili appena menzionate l’impianto si troverà in una condizione operativa piuttosto che
in un altra. Ebbene tra tutte le possibili condizioni operative in cui l’impianto può trovarsi si può prendere in considerazione il punto di lavoro riferito ai
15
valori nominali delle variabili in gioco, ossia:
descrizione dell’impianto
• ϑlu = 98◦ C
• ϑv = 120◦ C
• ϑli = 65◦ C
tale punto viene denominato punto di lavoro nominale e verrà considerato come il punto di lavoro di riferimento per condurre la maggior
parte delle analisi sullo scambiatore di calore.
16
2
A
M O D E L LO M AT E M AT I C O D E L LO
S C A M B I ATO R E D I C A LO R E
l fine di poter effettuare un analisi numerica dello scambiatore di
calore è necessario derivare per tale sistema un modello matematico sufficientemente affidabile. Prima di passare alla scrittura del modello
matematico del processo in esame è bene fare alcune assunzioni semplificative
[3].
• La temperatura del vapore saturo è uguale in ogni punto della camera
• Il fluido è assunto essere incompressibile
• Il fluido è assunto essere in regime laminare
• La serpentina è considerata come un tubo rettilineo
Si noti che l’aver assunto il fluido in regime laminare impone un
vincolo sulla velocità v dello stesso, la quale ovviamente non può
assumere valori qualunque. In particolare tale assunzione implica che:
v ∈ [0.1, 0.7] m s−1
Inoltre si supporrà nel seguito della trattazione che:
ϑli ≤ 80◦ C
(2.1)
◦
ϑli ≥ 100 C
(2.2)
A fronte delle ipotesi semplificative appena esposte, al fine di formulare un
modello matematico del processo in esame, si può effettuare un bilancio delle
energie scambiate da una sezione di lunghezza infinitesima della serpentina
in un intervallo temporale infinitesimo. A tal scopo è conveniente dapprima
fissare una ascissa di riferimento, che indicheremo con ξ, in modo da poter
individuare in maniera univoca ogni punto della serpentina. Detto ciò si indichino ora con ϑl (ξ, t) e con ϑr (ξ, t) rispettivamente, la temperatura del generico
punto del liquido e del rame nell’istante di tempo t e con ϑv (t) la temperatura
del vapore nella camera nell’istante di tempo t, dunque il bilancio energetico può essere effettuato considerando dapprima i flussi di energia entranti ed
uscenti, dalla porzione di liquido contenuta nel tratto di condotta di lunghezza
infinitesima in esame e successivamente prendendo in considerazione i medesimi flussi relativamente alle pareti della condotta (Figura 2). In particolare
relativamente ai flussi di energia riguardanti il liquido, si possono mettere in
evidenza principalmente tre forme di scambio termico:
• Energia che il liquido scambia al suo interno
• Energia scambiata con la parete della serpentina
• Energia che il liquido trasporta a causa del suo moto nella serpentina
17
modello matematico dello scambiatore di calore
Figura 2.: Tratto di condotta di lunghezza infinitesima considerato per il bilancio
energetico
Ebbene indicando con:
• αlm coefficiente di adduzione acqua-rame [ J m−2 s−1 C −1 ]
• cl calore specifico dell’acqua [ J Kg−1 C −1 ]
• Dl diametro interno della condotta [m]
• µl densità lineare dell’acqua [Kg m−1 ]
• v velocità del fluido nella condotta [m s−1 ]
gli scambi termici sopra elencati in maniera qualitativa possono
essere espressi in modo quantitativo come segue:
• Energia che il liquido scambia al suo interno
−cl µl dξ (ϑl (ξ, t + dt) − ϑl (ξ, t))
• Energia scambiata con la parete della serpentina
αlm dξ (ϑr (ξ, t) − ϑl (ξ, t))π Dl dξ dt
• Energia che il liquido trasporta a causa del suo moto nella serpentina
µl v dt cl (ϑl (t, ξ ) − ϑl (t, ξ + dξ ))
Viceversa relativamente ai flussi di energia riguardanti le pareti della condotta si possono mettere in evidenza, anche in questo caso, tre forme di scambio
termico:
• Energia che la parete scambia con se stessa
• Energia che la parete scambia con il liquido
• Energia che la parete scambia con il vapore
In questo caso indicando con:
18
modello matematico dello scambiatore di calore
• αv coefficiente di adduzione vapore-rame [ J m−2 s−1 C −1 ]
• cr calore specifico del rame [ J Kg−1 C −1 ]
• Dv diametro interno della condotta [m]
• µr densità lineare del rame [Kg m−1 ]
tali termini possono essere quantificati come segue:
• Energia che la parete scambia con se stessa
−µr dξ cr (ϑr (ξ, t + dt) − ϑr (ξ, t))
• Energia che la parete scambia con il liquido
−αlm π Dl dξ dt (ϑr (ξ, t) − ϑl (ξ, t))
• Energia che la parete scambia con il vapore
αv π Dv ξ dt (ϑr (ξ, t) − ϑv (t))
A questo punto imponendo la conservazione dell’energia si ottengono due bilanci energetici uno per il liquido ed uno per il rame della condotta.
Bilancio energetico del liquido:
− cl µl dξ (ϑl (ξ, t + dt) − ϑl (ξ, t)) + αlm dξ (ϑr (ξ, t) − ϑl (ξ, t))π Dl dξ dt+
+ µl v dt c1 (ϑl (t, ξ ) − ϑl (t, ξ + dξ )) = 0
(2.3)
Bilancio energetico del rame:
µr dξ cr (ϑr (ξ, t + dt) − ϑr (ξ, t)) + αlm π Dl dξ dt (ϑr (ξ, t) − ϑl (ξ, t))+
+ αv π Dv ξ dt (ϑr (ξ, t) − ϑv (t)) = 0
(2.4)
Pertanto al fine di ottenere delle equazioni di validità locale si può effettuare
un passaggio al limite e confondere gli incrementi parziali delle varie temperature con le relative derivate prime. In tal modo si ottiene il seguente sistema di
equazioni alle derivate parziali (E.D.P.):

∂ϑr (ξ, t)


αv (ϑv (t) − ϑr (ξ, t)) π Dv − αlm π Dl (ϑr (ξ, t) − ϑl (ξ, t)) − µr cr
=0


dt









∂(ϑl (ξ, t))
∂ϑ (ξ, t)


− cl µl
+ αlm π Dl (ϑr (ξ, t) − ϑl (ξ, t)) + µl v c1 l
=0



dt
dξ





ϑl (0, t) = ϑli (t) ∀t ≥ 0










ϑr (ξ, 0) = ϑr ( x ) ∀ ξ ∈ [0 , L]









ϑl (ξ, 0) = ϑl ( x ) ∀ ξ ∈ (0 , L]
19
modello matematico dello scambiatore di calore
(2.5)
Il significato delle boundary-conditions alle quali è soggetto il sistema di (E.D.P.)
è presto detto: La temperatura del liquido in ingresso deve poter essere specificata istante per istante mentre i profili spaziali delle temperature del rame e
del liquido nell’istante iniziale devono essere specificate al fine di ottenere un
problema ben posto. Fatte le necessarie precisazioni si può riscrivere il sistema
di (E.D.P.) (2.5) in maniera più compatta come segue:

∂ϑr (ξ, t)


=0
(ϑv (t) − ϑr (ξ, t)) − β (ϑr (ξ, t) − ϑl (ξ, t)) − τ2


dt









∂ϑ (ξ, t)
∂(ϑl (ξ, t))


+ (ϑr (ξ, t) − ϑl (ξ, t)) + v τ1 l
=0
− τ1



dt
dξ





ϑl (0, t) = ϑli (t) ∀t ≥ 0










ϑr (ξ, 0) = ϑr ( x ) ∀ ξ ∈ [0 , L]









ϑl (ξ, 0) = ϑl ( x ) ∀ ξ ∈ (0 , L]
(2.6)
Dove le variabili τ1 , τ2 , β sono espresse come segue:
β :=
Dl αlm
Dv α v
(2.7)
(2.8)
τ1 :=
cl µl
[s]
π Dl αlm
(2.9)
(2.10)
cr µr
τ2 :=
[s]
π Dv α v
(2.11)
Ottenuta una scrittura più snella per il sistema (E.D.P.) (2.5) si può passare ora
ad ottenere una soluzione analitica di tipo stazionario di tale sistema al fine
di poter valutare il profilo di temperatura che si stabilisce lungo la condotta a
transitorio esaurito a fronte di una una temperatura costante per il liquido in
ingresso
ϑli (t) ≡ ϑli0
Ebbene indicando con (ϑr (ξ ) , ϑl (ξ ) la soluzione cercata, allora essa certamente
deve soddisfare il seguente sistema di equazioni misto Algebrico-E.D.O:
20
modello matematico dello scambiatore di calore

(ϑv (t) − ϑr (ξ )) − β (ϑr (ξ ) − ϑl (ξ )) = 0











 (ϑ (ξ ) − ϑ (ξ )) + v τ ∂ϑl (ξ ) = 0
r
1
l
dξ







0


 ϑl (0, t) = ϑli ∀t ≥ 0


(2.12)
A fronte di opportune manipolazione e di una successiva integrazione rispetto
alla variabile spaziale ξ si ottengono le seguenti espressioni per la soluzione
cercata:

− τ v (xβ+1)

1
ϑ
(
ξ
)
=
e
(ϑli0 − ϑv ) + ϑv

l









 ϑr (ξ ) =
β
β+1
e
x
1 v ( β +1)
−τ
(2.13)
!
(ϑli0 − ϑv ) + ϑv +ϑv
Attraverso l’uso delle relazioni ottenute è possibile graficare i profili di temperatura del liquido e del rame che ci stabiliscono in condizioni stazionarie,
mostrati in (Figura 3)
120
115
110
105
100
θ(ξ ) ◦ C
95
90
85
80
75
θ r (ξ )
70
θ l (ξ )
65
60
0
0.5
1
1.5
2
2.5
ξ [m]
Figura 3.: Profili di temperatura in condizioni nominali
In particolare essendo la temperatura del fluido in uscita la grandezza di
maggior interesse è conveniente disporre di un espressione di tale valore dipendente dagli altri parametri di processo, tale espressione può essere ottenuta
a partire dalla (2.13) valutando ϑl (ξ ) nel punto ξ = L, l’espressione che si
ottiene è la seguente:
ϑ ( L) = ϑlu = e
L
1 v ( β +1)
−τ
(ϑli0 − ϑv ) + ϑv
21
modello matematico dello scambiatore di calore
Dal legame appena mostrato è possibile derivare il valore della velocità del
fluido necessario ad ottenere una data temperatura per il liquido in uscita dalla
condotta, in particolare la relazione suddetta è la seguente:
v=
L
(2.14)
ϑv
τ1 ( β + 1) ln( ϑϑluli −
− ϑv )
il legame appena pervenuto è di fondamentale importanza in quanto permette
di stabilire, a partire dalle condizioni operative del processo, se sia possibile o
meno condurre l’impianto in un equilibrio assegnato, ossia permette di stabilire se il valore di velocità da imporre al fluido in una data condizione operativa
dell’impianto si ammissibile o meno. A titolo di esempio si consideri il caso
in cui si voglia mantenere la temperatura del fluido in uscita al suo valore nominale ossia 98◦ C, anche a fronte di variazioni della temperatura del fluido in
ingresso e della temperatura del vapore, ebbene è evidente che esisteranno delle condizioni operative per le quali il valore della velocità richiesto non rientra
nei vincoli imposti dalla assunzione di regime laminare. Per ragioni che saranno chiare più avanti è conveniente definire l’insieme ammissibile dei punti di
equilibrio ammissibili in maniera precisa.
Definizione 2.1: Si definisce insieme ammissibile delle condizioni operative
Θϑlu l’insieme definito come segue:
(
)
L
Θϑlu : (ϑli , ϑv ) ∈ R2 : 0.1 ≤
≤ 0.7
ϑv
τ1 ( β + 1) ln( ϑϑluli −
− ϑv )
A supporto di quanto appena illustrato, in (Figura 4), viene riportato il grafico della relazione (2.14) ottenuto facendo variare ϑv e ϑli in un intervallo opportuno e fissando ϑlu al suo valore nominale ossia 98 [◦ C ] Inoltre per rendersi con-
0.7
0.8
0.7
0.6
0.6
v [m s−1]
0.5
0.5
0.4
0.4
0.3
0.2
0.3
0.1
0
120
120
100
100
60
60
θli[°C]
0.2
80
80
40
40
20
20
0
0
θv[°C]
0.1
Figura 4.: Grafico della mappa della velocità all’equilibrio
to di quale sia regione ammissibile, nel caso di temperatura del fluido di uscita
nominale, in (Figura 5) vengono discriminati, tramite un diversa colorazione, i
punti appartenenti alla regione ammissibile Θ98 dai restanti.
Prima di concludere questa sezione, si osservi che la relazione appena discussa può essere utilizzata per derivare il valore della velocità del fluido in condizioni nominali v̄ (velocità nominale del fluido), il quale a priori risulta non noto,
22
2.1 modellazione degli organi di misura e di attuazione
140
135
130
θv [ ◦C ]
125
120
115
110
105
100
50
55
60
65
70
75
80
θl i [ ◦C ]
Figura 5.: Piano (θli θv ) dei punti di equilibrio. In blu vengono indicati i punti al di
fuori della regione ammissibile mentre in rosso i punti appartenenti alla
regione ammissibile
ebbene a riguardo la ( 2.14 nella pagina precedente) porge il seguente valore
per la velocità nominale:
v̄ ≈ 0,31 m s−1
2.1
modellazione degli organi di misura e di attuazione
Come già indicato espressamente nella situazione in esame sono presenti degli
organi di misura e degli organi di attuazione, in particolare sono presenti una
termocoppia per la misura della temperatura del fluido in uscita e una servovalvola per la regolazione della velocità del fluido nella serpentina. Tali dispositivi
debbono essere modellati oculatamente al fine di ottenere un modello accurato
del processo sotto controllo. In particolare per quanto riguarda il dispositivo di
attuazione della valvola è possibile, data la sua elevata rapidità, trascurare la
sua presenza e tenere in conto solamente dell’intervallo delle velocità ammissibili per il fluido modellando la servovalvola come un dispositivo non lineare
statico descritto dalla caratteristica in (Figura 6). Relativamente alla termocop-
Figura 6.: Caratteristica statica delle servovalvola
pia, è possibile modellare tale dispositivo come un sistema a singola costante
m la misura della temperatura del liquido
di tempo. Con ciò indicando con ϑlu
23
2.2 analisi f.e.m
in ingresso fornita dalla termocoppia, con z lo stato associato a tale dispositivo
e scegliendo ragionevolmente per la costante di tempo del dispositivo il valore
di 1 [s], si ottiene il modello del dispositivo suddetto:
ż = −z + ϑlu
m
ϑlu
2.2
(2.15)
=z
(2.16)
analisi f.e.m
Il modello ottenuto, essendo descritto da un set di (E.D.P.), è di fatto un modello a dimensione infinita, dunque per poter effettuare una simulazione numerica del comportamento dello scambiatore di calore si deve realizzare una
analisi di tipo (F.E.M), ossia si effettua una discretizzazione spaziale (mesh) sufficientemente fitta della condotta, ottenendo in tale modo un set di equazioni
di (E.D.O.) ognuna delle quali rappresenta l’evoluzione della temperatura di
ognuno dei campioni spaziali (elementi finiti) ottenuti a fronte dell’operazione
di discretizzazione. In letteratura esistono molte metodologie per effettuare l’operazione di discretizzazione del dominio spaziale per il quale si sta cercando
la soluzione del sistema di (E.D.P.). Tuttavia in questa sede si utilizzerà il metodo più semplice di essi il quale, nel caso in esame, prevede di suddividere il
dominio spaziale in N elementi uguali fra loro di lunghezza
l=
L
N
e di considerare come stato la temperatura del liquido e del rame per ognuno
dei segmenti di condotta considerati. In tal modo si è implicitamente assunto
che la temperatura all’interno del generico elemento finito sia omogenea, tale
assunzione ovviamente introduce un errore che diminuisce all’aumentare del
numero di elementi che compongono la mesh. Al fine di ottenere il modello
(F.E.M) è necessario definire la mesh, ossia individuare gli elementi che comporranno la medesima. Tra tutte le scelte, quella che consente di ottenere un
modello ad elementi finiti per il quale sia possibile specificare istante per istante la temperatura del liquido ingresso è la seguente:
Elementi della Mesh relativa al liquido:


ϑl1 := ϑl (l )




 ϑl2 := ϑl (2 l )
 ...




 ϑ N := ϑ ( N l )
l
l
Elementi della Mesh relativa al rame:


ϑr1 := ϑr (l )




 ϑr2 := ϑr (2 l )
 ...




 ϑ N := ϑ ( N l )
r
r
(2.17)
(2.18)
24
2.2 analisi f.e.m
25
A fronte della scelta fatta per la mesh, approssimando le derivate parziali rispetto alla variabile spaziale ξ con i relativi rapporti incrementali all’indietro1 ,
presenti nelle (2.6) si ottiene il seguente sistema (2.19) di 2 N (O.D.E) seguente,
il quale rappresenta l’approssimazione (F.E.M) del sistema (2.6):

dϑri (t)

i
i
i

=0
(
ϑ
(
t
)
−
ϑ
(
t
)
−
β
(
ϑ
(
t
)
−
ϑ
(
t
))
−
τ

v
2
r
r
l

dt










ϑli (t) − ϑli−1 (t)
d(ϑli (t))

i
i

−
τ
+
(
ϑ
=0
(
t
)
−
ϑ
(
t
))
+
v
τ

1
1
r
l

dt
l



(2.19)
i


ϑri (0) = ϑr0









i

ϑli (0) = ϑl0









∀ i ∈ {1, . . . N }
Tale sistema può essere riscritto in forma matriciale come segue:
 ˙
ϑl1 (t)

..

.

 ˙N
 ϑl ( t )


 ˙1
 ϑr (t)


..

.
˙
N
ϑr (t)


 
 
 
 
 
 
=
 
 
 
 

−( vl +
v
l
0
β
τ2
0
0

v
l





+




0
0
0
..
.
0
1
τ1 )
0
0
..
.
0
..
v
. −( l +
0
..
.
0
0
0
β
τ2
1
τ1
1
τ1 )
0
0
0
..
.
0
0
1
τ1
0
..
.
0
−( τβ2 +
0
0
1
τ2 )
0
0
−( τβ2
0
+

 1
ϑl ( t )

..



.



  ϑlN (t) 




+

  ϑ1 (t ) 

 r


..



.
1
N
ϑr (t)
τ2 )

0
.. 
. 

0 
 ϑli
1 
ϑv
τ2 
.. 

. 
1
τ2
(2.20)
Per concludere, al fine di ottenere un modello che possieda come uscita la temm , è necessario introdurre la dinamica
peratura misurata del liquido in uscita ϑlu
della termocoppia ( 2.15 nella pagina precedente) caratterizzata dalla variabile
di stato z, ottenendo il modello dinamico seguente:
1 Utilizzando un approssimazione basata su rapporto incrementale in avanti, scrivendo l’equazione differenziale che descrive l’evoluzione della temperatura del liquido nel primo tratto di
condotta, si introdurrebbe un vincolo sulla derivata della temperatura del liquido in ingresso,
contraddicendo l’ipotesi di arbitrarietà di tale variabile.
2.2 analisi f.e.m
 ˙1
ϑl ( t )

..

.

 ϑ˙N (t)
 l


 ˙1
 ϑr (t)

..


.

 ϑ˙rN (t)
ż


 
 
 
 
 
 
 
=
 
 
 
 
 

−( vl +
v
l
0
β
τ2
0
0
0

v
l





+




0
0
0
..
.
0
1
τ1 )
0
..
.
..
.
0
..
.
0
...
1
τ1
0
0
−( vl
+
0
1
τ1 )
0
0
..
.
0
0
0
1
τ1
0
..
.
0
−( τβ2 +
0
β
τ2
1
0
0
0
1
τ2 )
0
...
0
0
+
−1
−( τβ2













1

)

τ2
26

ϑl1 (t)

..

.

ϑlN (t) 



+
1
ϑr (t) 

..

.


ϑrN (t) 
z

0
.. 
. 

0 
 ϑli
1 
ϑv
τ2 
.. 

. 
1
τ2
m
ϑlu
=z
(2.21)
Stabilita la struttura del sistema ad elementi finiti, al fine di poter utilizzare
il modello per le relative simulazioni si deve stabilire quale sia il numero di
elementi finiti opportuno nel caso in esame. Con lo scopo di ottenere un buon
modello è conveniente fare fede ad un opportuno criterio di ottimo minimizzando una particolare cifra di merito che tenga conto dell’errore commesso nella
discretizzazione del problema. Ebbene a fronte di numerose analisi statistiche è
stato osservato che una buona cifra di merito è data dalla seguente espressione:
J(N) =
1
2N
N
∑
k =1
c (ϑlk − ϑl (k l ))2 + (1 − c)(ϑrk − ϑr (k l ))2
(2.22)
In particolare scegliendo c = 0.8 si attribuisce un peso maggiore all’errore
commesso sull’approssimazione del profilo ti temperatura del liquido. Ebbene osservando l’andamento dell’indice J, al variare del numero di elementi
N, mostrato in (Figura 7) si nota che scegliendo N = 35 si ottiene un buon
compromesso tra precisione del modello (almeno in condizioni statiche) e efficienza computazionale. A posteriori delle scelta fatta per N si ottiene il modello
desiderato nella forma (2.21), nello specifico ponendo
x := [ϑl1 , . . . , ϑl35 , . . . . . . , ϑr1 , . . . , ϑr35 , z] T
(2.23)
w := [ϑli , ϑv ] T
2.2 analisi f.e.m
100
90
80
70
J (N )
60
50
40
30
20
10
0
0
10
20
30
40
50
60
70
80
90
100
N
Figura 7.: Andamento dell’indice J ( N )
il modello ottenuto può essere riscritto nella forma:
ẋ = A(v) x + B(v) w
(2.24)
y = Cx
(2.25)
avendo cura di indicare rispettivamente:

−( vl +
1
τ1 )
v
l






A(v) = 






0
β
τ2
0
..
.
0
0
0

v
l





B(v) = 




0
0
0
..
.
0
0
..
.
..
.
0
...

0
.. 
. 

0 

1 
τ2 
.. 

. 
0
1
τ1
0
0
−( vl +
0
0
β
τ2
1
1
τ1 )
0
−( τβ2
+
0
0
0
0
..
.
0





1

0
τ1


1
)
0
0

τ2

..

.
0


β
1
0 −( τ2 + τ2 ) 
...
−1
0
(2.26)
1
τ2
C = 0 ... 0 0 ... 1
Dunque al fine di valutare l’accuratezza del modello ad elementi finiti, si può
determinare quale sia il profilo di temperatura che si stabilisce in regime stazionario ottenuto da tale modello. A tal scopo si può procedere come segue:
1. Fissate ϑli ,ϑv e ϑlu si determina il valore v̄ della velocità del fluido all’
equilibrio tramite la relazione (2.14)
27
2.2 analisi f.e.m
2. Si determina il valore dello stato di equilibrio x̄ risolvendo la seguente
equazione lineare nell’incognita x̄
A(v̄) x̄ + B(v̄) w = 0
La procedura appena descritta può essere utilizzata nelle condizioni nominali,
ossia:
ϑli = 65◦ C
ϑlu = 98◦ C
(2.27)
◦
ϑv = 120 C
ottenendo il profilo mostrato in (Figura 8). Quanto appena mostrato di fatto
conclude la ricerca del modello ad elementi finiti necessario alla simulazione
dello scambiatore di calore.
100
θ l (k l )
θ l (ξ )
95
θ l (ξ ) ◦ C
90
85
80
75
70
65
0
0.5
1
1.5
2
2.5
ξ [m]
(a) Profilo di temperatura del liquido
117
θ r (ξ )
θ r (k l )
116
θ r (ξ ) ◦ C
115
114
113
112
111
0
0.5
1
1.5
2
2.5
ξ [m]
(b) Profilo di temperatura del rame
Figura 8.: Profili di temperatura del liquido e del rame in condizioni nominali confrontati con le relative approssimazione ottenute dal modello ad elementi
finiti
28
3
N
ANALISI NUMERICA DELLO
S C A M B I ATO R E D I C A LO R E
el capitolo precedente è stato ottenuto un modello ad elementi finiti
in grado di descrivere in maniera accurata lo scambiatore di calore. In questo capitolo si effettueranno delle analisi sul modello matematico ad
elementi finiti dello scambiatore di calore ottenuto in precedenza, atte a mettere
in luce quali siano le caratteristiche dello scambiatore stesso dal punto di vista
ingresso uscita. Tale caratterizzazione di fatto è necessaria ogni qual volta si
voglia effettuare il controllo di un determinato processo.
3.1
caratterizzazione input-output dello scambiatore di calore
Al fine di poter ottenere una caratterizzazione del processo in grado di fornire
delle indicazioni sulla sintesi di una legge di controllo è necessario effettuare
sul processo stesso delle analisi atte a mettere in evidenza quale sia la natura
dei legami input-output in questione. Prima di ogni altra cosa è bene fare chiarezza su quale sia il ruolo delle diverse variabili presenti nel modello; pertanto
considerando che in precedenza è stato indicato espressamente che l’unica variabile direttamente manipolabile è la velocità del fluido v, è evidente che tale
variabile è la variabile di controllo del processo, a questo punto non resta allora
che mettere in evidenza quale sia il ruolo delle restanti variabili di ingresso.
Ebbene è conveniente attribuire il ruolo di variabili di disturbo alle grandezze
ϑv , ϑli , infatti per esse non è possibile fissare arbitrariamente un valore.
Avendo chiarito quale sia il ruolo delle diverse variabili si può passare ora ad
analizzare le risposte del processo a fronte delle varie sollecitazioni rispetto ai
diversi ingressi, il tutto sempre tenendo presente che in alcuni casi si analizzerà
la risposta all’ingresso di controllo mentre in altri quello che si otterrà è la risposta a delle azioni disturbanti; tale circostanza è necessaria ad attribuire un peso
diverso ai risultati che si otterranno nelle diverse analisi. Prima di procedere
con le suddette analisi si tenga presente che ognuna di esse verrà effettuata a
partire da diverse condizioni operative dell’impianto, condizioni che si ottengono partendo da una opportuno stato iniziale e mantenendo costanti di volta in
volta alcune delle variabili di ingresso.
m
3.1.1 Analisi del legame v − ϑlu
Il primo legame che vuole studiare è quello tra la variabile di uscita e la variabile di controllo, certamente tale legame è molto importante in quanto esso
influenza la procedura di sintesi della legge di controllo che verrà mostrata più
avanti.
29
3.1 caratterizzazione input-output dello scambiatore di calore
3.1.1.1
Conduzione dell’impianto nello stato nominale
La prima delle analisi che si vogliono effettuare consiste nel condurre l’impianto da una condizione operativa iniziale al punto di lavoro nominale. Per effettuare tale analisi si sollecita l’impianto con un ingresso di velocità del liquido
v pari al valore nominale ossia 0.31 m s−1 a partire dalle seguenti condizioni
iniziali:
• Temperatura del rame: 90◦ C
• Temperatura del liquido: 65◦ C
• Temperatura del vapore: 120◦ C
tali condizioni operative, significative dal punto di vista fisico, si traducono in
dei valori opportuni dello stato iniziale del processo e dei valori temporali dei
vari ingressi. In particolare si ha:
• xi (0) = 90◦ C ∀ i ∈ {36, . . . , 70} elementi della mesh del rame
• xi (0) = 65◦ C ∀ i ∈ {1, . . . , 35} elementi della mesh le liquido
• x71 (0) = 65◦ C variabile di stato della termocoppia
• ϑli (t) = 65◦ C ∀ t ≥ 0
• ϑv (t) = 120◦ C ∀ t ≥ 0
m ( t ) a fronte dell’ingresso di velocità del valore
L’andamento della variabile ϑlu
pari al valore nominale, è riportato in (Figura 9). Dall’andamento in (Figura 9)
98
94
90
θ lmu [ ◦ C ]
86
82
78
74
70
66
62
0
5
10
15
time [s]
Figura 9.: Risposta ad un gradino di velocità pari a 0,31 m s−1 a partire dalla condizione
indicata in 3.1.1.1
si può notare che la temperatura del fluido in uscita non raggiunge esattamente
il valore di 98 ◦C in quanto è presente l’errore introdotto dalla approssimazione
agli elementi finiti, si può inoltre notare che l’andamento non presenta né ritardi né oscillazioni.
30
3.1 caratterizzazione input-output dello scambiatore di calore
Giunti a questo punto si vogliono effettuare ora delle analisi volte a mettere in
luce quale sia il comportamento del sistema a fronte di piccole variazioni, del
valore della velocità del liquido in ingresso, intorno al proprio valore nominale.
Ebbene questo tipo di analisi è molto importante dal punto di vista controllistico in quanto, come ben noto, sollecitando un qualunque sistema dinamico
non lineare con degli ingressi di ampiezza opportunamente piccola è possibile
ottenere in uscita una risposta principalmente caratterizzata dalla componente
lineare del sistema stesso, ottenendo in tal modo seppur in maniera qualitativa,
delle indicazioni sulle difficoltà che si riscontreranno nella sintesi di una opportuna legge di controllo. Quanto appena detto di fatto può essere osservato conducendo due tipi di analisi. La prima analisi consiste nel valutare la distorsione
armonica percentuale (T.H.D. %) della risposta del sistema, a fronte di una variazione sinusoidale di ampiezza opportuna della variabile v, se il sistema fosse
lineare la distorsione armonica della risposta in uscita sarebbe chiaramente nulla, ovvero conducendo questa analisi si può ottenere una misura quantitativa
dell’influenza dei termini non lineari, nella risposta del sistema. Per meglio
comprendere quale sia il comportamento del sistema, l’analisi suddetta viene
condotta per diversi valori di pulsazione dell’ingresso e per diversi valori di
ampiezza del medesimo. Si riportano i risultati in (Figura 10) e (Figura 11).
Osservando le figure di cui sopra si nota che all’aumentare dell’ampiezze dell’ingresso la distorsione armonica dell’uscita aumenta, tuttavia la distorsione
risulta comunque contenuta.
31
3.1 caratterizzazione input-output dello scambiatore di calore
∆θ lmu(t ) [ ◦ C ]
1
0
−1
−2
10
15
20
25
30
35
40
30
35
40
time [s]
v (t ) [ m s −1]
0.4
0.35
0.3
0.25
0.2
0.15
0.1
10
15
20
25
time [s]
(a) Andamento della velocità del fluido e della relativa variazione di temperatura del
liquido in uscita a partire da condizioni nominali
T H D = 9.00%
1.6
1.4
1.2
Mag
1
0.8
0.6
0.4
0.2
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Harmonic order
(b) Analisi armonica dell’andamento della variazione di temperatura del liquido in
uscita
Figura 10.: Sollecitazione dell’impianto con un ingresso di velocità
v(t) = v̄ + 0.1 sin(0.6 t) [m s−1 ]
In virtù della scarsa influenza dei termini non lineari nell’uscita dell’impianto, mostrata attraverso l’analisi appena condotta, risulta significativo ottenere
una buona approssimazione lineare dell’impianto in grado di predire in maniera sufficientemente accurata l’evoluzione temporale della variazione rispetm ( t ) a fronte di piccole
to al valore nominale, della temperatura in uscita ∆ϑlu
variazioni della velocità del liquido in ingresso. Una prima metodologia che
può essere utilizzata a tale scopo si basa sullo studio delle risposta indiciale
dell’impianto ottenuta per variazioni a gradino della velocità del fluido di ampiezza modesta, a titolo d’esempio si osservi la risposta mostrata in (Figura 12).
Più precisamente a partire dalla misura della risposta indiciale si identificano,
attraverso l’utilizzo di metodi numerici iterativi basati sulla minimizzazione
dello scarto quadratico medio, i parametri di opportuni modelli semplificati,
in modo tale da approssimare la risposta indiciale del processo stesso; tale
procedura di fatto rappresenta un metodo per descrivere il comportamento dinamico del processo, per piccole variazioni dell’ingresso, attraverso un numero
sufficientemente ridotto di parametri. Ebbene un primo modello che può es-
32
3.1 caratterizzazione input-output dello scambiatore di calore
∆θ lmu(t) [ ◦ C ]
4
2
0
−2
−4
10
15
20
25
30
35
40
30
35
40
time [s]
0.5
v [m s −1]
0.4
0.3
0.2
0.1
0
−0.1
10
15
20
25
time [s]
(a) Andamento della velocità del fluido e della relativa variazione di temperatura del liquido in uscita a partire da condizioni
nominali
T H D = 16.18%
3.5
3
2.5
Mag
2
1.5
1
0.5
0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Harmonic order
(b) Analisi armonica dell’andamento della variazione di temperatura del liquido in uscita
Figura 11.: Sollecitazione dell’impianto con un ingresso di velocità
v(t) = v̄ + 0.2 sin(0.6 t) [m s−1 ]
1.4
1.2
∆θ lmu(t ) [ ◦ C ]
1
0.8
0.6
0.4
0.2
0
0
10
20
30
40
50
60
70
80
Ti me [ s]
m ( t ) a fronte della variazioFigura 12.: Andamento della variazione di temperatura ∆ ϑlu
ne della velocità del fluido rispetto al valore nominale:
∆ v(t) = −0.0203 δ−1 (t − 20) [m s−1 ]
33
3.1 caratterizzazione input-output dello scambiatore di calore
sere considerato è conosciuto in letteratura come First Order Plus Dead Time
(F.O.P.D.T.) ed è descritto dalla funzione di trasferimento seguente:
P(s) =
Kp
e−Td s
1+τs
(F.O.P.D.T.)
(3.1)
La struttura del modello (3.1) è motivata da considerazioni fisiche, in particolare il modello (F.O.P.D.T. cerca di approssimare la risposta indiciale del
processo caratterizzando il comportamento statico del medesimo attraverso il
guadagno statico K p e caratterizzando invece il comportamento dinamico, tenendo conto della “inerzia” del processo modellata attraverso la costante di
tempo τ ed il ritardo Td . Pertanto avendo chiarito quale sia la giustificazione
della scelta del modello (F.O.P.D.T.) ed avendo messo in luce quale sia la procedura di determinazione dei parametri dei medesimi, al fine di valutare la bontà
del modello approssimante che è stato determinato, si riportano i parametri del
modello di tipo (F.O.P.D.T.) identificati dalla risposta (Figura 12), ed i relativi andamenti delle risposte indiciali confrontati con quella dell’impianto (Figura 13).
Modello (F.O.P.D.T.)
• K p = −66.385 [◦ C m−1 s]
• τ = 3.4801 [s]
• Td = 2.0473 [s]
P(s) =
−66.385 −2.0473 s
e
1 + 3.4801 s
(3.2)
Osservando la (Figura 13) si nota che il modello approssimato (3.2) commette
un errore nella rappresentazione della risposta indiciale dell’impianto abbastanza evidente. Al fine di ridurre tale errore è possibile fare ricorso ad un ulteriore
modello approssimato conosciuto in letteratura come Second Order Real Poles
Plus Dead Time (S.O.R.P.P.D.T.) e descritto della seguente funzione di trasferimento:
P(s) =
Kp
e−Td s
(1 + τ1 s) (1 + τ2 s)
(S.O.R.P.P.D.T.)
(3.3)
34
3.1 caratterizzazione input-output dello scambiatore di calore
v → θ lmu [ ◦ C ]
1.4
1.2
∆θ lmu(t ) [ ◦ C ]
1
0.8
0.6
0.4
0.2
F. O. P. D. T.
R eal P l ant
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 13.: Andamento reale ed approssimato con il modello (F.O.P.D.T.) della variam ( t ) rispetto al valore nominale a fronte della
zione di temperatura ∆ ϑlu
variazione della velocità del fluido rispetto al valore nominale ∆ v(t) =
−0.0203 δ−1 (t − 20) [m s−1 ]
Si noti che il modello (S.O.R.P.P.D.T.) di fatto prevede l’aggiunta di una ulteriore costante di tempo al modello (F.O.P.D.T.. In tal modo si cerca di ottenere
una maggiore aderenza tra la risposta dell’impianto e quella del modello approssimante. Si riportano di seguito i parametri identificati del modello di tipo
(S.OR.P.P.D.T.) identificati dalla risposta (Figura 12), ed i relativi andamenti delle risposte indiciali confrontati con quella dell’impianto (Figura 14).
Modello (S.O.R.P.P.D.T.) I parametri del modello (S.O.R.P.P.D.T.) identificati a
partire dalla risposta indiciale dell’impianto in (Figura 12) sono i seguenti:
• K p = −67.38 [◦ C m−1 s]
• τ1 = 2.2427 [s]
• τ2 = 2.2427 [s]
• Td = 0.90778 [s]
tali parametri conducono alla funzione di trasferimento seguente per il modello
in esame:
P(s) =
−67.38
e−0.90778 s
(1 + 2.2427 s) (1 + 2.2427 s)
(3.4)
In conclusione, al fine di confrontare i risultati ottenuti nella loro totalità, si
riportano in (Figura 15) le risposte indiciali dell’impianto reale, del modello
(F.O.P.D.T. e del modello (S.O.R.P.P.D.T.)
Osservando i risultati mostrati in (Figura 15) si nota che la migliore approssimazione dell’impianto è ottenuta adottando il modello (S.O.R.P.P.D.T.).
35
3.1 caratterizzazione input-output dello scambiatore di calore
v → θ lmu [ ◦ C ]
1.4
1.2
∆θ m l u(t ) [ ◦ C ]
1
0.8
0.6
0.4
0.2
S . O. R . P. P. D. T.
R eal P l ant
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 14.: Andamento reale ed approssimato con il modello (S.O.R.P.P.D.T.) delm ( t ) rispetto al valore nominale a fronla variazione di temperatura ∆ ϑlu
te della variazione della velocità del fluido rispetto al valore nominale
∆ v(t) = −0.0203 δ−1 (t − 20) [m s−1 ]
v → θ lmu [ ◦ C ]
1.4
1.2
∆θ m l u(t ) [ ◦ C ]
1
0.8
0.6
0.4
R eal P l ant
F. O. P. D. T.
0.2
S . O. R . P. P. D. T.
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 15.: Andamento reale ed approssimato con i modelli (S.O.R.P.P.D.T.) e
m ( t ) rispetto al valore no(F.O.P.D.T. della variazione di temperatura ∆ ϑlu
minale a fronte della variazione della velocità del fluido rispetto al valore
nominale ∆ v(t) = −0.0203 δ−1 (t − 20) [m s−1 ]
m
3.1.2 Analisi del legame ϑli − ϑlu
In questa sezione si vogliono effettuare delle analisi sul legame tra la temperatura misurata del fluido in uscita e la temperatura del fluido in ingresso, tale
legame, fissata la velocità del fluido in ingresso v, risulta di tipo lineare come
36
3.1 caratterizzazione input-output dello scambiatore di calore
m ◦
θ l i → θ lu
[ C]
2
∆θ lmu(t ) [ ◦ C ]
1.5
1
0.5
0
20
30
40
50
60
70
80
Ti me [ s]
m ( t ) a fronte della variaFigura 16.: Andamento della variazione di temperatura ∆ ϑlu
zione della temperatura del fluido in ingresso rispetto al valore nominale
∆ θli = 5 δ−1 (t − 20) [◦ C ]
si può osservare dalla ( 2.24 a pagina 27), pertanto nel caso in esame non verrà effettuata nessuna analisi di distorsione armonica, come mostrato nel caso
di cui sopra, viceversa la procedura di identificazione della risposta indiciale
illustrata nel paragrafo precedente verrà condotta nuovamente in questa sede,
considerando la risposta indiciale dell’impianto rispetto alla temperatura del
liquido in ingresso θli , a partire dalle condizioni nominali, mostrata in (Figura 16). Adottando la medesima procedura di identificazione illustrata nella
sezione precedente si determinano i parametri del modello in questione riportati di seguito:
Modello (F.O.P.D.T.)
• K p = 0.40474
• τ = 1.6055 [s]
• Td = 7.1951 [s]
I quali conducono alla seguente funzione di trasferimento:
P(s) =
0.40474
e−7.1951 s
(1 + 1.6055 s)
(3.5)
mentre in (Figura 17) si riporta il relativo andamento della risposta indiciale
confrontato con quella dell’impianto.
Osservando la (Figura 17) si può evincere che l’approssimazione fornita dal
modello (F.O.P.D.T.) è accettabile.
37
3.1 caratterizzazione input-output dello scambiatore di calore
m ◦
θ l i → θ lu
[ C]
2
1.5
1
0.5
F. O. P. D. T.
R eal P l ant
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 17.: Andamento reale ed approssimato con il modello (F.O.P.D.T.) della variam ( t ) rispetto al valore nominale a fronte della vazione di temperatura ∆ ϑlu
riazione della temperatura del fluido in ingresso rispetto al valore nominale
∆ θli (t) = 5 δ−1 (t − 20) [◦ C ]
m
3.1.3 Analisi del legame ϑv − ϑlu
Per concludere la caratterizzazione del comportamento Input-Output del processo sotto analisi si effettua una analisi sul legame tra la temperatura del vapore e la temperatura misurata del fluido in uscita. Le considerazioni sopra
esposte relativamente alla distorsione armonica restano valide anche in questo
m . Ebbene sollecitando anche in
caso in virtù della linearità del legame ϑv − ϑlu
questo caso l’impianto a partire dalla condizione operativa nominale, con una
variazione a gradino della temperatura del vapore, si ottiene per la risposta in
uscita quanto mostrato in (Figura 18). Ebbene anche in questo caso, al fine di
descrivere il legame in questione attraverso un numero sufficientemente ridotto
di parametri, si ripete la procedura di identificazione illustrata in precedenza.
In particolare considerando un modello di tipo (F.O.P.D.T.) si ottengono i seguenti valori per i parametri identificati:
Modello (F.O.P.D.T.)
• K p = 0.59525
• τ = 3.2133 [s]
• Td = 1.68 [s]
I quali conducono alla seguente funzione di trasferimento:
P(s) =
0.59525
e−1.68 s
(1 + 3.2133 s)
(3.6)
viceversa in (Figura 19) si riporta il relativo andamento della risposta indiciale
confrontato con quella dell’impianto.
38
3.1 caratterizzazione input-output dello scambiatore di calore
θ v → θ lmu [ ◦ C ]
3
∆θ lmu(t ) [ ◦ C ]
2.5
2
1.5
1
0.5
0
20
30
40
50
60
70
80
Ti me [ s]
m ( t ) a fronte della variaFigura 18.: Andamento della variazione di temperatura ∆ ϑlu
zione della temperatura del fluido in ingresso rispetto al valore nominale
∆ θv = 5 δ−1 (t − 20) [◦ C ]
θ v → θ lmu [ ◦ C ]
3
∆ θ lmu(t ) [ ◦ C ]
2.5
2
1.5
1
0.5
F. O. P. D. T.
R eal P l ant
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 19.: Andamento reale ed approssimato con il modello (F.O.P.D.T.) della variam ( t ) rispetto al valore nominale a fronte della vazione di temperatura ∆ ϑlu
riazione della temperatura del fluido in ingresso rispetto al valore nominale
∆ θv (t) = 5 δ−1 (t − 20) [◦ C ]
Anche in questo caso essendo l’approssimazione ottenuta tramite il modello
(F.O.P.D.T.) poco soddisfacente si ricorre ad un modello di tipo (S.O.R.P.P.D.T.)
per il quale la procedura di identificazione, basata sulla risposta indiciale in
(Figura 18), fornisce i seguenti valori dei parametri in questione:
39
3.1 caratterizzazione input-output dello scambiatore di calore
Modello (S.O.R.P.P.D.T.)
• K p = 0.59525
• τ1 = 2.0403 [s]
• τ2 = 2.0403 [s]
• Td = 0.69411 [s]
I quali conducono alla seguente funzione di trasferimento:
P(s) =
0.59525
e−0.69 s
(1 + 2.0403 s) (1 + 2.0403 s)
(3.7)
mentre la relativa risposta indiciale, confrontata con quella dell’impianto reale,
è mostrata in (Figura 20)
θ v → θ lmu [ ◦ C ]
3
∆ θ lmu(t ) [ ◦ C ]
2.5
2
1.5
1
0.5
S . O. R . P. P. D. T.
R eal P l ant
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 20.: Andamento reale ed approssimato con il modello (S.O.R.P.D.T.) della vam ( t ) rispetto al valore nominale a fronte delriazione di temperatura ∆ ϑlu
la variazione della temperatura del vapore rispetto al valore nominale
∆ θv (t) = 5 δ−1 (t − 20) [◦ C ]
40
3.1 caratterizzazione input-output dello scambiatore di calore
Anche in questo caso al fine di avere un quadro generale dei risultati ottenuti si mostrano in (Figura 21) le risposte indiciali dei modelli (F.O.P.D.T.),
(S.O.R.P.D.T.) e dell’impianto reale.
θ v → θ lmu [ ◦ C ]
3
∆ θ lmu(t ) [ ◦ C ]
2.5
2
1.5
1
R eal P l ant
F. O. P. D. T.
0.5
S . O. R . P. P. D. T.
0
20
30
40
50
60
70
80
Ti me [ s]
Figura 21.: Andamento reale ed approssimato con il modello (F.O.P.D.T.) e con il mom ( t ) rispetto al
dello (S.O.R.P.P.D.T.) della variazione di temperatura ∆ ϑlu
valore nominale a fronte della variazione della temperatura del vapore
rispetto al valore nominale ∆ θv (t) = 5 δ−1 (t − 20) [◦ C ]
41
3.2
3.2 analisi dell’impianto basata sul modello
analisi dell’impianto basata sul modello
Nella sezione precedente sono state condotte delle analisi, volte a mettere in
risalto la componente lineare dell’impianto, basate sullo studio della risposta
indiciale. Ebbene risulta evidente che tale tipo di analisi può essere condotta
anche nei casi in cui non si disponga di un modello dell’impianto, tuttavia nel
caso in esame disponendo di un modello del processo, ci si può giovare della
conoscenza di tale modello al fine di condurre delle analisi del comportamento dinamico del medesimo in maniera più sistematica; ad esempio effettuando
una linearizzazione dell’impianto attorno a diversi punti di equilibrio di interesse. A tal scopo si richiama la procedura generale di linearizzazione di un
generico sistema non lineare forzato, attorno ad un generico equilibrio.
Linearizzazione di un sistema non lineare
Si consideri un generico sistema non lineare
ẋ = f ( x, u)
(3.8)
y = h( x, u)
e sia ( x̄ , ū) un qualunque punto di equilibrio per il sistema (3.8) ossia tale che:
f ( x̄ , ū) = 0
si ponga allora:
x̃ := x − x̄
ũ := u − ū
ỹ := y − h( x̄, ū)
pertanto assumendo f ( x , u) e h( x, u) sufficientemente regolari si ha:
f ( x, u) =
∂ f ( x , u) ∂ f ( x , u) x̃ +
ũ + h.o.t.( x, u)
∂x
∂u
( x ,u)=( x̄ ,ū)
( x ,u)=( x̄ ,ū)
h( x, u) = h( x̄, ū) +
∂h( x , u) ∂h( x , u) x̃ +
ũ + h.o.t.( x, u)
∂x
∂u
( x ,u)=( x̄ ,ū)
( x ,u)=( x̄ ,ū)
Avendo indicato con h.o.t. i termini di ordine superiore al primo. Ebbene
trascurando i termini di ordine superiore e osservando che x̃˙ = ẋ si ha:
∂ f ( x , u) ∂ f ( x , u) x̃˙ =
x̃ +
ũ
∂x
∂u
( x ,u)=( x̄ ,ū)
( x ,u)=( x̄ ,ū)
ỹ =
∂h( x , u) ∂h( x , u) x̃ +
ũ
∂x
∂u
( x ,u)=( x̄ ,ū)
( x ,u)=( x̄ ,ū)
42
3.2 analisi dell’impianto basata sul modello
Il quale rappresenta l’approssimazione lineare del sistema di partenza, nell’intorno del punto ( x̄, ū).
Per comodità di scrittura di può porre:
∂ f ( x , u) ∂x
( x ,u)=( x̄ ,ū)
∂ f ( x , u) B :=
∂u
( x ,u)=( x̄ ,ū)
∂h( x , u) C :=
∂x
( x ,u)=( x̄ ,ū)
∂h( x , u) D :=
∂u
( x ,u)=( x̄ ,ū)
A :=
(3.9)
(3.10)
(3.11)
(3.12)
(3.13)
Potendo riscrivere in definitiva l’approssimazione lineare del sistema non
lineare originario come segue:
x̃˙ = A x̃ + B ũ
ỹ = C x̃ + D ũ
La procedura appena mostrata può essere condotta nel caso in esame prendendo in considerazione il modello dell’impianto (2.24), ossia nella forma:
ẋ = A(v) x + B(v) w
(3.14)
y = Cx
(3.15)
Ebbene al fine ottenere una maggiore aderenza con la notazione indicata nella
procedura di linearizzazione di cui sopra, si ponga:
u := [v , w T ] T
(3.16)
in tal modo, indicando con:
v̄, w̄, x̄
rispettivamente i valori della velocità del fluido, del vettore dei disturbi e dello
stato, assunti nel generico punto di lavoro, si ottiene per l’impianto, un modello
linearizzato, attorno al generico punto di lavoro, nella forma:
x̃˙ = A x̃ + B ũ
(3.17)
ỹ = C x̃
(3.18)
ove le matrici sopra indicate sono definite come segue:
A := A(v̄)
"
∂ A ∂ B x̄ +
v̄
B :=
∂v x̄, v̄, w̄
∂v x̄, v̄, w̄
(3.19)
#
B(v̄)
(3.20)
Pertanto resta solo da chiarire quale o quali siano i punti di lavoro scelti per
condurre l’analisi basata sulla linearizzazione dell’impianto. A tal scopo al
43
3.2 analisi dell’impianto basata sul modello
fine di effettuare una analisi, il più possibile generale, del comportamento dell’impianto vengono presi in considerazione cinque punti di lavoro, compreso il
nominale, scelti casualmente nella regione ammissibile in (Figura 5 a pagina 23).
In (Tabella 1) sono riportati i vari punti operativi presi in considerazione. PerTabella 1.: Valore della temperatura del liquido in ingresso e del vapore per i punti di
lavoro scelti. Il punto di lavoro nominale viene indicato con il numero 0.
Punto di lavoro
#
0
1
2
3
4
ϑli
◦
( C)
65
75
50
75
50
ϑv
◦
( C)
120
140
140
105
105
ciò, come già indicato in precedenza, al fine di determinare i diversi stati di
equilibrio dell’impianto, dapprima si fissa la temperatura di uscita del liquido
al valore nominale ossia, 98 ◦ C, successivamente si impone al vettore dei segnali esogeni w uno dei valori riportati in (Tabella 1), di seguito si determina la
velocità del fluido all’equilibrio v̄, per i diversi punti di lavoro, tramite la ( 2.14
a pagina 22) ed infine si calcola lo stato di equilibrio x̄ risolvendo la seguente
equazione algebrica nell’incognita x̄
A(v̄) x̄ + B(v̄) w = 0
(3.21)
La soluzione della quale porge la seguente espressione generale per il valore
dello stato di equilibrio x̄
x̄ = −A(v̄)−1 B(v̄) w̄
(3.22)
A vantaggio di quanto appena esposto, si riportano nei paragrafi seguenti le
analisi condotte, attraverso la procedura di linearizzazione effettuata attorno ai
punti di lavoro indicati in (Tabella 1).
3.2.1 Analisi del comportamento dinamico del processo attorno al punto di
lavoro nominale
Lo studio che ci si appresta a mostrare in questo paragrafo è relativo al comportamento dinamico del processo in un intorno del punto di lavoro nominale. Ebbene tale punto di lavoro come già detto più volte è caratterizzato dai seguenti
valori per le variabili di processo:
• ϑli = 98◦ C
• ϑv = 120◦ C
• v = 0.31m.s−1
44
3.2 analisi dell’impianto basata sul modello
a partire da tali informazioni e seguendo la procedura di linearizzazione sopra
indicata si ottiene il modello linearizzato dell’impianto attorno al punto di lavoro nominale.
Ottenuto il modello linearizzato del processo, si può passare ad effettuare un’analisi in frequenza del medesimo. In virtù di ciò è conveniente mettere in
evidenza l’espressione della risposta forzata del sistema, nel dominio della variabile complessa s, ebbene indicando con U (s) e Y (s) rispettivamente, la trasformata di Laplace dell’ingresso vettoriale u(t) definito nella (3.16) e dell’uscita
considerata, si ottiene, per la risposta forzata l’espressione seguente:
Y ( s ) = C ( s I − A ) −1 B U ( s )
(3.23)
In prima battuta, essendo il sistema a più ingressi, sembrerebbe opportuno effettuare un’analisi del massimo valore singolare della matrice di trasferimento
tra gli ingressi e l’uscita considerata; tuttavia avendo considerato la temperatura del liquido di ingresso e la temperatura del vapore, come due disturbi è
conveniente analizzare i singoli legami separatamente, considerando l’azione
di un ingresso alla volta. Per questo è necessario riscrivere la risposta forzata
(3.23), mettendo in evidenza l’azione distinta dei tre ingressi, in tal modo si
ottiene l’espressione seguente:
Y (s) = C (s I − A)−1 (b1 v(s) + b2 w1 (s) b3 w2 (s))
(3.24)
Ove w1 (s) , w2 (s) sono le due componenti del vettore dei disturbi definito in
( 2.23 a pagina 26), mentre le tre matrici degli ingressi sopra indicate sono
definite come segue:
"
#
∂ A ∂ B b1 :=
x̄ +
v̄
∂v x̄, v̄, w̄
∂v x̄, v̄, w̄
(3.25)
b2 := B1 (v̄)
b3 := B2 (v̄)
Avendo indicato con B1 e B2 rispettivamente la prima e la seconda colonna
della matrice B(v̄). La relazione (3.24), al fine di snellire la notazione, può
essere riscritta come segue:
Y (s) = Fy v (s) v(s) + Fy w1 (s) w1 (s) + Fy w2 (s) w2 (s)
dove si sono indicate rispettivamente:
Fy v (s) = C (s I − A)−1 b1
Fy w1 (s) = C (s I − A)−1 b2
Fy w2 (s) = C (s I − A)
−1
b3
(3.26)
(3.27)
(3.28)
Avendo distinto le tre azioni è possibile ottenere i diagrammi di Bode dei tre
legami, i quali vengono mostrati di seguito in (Figura 22), (Figura 23) e (Figura 24). Dalla analisi appena condotta si può stabilire, osservando il diagramma
mostrato in (Figura 23), che il motivo per il quale è stato osservato nell’analisi
nel tempo dell’impianto, un notevole ritardo nella risposta indiciale nel legame
m , è giustificato dall’elevato grado relativo della risposta in questione,
ϑli − ϑlu
tale caratteristica del legame in esame di fatto discende dalla presenza di un
45
3.2 analisi dell’impianto basata sul modello
Bode Diagram
60
Magnitude (dB)
40
20
0
−20
−40
−60
Phase (deg)
180
135
90
45
0
−3
10
−2
−1
10
0
1
10
10
Frequency (rad/sec)
10
2
10
Figura 22.: Diagramma di Bode della funzione di trasferimento Fy v (s)
Bode Diagram
Magnitude (dB)
0
−200
−400
−600
−800
Phase (deg)
−1000
0
−360
−720
−1080
−1440
−1800
−2160
−2520
−2880
−3240
−3600
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 23.: Diagramma di Bode della funzione di trasferimento Fy w1 (s)
ritardo di trasporto nell’impianto reale, il quale si manifesta nell’approssimazione ad elementi finiti attraverso un alto grado relativo.
Una ulteriore analisi che può essere condotta, disponendo delle funzioni di
trasferimento dei tre legami in questione, è la valutazione della bontà delle
46
3.2 analisi dell’impianto basata sul modello
Bode Diagram
Magnitude (dB)
0
−50
−100
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
0
10
1
10
Frequency (rad/sec)
10
Figura 24.: Diagramma di Bode della funzione di trasferimento Fy w2 (s)
approssimazioni, dei legami dell’impianto, fornite dai modelli (F.O.P.D.T.) e
(S.O.R.P.D.T.) e mostrate in (Capitolo 3.1 a pagina 29). Pertanto in (Figura 25),
(Figura 26) e (Figura 27) vengono mostrati dei confronti tra le tre funzioni di risposta armonica e le relative migliori approssimazioni ottenute nel (Capitolo 3.1
a pagina 29).
Bode Diagram
Re al P l ant
Magnitude (dB)
40
S. O. R. P. P. D . T.
20
0
−20
−40
Phase (deg)
180
135
90
45
0
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
Figura 25.: Diagramma di Bode della funzione di trasferimento Fy v (s) e della sua
approssimante (3.3)
47
3.2 analisi dell’impianto basata sul modello
Bode Diagram
10
F. O. P. D . T
Magnitude (dB)
0
Re al P l ant
−10
−20
−30
−40
−50
360
Phase (deg)
0
−360
−720
−1080
−1440
−1800
−2160
−2
10
−1
0
10
1
10
10
Frequency (rad/sec)
Figura 26.: Diagramma di Bode della funzione di trasferimento Fy w1 (s) e della sua
approssimante (3.5)
Bode Diagram
20
Magnitude (dB)
Re al P l ant
0
S. O. R. P. P. D . T.
−20
−40
−60
90
Phase (deg)
0
−90
−180
−270
−360
−450
−540
−2
10
−1
0
10
10
1
10
Frequency (rad/sec)
Figura 27.: Diagramma di Bode della funzione di trasferimento Fy w2 (s) e della sua
approssimante (3.7)
Osservazione 3.1:
Osservando le (Figura 25), (Figura 26) e (Figura 27) si nota che le approssimazioni
ottenute tramite identificazione dei parametri della risposta indiciale sono soddisfacenti
solamente a bassa frequenza ciò è deriva in parte della natura semplificata dei modelli
48
3.2 analisi dell’impianto basata sul modello
approssimanti considerati.
3.2.2 Analisi del comportamento dinamico del processo attorno ad alcuni punti
di lavoro ammissibili
In questo paragrafo si mostreranno solamente i diagrammi di Bode delle tre
funzioni di trasferimento Fy v (s), Fy w1 (s), Fy w2 (s) ottenute linearizzando il modello dello scambiatore di calore attorno ai punti di lavoro indicati in (Tabella 1).
Punto di lavoro 1 : θli = 75 ◦C, θv = 140 ◦C
In (Figura 28),(Figura 29) e (Figura 30) si mostrano rispettivamente i diagrammi
di Bode delle funzioni di trasferimento Fy v (s), Fy w1 (s), Fy w2 (s) relativamente al
punto di lavoro di cui sopra.
Bode Diagram
40
Magnitude (dB)
20
0
−20
−40
−60
−80
Phase (deg)
180
135
90
45
0
−45
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 28.: Diagramma di Bode della funzione di trasferimento Fy v (s) per il punto di
lavoro θli = 75 ◦C θv = 140 ◦C
49
3.2 analisi dell’impianto basata sul modello
Bode Diagram
200
Magnitude (dB)
0
−200
−400
−600
−800
−1000
Phase (deg)
0
−720
−1440
−2160
−2880
−3600
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 29.: Diagramma di Bode della funzione di trasferimento Fy w1 (s) per il punto di
lavoro θli = 75 ◦C θv = 140 ◦C
Bode Diagram
Magnitude (dB)
0
−50
−100
−150
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 30.: Diagramma di Bode della funzione di trasferimento Fy w2 (s) per il punto di
lavoro θli = 75 ◦C θv = 140 ◦C
50
3.2 analisi dell’impianto basata sul modello
Punto di lavoro 2 : θli = 50 ◦C θv = 140 ◦C
In (Figura 31), (Figura 32) e (Figura 33) si mostrano rispettivamente i diagrammi
di Bode delle funzioni di trasferimento Fy v (s), Fy w1 (s), Fy w2 (s) relativamente al
punto di lavoro di cui sopra.
Bode Diagram
60
Magnitude (dB)
40
20
0
−20
−40
−60
Phase (deg)
180
135
90
45
0
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 31.: Diagramma di Bode della funzione di trasferimento Fy v (s) per il punto di
lavoro θli = 50 ◦C θv = 140 ◦C
51
3.2 analisi dell’impianto basata sul modello
Bode Diagram
Magnitude (dB)
0
−200
−400
−600
−800
−1000
Phase (deg)
0
−720
−1440
−2160
−2880
−3600
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 32.: Diagramma di Bode della funzione di trasferimento Fy w1 (s) per il punto di
lavoro θli = 50 ◦C θv = 140 ◦C
Bode Diagram
Magnitude (dB)
0
−50
−100
−150
45
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 33.: Diagramma di Bode della funzione di trasferimento Fy w2 (s) per il punto di
lavoro θli = 50 ◦C θv = 140 ◦C
52
3.2 analisi dell’impianto basata sul modello
Punto di lavoro 3 : θli = 75 ◦C θv = 105 ◦C
In (Figura 34), (Figura 35) e (Figura 36) si mostrano rispettivamente i diagrammi
di Bode delle funzioni di trasferimento Fy v (s), Fy w1 (s), Fy w2 (s) relativamente al
punto di lavoro di cui sopra.
Bode Diagram
Magnitude (dB)
40
20
0
−20
−40
−60
−80
Phase (deg)
180
135
90
45
0
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 34.: Diagramma di Bode della funzione di trasferimento Fy v (s) per il punto di
lavoro θli = 75 ◦C θv = 105 ◦C
53
3.2 analisi dell’impianto basata sul modello
Bode Diagram
200
Magnitude (dB)
0
−200
−400
−600
−800
−1000
Phase (deg)
0
−720
−1440
−2160
−2880
−3600
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 35.: Diagramma di Bode della funzione di trasferimento Fy w1 (s) per il punto di
lavoro θli = 75 ◦C θv = 105 ◦C
Bode Diagram
Magnitude (dB)
0
−50
−100
−150
45
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 36.: Diagramma di Bode della funzione di trasferimento Fy w2 (s) per il punto di
lavoro θli = 75 ◦C θv = 105 ◦C
54
3.2 analisi dell’impianto basata sul modello
Punto di lavoro 4 : θli = 50 ◦C θv = 105 ◦C
In (Figura 37), (Figura 38) e (Figura 39) si mostrano rispettivamente i diagrammi
di Bode delle funzioni di trasferimento Fy v (s), Fy w1 (s), Fy w2 (s) relativamente al
punto di lavoro di cui sopra.
Bode Diagram
Magnitude (dB)
50
0
−50
−100
Phase (deg)
180
135
90
45
0
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 37.: Diagramma di Bode della funzione di trasferimento Fy v (s) per il punto di
lavoro θli = 50 ◦C θv = 105 ◦C
Infine di riportano in (Figura 40), (Figura 41) e (Figura 42) i diagrammi di
Bode determinati linearizzando l’impianto intorno a tutti i punti di lavoro
considerati.
55
3.2 analisi dell’impianto basata sul modello
Bode Diagram
200
Magnitude (dB)
0
−200
−400
−600
−800
−1000
Phase (deg)
0
−720
−1440
−2160
−2880
−3600
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 38.: Diagramma di Bode della funzione di trasferimento Fy w1 (s) per il punto di
lavoro θli = 50 ◦C θv = 105 ◦C
Bode Diagram
Magnitude (dB)
0
−50
−100
−150
45
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 39.: Diagramma di Bode della funzione di trasferimento Fy w2 (s) per il punto di
lavoro θli = 50 ◦C θv = 105 ◦C
56
3.2 analisi dell’impianto basata sul modello
Bode Diagram
F0
Magnitude (dB)
50
F1
F2
0
F3
F4
−50
Phase (deg)
180
135
90
45
0
−45
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 40.: Diagrammi di Bode della funzione di trasferimento Fy v (s) per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma
sono indicate con Fi (s) con i ∈ {0, 1, 2, 3, 4} a seconda del punto di
lavoro considerato, in maniera consistente alla numerazione riportata in
(Tabella 1)
Bode Diagram
200
F 0 (s)
Magnitude (dB)
0
F 1 (s)
−200
F 2 (s)
−400
F 3 (s)
−600
F 4 (s)
−800
−1000
Phase (deg)
0
−720
−1440
−2160
−2880
−3600
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 41.: Diagrammi di Bode della funzione di trasferimento Fy w1 (s) per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma
sono indicate con Fi (s) con i ∈ {0, 1, 2, 3, 4} a seconda del punto di
lavoro considerato, in maniera consistente alla numerazione riportata in
(Tabella 1)
57
3.2 analisi dell’impianto basata sul modello
Bode Diagram
50
F 0 (s)
Magnitude (dB)
0
F 1 (s)
F 2 (s)
−50
F 3 (s)
−100
F 4 (s)
−150
45
Phase (deg)
0
−45
−90
−135
−180
−225
−270
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 42.: Diagrammi di Bode della funzione di trasferimento Fy w2 (s) per i diversi punti di lavoro. Le diverse funzioni di trasferimento nel diagramma
sono indicate con Fi (s) con i ∈ {0, 1, 2, 3, 4} a seconda del punto di
lavoro considerato, in maniera consistente alla numerazione riportata in
(Tabella 1)
58
Parte II.
Controllo dello scambiatore di
calore
59
4
I
4.1
C O N T R O L LO D E L LO S C A M B I ATO R E
DI CALORE
descrizione del problema di controllo
l problema che si vuole affrontare in questa parte è il problema della
regolazione asintotica della temperatura misurata del fluido in uscita
m , anche a fronte di variazioni della condizione operativa
dallo scambiatore ϑlu
del processo, ossia a fronte di variazioni ammissibili delle temperature del liquido in ingresso ϑli e del vapore ϑv . In particolare il problema appena esposto,
con riferimento al modello dello scambiatore di calore indicato in ( 2.24 a pagina 27), può essere formalizzato come segue:
Problema 4.1 (Problema di regolazione asintotica della temperatura):
Sia assegnato il modello del processo P nella forma:



 ẋ (t) = A(v(t)) x (t) + B(v(t)) w(t)
y(t) = C x (t)


 x (t ) = x
0
(4.1)
0
dove:
x0 ∈ R71
t 0 ∈ R+
x (t) ∈ R71
C ∈ R1×71
u(·) : R+ → R
w(·) : R+ → R2
Si supponga inoltre che x (t) e d(t) non siano misurabili. Trovare se esiste un controllore
u(·) tale che:
1.
lim ||y(t) − y∗ || = 0 ∀ y∗ ∈ R, ∀ w∗ ∈ Θy∗ , ∀ x0 ∈ R71
t→+∞
2.
x (·) , u(·) ∈ L∞
3. Il controllo u(·) è tale che se:
x0 : y ( t0 ) = y ∗
60
(4.2)
allora deve risultare:
4.2 problema di regolazione asintotica
y ( t ) = y ∗ ∀ t ≥ t0
Osservando il problema appena esposto si può osservare che la specifica
(1) riassume la proprietà di regolazione asintotica in maniera più formale, la
specifica (2) impone che sia lo stato che il controllo siano limitati ed infine
la (3) impone che se l’uscita all’istante iniziale è uguale alla temperatura di
riferimento y∗ allora essa deve rimanere tale identicamente. E’ evidente che la
struttura del problema esposto è del tutto generale e prescinde dalla realtà fisica
in cui il problema centrale della presente trattazione di colloca. Pertanto, al fine
di mostrare come sia possibile risolvere il problema (4.1), si farà riferimento ad
un problema del tutto generale. Successivamente, una volta chiarito quali siano
le diverse strategie per la soluzione del problema di regolazione, si mostreranno
i risultati ottenuti applicando tali metodologie al caso in esame.
4.2
problema di regolazione asintotica
Il problema di controllo che si vuole affrontare si colloca nell’ambito del problema di regolazione [2]. Tale problema rappresenta un caso particolare del
problema di inseguimento di traiettoria. In particolare questo problema può
essere, in generale, formalizzato come segue:
Problema 4.2 (Problema di regolazione perfetta):
Sia P un sistema dinamico assegnato del tipo:



 ẋ (t) = f ( x (t), u(t), d(t))
y(t) = h( x (t))


 x (t ) = x
0
(4.3)
0
dove:
x 0 ∈ Rn
t 0 ∈ R+
x ( t ) ∈ Rn
y(t) ∈ R
u(·) : R+ → R
d(·) : R+ → R p
Si supponga inoltre che:
f (·) ∈ C 1 (Rn × R × R p ), h(·) ∈ C 1 (Rn )
x (t) e d(t) sono non noti.
Trovare se esiste un controllo u(·) tale che:
(4.4)
(4.5)
61
4.2 problema di regolazione asintotica
1.
lim ||y(t) − y∗ || = 0 ∀ y∗ ∈ R ∀ d∗ ∈ R p , ∀ x0 ∈ Rn
t→+∞
(4.6)
2.
x (·) , u(·) ∈ L∞
3. Il controllo u(·) garantisce che se:
x0 : h ( x0 ) = y ∗
allora:
h( x (t)) = y∗ ∀t ≥ t0
Per comodità di esposizione è conveniente introdurre la seguente grandezza
ỹ(t) = y ∗ −y(t)
detta errore di regolazione in tal modo si può riscrivere la (1) nella seguente
forma equivalente:
lim ||ỹ(t)|| = 0 ∀ y∗ ∈ R, ∀ d∗ ∈ R p , ∀ x0 ∈ Rn
t→+∞
Il problema esposto di fatto consiste di determinare, se esiste, un controllo u(t) limitato in grado di annullare asintoticamente l’errore di regolazione
a partire da qualunque condizione iniziale dell’impianto con l’aggiunta che la
traiettoria nello stato, e quindi sotto le ipotesi di regolarità della funzione h( x )
indicate in (4.4) dell’uscita, sia limitata. Infine si osservi che la specifica (3) garantisce che se l’errore di regolazione ỹ è nullo nell’istante iniziale allora esso si
mantiene nullo identicamente, l’aggiunta di tale specifica conferisce l’attributo
perfetto al problema di regolazione (Problema 4.2).
Al fine di proseguire la trattazione in maniera più agevole è conveniente indicare con Σ il sistema ottenuto controllando l’impianto P attraverso il generico controllo u(·), a tale sistema ci si riferirà nel corso della trattazione con
l’accezione Sistema di controllo Σ.
4.2.1 Primi passi verso la soluzione del problema di regolazione e assunzioni
sul problema
Il problema di regolazione appena esposto equivale ad un problema di stabilizzazione di un particolare equilibrio, infatti si osservi che le traiettorie di un
qualunque sistema dinamico, purché limitate, convergono sempre su insiemi
invarianti, perciò affinché l’errore di regolazione possa annullarsi asintoticamente, risulta necessario che, fissati y∗ e d∗ , il sistema di controllo Σ possieda
un punto di equilibrio per il quale si abbia y ≡ y∗ . Ebbene in virtù della ipo-
62
4.2 problema di regolazione asintotica
tesi di stazionarietà dell’impianto P e di quanto appena asserito, si osservi che
affinché il problema di regolazione ammetta soluzione è necessario che:
∃ ( x̄ , ū) ∈ Rn+1 : f ( x̄ , ū, d∗ ) = 0 , h( x̄ ) = y∗
. Tuttavia nella presente trattazione per semplicità si assumerà non solo che
la coppia ( x̄, ū), dalla sopracitata proprietà esista, ma anche che essa sia unica.
Ossia si farà la seguente assunzione:
Assunzione 4.1:
∀ (y∗ , d∗ ) ∈ R p+1 ∃! ( x̄ , ū) ∈ Rn+1 : f ( x̄ , ū, d∗ ) = 0 , h( x̄ ) = y∗
Si noti che il fatto di aver assunto l’unicità di tale coppia non è fortemente
restrittivo, ai fini primari della seguente trattazione, in quanto tale ipotesi può
essere ragionevolmente ritenuta soddisfatta nell’ambito del controllo dei processi.1
Per comprendere quale sia la struttura della legge di controllo in grado di risolvere il problema di regolazione, si osservi che il valore del controllo ū è non
noto e dipendente da y∗ e da d∗ , per cui sotto le ipotesi di non misurabilità di d∗
è immediato comprendere che la legge di controllo u(·) sia necessariamente in
retroazione dall’uscita dell’impianto P , ma non solo infatti dovendo il sistema
di controllo ammettere un equilibrio z̄ per il quale
y ≡ y∗
sembra ragionevole utilizzare un controllore in retroazione dall’errore di regolazione ỹ, infatti in tal modo l’equilibrio z̄ è indipendente dai singoli valori di
y∗ e y. In definitiva dunque il controllore K candidato per la soluzione del
problema di regolazione è espresso nella forma seguente:
(
ż = ϕ(z, ỹ)
K:
(4.7)
u = ε(z, ỹ)
E’ evidente che il problema di regolazione appena esposto, data la natura non
lineare del sistema considerato non è di semplice soluzione. Tuttavia rilassando
alcune specifiche è possibile determinare una legge di controllo lineare in grado
di assicurare una regolazione di tipo locale. In altre parole attraverso l’uso
di un controllore LTI in retroazione da ỹ è possibile, sotto (l’Assunzione 4.1),
garantire che:
lim ||ỹ(t)|| = 0 ∀ y∗ ∈ R, ∀ d∗ ∈ R p , ∀ x0 ∈ I( x̄ ) ⊂ Rn
t→+∞
. Quando detto induce a formulare un ulteriore problema di regolazione semplificato. In particolare tale problema verrà denominato problema di regolazione locale e verrà formulato come segue:
1 Si noti che nel caso dello scambiatore di calore, questa ipotesi è certamente verifica in quanto
m ,ϑ , ϑ esiste in corrispondenza un solo valore della velocità del fluido v
fissato un valore per ϑlu
v li
e un solo stato di equilibrio x
63
4.2 problema di regolazione asintotica
Problema 4.3 (Problema di regolazione perfetta locale):
Sia P un sistema dinamico assegnato del tipo:



 ẋ (t) = f ( x (t), u(t), d(t))
x ( t0 ) = x0


 y(t) = h( x (t))
(4.8)
dove:
x ( t ) ∈ Rn
y(t) ∈ R
u(·) : R+ → R
d(·) : R+ → R p
Si supponga inoltre che:
• f (·) ∈ C 1 (Rn × R × R p ) h(·) ∈ C 1 (Rn )
•
• x (t) e d(t) sono non noti.
Trovare se esiste un controllo u(·) tale che:
1.
lim ||ỹ(t)|| = 0 ∀ y∗ ∈ R, ∀d∗ ∈ R p , ∀ x0 ∈ I( x̄ (y∗ , d∗ )) ⊂ Rn (4.9)
t→+∞
2.
x (·) , u(·) ∈ L∞
3. Il controllo u(·) garantisce che se:
x0 : h ( x0 ) = y ∗
allora:
h( x (t)) = y∗ ∀t ≥ t0
Fatte le relative precisazioni si può passare a mostrare quale sia la soluzione
del problema 4.3. In particolare la soluzione che si vuole adottare consiste
nell’utilizzo di un controllore LTI in retroazione dall’errore di regolazione ỹ del
tipo:
(
ϕ̇ = Ψ ϕ + Λ ỹ
K:
(4.10)
u = K1 ϕ + K2 ỹ
64
4.2 problema di regolazione asintotica
con:
ϕ ∈ R`
K1 ∈ R
K2 ∈ R
(4.11)
1×`
Ψ ∈ R
(4.12)
(4.13)
`×`
(4.14)
Λ ∈ R`
(4.15)
L’adozione della legge di controllo (4.11) per il controllo dell’impianto P da
luogo allo schema di controllo mostrato in (Figura 43). Ebbene al fine di stabili-
Figura 43.: Schema di controllo adottato
re sotto quali condizioni esista un controllore K che sia soluzione del problema
di regolazione locale (4.3), è conveniente enunciare alcuni risultati.
65
4.2 problema di regolazione asintotica
Lemma 4.1:
In merito al sistema in (4.3), si fissino y∗ ∈ R, d∗ ∈ R p in modo qualunque ma a patto
che valga la Assunzione 4.1 con ū 6= 0. Si osservi ora che il controllore (4.10) fornisce
la seguente dinamica per il sistema di controllo:
ẋ = f ( x, ϕ, d∗ , y∗ )
(4.16)
ϕ̇ = Ψ ϕ + Λ (y∗ − h( x ))
Ebbene se sotto l’azione di tale controllo esiste, per il sistema di controllo (A.11), uno
stato di equilibrio ( x̂ , ϕ̂) ∈ R`+n tale che:
h( x̂ ) = y∗
allora valgono le seguenti affermazioni:
• ū = K1 ϕ̂.
• la matrice Ψ ammette λ = 0 come autovalore ed inoltre ϕ̂ ∈ Ker (Ψ)
• K1T ∈
/ Im Ψ
Dimostrazione. La prova del lemma è riportata nell’Appendice A
Il lemma appena enunciato di fatto afferma che affinché il K controllore sia
soluzione del problema di regolazione (4.3) esso deve necessariamente comprendere una azione integrale.
Teorema 4.1:
In merito al sistema in (4.3), si fissino y∗ ∈ R, d∗ ∈ R p arbitrariamente a patto che
valga la Assunzione 4.1. Ebbene in si indichino ora rispettivamente:
∂ f ( x , u, d∗ ) A :=
∂x
( x ,u)=( x̄ ,ū)
∂ f ( x , u , d∗ ) B :=
∂u
( x ,u)=( x̄ ,ū)
∂h( x ) C :=
∂x (x ,u)=( x̄ ,ū)
(4.17)
dunque se:
• la terna ( A, B, C ) è stabilizzabile e rilevabile
A B
• rank
= n+1
C 0
allora esiste un controllore K nella forma (4.10) che sia soluzione del problema (4.3) ed
inoltre esso è espresso nella forma seguente:


 ϕ̇1 = λ1 ỹ
K : ϕ̇2 = Ψ̂ ϕ2 + λ2 ỹ


u = k11 ϕ1 + K12 ϕ2 + K2 ỹ
(4.18)
66
4.2 problema di regolazione asintotica
dove:
ϕ1 ∈ R
(4.19)
λ1 ∈ R
(4.21)
k11 ∈ R
(4.23)
ϕ2 ∈ R`−1
λ2 ∈ R`−1
K12 ∈ R1×`−1
Dimostrazione. La prova del teorema è riportata nell’Appendice A
(4.20)
(4.22)
(4.24)
(4.25)
Osservazione 4.1:
La trattazione sviluppata fino ad ora di per se garantisce che fissati d∗ e y∗ è possibile
progettare un regolatore in grado di annullare asintoticamente l’errore di regolazione
a patto che lo stato iniziale dell’impianto sia sufficientemente prossimo allo stato di
equilibrio che garantisce l’errore nullo di regolazione. Ciononostante nella realtà sia il
riferimento y∗ sia il vettore dei disturbi d∗ non sono noti, ciò implica che di fatto il
modello linearizzato dell’impianto non risulta univocamente determinato! Si noti che
il problema che si presenta al variare delle condizioni operative è legato solamente alla
stabilizzazione del processo allargato. Infatti la presenza dell’azione integrale conferisce
al sistema di controllo di per se un elevata robustezza in quanto tale azione rende il
punto di lavoro desiderato un equilibrio per il sistema di controllo, qualunque siano d∗
e y∗ , per cui se si riesce determinare una legge di controllo in grado di rendere il generico
equilibrio asintoticamente stabile, il problema è risolto indipendentemente dal punto di
lavoro considerato. E’ evidente però che utilizzando tecniche di stabilizzazione lineare
la dipendenza dal punto di lavoro è un fattore critico in quanto tali tecniche sono basate
sul modello linearizzato dell’impianto il quale dipende fortemente dal punto di lavoro
considerato. Fortunatamente, sopratutto nell’ambito del controllo dei processi, è noto
almeno un punto di lavoro nominale dell’impianto, perciò in virtù di tale conoscenza
è possibile determinare un modello nominale sul quale basare la procedura di sintesi
del regolatore. E’ evidente che il modello di cui si dispone in questo modo fornisce
una conoscenza parziale dell’impianto per cui il progetto delle sintesi di controllo va
effettuato garantendo che il controllore stabilizzante sia in grado di stabilizzare una
famiglia di impianti e non il singolo impianto nominale.
Pertanto attraverso l’introduzione dei due risultati di cui sopra, e partizionando il controllore K come segue:

(

ϕ̇1 = λ1 ỹ



I:


u1 = k11 ϕ1







(
(4.26)
K:
ϕ̇2 = Ψ̂ ϕ2 + λ2 ỹ

K
:
s



u2 = K12 ϕ2 + K2 ỹ








 u = u +u
2
1
si perviene ad una legge di controllo che può essere rappresentata attraverso lo
schema in (Figura 44)
67
4.2 problema di regolazione asintotica
Figura 44.: Schema di controllo con azione integrale
4.2.2 Regolazione con azione feed-back+feed-forward
Il problema di regolazione è stato risolto, almeno nella sua versione locale, con
successo. Tuttavia si noti che la soluzione ottenuta, utilizzando il compensatore K in (4.10), potrebbe causare dei problemi in caso di variazioni a gradino
del riferimento y∗, infatti in tali situazioni l’andamento della variabile di controllo u(t), a causa del legame diretto ingresso uscita presente nel controllore
K , assumerebbe un andamento discontinuo provocando delle eccessive sollecitazioni dell’attuatore. Per risolvere tale problema è possibile effettuare una
sistema di controllo comprensivo anche di un azione in avanti, ossia invece di
utilizzare come riferimento il segnale y∗ si utilizza come riferimento l’uscita
di un sistema LTI avente come ingresso il riferimento y∗ (azione in avanti), lo
schema di controllo che deriva da quanto detto è mostrato in (Figura 45).
Figura 45.: Schema di controllo con azione avanti
Per comprendere se in tali condizioni il controllore K in (A) è ancora soluzione del problema di regolazione proposto, è conveniente enunciare il seguente
risultato:
68
4.2 problema di regolazione asintotica
Proposizione 4.1:
Sia H il seguente sistema espresso nella forma:
ε̇ = H ε + η y∗
m = Lε
(4.27)
e sia ΣH il sistema di controllo derivante dall’utilizzo del sistema H come azione in
avanti per il sistema di controllo Σ, ossia:
ẋ = f ( x, u , d∗ )
ε̇ = H ε + η y∗
ϕ̇ = Ψ ϕ + Λ ( L ε − h( x ))
(4.28)
u = K1 ϕ + K2 ( L ε − h( x ))
e sia il controllore K (4.10) soluzione del problema di regolazione allora esso risolve
anche il problema di regolazione con azione in avanti se e soltanto se:
• H ≺ 0
• L H −1 η = − 1
Dimostrazione. Se il controllore (4.10) risolve il problema di regolazione allora
∃( x̄ , ϕ̄) ∈ Rn+` tale che:
f ( x̄ , ϕ̄) = 0
(4.29)
Λ ϕ̄ = 0
(4.30)
y = y∗
(4.31)
Tuttavia per l’ipotesi poste sul sistema H nell’enunciato si ha che ∃ ε̄ ∈ R tale
che:
ε̄ = y∗
Si osservi quindi che il punto di equilibrio stato ( x̄, ε̄, ϕ̄) è tale che:
f ( x, z , d∗ ) = 0
Ψ ϕ̄ = 0
(4.32)
Ossia l’azione in avanti non altera le posizione del punto di equilibrio necessario alla regolazione. Per concludere l’assunto occorre e basta osservare che
linearizzando il sistema si controllo ΣH attorno al punto ( x̄, ε̄, ϕ̄) e utilizzando
la notazione introdotta nel Teorema A.1 si ottiene il seguente sistema LTI:
  
 
x̃˙
A − B K2 C B K1 B K2 L
x̃
 ϕ̃˙  =  −Λ C
(4.33)
Ψ
Λ L   ϕ̃
˙ε̃
0
0
H
ε̃
Ebbene è evidente che tale sistema è asintoticamente stabile in quanto gli autovalori della sua matrice dinamica sono l’unione degli autovalori del sistema di
controllo Σ, i quali sono certamente a parte reale negativa dato che K è soluzione del problema di regolazione per ipotesi, e degli autovalori del sistema H i
quali sono a parte reale negativa per ipotesi. Resta ora da mostrare la necessità
delle condizioni sopra enunciate. A tal scopo si osservi che se H fosse singolare allora non esisterebbe nessun punto di equilibrio per il sistema di controllo
69
4.3 metodi basati sulla risposta armonica per la regolazione
ΣH , viceversa se fosse H 0 allora il sistema di controllo ΣH ammetterebbe
il punto di equilibrio desiderato ma esso in questo caso sarebbe instabile. Per
quanto riguarda la necessità della seconda delle condizioni, si noti che se essa
non fosse vera allora, come conseguenza del Lemma A.1, il sistema di controllo
ammetterebbe come unico punto di equilibrio un punto per i quale si avrebbe:
y = − L H −1 η y ∗ 6 = y ∗
Concludendo in tal modo la dimostrazione dell’assunto.
4.2.3 Controllo Observer-Based
Per quanto è stato mostrato nelle sezione precedente, al fine di risolvere il
problema di regolazione è sufficiente progettare un compensatore in grado di
stabilizzare il processo allargato. Detto ciò osservando la struttura della legge
di controllo stabilizzante, è immediato constatare che la strategia di controllo
che può essere utilizzata allo scopo di stabilizzare il processo allargato può
essere ad esempio di tipo observer-based. Infatti, come è noto, dal principio di
separazione, si può progettare un compensatore stabilizzante in due passi. In
particolare indicando con η lo stato del processo allargato si può sintetizzare
una legge di controllo stabilizzante nel modo seguente:
1. si determina un controllo stabilizzante in retroazione dallo stato η
2. si sintetizza un osservatore asintotico dello stato del processo allargato e
si utilizza la stima dello stato η̂ fornita dall’osservatore per effettuare la
retroazione dallo stato determinata al passo precedente.
La procedura observer-based indicata tuttavia richiede la conoscenza dei parametri dell’impianto da controllare. Come già detto una conoscenza parziale
dell’impianto può essere ottenuta considerando un opportuno modello nominale. In questo modo è evidente che il progetto del compensatore deve essere
effettuato utilizzando delle metodologie che garantiscono una certa robustezza. Ebbene una tecnica di controllo observer-based in grado di garantire una
buona robustezza è il ben noto controllo LQG. Chiaramente l’utilizzo di una
metodologia di tipo observer-based essendo basata sul modello del processo è
intrinsecamente poco robusta. Per fare fronte a ciò nella sezione seguente si
mostrerà un ulteriore metodologia di sintesi per il compensatore, attraverso la
quale si può garantire un maggiore grado di robustezza.
4.3
metodi basati sulla risposta armonica per la
regolazione
In questo capitolo si vuole illustrare come sia possibile ricorrere a tecniche di
sintesi basate sulla risposta armonica per sintetizzare una legge di controllo
L.T.I. in grado di risolvere il problema delle regolazione locale di un sistema
non lineare.
Nello studio dei sistemi di controllo lineari, il criterio di Nyquist rappresenta uno strumento di notevole importanza soprattutto in fase di sintesi di un
70
4.3 metodi basati sulla risposta armonica per la regolazione
controllore, in quanto attraverso l’utilizzo di tale criterio è possibile progettare un controllore che assicuri la stabilità asintotica del sistema di controllo
facendo riferimento alla risposta armonica del sistema di controllo ad anello
aperto, sfortunatamente,nel controllo dei sistemi non lineari non è possibile fare ricorso, di per se al criterio di Nyquist. Tuttavia nell’ambito della presente
trattazione, essendo il progetto del regolatore fondato su tecniche lineari basate quest’ultime, sul teorema di approssimazione lineare di Lyapunov, sembra
ragionevole pensare che si possa fare anche in questo, caso uso del criterio di
Nyquist per garantire le proprietà di stabilità asintotica di un dato stato di equilibrio. Pertanto si osservi che sotto le ipotesi di utilizzo della legge di controllo
(controllore lineare), come già detto, al fine di studiare le proprietà di stabilità
esponenziale locale del punto di equilibrio ( x̄, ϕ¯1 ) per il sistema di controllo,
è sufficiente studiare la stabilità asintotica dell’approssimazione lineare di esso
(controllo) ottenuta linearizzando quest’ultimo attorno al punto ( x̄, ϕ¯1 ). Tale
approssimazione lineare di fatto può essere pensata come il sistema LTI che
si otterrebbe interconnettendo opportunamente il controllore (4.10), di per se
già lineare, e l’approssimazione lineare del processo. Attraverso l’ausilio di
questa osservazione è possibile utilizzare tutti gli strumenti noti della teoria lineare del controllo, basati sullo studio della risposta armonica considerando le
funzioni di trasferimento del controllore, la quale risulta ben definita essendo
quest’ultimo lineare, e la funzione di trasferimento dell’approssimazione lineare dell’impianto. Inoltre per quanto è stato detto in precedenza essendo a priori
non noti y∗ e d∗ non è nota univocamente la funzione di trasferimento del processo allargato.
4.3.1 Controllo PID based
Una delle leggi di controllo utilizzate soprattutto nell’ambito del controllo dei
processi è certamente la legge di controllo P.I.D.. Ebbene per quanto esposto
nei risultati presentati fin ora, sotto opportune condizioni sull’impianto da controllare, è sempre possibile determinare una legge di controllo L.T.I. in grado
di risolvere il problema delle regolazione, nella forma:
ϕ̇1 = λ1 ỹ
ϕ̇2 = Ψ̂ ϕ2 + λ2 ỹ
(4.34)
u = k11 ϕ1 + K12 ϕ2 + K2 ỹ
Pertanto per mostrare che un regolatore P.I.D. sia in grado di risolvere il problema di regolazione è sufficiente mostrare che tale controllore è un caso particolare di quella sopra indicata. A riguardo si riporta l’espressione della funzione
di trasferimento di un regolatore P.I.D. causale:
1
s
+ Kd
s
1+τs
K p , Ki , K d ∈ R τ ∈ R+
G ( s ) = K p + Ki
(4.35)
Per mostrare l’equivalenza di tale legge di controllo con quella più generale, è
conveniente determinare una rappresentazione con lo spazio di stato di tale fun-
71
4.3 metodi basati sulla risposta armonica per la regolazione
zione di trasferimento. A tal fine, tramite opportune manipolazioni algebriche,
si riscrive la (4.35) come segue:
(Ki − τ −1 Kd0 ) s + Ki τ −1
s 2 + τ −1 s
ove è stato posto per comodità
G (s) = K p + Kd0 +
(4.36)
Kd0 = Kd τ −1
A questo punto, indicando rispettivamente con ξ,v, w, lo stato, l’ingresso e l’uscita del controllore P.I.D., si può determinare una realizzazione minima della
funzione di trasferimento (4.36), ad esempio in forma canonica raggiungibile,
ottenendo:
0
1
0
ξ+
v
ξ̇ =
−
1
0 −τ
1
w = Ki τ −1 Ki − τ −1 Kd0 ξ + (K p + Kd0 ) v
(4.37)
ξ ∈ R2 , v ∈ R, w ∈ R
A questo punto di noti che la matrice dinamica del controllore (4.37) è banalmente diagonalizzabile, per cui per porre il controllore (4.37) nella forma
desiderata è sufficiente esprimere le equazioni (4.37) nella base degli autovettori della matrice in questione. Tale trasformazione, come è noto, può essere
effettuata ponendo:
ϕ = Tξ
(4.38)
dove la matrice T è quadrata di dimensione 2, non singolare e tale che la sua
inversa T −1 possieda come colonne una base di autovettori della matrice dinamica che si vuole rendere diagonale. In particolare nel caso in questione si ha:
T
−1
1
1
=
− τ −1 0
(4.39)
Perciò effettuando la trasformazione suddetta si ottiene:
− τ −1 0
−τ
ϕ+
v
0
0
τ
h 0
i
w = K2d Kτi ϕ + (K p + Kd0 ) v
τ
ϕ̇ =
(4.40)
ϕ ∈ R2 , v ∈ R, w ∈ R
Infine esplicitando le equazioni differenziali che descrivono l’evoluzione delle
due componenti del vettore ϕ si ottiene:
ϕ̇1 = −τ −1 ϕ1 − τ v
ϕ̇2 = τ v
K0
Ki
w = 2d ϕ1 +
ϕ2 + (K p + Kd0 ) v
τ
τ
(4.41)
72
4.4 valutazione delle prestazioni del sistema di controllo in piccolo
Confrontando le espressioni appena ottenute con le (4.34) si hanno le seguenti
corrispondenze tra le due diverse parametrizzazioni del controllore P.I.D. (4.35):
k11 = Ki τ −1
K12 = Kd τ −1
λ1 = τ
(4.42)
λ2 = − τ
K2 = K p + K d τ −1
E’ evidente che la trattazione appena svolta può essere ripetuta utilizzando
un controllore PID affiancato ad un altro compensatore lineare, giungendo alle medesime conclusioni, pertanto un generico compensatore comprensivo di
un azione integrale soddisfa in pieno i requisiti necessari per la soluzione del
problema di regolazione locale.
Osservazione 4.2:
Quanto mostrato mette bene in evidenza quali siano i ruoli delle tre azioni di controllo
del regolatore P.I.D. [1]. Infatti l’azione integrale rappresenta il modello interno, mentre l’azione proporzionale e l’azione derivativa rappresentano la parte stabilizzante del
controllore. Si noti che il complesso dell’azione derivativa e dell’azione proporzionale di
fatto può essere interpretata come l’azione di controllo ottenuta tramite l’utilizzo di un
osservatore di ordine ridotto. In particolare se il modello del processo allargato fosse del
secondo ordine allora l’azione proporzionale congiunta alla azione derivativa rappresenterebbe esattamente l’azione di un osservatore asintotico di ordine ridotto dello stato di
tale processo. L’aggiunta di ulteriori azioni dinamiche nel controllore fa si che l’azione
di stima dell’osservatore sia ricalcata il più possibile.
4.4
valutazione delle prestazioni del sistema di
controllo in piccolo
La trattazione affrontata fin ora, è stata condotta prevalentemente utilizzando
la procedura di linearizzazione come metodo per garantire la stabilizzazione
esponenziale locale dei punti di equilibrio di interesse. Tuttavia la procedura di
linearizzazione può essere utilizzata anche come strumento per la valutazione
delle prestazioni del sistema di controllo, in termini di risposta forzata, almeno
ingressi di ampiezza modesta, ossia attraverso l’utilizzo della linearizzazione
è possibile condurre una analisi sul comportamento del sistema di controllo
per piccoli segnali. Per far ciò si indichi con r l’ingresso di riferimento e si
consideri il sistema di controllo ΣH in (4.28), il quale può essere espresso nella
forma seguente dopo aver sostituito nella dinamica del processo l’espressione
del controllo u(t).
ẋ = f ( x, ϕ , d)
ε̇ = H ε + η r
ϕ̇ = Ψ ϕ + Λ ( L ε − h( x ))
(4.43)
73
4.4 valutazione delle prestazioni del sistema di controllo in piccolo
ebbene, indicando con ( x̂, ϕ̂, ε̂ ) il punto di equilibrio per il sistema di controllo
in corrispondenza ai valori costanti y∗ , d∗ per gli ingressi r, d e ponendo:
x̃ = x − x̂
(4.44)
ε̃ = ε − ε̂
(4.46)
ϕ̃ = ϕ − ϕ̂
(4.45)
∗
r̃ = r − y
d˜ = d − d∗
(4.47)
(4.48)
ỹ = y − y∗
(4.49)
allora per piccoli scostamenti dalla condizione di equilibrio la dinamica del
sistema di controllo ΣH è bene approssimata dal seguente sistema LTI2 :
x̃˙ = ( A − B K2 C ) x̃ + B K1 ϕ̃ + B K2 Lε̃ + Bd d˜
ε̃˙ = H ε̃ + η r̃
ϕ̃˙ = Ψ ϕ̃ − Λ C x̃ + Λ L ε̃
(4.50)
ỹ = −C x̃
dove le matrici A, B, C sono definite come di consueto mentre la matrice Bd è
definita come segue: pertanto ponendo:
Bd :=
∂ f ∂d x̂, ϕ̂, dˆ
(4.51)
al fine di ottenere una rappresentazione più compatta è possibile descrivere il
sistema di controllo ΣH linearizzato nel dominio di Laplace, ossia attraverso
una rappresentazione con matrici di trasferimento. In particolare indicando
con Ỹ (s), D̃ (s), R̃(s) rispettivamente le trasformate di Laplace delle funzioni
ỹ(t), d˜(t), r̃ (t)
e ponendo:
P ( s ) = C ( s I − A ) −1 B
G ( s ) = K1 ( s I − Ψ )
−1
Λ + K2
Pyd (s) = C (s I − A)−1 Bd
H(s) = (s I − H )−1 η
(4.52)
(4.53)
(4.54)
(4.55)
(4.56)
si ottiene il seguente legame tra le trasformate di Laplace dell’uscita controllata e le trasformate di Laplace dei segnali esogeni ỹ(t) e r̃ (t):
Ỹ (s) =
P(s) G (s)
1
H(s) R̃(s) +
P (s) D̃ (s)
1 + P(s) G (s)
1 + P(s) G (s) yd
tali relazioni possono essere espresse in maniera più compatta introducendo
due particolari funzioni di trasferimento, ossia la funzione di sensitività S(s)
2 Per giungere alla seguente rappresentazione si faccia riferimento al sistema (4.28) nel quale
compare in maniera esplicita l’espressione della legge di controllo
74
4.4 valutazione delle prestazioni del sistema di controllo in piccolo
e la funzione di sensitività complementare T (s), definite rispettivamente come
segue:
S(s) =
1
1 + P(s) G (s)
T (s) =
P(s) G (s)
1 + P(s) G (s)
(4.57)
ottenendo in definitiva il seguente legame tra le suddette grandezze:
Ỹ (s) = T (s) H(s) R̃(s) + S(s) Pyd (s) D̃ (s)
(4.58)
in tal modo il comportamento del sistema di controllo, per piccoli segnali, può
essere dedotto facendo riferimento al classico schema di controllo. Ottenuta
questa descrizione del sistema di controllo, è possibile fare riscorso ad alcuni
indicatori utilizzati nella teoria classica del controllo, per stabilire quali sia il
comportamento in piccolo del sistema di controllo. In particolare, al fine di
valutare la rapidità di risposta del sistema di controllo ad una variazione a
gradino della temperatura di riferimento, si può fare riferimento alla banda
passante delle funzione di risposta armonica H( j ω ) T ( j ω ), mentre si può fare
ricorso allo studio del modulo della funzione di risposta armonica della funzione dei sensitività, ossia |S( j ω )| per valutare quale sia l’effetto delle variazioni
della funzione di trasferimento P(s) sul comportamento del sistema di controllo. Questo ultimo fatto risulta molto importante nel caso in esame, in quanto a
fronte di una variazione del punto di lavoro del processo, certamente corrisponderà una variazione della funzione di trasferimento P(s). Inoltre osservando la
(4.58), si può evincere che lo studio della funzione di sensitività fornisce un
indicazione, almeno in prima approssimazione, dell’effetto delle variazioni del
disturbo d˜(t) sull’errore di regolazione ỹ. Si tenga bene presente infine che
le indicazioni ottenute tramite lo studio dell’approssimazione lineare del sistema di controllo, rappresentano solamente una stima del comportamento reale
del sistema di controllo, tuttavia sembra ragionevole cercare di ottenere delle
buone prestazioni in piccolo al fine di ottenerle per il reale sistema di controllo.
75
5
N
S I N T E S I D E I R E G O L ATO R I P E R LO
S C A M B I ATO R E D I C A LO R E
el capitolo precedente sono stati forniti alcuni strumenti di controllo necessari per la risoluzione del problema di regolazione
asintotica locale di un sistema non lineare, principalmente tali strumenti sono
incentrati attorno al Teorema A.1. Tuttavia questi strumenti possono essere
utilizzati per risolvere il Problema 4.1, a patto di rilassare alcune specifiche e
considerare i vincoli di ammissibilità imposti in tale problema. Infatti il Problema 4.1 differisce per alcuni aspetti dal Problema 4.3. Più precisamente le
differenze tra i due problemi sono principalmente tre:
1. Il Problema 4.3 è di tipo locale, mentre il Problema 4.1 di per se è di tipo
globale.
2. Per risolvere il Problema 4.3 è stata imposta l’Assunzione 4.1 che non è
detto sia verificata nel Problema 4.1
3. Nel Problema 4.3 non ci sono vincoli sulla variabile di controllo mentre
nel Problema 4.1 tali vincoli sono invece presenti
A fronte delle constatazioni appena riportate è necessario fare delle
opportune osservazioni. In particolare riguardo alla prima delle
tre constatazioni, si noti che essa impone che la soluzione determinata a partire
dagli strumenti di controllo presentati, sarà necessariamente di tipo locale. Ciononostante però questo fatto non è un grosso problema in quanto nella ( 3.1.1 a
pagina 29) è stato mostrato che il comportamento dell’impianto non è fortemente non lineare, ovvero ciò suggerisce intuitivamente che: anche se la soluzione
ottenuta per il problema di regolazione della temperatura è di tipo locale, la
regione di attrazione del generico equilibrio può essere ritenuta sufficientemente ampia da ritenere questo tipo di soluzione ingegneristicamente accettabile.
Anche il secondo dei punti indicati di fatto non rappresenta un problema in
quanto nel caso in esame l’Assunzione 4.1 è banalmente soddisfatta dato che
∗ e una coppia ( ϑ , ϑ ), purchè tale
fissata una temperatura di riferimento ϑlu
v
li
scelta dia luogo ad un valore della velocità del liquido v che sia ammissibile,
è sempre possibile trovare in corrispondenza una e una sola coppia ( x, v) che
sia un equilibrio dell’impianto e che fornisca una temperatura del liquido in
m uguale a quella di riferimento ϑ ∗ . Viceversa la terza delle consideuscita ϑlu
lu
razioni suddette, di per se rappresenta un problema non banale in quanto la
generica legge di controllo, sintetizzata tramite gli strumenti presentati nel capitolo precedente, garantisce che il valore del controllo sia ammissibile solamente
∗ , ϑ , ϑ in modo tale che in
in senso asintotico, ossia garantisce che fissati ϑlu
v
li
corrispondenza vi sia un valore di velocità v̂ ammissibile allora:
lim v(t) = v̂
t→ +∞
76
5.1 verifica delle condizioni per la sintesi del controllore
viceversa non è noto quello che accade nel generico istante di tempo finito.
Questo fatto rappresenta un grosso problema in quanto l’effetto della saturazione introdotta dalla servovalvola, negli istanti in cui il controllo assume dei
valori non ammissibili, non è stato in nessun modo preso in considerazione,
ciò implica che possano instaurarsi dei fenomeni oscillatori a priori non previsti! Tuttavia se il punto di lavoro che si desidera raggiungere è ammissibile, è
evidente che il fenomeno di saturazione può solamente presentarsi durante il
transitorio del sistema di controllo, perciò se per avventura, non esiste nessun
istante della risposta transitoria in cui il controllo assume dei valori non ammissibili, allora l’effetto della saturazione può essere del tutto ignorato garantendo
che la strategia di controllo, purché lo stato iniziale del sistema di controllo sia
sufficientemente prossimo al punto di lavoro desiderato, permetta di ottenere
la regolazione asintotica. E’ evidente che analiticamente, con gli strumenti presentati, non si può garantire che il controllo sia ammissibile per ogni istante
di tempo, d’altra parte è anche evidente che se lo stato iniziale del sistema di
controllo è sufficientemente prossimo al punto di lavoro desiderato e se la legge di controllo è progettata in modo da non generare grosse sovraelongazioni
della variabile di controllo, allora ingegneristicamente si può ritenere alquanto
improbabile l’attivazione della saturazione. Quanto detto riguardo alla mancata conoscenza a priori del comportamento del sistema di controllo, mette in
evidenza un fatto molto importante e di “stampo” prettamente ingengeristico,
ossia la legge di controllo deve necessariamente essere validata attraverso numerose simulazioni, in modo da poter stabilire quale sia il comportamento del
sistema di controllo nel maggior numero di situazioni. Avendo messo bene
in chiaro quali siano le limitazioni derivanti dall’utilizzo delle metodologie di
controllo presentate nel capitolo precedente, si può passare a l’utilizzo di tali
metodologie per la soluzione del Problema 4.1.
5.1
verifica delle condizioni per la sintesi del
controllore
Si osservi prima di ogni altra cosa che la trattazione svolta fino ad ora è applicabile al processo in esame a patto che siano verificate sia l’assunzione 4.1 sia le
condizioni dettate nel Teorema A.1. A tal scopo si noti che per il modello ( 2.5
a pagina 20), è stato mostrato che fissata una terna (ϑlu , ϑli , ϑv ) ammissibile esiste sempre un valore del velocità del fluido v tale da mantenere la temperatura
del fluido in uscita ad un valore costante e pari a ϑlu , ed inoltre è stato mostrato
che in corrispondenza a tale valore di velocità v esiste sempre un solo equilibrio
per il sistema a dimensione finita ( 2.20 a pagina 25), pertanto assumendo che
il modello ad elementi finiti approssimi in maniera sufficientemente accurata il
modello ( 2.5 a pagina 20) si può ritenere ragionevolmente che il suddetto valore della velocità del fluido sia unico, ossia l’assunzione 4.1 può essere ritenuta
ingegneristicamente verificata nel caso in esame. Per quanto riguarda le ipotesi
dettate nel Teorema A.1 si tenga presente che esse dovrebbero essere verificate
per infiniti punti di lavoro il che è chiaramente impossibile, tuttavia si osservi
che le condizioni imposte dal teorema in questione, sono delle condizioni per
la stabilizzazione, pertanto se esse sono verificate almeno dal punto di lavoro
77
5.2 controllo pid dello scambiatore di calore
nominale è possibile sintetizzare una legge di controllo basandosi sul modello
nominale, cercando di garantire il più possibile il mantenimento della proprietà
della stabilità anche a fronte di variazioni del punto di lavoro, a posteriori delle
sintesi del controllore sarà sufficiente effettuare delle analisi volte a valutare il
margine di robustezza del sistema di controllo. Per quanto detto si consideri il
ben noto modello dinamico dello scambiatore di calore:

































ẋ (t) = A(v(t)) x (t) + B(v(t)) w(t)
x ( t0 ) = x0
y(t) = C x (t)
x (t) ∈ R71
v(t) ∈ R
(5.1)
w(t) ∈ R2
A(v) ∈ R71×71
B(v) ∈ R71×2
con ciò linearizzando tale modello attorno al punto di lavoro nominale, ( x̄, v̄ , w̄),
si ottiene il seguente modello LTI:
x̃˙ = A x̃ + B ṽ
(5.2)
ỹ = C x̃
(5.3)
ove sono state effettuare per comodità le seguenti posizioni:
A := A(v̄)
∂ A B :=
∂v x̄, v̄, w̄
(5.4)
Pertanto le condizioni per verificare se le ipotesi del Teorema A.1 siano soddisfatte si deve verificare che:
1. La tripla ( A, B, C ) è stabilizzabile e rilevabile
A B
2. rank
= n+1
C D
Ebbene verificando tali condizioni attraverso un programma di calcolo si evince
che esse sono entrambi valide e con ciò è possibile sintetizzare un controllore
K LTI del tipo indicato in (A.23) in grado di risolvere il problema di regolazione
locale. Prima di procedere, per quanto detto in precedenza in merito all’analisi
del processo per piccoli segnali, è conveniente anche linearizzare l’impianto
attorno al punto di lavoro nominale considerando le relative variazioni delle
temperature del liquido in ingresso e della temperatura del vapore rispetto ai
loro valori nominali.
5.2
controllo pid dello scambiatore di calore
Il primo dei controllori sintetizzati per la soluzione del problema di regolazione locale è un controllore PID, come già mostrato in precedenza, tale controllore possiede tutti i requisiti per essere un compensatore candidato per la
78
5.2 controllo pid dello scambiatore di calore
soluzione del problema in questione. In particolare la sintesi di tale compensatore consiste nel determinare i quattro parametri della seguente funzione di
trasferimento:
Ki
Kd
C (s) = K p +
+
s
1+τs
in modo tale da garantire che il modello nominale (5.2) sia asintoticamente
stabile e in modo che si ottengano dei buoni margini di robustezza. In particolare la sintesi di tale controllore, come già detto in precedenza può essere fatta
utilizzando le tecniche di sintesi in frequenza, basate sulla funzione di trasferimento del processo nominale (5.2), che verrà indicata con P(s). Dunque al fine
di ottenere almeno una stima di primo tentativo dei parametri da determinare
è possibile utilizzare per il modello nominale uno dei modelli approssimati determinati nel Capitolo 3, ad esempio il modello ( 3.3 a pagina 34). Effettuando
vari tentativi e utilizzando gli strumenti messi a disposizione dai programmi
di calcolo si perviene alle seguente funzione di trasferimento per il controllore
PID:
C (s) =
−0.4339 s2 − 0.267 s − 0.03998
s2 + 15.15 s
(5.5)
Alla quale corrispondo i seguenti valori per i parametri sopra citati:
K p = −0.017449
(5.6)
Ki = −0.0026388
(5.7)
τ = 0.066007
(5.9)
Kd = −0.02749
(5.8)
L’adozione di tale controllore nel controllo dell’impianto (5.2) consente di ottenere i seguenti valori per il margine di fase m ϕ e per il margine di guadagno
mg :
m ϕ = 99.2◦
(5.10)
m g = +∞
(5.11)
Per comprende quale sia l’azione espletata dal controllore PID sintetizzato sul
processo è conveniente osservare il diagramma di Bode della risposta armonica
del medesimo, mostrata in (46).
Per valutare l’effetto del controllore PID sulla dinamica del sistema linearizzato, in (47) viene mostrato il diagramma di Bode della funzione di trasferimento
di anello L(s) = C (s) P(s), mentre in (48) è mostrato il diagramma di Nyquist
della funzione L(s). Per concludere la sintesi del sistema di controllo, è necessario sintetizzare il sistema H al fine di ridurre le sollecitazioni trasmesse
all’attuatore a fronte di variazioni a gradino della temperatura di riferimento
y∗ . Inoltre, come detto in precedenza, è possibile attraverso una opportuna
scelta del sistema H ottenere anche dei benefici in termini di risposta transitoria dell’errore di regolazione, pertanto la sintesi di tale sistema verrà effettuata
anche prendendo in considerazione tale obiettivo. Attraverso tecniche di tipo
trial-&-error basate su ripetute simulazioni numeriche, si perviene alla seguente
funzione di trasferimento per il sistema H.
H(s) =
s2
0.25
+ 0.23 s + 0.25
(5.12)
79
5.2 controllo pid dello scambiatore di calore
Bode Diagram
Magnitude (dB)
0
−10
−20
−30
Phase (deg)
−40
270
225
180
135
90
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
3
10
10
Figura 46.: Diagramma di Bode della risposta armonica del controllore PID
Bode Diagram
100
L (s)
Magnitude (dB)
50
P (s)
0
−50
−100
−150
Phase (deg)
−200
180
0
−180
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 47.: Diagramma di Bode della funzione d’anello L(s)
Al fine di valutare quale sia l’effetto dell’azione in avanti nel sistema di controllo si osservi in (Figura 49) che l’introduzione dell’azione in avanti produce un
aumento della banda passante del legame tra r̃ e ỹ, in particolare con l’introduzione dell’azione in avanti, la banda passante passa da circa 0.17 rad s−1 a circa
0.58 rad s−1 . Per quanto detto prima al fine di valutare, in prima approssimazione l’effetto delle variazioni della temperatura del fluido in ingresso e della
temperatura del vapore sull’errore di regolazione e l’effetto delle variazioni del
80
5.2 controllo pid dello scambiatore di calore
Nyquist Diagram
1.5
1
Imaginary Axis
0.5
0
−0.5
−1
−1.5
−1.5
−1
−0.5
0
Real Axis
0.5
1
1.5
Figura 48.: Diagramma di Nyquist della funzione d’anello L(s)
punto di lavoro sul comportamento del sistema di controllo, si può prendere
in considerazione il modulo della funzione di risposta armonica associato alla
funzione di sensitività S(s), il quale è mostrato in (50).
Per valutare la robustezza del sistema di controllo sintetizzato si riportano
in (Figura 51) i diagrammi di Bode delle funzioni di anello Li (s) ottenute considerando i modelli linearizzati attorno ai punti di equilibrio indicati nel corso
della trattazione. Si riportano inoltre in (Tabella 2) i margini di guadagno, i
margini di fase e le pulsazioni di attraversamento, ottenute per ogni punti di
lavoro, considerato.
81
5.2 controllo pid dello scambiatore di calore
Bode Diagram
T (s)
Magnitude (dB)
50
T (s) H (s)
0
−50
−100
−150
−200
Phase (deg)
0
−90
−180
−270
−360
−2
10
−1
10
0
10
Frequency (rad/sec)
1
2
10
10
Figura 49.: Diagrammi di Bode della funzione di sensitività complementare T (s) e
della funzione H(s) T (s)
Bode Diagram
10
0
Magnitude (dB)
−10
−20
−30
−40
−50
−3
10
−2
10
−1
0
10
10
Frequency (rad/sec)
1
10
2
10
Figura 50.: Diagramma di Bode del modulo della funzione di sensitività S(s)
82
5.2 controllo pid dello scambiatore di calore
Bode Diagram
60
L 1 (s)
Magnitude (dB)
40
20
L 2 (s)
0
L 3 (s)
−20
L 4 (s)
−40
L 0 (s)
−60
Phase (deg)
−80
−45
−90
−135
−180
−2
10
−1
0
10
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 51.: Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di
lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s)
con i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera
consistente alla numerazione riportata in (Tabella 1)
Tabella 2.: Valore dei margini di guadagno e di fase ottenuti per i punti di lavoro
scelti con il controllo PID. Il punto di lavoro nominale viene indicato con il
numero 0.
Punto di lavoro
#
mg
db
0
1
2
3
4
∞
∞
∞
∞
∞
mϕ
(deg)
99
47
106
96
81
ωt
(rad s−1 )
0,221
0,208
0,078
0,414
0,143
83
5.3
5.3 controllo pid con rete compensatrice
controllo pid con rete compensatrice
I risultati ottenuti con il controllore PID (5.5) rappresentano un buon compromesso tra rapidità di risposta e robustezza alle variazioni parametriche indotte
dalle variazioni del punto di lavoro. Tuttavia si nota in (Figura 51) che attraverso l’introduzione di una rete compensatrice si può aumentare lievemente la
pulsazione di attraversamento della funzione L(s), ottenendo in tal modo un
miglioramento delle prestazioni del sistema di controllo, anche in situazioni diverse dalla nominale. In particolare il miglioramento desiderato è in termini di
attenuazione degli effetti indotti, sull’errore di regolazione, da parte dalle variazioni della temperatura del liquido in ingresso e della temperatura del vapore,
tali benefici ovviamente devono essere ottenuti senza ridurre drasticamente i
margini di stabilità del sistema di controllo ottenuti attraverso il controllore
PID (5.5). In definitiva la funzione di compensazione sintetizzata è la seguente:
G (s) =
0.7979 s + 0.1535
s + 0.1094
(5.13)
in altre parole la funzione di compensazione è una funzione ritardatrice e viene progettata con lo scopo di ottenere un aumento del modulo della funzione
d’anello L(s) solamente a bassa frequenza, in modo da lasciare praticamente inalterate le prestazioni nel caso nominale, in quanto esse già abbastanza
soddisfacenti. Per capire quale sia l’effetto sulla funzione di anello L(s) del
compensatore G (s) in (Figura 52) è mostrato in diagramma di Bode del medesimo. Mentre in (Figura 53) è mostrato il compensatore risultante dall’azione
combinata dal controllore PID (5.5) e della rete ritardatrice (5.13).
Per quanto riguarda l’azione di feed-forward H(s) si sceglie di utilizzare la
medesima azione utilizzata nel controllo PID, ossia:
H(s) =
s2
0.25
+ 0.23 s + 0.25
(5.14)
Conclusa la sintesi del sistema di controllo, al fine di valutare i risultati ottenuti
in (Figura 54) e in (Figura 55) si mostrano, rispettivamente i diagrammi di Bode
e di Nyquist della funzione d’anello L(s) = G (s) C (s) P(s) nel caso nominale.
Mentre al fine di valutare, in prima approssimazione l’effetto delle variazioni
della temperatura del fluido in ingresso e della temperatura del vapore sull’errore di regolazione e l’effetto delle variazioni del punto di lavoro sul comportamento del sistema di controllo, si riporta in (Figura 56) il diagramma di Bode
del modulo della funzione di sensitività nel caso nominale.
Per quanto riguarda l’analisi della risposta nell’errore di regolazione del sistema di controllo a fronte di piccole variazioni della temperatura di riferimento
si mostra in (Figura 57) il diagramma di Bode della funzione di trasferimento
T (s) H(s) nel caso nominale. In particolare osservando la (Figura 57) si può
constatare che la banda passante del sistema di controllo nominale, per quanto
riguarda il legame tra la variazione delle temperatura di riferimento e l’errore
di regolazione è pari a:
ωb = 0.5550 rad s−1
Si riportano inoltre in (Tabella 3) i margini di guadagno, i margini di fase e le
pulsazioni di attraversamento, ottenute per ogni punti di lavoro, considerato.
84
5.3 controllo pid con rete compensatrice
Bode Diagram
Magnitude (dB)
3
2
1
0
−1
−2
−3
5
Phase (deg)
0
−5
−10
−15
−20
−3
10
−2
10
−1
10
Frequency (rad/sec)
0
1
10
10
Figura 52.: Diagrammi di Bode della funzione di compensazione G (s)
Magnitude (dB)
Bode Diagram
20
G(s) C (s)
10
C (s)
0
−10
−20
−30
Phase (deg)
−40
270
225
180
135
90
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 53.: Diagrammi di Bode del controllore C (s) G (s)
Infine per valutare la robustezza del sistema di controllo sintetizzato si ripor-
85
5.3 controllo pid con rete compensatrice
Bode Diagram
50
Magnitude (dB)
L (s)
0
P (s)
−50
−100
Phase (deg)
−150
180
90
0
−90
−180
−2
10
−1
10
0
1
2
10
10
Frequency (rad/sec)
3
10
10
Figura 54.: Diagrammi di Bode della funzione di trasferimento P(s) a confronto con la
funzione d’anello L(s) = P(s) C (s) G (s)
Nyquist Diagram
1.5
Imaginary Axis
1
0.5
0
−0.5
−1
−1.5
−1.5
−1
−0.5
0
Real Axis
0.5
1
1.5
Figura 55.: Diagramma di Nyquist della funzione di anello L(s) nel caso nominale
86
5.3 controllo pid con rete compensatrice
Bode Diagram
5
0
−5
−10
Magnitude (dB)
−15
−20
−25
−30
−35
−40
−45
−50
−3
10
−2
−1
10
0
1
10
10
Frequency (rad/sec)
10
Figura 56.: Diagramma di Bode della funzione sensitività nel caso nominale
Bode Diagram
Magnitude (dB)
50
0
−50
−100
−150
−200
Phase (deg)
0
−90
−180
−270
−360
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 57.: Diagramma di Bode della funzione T (s) H(s) nel caso nominale
tano in (Figura 58) i diagrammi di Bode delle funzioni di anello Li (s) ottenute
87
5.3 controllo pid con rete compensatrice
Tabella 3.: Valore dei margini di guadagno e di fase ottenuti per i punti di lavoro scelti
con il controllo PID+ Rete compensatrice. Il punto di lavoro nominale viene
indicato con il numero 0.
Punto di lavoro
#
mg
db
mϕ
(deg)
0
1
2
3
4
∞
∞
∞
∞
∞
86,4
33,4
94,2
91,7
65,1
ωt
(rad s−1 )
0,208
0,204
0,095
0,326
0,15
considerando i modelli linearizzati attorno ai punti di equilibrio indicati nel
corso della trattazione.
Bode Diagram
100
Magnitude (dB)
L 1 (s)
L 2 (s)
50
L 3 (s)
0
L 4 (s)
L 0 (s)
−50
Phase (deg)
−45
−90
−135
−180
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 58.: Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di
lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s)
con i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera
consistente alla numerazione riportata in (Tabella 1)
88
5.4
controllo lqg
5.4 controllo lqg
L’ultimo dei controllori che viene proposto è un controllore di tipo observer
based, in particolare un controllore LQG, progettato sul modello del processo
linearizzato attorno al punto di lavoro nominale. Tale metodologia di controllo,
come illustrato precedentemente, consiste nel progettare un controllore, basato
sulla metodologia ad assegnazione degli autovalori con osservatore asintotico
dello stato per il processo allargato definito come segue:
x̃˙ = A x̃ + B k11 ϕ̃1 + B v
ϕ̃˙ 1 = −λ1 C x̃
(5.15)
ỹ = −C x̃
Tale sistema avendo verificato in precedenza che le ipotesi del Teorema A.1 sono
soddisfatte, è certamente stabilizzabile e rilevabile, perciò certamente esiste un
controllore stabilizzante di tipo LQG. Ebbene, si osservi che come bene noto la
sintesi di un controllore LQG avviene utilizzando degli algoritmi numerici che a
partire da alcuni parametri di progetto e dal modello del processo sotto controllo, forniscono il compensatore cercato, pertanto dato che il processo allargato è
di dimensione considerevole è necessario effettuare, se possibile, una riduzione
dell’ordine di tale modello e dato che la dimensione elevata del processo allargato deriva naturalmente dalla dimensione elevata del processo soggetto al
controllo, è evidente che la riduzione dell’ordine andrà effettuata relativamente
al solo processo da controllare linearizzato attorno al punto di lavoro nominale.
A tal scopo è possibile effettuare una analisi ai valori singolari del processo
in esame, il quale esito è mostrato in (Figura 59). Osservando la (Figura 59)
si nota che l’ordine del modello, data la presenza di numerosi modi naturali
aventi scarsa influenza sul comportamento ingresso uscita del processo, può
essere ridotto notevolmente attraverso una realizzazione in forma bilanciata. In
particolare una scelta opportuna prevede di ridurre il modello del processo da
dimensione 71 a dimensione 6, a seguito di tale scelta si può osservare in (Figura 60) come la risposta in frequenza del modello originario sia sovrapponibile
a quella modello di ordine ridotto, la cui funzione di trasferimento viene indicata con Pr (s). Ottenuto il modello di ordine ridotto, è possibile sintetizzare
attraverso una procedura trial & error e attraverso l’uso di opportuni algoritmi numerici il controllore LQG sopra citato, pertanto includendo in fine nel
compensatore LQG sintetizzato la parte integrale del controllore come indicato
in (Figura 44), ottenendo un controllore descritto dalla seguente funzione di
trasferimento:
C (s) =
10.4126(s + 1)(s + 0.4477)(s2 + 0.6461s + 0.3566)(s2 + 0.4299s + 1.169)
s(s2 + 0.1777s + 0.6906)(s2 + 1.071s + 1.553)(s2 + 36.29s + 657.4)
(5.16)
Per comprendere l’azione del controllore progettato, si riporta in (Figura 61) il
diagramma di Bode della funzione di trasferimento C (s) di cui sopra. Osservando il diagramma di Bode in (Figura 61) si nota che il controllore progettato è
fortemente legato al modello nominale del processo, in quanto tale controllore
tende a imporre una dinamica al sistema di controllo compensando alcuni fenomeni presenti nella dinamica del processo nominale, pertanto risulta necessario
89
5.4 controllo lqg
Hankel Singular Values
45
40
35
30
abs
25
20
15
10
5
0
0
10
20
30
40
Order
50
60
70
80
Figura 59.: Valori singolari della funzione di trasferimento del processo linearizzato
attorno al punto di lavoro nominale
operare ad una riduzione dell’ordine di tale controllore in modo da evitare di
ottenere un controllore che sia poco robusto. A tal scopo si può operare una
riduzione dell’ordine basata sulle realizzazioni in forma bilanciata, come quanto è stato fatto in precedenza per il processo, in particolare in (Figura 62) si
mostrano i valori singolari del controllore.
osservando la (Figura 62) si nota che si può ad esempio ottenere un modello
del controllore ridotto del terzo ordine, in tal modo si ottiene il controllore
descritto dalla seguente funzione di trasferimento.
Cr (s) =
10.42s2 + 12.47s + 1.802
s3 + 36.28s2 + 653.5s
(5.17)
Per comprendere l’effetto della riduzione dell’ordine in (Figura 63) vengono
mostrati i diagrammi di Bode della funzione di trasferimento C (s) confrontata
con la funzione di traferimento Cr (s).
La sintesi del controllore LQG si può ritenere ultimata, tuttavia si nota attraverso numerose simulazioni che il controllore progettato, data l’assenza di
un legame diretto ingresso uscita risulta poco performante rispetto ai controllori mostrati in precedenza. Per risolvere tale problema si può aggiungere un
legame diretto ingresso-uscita al controllore, tuttavia si noti che tale modifica
potrebbe comportare una perdita di stabilità del sistema di controllo in quanto
come è stato mostrato la matrice dinamica del sistema di controllo linearizzato attorno ad una generico punto di lavoro, dipende anche dal legame diretto
ingresso-uscita K2 del controllore. Tuttavia quanto detto non presenta un problema a livello progettuale, in quanto attraverso delle simulazioni è possibile
90
5.4 controllo lqg
Bode Diagram
40
P r (s)
Magnitude (dB)
20
P (s)
0
−20
−40
−60
Phase (deg)
−80
180
135
90
45
0
−2
10
−1
10
0
10
Frequency (rad/sec)
1
10
2
10
Figura 60.: Diagrammi di Bode della funzione di trasferimento del processo nominale
e della sua approssimazione di ordine 6
determinare un valore del suddetto legame diretto ingresso-uscita che non pregiudichi la stabilità del sistema di controllo ma che ne migliori le prestazioni.
La procedura appena descritta conduce alla seguente scelta:
K2 = 0.01
in tal modo la funzione di trasferimento del controllore LQG con legame diretto
è la seguente:
Cd (s) =
0.01s3 + 10.78s2 + 19.01s + 1.802
s3 + 36.28s2 + 653.5s
(5.18)
Per ultimare la sintesi del sistema di controllo non resta che progettare l’azione
in avanti H(s), anche in questo caso tale azione viene progettata attraverso una
procedura trial & error, giungendo alla seguente funzione di trasferimento:
H(s) =
s2
0.16
+ 0.16s + 0.16
Completata la sintesi del sistema di controllo si mostrano, come di consueto i
diagrammi di Bode del compensatore Cd (s) (Figura 64), della funzione d’anello
L(s) (Figura 65),della funzione di sensitività S(s) (Figura 66), della funzione
di trasferimento tra la variazione di riferimento e la variazione dell’errore di
regolazione ossia H(s) T (s) (Figura 67) mentre in (Figura 68) si mostra in diagramma di Nyquist della funzione d’anello L(s) nel caso nominale.
In particolare osservando la (Figura 67) si può constatare che la banda passante del
91
5.4 controllo lqg
Bode Diagram
Magnitude (dB)
0
−10
−20
−30
Phase (deg)
−40
270
225
180
135
90
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 61.: Diagrammi di Bode della funzione di trasferimento del compensatore LQG
con azione integrale
sistema di controllo nominale, per quanto riguarda il legame tra la variazione
delle temperatura di riferimento e l’errore di regolazione, è pari a:
ωb = 0.5674 rad s−1
Ai fini di valutare la robustezza e in una qualche misura le perfomance del
sistema di controllo, si riportano inoltre in (Tabella 4) i margini di guadagno, i
margini di fase e le pulsazioni di attraversamento, ottenute per ogni punti di
lavoro, considerato.
Tabella 4.: Valore dei margini di guadagno e di fase ottenuti per i punti di lavoro
scelti con il controllo LQG. Il punto di lavoro nominale viene indicato con
il numero 0.
Punto di lavoro
#
mg
db
0
1
2
3
4
48,3
50,8
49,2
44,3
54,1
mϕ
(deg)
65,2
33
123
52,3
69,6
ωt
(rad s−1 )
0,416
0,26
0,116
0,582
0,234
Infine per valutare maggiormente la robustezza del sistema di controllo sintetizzato si riportano in (Figura 69) i diagrammi di Bode delle funzioni di anello
Li (s) ottenute considerando i modelli linearizzati attorno ai punti di equilibrio
indicati nel corso della trattazione.
92
5.4 controllo lqg
Hankel Singular Values (1 of jw−axis pole(s) = inf, not shown)
0.16
0.14
0.12
abs
0.1
0.08
0.06
0.04
0.02
0
1
2
3
4
Order
5
6
7
Figura 62.: Valori singolari associati ai modi naturali presenti nel legame ingressouscita del controllore LQG
Bode Diagram
20
C r (s)
Magnitude (dB)
10
C (s)
0
−10
−20
−30
Phase (deg)
−40
270
225
180
135
90
−3
10
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
Figura 63.: Diagramma di Bode di Cr (s) e di C (s) a confronto
3
10
93
5.4 controllo lqg
Bode Diagram
20
Magnitude (dB)
10
0
−10
−20
−30
−40
Phase (deg)
−50
270
225
180
135
90
−2
10
0
2
10
10
Frequency (rad/sec)
4
10
Figura 64.: Diagramma di Bode della funzione di trasferimento Cd (s)
Bode Diagram
Magnitude (dB)
100
0
−100
−200
−300
−45
Phase (deg)
−90
−135
−180
−225
−270
−2
10
0
2
10
10
Frequency (rad/sec)
4
10
Figura 65.: Diagramma di Bode della funzione di trasferimento L(s)
94
5.4 controllo lqg
Bode Diagram
10
0
Magnitude (dB)
−10
−20
−30
−40
−50
−3
10
−2
10
−1
0
1
10
10
Frequency (rad/sec)
2
10
10
Figura 66.: Diagramma di Bode del modulo della funzione di sensitività S(s)
Bode Diagram
100
Magnitude (dB)
0
−100
−200
−300
−400
−500
0
Phase (deg)
−90
−180
−270
−360
−450
−2
10
0
2
10
10
Frequency (rad/sec)
4
10
Figura 67.: Diagramma di Bode della funzione di trasferimento H(s) T (s)
95
5.4 controllo lqg
Nyquist Diagram
1.5
Imaginary Axis
1
0.5
0
−0.5
−1
−1.5
−2
−1.5
−1
−0.5
0
Real Axis
0.5
1
1.5
2
Figura 68.: Diagramma di Nyquist della funzione d’anello L(s)
Bode Diagram
100
L 1 (s)
Magnitude (dB)
50
L 2 (s)
0
L 3 (s)
−50
L 4 (s)
−100
L 0 (s)
−150
−200
−45
Phase (deg)
−90
−135
−180
−225
−270
−4
10
−2
10
0
2
10
10
Frequency (rad/sec)
4
10
6
10
Figura 69.: Diagrammi di Bode delle funzioni anello ottenute per i diversi punti di
lavoro. Le funzioni di trasferimento nel diagramma sono indicate con Li (s)
con i ∈ {0, 1, 2, 3, 4} a seconda del punto di lavoro considerato, in maniera
consistente alla numerazione riportata in (Tabella 1)
96
5.5
5.5 scelta del controllore più soddisfacente
scelta del controllore più soddisfacente
In questa sezione si vuole scegliere tra tutti i controllori progettati, quello che si
ritiene possa essere il più soddisfacente. Per tale valutazione si prenderanno in
considerazione, per i tre controllori, gli andamenti delle funzioni di sensitività
S(s), delle funzioni di trasferimento tra la variazione del riferimento e l’errore
di regolazione H(s) T (s), i margini di stabilità e la banda passante nel caso nominale.
97
5.5 scelta del controllore più soddisfacente
Margini di stabilità e banda passante
Tabella 5.: Valore dei margini di guadagno e di fase ottenuti per i punti di lavoro
scelti con il controllo LQG. Il punto di lavoro nominale viene indicato con
il numero 0.
Controlore
mg
db
mϕ
(deg)
PID
LQG
PID+Rete
∞
48,3
∞
99,2
65,2
84,6
ωb
(rad s−1 )
0,5721
0,5674
0,5550
Funzione di sensitività
Bode Diagram
20
LQG
PID
PID+Rete
10
0
Magnitude (dB)
−10
−20
−30
−40
−50
−60
−4
10
−3
10
−2
10
−1
0
10
10
Frequency (rad/sec)
1
10
2
10
3
10
Figura 70.: Diagrammi di Bode del modulo della funzione di sensitività per i tre
controllori nel caso nominale
Funzione di trasferimento H(s) T (s)
Osservando la (Figura 70), la (Figura 71) e la (Tabella 5) , si può evincere che
il controllore PID con rete compensatrice, rappresenta tra tutte la scelta per la
quale si ottiene il miglior compromesso tra robustezza e performance, in particolare tale controllore garantisce delle buone prestazioni, dei buoni margini di
robustezza e una buona reiezione dei disturbi a bassa frequenza. Pertanto al
fine di disporre di una realizzazione del controllore candidato che sia pronta
all’uso, se ne deriva una versione a tempo discreto che può essere facilemente
98
5.5 scelta del controllore più soddisfacente
Bode Diagram
Magnitude (dB)
50
LQG
PID
PID+Rete
0
−50
−100
Phase (deg)
−150
0
−90
−180
−270
−360
−3
10
−2
10
−1
10
Frequency (rad/sec)
0
10
1
10
Figura 71.: Diagrammi di Bode della funzione di trasferimento tra la variazione del
riferimento di temperatura e l’errore di regolazione H(s) T (s) per i tre
controllori nel caso nominale
implementata con un dispositivo digitale. In particolare scegliendo un tempo
di campionamento Ts = 0.01 s (un tempo certamente ragionevole considerando
un dispositivo di elaborazione di costo medio-basso) e discretizzando il regolatore in questione con il metodo di Tustin, si ottiene il regolatore a tempo
discreto descritto dalla seguente funzione di trasferimento.
Cd (z) =
0.323z3 − 0.9663z2 + 0.9637z − 0.3204
z3 − 2.858z2 + 2.716z − 0.8582
Al fine di stabilire se le performance del sistema di controllo a ottenuto con il
regolatore a tempo continuo siano sufficientemente recuperate in (72), si riportano i diagrammi di Bode delle funzioni di trasferimento Cd (z) e del regolatore
originario C (s)
99
5.5 scelta del controllore più soddisfacente
Bode Diagram
0
C (s)
C d( z )
−5
−10
Magnitude (dB)
−15
−20
−25
−30
−35
−40
−2
10
−1
10
0
1
10
10
Frequency (rad/sec)
2
10
3
10
Figura 72.: Diagrammi di Bode delle funzioni di trasferimento C (s) e Cd (z) a confronto
100
6
A L C U N I R I S U LT A T I N U M E R I C I
In questo capitolo vengono mostrati alcuni risultati ottenuti utilizzando i controllori sintetizzati nel corso della presente trattazione. In particolare vengono
valutate le risposte del sistema di controllo a variazioni a gradino della temperatura di riferimento e a variazioni a gradino delle temperatura del liquido in
ingresso e della temperatura del vapore. Nello specifico vengono imposte le
variazioni illustrate in (Figura 73) e in (Figura 74).
102
100
◦
C
98
96
94
92
90
0
100
200
300
400
500
600
700
800
900
1000
ti me [ s]
Figura 73.: Andamento del riferimento y∗ (t)
160
ϑ v (t )
ϑ l i(t )
140
◦
C
120
100
80
60
40
0
100
200
300
400
500
600
700
800
900
1000
ti me [ s]
Figura 74.: Andamento delle variabili ϑli (t) e ϑv (t)
A fronte degli andamenti imposti per le grandezze esogene, controllando il
processo con i controllori ottenuti nel corso della trattazione, si ottengono i
risultati mostrati in (Figura 75) e in (Figura 76).
I risultati mostrati di cui sopra non evidenziano molto bene l’andamento della risposta transitoria del sistema di controllo, pertanto nei paragrafi seguenti
101
6.1 controllo pid con azione in avanti
115
PID
PID+RETE
LQG
PID+RETE discreto
110
◦
C
105
100
95
90
85
0
100
200
300
400
500
600
700
800
900
1000
ti me [ s]
Figura 75.: Andamento dell’uscita controllata y(t)
0.8
PID
PID+RETE
LQG
PID+RETE discreto
0.7
m s −1
0.6
0.5
0.4
0.3
0.2
0.1
0
100
200
300
400
500
600
700
800
900
1000
ti me [ s]
Figura 76.: Andamento del controllo u(t)
verrano mostrati singolarmente, gli andamenti delle risposte transitorie di cui
sopra, ottenute con tutti i controllori sintetizzati.
6.1
controllo pid con azione in avanti
Risultati ottenuti con il controllo PID.
102
6.1 controllo pid con azione in avanti
99
P l a nt Ou tp u t
R eferen ce
98
[ ◦C ]
97
96
95
94
93
60
80
100
120
140
160
180
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.39
0.38
v (t )[ m s −1]
0.37
0.36
0.35
0.34
0.33
0.32
0.31
0.3
90
100
110
120
130
140
150
160
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 77.: Variazione delle temperatura di riferimento da 98 ◦ a 94 ◦
103
6.1 controllo pid con azione in avanti
100
P l a nt Ou tp u t
R eferen ce
99
98
[ ◦C ]
97
96
95
94
93
190
200
210
220
230
240
250
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
v (t )[ m s −1]
0.35
0.3
0.25
200
210
220
230
240
250
260
270
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 78.: Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
104
6.1 controllo pid con azione in avanti
100
P l a nt Ou tp u t
R eferen ce
99
98
[ ◦C ]
97
96
95
94
93
310
320
330
340
350
360
370
380
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.3
0.28
0.26
v (t )[ m s −1]
0.24
0.22
0.2
0.18
0.16
0.14
0.12
0.1
320
325
330
335
340
345
350
355
360
365
370
375
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 79.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
105
6.1 controllo pid con azione in avanti
108
P l a nt Ou tp u t
R eferen ce
106
[ ◦C ]
104
102
100
98
96
510
520
530
540
550
560
570
580
590
600
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.7
0.6
v (t )[ m s −1]
0.5
0.4
0.3
0.2
0.1
520
530
540
550
560
570
580
590
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 80.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
106
6.1 controllo pid con azione in avanti
P l a nt Ou tp u t
100
R eferen ce
98
[ ◦C ]
96
94
92
90
88
690
700
710
720
730
740
750
760
770
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.65
v (t )[ m s −1]
0.6
0.55
0.5
0.45
0.4
0.35
0.3
700
710
720
730
740
750
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 81.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
107
6.1 controllo pid con azione in avanti
104
P l a nt Ou tp u t
102
R eferen ce
100
[ ◦C ]
98
96
94
92
90
88
835
840
845
850
855
860
865
870
875
880
885
890
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.6
0.55
0.5
v (t )[ m s −1]
0.45
0.4
0.35
0.3
0.25
0.2
0.15
0.1
835
840
845
850
855
860
865
870
875
880
885
890
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 82.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
108
controllo pid con rete compensatrice con
azione in avanti
99
P l a nt Ou tp u t
R eferen ce
98
97
96
95
94
93
92
90
100
110
120
130
140
150
160
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.39
0.38
0.37
v (t )[ m s −1]
6.2
6.2 controllo pid con rete compensatrice con azione in avanti
0.36
0.35
0.34
0.33
0.32
0.31
0.3
90
100
110
120
130
140
150
160
170
180
190
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 83.: Variazione delle temperatura di riferimento da 98 ◦ a 94 ◦
109
6.2 controllo pid con rete compensatrice con azione in avanti
100
P l a nt Ou tp u t
99
R eferen ce
98
97
96
95
94
93
92
190
200
210
220
230
240
250
260
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.38
0.36
v (t )[ m s −1]
0.34
0.32
0.3
0.28
0.26
0.24
0.22
0.2
190
200
210
220
230
240
250
260
270
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 84.: Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
110
6.2 controllo pid con rete compensatrice con azione in avanti
101
P l a nt Ou tp u t
100
R eferen ce
99
98
97
96
95
94
93
320
330
340
350
360
370
380
390
400
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.3
v (t )[ m s −1]
0.25
0.2
0.15
0.1
310
320
330
340
350
360
370
380
390
400
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 85.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
111
6.2 controllo pid con rete compensatrice con azione in avanti
110
P l a nt Ou tp u t
108
R eferen ce
106
104
102
100
98
96
94
520
530
540
550
560
570
580
590
600
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.7
0.6
v (t )[ m s −1]
0.5
0.4
0.3
0.2
0.1
510
520
530
540
550
560
570
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 86.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
112
6.2 controllo pid con rete compensatrice con azione in avanti
102
P l a nt Ou tp u t
R eferen ce
100
98
96
94
92
90
88
690
700
710
720
730
740
750
760
770
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.6
v (t )[ m s −1]
0.55
0.5
0.45
0.4
0.35
0.3
700
710
720
730
740
750
760
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 87.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
113
6.2 controllo pid con rete compensatrice con azione in avanti
102
P l a nt Ou tp u t
R eferen ce
100
98
96
94
92
90
840
850
860
870
880
890
900
910
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.35
v (t )[ m s −1]
0.3
0.25
0.2
0.15
0.1
840
850
860
870
880
890
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 88.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
114
controllo lqg con azione in avanti
99
P l a nt Ou tp u t
R eferen ce
98
97
96
95
94
93
90
100
110
120
130
140
150
160
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.38
0.36
0.34
v (t )[ m s −1]
6.3
6.3 controllo lqg con azione in avanti
0.32
0.3
0.28
0.26
0.24
0.22
0.2
90
100
110
120
130
140
150
160
170
180
190
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 89.: Variazione delle temperatura di riferimento da 98 ◦ a 94 ◦
115
6.3 controllo lqg con azione in avanti
100
P l a nt Ou tp u t
R eferen ce
99
98
97
96
95
94
93
200
210
220
230
240
250
260
270
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
v (t )[ m s −1]
0.35
0.3
0.25
200
210
220
230
240
250
260
270
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 90.: Variazione delle temperatura di riferimento da 94 ◦ a 98 ◦
116
6.3 controllo lqg con azione in avanti
105
P l a nt Ou tp u t
R eferen ce
100
95
90
320
340
360
380
400
420
440
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.4
0.35
v (t )[ m s −1]
0.3
0.25
0.2
0.15
0.1
320
340
360
380
400
420
440
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 91.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (65◦ , 120◦ ) a (75◦ , 140◦ )
117
6.3 controllo lqg con azione in avanti
P l a nt Ou tp u t
R eferen ce
110
105
100
95
90
510
520
530
540
550
560
570
580
590
600
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.6
v (t )[ m s −1]
0.5
0.4
0.3
0.2
0.1
510
520
530
540
550
560
570
580
590
600
610
620
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 92.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 140◦ ) a (50◦ , 140◦ )
118
6.3 controllo lqg con azione in avanti
P l a nt Ou tp u t
104
R eferen ce
102
100
98
96
94
92
90
88
690
700
710
720
730
740
750
760
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.65
0.6
0.55
v (t )[ m s −1]
0.5
0.45
0.4
0.35
0.3
0.25
0.2
690
700
710
720
730
740
750
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 93.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (50◦ , 140◦ ) a (75◦ , 105◦ )
119
6.3 controllo lqg con azione in avanti
P l a nt Ou tp u t
104
R eferen ce
102
100
98
96
94
92
90
88
840
850
860
870
880
890
900
910
Ti me [ s]
(a) Andamento della temperatura del fluido in uscita misurata ϑlu
0.5
0.45
v (t )[ m s −1]
0.4
0.35
0.3
0.25
0.2
0.15
0.1
840
850
860
870
880
890
Ti me [ s]
(b) Andamento della variabile di controllo v
Figura 94.: Variazione della temperatura del fluido in ingresso e della temperatura del
vapore da (75◦ , 105◦ ) a (50◦ , 105◦ )
120
7
L
CONCLUSIONI
’esito del progetto si può ritenere positivo in quanto la regolazione
del processo avviene con successo anche in condizioni operative
molto diverse dalla nominale. Tuttavia l’approccio utilizzato non consente di
certificare quale sia l’insieme delle condizioni operative per le quali l’obiettivo
della regolazione venga conseguito. Un altro grosso problema è quello della
saturazione dell’attuatore che non è stata presa in considerazione in fase di
sintesi del sistema di controllo. Pertanto come sviluppi futuri si potrebbe pensare di studiare il problema della stabilizzazione robusta del processo con un
approccio basato sui dei criteri di stabilità robusta in presenza di incertezze parametriche, mentre per quanto riguarda la saturazione si può sintetizzare uno
schema di controllo con azione di anti-windup.
121
A
STRUMENTI NECESSARI PER LA
S I N T E S I D E I R E G O L ATO R I
In questa appendice verranno illustrati gli strumenti concettuali utilizzati nella
presente trattazione, per la sintesi dei regolatori mostrati.
Lemma A.1:
Sia K un generico sistema LTI nella forma:
ż = A z + B u
(A.1)
z ∈ Rn
(A.2)
con:
u∈R
(A.3)
dunque se A è singolare e se ( A, B) stabilizzabile, allora preso ( x̄, ū) tale che:
A x̄ + B ū = 0
necessariamente deve risultate:
ū = 0
Dimostrazione. Prima di ogni cosa si noti che essendo ( A, B) stabilizzabile si ha
necessariamente che:
rank[ A B] = n
e quindi segue anche che:
dim{Ker A} = 1
pertanto il seguente sistema equazioni ammette ∞1 soluzioni. Ossia:
x̄
=0
A B
ū
(A.4)
ebbene in virtù del fatto che la matrice A è singolare si ha che:
∗
x
∗
∀ x ∈ Ker A : [ A B]
=0
0
tuttavia dato che:
dim{Ker A} = 1
si ha che tutte e sole le soluzioni del sistema (A.4) sono del tipo ( x ∗ , 0) con
x ∗ ∈ Ker A, infatti se per assurdo esistesse un soluzione del sistema (A.4) del
tipo: ( x̂ û) ∈ Rn con û 6= 0 allora vuol dire il sottospazio delle soluzioni per il
sistema (A.4) ammetterebbe come base B il seguente insieme dei vettori:
( )
x∗
0
B=
,
0
û
122
strumenti necessari per la sintesi dei regolatori
ma ciò contraddirebbe il fatto che il sottospazio delle soluzioni del sistema (A.4)
sia di dimensione unitaria, contraddicendo in tal modo l’ipotesi di stabilizzabilità della coppia ( A, B).
Lemma A.2:
In merito al sistema in (4.3), si fissino y∗ ∈ R, d∗ ∈ R p in modo qualunque ma a patto
che valga la Assunzione 4.1 con ū 6= 0. Si consideri ora la seguente legge di controllo
L.T.I. in retroazione dall’errore di regolazione ỹ = y∗ − y:
ϕ̇ = Ψ ϕ + Λ ỹ
(A.5)
u = K1 ϕ + K2 ỹ
Con:
ϕ ∈ R`
K1 ∈ R
(A.6)
1×`
(A.7)
K2 ∈ R
(A.8)
Λ ∈ R`
(A.10)
Ψ ∈ R`×`
(A.9)
La quale fornisce la seguente dinamica per il sistema di controllo:
ẋ = f ( x, ϕ, d∗ y∗ )
(A.11)
ϕ̇ = Ψ ϕ + Λ (y∗ − h( x ))
Ebbene se sotto l’azione di tale controllo esiste, per il sistema di controllo (A.11), uno
stato di equilibrio ( x̂ , ϕ̂) ∈ R`+n tale che:
h( x̂ ) = y∗
allora valgono le seguenti affermazioni:
• ū = K1 ϕ̂.
• la matrice Ψ ammette λ = 0 come autovalore ed inoltre ϕ̂ ∈ Ker (Ψ)
•
K1T ∈
/ Im Ψ
Dimostrazione. Se lo stato ( x̂ , ϕ̂) è di equilibrio per il sistema (A.11) allora si ha:
f ( x̂, ϕ̂, d∗ y∗ ) = 0
Ψ ϕ̂ + Λ (y∗ − h( x̂ )) = 0
(A.12)
Tuttavia in vista dell’Assunzione 4.1 si ha necessariamente:
x̄ = x̂
ū = K1 ϕ̂
(A.13)
123
strumenti necessari per la sintesi dei regolatori
Provando in tal modo in primo e il secondo punto dell’enunciato.
Allo scopo di completare la dimostrazione del lemma si osservi che le espressioni in (A.12) si modificano come segue:
f ( x̂, ϕ̂, d∗ y∗ ) = 0
(A.14)
Ψ ϕ̂ = 0
(A.15)
ossia il vettore ϕ̂ deve soddisfare il seguente sistema di equazioni lineari:
Ψ
0
ϕ̂ = `×1
K1
ū
(A.16)
Per cui affinché esista almeno una soluzione di tale sistema di equazioni, occorre e basta, per il teorema di Rouchè-Capelli che:
Ψ
rank
K1
Ψ 0`×1
= rank
K1 ū
(A.17)
ebbene si osservi ora che dato che per ipotesi ū allora:
Ψ 0`×1
= rank Ψ + 1
rank
K1 ū
(A.18)
pertanto dalla (A.17) si ha che le colonne della matrice Ψ sono linearmente
indipendenti dal vettore K1T concludendo la dimostrazione del lemma.
Osservazione A.1:
Il lemma appena enunciato e provato, di fatto indica che un qualunque controllore L.T.I.
che risolva il problema di regolazione contiene necessariamente una azione integrale.
Tale proprietà è ben nota ed in particolare essa è denominata in alcuni testi con in nome
di Principio del modello interno, facendo riferimento al fatto che il sistema di controllo
contiene al suo interno un modello degli ingressi esogeni da compensare, siano essi
riferimenti o disturbi, i quali nel caso del problema di regolazione sono di tipo costante,
ossia descritti da un sistema dinamico di tipo integrale non forzato. Infine si noti che
nelle ipotesi del lemma si assume che ū 6= 0, in caso contrario il lemma perde la
sua validità di condizione necessaria, tuttavia essendo in generale ū qualunque, tale
condizione rimane necessaria per gli scopi preposti. Si noti infine che la condizione sul
rango enunciata nel lemma implica che se il sistema (A.5) è rilevabile allora può esistere
uno e un solo valore di ϕ̂ per il quale vale la proprietà di regolazione.
Prima di proseguire è conveniente enunciare un ulteriore risultato.
Proposizione A.1:
Il sistema (A.5) è rilevabile soltanto se:
rank [Ψ − λ I] ≥ ` − 1 ∀ λ ∈ C+
Dimostrazione. La dimostrazione dell’assunto è immediata e discende direttamente dal P.B.H. test di rilevabilità. Infatti il sistema (A.5) è rilevabile se e solo
se:
Ψ−λI
rank
= ` ∀ λ ∈ C+
K1
124
strumenti necessari per la sintesi dei regolatori
ma ciò è possibile solo se:
rank Ψ − λ I ≥ ` − 1 ∀ λ ∈ C+
dimostrando l’assunto.
Proposizione A.2:
Sia data la legge di controllo (A.5), allora se sotto l’azione di tale controllo esiste, per il
sistema di controllo (A.11), un solo stato di equilibrio ( x̂ , ϕ̂) ∈ R`+n tale che:
h( x̂ ) = y∗
e se la coppia (Ψ, K1 ) è rilevabile allora esiste una trasformazione di coordinate invertibile ζ = T ϕ tramite la quale si può esprimere il sistema (A.5) nella seguente forma:
01×1 0
ζ + Λ J ỹ
ζ̇ =
0 Ψ̂
u=
K1J
(A.19)
ζ + K2 ỹ
Dimostrazione. Si indichi con JΨ la forma di Jordan della matrice Ψ, ebbene è
noto che certamente esiste una matrice quadrata e non singolare TJ tale che:
TJ Ψ TJ−1 = JΨ
Tuttavia dalle ipotesi per il Lemma 4.1 e per la Proposizione A.1, λ = 0 è
un autovalore della matrice Ψ con molteplicità geometrica unitaria, quindi la
forma di Jordan della matrice Ψ può essere scritta nel modo seguente:


0 1 0 ... 0
 .. . .

. . ..
 .

. 1
. .




..
O
. 0
 0 0 0



JΨ =  .. . . . . . .

.
.
.
 .

1


 0 0 0

0
0




O
ĴΨ
dove il blocco di Jordan relativo all’autovalore nullo è di dimensione pari alla molteplicità algebrica di tale autovalore. A fronte di quanto mostrato si
consideri il seguente cambio di coordinate:
ζ̄ = TJ x
il quale fornisce la seguente rappresentazione equivalente per il sistema (A.5):
0 1
 .. . .
 .
.


 0 0

ζ̄˙ =  .. . .
.
 .

 0 0



0
...
..
.
1
..
.
0
.. ..
.
.
0
0
O
u = K̄1J ζ̄ + K2 ỹ
0
..
.
0
1
0




O 


 ζ̄ + Λ̄ J ỹ





ĴΨ
(A.20)
125
strumenti necessari per la sintesi dei regolatori
Infine attraverso l’ulteriore trasformazione:
ζ = T̂ ζ̄
con:

0







T̂ = 







0
0
0
1
..
.
0
.. . .
.
.
0
1
..
.
0
1
0
0
0
0
0
..
.
1
..
.
..
0
.
..
. 0
0 1
... 0
O






O 









I
Si ottiene per il sistema (A.5) la seguente rappresentazione:


0 0 0
0 0
 . .

 .. . . 1 . . . ...





.
.
 0 0 0

. 0

O 


.. .. ..
 ζ + Λ J ỹ
.
.
. 0
ζ̇ = 
 1

 0 0 0

0 0




0
1
0
.
.
.
0






O
ĴΨ
(A.21)
u = K1J ζ + K2 ỹ
In definitiva dato che la trasformazione:
ζ = T̂ Tj x
conduce il sistema (A.5) nella forma (A.21) la dimostrazione dell’enunciato è
conclusa.
Osservazione A.2:
Si osservi che nella rappresentazione (A.21) ottenuta per il sistema (A.5) la matrice ĴΨ
è in forma di Jordan, tuttavia ciò non da luogo a nessuna perdita di generalità in quanto
è possibile, attraverso un ulteriore cambio di coordinate, portare la matrice ĴΨ in una
forma del tutto generale.
Teorema A.1:
In merito al sistema in (4.3), si fissino y∗ ∈ R, d∗ ∈ R p arbitrariamente a patto che
valga la Assunzione 4.1. Ebbene in si indichino ora rispettivamente:
∂ f ( x , u, d∗ ) A :=
∂x
( x ,u)=( x̄ ,ū)
∗
∂ f (x , u , d ) (A.22)
B :=
∂u
( x ,u)=( x̄ ,ū)
∂h( x ) C :=
∂x (x ,u)=( x̄ ,ū)
dunque se:
126
strumenti necessari per la sintesi dei regolatori
• la terna ( A, B, C ) è stabilizzabile e rilevabile
• rank
A B
= n+1
C 0
allora esiste un controllore K nella forma (4.10) che sia soluzione del problema (4.3) ed
inoltre esso è espresso nella forma seguente:


 ϕ̇1 = λ1 ỹ
(A.23)
K : ϕ̇2 = Ψ̂ ϕ2 + λ2 ỹ


u = k11 ϕ1 + K12 ϕ2 + K2 ỹ
con:
ϕ1 ∈ R
ϕ2 ∈ R`−1
λ1 ∈ R
(A.24)
(A.25)
(A.26)
λ2 ∈ R`−1
(A.27)
K12 ∈ R1×`−1
(A.29)
k11 ∈ R
(A.28)
(A.30)
Dimostrazione. La dimostrazione del risultato appena enunciato è di tipo costruttivo e conduce alla costruzione del controllore che risolve il problema di
regolazione in questione.
Si consideri l’espressione della dinamica ad anello chiuso, ottenuta controllando il sistema (4.3) con il controllore (A.23) :
ẋ = f ( x, ϕ, d∗ y∗ )
ϕ̇ = Ψ ϕ + Λ (y∗ − h( x ))
(A.31)
Ebbene affinché l’errore di regolazione possa annullarsi asintoticamente come
richiesto, è necessario che il sistema di controllo ammetta uno stato di equilibrio
( x̂ , ϕ̂) ∈ R`+n tale che:
h( x̂ ) = y∗
Dunque a fronte del Lemma A.2 è necessario per la sintesi del controllore che:
0 ∈ σ (Ψ) Ebbene a questo punto essendo arbitraria la molteplicità algebrica
dell’autovalore nullo della matrice Ψ si può scegliere essa unitaria e quindi
imponendo evidentemente che il compensatore sia rilevabile, dal lemma A.2,
partizionando i vettori ϕ e Λ come segue:
ϕ := [ ϕ1 , ϕ2 T ] T
(A.32)
Λ := [λ1 , λ2 T ] T
(A.33)
K1 := [k11 , K12 ]
(A.34)
si può scrivere il controllore (A.23) in questo modo:
ϕ̇1 = λ1 ỹ
ϕ̇2 = Ψ̂ ϕ2 + λ2 ỹ
u = k11 ϕ1 + K12 ϕ2 + K2 ỹ
(A.35)
127
strumenti necessari per la sintesi dei regolatori
dove:
ϕ1 ∈ R
ϕ2 ∈ R
λ1 ∈ R
λ2 ∈ R
(A.36)
`−1
(A.37)
(A.38)
`−1
k11 ∈ R
K12 ∈ R
(A.39)
(A.40)
1×`−1
(A.41)
(A.42)
In tal modo la dinamica del sistema di controllo può essere riscritta in tal modo:
Dinamica del sistema di controllo
ẋ = f ( x, ϕ1 , ϕ2 , d∗ , y∗ )
ϕ̇1 = λ1 (y∗ − h( x ))
(A.43)
∗
ϕ̇2 = Ψ̂ ϕ2 + λ2 (y − h( x ))
Prima di procedere si noti che la scelta fatta per la struttura del controllore
impone che se esso è stabilizzabile e rilevabile, allora k11 6= 0 e λ1 6= 0. Si
osservi ora che a fronte dell’Assunzione 4.1 e del Lemma 4.1
∃( x̄ , ϕ̄1 , ϕ̄2 )
tale che:
f ( x̄, ϕ̄1 , ϕ̄2 , d∗ , y∗ ) = 0
Ψ̂ ϕ̄2 = 0
ū = k11 ϕ̄1 + K12 ϕ̄2
(A.44)
h( x̄ ) = y∗
Ossia il sistema di controllo ammette come punto di equilibrio lo stato che garantisce un errore di regolazione identicamente nullo, pertanto per conseguire
l’obiettivo della regolazione locale, occorre e basta provvedere alla stabilizzazione asintotica, almeno locale, di tale equilibrio. Pertanto come noto dalla teoria
dei sistemi dinamici non lineari, al fine di conseguire l’obiettivo della stabilità
esponenziale locale di un dato equilibrio, è sufficiente rendere asintoticamente
stabile l’approssimazione lineare del sistema in esame attorno all’equilibrio considerato. A tal scopo si consideri la linearizzazione, attorno al punto ( x̄ , ϕ̄1 , ϕ̄2 ),
dei termini a secondo membro della (A), la quale porge le seguenti espressioni:
x̃˙ = ( A − B K2 C ) x̃ + B (k11 ϕ̃1 + K12 ϕ̃2 )
ϕ̃˙ 1 = −λ1 C x̃
ϕ̃˙ 2 = Ψ̂ ϕ̃2 − λ2 C x̃
(A.45)
128
strumenti necessari per la sintesi dei regolatori
dove:
∂ f ( x , u d) ∂x
( x ,u)=( x̄ ,ū)
∂ f ( x , u) B :=
∂u
( x ,u)=( x̄ ,ū)
∂h( x ) C :=
∂x (x ,u)=(x̄ ,ū)
A :=
(A.46)
(A.47)
(A.48)
x̃ = x − x̄
(A.49)
ϕ̃1 = ϕ1 − ϕ̄1
(A.51)
ũ = u − ū
(A.50)
ϕ̃2 = ϕ2 − ϕ̄2
(A.52)
(A.53)
Al fine di stabilizzare il sistema ottenuto, ci si può giovare della arbitrarietà del
sottosistema caratterizzato dallo stato ϕ2 . Infatti si può osservare che il sistema (A.45) può essere riscritto come l’interconnessione di due sistemi L.T.I., di
cui uno dei due, il quale verrà denominato processo allargato, è composto dall’impianto e dalla componente della legge di controllo caratterizzata dallo stato
ϕ1 ; mentre l’altro,che verrà indicato come controllore stabilizzante, è composto
dalla componente della legge di controllo caratterizzata dallo stato ϕ2 e verrà
sintetizzato arbitrariamente con il solo scopo di conseguire l’obiettivo della stabilità asintotica dell’intero sistema di controllo (A.45). Quanto detto conduce
alla seguente scrittura per il sistema (A.45).
Processo Allargato
x̃˙ = A x̃ + B k11 ϕ̃1 + B v
ϕ̃˙ 1 = −λ1 C x̃
(A.54)
ỹ = −C x̃
Controllore stabilizzante
ϕ̃˙ 2 = Ψ̂ ϕ̃2 + λ2 ỹ
(A.55)
Per concludere quindi la prova del teorema occorre e basta mostrare che il
processo allargato (A.54) sia stabilizzabile e rilevabile. Per cui si consideri il
processo allargato riscritto in forma compatta come segue:
x̃˙
A
B k11
x̃
B
=
+
v
ϕ̃1
0
ϕ̃˙ 1
− λ1 C
0
x̃
ỹ = −C 0
ϕ̃1
(A.56)
a questo punto il processo allargato è stabilizzabile se e soltanto se:
A − α I B k11 B
rank
= n + 1 ∀ α ∈ C+
− λ1 C − α 0
(A.57)
129
strumenti necessari per la sintesi dei regolatori
come noto, la condizione sopra indicata va verificata solamente per i valori di
α che sono autovalori dell’impianto allargato aventi parte reale non negativa,
ossia in questo caso, per α = 0 e per tutti gli α ∈ C+ che sono autovalori di A.
Perciò per α = 0 la condizione (A.57) diviene ha:
A
B k11 B
rank
= n+1
(A.58)
− λ1 C
0
0
ma sotto l’ipotesi che:
A
B
rank
= n+1
− λ1 C 0
(A.59)
ebbene si vede subito che tale condizione è certamente verificata in quanto:
A
B k11 B
A B
rank
= rank
= n+1
(A.60)
− λ1 C
0
0
−C 0
Resta dunque da mostrare che la condizione di cui sopra sia verificata anche
per tutti gli α 6= 0 ∈ C+ che sono autovalori di A. Ma ciò è immediato da
constatare in quanto data l’ipotesi di stabilizzabilità della coppia ( A, B) tutte le
righe della matrice:
A − α I B B k11
(A.61)
rank
− λ1 C − α
0
sono linearmente indipendenti gli α 6= 0 ∈ C+ che sono autovalori di A. Infatti le prime n sono linearmente indipendenti, per l’ipotesi di stabilizzabilità
della coppia ( A, B) e l’ultima riga di tale matrice è linearmente indipendente
da quelle soprastanti dato che l’elemento di posto (n + 1 , n + 2) è nullo, dimostrando di fatto che il rango della matrice considerata sia proprio n + 1.
Rimane dunque da dimostrare che l’impianto allargato sia rilevabile dall’errore
di regolazione ỹ. Ossia bisogna verificare che:

A − α I B k11
rank  −λ1 C −α  = n + 1 ∀ α ∈ C+
−C
0

(A.62)
La verifica di tale proprietà va condotta per i valori di α indicati nel caso precedente. In particolare per α = 0 si ha che la condizione (A.62) è banalmente
verificata sotto le ipotesi che:
A
B
rank
= n+1
(A.63)
− λ1 C 0
in quanto:


A
B k11
A
B
rank −λ1 C
= n+1
0  = rank
−C 0
−C
0
(A.64)
Mentre ragionando in maniera duale, rispetto al caso della proprietà di stabilizzabilità analizzata in precedenza, si giunge anche in questo caso, sotto l’ipotesi
di rilevabilità della coppia ( A, C ) alla conclusione che il processo allargato sia
130
strumenti necessari per la sintesi dei regolatori
rilevabile, pertanto è possibile sintetizzare un controllore stabilizzante in retroazione da ỹ, ottenendo in tal modo la stabilità esponenziale locale dell’equilibrio
di regolazione. Per concludere la dimostrazione resta da provare che il controllo u(t) sia limitato, ma ciò è immediato in quanto tale controllo per struttura
è combinazione lineare di variabili di stato che sono limitate per condizioni
iniziali contenute nel bacino di attrazione del punto di equilibrio in questione, concludendo di fatto la dimostrazione dell’enunciato. Prima di concludere
la dimostrazione sembra opportuno ribadire che le condizioni di rilevabilità e
stabilizzabilità in questo caso sono solamente sufficienti per la soluzione del
problema in quanto, nulla si può dire se la matrice dinamica del sistema di
controllo contiene uno o più autovalori a parte reale nulla.
131
BIBLIOGRAFIA
[1]
Francesco Ferrante. “I regolatori standard PID: Il punto di vista delle teoria
del controllo”. In: (2011).
[2]
Hassan Khalil. Nonlinear Systems, Third Edition. Prentice-Hall, 2002.
[3]
Riccardo Marino. Appunti del corso di Controllo dei Processi. 2011.
132
Scarica

Analisi e controllo di uno scambiatore di calore