Příklady k přednášce 25 – Dopravní zpoždění
Michael Šebek Automatické řízení 2013 21-4-13
Dopravní zpoždění v Laplaceově transformaci Automatické řízení - Kybernetika a robotika
f (t ) : f (t )= 0 ∀t < 0, L { f (t )}= f (t − τ ) : τ > 0, L { f (t − τ )} = ? t − τ )} L { f (=
− st e ∫ − f (t )dt= f (s) 0
∞
− st e t −τ v = ∫ − f (t − τ )dt ← 0
=∫
∞
=e
− sτ
=e
−s
−τ
− s ( v +τ ) e f (v)dv − ∞
∫τ τ ∫ e −
∞
− sv e f (v)dt − − sv
0−
= e − sτ f ( s )
Michael Šebek
∞
f (v)dt − sτ e L {δ (t= − τ )} e − sτ , L {1(t= − τ )} , s ARI-25-2012
2
Příklad: Válcovací stolice Automatické řízení - Kybernetika a robotika
• dopravní zpoždění na výstupu • přenos
H ( s, e −τ s ) = G ( s )e −τ s • tedy obsahuje dynamiku bez zpoždění + zpoždění
Michael Šebek
ARI-Pr-25-2012
3
Příklad: Pásový dopravník Automatické řízení - Kybernetika a robotika
Těžba fosfátu v Bou Craa, Západní Sahara: 100 km systém dopravníků
Bangladéš: 1 pás 17 km Michael Šebek
ARI-Pr-25-2012
4
Příklad: Regulace v buňce Automatické řízení - Kybernetika a robotika
x1 (t ) = −λ1 x1 (t ) + c1 x2 (t − τ 1 )
x2 (t ) = −λ2 x2 (t ) + g ( x1 (t − τ 1 ) )
Michael Šebek
ARI-Pr-25-2012
5
Příklad: Hematologie Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
6
Příklad: Mechanismus aktivace enzymu Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
7
Příklady: Operační výzkum Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
8
Příklad: Tepelný systém Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
9
Příklad: Síťové řídicí systémy Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
10
Příklad: Router Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
11
Příklad: Automatické řízení - Kybernetika a robotika
Michael Šebek
ARI-Pr-25-2012
12
Čisté zpoždění v Matlabu Automatické řízení - Kybernetika a robotika
D( s) = e− s
>> D = tf(1,1,'InputDelay',1) Transfer function: exp(-1*s) * (1) >> s = tf('s'); D = exp(-s) Transfer function: exp(-s) * (1)
Step Response 1.5
Amplitude
1
0.5
>> step(tf(1),D) >> bode(D) >> nyquist(D)
0
1.5
1
0.5
0
2
2.5
3
5
4.5
4
3.5
Time (seconds)
Nyquist Diagram 1 0.8
Bode Diagram 1
0.6 0.4 0
Imaginary Axis
Magnitude (dB)
0.5
-0.5
-1 0
0.2 0 -0.2
-90
-0.4 Phase (deg)
-180 -270
-0.6
-360 -450
-0.8
-540
-1
-630 0 10
1
10
-1
Frequency (rad/s)
Michael Šebek
ARI-Pr-25-2012
-0.8
-0.6
-0.4
-0.2
0 Real Axis
0.2
0.4
0.6
0.8
1
13
Systém se zpožděním na vstupu/výstupu v Matlabu Automatické řízení - Kybernetika a robotika
0.12 Amplitude
>> G = tf(1,[1 10],'InputDelay',2.1) Transfer function: 1 exp(-2.1*s) * -----s + 10 >> G = tf(1,[1 10],'OutputDelay',2.1); >> s = tf('s');GG= 1/(10+s); G = exp(-2.1*s)*GG;
1 G ( s) = e −2.1s 10 + s 0.1 0.08 0.06 0.04 0.02 0
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Time (seconds)
Nyquist Diagram 0.1
Bode Diagram -20
0.08 0.06
-21 -21.5
0.04
-22 Imaginary Axis
Magnitude (dB)
-20.5
-22.5 -23 -23.5 0
0 -0.02 -0.04
-180
Phase (deg)
0.02
-360
-0.06
-540
-0.08
-720
-0.1 -0.2
-900 -1
10
Michael Šebek
0
10
Frequency (rad/s)
ARI-Pr-25-2012 1
10
-0.15
-0.1
-0.05
0 Real Axis
0.05
0.1
0.15
0.2
14
Systém s vnitřním zpožděním v Matlabu Automatické řízení - Kybernetika a robotika Step Response 0.12
1 e −2.1s T ( s ) = 10 + s 1 e −2.1s 1+ 10 + s e −2.1s = 10 + s + e −2.1s
0.1
0.08
p tude
>> s = tf('s');GG= 1/(10+s); >> G = exp(-2.1*s)*GG; >> T=G/(1+G) ... Output delays (seconds): 2.1 Internal delays(seconds):2.1 Continuous-time model. >> step(T)
0.06
0.04
0.02
0
0
1
2
3
4
5
6
7
Time (seconds)
Nyquist Diagram
Bode Diagram -19 -19.5
0.15
-20.5 0.1
-21 -21.5 Imaginary Axis
Magnitude (dB)
-20
-22 -22.5 -23
0.05
0
-23.5 0 -180
-0.05
Phase (deg)
-360 -540 -0.1
-720
-0.1
-900
-1260 -1 10
-0.08
-0.06
-0.04
-0.02
0
0.02
0.04
0.06
0.08
Real Axis
-1080 0
10
1
10
Frequency (rad/s)
Michael Šebek
ARI-Pr-25-2012
15
0.1
Příklad: Složitější systémy Automatické řízení - Kybernetika a robotika
H ( s, e − s ) =
>> delay=tf(1);set(delay,'ioDelay',1) >> E=tf(delay),S=tf(s); Transfer function: exp(-1*s) * (1) >> H=E/(S+a+S*b*E),step(H,20)
−s
e s + 1 + se − s
Step Response 1.4
Bode Diagram
1.2
50 40 30
1
Magnitude (dB)
20
Amplitude
0.8
Nyquist Diagram 100 90
-30 -40
80 0.4
0 -10 -20
0 dB 0.6
10
-50 0
70
-720
0
0
2
4
6
8
10 Time (seconds)
12
-1440
Phase (deg)
0.2
Imaginary Axis
60 50 40 30
14
16
20
18
-2160 -2880 -3600 -4320 -5040
20
-5760
10
-6480 -2 10
-1
10
0
10
1
2
10
10
2 dB-2 dB 0 -10 -50
-40
-30
-20
-10
0
10
20
30
40
50
Real Axis
Michael Šebek
ARI-Pr-25-2012
16
Příklad: Složitější systémy Automatické řízení - Kybernetika a robotika
H ( s, e
−0.2 s
,e
−0.3 s
)=
e
>> delay2=tf(1);set(delay2,'ioDelay',.2) >> delay3=tf(1);set(delay3,'ioDelay',.3) >> H=E2/(S+E3+2*E2);step(H,20) >> E2=tf(delay2),E3=delay3,S=tf(s); Transfer function: exp(-0.2*s) * (1) Transfer function: exp(-0.3*s) * (1) >> H=E2*S/(S+E3+2*E2);step(H,2)
−0.2 s
s + e −0.3 s + 2e −0.2 s
Nyquist Diagram 2
1.5
1
Bode Diagram 10
Imaginary Axis
1.2
1
Magnitude (dB)
Step Response 0.5
0
-0.5
0.8
Amplitude
-30 x 104 4.608
-1.5
-2 -1.5
0.4
-10
-20
-1
0.6
0
-1
-0.5
0
0.5
1
1.5
2
0
Real Axis
0.2
-4.608 -9.216
0
-0.2
-13.824 -1 10 0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
1.8
2
0
10
1
10
2
10
3
4
10
10
Frequency (rad/s)
Time (seconds)
Michael Šebek
ARI-Pr-25-2012
17
Pro zajímavost: Lambertova funkce Automatické řízení - Kybernetika a robotika
Lambertova W-funkce (také omega funkce) • je inverzní funkce k f (W ) = WeW • v reálném oboru rozumná, ale v komplexním divoká (nekonečně mnoho větví)
reálná
imaginární část Lambertovy funkce (analytického prodloužení) • http://mathworld.wolfram.com/LambertW-Function.html Michael Šebek
a
ARI-Pr-25-2012
18
Příklad: Retardovaný a neutrální systém Automatické řízení - Kybernetika a robotika
• Retardovaný systém
>> solve('exp(-tau*x)+x=0') ans = lambertw(0, -tau)/tau >> tau=1, r=lambertw(-10:10,-tau)./tau; >> plot(real(r),imag(r),'+r')
cCL ( s )= s + e − s 80 60 40 20 0
• Neutrální systém −s
(s) e s + 1 cCL =
-20 -40 -60 -80 -5
-4
-3
-2
-1
0
1
2
3
4
5
>> solve('1+exp(-tau*x)*x=0') ans = -lambertw(0, tau)/tau >> tau=1, r=-lambertw(-10:10,tau)./tau; >> plot(real(r),imag(r),'+b') Michael Šebek
ARI-Pr-25-2012
19
Podmínky stability jednoduchého kvazipolynomu Automatické řízení - Kybernetika a robotika
(Kharitonov et al., 2004, p. 40) • Systém s charakteristickým kvazipolynomem a ( s, e −τ s ) = s + a + be −τ s
kde a + b > 0 (pokud ne, pak není stabilní ano bez zpoždění) je • Stabilní nezávisle na velikosti zpoždění („iod“) právě když a ≥ b . Jinak b • Pokud je a > 0, b > 0 , je stabilní π − arccos ( a b ) ∀τ < 2 2 3
2
1
b −a • Pokud je a < 0, b > 0 , je stabilní
∀τ < Michael Šebek
arccos ( a b ) b −a 2
0
-1
-2
2
-3 -3
ARI-Pr-25-2012
-2
-1
0
1
a
2
3
20
Příklad: Destabilizující efekt zpoždění Automatické řízení - Kybernetika a robotika
1 τ 0= : s1 − K = G (s) = , C ( s ) Ke −τ s= s τ = 0+ : max {Re si } <<< 0 cCL ( s )= s + Ke −τ s τ =↑: max {Re si } ↑ 0
4
8
x 10
τ = 0, τ = 0++ , τ = 0+
6
4
2
0
-2
-4
>> solve('x+k*exp(-tau*x)=0') ans = lambertw(0, -k*tau)/tau >> k=1,tau=0.5 >> r1=lambertw(-10:10,-k*tau)/tau; >> plot(real(r),imag(r),'+r')
-6 -12000
-10000
-8000
-6000
-4000
-2000
0
40
τ c = 1.5708
30
0.5 max Re {si }
20
0 10
-0.5 0
-1 -10
-1.5
τ =1
-20
-2 -2.5
τ = 0.5
0
τc
2
Michael Šebek
4
6
8
10
τ
-30
12
ARI-Pr-25-2012
-40 -3.5
-3
-2.5
-2
-1.5
-1
-0.5
0
21
0.5
Příklad: Stabilizující efekt zpoždění Automatické řízení - Kybernetika a robotika
cCL ( s ) = s 2 + 9 + 1.5e −τ s
• je nestabilní pro τ = 0, • ale stabilní pro malá nenulová zpoždění
• Srovnej jeho odezvu s použitím PD regulátoru Michael Šebek
ARI-Pr-25-2012
22
Příklad: Zpoždění jako derivační ZV Automatické řízení - Kybernetika a robotika
u (t ) je nestabilní. • Soustava s rovnicí y (t ) − 0.1y (t ) + y (t ) =
• Můžeme ji stabilizovat derivační ZV s zesílením k > 0.1 u (t ) = ky (t )
• Alternativně ji můžeme stabilizovat „zpožděnou ZV“ u (t ) = y (t − τ ) − y (t )
Step Response
k =1
• Což můžeme interpretovat jako ZV s konečnou diferencí y (t ) − y (t − τ ) u (t ) = −τ τ
τ = 0.2 τ = 0.5
• Aproximujícím derivaci s
k =τ
Michael Šebek
1
ARI-Pr-25-2012
2
3
4
5
6
7
8
10
9
23
Padého Aproximace Automatické řízení - Kybernetika a robotika
• Henri Eugene Padé - francouzský matematik (1863-1953) • dnes znám hlavně jako autor aproximace obecné funkce pomocí racionální funkce, • která je často lepší než Taylorova Padého aproximace • pro danou funkci f a přirozená čísla m,n je Padého aproximant řádu (m,n) p0 + p1 x + p2 x 2 + + pm x m R( x) = 1 + q1 x + q2 x 2 + + qn x n
kde
f (0) = R (0) f ′(0) = R′(0) f ′′(0) = R′′(0) f ( m + n ) (0) = R ( m + n ) (0)
• součet prvních m+n+1 členů Taylorových řad f a R je stejný Michael Šebek
ARI-Pr-25-2012
24
Příklad: Padého aproximace Automatické řízení - Kybernetika a robotika
>> del=tf(1); >> set(del,'ioDelay',5); del exp(st) Transfer function: exp(-5*s) * 1 >> pade1=pade(del,1) Transfer function: (1,1) -s + 0.4 -------s + 0.4 >> pade2=pade(del,2) Transfer function: (2,2) s^2 - 1.2 s + 0.48 -----------------s^2 + 1.2 s + 0.48 (3,3) >> pade3=pade(del,3) Transfer function: -s^3 + 2.4 s^2 - 2.4 s + 0.96 ----------------------------s^3 + 2.4 s^2 + 2.4 s + 0.96 >> step(del,pade1,pade2,pade3) >> bode(del,pade1,pade2,pade3)
Michael Šebek
ARI-Pr-25-2012
25
Příklad na „přesný“ návrh: P regulátor Automatické řízení - Kybernetika a robotika
•= pro G (s)
K −τ s e , C (s) K p = Ts + 1
je
• a pro hodnoty K= T= τ= 1, K P= 2 • CL charakteristický kvazipolynom je
T (s) =
je
KK p e −τ s (Ts + 1) + KK p e −τ s
2e − s T (s) = s + 1 + 2e − s
cCL ( s ) = s + 1 + 2e − s
• Pro K P = 5
je
cCL ( s ) = s + 1 + 5e − s
>> solve('x+1+5*exp(-x)=0') ans = -1+lambertw(-5*exp(1)) >> r=-1+lambertw(-10:10,-5*exp(1)); >> plot(real(r),imag(r),'*') >> solve('x+1+2*exp(-x)=0') ans = -1+lambertw(-2*exp(1)) >> r=-1+lambertw(-10:10,-2*exp(1)); >> plot(real(r),imag(r),'*') Michael Šebek
ARI-Pr-25-2012
26
Příklad: CL stabilita Automatické řízení - Kybernetika a robotika
b −τ s e s+a C (s) = k G ( s )C ( s ) T (s) = 1 + G ( s )C ( s ) kb −τ s e = s+a kb −τ s e 1+ s+a kbe −τ s = s + a + kbe −τ s
G (s) =
Michael Šebek
ARI-25-2012
27
Příklad na přesný návrh s I regulátorem Automatické řízení - Kybernetika a robotika
• Pro
K 1 = G (s) = , C (s) Ts + 1 TI s
• je
Ke −τ s T (s) = TI s (Ts + 1) + Ke −τ s
• a pro hodnoty K= TI= T= τ= 1
e− s T (s) = s ( s + 1) + e − s
1.0e+002 * -0.09225032465054+1.00357526290690i -0.09096327458520+0.94065517770892i -0.08958771536361+0.87772450636321i -0.08811056322665+0.81478112589439i -0.08651559672926+0.75182228804814i -0.08478236381840+0.68884436592290i -0.08288456611943+0.62584246447846i -0.08078758735029+0.56280979999173i -0.07844455447831+0.49973666852271i -0.07578973834359+0.43660864096520i -0.07272678028205+0.37340319845989i -0.06910590713751+0.31008293641609i -0.06467468145853+0.24658031170133i -0.05895295773482+0.18275807028892i -0.05081944749196+0.11828269472765i -0.03610143271894+0.05213342976426i -0.00037250765679+0.00819937714011i
• CL charakteristický kvazipolynom
cCL ( s ) = s 2 + s + e − s
• má nekonečně mnoho kořenů část kořenů nad reálnou osou Michael Šebek
ARI-Pr-25-2012
28
Příklad: Smithův prediktor Automatické řízení - Kybernetika a robotika
• P regulátor se zesílením 2
2 −s T (s) = e s + 1 + 2e − s
2 −s T (s) = e s+3 Michael Šebek
ARI-Pr-25-2012
29
Příklad: Smithův prediktor Automatické řízení - Kybernetika a robotika
• P regulátor se zesílením 5
5 −s T (s) = e s + 1 + 5e− s
5
5
T (s) =
5 −s e s+6
• po odpojení větve bez prediktoru Michael Šebek
ARI-Pr-25-2012
30
Příklad: Smithův prediktor Automatické řízení - Kybernetika a robotika
• I regulátor a Smithův prediktor 1 s
1 s
Michael Šebek
ARI-Pr-25-2012
31
Příklad: Nestabilní soustava Automatické řízení - Kybernetika a robotika
• Smithův prediktor nefunguje pro nestabilní soustavu! Proč? −0.5 s b s e ( ) • Ale např. G= přesto můžeme stabilizovat (s) = a(s)
s −1
• Pro toto konkrétní zpoždění, „iod“ to nejde • P-regulátor se zesílením k = 1.5 dává stabilní cCL ( s ) = s − 1 + 1.5e −0.5 s 1.5e −0.5 s T (s) = s − 1 + 1.5e −0.5 s
>> Gss=ss(1/(s-1)); >> set(Gss,'ioDelay',.5) >> Tss=feedback(1.5*Gss,1) >> step(Gss,Tss,9) Michael Šebek
>> solve('x + 3/(2*exp(x/2)) - 1 = 0') ans = 2*lambertw(0, -3/(4*exp(1/2))) + 1 >> pol=2*lambertw(-10:10,-3/(4*exp(1/2)))+1; >> plot(real(pol),imag(pol),'+r'),grid ARI-Pr-25-2012
32
Příklad: Přiřazení konečného počtu pólů Automatické řízení - Kybernetika a robotika
• Nestabilní soustava
>> solve('1+exp(-x)=0') ans = pi*i >> zer=(pi.*(-10:10)); >> solve('x-exp(-x)=0') ans = lambertw(0, 1) >> pol=lambertw(-10:10,1); >> plot(real(zer),imag(zer),… 'ob',real(pol),imag(pol),'+r')
b( s ) 1 + e − s G= (s) = a(s) s − e− s
• Char. kvazipolynom
( s − e ) p(s) + (1 + e ) q(s) =s + 1 −s
−s
• Regulátor a CL q( s) 1 = → p( s) 1
1 + e− s T (s) = s +1
• Pro simulace s pade(3,3) e− s =
>> >> >> >> >>
120 − 60 s + 12 s − s 120 + 60 s + 12 s 2 + s 3
Michael Šebek
2
3
ARI-Pr-25-2011
del=tf(1);set(del,'ioDelay',1); pade3=sdf(pade(del,3)); G=(1+pade3)./(s-pade3); T=(1+pade3)./(s+1); step(tf(G),tf(L),0:.01:8) 33
Příklady k přednášce 25 – Systémy proměnné v čase
Michael Šebek Automatické řízení 2012 21-4-13
Příklad: Houpačka Automatické řízení - Kybernetika a robotika
• dětská houpačka je kyvadlo s rovnicí (standardní předpoklady): d (ml 2ϕ ) + mgl sin ϕ = 0 dt l (t ) ∈ [l − , l + ], L = ½ (l − + l + ) • ale délka je zde proměnná l = • tedy celkem (pozor při derivaci!) d swing.mdl (l (t ) 2 ϕ ) + gl (t )sin ϕ = 0 dt 2lϕ g sin ϕ ϕ + 0 + = l l
Michael Šebek
ARI-Pr-25-2012
35
Pokračování: Parametrická rezonance Automatické řízení - Kybernetika a robotika
• Pro teoretické zkoumání označ. ν =lϕ ,ν =lϕ + lϕ ,ν = l ϕ + 2lϕ + lϕ a zjednodušíme rovnici sin ϕ ≅ ϕ 1 2lϕ g sin ϕ ν + ( g − l )ν = 0 ϕ + + = 0 l
l
l
= ω x x′ dx dt • Dosadíme l (t ) ≅ L(1 + ε cos ωt ) , označíme t ≅ ωt →= 2 2 ′′ = ≅ ω δ ω x x , g ( L ) δ + ε cos t a dostaneme
ν ′′ +
ν= 0 1 + ε cos t 2 2 cos t (1 ) cos t δ + ε − δ ε • Použijeme aproximaci 1. řádu = δ + (1 − δ )ε cos t + 1 + ε cos t 1 + ε cos t v ε a dostaneme tzv.
• Mathieuovu rovnici
swingMat.mdl
ν ′′ + (δ + ε (1 − δ ) cos t )ν = 0
• Ta má neomezené řešení pokud ε = 0.1, δ = ¼ → ω = 2 g L = 2x přirozená frekvence kyvadla Michael Šebek
ARI-Pr-25-2012
36
Stabilita LTV Automatické řízení - Kybernetika a robotika
• Pro lineární systém proměnný v čase = x A(t ) x + B(t )u ,
x(t0 )
t ≥ t0
= y C(t )x + D(t )u
• Je řešení je dáno stavovou maticí přechodu s počáteční hodnotou x(t ) = Φ(t , t )x(t ) 0
0
Φ(t0 , t0 ) = I
• Definice stability je podobná jako u LTI, přesněji • ekvilibrium v počátku je globálně stejnoměrně asymptoticky stabilní právě když − γ ( t −t ) Φ(t , t0 ) ≤ ke
0
, ∀t ≥ t0 ≥ 0
• Stabilitu ale nelze charakterizovat vlastními čísly matice A ani v případě, že jsou tato čísla konstantní ! Michael Šebek
ARI-Pr-25-2012
37
Příklad: Stabilita LTV systému Automatické řízení - Kybernetika a robotika
• LTV systém 2. řádu s A (ostatní matice jsou nulové) −1 + 1.5cos 2 t 1 − 1.5sin t cos t A(t ) = 2 −1 − 1.5sin t cos t −1 + 1.5sin t
• má vlastní čísla nezávislá na t a ležící v levé polorovině » syms t s » A=[-1+1.5*cos(t)^2,1-1.5*sin(t)*cos(t);-11.5*sin(t)*cos(t),-1+1.5*sin(t)^2] A = [ -1+3/2*cos(t)^2, 1-3/2*sin(t)*cos(t)] [ -1-3/2*sin(t)*cos(t), -1+3/2*sin(t)^2] » eig(A) ans = [ -1/4+1/4*i*7^(1/2)] [ -1/4-1/4*i*7^(1/2)]
• Tedy by se zdálo, že je systém stabilní? Michael Šebek
ARI-Pr-25-2012
38
Příklad: Stabilita LTV systému Automatické řízení - Kybernetika a robotika
• Přitom ale je neboť
e0.5t cos t e − t sin t Φ(t , 0) = 0.5t −t e t e t sin cos −
» PHI=[exp(t/2)*cos(t),exp(-t)*sin(t);-exp(t/2)*sin(t),exp(-t)*cos(t)]
PHI =[ exp(1/2*t)*cos(t), [ -exp(1/2*t)*sin(t),
exp(-t)*sin(t)] exp(-t)*cos(t)]
» [simplify(A*PHI(:,1)-diff(PHI(:,1),t)), simplify(A*PHI(:,2)-diff(PHI(:,2),t))]
ans = [ 0, 0] [ 0, 0]
x(t ) = Φ(t , 0)x(0) • Jelikož • Tak zřejmě pp. libovolně blízko počátku, pro které řešení uteče do nekonečna - systém je tedy nestabilní • Pro časově proměnné systémy vlastní čísla nefungují! Michael Šebek
ARI-Pr-25-2012
39
