Matlab
Identificarea sistemelor folosind pachetul de programe matlab. obiecte si proceduri utilizate de subpachetul system identification toolboxIdentificarea sistemelor folosind pachetul de programe Matlab. obiecte Si proceduri utilizate de subpachetul system identification toolbox 1. OBIECTIVELE LUCRARII Familiarizarea cu facilitatile si procedurile de identificare oferite de subpachetul System Identification Toolbox din pachetul de programe Matlab. In etapa acoperita de lucrarea prezenta se exploreaza obiectele principale din subpachet asociate cu sistemele dinamice in reprezentari variate si procedurile ar, arx, armax 2. BREVIAR TEORETIC Procedura de identificare experimentala a unui sistem se realizeaza in contextul unei lipse totale a informatiei privitoare la forma modelului matematic, avand la dispozitie sistemul a carei dinamica trebuie modelata. Sistemul este tratat ca „black box”, pentru care modelul se va sintetiza din informatiile de tip I/E. System Identification Toolbox – SIT – se refera la problema construirii modelelor matematice ale sistemelor dinamice bazate pe date experimentale. SIT este alcatuit dintr-o colectie de tehnici de baza, bine intelese si verificate in aplicatii practice. Identificarea modelelor folosind datele experimentale presupune decizia factorului uman implicat in cautarea modelelor, ca si prelucrarea acestor date in scopul furnizarii bazei acestor decizii. De regula, trebuie sa fie parcurse mai multe bucle, revizuind deciziile anterioare. In acest context, softul interactiv reprezinta modalitatealogica de abordare practica a identificarii sistemelor, fiind totodata si un mijloc elegant de a „compacta” notiunile teoretice extensive, facandu-le accesibile utilizatorului. Matlab este un mediu excelent pentru asemenea calcule interactive, date fiind conceptul sau de workspace, facilitatile grafice si modul de lucru cu datele. SIT este o colectie de fisiere tip .m care implementeaza cele mai uzuale tehnici parametrice/neparametrice disponibile. Este, de asemenea, dedicat pentru utilizator, nu numai din punct de vedere al calculului, ci si al evaluarii modelelor. Pentru automatisti, SIT este de foarte mare utilitate in problema identificarii sistemelor. SIT se afla in directorul ident.off al directorului principal Matlab, fisierul demonstrativ putand fi accesat cu comanda iddemo O lista completa a functiilor pe care subpachetul de identificare le contine se poate vizualiza, dupa o prealabila invocare a platformei Matlab, prin comanda help toolbox/ident care are ca efect afisarea listei respective cu numele functiilor pe categorii, insotite de scurte explicatii in limba engleza. Detaliile asupra functiilor si functionarii acestor programe se obtin prin comanda help urmata de numele programului/functiei type, urmata de numele programului respectiv, sau dbtype urmata de numele programului si de o pereche de intregi separatii de care indica linia de inceput si linia de sfarsit a partii din program afisate pe monitor. 2.1. Obiecte din subpachetul System Identification Toolbox Prima sectiune din lucrare are ca scop familiarizarea cu functiile Matlab care creeaza obiecte asociate cu sisteme dinamice diverse in variatele forme in care acestea se pot prezenta. Se verifica pe rand si apoi in relatia lor logica functiile din subpachetul System Identification Toolbox si din alte subpachete cu care programele de identificare au legatura. Functia tf cu doua argumente se studiaza dupa cum urmeaza num=[1 2];% se initializeaza coeficientii numaratorului functiei de transfer. den=[1 2 5]; % se initializeaza coeficientii numitorului functiei de transfer. sis=tf(num,den) % se genereaza un obiect de tipul tf care contine informatia structurala si parametrica a sistemului dinamic modelat prin functia de transfer afisata. Raspunsul la comanda ultima (fara caracterul „;” la final) este: transfer function s + 2 s^2 + 2s +5 Se poate crea astfel un obiect a carui structura se poate citi prin comanda get(sis),cu rezultatul in mare parte usor de interpretat: num = den= Variable = ’s’Ts = 0 Td = 0InputName =OutputName = Notes =UserData = [] Aceasta este convenabila in lucrul cu sisteme dinamice, cu structuri si parametrizari diverse. De exemplu, prin comanda save nume fisier sis se poate depune in memorie in fisierul nume fisier.mat informatia din sis care poate fi oricand recuperata prin comanda load nume fisier sis. In plus, prin apelurile simple la alte functii, cum sunt impulse (sis) sau step(sis) se pot obtine raspunsurile sistemului la intrarile tipice - impulsul unitar si treapta unitara. Recuperarea „pe bucati” a informatiei din obiectul cu numele sis, de tipul tf se obtine prin comenzile [a,b] = tfdata (sis) sau [a,b] = tfdata(sis,’v’) cu rezultatul depunerii in a si b a numaratorului si numitorului functiei de transfer inclusa in sis, in varianta a doua sub forma vectoriala explicita. Modificarea unei parti a unui obiect tf se poate face prin intermediul comenzilor de genul set (sis,’numeproprietate’,valoare).
Un obiect din categoria tf poate servi ca argument multor altor functii cum sunt bode, freqrest, nyquist, ltiview, nichols etc. si se recomanda incercarea lor si urmarirea efectelor. Aceeasi functie tf poate crea obiecte asemanatoare, dar care se refera la sisteme discrete in timp. Pentru aceasta se include obligatoriu un al treilea argument care reprezinta intervalul de esantionare: sisd = tf (num,den,0.1); Prin invocarea numelui rezultatului sisd se obtine continutul obiectului, adica forma functiei de transfer si intervalul de discretizare. Variabila care apare in scrierea simbolica a functiei de transfer este z. Daca este vorba de modele continue in timp sau discrete in timp, exista posibilitatea ca functia sa fie scrisa in alte variabile cum ar fi p pentru reprezentari continue, z -1 sau q pentru reprezentari discrete. Un alt gen de obiect de un tip diferit, tipul ss, se refera la modelele ecuatie-de-stare-ecuatie-de-observare. Secventa care urmeaza duce la crearea unui obiect de acest tip. a = [-2 1;1 -3]; b = [1 1] ’; c = [1 0;0 1]; d = [0 0]’; sistem=ss (a,b,c,d); Ultima comanda, fara caracterul suprimat al afisarii ”;” expune pe monitor cele patru matrici care intervin in forma modelului dx = Ax(t)+Bu(t) (1.1) dt y(t) = Cx(t)+Du(t) care ia in considerare evolutia starii sistemului in forma continua sau in forma discreta x(t+Δt) = Ax(t)+Bu(t) (1.2) y(t) = Cx(t)+Du(t) daca, din nou, se mai adauga functiei ss inca un argument care este intervalul de esantionare. Se intelege ca in varianta din urma timpul nu poate fi decat multiplu intreg al intervalului Δt. Cu un obiect de acest tip, cum este obiectul sistem, se pot evalua raspunsurile sistemului la intrari variate, cu aceleasi comenzi utilizate si cu obiectul sis, de un alt tip, de tipul tf. Matricile se pot recupera pe rand si separat prin comenzi de genulal = sistem.a; Foarte populara printre automatisti este descrierea unui sistem prin zerourile, polii si amplificarea/amplificarile lui. Pentru aceasta, System Identification Toolbox include o functie zpk destinata crearii de obiecte de un tip diferit de cele descrise sumar mai devreme. Sintaxa acestei functii este sys=zpk(z,p,k) sau sys=zpk(z,p,k,t) din nou dupa cum este vorba de reprezentari continue sau discrete in timp. Daca sistemul este de tipul o intrare-o iesire (SISO- Single Input-Single Output), atunci zerourile si polii se introduc ca vectori (vectorul vid uneori!), iar amplificarea ca un scalar. Daca sistemul este de tip MIMO (Multiple Inputs-multiple Outputs), cu mai multe intrari si mai multe iesiri, atunci argumentele z si p sunt matrici celulare cu cate o celula de vectori z[i, j],p[i, j] pentru fiecare intrare (j)-iesire (i), iar k este o matrice rectangulara care contine amplificarile, de asemenea pentru fiecare canal intrare-iesire. ExempluComanda sist=zpk , ,[-5;1] este pentru crearea obiectului /sistemului sist cu o intrare si doua iesiri [ -5/(s+1) ] [(s-2) (s-3)/s(s+1)] Deoarece in situatii diferite pot fi necesare modele ale sistemelor in forme diferite, exista un numar de functii care permit transformarea modelelor de un anumit tip in modele de un alt tip. Astfel:
Desigur, aici nu se pot da toate detaliile referitoare la functiile din System Identification Toolbox.. Prin comenzile specifice help se pot observa si studia mai in amanunt functiile invocate. Exemple simple tratate cu functiile subpachetului sunt, de asemenea, recomandate. O alta forma in care se stocheaza si se manevreaza structuri si parametri de modele este forma/obiectul theta. O comanda help theta dezvaluie structura acestui obiect care este o matrice. O matrice theta contine informatii de structura, parametrii (estimati) pe structura data si acuratetea lor (de asemenea estimata). Pentru modelele structurilor cu iesiri multiple si cu evidentierea spatiului starilor exista reprezentarea thss, analoga in buna masura cu modul de reprezentare theta. Comanda help thss lamureste proprietatile matricilor thss si diferentele dintre cele doua reprezentari din aceeasi familie. Nu sunt de ignorat elementele de datare calendaristica si orara continute in aceste reprezentari, utile in identificarea in sens strict a unor rezultate stocate in aceste forme. Matricile theta au in vedere structura de model intrare-iesire polinomiala foarte generala: A(q)y(t) = [B(q) / F(q)] u(t-nk) + [C(q) / D(q)] e(t) (1.3) cu A, B, C, D si F polinoame in operatorul de intarziere q, de ordinul na, nb, nc, nd, df, respectiv. Matricile theta sunt create de comenzile poly2th, pem, iv, iv4, arx, armax, oe, bj, ar, ivar. Ele pot fi transformate in alte reprezentari cu comenzile trf, zp, th2poly, th2th, th2ss si th2par. Matricile thss sunt structuri care definesc in general modele liniare in spatiul starilor. Sunt folosite de comenzile/functiile pem, th2arx, th2par,th2ss, idsim, present, thinit, fixpar, unfixpar, th2th, sunt create prin comenzile ms2th, mf2th, arx, iv4, arx2th si pot fi modificate prin comenzile pem, thinit, fixpar si unfixpar. Detalierea continutului matricilor theta se poate obtine cu ajutorul functiei/comenzii present. Apelul la comanda help urmata de un nume din lista clarifica actiunea fiecarei functii din cele mentionate. 2.2. Proceduri din subpachetul System Identification Toolbox Forma generala a modelelor intrare-iesire pentru sistemele cu iesire unica este (1.4) cu ut intrarea numarului i din cele nu intrari, cu y iesirea, cu A, Bt, C, D si Ft polinoame in operatorul de deplasare q –1 (sau z -1). Structura generala este definita prin intarzierile nki si prin ordinele polinoamelor implicate care coincid cu numarul de poli si de zerouri ale modelului dinamic de la intrari la iesire si ale modelului perturbatiilor de la e la y. Cateva cazuri speciale, particulare, ale modelului de mai sus sunt urmatoarele: ARX:A(q)y(t) = B(q) u(t-nk)+e(t) ARMAX: A(q)y(t) = B(q) u(t-nk)+C(q) e(t) OE: y(t) = B(q) u(t-nk)+e(t) (1.5) F(q) BJ: y(t) = B(q) u(t-nk)+C(q) e(t) (Box-Jenkins) F(q) D(q) Procedura arx este utilizata pentru identificarea sistemelor de tip ARX. Calculele se fac pe baza unor date (y,u) observate experimental. Secventa care urmeaza executa un exemplu de identificare pe date „experimentale” generate prin simularea unui sistem dinamic de forma cunoscuta. clear sisc=tf % se creeaza un sistem continuu in timp step (sisc) % se reprezinta grafic raspunsul la treapta unitara sisd=c2d (sisc,0.01)% se converteste sistemul la timp discret hold on % pe acelasi grafica step (sisd) % se reprezinta raspunsul la treapta unitara [y,t]= step (sisd); % se evalueaza raspunsul la treapta unitara u1=sign ( randn (size (t) ) ); % se creeaza o intrare binar- aleatoare [y1,t]=lsim (sisd, u1, t);% se evalueaza raspunsul prin simulare z=[y1 u1]; % se creeaza matricea cu coloanele y1,u1 thd=arx (z, [2 2 1])% se aplica procedura de identificare pe z,pe structura cunoscuta [na,nb,nk]=[2 2 1] e=pe(z,thd); % se calculeaza erorile de reprezentare figure% pe un grafic nou h2=gcf;% cu handler-ul h2 set(h2,’Position’,[150 90 560 420]); % dupa pozitionarea graficului plot(t,e)% se reprezinta erorile de modelare y1=y1+0.0001*randn(size(t)); % se altereaza iesirile observate z1=[y1 u1];% se creeaza noua matrice z sub numele de z1 thd1=arx(z1,[2 2 1])% se aplica din nou procedura de identificare e1=pe(z1,thd1);% se calculeaza erorile de reprezentare title (’Erori de modelare’) ylabel (’Date curate’) xlabel (’Timp (s) ’) subplot(2,1,2) plot(t,e1) % se reprezinta erorile de modelare ylabel (’Date cu erori’) xlabel (’Timp (s) ’) Procedura armax este utilizata pentru identificarea sistemelor de tipul ARMAX. O simpla inlocuire in secventa de mai sus a apelului la functia arx prin apelul la functia armax, inlocuire insotita de modificarea argumentului relativ la structura, [2 2 0 1] in loc de [ 2 2 1], face secventa utilizabila pentru ilustrarea identificarii unui sistem ARMAX. Procedura oe este utilizata pentru identificarea sistemelor modelate prin relatii de tipul OE. Si de data aceasta, o simpla inlocuire in secventa de mai sus a apelului la functia arx prin apelul la functia oe face secventa utilizabila identificarii unui sistem OE. Desi semnificatia componentelor vectorului de structura se modifica, acesta poate ramane neschimbat: [2 2 1]. Pentru modele de tipul BJ (Box-Jenkins) este utilizata procedura bj. Si de aceasta data, in secventa de mai sus a apelului la functia arx prin apelul la functia bj face secventa utilizabila pentru ilustrarea identificarii unui sistem BJ. Vectorul de structura se modifica din [2 2 1] in [2 0 0 2 1 ]. Invocarea comenzii help urmata de numele functiei particulare utilizate lamureste semnificatia vectorului de structura pentru fiecare caz in parte. In secventa Matlab de mai sus apare si functia/comanda pe. Aceasta comanda, in varianta de apelare e = pe (z, th); produce diferentele (erorile) dintre datele experimentale z si modelul th dat in forma theta. O invocare cu o lista mai completa de restituiri [ e, v, w, a, b, c, d, f] = pe (z, th) genereaza nu numai erorile de modelare ci si alte informatii asupra modelului. In a, b, c, d si f sunt depuse polinoamele care apar in forma generala a modelului, iar in ceilalti doi vectori cantitatile w=(b/f) u si v = a [y -w]. O comanda/ functie care poate face tot ce pot face functiile referite mai devreme este functia pem, cu urmatoarele posibilitati de apelare: th=pem (z, thstruc); th=pem (z, thstruc, index); si cu alte posibilitati sondabile prin comanda help pem. In aproape toate cazurile, thstruc este o matrice cu structura modelului si cu o parametrizare, fie ea si provizorie, in genul celor generate de functia poly2th. Se poate insa apela functia pem si cu thstrc in forma vectoriala, simpla [na nb nc nd nf nk]. Forma de apelare fara index face operatia de identificare cu initializarea libera sau inspirata de parametrii din matricea thstruc, daca aceasta este completa. Forma cu argumentul index face aceeasi operatie, dar numai pentru parametrii din linia specificata prin acel index. Se intelege ca, in acest din urma caz, matricea thstruc trebuie sa fie completa. In ambele cazuri, introducerea unui ultim argument de tip sir cu valoarea ’trace’ produce afisarea pe monitor a etapelor calculului de estimare. 3. MODUL DE LUCRU Daca nu este deja creat, se creeaza un director/folder de lucru. Se activeaza platforma Matlab si se introduce comanda Cd (calea spre directorul/folderul de lucru) Se apeleaza functiile care creeaza sisteme cu obiecte de toate tipurile, continue sau discrete in timp. Se memoreaza sistemele si apoi se recupereaza, dupa o comanda intermediara clear ( care elibereaza memoria din spatiul Matlab de lucru curent). Se verifica continutul spatiului de lucru curent cu comenzile who si/sau whos. Din imaginatie sau din lucrarile anterioare se creeaza modele de tipuri variate ale unor sisteme dinamice. Se aplica sistemelor intrarea impulse sau intrarea step, prin comenzile adecvate descrise sumar mai sus, si se urmaresc raspunsurile sistemelor. Se exploreaza conversiile posibile intre diferite forme de modele si intre diferite forme sintetice de stocare a caracteristicilor lor. Se extrag caracteristicile numerice ale raspunsurilor - lista de valori ale iesirilor si vectorul momentelor la care raspunsurile sunt observate. Se creeaza un script dupa modelul dat in sectiunea Breviar teoretic. Scriptul va fi executat repetat cu schimbarile semnalate in acea sectiune. Se exploreaza identificarea in conditii de zgomot variate. Exemplul dat introduce un zgomot cu abaterea medie patratica de 0,0001, coeficientul numeric din linia y1=y1+0.0001*randn (size (t)); Se modifica structura propusa si se repeta operatia de estimare de parametri. Din imaginatie sau din lucrari anterioare se creeaza sisteme dinamice si modele de structuri diverse. Se introduc in script si se face succesiv identificarea lor. Se supun observarii critice in identificare si legatura lor cu structurile alternative propuse si cu acuratetea datelor.
|