II. Folytonos lineáris rendszermodellek időtartományi analízise MATLAB környezetben

II.1. Elméleti alapfogalmak

Egy folytonos lineáris egyváltozós (egy bemenettel és egy kimenettel rendelkező) rendszer matematikai modelljét az időtartományban egy n-ed rendű lineáris (valós, állandó együtthatójú) differenciál egyenlettel fejezzük ki. Ennek általános alakja:

     (II.1)

ahol:   - u(t) a rendszer bemenőjele (az ok, a gerjesztés),

- y(t) a rendszer kimenőjele (az okozat, a hatás)

- , a differenciál egyenlet együtthatói,

- , a fizikai megvalósíthatóság feltétele

- , a k-ad rendű (időszerinti) deriválási művelet jele.

Ugyanazt a lineáris rendszert az operátortartományban a W(s) átviteli függvény segítségével jellemezzük. Ez utóbbit a lineáris modell diferenciálegyenletéből kapjuk meg úgy, hogy alkalmazzuk rá a Laplace transzformációt, zérus kezdeti feltételek mellett, majd a kimenőjel és a bemenőjel Laplace transzformáltjainak hányadosát fejezzük ki.. Ennek eredményeként következik, hogy:

                                                 (II.2)

vagyis, az átviteli függvény egy polinomiális törtfüggvény s-ben, ahol s a Laplace változó. Az átviteli függvény alakját úgy is megkaphatjuk, ha a dinamikus viselkedést leíró differenciál egyenletben a d/dt differenciáló műveletet egy s operátorral helyettesítjük, majd a ki- és bemenő jelek arányát képezzük.

Az átviteli függvény a lineáris dinamikus rendszerek legelterjedtebb, legáltalánosabb kifejezési módja. A közvetlen kapcsolat miatt, a lineáris rendszermodellek analízise során általában ebből a kifejezésmódból származtatják az összes többi leírási módot, mint pl. az időtartományban a diferenciálegyenletet, a frekvenciatartományban a frekvenciafüggvényt, stb. Az átviteli függvény nevezőjét karakterisztikus polinom-nak, annak gyökeit pedig pólusok-nak nevezzük. A számláló polinom gyökeit zérusok-nak nevezzük.

Az időtartományi rendszerleírásban a rendszer válaszát az y(t) kimenőjel (kimeneti változó) jelenti, az u(t) függvény (bemenőjel vagy bemeneti változó) ismeretében. Bár a választ általában meg lehet határozni analitikusan is, mégis a grafikus ábrázolás áll a legtöbb esetben a rendelkezésünkre. A folytonos lineáris rendszermodellek tanulmányozására, az általános eset kivételével, amikor u(t) egy tetszőleges alakú jel (időfüggvény), jellegzetes, ún. típikus vizsgálójeleket alkalmazunk. Ilyen vizsgálójelek az egységugrás jel, az egységimpulzus jel illetve a harmonikus függvény jel. Ezekre a jelekre a rendszermodell válaszjele egyike a következőknek:

1.    Átmeneti függvény: a rendszer válasza az időtartományban egységugrás bemenő jelre, vagyis ha , zérus kezdeti feltételek mellett;

2.    Súlyfüggvény: a rendszer válasza az időtartományban egységimpulzus bemenő jelre, vagyis ha , zérus kezdeti feltételek mellett;

3.    Válasz szinuszos gerjesztésre: a rendszer válasza az időtartományban szinuszos bemenő jelre, vagyis amikor  , .

Elemzéseink elvégzéséhez a fenti tipikus válaszjelek grafikus ábrázolásából és esetleges további tanulmányozásából általában elégséges információhoz jutunk.

Arányos jellegű rendszerek esetén a kimenő jel a tranziens folyamat idején is arányos a bemenő jellel. Az ilyen rendszermodellek átviteli függvénye egy állandó (az arányossági tényező) érték, azaz:

                                                                                                                     (II.3)

Az átmeneti függvénykép egy, a vízszintes tengellyel Ka értéknél húzott párhuzamos egyenes.

Egy elsőfokú differenciál egyenlet által leírt rendszermodell (egytárolós arányos vagy PT1 tag) esetén az átviteli függvény analitikus kifejezése:

                                                                                                             (II.4)

ahol:   - Ka, az arányosságot kifejező arányossági tényező,

- T1, a tárolós tulajdonságot (tehetetlenséget) kifejező időállandó.

Ebben az esetben az átmeneti függvény egy, az arányossági tényezőnek megfelelő vízszintes aszimptótához közelítő exponenciális görbe. Az origóban az exponenciális görbéhez húzott érintő az aszimptótán a tárolós tulajdonság (tehetetlenség) időállandóját határozza meg.

Egy másodfokú differenciál egyenlet által leírt rendszermodell (kéttárolós arányos vagy PT2 tag) esetén az átmeneti függvény szintén az arányossági tényezőnek megfelelő vízszintes aszimptótához tart. A tranziens átmenetet jellemezhetjük egy paraméter, a ζ csillapítási tényező, segítségével. Ennek értékeitől függően az átmeneti függvény időbeni változása:

-           esetén, csillapított, aperiodikus felfutású,

-     esetén, csillapított, periodikus lengésekkel rendelkező,

-          esetén, csillapítatlan, állandó amplitudójú lengésekkel rendelkező,

- esetén, csillapítatlan, növekvő amplitudójú lengésekkel rendelkező,

-        esetén, csillapítatlan, aperiodikusan növekvő amplitudóval rendelkező.

A valós rendszerekre jellemző csillapított tranziens viselkedés matematiikai leírására használt átviteli függvény analitikus kifejezésére két különböző, de egyenértékű változattal találkozunk, attól függően, hogy a ζ csillapítási tényező értéke nagyobb vagy kisebb egynél. A két kifejezésmód:

,  ha                                                                    (II.5.a)

,          ha                                                                    (II.5.b)

Az első esetben a karakterisztikus polinom két gyöke (pólusok) negatív valós értékek. Tehát a kéttárolós arányos tulajdonság két egytárolós arányos tag soros kapcsolatának köszönhető. A két tárolós tulajdonság (tehetetlenség) jellemzésére szolgál a két időállandó (T1 és T2).

A második képletben a karakterisztikus polinom gyökei (pólusok) negatív valós részű komplex konjugált értékek. Ebben az esetben a kéttárolós jelleget két másik paraméter segítségével szoktuk kifejezni. Az egyik a már említett ζ csillapítási tényező, amely a csillapodó tranziens mozgás során fellépő energiapazarló (energiaveszteséget okozó) erők hatásával hozható kapcsolatba, illetve a csillapodó lengések amplitudóit összekötő exponenciális burkológörbe analitikus kifejezésében jelenik meg. A másik paraméter a T0, és a saját periodusidő nevet viseli. Ezzel a paraméterrel jellemzzük azoknak az ideális helyzetben fellépő lengéseknek periodusidejét, amikor a rendszerünk mozgásában nincs energiaveszteség, azaz csillapítási tényezőnk nulla. Ilyenkor tipikusan a magára hagyott másodfokú rendszer állandó amplitudóval és a saját periodusidővel jellemzett lengő (rezgő) mozgást végez. Valós rendszerek esetén a csillapítás csökkenti nem csak a harmonikus lengések amplitudóját, hanem azok periodusidejét is. Természetesen, a valós rendszer jellege határozza meg a lengéseknek fizikai (mechanikai, villamos, termikus, stb) mennyiségekben való megnyilvánulását is.

Az egyszerű integráló jellegű rendszerek esetén, a kimenő jel lineáris függvénye a bemenő jelnek. Az ilyen rendszereket leíró matematikai modell differenciálegyenlete:

        vagy                                                                     (II.6.a)

A fenti rendszermodellhez tartozó átviteli függvény

 ,                                                                                                        (II.6.b)

ahol Ki a bemenő jelváltozás és a kimenő jel közötti arányosságot kifejező integrálási tényező.

Az átmeneti függvény ebben az esetben egy sebességugrás (egyenletesen növekvő vagy csökkenő) jelnek felel meg, amelynek meredekségét a Ki tényező értéke adja. A gyakorlatban a Ki tényező fordított arányát, a Ti integrálási időállandó jelölést is szoktuk használni (lásd fenti képlet).

Ha az integráló jellegű rendszer egy tárolós tulajdonsággal (tehetetlenséggel) is bír, akkor az átmeneti függvény exponenciális jelleggel tart a sebességugrás jelhez, mint aszimptótához. Ha az integráló jellegű rendszer két tárolós tulajdonsággal is rendelkezik, akkor az átmeneti függvény szintén a sebességugrás jelhez tart, de a tranziens átmenetet a már bevezetett ζ csillapítási tényező segítségével, a másodfokú arányos tagokhoz hasonlóan jellemezhetjük.

Összetett rendszerek matematikai modelljét általában a fent említett alapvető rendszermodell típusok soros, párhuzamos és visszacsatolásos kapcsolásával nyerjük. Az összetett rendszer átviteli függvényét könnyen, a kapcsolási tipusokra jellemző képletekkel lehet meghatározni. Az így kapott rendszermodell átviteli függvénye vagy átmeneti függvénye általában magán hordozza az összetevő részrendszerek átviteli vagy átmeneti függvényeinek jellegét.

Az időtartományi analízis eredménye rendszerint egy válaszfüggvény grafikus ábrázolása. Ezek közül is kiemelkedik az átmeneti függvény képe, amelynek segítségével a lineáris (nyított vagy zárt) valós rendszerekre a következő, különösen fontos minőségi jellemzők határozhatók meg:

1.    A relatív maximális túllendülés (Δv): egyenlő az átmeneti függvény maximális értéke és a végértéke közötti különbségének a végértékhez viszonyított arányával. Gyakran százalékban is kifejezik. Csak olyan rendszerekre definiálható, amelyeknél a kimenő jel rendelkezik végértékkel. Így, integráló jellegű rendszerek esetén a maximális túllendülés nem határozható meg.

Összetett rendszerek esetén a maximális túllendülés kézi számítása általában nagyon nehéz feladat, főleg ha nem állnak rendelkezésünkre számítógépes szimulációs módszerek. A gyakorlatban sokszor megelégszünk egy kevésbé pontos értékkel, amelyet a domináns pólus vagy a domináns komplex póluspár figyelembevételével, egyszerűbb képletek alapján kapunk meg.

A domináns komplex póluspár módszere alapján az összetett rendszert egy másodfokú rendszer leírásával közelítjük meg, úgy, hogy a magasabb fokú rendszernek csak azt a komplex póluspárját vesszük figyelembe, amelyik a legközelebb esik a képzetes tengelyhez. Így a magasabbfokú rendszerek relatív maximális túllendülését a csillapított lengésekkel rendelkező másodfokú arányos tagokra jellemző képlettel, az alábbi összefüggés szerint határozunk meg:

                                                                                (II.7)

ahol:   -  az átmeneti függvény maximális értéke,

-  az átmeneti függvény állandósult (stacionárius) végértéke.

A fenti képlet alapján különböző csillapítási tényező értékekre a relatív maximális túllendülések értékeit könnyen ki lehet számítani. Néhányat ezek közül az alábbi táblázat tartalmaz.

ζ

Δv

ζ

Δv

ζ

Δv

0.1

0.7292

0.4

0.2538

0.7

0.0460

0.2

0.5266

0.5

0.1630

0.8

0.0152

0.3

0.3723

0.6

0.0948

0.9

0.0015

2.     A felfutási (növekedési, emelkedési) idő (): jelenti azt az időt, amely alatt az átmeneti függvény először közelíti meg állandósult értékét. A megközelítés mértéke általában a végérték 90÷95 százaléka tehát , ahol .

3.    A beállási (szabályozási) idő ():jelenti azt az időt, amelynek eltelte után az átmeneti függvény a végértékéhez közelít és attól egy általunk előírt értéknél kisebbel már nem tér el, vagyis , ahol általában  vagy akár ennél kisebb érték. A beállási időt a gyakorlatban általában a tranziens folyamat befejezésének tekintjük.

A MATLAB egyik fontos kiegészítője, a CONTROL SYSTEM Toolbox, amely lényegében tartalmazza mindazokat az M-file-okat, mint beépített függvényeket, amelyek az időtartományi szimulációkhoz szükségesek. A következőkben röviden ismertetetjük ezeket.

II.2. Lineáris rendszermodellek időtartományi szimulációja

II.2.1. A lineáris rendszermodell paramétereinek megadása

Mint láttuk, a folytonos lineáris egyváltozós (egy bemenettel és egy kimenettel rendelkező) rendszereket a legegyszerűbb az átviteli függvény segítségével leírni. A MATLAB a W(s) átviteli függvényt két sorvektor alakjában értelmezi, amelyeknek rendezői a B(s) és az A(s) polinomoknak az együtthatói, az s hatványainak csökkenő sorrendjében. Ennek megfelelően, munkánk megkezdésekor definiálnunk kell a "" és az "" sorvektorokat.

II.2.2. A szimuláció időértékeinek kiválasztása

Az összes időtartománybeli szimulációt szolgáló MATLAB funkciók számára a rendszer paraméterein kívül definiálni kell egy idővektort is. Legyen ez egy "t" sorvektor, amelynek rendezői általában egyenlő távolságra lévő időpontok. Ezek a pontok lesznek azok az időértékek, amelyekre az időbeni válasz kiszámításra kerül. Általában a szimuláció kezdeti időpontját nullának vesszük.

Az idővektort Shanon tételének a figyelembevételével kell meghatározni. Gyakorlatilag elfogadható eredményeket kapunk, ha a rendezők értékei közötti "dt" (mintavételezési periódusidő) időlépés értékeként a rendszer legkisebb időállandójánál egy nagyságrenddel kisebb értéket veszünk.

Az időintervallum "tend" nagysága, vagyis az időtartományi szimuláció kiterjedése elegendő kell legyen ahhoz, hogy a kimenő jel elérje a beállási idő által meghatározott állandósult állapotot. Gyakorlatilag elfogadható eredményt kapunk, ha a modell legnagyobb időállandójának 3÷5-szörösét vesszük. Ezek után a "t" idővektort könnyen generálni tudjuk a ”t = [0:dt:tend]” parancssorral.

II.2.3. A "step" függvény

A CONTROL SYSTEM Toolbox beépített "step" függvénye a rendszermodell átmeneti függvényének értékeit (a rendszer egységugrás bemenő jelre adott válaszát) számítja ki. Szintaxisa:

v=step(m,n,t)

amelyben "m" és "n" rendre az átviteli függvény számláló és nevező polinomjainak sorvektorai, "t" az idő sorvektora, "v" pedig a válasz (átmeneti függvény) oszlopvektora.

Az átmeneti függvény állandósult értéke az átviteli függvénybe történő, s=0 behelyettesítéssel érhető el. Ezt az értéket jól megközelítheti a helyes szimulációnk utolsó időpontjában meghatározott függvényérték, amelyet az alábbi utasítások egyikével kaphatunk meg:

vst=v(size(v,1))

vst=v(length(v))

A "size" függvény egy mátrix dimenzióját adja, a "length" függvény egy vektor "hosszát", a "vst" (egy elemből álló) mátrix pedig az átmeneti függvény (megközelítő) végértéket tartalmazza.

A maximális relatív túllendülés meghatározásához a "max" függvényt kell felhasználni:

dv=(max(v)-vst)/vst

ahol a "max" függvény egy mátrix elemei közül a legnagyobbat határozza meg, míg a "dv" változó a maximális relatív túllendülés értékét tartalmazza.

A felfutási idő meghatározásához megkeressük azt az időpontot, amikor az átmeneti függvény eléri a , értéket. Ezt például a következő utasításokkal végezhetjük el:

k=1;beta=0.9;

while v(k)<beta*vst,k=k+1;end

tf=t(k)

A "tf" (egy elemből álló) mátrix a felfutási idő értékét fogja tartalmazni.

A beállási (szabályozási) idő meghatározásához megkeressük azt az időpontot, amelyiktől a "v" oszlopvektor értékei a  végérték és a paraméter által meghatározott  határok közé esnek, a következő lehetséges módon:

delta=0.05;kmax=length(t);k=kmax;i=kmax;

while (1-delta)*vst<v(k),k=k-1;end

[k,v(k),t(k)]

while v(i)<(1+delta)*vst,i=i-1;end

[i,v(i),t(i)]

ts=max(t(k),t(i))

A "ts" mátrix a beállási (szabályozási) idő értékét fogja tartalmazni.

II.2.4. Az "impulse" függvény

A CONTROL SYSTEM Toolbox beépített "impulse" függvénye a rendszermodell súlyfügg-vényének értékeit (a rendszer egységimpulzus bemenő jelre adott válaszát) számítja ki. Szintaxisa:

w=impulse(m,n,t)

amelyben a "w" a válaszfüggvény (súlyfüggvény) oszlopvektora.

II.2.5. Az "lsim" függvény

A CONTROL SYSTEM Toolbox beépített "lsim" függvénye meghatározza egy lineáris rendszermodell válaszát tetszőleges bemenő jelre. Szintaxisa:

y=lsim(m,n,u,t)

amelyben "u" a bemenő jel (gerjesztő függvény) sorvektora, míg "y" a válaszfüggvény (kimenő jel) oszlopvektora. Az "u" sorvektort a parancsfüggvény meghívása előtt meg kell határozni, úgy, hogy mérete (lehetőleg) megegyezzen a "t" idővektor dimenziójával.

Ha a "step", az "impulse" vagy az "lsim" parancsfüggvényeket az egyenlőség baloldala nélkül írjuk fel és úgy használjuk, akkor a számítás után a képernyőre a megfelelő válaszgörbék automatikus kirajzolása is megtörténik.

II.2.6. A "damp" függvény

Az "damp" függvény kiszámítja egy adott lineáris rendszermodell  (saját, természetes) lengési körfrekvenciáit és ζ csillapítási tényezőit. Hívási szintaxisa:

[omega0,zeta]=damp(a)

ahol "omega0" a lengési körfrekvenciák, illetve "zeta" a csillapítási tényezők oszlopvektora, "a" pedig az átviteli függvény (már említett) nevező polinomjának együtthatóit tartalmazó sorvektor.

II.2.7. A "pade" függvény

Ha a folytonos lineáris rendszer holtidővel is rendelkezik, akkor az átviteli függvény egy exponenciális tagot is tartalmaz, amelynek kitevőjében a Th holtidő értéke megjelenik

                                                                      (II.8)

A "step", az "impulse" illetve az "lsim" parancsfüggvények a holtidős rendszereket nem tudják közvetlenül kezelni, mivel polinomiális törtfüggvényekre vannak definiálva. Ezért is, meg a könnyebb kezelhetőség miatt is, a gyakorlatban a holtidős rendszerek átviteli függvényében a holtidőt tartalmazó exponenciális tényezőt különböző polinomiális megközelítéssel szoktuk helyettesíteni. A leggyakoribb helyettesítések a következők:

-  a hatványfüggvény (Strejc) approximáció, amelynek képletei a Taylor sorba való fejlesztésből következnek az alábbi összefüggés szerint. Pontosságuk az n értékével nő:

  , ahol  n=1, 2, 3, ...                                                                           (II.9)

-  az n-ed fokú Padé approximáció, amelynek képletei polinomiális törtfüggvénnyel közelítik meg az exponenciális kifejezést. Ezt a megközelítést építették be a CONTROL SYSTEM Toolbox pade.m állományába.

A "pade" parancsfüggvény a lineáris rendszermodell holtidős részének az n-ed rendű Padé approximációs modelljét állítja elő. Szintaxisa:

[num,den]=pade(th,n)

amelyben "th" a holtidő értéke (egy elemből álló mátrix), "n" a megközelítés foka (egy elemből álló mátrix), "num" és "den" pedig a megközelítő törtfüggvény számláló illetve nevező polinomjainak sorvektorai. Természetesen, az eredő megközelítő átviteli függvény kiszámításánál a holtidőtől független részt (az átviteli függvény holtidő nélküli törtfüggvény részét) is figyelembe kell venni.

II.3. Jellegzetes lineáris rendszerek időtartományi analízise

II.3.1. Elsőfokú időkésleltetéses rendszer

Legyen feladatunk az, hogy határozzuk meg és ábrázoljuk grafikusan a:

                                                                                                           (II.10)

differenciál egyenlettel megadott elsőfokú időkésleltetéses (egytárolós arányos) rendszer átmeneti függvényét és súlyfüggvényét valamint válaszát az  bemenő jelre. Határozzuk meg a felfutási idő értékét, amely alatt az átmeneti függvény eléri végértékének 95%-át. A megoldás lépéseit könnyen követni tudjuk az alábbiakban.

II.3.1.1. A rendszer leírása

A differenciál egyenletből könnyen meghatározható a rendszermodell átviteli függvénye:

                                                                                                                 (II.11)

Az átviteli függvényhez tartozó számláló illetve nevező polinomok:

b1=[0 1];a1=[4 1];

A számláló polinomot ugyanolyan dimenziójúnak, egy kiegészítő 0 segítségével adtuk meg.

II.3.1.2. A szimuláció idővektorának meghatározása

A rendszert jellemző időállandó értéke . Válasszuk a tanulmányozandó időintervallumot 25-nek és a lépésközt 0.2-nek. Ennek megfelelően:

t1=[0:0.2:25];

II.3.1.3. Az átmeneti függvény meghatározása

Az átmeneti függvény "v1" oszlopvektorát az alábbi utasításokkal határozzuk meg és ábrázoljuk grafikusan (lásd alábbi ábra):

v1=step(b1,a1,t1);step(b1,a1,t1),grid,title('Átmeneti függvény 1')

Megjegyzés: Ugyanazt a képet kapnánk, de a grafikus ábrázolási lehetőségekben szélesebb palettát alkalmazhatnánk, ha a másodszor meghívott "step" függvény helyett a "plot" függvényt használnánk. A tengelyek nevesítéséről nekünk kell gondoskodnunk mindkét esetben!

 

II.3.1.4. A súlyfüggvény meghatározása

A súlyfüggvény "w1" oszlopvektorát az alábbi utasítások segítségével határozzuk meg és ábrázoljuk grafikusan (lásd alábbi ábra):

w1=impulse(b1,a1,t1);impulse(b1,a1,t1),grid,title('Súlyfüggvény 1')

Megjegyzés: Ugyanazt a képet kapnánk, de a grafikus ábrázolási lehetőségekben szélesebb palettát alkalmazhatnánk, ha a második "impulse" függvény helyett a "plot" függvényt használnánk. A tengelyek nevesítéséről nekünk kell gondoskodnunk mindkét esetben!

  

II.3.1.5. Az átmeneti függvény és a súlyfüggvény összehasonlítása

Hasonlítsuk össze a két grafikus ábrázolást. Ezt a legkönnyebben úgy tehetjük meg, hogy ugyanabban a grafikus ablakban ábrázoljuk mindkét függvényt. A szükséges utasítások:

subplot(2,1,1)

step(b1,a1,t1),grid

subplot(2,1,2)

impulse(b1,a1,t1),grid

  

Vegyük észre a különbséget a kezdeti értékek és végértékek között. Az átmeneti függvény végértékét az alábbi utasítások begépelése után kapjuk meg:

vst1=v1(length(v1));disp(vst1)

A kapott értéket a "disp" függvényparancs a képernyőre íratja

 0.9981

azaz az időintervallum végére az elméleti végértéket 0.2%-nál kisebb hibával közelítjük meg.

II.3.1.6. A szinuszos bemenő jelre adott válasz meghatározása

A szimuláció céljából definiáljuk a színuszos bemenő jel "u1" sorvektorát, majd határozzuk meg a válaszjel "y1" oszlopvektorát és rajzoltassuk azt is ki. A következő parancsok szükségesek:

u1=sin(1.8*t1);y1=lsim(b1,a1,u1,t1);

lsim(b1,a1,u1,t1),grid,title('Szinuszos jelre adott válasz 1')

amelyek az alábbi grafikus ábrához vezetnek:

  

Értékeljük a kapott eredményt. Vegyük észre a tárolós tulajdonság időállandójának hatását a tranziens válaszban! Ezután rajzoljuk ki ugyanabba a grafikus ablakba mindkét jel időfüggvényét:

lsim(b1,a1,u1,t1),grid,title('Szinuszos jelre adott válasz 1')

hold on,plot(t1,u1),hold off  

  

Hasonlítsuk össze a két jelet! Nyomban következtetni lehet a tanulmányozott tag csillapító (szűrő) és fáziskésleltető hatására. Ennek további tanulmányozására végezzünk szimulációkat a színuszos bemenő jel más frekvencia értékeire is! Milyen következtetéseket lehet levonni?

II.3.1.6. A felfutási idő illetve a beállási idő meghatározása

Hasonlítsuk össze úgy az átmeneti függvény mint a súlyfüggvény esetében a felfutási időket és a beállási időket. A grafikus válaszokból kitűnik, hogy mindkét függvény esetén ezek az idők megegyeznek ugyanazokra a β illetve δ értékekre. Ha β = 0.95 illetve δ = 0.05 gyakorlati értékeket vesszük, akkor a beállási idők megegyeznek a felfutási időkkel. Ezt a közös értéket a már felvázolt utasítások illetve az alábbiak szerint határozzuk meg:

k=1;beta=0.95;

while v1(k)<beta*vst1, k=k+1; end

tf1=t1(k);disp(tf1)

A kapott érték

   12 

azaz 3-szor a  tárolós időállandó értéke. Jegyezzük meg jól ezt az összefüggést!

II.3.2. Másodfokú időkésleltetéses rendszer

Legyen feladatunk az, hogy határozzuk meg és ábrázoljuk grafikusan a:

                                                                               (II.12)

differenciál egyenlettel megadott másodfokú, csillapított lengésekkel rendelkező időkésleltetéses rendszer átmeneti függvényét és súlyfüggvényét. Határozzuk meg a rendszer pólusait, csillapítási tényezőjét, lengési körfrekvenciáját, valamint a maximális túllendülést, a felfutási időt és a beállási időt, amely alatt az átmeneti függvény eléri végértékének 95%-át. Tanulmányozzuk az adott másodfokú rendszer viselkedését más csillapítási tényező értékekre is.

II.3.2.1. A rendszer leírása

A differenciál egyenletből könnyen meghatározható a rendszermodell átviteli függvénye:

                                                                                                      (II.13)

Az átviteli függvényhez tartozó számláló illetve nevező polinomok:

b2=[0 0 4];a2=[1 0.8 4];

II.3.2.2. A szimuláció idővektorának meghatározása

A rendszer átviteli függvényét felírhatjuk még az alábbi alakba is:

 ,                                                            (II.14)

amelyből könnyen meghatározható az adott rendszer  saját (természetes) körfrekvenciája és a  csillapítási tényező érték. Válasszuk ezeknek az értékeknek megfelelően a tanulmányozandó időintervallumot 10 egységnek és a lépésközt 0.02 egységnek, vagyis:

t2=[0:0.01:10];

II.3.2.3. Az átmeneti függvény és a súlyfüggvény meghatározása

Az átmeneti függvény "v2" és a súlyfüggvény "w2" oszlopvektorait a már ismert "step" illetve "impulse" parancsfüggvényekkel határozzuk meg. A jobb áttekinthetőség végett, ábrázoljuk ezúttal a két függvényt ugyanabban a grafikus ablakban (lásd alábbi ábra). A szükséges utasításokból álló program ezúttal:

v2=step(b2,a2,t2);step(b2,a2,t2),title('Átmeneti függvény és súlyfüggvény 2')

hold on

w2=impulse(b2,a2,t2);impulse(b2,a2,t2),grid

hold off 

  

Az átmeneti függvény végértékét az alábbi utasítások begépelése után kapjuk meg:

vst2=v2(length(v2));disp(vst2)

A kapott érték

0.9840

azaz az időintervallum végére az elméleti végértéket most is 0.2%-nál kisebb hibával közelítjük meg.

II.3.2.4. A rendszerpólusok meghatározása

Egy lineáris rendszert jellemző átviteli függvény nevező (karakterisztikus) polinomjának gyökei a pólusok. Ezek hasznos információkat szolgáltatnak az illető rendszer dinamikus viselkedéséről, belső struktúrájáról (lásd szakirodalom!). Adott rendszerünk pólusait könnyen meghatározhatjuk a "roots" beépített függvényparancs segítségével:

p2=roots(a2);disp(p2)

A kapott értékek

  -0.4000 + 1.9596i

  -0.4000 - 1.9596i

A pólusok komplex (konjugált) alakja igazolja az átmeneti illetve a súlyfüggvény lengő voltát.

II.3.2.5. A csillapítási tényező és a lengési körfrekvencia meghatározása

A csillapítási tényező és a lengési (saját) körfrekvencia értékeit a "damp" függvény adja:

[om2,ze2]=damp(a2)

A képernyőre kapott válasz

om2 =

     2

     2

ze2 =

    0.2000

    0.2000

Ha a "damp" függvényt csak magára hívjuk meg:

damp(a2)

úgy az alábbiak szerinti választ kapjuk:

    Eigenvalue        Damping        Freq. (rad/sec)

  -0.4000 + 1.9596i   0.2000             2.0000         

  -0.4000 - 1.9596i   0.2000             2.0000            

Figyeljük meg a kétféle utasítás közötti különbséget és a kapott eredményeket! A második esetben a pólusokat mint sajátértékeket tünteti fel a számítógép. Ez annak köszönhető, hogy a pólusok megegyeznek az állapottérben való leírásban az állapotmátrix sajátértékeivel (lásd később). Az eredmények természetesen találnak az előbb meghatározott értékekkel!

II.3.2.6. A felfutási idő, a maximális relatív túllendülés és a beállási idő meghatározása

Ha most is a β = 0.95 viszonyítási értéket vesszük, akkor a felfutási idő értéke az alábbiakban felvázolt utasítások szerint kerül kiszámításra:

k=1;beta=0.95;

while v2(k)<beta*vst2, k=k+1; end

tf2=t2(k);disp(tf2)

A kapott eredmény:

    0.8600

amit a grafikus ábrázolásból is meghatározhatunk (pontatlanabbul!)

A maximális relatív túllendülés értékét az alábbiak begépelésével kapjuk:

dv2=(max(v2)-vst2)/vst2;disp(dv2)

Az eredmény:

    0.5514

ami kb. 55%-os relatív túllendülést jelent.

A beállási (szabályozási) idő értékét δ = 0.05 tűréshatárral az alábbi utasításokkal kapjuk meg:

delta=0.05;kmax=length(t2);k=kmax;i=kmax;

while (1-delta)*vst2<v2(k),k=k-1;end

while v2(i)<(1+delta)*vst2,i=i-1;end

ts2=max(t2(k),t2(i));disp(ts2)

A program által szolgáltatott eredmény:

    8.3200

vagyis kb. 8 időegység elteltével állandósult állapotot tekinthetünk.

II.3.2.7. Az átmeneti függvény és a súlyfüggvény más csillapítási tényező értékekre

Az adott másodfokú rendszer viselkedését az átmeneti függvény és a súlyfüggvény grafikus ábrázolásának segítségével más csillapítási tényező értékekre is tanulmányozhatjuk. Észrevesszük, hogy a csillapítási tényező értéke az átviteli függvény nevezőjének második tagjában szerepel.

Legyenek a csillapítási tényező értékei rendre: 0, 0.2, 0.4, 0.6, ..., 1.8, 2.0. Mivel a lépésköz állandó, a nevező vektorának megadását egy "for" ciklusutasítással végezhetjük el. A ciklusban szerepelnie kell rendre az átmeneti illetve a súlyfüggvény kiszámításának és grafikus ábrázolásának. A viselkedés jobb megértése végett ábrázoljuk grafikusan egy-egy diagrammon a kapott átmeneti- és súlyfüggvény seregeket. Az átmeneti függvények esetén a lehetséges utasítások:

step(b2,[1 0 4],t2),grid, title('Átmeneti függv. 2 - változó csillapítás')

hold on,for k=1:10,step(b2,[1 4*k*0.2 4],t2),end,hold off 

  

Végezzük el ugyanazt a súlyfüggvény esetében is, azaz:

impulse(b2,[1 0 4],t2),grid, title('Súlyfüggvény 2 - változó csillapítás')

hold on,for k=1:10,impulse(b2,[1 4*k*0.2 4],t2),end,hold off 

  

Az ábrákból kitűnik, hogy minél kisebb a csillapítási tényező, annál nagyobb amplitudójú és lassúbb lecsengésű válaszokat kapunk, illetve ha a csillapítási tényező értéke meghaladja az egységet, a válaszok már nem tartalmaznak lengő összetevőket. Ha a csillapítási tényező értéke nulla, a lengések állandósultak a saját frekvenciaértéknek megfelelő periódusidővel.

II.3.3. Összetett folytonos lineáris rendszer

Feltételezzük, hogy adott az alábbi átviteli függvénnyel rendelkező rendszermodell:

                                                                             (II.15)

Határozzuk meg a rendszer átmeneti függvényét; pólusait, a maximális túllendülést, a felfutási (emelkedési) időt és a beállási (szabályozási) időt, valamint határozzuk meg a rendszer válaszát egy  egység periódusidejű u(t) fűrészfog alakú bemenő jelre, amelynek maximális u(5)=4 értéke T1=5 időegységnél van. Természetesen, ahol szükséges, ábrázoljuk grafikusan is az eredményeket.

II.3.3.1. A rendszer leírása

MATLAB parancsok segítségével meghatározzuk az átviteli függvényhez tartozó számláló illetve nevező polinomokat az alábbiak szerint:

b3=12.5*[0,0,0,1,0.8];

a3=conv([1,1],[1,4]);a3=conv(a3,[1,3,6.25]);

Vegyük észre a számláló és nevező polinomok megadási módját. Elvileg, a számláló polinom sorvektora ugyanolyan dimenziójú lehet, mint a nevezőé. A szükséges kiegészítést 0-k végzik.

II.3.3.2. A szimuláció idővektorának meghatározása

A rendszer átviteli függvényéből meghatározhatjuk a legnagyobb és legkisebb időállandó illetve a saját periódusidő értékét (végezzük el!) Ezek alapján írhatjuk fel az alábbi idővektort:

t3=[0:0.01:6];

II.3.3.3. Az átmeneti függvény és a súlyfüggvény meghatározása

A "v3" átmeneti függvény és a "w3" súlyfüggvény oszlopvektorait az előbbiekhez hasonlóan határozzuk meg. Ábrázoljuk most is a két függvényt ugyanabban a grafikus ablakban (lásd program):

v3=step(b3,a3,t3);step(b3,a3,t3),title('Átmeneti függvény és súlyfüggvény 3')

hold on

w3=impulse(b3,a3,t3);impulse(b3,a3,t3),grid

hold off

  

Az átmeneti függvény végértékét az alábbi utasítások begépelése után kapjuk meg:

vst2=v2(length(v2));disp(vst2)

A kapott érték

0.4005

azaz az általunk választott időintervallum végére az elméleti végértéket (számítsuk ki!) most is igen kis hibával közelítjük meg. Mekkora ez a hiba %-ban?

II.3.3.4. A pólusok, a csillapítási tényezők és a lengési körfrekvenciák meghatározása

A legegyszerűbb, ha a "damp" függvényt magára meghívjuk, azaz:

damp(a3)

A kapott válasz:

Eigenvalue            Damping             Freq. (rad/sec)

  -1.0000             1.0000             1.0000         

  -1.5000 + 2.0000i   0.6000             2.5000         

  -1.5000 - 2.0000i   0.6000             2.5000         

  -4.0000             1.0000             4.0000         

amelyből nyomban kitűnik, hogy két valós és két komplex konjugált pólus jellemzi modellünket. A komplex póluspárhoz tartozó csillapítási tényező  illetve a saját frekvencia .

II.3.3.5. A felfutási idő, a maximális relatív túllendülés és a beállási idő meghatározása

Ha most a β = 0.9 viszonyítási értéket vesszük, akkor a felfutási idő értéke az alábbiakban felvázolt utasítások szerint kerül kiszámításra:

k=1;beta=0.9;

while v3(k)<beta*vst3, k=k+1; end

tf3=t3(k);disp(tf3)

A kapott eredmény:

    1.0600

amit a grafikus ábrázolásból is meghatározhatunk (pontatlanabbul!)

A maximális relatív túllendülés értékét az alábbiak begépelésével kapjuk:

dv3=(max(v3)-vst3)/vst3;disp(vh3)

A kapott eredmény:

    0.1769

ami közel 18%-os relatív túllendülést jelent.

A beállási (szabályozási) idő értékét δ = 0.02 tűréshatárral az alábbi utasításokkal számítjuk ki:

delta=0.02;kmax=length(t3);k=kmax;i=kmax;

while (1-delta)*vst3<v3(k),k=k-1;end

while v3(i)<(1+delta)*vst3,i=i-1;end

ts3=max(t3(k),t3(i));disp(ts3)

A program által szolgáltatott eredmény:

    3.0200

vagyis kb. 3 időegység elteltével állandósult állapotot tekinthetünk.

II.3.3.6. A fűrészfog bemenő jelre adott válasz meghatározása

A fűrészfog bemenő jel előállításához előbb egy megfelelő új idővektort generálunk:

ujt3=0:0.05:10;

majd a bemenő jel vektorát a következőképpen építjük fel:

u3=(4/5)*ujt3.*[ones(1,101),zeros(1,100)];

Ezek után a válaszjelet az "lsim" parancsfüggvény segítségével határozzuk meg:

y3=lsim(m3,n3,u3,ujt3);

majd a jobb megértés végett ábrázoljuk egy grafikus ablakban a bemenő jelet és a rendszer válaszát. A szükséges programot gyakorlatként kérjük megírni. Az alábbi ábra munkánk ellenőrzésére szolgál.

  

Figyeljük meg az u(t) diszkontinuitásának a hatását az y(t) kimenetre. Adjunk magyarázatot a rendszermodell válaszjelének időbeni változására!

II.3.4. Holtidős rendszer

Tételezzük fel, hogy a (II.16) rendszer holtidővel is rendelkezik, azaz átviteli függvénye:

                                                                     (II.16)

Legyen feladatunk az, hogy határozzuk meg és ábrázoljuk grafikusan az átmeneti függvényt, alkalmazva a holtidős részre először a hatodfokú hatványfüggvény megközelítést, majd másodszor az első- illetve másodfokú Padč megközelítést.

II.3.4.1. A szimuláció idővektorának meghatározása

A holtidős rendszermodell szimulációjához gyakorlatilag az előbbi szimuláció időtartományát meg kell hosszabbítanunk a holtidő nagyságával, azaz 2 egységgel. Az új idővektor:

t4=[0:0.01:8];

II.3.4.2. A rendszer leírása hatványfüggvény megközelítéssel

Az exponenciális (holtidős) tényezőt az alábbi hatodfokú hatványfüggvénnyel közelítjük:

, ahol Th = 2                                                           (II.17)

A rendszermodell eredő átviteli függvényének nevező polinomja tizedfokú lesz. A számláló és a nevező megfelelő sorvektorait a következő utasításokkal határozhatjuk meg:

b4=[0 0 0 0 0 0 0 12.5 10];

a41=[1/3 1];a42=conv(a41,a41);a44=conv(a42,a42);a46=conv(a44,a42);

a40=conv([1,1],a46);a4=conv([1,4],a40);

Figyeljük meg a MATLAB programnyelv adta lehetőségek egyszerűségét, azok kihasználását.

II.3.4.3. Az átmeneti függvény meghatározása hatványfüggvény megközelítéssel

Az átmeneti függvényt az előbbi esetekhez hasonlóan határozzuk meg és ábrázoljuk:

v4=step(b4,a4,t4);step(b4,a4,t4)

title('Átmeneti függvény 4 - hatványfüggv.')

  

II.3.4.4. A rendszer leírása elsőfokú Padé megközelítéssel

Az elsőfokú Padé megközelítéssel az exponenciális tényezőt az alábbi képlettel helyettesítjük:

,    vagyis                                                                          (II.18)

Így a rendszer eredő átviteli függvényének nevező polinomja ötödfokú lesz. A számláló és nevező polinomoknak megfelelő sorvektorokat a következő utasításokkal határozhatjuk meg:

b4p1=conv(b3,[-1 1]);a4p1=conv(a3,[1 1]);

  

II.3.4.5. Az átmeneti függvény meghatározása elsőfokú Padé megközelítéssel

Az átmeneti függvényt hasonlóan határozzuk meg és ábrázoljuk (lásd előző ábra):

h4p1=step(b4p1,a4p1,t4);step(b4p1,a4p1,t4)

title('Átmeneti függvény 4 - Pade 1 ')

Próbáljunk magyarázatot találni a fenti válaszgörbe alakjára!

II.3.4.6. A rendszer leírása másodfokú Padé megközelítéssel

A másodfokú Padé megközelítéssel az exponenciális tényezőt az alábbi képlettel helyettesítjük:

,    vagyis                                               (II.19)

Így a rendszer eredő átviteli függvényének nevező polinomja hatodfokú lesz. A számláló és nevező polinomoknak megfelelő sorvektorokat a következő utasításokkal határozhatjuk meg:

b4p2=conv(b3,[0.333 -1 1]);a4p2=conv(a3,[0.333 1 1]);

II.3.4.7. Az átmeneti függvény meghatározása másodfokú Padé megközelítéssel

Az átmeneti függvényt az ismert módon határozzuk meg és ábrázoljuk grafikusan:

h4p2=step(b4p2,a4p2,t4);step(b4p2,a4p2,t4)

title('Átmeneti függvény 4 - Pade 2 ') 

  

Próbáljunk magyarázatot találni a fenti válaszgörbe alakjára is! Milyen következtetésekre jutunk?

Hasonlítsuk össze a holtidős rendszer különböző megközelítésekkel meghatározott átviteli függvényeit az elméleti válasszal, kirajzolva azokat egyetlen egy grafikus ablakba. A szükséges MATLAB programot gyakorlatként kérjük megírni, kihasználva a színes képernyő által adott színes grafikus lehetőségeket.

További gyakorlatként oldjuk meg a jelen fejezethez kapcsolódó feladatokat (II. Gyakorlat).