Příklady k přednášce
21 - Diskrétní modely
spojitých systémů
Michael Šebek
Automatické řízení 2013
3-5-14
Na úvod: CL stabilita při spojitém a diskrétním řízení
Automatické řízení - Kybernetika a robotika
• Pozor při návrhu emulací: CL stabilita spojitého řízení nezaručuje CL
stabilitu diskrétního řízení! Musíme testovat „diskrétní stabilitu“!
a
• Pro soustavu a P regulátor
je
P( s) =
s+a
, a > 0, C ( s ) = k p
• „spojitý“ CL char. polynom cCL ( s ) = s + a + ak P stabilní ∀k P ∈ ( −1, ∞ )
Root Locus
0.4
C=
( z ) C=
(s) k p
a =1
0.8
0.6
0.2
kp = 0
0.1
0
0.4
Imaginary Axis
Imaginary Axis (seconds -1)
• Zde diskrétní regulátor =
= spojitý regulátor
0.3
Root Locus
1
k p = −1
-0.1
0.2
k P = 20
0
k P = −1
-0.2
-0.4
-0.2
-0.6
-0.3
-0.8
-0.4
-6
-5
-4
-3
-2
-1
0
Real Axis (seconds-1)
1
-1
-1.5
-1
-0.5
0
0.5
1
1.5
Real Axis
• Po připojení k „diskrétní soustavě“ (spojitá + vzorkovač + ZOH) je
1 − e − ah
P( z ) =
z − e − ah
→ cCL ( z ) = z − e − ah + (1 − e − ah )k P
• Což je nestabilní nejen pro
Michael Šebek
1 + e − ah
k P ≤ −1 , ale i pro k P ≥
1 − e − ah
Pr-ARI-21-2012
2
CL stabilita při spojitém a diskrétním řízení
Automatické řízení - Kybernetika a robotika
1 − e − ah
a
P( s ) = , a > 0, C ( s ) =
k p → P( z ) = − ah
s+a
z −e
Nyquist Diagram
Nyquist Diagram
15
1.5
Imaginary Axis
1
10
kp ↑
0.5
kp ↑
5
0
0
-0.5
-5
-1
-1.5
-1
-0.5
k p ∈ ( −1, 0]
0
0.5
1
1.5
2
2.5
kp > 0
3
10
15
-5
0
5
10
15
20
• Porovnej GM !
• Proč? V čem je rozdíl?
Michael Šebek
Pr-ARI-21-2012
3
CL stabilita při spojitém a diskrétním řízení
Automatické řízení - Kybernetika a robotika
• Rozdíl je ve vzorkování + ZOH !
− sh 2
e
• ZOH vnáší, zhruba řečeno, dopravní zpoždění
• Srovnej
a =1
a − sh 2
Pse ZOH ( s ) =
e
s+a
Nyquist Diagram
0.6
0.4
0.2
• Což má GM = 30!
Tedy konečné
• Není to sice těch 20,
ale je to jen zhruba
• Proto je lépe s tím počítat
už při spojitém návrhu
Michael Šebek
System: Pzoh
Gain Margin (dB): 30.1
At frequency (rad/s): 32
Closed loop stable? Yes
0
-0.2
-0.4
-0.6
-0.8
-1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
Real Axis
Pr-ARI-21-2012
4
Odvození metod aproximace
Automatické řízení - Kybernetika a robotika
• Protože obecný regulátor (systém) můžeme realizovat sestavou s integrátory
u ( s ) an s n +  + a1s + a0
C=
(s) =
e( s ) bn s n +  + b1s + b0
u (s)
b0
u (s) 1
C=
(s) =
e( s ) s
− a0
− a1
1
s
y(s)
−an −1
t
u=
(t ) u (0) + ∫ e(τ )dτ
0
• Výstup za jednu vzorkovací periodu je
u (kh + h=
) u (kh) + ∫
kh + h
kh
e ( kh )
u (kh + h)
− u (kh)
kh
e ( kh + h )
kh + h
t
Pr-ARI-21-2012
kh + h
kh
∫
e(t )
Michael Šebek
bn
1
s
1
s
• odvodíme diskrétní aproximaci
pro jeden (každý) integrátor
bn −1
b1
e(τ )dτ
e(τ )dτ= u (kh + h) − u (kh)
Různé populární metody různě
aproximují integrál použitím
hodnot v diskrétních okamžicích
vzorkování
5
Metoda přímé diference
Automatické řízení - Kybernetika a robotika
• Náhradě derivace přímou diferencí
odpovídá v tomto grafu náhrada
červené plochy zeleným obdélníkem
∫
kh + h
kh
e(t )
he(kh)
e ( kh )
e(τ )dτ ≈ he(kh)
) u (kh) + ∫
u (kh + h=
kh + h
kh
e(τ )dτ
e ( kh + h )
h
t
kh + h
kh
= u (kh) + he(kh)
• Použitím z-transformace
zu=
( z ) u ( z ) + he( z )
u( z)
h
=
e( z ) z − 1
1
h
≅
s z −1
s≅
z −1
h
• Metoda se také nazývá Eulerova aproximace
Michael Šebek
Pr-ARI-21-2012
6
Metoda zpětné diference
Automatické řízení - Kybernetika a robotika
• Náhradě derivace zpětnou diferencí
odpovídá v tomto grafu náhrada
červené plochy zeleným obdélníkem
∫
kh + h
kh
e(t )
u (kh + h) = u (kh) + ∫
kh
e ( kh + h )
e ( kh )
e(τ )dτ ≈ he(kh + h)
kh + h
he(kh + k )
e(τ )dτ
h
t
kh + h
kh
= u (kh) + he(kh + h)
• Použitím z-transformace
zu=
( z ) u ( z ) + zhe( z )
u( z)
zh
=
e( z ) z − 1
Michael Šebek
1
zh
≅
s z −1
Pr-ARI-21-2012
s≅
z −1
zh
7
Metoda Tustinova neboli bilineární
Automatické řízení - Kybernetika a robotika
• Bilineární Tustinově transformaci
odpovídá v tomto grafu náhrada
červené plochy zeleným lichoběžníkem
kh + h
h
τ
τ
≈
(
)
e
d
[e(hk ) + e(hk + h)]
∫kh
2
u (kh + h=
) u (kh) + ∫
kh + h
kh
= u (kh) +
he(kh)
e ( kh )
e ( kh + h )
h
kh
t
kh + h
obsah zeleného
lichoběžníku
e(τ )dτ
h
[e(hk ) + e(hk + h)]
2
• Použitím z-transformace
h
h
zu ( z ) =u ( z ) + e( z ) + ze( z )
2
2
u( z) h z + 1
=
e( z ) 2 z − 1
Michael Šebek
e(t )
1 h z +1
≅
s 2 z −1
Pr-ARI-21-2012
h
[e(hk ) + e(hk + h)]
2
2 z − 1 2 1 − z −1
=−1
s≅
h z +1 h 1+ z
8
Vlastnosti aproximací
Automatické řízení - Kybernetika a robotika
Řád
• Všechny tyto transformace zachovávají řád systému a tedy i počet pólů
• Aproximace vyšších řádů se neužívají proto, že by řád zvyšovaly
Aliasing
• Pozor na stroboskopický efekt: regulátor nesprávně reaguje na nesprávně
vzorkovaný (alised) signál: poruchu nebo referenci
• Může pomoci „anti-aliasing filter“: Pak VF signály nepůsobí chybně, jsou
neviditelné
• Ale filtr přidává fázové zpoždění a potenciálně destabilizuje CL
Stabilita OL
• Je-li C(s) stabilní, je stabilní i C(z)? Porovnání dále. Podobně minimální fáze
Stabilita CL:
• I když je při předběžném návrhu spojitý CL systém navržen jako stabilní,
po připojení aproximovaného diskrétního regulátoru stabilní být nemusí
• Musíme vypočítat diskrétní model soustavy se vzorkovačem a ZOH, spojit
ho s diskrétním regulátorem a testovat toto CL spojení na diskrétní stabilitu!
Michael Šebek
Pr-ARI-21-2012
9
Stabilita spojitého regulátoru a jeho aproximace (OL)
Automatické řízení - Kybernetika a robotika
Přímá diference
spojitý
diskrétní
z= e sh ≈ 1 + sh
Re s < 0 → Re s < 1
z −1 z 1
=
−
h
h h
Re z < 1 → Re 1 + sh < 1
s=
0
−1 h
0
1
1
• Stabilní spojitý regulátor s módy (póly) VF nebo slabě tlumenými
má nestabilní diskrétní aproximaci
Michael Šebek
Pr-ARI-21-2012
10
Stabilita spojitého regulátoru a jeho aproximace (OL)
Automatické řízení - Kybernetika a robotika
Zpětná diference
1
0
z=
1 − sh
• Stabilita zachována
• I když má spojitý málo tlumené módy, diskrétní je nemá
1 1
−
h hz
Tustinova metoda
1 + sh 2
2 z −1
=
z =
,s
1 − sh 2
h z +1
s=
1
0
• Nejen že zachová stabilitu
(a minimální fázovost)
• Zobrazení stabilní oblasti je
one-to-one, proto se užívá nejčastěji
Michael Šebek
1
Pr-ARI-21-2012
1h
0
1
11
Tustinova metoda: Numerický příklad
Automatické řízení - Kybernetika a robotika
• Ručně
a
C (s) =
a+s
(1 + z −1 ) ah
a
CTustin ( z ) =
=
2  1 − z −1  ah + 2 (1 − z −1 )
a+ 
−1 
h
1
z

+

−1
2  1− z 
s= 

h  1 + z −1 
>> a=2;h=4;C=a./(a+s)
C =
• V Matlabu – CSTbx:
Michael Šebek
2
----2 + s
>> CTustin=c2d(tf(C),h,'tustin')
Transfer function:
0.8 z + 0.8
----------z + 0.6
Sampling time: 4
Pr-ARI-21-2012
12
Příklady
Automatické řízení - Kybernetika a robotika
Spojitý regulátor s přenosem
C (s) =
a
a+s
• Aproximujeme přímou diferencí
( z)
Cforward
=
• Zpětnou diferencí
( z)
Cbackward
=
• A Tustinovou metodou
a
ah
=
z − 1 z + ah − 1
a+
h
a
ahz
=
z − 1 z (ah + 1) − 1
a+
zh
a
=
CTustin ( z ) =
2 z −1 
a + 

h  z +1 
Michael Šebek
Pr-ARI-21-2012
ah ( z + 1)
( ah + 2 ) z + ah − 2
13
Příklady
Automatické řízení - Kybernetika a robotika
Spojitý regulátor s přenosem
Aproximujeme:
• přímá diference
• zpětná diference
• tustin
35
30
25
Magnitude (dB)
s +1
C (s) =
( 0.1s + 1)( 0.01 + s )
Bode Diagram
40
20
15
10
5
0
-5
-10
90
Phase (deg)
45
>> C=(s+1)/((0.1*s+1)*(0.01*s+1));h=.05;
>> fs=1/h,fN=fs/2,oms=2*pi*fs,omN=2*pi*fN 0
h = 0.0500, fs = 20, fN = 10, oms = 12^, omN = 63
>> S=(z-1)/h;
-45
>> Cpd=(S+1)/((0.1*S+1)*(0.01*S+1)),props(Cpd,h);
10
Cpd = 50(z-0.9500) / (z+4.0000)(z-0.5000)
>> S=(z-1)./h*z;
>> Czd=(S+1)/((0.1*S+1)*(0.01*S+1)),props(Czd,h);
>> Ctu=c2d(tf(C),h,'tustin')
Transfer function:
5.857 z^2 + 0.2857 z - 5.571
---------------------------z^2 - 0.1714 z - 0.2571
>> bode(tf(C),tf(Cpd),tf(Czd),Ctu)
-1
Michael Šebek
Pr-ARI-21-2012
2
1
0
10
10
10
Frequency (rad/s)
ωS 2 =
ωN = 63rad s
Tustinova aproximace: nejlepší,
až na okolí Nyquistovy frekvence!
14
Zkreslení frekvence při Tustinově transformaci
Automatické řízení - Kybernetika a robotika
• Tustinova trans. zobrazuje imaginární osu na jednot. kružnici jedenkrát!
• Všechny „spojité frekvence“ s =jω ∈ j [ 0, ∞ ] zobrazí na „diskrétní
z e jΩh , Ω ∈ [ 0, π .h ] . Např. s = j∞ → z = e jπ = −1
frekvence“=
• Naproti tomu vzorkování zobrazuje frekvence přes ω → Ω , takže
s =j π h → z =e jπ =−1
2 z −1
s
=
=
s j=
ω , z e jΩh
• Dosadíme-li do transformačního vztahu
za
h z +1
dostaneme
2 e jΩh − 1 2 e jΩh 2 − e − jΩh 2 2 j sin(Ωh 2) 2 j
 Ωh 
=
jω =
=
=
tan 

jΩh
jΩh 2
− jΩh 2
h e +1 h e
+e
h cos(Ωh 2)
h
 2 
• Z toho vychází výsledné zkreslení
jω
2
 ωh 
Ω = arctan 

h
 2 
− jω
s
z
e jω
přesně
e jΩ aproximace
e − jΩ
e − jω
Michael Šebek
Pr-ARI-21-2012
15
Zkreslení frekvence při Tustinově transformaci
Automatické řízení - Kybernetika a robotika
2
 ωh 
• Ω = arctan  
h
 2 
 (ω h) 2 
je přibližně Ω ≈ ω 1 −

12


Diskrétní frekvence Ω [ rad s ]
Tustinova
transformace
(−1) n x 2 n +1
arctan x = ∑ 0
2n + 1
3
x
=x − +  ,
3
x ≤ 1, x =
i, −i
∞
Spojitá frekvence ω [ rad s ]
Ω =ω h vzorkování
• Toto zkreslení měřítka frekvencí - frequency warping - je nepříjemné.
• Např. digitální realizace spojitého filtru typu pásmová propusť (zádrž)
nemusí správně zachovat požadované pásmo
Michael Šebek
Pr-ARI-21-2012
16
Prewarping
Automatické řízení - Kybernetika a robotika
• Pro jednu vybranou frekvenci můžeme zkreslení napravit
• Úpravou transformace ne s = α z − 1 , které také Re s < 0 → z < 1
z +1
• a stupeň volnosti α můžeme využít k modifikaci zkreslení:
Prewarping
jω h
)
• Vybereme α , aby pro jednu frekvenci ωpw bylo C ( jω pw ) = Cz (e
• Pak bude C ( s ) = Cz ( z ) přesně=
platit pro s j=
0, (DC) a s jω pw
pw
jω pw
ω pw
e pw − 1
= α jω pwh
= jα → α =
tan (ω pw h 2 )
e
+1
jω h
• Frekvence pro prewarping musí být v rozsahu ω pw ∈ [ 0, π h )
• Pro ω pw = 0 dostáváme klasikou transformaci s α = 2 h
• Vybíráme např. ω pw = ωC (pomůže zachovat PM), nebo frekvenci
kritického vrubu nebo frekvenci kritického oscilačního módu
Michael Šebek
Pr-ARI-21-2012
17
Prewarping
Automatické řízení - Kybernetika a robotika
2
 ωh 
• Klasická Tustinova
Ω = arctan 

h
transformace: zkreslení
 2 
• Tustinova transformace
ω pw
2
ω 
n  , α
=
Ω pw
=
arcta
s prewarpingem: zkreslení
h
tan (ω pw h 2 )
α 
prewarping
Diskrétní frekvence Ω [ rad s ]
Tustinova
transformace
ω pw
Spojitá frekvence ω [ rad s ]
Ω =ω h vzorkování
Michael Šebek
Pr-ARI-21-2012
18
Příklad
Automatické řízení - Kybernetika a robotika
Spojitý regulátor s přenosem
16
14
Magnitude (dB)
Aproximujeme:
• Tustin
• Tustin s prewarpingem
ω pw = 50 rad s
18
12
10
8
ω pw = 50 rad s
6
4
2
0
90
45
Phase (deg)
s +1
C (s) =
( 0.1s + 1)( 0.01 + s )
Bode Diagram
20
0
>> C=(s+1)/((0.1*s+1)*(0.01*s+1));h=.05;
-45
>> Ctustin=c2d(tf(C),h,'tustin')
Transfer function:
-90
10
5.857 z^2 + 0.2857 z - 5.571
---------------------------z^2 - 0.1714 z - 0.2571
Sampling time (seconds): 0.05
>> opt = c2dOptions('Method','tustin','PrewarpFrequency',50);
>> Cpw = c2d(tf(C),h,opt)
Transfer function:
5.675 z^2 + 0.6444 z - 5.031
---------------------------z^2 + 0.4666 z - 0.1777
>> bode(tf(C),Ctu,Cpw)
1
Michael Šebek
Pr-ARI-21-2012
Frequency (rad/s)
ωS 2 =
ωN = 63rad s
19
2
10
Příklad: Prewarping
Automatické řízení - Kybernetika a robotika
• Soustavu s přenosem
G(s) =
1
s ( s + 1)
řídíme diskrétně s h = 0.1 s
1
GZOH ( s ) =
0.05s + 1
• Pro spojitý návrh započteme vliv ZOH
takže během návrhu pracujme s přenosem soustavy
=
Gdes ( s ) G=
( s ) GZOH ( s )
1
s ( s + 1)( 0.05s + 1)
• Při spojitém návrhu (který tu přeskočíme) jsme navrhli spojitý
regulátor
0.593s + 1
C ( s ) = 25.7
0.0102 s + 1
• Kterým jsme dosáhli hodnot
=
=
PM 46.9°
ωC 12.8,
Michael Šebek
Pr-ARI-21-2012
20
Příklad: Prewarping
Automatické řízení - Kybernetika a robotika
• Spojitý regulátor
0.593s + 1
C ( s ) = 25.7
0.0102 s + 1
• Nyní aproximujeme Tustinovou metodou
s
z −1
z − 0.844
→ Ctustin ( z ) = 274.2
20
z +1
z + 0.659
>>G=tf(1./(s*(s+1)));h=.1;C=tf(25.7
*(.593*s+1)./(.0102*s+1));Gdes=G*tf
(1/(.05*s+1));omegaC=12.8;
>> Ctustin=c2d(C,h,'tustin')
Transfer function:
274.5 z - 231.8
--------------z + 0.6611
>> Cpw=c2d(C,h,'prewarp',omegaC)
Transfer function:
244.8 z - 201
------------z + 0.7016
• A také použitím prewarping na frekvenci =
ωC 12.8
= ω pw
s
ω pw
z −1
12.8 z − 1
z −1
z − 0.822
=
= 17.2
→ C pw ( z ) = 244.5
z +1
z + 0.700
 ω pwT  z + 1 tan ( 0.64 ) z + 1
tan 

2


Michael Šebek
Pr-ARI-21-2012
21
Příklad: Bodeho grafy regulátorů
Automatické řízení - Kybernetika a robotika
>> bode(C,Ctustin,Cpw)
Michael Šebek
Pr-ARI-21-2012
ω pw = 12.8 rad s
22
Příklad: Bodeho grafy přenosu OL
Automatické řízení - Kybernetika a robotika
>> Gdt=c2d(G,h,'zoh');
Transfer function:
0.004837 z + 0.004679
---------------------z^2 - 1.905 z + 0.9048
>>
>>
>>
>>
Ldes=Gdes*C
Ltustin=Gdt*Ctustin
Lpw=Gdt*Cpw,
bode(Ldes,Ltustin,Lpw)
Michael Šebek
ω pw = 12.8 rad s
Pr-ARI-21-2012
23
Příklad: Odezvy na skok
Automatické řízení - Kybernetika a robotika
>>
>>
>>
>>
Michael Šebek
ARI-21-2012
Tdes=Ldes/(1+Ldes);
Ttustin=Ltustin/(1+Ltustin);
Tpw=Lpw/(1+Lpw);
step(Tdes,Ttustin,Tpw,1)
24
Metoda MPZ
Automatické řízení - Kybernetika a robotika
Metoda sladěných nul a pólů - MPZ (Matched pole-zero)
sh
• vychází ze vztahu mezi póly/nulami zi = e i
spojitého a vzorkovaného signálu
• použitého zde na impulsní charakteristiku regulátoru
• Navíc, je-li to možné, přidáváme nuly v z −1 = −1 , tj. členy ( z −1 + 1)
do čitatele, což vede k zprůměrování současné a předchozí hodnoty
• Metoda je jednoduchá a praktická, i když ne moc podložená
Postup MPZ
1. vypočti nuly a póly spojitého regulátoru C(s)
sh
2. sestav C(z) tak, aby pro jeho nuly a póly platilo zi = e
3. je-li to možné, přidej do čitatele členy ( z + 1) tak, aby se
stupeň čitatele = stupeň jmenovatele
4. nastav zesílení C(z) pro nulové nebo nízké frekvence stejné jako
bylo zesílení v C(s)
i
Michael Šebek
Pr-ARI-21-2012
25
Metoda MPZ
Automatické řízení - Kybernetika a robotika
MPZ pro
s+a
C (s) = KC
s+b
z − e − ah
CMPZ ( z ) = K D
z − e − bh
zi = e si h
a !
1 − e − ah
C=
(0) K=
CMPZ=
(1) K D
C
b
1 − e − bh
a 1 − e − bh
K D = KC
b 1 − e − ah
MPZ pro
s+a
z − e − ah
( z + 1)( z − e − ah )
C ( s=
→ CMPZ ,1 ( z=
→ CMPZ ( z=
) KC
) KD
) KD
s ( s + b)
( z − 1)( z − e − bh )
( z − 1)( z − e − bh )
a 1 − e − bh
K D = KC
2b 1 − e − ah
Michael Šebek
Pr-ARI-21-2012
26
Příklad: ZOH
Automatické řízení - Kybernetika a robotika
r (s)
C ( z)
ZOH
y(s)
G (s)
r ( z)
{
C ( z)
G( z)
y( z)
}
a  −1 1 1
−1 
1 − e − akh

= 1 − e− at

 =  −
pro spojitý přenos
s s+a
 s( s + a) 
a
z
z
G ( s) =

−
1 − e − akh } =
{
s+a
z − 1 z − e − ah
je diskrétní přenos
1
1
=
−
−1
− ah −1
−
−
1
1
z
e
z
a 

−
1
(1 − z ) 
− ah
−1
G( z) =

=
−
(1
)
e
z
 s( s + a) 
=
(1 − z −1 )(1 − e − ah z −1 )
− ah
−1
(1 − e ) z
= (1 − z −1 )
(1 − z −1 ) (1 − e− ah z −1 )
>> sdf(c2d(tf(1/(s+1)),1,'zoh'))
1−α
ans = 0.6321
=
α = e − ah
,
---------reduced
(z-0.3679)
z −α
Michael Šebek
Pr-ARI-21-2012
27
Příklad: Dopravní zpoždění
Automatické řízení - Kybernetika a robotika
teplá
• spojitý přenos části mixéru je
Te ( s )
s
G (s) =
eτ d s=
F ( s ) eτ d s
=
Tec ( s )
s+a
Tec zanedbatelné ztráty tepla
studená
Te
1.5 najdeme diskrétní přenos
a 1,=
h 1,τ=
• Pro hodnoty =
d
• Protože dopravní zpoždění τd není celočíselným násobkem periody
vzorkování h, rozdělíme ho na
τ d =lh − mh, l ∈ Z , m < 1, m ∈ R
1.5 =
2 − 0.5, h ==
1, l 2, λ =
0.5
• tak dostaneme
− mhs
G (s)
F (s)
− lhs e
=e
s
s
• přičemž člen e
Michael Šebek
− lhs
nakonec přejde na z
Pr-ARI-21-2012
−l
28
Příklad: Dopravní zpoždění
Automatické řízení - Kybernetika a robotika
• po dosazení a rozkladu na parciální zlomky dostaneme
{ }
G (s)
G ( z )= (1 − z ) 
=
s
−1
 e − mhs e − mhs 
−
(1 − z ) z 


s+a
 s
−1
−l
o mh sekund posunutý jednotkový skok a stejně posunutá exponenciála
• protože posunutí jsou o méně než celou periodu (m<1),
nebere se vzorek pro t < 0
• vzorky jsou
• tedy
G( z)
Michael Šebek
1(kh)
z ( z − 1)
e − ah ( k + m )1(kh)
ze − amh ( z − e − ah )
z −1 1  z
z +α
ze − amh 
− amh
=
(
)
−

 1− e
z z l  z − 1 z − e − ah 
z l ( z − e − ah )
Pr-ARI-21-2012
29
Příklad: Dopravní zpoždění
Automatické řízení - Kybernetika a robotika
• pro hodnoty =
1.5
a 1,=
h 1,τ=
d
>> G=tf([1],[1 1],'iodelay',1.5)
Transfer function:
1
exp(-1.5*s) * ----s + 1
>> Gd=c2d(G,1,'zoh')
Transfer function:
0.3935 z + 0.2387
z^(-1) * ----------------z^2 - 0.3679 z
Sampling time: 1
>> Gd/.3935
Transfer function:
z + 0.6065
z^(-1) * -------------z^2 - 0.3679 z
Sampling time: 1
• spojitý přenos nemá nulu
• diskrétní přenos má nulu v
e − amh − e − ah
−α =−
1 − e − amh
α
m
>> m=0.1:0.01:1;alpha=(exp(-m)-exp(-1))./(1-exp(-m)); plot(m,-alpha);
>> syms m; m_sb=solve('(exp(-m)-exp(-1))/(1-exp(-m))=1')
m_sb = -log(1/2*exp(-1)+1/2)
>> vpa(m_sb,3)
ans = .380
Michael Šebek
Pr-ARI-21-2012
30
Příklad
Automatické řízení - Kybernetika a robotika
r (s)
C ( z)
• pro spojitý přenos
FOH
y(s)
G (s) =
{ } {}
r ( z)
G (s)
C ( z)
G( z)
y( z)
1 napřed vypočteme
s2
{
}
t3
h3
h3 z ( z 2 + 4 z + 1)
(kh)3
G (s)
1
1
3
3
 2 =
4 
=→
→  (kh) = 
k }=
{
3!
6
6
6
( z − 1) 4
s
s
6
• Pak je diskrétní přenos
>> Gz=c2d(tf(1/s^2),1,'foh')
Transfer function:
0.1667 z^2 + 0.6667 z + 0.1667
-----------------------------z^2 - 2 z + 1
Sampling time: 1
>> Gzp=sdf(Gz)
Gzp =
0.1667(z+3.7321)(z+0.2679)
-------------------------(z-1)(z-1)
( z − 1) 2 h3 z ( z 2 + 4 z + 1)
G( z) =
hz 6
( z − 1) 4
h 2 ( z 2 + 4 z + 1)
=
6 ( z − 1) 2
Michael Šebek
Pr-ARI-21-2012
31
Příklad
Automatické řízení - Kybernetika a robotika
>> D=5/(s+5);DD=zpk(D); T=1/15;
>> DDtustin=c2d(DD,T,'tustin'),
DDmpz=c2d(DD,T,'matched'), ...
DDzoh=c2d(DD,T,'zoh'),DDfoh=c2d(DD,T,'foh'), ...
Zero/pole/gain:
0.14286 (z+1)
------------Sampling time: 0.066667
(z-0.7143)
Zero/pole/gain:
0.28347
---------Sampling time: 0.066667
(z-0.7165)
Zero/pole/gain:
0.28347
---------Sampling time: 0.066667
(z-0.7165)
Zero/pole/gain:
0.14959 (z+0.8949)
-----------------Sampling time: 0.066667
(z-0.7165)
>> bode(DD,DDtustin,DDmpz,DDzoh,DDfoh)
>> omegas=2*pi/T
omegas = 94.2478
Michael Šebek
Pr-ARI-21-2012
32
Porovnání frekvenčních charakteristik
Automatické řízení - Kybernetika a robotika
MPZ
a
ZOH
5
s+5
h = 1 15 s
ωs ≈ 94 rad
G (s) =
ωs 4
Tustin
a
FOH
• všechny
OK asi do
spojitý
Tustin
a
FOH
ω N 2 = ωs 4
ω N 2 MPZ
a
= ωs 4
spojitý
ZOH
ω N = ωs 2
Michael Šebek
Pr-ARI-21-2012
33
33
Diskretizace stavového modelu: odvození
Automatické řízení - Kybernetika a robotika
Diskrétní stavový model soustavy + tvarovacího členu 0. řádu
u (kh)
u (t )
ZOH
y (t )
=
x Ax + Bu
=
y Cx + Du
y (kh)
Y (s)
u (kh)
x=
Fx k + Guk
k +1
yk Cx k + Duk
=
y (kh)
• Při odvození diskrétního modelu vyjdeme z modelu spojitého
=
x Ax + Bu
=
y Cx + Du
• Je-li tento systém v čase t0 ve stavu x(t0 ), pak v čase t ≥ t0
je ve stavu
=
x(t ) e
A ( t − t0 )
t
x(t0 ) + ∫ e A (t −τ ) Bu (τ )dτ
t0
• kde pro výpočet musíme znát vstup v celém intervalu [t0 , t )
Michael Šebek
ARI-21-2012
34
Diskretizace stavového modelu: odvození
Automatické řízení - Kybernetika a robotika
• Při hledání diskrétního modelu nás zajímá, jaký závisí stav v čase tk +1
na stavu v čase tk za předpokladu ZOH, tj. konstantního vstupu
=
uk u (τ ),τ ∈ [tk , tk +1 ) během celého intervalu mezi vzorkováními
h tk +1 − tk a použijeme předchozí vzoreček, dostaneme
• Označíme-li =
x(tk +1 ) e
=
A ( tk +1 −tk )
x(tk ) + ∫
tk +1
tk
= e x(tk ) + ∫
Ah
tk +1
tk
= e x(tk ) +
Ah
(∫
h
0
e A (tk +1 −τ ) Bu (τ )dτ
e A (tk +1 −τ ) dτ Bu (tk )
)
e Aν dν Bu (tk )
• Protože výstupní rovnici vzorkování nezmění,
dostaneme diskrétní model
x=
(tk +1 ) Fx(tk ) + Gu (tk )
=
y (tk ) Cx(tk ) + Du (tk )
Michael Šebek
=
ν tk +1 − τ
ARI-21-2012
F = e Ah
G=
(∫
h
0
)
e Aν dν B
35
Výpočet maticové exponenciály
Automatické řízení - Kybernetika a robotika
Existuje mnoho metod pro výpočet
maticové exponenciály a jejího integrálu
Ah
, G
=
F e=
(∫
h
0
)
e Aν dν B
• Rozklad exponenciály v Taylorovu řadu
2
2 3
i i +1
h
h
h
A
A
A
V = ∫ e Aν dν = Ih +
+
++
+
0
2!
3!
(i + 1)!
F=
I + AV, G =
VB
h
• Přes Jordanův tvar (vlastní čísla)
A = V diag {λi }V −1
>> expmdemo2
{ }
e Ah = V diag ehλi V −1
• Caylay-Hamiltonův teorém
>> expmdemo3
• Matlabská funkce expm – Padého aproximace
>> expmdemo1
Michael Šebek
ARI-21-2012
36
Příklad
Automatické řízení - Kybernetika a robotika
• Pro spojitý systém 1. řádu
=
x α x + β u
• a periodu vzorkování h je
Ah
F e=
=
eα h
(
)
β αν
G ∫ e d=
=
ν β
e − 1)
(
0
α
h
αν
• Takže diskrétní popis vzorkovaného systému se ZOH je
β αν
x(k =
+ 1) e x(k ) + ( e − 1) u (k )
α
αh
Michael Šebek
ARI-21-2012
37
Příklad
Automatické řízení - Kybernetika a robotika
• Pro dvojitý integrátor
se stavovým popisem
• vzorkovaný s periodou h
• je
0 1 
0 

=
x(t ) 
 x(t ) + 1  u (t )
0
0


 
y (t ) = [1 0] x(t )
1 0   0 h 
1 h 
0
+
+
+
=
F =e Ah =I + Ah + A 2 h 2 2 +  =

 

0 1 
0 1  0 0 


h
h ν 
h2 2
Aν
=
dν 
G ∫ e=
Bdν ∫=

0
0 1 
h
 


(
)
• Takže diskrétní popis vzorkovaného systému (se ZOH) je
h2 2
1 h 
x(k + 1) 
x( k ) + 
=
 u (k )

0 1 
 h 
y (k ) = [1 0] x(k )
Michael Šebek
ARI-21-2012
38
Problém okamžitého výpočtu
Automatické řízení - Kybernetika a robotika
• v předchozích příkladech: stupeň čitatele v z = stupeň jmenovatele v z
• tedy diferenční rovnice regulátoru je
u(k) + členy s k-1,… = ce(k) + členy s k-1, k-2,…
• takový číslicový regulátor musí počítat okamžitě, tedy zpoždění
plynoucí z nenulové doby výpočtu je zanedbáno
• to je prakticky přijatelné jen pokud výpočetní čas < 1/10 h, jinak
musí mít diskrétní regulátor aspoň zpoždění 1 krok
• stupeň čitatele v z < stupeň jmenovatele v z
(dostaneme ho třeba tak, že v
MPZ metodě stupeň nedorovnáme)
• anebo to zpoždění musíme „přidat“ do soustavy
• nejde o dopravní zpoždění, je to
jen způsob indexování
Michael Šebek
Pr-ARI-21-2012
39
Frekvenční model ZOH: Značení
Automatické řízení - Kybernetika a robotika
• Jednotkové pulzy
1, t ≤ 1 2
rect(t ) = 
0, t > 1 2
t − h 2 1, t ∈ [ 0, h ]
rect(
)=
h
0, t ∉ [ 0, h ]
1
1
t
t
−1 2 1 2
• Funkce sinc
sinc( x) =
0
h
sin x
x
• Normalizovaná funkce sinc
sinc( x) =
sin(π x)
πx
• používá se ve zpracování signálů
a má integrál (přes všechny x) = 1
Michael Šebek
ARI-21-2012
40
Frekvenční model ZOH?
Automatické řízení - Kybernetika a robotika
∞
• ZOH rekonstruuje spojitý signál
t − h 2 − kh
u (t ) = ∑ u (k ) rect(
)
z diskrétního podle vztahu
h
k =0
• Můžeme ho také modelovat LTI filtrem s obdélníkovou impulzní
odezvou, ne jehož vstup přichází série Diracových pulzů
∞
t − kh
=
us (t ) ∑ u (k=
)δ (
) h∑ u (k ) δ (t − kh)
h
=
k 0=
k 0
∞
• Škálování h vzniká z časového škálování Diraců. Díky jemu je střední
hodnota us (t ) rovna středná hodnotě u (k ) a výsledný filtr má DC zesílení
= 1. Někteří autoři ho vypouštějí, a pak má výsledný filtr DC zesílení = h.
• ZOH tedy můžeme považovat za
ZOH
hypotetický LTI filtr (systém), který mění
výše uvedenou vstupní posloupnost na výše uvedenou
výstupní funkci.
1
t 1 1 h , t ∈ [ 0, T ]
• Jeho impulzní odezva je obdélník g =
(
t
)
rect(
) 
−=
ZOH
h
h 2  0, t ∉ [ 0, T ]
délky h ,výšky 1/h, s plochou =1
Michael Šebek
ARI-21-2012
41
Frekvenční model ZOH?
Automatické řízení - Kybernetika a robotika
• Použitím Fourierovy transformace vypočteme jeho frekvenční přenos
jπ fh
1 − e − j 2π fh
− e − jπ fh
− jπ fh e
− jπ fh j 2sin(π fh)
− jπ fh
sinc( fh)
G
f
e
e
e
=
=
=
=
(
)
ZOH
j 2π fh
j 2π fh
j 2π fh
1 − e − jω h
GZOH
=
( jω ) = e− jω h 2sinc(ω h 2π )
jω h
s j=
ω j 2π f dostaneme přenos v Laplaceově transformaci
• Dosazením=
1 − e − sh
GZOH ( jω ) =
sh
u (t )
h
u (k )
u (t )
ZOH
• Tento přenos má DC zesílení = 1.
− sh
1
−
e
• Někteří požívají jako přenos ZOH GZOH ( jω ) =
s
Michael Šebek
ARI-21-2012
s DC zesílením = h
42
Ideální rekonstrukce
Automatické řízení - Kybernetika a robotika
• Inverzní proces ke vzorkování: interpolační vzorec, Whittaker–Shannon
u (t )
h
• Ideální dolní propusť: Její
• impulsní odezva je
u (k )
sin(π x)
sinc( x) =
πx
t
sinc( )
h
u (t )
a Bodeho graf je
1, ω < π h
G ( jω ) = 
0, ω ≥ π h
• Nekauzální, nelze implementovat - Lze aproximovat, ale lepší
aproximace mají větší zpoždění, což je špatné pro CL stabilitu
Michael Šebek
ARI-21-2012
43
Jak dobře rekonstruuje ZOH?
Automatické řízení - Kybernetika a robotika
• ZOH je nejjednodušší a nejlevnější rekonstrukční filtr. Jak je dobrý?
− jω h 2
1 − e − sh
1 − e − jω h
2 j sin (ω h 2)
− e − jω h 2
− jω h 2 e
GZOH ( s ) =
e
e − jω h 2
→ GZOH ( jω ) =
=
=
hs
jhω
jhω
jω h
− jω h 2 sin (ω h 2)
e=
e − jω h 2sinc(ω h 2)
ωh 2
• Je tedy něco jako dolní propusť
s dopravním zpožděním h/2
• Moc se nepodobá ideální
dolní propusti a má
dost velké fázové zpoždění
• Zejména poblíž ωN = π h
Michael Šebek
ARI-21-2012
π
h
2π
h
44
Download

Příklady k přednášce 21 - Diskrétní modely spojitých systémů