1. Bevezetés a Matlab Használatába

Bevezetés a Matlab Használatába

Hetthéssy JenĘ, BarsRuth, Barta András, 2004



» A=[7, 8, 9; 5, 6, 7] % A egy 2x3-as mátrix, » MATRIX=[row1; row2; . . . ; rowN];

1. Bevezetés a Matlab Használatába A Matlab egy interaktív programcsomag tudományos és mérnöki számítások, szimulációk és grafikus adatmegjelenítés elvégzésére. A Matlab erĘs hátteret biztosít mátrix algebra, differenciálegyenletek és más matematikai problémák megoldására. A Matlab elterjedtsége és széleskörĦ alkalmazása abból adódik, hogy a Matlab utasításkészlete kiterjeszthetĘ toolbox-ok segítségével A toolbox tulajdonképpen egy függvénykönyvtár, amit különféle szakterületek támogatására fejlesztettek ki. Ilyen szakterületek például a jelfeldolgozás (Signal Processing Toolbox), a szabályozástechnika (Control System Toolbox), a képfeldolgozás (Image Processing Toolbox), a neurális hálózatok alkalmazása (Neural Network Toolbox), stb. Ennek a bevezetésnek a célja, hogy minél gyorsabban el lehessen jutni egy olyan szintre, ami szükséges az alapproblémák megoldására. Részletesebb információt a Matlab felhasználói segédletek tartalmaznak. Elektronikus formában ezek a ‘matlab/help’ könyvtárban találhatók meg. Egy on-line súgó is a felhasználó rendelkezésére áll. » helpdesk A help utasítás nyújt segitséget az egyes utasítások használatárol, például: » help sqrt A Matlab utasításokból egy szkript fájlt is létre lehet hozni. Ez egy text fájl .m kiterjesztéssel. Ez a szkript új utasításként használható (kiterjesztés nélkül).

Változó nevek: A változók nevének mérete maximálisan 31 karakter lehet (betĦk, számok és a ’_’ karakter). Az elsĘ karakter nem lehet szám. A kis és nagybetĦk nem azonosak. A casesen utasítás megváltoztatja ezt a mĦködési módot (azaz azonossá teheti a kis és nagybetĦket, de ennek használata nem ajánlott). Minden vátozó egy mátrix. Egy skaláris változó egy egyszer egyes mátrixnak tekinthetĘ. Adatbevitel: A Matlab számos adattípust használ, de ezeket nem kell megadni, automatikusan deklarálja Ęket. Integer: » k=2 Ha az utasítás pontosvesszĘvel végzĘdik, akkor az eredmény nem jelenik meg. » J=-4; Valós: » s=3.6 » F2=-12.6e-5 Komplex: » z=3+4*i » r=5*exp(i*pi/3) Az i=sqrt(-1) változó, az egység képzetes vektor elĘre meghatározott. Ehelyett azonban más változó is használható, például: » j=sqrt(-1) Vektor: » x=[1, 2, 3] % sorvektor, az elemeket vesszĘvel vagy üres karakterrel (space) kell elválasztani » q=[4; 5; 6] % oszlopvektor, az elemeket pontosvesszĘvel kell elválasztani A sorvektorból a transzponálás mĦveletével is képezhetĘ oszlopvektor: » v=[4, 5, 6]’ % ugyanaz mint q Megjegyzés: A transzponálás komplex változók esetén a transzponált vektor komplex konjugáltját képezi. » i’ 0 - 1.000i Mátrix:

1

Speciális vectorok és mátrixok: » u=1:3; » w=1:2:10; » E=eye(4)

% létrehozza a v=[1 2 3] vektort, » u=start:stop % létrehozza a w=[1 3 5 7 9] vektort, » =start:increment:stop

E= 1 0 0 0

0 1 0 0

0 0 1 0

0 0 0 1

1 0 0

0 1 0

0 0 1

0 0 0

0 0

0 0

0 0

0 0

1 1 1

1 1 1

1 1 1

1 1 1

» B=eye(3,4) B=

» C=zeros(2,4) C= » D=ones(3,5) D= 1 1 1

Változó értékek: A változó nevének beírása kijelzi a változó értékét. » A A= 7 5

8 6

9 7

» A(2, 3) ans = 7 Az elsĘ index a sor számát jelzi, a második pedig az oszlop sorszámát adja meg. A választ az ans változó veszi fel. Egy mátrix vagy vektor egy elemének megváltoztatása után a változó összes eleme látható a képernyĘn, kivéve ha az utasítás pontosvesszĘvel végzĘdik. » v(2)= -6 v= 4 -6 6 Indexelés (Subscripting): A kettĘspont (:) segítségével a mátrix több eleme is hozzáférhetĘ egyszerre. Többféleképpen is használható ez az operátor. Kezdeti index : végsĘ index - egy mátrix egy részéhez enged hozzáférni : - a kettĘspont az összes elemet jelöli ki egy adott sorban vagy oszlopban. Vektorok esetén: v=[v(1) v(2) . . . v(N)] Mátrixok esetén: M=[M(1,1)...M(1,m); M(2,1)...M(2,m); ... ; M(n,1)...M(n,m)] Ha B egy 8x8 mátrix, akkor

2

Bevezetés a Matlab Használatába B(1:5,3) B(2:3,4:5) B(:,3) B(2,:) B(1:3,:) » A(2,1:2) » A(:,2)


egy oszlopvektor, [B(1,3); B(2,3); B(3,3); B(4,3); B(5,3)] egy mátrix [B(2,4) B(2,5); B(3,4) B(3,5)] a B mátrix mindegyik sorát és harmadik oszlopát, azaz a harmadik oszlopot jelenti a B mátrix második sorát jelenti a B mátrix elsĘ három sorát jelenti. % az A mátrix második sorának elsĘ és második elemét jelöli ki % az A mátrix második oszlopát jelöli ki

Munkaterület (Workspace): A használt változók a workspace-nek nevezett memória területen helyezkednek el. A workspace a következĘ utasításokkal jeleníthetĘ meg: » who » whos % kiírja a változó méreteit is



» 2*x ans= -2 0 4 Mátrix szorzása skalárral: » 3*A ans= 3 9

6 12

Skaláris szorzat: » s=x’*y s= 4 » y’*x

A változók méreteit a length és a size utasítással lehet megadni. Vektorok esetén: » lng=length(v) lng= 3 Vektorok és mátrixok esetén: » [m,n]=size(A) m= 4 % a sorok száma n= 4 % az oszlopok száma A workspace-t el lehet menteni, betölteni és kitörölni: » save % elmenti a workspace-t a default matlab.mat fájlba. » save filename.mat % elmenti a workspace-t a filename.mat fájlba. » clear % kitörli az összes változót. » load % betölti a workspace-be a matlab.mat fájl tartalmát » load filename.mat % betölti a workspace-be a filename.mat fájl tartalmát

ans= 4 Diadikus szorzat( mátrixot eredményez): » M=x*y’ M= 2 1 0 0 -4 -2 » y*x’ ans= 2 0 1 0 -1 0

-4 0 8

4];

Osztás: » C/2 Mátrixok esetén: B/A A\B

C= 2 5

5 8

megfelel BA-1-nek megfelel A-1B-nek

Hatványozás:

» D=A-B D= » x=[-1 » y=x-1

-4 -2 2

Mátrix szorzása vektorral: » b=M*x b=

Aritmetikai operátorok: Összeadás és kivonás: » A=[1 2; 3 » B=A’; » C=A+B;

-1 0 2

0

0 -1 1 0 2]’; (Vegyük észre, hogy a vektor mindegyik eleme megváltozott!) y= -2 -1 1

Szorzás: Vektor szorzása skalárral:

3

A^p , ahol A egy négyzetes mátrix és p egy valós konstans például A inverze : A^(-1) , vagy ezzel egyenértékĦ inv(A) MĦveletek komplex számokon: » r r= 2.5000 + 4.3301i » real(r) ans= 2.5000 » imag(r) ans= 4.3301

4



» abs(r)

rat expm logm sqrtm

ans= 5 » angle(r) ans=


rational approximation matrix exponential matrix logarithm matrix square root.

Például: » help sqrt » g=sqrt(2)

1.0472 A szög fokokban » 180/pi*angle(r) ans= 60.0000

Grafikus kimenet:

Tömbök elemein végzett mĦveletek: A tömbök elemein végzett aritmetikai mĦveletek a mĦveletek fontos csoportját alkotják. Ha az aritmetikai mĦveletet megelĘzi egy pont, akkor ez a mĦvelet a változók elemeire vonatkozik. A .* mĦvelet összeszorozza a változók megfelelĘ indexĦ elemeit. Például sorvektorok esetén: a.*b=[a(1)*b(1), a(2)*b(2), ..., a(n)*b(n)]. A változók méretének meg kell egyeznie. Például: » a=[2 4 6] » b=[5 3 1] » a.*b ans= 10 12 6 Az elemenként végzendĘ mĦveletek értelemszerĦen más operátorokra is alkalmazhatók, például osztás és hatványozás esetén: ./, .^

AlapvetĘ matematikai függvények: (Az egyes függvények használatát a help utasítás mutatja meg) abs sqrt real imag conj round fix floor ceil sign rem sin cos tan asin acos atan atan2 sinh cosh tanh exp log log10 bessel


absolute value or magnitude of complex numbers square root real part imaginary part comlplex conjugate round to nearest integer round towards zero round towards -infinity round towards +infinity signum function remainder sine cosine tangent arcsine arccosine arctangent four quadrant arctangent hyperbolic sine hyperbolic cosine hyperbolic tangent exponential base e natural logarithm log base 10 Bessel function

5

A legegyszerĦbb grafikus utasítás a plot. » plot(2,3) % kirajzolja az x=2, y=3 pontot Egyszerre több pontot is ki lehet rajzolni, ha a koordináta értékeket egy vektorban tároljuk el. » x=[1,2,3] » y=[0,2,1] » plot(x,y) % a pontok össze vannak kötve egyenesekkel » plot(x,y,’*’) % így az utasítás csak a pontokat rajzolja ki Ezzel a módszerrel bonyolultabb görbék is kirajzolhatók. Például: » t=0:0.05:4*pi; » y=sin(t); » plot(t,y) » title(‘Sine function’),xlabel(‘Time’),ylabel(‘sin(t)’),grid; ahol a title , xlabel , ylabel és grid utasítások opcionálisak. Egyszerre több görbe is kirajzolható: » y1=3*sin(2*t); » plot(t,y,’r’,t,y1,’b’); % r-red (piros), b-blue (kék) A kirajzolás színét és típusát a következĘ formátumban lehet megadni: » plot(t,y,’@#’), ahol ‘@’ jelenti a vonal típusát: - folytonos -- szaggatott : pontvonal . pontok + plusz * csillag o kör x x-jel és ‘#’ jelenti a színt:

r g b w y

red, piros green, zöld blue, kék white, fehér yellow, sárga

Polinomok: A Matlab-ban a polinomokat együtthatóikkal adhatjuk meg. Például a P(x)=4x4-6x3+9x2-5 polinomot megadhatjuk egy vektorral, aminek az elemei a polinom együtthatói. » v=[4 -6 9 0 -5] A polinom gyökeit, azaz a P(x)=0 egyenlet megoldásait a roots utasítás számolja ki. » p=roots(v) p = 0.6198 + 1.4337i 6



0.6198 - 1.4337i 0.8577 -0.5974 A polinom p1, p2, ... , pn gyökeibĘl kiszámolhatjuk a P(x)=(x- p1)(x- p2)...(x- pn) polinom együtthatóit. » v1=poly(p) v1 = 1.0000 -1.5000 2.2500 0 -1.2500 A poly utasítás a legnagyobb hatványú együtthatót egynek veszi (normálja). Ha vissza akarjuk kapni az eredeti polinomot, akkor meg kell szorozni az összes együtthatót az elsĘ együtthatóval, a mi esetünkben 4-gyel. » v2=4*poly(p) v2 = 4.0000 -6.0000 9.0000 0 -5.0000 Vegyük fel a következĘ mátrixot: » M=[3 5 ; 7 -1] Az M mátrix sajátértékeit az eig utasítás határozza meg. » e=eig(M) e = 7.2450 -5.2450 A sajátértékeket egy polinom gyökeinek tekintve a poly utasítással egy polinom meghatározható: » poly(e) ans = 1 -2 -38 A fenti polinom, azaz x2-2x-38 az M mátrix karakterisztikus polinomja, ami a det(xI-M) egyenlettel határozható meg. A karakterisztikus polinom az M mátrixból közvetlenül is kiszámítható: » poly(M) ans = 1 -2 -38 Látható, hogy sokszor egy utasítás többféleképpen is meghívható. A Matlab help egy utasítás valamennyi meghívási lehetĘségét megadja.

7

Bevezetés a Matlab Control System Toolbox Használatába

Hetthéssy JenĘ, Bars Ruth, Barta András, 2004

2. Bevezetés a Matlab Control System Toolbox Használatába Tekintsünk egy egy-bemenetĦ, egy-kimenetĦ folytonos lineáris rendszert, amely átviteli függvényével adott:

u(t) U(s)

H ( s)

Állítsuk elĘ a Matlab segítségével a H ( s)

Y (S )

num

U ( s)

den

y(t) Y(s)

2 rendszer átmeneti függvényét, azaz az y(t) s 2 2s 4

kimenetet egységugrás bemenet esetén. Az átviteli függvényt a számláló és a nevezĘ polinom felvételével határozzuk meg. Egy polinomot az együtthatóival adhatunk meg, így a fenti példában num 2 , den s 2 2s 4 . » num=2 Step response 0.7 » den=[1 2 4] Általános esetben az együtthatókat csökkenĘ fokszám szerinti 0.6 sorrendben kell egy vektorba gyĦjteni. Az átmeneti függvényt a 0.5 Matlab step utasításával lehet megjeleníteni: 0.4 » step(num,den); 0.3 Ugyanezt közvetlenül is meg lehet adni: » step(2,[1 2 4]); 0.2 Vegyük észre, hogy az idĘskálát a Matlab automatikusan 0.1 választotta. 0 Amennyiben tárolni kívánjuk a kiszámított átmeneti függvény 0 1 2 3 4 5 6 értékeit, akkor a step függvényt bal oldali argumentummal annak a változó nevének megadásával, amelyben tárolni kívánjuk a kiszámított függvényértékeket - kell meghívnunk: » [y,x,t]=step(num,den) vagy egyszerĦbben » y=step(num,den) Látható és érdemes megjegyezni, hogy a Matlab függvényhívásoknak jellemzĘje, hogy az argumentumok száma változó lehet mindkét oldalon. A fenti három bal oldali változót megjelölĘ függvényhívásban az y a kimenĘjel értékeit, az x a rendszer belsĘ állapotváltozóinak értékeit és a t a számítások idĘpontjainak értékeit tartalmazza. Vegyük észre, hogy ebben az esetben grafikus ábrázolás nem történt. Ha pontosvesszĘt teszünk a sor végére, akkor az adatokat sem írja ki a Matlab. Az értékek kiírhatók a változó nevének beírásával. » y Látható, hogy egy oszlopvektort kaptunk, amelynek elemei a mintavételezési pontokban felvett értékek. A mintavétel természetes, hiszen nem analóg számítógéppel dolgozunk, így csakis diszkrét idĘpontokban várhatjuk a szimuláció eredményét. Kérdés, hogy milyen idĘfelbontással dolgoztunk. Az ábrából látható, hogy a Matlab automatikusan a 0 d t d 6 idĘtartományt választotta ki az átmeneti függvény megjelenítésére. A Matlab a rendszer dinamikai tulajdonságai (zérusai, pólusai) alapján választja meg az idĘvektort. Az y vektor elemeinek számából így már vissza tudunk következtetni a mintavételi idĘre. » n=length(y) n = 109 » T=6/n T = 0.055 A számított mintavételi idĘ: T=6/109=0.055 sec. A help utasítás megmutatja az utasítás részletes használati módjait: » help step



Természetesen számtalan módon variálható a fenti szimulációs futtatás. Az utasítás egy fontos alakja, amikor nem kívánjuk a Matlab automatikus szimulációs idĘközét és idĘtartamát választani, hanem ezt rögzíteni kívánjuk a szimuláció végrehajtása elĘtt, pl. 10 sec-ig szeretnénk látni az átmeneti függvényt, mégpedig 0.1 sec felbontási idĘvel. Adjuk meg elĘször a 0-tól 10-ig terjedĘ intervallumot 0.1-es lépésközökkel »t=0:0.1:10 majd jelenítsük meg az átmeneti függvényt: » y=step(num,den,t) A kimenĘjelet grafikusan megjeleníthetjük a plot utasítással: » plot(t,y); A könnyebb értelmezhetĘség érdekében egy koordináta hálót is rajzolhatunk az ábrára a grid utasítással. » plot(t,y),grid; A plot utasítás a pontok között lineáris interpolációt alkalmaz, azaz összeköti a pontokat egy egyenessel. Az egyes pontok összekötése csak a vizuális megjelenítést (végül is folytonos rendszer folytonos kimenetét szeretnénk látni) szolgálja, valójában nincs pontosan kiszámított információnk a kimenĘjel két szimulációs pont közötti viselkedését illetĘen. Ezt az interpolációt elhagyhatjuk, ha pontokkal ábrázoljuk a mintavételi értékeket. » plot(t,y,'.'); Ezeket az y értékeket további számításokra is fel lehet használni. Az átmeneti függvény maximális értékét megkaphatjuk a max függvény hívással. » ym=max(y) ym = 0.5815 Stacionárius értékét a dcgain utasítás adja meg. » ys=dcgain(num,den) ys = 0.5 A százalékos túllüvést ezekbĘl ki tudjuk számolni. » yovrsht=(ym-ys)/ys*100 yovrsht = 16.2971


» [r,p,k]=residue(num,den) r = 1.0000 -4.0000 2.0000 p = -3.0000 -3.0000 -2.0000 k = []

A eredmény a Laplace tartományban (komplex frekvencia tartományban): Y (s)

r (1) r (2) r (3) k s p (1) ( s p (2)) 2 s p (3) y (t )

e 3t 4te3t 2e2t , t t 0

Figyeljük meg a kettĘs pólusnak megfelelĘ struktúrát mind az r és p vektorban, mind a részlettörtes Laplace transzformált alakban, mind pedig az idĘtartománybeli kifejezésben. Az analitikus kifejezés alapján természetesen elĘ tudjuk állítani az idĘfüggvényt: » t=0:0.05:6; » y=r(1)*exp(p(1)*t)+r(2)*t.*exp(p(2)*t)+r(3)*exp(p(3)*t); (A második kifejezésben a t melletti pont azt jelöli, hogy a mĦveletet nem a vektorra, hanem annak elemeire értelmezzük). Ugyanezt egyszerĦbben is meg lehet határozni numerikusan. » yi=impulse(num,den,t); » plot(t,[y,yi]),grid;

LTI modell struktúra (sys): Az egyszerĦbb használat érdekében a Control System Toolbox bevezeti az LTI adatstruktúrát (Linear Time Invariant system). Egy LTI struktúrában egy lineáris rendszert háromféle alakban lehet megadni.

r s p r

s m bm 1 s m 1 .... b2 s 2 b1 s b0 an s n an 1 s n 1 .... a2 s 2 a1 s a0 ( s z1 )( s z2 )...( s zm ) ( s p1 )( s p2 )...( s pn )

2 s 2 3s 2

2 ( s 1)( s 2)

1

L o re pt L

1

x Állapotteres alak: pt

o rte (s p)2 Ezt az alakot részlettörtekre bontással lehet meghatározni, amit a residue utasítással lehet végrehajtani. Legyen a feladat az Y(s)=

x Átviteli függvény alak: H tf ( s )

x Zérus-pólus-erĘsítés alak: H zpk ( s) k

1

L k o k 1(t )

1 4 2 s 3 ( s 3) 2 s 2

és az idĘtartományban:

Inverz Laplace Transzformáció: Analitikus vizsgálatokhoz szükséges az inverz Laplace transzformáció számítása. A módszer lényege, hogy a függvény Laplace transzformáltját olyan tagok összegére bontjuk fel, amelyeknek inverz transzformáltját már ismerjük. A leggyakrabban elĘforduló tagok:


3s 2 +13s+16 (s+2)(s+3) 2

valós törttel a Laplace tartományban leírt jel visszatranszformálása az idĘtartományba, azaz y(t) analitikus meghatározása. ElĘször adjuk meg a jel Laplace transzformáltját számláló és nevezĘ polinomjával: » num=[3 13 16]; » den=poly([-3 -3 -2]); A részlettörtekre való bontás:

x = Ax + Bu y = Cx + Du

ª 3 1º A=« »,B ¬2 0¼

ª1 º «0 » , C ¬ ¼

>0 1@ , D

0. .

A Matlab-ban a tf, zpk és ss utasításokkal lehet megadni egy lineáris rendszert LTI sys struktúrában. Legyen a rendszer átviteli függvénye H(s). » num=2 » den=[1, 3, 2] » H=tf(num,den) Transfer function: 2 ------------s^2 + 3 s + 2 vagy közvetlenül: » Htf=tf(2,[1, 3, 2]) A többi modell hasonlóan adható meg:

» Hzpk=zpk([],[-1, -2],2) Zero/pole/gain: 2 ----------(s+1)(s+2) » A=[-3, -1; 2, 0]; B=[1; 0]; C=[0, 1]; D=0; » Hss=ss(A,B,C,D); A modellek átkonvertálhatók egymásba. » Hzpk1=zpk(Htf) » Hss1=ss(Htf) » Htf1=tf(Hss) A modellek számos paramétert (properties) foglalnak magukban. Ezeket a get utasítással nézhetjük meg. » get(Htf) » get(Hzpk) » get(Hss) Az LTI modell paramétereket megkaphatjuk a tfdata, zpkdata, ssdata utasításokkal. A 'v' flag szükséges, hogy vektor alakban kapjuk meg a paramétereket. » [num1,den1]=tfdata(Htf,'v') num1 = 0 0 2 den1 = 1 3 2 » [z,p,k]=zpkdata(Hzpk,'v') A modell paraméterek közvetlenül is elérhetĘk a nevükön keresztül: » num2=Htf.num{1} num2 = 0 0 2 Az LTI adatstruktúra néhány paramétert ún. cell array típusban tárol, hogy több bemenetĦ, több kimenetĦ (MIMO, multi- input multi-output) rendszerek esetén is használható legyen. A cell array egy olyan mátrix, melynek elemei mátrixok. Egy cell array például a következĘképpen adható meg. » ca={1, [1,2],[1,2,3]} ca = [1] [1x2 double] [1x3 double] » ca{2} ans = 1 2 A {1} jelöléssel hivatkozunk a cell array elsĘ elemére. Egy bemenetĦ, egy kimenetĦ (SISO single-input single- output) rendszerre mindig {1} jelölést kell használni, mivel csak az elsĘ elem létezik. A num, den, z, p paraméterek cell array formátumban adottak. A get(…) utasítással lehet megnézni a paraméterek típusait. A többi modell paraméter hasonló módon érhetĘ el: » Hzpk.k=2; Hzpk.z=[]; Hzpk.p=[-1, -2]; » Hzpk k: 2 z: [] p: [-1 -2] » Hss.a ans = -3 -1 2 0 . Szimbólikus adatbevitellel tovább egyszerĦsíthetĘ az adat megadás. Az s a Laplace transzformációs változó. » s=tf('s') » H=1/((s+1)*(s+2))



Aritmetikai mĦveletek is alkalmazhatók az LTI struktúrákra. A leggyakrabban használatos mĦveletek: +, - ,*, / , \ , ' , inv, ^ . Például egy zárt rendszer eredĘ átviteli függvénye kiszámítható közvetlenül az alábbi szimbólikus összefüggéssel: » Hcl=H/(1+H) Az LTI sys struktúrák egy hierarchikus sorrendben helyezkednek el: tf -> zpk -> ss. Ha egy utasításban vagy számításban különbözĘ modellek szerepelnek, akkor az eredmény mindig a magasabb hierarchiájú modell formájában keletkezik. Például a » Htf*Hzpk utasítás eredménye zpk alakban jelenik meg.

IdĘtartománybeli vizsgálat: A Control System Toolbox több utasítást tartalmaz, amelyek segítségével a lineáris rendszerek idĘtartománybeli viselkedését lehet vizsgálni. Vizsgáljuk ismét a következĘ rendszert: H ( s)

2 s 2 2s 4

» H=2/(s^2+2*s+4) Vegyünk fel egy idĘtartományt: » t=0:0.1:10; x Átmeneti függvény (step response): Az átmeneti függvény a rendszer kimenĘjele egységugrás (step) bemenĘjel esetén. A step utasítás elĘzĘekben tárgyalt alakjai mellett az alábbi alak is hívható: » step(H); x Súlyfüggvény (impulse response): A súlyfüggvény a rendszer kimenete Dirac-delta bemenet esetén. » impulse(H); » yi=impulse(H,t); » plot(t,yi) x Kezdeti feltétel (initial condition): A rendszer viselkedése kezdeti feltétel esetén az initial utasítás segítségével vizsgálható. Csak állapotteres reprezentációra alkalmazható ez az utasítás. Bode Diagrams » H=ss(H) From: U(1) » x0=[1, -2] 0 » [y,t,x]=initial(H,x0); -10 » plot(t,y),grid -20 Vegyük észre, hogy az x változó egy mátrix, amelynek annyi -30 -40 oszlopa van, ahány állapotváltozója van a rendszernek 0 (ebben az esetben 2). Sorainak száma megegyezik az -50 idĘvektor elemeinek számával. EllenĘrizzük: -100 » size(x) -150 ans = 109 2 -200 10 10 10 A rendszer állapottrajektóriáját is tudjuk vizsgálni. Az x mátrix elsĘ oszlopa tartalmazza az elsĘ állapotváltozót, a Frequency (rad/sec) második pedig a másodikat. A : kiválasztja az adott vektor minden egyes elemét. Az állapottrajektória az egyik állapotváltozót a másik állapotváltozó függvényében ábrázolja. » x1=x(:,1); x2=x(:,2); » plot(x1,x2) Phase (deg); Magnitude (dB)


To: Y(1)


-1

x

0

TetszĘleges bemenet: A rendszer kimenetét tetszĘleges bemenet esetén is ki lehet számítani. Határozzuk meg a kimenetet u(t)= 2sin(3t) bemenet esetén. » usin=2*sin(3*t); » ysin=lsim(H,usin,t); Ábrázoljuk a bemenĘjelet piros, a kimenĘjelet pedig kék görbével. » plot(t,usin,'r',t,ysin,'b'), grid;

1



Frekvenciatartománybeli vizsgálat: A rendszer viselkedése a frekvenciatartományban is vizsgálható. x Bode diagram A rendszer Bode diagramját a bode utasítással számíthatjuk. Számos módon alkalmazhatjuk ezt az utasítást. Egy megadott frekvencián kiszámolhatjuk a rendszer erĘsítését és fázisszögét. Számítsuk ki az erĘsítést az Z =5 frekvencián. » w=5; » [gain,phase]=bode(H,w); Az eredmény: gain = 0.0860, phase = -154.5367 A számítás megismételhetĘ különbözĘ frekvenciákon. A bode utasítás azonban meghívható úgy is, hogy egyidejĦleg több frekvencián is kiszámítja a frekvenciafüggvény abszolút értékét és fázisszögét. A Bode diagram közvetlenül grafikusan is megjeleníthetĘ. Ebben az esetben a Matlab automatikusan felvesz egy frekvenciatartományt a rendszer dinamikus tulajdonságai alapján (a step utasításhoz hasonlóan). » bode(H),grid A frekvencia pontok logaritmikus skálán helyezkednek el, mivel ez szemléletesebb, nagyobb frekvenciatartományt átfogó megjelenítést biztosít. A frekvencia vektor közvetlenül is meghatározható: » w=logspace(-1,1,200); A logspace utasítás létrehoz egy 200 pontból álló logaritmikusan egyenlĘ távolságra elhelyezkedĘ pontokból álló frekvencia vektort a 0.01 and 100 közti frekvenciatartományban. Ezzel a Bode diagram értékei már kiszámíthatók és ábrázolhatók. » [gain,phase]=bode(H,w); x

Nyquist diagram : A rendszer erĘsítését és fázistolását ábrázoljuk a komplex síkon úgy, hogy a frekvencia változik. » nyquist(H); x Margin : ez az utasítás kiértékeli a frekvenciafüggvény jellemzĘit, kiszámítja a stabilitásra jellemzĘ mutatókat. » margin(H); Nyquist Diagrams x Zérusok, pólusok: A karakterisztikus polinom (az átviteli függvény nevezĘje) gyökei a rendszer pólusai. » [num,den]=tfdata(H,'v'); » poles=roots(den); A zérusok az átviteli függvény számlálójának a gyökei. » zeros=roots(num); A zérusok és pólusok megkaphatók közvetlenül a zpk modellbĘl: Real Axis » [z,p,k]=zpkdata(H,'v'); A zérusok és pólusok ábrázolhatók a komplex síkon. » subplot(111); » pzmap(H); A damp utasítás megadja az összes pólust és komplex póluspár esetén a csillapítási tényezĘt és a sajátfrekvenciát is: » damp(H); A rendszer DC (egyenáramú) erĘsítése meghatározható: » K=dcgain(H); From: U(1)

0.6

0.4



Egy lineáris rendszer egyszerĦen vizsgálható az LTI Viewer segítségével. Az LTI Viewer egy grafikus felület egy rendszer idĘ- és frekvenciatartománybeli tulajdonságainak vizsgálatára. Több rendszer egyszerre is vizsgálható a különbözĘ menĦpontok és a jobb egér billentyĦ segítségével. » ltiview vagy konkrétabban, » ltiview('bode',H);

Simulink: A Simulink egy grafikus programcsomag, amely blokk diagramon alapuló szimuláció végrehajtására képes. Egy Simulink szimuláció két lépésbĘl áll. ElĘször létrehozzuk a modellt, majd futtatásokat végzünk a kívánt paraméterbeállítások mellett. A Simulink-ben egy dinamikus rendszert blokk diagrammal lehet megadni. A rendszer definiálása a blokk diagram megrajzolását jelenti. A blokkokat nem kell külön elĘállítani, egy könyvtárból ki lehet Ęket másolni. A könyvtár blokkjai mĦködésük szerint külön alkönyvtárakban helyezkednek el. A Simulink blokk könyvtár megnyitható a Matlabból a simulink utasítással, vagy a Simulink ikonnal. Egy új modellt létre lehet hozni –File -> New menĦ aktíválásával. Ez létrehoz egy új ablakot, amelybe a blokkok már bemásolhatóak a könyvtárból. Az egér bal gombjának lenyomva tartásával lehet a könyvtári blokkokat áthúzni a fájlunkba. A blokkokat össze tudjuk kötni az egér bal gombjával, ráklikkelünk a blokkok bemenetére vagy kimenetére, majd a vonalat odahúzzuk a kívánt helyre. Egy blokk duplikálható a jobb egér billentyĦ klikkelésével. A blokkokat méretezni és elforgatni is lehet. A blokk paramétereit dupla klikkeléssel lehet megváltoztatni. Az új fájl elmenthetĘ a File->Save menĦponttal. A szimuláció elindítható a Start -> Simulation menĦvel, vagy a Run ikonnal (Ź). A szimuláció beállítási paraméterei megváltoztathatók.A futtatási eredményeket a Scope blokkal lehet megvizsgálni vagy a To Workspace blokkal visszaküldeni a Matlab felületre. A Simulink a modellek szimulációját közönséges differenciálegyenletek integrálásával végzi. A Simulink számos integrálási módszer használatát teszi lehetĘvé. Pontos eredmények elérésére a megfelelĘ módszer kiválasztása is egy fontos szempont lehet. A Simulink mĦködésének megismerése céljából vegyünk egy egyszerĦ példát. Hozzunk létre egy új fájlt és másoljuk be a különféle blokkokat. A blokkok paramétereit változtassuk meg a kívánt értékre. Állítsuk be a menĦbĘl a Simulation–>Parameters–>Stop time paramétert 50-re. A Simulink a Matlab-ban definiált változókat használja.

To: Y(1)

Imaginary Axis

0.2

0

t

-0.2

-0.4

Clock

y

time

output

-0.6

-0.8 -1

LTI Viewer:

-0.5

0

0.5

1.5 Step Input

Gain

H LTI SystemPs

Transport Delay

Scope1

H(s): Control System Toolbox –>LTI system : H KülönbségképzĘ: Simulink–>Math–>Sum: +Késleltetés, holtidĘ: Simulink–>Continuous–>Transport Delay: 1 ErĘsítés: Simulink–>Math–>Gain: 1.5 Step input: Simulink–>Sources–>Step Scope: Simulink–>Sinks–>Scope Clock: Simulink–>Sources–>ClockOutput, time: Simulink–>Sinks–>To Workspace: y,t



A frekvenciafüggvény


Változtassuk a Gain paraméter értékét 0.5 és 2 között. Állapítsuk meg, a paraméter milyen értékénél lépnek fel a rendszerben állandósult lengések.

3. A frekvenciafüggvény Stabilis lineáris rendszerek egy alapvetĘ tulajdonsága, hogy szinuszos bemenĘjelekre állandósult állapotban, a tranziensek lecsengése után a bemenĘjel frekvenciájával megegyezĘ frekvenciájú szinuszos jelekkel válaszolnak. Legyen

Au sin(Zt M u )

u( t )

t t 0 esetén a rendszerünk bemenete. A kimenĘjel

y (t )

yállandósult (t ) ytranziens (t )

A kimenĘjel állandósult (kvázistacionárius) állapotban

Ay sin(Zt M y )

y állandósult (t )

u (t )

Au sin(Zt Mu )

y (t )

H ( s)

Ay sin(Z t M y ) ytranziens

A frekvenciafüggvény az Ay/Au amplitúdó arány és a (My-Mu) fáziskülönbség frekvenciafüggését leíró függvény. Mivel a frekvenciafüggvény szimultán módon két rendszerjellemzĘ tulajdonság frekvenciafüggését kívánja reprezentálni, nem meglepõ, hogy a frekvenciafüggvény komplex függvény. Bebizonyítható, hogy formailag a frekvenciafüggvényt az átviteli függvénybĘl az s=jZ helyettesítéssel lehet származtatni.

H ( jZ )

H ( s) s

jZ

M (Z )e jM (Z )

A frekvenciafüggvény kifejezésében M(Z) az amplitúdófüggvény (a frekvenciafüggvény abszolút értéke), M(Z) pedig a fázisfüggvény:

M (Z )

H ( jZ )

Ay (Z )

M (Z ) arg^H ( jZ )` M y (Z ) M u (Z )

Au (Z )

A frekvenciafüggvény többféleképpen is megjeleníthetĘ. Az M(Z) és M(Z) függvényeket külön-külön ábrázolhatjuk egy kijelölt frekvenciatartományban. Ez a technika a Bode diagramot eredményezi. A frekvenciaskála léptéke rendszerint logaritmikus. A Nyquist diagram a frekvenciafüggvényt a komplex számsíkon ábrázolja. A kiválasztott frekvenciatartomány minden egyes Z értékére a komplex síkban az M(Z) és M(Z) értékpárnak megfelelĘ pontot adhatunk meg. E pontok kontúrral való összekötése eredményezi a yNquist diagramot. A N yquist diag ram ábrázolásakor a frekvenciát rendszerint nulla és végtelen között változtatjuk. Példa: Legyen egy rendszer átviteli függvénye H ( s)

10 s 2 2 s 10

Határozzuk meg a rendszer kimenĘjelét, ha bemenĘjele u(t)= 2*sin(3*t). » num=10; » den=[1, 2, 10]; » H=tf(num,den); » t=0:0.05:10; » u=2*sin(3*t); » yl=lsim(H,u,t); Ábrázoljuk a rendszer bemenĘjelét és kimenĘjelét egy diagramban (bemenĘjel: piros (red), kimenĘjel: kék (blue)):

1



» plot(t,u,'r',t,yl,'b'), grid; Állandósult állapotban, a tranziens lecsengése után a kimenĘjel szinuszos, frekvenciája megegyezik a bemenĘjel frekvenciájával, amplitudója és fázisszöge függ a bemenĘjel frekvenciájától. Az amplitudó és a fázisszög kiszámítható a H ( s jZ ) frekvenciafüggvénybĘl. A Matlab-ban a bode utasítás alkalmazható ezen értékek kiszámítására egy adott frekvencián vagy egy adott frekvenciatartományban. Például az Z 3 frekvencián: » [gain,phase]=bode(H,3); Vizsgáljuk a rendszer válaszát az alábbi bemenĘjelekre: u(t)= sin( Z t), Z1 1, Z 2 2, Z3 5, Z4 10 . » w1=1; w2=2; w3=5; w4=10; » u1=sin(w1*t); u2=sin(w2*t); u3=sin(w3*t); u4=sin(w4*t); » y1=lsim(H,u1,t); y2=lsim(H,u2,t); y3=lsim(H,u3,t); » y4=lsim(H,u4,t); » plot(t,y1,'r',t,y2,'g',t,y3,'b',t,y4,'m'), grid; Számítsuk ki az amplitudót és a fázisszöget a megadott frekvenciákra: » [gain1,phase1]=bode(H,w1); » [gain2,phase2]=bode(H,w2); » [gain3,phase3]=bode(H,w3); » [gain4,phase4]=bode(H,w4); Készítsünk táblázatot a fenti frekvenciák, amplitudók és fázisszögek értékeivel. Az amplitudó és fázisértékek egyszerĦbben számíthatók, ha a frekvenciákat vektorban adjuk meg. Hozzunk létre egy lineáris frekvencia vektort. » w=1:0.1:10; » [gain,phase]=bode(num,den,w); Az amplitudó és a fázisdiagramok egyidejĦleg a képernyĘ különbözĘ ablakaiban jeleníthetĘk meg a subplot parancs segítségével. (A subplot(211) parancs 2 x 1 ablakot hoz létre a képernyĘn és ezek közül az elsĘre hivatkozik.) » subplot(211), plot(w,gain) » subplot(212), plot(w,phase) E két görbe alkotja a rendszer Bode diagramját. » subplot(111) % visszaállítjuk az egy ablakos ábrázolást. AN yquist diagram a komplex számsíkon ábrázolja a frekvenciafüggvény pontjait. » clf % clear figure, törli az elĘzĘ ábrát. » re=real(gain.*exp(j*phase*pi/180)); » im=imag(gain.*exp(j*phase*pi/180)); » plot(re,im) Az adott frekvenciavektor pontjaihoz a frekvenciafüggvény valós és képzetes részei közvetlenül is számíthatók a nyquist paranccsal. » [re,im]=nyquist(num,den,w) » plot(re,im) A frekvenciafüggvényt közvetlenül is ábrázolhatjuk. Megjelenítés: Ha a bode parancsot bal oldali változók megadása nélkül hívjuk meg, a diagram megjelenik a képernyĘn, de a számított adatok nem ĘrzĘdnek meg. » clf % törli az ábra ablakot. » bode(H),grid Figyeljük meg, hogy mind a frekvencia, mind pedig az amplitúdó léptéke logaritmikus. Hasonlóan a N yquist diagram az alábbi utasítás meghívásával jeleníthet Ę meg: » nyquist(H)

2



A frekvenciafüggvény pontjait elĘzetesen kiszámíthatjuk, majd azután ábrázolhatjuk. Számítás: Az amplitudó és fázisértékeket közvetlenül kiszámíthatjuk, ha az utasításokat bal oldali változók megadásával hívjuk meg. A számításhoz az átviteli függvényt számláló és nevezĘ polinomjaival adjuk meg (num,den), mivel az LTI struktúra egy 3dimenziós elrendezést (array) hoz létre. » [gain,phase,w]=bode(num,den) Ebben az esetben a Matlab számítja a frekvenciavektort is. (Egy Matlab utasításban a jobb oldali változók bemeneti, a bal oldali változók kimeneti változók.) Ábrázoljuk az amplitúdót kétszeresen logaritmikus léptékben a Matlab loglog utasításának meghívásával. Ez az utasítás hasonló a plot utasításhoz, de a megjelenítéshez mind az x, mind az y tengelyen logaritmikus skálát használ. » loglog(w,gain),grid Alkalmazzuk ismét a subplot utasítást az amplitúdó és a fázisszög megjelenítéséhez. » subplot(211) » loglog(w,gain) » subplot(212) » semilogx(w,phase) Logaritmikus frequencia vektort közvetlenül a logspace utasítással hozhatunk létre. » w=logspace(-2,2,100) Ez az utasítás 100 logaritmikusan ekvidisztáns pontot számol a 10-2=0.01 és 102=100 alsó és felsĘ frekvenciaértékek között. (Alapértelmezésben, a harmadik változó megadása nélkül az utasítás 50 pontot számol.) Ezzel a frekvenciavektorral számítsuk újra a Bode diagramot. » [gain,phase]=bode(num,den,w) » subplot(111) » loglog(w,gain) Vegyük észre, hogy a frekvenciatartomány eltér az elĘbbitĘl.

Feladat: A rendszer átviteli függvénye: H(s)=10/(s^2+ 2s+ 5) A bemenĘjel: u(t)=A*sin(2t+ phase), t t 0 A kimenĘjel kvázistacionárius állapotban: y(t)=5*sin(2t+ phase) Határozzuk meg az A és phase paraméterek értékeit (alkalmazzuk a bode utasítást). A N yquist diagram alakja jellemzi a rendszert. An alizálva a yNquist diagramot a rendszer fontos tulajdonságairól kaphatunk képet. A Bode diagram egyik elĘnye rendszerek soros eredĘi frekvenciafüggvényeinek számításakor mutatkozik. A soros eredĘ a frekvenciafüggvények összeszorzásával kapható meg. A logaritmikus lépték miatt a Bode diagramok egyszerĦen összeadódnak. A Bode diagram másik elĘnye, hogy rendszerint jól közelíthetĘ aszimptotáival. A közelítĘ diagramok jellegébĘl és töréspontjaiból a rendszer tulajdonságairól gyors értékelést adhatunk. Az alábbi példa szemlélteti, hogy a közelítĘ Bode amplitúdó diagram széles frekvenciatartományban valóban jól közelíti a pontos görbét.

Pontos és közelítĘ Bode diagram: Határozzuk meg az alábbi átviteli függvénnyel adott rendszer pontos és közelítĘ Bode amplitúdó diagramját.

H ( s)

10 s (1 s )(1 10s ) 2

3



Lineáris Rendszer Elemei


A rendszer H(s) átviteli függvénye ún. idĘállandós formában adott. Határozzuk meg erĘsítés-zérus-pólus alakban is.

H ( s ) 0.1

1 1 1 s (0.1 s ) 2 1 s

Egy lineáris rendszer ún. idĘállandós alakban általánosan az alábbi átviteli függvénnyel adható meg: c

» s=zpk(’s’) » H=10/(s*(1+s)*(1+10*s)*(1+10*s)) » [num,den]=tfdata(H,’v’) » p=roots(den) A rendszer pólusai: 0, - 0.1, - 0.1, -1. A pólusok (a negatív elĘjel nélkül) adják a közelítĘ Bode diagram töréspontjait. Vezessük be az alábbi frekvenciatartományokat: w 1:0.02-0.1; w 2:0.1-1;w 3:1-4 » w1= logspace (log10(0.02),-1,10); » w2=logspace(-1,0,10); » w3=logspace(0,log10(4),10); A teljes frekvenciatartományt a résztartományok egyesítésével képezzük: » w= [w1 w2 w3]; 1. Számítsuk ki a közelítĘ amplitúdó értékeket az

4. Lineáris Rendszer Elemei

1 tényezĘre. s

H (s)

j 0j

e sTD

1 f

j 0j

2 0j

ss T )

1

ahol k az átviteli vagy erĘsítési tényezĘ, i az integrátorok száma, TD a holtidĘ, Ĳ és T idĘállandók, a Ȣ paraméterek a csillapítási tényezĘk. Az alábbiakban a legfontosabb elemek idĘ- és frekvencia tartománybeli tulajdonságait vizsgáljuk. Ezek az elemek az arányos, integráló, differenciáló, holtidĘs és tárolós jellegĦ tagok, illetve ezek soros kapcsolásával adódó eredĘk.

1. Arányos (P) tag: H ( s)

k

Az erĘsítési tényezĘ k , a fázisszög pedig zérus minden frekvencián.

1 tényezĘre. A törésponti frekvenciánál (0.1 s ) 2 1 konstanssal közelíthetĘ. A töréspont felett kisebb frekvenciákra (kisfrekvenciás közelítés) a tag az 0.12 1 (nagyfrekvenciás közelítés) a tag az 2 átviteli s 4

10

2. Integráló (I) tag: H ( s )

0

10

-2

10

» w3L=[w1 w2]; w3H=w3; 10 » gain3L=0*w3L+1; gain3H=1./w3H; 10 10 10 10 » gain3=[gain3L gain3H]; A teljes rendszer erĘsítése: » gainappr=0.1*gain1.*gain2.*gain3; Számítsuk ki a pontos Bode diagram amplitúdó értékeit. » gainexact=bode(num,den,w); Ábrázoljuk a két görbét egy diagramban. » loglog(w, gainexact,’b’,w,gainappr,’r’),grid Látható, hogy a közelítés elfogadható, az aszimptoták jól közelítik a pontos amplitúdó görbét. A legnagyobb az eltérés az Z 0.1 frekvencián, ami egy kettĘs pólushoz tartozó töréspont. Megjegyezzük, hogy az eltérés jóval nagyobb lehet gyengén csillapított konjugált komplex pólusok esetén ( ] 0.3 ). -4

-2

-1

0

1

k s

i=1, a szakasz egyetlen integráló hatást tartalmaz. További idĘállandós tagok nincsenek. Az integráló tag memória tulajdonságú, kimenetén a jel akkor lehet konstans, ha bemenetén a jel értéke zérus. A kimenĘjel aktuális értéke a múltbeli bemenĘjel értékeitĘl függ. Vizsgáljuk a

k átviteli függvénnyel adott tiszta integráló tag viselkedését k 1 és k s 1 , s

H 2 (s)

5 s

» clear % törlünk minden korábban definiált változót. » s=zpk(‘s’) % definiáljuk az s szimbolikus változót zpk alakban. » H1=1/s » H2=5/s Az átmeneti függvény (egységugrásra adott válasz): » figure(1), step(H1,’r’,H2,’g’), grid Bode diagram: » figure(2), bode(H1,’r’,H2,’g’), grid N yquist diagram: » figure(3), nyquist (H1,’r’,H2,’g’),grid Az átmeneti függvény végértékei: » y1inf=dcgain(H1) » y2inf=dcgain(H2) Látható, hogy az átmeneti függvény végértéke végtelen. A frekvenciafüggvény amplitúdója kisfrekvencián végtelen, míg fázisszöge minden frekvencián –90º.

3. Egytárolós arányos (PT1) tag: H ( s )

k 1 Ts

Adjuk meg az alábbi tag átmeneti függvényét, Bode és yNquist diagramját.

4

5 erĘsítési

tényezĘk mellett.

H1 ( s )

2

10

1 tényezĘre. 1 s

(1 sT ) (1 2] T

2

1 e

j

2. Határozzuk meg a közelítĘ amplitúdó értékeket az

értékeket az

s s 2W 02 j )

j

1

» gain1=1./w;

függvénnyel közelíthetĘ. » w2L=w1; w2H=[w2 w3]; » gain2L=1./((0.1+0*w2L).^2); » gain2H=1./(w2H.^2); » gain2=[gain2L gain2H]; 3.Hasonlóan adjuk meg a közelít Ę amplitúdó

k si

d

(1 sW ) (1 2] W



H1 (s)

k 1 Ts

2 1 10 s



(1-piros, 2-zöld, 3- kék) » H2=2/((1+10*s)*(1+2*s)) » H3=2/((1+10*s)*(1+2*s)*(1+0.1*s))

D efiniáljuk a rendszert. » H1=2/(1+10*s) illetve » H1=tf(2,[10, 1]) Az átmeneti függvény: » t=0:0.1:50; » y1=step(H1,t); » plot(t,y1),grid vagy egyszerĦen » step(H1)

Átmeneti függvények: » figure(1), step(H1,’r’,H2,’g’,H3,’b’),grid Bode diagramok: » figure(2), bode(H1,’r’,H2,’g’,H3,’b’), grid N yquist diagramok: » figure(3), nyquist(H1,’r’,H2,’g’,H3,’b’),grid

Vizsgáljuk a k és T paraméterek jelentését, hatását a rendszer válaszaira. A kimenĘjel állandósult értéke, y (t o f) a Laplace transzformáció végértéktételével számítható.

4. Integráló tárolós (IT) tag: H ( s )

Megjegyezzük, hogy a görbék kijelölt része megnagyítható a figure ablakban aktivizált zoom paranccsal.

y (t o f) lim sY ( s ) lim sH ( s )U ( s ) s o0

ahol U ( s ) a bemenĘjel Laplace transzformáltja. Egységugrás bemenĘjel esetén s o0

Ha H ( s )

s o0

1 s

lim H ( s ) s o0

H1 ( s ) , y (t o f)

2 1 10 s

2 s 0

Matlab segítségével: » y1inf=dcgain(H1) Vizsgáljuk a T idĘállandó hatását a tranziens viselkedésre. Ismételjük meg az » y1=step(2,[T 1],t); Matlab utasítást különbözĘ T értékekre. A szakasz átviteli függvénye megadható zérus-pólus alakban is:

H1 (s)=

kp

2 10(0.1 s )

s p

2 2 2 , H 2 (s) , H 3 (s) s s (1 10 s) s(1 10s )(1 2s ) Az integráló taghoz sorbakapcsolt tárolós tagok csatlakoznak. » H1=2/s » H2=2/(s*(1+10*s)) » H3=2/(s*(1+10*s)*(1+2*s)) Átmeneti függvények: » figure(1), step(H1,’r’,H2,’g’,H3,’b’),grid Bode diagramok: » figure(2), bode(H1,’r’,H2,’g’,H3,’b’),grid, N yquist diagramok: » figure(3), nyquist (H1,’r’,H2,’g’,H3,’b’),grid H1 ( s )

s o0

y (t o f) lim sH ( s )U ( s ) lim sH ( s )

0.2 . s 0.1

1 , T

kp

k T

A pólus abszolút értéke megadja a közelítĘ amplitúdó-frekvencia diagram töréspontját (sarokfrekvencia). Az átmeneti függvény állandósult értéke

y (t o f) lim H ( s ) s o0

megegyezik az amplitúdó-frekvencia függvény kisfrekvenciás értékével. A tag Bode és Nyquist diagramja: » bode(H1); » nyquist(H1); Hasonlóan vizsgáljuk az alábbi egytárolós, kéttárolós és háromtárolós arányos (PT1, PT2, PT3) tagok átmeneti és frekvenciafüggvényeit: 2 2 2 H1 ; H2 ; H3 1 10 s (1 10 s )(1 2 s ) (1 10 s )(1 2 s )(1 0.1s )

1 s 2T02 2] T0 s 1

5. Kéttárolós lengĘ tag: H ( s ) Vizsgáljuk az alábbi rendszert:

H ( s)

ahol

p

k s (1 Ts )

ahol Z 0

1 9s 2 2 s 1

1 s 2T02 2 ] T0 s 1

1 a természetes körfrekvencia és ] a csillapítási tényezĘ ( T0 T0

» num=1; » den=[9, 2, 1] » H=tf(num,den) A rendszer pólusai a roots utasítással számolva » roots(den) illetve a damp utasítással: » damp(H) A két konjugált komplex pólus megadható az alábbi alakban: p1 Az átmeneti függvény vt túllövése a következĘképpen számítható: ]S Sa a 1] 2 2 2 2 e b Z 0 a b , ] , vt e

b

1

Z0

a jb, p2

3 , ] =1/3).

a jb



A lengési körfrekvencia:

Zp

b Z0 1 9 2

Az átmeneti függvény elsĘ maximumának idĘpontja (csúcsidĘ): Tp

S Zp

S b

.

» zeta=1/3 » vt=exp(-zeta *pi/sqrt(1- zeta * zeta)) Az átmeneti függvény: » [y,t]=step(H); » plot(t,y), grid Az átmeneti függvény maximális értéke: » ym=max(y) Állandósult értéke: » ys=dcgain(H) A túllövés értéke másképp számítva: » yo=(ym-ys)/ys


» H2=(2*s)/(1+10*s) » H3=(2*s )/((1+10*s)*(1+2*s)) » H4=(2*s )/((1+10*s)*(1+2*s)*(1+2*s)) H1 ( s ) átmeneti függvénye: » figure(1), step(H1) A Matlab nem tudja kiértékelni a rendszer viselkedését, mivel a tag nem realizálható, számlálója magasabbfokú a nevezĘjénél. (Az átmeneti függvény iDrac delta). Átmeneti függvények: » figure(1), step(H2,’r’,H3,’g’,H4,’b’),grid Bode diagramok: » figure(2), bode(H2,’r’,H3,’g’,H4,’b’),grid, N yquist diagramok: » figure(3), nyquist(H2,’r’,H3,’g’,H4,’b’),grid

7. A zérusok hatása: H ( s )

Vizsgáljuk az átmeneti függvényeket, a Bode és a yNquist diagramokat különböz Ę csillapítási tényezĘk mellett. ] =0.3,0.7,1 . » zeta1=0.3, zeta2=0.7, zeta3=1 » T0=3 » H1=1/(s*s*T0*T0+2*zeta1*T0*s+1) » H2=1/(s*s*T0*T0+2*zeta2*T0*s+1) » H3=1/(s*s*T0*T0+2*zeta3*T0*s+1) Az átmeneti függvények: » figure(1), step(H1,’r’,H2,’g’,H3,’b’),grid A Bode diagramok: » figure(2), bode(H1,’r’,H2,’g’,H3,’b’),grid AN yquist diagramok: » figure(3), nyquist(H1,’r’,H2,’g’,H3,’b’),grid A pólus-zérus konfigurációk: » pzmap(H1) » pzmap(H2) » pzmap(H3) Látható, hogy 0.3 csillapítási tényez Ę mellett az átmeneti függvény a leglengĘbb, a túllendülés a legnagyobb, és a Bode amplitudó diagram kiemelése is a legnagyobb. A N yquist diagram metszéke az imaginárius tengellyel ennél a csillapításnál a legnagyobb. A pólusok imaginárius értéke, ami a lengési körfrekvenciát meghatározza szintén itt a legnagyobb. A Bode amplitudó diagram kiemelése összefüggésben van az átmeneti függvény túllendülésével. Ha el akarjuk kerülni az idĘtartományban a nagy túllendülést, nem engedünk meg a frekvenciatartományban nagy amplitúdó kiemelést. A Ȣ=0.7 csillapítási tényezĘ kedvezĘ viselkedést, gyors beállást és kis túllendülést mutat. Szabályozási rendszerek tervezhetĘk ilyen dinamikus viselkedésre.

6. Differenciáló (D és DT) tag: H ( s ) sTd 2s 2s H1 ( s ) sTd 2s, H 2 ( s ) , H3 , H4 1 10s (1 10 s )(1 2 s )


k ( s z1 )( s z2 )...( s zM ) D( s)

Vizsgáljuk a zérus hatását az átmeneti függvényre és a frekvenciafüggvényre az alábbi átviteli függvény esetén: 1W s H (s)

1 s 1 10 s

W , a számláló idĘállandója (a zérus értéke 1/ W ) –4 és 4 értékek között változik. Pozitív (a komplex számsík jobb oldalán lévĘ) zérus esetén a tagot nem-minimumfázisúnak nevezzük. » s=tf('s') » tau=[-4 -2 0 2 4]; » D=50*s*s+11*s+1; » for i=1:5, » H(i)=(s*tau(i)+1)/D » end » t=0:0.1:60; » Y1=step(H(1),t); Y2=step(H(2),t); Y3=step(H(3),t); » Y4=step(H(4),t); Y5=step(H(5),t); » plot(t,[Y1,Y2,Y3,Y4,Y5]),grid,shg Megjegyezzük, hogy a nem-minimumfázisú rendszerek szokatlan viselkedést mutathatnak, például egy jobb oldali zérus esetén az átmeneti függvény kezdetben az állandósult értékével ellentétes irányban indul el, majd irányt váltva éri el végül állandósult értékét. A zérus beiktatása gyorsítja a rendszert. Bode diagramok: » figure(1) » bode(H(1),’r’,H(2),’g’,H(3),’k’,H(4),’y’,H(5),’b’),grid N yquist diagramok: » figure(2) » nyquist(H(1),’r’,H(2),’g’,H(3),’k’,H(4),’y’,H(5),’b’),grid É rtékeljük a zérus hatását a frekvenciatartományb an a Bode és a yNquist diagramok alakulására!

8. HoltidĘs tag: H ( s ) e sTD

2s (1 10 s )(1 2s )(1 0.1s )

Adjuk meg ismét a tag átmeneti függvényét, Bode és N yquist diagramját. » H1=2*s

A holtidĘs tag leírása az idĘ- és a Laplace operátor tartományban:

y (t ) o y (t TD )

Y ( s ) o Y ( s )e sTD Az átviteli függvény:


H D ( s)


H ( s )e sTD

A frekvenciatartományban:

e jZTD erĘsítés:

1, arg ê jZTD ` ZTD HD

H

fázisszög: arg( H D ) arg( H ) ZTD Vizsgáljuk az alábbi átviteli függvényĦ tagok frekvenciafüggvényeit:

H1 ( s )

1 , H 2 (s) 1 10s

1 e 2 s (1 10 s )

A holtidĘs tag amplitúdója és fázisszöge: H1 H 2 , gain2=gain1

arg( H 2 ) arg( H1 ) ZTD , phase2= phase1 ZTD » Td=2 » num1=1 » den1=[10, 1] A Bode diagram számításakor most a num, den számláló és nevezĘpolinommal adott alakot kell használni. ElĘször hozzuk létre a frekvenciavektort. » w=logspace(-2,0,100) » [gain1,phase1]=bode(num1,den1,w) » delay=180/pi*Td*w' % a holtidĘs tag ZTD radiánról fokra számítjuk át. Az amplitúdó és a fázisszög a holtidĘ figyelembe vételével: » gain2=gain1 » phase2=phase1-delay » subplot(211),loglog(w, gain1,’r’,w, gain2,’b’),grid; » subplot(212),semilogx(w, phase1,’r’,w, phase2,’b’),grid A fázisszög linearitása jól látható, ha ábrázolására a semilogx utasítás helyett a lineáris léptékĦ plot utasítást használjuk. » figure(2),subplot(111),plot(w, phase1,’r’,w, phase2,’b’),grid Ábrázoljuk a N yquist diagramot. Számítsuk ki el Ęször a valós és a képzetes értékeket. » h1= gain1.*exp(j*phase1*pi/180) » h2= gain2.*exp(j*phase2*pi/180) » plot(real(h1),imag(h1),’r’, real(h2),imag(h2),’b’) A nagyfrekvenciás viselkedés jobban vizsgálható, ha nagyobb frekvenciatartományt jelölünk ki a logspace utasítással. Az idĘtartománybeli viselkedés Simulink blokk-diagram megépítésével és a program futtatásával követhetĘ a legjobban, mivel a holtidĘs tag a Simulink program építĘeleme. A holtidĘs tag átviteli függvénye nem racionális törtfüggvény, a Matlab program csupán közelíteni tudja egy nem-minimumfázisú racionális törtfüggvénnyel, az ún. Pade közelítéssel. A holtidĘs átviteli függvény és a Pade közelítés Taylor sorának elsĘ elemei megegyeznek. Minél magasabbfokú a törtfüggvény, annál jobb a közelítés. A közelítĘ átmeneti függvény zérus helyett + 1 vagy –1 értékb Ęl indul. A Matlab-ban a pade utasítás szolgál a közelítés megadására. A Pade közelítés szemléltetésére alkalmazzunk ötödfokú közelítést.

H deadtime ( s ) e sTD



» H1=tf(num1,den1) » [numpade,denpade]=pade(Td,5) » Hdeadtime=tf(numpade,denpade) » H2=H1* Hdeadtime Az átmeneti függvény: » figure(1), step(H1,’r’,H2,’g’),grid A Bode diagram: » figure(2), bode(H1,’r’,H2,’g’),grid AN yquist diagram: » figure(3), nyquist(H1,’r’,H2,’g’),grid Feladat: Vizsgáljuk a Pade közelítés jóságát a fenti egytárolós holtidĘs tag esetén. É pítsünk egy Simulink programot a “transport delay”blokk, illetve a Pade közelít Ę ötödfokú racionális törtfüggvénnyel. Hasonlítsuk össze az átmeneti függvényeket.

9. KettĘs integrátor: H ( s )

H1 ( s )

4 s2

, H 2 (s)

k s2 4

s 2 (1 10 s )

, H3

4 s 2 (1 10 s )(1 2 s )

, H4

4 s 2 (1 10 s )(1 2 s )(1 0.1s )

Feladat: Határozzuk meg a megadott tagokra az átmeneti függvényeket valamint a Bode és a yNuqist diagramokat.

Összefoglalás: Az átmeneti függvény értéke állandósult állapotban ( t o f ) és a frekvenciafüggvény amplitúdója Z o 0 értéknél megegyezik. Az arányos jellegĦ (tiszta arányos tagot és tárolós tagokat tartalmazó) elemek yNuqist diagramja Z 0 értéknél a komplex számsík valós tengelyének egy pozitív pontjából indul, ez az érték a szakasz átviteli tényezĘje. Integráló jellegĦ elemek yNuqist diagramja az imaginárius tengely negatív végtelen irányából indul. A kétszeresen integráló (esetleg tárolós tagokat is tartalmazó) elem N yuqist diagramja a negatív valós tengely végtelen irányából indul. iDfferenciáló jelleg Ħ elemek yNuqist diagramja a komplex számsík origójából indul a pozitív imaginárius tengely irányába. Tárolós tagok esetén (zérusok nélkül) a yquist diagram annyi síknegyeden halad át, ame nnyi a tárolók száma. A fázisszög a frekvenciával N monoton módon változik. A zérusok eltorzítják a fázisszög monoton változását, meghatározott frekvenciatartományokban a görbén b„ ehor padásokat”okoznak. Arányos jelleg Ħ tagok Bode amplitudó diagramja a frekvencia tengellyel párhuzamosan indul, zérus fáziseltolással. Egy integrátort tartalmazó rendszer Bode amplitudó diagramja 20dB / dekád meredekséggel indul 09 0 fázisszöggel, míg a kettĘs integrátort tartalmazó rendszer Bode amplitúdó diagramja 40dB / dekád meredekséggel indul 18 0 0 fázisszöggel. A tárolós tagok a töréspontokban 20dB / dekád meredekséggel változtatják az aszimptotikus Bode amplitúdó görbe meredekségét. A zérusok megjelenése 20dB / dekád meredekségváltozást okoz.

Visszacsatolás Vizsgálata




1. Példa.

5. Visszacsatolás Vizsgálata

A felnyitott kör átviteli függvénye: A szabályozástechnika alapvetĘ struktúrális eleme a visszacsatolás. A zárt körben az egyes jelek közötti eredĘ átviteli függvények (pl. a szabályozott jellemzĘ és az alapjel közötti eredĘ átviteli függvény) könnyen kiszámíthatók. Az eredĘ frekvenciafüggvények egyszerĦen származtathatók az átviteli függvényekbĘl s=jZ helyettesítéssel. A MATLAB közvetlenül támogatja a tipikus blokkdiagram-algebrai struktúrákat, a sorosan kapcsolt illetve a visszacsatolt körök átviteli függvényeinek számítását. E függvényhívások rövid összefoglalására tekintsük az alábbi zárt szabályozási kört (C - szabályozó, P folyamat, F - visszacsatolás):

L( s )

K s(1 sT1 )

Egységnyi negatív (merev) visszacsatolást alkalmazunk. A zárt kör eredĘ átviteli függvénye:

T (s)

K s 2T1 s K

1 T 1 1 s s2 1 K K

Az eredĘ másodrendĦ lengĘ tag, amely általánosan az alábbi alakban írható fel:

r(t) R(s)

-

e(t) E(s)

C ( s)

u(t) U(s)

P( s)

y(t) Y(s)

1 1 29 T0 s T02 s 2

Az együtthatók összehasonlításával: T0

T1 K

, 9

1 2 KT1

F (s) Ábrázoljuk a felnyitott és a zárt rendszer Bode amplitúdó diagramjait egy ábrában T1=1, és K=0.1 valamint K= 4 esetén. Számítsuk ki a felnyitott kör frekvenciafüggvényének jellegzetes pontjait. A szabályozási körben az alábbi eredĘ átviteli függvények számíthatók: A felnyitott kör átviteli függvénye (hurokátviteli függvény): L( s ) C ( s ) P( s ) F ( s ) A zárt szabályozási kör eredĘ átviteli függvényei, a hibajel, a beavatkozójel, illetve a kimenĘjel (szabályozott jellemzĘ) alapjelre vonatkozó eredĘ átviteli függvényei: E (s) 1 U (s) C (s) Y ( s) C ( s) P( s) , Wu ( s ) , T (s) We ( s ) R ( s ) 1 L( s ) R( s ) 1 L( s ) R( s) 1 L( s) Az elĘrevezetĘ ág átviteli függvénye: E=CP: » E=series(C,P); vagy » E=C*P; A hurokátviteli függvény: L=CPF=EF » L=series(E,F); vagy » L=C*P*F; A zárt körban a kimenĘjelnek az alapjelre vonatkozó eredĘ átviteli függvénye: » T=feedback(E,F,-1); Egységnyi visszacsatolás esetén (F=1) a zárt kör eredĘ átviteli függvénye a cloop utasítással számítható: » T=cloop(E,-1); Az eredĘ átviteli függvény közvetlenül az LTI struktúrával szimbólikusan is számítható: » T=C*P/(1+C*P*F) A minreal utasítás meghívásával kiejthetĘk a megegyezĘ zérusok és pólusok. » T=minreal(T) Az egyszerĦsítéshez tolerancia is megadható. A tolerancia értéke alapértelmezésben sqrt(eps)=1.4901e-008. » T=minreal(T,0.001) Ezesetben az egyszerĦsítés akkor történik meg, ha a pólusok és zérusok értékének különbsége kisebb 0.001-nél. Az átviteli függvényekbĘl ezután a bode vagy a nyquist függvényhívással számítható illetve jeleníthetõ meg a frekvenciafüggvény.

A frekvencia függvény: Vizsgáljuk a felnyitott és a zárt rendszer frekvenciafüggvényeinek kapcsolatát, a zárt szabályozási kör átmeneti függvényét és annak kapcsolatát a frekvenciafüggvénnyel.

» » » » » »

s=tf('s') T1=1 K1=0.1 K2=4 L1=K1/(s*(s+1)) L2=K2/(s*(s+1))

A kimenĘjel és az alapjel közötti eredĘ átviteli függvény az alábbi összefüggéssel számítható: L( s ) num num L( s ) , T ( s) 1 L( s ) num den den » T1=L1/(1+L1) » T1=minreal(T1) % kiejti az azonos pólus-zérus párokat. vagy » T1=feedback(L1,1) » T2=feedback(L2,1) Adjuk meg az átmeneti függvényeket, határozzuk meg a felnyitott és a zárt kör frekvenciafüggvényeit. » » » »

figure(1) step(T1,'r',T2,'y') figure(2) bode(L1,'b',L2,'c',T1,'r',T2,'y')

Egységnyi merev visszacsatolás esetén a zárt rendszer Bode amplitudó-körfrekvencia diagramja a kisfrekvenciás tartományban közel egységnyi (nagy kisfrekvenciás erĘsítés esetén), a nagyfrekvenciás tartományban közelítĘen egybeesik a felnyitott rendszer amplitúdó menetével. A felnyitott rendszer diagramjából számított vágási körfrekvencia környékén nagy K értékeknél a zárt rendszer amplitúdó diagramjában nagy kiemelés lehet. Megfigyelhetjük, hogy minél nagyobb a zárt rendszer frekvenciafüggvényének kiemelése, annál nagyobb az átmeneti függvény túllendülése.



Stabilitásvizsgálat


Határozzuk meg a zárt rendszer csillapítási tényezĘjét.

6. Stabilitásvizsgálat

» damp(T1) Eigenvalue -1.13e-001 -8.87e-001

Damping 1.00e+000 1.00e+000

Egy rendszer stabilis, ha kitérítés után a magára hagyott rendszer visszatér nyugalmi állapotába, tranziensei lecsengenek. A stabilitás egy másik megfogalmazása szerint egy rendszer stabilis, ha korlátos bemenet esetén a kimenet is korlátos.

Freq. (rad/s) 1.13e-001 8.87e-001

» damp(T2) Eigenvalue -5.00e-001 + 1.94e+000i -5.00e-001 - 1.94e+000i

Damping 2.50e-001 2.50e-001

Freq. (rad/s) 2.00e+000 2.00e+00

Tárolhatjuk a a csillapítási tényezĘt és a frekvencia értékét változókban is. » [w01,zeta1]=damp(T1) » [w02,zeta2]=damp(T2) A különbözĘ erĘsítési tényezĘk lényegesen eltérĘ viselkedést eredményeznek. Az elsĘ esetben a pólusok valósak, míg a második esetben konjugált komlexek, amelyek az idĘtartományban lengĘ viselkedést eredményeznek.

Tekintsük az alábbi hatásvázlattal megadott szabályozási kört. A felnyitott rendszer átviteli függvénye L(s) . A rendszerben merev (egységnyi) negatív visszacsatolást alkalmazunk:

r(t)

-

e(t) E(s)

u(t) U(s)

L( s)

T ( s)

k (1 s )(1 5s )

10 (1 s )(1 5s )

1 s

A felnyitott kör és a zárt kör átmeneti függvényeinek állandósult értéke: s o0

yzárt (t o f) lim T ( s ) s o0

k 10

ynyitott (t o f) 1 ynyitott (t o f)

A hiba értéke állandósult állapotban: e(t o f) 1 ynyitott (t o f)

T (s)

y(t) Y(s)

L( s ) 1 L( s )

A rendszer stabilitása meghatározható a zárt rendszer pólusai alapján. A rendszer stabilis, ha pólusai a komplex bal félsíkra esnek, azaz összes pólusának a valós része negatív.

1. Példa Egy zárt rendszer eredĘ átviteli függvénye:

T (s)

1 y (t o f) lim sR( s ) H ( s ) lim s H ( s ) lim H ( s ) s o0 s o0 s o0 s ynyitott (t o f) lim L( s )

r(t)

1. Stabilitásvizsgálat a zárt rendszer pólusai alapján:

Vizsgáljuk a felnyitott és a zárt rendszer viselkedését. » s=zpk('s') » L=10/((1+s)*(5*s+1)) » T=feedback(L,1) Határozzuk meg a felnyitott és a zárt kör átmeneti függvényeinek állandósult értékeit. A Laplace transzformáció végértéktétele szerint egységugrás alapjel esetén

r (t ) 1(t ), R( s )

{

A zárt rendszer eredĘ átviteli függvényét a következĘképpen számolhatjuk:

2. Példa. Egy rendszer hurokátviteli függvénye:

L( s )

y(t) Y(s)

k 1 k

10 1 10

1 1 k

1 11

Ábrázoljuk a felnyitott és a zárt kör átmeneti függvényeit: » step(L,’r’,T,’b’) Az állandósult értékek leolvashatók a görbékrĘl vagy kiszámíthatók a dcgain utasítással. » yos=dcgain(L) » ycs=dcgain(T) Jelenítsük meg a felnyitott és a zárt kör Bode diagramjait egy ábrán. » bode(L,’r’,T,’b’) Határozzuk meg az állandósult hiba értékét k=1, 20, 100 esetére.

s5 . s 5 3s 4 4 s 3 10s 2 5s 10

Állapítsuk meg, hogy a rendszer stabilis-e. Megoldás: Vizsgáljuk meg a pólusok elhelyezkedését. » num=[1, 5] » den=[1, -3, 4, 10, 5, -10] » T=tf(num,den) » poles=roots(den) poles = 2.1150 + 2.1652i 2.1150 - 2.1652i -0.9824 + 0.7214i -0.9824 - 0.7214i 0.7348 vagy másképpen » [z,p,k]=zpkdata(T,'v') A rendszer labilis, mivel vannak pozitív valós részĦ pólusai. A pzmap utasítás grafikusan jeleníti meg a pólusok elhelyezkedését: » pzmap(T) Látható, hogy a rendszer labilis, mivel vannak pólusai a komplex számsík jobb oldali részén.



2. A Nyquist stabilitási kritérium alkalmazása: A zárt rendszer stabilitása meghatározható a felnyitott rendszer frekvenciafüggvénye alapján is. a. Az egyszerĦsített Nyquist kritérium akkor használható, ha a felnyitott rendszernek nincs labilis (pozitív valós részĦ) pólusa. A zárt rendszer stabilis, ha a felnyitott rendszer Nyquist diagramja nem veszi körül a (–1+0j) pontot. b. Az általánosított Nyquist kritériumot kell használni akkor, ha a felnyitott rendszernek van labilis pólusa. A zárt rendszer stabilis, ha a felnyitott rendszer Nyquist diagramja annyiszor veszi körül a (–1+0j) pontot az óramutató járásával ellentétes (pozitív) irányban, ahány labilis pólusa van a felnyitott rendszernek.



» T=feedback(L,1) » step(T) » [z,p,k]=zpkdata(T,'v') » pzmap(T) Vizsgáljuk meg a rendszer stabilitását, ha a pólusok elĘjelét megváltoztatjuk.

L( s )

5 (1 10s )(1 0.1s )

Milyen irányban kerüli meg a Nyquist diagram a -1 pontot? Stabilis-e a zárt rendszer?

2. Példa A felnyitott kör átviteli függvénye

L( s )

10 . (1 10 s)(1 s )

Negatív merev visszacsatolást alkalmazunk. A Nyquist kritérium alapján határozzuk meg, hogy a zárt rendszer stabilis-e. » s=tf('s') » L=10/((1+10*s)*(1+s)) » [z,p,k]=zpkdata(L,'v')

Van-e labilis (pozitív valós részĦ) pólusa a nyílt rendszernek? (Nincs, tehát az egyszerĦsített Nyquist kritérium használható.) » nyquist(L),grid Körülveszi-e a Nyquist diagram a (–1+0j) pontot? (Nem, tehát a rendszer stabilis.) Stabilis marad-e a rendszer, ha a 10-es erĘsítést a számlálóban megnöveljük? EllenĘrizzük a stabilitást a zárt rendszer pólusainak vizsgálatával. A feedback utasításban az 1 a merev (egységnyi), a -1 pedig a negatív visszacsatolást jelenti. (A -1 elhagyható, mert az alapértelmezés a negatív visszacsatolás)

» T=feedback(L,1,-1) vagy » T=L/(1+L); T=minreal(T) A minreal utasítás elvégzi a lehetséges egyszerĦsítéseket, kiejti a közös zérusokat és pólusokat. » [z,p,k]=zpkdata(T,'v') Ábrázoljuk a zárt rendszer pólusait: » pzmap(T) Láthatjuk, hogy a rendszer struktúrálisan stabilis. A Nyquist diagram nem veszi körül a (–1+0j) pontot még megnövelt erĘsítés esetén sem.

3. Példa A felnyitott rendszer átviteli függvénye

L( s )

5 . (1 10s )(1 0.1s )

Negatív merev visszacsatolást alkalmazunk. A Nyquist kritérium alapján határozzuk meg, hogy a zárt rendszer stabilis-e? » L=-5/((1-10*s)*(1+0.1*s)) » [z,p,k]=zpkdata(L,'v') Van-e labilis pólusa a felnyitott rendszernek? (Igen, p2 0.1 ) » nyquist(L) Milyen irányban kerüli meg a görbe a (–1+0j) pontot? (Az óramutató járásával ellentétes irányban) Stabilis-e a zárt rendszer? (Igen) EllenĘrizzük az eredményt a zárt rendszer pólusainak kiszámításával.

4. Példa Egy felnyitott rendszer átviteli függvénye

L( s )

k

1 s . (1 s)(1 0.5s)

Negatív merev visszacsatolást alkalmazunk. Határozzuk meg a k paraméter azon értékeit, amelyekre a zárt rendszer stabilis lesz. a. k=1 értéket feltételezünk: » L=(1-s)/((1+s)*(1+0.5*s)) » [z,p,k]=zpkdata(L,'v') Van-e labilis pólusa a felnyitott rendszernek? (Nincs, tehát az egyszerĦsített Nyquist kritérium használható) » nyquist(L),grid Határozzuk meg azt a pontot, ahol a görbe metszi a valós tengelyt (-0.666). A zoom utasítás vagy a zoom menĦpont és az egér segítségével a Nyquist diagramról leolvashatjuk ezt az értéket. Ha egy rendszer erĘsítését k-val megszorozzuk, akkor a Nyquist diagramját k-szorosára nagyítjuk. A rendszer akkor van a stabilitás határhelyzetében, ha a metszéspont a -1 értékre esik. Tehát a rendszer erĘsítését még megnövelhetjük k=1/0.666=1.5 értékkel. A rendszer stabilis, ha teljesül, hogy 0< k<1.5 (pozitív k-ra) b. k=-1 értéket feltételezünk: » nyquist(-L), grid Határozzuk meg ismét azt a pontot, ahol a görbe metszi a valós tengelyt (-1). Tehát a rendszer stabilis, ha k>-1>0. Tegyük össze a két feltételt. A zárt rendszer stabilis, ha –1
Fázistartalék (fázistöbblet), erĘsítési tartalék: Egy rendszer stabilitásvizsgálatánál nemcsak az a fontos, hogy a rendszer stabilis-e, hanem az is, hogy milyen messze van a stabilitás határhelyzetétĘl. A fázistartalék és az erĘsítési tartalék fejezik ki ezeket a jellemzĘket. A fázistartalék azt adja meg, hogy a vágási körfrekvenciánál a rendszer fázisszögét még mennyivel lehet megnövelni, hogy a rendszer a stabilitás határhelyzetébe ( -180°) kerüljön. A fázistartalékot meghatározhatjuk az Z c vágási körfrekvencián felvett fázisszög értékébĘl.

Mm

M (Zc ) 180q



c

Ha a fázistartalék pozitív, a rendszer stabilis. Például, ha a rendszer fázisszöge a vágási körfrekvencián M 120q , akkor a fázistartalék Mm M 180q 120q 180q 60q , azaz a rendszer stabilis. Az erĘsítési tartalékkal a rendszer erĘsítését megszorozva a rendszer a stabilitás határhelyzetébe kerül.

1 L (ZS )

,

ahol

ZS : M (Z )Z ZS 180 ZS az a frekvencia, ahol a fázisszög -180˚. Ha az erĘsítési tartalék nagyobb mint egy, akkor a rendszer stabilis. A fázistartalék szemléltethetĘ a Nyquist és a Bode diagramon is. A margin utasítás segítségével egyszerĦen kiértékelhetĘk és ábrázolhatók.

5. Példa Egy felnyitott rendszer átviteli függvénye

L( s )

1 . (0.5 s )( s 2 2 s 1)

Negatív merev visszacsatolást alkalmazunk. Határozzuk meg a fázistartalékot, az erĘsítési tartalékot és a felnyitott rendszer vágási körfrekvenciáját. » L=1/((0.5+s)*(s^2+2*s+1)) » [gm,pm,wg,wc]=margin(L) A gm (gain margin) az erĘsítési tartalék, pm (phase margin) a fázistartalék , wc (cut-off frequency) a vágási körfrekvencia és wg az a körfrekvencia, ahol a fázisszög értéke -180º. Grafikusan is meg lehet jeleníteni ezen értékek elhelyezkedését: » margin(L); A grafikus megjelenítés esetén a Gm érték decibelben adott, Gm=20*log10(gm) Az erĘsítés, fázis és frekvencia értékek táblázatban is megjeleníthetĘk. » w=logspace(-1,1,100); » [num,den]=tfdata(L,'v') » [mag,phase]=bode(num,den,w); » Tabl=[mag, phase, w’] Bode Diagrams Gm=13.064 dB (at 1.4142 rad/sec), Pm=72.227 deg. (at 0.56754 rad/sec)

phase

20

w

1.1123 -99.5242 0.5094 1.0643 -103.0406 0.5337 1.0158 -106.6104 0.5591§wc 0.9669 -110.2286 0.5857 ……. 0.2449 -176.7848 1.3530 0.2211 -180.1658 1.4175§wg 0.1991 -183.4774 1.4850 0.1789 -186.7160 1.5557 …….

A fázistartalék kiolvasható a táblázatból. Az 1.0 erĘsítés értékhez tartozó fázisszög

0

Phase (deg); Magnitude (dB)

mag

-20 -40 -60 -80 0

-100

-200

-300 10 -1

10 0

Frequency (rad/sec)


értékébĘl lehet kiszámolni. pm=180-106.6=73.4. A vágási körfrekvencia w=wc=0.5591. Az erĘsítési tartalék a -180°-hoz tartozó erĘsítés reciproka lesz, gm=1/0.221=4.52. Az így számolt értékek pontosságát a táblázat felbontása korlátozza.

Z c az a frekvencia, ahol a felnyitott rendszer frekvenciafüggvényének abszolút értéke 1. Z c : L( jZ ) Z Z 1

gm


10 1

Soros PID Kompenzáció


7. Soros PID Kompenzáció Tekintsük az alábbi zárt szabályozási kört. P( s ) a folyamat (szabályozott szakasz), C ( s ) a szabályozó átviteli függvénye. r(t)

e(t)

R(s)

E(s)

-

C (s)

u(t) U(s)

P(s )

y(t) Y(s)

Az adott szabályozandó folyamathoz olyan sorbaiktatott szabályozót kívánunk megadni, amellyel a zárt szabályozási rendszer stabilis mĦködésĦ és megfelel a minĘségi elĘírásoknak. ElĘször megfogalmazzuk a leglényegesebb minĘségi elĘírásokat, majd néhány szabályozó tervezési eljárást mutatunk be.

MinĘségi elĘírások: Zárt szabályozási körökkel szemben általában az alábbi követelmények, gyakorlati igények merülnek fel: Stabilitás: AlapvetĘ követelmény a szabályozási rendszer stabilis mĦködése. A stabilitást különbözĘ módon fogalmazhatjuk meg. A BIBO (Bounded Input - Bounded Output) stabilitás azt jelenti, hogy a rendszer korlátos bemenĘjelre korlátos kimenĘjelet ad. Aszimptotikusan stabilis a rendszer, ha tranziensei lecsengenek. Robusztusság: a zárt szabályozási kör viselkedése ne legyen érzékeny a folyamatról rendelkezésünkre álló információ pontatlanságára. A stabilitást akkor is biztosítani kell, ha a rendszer paraméterei nem pontosan ismertek, a névleges értékük körül adott tartományban változnak. Statikus viselkedés: Egy másik fontos követelmény a szabályozási rendszer statikus vagy állandósult állapotbeli pontossága. Az állandósult állapotra vonatkozó statikus elĘírások: -Alapjelkövetés: a szabályozott jellemzĘ kövesse az alapjelet (referencia jelet). A követési hiba az elĘírt értéken belül legyen. -Zavarelhárítás: a szabályozás állandósult állapotban küszöbölje ki a fellépĘ zavarások hatását. A statikus pontosság függ a rendszer struktúrájától és a bemenĘjelektĘl is. Tranziens viselkedés: A rendszer tranziens válasza a szabályozási kör lényeges dinamikus tulajdonsága. A tranziens viselkedés tulajdonságai megfogalmazhatók az idĘtartományban a rendszer egységugrás alapjelre vagy zavarásra adott válaszának tranziens tulajdonságaival. A tranziens beállásra vonatkozó elĘírások: -A kimenĘjel túllövése, - beállás ideje (szabályozási idĘ). A beállási idĘ alatt a szabályozott jellemzĘ 1-2%-on belül megközelíti állandósult értékét. Általában egy kis, 5-10%-on belüli túllövés tolerálható, de vannak olyan folyamatok, ahol aperiodikus tranzienseket írnak elĘ. A minĘségi elĘírások a frekvenciatartományban is megfogalmazhatók. Az átmeneti függvény túllendülése összefüggésben van a zárt kör amplitúdó - frekvencia függvényének maximális értékével illetve a felnyitott rendszer frekvenciafüggvényébĘl számítható fázistartalékkal. Ha a fázistartalék körülbelül 60q, akkor a zárt rendszer átmeneti függvényének túllövése általában 5-10%-on belül lesz.



A beállási idĘ a vágási körfrekvenciától függ és az alábbi összefüggéssel becsülhetĘ:

3

Zc

d ts d

10

Zc

.

Beavatkozó jel nagysága: A gyakorlatban a beavatkozójel rendszerint korlátos. Követelmény, hogy a beavatkozójel maximuma ne haladja meg az elĘírt korlátot. A szabályozó paramétereit a tervezési követelmények alapján kell megválasztani. Bonyolultabb szabályozási problémák esetén a kimenĘjel és a beavatkozójel teljes lefolyására elĘírhatók korlátozások. Például a szabályozási hiba és a beavatkozójel négyzetintegrálját minimalizáljuk. Ezek a korlátozások nemlineárisak és egymásnak ellentmondóak is lehetnek. Általános esetben csak valamilyen optimalizálási eljárással lehet a szabályozó paramétereit meghatározni. Itt egy egyszerĦ és a gyakorlatban is gyakran alkalmazott módszert mutatunk be, amikor a rendszer fázistartalékát írjuk elĘ. A tervezés a felnyitott kör frekvenciafüggvénye alapján történik, azaz a felnyitott kör viselkedésébĘl a zárt rendszer tulajdonságaira következtetünk. A felnyitott kör átviteli függvénye L( s ) C ( s ) P( s ) . A zárt szabályozási kör eredĘ átviteli függvényei: E (s) 1 U ( s) C ( s) Y (s) C ( s) P( s) , Wu ( s ) , T (s) We ( s ) R(s) 1 C (s) P(s) R( s) 1 C ( s) P( s) R( s) 1 C ( s) P( s) Itt r(t) az alapjel, u(t) a beavatkozójel és y(t) a kimenĘjel (szabályozott jellemzĘ).

Soros, a szakasszal sorba beiktatott P (Proporcionális), PI (Proporcionális+Integráló), PD (Proporcionális+Differenciáló) és PID (Proporcionális+Integráló+Differenciáló) jellegĦ szabályozókat alkalmaznak a leggyakrabban, amelyeknek átviteli függvényei: CP ( s )

A

kc

1 Ti s 1 ) kc sTi s sW 1 Td s ) kc CPD ( s ) kc (1 1 T1s 1 T1s 1 Ti s 1 Td s 1 sW ) | kc CPID ( s ) CPI ( s )CPD ( s ) A(1 1 T1s sTi 1 T1s s

C PI ( s )

A(1

A PID tag PI és PD tagok szorzatával megadott közelítése akkor alkalmazható, ha Ti ! W ! T1 . A szabályozó átviteli függvénzét megadhatjuk összeg illetve szorzat alakban. PI szabályozót akkor alkalmazunk, ha állandósult állapotban pontos beállást követelünk meg egységugrás alapjelre. Az integráló hatás fogja biztosítani a hibamentes beállást. A PD szabályozó beiktatása gyorsítja a rendszert. Ideális PD szabályozóban a T1 paraméter értéke zérus. Ez a szabályozó azonban nem realizálható. PID szabályozót akkor alkalmazunk, ha a szabályozás pontosságát és gyorsaságát is növelni kívánjuk. A PD tag a beavatkozójel jelentĘs kezdeti megnövekedését eredményezi, amely a rendszer gyorsításáért felelĘs. A szabályozót a folyamathoz (illetve annak modelljéhez) tervezzük a minĘségi követelmények kielégítésére. Gyakori kompenzálási technika a póluskiejtés, amikor a szabályozó átviteli függvényének zérusaival kiejtjük a szakasz kedvezĘtlen pólusait, és a szabályozó pólusaival kedvezĘbb dinamikát



biztosítunk a szabályozási körben. A póluskiejtéses technika alkalmazásához célszerĦ a szabályozó átviteli függvényét szorzat alakban megadni. Egy PID szabályozóban négy modell paraméter van: kc , TI , TD , T1 . Póluskiejtéses szabályozó tervezés esetén ezeket a következĘképpen választjuk meg: a Ti paramétert a folyamat legnagyobb idĘállandójával megegyezĘnek választjuk (legkisebb frekvenciájú pólus), a Td paramétert pedig a második legnagyobb idĘállandóval tesszük egyenlĘvé. Így elérhetjük, hogy a behozott zérusok kiejtsék a folyamat pólusait. A Td T1 paramétert T1 alakban adjuk meg, ahol n p a pólus áthelyezési arány, amely megadja, hogy a PD np tag hányszoros frekvenciára tolja el a kompenzált pólust. Jó gyakorlati szabály, hogy n p -t az 2-10 tartományba vesszük fel. Ha nagyobbra választjuk, akkor a rendszer gyorsabb lesz, de ennek az ára az, hogy a beavatkozó jel maximuma is nagyobb lesz. Mivel a kc paraméter értéke nem befolyásolja a felnyitott kör fázismenetét, ezért ezt arra használhatjuk, hogy beállítsuk vele az elĘírt fázistartalék értéket. A soros P, PI, PD, PID szabályozó tervezés lépéseit az alábbi folyamat kompenzálására mutatjuk be. 1 P(s) (1 10 s )(1 s )(1 0.2 s ) Tervezzünk soros P, PI, PD, PID típusú szabályozókat 60q fázistöbblet biztosítására. Adjuk meg a minĘségi jellemzĘk értékeit. Számítsuk ki és ábrázoljuk a kimenĘjelet és a beavatkozójelet egységugrás alapjel esetén..

P szabályozó tervezése: A szabályozót C ( s ) kc alakban adjuk meg. Tehát csak a kc paraméter értékét kell meghatározni. Adjuk meg a szakasz átviteli függvényét. » s=tf('s') » P=1/((1+10*s)*(1+s)*(1+0.2*s)) » P=zpk(P) ElsĘ lépésben vegyük fel a kc=1 értéket. » kc=1 » C=kc » L=C*P Ábrázoljuk a felnyitott rendszer Bode diagramját és határozzuk meg annak jellemzĘ paramétereit (fázistartalék, erĘsítési tartalék, vágási körfrekvencia) a margin utasítás segítségével. » margin(L) Látható, hogy a rendszernek jelentĘs fázis- illetve erĘsítési tartaléka van. A fázistolás monoton csökken nullától -270q-ig, tehát az erĘsítés megváltoztatásával be lehet állítani a kívánt fázistartalékot. A M M m 180q -hoz tartozó erĘsítés reciproka lesz a szabályozó erĘsítési tényezĘje. Ezt le lehet olvasni a rendszer Bode diagramjárol, vagy ki lehet olvasni egy táblázatból vagy a margin utasítás segítségével lehet számolni. A margin utasítás az erĘsítési tartalék (gm) értékét a -180q-hoz tartozó erĘsítési értékbĘl számolja. Ha a rendszer fázis értékeit lecsökkentjük a kívánt fázistartalék értékkel, akkor a margin utasítás éppen az adott fázistartalékhoz tartozó erĘsítési érték reciprokát adja vissza gm paraméterként. » [mag,phase,w]=bode(L); » gm=margin(mag,phase-60,w)



» kc=gm Tehát a szabályozó: CP ( s ) kc 7.51 Ezzel egyenértékĦen a kc parameter meghatározható a frekvenciafüggvény adatait tartalmazó táblázatból is az alábbiak szerint: » w=logspace(-1,1,100); » [num,den]=tfdata(L,’v’) » [mag,phase]=bode(num,den,w); » Table=[mag, phase, w’]

mag phase w 0.1527 -115.4479 0.5591 0.1442 -117.3498 0.5857 0.1361 -119.2727 0.6136 0.1283 -121.2166 0.6428 A -120q- hoz tartozó erĘsítés reciproka lesz a kc, a hozzátartozó frekvencia pedig a vágási körfrekvencia: kc=1/0.1361=7.3475, wc=0.6136 A w frekvenciavektor felbontását finomítva mindkét módszer azonos eredményt ad. Fontos, hogy ellenĘrizzük a rendszer viselkedését. » C=kc » L=kc*L EllenĘrizzük a stabilitási tartalékra jellemzĘ paramétereket (erĘsítési tartalék, fázistöbblet). » margin(L) » [gm,pm,wg,wc]=margin(L) A fázistöbblet valóban 60q. A vágási körfrekvencia Z c =0.6245. Számoljuk ki a zárt rendszer eredĘ átviteli függvényét. » T=feedback(L,1) Adjuk meg a felnyitott és a zárt kör Bode diagramját. » bode(L,'r',T,'b') Ábrázoljuk a zárt rendszer átmeneti függvényét. » step(T) Számítsuk is ki az átmeneti függvényt. » t=0:0.05:10; » y=step(T,t); Az átmeneti függvény maximális értéke: » ym=max(y) Az átmeneti függvény állandósult értéke: » ys=dcgain(T) EzekbĘl az értékekbĘl a százalékos túllendülés értéke kiszámítható. » yt=(ym-ys)/ys Az állandósult hiba: » es=1-ys Még meg kell vizsgálni az u(t) beavatkozó jel viselkedését is. Ez azért fontos, mert ez a valóságos folyamat bemenĘjele és itt lépnek fel a fizikai korlátozások. Számoljuk ki a beavatkozójel és az alapjel közötti eredĘ átviteli függvényt. » U=feedback(C,P) vagy » U=C/(1+C*P) Egységugrás alapjel esetén: » u=step(U,t); » plot(t,u) A leggyakrabban korlátozást a beavatkozó jel maximumára adnak meg, számoljuk ezt ki.





» um=max(u)

P,PI,PD,PID szabályozók tervezése: A PI,PD,PID szabályozókat hasonló módon lehet megtervezni. Hozzunk létre egy m file-t a szabályozó tervezésre a szabályozás jellemzĘ paramétereinek kiszámítására. Az alábbi táblázat összefoglalja az eredményeket.

C ( s)

L( s )

Szakasz

C ( S ) P( s )

Zc

kc

yt

es

um

~ ts

0

0

0.5

2

12

7.51

0.62

0.153

0.117

7.51

8

0.504

0.46

0.078

0

5.16

9

21.69

1.98

0.103

0.044

217

2

19.34

1.79

0.076

0

193

2

1

1 10s 1 s 1 0.2s

P

kc

kc

1 10s 1 s 1 0.2s PI

1 10s

kc

s

s 1 s 1 0.2 s

1 s 1 0.1s

1 10 s 1 0.1s 1 0.2s

1 10 s 1 s s 1 0.1s

s 1 0.1s 1 0.2 s

kc

PD

kc

PID

kc

kc

kc

Step response

Control signal, u(t)

8

1.4

PI 6 1.2 PI

P

4

PID 1

2

PD P

0.8

0

0.6

-2

0.4

250

0.2

200

0

2

4

6

8

10

8

10

Control signal,

150 0

0

1

2

3

4

5

6

7

8

9

10

100 50 PID 0 -50

PD 0

2

4

6

A szabályozott rendszerrel szemben támasztott követelmény a gyors beállás és az alapjel követés volt. Látható, hogy a P kompenzáció egyik feltételt sem teljesíti. A beállás elég lassú és nem is éri el a várt y=1 értéket. A PI szabályozás hatására a stacionárius hiba nullára csökkent, de a gyorsaság nem változott. A PD kompenzáció felgyorsítja a rendszert, de csak maradó hibával követi az egységugrást. Ezt a gyorsítást úgy érte el, hogy jelentĘsen megnövekedett az u(t) beavatkozó jel értéke. A PID szabályozó hatására a rendszer gyorsabb is lett és a stacionárius hiba is lecsökkent nullára





MásodrendĦ lengĘ tag kompenzálása Írjunk Matlab programot (hozzunk létre egy m fájlt, ez egy szövegfájl m kiterjesztéssel). Nyissunk meg egy új m-fajlt a Matlab file menĦjébĘl. Írjuk a Matlab utasításokat az üres fájlba. Mentsük el a fájlt: Save As, C:/Matlab/work/myfile.m A program meghívásához egyszerĦen írjuk be a program nevét kiterjesztés nélkül a Matlab parancs ablakába. » myfile

A szakasz az alábbi átviteli függvénnyel adható meg:

Az alábbi Matlab program szolgálja a szabályozó tervezést.

diagram meredeksége 0-ról –40dB/dekádra változik. Ezesetben egy lehetséges PID póluskiejtéses kompenzálási technika, ha mind a PI, mind pedig a PD tag idĘállandóját a természetes körfrekvencia 1 T0 s 1 T0 s . reciprokára választjuk, vagyis Ti T0 és Td T0 , C ( s ) kc s 1 T1s Egy másik lehetĘség egy tiszta integráló tag alkalmazása kompenzáló tagként, amelynek erĘsítését úgy választjuk meg, hogy a fázistöbblet 60º körüli legyen.

clear s=tf('s'); P=1/((1+10*s)*(1+s)*(1+0.2*s )) P=zpk(P) Cp=1 Cpi=(1+10*s)/(10*s) Cpd=(1+s)/(1+0.1*s) Cpid=Cpi*Cpd Cpid=zpk(Cpid) [mag,phase,w]=bode(Cp*P); kp=margin(mag,phase-60,w) Cp=kp*Cp; [mag,phase,w]=bode(Cpi*P); kpi=margin(mag,phase-60,w); Cpi=kpi*Cpi; [mag,phase,w]=bode(Cpd*P); kpd=margin(mag,phase-60,w); Cpd=kpd*Cpd; [mag,phase,w]=bode(Cpid*P); kpid=margin(mag,phase-60,w); Cpid=kpid*Cpid; Tp=feedback(Cp*P,1); Tpi=feedback(Cpi*P,1); Tpd=feedback(Cpd*P,1); Tpid=feedback(Cpid*P,1); Up=feedback(Cp,P); Upi=feedback(Cpi,P); U d f db k(C d P)

figure(1),step(Tp,'r',Tpi,'b' ,Tpd,'g',Tpid,'y') figure(2),step(Up,'r',Upi,'b' ) figure(3),step(Upd,'g',Upid,' y') t=0:0.05:10; yp=step(Tp,t); ypi=step(Tpi,t); ypd=step(Tpd,t); ypid=step(Tpid,t); ysp=dcgain(Tp) yspi=dcgain(Tpi) yspd=dcgain(Tpd) yspid=dcgain(Tpid) ep=1-ysp epi=1-yspi epd=1-yspd epid=1-yspid ytp=(max(yp)-ysp)/ysp ytpi=(max(ypi)-yspi)/yspi ytpd=(max(ypd)-yspd)/yspd ytpid=(max(ypid)-yspid)/yspid up=step(Up,t); upi=step(Upi,t); upd=step(Upd,t); upid=step(Upid,t); upim=max(upi) updm=max(upd) upidm=max(upid)

A szabályozási kör viselkedése vizsgálható Simulink blokk-diagram megépítésével és futtatásával az adott szakasszal és a megtervezett szabályozókkal.

P( s)

A ( s p1 )( s p2 )

A , ahol p1,2 s 2T02 2] T0 s 1

a r jb

A pólusok konjugált komplexek. A Bode diagram törésponti frekvenciája Z 0

1 Itt a közelítĘ Bode T0

Érzékenység vizsgálata: A szakasz paraméterei általában nem pontos értékek, hiszen méréssel, identifikációval határozzuk meg Ęket. Az átmeneti függvénybĘl például közelítĘen meghatározhatjuk a szakasz jellegét (pl. egytárolós holtidĘs tag, kéttárolós tag, stb.) és megbecsülhetjük paramétereit. A paraméterek értéke egy minimális és egy maximális érték között változhat. A szabályozási rendszer viselkedésének vizsgálatakor fontos szempont a paraméter bizonytalanságok hatásának figyelembe vétele. A szabályozó tervezésekor figyelembe kell vennünk a paraméterek értékeinek bizonytalanságát. A szabályozásnak a paraméterek névleges értékeinél az elĘírásoknak megfelelĘen kell viselkednie, de elfogadható viselkedést kell kapnunk, ha a paraméterek a névleges értéküktĘl eltérnek, de még a megadott sávban vannak. KülönbözĘ szabályozó tervezési eljárások vannak a robusztus viselkedés biztosítására. A paraméter bizonytalanságok hatásának szemléltetésére tervezzünk egy PID szabályozót az alábbi lengĘ taghoz, és vizsgáljuk a szabályozási rendszer viselkedését, ha a szakasz csillapítási tényezĘje változik.

P( s)

1 , a 9 csillapítási tényezĘ 0.2 és 1 között változik. 4 s 2 49 s 1

Tervezzük a szabályozót 60º fázistöbbletre. A Matlab utasítások: » clear » s=tf('s'); Adjuk meg a szakaszt: » P1=1/(4*s*s+0.8*s+1) %zeta=0.2 » P2=1/(4*s*s+2.8*s+1) %zeta=0.7, névleges » P3=1/(4*s*s+4*s+1) %zeta=1 » P=zpk(P); P1=zpk(P1); P2=zpk(P2); P3=zpk(P3); A szabályozó kp=1 erĘsítéssel: » Cpid=(1+2*s)*(1+2*s)/((2*s)*(1+0.2*s)) » Cpid=zpk(Cpid) A Bode diagram a névleges P2 szakaszra: » [mag,phase,w]=bode(Cpid*P2); Határozzuk meg a szabályozó átviteli tényezĘjét a 60º fázistöbblet biztosításához: » kp=margin(mag,phase-60,w)



» Cpid=kp*Cpid; » figure(1); bode(Cpid*P2) % a Bode diagram ábrázolása » Tpid=feedback(Cpid*P2,1); % az eredĘ átviteli függvények számítása » Upid=feedback(Cpid,P2); Számítsuk ki és ábrázoljuk a kimenĘjelet és a beavatkozójelet. » t=0:0.05:10; » ypid=step(Tpid,t) ; » figure(2); plot(t,ypid),grid,shg » upid=step(Upid,t); » figure(3); plot(t,upid),grid,shg » Tpid1=feedback(Cpid*P1,1) » Tpid3=feedback(Cpid*P3,1) Számítsuk ki és ábrázoljuk a kimenĘjelet a névleges szakasszal és a megváltozott csillapítási tényezĘkkel. » ypid1=step(Tpid1,t); » ypid3=step(Tpid3,t); » figure(4); plot(t,[ypid,ypid1,ypid3]),grid

Az alábbi ábra mutatja a szabályozás kimenĘjelét a névleges szakaszhoz ( 9 0.7 ) tervezett szabályozóval a névleges szakasszal és a módosított csillapítási tényezĘjĦ szakasszal. Látható, hogy a túllendülés lényegesen nagyobb a kisebb csillapítási tényezĘvel. A szabályozás érzékeny a csillapítási tényezĘ megváltozására. 1.4

HoltidĘs Rendszer Soros Kompenzációja


8. HoltidĘs Rendszer Soros Kompenzációja HoltidĘs tagot tartalmazó rendszer kompenzálása bonyolultabb, mint egy holtidĘ mentes rendszeré, mivel a holtidĘs tag nem reprezentálható pontosan egy racionális törtfüggvénnyel. A holtidĘs tag által létrehozott fázistolást a tervezés során külön kell figyelembe venni. Tekintsük az alábbi szabályozási kört:

r(t) R(s)

-

e(t) E(s)

L( s ) C ( s ) Pd ( s ) a felnyitott kör átviteli függvénye. Vegyük a következĘ példát: Pd (s) P(s)e sTd

e s . 1 20s

A C ( s ) szabályozót úgy kell megválasztani, hogy az elĘírt minĘségi jellemzĘk teljesüljenek. ElĘírások: Egységugrás alapjel esetén a szabályozott rendszer statikus hiba nélkül kövesse az alapjelet és a kimenĘjel túllövése körülbelül 10% legyen. Ezek alapján egy PI kompenzáló tagot választunk: C (s)

0.7

y(t) Y(s)

Pd ( s )

ahol P ( s ) a folyamat átviteli függvénye holtidĘ nélkül, C ( s ) a szabályozó átviteli függvénye és

1.2

kc

0.2

1 20 s s

A zárt rendszer eredĘ átviteli függvénye:

1 1

1 20 s e s e s kc s 1 20 s s A kc konstanst úgy választjuk meg, hogy a fázistartalék 60q legyen. L( s )

0.8

0.6

C (s) P( s)

kc

A folyamat frekvenciafüggvényének amplitúdóját meghatározhatjuk a késleltetés nélküli rendszerbĘl:

0.4

Pd ( jZ )

0.2

0

u(t) U(s)

C ( s)

0

1

2

3

4

5

6

7

8

9

10

P( jZ )e

jZTd

P( jZ ) , mivel e

jZTd

1;

A fázisszög pedig Arg ^Pd ( jZ )` Arg ^P( jZ )` Arg ê

jZTd

`

Arg ^P( jZ )` Z Td

» s=zpk(‘s’) » P=1/(1+20*s) » Td=1 » kc=1 » C=kc*(1+20*s)/s A felnyitott kör átviteli függvénye: » L=C*P » L=minreal(L) A felnyitott kör frekvenciafüggvényének amplitúdója megegyezik a késleltetés nélküli rendszer amplitúdójával, a fázisszöge pedig egy lineáris taggal módosul: » [num,den]= tfdata(L,’v’) » [mag,phase,w]=bode(num,den);



» magd=mag; » phased=phase-w*Td*180/pi; A kc erĘsítést kétféleképpen számolhatjuk: 1. módszer: A margin utasítás segítségével: Az erĘsítési tartalék -120° -os fázisra: » gm=margin(magd,phased-60,w) Ez lesz a kc tényezĘ értéke » kc=gm 0.5229

Save data to workspace Variable name: ty (tu for the control signal) Matrix format Így a t idĘ és az y kimenĘjel vektorai egyszerĦen kinyerhetĘek a szimuláció után. EzekbĘl pedig a minĘségi jellemzĘk meghatározhatók (túllövés, beállási idĘ, maximális beavatkozó jel, stb.). » t=ty(:,1) » y=ty(:,2) » plot(t,y),grid

Pade közelítés: HoltidĘs rendszer kezelése Pade közelítés segítségével is elvégezhetĘ. A késleltetĘ tag sT közelíthetĘ egy racionális törtfügvénnyel, PPADE (s) | e d . (A Pade törtfüggvény Taylor sorának elsĘ néhány tagja megegyezik a holtidĘs tag átviteli függvénye Taylor sorának elsĘ néhány tagjával.)

r(t) R(s)

A rendszer viselkedését a következĘ Simulink modellel lehet vizsgálni. A Simulink lehetĘvé teszi, hogy a holtidĘt egyszerĦen szimuláljuk. 1.4 1.2

Scope1

1 0.8

C L T I System

P L T I System1

0.6

Transport Delay

Scope2

0.4 0.2 0

0

5

10

15

20

25

30

A Scope block segítségével egyszerĦen visszaküldhetjük a szimuláció eredményét a Matlab felületre. Ezeket az adatokat felhasználhatjuk további vizsgálatokra vagy grafikus megjelenítésre. A Scope grafikus ablak paramétereit állítsuk be a következĘképpen: A ‘properties’ menĦ alatt.

-

e(t) E(s)

C (s)

u(t) U(s)

P( s)

PPADE ( s )

y(t) Y(s)

Matlab-ban a pade utasítás számítja ki a közelítést a kívánt fokszámra. Például használjunk 5-ödrendĦ közelítést. » [numpade,denpade]=pade(Td,5) » Ppade=tf(numpade,denpade) » Pd=P*Ppade A kc erĘsítés ezzel már meghatározható. » [mag,phase,w]=bode(Pd); » kc= margin(magd,phased-60,w) 0.5229 Ki kell hangsúlyozni, hogy az elsĘ részben megadott módszer pontosabb, mint a Pade közelítés alkalmazása és a kapott felnyitott és zárt átviteli függvények jóval bonyolultabbak a magasabb fokú közelítés miatt. A Pade közelítés elĘnye, hogy a tervezés menete nagyon hasonló a késleltetés nélküli rendszerek tervezési módszeréhez.

<=

A mag értéke phased 120q -nál 1.8738 . EbbĘl kc=1/1.8738 » kc=1/1.8738 0.5337 Számoljuk újra a szabályozót: » C=kc*(1+20*s)/s EllenĘrizzük, hogy a fázistartalék valóban 60° lett-e.

Step


Data history:

2. módszer: Táblázat segítségével: » w=logspace(-2,1,100); » [mag,phase]=bode(num,den,w); » magd=mag » phased=phase-w’*Td*180/pi Hozzunk létre egy táblázatot » T=[ phased, magd, w’]

-116.5943 2.1544 0.4642 -118.5162 2.0092 0.4977 -120.5770 1.8738 0.5337 -122.7868 1.7475 0.5722 -125.1562 1.6298 0.6136



Bode Diagrams


P(s)

From: U(1) -20

1 ( s 2)( s 5)

-40 -60 -80

To: Y(1)

Stabilizálható-e a folyamat egy arányos (P) -100 C ( s ) kc szabályozóval? 180 » s=zpk('s') 170 » P=1/((s+2)*(s-5)) 160 » figure(1); nyquist(P); 150 » figure(2); 10 10 10 10 » figure(3);rlocus(P); Frequency (rad/sec) Látható, hogy a stabilitást nem lehet elérni, mivel a Nyquist diagram nem veheti körbe a (-1+j0) pontot az óramutató járásával ellentétes irányban, és a fázistartalék is mindig kisebb mint nulla. Ugyanez az eredmény kapható meg a karakterisztikus egyenlet és a pólusok alapján. A pólusokat az s2-3s-10+kc=0 karakterisztikus egyenletbĘl határozhatjuk meg. A stabilitás szükséges feltétele, hogy az együtthatóknak azonos elĘjelĦeknek kell lenniük. Ez nem teljesülhet, mert a kc értéke nem befolyásolja a -3 együtthatót. A gyökhelygörbe (0 d kc < f) alapján is ugyanezt kapjuk. -1

0

1

2

2. Módszer: A maximális értéket meghatározhatjuk a max utasítással is. » [maxphase,index]=max(phase) » kc=1/mag(index) Számoljuk újra a szabályozót: » C=kc; » L=C*P; A margin utasítással grafikusan is ellenĘrizhetjük a fázistartalékot. » margin(L); A rendszernek elég kicsi a fázistartaléka, Pm 25.4q (60º-ot szeretnénk elérni). Nézzük meg a zárt rendszer viselkedését: » T=feedback(L,1) » step(T) A stabilitást sikerült elérni, de a viselkedés nem megfelelĘ, a rendszernek közel 100%-os statikus hibája van. PID szabályozó alkalmazásával a szabályozás minĘségi viselkedése javítható.

3. Példa. Tervezzünk PID szabályozót a statikus hiba lecsökkentése érdekében. C (s)

2. Példa. Vizsgáljuk meg most a következĘ folyamatot: P(s)


A Bode fázisgörbe a M = -154.62 mag=0.045 Z =3.16 -nél éri el a maximumát. Az arányos szabályozóval elérhetĘ maximális fázistartalék tehát pm=180-154.6=25.8°. Ebben az esetben az erĘsítést kc=1/0.045=22.22 értékre kell választani. » kc=22.22

9. Labilis Rendszer Soros Kompenzációja 1. Példa. Vizsgáljuk a következĘ labilis folyamatot:

Labilis Rendszer Soros Kompenzációja

1 ( s 2)( s 5)

To: Y(1)


Stabilizálható -e a folyamat egy arányos (P) C ( s) kc szabályozóval? » clear » s=zpk('s') » P=1/((s-2)*(s+5)) » figure(1); nyquist(P); » figure(2); bode(P); » figure(3); rlocus(P); Ebben az esetben már stabilizálni lehet a folyamatot. Az erĘsítés növelésével teljesíthetĘ az általános Nyquist kritérium, kc >10 esetben a zárt rendszer stabilis lesz. A Bode diagramról is leolvasható, hogy a fázistartalék pozitív lehet. Válasszuk meg a kc értékét úgy, hogy a vágási körfrekvencia a fázisszög maximumára essen. ElĘször vizsgáljuk a felnyitott kört kc=1 esetén. Bode Diagrams » C=1 From: U(1) -20 » L=C*P -40 » [num,den]=tfdata(L,’v’) -60 » [mag,phase,w]=bode(num,den); -80 1. Módszer: Táblázatból olvassuk ki az értékeket. -100 » T=[ phase, mag, w] -150 -155.9825 0.0613 2.2122 -160 -154.7797 0.0506 2.8072 -170 -154.6231 0.0452 3.1623 -154.7797 0.0399 3.5622 -180 10 10 10 10 -155.9825 0.0300 4.5204 -1

0

Frequency (rad/sec)

1

2

kc

( s 2) ( s 5) . s s 50

A labilis p1 2 pólust nem ejthetjük ki közvetlenül egy labilis zérussal, mivel a paramétereket rendszerint mérési eredmények alapján határozzuk meg, és a pólus és az Ęt kiejtĘ zérus kis eltérése esetén is már a rendszer labilissá válik. A labilis pólust ehelyett egy stabilis PI taggal kompenzáljuk. A rendszer gyorsítása érdekében a stabilis p1 5 polust nagyobb frekvenciára toljuk el egy PD taggal (p=-50, a póluseltolási arány 10). A kc konstanst ismét úgy választjuk meg, hogy a maximális fázistartalékot kapjuk. Bode Diagrams Zárjuk be a grafikus ablakokat. Gm=-17.25 dB (at 2.0851 rad/sec), Pm=58.122 deg. (at 14.596 rad/sec) » clear 100 » s=zpk('s') 50 » P=1/((s-2)*(s+5)) 0 » C=((s+2)*(s+5))/(s*(s+50)) -50 » L=C*P -100 -100 » L=minreal(L) » bode(L) -150 » [mag,phase,w]=bode(L); -200 Határozzuk meg az erĘsítést a maximális fázishoz: -250 » [maxphase,index]=max(phase) -300 10 10 10 10 10 10 » kc=1/mag(index) Frequency (rad/sec) Az erĘsítési tényezĘ kc=760.26, a fázistöbblet pedig pm=180+maxphase =59. EllenĘrizzük a rendszer viselkedését. » kc=760.26 » C=kc*((s+2)*(s+5))/(s*(s+50)) » L=C*P, L=minreal(L) » margin(L) A zárt rendszer átmeneti függvényét az alábbi ábra mutatja. » T=feedback(L,1) Phase (deg); Magnitude (dB)


-2

-1

0

1

2

3



A statikus hiba lecsökkent nullára. » es=1-dcgain(T) EllenĘrizzük az u beavatkozó jelet: » U=feedback(C,P); » step(U) » u=step(U) » um=max(u)

Mintavételes Rendszerek


10. Mintavételes Rendszerek z-transzformáció:

A beavatkozó jel um maximuma elég nagy. Ez az érték a póluseltolási arány csökkentésével mérsékelhetĘ.

Digitális szabályozási rendszerekben a jeleket digitalis alakban tároljuk és dolgozzuk fel. Ezekben a rendszerekben mĦködési elvük következtében a jelfeldolgozás csak diszkrét mintavételi idĘpontokban történik. EbbĘl adódóan a szabályozási körökben megjelenĘ analóg jeleket mintavételezni és digitalizálni kell. Egy folytonos y (t ) jel digitális alakját eltólt impulzusok összegével lehet leírni.

yd (t )

^ y (nTs )`

y (0)G (t ) y (Ts )G (t Ts ) y (2Ts )G (t 2Ts ) y (3Ts )G (t 3Ts ) ...

Az yd (t ) jel Laplace transzformáltja:

Step Response From: U(1)

yd ( s )

1.5

y (0) y (Ts )e sTs y (2Ts )e 2 sTs y (3Ts )e3 sTs ...

Vezessük be a következĘ jelölést:

To: Y(1)

Amplitude

z

e sTs

e sTs , z 1

1

yd ( z )

y (0) y (Ts ) z 1 y (2Ts ) z 2 y (3Ts ) z 3 ...

f

¦ y(kT ) z

k

s

k 0

Ez az átalakítás a mintavételezett jel z-transzformációja, ahol a z 1 -et az eltolási operátorként lehet értelmezni.

0.5

0 0

0.2

0.4

0.6

Time (sec )

0.8

1

1.2

1. Példa. Egy jel z-transzformáltja

2z2 z , Ts 0.5 z z 0.24 Számoljuk ki a jel inverz z-transzformációjának elsĘ N s 5 mintavételi pontban felvett értékét. yd

2

a/ Matlab utasítások: » num=[2, -1, 0] » den=[1, -1, 0.24] » Ns=5 » yd=dimpulse(num,den,Ns); » plot(1:Ns,yd,'*') Szimbólikus adatbevitellel és LTI (sys) struktúrával is elvégezhetĘ a számítás. Hasonlóan az s változóhoz egy z változó is definiálható, de itt a Ts mintavételi idĘt is meg kell adni. » Ts=0.5 » z=tf(‘z’,Ts) » Y=(2*z^2-z)/(z*z-z+0.24) Az LTI (sys) struktúra használatának egyik elĘnye, hogy ugyanazt az utasítást lehet alkalmazni folytonos és diszkrét rendszerekre. Ha a sys.Ts mintavételi idĘ (sys.Ts) paraméter nulla, akkor a Matlab a rendszert folytonosnak, ha pedig nullától eltérĘ, akkor diszkrétnek feltételezi. A help utasítás megadja az impulse utasítás pontos használatát. » help impulse IMPULSE(SYS,TFINAL) simulates the impulse response from t=0 to the final time t=TFINAL. A tfinal 2 idĘ szükséges Ts 0.5 mintavételi idĘ esetén ahhoz, hogy 5 mintavételi pontot kapjunk (Az elsĘ érték t=0-hoz tartozik) » tfinal=2; » Yd=impulse(Y,tfinal); » plot(Yd,’*’); 1





Határozzuk meg az alábbi folytonos folyamat impulzusátviteli függvényét. A mintavételezési idĘ Ts

1 5s , (1 10 s)(1 8s)(1 4 s )(1 2s )

P( s) b/ Részlettörtekre való bontás Egy exponenciális tag z-transzformáltja: Z o y (t ) e at , yd (nTs ) e anTs

num den

z

A folytonos folyamat átviteli függvénye zerus-pólus alakban:

P( s)

§ r r · z ¨ k0 1 2 ¸ z p1 z p2 ¹ ©

num1 den

s 0.2 ( s 0.5)( s 0.25)( s 0.125)( s 0.1)

A Matlab-ban egy folytonos folyamat diszkrét imulzusátviteli függvénye a c2d utasítással számolható: sysd = c2d(sysc,ts,method), ahol ts a mintavételi idĘ és method pedig a tartószerv típusát adja meg. A default értéke ‘zoh’, ami nulladrendĦ tartószervet jelent. » Ts=1; » Pd=c2d(P,Ts,’zoh’)

Pd ( z )

k0 =

§ r r · 1 · z z § 1 z¨ 1 2 ¸ z¨ ¸ z p z p z 0.6 z 0.4 z 0.6 z 0.4 © ¹ 1 2 ¹ © anTs bnTs aTs bTs yd (t ) ^ y (nTs )` e e , where e 0.6, e 0.4

yd ( z )

^

`

» a=log(0.6)/Ts » b=log(0.4)/Ts

1.02, b

ln(0.4) Ts

y (t ) e 1.02 t e1.83 t ,

0.0011

( z 3.088)( z 0.8187)( z 0.2198) ( z 0.9048)( z 0.8825)( z 0.7788)( z 0.6065)

A diszkrét zérusok és pólusok: » [zd,pd,kd]= zpkdata(Pd,’v’) A diszkrét pólusok: 0.9048, 0.8825, 0.7788, 0.6065 A diszkrét zérusok: -3.088, 0.8187, -0.2198 sT Vegyük észre, hogy a diszkrét pólusokat a z e s transzformációval kapjuk a folytonos pólusokból. » exp(-0.5) 0.6065 Hasonlóan » exp(p) Ez a transzformáció a folytonos zérusok és a diszkrét zérusok közt nem érvényes, az összefüggés ennél bonyolultabb.

[] » Ts=0.5

ln(0.6) Ts

0.0078

» [z,p,k]=zpkdata(P,’v’) A folytonos rendszer pólusai: -0.5, -0.25, -0.125, -0.1 A folytonos rendszer zérusa: -0.2

» num1=[2, -1] » den=[1, -1, 0.24] » [r,p,k0]=residue(num1,den) r = 1 1 p = 0.6000 0.4000

a

1

» s=zpk(‘s’) » P=(1+5*s)/((1+10*s)*(1+8*s)*(1+4*s)*(1+2*s))

z z e aTs

A jelet részlettörtekre bonthatjuk a Matlab residue utasításával.

yd ( z )

Ts

1.83

tt0

Végérték tételek Egy diszkrét jel kezdeti értékét meghatározhatjuk a következĘ összefüggéssel: lim y (nTs ) lim y ( z ) ,

» t=[0:Ts:2]’ » yd=exp(a*t)+exp(b*t) » plot(t,yd,’*’);

n o0

z of

végértékét pedig az alábbi összefüggéssel adhatjuk meg: lim y (nTs ) lim(1 z 1 ) y ( z ) . n of

Impulzusátviteli függvény (Diszkrét átviteli függvény): Egy folytonos folyamatot a D/A átalakítóval együtt az impulzusátviteli fügvénnyel irhatjuk le. Az impulzusátviteli függvény a kimenĘjel és a bemenĘjel z-transzformáltjainak hányadosát adja meg nulladrendĦ tartószervet feltételezve a bemeneten.

z o1

3. Példa. Az elĘzĘ példában megadott rendszerre egységugrás bemenetet alkalmazunk. Az egységugrás diszkrét megfelelĘje egy egységnyi amplitúdójú impulzus sorozat, amelynek a z-transzformáltja

u[k] U(z)

D / A u(t) ZOH U(s)

P (s )

1.

y(t)

y[k]

Y(s)

Y(z)

{

U(z)

Pd (z)

Y(z) U(z)

Y(z)

Határozzuk meg a kimenĘjel kezdeti- és végértékét.

y( z) 2. Példa.

2

z Pd ( z ) z 1

3

z . z 1



A kezdeti érték nulla, mivel a számláló fokszáma magasabb, mint a nevezĘ fokszáma. A végértéket meghatározhatjuk a Pd ( z ) impulzusátviteli függvénybĘl a z=1 helyettesítéssel. A Matlab-ban a végértéket a dcgain utasítással számolhatjuk folytonos és diszkrét esetre is (diszkrét esetre a ddcgain(numd,dend) utasítás szintén használható). Hasonlóan itt is a mintavételi idĘ paraméter határozza meg, hogy a rendszer folytonos vagy diszkrét. » P.Ts » A=dcgain(P) » Ad=dcgain(Pd) Mindkét, A (folytonos) és Ad (diszkrét) értékre 1 értéket kaptunk, ahogy ezt vártuk.

Stabilitás Egy diszkrét rendszer stabilis, ha a pólusai (karakterisztikus egyenletének, azaz a nevezĘjének a gyökei) a complex sík egységsugarú körének a belsejébe esnek.

Mintavételes PID Szabályozó Tervezése Kisfrekvenciás Közelítés Alapján


11. Mintavételes PID Szabályozó Tervezése Kisfrekvenciás Közelítés Alapján Mintavételes szabályozási körben a szabályozó tervezhetĘ a frekvenciatartományban annak szem elĘtt tartásával, hogy a mintavételezés és zérusrendĦ tartószerv alkalmazása úgy tekinthetĘ, mintha járulékos holtidĘ lépne fel a rendszerben. Egy egytárolós arányos tagot szinuszos jellel gerjesztünk. A tag kimenĘjelét 1 sec mintavételi idĘvel mintavételezzük, majd zérusrendĦ tartószervet alkalmazunk. Ábrázoljuk a bemenĘjelet, valamint a folytonos és a mintavételezett szakasz kimenĘjelét egy diagramban az alábbi ábra szerint.. 0.8

A folytonos kimenĘjel kvázistacionárius állapotban a bemenĘjelhez képest amplitudóban és fázisban eltér, de a bemenĘjellel megegyezĘ frekvenciájú szinuszos jel. Szemmel láthatóan a mintavételezett jel késik a folytonos kimenĘjelhez képest, és annak alapharmonikusa kb. fél mintavételnyi idejĦ holtidĘs késleltetést mutat a folytonos kimenĘjelhez képest. Elemezzük a jelenséget a frekvenciatartományban. Hasonlítsuk össze egy folytonos egytárolós arányos tag frekvenciafüggvényét a mintavételezĘvel és zérusrendĦ tartószervvel ellátott egytárolós arányos tag frekvenciafüggvényével. Ts Ts A A zérusrendû tartószerv 1 sT1 1 sT1

4. Példa.

0.6

Stabilis-e az alábbi diszkrét rendszer?

Pd (z)

0.4

z 2 - 0.3 z - 0.1 , Ts z 3 + 3 z 2 + 2.5 z + 1

0.2

1

0

-0.2

Pole-zero map

-0.4 -0.6 -0.8

4

To: Y(1)

Amplitude

Imag Axis

1 Definiáljuk a diszkrét rendszert: 0.8 » Ts=1 0.6 » z=tf(’z’,Ts) 0.4 » Pd=(z*z-0.3*z0.2 1)/(z^3+3*z^2+2.5*z+1) 0 Határozzuk meg a pólusokat: -0.2 -0.4 » [z,p,k]=zpkdata(Pd,’v’) -0.6 A pólusok abszolút értékei: -0.8 » abs(p) -1 -2.5 -2 ans = 2.0000 0.7071 0.7071 A rendszer labilis, mert egyik pólusának abszolút értéke nagyobb mint 1. A stabilitás meghatározható a pólusok complex síkon történĘ ábrázolásával. » pzmap(Pd) » zgrid A rendszer labilis, mert a -2 pólus az egységsugarú körön kívül helyezkedik el. A zgrid utasítás kirajzolja az egységsugarú kört és a konstans csillapítási tényezĘhöz és természetes frekvenciához tartozó vonalakat. A stabilitás az idĘtartományban is meghatározható az átmeneti függvénybĘl. » step(Pd) A kimenet nem korlátos, tehát a rendszer labilis.

-1.5

-1

-0.5

0

0.5

1

Real Axis

0

2

4

6

8

10

12

14

16

18

20

Folytonos rendszer

A mintavételes rendszer impulzusátviteli függvénye:

Step Response From: U(1) 20

P( z )

15 10

1 e Ts / T1 z e Ts / T1

A frekvenciafüggvényt z Taylor sorukkal.

5 0

e jZTs helyettesítéssel nyerjük. Az exponenciális tagokat közelítsük 2

-5 -10

Mintavételes rendszer

0

1

2

Time (sec )

3

4

5

P( z

e jZTs )

1 e Ts / T1 e jZTs e Ts / T1

Ha Ts
Ts 1 § Ts · ¨ ¸ ...) T1 2 © T1 ¹ | 2 T 1§T · (ZTs ) 2 1 jZTs ... (1 s ¨ s ¸ ...) 2 T1 2 © T1 ¹ 1 (1

Ts/T1 valamint XTs magasabb hatványai elhanyagolhatók, és a frekvenciafüggvénye közelítĘen megegyezik a folytonos szakasz

1




Ts / T1 1 jZTs Ts / T1 1 jZT1 Az elhanyagolásokat egy további T j járulékos holtidĘvel vehetjük figyelembe, amelynek értéke P( z

e jZTs ) |

Ts/2 és Ts között becsülhetĘ. P ( z

e jZTs ) |

r[k]

e[k]

C (z )

-

1 jZ T e j. 1 jZT1

1. Példa Hasonlítsuk össze a folytonos rendszer és a mintavételes rendszer frekvenciafüggvényeit. Legyen A=1, T1=0.1 és a mintavételezési idĘ Ts=0.1. » numPs=1; » denPs=[0.1 1]; » Ts=0.1; » [numPz,denPz]=c2dm(numPs,denPs,Ts,'zoh') » omega=logspace(-1,2,200); » [Mags,Phases]=bode(numPs,denPs,omega); » [Magz,Phasez]=dbode(numPz,denPz,Ts,omega); » loglog(omega,Mags,'r',omega,Magz,'b'); » title('folytonos (piros), és z-transzformációs (kék) amplitudó') » pause » semilogx(omega,Phases,'r',omega,Phasez,'b'),grid » title('folytonos (piros) és z-transzformációs (kék) fázisszög') » pause » PhasesWithDel=Phases-omega'*(Ts/2)*180/pi; » semilogx(omega,Phases,'r',omega,Phasez,'b',omega, PhasesWithDel,'g') » title('folytonos (piros), z-transzf. (kék), jár.holtidĘs(zöld) fázis') Látható, hogy a z-transzformációnak megfelelĘ (kék) amplitudó-frekvencia görbék X=1/Ts=10-ig jól követik a folytonos rendszer (piros) amplitudó-frekvencia görbéjét. A nagyobb frekvenciák tartományában a közelítés nem elfogadható, az eltérések igen nagyok. A fázisszögben az eltérés a folytonos (piros) és a mintavételes (kék) rendszer között már a kisfrekvenciás tartományban is észrevehetĘ, X=1/Ts=10-nél már kb. 0.5 radián, a járulékos holtidĘs fázisszöggel kiegészített folytonos fázisszög görbe (zöld) azonban jó közelítést ad. A mintavételezésbĘl eredĘ járulékos holtidĘ (bármilyen kicsi is legyen az) megváltoztatja az eredeti rendszer struktúrális tulajdonságait, pl. egy struktúrálisan stabilis rendszer mintavételes szabályozás megvalósításakor nem lesz struktúrálisan stabilis.

u[k] U(z)

D/ A


u(t)

y(t)

P(s )

Y(s)

y[k] Y(z)

A/ D P( z )

Y ( z) U ( z)

A fenti ábrán C(z) a tervezendĘ diszkrét idejĦ szabályozó, P(s) az adott folytonos szakasz. A szabályozással szemben minĘségi követelményeket támasztunk, amelyek egyrészt a statikus alapjelkövetési illetve zavarelhárítási tulajdonságokat, másrészt a rendszer dinamikus tulajdonságait írják elĘ. Tervezhetünk PID jellegĦ folytonos szabályozót a folytonos szakaszhoz, figyelembe véve, hogy a mintavételezés járulékos holtidĘt jelent. Ezután a folytonos szabályozót diszkrét algoritmussá transzformáljuk. CélszerĦ azonban közvetlenül a D/A nulladrendĦ tartószerv és a P(s) folytonos szakasz együttes mintavételes alakjából, a szakasz P(z) impulzusátviteli függvényébĘl kiindulva közvetlenül diszkrét PID jellegĦ szabályozó algoritmust tervezni. A tervezés póluskiejtéses technikán alapulhat. A diszkrét felnyitott rendszer Bode diagramját úgy módosítjuk, hogy kiejtjük a rendszer kedvezĘtlen pólusait, és helyettük megfelelĘ pólusokat hozunk be. Az eljárás lényege a következĘ: A tárolós jellegĦ szakaszok impulzusátviteli függvényének nevezĘje ( z e Ts / T1 )( z e Ts / T2 )... alakú tényezĘket tartalmaz. A diszkrét P, PI, PD, PID jellegĦ algoritmusok az alábbi impulzusátviteli függvényekkel adhatók meg: P szabályozó:

C(z)=A

PI szabályozó:

C ( z)

A

z e Ts / T1 z 1

(KiejthetĘ a szakasz legnagyobb idĘállandója, helyette integráló hatást hozunk be.) Differenciaegyenlet: u[k ] Ae[ k ] A exp(Ts / T1 )e[k 1] u[ k 1] , ahol u[k] a beavatkozójel, e[k] pedig a hibajel aktuális értékét jelöli. PD szabályozó:

C ( z)

A

z e Ts / T2 *

z e Ts / T2

, ahol T2* T2 . (KiejthetĘ a szakasz egy kedvezĘtlen

idĘállandója, amely helyett egy kisebb idĘállandót hozunk be.) Ideális PD szabályozó:

Mintavételes PID Szabályozó Tervezése A mintavételes zárt szabályozási rendszer alapvetĘ struktúráját az alábbi blokkvázlat adja meg. Ez a rendszer hibrid rendszer abban az értelemben, hogy folytonos és diszkrét idejĦ jelek egyaránt szerepelnek benne.

2

z e Ts / T2 z Differenciaegyenlete: u[k ] C ( z)

A

Ae[k ] A exp(Ts / T1 )e[k 1] . (A folytonos PD algoritmustól eltérĘen diszkrétben az ideális PD hatás realizálható, mivel nem eredményez végtelen túlvezérlést.) 3


PID szabályozó:

C ( z)

A


( z e Ts / T1 )( z e Ts / T2 ) ( z 1) z

Differenciaegyenlete: u[ k ]


Ae[ k ] A(exp( Ts / T1 ) exp( Ts / T2 ))e[ k 1] A(exp( Ts / T1 ) exp(Ts / T2 ))e[ k 2] u[ k 1]

Megjegyzés: Többszörös PI, illetve PD hatás is alkalmazható, amennyiben az u beavatkozójelre megadott korlátot nem lépjük túl. A kompenzáló algoritmusoknak megfelelĘ differenciaegyenlet rekurzív összefüggés, amely valós idĘben realizálható.

P s

A szabályozók zérusaival kiejtjük a szakasz nemkívánatos pólusait. A szabályozó A átviteli tényezĘjét úgy választjuk meg, hogy a diszkrét rendszer Bode diagramja alapján a fázistöbblet az elĘírt (rendszerint ~60q) legyen. A Bode diagram a dbode (ill. a rendszer LTI megadásakor a bode Matlab utasítással számítható. A frekvenciatartomány kb. az 1/Ts értékig veendĘ fel. A vágási körfrekvencia nagyjából az Z c | 1/ 2(TholtidĘ Ts ) értékre fog adódni. (A szakasz mintavételezésébĘl Ts/2, a szabályozó algoritmusból további Ts/2 járulékos holtidĘ vehetĘ figyelembe.) A frekvenciatartománynak tehát csupán a kisfrekvenciás része érdekes a szabályozás tervezése szempontjából, ahol a tárgyalt közelítés érvényes. Mivel a fenti szabályozási algoritmusok a szakasz zérusait nem kompenzálják, a mintavételi pontok között nem várhatók lengések. A beavatkozójelre rendszerint korlátok érvényesek. EllenĘrizni kell, hogy a beavatkozójel a korlátokon belül marad-e.

e s 1 10s 1 5s

P1 ( s )e s

A rendszerben egy holtidĘs tag is szerepel 1 késleltetéssel. A mintavételi idĘ Ts=1 . Tervezzünk diszkrét soros szabályozót, amely teljesíti a következĘ minĘségi elĘírásokat. - Fázistartalék#60q - A beállási idĘ minél kisebb legyen - Egységugrás alapjel esetén a statikus hiba legyen nulla. ElĘször vegyük fel a folytonos rendszert késleltetés nélkül: » s=zpk(’s’); » P1s=1/((1+10*s)*(1+5*s)); Határozzuk meg a diszkrét megfelelĘjét nulladrendĦ tartószerv alkalmazásával.

2. Példa Vizsgáljuk az egyes szabályozók egységugrás válaszát. »T1=10; »T2=4; »Ts=1; »A=2; »numPI=A*[1 -exp(-Ts/T1)]; »denPI=[1 -1]; »y=dstep(numPI,denPI,20); »stairs(y),grid, title('PI szabályozó'); »pause »numPD=A*[1 -exp(-Ts/T2)]; »denPD=[1 -exp(-Ts/(T2/5))]; »y=dstep(numPD,denPD,20); »stairs(y),grid, title('PD szabályozó'); »pause »numPDid=A*[1 -exp(-Ts/T2)];»denPDid=[1 0]; »y=dstep(numPDid,denPDid,20);»stairs(y),grid, title('ideális PD szabályozó'); »pause »numPID=A*[1 -exp(-Ts/T1)-exp(-Ts/T2) exp(-Ts/T1)*exp(Ts/T2)]; »denPID=[1 -1 0]; »y=dstep(numPID,denPID,20);»stairs(y),grid, title('PID szab'); »pause


P ( s) ½ P1 z (1 z 1 ) Z ® 1 ¾ ¯ s ¿ » Ts=1; » P1z=c2d(P1s,Ts,’zoh’); (A default ‘zoh’ method megadása nem kötelezĘ). A késleltetést hozzáadjuk a rendszerhez, ha a P1 ( z ) impulzusátviteli függvényt megszorozzuk

z 1 -gyel, mivel

Z ê s `

z 1 .

A z változót definiáljuk hasonlóan mint a folytonos rendszerekre az s változót. Itt a mintavételi idĘt is meg kell adni. » z=zpk(‘z’,Ts) » Pz=P1z/z

P z

z 1 P1 z

z 0.9048 1 P1 z 0.00905 z z 0.9048 z 0.8187 z

Határozzuk meg a diszkrét folyamat zérusait és pólusait: » [zd,pd,kd]=zpkdata(Pz,’v’) PI kompenzáló tag szükséges a nulla statikus hiba eléréséhez, a PD kompenzáló tagot pedig a rendszer gyorsítása érdekében alkalmazzuk. A PI és PD tagok törési frekvenciáját a folytonos folyamat idĘállandói (pólusai) alapján választjuk meg: A Ti integrálási idĘállandót a folyamat legnagyobb idĘállandójával megegyezĘ értékre választjuk, a Td differenciálási idĘállandót pedig a második legnagyobb idĘállandóval tesszük egyenlĘvé. Ezen sarokfrekvenciák diszkrét megfelelĘit a z számíthatjuk:

PI sarokfrekvencia: 0.1 e PD sarokfrekvencia: 0.2 e

Ts Ti

T s Td

e

e

1 10

1 0.5

e sTs transzformáció alapján e 0.1 = 0.9048 e 0.2 = 0.8187

A diszkrét szabályozó tehát:

C ( z)

kc

z 0.9048 z 0.8187 z 1 z

A kc paramétert úgy választjuk meg, hogy körülbelül 60° fázistartalékot érjünk el.

3. Példa A folytonos rendszer adott az átviteli függvényével:

4

5



ElĘször tételezzük fel, hogy kc =1 : » kc=1; » Cz=((z-0.9048)*(z-0.8187))/((z-1)*z) vagy közvetlenül a diszkrét átviteli függvény pólusaiból: » Cz=((z-pd(1))*(z- pd(2)))/((z-1)*z) A felnyitott kör diszkrét átviteli függvénye: L( z ) C ( z ) P ( z ) . » Lz=Cz*Pz; % ejtsük ki a azonos zérus-pólus párt » Lz=minreal(Lz, 0.001); A kc paramétert a margin utasítással számolhatjuk vagy kiolvashatjuk egy táblázatból.

a. A margin utasítás használata. A bode utasítás a Ts mintavételi idĘ alapján dönti el, hogy a folyamat folytonos vagy diszkrét. » [mag,phase,w]=bode(Lz); » kc=margin(mag,phase-60,w); b. Táblázat alkalmazásával » [numz,denz]=tfdata(Lz,’v’); » [mag,phase,w]=dbode(numz,denz,Ts); » T=[mag, phase, w] % táblázat képzése T= 0.0781 -114.8936 0.2200 0.0697 -117.8611 0.2462 0.0622 -121.1822 0.2756 0.0555 -124.8991 0.3084 0.0495 -129.0589 0.3452 Látható, hogy körülbelül -120q fázisszöget 0.062 erĘsítésnél érhetünkk el, ezért a kc -t 1/0.0622=16-re kell választani. » kc=1/0.0622 Megjegyezzük, hogy a margin utasítás kétféleképpen hívható meg. A bemenetei lehetnek a Bode amplitúdó, fazis és körfrekvencia értékei vagy az átviteli függvény LTI (sys) alakban. EllenĘrizzük a viselkedést. Számoljuk ki a fázistartalék értékét a margin utasítással. » Cz=kc*((z-pd(1))*(z- pd(2)))/((z-1)*z); » Lz=Cz*Pz; » Lz=minreal(Lz); » margin(Lz); Használjuk a Simulink-et a zárt rendszer viselkedésének vizsgálatára. » simulink % elindítja a Simulink-et



Data history:

Save data to workspace Variable name: My for ScopeY and Mu for ScopeU Matrix format Változtassuk meg a Simulation–>Parameters–>Stop Time parametert 20-ra.

ScopeU

Cz Step

P1s

LTI System Zero-Order LTI System1 Hold

Transport Delay

Így a t idĘ és az y kimenĘjel vektorai egyszerĦen kinyerhetĘek a szimuláció után. EzekbĘl pedig a minĘségi jellemzĘk meghatározhatók (túllövés, beállási idĘ, maximális beavatkozó jel, stb.). » ty=My(:,1), y=My(:,2) » tu=Mu(:,1), u=Mu(:,2) Ábrázoljuk az y(t) kimenĘjelet » subplot(211), plot(ty,y), grid és az u(t) beavatkozó jelet, amely a nulladrendĦ tartószerv kimenete. » subplot(212), stairs(tu,u) , grid A mintavételezett u[k] jelet is ki lehet rajzolni. » hold on, plot(tu,u,’*’)

1.5

1

0.5

0

0

2

4

6

8

10

12

14

16

18

20

0

2

4

6

8

10

12

14

16

18

20

20 15 10

A rendszer viselkedését a következĘ Simulink modellel lehet vizsgálni. A Simulink lehetĘvé teszi, hogy a folytonos folyamatot, a diszkrét szabályozót és a holtidĘt egyszerre szimuláljuk. Hozzunk létre egy új fájlt és másoljuk be az egyes blokkokat. A blokkok paramétereit változtassuk meg a kívánt értékre: - C(z), P1(s): Control System Toolbox –>LTI system :Cz, P1s - Sum: Simulink–>Math–>Sum: +- Dead time, delay: Simulink–>Continuous–>Transport Delay: 1 - Step input: Simulink–>Sources–>Step, Változtassuk meg a Step time paramétert nullára - Zero-Order-Hold: Simulink–>Discrete–> Zero-Order-Hold, Sampling time: Ts - Scope: Simulink–>Sinks–>Scope A Scope blokk segítségével egyszerĦen visszaküldhetjük a szimuláció eredményét a Matlab felületre. Változtassuk meg a paramétereket a ‘properties’ menĦ alatt a Scope grafikus ablakában. 6

ScopeY

5 0 -5

7

Mintavételes Szabályozó Tervezése w-Transzformációval


12. Mintavételes Szabályozó Tervezése Tustin-transzformációval A Tustin-transzformáción (w- transzformáción) alapuló mintavételes szabályozótervezés motivációja az, hogy a folytonos rendszerekre kialakított tervezési módszereket a mintavételes környezetben is alkalmazni tudjuk. A módszer alapja egy gyakran használt numerikus integrálási módszer, a trapéz szabály szerinti integrálás alkalmazása. Egy folytonos függvény menetét két mintavételi pont között egy egyenes szakasszal helyettesítve a görbe alatti területet egy trapézzal közelíthetjük: kTs

Ts

³ f (t )dt # 2 > f ( kT ) f ( kT T )@ s

s

s

( k 1)*Ts

Ez a közelítés az kTs

³ f (t )dt

i( kTs )

0

integrál növekményére a kTs

Ts

³ f (t )dt # 2 > f ( kT ) f ( kT T )@

i( kTs ) i( kTs Ts )

s

s

s


A w-transzformáció egy fontos tulajdonságának megismeréséhez tekintsük az alábbi mintapéldát. Adott egy folytonos integrátor a 1 P(s) s átviteli függvénnyel. Ts=0.1 másodperces mintavételi idõ mellett a folytonos szakasz zérusrendĦ tartószervvel együtt vett mintavételes alakja (azaz impulzusátviteli függvénye) analitikusan 1 e sT

Z®

Zî( kTs )` Z^ f ( kTs )`

Ts 1 z 1 2 1 z 1

amely átviteli függvénynek a Laplace operátoros tartományban az s-sel való osztás felel meg:

1 s

Ts 1 z 1 2 1 z 1

s

2 1 z 1 Ts 1 z 1

avagy

Ezt az összefüggést bilineáris vagy Tustin transzformációnak hívják, fordított irányban a

z

sT 1 s 2 sTs 1 2

P ( w)

w

2 1 z Ts 1 z 1

wTs 2 wTs 1 2

1

1

z

1 e 0.1s 1 ½ ¾ s¿ ¯ s

0.1 z 1

Z®

-0.05

w-20 w

Ez a folytonos átviteli függvény szemmel láthatóan különbözik attól a P(s) átviteli függvénytĘl, amibĘl kiindultunk. A különbség a P(s) és P(w) átviteli függvények frekvenciafüggvényének összehasonlításával elemezhetĘ: » w= logspace(-1,2,200); » bode(Ps,'r',Pw,'b',Pz,’g’,w),grid; Jól látható, hogy P(w) frekvenciafüggvénye a kisfrekvenciás tartományban, az Z 1/ Ts körfrekvenciáig nagyon jól közelíti P(s) és P(z) amplitúdó-körfrekvencia függvényét (-20dB/dekád meredekségĦ egyenes), fázisszöge pedig ebben a tartományban egybeesik a diszkrét rendszer fázisszögével, amely a folytonos rendszer -90˚-os fázisszögétĘl ebben a tartományban is kissé eltér, mivel a mintavételezés kb. Ts/2 járulékos holtidĘt hoz be. A w transzformált folytonos frekvenciafüggvény fázisszöge tehát a kisfrekvenciás tartományban jól követi a diszkrét rendszer fázisszögét, így a folytonos szabályozótervezésben a mintavételezésbĘl adódó járulékos holtidĘ hatását figyelembe veszi. A pontos vizsgálatokhoz tekintsük a w-transzformált frekvenciafüggvényét. Jelöljük a w-transzformált komplex változójának képzetes részét Q-vel (ez felel meg az s komplex frekvencia w komponensének), nevezzük ezt fiktív frekvenciának. Ahogy az s tartományban formálisan az s=jZ helyettesítéssel származtatjuk a frekvenciafüggvényt, úgy a w tartományban a w=jQ helyettesítéssel élünk: ww

összefüggést jelenti. Abban az esetben, ha egy feladat kapcsán egyrészt folytonos rendszer átviteli függvényérĘl, másrészt diszkrét átviteli (más szóval impulzusátviteli) függvényrĘl, továbbá a fenti transzformációval származtatott átviteli függvényrĘl van szó egyidejĦleg, a bilineáris transzformációval származtatott Laplace operátoros átviteli függvény változóját s helyett w-vel szokásos jelölni, ebben az értelemben beszélhetünk tehát w-transzformációról:

¿

a Matlab segítségével számított módon pedig » Ts=0.1 » Ps=tf(1,[1, 0]) » Ps=zpk(Ps) » Pz=c2d(Ps,Ts,’zoh’) A Control System Toolbox a diszkrét és folytonos átviteli függvényeket hasonlóan kezeli, csak a Ts mintavételi idĘ különbözteti Ęket meg. Folytonos rendszerek esetén a mintavételi idĘ mindig nulla. Az impulzusátviteli függvénybõl a w-transzformált egy from-discrete-to-continuous jellegĦ utasítással számítható: » Pw=d2c(Pz,'tustin')

értéket adja. Figyelembe véve a z operátor idĘben való egylépéses késleltetését, írhatjuk, hogy

tehát a diszkrét idejĦ integrátor (digrátor) átviteli függvénye

P( s) ¾

s

¯

-1

Ts (1 z 1 ) f ( kTs ) 2

½

s

P( z)

kTs Ts

(1 z 1 )i( kTs ) #


jQ

jQ

2 z 1 Ts z 1 z

j e jZTs

2 Ts

tg (

Z Ts 2

)

azaz Q

tehát kis Z értékek esetén

2 Ts

tg (

ZTs 2

)

Q #Z

a növekvõ Z értékekkel viszont egyre növekvõ frekvencia torzulás lép fel. Mindez azt jelenti, hogy a wtranszformációt alkalmazó szabályozótervezés eredményét fokozott körültekintéssel kell ellenĘrizni.



Nézzük ezek után a fenti transzformációt felhasználó tervezési módszer alapgondolatát. Induljunk ki a mintavételes zárt szabályozási rendszerek alapvetĘ struktúráját leíró blokkvázlatából. Ez a rendszer egy hibrid rendszer abban az értelemben, hogy folytonos és diszkrét idejĦ jelek egyaránt szerepelnek benne. r[k]

e[k]

u[k]

C (z )

-

U(z)

D/ A

u(t)

y[k] A/ D Y(z)

y(t)

P (s )

Y(s)

P( z )

Y ( z) U ( z)

A fenti ábrán C(z) a tervezendĘ diszkrét idejĦ szabályozó, P(s) az adott folytonos szakasz. ElsĘ lépésben a D/A nulladrendĦ tartószerv és a P(s) folytonos szakasz együttes mintavételes alakját, a P(z) impulzusátviteli függvényt határozzuk meg. Ehhez a P(s) átviteli függvényen kívül a Ts mintavételi idĘt is meg kell majd adnunk. A tervezésnek ezt a lépését az alábbi blokkvázlat tükrözi. r[k]

e[k]

-

C (z )

u[k]

P( z )

y[k]

Most áttérünk a folytonos frekvenciatartományra, P(z) => P(w) transzformációt hajtunk végre, azaz P(z)ben formálisan z helyére a korábban megadott összefüggés szerint a w változót hozzuk be, C(z) helyett pedig egyszerĦen C(w)-t írhatunk, hiszen e pillanatban sem C(z), sem pedig C(w) nem ismert. r(t)

e(t)

-

C ( w)

u(t)

P ( w)

y(t)

Ez egy folytonos szabályozótervezési feladat, a w változó csupán azt jelzi, hogy transzformált modellrĘl van szó. C(w) megtervezését követĘen magával a tervezéssel végeztünk is, azonban a szabályozási rendszer megvalósítása érdekében egy diszkrét idejĦ C(z) átviteli függvényre van szükségünk, már csak formai okokból is, hiszen egy D/A átalakítót kell meghajtania a szabályozónak. A C(w) => C(z) átalakítás egy inverz w-transzformáció. Természetesen a teljes tervezési folyamat a zárt kör viselkedésének ellenĘrzésével egészül ki, és abban az esetben, ha a zárt kör viselkedése nem felel meg a tervezési specifikációnak, vissza kell térni a C(w) áttervezéséhez. Egy konkrét feladat kapcsán tekintsünk végig egy a fenti tervezési elvet követĘ tervezési folyamatot. Legyen a szabályozandó szakasz

P( s )

1 s(10s 1)

» Ps=tf(1,[10 1 0]) » Ps=zpk(Ps) Legyen a mintavételes szabályozás mintavételi ideje Ts=1 sec: » Ts=1 Származtassuk a P(z) és a P(w) átviteli függvényeket: » Pz=c2d(Ps,Ts,'zoh') » Pw=d2c(Pz,'tustin')


P ( w)

-0.00041625


(w+120) (w-2) w (w+0.09992)

Hasonlítsuk össze a P(s), P(z) illetve a P(w) átviteli függvényĦ folytonos modellek frekvenciafüggvényének amplitúdó és fázis függvényét: » bode(Ps,'r',Pw,'b'); Most látjuk, hogy a számunkra érdekes frekvencia tartomány 10-3 és valamennyivel a közelítés jóságának határa (|1) fölött van, válasszuk a [10-3...1] tartományt: » w=logspace(-3,0,200); » bode(Ps,'r',Pw,'b',Pz,’g’,w); C(w)=1 arányos szabályozó mellett vizsgáljuk meg a viszonyokat: » Cw=1 Számítsuk ki az L ( w) C ( w) P ( w) hurokátviteli függvényt: » Lw=Cw*Pw A fázistöbblet és az erĘsítési tartalék számításához számítsuk ki a felnyitott kör frekvenciafüggvényét, majd használjuk a margin függvényt: » margin(Lw) » [gm,pm,wpi,wc]=margin(Lw) gm= 2.0339 pm= 9.1835 wpi= 0.4510 wc= 0.3102 Látható, hogy a fázistartalék és a vágási körfrekvencia növelése érdekében az arányosnál dinamikusabb szabályozásra van szükség. Legyen a szabályozó PD jellegĦ. Mivel a folyamat már tartalmaz egy integrátort, nem szükséges még egyet betenni a szabályozóba. w 0.09992 C ( w) kcw w6 Vizsgáljuk a szabályozási kör jellemzĘit elĘször kcw =1 feltételezésével. » Cw=tf([1, 0.09992 ],[1, 6]) Számítsuk ki az L( w) C ( w) P ( w) hurokátviteli függvényt: » Lw=Cw*Pw Végeztessük el a Matlab-bal a lehetséges egyszerĦsítéseket 10-2-os pontossággal: » Lw=minreal(Lw,1e-2) A fázistöbblet és az erĘsítési tartalék számításához számítsuk ki a felnyitott kör frekvenciafüggvényét, majd használjuk a margin függvényt: » margin(Lw) Válasszuk meg most a szabályozó kcw átviteli tényezĘjét úgy, hogy a fázistöbblet kb. 60º legyen. » [mag,phase,w]=bode(Lw); » gm=margin(mag,phase-60,w) » kcw=gm Számoljuk ki a szabályozót az új erĘsítéssel. » Cw=kcw*Cw » Lw=kcw*Lw; A rendszer jellemzĘinek ellenĘrzése: » margin(Lw) Transzformáljuk át a szabályozót a ‘z’ tartományba: » Cz=c2d(Cw,Ts,'tustin') A Cz zérus-pólus-erĘsítés alakja: » Cz=zpk(Cz)



z-0.9048 z+0.5 Ha pontos képet akarunk kapni a szabályozott rendszer viselkedésérĘl, akkor a diszkrét szabályozó és a folytonos folyamat viselkedését kell vizsgálni. Ennek a szimulációja a Matlab környezetben nehézkes, azonban a Simulink segítségével könnyen elvégezhetĘ. Indítsuk el Simulink-et és hozzuk létre az alábbi modellt és vizsgáljuk a viselkedést. » simulink Hozzunk létre egy új file-t és másoljuk be a különféle blokkokat. A blokkok paramétereit változtassuk meg a kívánt értékre. A Simulink a Matlab-ban definiált változókat használja. C ( z ) , P ( s ) : Control System Toolbox –>LTI system NulladrendĦ tartószerv: Simulink–>Discrete–>Zero-Order-Hold KülönbségképzĘ: Simulink–>Math–>Sum Step input: Simulink–>Sources–>Step Clock: Simulink–>Sources–>Clock Scope: Simulink–>Sinks–>Scope Output: Simulink–>Sinks–>To Workspace C ( z ) 12.1021

Clock

t

u

time

control

Step Input

y Scope

Cz

Ps

LTI System

LTI SystemPs

Zero-Order Hold

output

Scope1

Véges beállású Mintavételes Szabályozó Tervezésének Alapjai


13. Véges Beállású Mintavételes Szabályozó Tervezésének Alapjai A mintavételes szabályozásokban a szabályozó szerepét egy intelligens eszköz, tipikusan egy mikroprocesszoros PLC (Programmable Logic Controller) tölti be. A PLC feladata egy szabályozási algoritmus megvalósítása, valamint a szabályozó mĦködtetéséhez kapcsolódó jelek kezelése (szĦrés, ADC és DAC átalakítók, jelkondícionáló interfészek). A szabályozási algoritmusokat illetĘen a szoftverrel történĘ megvalósítás számos különleges tervezési módszer kialakítására nyújt lehetĘséget. Egy ilyen szabályozási algoritmus a kimenĘjel véges beállási idejét biztosító, az irodalomban többnyire dead-beat szabályozónak nevezett algoritmus. A véges beállású szabályozó a kimenet hibamentes beállását a mintavételi idĘ egész számú többszöröse alatt kívánja biztosítani. A módszer lényegének megismeréséhez az egyszerĦség kedvéért vizsgálatainkat holtidĘ mentes, stabilis szakaszra végezzük el, az alapjelrĘl pedig feltételezzük, hogy egységugrás. Megjegyezzük, hogy a véges beállású szabályozó akkor is megtervezhetĘ - a bemutatottnál kissé bonyolultabb módon - ha a felsorolt feltételek nem teljesülnek. A következĘkben három lépésben mutatunk be egy gyakorlati igényeket is kielégítĘ megoldást. ElsĘ lépésben a lehetĘ leggyorsabb beállást, azaz az egyetlen mintavétel alatti beállást biztosító szabályozót tervezzük meg. Látni fogjuk, hogy az egylépéses beállás biztosítható, de tipikus mintavételezési frekvenciák esetén hatalmas bemenĘjeleket igényel, ráadásul lengéseket eredményez a mintavételi pontok között. A lengések elkerülésére a tervezési módszerben szét fogjuk választani a kiejthetĘ és ki nem ejthetĘ folyamat-zérusokat. Látni fogjuk, hogy a lengések elkerülése a beállási idĘ megnövelésével válik lehetségessé. Amennyiben a lengésmentesítés után is túlságosan nagy bemenĘjeleket eredményezne a tervezésünk, úgy egy tervezési polinommal tovább finomítjuk a megoldást, növeljük (de továbbra is véges értéken tartjuk) a beállási idĘt. Ami a tervezési metodikát illeti, a tervezés végig a ‘z’ tartományban történik. A módszer érdekessége, hogy bizonyos nem kívánt idĘtartománybeli tulajdonságokat (lengések, túlságosan nagy bemenĘjelek) ‘z’ tartománybeli megfontolásokkal küszöbölünk ki. Az alapfeladat egy mintavételes szabályozó tervezése:

Az eredményt közvetlenül a Scope blokkokkal vizsgálhatjuk vagy a To Workspace output blokkok segítségével visszaküldjük a Matlab felületre és ott ábrázoljuk. »subplot(211); plot(t,y); grid; »subplot(212); plot(t,u); grid;

r[k]

e[k]

-

C (z )

u[k] U(z)

D/ A

u(t)

y(t)

P(s )

Y(s)

y[k] 1.5

Y(z)

A/ D P( z )

1 0.5 0

0

2

4

6

8

10

20

0

-20

0

2

4

6

8

10

Y ( z) U ( z)

A fenti hibrid (folytonos-diszkrét) problémát elĘször tisztán diszkrét idejĦ feladattá konvertáljuk azáltal, hogy meghatározzuk a D/A tartószervet és a P(s) folytonos szakasz együttes mintavételes alakját, a P(z) impulzusátviteli függvényt. Ehhez a P(s) átviteli függvényen kívül a Ts mintavételi idĘt is meg kell majd adnunk. Ezt követĘen megtervezzük a C(z) szabályozót, majd ellenĘrizzük a zárt kör mĦködését. A tervezés lényege, hogy elĘírjuk, hogy a szabályozott rendszer hogyan viselkedjen, azaz elĘírjuk a zárt rendszer T ( z ) eredĘ átviteli függvényét. Véges beállású szabályozó esetén olyan eredĘ átviteli függvényt írunk elĘ, amely biztosítja a véges beállást. C ( z ) P( z ) 1 C ( z ) P( z )

T ( z)

EbbĘl az egyenletbĘl a szabályozó C ( z ) impulzusátviteli föggvénye kifejezhetĘ:

1


C ( z)




Hozzunk létre egy új fájlt és másoljuk be a különféle blokkokat. A blokkok paramétereit változtassuk meg a kívánt értékre. Állítsuk be a menĦbĘl a Simulation–>Parameters–>Stop time paramétert 20-ra. A Simulink a Matlab-ban definiált változókat használja.

T ( z) P( z )(1 T ( z ))

Mintapéldaként tekintsük a

1 (1 5s )(1 10 s ) folytonos folyamatot, a mintavételi idĘ legyen Ts=1 sec. ElĘször is definiáljuk az s és z változókat a szimbólikus adatbevitelhez. » Ts=1; » s=tf('s'); » z=tf('z',Ts); » Ps=1/((1+5*s)*(1+10*s)) Alakítsuk át zérus-pólus alakra. » Ps=zpk(Ps) Annak érdekében, hogy legyen benyomásunk arról, hogy mit jelentenek a minimális beállással kapcsolatos követelmények, nézzük meg a szabályozandó folyamat átmeneti függvényét » step(Ps); A folytonos szakasz és a nulladrendĦ tartó együttes mintavételes modellje: » Pz=c2d(Ps,Ts) P( s)

P( z )

B( z ) A( z )

0.0090559

( z 0.9048) ( z - 0.8187) ( z - 0.9048)

A diszkrét pólusok közvetlenül transzformálódtak a folytonos pólusokból a zi alapján, de egy járulékos zérus is megjelent.

e piTS

e

T S Ti

összefüggés

ElĘször tervezzük meg az egylépéses beállást biztosító szabályozót. Az egylépéses beállás feltétele az, hogy a zárt kör eredĘ átvitele az egységugrásnak feltételezett alapjel és a kimenĘjel között éppen az egylépéses késleltetést leíró z-1 érték legyen: T ( z ) z 1 Rendezés után a szabályozó egyenlete: 1 C ( z) P ( z )( z 1) » Tz=1/z » Cz=Tz/(Pz*(1-Tz)) » Cz=minreal(Cz) C ( z ) 110.425

( z - 0.8187) ( z - 0.9048) ( z 0.9048) ( z -1)

Ha az így kapott szabályozóval megvizsgáljuk a rendszer viselkedését diszkrét idĘben, akkor meg is kapjuk a várt eredményt. Az átmeneti függvény egylépéses késleltetést mutat. » step(Cz*Pz/(1+Cz*Pz)) Azonban ha pontos képet akarunk kapni, akkor a folytonos folyamat viselkedését kell vizsgálni a mintavételi pontok közötti idĘbeli lefolyást is tekintve. Ennek a szimulációja nehézkés Matlab környezetben, azonban a Simulink segítségével ez könnyen elvégezhetĘ. Indítsuk el a Simulink -et és hozzuk létre az alábbi modellt. » simulink

2

Clock

t

u

time

control

Cz Step Input

LTI System

y Scope

output

Ps Zero-Order Hold

LTI SystemPs

Scope1

C ( z ) , P ( s ) : Control System Toolbox –>LTI system NulladrendĦ tartószerv: Simulink–>Discrete–>Zero-Order-Hold KülönbségképzĘ: Simulink–>Math–>Sum Step input: Simulink–>Sources–>Step Clock: Simulink–>Sources–>Clock Scope: Simulink–>Sinks–>Scope Output: Simulink–>Sinks–>To Workspace

Az eredményt közvetlenül a Scope blokkokkal vizsgálhatjuk vagy a To Workspace output blokkok segítségével visszaküldjük a Matlab felületre és ott ábrázoljuk. »subplot(211); plot(t,y); grid; »subplot(212); plot(t,u); grid; A szimuláció jól mutatja, hogy a mintavételi pontok között lengések jelentek meg. E lengések oka az, hogy a szabályozó a szakasz B(z) polinomját kompenzálta, és ennek 1.5 következtében van olyan pólusa, amely az u[k] jelben oszcillációt eredményez. Ha a folyamatot a 1 B( z ) 0.5 P( z ) A( z ) 0 0 5 10 15 20 alakban írjuk fel, akkor a szabályozó 200 A( z ) C( z) , B( z )( z 1) 0 és a nyílt kör -200 A( z ) B( z ) 0 5 10 15 20 L( z ) C ( z ) P ( z ) B( z )( z 1) A( z ) Jól látható, hogy a folyamat zérusai megjelennek a szabályozóban pólusként. A P( z ) impulzusátviteli függvénynek egy zérusa van, z1 = - 0.9048 . Vizsgáljuk meg most ezt közelebbrĘl. A szabályozóban ez pólusként jelentkezik. » C1z=1/(z+0.9048) » step(C1z) Látható, hogy ez a komponens okozza a lengéseket. Vizsgáljuk meg, hogyan néz ki ez a tag a folytonos tartományban. Mivel z e sTS , s ln( z ) / TS » p1=log(-0.9048) p1 = -0.1000 + 3.1416i

3



Látható hogy ez egy kis csillapítási tényezĘjĦ komplex pólus-párnak felel meg, amelynek a lengési S 3.14 , épp a mintavételezési körfrekvencia kétszerese. Nyilvánvalóan ez nem egy valós frekvenciája Ts

zérus, csak a mintavételezés miatt jelent meg, ezért ezt nem kell kompenzálni. A lejjebb található függelék megmutatja, hogy a komplex pólus-párok hogyan transzformálódnak a z tartományba. Válasszuk szét a kiejtéssel kompenzálható illetve nem kompenzálható folyamat-zérusokat: B( z ) B ( z ) B ( z )


eredĘ átviteli függvényében.

T ( z)

ahol

k

Bn ( z ) , zk

1 deg ( B ) , Bn ( z )

B ( z) . B (1)

A z k komponens szükséges, hogy a rendszer kauzális, megvalósítható legyen (a nevezĘ fokszáma nagyobb legyen mint a számlálóé). A szabályozóval szemben támasztott másik követelmény, hogy az egységugrás alapjelet hiba nélkül kövesse. Ezt elérhetjük, ha B ( z ) -t leosztjuk, normáljuk a stacionárius (z=1 helyettesítéssel felvett) értékével. B( z ) B ( z ) B (1)

B ( z) B (1)

T ( z)

z 0.9048 1.9048

0.525 z + 0.475 és Bn ( z )

1.5

Bn ( z ) zk

1 0.5

A( z ) F ( z ) Bn ( z )[ z k Bn ( z ) F ( z )]

0

0

2

4

6

8

0

2

4

6

8

20

Figyeljük meg, hogy F (1) 1 , azaz a tervezési polinom nem befolyásolja a zérus statikus hibát. » Tz=Fz*Bmn/(z^2) » Cz=Tz/(Pz*(1-Tz)) » Cz=minreal(Cz)

0

-20

A korábbiakban említettük, hogy a véges beállást biztosító szabályozó akkor is megtervezhetĘ, ha a szabályozandó folyamatnak holtideje is van. Tegyük fel, hogy a folytonos folyamat Td holtideje egész számú többszöröse a Ts mintavételi idĘnek, és legyen ez az arány d=Td/Ts (dt0): P ( s )e sTd

==>

B( z ) d z A( z )

P( z )

Továbbra is feltételezve, hogy a folytonos folyamat önmagában stabilis, valamint azt, hogy az alapjel egységugrás, a zárt kör impulzusátviteli függvényére vonatkozó feltételünk továbbra is k pz B (1)=0.0173

Bn ( z ) zk

1.5

T ( z) 1

EbbĘl a kiadódó szabályozó C ( z)

F ( z)

így ebbĘl következĘen

Bn ( z ) Bn ( z ) .

Az így kapott értékek: B (1) 1 0.9048 , Bn ( z )

z2 z 1 3z 2

tervezési polinomot választva, annak simító hatása az átmeneti függvénye alapján érzékelhetĘ » plot(0:5,[0 1/3 2/3 1 1 1],0:5,[0 1/3 2/3 1 1 1],'*r'),grid; vagy diszkrét idĘben » Fz=(z^2+z+1)/(3*z^2) » step(Fz) A tervezési polinomot felhasználva a szabályozási kritérium

C ( z)

A következĘ lépésben olyan szabályozót tervezünk, amely nem kompenzálja ezt az oszcillációt okozó zérust. A kimenetet tehát a következĘ alakban írjuk elĘ:

1 z 1 z 2 3 3 3

F ( z)

B ( z ) tartalmazza a kompenzálható gyököket és B ( z ) a nem kompenzálhatókat. Ha nem kompenzáljuk a B ( z ) valamelyik zérusát, akkor az közvetlenül meg fog jelenni a kimenĘjelnek az alapjelre vonatkozó


azonban k értékét a korábbiakhoz képest meg kell növelni a holtidĘnek megfelelĘ d lépésszámmal:

0.5

A( z ) B ( z )[ z k Bn ( z )] n

0

k 0

2

4

6

8

Ezt a következĘképpen számíthatjuk: 100 % vagy » Bm =(z-Pz.z{1}) 50 » Bm =(z+0.9048) 0 » Bpn=Pz.k*dcgain(Bm) -50 » Bmn=Bm/dcgain(Bm) 0 2 4 6 8 » Tz=Bmn/(z^2) » Cz=Tz/(Pz*(1-Tz)) % vagy »(Cz=minreal(Cz,0.001)) » Cz=minreal(Cz) Szimuláljuk ismét a viselkedést a Simulink modellel. Látható, hogy sikerült kiküszöbölnünk a lengéseket, és egyúttal a bemenĘjel “forszírozása” is korlátozottabbá vált. A rendszer lassúbbá vált abban az értelemben, hogy most két mintavételi pont után éri el a kimenĘjel mintavett értéke az állandósult állapotot. Mivel a mintegy ötvenszeres túlvezérlés még így is jóval meghaladhatja egy átlagos beavatkozó szerv lehetĘségeit, módszert kell keresnünk arra az esetre, ha a véges beállási idĘ biztosítása mellett tovább kívánjuk csökkenteni a túlvezérlést. Egészítsük ki a szabályozási algoritmust egy tervezési, szĦrĘ polinommal, amelynek segítségével mintegy “megvezetjük” a véges idejĦ felfutást. Például az

4

1 deg ( B ) d

A fentiek következtében a keresett szabályozó mellett a hurokátviteli függvény C ( z ) P( z )

z k Bn ( z ) 1 z k Bn ( z )

, amibĘl

C ( z)

z d A( z ) Bn ( z )[ z k Bn ( z )]

.

A mintapéldában megadott folytonos folyamatot egészítsük ki egy Td=1 másodperces holtidĘvel (d=1), és szimuláljuk a mintavételes szabályozás mĦködését. » Td=1 » d=Td/Ts » Tz=Bmn/(z^(2+d)) » Cz=Tz/(Pz*(1-Tz)) » Cz=minreal(Cz) Tegyünk be egy delay blokkot (Simulink–>Continuous–>Transport Delay) a Simulink modellbe és állítsuk be a paraméterét Td-re. Az elĘzĘkhöz hasonlóan a szabályozó egy további F(z) tervezési polinommal rugalmasabb specifikációknak is eleget tud tenni. Függelék: Nézzük meg, hogy pontosan hogyan néz ki egy adott csillapítású komplex konjugált póluspár kontúrja a ztartományban. Az s-tartományban a konstans ] vonalak az origón áthaladó s=V+jZ egyenesek, ahol egy

5



V 1 9 2 tartozik, amely komplex s értékeket a z e sTs leképzés szerint szív alakú ] görbékbe képez le. Ezt szemléltetendĘ legyen például » szigma=0:0.01:1.6; » zeta=0.4; » Ts=1; » z=exp(Ts*(-szigma+j*sqrt(1-zeta*zeta)*szigma/zeta)); » plot(real(z),imag(z),real(z),-imag(z)),grid; A zárt görbén belüli területre (ahol a csillapítás nagyobb, mint a kontúr mentén) esnek a B( z ) polinom

HoltidĘs Rendszer Szabályozása Smith Prediktorral

14. HoltidĘs Rendszer Szabályozása Smith Prediktorral

adott V értékhez Z

gyökei közül azok, amelyek a B ( z ) polinomot alkotják, és a zárt görbére vagy azon kívüli területre esnek a B( z ) polinom gyökei közül azok, amelyek a B ( z ) polinomot alkotják.


HoltidĘs rendszerek szabályozására az 1950-es években kifejlesztett Smith prediktor jól használható. A Smith prediktoros tervezés alapötlete, hogy a holtidĘs rendszert hatásvázlat átalakítással egy olyan struktúrává alakítjuk át, ahol a holtidĘ a zárt körön kívül jelenik meg. r(t)

e(t)

Csp (s )

u(t)

P( s)e sTD

y(t)

r(t)

-

e(t)

C (s)

u(t)

P( s)

-

e sTD

y(t)

A holtidĘmentes szakaszhoz tervezzük a C(s) szabályozót, majd a holtidĘs szakaszhoz alkalmazott Csp(s) szabályozót az alábbi számítással határozzuk meg: C (s) . Csp ( s ) 1 (1 e sT )C ( s ) P( s ) D

Az egyenlĘség könnyen igazolható, ha felírjuk a zárt rendszer átviteli függvényét a két esetre: sT Csp ( s ) P( s )e D C ( s) P( s ) sTD

1 Csp( s ) P( s )e sTD Csp ( s ) 1 C ( s ) P( s ) Csp ( s )

1 C ( s) P( s)

e

C ( s ) 1 Csp( s) P( s )e sTD

C ( s) 1 1 e sTD C ( s ) P( s )

Elméletileg a Smith prediktort folytonos és mintavételes rendszerekre is lehet alkalmazni, de mivel a késleltetés megvalósítása folytonos rendszerek esetén nehézkés, ezért általában csak mintavételes szabályozóként alkalmaznak Smith prediktort. A tervezés két lépésben oldható meg. Az elsĘ lépésben a C ( s ) szabályozót kell megtervezni a holtidĘmentes P( s ) folyamathoz, a második lépésben pedig a Csp ( s ) Smith prediktoros szabályozó a fenti összefüggéssel kiszámítható. Példaként tervezzünk Smith prediktoros szabályozót egy három idĘállandós holtidĘs rendszerhez. A folytonos folyamat átviteli függvénye:

PD s K

5, T1 10, T2

P s e sTD

4, T3 1, TD

Ke sTD , 1 sT1 1 sT2 1 sT3

10 paraméterekkel.

Tervezzük meg a C ( s ) szabályozót úgy, hogy a fázistartalék 60q legyen és a beavatkozójel maximuma ne haladja meg az umax=10 értéket. »K=5; T1=10; T2=4; T3=1; Td=10; ft=60; Definiáljuk a 's' változót szimbólikus adatbevitelhez. »s=tf('s'); Használjuk a Control System Toolbox LTI struktúráit: »Ps=1/((1+s*T1)*(1+s*T2)*(1+s*T3)) Átalakíthatjuk az átviteli függvényt zérus-pólus alakra. »Ps=zpk(Ps)

6

1



ElsĘ lépés: A C ( z ) szabályozó tervezése A mintavételezési idĘ megválasztására jó gyakorlati szabály, hogy legyen kisebb a legkisebb idĘállandónál. (Szimulációs feladatokhoz választhatjuk a legkisebb idĘállandó tizedére, szabályozáshoz megfelelĘ, ha a legkisebb idĘállandó felére, harmadára vesszük fel. Ezenkívül célszerĦ úgy felvenni, hogy a holtidĘ a mintavételezési idĘ egész számú többszöröse legyen. »Ts=0.5; A C ( z ) szabályozót a kisfrekvenciás közelítés módszerével tervezzük. A folytonos szakasz mintavételi tartószervvel együtt vett impulzusátviteli függvénye (z tartomány): » Pz=c2d(Ps,Ts) A kapott impulzusátviteli függvény:

P( z )

0.0004416

(z+0.2254) (z+3.167) (z-0.9512) (z-0.8825) (z-0.6065)



»margin(Lz) Vizsgáljuk meg diszkrét rendszer kimenĘjelét és beavatkozójelét. »Hz=Cz*Pz/(1+Cz*Pz) »Uz=Cz/(1+Cz*Pz) »subplot(211); step(Hz) »subplot(212); step(Uz) A hibrid (diszkrét - folytonos) rendszer viselkedését Simulink környezetben vizsgálhatjuk. Hozzuk létre az alábbi rendszert: A modell a Cz, Ps, paramétereket a Matlab felületrĘl veszi, majd a szimuláció eredményét a grafikus megjelenítés mellett az u, y változókban tárolja. A holtidĘt most állítsuk nulla értékre. Állapítsuk meg a kimenĘjel túllövését és beállási idejét és ellenĘrizzük, hogy a beavatkozójel nem haladja-e meg az elĘírt maximális értéket.

Clock

t

u

time

control

y Scope

Póluskiejtéses PID szabályozót tervezünk.

A P ( z ) legnagyobb idĘállandóhoz tartozó pólusa, p1

e

Ts T1

0.9512 helyett behozunk egy integrátort (PI), és a második legnagyobb idĘállandóhoz tartozó pólusa, p2 0.8825 helyett

Step Input

output

Cspz

Ps

LTI System

LTI SystemPs

Zero-Order Hold

Transport Delay

Scope1

egy ideális deriváló tagot használunk (PD). Ezek alapján a szabályozó impulzusátviteli függvénye:

C z

kc

(z-0.9512) (z-0.8825) . z 1 z

A szabályozóban szereplĘ kc konstanst úgy állítjuk be, hogy a fázistartalék a kívánt érték legyen. Az elsĘ lépésben egynek vesszük fel. »p1=exp(-Ts/T1) »p2=exp(-Ts/T2) »kc=1; »Cz=zpk([p1 p2],[0 1],kc,Ts) Ezek alapján a szabályozó impulzusátviteli függvénye:

C z

kc

(z-0.9512) (z-0.8825) . z 1 z

A nyílt kör impulzusátviteli függvénye: »Lz=Cz*Pz Végezzük el a Matlab-bal a lehetséges egyszerĦsítéseket: »Lz=minreal(Lz) A kc értéket leolvashatjuk az ltiview grafikus megjelenítĘrĘl. A -120˸-hoz tartozó erĘsítés reciproka lesz az értéke. »ltiview('bode',Lz); Közvetlenül a margin utasítással is számolható, »[mag,phase,w]=bode(Lz); »kc=margin(mag,phase-60,w); A kapott konstans kc=33.46. Számoljuk át újra az átviteli függvényeket. »Cz=kc*Cz; »Lz=kc*Lz; EllenĘrizzük az így kapott rendszert. EllenĘrizzük a fázistartalékot.

2

Második lépés: A Csp( z ) Smith prediktor kiszámítása: Csp ( z )

C ( z)

1 (1 z d )C ( z ) P ( z )

, d

»z=tf('z',Ts); »d=Td/Ts »Cspz=Cz/(1+(1-z^(-d))*Pz*Cz)

TD Ts 1.5 1 0.5

A rendszert a fenti Simulink modellel szimuláljuk. Állítsuk be a holtidĘt a Td változóra. A Simulation-Parameters-Stop time paramétert állítsuk be 30ra. A kimenĘjel és a beavatkozójel idĘbeli lefolyását az alábbi ábra mutatja. A szabályozás sokkal gyorsabb mĦködésĦ, mint a holtidĘs szakaszhoz tervezett soros PID szabályozóval lenne.

3

0 0

10

20

30

0

10

20

30

40 20 0 -20

Bevezetés a Két Szabadságfokú Polinomiális SzabályozóTervezésbe


15. Bevezetés a Két Szabadságfokú (2DF) Polinomiális Szabályozó Tervezésbe A két szabadságfokú (2DF – two degree of freedom) szabályozási struktúra különbözik a szokásosnak nevezhetĘ, egyetlen soros kompenzáló elemet tartalmazó zárt szabályozási körök struktúrájától, elnevezésével összhangban két szabályozót tartalmaz (elĘszĦrĘ és visszacsatolás) és az alábbiak szerint épül fel : zavarás

r(t)

y(t)

u(t)


alapjelkövetésre is vonatkozzék, miután a nem kielégítĘ dinamikát az elĘszĦrĘvel még módunkban áll módosítani. A fenti elvek szisztematikus alkalmazása a robusztus szabályozáselmélet tárgya, jelen bevezetĘ jellegĦ anyag célja csupán az alapvetĘ összefüggések felvázolása és az alkalmazott technika bemutatása.

A szabályozó alapegyenlete Folytonos, egy bemenetĦ-egy kimenetĦ lineáris folyamatokat vizsgálunk, amelyeket mintavételesen, zárt körben szabályozunk. Ts mintavételi idĘt feltételezve legyen a nulladrendĦ tartószerv és a P(s) átviteli függvényĦ folytonos folyamat együttesen vett impulzusátviteli függvénye P(z), ahol

folyamat

elõszûrõ

P ( z ) (1 z 1 )Z{

-


P( s) } s

Tegyük fel továbbá, hogy P(z) racionális törtfüggvény, melynek alakja az A(z) és B(z) polinomok bevezetésével

B( z ) A( z )

P( z )

visszacsatolás

A bevezetett polinomok a z-transzformáció z operátorának változói:

A( z ) mérési zaj

Zárt szabályozási körökkel szemben általában az alábbi követelmények, gyakorlati igények merülnek fel: - Alapjelkövetés: a szabályozott jellemzĘ kövesse az alapjelet (referencia jelet) - Zavarelhárítás: a szabályozott jellemzĘ viselkedését ne befolyásolják a fellépĘ zavarások - Robusztusság: a folyamatról rendelkezésünkre álló információ pontatlanságára ne legyen érzékeny a zárt szabályozási kör viselkedése. Az egy szabadságfokú struktúrák egyetlen szabályozót tartalmazván a fenti három követelménynek azonos módon tesznek eleget. A felsorolt szabályozási funkciók a fenti 2DF struktúrával kettéválaszthatók, a zavarelhárítás és a robusztusság (és természetesen a zárt kör stabilitásának biztosítása is!) a visszacsatolásban lévĘ kompenzálással, míg az alapjelkövetés (szervó) az elĘszĦrĘvel tartható kézben. A logikus tervezés tehát a visszacsatolás tervezésével kezdĘdik és ennek eredményét figyelembe véve az elĘszĦrĘ tervezésével zárul. A 2DF struktúra alkalmazásával a funkcionális szétválasztás eredményeképpen az alapjelkövetést és a zavarelhárítást illetĘen különbözĘ követelményeknek tehetünk eleget. Ezen nyilvánvaló tervezési rugalmasságon túlmenĘen azonban a 2DF struktúra rendelkezik egy további elĘnnyel is. Abban az esetben ugyanis, ha a szabályozott jellemzĘ érzékelésekor jelentĘs mérési zajokra kell számítanunk, akkor a visszacsatolás megtervezéséhez egymásnak ellentmondó feltételeket kell kielégítenünk, ugyanis a hatékony zavarelhárítás egyúttal a mérési zajok kiemelését is eredményezi, ami nyilvánvalóan elkerülendĘ. A visszacsatolás csak kompromisszumok árán tervezhetĘ meg, a jó tervezés viszonylag hatékony zavarelhárítást és még elviselhetĘ zajkiemelést fog megvalósítani. A mérési zaj rendszerint nagyobb frekvenciájú összetevĘket tartalmaz, mint a zavarás, így a megfelelĘ zavarelhárítási és zajelnyomási tulajdonságok más-más frekvenciatartományban biztosíthatók. E kompromisszum ugyanakkor nem kell, hogy az

1

ao z n a1 z n1 ... an

B ( z ) bo z m b1 z m1 ... bm

a polinomok fokszámára pedig alkalmazzuk a továbbiakban az alábbi jelölést:

degreeÂ( z )` wA n

degree^B( z )` wB

m

Feltételezzük, hogy az A(z) és B(z) polinomoknak nincsenek közös gyökei, továbbá azt is, hogy a gyakorlati eseteknek megfelelĘen mdn (kauzális rendszer). A mintavételezĘ és a zérusrendĦ

r[k]

y[k]

u[k]

T (z) R( z)

P(z) -

S ( z) R( z )

tartószerv beiktatása után a szabályozási struktúra a következĘ: A szabályozó egyenlete:

u[k ]

T ( z )r[k ] S ( z ) y[k ] R( z )

A szabályozó egyenletét a folyamatot és a nulladrendĦ tartószervet együttesen leíró

y[k ]

B( z ) u[k ] A( z )

egyenletbe behelyettesítve

2



B( z )T ( z ) r[k ] A( z ) R( z ) B( z ) S ( z )

y[k ]

adódik a zárt rendszer átvitelére az r[k] alapjel és az y[k] szabályozott jellemzĘ között. A szabályozótervezés célját nagyon egyszerĦen meg tudjuk fogalmazni: azt szeretnénk elérni, hogy a zárt rendszer átviteli tulajdonságai egy általunk megadott modell szerint alakuljanak, nevezetesen az r[k] alapjel és a y[k] szabályozott jellemzĘ között egy követni kívánt referenciamodell átviteli függvényét írjuk elĘ:

y[k ]

Bm ( z ) r[k ] Am ( z ) Bm ( z ) Am ( z )

egyenlet adja az {R(z),S(z),T(z)} polinomhármassal jellemzett szabályozó feltételi egyenletét, amelyben a B ( z ), A( z ), Bm ( z ), Am ( z ) polinomok ismertek.

felbontással legyen

Bm ( z )

Bm ( z ) Ao ( z ) Am ( z ) Ao ( z )

ElĘfordulhat, hogy az Am és Bm modell polinomok közös gyököket tartalmaznak. Mivel a legegyszerĦbb (legalacsonyabb fokszámú) megoldást keressük, a közös gyököket is az A0(z) polinomba emeljük ki, így a továbbiakban Am(z) és Bm(z) már relatív prím (közös gyök nélküli) polinomokat jelölnek. A továbbiakban a fenti egyenletet kívánjuk megoldani úgy, hogy a keresett polinomok a lehetĘ legkisebb fokszámúak legyenek, az átviteli függvények kauzális szĦrĘket eredményezzenek, valamint a mintavételi idĘpontok között ne legyen lengés.

ahol B (z ) tartalmazza a nem kompenzálandó zérusokat. Az R( z ) polinommal viszont egyszerĦsíthetjük a kompenzálható B ( z ) polinomot.

B ( z ) Rc( z )

R( z ) Tehát Ennek megfelelĘen

B ( z ) B ( z )T ( z ) A( z ) R( z ) B( z ) S ( z )

Bm ( z ) Ao ( z ) Am ( z ) Ao ( z )

B ( z ) Bmc ( z ) Ao ( z ) Am ( z ) Ao ( z )

B ( z )T ( z ) A( z ) R ( z ) B ( z ) S ( z )

Bmc ( z ) Ao ( z ) Am ( z ) Ao ( z )

a szabályozó polinomok feltételi egyenlete. Vegyük észre, hogy mivel az R ( z ) polinom tartalmazza a B ( z ) tényezĘt, B (z ) az A( z ) R ( z ) B ( z ) S ( z ) polinomnak osztója.

B ( z )[ A( z ) Rc( z ) B ( z ) S ( z )]

így

T ( z) A( z ) Rc( z ) B ( z ) S ( z )

Bmc ( z ) Ao ( z ) Am ( z ) Ao ( z )

ahonnan

T ( z)

Bmc ( z ) Ao ( z )

és

A( z ) Rc( z ) B ( z ) S ( z )

Am ( z ) Ao ( z )

adódik. A szabályozó meghatározásához a fenti egyenletet kielégítĘ módon kell megkeresnünk az R’(z) és S(z) polinomokat. A szabályozási kör blokkvázlata alapján láthatjuk, hogy a soros kompenzációhoz hasonlóan a hurokba tehetünk integráló hatást az R'(z) polinomon keresztül. Ha az integrátorok száma l, akkor

R ' ( z ) ( z 1)l R1 ( z ) Vagyis az R(z) polinomnak z=1 annyiszoros gyöke, ahány integrátort a szabályozó tartalmaz.

Rendszertechnikai megfontolások A kauzalitás feltételei az ábrán megadott visszacsatoló és elĘrevezetĘ szabályozókra:

wS d wR wT d wR A véges beállást biztosító szabályozó tervezésekor a korábbiakban már láttuk, hogy az eredĘ átviteli függvényre tett elvárásaink megfogalmazásakor nem csupán a kimenĘjel mintavételi idĘpontokban felvett értékeire kell figyelnünk, hanem a mintavételi pontok közötti esetleges lengések elkerülésére is. A különbözĘ tartományokban végzett vizsgálatok eredményének kiértékelésébĘl tudjuk, hogy az idĘtartománybeli nem kívánt viselkedés elkerülése ztartománybeli megfontolások alapján is elvégezhetĘ. A véges beállást biztosító szabályozó tervezéséhez hasonlóan most is úgy intézkedünk, hogy a B(z) polinom lengéseket okozó gyökei jelenjenek meg a zárt kör átviteli függvényében, más szóval a modell tartalmazza a folyamat ún. inverz labilis és gyengén csillapított (“szívgörbén kívüli”) zérusait, tehát a

B( z )

B ( z ) Bmc ( z )

A( z ) R( z ) B( z ) S ( z ) Látható, hogy a fenti egyenletben a megoldás a számláló és a nevezĘ polinomok egyenlĘségére vezethetĘ vissza. Két polinom viszont akkor és csak akkor egyenlĘ, ha fokszámuk és összes együtthatójuk egyenlĘ. Ez igen szigorú feltételt ad a baloldal polinomjainak fokszámára. Ezért a jobboldalt lényegileg nem változtatva bevezetünk egy A0(z) polinomot, amely lehetĘvé teszi, hogy a fokszámok megfelelĘek legyenek mindkét oldalon.

B( z )T ( z ) A( z ) R ( z ) B ( z ) S ( z )


majd az egyszerĦsítés után

amibĘl a

B( z )T ( z ) A( z ) R ( z ) B ( z ) S ( z )


B ( z) B ( z)

3

A szabályozó feltételi egyenletének megoldása A szabályozót meghatározó

A( z ) Rc( z ) B ( z ) S ( z )

Am ( z ) Ao ( z )

egyenlet megoldása további megfontolásokat igényel, hiszen gondoljunk csak arra, hogy egy egyenletbĘl két ismeretlent kell meghatároznunk, ráadásul mindkét ismeretlenünk egy-egy polinom. Általában, ha egy

AX BY

C

alakú (valós együtthatós) polinomiális egyenletet adott A, B és C polinomok mellett X és Y polinomokra kell megoldanunk, ún. diophantoszi egyenlet megoldásáról beszélünk. Ha arra gondolunk, hogy egy ilyen jellegĦ egyenletet a két oldal együtthatóinak összehasonlítása révén oldhatunk meg, nyilvánvaló, hogy elĘször a megoldhatóság feltételeit kell tisztáznunk. A megoldhatóság feltétele ismert algebrai eredmény: A és B legnagyobb közös osztója osztója kell, hogy legyen C-nek is. Könnyen belátható, ha találunk egy {Xo , Yo} megoldást, akkor

X

X o QB

Y

4

Yo QA





wS d wA l 1

is megoldás, mivel

AX BY

A(Xo QB) B(Yo QA)

AXo BYo

C

Egyetlen megoldás akkor létezik, ha a fokszámokra vonatkozóan teljesülnek a következĘ feltételek: wX wB vagy wY wA . Ha ugyanis a fenti feltételek bármelyike teljesül, a fokszám kritériumok úgy alakulnak, hogy a diophantoszi egyenlet egyértelmĦen meghatározza X vagy Y fokszámát, az együtthatók pedig nem engednek szabadsági fokot.

EgyszerĦsítési lehetĘségek, fokszámfeltételek és statikus hiba Foglaljuk tehát össze az eddigi eredményeinket: Kauzalitás

wS d wR wT d wR wB d wA Polinomok szabad paramétereinek megválasztása. Két polinom hányadosában mindig elérhetĘ az, hogy az egyik polinom vezetĘ együtthatója (a legmagasabb fokszámú tag együtthatója) 1 legyen (elosztjuk a törtfüggvényt a nem 1 vezetĘ együtthatóval). Így az alábbi monic (1 vezetĘ együtthatójú) polinomok lesznek: A, R 1' , Am, A0.

6. Megoldjuk a diophantoszi egyenletet. Általában a modellt is igyekszünk minél egyszerĦbbre választani, így sokszor a Bm(z)=konstans választással élünk.

Megvalósítás Fontos gyakorlati szempont, hogy amennyiben R(z) integrátor(oka)t tartalmaz, az integrátor(oka)t ki kell emelni az elĘrevezetĘ és a visszavezetĘ ágból is, és egyetlen példányban el kell helyezni a hurokban, az összegzĘ pontot követĘen a szabályozott szakasz elé az alábbi ábra szerint. EllenkezĘ esetben a beavatkozójelet két labilis tag kimenete különbségeként kellene elĘállítani, ami nyilvánvaló gyakorlati nehézségekbe ütközne. Az integráló hatást tehát „bevisszük” a különbségképzĘ mögé, amit az átviteli függvények szétválasztásával érünk el. Az átviteli függvényeket a z 1 eltolási operátor függvényében adjuk meg. r(z)

T(z-1)

e(z)

rF (z)

1 R(z-1)

yF (z)

és

max(wA wR',wB wS) wA m wA 0

Viszont akkor járunk legjobban (akkor lesz a legegyszerĦbb R’és S), ha az elsĘ tagot választjuk maximálisra, vagyis:

wA wR' wA m wA 0 wR' wA m wA 0 wA A wX wB vagy wY wA alapján az egyértelmĦ megoldás feltétele:

wS d wA l

B m (1) 1 A m (1) Ezt egy konstanssal való szorzással könnyen biztosíthatjuk. A fentiek alapján a szükséges lépések: 1. A modellt úgy választjuk meg, hogy átviteli tényezĘje egységnyi legyen. 2. Megválasztjuk az integrátorok számát (l). 3. A wR' wA m wA 0 wA összefüggés alapján meghatározzuk R’ fokszámát, szükség esetén az A0-t is felhasználjuk, ha nincs szükség rá, akkor A0 konstans. 4. Meghatározzuk T fokszámát a wT wB'm wA 0 összefüggés alapján. 5. S fokszámát a lehetĘ legkisebbre választjuk. Általában

5

S(z-1)

Példa. Kéttárolós lengĘ szakasz szabályozása A fentiek szemléltetésére vizsgáljunk meg egy példát. A folytonos szakasz legyen egy kéttárolós lengĘ tag, melynek idĘállandója Tp=1 illetve csillapítási tényezĘje ȗp= 0.25 adott. A szakasz átviteli függvénye:

P(s)= A statikus hiba elkerüléséhez a zárt kör átvitelének állandósult állapotban egységnyinek kell lennie. Ez a kezdeti és végérték tétel értelmében az eredĘ átviteli függvény vizsgálatát jelenti z=1 esetén. Ez lényegében a modellen múlik:

y(z)

P(z)

-

Fokszámfeltételek.

wT wB 'm wA 0

u(z)

1 1 = 2 2 Tp s 2 +2Tp Ȣ ps+1 s +0.5s+1

A Matlab programrész a szakasz megadására: » Tp=1; % idoallando » wp=1/Tp; % sajatfrekvencia » ZetaP=0.25; » numP=1; » denP=[Tp^2 2*ZetaP*Tp 1]; » Hp=tf(numP,denP) A szakasz természetes frekvenciája (Ȧp=1/Tp ), Ȧp= 1 lesz, átmeneti függvénye pedig: » step(Hp); A lengéseket szeretnénk elsĘsorban csillapítani, illetve a folyamat lefolyását kissé felgyorsítani. Mintavételezve Ts=0.2 sec. mintavételi idĘvel, a folyamat és a zérusrendĦ tartószerv együttes impulzusátviteli függvénye a következĘ:

6


P(z)=(1-z -1 )Z{


P(s) B(z) }= s A(z)

Pm (z)=

A Matlab programrészlet: » Ts=0.2; % mintaveteli ido » Hpz=c2d(Hp,Ts,'zoh') » [numPz,denPz]=tfdata(Hpz,'v');

B (z )= B + (z )B - (z ) szerinti felbontásában kövessük a véges beállású szabályozónál alkalmazott konvenciót ( B (1) 1 ): és

B m (z) A m (z)

» Tm=0.75; » ZetaM=0.707; » numM=1; » denM=[Tm^2 2*ZetaM*Tm 1]; » Hm=tf(numM,denM); A referenciamodell impulzusátviteli függvénye: » Hmz=c2d(Hm,Ts,'zoh'); » [numMz,denMz]=tfdata(Hmz,'v'); A

diszkrét

átviteli

függvénybĘl

csak

az

Am

nevezĘt

vesszük

át,

a

számláló

A polinomok, melyeknek vezetĘ együtthatójuk 1 (monic): A, R1, Am, Ao; a polinomok melyeknek vezetĘ együtthatójuk nem 1 (non-monic): B-,B+, Bm’, S.

B m (z)=B cm (z)B - (z)=const B(z)/B(1) A Matlab programrészlet amely a felbontást realizálja:

Az eddig ismert fokszámok a következĘk:

Bz=numPz; Az=denPz; Bmz=Bz/(Bz*ones(length(Bz),1)); Bmz=[Bmz(2) Bmz(3)]; % B-, a nem kiejtheto zerusok Bpz=Bz*ones(length(Bz),1); % B+, a kiejtheto zerusok

w A=2, w B=1, w B + =0, w B - =1 . A további fokszámfeltételek a következĘk:

wA m =1+wB - =2

Legyen a referenciamodell is másodrendĦ (lásd fokszámfeltételek), melynek folytonos paraméterei legyenek Ȧm és Ȣm (Ȣm=cosĳ). A pólusok elhelyezkedését a komplex számsíkon az alábbi ábra mutatja. Im{s}

so1

(szükséges egy másodfokú megfigyelĘ polinom)

R(z)=(z-1)(z+r)=(z-1)R1' (z) 2

S(z)=s0 z +s1z+s 2 , ahol l=1 az integrátorok számát jelöli. wS=wA+l-1=2 Mivel T(z)=Bcm (z)A o (z)=const A o (z)=A m (1)A o (z) következik, hogy wT=2 A fokszámfeltételeket felírva következik, hogy szükséges egy másodfokú Ao megfigyelĘ polinom. Legyen ez is, a referenciamodellhez hasonlóan, egy kéttárolós lengĘ tag diszkrét átviteli függvényének a nevezĘje. Folytonos idĘben az átviteli függvény:

Zo sm1 Zm

M

(tehát másodfokú referenciamodell)

wA o =wA+l-1-0=2

wR =2

Po (s)=

Re{s}

Ȧo2 s +2ȗ o Ȧo s+Ȧo2 2

A csillapítás ebben az esetben is legyen Ȣo= 0.707, a sajátfrekvencia pedig legyen Ȧo=2 . sm1

so1

A csillapítás (Ȣm= cosĳ, értékét Ȣm= 0.707 (Ȣm = cos45˚) választjuk, mivel ezzel 10%-on belüli túllendülés érhetĘ el. A gyorsítás szempontjából az Ȧm értéke nagyobb kell legyen (tehát a pólus “távolabb” kell essen a képzetes tengelytĘl). Legyen Ȧm=1.333 (aminek megfelel a Tm=1/Ȧm=0.75 idĘállandó). A folytonos átviteli függvény a következĘ formát ölti:

Pm (s)=

a

B m (z)=B cm (z)B - (z)=const B(z)/B(1) egyenletnek megfelelĘen alakul.

B - (z)=B(z)/B(1)

Ennek megfelelĘen

» » » » »


A Matlab-ban a a referencia modellt az alábbi parancsokkal számítjuk ki:

A folyamat egyetlen zérusa nem kompenzálható, mivel ugyan az egységkörön belül található, viszont ennek a bal félsíkjában, így kompenzálása lengéseket eredményezne a mintavételi pontok között. A mintavételi pontok közötti lengések elkerülésére és az egységnyi átvitel a B(z) polinomot felbontjuk a B ( z ) kompenzálható és a B ( z ) nem kompenzálható tényezĘkre.

B + (z)=B(1)


Ȧ 2m s +2ȗ m Ȧ m s+Ȧ 2m

» » » » » »

wo=2; % sajatfrekvencia To=1/wo; % ennek megfelelo idoallando ZetaO=0.707; numO=1; denO=[To^2 2*ZetaO*To 1]; Ho=tf(numO,denO);

A diszkrét idejĦ átviteli függvény:

Po (z)=

2

illetve az impulzusátviteli függvény:

7

B o (z) A o (z)

8




melynek nevezĘjébĘl származtatható az Ao(z) megfigyelĘ polinom.

ª 1 « « a 1 -1 « a 2 -a 1 « ¬« -a 2

» Hoz=c2d(Ho,Ts,'zoh'); » [numOz,denOz]=tfdata(Hoz,'v'); A szabályozott szakasz (p), a modell (m), és a megfigyelĘ polinom (o) paraméterei:

b -o b1-

0 b o-

0

b1-

0

0


0 º ª r º ª Ȗ1 +1-a 1 º »« » « » 0 » «s o » « Ȗ 2 +a 1 -a 2 » = - » b o « s1 » « Ȗ 3 +a 2 » »« » « » b1- ¼» ¬« s 2 ¼» ¬« Ȗ 4 ¼»

amibĘl Természetes frekvencia

Csillapítás

Zp=1,0 Zm=1,333 Zo=2,0

]p =0.25 ]m =0,707 ]o =0,707

Szabályozott szakasz Modell MegfigyelĘ

m

(z )A o (z )

diophantoszi egyenletet kell megoldani. A polinomiális egyenlet felírható a következĘ formában:

C(z)=A m (z)A o (z)=z 4 +Ȗ1z3 +Ȗ 2 z 2 +Ȗ 3 z+Ȗ 4 » Cz=conv(denMz,denOz); Mindkét oldalt a megfelelĘ fokszámú tagokkal kifejezve következik

(z 2 +a1z+a 2 )(z-1)(z+r)+(bo- z+b1- )(so z 2 +s1z+s2 )=z 4 +Ȗ1z 3 +Ȗ 2 z 2 +Ȗ 3 z+Ȗ 4

z4:

1=1

z3:

r-1+ a 1 + s o b o- = Ȗ 1

z2:

-r+a1 (r-1)+a 2 +s1bo- +s o b1- =Ȗ 2

z1:

-a1r+a 2 (r-1)+s 2 bo- +s1b1- =Ȗ3

z0:

-a 2 r+s 2 b1- =Ȗ 4

Vektoros formában a fenti egyenletek:

9

b 10 0

b ob 10

0º » 0» - » bo » b 1- »¼

-1

ª Ȗ 1 +1-a 1 º « Ȗ +a -a » « 2 1 2» « Ȗ 3 +a 2 » « » Ȗ4 «¬ »¼

annak érdekében, hogy statikus hiba egységugrás alakú alapjel esetén zérus legyen.

Amennyiben 1-es típusú kört kívánunk létrehozni, az

majd az együtthatókat összehasonlítva az alábbi egyenleteket kapjuk:

0

T (z)= B cm (z)A o (z)= co n st A o (z)= A m (1 )A o (z)

subplot(3,1,1); step(Hp,t),grid,title('Folyamat'); subplot(3,1,2); step(Hm,t),grid,title('Referenciamodell'); subplot(3,1,3); step(Ho,t),grid,title('Megfigyelo');

A ( z ) ( z -1 ) ( z + r ) + B - ( z ) ( s o z 2 + s 1 z + s 2 ) = A

b -o

végül a T(z) polinom

Nézzük meg mindhárom komponens átmeneti függvényét: » » » » » »

ªrº ª 1 «s » « « o » = « a 1 -1 « s1 » « a 2 -a 1 « » « «¬ s 2 »¼ «¬ -a 2

A másodfokú szakasz és a példában létezĘ fokszámfeltételek mellett a diophantoszi egyenlet a findSzRz program segítségével oldható meg: » findSzRz; A findSzRz program kódja a következĘ: » M=zeros(4,4); » M(1,1)=1; » M(1,2)=Bmz(1); » M(2,1)=Az(2)-1; » M(2,2)= Bmz(2); » M(2,3)= Bmz(1); » M(3,1)=Az(3)-Az(2); » M(3,3)= Bmz(2);; » M(3,4)= Bmz(1); » M(4,1)=-Az(3); » M(4,4)= Bmz(2); » v(1)=Cz(2)+1-Az(2); » v(2)=Cz(3)+Az(2)-Az(3); » v(3)=Cz(4)+Az(3); » v(4)=Cz(5); » x=inv(M)*v'; » r=x(1); » Sz(1:3)=x(2:4); » Rz=conv([1 r],[1 -1]); Tovább haladva a megadott egyenletek szerint: » Tz=denOz*(denMz*ones(length(denMz),1)); » RA=conv(Rz,Az); EllenĘrizzük a diophantoszi egyenlet megoldását, a baloldal A(z)R c(z)+B- (z)S(z) : » RA+[0 conv(Bmz,Sz)] és a jobboldal A o (z)A m (z) : » Cz összehasonlításával. Látható, hogy a két polinom azonos.

10



A szabályozó polinomjai az Rz, Sz, Tz változókban vannak tárolva (Matlab). A szimulációhoz Simulink diagramot használunk az alábbi ábra szerint.

Step1

Tz(z) Step

1 Rz(z) Discrete Filter

Transfer Fcn

Zero-Order Hold


Vizsgáljuk meg a megfigyelĘ polinom hatását a zárt kör alapjelkövetésére és zavarelhárítására. Ismételjük meg a számítást Z 0 8 mellett. Látni fogjuk, hogy a megfigyelĘ polinom módosítása nem változtatja meg a zárt kör alapjelkövetését, viszont hatással van a zavarelhárítási tulajdonságokra. Vizsgáljuk meg, hogyan viselkedik a szabályozási kör, ha a szakasz paramétereit r25% -kal megváltoztatjuk.

1 s2 +0.5s+1

1 Discrete Filter1


Scope

Discrete Filter2 Sz(z)

Kapcsolat a 2DF és a véges beállású tervezési módszer között

y_out

1

To Workspace1

A véges beállású tervezésnél a tervezés célspecifikációját

t Clock

Y ( z) R( z )

To Workspace

Implementálhatósági szempontokat figyelembe véve az Rz polinomot (mely tartalmazza az integrátort) csak egyszer szabad a szabályozási körbe elhelyezni, éspedig közvetlenül a szakasz elé. Rz, Sz és Tz implementálása diszkrét filterekkel történik (Simulinkben a Simulink / Discrete / Discrete Filter blokk segítségével, amely z-1 hatványaiban realizálja az átviteli függvényeket), ahogy a következĘ ábrán látható. r

e

rF

yF

y

u

1 1 (r 1) z 1 rz 2

t 0 t1 z 1 t 2 z 2

z k B ( z)

formában fogalmaztuk meg, ahol

k 1 wB ( z ) Nyilvánvaló, hogy a 2DF terminológia szerint és Am ( z ) z k Ezzel a modellel a beállás a korábbi szerényebb specifikációt kielégítĘ tervezésnél lényegesen dinamikusabb beavatkozást igényel.

Bmc ( z ) 1

Szakasz

s 0 s1 z 1 s 2 z 2

2 DOF

1.5

1

0.5

0

0

2

4

6

8

10

12

14

16

18

20

0

2

4

6

8

10

12

14

16

18

20

2

1.5

1

0.5

A szabályozott jellemzĘ és a beavatkozójel lefolyását egységugrás alapjelre a fenti ábra mutatja.

11

1 B ( z) zk

12

Állapotteres Leírás, Irányíthatóság, MegfigyelhetĘség


16. Állapotteres Leírás, Irányíthatóság, MegfigyelhetĘség Egy lineáris rendszer jellemezhetĘ bemenĘjeleivel, kimenĘjeleivel és állapotváltozóival. Az állapotváltozók olyan mennyiségek, amelyek a bemenĘjel hirtelen megváltozására nem változnak meg hirtelen, pillanatnyi értéküket a rendszert a múltban ért hatások alakítják ki. Másképpen a rendszer azon változói, amelyek tároló funkciót látnak el. Tekintsünk egy egy bemenetĦ, egy kimenetĦ (SISO - single input - single output) rendszert. A rendszer u bemenĘjele, x állapotváltozói és y kimenĘjele között az alábbi állapotegyenlet adja meg az összefüggést:

x

Ax bu

y

cx du

ahol az ^ A, b, c, d ` mátrixok a rendszert jellemzĘ paramétermátrixok. Az állapotegyenlettel egy n-edrendĦ differenciálegyenlettel megadható rendszert n számú elsĘrendĦ differenciálegyenlettel írunk le. Az állapotegyenletnek végtelen sok reprezentációja létezik, mivel az állapotváltozók bármilyen lineáris kombinációja új állapotváltozókat eredményez. A különbözĘ reprezentációk a bemenet és a kimenet között ugyanazt az átviteli kapcsolatot adják. A MATLAB az átviteli függvénybĘl kiindulva az állapotegyenlet egy lehetséges reprezentációját adja meg a tf2ss (transfer function to state space) utasítással.


Példa: » A=[-1 –0.5 0.5; 2 –3 0; 2 –1 –2] » b=[2;3;1] » c=[0 0 1] » d=0 » H=ss(A,b,c,d) Határozzuk meg a sajátvektorokat (V) és a sajátértékeket (ev). » [V,ev]=eig(A) Transzformáljuk a rendszert kanonikus alakra. Ekkor a sajátvektorokat megadó V mátrix inverze a transzformációs mátrix. » Pi=V; P=inv(V) » Ap=P*a*Pi » bp=P*b » cp=c*Pi » dp=d Figyeljük meg, hogy Ap valóban diagonális. A fenti számítások egyszerĦbben is elvégezhetĘk: » [Ap,bp,cp,dp]=ss2ss(A,b,c,d,inv(V)) vagy » [Ap,bp,cp,dp]=canon(A,b,c,d,’modal’) ª x1 º ª 1 0 0 º ª x1 º ª 1.73 º

« x » « 2» ¬« x3 ¼»

Példaként tekintsük az alábbi átviteli függvénnyel adott másodrendĦ rendszert:

H ( s)


«0 « ¬« 0

3

0 » « x2 » « 0 » u

0

2 ¼» ¬« x3 ¼» ¬« 2.23¼»

» « » « ª x1 º

1 s2 s 1

y

Adjuk meg a rendszer egy állapotteres reprezentációját. » num=1; » den=[1 1 1]; » [A,b,c,d]=tf2ss(num,den) % transzformálás állapotteres alakra » [num1,den1]=ss2tf(A,b,c,d) % ellenĘrzésképpen visszatranszformálás átviteli függvény alakra A kétféle megadási forma egyenértékĦségét nézzük meg az átmeneti függvények felrajzolásával. » step(num,den); » step(A,b,c,d); A rendszer egy másik állapotreprezentációját koordináta-transzformációval (hasonlósági transzformáció) kaphatjuk meg. Az új x állapotváltozók és az eredeti x állapotváltozók között a P transzformációs mátrix adja meg a kapcsolatot.

> 0.57

0@ « x2 » 0 u

0.707

« » ¬« x3 ¼»

Rajzoljuk fel a rendszer blokk-diagramját. Ez a rendszer párhuzamos reprezentációja. Figyeljük meg, hogy az átviteli tényezĘ értéke b1 c1 1.73 0.57 1 . Bármilyen b1 és c1 felbontással, amelyre

b1 c1 1 a kimenĘjel és a bemenĘjel között azonos összefüggést kapunk, miközben az állapotegyenlet paraméterei különböznek. 1.73

x1

³

x1

x2

u

x

Ax bu

³

y

cx du

-3

A

Px => x -1

PAP , b

Pb, c

P -1 x 1

cP , d

-2.23

d

Ha a P 1 mátrix oszlopvektorai az A mátrix sajátvektorai, az A mátrix diagonális lesz (kanonikus transzformáció).

0.57

-1

x2

ahol

x

»

x3

³ -2

x3

-0.707

y



A fenti párhuzamos struktúrából (az állapotegyenlet kanonikus alakjából) közvetlenül megadhatjuk a rendszer irányíthatósági és megfigyelhetĘségi tulajdonságait. Látható, hogy az u bemenĘjel nem befolyásolja az x2 állapotváltozót (az x2 állapotváltozó nem irányítható), az y kimenĘjel pedig nem tartalmaz információt az x3 állapotváltozóról (az x3 állapotváltozó nem megfigyelhetĘ). Az y kimenĘjel irányítható, mivel az u bemenĘjel az x1 állapotváltozón keresztül befolyásolja. Ezeket a tulajdonságokat az eredeti alakból kiindulva is megállapíthatjuk a Kalman-féle irányíthatósági és megfigyelhetĘségi mátrixok rangjának vizsgálatával.

Állapotvisszacsatolás


17. Állapotvisszacsatolás Tekintsük elĘször a folytonos rendszereket. Egy folytonos, lineáris, idĘinvariáns, több bemenetĦ/több kimenetĦ rendszer állapotteres reprezentációja az A, B, C, D paraméter mátrixokkal írható le:

x = Ax + Bu y = Cx + Du

A rendszer irányíthatósági mátrixa Co ª¬b Ab ... An 1 b ¼º számítható a ctrb utasítással. » Co=ctrb(a,b) vagy az LTI struktúrával: » Co =ctrb(H) A rendszer teljesen állapotirányítható, ha Co rangja megegyezik az állapotváltozók n számával.

Ismeretes, hogy a fenti egyenletek (nevezetesen az állapotegyenlet illetve a kimeneti egyenlet) az u bemenet és az y kimenet között a

» rank(Co) Esetünkben a rendszer nem teljesen állapotirányítható, mivel rank (Co ) 2 n . A rendszer kimeneti irányítható, ha a coy c * Co mátrix rangja megegyezik a kimenĘjelek számával.

névvel is illetik, amely modell pólusait a

» rank(c*Co) Mivel ez az érték 1, megegyezik a kimenĘjelek számával, a rendszer kimeneti irányítható.

A rendszer állapotváltozóiról megvalósított, az alábbi ábrán bemutatott visszacsatolás az r referencia jel (alapjel) és az y kimenĘjel közötti átvitelt, mind statikus, mind pedig dinamikus értelemben módosítja:

A rendszer megfigyelhetĘségi mátrixa

Ob

ª c º « cA » « » « ... » « n 1 » ¬ cA ¼

C ( s I A ) 1 B D átviteli függvényt valósítja meg. Az Â, B, C, D` négyessel adott rendszerleírást állapotmodell

det( sI A) 0 karakterisztikus egyenlet gyökei határozzák meg.

D

r

A MATLAB obsv utasításával számolva: » Ob=obsv(A,c) » Ob =obsv(H) » rank(Ob) A rendszer megfigyelhetĘ, ha a megfigyelhetĘségi mátrix rangja n. Mivel rank (Ob ) 2 n , a rendszer nem teljesen megfigyelhetĘ.

G

u

x

B

³

-

x

C

y

A K

Határozzuk meg a rendszer átviteli függvényét. » [num,den]=ss2tf(A,b,c,d) Az átviteli függvény zérus-pólus alakban: » Hzpk=zpk(H) Figyeljük meg, hogy két zérus kiejt két pólust. » Hzpk=minreal(Hzpk) Az átviteli függvény alakban információt vesztünk a rendszer állapotairól. Az átviteli függvény csak az irányítható és megfigyelhetĘ állapotváltozókról tartalmaz információt.

A visszacsatolásban K visszacsatolási mátrixot, az elĘrevezetĘ ágban G erĘsítési mátrixot feltételezve a visszacsatolt (zárt) rendszer egyenletei u = Gr - Kx következtében az alábbiak szerint írhatók fel:

A párhuzamos reprezentációt megkaphatjuk a részlettörtes felbontással is: » [num,den]=tfdata (H,’v’) » [r,p,k]=residue(num,den)

adódik a zárt rendszerre, amelynek karakterisztikus egyenlete:

H ( s) Megjegyzés: Látható, hogy az

0 0 1 s 3 s 2 s 1

1 átviteli függvény az állapotegyenlettel megadott teljes rendszernek s 1

Bevezetve az A F

A BK , B F

x = (A BK ) x + GBr . y = (C DK ) x + DGr GB, C F C DK , D F DG jelöléseket x = A F x + B F r y = CF x + DFr

det( sI A F ) det(sI - A + BK ) 0 . A nyitott és zárt rendszer karakterisztikus egyenletét összevetve látható, a nyitott rendszer pólusai csak A értékétĘl, míg a zárt rendszer pólusai az Â, B, K ` hármastól függenek. A zárt rendszer tulajdonságait úgy tervezzük meg, hogy elĘírjuk, hová kerüljenek a zárt rendszer pólusai a komplex síkban. Matematikailag olyan K visszacsatolási mátrixot keresünk, amely a karakterisztikus egyenlet gyökeit éppen e kívánt értékekre helyezi el.

csak egy részrendszerét írja le. Ez a részrendszer a teljes rendszer irányítható és megfigyelhetĘ része.

1



Az állapotvisszacsatolásos rendszer tervezése az alábbi lépésekben történik: - válasszuk meg a zárt rendszer pólusait - egy bemenetĦ/egy kimenetĦ rendszerek esetén a K visszacsatolási mátrixot az Ackermann formula segítségével határozhatjuk meg (a Matlab acker parancsának közvetlen alkalmazásával) - az állandósult állapotbeli átvitelre megfogalmazott követelmények alapján határozzuk meg a G erĘsítési mátrixot.



» GainOL=dcgain(A,b,c,d) GainOL = 1 » GainCL=dcgain(A-b*k,b,c,d) GainCL = 0.04 » G=GainOL/ GainCL » subplot(212),step(A-b*k,G*b,c,d) Az alábbi ábra a nyitott és zárt rendszer átmeneti függvényét mutatja be, a zárt rendszer viselkedését G 1 , valamint G GainOL / GainOC 25 értékre is láthatjuk:

Példa. Tekintsük az alábbi átviteli függvénnyel adott rendszert:

1 0.8

6 H (s) ( s 1)( s 2)( s 3) A választott rendszer egy bemenetĦ/egy kimenetĦ, így a K visszacsatolási mátrix ebben az esetben visszacsatolási vektort, a G erĘsítési mátrix pedig erĘsítési vektort fog jelenteni. Vegyük észre, hogy a rendszer statikus erĘsítése egységnyi. A rendszert Matlab szinten adjuk meg p o pólusaival, majd

0.4 0.2 0

A

11 6 º 0 0 »» , B 1 0 »¼

>6

3 4 j

1

2

3

4

0.5

ª1 º «0» , C « » «¬ 0 »¼

output without gain compensation

>0

0

6@, D

0

0

po

>1

2

1

2

3

4

k [0,..., 0,1] C D(c A ) C x a nyitott rendszer irányíthatósági

mátrixa, D c ( s ) a zárt rendszer karakterisztikus egyenlete (amelyet az elĘírt pólusok határoznak meg!),

D c ( A) pedig e polinom, mint mátrix polinom helyettesítési értéke az A helyen. Matlab szinten mindez egyetlen utasítás:

6

Integrátor beiktatása a visszacsatolási hurokba Ami a zárt kör dinamikus tulajdonságait illeti, az állapotvisszacsatolás önmagában a korábbiakban frekvenciatartománybeli megfontolások alapján származtatott PD kompenzáció alkalmazásánál tapasztalt hatást idézi elĘ. A statikus viszonyokat ugyanakkor egy a visszacsatolási hurkon kívül elhelyezett erĘsítés határozza meg, amelyet a rendszerparaméterek ismeretében határoztunk meg. Ez azt jelenti, hogy az erĘsítés meghatározása érzékeny a rendszerparaméterek ismeretének pontosságára. Ezen túlmenĘen, a kimenetet érĘ zavaró hatásokat hurkon kívül elhelyezett elemekkel nem lehet kompenzálni. Következésképpen a statikus hiba elhárítására - a frekvenciatartománybeli megfontolásokhoz hasonló módon - integrátort iktatunk be a hurokba. A beiktatott integrátort is tartalmazó megoldást az alábbi ábra részletezi:

» k=acker(A,b,pc)

50 144 @

» t=0:0.1:6; » subplot(211),step(A,b,c,d,1,t); » subplot(212),step(A-b*k,b,c,d,1,t);

D

xi

r

-

³

xi

ki

u

B

x

³

-

A

Látható, hogy a pólusok balra tolásával az átmeneti függvény tranziense gyorsabb lefutásúvá vált, ugyanakkor – feltételezve, hogy a zárt körben is egységnyi statikus átvitelt szeretnénk biztosítani - a zárt rendszer állandósult állapotbeli viselkedése nem elfogadható. A megfelelĘ állandósult állapot bizosításának egy lehetséges módja a G erĘsítési vektor beállítása. Az eljárás részletei jól követhetĘk az alábbi Matlab utasítás sorral, ahol GainOL a nyitott kör, GainCL pedig a zárt kör statikus erĘsítését jelöli: 2

5

3@

3 4 j @ pólusaiba helyezi át. Az Ackermann

formában állítja elĘ a keresett visszacsatolási vektort, ahol

0

Az ábrán jól látható a dinamikus tulajdonságoknak a pólusáthelyezés következtében elĘálló javulása.

1 x

>6

6

output with gain compensation

formula analitikus alakja

k

5

1

Jelöljük ki a zárt rendszer p c pólusait: » pc=[-6;-3+i*4;-3-i*4]; majd határozzuk meg azt a visszacsatolási vektort, amely a nyitott rendszer pólusait a zárt rendszer pc

0

1.5

transzformáljuk át a rendszert állapotteres formába: » po=[-1;-2;-3]; » [A,b,c,d]=tf2ss(6,poly(po))

ª 6 «1 « «¬ 0

plant output

0.6

K

3

x

C

y



Látható, hogy az integrátor xi kimenetét állapotváltozóként figyelembe véve az állapotváltozók száma eggyel növekedett. Egy bemenetĦ/egy kimenetĦ rendszert és D = 0 közvetlen input-output átvitelt feltételezve a zárt rendszer egyenletei egyetlen vektor-mátrix egyenletbe foglalhatók:

ª x º « x » ¬ i¼

x e

bki º ª x º ª 0 º r 0 »¼ «¬ xi »¼ «¬ 1 »¼

ª A - bk « -c ¬

A cxe bcr ,

ahol A c tovább részletezhetĘ az alábbiak szerint:

Ac

ª A - bk « -c ¬

bki º 0 »¼

ªA « -c ¬

0 º ªb º >k 0 »¼ «¬ 0 »¼

ki @

A e b ek e

bc = 0 0 0 1 » cc=[c 0] cc = 0 0 6 0 » dc=0; » subplot(111),step(Ac,bc,cc,dc,1,t)

0

1.4

karakterisztikus egyenlet gyökei határozzák meg, így a kibĘvített rendszer megválaszott p ce pólusait

1.2

biztosító k e visszacsatolási vektort ismét az Ackermann formulával számíthatjuk ki. Ehhez azonban elĘször a kibĘvített rendszer paramétermátrixait kell elĘállítani:

>15

1

0.8

0.6

» z=[0;0;0]; » Ae=[A z;c 0] Ae = -6 -11 -6 0 1 0 0 0 0 1 0 0 0 0 6 0 » be=[b;0] be = 1 0 0 0 Jelöljük ki a zárt rendszer pólusait: » pce=[-9 -6 -3+i*4 -3-i*4]; » ke=acker(Ae,be,pce)

ke


A zárt kör átmeneti függvénye:

A zárt rendszer pólusait a

det( s I A e b e k e )


0.4

0.2

0

0

1

2

3

4

5

6

Látható, hogy a zárt kör dinamikus és statikus viselkedése is megfelelĘ.

Állapotbecslés

158 693

225

@

A kibĘvített erĘsítési vektor elsĘ három eleme az eredeti állapotváltozóról történĘ visszacsatolást valósítja meg, míg a negyedik k i elem a mesterségesen bevezetett integrátorhoz tartozik: » k=ke(1:3) k = 15 158 693 » ki=ke(4) ki = 225 » Ac=[A-b*k b*ki;-c 0] Ac = -21 -169 -699 225 1 0 0 0 0 1 0 0 0 0 -6 0 » bc=[z;1]

4

A gyakorlati alkalmazások döntĘ részében a folyamatok mĦszerezése magába foglalja a rendszer kimenetének mérését, azonban az állapotváltozók közvetlen méréssel nem hozzáférhetĘk. Ilyen esetekben az állapotvisszacsatolásos szabályozási technika alkalmazhatóságát a nem mért állapotváltozók becslése kell, hogy kiegészítse. Az állapotváltozók becslésére szolgáló ún. Kalman szĦrĘ blokkvázlata az alábbi ábrán látható. A Kalman szĦrĘ részét képezĘ modell felépítése követi a rendszer felépítését, a rendszer és a modell kimenetének különbsége pedig egy hibajelet állít elĘ. E hibának a modellbe történĘ visszacsatolása biztosítja azt, hogy a modell állapotváltozói minél kisebb hibával kövessék (másolják le) a rendszer állapotváltozóit. A becslĘ hálózat Simulink-ben könnyen realizálható, ugyanis a becslĘ hálózat elĘírt pólusai alapján az L erĘsítés ismét az Ackermann formulával számítható. A becslĘ hálózat pólusait úgy szokásos megválasztani, hogy a zárt rendszer tranziens viselkedését az állapotbecslés túlságosan ne lassítsa le, ugyanis a zárt rendszer állapotvisszacsatolása a becsült állapotváltozó felhasználásával történik.

5



u

x

B

x

³

y

C

A

Model

D x

B

x

³

y

C


A szakasz legyen a szabályozási példánál is vizsgált háromtárolós arányos szakasz (az integráló állapotváltozóval való kiterjesztés nélkül). Legyen mindhárom állapotváltozó kezdeti értéke 1. Az alapjel és a zavarás értéke legyen zérus. A becslĘ hálózat elĘírt pólusai legyenek » Se=[-7 -7 -7] Ezzel L’=[-17.3333 7.66672.5000] érték adódik. Szimuláljuk az állapotbecslés lefolyását, majd a Matlab felületrĘl ábrázoljuk az elsĘ állapotváltozó és becsült értékének idĘbeli lefolyását. » plot(t,[X(:,1),Xe(:,1)]),grid A szimuláció az állapotváltozók gyors beállását mutatja.

System

D


1

-

e

0

A

X1

-1

-2

L -3 Xe1

A becslĘ hálózat Simulink modelljét az alábbi ábra adja meg.

-4

-5

K uz1

uz

X

x' = Ax+Bu y = Cx+Du Sum 2

Bp

K Cp2

Ap,I,I,zeros(3,3 )

Sum 3

Scop e

K L x' = Ax+Bu y = Cx+Du Sum 4 t Clock

0

0.5

1

1.5

2

2.5

3

3.5

4

y y

K ua

-6

X

Bz

Su m

Feladat: Állítsuk a kezdeti feltételek értékét zérusra, és az uz zavarás értékét 1-re. A szimulációt lefuttatva látható, hogy statikus eltérés lép fel a tényleges és becsült állapotváltozók között. Ennek kiküszöböléséhez a zavarás állapotváltozóival célszerĦ kibĘvíteni a szakasz állapotegyenletét és az állapotbecslést a kibĘvített rendszerre ajánlatos elvégezni. (Erre azonban a továbbiakban nem térünk ki.)

K

Ak,I,I,zeros(3,3 )

Cp

Xe Xe

t

A szakaszt és modelljét a Simulink Continuous könyvtárának State-Space valamint Math könyvtárának Matrix Gain blokkjaiból építjük fel. A C paraméter különválasztására azért van szükség, mert nemcsak a kimenĘjelhez, hanem az állapotváltozókhoz is hozzá akarunk férni. A B paramétert is különválasztjuk az állapotmodell blokktól, mivel a modell állapotváltozóinak deriváljait módosítjuk, így a deriváltakhoz is hozzá kell férnünk.

Állapotvisszacsatolás a becsült állapotváltozókról Az állapotbecslés és állapotvisszacsatolás egymástól függetlenül elvégezhetĘ feladatok (szeparációs elv). Ha az állapotváltozók nem hozzáférhetĘk, az állapotvisszacsatolásos szabályozást a becsült állapotváltozók visszacsatolásával valósíthatjuk meg azzal a K visszacsatoló mátrix-szal, amellyel a hozzáférhetĘ állapotváltozókat csatolnánk vissza. A Simulink vázlatot az alábbi ábra mutatja. Kezdeti feltételekre és alapjelre az alapjel statikus kompenzációjával a szabályozás jól mĦködik. A zavarások hatását azonban csak statikus hibával tudja kiküszöbölni.

A becslĘ hálózat paraméterei (az L vektor elemei) a zárt becslĘ kör karakterisztikus egyenlete gyökeinek elĘírásával az Ackermann képlettel, a Matlab » L=acker(A’,C’,Se)’ utasításával határozhatók meg. Se a becslĘ kör elĘírt pólusainak vektora. A becslĘ kör legyen gyorsabb, mint a szakasz, és legyen gyorsabb, mint a szabályozás az állapotvisszacsatolással. (A transzponálásra a mátrixok, vektorok méreteinek egyeztetése miatt van szükség.)

6

7



x' = Ax+Bu y = Cx+Du

gain K ua

Gain Sum 5

uz2

To Workspace1

Bz

K Su m

e AT(x iT ) A(1 e AT )I Bu ( iT )

x([ i 1) T] xi 1

Xe

-Kb

állandó u (iT ) , és így kiemelhetĘ az

A mintavételezett rendszer állapotegyenlete:

yi

To Workspace

t

u (W )

vége ekkor t (i 1)T . Ebben a tartományban integrálás elé. Az integrálást elvégezve

K Cp

Ap,I,I,zeros(3,3 ) K

Clock

t0

A mintavételezés során alkalmazzunk zérusrendĦ tartószervet. Legyen t0 egy mintavételi intervallum kezdete, t0 iT , ahol T a mintavételezési idĘ, az intervallum

Sum 3

Cp1


x(t ) e A (t t0 ) x(0) ³ e A ( t W ) Bu (W )dW

outpu t

L

Sum 4

t

Scope 1 y

K

Ap,I,I,zeros(3,3) ,

Sum 2

Bp 1


Induljunk ki a folytonos állapotegyenlet általános megoldásából:

X

K uz1


Ad xi Bd ui Cd xi Dd ui

ahol

tim e

Ad

Scop e

Futtassuk le azonos kezdeti értékekre 1[ 0 0] és egységugrás alapjelre a Simulink programot a visszacsatolást a tényleges, illetve a becsült állapotváltozókról megvalósítva. A kimenĘjeleket (illetve az állapotváltozókat) ábrázoljuk egy diagramban. Látható, hogy a becsült állapotváltozókról való visszacsatolással a túllendülés kissé nagyobb. Hasonlítsuk össze a beavatkozójeleket is.

A1( e AT I) B;

e AT ; Bd

C ;

D=d .D

det( zI A d ) 0 karakterisztikus egyenlete a z tartományban adja meg a rendszer pólusait. Az Ackermann formula a megoldás egyik kulcseleme volt a folytonos rendszereknél, vizsgáljuk meg ezt az összefüggést mintavételes rendszerekre. A nyitott rendszer diszkrét állapotegyenlete:

1.4

xi 1

1.2

Cd

Ismeretes, hogy a gyakorlati rendszerekre tipikusan Dd=0 áll fenn, továbbá a nyitott rendszer

Ad xi Bd u i

y i C d xi Állapotvisszacsatolás: u i Kxi A zárt rendszer állapotegyenlete: xi 1 ( Ad Bd K ) xi

1

0.8

A zárt rendszer elĘírt karakterisztikus polinomja:

0.6

D c ( z ) det( zI ( A d B d K )) ( z pc1 )( z pc 2 )...( z pcn ) pc1 , pc 2 ,..., pcn a zárt rendszernek a z tartományban elĘírt pólusai. Legyen a szakasz irányíthatósági mátrixa: M c

0.4

ahol

0.2

0

0

1

2

3

4

5

6

Mc

ª Bd ¬

Ad Bd

... Ad

n 1

Bd º¼ ,

n dim x

Az Ackermann képlet szerint a szakasz Ad, Bd mátrixaiból és a zárt rendszer elĘírt pólusaihoz tartozó karakterisztikus polinomból a K visszacsatoló vektor az alábbi összefüggéssel számítható:

Feladat: É pítsünk Simulink vázlatot az integráló állapotváltozóval kib Ęvítve a rendszert, állapotbecsléssel és a becsült állapotváltozókról történĘ visszacsatolással.

Állapotvisszacsatolás mintavételes rendszerekre Vizsgáljuk meg az állapotvisszacsatolást mintavételes rendszerekre. A folytonos rendszerek körében folytatott vizsgálataink azt mutatták, hogy a szabályozási kör pólusai a folyamat állapotváltozóinak a folyamat bemenetére való visszacsatolásával elĘírt értékekre állíthatók. Most az állapotvisszacsatolást mintavételes rendszerre kívánjuk alkalmazni.

8

K

>0

... 0 1@ M c1D c ( Ad )

ahol D c ( Ad ) a zárt rendszer karakterisztikus polinomja z A d helyettesítéssel. Az állapotokat visszacsatoló K vektor kiszámítását a Matlab az acker utasítással támogatja: K=acker(Ad,Bd,Sd) Sd a visszacsatolt rendszer elĘírt pólusait (az D c ( z ) 0 ) egyenlet elĘírt gyökeit) tartalmazó vektor. Ha az elĘírt pólusokat folytonos idĘben adjuk meg, ezek diszkrét megfelelĘit a z=esT transzformációval határozhatjuk meg, ahol T a mintavételezési idĘ.

Példa Legyen a folytonos szabályozott szakasz az alábbi háromtárolós arányos tag:

9



u

A 1 sT3

x3

x2

1 1 sT2

1 1 sT1

x1=y

A paraméterek értékei: A=1, T1=0.2, T2=2, T3=4. Az állapotváltozókat vegyük fel az ábra szerint. Adjuk meg a szakasz folytonos állapotváltozós egyenleteit, majd T=0.2 mintavételezési idĘ mellett határozzuk meg a diszkrét állapotegyenletet zérusrendĦ tartószerv feltételezésével.


Alkalmazzuk az Ackermann formulát az állapotvisszacsatoló vektor meghatározásához. » K=acker(Ad,Bd,Sd) Az állapotvisszacsatolásos kompenzálás Simulink blokk diagramja az alábbi ábrán látható. A kapcsolást kiegészítettük a folytonos szakasz modelljével is, hogy a mintavételi pontok közötti viselkedés is vizsgálható legyen. Vizsgáljuk a szabályozás mĦködését [1 1 1]’ kezdeti feltételek lecsengetésére illetve egységugrás alapjel követésére gain=1 beállítás mellett. A zavarás a második állapotváltozó deriváltjára hat, Bz=[0 1 0]’. Lá tható, hogy a szabályozás a kezdeti feltételeket az elĘírt dinamika szerint lecsengeti, azonban az alapjel követésében jelentĘs statikus hiba mutatkozik. A zavarás hatását sem tudja elhárítani.

a./ Tervezzünk állapotvisszacsatolást a diszkrét zárt rendszer pólusainak elĘírásával. A pólusokat folytonos idĘben specifikáljuk, majd a z e sT transzformációval térjünk át diszkrét idĘre. Bevezetjük a Tsum T1 T2 T3 összeg idĘállandót (ennek az A átviteli tényezĘvel megszorzott értéke adja a szakasz átmeneti függvénye és annak állandósult értéke eltérésének integrálját, az ún. lassító területet). Ehhez képest gyorsítjuk fel a rendszert. A zárt rendszer elĘírt pólusai legyenek:

s1, 2 egy kéttárolós lengĘ tag pólusai, Z 0 s3

5 Tsum

; ]

0.7 .


K B

z gain

Bz

+

u control

max{Z 0 ,1 / min{T1 , T2 , T3 )}

+

K

+ Sum5


Sum6

+

Gain

r

+

K

+

Bd

Sum4

control

Megoldás

ª x1 « x « 2 «¬ x 3 y

>1

Discrete State-Space

Cd

K

states

yd yd

output

x x

Statikus kompenzáció az alapjel követéséhez

0 0@X 0 u

A zárt rendszerre elĘírt pólusok: s1, 2

foutput K

td time

Clock

5 0 º ª x1 º ª 0 º ª 5 « 0 0 .5 0 .5 »» «« x 2 »» «« 0 »»u « «¬ 0 0 0.25»¼ «¬ x3 »¼ «¬0.25»¼

º » » »¼

C

x(n+1)=Ax(n)+Bu(n) y(n)=Cx(n)+Du(n)

Kd

Vizsgáljuk a rendszer mĦködését kezdeti feltételekre, illetve alapjelkövetésre és zavarelhárításra.

A folytonos szakasz állapotegyenlete:

A,I,I,zeros(3,3)

yf yf

K

]Z 0 r jZ 0 1 ] 2

0.5645 r j 0.57 6

s3 5 . » S=[-0.5645+i*0.576,-0.5645-i*0.576,-5] Diszkrét megfelelĘje Ts=0.2 mintavételezési idĘvel: » Ts=0.2 » Sd=exp(S*Ts) A kapott értékek: 0.8873 + 0.1027i, 0.8873 - 0.1027i , 0.3679.

Az állapotvisszacsatolásos kompenzációval r ugrásalakú alapjel esetén szeretnénk elérni, hogy állandósult állapotban az y kimenĘjel megegyezzen az alapjel konstans értékével. Ekkor az állapotváltozók deriváltjai zérussal egyenlĘk. Az alapjel egy k1 korrekciós erĘsítési tényezĘn keresztül hat a szabályozás bemenetére. Az összefüggést SISO (egy-bemenetĦ, egy-kimenetĦ esetre írjuk fel. Állandósult állapotban az állapotváltozók n+1-edik értéke megegyezik az n-edikkel. A kimenĘjel állandósult értéke pedig megegyezik az alapjellel.

rz0 xf

Ad xf Bd u f

uf

k1 r Kxf

yf

C d xf

xf

( I Ad Bd K ) 1 Bd k1 r

yf

C d ( I Ad Bd K ) 1 Bd k1 r

r

ahonnan

Feladat: A várható gyorsítás szemléltetésére ábrázoljuk a szakasz és az elĘírt pólusokkal rendelkezĘ rendszer átmeneti függvényeit. Alkalmazzuk a step ill. dstep MATLAB utasításokat. Határozzuk meg a diszkrét állapotegyenletet: » [Ad,Bd,Cd,Dd]=c2dm(A,B,C,D,Ts,'zoh')

10

r

és a korrekciós tényezĘ:

k1

1 /(C d ( I Ad Bd K ) 1 Bd )

11



Példánkban »k1=1/(Cd*inv(eye(3)-Ad+Bd*K)*Bd, értéke gain=k1=5.007. Ezzel a statikus kompenzációval újra futtatva a programot a kimenĘjel pontosan beáll az alapjel értékére. A zavarás elhárításánál azonban megmarad a statikus hiba.



Integráló szabályozás állapotbĘvítéssel

Az integráló állapotváltozó visszacsatolása: » Ki=kk(4) Szimuláljuk a rendszer mĦködését a Simulink program futtatásával kezdeti feltételekre, illetve alapjelre és zavarásra. Látható, hogy a szabályozás az alapjelet statikus állapotban hiba nélkül követi (statikus kompenzáció nélkül is), és a zavarás hatását is kiküszöböli.

Az ugrásalakú alapjel követésére, a zavarás hatásának csökkentésére célszerĦ a szabályozóban integrátort elhelyezni. BĘvítsük a rendszert egy új állapotváltozóval, amely a kimenĘjel integrálja.

Állapotbecslés

y

xI ,i Tyi

xI ,i 1

Tz 1 1 z 1

xI

Ha az állapotváltozók nem hozzáférhetĘk, nem mérhetĘk, azokat becsülnünk kell. A becsléshez a Kalman szĦrĘnek megfelelĘ kapcsolást alkalmazhatjuk az alábbi ábrán látható Simulink program szerint.

xI ,i TCxi

X X

A bĘvített állapotegyenlet:

ª xi 1 º «x » ¬ I ,i 1 ¼

ª Ad 0 º ª xi º ª Bd º «TC I » « x » « 0 » u ¬ ¼ ¬ I ,i ¼ ¬ ¼

r

Ts

-

1-z -1

sum1

Filter

Kint

+ + sum2

x(n+1)=Ax(n)+Bu(n) y(n)=Cx(n)+Du(n) K Bd

+

Discrete State-Space1

Bp

Discrete State-Space1

+

Cd

Sum2

y

Mux + K

+ +

t

output Clock

time

-

L x(n+1)=Ax(n)+Bu(n) y(n)=Cx(n)+Du(n)

Sum3 Discrete State-Space2

K Cd

K

Mux

Graph

Sum K Cd

Xe Xe

sum3 K

u control

Clock

ua

y y

+

Bz Ki

x(n+1)=Ax(n)+Bu(n) y(n)=Cx(n)+Du(n)

K

K

+

y +

Az új állapotváltozóval kibĘvített rendszer állapotvisszacsatolásának konstansai az Ackermann formulával kiszámíthatók. A módosított irányítási struktúrának megfelelĘ Simulink vázlat az alábbi ábrán látható.

z

uz1

td time

Példánkban a módosított állapotmátrixok: » A1=[Ad zeros(3,1);Ts*Cd 0] » B1=[Bd;0] Legyen a zárt rendszer elĘírt pólusainak vektora: » S1=[Sd;0.3679] Alkalmazzuk az Ackermann formulát a kibĘvített rendszerre: » kk=acker(A1,B1,S1) A visszacsatoló mátrix elemei: 5[ .0070 17.5259 13.7103 15.8251] Az eredeti állapotváltozók visszacsatolása: » Kc=kk(1:3)

12

-Kc

states x states

Ha a szakasz ismert, modelljét megépítjük. A szakaszt és modelljét ugyanazzal a bemenĘjellel gerjesztjük. A szakasz és a modell kimenĘjeleit összehasonlítva a hibajelet használjuk fel a modell állapotváltozóinak beállítására az L paramétervektoron keresztül. Akkor várható, hogy a modell állapotváltozóinak értéke gyorsan megközelíti és követi a tényleges állapotváltozók értékét, ha a becslĘ hálózat jóval gyorsabb a szakasz tranzienseinél. ElĘírhatjuk a becslĘ hálózat pólusait, és az Ackermann formula alkalmazásával meghatározhatjuk az L vektor paramétereit. A szakasz legyen a szabályozási példánál is vizsgált háromtárolós arányos szakasz. Legyen mindhárom állapotváltozójának kezdeti értéke 1. Az alapjel és a zavarás értéke legyen zérus. A becslĘ hálózat elĘírt pólusai legyenek valósak és azonos értékĦek, és az állapotváltozókról visszacsatolt hálózat leggyorsabb idĘállandójánál 3-szor gyorsabb viselkedést írjanak elĘ (konjugált komplex pólusok esetén tekintsük az abszolút érték reciprokát). A folytonos zárt rendszer elĘírt pólusai voltak: » S=[-0.5645+i*0.576,-0.5645-i*0.576,-5] » abs(S)=[0.8065 0.8065 5] Legyen tehát » Se=[-15 -15 -15]. Diszkrétben » Sed=exp(Ts*Se) A becslĘ hálózat paraméterei (az L vektor elemei) az » L=acker(Ad’,Cd’,Sed)’ utasítással határozhatók meg. A kapott értékek: 2[ .0746 3.1280 12.8569]'

13



A szimuláció az állapotváltozók gyors beállását mutatja. Ábrázoljuk a tényleges és becsült állapotváltozók lefolyását egy diagramban. » plot(t,X,t,Xe),grid Állítsuk a kezdeti feltételek értékét zérusra, és az uz1 zavarás értékét 1-re. A szimulációt lefuttatva látható, hogy statikus eltérés lép fel a tényleges és becsült állapotváltozók között. Ennek kiküszöböléséhez a zavarás állapotváltozóival célszerĦ kibĘvíteni a szakasz állapotegyenletét és az állapotbecslést a kibĘvített rendszerre ajánlatos elvégezni (Erre azonban itt nem térünk ki.) A szabályozást a becsült állapotváltozók visszacsatolásával valósíthatjuk meg.

14

1. Bevezetés a Matlab Használatába

Recommend Documents