Proceedings of International Scientific Conference of FME Session 4: Automation Control and Applied Informatics
Paper 41
Identifikace dynamických vlastností soustavy s ruční zpětnou vazbou TŮMA, Jiří1 1
Doc.Ing.CSc., VŠB - TU Ostrava, Fakulta strojní, Katedra ATŘ, tř. 17. listopadu, 708 33 Ostrava - Poruba,
[email protected]
Abstrakt: Referát se zabývá identifikací dynamického modelu soustavy ve tvaru diferenční rovnice. Soustava je během identifikace řízena ručně, což způsobí vychýlení identifikovaných parametrů. Na soustavu působí náhodné poruchy, které znemožní využít klasické postupy identifikace testovacími signály, a ruční řízení nelze zastavit. V příspěvku je navržen postup identifikace včetně způsobu korekce vychýlení a identifikace ruční zpětné vazby. Klíčová slova: identifikace, ruční zpětná vazba, korelační funkce
1. Úvod Před zautomatizováním řízení jsou mnohé technologické procesy řízeny ručně. Účinnost řízení, která je hodnocena například kvadratickým kritériem, dosahuje jisté úrovně. Při identifikaci proto není žádoucí ruční řízení zastavit. Rovněž je nepřijatelné použít klasické způsoby identifikace, které jsou založeny na deterministických nebo náhodných testovacích signálech. Při průběžné identifikaci byla teoretiky zavedena randomizace akční veličiny, která nahrazovala testovací náhodný signál [Peterka, 1975]. Při ručním řízení lze za randomizující prvek považovat rozhodování obsluhy. Popis postupu identifikace se bude opírat o příklad dat pro výpočet dílčího modelu vysoké pece, jmenovitě závislosti obsahu křemíku (Si) v surovém železe na přídavku vodní páry do větru (w). Složení surového železa je určováno ze vzorků odebraných v průběhu odpichu a jejich frekvenci je přizpůsobeno průměrování přídavku páry. Zmíněný obsah křemíku je regulovaná veličina, která se ručně ovládá změnami přídavku vodní páry do větru. Statistické vlastnosti dat lze posoudit podle korelačních funkcí na obr. 1. Průběh autokorelační funkce Si x Si v okolí nulového posunutí ukazuje nejen na přítomnost aditivního bílého šumu ve významu chyby měření nebo nereprezentativnosti odebraného vzorku surového železa k analýze na kvantometru, ale také malou „setrvačnost“ změn složení surového železa. Naproti tomu změny přídavku vodní páry do větru vykazují podle autokorelační funkce w x w velmi vysokou setrvačnost, tj. jen pomalé změny. Experimentálním důkazem významu zpětné vazby je průběh vzájemné korelační funkce w x Si mezi centrovanými řadami obsahů Si v surovém železe a centrovanou řadou přídavků vodní páry do větru. Teoreticky je vzájemná korelační funkce mezi výstupem a vstupem soustavy pro Diracův impuls na vstupu roven impulsní charakteristice, která je z důvodu kausality rovna nule pro záporná posunutí. Toto tvrzení platí pro jakýkoliv vstupní signál za podmínky, že soustava není součástí uzavřené smyčky. Vysoký stupeň korelace (v obr. 1 označeno kroužkem) pro záporné posunutí mezi oběmi řadami, jako zdánlivé porušení kausality mezi vstupem a výstupem soustavy, signalizuje závislost vstupu soustavy na jeho výstupu. Tato závislost vyplývá ze přítomnosti zpětné vazby.
1 w x 0.5 w 0 -0.5 -100
-80
-60
-40
-20
0
20
40
60
80
100
-60
-40
-20
0
20
40
60
80
100
-60
-40
40
60
80
100
1
Si x Si
0.5 0 -0.5 -100 -80 4 x 10 1 0.5
wx Si
0 -0.5 -1 -100
-80
-20 0 20 posunuti v poctu odpichu
Obr. 1. Autokorelační funkce přídavku páry a obsahu křemíku v surovém železe a jejich vzájemná korelační funkce
2. Vliv zpětné vazby na zkreslení parametrů modelu Nechť model identifikované soustavy popisuje následující výchozí diferenční rovnice y k = a1 y k −1 + a 2 y k − 2 + ... + b0 u k ,
(1)
kde y označuje regulovanou veličinu a u akční veličinu, které jsou vzorkovány v čase k = 0, 1, 2,... . Parametry modelu jsou a1 , a 2 ,... a b0 , který má význam zesílení soustavy, tj. vlivnost x na velikost y. Obraz diferenční rovnice (1) v Z-transformaci je následující b0 Y (z ) = . −1 U (z ) 1 − a1 z − a 2 z − 2 − ...
(2)
Ruční zpětnou vazbu je vhodné popsat přenosem mezi regulační odchylkou ∆y a akční veličinou u U (z ) = k 0 + k1 z −1 + k 2 z − 2 + ... ., ∆Y (z )
(3)
kde k 0 , k1 , k 2 ,... jsou parametry ruční zpětné vazby. Součástí uzavřené smyčky je jeden krok dopravního zpoždění. Tato forma je zvolena proto, aby bylo možné posoudit váhu zpětných informací o analýzách křemíku, které obsluha velínu využije k rozhodování o velikosti přídavku páry do větru. Posloupnost parametrů přenosu (3) představuje impulsní odezvu ruční zpětné vazby. Blokové schéma identifikované soustavy se zpětnou vazbou je znázorněno na obr. 2. Podle tohoto schématu ovlivňuje regulační odchylka vstup identifikované soustavy.
Plant
u
y
b0 −1 1 − a1 z − a2 z − 2 − ... Delay z
−1
Manual Feed Back −1
ySP
−2
k0 + k1 z + k 2 z + ... ∆y
Obr. 2. Blokové schéma regulačního obvodu Působením zpětné vazby se přenos identifikované soustavy změní následujícím způsobem (4) b0 Y (z ) = . U (z ) −1 −2 T 1 − (a1 − b0 k 0 ) z − (a 2 − b0 k1 ) z − ... Místo výchozího modelu se tedy identifikuje jiný model. Diferenční rovnici tohoto modelu lze odvodit z předchozího přenosu (4) y k = (a1 − b0 k 0 ) y k −1 + (a 2 − b0 k1 ) y k − 2 + ... + b0 u k
(5)
a po úpravě značení y k = a1* y k −1 + a 2* y k − 2 + ... + b0 u k .
(6)
Poslední rovnice ukazuje, že identifikovaný model má zkreslené (vychýlené) parametry vlivu zpožděných hodnot výstupní veličiny (y) na její aktuální velikost. Pouze vlivnost b0 vstupní veličiny (u – obsah vodní páry ve větru) na výstupní veličinu (y – obsah Si v surovém železe ve výše uvedeném příkladě) je identifikována správně. Zkreslení parametrů dynamického modelu je závislé na přenosu zpětné vazby. Neznámé koeficienty a1* , a 2* ,... v poslední rovnici lze určit z experimentálních dat například známou metodou nejmenších čtverců. Pro konkrétní výpočet je třeba zvolit počet zpožděných hodnot regulované veličiny, které budou do tohoto modelu typu ARX zahrnuty. Kritérium pro rozhodnutí může být koeficient determinace nebo směrodatná odchylka chyby modelu (6). Nechť je identifikován model řádu N, tj. poslední nenulový parametr je a *N .
3. Kompenzace vlivu zpětné vazby na výpočet statistického modelu Jak je zřejmé z blokového schématu na obr. 2, přenos regulátoru (3) při působení zpětné vazby se změní následujícím způsobem U (z ) U (z ) ∆Y (z ) . ∆Y (z ) = T 1 + z −1 U (z ) Y (z ) ∆Y (z ) U (z )
(7)
Po dosazení detailních výrazů pro přenosové funkce (2) a (3) a přepisu parametrů podle vzorce (6) je výsledný tvar přenosu regulátoru
(
)(
)
k 0 + k1 z −1 + k 2 z −2 + k 3 z −3 + ... 1 − a1 z −1 − a 2 z −2 − a3 z −3 + ... U (z ) = . ∆Y (z ) 1 − a1* z −1 − a 2* z − 2 − a3* z −3 + ... T Přepis přenosu regulátoru (8) na diferenční rovnici je následující
(8)
u k − a1* u k −1 − a 2* u k − 2 − ... = k 0 ∆ y k + (k1 − a1 k 0 )∆ y k −1 + (k 2 − a1 k1 − a 0 k 2 )∆ y k − 2 + ...
(9)
Parametry na levé straně diferenční rovnice (9) jsou shodné se zkreslenými parametry pravé straně rovnice modelu (6). Mohou být proto použity k výpočtu pomocné veličinu uk* , (10)
u k* = u k − a1* u k −1 − a 2* u k − 2 − ... , pomocí které lze vytvořit lineární model závislosti této pomocné akční veličiny na zpožděných regulačních odchylkách ∆ y , tj.
(11)
u k* = k 0 ∆y k + k1* ∆y k −1 + k 2* ∆y k − 2 + ...
Neznámé parametry k 0 , k1* , k 2* ,... lze stejně jako u rovnice (6) určit z přepočítaných experimentálních dat metodou nejmenších čtverců. Součinitel u proměnné ∆ y k je roven parametru k 0 . Znalost velikosti tohoto parametru
k0
(a samozřejmě také
a 1* ) umožní
vypočítat nevychýlený parametr a1 ze vztahu a1* = a1 − b0 k 0 . Znalost velikosti a1 a hodnoty výrazu k1* = k1 − a1 k 0 z modelu (11) umožní vypočítat parametr k 1 . Tímto postupem lze postupně určit všechny neznámé parametry přenosů (2) a (3). Pořadí výpočtů je následující: 1)
a1 = a1* + b0 k 0
(12)
2)
k1 = k1* + a1 k 0
(13)
3)
a 2 = a 2* + b0 k1
(14)
4)
k 2 = k 2* + a1 k1 + a 2 k 0
(15)
5)
a3 = a3* + b0 k 2
(16)
6)
k 3 = k 3* + a1 k 2 + a 2 k1 + a3 k 0
(17)
obecně ai = ai* + b0 k i −1 i
k i = k + ∑ a j k i− j * i
(18) (19)
j =1
atd. Stejně jako u rovnice (6) vzniká problém volby řádu diferenční rovnice (11). Nechť je identifikován model řádu M, tj. poslední nenulový parametr je k M* . Podle vzorců (12) až (17) může teoreticky pokračovat výpočet k i , a i do nekonečna. Pro i > N je a i = b0 k i −1 a pro i
i > M je k i = ∑ a j k i − j , kde M a N jsou řády rovnic (6) a (11). j =1
Jestliže pro i > N je předpokládáno, že a i = 0 , pak počet rovnic (18) je jen N a pro stabilní přenos (2) parametry k i , i = 0, 1, 2,... konvergují k nule.
4. Příklad Příklad navazuje na korelační funkce v obr. 1. Řády modelů (6) a (11) byly zvoleny 8. Vychýlené a nevychýlené parametry modelů (6) a (1) jsou v tabulce 1. Na obr. 3 je impulsní odezva příslušná modelu (6), tj. včetně ruční zpětné vazby, a na obr. 5 impulsní odezva samotné identifikované soustavy. Obě impulsní charakteristiky neobsahují jeden krok dopravního zpoždění, který k uzavřené smyčce logicky patří. Obr. 4 obsahuje impulsní odezvu ruční zpětné vazby. Všechny impulsní odezvy jsou trvají teoreticky nekonečnou dobu. Tab. 1. Vlivnost normovaných nezávisle proměnných na závisle proměnnou a1* = 0,5915
a1 = 0,6230
a 2* = 0,0058
a 2 = −0,0350
a3* = −0,2093
a3 = −0,3395
a = 0,1619
a 4 = 0,0415
a = −0,0898
a 5 = −0,1438
a 6* = −0,0636
a 6 = −0,0748
a 7* = 0,1223
a7 = 0,1639
a = 0,0817
a8 = 0,1552
* 4 * 5
* 8
0,04
h i*
0,03 0,02 0,01 0,00
-0,01 0
5
10
15
20
počet kroků i
Obr. 3. Impulsní odezva modelu podle rovnice (6) 0,3
ki
0,2 0,1 0 -0,1 -0,2 0
5
10
15
počet kroků i
Obr. 4. Impulsní odezva ruční zpětné vazby
20
0,04
hi
0,03 0,02 0,01 0,00 -0,01 -0,02 0
50
100
150
200
počet kroků i
Obr. 5. Impulsní odezva modelu podle rovnice (1) Na obr. 6 je znázorněna poloha pólů přenosových funkcí (2) a (4) v komplexní rovině, ve které je vyznačena oblast stability pro póly diskrétních přenosových funkcí. Zajímavý je přesun pólů s nejnižší frekvencí kmitů k hranici stability (označeno kroužky). Kmity regulované veličiny s periodou asi 9 vzorkovacích intervalů jsou zjevné v autokorelační funkci Si x Si na obr. 1. Při návrhu řízení bude třeba tuto skutečnost respektovat. Řešením je suboptimální regulace [Tůma, 1998]. &'()*+,-.'*/,01+)*2,&3045 "
"
#
%$#
%$#
#
!#
%$!#
%$!" !"
!#
%$#
#
%$"
!" !"
4(0,01+)*2,&3045
!#
%$#
#
%$"
Obr. 6. Poloha pólů přenosových funkcí (2) a (4) v komplexní rovině
5. Závěr V referátu je popsán postup identifikace soustavy řízené ruční zpětnou vazbou. Kromě identifikace vlastní řízené soustavy je zjištěna impulsní odezva ruční zpětné vazby. Blízkost kořenů jmenovatele přenosu řízené soustavy je třeba zohlednit při návrhu regulátoru.
6. Použita literatura PETERKA, V. 1975. Číslicové řízení procesů s náhodnými poruchami a neurčitými charakteristikami. Doktorská disertační práce, Praha 1975. 305 s. TŮMA, J. 1999. Složité systémy řízení. 2. vyd. Ostrava: VŠB-TU Ostrava, 1998. 151 s.