+ All Categories
Home > Documents > Indrumar de Laborator Matlab

Indrumar de Laborator Matlab

Date post: 25-Jul-2015
Category:
Upload: roxana-munteanu
View: 234 times
Download: 17 times
Share this document with a friend
74
Corneliu Octavian TURCU Cristina Elena TURCU Facultatea de Inginerie Electrică Universitatea “Ştefan cel Mare” SUCEAVA pentru uzul studenţilor Suceava, 1996 Matlab pentru introducerea în automatică
Transcript
Page 1: Indrumar de Laborator Matlab

Corneliu Octavian TURCU Cristina Elena TURCU

Facultatea de Inginerie Electrică

Universitatea “Ştefan cel Mare” SUCEAVA

pentru uzul studenţilor

Suceava, 1996

Matlabpentru introducerea în automatică

Page 2: Indrumar de Laborator Matlab

Capitolul I

1. Elemente de bază ale mediului MATLAB 1.1 Introducere

Matlab este un mediu interativ, de înaltă performanţă, utilizat pentru calcule ştiinţifice şi inginereşti.

MATLAB înglobează diverse posibilitaţi de lucru, cum ar fi analiza numerică, calculul matricial, procesare de semnale şi trasarea grafică, într-un mediu uşor de utilizat, enunţurile şi soluţiile fiind exprimate exact cum sunt scrise matematic, nefiind necesară o programare tradiţională.

Numele acestui mediu provine de la matrix laboratory. A fost dezvoltat pornind de la doua pachete de programe, LINPACK şi EISPACK, reprezentând bazele soft-ului pentru calculul matricial.

Elementele de baza ale MATLAB-ului sunt matriciale, nefiind necesară o dimensionare a acestora. Acest mediu de programare este util în rezolvarea problemelor numerice.

1.2 Utilizarea mediului MATLAB

1.2.1 Lansarea mediului

Pentru lansare se deschide, în mediul de programare Windows, grupul Matlab for Windows şi se da un click pe icon-ul Matlab. La iniîializare, se afişează prompterul mediului, constând în simbolul >>. La prompter se pot scrie liniile de comandă.

1.2.2 Liniile de comandă

Pentru scrierea liniilor de comandă se utilizează tastatura calculatorului. Tastele săgeţi pot fi utilizate pentru editarea liniilor sau pentru rechemarea unei comenzi anterioare (tasta săgeată sus).

O linie de comandă arată în modul următor: log(sqrt(atan2(3,4)))

Tastele săgeţi acţionează asupra unor copii al liniilor de comandă anterioare, aceste copii fiind salvate într-un buffer. Funcţiile tastelor sunt prezentate în continuare:

Săgeată sus Recheamă comanda dată anterior Săgeată jos Recheamă linia de comandă următoare Săgeată stânga Deplasare la stânga cu un caracter Săgeată dreapta Deplasare la dreapta cu un caracter Ctrl+săgeată stânga Deplasare la stânga cu un cuvânt Ctrl+săgeată dreapta Deplasare la dreapta cu un cuvânt Home Deplasare la începutul liniei de comandă End Deplasare la sfârşitul liniei de comandă Esc Abandonează linia de comandă curentă Ins Trece între modurile Insert şi Overtype Backspace Şterge un caracter la stânga cursorului

1.3 Fundamente MATLAB MATLAB lucrează cu un singur tip de elemente, şi anume matrici rectangulare, elementele acestora

putând fi numere reale sau chiar complexe. În unele situaţii, acest mediu lucrează şi cu scalari (reprezentaţi ca matrici cu o linie şi o coloană) sau cu matrici care au o singură linie sau coloană, acestea din urmă fiind vectori.

Page 3: Indrumar de Laborator Matlab

1.3.1 Introducerea matriciilor

Matriciile pot fi introduse în MATLAB în diferite moduri: • explicit; • ca rezultat al unor calcule; • create în fişiere de tip M; • incărcate din fişiere externe. Cea mai simlpă metodă de a introduce o matrice este scrierea explicită. De exemplu, introduând

linia de comandă: A=[1 2 3 ;4 5 6;7 8 9] Se obţine rezultatul: A=

1 2 3 4 5 6 7 8 9

Matricea A este salvată în vederea unor prelucrări ulterioare. Matriciile de dimensiuni mari pot fi introduse pe mai multe linii, utilizand RETURN în locul “,”. Pentru exemplificare, o matrice poate fi introdusă astfel:

A= [ 1 2 3

4 5 6 7 8 9 ]

1.3.2. Elementele Matricilor

Elementele matriciilo pot fi orice expresii MATLAB. De exemplu comanda x=[-1.3 sqrt(3) (1+2+3)*4/5] conduce la rezultatul: x= -1.3000 1.7321 4.8000

Elementele matricilor pot fi referite şi individual, cu indicele cuprins între paranteze rotunde. De exemplu,

x(5)=abs(x(1)) va produce x= -1.3000 1.7321 4.8000 0 1.3000

Matricile de dimensiuni mari pot fi construite utilizând ca elemente matrici de dimensiuni mai mici. De exemplu, se poate ataşa un nou rând matricii A prin comenzile:

r=[10 11 12] A=[A;r] rezultatul fiind: A=

1 2 3 4 5 6 7 8 9

10 11 12

Page 4: Indrumar de Laborator Matlab

O matrice poate fi extrasă dintr-o matrice de dimensiune mai mare. De exemplu, comanda: A=A(1:3,:) va avea ca efect extragerea din matricea A a unei matrici formată din primele trei rânduri şi toate coloanele matricii curente.

1.3.3. Variabile şi declaraţii

Expresiile introduse de utilizator sunt interpretate şi evaluate. În MATLAB declaraţiile sunt, de cele mai multe ori de tipul:

varabilă = expresie.

sau simplu:

expresie.

Evaluarea unei expresii produce o matrice afişată pe ecran şi alocată unei variabile, în vederea

utilizării ulterioare. Dacă numele variabilei şi semnul “=” sunt omise, este creată automat o variabilă cu numele ans (de la answer).

Dacă o linie de comandă se termină cu caracterul “,” afişarea rezultatului este inhibată. Dacă linia de comandă este mai lungă decât o linie de ecran, se poate utiliza secvenţa “…” pentru indicarea continuării scrierii liniei de comandă pe următoarea linie de ecran.

Numele variabilelor sunt formate din o literă, urmată de litere şi cifre. Sunt memorate doar primele 19 caractere ale numelor variabilelor.

MATLAB este sensibil la tipul de caractere utilizate, mici sau mari. Toate numele funcţiilor se scriu cu litere mici.

1.3.4. Informaţii asupra spaţiilor de lucru

Variabilele create anterior sunt memorate în spaţiul de lucru al MATLAB-ului.

Execuţia comenzii:

who

va avea ca rezultat afişarea tuturor variabilelor din acest spaţiu.

1.3.5. Expresii aritmetice

Pentru numere este utilizată reprezentarea cu punct zecimal, cu semnul “_” în cazul numerelor negative, sau în modul bază exponent. Exemplu de numere declarate corect sunt: 3 3.83745683 -77 1.896E-10 0.002 7.231E21 Intervalul în care pot fi introduse numere este de 10-308 la 10308. MATLAB are implementate multe funcţii in fişiere de tip m. De exemplu,

abs,sqrt, log, sin. De asemenea există funcţii ce returnează valori speciale, cum ar fi pi, ce returnează valoarea lui π, sau Inf, ce returnează valoarea ∞.

Page 5: Indrumar de Laborator Matlab

1.3.5. Numere complexe Numere complexe sunt introduse in MATLAB utilzându-se funcţiile speciale I sau j. Exemple:

z=3+4*I; z=3+4j; w=r*exp(I*theta).

Pentru introducerea matricilor complexe se poate recurge la:

A=[1 2 ; 3 4]+I*[5 6 ; 7 8] sau

A=[1+5*i 2+6*i ; 3+7*I 4+8*i]

comenzi care vor conduce la acelaşi rezultat.

1.3.6. Formatul de ieşire Rezultatul declaraţiilor şi evaluărilor executate sunt afişate pe ecran. Pentru a explicita formatele

posibile, se consideră următorul exemplu: dându-se

x=[4/3 1.2345E-6]

atunci formatele ce pot fi utilizate sunt:

short 1.3333 0.0000 short e 1.3333E+000 1.2345E-006 long 1.333333333333338 0.000001234500000 long e 1.333333333333338E+000 0.000001234500000E-0006 hex 3FF5555555555555 3EB4B6231ABFD271 format + + +

1.3.7. Facilităţi oferite de HELP Pentru a obţine informaţii ajutătoare se va tasta:

help Dacă se doreste o informaţie specifică, atunci se va tasta:

help eig Comanda anterioară va vea ca efect afişarea informaţiilor referitoare la valorile proprii ale unei matrici (eig=eigenvalues), iar comanda:

help [ va furniza informaţii referitoare la utilizarea parantezelor pătrate.

1.3.9. Salvarea spaţiului de lucru Pentru părăsirea mediului se va tasta quit sau exit. Terminarea sesiunii de lucru are ca efect

pierderea valorilor variabilelor din spaţiul de lucru. Înainte de terminare, se poate salva conţinutul spaţiului de lucru utilizând comanda:

save

Ca urmare a acestei comenzi, se vor salva toate variabilele din spaţiul de lucru într-un fişier numit

matlab.mat. La următoarea sesiune de lucru, se poate restaura spaţiul de lucru cu comanda:

load

Page 6: Indrumar de Laborator Matlab

Comenzile anterioare pot fi folosite în joncţiune cu nume de fişiere sau numai pentru variabile sectate.

save temp X

save temp x z y

vor avea ca efect salvarea în fişierul temp.mat a variabilelor x si respectiv x,y,z. 1.3.9.1. Funcţii MATLAB pune la dispoziţia utilizatorului o multitudine de funcţii. Unele dintre acestea sunt intrinsece,

iar altele sunt construite. De asemenea, există funcţii specifice domeniilor de lucru ale MATLAB-uli, disponibile în modulele de istalare (MATLAB TOOLBOX). O caracteristică importantă a acestui mediu este că utilizatorul îşi poate construi propriile funcţii, ele putând fi apelate exact ca funcţiile intrisece ale mediului.

1.4 Operaţii cu matrici

1.4.1. Transpusa unei matrici Caracterul special ’ (apostrof) semnifică transpusa unei matrici. Declaraţiile:

A=[1 2 3 ; 4 5 6 ; 7 8 0], B=A’

conduc la: B= 1 4 7 2 5 8 3 6 0 iar:

x=[-1 0 2]’ produce: x= -1 0 2

În principal, utiliându-se apostrof, are loc o transpunere formală a matricii. Dacă matricea este complexă, atunci A’ este conjugata complexă transpusă. Acest lucru poate conduce uneori la rezultate eronate, mai ales în cazul lucrului cu numere complexe. Pentru transpusa neconjugată se va utiliza A.’ sau conj(A’).

1.4.2. Adunarea şi scăderea matricilor Operaţiile de adunare şi scădere a matricilor sunt simboloizate prin + şi -. Operatorii sunt definiţi

numai dacă matricile au aceaşi dimensiune. Operaţia de adunare şi scădere mai este definită şi atunci când unul din operanţi este scalar (matrice cu o linie şi o coloană). De exemplu:

y=x-1 va produce:

y= -2 -1 1

Page 7: Indrumar de Laborator Matlab

1.4.3. Înmulţirea matricilor Operaţia este simbolizată prin *. Operaţia este validă numai în ipoteza concordanţei dintre

dimensiunile matricilor (indicele de interior). De asemenea se poate vorbi de produsul scalar a doi vectori, definit ca:

x’*y sau y’*x rezultatul fiind acelaşi. Sunt posibile şi produsele:

x*y’ sau y*x’ Alte operaţii de înmulţire pot fi executate în cazul operanzilor de tip matrice şi vector, şi scalar şi vector, scalar şi matrice.

1.4.4. Împarţirea matricilor În MATLAB sunt doua tipuri de împărţiri, simbolizate prin \ şi /. Dacă A este o matrice nesingulară,

atunci A\B şi B/A corespund formal înmulţirii la stânga şi la dreapta a matricii B cu inversa lui A, adică inv(A)*B şi B*inv(A), rezultatul fiind obţinut direct, fară calculul inversei matricii A. În general,

X=A\B este soluţia ecuaţiei A*X=B X=A/B este soluţia ecuaţiei X*A=B 1.4.5. Ridicarea la putere Expresia A^p denotă matricea A ridicată la puterea p. Operaţia este validă dacă A este matrice

pătratică si p este scalar.

1.5 Operaţii cu tablouri Utilizarea expresiei operaţii cu tablouri se referă la operaţii aritmetice executate element cu element

(spre deosebire de operaţiile liniare cu matrici, simbolizate prin * / \ ‘ ). Dacă un operator este precedat de un punct . atunci înseamnă că operaţia este executată element cu element.

1.5.1. Adunarea şi scăderea Pentru aceste operaţii, situaţia este identică cu cea din cazul matricilor. Se utilizează aceleaşi

simboluri, + şi - . 1.5.2. Înmulţirea şi împărţirea Înmulţirea element cu element este simbolizată prin .*. Dacă A şi B au aceeaşi dimensiune, atunci

A.*B denotă produsul element cu element al celor două matrici. De exemplu, dacă:

x=[1 2 3];y=[4 5 6]; z=x.*y

z=[4 10 18] Expresiile A./B şi A.\B dau ca rezultat câtul împarţirii element cu element. z=x.\y z=[4.0000 2.5000 2.0000]

Page 8: Indrumar de Laborator Matlab

1.5.3. Ridicarea la putere Ridicarea la putere element cu element este simbolizată prin .^.

Exemle: z=x.^y z=[1 32 729] z=x.^2 z=[1 4 9] z=2.^[x y] z=[2 4 8 16 32 64] Important: ultimul exemplu ilustrează una din subtilităţile MATLAB-ului. Este dificil de văzut că exista un spaţiu între 2 şi punct. Dacă acest spaţiu nu ar fi existat, punctul ar fi fost interpretat ca punct zecimal asociat cu 2. MATLAB ar fi sesizat numai ^ şi ar fi înceracat să calculeze puterea unei matrici, lucru care ar fi condus la o eroare, deoarece exponentul nu este pătratic. O alternativă pentru această operaţie ar fi utilizarea parantezelor.

1.5.4. Operatori relaţionali MATLAB pune la dispoziţia utilizatorului un număr de 6 operatori releţionali: < mai mic <= mai mic sau egal

> mai mare >= mai mare sau egal == egal = non egal Comparaţia este făcută între elemente corespondente. Rezultatul este o matrice formată din elemente 0 şi 1, reprezentând FALSE şi TRUE. De exemplu,

2+2=4 va fi 0. Operatorii relaţionali pot fi utilizaţi pentru localizarea elementelor unei matrici, elemente ce satisfac anumite condiţii. De exemplu, se consideră o matrice magică (matrice având elemente cuprinse între 1 şi n2, având sumele pe linii şi coloane identice).

A=magic(6) A= 35 1 6 26 19 24 3 32 7 21 23 25

31 9 2 22 27 20 8 28 33 17 10 15 30 5 34 12 14 16 4 36 29 13 18 11

Se poate observa că elementele divizibile cu 3 sunt dispuse la fiecare a treia diagonală a matricii. Pentru

a vedea acest lucru, se tastează: P=(rem(A,3)==0)

ceea ce semnifică un test al egalitaţii cu 0 al resturilor împărţirii elementelor matricii A la 3.P va deveni o matrice cu elemente 0 şi 1.

P= 0 0 1 0 0 1 1 0 0 1 0 0

0 1 0 0 1 0 0 0 1 0 0 1 1 0 0 1 0 0 0 1 0 0 1 0

Page 9: Indrumar de Laborator Matlab

Funcţia find este utilă în lucrul cu operatori relaţionali, servind la găsirea elementelor neutre într-o matrice cu elemente 1 şi 0, si chiar a elementelor care satisfac o anumităm condiţie. De exemplu, dacă Y este un vector, find(Y 3.0) va returna un vector conţinând indicii elemnetelor lui Y care sunt mai mici decât 3.

Declaraţiile: i=find(y>3.0) Y(I)=10*ones(i)

vor conduce la înlocuirea cu 10.0 a elementelor mai mari decât 3. Aceste operaţii sunt valide şi la lucrul cu matrici, deoarece o matrice poate fi referită ca un vector coloană cu un singur indice.

1.5.5. Operaţii logice Operatorii logici sunt: & AND | OR ~ NOT Primul dintre ei compară doi scalari sau două matrici de dimensiuni egale. Pentru matrici, aceşti operatori lucrează asupra elementelor corespondente. Dacă A şi B sunt matrici formate din elemente 1 şi 0, atunci A&B va fi o altă matrice 0-1 reprezentând produsul logic al elementelor corespondente. Acest operator consideră orice elemnet diferit de 0 ca fiind 1. Operatorul NOT este un operator unar. Expresia A returnează 0 acolo unde A este non zero şi 1 unde A este 0. Expresiile: P&(~P) P|(~P) Vor returna zerouri si respectiv elemente 1. Funcţiile any şi all sunt utile în confuncţie cu operatori logici. Dacă x este un vector de elemente 1 şi 0, atunci any(x) va returna 1 dacă cel puţin un element al lui x este diferit de 0, şi 0 în restul cazurilor. Funcţia all(x) va returna 1 dacă toate elementele lui x sunt nenule. Aceste funcţii sunt utile în declaratii de tipul.

If all(A<.5) execută secvenţă de operaţii end Unele din funcţiie logice disponibile în MATLAB sunt prezentate în tabelul următor: any condiţie logică all condiţie logică find găseşte indicii valorilor logice exit caută dacă variabilele există isnan detectează elemente care nu sunt numere finite detectează valoarea infinit isempty detectează o matrice fară elemente isstr detectează şiruri de caractere strcmp compară variabile de tipul şirurilor de caractere

1.5.5. Funcţii matematice elementare

MATLAB pune la dispoziţie o gamă largă de funcţii elementare şi trigonometrice, prezentate în tabelele următoare:

sin sinus cos cosinus tan tangent asin arc sinus acos arc cosinus atan arc tangent sinh sinus hiperbolic cosh cosinus hiperbolic tanh tangent hiperbolic asinh arc sinus hiperbolic acosh arc cosinus hiperbolic atanh arc tangent hiperbolic

Page 10: Indrumar de Laborator Matlab

Funcţii elementare abs valoarea absolută sau modul angle argumentul sqrt rădăcină pătrată real partea reală imag partea imaginară conj conjugatul unui număr complex round rotunjire către cel mai apropiat întreg sign funcţia signum rem restul împărţirii exp exponenţial, baza e log logaritm natural log 10 logaritm zecimal

1.6. Lucrul cu matrici şi vectori

1.6.1. Generarea vectorilor Simbolul : este foarte improtant în MATLAB. Declaraţia:

X=1:5 generează un vector conţinând elemente de la 1 la 5, cu incrementul I:

x=[1 2 3 4 5] Dacă se doreşte un alt increment (pozitiv sau negativ), atunci declaraţia va fi de tipul:

y=0:pi/4:pi O altă posibilitate de a genera vectori este utilizarea funcţiilor linspace sau logspace, care determină generarea unui număr de puncte într-un interval considerat. De exemplu,

k=linspace(-pi,pi,4)

conduce la: k= -3.1416 -1.0472 1.0472 3.1416

1.6.2. Utilizarea indicilor Elementele unei matrici pot fi accesate prin utilizarea indicilor lor, plasaţi între parenteze. O expresie utilizată ca indice este rotunjită către cel mai întreg. Exemplu:

A=(3,3)=A(1,3)+A(3,1)

În l9oc de indice se poate folosi şi un vector. Dacă, de exemplu, X şi V sunt vectori, atunci X(V) este [X(V(1)), X(V(2))……….., X(V(n))]. Pentru matrici, indicii se pot utiliza pentru a referi submatrici. Dacă A este matrice 10x10, atunci :

A(1:5,3)

va referi o submatrice cu 5 linii şi o coloană, constând din primele 5 elemente ale coloanei 3 a matricii A. În mod identic,

A(1:5,7:410)

este submatricea formată din elementele primelor 5 rânduri ale ultimelor 4 coloane. Dacă se utilizează doar cele două puncte, fară indici, atunci acestea se vor referi la toate elementele liniilor sau coloanelor. De exemplu,

A(:,3) A(1:5,:) reprezintă a treia coloană a matricii şi respectiv primele 5 rânduri ale matricii A. Declaraţia:

A(:,[3 5 10])=B(:,1:3)

va conduce la înlocuirea coloanelor 3,5 şi 10 ale matricii A cu primele 3 coloane ale matricii B. Funcţiei de plasamentul celor două puncte, se pot întâlni situaţiile:

Page 11: Indrumar de Laborator Matlab

a) A=[1 2; 3 4; 5 6] b=A(:)

va conduce la obţinerea unui vector coloană cu toate elementele nmatricii A. b) dacă acestea sunt plasate în partea dreaptă, ele pot fi utilizate pentru redimensionarea unei matrici. Pentru aceasta, matricea A trebuie să existe. În aceste condiţii, A(:) denotă o matrice de aceeaşi dimensiune cu A, dar cu conţinut schimbat în partea dreaptă. Dacă A are 3 rânduri şi 2 coloane, atunci:

A(:)=11:16 Aranjează un vector de 6 elemente într-o matrice 3x2.

1.6.2. Indexarea utilizând vectori de elemente 1 şi 0 Vectorii de acest gen, creaţi,de obicei, prin operaţii relaţionale, pot fi utulizaţi pentru a referi

submatrici. Dacă, de exemplu, A este o matrice m x n şi L este un vector de lungime m, format din elemente 0 şi 1, atunci:

A(L,:) specifică rândurile matricii A unde elemntele lui L sunt nenule. Pentru eliminarea unor elemente ale unui vector, se pot utiliza secvenţe de tipul:

x=x(x<=3*std(x)) sau

L=X(:,3)>100; X=X(L,:);

ultima conducând la înlocuirea lui x cu acele rânduri ale lui X ale căror a treia coloană este mai mică decât 100.

1.6.4. Matrici vide Declaraţia

x=[ ] determină generarea unei matrici de dimensiune 0 x 0. Utilizarea acestei comenzi nu va conduce la o eroare. Matricile vide există în spaţiul de lucru (deci declaraţia anterioară este diferită de declaraţia clear x, care şterge pe x din lista variabilelor curente) dar au dimensiunea 0. Funcţia exit poate fi utilizată pentru verificarea existenţei unei matrici, iar isempty indică dacă o matrice este vidă. O posibilitate de a şterge rânduri sau coloane ale unei matrici este de a le declara vide. Astfel,

A(:,[2 4])=[ ] va şterge coloanele 2 şi 4 ale matricii A.

1.6.5. Declararea matricilor mari

Matricile mari pot fi create din matrici mici. De exemplu,

C=[A eye(4);ones(A) A^2]

creează o matrice, presupunând că A are 4 rânduri.

1.6.6. Calculul valorilor proprii ale unei matrici Dacă A este o matrice de dimensiune n x n, cele n numere ce satisfac relaţia Av=λv poartă

denumirea de valori proprii ale matricii A. Ele sunt calculate utilizând funcţia eig(A) ce returnează valorile proprii într-un vector coloană. Valorile şi vectorii proprii se pot obţine utilizând comanda [V,D]=eig(A). V este matricea vectorilor proprii. Exemplu: A=[0 1 -1 ; -6 -11 6; -6 -11 5]; [V,D]=eig(A)

Page 12: Indrumar de Laborator Matlab

Rezultatul obţinut este: V= -0.7071 0.2182 -0.0921 0.0000 0.4364 -0.5523 -0.7071 0.8729 -0.8285 D= -1.0000 0 0 0 -2.0000 0 0 0 -3.0000 1.7. Controlul fluxului de operaţii

MATLAB are posibilitatea controlului fluxului de operaţii, exact ca orice alt limbaj de programare.

Structurile repetitive şi decizionare puse la dispoziţia utilizatorului sunt prezentate în continuare.

1.7.1 Bucla FOR MATLAB posedă o versiune mai veche a buclelor DO şi FOR. Acestea presupun executarea unei

instrucţiuni (sau secvenţă de instrucţiuni) de un număr determinat de ori. Structura generală a unei bucle FOR este: for v=expresie secvenţă de instrucţiuni end

Exemple de bucle For: for i=1 : m for j=1 : n A(i,j)=1/(i+j-1); end end A n=max(size(t));

for j=1 : n for i=1 : n A(i,j)=t(i)^(n-j); end end

A(:,n)=onese(n,1) for j=n-1:-1:1 A(:,j)=t.*A(:j+1) end 1.7.2. Bucla WHILE Această buclă este utilizată atunci când se doreşte executarea unei instrucţiuni(sau secvenţă de instrucţiuni) de un număr indefinit de ori, sub controlul unei condiţii logice. Forma generală este următoarea: while expresie secvenţă de instrucţiuni end

Page 13: Indrumar de Laborator Matlab

Exemple de bucle While: n=1 while prod(1:n)<1.e100, n=n+1; end n E=zeros(A); F=eye(A); k=1; while norm(E+F-E,1)>0 E=E+F; F=A*F/k; k=k+1; end

1.7.2 Comenzile IF şi BREAK Pentru a ilustra modul de folosire a acestor comenzi, se consideră două exemple. În primul este

efectuat un calcul dependent de semnul şi paritatea lui n, iar în cel de-al doilea este prezentată o problemă din teoria numerelor. Exemplul 1. If n<0 A=negative(n) elseif rem(n,2)==0 A=event(n) %par else A=odd(n) %impar End Exemplu 2. Se consideră un întreg pozitiv, oarecare. Dacă este par, se va divide cu 2, iar dacă este

impar, se va multiplica cu 3 şi se va aduna 1 la produs. Se va repeta acest algoritm pâna când se va ajunge la valoarea 1. În teoria numerelor se pune problema dacă există un întreg pentru care programul rulează la infinit. Programul ilustrează modul de folosire a comenzilor while şi if. De asemenea, sunt utilizate funcţiile input (acceptă introducerea unei secvenţe de la tastaură) şi break.

while 1

n=input(‘Introduceţi un intreg pozitiv: ‘); if n<=0, break, end while n>1 if rem(n,2)==0 n=n/2 else n=3*n+1 end; end end

1.8. Trasarea graficelor in MATLAB 1.8.1. Reprezentari bidemensionale MATLAB poate crea diverse tipuri de grafice, cum ar fi: 2D, 3D, liniare, logaritmice, semilogaritmice,

reprezentări în coordonate polare, etc. Câteva dintre modalitaţile de reprezentare a graficelor 2D sunt: plot, loglog, semilogx, semilogy, polar şi bar. Comanda grid adaugă linii ajutătoare graficelor, iar comenzile title(‘ text’), xlabel (‘text’), ylabel(‘text’) şi text (x,y,’text’) pot fi utilizate pentru plasarea de etichete sau texte pe grafice. Sintaxa comenzii pentru trasarea graficelor include şi o serie de simboluri opţionale ( . , +, * , o, x) precum şi câteva litere asignate culorilor (r,b,g,w) cu care lucrează mediul de programare. De exemplu, comanda:

plot(t, y1, ‘r’, t, y2, ‘+b’) va determina trasarea unor curbe de culoare roşie (y1) şi respectiv albastră (y2). Scalarea graficelor este efectuată automat. Comana: axis([xmin xmax ymin ymax])

Page 14: Indrumar de Laborator Matlab

forţează o scalare manuală a graficului. De exemplu, comanda: axis([-10 40 –60 60]) va produce o reprezentare cuprinsă intre -10 şi 40 pe axa x şi respectiv –60 şi 60 pe axa y. ăn continuare se vor prezenta câteva exemple ce ilustrează posibilitaţile MATLAB-ului de a reprezenta diverse funcţii. Exemplu 1. Să se traseze în coordonate x-y evoluţia:

t 0 1 2 3 4 5 6 7 8 9 10 11 12 y 0 0,55 1 2 4 7,2 11 14 15,1 16 16 16 16

Pentru un numar mic de date, acestea pot fi introduse în MATLAB utilizând paranteze pătrate: t=[0 1 2 3 4 5 6 7 8 9 10 11 12]; y=[0 0.55 1 2 4 7.2 11 14 15.1 16 16 16 16]; plot(t,y) grid Rezultatul execuţiei secvenţei de comenzi anterioare este afişarea următorului grafic:

Fig. 1.1. Evoluţia semnului din exemplu 1.

Pentru un numâr mare de date, se poate utiliza un editor de texte pentru a crea un fişier cu extensia .m. Tastarea comenzii nume_fişier va produce aducerea acelor date din fişier în spaţiul de lucru. Exemplul 2. Să se reprezinte grafic funcţia y=sin(x) pentru -4π x 4π. x=-4*pi: 0.05:4*pi; y=sin(x)/x; plot(x,y), title(‘Exemplul 2- Funcţiasin(x)/x’) xlabel(‘Radian’) ylabel(‘sin(x)/x) grid

Page 15: Indrumar de Laborator Matlab

Fig. 1.2. Evoluţia semnalului din exemplu 2. subplot va împărţi spaţiul pentru trasarea graficelor în mai multe subspaţii, această posibilitate de lucru fiind utilă atunci când se doreşte trasarea mai multor grafice în acelaşi timp. Comanda subplot(mnp) împarte spaţiul în m x n subspaţii, utilizându-l pe cel cu numărul de ordine p pentru reprezentarea grafică. Acest lucru este demonstrat în exemplul următor. Exemplul 3. a) să se traseze grafic funcţia v=120sinωt şi i=100sin(ωt-π/4) în raport cu ωt (în colţul din stânga sus). b) să se traseze p=vi (dreapta sus) c) dacă Fm=3.0 , să se traseze fa=Fmsinωt, fh=Fmsin(ωt-2π/3) şi fc=Fmsin(ωt-4π/3) în funcţie de ωt

(stânga jos). d) dacă fR=3Fm să se construiască o reprezentare în coordonate polare (dreapta jos).

ctg %şterge fereastră wt=0:0.05:3*pi; i=100*sin(ωt-pi/4); v=120*sin(ωt); p=v.*I; %puterea instantanee subplot(221),plot(wt,v,wt,I) %trasare v şi I în raport cu wt title(‘tensiune ţi curent’) xlabel(‘wt-radieni’); subplot(222),plot(wt,p) %puterea instantanee title(‘Putere’) xlabel(‘wt-radieni’); ylabel(‘Wa tt’); Fm=3.0 fa=Fm*sin(ωt); fb=Fm*sin(ωt-2*pi/3); fc= Fm*sin(ωt-4*pi/3); subplot(223),plot(wt,fa,wt,fb,wt,fc) title(‘Faze’) xlabel(‘wt-radieni’); fR=3/2*Fm(wt+1)/(wt+1); subplot(224),polar(wt,fR) title(‘Coord polare’)

Page 16: Indrumar de Laborator Matlab

Rezultatul execuţiei comenzilor de mai sus este prezentat în figura următoare:

Fig.1.3. Evoluţia semnalelor de la exemplu 3.

1.8.2. Reprezentări tridimensionale Comanda mesh(Z) crează o reprezentare tridimensională a elementelor matricii Z. O suprafaţă

tridimensională este definită de coordonatele x, y, z ale punctelor. Reprezentarea este formată prin unirea punctelor adiacente cu linii drepte. Funcţia meshdom transformă un domeniu specificat prin vectorii x şi y într-un tablou X şi Y. Un exemplu des întâlnit este reprezentarea funcţiei Bessel.

J022 yx +

in domeniul –12<x<12, -12<y<12. Următoarele comenzi: clg [x,y]=meshdom(-12:.6:12, -12:.6:12); r=sqrt(x.^2+y.^2); z=bessel(0,r); m=[-45,60]; mesh(z,m) vor conduce la afişarea următorului grafic:

Fig.1.4. Reprezentare tridimensională

Page 17: Indrumar de Laborator Matlab

1.9. Fişiere MATLAB este utilizat, cel mai adesea, în modul comandă. Atunci când la prompterul mediului este

procesată şi sunt afşate reultatele. Dar, acest mediu este capabil de a executa secvenţe de comenzi ce au fost scrise anterior în fişiere.

Fişierele ce conţin comenzi MATLAB sunt denumite fişiere de tip m, deoarece ele această extensie. Un fişier de tip m constă într-o secvenţă de instrucţiuni normale, fiind posibilă apariţa unor apeluri ale altor fişiere de tip m. Un astfel de fişier se poate apela el însuşi.

Una din unitaţile fişierelor m este înglobarea unor secvenţe lungi de comenzi. Astfel de fişiere se numesc fişiere script (text). O altă categorie de fişiere sunt fişierele funcţii. Ele permit adăugarea unor noi funcţii la cele existente. Capacitatea de a se crea noi funcţii, utile pentru rezolvarea unor probleme specifice, fac din MATLAB un mediu puternic şi flexibil.

Amândouă tipurile de fişiere prezentate sunt, în esenţă, foşiere ASCII şi pot fi create utilizând un editor de texte.

1.9.1. Fisiere script(script files)

Când este apelat un astfel de fişier, MATLAB execută comenzile pe care le găseşte în aceasta, fără a

aştepta introducerea lor de la tastaură. Aceste fişiere sunt utile în rezolvarea unor probleme ce implică un număr mare de instruncţiuni. Aceste instrucţiuni operează asupra variabilelor din spaţiul de lucru. După terminarea execuţiei, variabilele utilizate ramân în sapţiul de lucru.

De exemplu, secvenţa de instrucţiuni următoare poate fi scrisă într-un fişier f.m %Fisier m pentru calcularea numerelor lui Fibonacci f=[1 1]; i=1;

while f(i)+f(i+1)<100

f(i+2)=f(i)+f(i+1);

i=i+1;

end

plot(f)

Dacă la promterul mediului se tastează f, atunci se va executa secvenţa anterioară ce calculează o secventă a numerelor Fibonacci.

1.9.2. Fişiere funcţii (function files) Dacă prima linie a unui fişier conţine cuvântul function, atunci fişierul este de tip funcţie. O funcţie

diferită de un fişier text prin aceea că argumentele pot fi transmise, varibilele definite şi manipulate în interiorului fişierului fiind locale. De exemplu, o funcţie pentru calcularea valorii medii este următoare:

Function y=medie(x) %Calculeaza media. %Pentru vectori returneaza media elemnetelor. %Pentru matrici returneaza un vector linie ce contine %media elementelor fiecarei coloane a matricii. [m,n]=size(x); if m==1 m=n end y=sum(x)/m; Existenţa acestui fişier defineşte o nouă funcţie MATLAB. Presupunând că numele fişierului este

medie.m, atunci această funcţie poate fi utilizată ca orice alta. Elementele ce apar la un fişier funcţie sunt:

• prima linie conţine numele funcţiei, argumenetele ei precum şi argumentele ce vor fi returnateş • simbolul % indică faptul că restul liniei este un comentariu;

Page 18: Indrumar de Laborator Matlab

• liniile de comentariu de la începutul fişierului vor fi afişate ca informaţii ajutătoare, la tastarea comenzii help medie;

• variabilele m, y şi n sunt variabile locale; • o funcţie poate returna şi un număr mai mare de argumente sau poate primi un număr mai mare de

argumente, prima linie punând în evidenţă acest lucru. De exemplu, function[x, y] = nume functie(z) sau function x = nume functie(z, y, m). Se poate observa că numărul de argumente de intrare şi ieşire sunt elemente importante în construcţia funcţiilor. MATLAB este prevăzut cu două funcţii pentru aflarea acestor numere, şi anume nargin şi narg.

Page 19: Indrumar de Laborator Matlab

Capitolul II

2. Reprezentarea sistemelor automate 2.1. Introducere MATLAB lucreaza cu un singur tip de obiecte, si anume matrici numerice, eventual avand ca elemente numere complexe. Modulul Control System Toolbox poate fi utilizat pentru sisteme liniare, invariante in timp. 2.2. Modalitati de reprezentare asistemelor in MATLAB 2.2.1. Reprezentarea cu ajutorul ecuatiilor de stare Un sistem liniar, invariant in timp, poate fi reprezentat cu ajutorul unui set de ecuatii diferentiale de ordinul I. In spatiul starilor, aceste ecuatii pot fi scrise dupa cum urmeaza:

( ) ( ) ( )tButAxtx +=.

( ) ( ) ( )tDutCxty += unde u este vectorul marimilor de intrare, x este vectorul variabilelor de stare si y este vectorul variabilelor de iesire. Sistemele astfel reprezentate pot fi usor manipulate in MATLAB. De exemplu daca un sistem de intirziere de ordinul II, are polii caracterizati de pulsatia naturala 5.1=nω , iar factorul de amortizare este

2.0=ξ , ecuatiile de stare corespunzatoare sistemului pot fi reprezentate in MATLAB astfel: wn= 1.5; z=0.2; a=[0 1; -wn^2 –2*z*wn]; b=[0;wn^2]; c=[1 0]; d=0; Aceasta reprezentare este cea mai utilizata in MATLAB.

2.2.2 Reprezentarea sistemelor prin functii de transfer

O reprezentare echivalenta a sistemelor este functia de transfer. Aceasta este defenita ca fiind:

( ) ( ) ( )sUsGsY = unde ( ) ( ) DBAsICsG +−= −1 Se va exemplifica ulterior, modalitatea de trecere de la reprezentare prin ecuatii de stare la functie de transfer. 2.2.3.Factorizarea poli zerouri a functiei de transfer O functie de transfer poate fi factorizata sub forma:

( ) ( )( )

( )( ) ( )( )( ) ( )n

m

psKpspszsKzszs

ksNsMsG

−−−−−−

==21

21

Page 20: Indrumar de Laborator Matlab

unde ( )sN si ( )sM sunt polinoame in s , avind gradele n respectiv m. Conditia de realizabilitate fizica a sistemului se poate exprima prin mn ≥ . Radacinile celor doua polinoame pot fi calculate utilizand functia roots. Functia poly returneaza coeficientii polinomului ale carui radacini sunt elementele unui vector. Exemplu: Fie p=[1 3 5 2]. Atunci: r=roots(p) r= -1.2267 +1.4677i -1.2267 -1.4677I -0.5466 pp=poly(r) pp=

1.0 3.0000 5.0000 2.0000

2.2.4. Descompunerea in fractii simple

O functie de transfer poate fi descompusa in fractii simple. Pentru un sistem cu o singura intrare si o singura iesire:

( ) ( )( ) ( ) ( ) ( ) )(...

2

2

1

1 skps

rps

rps

rsNsMsG

n

n +−

++−

+−

==

Un vector p contine polii functiei de transfer, un vector r va contine reziduurile corespondente polilor din p. Functia de transfer poate fi convertita la aceasta reprezentare utilizand functia residue. Exemplu: Sa se descompuna in fractii simple:

( )44

19223

3

+++++=sss

sssG

b=[2 0 9 1]; a=[1 1 4 4]; [r, p, k]=residue(b,a) r=

0.0000-0.2500i 0.0000+0.2500i -2.0000

p= 0.0000+2.0000i 0.0000-2.0000i -1.0000 k= 2 Prin urmare, descompunerea in fractii simple este:

( )4

11

22225.0

225.0

122 2 +

++

−+=−

−++

++

−+=ssjs

jjs

js

sG

Functia [b,a]=residue(r, p, k) converteste dezvoltarea in fractii simple in reprezentarea polinomiala.

Page 21: Indrumar de Laborator Matlab

2.3. Reprezentarea sistemelor prin ecuatii de stare Ecuatia diferentiala a unui sistem liniar poate fi scrisa sub forma :

( ) ( ) ( )tButAxtx +=.

Aceasta ecuatie diferentiala de ordinul intai a unui sistem reprezinta ecuatia de stare a sistemului, unde x este vectorul de stare. Un avantaj al metodei reprezentarii in spatiul starilor este ca aceasta forma duce usor la determinarea solutiei ecuatiei prin metode de calcul discrete si/sau analogice.In plus, metoda reprezentarii in spatiul starilor poate fi usor extinsa la analiza sistemelor neliniare. Ecuatia de stare poate fi obtinuta dintr-o ecuatie diferentiala de ordin n sau direct din modelul sistemului, prin identificarea variabilelor de stare proprii. 2.3.1. Variabile de stare Pentru a ilustra modul de selectie a unui set de variabile de stare, consideram un model al unui sistem liniar de ordin n descris de urmatoarea ecuatie diferentiala:

( )tuyadtdya

dtyda

dtyd

n

n

nn

n

=++++ −

− 011

1

1 ....

unde ( )ty este iesirea sistemului si ( )tu este intrarea sa. Un model de stare pentru acest sistem nu este unic, fiind dependent de alegerea unui set de variabile de stare. O modalitate de alegere a unui set de variabile de stare, este urmatoarea:

1

.

2

1

........−=

=

=

nn yx

yx

yx

Vom exprima : 1

.

+= kk xx pentru 1,.....,2,1 −= nk

Apoi vom rezolva pentru n

n

dtyd

si vom inlocui y si derivatele sale prin variabilele de stare corespunzatoare,

rezultand:

( )tuxaxaxax

xx

xx

xx

nnn

nn

+−−−−−=

=

=

=

1121110

.

1

.

32

.

21

.

....

........

sau in forma matriciala:

( )tu

xx

xx

aaaaxx

xx

n

n

nn

n

+

−−−−

=

10:00

:

...1...000:...:::0...1000...010

:1

2

1

1210.

1

.

2

.1

.

ecuatia iesirii fiind: [ ]xy 0...001=

Page 22: Indrumar de Laborator Matlab

Exemplul 1. Functia ode2phv.m converteste o ecuatie diferentiala ce caracterizeaza un sistem in ecuatiile de stare ale sistemului. Pentru sistemul descris de:

( )tuydtdy

dtyd

dtyd 108642 2

2

3

3

=+++

aceasta functie va returna matricile A,B,C. Apelarea ei se poate face in urmatorul mod: ai=[2 4 6 8]; k=10: [A,B,C]=ode2phv(ai,k) rezultatul fiind:

−−−=

234100010

A

=

500

B [ ]001=C

Exsemplul 2. In cazul sistemelor electrice, variabilele de stare vor fi legate de elementele ce inmagazineaza energie. Pentru exemplificare, se considera sistemul din figura urmatoare:

2 H 4 Ω

0.25 F 0.5 F

iR iL

V1

+ -

iS

+ +

Fig 2.1.Schema electrica a circuitului Ecuatiile de functionare ale sistemului sunt:

02

05.0

4

025.0

12

22

1

1

=−+

=−+−

−=

=−+

ccL

Sc

Lc

ciR

RLc

vvdtdi

il

vi

dtdv

vvi

iidt

dv

Aceste ecuatii pot fi scrise sub forma intrare-stare-iesire, dupa cum urmeaza:

+

−−

−−=

S

i

L

c

c

L

c

c

iv

ivv

ivv

002001

05.05.0220401

2

1

2

1

Page 23: Indrumar de Laborator Matlab

2.3.2. Conversia: functie de transfer-spatiul starilor Documentatia “Control System Toolbox” contine un set de functii pentru conversia modelului dintr-o reprezentare in alta. Functia:

[A, B, C, D]=tf2ss(num, den) converteste sistemul dat sub forma functiei de transfer in reprezentare in spatiul starilor. Exemplu: Pentru un sistem descris de functia de transfer:

( )24269

2723

2

+++++=

ssssssG

apelurile urmatoare duc la determinarea matricilor ce definesc ecuatiile de stare ale sistemului. Mentionam ca aceasta reprezentare nu este unica. num=[1 7 2]; den=[1 9 26 24]; [A, B, C, D]=tf2ss(num, den) conduc la: A=[-9 –26 –24; 1 0 0; 0 1 0] B=[1; 0; 0] C=[1 7 2] D=0 2.3.3. Conversia: spatiul starilor-funstia de transfer Considerand ecuatiile de stare si ale iesirii:

( ) ( ) ( )tButAxtx +=.

( ) ( ) ( )tDutCxty += aplicand transformata Laplaces si rearanjand, se obtine: ( ) ( ) ( ) ( )sUDsUBAsICsY ⋅+⋅⋅−= −1 sau:

( ) ( )( ) ( ) DBSsICsUsYsG +⋅−== −1

[num,den]=ss2tf(A, B, C, D, i) converteste ecuatia de stare pentru functia de transfer pentru intrarea i. Exemplu: Se considera sistemul descris de urmatoarele ecuatii de stare:

uxxx

x

x

x

+

−=

00

10

321110010

3

2

1

.

3

.

2

.

1

Sa se determine functia de transfer a sistemului: ( ) ( ) ( )sUsYsG = A=[0 1 0; 0 0 1; -1 –2 -3]; B=[10; 0; 0]; C=[1 0 0]; D=[0]; [n,m]=ss2tf(A,B,C,D,1) n=

0 10.0000 30.0000 20.0000 m=

1.0 3.0000 2.0000 1.0000 printsys(n,m) num/den= 10 s^2 + 30 s + 20 ------------------------------- s^3 + 3 s^2 +2 s + 1

Page 24: Indrumar de Laborator Matlab

De asemenea, [z,p]=ss2tf(A,B,C,D,1) converteste ecuatia de stare in functia de transfer sub forma factorizata(pune in evidenta polii si zerourile functiei de transfer) [z,p]=ss2zp(A,B,C,D,1) z= -1 -2 p= -0.3376 + 0.5623i -0.3376 - 0.5623i -2.3247 2.3.4. Transformari ale starii 2.3.4.1. Diagonalizarea matricii A Unul din motivele diagonalizarii matricii A, presupunand ca avem valori proprii distincte, este ca acestea sunt toate localizate pe diagonala principala. Rezulta deci ca matricea de tranzitie a starilor este de asemenea diagonala,cu elementele de pe diagonala ttt nxKxx λλλ ,,, 21 .

Dandu-se sistemul liniar ( )tBuAxx +=.

, unde A are valori proprii distincte, este de dorit a se gasi matricea nesingulara P astfel incat aplicand transformarea: ( ) ( )tPxtx =&~ sa transforme ecuatia de stare de mai sus in forma canonica ( ) ( ) ( )tuBtxAtx ~~~ +=& cu a dat de matricea diagonala: APPA 1~ −= si BPB 1~ −= In general, sunt cateva metode de a gasi matricea P. Astfel, matricea P poate fi formata prin utilizarea valorilor proprii ale matricii A. Exemplu: Fie sistemul descris de ecuatiile de stare:

[ ]xy

uxxx

x

x

x

001

00

10

321110010

3

2

1

.

3

.

2

.

1

=

+

−=

[P,L]=eig(A);L este matricea valorilor proprii P este matricea corespunzatoare vectorilor proprii a=inv(P)*A*P; b=inv(P)*B; P= -0.7071 0.2182 -0.0921 0.0000 0.4364 -0.5523 -0.7071 0.8729 -0.8285 a= -1.0000 0.0000 0.0000 0.0000 -2.0000 0.0000 0.0000 0.0000 -3.0000

Page 25: Indrumar de Laborator Matlab

b= 2.8284 13.7477 10.8628 2.3.4.2. Transformari nesingulare ale starii In ecuatia de stare a unui sistem liniar invariant in timp:

BuAxx +=&

daca matricea

[ ]BABAABBS n 12 ... −= este nesingulara, atunci exista o transformare nesingulara:

( ) ( )tPxtx =&~ sau ( ) xPtx &~1−=

ca transforma ecuatia de stare de mai sus in forma:

( ) ( ) ( )tuBtxAtx ~~~ +=&

unde:

−−−−

=

naaaa

A

K

K

MMMMM

K

K

321

1000

01000010

~ si

=

10

00

~MB

1~ −= PAPA si PBB =~

Matricea P este data de:

=

−11

1

1

nAP

APP

PM

unde:

[ ][ ] 1121 100 −−= BAABABBP nKK

Exemplu. Fisierul ss2phv.m este dezvoltat pentru a indeplini transformarea de ami sus [a,b]=ss2phv(A,B) va returna a si b, reprezentand matricile ecuatiei de stare, dupa transformarea nesingulara a starii. A=[0 1 0; 3 0 2; -12 -7 -6]; B=[-1; 2; 3]; [a,b]=ss2phv(A,B)

Page 26: Indrumar de Laborator Matlab

a= 0 1.0000 0 0 0.0000 1.0000 -6.0000 -11.0000 -6.0000 b= 0 0 1.0000 2.3.5. Solutia ecuatiei de satre Solutia ecuatiei de stare liniara, neomogena:

( ) ( ) ( )tButAxtx +=&

poate fi obtinuta prin abordarea transformatei Laplace: ( ) ( ) ( ) ( )sBUsAXxssX +=− 0

sau: ( ) ( ) ( ) ( ) ( )sBUsxssX ⋅+⋅= φφ 0

unde: ( ) ( ) 1−−=Φ AsIs

( ) ( )[ ]sLt Φ= −1φ este cunoscuta ca matricea de tranzitie a starilor. Astfel, solutia ecuatiei de stare este:

( ) ( )[ ] ( ) ( ) ( )[ ]sUBsLxsLtx ⋅⋅Φ+Φ= −− 11 0

Putem, de asemenea, exprima ecuatia de mai sus in functie de ( )tφ si integrala de convolutie.

( ) ( ) ( ) ( ) ( ) τττφφτ

dtBuxttx −⋅+⋅= ∫0

0

Daca a este nesingulara, atunci ecuatia de mai sus poate fi simpliicata pentru a da raspunsurile la impuls, treapta si rampa. Pentru intrare de tip impuls, raspunsul este:

( ) ( ) ( ) ( )BKtxttx φφ +⋅= 0

Pentru intrarea treapta ( ) ,Ktu = raspunsul este:

( ) ( ) ( ) ( ) KBItAxttx ⋅⋅−+⋅= − ][0 1 φφ

Pentru intrarea rampa ( ) ,Kttu = raspunsul este:

( ) ( ) ( ) ( ) ( ) KBAtItAxttx ⋅−−+⋅=− ][0 12 φφ

Page 27: Indrumar de Laborator Matlab

2.3.6. Transformata Laplace a matricii de tranzitie a starii, ( )sΦ

( )sΦ este obtinut din algoritmul Fadeeva, dat de:

( ) ( )01

11

022

11

1

asaKsasaEKEsEs

AsIs nn

nn

nn

nn

++++++⋅+⋅

=−=Φ −−

−−

−−

unde matricile E sunt:

IEn =−1 IaAEE knknkn −−−− +=1 , k=1,2,3,…,n-1

Fisierul Itstm este dezvoltat pentru a calcula ( )sΦ in concordanta cu algoritmul de mai sus. Exemplu: A=[-2 -1; 2 -5]; Itstm(A) Matricile E, in ordinea descrescatoare a puterilor lui s sunt: E=

1 0 0 1 E= 5 -1 2 2 a= 1 7 12 Deci:

( ) ( )( )4322

15

1272215

1001

2 ++

+−+

=++

−+

=Φsss

s

ss

ss

2.3.7. Evaluarea lui ( )sΦ din valorile caracteristice ale lui A. Metode Cayley-Hamilton. Teorema Cayley-Hamilton stabileste ca daca ecuatia caracteristica a unei matrici patratice A este:

011 =−−− −

nnn aa Kλλ

atunci A satisface ecuatia matriceala:

01

1 =−−− − IAA nnn αα K

Deci, fiecare matrice patrtica satisface propria sa ecuatie caracteristica. Fie:

( ) ( ) ( ) 121

−−−−= nn

At AtKAtKItKe K

Page 28: Indrumar de Laborator Matlab

Se poate arata ca o ecuatie scalar, echivalenta cu ecuatia de mai sus este satisfacuta cand A este Inlocuit de ,λ astfel:

Ke

nnnn

n

n

n

t

=

12

13

233

12

2212

11

211

1

111

λλλ

λλλλλλλλλ

λ

K

MMMMM

K

K

K

unde valorile λ sunt valori proprii distincte ale matricii A. Cand soua valori proprii sunt egale, de exemplu cand ,32 λλ = atunci al treilea rand al matricii de mai sus este inlocuit de:

( )nn

t

kkkdtd

dde 1

11 ... −+++= λλλ

λ

Fisierul strm.m este dezvoltat bazandu-se pe metoda Cayley-Hamilton. Aceasta functie evalueaza matricea de tranzitie a starii in forma inchisa, Valorile proprii multiple sunt de multiplicitate 2. exemplul urmator demonstreaza aceste functii. Exemplu. Sa se determine matricea de tranzitie a starilor pentru sistemul de la exemplul anterior. A=[-2 -1; 2 -5]; Itstm(A)

( ) ( ) qPASIinvs /=−=φ where, P=s**(n-1)E(n-1)+S**(n-2)E(n-2)+…+E(0) Q=a(n)s**n+a(n-1)s**n-1+a(1)s+…+a(0) The E matrices in descending power of s are: E= 1 0 0 1 E= 5 -1 2 2 a= 1 7 12 Deci,

( ) ( )( )4322

15

1272215

1001

2 ++

+−+

=++

−+

=Φsss

s

ss

ss

Rezultatul este urmatorul:

( )

+−−+−−

=

−−

+

−−

= −−−−

−−−−−−

tttt

tttttt

eeeeeeee

eet 4343

434343

22

2211

1212

φ

Page 29: Indrumar de Laborator Matlab

2.3.8. Determinarea solutiei ecuatiei de stare Procedura practica pentru gasirea raspunsului in timp al unui sistem este utilizarea simularii digitale. Reprezentarea in spatial starilor ne perminte sa simulam un sistem pe calculator. O prima posibilitate de aflare a solutiilor este utilizarea functiilor ode23 si ode45. Pentru sistemele liniare continue in timp, documentatia aferenta pentru MATLAB, Control system toolbox, furnizeaza functiile [y,x]=impulse(A, B, C, D, iu, t) si [y,x]=step(A, B. C, D, iu, t) care obtine raspunsul la impuls si raspunsul la treapta utilizand ecuatia de stare. Functia [y,x]=lsim[A, B, C, u, t] simuleaza ecuatia de stare cu intrare arbitrara. Exemplu. Pentru un sistem descris de:

uxxx

xxx

+

−−−=

111

6116100010

3

2

1

3

2

1

&

&

&

[ ]xy 011=

se da:

( )

−=

5.05.0

10x

Se cere sa sedetermine x(t) si y(t), r(t) fiind intrare de tip treapta. A=[0 1 0; 0 0 1; -6 -11 -6]; B=[1;1;1]; C=[1 1 0]; D=0; x0=[1 0.5 -0.5]; t=0:.05:4; u=ones(1, length(t)); [y,x]=lsim(A,B,C,D,u,t,x0); plot(t,x,t,y) title(’Solutiile ecuatiei de stare’) xlabel(‘Timp-sec’) text(3.8,1.8,’y’),text(3.8,2.6,’x1’), text(3.8,-0.8,’x2’), text(3.8,-1.4,’x3’)

Exeplu. Pentru sistemul anterior sa se traseze y(t) si x(t), intrarea semnalului fiind r(t)=sin(2π t). A=[0 1 0; 0 0 1; -6 -11 -6]; B=[1;1;1]; C=[1 1 0]; D=0; x0=[1 0.5 -0.5]; t=0:.05:4; u=sin(2*pi*t); %genereaza un vector u [y,x]=lsim(A,B,C,D,u,t,x0); plot(t,x,t,y) title(’Solutiile ecuatiei de stare’) xlabel(‘Timp-sec’) text(.1,1.7,’y’),text(.1,1.25,’x1’), text(.1,.55,’x2’), text(.1,-1,’x3’)

Page 30: Indrumar de Laborator Matlab

2.3.9. Transfigurarea schemelor bloc Fisierul blkbuild (MATLAB Control system toolbox) si functia connect convertesc diagrame bloc la modele in spatial starilor. Blocurile functiei de transfer sunt numerotate secvential de la 1 la nr. De blocuri. Functia nblocks defineste numarul total de blocuri si blbblock converteste fiecare bloc la o reprezentare in spatial starilor nelegata (neconectata). Expresia [A,B,C,D]=connect(a,b,c,d,q,iu,iy) conecteaza blocurile in concordanta cu o matrice predefinita q care specifica interconectarile. Primul element al fiecarui rand al matricei q este numarul blocului. Restul elementelor indica sursa blocurilor ce intra in sumator. Cand intrarea in sumator este negative numarul blocului este introdus cu semnul minus. Elementele in si iy sunt doi vectori linie, memorand blocurile de intrare si iesire. In final, pentru a obtine functia de transfer globala, [num,den]=ss2tf(A,B,C,D,iu) calculeaza functia de transfer de la intrarea iu. Exemplu. Sa se determine reprezentarea intrare-stare-iesire a sistemului din figura urmatoare:

R(s)

1

1

5

2

0.5 4

4+s 2

1+s 3

1+s

C(s) + + - - -

1 2 3 4 5

6

7

8 Fig. 2.4. Schema bloc functionala a sistemului

n1=1; d1=1; % coeficientii fiecarei functii de transfer (numerator si numitor) n2=5; d2=1; n3=4; d3=[1,4]; n4=1; d4=[1 2]; n5=1; d5=[1 3]; n6=2; d6=1; n7=5; d7=1; n8=1; d8=1; nblocks=8; % numarul be blocuri blkbuild q=[1 0 0 0 0 % q=matricea configuratiei 2 1 -6 -7 -8 3 2 0 0 0 4 3 0 0 0 5 4 0 0 0 6 3 0 0 0 7 4 0 0 0 8 5 0 0 0]; iu=[1]; % intrarea sistemului iy=[8]; % iesirea sistemului [A,B,C,D]=connect(a,b,c,d,q,iu,iy) % realizeaza conexiunea intre blocuri si determina reprezentarea

% intrare-stare-iesire

Page 31: Indrumar de Laborator Matlab

[num,den]=ss2tf(A,B,C,D,1) % determina functia de transfer a sistemului A= -8.0000 -2.0000 -0.5000 4.0000 -2.0000 0 0 1.0000 -3.0000 B= 0.5000 0 0 C= 0 0 1 D= 0 num= 0 0.0000 0.0000 2.0000 den= 1.0000 13.0000 56.0000 80.0000 Deci functia de transfer a sistemului este data de:

( )805613

223 +++

=sss

sG

Control system toolbox contine 4 functii care sunt utilizate in constructia unui model. append combina functionalitatea a doua sisteme in spatial starilor formand un model mai complex. Functiile parallel si series conecteaza doua sisteme in spatiul starilor inparalel, respectiv in serie. In final, functia ode genereaza matricile A,B,C,D pentru un sistem de ordin doi.

Page 32: Indrumar de Laborator Matlab

Capitolul III

3.Răspunsul sistemelor automate

3.1. Răspunsul în timp al sistemelor automate

Evaluarea performanţelor în domeniul timp ale modelelor este importantă deoarece sistemele de control sunt studiate din două puncte de vedere : performanţe în domeniul timp şi în domeniul frecvenţei. Performanţele sistemelor în domeniul timp pot fi definite în termeni ai răspunsului în timp la intrări standard de test. Una din cele mai frecvent întâlnite intrări într-un sistem automat este funcţia treaptă (funcţia Heaviside). Dacă răspunsul la o intrare treaptă este cunoscut, matematic este posibil să calculăm răspunsul la orice intrare. Un alt tip de intrare, de o importanţă majoră, este funcţia sinusoidală. O ieşire stabilă sinusoidală este obţinută când un sistem liniar asimptotic stabil este supus unei intrări sinusoidale. De aceea, dacă cunoaştem răspunsul unui sistem liniar invariant în timp la intrări sinusoidale pentru toate frecvenţele, avem o descriere completă a sistemului.

3.1.1.Răspunsul sistemelor de întârziere de ordinul doi Forma standard a funcţiei de transfer de ordinul doi este dată de:

22

2

2)(

nn

n

sssG

ωξωω

++=

unde ωn este pulsaţia naturală. Pulsaţia naturală este pulsaţia de excitaţie, dacă toate amortizările sunt şterse. Valoarea sa ne dă o indicaţie asupra vitezei răspunsului. ξ este factorul de amortizare. Acesta oferă informaţii despre natura răspunsului tranzitoriu. Răspunsul în timp al unui sistem de întârziere de ordinul 2 este dat de:

)1sin(1

1)( 2

2ϕξω

ξ

ξω

+−−

−=−

tety n

ntn

unde δ este un factor care depinde numai de ξ.

3.2.Performanţele în domeniul timp

Criteriile de performanţă utilizate pentru a caracteriza răspunsul tranzitoriu la o intrare treaptă sunt: suprareglajul, timpul de creştere, timpul de suprareglare, timpul de răspuns. Vom defini timpul de creştere tc ca fiind timpul în care răspunsul sistemului creşte de la 10% din valoarea de regim staţionar la 90% din aceasta (dacă răspunsul nu atinge valoarea de regim staţionar) sau timpul după care răspunsul sistemului atinge prima dată această valoare. Timpul de suprareglare ts este timpul la care se obţine valoarea maximă a răspunsului sistemului. Viteza răspunsului este măsurată de tc şi ts. Timpul de răspuns reprezintă timpul după care răspunsul sistemului intră definitiv într-o bandă de ±5% din valoarea de regim staţionar. Suprareglajul sistemelor este definit ca:

(%)100.

.max.finalvaloare

finalvaloarevaloare −=σ

Timpul de suprareglare este obţinut prin anularea derivatei expresiei lui y(t).

Page 33: Indrumar de Laborator Matlab

21 ξω

π−

=n

st

Corespunzător acestui timp se obţine:

211)( ξ

πξ

−−

+= ety s Deci, suprareglajul pentru un sistem de ordinul 2 este:

(%)10021 ξ

πξ

σ −−

= e Timpul de răspuns este considerat ca fiind de ordinul a 4 constante de timp, ceea ce înseamnă că:

ts=4/ξωn

Funcţia timespec ( num, den ) returnează performanţele în domeniul timp al sistemelor automate: suprareglajul (overshoot), timpul de suprareglare (peak time), timpul de creştere (rise time) şi timpul de răspuns (settling time).

3.3.Efectele adăugării polilor şi zerourilor 3.3.1.introducerea unui zerou suplimentar în funcţia de transfer

Zerourile unei funcţii de transfer afectează amplitudinea răspunsului, dar nu şi natura acestuia. Dacă zeroul este apropiat de polii dominanţi ai sistemului, acesta va avea influenţe asupra răspunsului tranzitoriu al sistemului. Timpul de creştere şi timpul suprareglare vor avea valori mai mici în timp ce suprareglajul va creşte. Dacă poziţia zeroului este depărtată de polii dominanţi, atunci răspunsul este aproximativ identic cu cel al unui sistem de ordinul 2.

Adăugarea unui pol Deoarece polii funcţiei de transfer în buclă închisă sunt rădăcinile ecuaţiei caracteristice, ei controlează direct răspunsul tranzitoriu al sistemului . Timpul de creştere şi timpul de suprareglare vor avea valori mărite, suprareglajul va fi mai mic, acestea ducând la obţinerea unui răspuns mai greoi al sistemului. Pe măsură ce polul introdus se deplasează la stânga polilor dominanţi, efectul acestuia scade, deoarece termenul exponenţial introdus de el se anulează cu atât mai repede cu cât este mai depărtat de polii dominanţi . Dacă este suficient de depărtat , atunci sistemul poate fi aproximat cu un sistem de întârziere de ordinul 2. Răspunsul unui sistem de ordinul 3 cu un zerou suplimentar poate fi exprimat prin:

sK

ssTssa

sYn

n

)12)(1()1(

)( 2

2

+++++

=ξω

ω

Aplicând transformata Laplace inversă se obţine:

T

nn

nn

t

nn

nn eTT

TaTKte

TTaaKKty n

1

22

22

22

22

2 21)(

)1sin(2121

1)(

+−−

++−+−+−

−+=

ωωξωϕξω

ωωξωωξ

ξξω

unde:

ξξ

ξωωξ

ξωωξϕ

−−

−−−

−−

−=

222 111

11

arctgT

Tarctg

aa

arctgn

n

n

n

Page 34: Indrumar de Laborator Matlab

Funcţia y=stepzwn ( z , ωn , R , a , T , t ) determină răspunsul pentru sistemul de ordinul 3 descris anterior, unde z este factorul de amortizare, ωn pulsaţia naturală, R este amplitudinea treptei de intrare . Pentru un sistem de ordinul 2, a şi T sunt setate la 0. Control System Toolbox oferă utilizatorului funcţiile y=impuls ( num , den , t ), y=step ( num , den , u , t ) şi y=Isim ( num , den , u , t ) pentru determinarea răspunsului în timp al sistemului. Exemplul 1. Să se determine răspunsul în timp al sistemului de ordinul 2, cu ωn=5 şi ξ=0.6. num = [ 10 , 25] ; den = [ 0.16 1.96 10 25 ] ; t = 0 : 0.02 : 2 ; y = step(num, den, t); plot(t, y); xlabel(‘t-sec’) ylabel(‘y(t)’), grid, pause timespec(num, den) Rezultatele obţinute sunt: Peak time = 0.553333 Percent overshoot = 37.9675 Rise time = 0.206667

Settling time = 1.59

Exemplul 2. Să se determine răspunsul indicial, precum şi performanţele în domeniul timp ale sistemului descris de funcţia de transfer:

)256)(16.01(

)4.01(25)()()( 2 +++

+==sss

ssUsYsG

Programul scris în Matlab este următorul: num = [ 10, 25]; den = [0.16 1.96 10.25]; t = 0:0.02:2; c = step(num, den, t); plot(t, c), xlabel(‘t-sec’), ylabel(‘y(t)’), grid, pause timespec(num, den)

Page 35: Indrumar de Laborator Matlab

iar reprezentarea răspunsului indicial este prezentată în figura următoare. Performanţele în domeniul timp ale sistemului sunt date de: Peak time = 0.553333 Percent overshoot = 37.9675 Rise time = 0.206667 Settling time = 1.59

Exemplul 3. Schema bloc a unui servomecanism este reprezentată în figura următoare G1(s) R(s) + E(s) C(s) - G2(s)

Fig. 3.3. Schema bloc funcţională

d s(s + 1)

1 + es

Page 36: Indrumar de Laborator Matlab

Să se determine valorile constantelor d şi e astfel încât suprareglajul maxim al sistemului să fie de 40% şi timpul de suprareglare să fie 0.8 secunde. Programul scris în Matlab constă în următoarele linii: sigma = 40; tmax = .80; z = log (100/sigma)/sqrt(pi^2 + (log(100/sigma))^2) wn = pi/(tmax*sqrt(1 – z^2)) num = wn^2; den = [1 2*z*wn wn^2]; t = 0:0.02:4; timespec(num, den) Rezultatul execuţiei programului este: z = 0.2800 wn = 4.0906 Peak time = 0.803239 Percent overshoot = 39.9965 Rise time = 0.314311 Settling time = 3.37011 Din schema bloc rezultă funcţia de transfer a sistemului: G(s) = d/ [s2 + (de + 1)s + d] Ecuaţia caracteristică este: s2 + (de + 1)s + d = s2 + 2ξωns + ωn

2 = 0 Din valorile determinate anterior se deduce: d= ωn

2=4.09062=16.733 de +1 = 2(0.28)(4.0906) deci e=0.077

3.4.Răspunsul în frecvenţă al sistemelor Răspunsul în frecvenţă al sistemului este definit ca răspunsul sistemului atunci când la intrarea acestuia se aplică un semnal de tip sinusoidal. Să considerăm un sistem cu funcţia de transfer G(s) şi o intrare sinusoidală u(t)=Asinωn. Utilizând transformata lui r(t), transformata Y(s) a ieşirii sistemului este:

22

)()(ω

ω+

=s

sGAsY

Dezvoltarea în fracţii simple este:

∑++

+−

= G(s) lui polii de generati termeni'

)( 11

jsK

jsK

sY

Polii lui G(s), pentru situaţia în care sistemul este stabil, au valori reale negative sau, dacă sunt

complecşi conjugaţi au partea reală negativă. Ca urmare ei nu au aport la răspunsul sistemului, y(t).Deci, răspunsul în spaţiul sistemului este dat de transformata Laplace inversă a primilor doi termeni ai lui Y(s). y(t) = A/G(jω)/sin(ωt + θ ) Din această ecuaţie se poate vedea că ieşirea sistemului are aceeaşi frecvenţă ca şi intrarea, şi poate fi obţinută prin multiplicarea mărimii de intrare cu / G(jω) / şi modificarea fazei cu argumentul lui G(jω) şi unghiul θ pentru toate ω, constituie răspunsul în frecvenţă al sistemului. Corelarea între răspunsul în frecvenţă şi răspunsul tranzitoriu este indirectă, exceptând cazul sistemelor de ordin doi. În practică

Page 37: Indrumar de Laborator Matlab

caracteristica răspunsului în frecvenţă este pusă la punct prin utilizarea unor criterii diferite de proiectare, care vor avea ca rezultat obţinerea unui răspuns tranzitoriu acceptabil. Să considerăm răspunsul în frecvenţă al unui sistem de ordin întâi cu următoarea funcţie de transfer:

1

1)(+

=s

sGτ

Pentru s = jω se poate scrie:

)(

2211)( ω

ωτω Φ

+= jejG

unde Φ(ω) = -tan-1τω. Caracteristicile Bode pentru sistemul de mai sus ( τ=10 ) sunt date în figura următoare:

Lăţimea de bandă a sistemului este definită ca valoarea pulsaţiei la care mărimea răspunsului în frecvenţă este redus la 1/√2 din cea mai mică valoare a frecvenţei . Această pulsaţie este notată cu ωB. Pentru sistemul de ordinul întâi, lăţimea de bandă este dată de ωB = 1/τ. Acum se consideră funcţia de transfer standard a unui sistem de ordin doi

2n

2

2 s2)(

n

n

ssG

ωξωω

++=

Răspunsul în frecvenţă este dat de :

Page 38: Indrumar de Laborator Matlab

)(2/12

n22 /(2 /-((1

1)( ω

ωξωωωω Φ

))+)= j

nn

ejG

Reprezentarea răspunsurilor în frecvenţă pentru ξ =0.5 şi ωn = 2 daţi este dată în figura următoare:

Pentru ξ constant, creşterea lui ωn determină creşterea lăţimii de bandă cu acelaşi factor. Acesta corespunde la o descreştere a timpului de suprareglare şi a timpului de creştere a răspunsului tranzitoriu. Deci pentru creşterea vitezei răspunsului unui sistem este necesară creşterea lăţimii de bandă a sistemului. Pentru un sistem particular o legătură aproximativă este dată de : ωBtc ≈ constant unde această constantă are valoare aproximativ egală cu 2. Pulsaţia la care se obţine suprareglajul este obţinută prin anularea derivatei expresiei:

22

2

2)(

nn

n

sssG

ωξωω

++=

Pentru ξ < 0.707, frecvenţa de rezonanţă ωr este dată de:

221 ξωω −= nr Valoarea maximă a răspunsului la treaptă, notată prin Mpω este:

2p121M

ξξω

−=

Răspunsul în frecvenţă este obţinut utilizând funcţia MATLAB g=freqs( num, den, w). Această funcţie

evaluează g pentru domeniul specificat al frecvenţei ω. num şi den sunt vectori linie conţinând coeficienţii

Page 39: Indrumar de Laborator Matlab

numărătorului şi numitorului funcţiei de transfer. Ultima versiune de MATLAB pune la dispoziţia utilizatorului funcţia bode(num, den) care permite trasarea caracteristicilor Bode.

Funcţia frqspec(w, mag) este dezvoltată pentru returnarea lui ωr, Mpω şi ωB bazată pe argumentele w, mag.

Exemplu . Se dă sistemul descris de funcţia de transfer de ordinul 3:

75020536750)( 23 +++

=sss

sG

1.Să se determine polii dominanţi ai sistemului. 2.Să se găsească un model de ordin redus al sistemului. Să se determine timpul de creştere, de

suprareglare şi suprareglarea în procente pentru răspunsul indicial. De asemenea, să se determine pulsaţia de rezonanţă, modulul maxim şi lărgimea de bandă pentru sistem.

3.Să se determine valoarea exactă a parametrilor şi să se compare cu valorile de la punctul 2. Comenzile: a = [1 36 205 750] r = roots(a) vor determina polii funcţiei de transfer, care sunt: s1= -3 + j4, s2= -3 –j4, s3= -30. Deci funcţia de

transfer poate fi scrisă ca: G(s) = 750[(s - 30)(s2 + 6s + 25)] = 25/[(1 + 0.333s)(s2 + 6s + 25)] Deoarece polul s3 este destul de depărtat de polii dominanţi s1 şi s2, efectul lui poate fi neglijat. Deci

un sistem echivalent pentru sistemul iniţial este cel descris de funcţia de transfer: G(s) = 25/(s2 + 6s + 25) num1 = 25; den1 = [1 6 25]; %sistemul de ordinul 2 t = 0:.02:2; y1 =step(num1, den1, t) timespec(num1, den1) w = 0: .02:8; g = freqs(num1, den1, w) mag1 = abs(g) frqspec(w, mag1) num2 = 750; den2 = [1 36 205 750]; % sistemul iniţial y2 = step(num2, den2, t); timespec(num2, den2); g =freqs(num2, den2, w); mag2 = abs(g); frqspec(w,mag2) subplot(121), plot(t, y1, t, y2), xlabel(‘t-sec’),grid title(‘răspunsurile indiciale ale celor 2 sisteme’) subplot(122),plot(w, mag1, w, mag2), xlabel(‘w-rad/s’),grid title(‘modulele funcţiilor de transfer ’) subplot(111) Răspunsurile indiciale şi modulele funcţiilor de transfer ale celor 2 sisteme sunt reprezentate în

fig.3.5. Comparând cele rezultate se poate concluziona că aproximarea făcută duce la schimbări

neesenţiale în evoluţia celor 2 sisteme. Specificaţiile pentru performanţele celor 2 sisteme sunt următoarele: Peak time = 0.786667 Percent overshoot = 9.47783

Page 40: Indrumar de Laborator Matlab

Rise time = 0.373333 Settling time = 1.18667 Peak Mag. = 1.04 wr = 2.64 Bandwidth = 5.75 Peak time = 0.823333 Percent overshoot = 9.32926 Rise time =0.376667 Settling time = 1.22333 Peak Mag. =1.04 wr = 2.58 Bandwidth = 5.67

După cum se poate observa din figura următoare, aproximările făcute sunt suficient de bune.

Page 41: Indrumar de Laborator Matlab

Capitolul IV

4. Caracteristicile unui sistem de automat

Obiectivul unui sistem de control este de a controla ieşirea y in funcţie de intrarea u prin elemente ale sistemului automat. Câteva din caracteristicile esenţiale ale sistemelor de reglare sunt descrise in următoarele subcapitole.

4.1. Stabilitatea

Pentru ca un sistem sa poată fi utilizat, trebuie sa fie stabil. Un sistem liniar invariant in timp este stabil daca fiecare variaţie a intrării produce o variaţie a ieşirii. Aceasta caracteristică se numeşte stabilitate. Stabilitatea unui sistem poate fi privită din doua puncte de vedere: stabilitate IMEM şi stabilitate internă. Toate abordările din acest capitol vor fi considerate numai din punctul de vedere al stabilităţii IMEM.

Un sistem este stabil IMEM dacă si numai dacă polii sistemului in buclă închisă sunt plasaţi în semiplanul complex stâng. De aceea, o condiţie necesară si suficientă ca un sistem de reglare să fie stabil IMEM este ca toţi polii funcţiei de transfer a sistemului să aibă partea reala negativă.

Stabilitatea unui sistem liniar invariant în timp poate fi verificată prin utilizarea funcţiei impulse (Control System Toolbox) pentru a obţine răspunsul la impuls al sistemului. Sistemul este stabil daca răspunsul său la impuls tinde la zero când timpul tinde la infinit. O altă modalitate de a determina stabilitatea unui sistem este prin simulare.

Funcţia Isim poate fi utilizată pentru a observa ieşirea pentru intrări tipice. Aceasta este utilizată în particular pentru sisteme neliniare.

Funcţia MATLAB roots poate fi utilizată pentru a obţine rădăcinile ecuaţiei caracteristice. În teoria clasică a controlului, câteva tehnici au fost dezvoltate pentru analiza stabilităţii. Una din aceste tehnici este criteriul Routh-Hurwitz. În acest capitol se va prezenta utilizarea unui program simplu bazat pe criteriul Routh-Hurwitz. Se va considera că gradul de stabilitate al unui sistem oferă informaţii despre comportarea lui. Se pune întrebarea când un sistem este stabil si când este relativ stabil. Uzual, stabilitatea relativă este corelata cu viteza răspunsului si cu suprareglajul. Alte metode des utilizate pentru studiul stabilităţii sunt: diagramele Bode, locul rădăcinilor, criteriul Nyquist si criteriul Liapunov.

4.2. Criteriul de stabilitate Routh-Hunvitz

Criteriul Routh-Hurwitz furnizează o metodă rapidă pentru determinarea stabilităţii absolute care poate fi aplicată unei ecuaţii caracteristice de ordin n, de forma:

ansn +an-1sn-1+ ... +a1s+a0 =0 Criteriul este aplicat prin utilizarea schemei Routh definită ca:

sn an an-2 an-4 …….. sn-1 an-1 an-3 an-5 …….. sn-2 b1 b2 b3 …….. sn-3 c1 c2 c3 …….. … … … … ……..

an, an-1,.... , a0 fiind coeficienţii ecuaţiei caracteristice si:

b1=1

321

−−− −

n

nnnn

aaaaa

b2=1

541

−−− −

n

nnnn

aaaaa

…..

Page 42: Indrumar de Laborator Matlab

c1=1

2131

bbaab nn −− −

c2=1

3151

bbaab nn −− −

…..

Calculele sunt continuate pentru fiecare rând până când rămân doar elemente zero. O condiţie necesară si suficientă pentru ca toate rădăcinile ecuaţiei caracteristice să fie situate in semiplanul complex stâng este ca elementele primei coloane a schemei Routh să aibă acelaşi semn. Daca sunt schimbări de semn la elementele primei coloane, numărul schimbărilor de semn indică numărul rădăcinilor cu partea reală pozitivă.

O funcţie numită routh(a) construieşte schema Routh si determină rădăcinile cu partea reală pozitivă. a este un vector linie ce conţine coeficienţii ecuaţiei caracteristice.

Exemplu. Să se determine dacă sistemul având ecuaţia caracteristică:

s4 + 10s3+35s2+50s+2+24=0

este stabil sau nu.

a=[1 10 35 50 24]; routh(a)

Routh-Hurwitz Array 1 35 24 10 50 0 30 24 0 42 0 0 24 0 0

System is stable

Exemplu. Să se aplice criteriul lui Routh pentru următoarea ecuaţie caracteristică şi să se determine câte rădăcini sunt plasate în semiplanul complex drept.

s4+4s3-7s2-22s+24=0

a=[1 4 -7 -22 24] routh(a)

Routh/Hurwitz Array

1.000 -7.000 24.000 4.000 -22.000 0 -1.500 24.000 0 42.000 0 0 24.000 0 0

There are 2 roots in the right half s-plane.

4.2.1. Cazuri speciale

Cazul 1: Daca primul element dintr-un rând este 0, atunci este înlocuit printr-un număr pozitiv foarte mic ε si calculul tabloului este complet.

Exemplu. Să se utilizeze criteriul Routh-Hurwitz pentru a determina numărul rădăcinilor ecuaţiei caracteristice:

s4-5s2+20s+24 =0

a=[1 0 -5 20 24]; routh(a) Zero in the first column is replaced by 0.00001

Page 43: Indrumar de Laborator Matlab

Routh-Hurwitz Array

1.0000e+000 -5.0000e+000 2.4000e+001

1.0000e-005 2.0000e+001 0

-2.0000e+006 2.4000e+001 0

2.000e+001 0 0

2.400e+001 0 0 There are 2 roots in the right half s-plane

Cazul 2: Dacă toate elementele unui rând sunt zero, sistemul are poli pe axa imaginară ( perechi de rădăcini complex conjugate simetrice faţă de originea planului complex ), sau perechi de rădăcini reale cu semne opuse. În acest caz, o ecuaţie auxiliară este formată din rândul precedent. Tot rândul zero este apoi înlocuit cu coeficienţii obţinuţi prin diferenţierea ecuaţiei auxiliare.

Exemplu. Să se construiască schema lui Routh pentru sistemul automat a cărui ecuaţie caracteristică este dată de:

a=[1 10 36 60 59 50 24]; routh(a)

Elements of row 6 are all zero.

They are replaced by the auxiliary Eq. coefficients. Routh-Hurwitz Array

1 36 59 24 10 60 50 0 30 54 24 0 42 42 0 0 4 24 0 0 48 0 0 0 24 0 0 0

Characteristic Equation include roots on jw-axis or

pairs of real or complex roots symmetrical about jw-axis

4.3. Sensibilitatea

Unul din obiectivele utilizării reacţiei într-un sistem de control este de a reduce sensibilitatea sistemului la variaţiile parametrilor si la zgomot. Un sistem de control bun trebuie să fie insensibil la modificările parametrilor sau la perturbaţiile externe. Sensibilitatea unui sistem poate fi măsurată ca raportul între procentul schimbării funcţiei de transfer si procentul schimbării valorii unui parametru. De exemplu, sensibilitatea funcţiei de transfer T(s) la variaţia parametrului b este definită ca:

STb =

bbsTsT

/)(/)(

∆∆

=bsT

∆∆ )(

.)(sT

b

Cu ∆b tinzând la 0, sensibilitatea lui T la variaţiile lui b este :

STb =

bsT

∂∂ )(

.)(sT

b

Page 44: Indrumar de Laborator Matlab

Sensibilitatea statică este valoarea lui S pentru s→0. Sensibilitatea dinamică este calculată uzual prin înlocuirea lui s prin jω si reprezentând S ca funcţie de frecvenţa ω. Mărimea lui S(jω) măsoară erorile sistemului. Deci, se impune minimizarea lui. Aceasta este relativ uşor de făcut la o frecvenţă scăzuta datorită lărgimii de bandă scăzute a echipamentului fizic. Orice sistem fizic are o lăţime de bandă finită. De aceea funcţia de transfer a unui sistem de control real tinde totdeauna la 0, şi sensibilitatea sa tinde la unitate o dată cu creşterea lui ω(S(jω)→1 dacă ω→∞). Aceasta condiţie conduce la erori mari şi nu pot fi rejectate perturbaţiile

Exemplu. Se consideră sistemul reprezentat prin schema bloc funcţională din figura 4.1. unde b = 4 si h = 0.5.

1. Să se determine sensibilitatea lui G(s) în raport cu b. Să se traseze modulul "funcţiei de sensibilitate" în funcţie de pulsaţie pentru K=1 si K=0.5.

2. Să se determine sensibilitatea lui G(s) în raport cu h. Să se traseze modulul "funcţiei de sensibilitate" în funcţie de pulsate pentru K = 2 si K = 0.5.

Fig. 4.1. Schema bloc funcţionala a sistemului.

Funcţia de transfer a sistemului este dată de G(s) = Kb/(s+l+Kbh)

Pentru b = 4 si h = 0.5 lărgimea de bandă a sistemului este 1+2K. Sensibilităţile funcţiei de transfer G(s) În raport cu b si h, evaluate pentru valorile anterioare sunt:

S Kb =

Kbhss+++

11

=Ks

s211

+++

S Kh =

KbhsKbh++

−1

=Ks 21

2++

Următorul program calculează si afişează modulele funcţiilor calculate anterior:

clg k1=2; k2=0.5; num=[1 1]; den1=[1 1+2*k1]; den2=[1 1+2*k2]; w=0:.02:15; SGb1 = abs(freqs(num,den1,w)); SGb2 = abs(freqs(num,den2,w)); y1=0:.02:.721;wB1=5; x1 =wB1 *ones(1 ,length(y1)); y2=0:.02:.7905;wB1=2; x2=wB1 *ones(1 ,length(y2)); subplot(211),plot(w,SGb1, w, SGb2, x1, y1, x2,y2) text(2.1 ,.06;’wB=2'),text(5.1, .06,'wB=5'),text(11 ,-.2,'w-Radian/s') text(3,.77,'K=0.5'),text(3,.44,'K=2') title('Sensibilitatea lui G in raport cu b') num=-2; SGb1 = abs(freqs(num,den1,w));

K b/(s+1)

h

Y(s)R(s)

+ -

Page 45: Indrumar de Laborator Matlab

SGb2 = abs(freqs(num,den2,w)); subplot(212),plot(w,SGb1, w, SGb2), text(11 ,-.2,w-Radian/s'),text(2.7,.20,’K=0.5'),text(2.7,.6,'K=2’) title(' Sensibilitatea lui G in raport cu h') subplot(111)

Sensibilitatea in raport cu b

Fig. 4.2. Sensibilitatea sistemului in raport cu b si h După cum se poate observa din figură, sensibilitatea sistemului in raport cu b descreşte odată cu

creşterea factorului de amplificare al buclei deschise, în timp ce sensibilitatea în raport cu h creşte odată cu creşterea lui K.

4.4. Erori staţionare Pe lângă condiţia de stabilitate, un sistem automat trebuie sa îndeplinească anumite performante când

referita sa se schimba, sau când este supus unei acţiuni perturbatoare externe. Performanţele unui sistem pot fi privite nu numai din punctul de vedere al răspunsului tranzitoriu, ci si din punctul de vedere al erorilor staţionare. Eroarea staţionară este eroarea sistemului la valoarea timpului infinită. Un câştig mare corelat cu o sensibilitate redusă va duce la reducerea erorii staţionare.

Eroarea staţionară pentru un sistem automat este clasificată în concordanţă cu răspunsul la o intrare polinomială. Un sistem poate să nu aibă eroare staţionară la o intrare treaptă, dar acelaşi sistem poate să aibă o eroare staţionară nenulă la o intrare de tip rampă. Aceasta depinde de tipul funcţiei de transfer in buclă deschisă.

Se consideră sistemul prezentat in figura 4.3.

Funcţia de transfer in buclă închisă este :

)(*)(1

)()()(

sHsGsG

sUs

+=Υ

Fig. 4.3. Schema unui sistem cu reacţie neunitară.

G(s)

H(s)

Y(s)U(s) +

-

Page 46: Indrumar de Laborator Matlab

Eroarea sistemului în buclă închisă este:

E(s)=U(s)-H(s)*Y(s)= )(*)(*)(1

1 sUsHsG+

Utilizând teorema valorii finale, vom avea

ess=)(*)(1

)(*lim0 sHsG

sUss +→

Pentru intrări polinomiale, ca treaptă, rampă şi de tip parabolă eroarea staţionară din ecuaţia de mai sus va fi :

Intrare treaptă :

ess=ps

KisHsG +=

+→

1)(*)(lim1

1

0

Intrarea rampă :

ess=vs

KsHsGs1

)(*)(*lim1

0

=→

Intrare parabolică :

ess=as

KsHsGs1

)(*)(*lim1

2

0

=→

În scopul definirii tipului sistemului funcţiei de transfer generală în buclă deschisă sub forma următoare:

G(s)*H(s)=)*1(**)*1)(*1(*

)*1(**)*1)(*1( 21

sTKsTsTssTKsTsTK

nbaj

m

++++++

Funcţie de numărul de integratoare j pe care le conţine sistemul în buclă deschisă se pot face aprecieri asupra erorilor staţionare ale sistemului în buclă deschisă.

Erorile staţionare pentru sisteme de tip cu j = 0, 1 si 2 sunt prezentate in tabelul următor:

j/U(s) 1/s 1/s2 1/s3 0 l/(l+Kp) inf inf Kp=lim G(s)H(s)1 0 1/Kv inf Kv= lim sG(s)H(s) 2 0 0 1/Ka Ka= lim s2 G(s)H(s)

Valorile Kp, Kv si Ka poartă numele de coeficienţi generalizaţi ai erorii. ( p, v si a provin de la poziţie, viteza şi acceleraţie ).

Funcţiile errorzp(z,p,k) şi errortf(num,den) permit calculul erorilor staţionare la intrări tipice, cum ar fi treaptă, rampă şi de tip parabolă ( prima este utilizată atunci când sistemul este dat sub forma unei factorizate cu z = zerouri, p = poli, k = factor de amplificare, iar a doua atunci când sistemul este dat sub forma unei funcţii de transfer nefactorizate ). Dacă gradul numărătorului m este mai mic decât gradul numitorului n, atunci sunt n-m zerouri la infinit si vectorul z trebuie completat cu inf până la o dimensiune egală cu cea a vectorului polilor.

Page 47: Indrumar de Laborator Matlab

Exemplu. Să se determine erorile staţionare de poziţie, de viteza şi de acceleraţie, pentru sistemul având funcţia de transfer de forma:

G(s)=)5)(2)(1(

)4(10+++

−ssss

s

k=10;

z=[-4; inf; inf; inf]; p=[0; -1;-2; -5]; errorzp(z,p,k) System type is 1

Error Constants: Kp Kv Ka Inf 4 0 Steady-state Errors: Step Ramp Parabolic 0 0.2500 inf Exemplu : Să se determine erorile staţionare de poziţie, de viteză şi de acceleraţie, pentru sistemul având funcţia de transfer de forma:

G(s)=5014

102 ++ ss

num=10;

den=[1 14 50]; errortf(num, den)

System type is 0 Error Constants: Kp Kv Ka 0.2000 0 0 Steady-state Errors:

Step Ramp Parabolic

0.8333 Inf Inf

Page 48: Indrumar de Laborator Matlab

Capitolul V

5. Locul rădăcinilor - analiză şi proiectare

Metoda locului rădăcinilor ne permite să determinăm polii sistemului în buclă închisă cunoscând polii în bucla deschisă, pentru toate valorile factorului de amplificare al funcţiei de transfer în buclă deschisă. Locul rădăcinilor unui sistem este definit ca fiind reprezentarea poziţiilor rădăcinilor ecuaţiei caracteristice a sistemului atunci când factorul K este variabil. Bazat pe aceste afirmaţii, proiectantul trebuie să selecteze un factor de amplificare adecvat pentru a se obţine performanţele dorite. Dacă performanţele cerute nu pot fi atinse, poate fi adăugat sistemului un controler pentru a modifica locul rădăcinilor în maniera cerută.

Funcţia (Control System Toolbox) rlocus (num, den, K) obţine locul rădăcinilor pentru un domeniu specificat al lui K. Dat fiind faptul că se lucrează cu rădăcinile ecuaţiei caracteristice a sistemului, se poate trage concluzia că locul rădăcinilor este util pentru studiul stabilităţii IMEM a respectivului sistem Regulile de bază ale trasării locului rădăcinilor sunt prezentate pe scurt în cele ce urmează.

5.1. Metoda locului rădăcinilor

Să considerăm sistemul de control cu reacţie dat în figura urmatoare.

Fig. 5.l.Structura sistemului automat pentru analiza

locului rădăcinilor

Funcţia de transfer în buclă închisă a sistemului este:

(5.1)

În general, funcţia de transfer în buclă deschisă este dată de:

(5.2)

unde m este numărul zerourilor finite şi n este numărul polilor finiţi ai funcţiei de transfer a buclei. Dacă n > m sunt (n-m) zerouri la infinit.

Ecuaţia caracteristică a funcţiei de transfer a buclei închise este:

(5.3)

Atunci:

(5.4)

Din (5.4) se poate vedea că pentru ca un punct din planul s să aparţină locului rădăcinilor, când 0<K<∞, trebuie să fie satisfăcute următoarele două condiţii:

(5.5)

K

H(s)

G(s)U(s) Y(s)

-

( ))()(1

)()()(

0 sHsGKsGK

sHsYsG

⋅⋅+⋅==

)())(()())(()()(

21

21

n

m

psKpspszsKzszsKsHsGK

++++++=⋅⋅

0)()(1 =⋅⋅+ sHsGK

KzsKzszspsKpsps

m

n −=++++++

)())(()())((

21

21

KA

A

m

jz

n

ip

=∏

=

=

1

1

α

β

Page 49: Indrumar de Laborator Matlab

şi (5.6) unde K,5,3,1 ±±±=r şi )exp()(),exp()( αααβββ θθ zzpp jAzsjAps =−=−

5.2. Sumar al regulilor generale pentru construirea locului rădăcinilor 1. Numărul ramurilor. Pentru n m, numărul ramurilor este egal cu numărul polilor funcţiei de transfer G(.'i) *H(s) Locukl rădăcinilor este simetric faţă de axa reală. 2. Puncte de plecare şi sosire. Pentru K luând valon de la 0 la infinit, locul pentru bucla închisă porneşte din polii funcţiei de transfer corespunzătoare buclei deschise (K=0) si se termină în zerourile funcţiei de transfer a circuitului in bucla deschisă. Se spune că zerourile tind să atragă locul rădăcinilor şi polii tind să îl respingă. 3. Porţiunile locului rădăcinilor pe axa reală. Pentru K>0, locul rădăcinilor pe axa reală indeplineşte următoarea proprietate: un punct al locului rădăcinilor aparţine axei reale dacă la dreapta lui se află un număr impar de poli şi zerouri (situaţi pe axa reală). 4. Intersecţia cu axa imaginară Utilizând criteriul Routh-Hurwitz se determină intersecţiile locului rădăcinilor cu axa imaginară. Valorile lui K şi 0) pot fi obţinute din schema lui Routh. 5. Asimptote (Pentru n ≠ m). Pentru cele mai multe sisteme practice n > m Pentru n > m există (n-m) zerouri la infinit şi pentru 0<K<∞ (n-m) ramuri ale locului rădăcinilor sfârşesc în zerourile de la infinit.

Unghiurile dintre asimptote şi axa reală pozitivă sunt date de relaţia:

mnr

−⋅=

0180θ K,5,3,1 ±±±=r (5.7)

Astfel, se poate construi următorul tabel:

n-m Unghiurile asimptotelor 0 nu există asimptote

1 180°

2 ±90°

3 180°,±60°

4 ±45°,±135°

Asimptotele pornesc dintr-un punct a! axei reale având abscisa dată de:

mnsHszerouriGsHspoliiG

a −−

= ∑ ∑ )()()()(σ (5.8)

6. Unghiuri de plecare şi de sosire

Considerând un punct arbitrar s0 în apropierea polilor (pentru plecare) sau a zerourilor (pentru sosire) şi aplicând relaţia fundamentală pentru unghiuri se obţine:

(5.9+5.10)

unde r= ±1, ±3....

∑ ∑ =− )180( 0rzp αβ θθ

)180( 0ri

zii

pia +−= ∑∑ θθθ

)180( 0ri

pii

zid +−= ∑∑ θθθ

Page 50: Indrumar de Laborator Matlab

7. Puncte de ramificare Acestea sunt puncte pe axa reală unde două sau mai multe ramuri ale locului rădăcinilor pleacă

sau sosesc pe axa reală. Aceste puncte pot fi determinate exprimând pe K funcţie de s

(5.11) şi rezolvând:

Rădăcinile reale ale acestei ecuaţii care satisfac regula 3 reprezintă abscisele punctelor de separare sau de reintrare.

Exempiu 1. Să se traseze locul rădăcinilor pentru sistemul reprezentând un servomotor din figura 52.

Fig. 5.2. Schema bloc funcţională a servomotorului. Funcţia de transfer în circuit deschis este:

Pentru trasarea locului rădăcinilor este utilizată funcţia rlocus(num, den). Se pot pune în evidenţă regulile de trasare a locului rădăcinilor, astfel:

• sunt n-m=2 zerouri la infinit; • două asimptote făcând cu axa reală unghiurile de ±900 , • centrul de greutate a asimptotelor este: • punctele de ramificare de pe axa reală sunt date de :

adică 2−=s

K 1/|s(s+1)|U(s) Y(s) -

)()(1

sHsGK

⋅=

0)( ==SBSds

sdK

ssk

ssksHsG

4)4()()( 2 +

=+

=

22

24 −=−−=aσ

0)4( 2 =+−= ssdsd

dsdk

Page 51: Indrumar de Laborator Matlab

Locul rădăcinilor este trasat în figura 5.3.

Fig. 5.3. Locul rădăcinilor pentru exemplul 1

Exemplul 2. Să se traseze locul rădăcinilor pentru sistemul având funcţia de transfer în circuit deschis de forma:

• sunt n-m=3 zerouri la infinit; • trei asimptote făcând cu axa reală unghiurile de ±60° şi 180°; • centrul de greutate al asimptotelor este : • punctele de ramificare de pe axa reală sunt date de:

0)13197( 23 =+++−= sssdsd

dsdk

rădăcinile acestei ecuaţii fiind

Acestea fiind complexe, rezultă că nu există puncte de ramificare pe axa reală. Unghiurile de plecare din polii complecşi conjugaţi sunt de -45° si +45°. Intersecţia locului rădăcinilor cu axa imaginară poate fi determinată cu ajutorul schemei lui Routh. Aceasta este:

Impunând condiţia ca al treilea element al primei coloane să fie pozitiv, rezultă k = 120 fi respecliv s = ± j4.36. Prin urmare, pentru valori ale lui k mai mari decât 120, sistemul devine instabil IMEM.

Locul rădăcinilor este trasat in figura 5.3.

6

2 4

0 -2

-4 -6

8

-8 8 0 4 6 -4 -2 2 -8 -6

Imag Axis

Real Axis

13197)23)(23)(1()()( 23 +++

=++−++

=sss

kjsjss

ksHsG

33.23

133 −=−−−=aσ

94.033.294.033.2

jsjs

−−=+−=

07

12013

1971

kk

−+

Page 52: Indrumar de Laborator Matlab

Fig. 5.4. Locul rădăcinilor pentru sistemul din exemplul 2. 5.3. Proiectarea locului rădăcinilor

Specificaţiile de proiectare considerate aici sunt limitate de acurateţea sistemului şi de performanţele lui în domeniul timpului. Aceste performanţe pot fi definite în concordanţă cu localizarea polilor dominanţi ai sistemului în buclă închisă. Locul rădăcinilor poate fi utilizat pentru determinarea valorii caştigului K, care să ducă la o comportare satisfăcătoare în bucla închisă. Acest K este echivalentul unui controler proporţional şi valoarea lui are importanţă în stabilitatea sistemului. Valori mari ale câştigului determină instabilitate şi există nişte limite practice care impun ca factorul K să ia valori într-un anumit domeniu. Dacă locul rădăcinilor este de aşa natură încât performanţele impuse nu pot fi asigurate prin ajustarea lui K, atunci este necesar ca locul rădăcinilor sa fie reconfigurat prin introducerea unui controler adiţional Gc(s) în funcţia de transfer în circuit deschis. Gc(s) trebuie să fie aleasă astfel încât locul rădăcinilor să fie plasat într-o zonă corespunzătoare a planului s.

Pentru un regulator proporţional acţiunile sale sunt determinate de valoarea prezentă a erorii. Dar, un controler poate face corecţii şi în concordanţă cu valorile trecute şi viitoare ale erorii. Astfel, pot exista controlere cu acţiune proporţional integrală (PI), proporţional derivativă (PD) sau proporţional - integral - derivativă (PID).

(5.12)

Controlerele de tip integral şi diferenţial necesită utilizarea amplificatoarelor operaţionale. Alte controlere care pot fi utilizate, realizate doar cu elemente pasive, sunt controlerele cu avans, cu intârziere şi cu avans-intârziere de fază.

Un controler de ordinul I având un zero şi un pol în funcţia lui de transfer este :

Pentru proiectarea locului rădăcinilor au fost dezvoltate mai multe funcţii:

Functia Regulator [numopen, denopen,denclsd]=pcomp(num,den,ξ) Proporţional

[numopen. denopen, denclsd]=phlead(num,den,sl) Phase-Lead [numopen, denopen. denclsd]=phlag(num,den,ξ) Phase-Lag [numopen. denopen, denclsd]=pdcomp(num,den,sl) PD [numopen, denopen, denclsd]=picomp(num,den,sl) PI [numopen, denopen, denclsd]=pidcomp(num,den,sl) PID

6

2

4

0

-2

-4

-6

0 4 6 -4 -2 2 -6 Real Axis

Imag Axis

sKs

KKsG dl

pc ++=)(

0

0 )()(ps

zsKsG cc +

+=

Page 53: Indrumar de Laborator Matlab

Alternativ, funcţia [numopen, denopen, denclsd] = rldesign(num, den, s1) permite utilizatorului să selecteze unul din tipurile de controlere de mai sus. s1 =σ + Jω reprezintă polul dorit pentru funcţia de transfer în buclă închisă, exceptând pcomp şi phlag unde ξ; factorul de amortizare al polilor dominanţi este substituit cu s1.

num şi den sunt vectorii linie ai coeficienţilor polinoamelor funcţiei de transfer a sistemului în buclă deschisă (necompensată).

Funcţia phlead(num, den, s1) poate fi utilizată pentru proiectarea controlerelor cu întârziere. Pentru a face asta, poziţia lui s1 trebuie aleasă la dreapta polilor din sistemul necompensat.

5.4. Proiectarea regulatorului de tip P Proiectarea regulatorului de tip P constă în alegerea valorii lui K0 astfel încât să se obţină o

comportare satisfăcătoare a sistemului ( a răspunsului său indicial ). Funcţia [numopen, denopen, denclsd]=rldesign(num,den,ξ) afişază şase opăţiuni referitoare la proiectarea cu ajutorul locului rădăcinilor. Dacă se alege opţiunea 1 este apelată funcţia [numopen, denopen, denclsd]=pcomp(num,den,ξ) care va calcula valoarea K0, în concordanţă cu factorul de amortizare dorit al polilor dominanţi.

Exemplu. Să se determine valoarea K0 a regulatorului proporţional pentru un sistem a cărui funcţie de transfer în buclă deschisă este: astfel încât factorul de amortizare să fie egal cu 0.6. Să se traseze locul rădăcinilor, răspunsul indicial şi să se determine performanţele în domeniul timp ale sistemului în buclă închisă. num=1; den=[1 5 4 0];zeta=0,6; [numopen, denopen, denclsd]=rldesign(num,den,zeta); %Returnează

%numărătorul fdt pentru sistemul % deschis şi numitorii pentru sistemele închis şi deschis

k=0:.02:2; r=rlocus(numopen, denopen, k);

% Locul răd. se trasează pentru sistemul în buclă deschisă t=0:.01:10; c=step(numopen, denclsd, t); timespec(numopen, denclsd); % Performanţele în timp în buclă închisă subplot(221), plot(r,'.'); title('Locul rădăcinilor'), grid subplot(222), plot(t, c); title('Răspunsul indicial'), grid subplot(111) Controller gain: KO = 2.05 Row vectors of polynomial coefficients of the compensated system: Open-loop num. 0 0 0 2.0500 Open-loop den. 1 5 4 0 Closed-loop den 1.0000 5.0000 4.0000 2.0500 Roots of the compensated characteristic equation: -4.1563 -0.4219 + 0.5615i -0.4219-0.5615i Peak time = 5.85495 Percent overshoot = 9.28566 Rise time = 2.70229 Settling time = 8.69947

)4)(1(1)(

++=

ssssGd

Page 54: Indrumar de Laborator Matlab

Fig.5.5.Locul rădăcinilor şi răspunsul indicial al sistemului.

După cum se poate observa K0=2.05. Pentru această valoare, coeficientul Kv=2.05/ 4=05125. Eroarea staţionară de viteză este ess=1/Kv=1.95. Funcţia de transfer a sistemului în circuit închis este:

5.5. Proiectare "phase-lead" În relaţia (5.13) controlerul are o comportare de filtru trece sus dacă p0>zo. Reţeaua care

îndeplineşte o asemenea condiţie contribuie cu un unghi pozitiv la criteriul (5.6) al locului rădăcinilor şi are tendinţa de a deplasa locul rădăcinilor spre partea stângă a planului s. Efectele acesteia se manifestă asupra comportării dinamice a sistemului, ducând la o mărire a vitezei de răspuns. Ca urmare, o reţea cu avans de fază aproximează un termen derivativ. Dacă p0< z0, controlerul are o comportare de tip trece jos ( phase-lag ). Acesta contribuie cu unghi negativ la criteriui (5.6) al locului rădăcinilor având tendinţa de a-l deplasa spre partea dreaptă a planuiui s. Unghiul introdus trebuie să fie destul de mic pentru a menţine stabilitatea sistemului. Acest tip de controler este un controler de tip integrator.

Pentru controler se poate determina factorul lui de ampiificare care este:

Pentru o poziţie dorită a polului sistemului în buclă deschisă si, determinarea funcţiei de transfer a controlerului se poate face în următorul mod: se selectează o valoare corespunzătoare a lui z0 şi se utilizează criteriul (5.6) pentru determinarea lui p0. După aceasta, se determină Kc aplicând criteriul (5.5). Dacă factorul de ampiificare al controlerului a0 este specificat, atunci, pentru o poziţie dorită a polului sistemului în circuit închis (5.15) Z0 si p0 sunt obţinuţi astfel încât ecuaţia:

0)()(1 =++ sGsG dc (5.16) să fie îndeplinită. Parametrii pot fi determinaţi cu ajutorul

relaţiilor:

şi (5.17)

unde:

(5.18)

000 )0( pzKGa cc ==

05.24505.2)( 230 +++

=sss

sG

)exp(11 βjss =

100 aaz = 10 1 bp = 000 zpaKc =

ψψψβ

ψψββ

sinsin)sin(

sin)sin(sin

1

01

1

01

⋅++=

⋅−+=

sMab

MsMaa

Page 55: Indrumar de Laborator Matlab

M şi Ψ fiind amplitudinea şi faza funcţiei de transfer în circuit deschis evaluate în s1, astfel încât: (5.19)

Pentru cazurile în care Ψ este 0 sau 180, valorile a1 şi b1 îndeplinesc relaţia:

(520)

unde semnul + este aplicat pentru valoarea 0 iar - pentru valoarea 180. În această situaţie, zeroul controlerului trebuie de asemenea specificat. Bazată pe aceste relaţii, funcţia:

[numopen, denopen, denclsd]=phlead(num,den,s1) returnează funcţia de transfer a sistemului în circuit închis şi rădăcinile ecuaţiei caracteristice. Utilizatorul trebuie să introducă factorul de amplificare al controlerului. De asemenea, funcţia returnează numărătorii sistemelor închis şi deschis şi numitorul funcţiei de transfer în circuit închis. mum şi den reprezintă numărătorul şi numitorul funcţiei de transfer a sistemului în circuit deschis, iar s1 polul dorit al sistemului în circuit închis. Exemplu. Sistemul din exemplul anterior este nacesar ca să aibă o viteză mai mare de răspuns şi o eroare.staţionară de viteză mai mică. Controlerul cu avans de fază va fi proiectat astfel încât să îndeplinească următoarele condiţii: 1. Constanta de timp τ=1/ξωn=0.6667 sec 2. Factorul de amortizare ξ=0.6 3 Eroarea stationara de viteza ess=0.5.

Să se traseze locul rădăcinilor, răspunsul indicial şi să se specifice performanţele în timp ale sistemului în circuit închis.

Din primele două specificaţii ξωn=1/τ=1.5 şi θ=cos-10.6=53.13. Deci, poziţia polului s1 al sistemului în circuit închis este s1= -1.5+2J. Din a treia specificaţie rezultă: ess=1/Kv=0.5 sau Kv=2 unde:

Din relaţia anterioară se poate deduce a0=4, Kv=8. Deci factorul de amplificare al controlerului este 8.

num=1;den=[1 5 4 0];j=sqrt(-1);s1=-1.5+j*2; [numopen, denopen, denc!sd]=rldesign(num,den,s1); k=0:.02:2; r=rlocus(numopen, denopen, k); t=0:.01:3; c=step(numopen, denclsd, t); timespec(numopen, dendsd); subplot(121), plot(r,'.'); tit!e('Locul rădăcinilor'), grid subplot(122), plot(t, c); title('Răspunsul indicial'), grid subplot(111)

Compensator type Enter

Gain compensation 1 Phase-lead (or phase-lag ) 2 Phase-lag (Approximate K ÷ KO/Kc) 3 PD Controller 4 PI, Controller 5

)exp()( 1 ψjMsGd =

01coscos 011

11 =+±⋅

±⋅ aMM

sbsa ββ

41

)4)(1()()(lim

0

0

0

0 ⋅=+++

+=pzK

ssspszsKsK cc

v

Page 56: Indrumar de Laborator Matlab

PID Controller 6 To quit. 0

Enter your choice -> 2

Enter the compensator DC Gain -> 8 Gc(0) =8, Gc= 82.2812(s + 1.11407)/(s + 11.4583) Row vectors of polynomial coefficients of the compensated system: Open-loop num. 0 0 0 82.2812 91.6667 Open-loop den. 1.0000 16.4583 61.2917 45.8333 0 Closed-loop den 10000 164583 612917 1281146 91.6667 Roots of the compensated characteristic equation:

-12.2623 -1.5000+2.0000i -1.5000-2.0000i -1.1961

Peak time = 1.6136 Percent overshoot = 12.0833 Rise time = 0.710653 Settling time = 2.53327

Rezultatele sunt prezentate in figura 5.

Locul radacinilor 4,

Fig. 5 6. Locul radacinilor pentru regulatorul "phase-lead"

5.6. Proiectarea unui controler cu intârziere de fază

În această situaţie, polii şi zerourile controlerului sunt foarte apropiaţi, iar combinaţia lor este relativ apropiată de originea planului s. Deci, locul rădăcinilor este mişcat relativ puţin faţă de poziţia lui iniţială. Acest tip de controler este utilizat atunci când răspunsul tranzitoriu este satisfăcător dar se doreşte o reducere a erorii staţionare. Funcţia [numopen, denopen, denclsd]=phlead(num, den, s1) poate fi folosită pentru acest tip de proiectare cu specificaţia că polul s1 trebuie plasat uşor la dreapta polului sistemului în circuit deschis. O altă modalitate de a aborda această proiectare este una care se bazează pe nişte aproximaţii, prezentate în cele ce urmează:

Se consideră factorul de amplificare al controlerului ca fiind 1.

(5.21)

Răspunsul indicial

1)0( 0000 === pzKGa c

Page 57: Indrumar de Laborator Matlab

Deci, 00 zpKc = şi acum 100 << cKzp (5.22)

Dacă K0 este factorul de amplificare pentru polu sistemului închis s1, atunci din (5.3)

(5 23)

Daca se plaseaza polul si zeroul controlerului foarte apropiati unul de celalalt, cu amplitudinea mai mica decat a lui S1, atunci:

(5.24) Pentru a plasa polul sistemului în circuit închis aproximativ în s1 factorul de amplificare K necesar este dat de:

(5.25)

Deoarece Kc<l, atunci K>K0. Următorul pas este alegerea zeroului z0, foarte mic. Din relatia (5.21)

polul controlerului este: (5.26)

Funcţia de transfer a sistemului cu controler este:

(5.27)

Functia [numopen, denopen, denclsdi] = phlag(dum, den, si) permite obţinerea unui controler "phase-lag" pe baza metodei aproximative. num şi den sunt numărătorul şi numitorul funcţiei de transfer a sistemului în buclă deschisă, iar ξ este factorul de amortizare dorit pentru polii funcţiei de transfer în circuit închis. Utilizatorul trebuie să introducă valoarea factorului de amplificare K pentru a realiza eroarea staţionară dorită şi zeroul z0 al controlerului.

Exemplu. Pentru sistemul de la exemplul anterior este necesar să se obţină aproximativ aceeaşi poli dominanţi, acelaşi factor de amortizare, dar să se utilizeze un controler de tip phase-lag iar eroarea staţionară de viteză să fie ess = 0.125 Să se traseze locul rădăcinilor, răspunsul la impuls şi să se determine performanţele în domeniul timp ale sistemului în circuit închis.

Locul radacinilor

Fig. 5.7. Locul rădăcinilor.

Raspunsul indicial

)(1 10 sGK d−=

cc

c Kps

ssKsG ≈+

+=0

0 )()(

cdcdc KK

sGKsGsGK 0

111 )(1

)()(1 ≈≈=

000 zKp =

cccd GpszsKKGKG

0

0

++=

Page 58: Indrumar de Laborator Matlab

Factorul de amplificare K care rezultă din ess este obţinut pe baza relaţiei:

De aici rezultă K=32. Apelând funcţia Irdesign cu opţiunea 3, factorul de amplificare al

controlerului fiind 32 şi o poziţie a zeroului z0= -0.1 se obţin următoarele rezultate:

Enter your choice -> 3 Enter gain K required for the steady-state error specification -> 32 Enter magnitude of the compensator zero -> 0.1 Gain for the desired closed-loop pole KO = 2.05 Gain for the desired steady-state response K = 32 Gc(0) = 1, Gc = 0.0640625(s + 0.1 )/(s + 0.00640625) Row vectors of polynomial coefficients of the compensated system: Open-loop num. 0 0 0 2.0500 0.2050 Open-loop den. 1.0000 5.0064 4.0320 0.0256 0 Closed-loop den 1.0000 5.0064 4.0320 2.0756 0.2050

Roots of the compensated characteristic equation: -4.1530 -0.3646 + 0.5142i -0.3646 - 0.5142i -0.1242 Peak time = 6.03676 Percent overshoot = 28.1944 Rise time = 2.33421 Settling time = 22.0543

Rezultatele sunt prezentate în figura 5.7. După cum se poate observa, polii complex conjugaţi sunt plasaţi aproximativ în aceeaşi poziţie ca şi la exemplul anterior. Eroarea staţionară de viteză este sensibil redusă dar suprareglajul sistemului creşte.

Indicaţii de lucru

În mod asemănător proiectării regulatoarelor de mai sus, se pot sintetiza regulatoare de tipul P, PI, PD, PID. Să se utilizeze funcţiile menţionate pentru acestea astfel încât sistemele obţinute să fie stabile IMEM.

)4)(1(lim81

++===

∞→ sssKs

eK

sss

v

Page 59: Indrumar de Laborator Matlab

Capitolul VI 6. Analiza răspunsului în frecvenţă. Proiectare în frecvenţă.

Răspunsul în frecvenţă al unui sistem este răspunsul acestuia la un semnal de intrare sinusoidal.

Metoda răspunsului în frecvenţă şi metoda locului rădăcinilor reprezintă două căi simple, diferite, de aplicare a câtorva principii de baza ale analizei. Aceste metode se completează una pe cealaltă şi, în multe probleme practice de proiectare, sunt utilizate ambele tehnici. Un avantaj al metodei răspunsului în frecvenţă este ca funcţia de transfer a sistemului poate fi determinată experimental prin teste asupra răspunsului în frecvenţă. Mai mult, proiectarea unui sistem în domeniu frecvenţei permite proiectantului un control asupra lărgimii de banda a sistemului, asupra efectului zgomotului şi perturbaţiilor în răspunsul sistemului

In acest capitol diagrama Bode a funcţiei de transfer în buclă deschisă, specificaţiile marginii de câştig şi marginii de fază sunt obţinute utilizând funcţii din MATLAB Control System Toolbox. Este examinată de asemenea şi stabilitatea relativă a sistemului în buclă deschisă, utilizând criteriul Nyquist. In plus, sunt obţinute răspunsul în frecvenţă în buclă deschisă, amplitudinea maximă şi lărgimea de bandă. De asemenea, se va demonstra modul de utilizare a câtorva funcţii pentru proiectarea unui sistem de reglare în domeniul frecvenţei.

6.1. Răspunsul în frecvenţă

Răspunsul unui sistem liniar invariant m timp la un semnal de intrare sinusoidal )sin()( tAtr ω=

este dat de: )](sin[|)(|)( ωϕωω += tjGAty (6.1)

unde funcţia de transfer )( ωjG este obţinută înlocuind s cu . jω în expresia lui G(s). Funcţia de transfer rezultată poate fi scrisă astfel:

)(|)(|)( ωϕωω jejGjG ∗= (6.2)

Dar, funcţia de transfer poate fi reprezentată în planul complex sub următoarea expresie: )()()(Im)(Re)( ωωωωω jjXjRjGjGjG +=+= (6.3)

Cea mai utilizată reprezentare grafică a răspunsului în frecvenţă este diagrama Bode. Alte reprezentări ale funcţiei de transfer sinusoidale sunt contururile Nyquist.

6.1.1. Diagrame Bode

Diagrama Bode constă în două grafice trasate cu axa verticală liniară şi cea pe orizontală

logaritmică. Primul grafic este o reprezentare a modulului răspunsului )( ωjG în decibeli dB funcţie de lgω . Al doilea grafic al diagonalei Bode prezintă funcţia de fază )(ωϕ în funcţie de ωlg . Reprezentarea logaritmică este folositoare deoarece arată ambele caracteristici în frecvenţă, low şi high, ale funcţiei de transfer, pe aceeaşi diagramă. Mai mult, răspunsul în frecvenţă al sistemului poate fi aproximat printr-o serie de segmente de dreaptă.

Funcţia din MATLAB Control System Toolbox utilizată pentru trasarea diagramelor Bode este:

Page 60: Indrumar de Laborator Matlab

(mag, phase) = bode(num,den,ω ))

unde mum şi den reprezintă vectorii linie corespunzători polinoamelor de la numărătorul şi numitorul funcţiei de transfer a sistemului.

Exemplu. Pentru sistemul descris de fdt:

sssK

sssKsGH

10052)50)(2()( 23 ++

=++

=

să se traseye diagramele Bode.

Fig. 6 1. Diagramele Bode pentru exemplul considerat.

6.1.2. Locul de transfer

Locul de transfer, numit şi diagrama Nyquist, este o reprezentare a lui ImG(jω ) funcţie de ReG(jω ) cu ω variind de la ∞− la ∞+

Funcţia din MATLAB Control System Toolbox

[Re, Im] = nyquist(num, den, ω ) returnează partea reală şi partea imaginară a funcţiei de transfer pentru un domeniu specificat de frecvenţe

6.1.3. Caracteristica amplificare - fază Caracteristica amplificare - fază arată evoluţia amplificării în decibeli, funcţie de fază, pentru un

domeniu de pulsaţi de interes, reprezentat pe graficul Nichols. Graficul Nichols conţine linii de amplitudine şi fază constante, evidenţiind relaţia dintre răspunsul sistemului în circuit deschis şi cel în circuit închis.

6.2. Stabilitatea relativă

Funcţia de transfer în buclă închisă a unui sistem automat este dată de:

)(1)(

)()()(0 sKGH

sKGsUsYsG

+== (6.4)

Pentru ca sistemul să fie stabil IMEM polii lui 0G (s) trebuie să fie situaţi în jumătatea stângă a

planului complex. Deoarece zerourile lui 1+KGH(s) sunt poli pentru 0G (s), sistemul este stabil IMEM dacă rădăcinile ecuaţiei caracteristice 1+KGH(s)=0 se găsesc în jumătatea stângă a planului s. Locul rădăcinilor, care reprezintă poziţia rădăcinilor ecuaţiei caracteristice pentru K variind de la zero la infinit, a tost dezbătut în capitolul anterior. Toate punctele de pe locul rădăcinilor satisfac următoarele condiţii:

1|)(| =sFGH şi 0180)(arg −=sGH (6.5)

Cu creşterea lui K, sistemul devine asimptotic stabil când părţile reale ale rădăcinilor complexe dominante sunt zero. Aceasta corespunde unei intersecţii a locului rădăcinilor cu axa imaginară, când s=jω . Pentru acest sistem, valoarea lui K pentru care sistemul este la iimita de stabilitate este K c = 5200

101−

100

101

102

Amplificare-pulsaţie 50 0 -50 -100

Page 61: Indrumar de Laborator Matlab

Intersecţia locului de transfer cu axa reală negativă are unghiul de fază de -180°. Frecvenţa pcω corespunzând acestui punct este cunoscută ca frecvenţă de tăiere. În plus, cu cât câştigul buclei este crescut, locul de transfer trecând prin punctul (-1,0) are proprietatea descrisă de: 1|)(| =pcc jGHK ω şi 0180)(arg −=pcjGH ω (6.6)

Răspunsul în buclă închisă devine asimptotic stabil când amplitudinea răspunsului în frecvenţă este 1 si unghiul de fază este -180 0 . Frecvenţa la care locul de transfer intersectează punctul (-1,0) este aceeaşi frecvenţă la care locul rădăcinilor intersectează axa jω . Pentru o valoare mai mare a lui K, locul de transfer va înconjura punctul (-1,0) şi sistemul este instabil.

Deci, sistemul este stabil dacă:

1|)(| ⟨ωjKGH şi 0180)(arg −=pcjGH ω (6.7)

Aprecierea reprezentării lui KGH(Jω ) în coordonate polare faţă de punctul (-1,0) dă o indicaţie asupra stabilităţii sistemului în buclă închisă.

6.2.1. Marginile de amplificare şi de fază.

Marginea de amplificare şi marginea de fază sunt două criterii obişnuite de proiectare, relative la

răspunsul în frecvenţă în buclă deschisă. Marginea de amplificare este valoarea cu care trebuie crescut factorul de amplificare unui sistem stabil pentru ca locul de transfer să treacă prin punctul (-1,0). Marginea de amplificare este definită ca:

margineaKKmaeamplificarde c=)(__ (6.8)

unde K c este câştigul critic al buclei pentru limita de stabilitate şi K este câştigul actual. Raportul anterior poate fi scris astfel:

dGHKjGHK

jGHKma

pcpc

pcc 1|)(|

1|)(||)(|

===ωω

ω (6.9)

În decibeli, marginea de amplificare este:

djKGHMGma pcdb 101010 log20)(log20.).(log20 −=−== ω (6.10) Marginea de amplificare este factorul cu care trebuie modificat K pentru a face sistemul instabil. În

cazul în care se studiază stabilitatea unui sistem, folosirea numai a marginii de amplificare nu este de ajuns, mai ales în situaţiile în care parametrii sistemului afectează faza. Deci mai este necesară o altă mărime pentru a indica gradul de stabilitate a sistemului. Dacă se notează cu pcω pulsaţia la care amplificarea sistemului este 0 ( deci modulul funcţiei de transfer este 1 ), atunci marginea de fază reprezintă unghiul cu care trebuie rotit locul de transfer în jurul originii, astfel încât acesta să treacă prin punctul critic. Marginea de fază este dată de relaţia:

))180()(arg( −−=Φ gcm jGH ω (6.11)

Pentru performanţe satisfăcătoare, marginea de fază trebuie să fie cuprinsă între 30 o şi 60°, iar marginea de amplificare trebuie să fie mai mare decât 6 B.

Funcţia margin(mag,phase,ω ) permite evaluarea marginilor de amplificare şi de fază, precum şi ale pulsaţiilor la care acestea sunt obţinute.

Page 62: Indrumar de Laborator Matlab

6.3. Proiectarea in domeniul frecventelor Răspunsul în frecventă al unui sistem poate fumiza informaţii despre răspunsul staţionar al sistemului,

despre marginea de stabilitate şi lărgimea de bandă, Răspunsul tranzitoriu al sistemului poate fi estimat indirect în concordanţă cu marginea de fază, marginea de amplificare şi amplificarea maximă a sistemului. Suprareglajul este redus odată cu creşterea marginii de fază, iar viteza de răspuns creşte odată cu creşterea lărgimii de bandă. Deci pulsaţia la care se obţine amplificare 0, pulsaţia de rezonanţă şi lărgimea de bandă oferă informaţii despre viteza răspunsului tranzitoriu.

La fel ca şi la proiectarea prin metoda locului rădăcinilor, o posibilitate de proiectare în domeniul frecvenţei constă în modificarea factorului de amplificare al sistemului în circuit deschis, astfel încât să se obţină performanţele dorite. Acesta este în esenţă un regulator de tip P. Dacă nu sunt îndeplinite condiţiile impuse asupra marginilor de amplificare şi fază, atunci este necesară introducerea unui regulator (controler) având funcţia de transfer G(s). Aceasta funcţie de transfer poate fi aleasă astfel încât sistemul să îndeplinească condiţiile impuse. Ca şi la proiectarea prin metoda locului rădăcinilor, elementul având funcţia de transfer G c (s) poate fi realizat cu componente active,

sKs

KKsG Dps ++= 1)( (6.17)

sau cu componente pasive:

0

0 )()(ps

zsKsG cc +

+= (6.18)

Funcţiile dezvoltate pentru proiectarea în domeniul frecvenţei sunt apelate opţional la apelul funcţiei [numopen, denopen, denclsd]=frdesign(num, den) unde mum si den au semnificaţia cunoscută. Opţiunile sunt: frqp ( regulator proporţional ), frqlead ( regulator cu avans de fază ), frqlag (regulator cu întârziere de fază ), frqpd, frqpi, frqpid.

6.4. Proiectarea in frecventa a regulatorului proporţional

Proiectarea unui astfel de regulator constă în alegerea valorii K p pentru sistemul fără regulator

astfel încât să fie îndeplinite condiţiile de eroare staţionară. Atunci când K p variază, faza nu va fi afectată

Caracteristica amplificare-pulsaţie va fi ridicată sau coborâtă dacă valoarea lui K p creşte sau scade. La apelul funcţiei frdesign cu opţiunea corespunzătoare, utilizatorul trebuie să introducă

valoarea dorită a lui K p .

Exemplu. Să se traseze diagramele Bode şi să se determine marginile de amplificare şi de fază pentru sistemul având funcţia de transfer în circuit deschis:

sssssssG

458

)4)(1(8)( 23 ++

=++

=

Să se determine valoarea lui K p astfel încât eroarea staţionară la o intrare de tip rampă să fie 0.25. Să se determine răspunsul şi marginile de amplificare şi de fază pentru sistemul compensat.

Din valoarea erorii staţionare se determină valoarea lui K p 25.0/1 == vss Ke deci 4=vK . Dar:

pp

v Ksss

KsK 2

)4)(1(8

lim =++

=

de unde 2=pK . Apelul funcţiei frdesign cu opţiunea 1 şi introducerea valorii lui K p va duce la obţinerea următoarelor rezultate: Enter your choice -> 1 Uncompensated control system Gain Margin = 2.5 Gain crossover w = 1.22 Phase Margin = 22.5 Phase crossover w = 2

Page 63: Indrumar de Laborator Matlab

Enter the desired gain factor Kp -> 2 Gain & Phase Margins with gain compensation, Kp = 2 Gain Margin = 1.25 Gain crossover w= 1.79 Phase Margin = 5.21 Phase crossover w = 2 . Peak Mag. = 11.6 wr= 1.8 Bandwidth = 2.75 Row vectors of polynomial confidents of the compensated system: Open-loop num. 0 0 0 16 Open-loop den. 1 5 4 0 Closed-loop den 1 5 4 16 Roots of the compensated characteristic equation: -4.8549 -0.0725+1.81391 -0.0725-1.8139i

Fig.6.2. Caracteristicile Bode pentru exemplul considerat.

Se poate observa ca prin mărirea factorului de amplificare se micşorează marginile de fază şi de amplificare. Aceasta va duce la o stabilitate proastă a sistemului. Este necesară o reproiectare a sistemului, cu un alt tip de controler pentru a modifica răspunsul în frecvenţă al sistemului.

6.5. Proiectarea regulatorului cu avans de fază

Acest tip de regulator este obţinut dacă. 00 zp > . El contribuie cu un unghi pozitiv la faza sistemului şi are tendinţa de a mări marginea de fază şi de a îmbunătăţi stabilitatea sistemului. De asemenea este nevoie să se mărească pulsaţia de tăiere. Această mărire va duce la mărirea lărgimii de bandă rezultând un răspuns tranzitoriu mai rapid.

Dacă 0ω şi pω sunt pulsaţiile marginale ale regulatorului cu avans de fază, atunci valoarea maximă

a fazei regulatorului se obţine pentru pulsaţia mω dată de media geometrică a celor două frecvenţe. Este necesar ca să se selecteze polul şi zeroul controlerului astfel încât avansul de fază să apară la noua frecvenţă de tăiere fără o deformare a caracteristicii amplitudine pulsaţie în apropierea acestei frecvenţe

Metodologia de proiectare include specificaţii asupra marginii de fază şi erorii staţionare necesare. Se poate scrie:

0

00 )0(

pzKGa c

c == (6.20)

01

01

0

0

)()()(

bsbasa

pszsKsG c

c ++=

++= (6.21)

Caracteristica amplificare-pulsaţie (db)

101−

100

101

50 0 -50

101−

100

101

Fază-pulsaţie (grade) 0 -100 -200 -300

Page 64: Indrumar de Laborator Matlab

unde

01 / pKa c= ; 01 /1 pb = ; 10 =b ;

a 0 poate fi determinat din condiţia de eroare staţionară. Pentru o eroare staţionară de viteză e ss , se poate scrie:

)(lim)()(lim10 ssGHasGHssG

eK c

ssv === (6.22)

Din relaţia de mai sus se determină a 0 . Ecuaţia caracteristică a sistemului în circuit închis este:

0)()(1 =+ sGHsGc (6.23)

Dacă mΦ este marginea de fază la noua frecvenţă de tăiere gcω , atunci din (6.23) se poate scrie:

)]180(exp[1)()( mgcgcc jjGHjG Φ+−=ωω (6.24)

Înlocuind )( gcc jG ω şi egalând părţile reale şi imaginare ale ecuaţiei de mai sus se obţin parametrii regulatorului astfel:

θωθ

sincos1

1 MMa

cg

−= (6.25)

θωθ

sincos 0

1cg

Mab −= (6.26)

unde θ este argumentul funcţiei de transfer a controlerului evaluat în gcω .

Ψ−Φ+−== mgcc jG 180)(arg ωθ (6.27)

unde M şi Ψ sunt modulul şi argumentul funcţiei de transfer a sistemului în circuit deschis evaluate în j gcω adică:

)exp()( Ψ= jMjGH gcω (6.28)

Pentru un controler cu avans de fază argumentul controlerului trebuie să fie pozitiv, deci din (6.27) rezulta:

mΦ−=Ψ 180 (6.29)

De asemenea, pentru un astfel de controler )(0 gcc jGa ω< , şi din (6,24) rezultă:

10 <Ma (6.30)

Pentru un controler stabil este necesar ca 1a şi 1b să fie mai mari ca 0.

Proiectarea constă în apelarea funcţiei frdesign cu opţiunea 2. Utilizatorul trebuie să introducă marginea de fază dorită şi factorul de amplificare al controlerului G c (0). Programul va determina un domeniu de pulsaţii pentru un controler stabil. Utilizatorul va specifica o pulsaţie din acest domeniu, şi programul va determina funcţia de transfer a sistemului, precum şi performanţele în domeniul frecvenţei.

Exemplu. Să se proiecteze un regulator pentru sistemul de la exemplul anterior, astfel ca sistemul să aibă o margine de fază de 45 grade şi o eroare staţionară de 0.25 la o intrare de tip rampă.

Page 65: Indrumar de Laborator Matlab

Apelarea funcţiilor sus menţionate va duce la obţinerea următoarelor rezultate:

Eter your choice -> 2 Enter the compensator DC Gain -> 2 Enter desired Phase Margin -> 45 For a stable controller select a compensated gain crossover frequency wgc between 2.32 and 4.82 Suggested wgc for max. phase lead is 3.57 Enter wgc -> 3.57 w1 =

3.5700 Uncompensated control system Gain Margin = 2.5 Gain crossover w = 1.22 Phase Margin = 22.5 Phase crossover w = 2 Controller transfer function

Gc(0) =2, Gc = 83.535(s + 0.821555)/(s + 34.3143) Row vectors of polynomial coefficients of the compensated system: Open-loop num. 000 668.2798 549.0286 Open-loop den. 1.0000 39.3143 175.5714 137.2572 0 Closed-loop den 1.0000 39.3143 175.5714 805.5369 549.0286 Gain Margin = 8.27 Gain crossover w = 3.57 Phase Margin = 45 Phase crossover w = 12 Peak Mag. = 1.31 wr = 3.6 Bandwidth = 5.95 Roots of the compensated characteristic equation: -34.9359 -1.7915+4.0680i -1.7915-4.0680i -0.7954 Peak time = 0.804629 Percent overshoot = 22.0277 Rise time = 0.352025 Settling time = 1.96128

Dacă se doreşte cunoaşterea performanţelor pentru sistemul necompensat, se obţin următoarele valori:

Peak Mag. = 2.63 wr=1.3 Bandwidth = 1.95

Peak time = 2.63 Percent overshoot = 22.8697

Rise time = 0.97227

Settling time = 15.198 După cum se poate observa din valorile obţinute, utilizarea unui astfel de controler pentru mărirea

marginii de fază, va avea ca rezultat creşterea pulsaţiei de tăiere. Pentru sistemul în circuit închis, pulsaţia la care se obţine amplificarea maximă este redusă iar lărgimea de bandă este mărită. Totuşi, pentru regulatorul proiectat, raportul între factorul de amplificare pentru frecvenţe înalte şi cel pentru frecvenţe joase este 83 535/2=41 76 Această valoare mare va duce la probleme cu zgomotele în domeniul frecvenţelor înalte. În acest caz, raportul trebuie să fie mai mic decât 10. Pentru a obţine acest lucru este necesară o reproiectare şi se alege o pulsaţie de tăiere mai mică sau se reduce marginea de fază.

Caracteristica amplificare-pulsaţie (db)

101−

100

101

50 0 -50

Page 66: Indrumar de Laborator Matlab

Fig.6.3. Diagramele Bode ale sistemului

Proiectarea unui regulator cu întârziere de fază se face în acelaşi mod ca mai sus, utilizând opţiunea 3 a funcţiei frdesign

6.6. Proiectarea unui regulator de tip PID

Funcţia de transfer a unui astfel de regulator este, pentru j gcω de forma

gcDgc

pgc jKjKKjG ωω

ω ++= 1)( (6.33)

Parametrii controlerului sunt obţinuţi prin înlocuirea lui (6.33) în (6.24) şi egalarea părţilor reale şi imaginare.

M

K pθcos= (6.34)

gcgc

D MKK

ωθ

ωsin

21 += (6.35)

unde θ , M şi Ψ sunt date de ecuaţiile (6.27) şi (6.28). Pentru un controler stabil, K p şi K D trebuiesc să fie pozitive. Deci gcω trebuie selectat astfel încât

θ în (6.27) să fie mai mic de 90 grade. Funcţia utilizată pentru proiectarea regulatoarelor PI, PD, PID este funcţia frdesign cu opţiunile 4, 5 şi 6. Utilizatorul trebuie să introducă marginea de fază. Funcţia va returna domeniul de pulsaţii pentru care controlerul este stabil, din care utilizatorul va alege o pulsaţie Pentru controlerul de tip PID factorul K1 trebuie specificat. Funcţia va returna toate elementele amintite la proiectarea. regulatorului cu avans de fază.

101−

100

101

Caracteristica fază-pulsaţie (grade) 0 -100 -200 -300

Page 67: Indrumar de Laborator Matlab

Capitolul VII

7. Proiectarea moderna a regulatoarelor Tehnicile clasice de proiectare descrise în capitolele anterioare sunt bazate pe locul rădăcinilor în frecvenţă , care utilizează doar ieşirea părţii fixate cuprinsă intr-o configuraţie tip reacţie cu regulator. În acest capitol vor fi folosite tehnici moderne de proiectare a regulatoarelor, care necesită utilizarea tuturor variabilelor de stare pentru a determina un regulator. Proiectarea modernă a regulatoarelor este utilă in special în cazul sistemelor multivariabile; totuşi , in acest capitol, tehnicile de proiectare în spaţiul stărilor sunt ilustrate utilizând sisteme cu o singură intrare, o singură ieşire .

O metodă de abordare în controlul modern al sistemelor, însoţit de utilizarea reacţiei , este cunoscută ca proiectarea prin alocarea polilor . Această metoda de proiectare permite ca toate rădăcinile ecuaţiei caracteristice a sistemului să fie plasate în locaţiile dorite . Astfel rezultă un regulator cu constanta câştig vectorul k . Conceptul de variabilă de stare în reacţie cere ca toate stările sa fie accesibile într-un sistem fizic , dar la majoritatea sistemelor nu este îndeplinită această cerinţă , de exemplu câteva stări fiind inaccesibile . Pentru sistemele în care toate stările nu sunt “disponibile” pentru reacţie , poate fi proiectat un estimator de stare(observer) . O altă modalitate de proiectare a sistemelor automate este problema controlului optimal , în care este minimizat un criteriu de performanţă specificat matematic .

7.1 Proiectarea prin alocarea polilor

Se consideră sistemul descris de ecuaţiile de stare :

( ) ( ) ( )tButAxtx +=.

( ) ( )tCxty = ( 7.1 ) Se consideră diagrama bloc a sistemului prezentată în figura 7.1 cu următoarea comandă : u(t)=-Kx(t) ( 7.2 )

unde K este un vector de dimensiune l x n format din coeficienţi plasaţi pe buclele de reacţie . Intrarea sistemului , r(t) , este presupusă a fi zero . Obiectivul acestui sistem este de a aduce toate variabilele de stare la valoarea zero când starea a fost perturbată . Înlocuind ( 7.2 ) în ( 7.1 ) , rezultă reprezentarea sistemului in buclă închisă în planul stărilor :

( ) ( ) ( ) ( )txAtxBKAtx f=−=.

Page 68: Indrumar de Laborator Matlab

Ecuaţia caracteristică a sistemului în bucla închisă este : | sI - A + BK | = 0 ( 7.4 )

Fig.7.1 Structura sistemului pentru proiectarea prin alocarea polilor Se presupune că sistemul este reprezentat în forma canonica astfel :

( )tu

xx

xx

aaaax

x

x

x

n

n

nn

n

+

−−−−

=

10

00

.... 1 .... 0 0 0

0 ... . 1 0 0 0 .... 0 1 0

.

.

.

.

1

.

.

.

.

2

1

1210

.

.

.

.

.

.

.

1

.

.

.

2

1

(7.5)

Înlocuind A si B în ecuaţia ( 7.4 ) este determinată ecuaţia caracteristică pentru sistemul în

buclă deschisă : | sI – A + BK | =sn + ( an-1 + kn )sn-1 + K + ( a1 + k2)s + k1 = 0 ( 7.6 )

Pentru poziţii specificate ale polilor in buclă închisa –λ1 . . . . –λn ecuaţia caracteristică dorită este :

αc(s) = (s+λ1) . . . (s+λn) = sn + αn-1 sn-1 + . . . + α1s + α0 = 0 ( 7.7)

Obiectivul proiectării este găsirea matricii K astfel încât ecuaţia caracteristică pentru sistem să fie

identică cu ecuaţia caracteristică dorită . Astfel , vectorul K este obţinut egalând coeficienţii ecuaţiilor (7.6) si

(7.7) :

k1 = αi-1 – ai-1 (7.8)

Daca modelul de stare nu este in forma canonică, vor putea fi utilizate tehnicile de transformare

prezentate în capitolul 2, pentru a transforma modelul de stare în forma canonică . Factorul de câştig este

kn

k2

k1

r(t)=0 u(t) y(t)

xn(t) x2(t) x1(t)

Parte fixată - -

Page 69: Indrumar de Laborator Matlab

obţinut pentru model şi apoi va fi transformat înapoi pentru a corespunde cu modelul original . Acest

procedeu este descris de formula următoare , cunoscută ca formula Ackermann :

k = [ 0 0 . . . 0 1 ] S-1 αc(A) (7.9)

unde matricea S este data de S = [ B AB A2B . . . An-1B] (7.10)

si αc(A) dat de :

αc(A)=An + αn-1An-1 + . . . + α1A + α0I (7.11)

Funcţia [K,Af]=placepol(A,B,C,p) este utilizată pentru proiectarea prin alocarea polilor . A,B,C sunt matricile sistemului si p este vector linie continuând polii doriţi in buclă închisă . Această funcţie returnează

vectorul K si matricea sistemului în buclă închisa Af . De asemenea MATLAB Control System Toolbox

conţine două funcţii pentru proiectarea prin plasarea polilor . Funcţia K = acker(A,B,p) este pentru sisteme

cu o singură intrare , şi K = place(A,B,p) , care utilizează un algoritm mai sigur , este pentru sisteme cu mai

multe intrări .

Condiţia care trebuie să existe pentru a plasa polii în buclă închisă în locurile dorite este de a putea transforma modelul de stare în formă canonică . Aceasta înseamnă că matricea de controlabilitate S , dată în

(7.10) trebuie sa aibă un determinant nenul . Această caracteristică este cunoscută drept controlabilitate .

7.2 Controlabilitate

Se spune că un sistem este controlabil când intrarea părţii fixate u poate fi utilizată pentru a

transfera partea fixată din orice stare iniţială într-o stare arbitrară în timp finit. Partea fixată descrisă de ( 7.1 ) cu matricea sistem având dimensiunea n x n este de stare complet controlabilă dacă şi numai dacă

matricea de controlabilitate S din (7.10) are rangul n . Funcţia S = cntrable(A,B) este dezvoltată pentru a

returna matricea de controlabilitate S si determină dacă sistemul este sau nu controlabil .

7.3 Proiectarea observer-lui

În proiectarea sistemelor de conducere prin metoda alocării polilor , se presupune că toate variabilele de stare sunt disponibile pentru reacţie . Totuşi , în practică este dificilă instalarea tuturor

traductoarelor necesare pentru măsurarea tuturor stărilor . Dacă variabilele de stare nu sunt disponibile

datorită configuraţiei sistemului sau costului , este necesar un observer sau un estimator . Observerul este

Page 70: Indrumar de Laborator Matlab

un algoritm estimativ bazat pe modelul matematic al sistemului. Observerul crează o estimare ( )tx^

a stării

din măsurarea ieşirii y(t) . Starea estimată intră apoi în controler , conform schemei din fig . 7.2

Se consideră un sistem reprezentat de ecuaţiile de stare şi de ieşire :

( ) ( ) ( )tButAxtx +=.

(7.12)

( ) ( )tCxty = (7.13)

Se presupune că starea x(t) este aproximata de starea ( )tx^

a modelului dinamic :

( ) ( ) ( ) ( ) ( )

−++= txtyGtButxAtx

^^^*

(7.14)

( ) ( )txCty^^

= (7.15)

Fig. 7.2 Estimarea stării

Scăzând (7.14) din (7.12) si (7.15) din (7.13) , rezulta :

( ) ( ) ( ) ( ) ( ) ( )

−−

−=− tytyGtxtxAtxtx

^^^.)(

*

(7.16)

( ) ( ) ( ) ( )

−=

− txtxCtyty

^^ (7.17)

unde ( ) ( )txtx^

− este eroarea dintre vectorul de stare actual şi vectorul estimat şi ( ) ( )tyty^

− este eroarea dintre ieşirea actuală si ieşirea estimată . Înlocuind ecuaţia ieşirii în ecuaţia de stare , se obţine ecuaţia erorii dintre vectorul de stare estimat şi vectorul de stare actual .

( ) ( ) ( ) ( )teAteGCAte e=−=.

(7.18) unde :

( ) ( ) ( )txtxte^

−= (7.19) Daca G este ales astfel încât valorile proprii ale matricii A – GC sa aibă toate partea reală

negativă, atunci valoarea staţionară a vectorului eroare a stării estimate e(t) pentru orice condiţii

iniţiale va tinde la zero . Deci, ( )tx^

converge la x(t) . Proiectarea observatorului este asemănătoare cu proiectarea controlerului . Totuşi , valorile proprii pentru A-GC trebuie sa fie selectate la stânga valorilor proprii ale lui A . Aceasta asigură că observerul este mai rapid decât controlerul pentru asigurarea unei rapide actualizări a vectorului de stare estimat .

r(t)=0 u(t) y(t) Parte fixată

Regulator Observer

Stare estimată-

Page 71: Indrumar de Laborator Matlab

Ecuaţia caracteristică a estimatorului este dată de : | sI – A + GC | = 0 ( 7.20)

Pentru o viteză specifică a răspunsului, ecuaţia caracteristică dorită pentru estimator este :

αe(s) = sn + αn-1sn-1

+ K + α 1 s + α 0 = 0 ( 7.21)

Astfel , caştigul estimatorului G este obţinut prin egalarea coeficienţilor din ( 7.20) si (7.21) . Această

metoda este identică cu cea de la alocarea polilor si G este găsit aplicând formula Ackermann :

( )

=

1

001

1M

CACAC

AGn

eα (7.22)

si notaţia αe(A) este data de : αe(A) = An + αn-1 An-1 + K + α1 A + α0I (7.23) Funcţia [ G,Ae] = observer (A , B, C,pe) este dezvoltată pentru determinarea estimatorului ; pe conţine valorile proprii ale estimatorului dorit . Această funcţie returnează vectorul G şi matricea sistemului în buclă închisă Af . Condiţia necesară pentru proiectarea unui observer este ca toate stările să poată fi observate din măsurarea ieşirii. Această caracteristică este cunoscută ca observabilitate. Exemplu. Pentru partea fixată

nxx

x

x

+

=

10

0 161 0

2

1

2

.

1

.

[ ]0 1=y să se proiecteze un observer de stare complet ştiind că observerul este amortizat critic cu valorile proprii la –8 şi –8. Estimator gain vector G 16 80 Open-loop Plant transfer function: Numerator 0 0 1 Deominator 1 0 -16 Error matrix A-G*C -16 1 -64 0

7.4 Observabilitate

Se spune că un sistem este observabil dacă vectorul iniţial x(t) poate fi găsit prin măsurarea lui u(t) si y(t) . Partea fixată descrisă de (7.12) este complet observabilă de stare dacă există inversa matricii în (7.22) .

Funcţia O = obsvable(A,C) returnează matricea de obsevabilitate C şi determină dacă sistemul este sau nu observabil.

Page 72: Indrumar de Laborator Matlab

7.5 Proiectarea combinata controler – observer Se considera sistemul reprezentat de ecuaţiile de stare si ieşire (7.12) si (7.13) cu controlul reacţiei

după stare bazat pe stare observată ( )tx^

data de :

( ) ( )txKtu^

−= (7.24) Înlocuind în (7.12) ecuaţia de stare devine :

( ) ( ) ( ) ( )tBKetxBKAtx +−=.

(7.25) Combinând ecuaţia anterioară cu ecuaţia erorii dată de (7.19), rezultă:

( )( )

( )( )

−=

tetxBKA

te

txGC-A 0

BK .

.

(7.26)

Funcţia [K,G,Ac] = placeobs (A,B,C,p,pe) este utilizată pentru proiectarea combinată controler +

observer . A este matricea sistemului , B este vectorul intrării si C este vector linie al ieşirii , p este vectorul linie conţinând polii în buclă închisă doriţi şi pe conţine valorile proprii ale observerului dorit .Funcţia afişează vectorii K şi G, funcţia de transfer a părţii fixate în buclă deschisă şi funcţia de transfer în buclă închisă a sistemului automat. De asemenea, funcţia returnează vectorul K şi matricea sistem combinată în (7.26).

7.6 Proiectarea unui regulator optimal

Scopul proiectării unui regulator optimal este determinarea comenzii optimale (optimal control)

( )txu ,.

care poate transfera sistemul din starea sa iniţială în starea finală (cu zero intrări sistem) astfel încât să fie minimizat un indice de performanţă dat. Indicele de performanţă care este mult utilizat în proiectarea optimală (optimal control design) este cunoscut ca indice pătratic de performanţă şi este bazat pe criteriul de minimizare a erorii şi de minimizare a energiei.

Se consideră partea fixată descrisă de:

( ) ( ) ( )tButAxtx −=.

(7.27)

Problema este de a găsi vectorul K(t) pentru expresia comenzii: ( ) ( ) ( )txtKtu −= (7.28) care minimizează valoarea indicelui pătratic de performanţă J de forma:

( )∫ +=1

0

''t

t

dtRuUQxJ (7.29)

supus ecuaţiei dinamice a părţii fixate în (7.27). În (7.29) Q este o matrice pozitiv semidefinită şi R este o matrice simetrică reală. Q este pozitiv semidefinită dacă toţi minorii săi principali sunt nenegativi. Alegerea elementelor pentru Q şi R permite orientarea relativă a variabilelor de stare individuale şi a intrărilor separate.

Pentru a obţine o soluţie formală, se poate utiliza metoda multiplicatorilor lui Lagrange. Problema este rezolvată prin adăugarea lui (7.27) la (7.29) utilizând un vector cu dimensiunea n, vectorul multiplicatorilor Lagrange. Problema se reduce la minimizarea următoarei funcţii:

Page 73: Indrumar de Laborator Matlab

( ) [ ]

−+++=

.''',,, xBuAxRuuQxxtux λλα (7.30)

Valorile optimale (notate cu steluţă *) se găsesc egalând derivatele parţiale cu zero.

****** 0 BuAxxxBuAxx

+=⇒=−+=∂∂α

(7.31)

BRuBRuu

'210'2 1** λλα −−=⇒=+=

∂∂

(7.32)

λλλλα '20'2 *..

* AQxAQxx

−−=⇒=++=∂∂

(7.33)

Se presupune că există o matrice p(t) simetrică pozitiv definită, variantă în timp, care satisface: ( ) *2 xtp=λ (7.34) Înlocuind (7.34) în (7.32) se ajunge la expresia comenzii optimale în buclă închisă: ( ) ( ) *1* ' xtpBRtu −−= (7.35) Calculând derivata în (7.34), rezultă:

+=

.**

..2 xpxpλ (7.36)

În final, egalând (7.33), se obţine:

( ) ( ) ( ) ( ) ( )tpBBRtpQtpAAtptp '' 1.

−+−−−= (7.37)

Ecuaţia anterioară este cunoscută ca ecuaţia matricială Riccati. Condiţia limită pentru (7.37) este p(tg) = 0. Deci, (7.37) trebuie să fie integrată.

Funcţia [τ,p,k,t,x] = riccati este dezvoltată pentru soluţia în domeniul timp a ecuaţiei Riccati. Funcţia returnează soluţia ecuaţiei matriciale Riccati p(τ), vectorul câştig al reacţiei optimale k(τ) şi răspunsul stării iniţiale x(t) . Pentru a utiliza această funcţie, utilizatorul trebuie să declare funcţia:

[A,B,Q,R,t0,tf,x0] = system(A,B,Q,R,t0,tf,x0)

conţinând matricile sistem şi matricile indicelui de performanţă într-un fişier de tip m, numit system.m. Pentru sisteme liniare invariante în timp, p = 0 când procesul are o durată infinită, deci ty = ∞, rezultând reducerea la ecuaţia algebrică Riccati: pA + A’p + Q - pBR-1B’p = 0 (7.38)

Funcţia din MATLAB Control System Toolbox [k,p] = lgr2 (A,B,Q,R) poate fi utilizată pentru soluţionarea ecuaţiei algebrice Riccati.

Page 74: Indrumar de Laborator Matlab

Bibliografie

• Clark, R.N., Introduction to Automatic Control System, John Wilez, 1962. • Dodescu, Gh., Simularea sistemelor, Editura Academiei, Bucureşti. • Dumitrache, I., Dumitriu, S., Mihu, I., ş.a., Automatizări electronice, Editura

Didactică şi Pedagogică, Bucureşti, 1993. • Nitu, C., Structuri de reglare analogică, note de curs 1992 – 1993. • Ogata, K., Modern Control Engineering, Prentice Hall, 1991. • Voicu, M., Tehnici de analiză a stabilităţii sistemelor automate, Editura Tehnică,

Bucureşti, 1986.


Recommended