Cruise Control:Metoda Spatiul Starilor pentru proiectarea unui controller

In acesta prezentare vom proiecta si observa un controller pentru sistemul cruise-control

folosind modelul reprezentarii de stare.

Comenzi cheie din MATLAB folosite in aceasta prezentare: ss , feedback

Continut:

Ecuatiile de stare

Cerintele de proiectare

Proiectarea controller-ului folosind pozitionarea cu poli pole placement

Referinta de intrare

Ecuatiile de stare

Ecuatiile de miscare in spatiul starilor sunt urmatoarele:

unde:

(m) masa vehicolului 1000 kg

(b) coeficientul de amortizare 50 N.s/m

(u) forta nominala de control 500 N

(v) viteza vehicolului unde y=v este iesirea din sistem

Cerintele de modelare

Timpul de raspuns < 5s

Abatere < 10%

Eroarea la stabilizare < 2%

Proiectarea controller-ului folosind pozitionarea cu poli pole placement

Schema de raspuns a sistemului este prezentata mai jos:

unde:

K = matricea de raspuns spatiului starilor

u = r Kv = controlul intrarii

Va amintiti din tutorialul Spatiul Starilor ca putem folosi tehnica de pozitionare a polilor

pentru a obtine rezultatul dorit. Polii unui sistem cu bucla inchisa pot fi gasiti din ecuatia

caracteristica: determinantul matricei [sI-(A-B*K). Daca polii sistemului pot fi plasati in

locul dorit prin gasisrea unei matrice potrivite K, atunci poate fi obtinuta iesirea dorita. In

aceasta prezentare, vor fi alesi polii, apoi vom folosi MATLAB pentru a gasi matricea

corespunzatoare K.

Acum, trebuie sa gasim locul de pozitionare al polilor pentru sistemul nostru. Deoarece, la

noi, [sI-(A-B*k)] este de forma 1x1, avem doar un pol de pozitionat. Pozitionam arbitrar

polul la -1,5. La fel ca in tutorialul Spatiul starilor, comanda MATLAB place va fi folosita

pentru a gasi matricea de control K. Creati un nou fisier m-file si introduceti comenzile

urmatoare. Rularea fisierului m-file in fereastra de comanda a programului MATLAB ar

trebui sa ne returneze matricea de control si raspunsul trapta de mai jos:

m = 1000;

b = 50;

t = 0:0.1:10;

u = 500*ones(size(t));

A = [-b/m];

B = [1/m];

C = [1];

D = [0];

sys = ss(A,B,C,D);

x0 = [0];

p1 = -1.5;

K = place(A,B,[p1])

sys_cl = ss(A-B*K,B,C,D);

lsim(sys_cl,u,t,x0);

axis([0 10 0 0.35])

K = 1450

Cum puteti vedea, timpul de raspuns e satisfacator, dar eroarea la stabilizare este prea mare.

Referinta de intrare

Un factor de scalare numit Nbar (schema este aratat mai jos) poate fi folosit pentru a elimina

eraorea la stabilizare. Putem folosi functia rscale pentru a calcula factorul de scalare. Intrarea

este deja multiplicata de 500 de ori, si vrem viteza de stabilizare sa fie de 10 m/s, deci trebuie

sa tinem cont si de acesti factori.

Copiati urmatoarele comenzi intr-un fisier m-file in fereastra de comanda MATLAB. Trebuie

sa obtineti raspunsul treapta aratat mai jos.

Nbar = rscale(sys,K)*10/500;

sys_cl = ss(A-B*K,B*Nbar,C,D);

lsim(sys_cl,u,t,x0);

axis([0 10 0 11])

Cum observati, eraoare la stabilizare a fost eliminata. Timpul de raspuns este mai mic de 5s

iar abaterea este de fapt 0. Toate cerintele de proiectare sunt satisfacute.

Cruise Control: Proiectare unui Controller digital

In versiunea digitala de control a problemei sistemului Cruise Control, vom folosi metoda de

proiectare root-locus pentru controllerul digital.

Comenzi cheie din MATLAB folosite: tf , c2d ,rlocus , zgrid , feedback , step

Continut:

Modelul sistemului

Parametri sistemului

Specificarea performantelor

Functia discreta de transfer

Root-locus in Planul-Z

Compensarea folosind Controller-ul digital

Modelul sistemului

Modelul functiei de transfer pentru problema sistemului Cruise control este dat mai jos.

Pentru derivata, va rugam sa vedeti pagina Cruise Control: Modelarea Sistemului.

Parametri sistemului

In acest exemplu sa presupunem ca parametri sistemului sunt:

(m) masa vehicolului 1000 kg

(b) coeficientul de amortizare 50 N.s/m

(r) viteza de referinta 10 m/s

(u) forta nominala de control 500 N

Specificarea performantelor

Avem nevoie de un controller care sa corespunda urmatoarelor criterii de proiectare:

timp de raspuns < 5s

abaterea < 10%

eroarea la stabilizare < 2%

Functia discreta de transfer

Primul pas in efectuarea unei analize discrete a sistemului este gasirea functiei discrete de

transfer echivalenta pentru partea continua. Vom transforma functia de transfer de mai sus

(Y(s)/U(s)) intr-o functie discreta de transfer folosind functia MATLAB c2d. Pentru a folosi

aceasta functie, trebuie specificate trei argumente: timpul de esantionare Ts, si metoda. Ar

trebui deja sa fiti familiari cu introducerea de matrici num si den. Timpul de esantionare Ts,

in unitatea de sec/esantion, ar trebui sa fie mai mic decat 1/(30*BW), unde BW este

frecventa latimii de banda a buclei inchise. Pentru metoda, vom folosi Zero-order hold.

(zoh).

Timpul de esantionare il lasam 1/50 sec; este suficient de rapid presupunand ca latimea de

banda este 1 rad/sec. Acum introduceti comenzile urmatoare intr-un fisier m-file in fereastra

de comanda din MATLAB.

m = 1000;

b = 50;

u = 500;

s = tf('s');

P_cruise = 1/(m*s+b);

Ts = 1/50;

dP_cruise = c2d(P_cruise,Ts,'zoh')

dP_cruise =

1.999e-05

---------

z - 0.999

Sample time: 0.02 seconds

Discrete-time transfer function.

Root-locus in Planul-Z

Va amintiti din tutorialul Control Digital ca functia MATLAB zgrind poate fi folosita

pentru a gasi o regiune acceptabila a root-locus discreta care na da un raspuns K dorit.

Comanda zgrid necesita doua argumente: frecventa naturala (Wn) si rata de amortizare

(zeta). Aceste doua argumente pot fi aflate din timpul de raspuns si cerintele de abatere si

urmatoarele 2 ecuatii:

unde:

zeta = rata de amortizare

Wn = frecventa naturala (rad/sec)

Tr = timp de raspuns

% OS = abaterea maxima

Deoarece cerintele pentru timpul nostru de raspuns si abatere sunt de 5 sec si 10% , putem

determina ca frecventa naturala (Wn) trebuie sa fie mai amre decat 0.36 rad/sec iar rata de

amortizare (zeta) trebuie sa fie mai mare decat 0.6.

Sa generam root-locus si sa folosim comanda zgrid pentru a gasi regiunile aceptabile ale

root-locus. Dra inainte de asta, argumentul frecventei naturale pentru zgrid trebuie sa fie

in unitati de masura rad/esantion, deci Wn = 0.36*Ts = 0.0072 rad/esantion. Acum

adaugati urmatoarele comenzi si rulati programul. Trebuie sa obtineti:

Wn = 0.0072;

zeta = 0.6;

rlocus(dP_cruise)

zgrid(zeta, Wn)

axis ([-1 1 -1 1])

Regiunea din planul complex care ne intereseaza este cea aproape de punctul (1,0), deci

ar trebui sa mariti aceasta zona. Rulati din nou programul folosind urmatoarele comenzi

pentru axe si figura ar trebui sa arate ca mai jos:

axis ([0.95 1 -.1 .1])

Linia punct din dreapta indica locatiile constantei frecventei naturale (Wn), si frecventa

naturala este mai mare decat 0.0072 in exteriorul liniei. Cealalta linie punct indica locatia

constantei ratei de amortizare (zeta), iar rata de amortizare este mai mare decat 0.6 in

interiorul liniei. Linia verticala franta este o portiune din cercul unitate care este calculat

la rezolutie mica.

In figura de mai sus, puteti observa ca partea din root-locus este in interiorul zonei dorite.

Haideti sa gasim un raspuns specific (K) folosind functia MATLAB rlocfind si apoi sa

obtinem raspunsul treapta corespunzator. Ruland comanda [K,poles] = rlocfind

(dP_cruise) in fereastra de comanda din MATLAB, ne va indruma sa selectam un punct

pe root-locus. Va amintiti ca daca selectati un pol prea in interiorul cercului unitate,

atunci rspunsul treapta va fi prea rapid indicand o acceleratie fizica imposibila. Asadar,

trebuie sa selectati un pol apropiat de intersectia constantei frecventei naturale si axa

reala. Selectati punctul aproape de 0.99 aratat in imaginea de mai jos:

Dupa selectarea unui punct, trebuie sa vedeti urmatoarele rezultate afisate in fereastra de

comanda din MATLAB indicand punctul ales, punctul pe root-locus cel mai apropiat de

pole, si raspunsul K care pozitioneaza polul buclei inchise.

Select a point in the graphics window

selected_point =

0.9900 - 0.0003i

K =

451.1104

poles =

0.9900

Pentru a raspunsul treapta in bucla inchisa, adaugati urmatorul cod in fisierul m-file:

K = 451.1104;

sys_cl = feedback(K*dP_cruise,1);

r = 10;

figure

step(r*sys_cl,10);

Acest raspuns satisface timpul de raspuns si abaterea cerute. Dar eraoare la stabilizare

este de 11%. Pentru a obtine abaterea la stabilizare dorita, vom modifica controllerul

digital.

Compensarea folosind Controller-ul digital

Un compensator cu intarziere a fost adaugat sistemului pentru a obtine raspunsul dorit. In

aceasta versiune de control digital al cruise controlului , vom modifica controllerul digital

existent prin adaugarea unui compensator cu intarziere de forma aratat mai jos:

Acesta este un ghid pentru proiectarea unui compensator digital de conducere si intarziere si

un chid pentru proiectarea unui compensator de timp continuu de conducere si intarziere.

Metoda discreta descrisa spune ca zero pentru compensator trebuie ales aproximativ sa

anuleze unul dintre poli cat timp este stabil. Astfel, vom alege ca zero, zo =0.999.

Pentru a reduce eroarea la stabilizare cu un factor de 5, alegem po = 0.9998. Pentru a avea un

raspuns de 1 la frecventa 0, numaratorul este multiplicat de Kd=(1 zp )/(1 zo ) =0.2 ori

inainte de a folosi root-locus. A se mentiona ca tot compensatorul este multiplicat cu rapunsul

in bucla determinat de root-locus.

Acum avem functia discreta a transfer a compensatorului. Sa generam root-locus si sa

obtinem raspunsul treapta. Prima data creati un nou fisier m-file si introduceti comenzile:

m = 1000;

b = 50;

u = 500;

s = tf('s');

P_cruise = 1/(m*s+b);

Ts = 1/50;

dP_cruise = c2d(P_cruise,Ts,'zoh');

z = tf('z',Ts);

C = 0.2*(z - 0.999)/(z - 0.9998);

Wn = 0.0072;

zeta = 0.6;

rlocus(C*dP_cruise)

zgrid(zeta, Wn)

axis([0.98 1 -0.01 0.01])

Apoi introduceti comenzile [K,poles] = rlocfind(C*dP_cruise) in linia de comanda si

alegeti iar o locatie a polului apropiata de 0.99 cum este aratat in figura de mai jos:

Dupa ce ati ales un punct, trebuie sa vedeti afisat pe ecran in fereastra de comanda

MATLAB indicand punctul ales, punctul de pe root-locus cel mai apropiat de pol si

raspunsul K ce pozitioneaza polul buclei inchise:

Select a point in the graphics window

selected_point =

0.9900 - 0.0000i

K =

2.4454e+03

poles =

0.9900

0.9900

Adaugati codul:

K = 2.4454e+03;

sys_cl = feedback(K*C*dP_cruise,1);

r = 10;

step(r*sys_cl,10);

Acest timp de raspuns la fel de rapid ca cel precedent, are eraoarea la stabilizare redusa la

2%. Sistemul satisface toate cerintele de proiectare in timp ce necesita un efort de control

rezonabil.

Nota: O problema de proiectare nu are neaparat un raspuns unic. Pentru exercitiu, puteti

folosi altfel de compensator pentru a obtine un raspuns mai bun decat cel aratat mai sus.

SIMULINK MODELARE

Continut:

Setare fizica si ecuatiile sistemului

Construirea modelului

Raspunsul in bucla deschisa

Setarea fizica si ecuatiile sistemului

Modelul sistemului Cruise Control este realtiv simplu. Daca se presupune ca rezistenta al

rostogolire si rezistenta aerului sunt proportionale cu viteza de deplasare a autovehicolului,

atunci problema se reduce la sistemul de masa si amortizare de mai jos:

Folosind a doua lege a lui Newton, ecuatia ce guverneaza sistemul devine:

unde u este forta generata intre raota/cauciuc si poate fi controlata direct. Pentru acest

exemplu, presupunem ca:

m = 1000 kg

b = 50 N.sec/m

u = 500 N

Construirea modelului

Acest sistem va fi modelat prin insumarea fortelor ce actioneaza asupra masei si integrarea

acceleratiilor pentru a rezulta viteza. Deschideti Simulink si o noua fereastra de model. Intai,

vom modela integrala acceleratiei:

Inserati un Bloc Integrator ( din libraria Continuous) si trasati lini spre si de al el din

terminalele de iesire si intrare.

Numiti linia de intrare vdot iar pe cea de iesire v. pentru a realiza aceasta operatie

dublu click in locul liber de deasupra liniei.

Deoarece acceleratia (dv/dt) este egala cu suma fortelor impartita de masa, vom imparti

semnalul de intrare cu masa.

Inserati un bloc gain (libraria Math Operations) conectat la blocul integrator la

terminalul de intrare si trasati o linie catre intrarea blocului Gain.

editati blocul Gain prin dublu click pe el si schimbati valoarea la 1/m

schimbati eticheta blocului Gain cu Inertia prin click pe cuvantul Gain de sub bloc.

Acum, vom adauga fortele din ecuatia (1).Prima, va fi forta de amortizare.

Adaugati un Bloc Sum (Math Operations) la linia ce duce al blocul Inertia Gain

Schimbati semnul blocului Sum cu +

Inserati un bloc Gain sub blocul Inertia, selectatil prin click simplu si selectati Flip

din meniul Format pentru o intoarcere de al stanga la dreapta.

Setati valoarea blocului la b si redenumitil damping

trasati (tineti apasat CTRL) o linie de la iesirea blocului integrator si conectati la

intrarea blocului Damping.

trasati o linie de la iesirea blocului Damping la intrarea negativa la blocului Sum

Cea de a doua forta ce actioneaza asupra masei este intrarea de control, u. Vom aplica o

intrarea treapta.

Inserati un Bloc Step (libraria Sources) si conectati-l cu o linie la intrarea pozitiva a

blocului Sum.

Pentru a vedea viteza de iesire, inserati un bloc Scope (libraria Sinks) conectat la

iesirea integratorului.

pentru a furniza un semnal treapta corespunzator de 500 al timp = 0, dublu click pe

blocul Step si steati Step Time la 0 si valoarea finala la u

Raspunsul in bucla deschisa

Pentru a simula sistemul, intai, trebuie setat un timp corespunzator pentru simulare.

Selectati Parameters din meniul Simulation si introduceti 120 in campul Stop Time.

120 de secunde este destul de mult pentru a vizualiza un raspuns in bucla deschisa.

Acum trebuie setati parametri fizici. Rulati urmatoarele comenzi in MATLAB:

m = 1000;

b = 50;

u = 500;

Rulati simularea (CTRL+T sau Start din Simulation Menu). Dupa ce se termina simularea,

dublu click pe blocul Scope si dati clik pe Autoscale. Veti vedea ceva asemanator cu

imaginea de mai jos:

Observand graficul, vrem sa imbunatatim raspunsul sistemului.