Zatvaranje povratne sprege, primer

za početak, da učitamo toolbox

In [1]:
pkg load control

da prebacimo slike na vektorsku grafiku, ovo nije Octave komanda, jupyter je!

In [2]:
%plot inline -f "svg"

ovo je objekt upravljanja, njegova funkcija prenosa

In [3]:
T0 = tf([1],[1, 3, 3, 1])

Transfer function 'T0' from input 'u1' to output ...

                1          
 y1:  ---------------------
      s^3 + 3 s^2 + 3 s + 1

Continuous-time model.
In [4]:
bode(T0)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -120 -100 -80 -60 -40 -20 0 10-2 10-1 100 101 102 Magnitude [dB] Bode Diagram gnuplot_plot_1a -250 -200 -150 -100 -50 0 10-2 10-1 100 101 102 Phase [deg] Frequency [rad/s] T0 T0

funkcija prenosa ima malo pojačanje, nema čak ni crossover frequency (ako izuzmemo DC), stalno je manja od 1 po amplitudi; moramo da podignemo pojačanje i postavimo crossover frequency; zvuči jednostavno, čak i jeste, ali je problem gde postaviti crossover frequency? možda se zbog toga vratimo u drugom prolazu na ovo mesto, ako velik propusni opseg suviše zahteva od sistema (napona, struja, duty ratio van opsega $0 \le d \le 1$, . . . ); postavimo crossover frequency na $\omega_C = 10$ radijana u sekundi, $k = 1000$

In [5]:
T1 = 1000 * T0

Transfer function 'T1' from input 'u1' to output ...

              1000         
 y1:  ---------------------
      s^3 + 3 s^2 + 3 s + 1

Continuous-time model.
In [6]:
bode(T1)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -60 -40 -20 0 20 40 60 10-2 10-1 100 101 102 Magnitude [dB] Bode Diagram gnuplot_plot_1a -250 -200 -150 -100 -50 0 10-2 10-1 100 101 102 Phase [deg] Frequency [rad/s] T1 T1
In [7]:
nyquist(T1)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -500 0 500 -200 0 200 400 600 800 1000 Imaginary Axis Real Axis Nyquist Diagram T1 T1 gnuplot_plot_2a gnuplot_plot_3a

liči na nestabilno, ali čik razaznajte šta je oko (-1, 0) sa ove slike? detaljnije je dole prikazano, zoom je česta komanda kada se radi nyquist()

In [8]:
nyquist(T1)
axis([-1, 1, -1, 1])
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 Imaginary Axis Real Axis Nyquist Diagram T1 T1 gnuplot_plot_2a gnuplot_plot_3a

vrlo nestabilno, ali treba da se tumači; da probamo sa mojom kompresiojom?

In [9]:
compressed(T1)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a

a sada najava nečeg što sledi, da dodamo obilaznice i brojanje obuhvata kritične take:

In [10]:
dodato_nestabilnih_polova = pypass(T1)

dodato_nestabilnih_polova = 2
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a

krajnje očigledno nestabilno, obuhvata kritičnu tačku; da napravimo lead kompenzator da popravimo faznu marginu? ovde je glavno da kompenzator ne pokvari $\omega_C$, tako se pravi funkcija prenosa lead kompenzatora, da na $\omega_C$ ima pojačanje 1, odnosno 0 dB

In [11]:
H1 = feedback(T1)

Transfer function 'H1' from input 'u1' to output ...

                1000          
 y1:  ------------------------
      s^3 + 3 s^2 + 3 s + 1001

Continuous-time model.
In [12]:
pzmap(H1)
legend('location', 'northwest')
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -10 -5 0 5 10 -15 -10 -5 0 5 Imaginary Axis Real Axis Pole-Zero Map H1 H1

tačno tako, dodata dva nestabilna pola, tj. dva pola premeštena u desnu kompleksnu poluravan; igorišite ružnu oznaku u legendi, tu nije pol već oznaka za pol

In [13]:
Hlead = 0.1 * tf([1, 1], [0.01, 1])

Transfer function 'Hlead' from input 'u1' to output ...

      0.1 s + 0.1
 y1:  -----------
      0.01 s + 1 

Continuous-time model.
In [14]:
bode(Hlead)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -20 -10 0 10 20 10-1 100 101 102 103 Magnitude [dB] Bode Diagram gnuplot_plot_1a 20 40 60 80 10-1 100 101 102 103 Phase [deg] Frequency [rad/s] Hlead Hlead

na $\omega_C$ je amplituda lead kompenzatora 0 dB; neće uticati na $\omega_C$

In [15]:
T2 = series(T1, Hlead) # overloaded operator? ranije rešavano množenjem

Transfer function 'T2' from input 'u1' to output ...

                      100 s + 100                
 y1:  -------------------------------------------
      0.01 s^4 + 1.03 s^3 + 3.03 s^2 + 3.01 s + 1

Continuous-time model.

ovde je popravljena faza; uočite kako red sistema raste! nažalost, to će nas pratiti do kraja

In [16]:
bode(T2)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -100 -80 -60 -40 -20 0 20 40 10-2 10-1 100 101 102 103 Magnitude [dB] Bode Diagram gnuplot_plot_1a -250 -200 -150 -100 -50 0 10-2 10-1 100 101 102 103 Phase [deg] Frequency [rad/s] T2 T2

i, šta kaže Nyquist? čik razaznajte odavde!

In [17]:
nyquist(T2)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -60 -40 -20 0 20 40 60 0 20 40 60 80 100 Imaginary Axis Real Axis Nyquist Diagram T2 T2 gnuplot_plot_2a gnuplot_plot_3a

jedan zoom . . .

In [18]:
nyquist(T2)
axis([-5 5 -5 5])
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -4 -2 0 2 4 -4 -2 0 2 4 Imaginary Axis Real Axis Nyquist Diagram T2 T2 gnuplot_plot_2a gnuplot_plot_3a

i još jedan . . .

In [19]:
nyquist(T2)
axis([-1 1 -1 1])
axis square
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 Imaginary Axis Real Axis Nyquist Diagram T2 T2 gnuplot_plot_2a gnuplot_plot_3a

uz zoom se vidi da je jedva stabilno; da probamo kompresiju?

In [20]:
compressed(T2)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a

lako i jasno se vidi da je stabilno, ali i mala margina faze; da probamo step response u zatvorenoj sprezi?

još jedan osvrt na ono što sledi . . .

In [21]:
dodato_nestabilnih_polova = pypass(T2)

dodato_nestabilnih_polova = 0
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a
In [22]:
# nekad minreal(T2 / (1 + T2)), dakle nije overload
H2 = feedback(T2)

Transfer function 'H2' from input 'u1' to output ...

                      100 s + 100                 
 y1:  --------------------------------------------
      0.01 s^4 + 1.03 s^3 + 3.03 s^2 + 103 s + 101

Continuous-time model.

sistem je ostao četvrtog reda, nije postao sistem sedmog reda zbog konačne dužine reči, tzv. numeričke greške! dakle, imamo funkciju feedback() koja lepo radi!

In [23]:
step(H2)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 0 0.5 1 1.5 0 2 4 6 8 10 y1 Time [s] Step Response H2 H2

uh, suvise oscilatorno! nema smisla sa ovim ići dalje, bolje da još popravimo marginu faze, još jedan lead kompenzator dodajemo na red!

odziv je isti! naravoučenije: zaboravite na minreal ako baš ne mora; biće još reči o ovome!

In [24]:
T22 = series(T2, Hlead)

Transfer function 'T22' from input 'u1' to output ...

                          10 s^2 + 20 s + 10                    
 y1:  ----------------------------------------------------------
      0.0001 s^5 + 0.0203 s^4 + 1.06 s^3 + 3.06 s^2 + 3.02 s + 1

Continuous-time model.

sada je sistem petog reda; takav nam je postupak, stalno povećavamo red sistema; rešavamo jedne probleme, a pravimo druge!

In [25]:
bode(T22)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -80 -60 -40 -20 0 20 10-2 10-1 100 101 102 103 Magnitude [dB] Bode Diagram gnuplot_plot_1a -250 -200 -150 -100 -50 0 10-2 10-1 100 101 102 103 Phase [deg] Frequency [rad/s] T22 T22

i tako dođosmo do sistema petog reda; šta kaže Nyquist?

In [26]:
nyquist(T22)
axis([-1 1 -1 1])
axis square
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 Imaginary Axis Real Axis Nyquist Diagram T22 T22 gnuplot_plot_2a gnuplot_plot_3a

ok, stabilno, čak prevelika margina faze, možda bude presporo, ali rešićemo kasnije, bolje je za sada da bude stabilnije; da probamo kompresiju?

In [27]:
compressed(T22)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a

vrlo jasno sta je, margina faze skoro $90^\circ$, nalik na sistem prvog reda; da probamo odziv na step pobudu u zatvorenoj sprezi?

još jedna digresija u budućnost:

In [28]:
dodato_nestabilnih_polova = pypass(T22)

dodato_nestabilnih_polova = 0
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a gnuplot_plot_4a gnuplot_plot_5a

usput, uočite kako je pojačanje na niskim frekvencijama palo, tracking neće biti dobar! rešili smo marginu faze, ali smo spustili pojačanje; u skladu sa najavom: rešavanje jednog problema je napravilo drugi!

In [29]:
H22 = minreal(T22 / (1 + T22))
H22 = feedback(T22)

Transfer function 'H22' from input 'u1' to output ...

                                1e+05 s^5 + 5e+05 s^4 + 1e+06 s^3 + 1e+06 s^2 + 5e+05 s + 1e+05                          
 y1:  -------------------------------------------------------------------------------------------------------------------
      s^8 + 206 s^7 + 1.122e+04 s^6 + 1.63e+05 s^5 + 6.54e+05 s^4 + 1.203e+06 s^3 + 1.151e+06 s^2 + 5.602e+05 s + 1.1e+05

Continuous-time model.

Transfer function 'H22' from input 'u1' to output ...

                           10 s^2 + 20 s + 10                      
 y1:  -------------------------------------------------------------
      0.0001 s^5 + 0.0203 s^4 + 1.06 s^3 + 13.06 s^2 + 23.02 s + 11

Continuous-time model.

NEKAD: peti red sistema posta osmi, optimizacija ili skracivanje (minreal) nije baš sjajno odradila posao SAD: ostao je peti red; drugi metod je zaobišao problem numeričke greške

In [30]:
step(H22)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 0 0.2 0.4 0.6 0.8 0 1 2 3 4 5 6 y1 Time [s] Step Response H22 H22

očekivano, smirene oscilacije, usporen sistem; samo, vidite li da je pojačanje na niskim frekvencijama malo, imate grešku ustaljenog stanja! Da bi to popravili dodajemo lag kompenzator; on fazu kvari, pa ga postavljamo tako da je ne pokvari puno; $\omega_C$ isto ne sme da pokvari!

In [31]:
Hlag = tf([1, 1], [1 0])

Transfer function 'Hlag' from input 'u1' to output ...

      s + 1
 y1:  -----
        s  

Continuous-time model.
In [32]:
bode(Hlag)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 0 10 20 30 40 10-2 10-1 100 101 102 Magnitude [dB] Bode Diagram gnuplot_plot_1a -80 -60 -40 -20 0 10-2 10-1 100 101 102 Phase [deg] Frequency [rad/s] Hlag Hlag

na $\omega_C$ pojačanje je blisko 1, ne kvari se crossover frequency

In [33]:
T3 = series(T22, Hlag)

Transfer function 'T3' from input 'u1' to output ...

                      10 s^3 + 30 s^2 + 30 s + 10                 
 y1:  ------------------------------------------------------------
      0.0001 s^6 + 0.0203 s^5 + 1.06 s^4 + 3.06 s^3 + 3.02 s^2 + s

Continuous-time model.

eto sistema šestog reda!

In [34]:
bode(T3)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -80 -60 -40 -20 0 20 40 60 10-2 10-1 100 101 102 103 Magnitude [dB] Bode Diagram gnuplot_plot_1a -250 -200 -150 -100 10-2 10-1 100 101 102 103 Phase [deg] Frequency [rad/s] T3 T3

šta kaže Nyquist? razaznajte iz slike, ako možete . . .

In [35]:
nyquist(T3)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -300 -200 -100 0 100 200 300 -1 -0.8 -0.6 -0.4 -0.2 0 Imaginary Axis Real Axis Nyquist Diagram T3 T3 gnuplot_plot_2a gnuplot_plot_3a

da probamo kompresiju? bolje, ali ne sjajno, nema obilaznice, ovo je samo komprimovan dijagram koji daje nyquist() funkcija iz paketa control za octave; naša funkcija prenosa u otvorenoj sprezi sada ima pol u koordinatnom početku (nuli, (0, 0))

In [36]:
compressed(T3)
Gnuplot Produced by GNUPLOT 5.2 patchlevel 8 -2 -1 0 1 2 -2 -1 0 1 2 gnuplot_plot_1a gnuplot_plot_2a gnuplot_plot_3a