ĆWICZENIE 1
Identyfikacja obiektów dynamicznych
1.1 Cel ćwiczenia
Celem ćwiczenia jest ilustracja częstotliwościowych i czasowych metod identyfikacji obiektów
dynamicznych.
1.2 Identyfikowane modele
Badane obiekty dynamiczne modelowane są odpowiednią operatorową transmitancją:
a) układ inercyjny pierwszego rzędu
p
p
sT
k
s
G
+
=
1
)
(
(1.1)
gdzie
p
k jest statycznym wzmocnieniem,
p
T oznacza stałą czasową.
b) układ inercyjny pierwszego rzędu z opóźnieniem transportowym
0
1
)
(
sT
p
p
e
sT
k
s
G
−
⋅
+
=
(1.2)
gdzie
p
k jest statycznym wzmocnieniem,
p
T oznacza stałą czasową, zaś
0
T reprezentuje
opóźnienie transportowe.
c) układ całkujący
i
sT
s
G
1
)
(
=
(1.3)
gdzie
i
T oznacza stałą całkowania.
d) układ drugiego rzędu
2
2
2
2
2
2
2
1
2
2
1
1
1
1
)
(
s
s
s
s
a
s
sa
s
G
n
n
n
+
ζω
+
ω
ω
=
τ
+
ζτ
+
=
+
+
=
(1.4)
gdzie
ζ jest współczynnikiem tłumienia, zaś
τ
=
ω
/
1
n
reprezentuje pulsację naturalną
(pulsację drgań nietłumionych).
e) układ nieminimalnofazowy
y
x
sT
sT
s
G
+
−
=
1
1
)
(
(1.5)
gdzie
x
T jest stałą czasową zera, zaś
y
T reprezentuje stałą czasową bieguna.
Identyfikacji podlegają odpowiednie parametry transmitancji (1.1)-(1.5).
1.2.1 Charakterystyki czasowe i częstotliwościowe obiektu inercyjnego pierwszego rzędu
Odpowiedź impulsową obiektu inercyjnego pierwszego rzędu (1.1) opisuje wzór
)
(
)]
(
[
)
(
/
1
t
e
T
k
s
G
L
t
g
p
T
t
p
p
1
⋅
=
=
−
−
(1.6)
Odpowiedź skokowa tego członu dana jest wzorem
)
(
)
1
(
]
/
)
(
[
)
(
/
1
t
e
k
s
s
G
L
t
h
p
T
t
p
1
⋅
−
=
=
−
−
(1.7)
zaś przykładowy przebieg tej odpowiedzi pokazano na rys. 1.1.
0
0.5
1
1.5
2
2.5
3
3.5
4
0
0.2
0.4
0.6
0.8
1
←
→
T
p
h(
t)
t [s]
Rys. 1.1. Odpowiedź skokowa członu dynamicznego pierwszego rzędu dla
1
=
p
k
i
s
T
p
1
= .
Ze wzoru (1.6) wynika, iż czas ustalania T
s
∆
odpowiedzi skokowej (1.7), definiowany jako
)}
1
(
)
(
:
{
∆
−
=
=
∆
p
s
k
t
h
t
T
,
0
1
< ≤
∆
(1.8)
wynosi
∆
⋅
−
=
∆
ln
p
s
T
T
(1.9)
Zachodzi ponadto
p
p
T
t
k
k
e
t
h
p
⋅
≅
−
=
−
=
6321
.
0
)
1
(
)
(
1
(1.10)
p
t
k
t
h
=
→∞
)
(
(1.11)
Widmowa charakterystyka członu pierwszego rzędu (1.1) dana jest wzorem
( )
p
T
j
p
p
j
j
s
e
T
k
e
M
s
G
ω
⋅
−
ω
φ
ω
=
ω
+
=
ω
=
arctg
2
2
)
(
1
)
(
)
(
(1.12)
Przykładowy przebieg charakterystyki amplitudowej
)
(
ω
M
oraz fazowej
)
(
ω
φ
tego członu
pokazano na rys. 1.2. Ze wzoru (1.12) wynika, iż pulsacja trzydecybelowego pasma
przenoszenia członu (1.1) wynosi
p
T
1
dB
3
=
ω
(1.13)
czemu odpowiadają następujące wartości charakterystyki amplitudowej oraz fazowej tego członu
2
)
(
dB
3
P
k
M
=
ω
ω
=
ω
(1.14)
o
45
)
(
dB
3
−
=
ω
φ
ω
=
ω
.
(1.15)
Wzory (1.9)-(1.11) oraz (1.13)-(1.15) mają istotne praktyczne znaczenie, stanowiąc podstawę
prostych reguł identyfikacji członu (1.1).
0
1
2
3
4
5
0
0.707
1
ω
3dB
M(
ω
)
ω
[rad/s]
0
1
2
3
4
5
-90
-45
0
ω
3dB
φ
(
ω
) [
°
]
ω
[rad/s]
Rys. 1.2. Częstotliwościowe charakterystyki układu pierwszego rzędu dla
1
=
p
k
i
s
1
=
p
T
.
1.2.2 Charakterystyki czasowe i częstotliwościowe obiektu inercyjnego pierwszego rzędu
z opóźnieniem transportowym
Skokowa odpowiedź obiektu opisanego transmitancją (1.2) dana jest wzorem
).
(
)
1
(
]
/
)
(
[
)
(
0
/
)
(
1
0
T
t
e
k
s
s
G
L
t
h
p
T
T
t
p
−
⋅
−
=
=
−
−
−
1
(1.16)
Zachodzą przeto następujące związki, które można wykorzystać do identyfikacji rozważanego
modelu:
0
)
(
=
t
h
dla
t T
≤
0
⇒
{
}
T
t h t
0
0
=
=
max : ( )
, (1.17)
p
t
k
t
h
=
∞
→
)
(
lim
,
(1.18)
p
p
T
T
t
k
e
k
t
h
p
⋅
≅
−
=
−
+
=
632
.
0
)
1
(
)
(
1
0
, (1.19)
p
p
T
t
T
k
t
h
t
=
=
0
)
(
d
d
.
(1.20)
Przebieg przykładowej skokowej odpowiedzi h t
( ) pokazano na rys. 1.3.
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
0
0.2
0.4
0.6
0.8
1
←
→
T
p
← →
T
0
h(
t)
t [s]
Rys. 1.3. Skokowa odpowiedź układu pierwszego rzędu z opóźnieniem transportowym
dla 1
=
p
k
, s
1
=
p
T
oraz
s
5
.
0
0
=
T
.
Transmitancji (1.2) przyporządkować można charakterystykę amplitudową
)
(
ω
M
oraz fazową
)
(
ω
φ
:
2
2
)
(
1
)
(
)
(
p
p
j
j
s
T
k
e
M
s
G
ω
+
=
⋅
ω
=
ω
φ
−
ω
=
,
(1.21)
( )
0
arctg
)
(
T
T
p
ω
−
ω
−
=
ω
φ
.
(1.22)
Zachodzą przy tym następujące związki, mogące stanowić podstawę prostych procedur
identyfikacji:
p
k
M
=
)
0
(
,
(1.23)
2
)
(
dB
3
P
k
M
=
ω
ω
=
ω
,
p
T
1
dB
3
=
ω
⇒
dB
3
1
ω
=
p
T
, (1.24)
p
T
T
0
4
)
(
dB
3
−
π
−
=
ω
φ
ω
=
ω
⇒
⎟
⎠
⎞
⎜
⎝
⎛
π
+
ω
ϕ
⋅
−
=
4
)
(
dB
3
2
0
p
T
T
. (1.25)
Przykładowe charakterystyki
)
(
ω
M
oraz
)
(
ω
φ
, dla
1
=
p
k
, s
1
=
p
T
oraz
s
5
.
0
0
=
T
, pokazano na
rys. 1.4 oraz rys. 1.5 charakterystykę Nyquista.
0
5
10
15
0
0.707
1
ω
3dB
M(
ω
)
ω
[rad/s]
0
1
5
10
15
-360
-45
0
↑
↓
∆φ
(
ω
)
ω
3dB
φ
(
ω
) [
°
]
ω
[rad/s]
Rys. 1.3. Charakterystyki częstotliwościowe identyfikowanego modelu.
-0.4
-0.2
0
0.2
0.4
0.6
0.8
1
-0.8
-0.6
-0.4
-0.2
0
0.2
0.4
ReG(j
ω
)
Im
G(
j
ω
)
ω
=0
ω
=
∞
Rys.1.4. Charakterystyka Nyquista identyfikowanego modelu
1.2.3 Charakterystyki czasowe i częstotliwościowe obiektu całkującego
Odpowiedź impulsową obiektu całkującego (1.3) opisuje wzór
)
(
1
)]
(
[
)
(
1
t
T
s
G
L
t
g
i
1
⋅
=
=
−
.
(1.26)
Odpowiedź skokowa tego członu dana jest wzorem
)
(
]
/
)
(
[
)
(
1
t
T
t
s
s
G
L
t
h
i
1
⋅
=
=
−
,
(1.27)
zaś przykładowy przebieg tej odpowiedzi pokazano na rys. 1.5.
0
1
2
3
4
0
1
2
3
4
tg
α
= 1 / T
i
h(
t)
t [s]
Rys. 1.5. Odpowiedź skokowa członu całkującego dla
s
T
i
1
= .
Charakterystyka częstotliwościowa członu całkującego (1.3) dana jest wzorem
2
)
(
1
)
(
)
(
π
−
ω
φ
ω
=
ω
=
ω
=
j
i
j
j
s
e
T
e
M
s
G
.
(1.28)
Ze wzoru (1.28) wynika, iż pulsacja odcięcia
gc
ω członu całkującego (1.3) wynosi
i
gc
T
1
=
ω
(1.29)
czemu odpowiada następująca wartości charakterystyki amplitudowej tego członu
1
dB
0
)
(
=
=
ω
ω
=
ω
gc
M
.
(1.30)
Przykładowy przebieg charakterystyki amplitudowej
)
(
ω
M
oraz fazowej
)
(
ω
φ
tego członu
pokazano na rys. 1.6.
10
-1
10
0
10
1
-20
0
20
ω
gc
M(
ω
) [
d
B
]
ω
[rad/s]
10
-1
10
0
10
1
-90
0
φ
(
ω
) [
°
]
ω
[rad/s]
Rys. 1.6. Charakterystyki Bodego układu całkującego dla
s
T
i
1
= .
1.2.4 Charakterystyki czasowe i częstotliwościowe członu drugiego rzędu
Odpowiedź impulsowa członu drugiego rzędu (1.4) wyraża się wzorem
)
(
sin
)
1
(
)]
(
[
)
(
0
2
1
t
t
e
s
G
L
t
g
t
n
n
1
⋅
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
ω
ζ
−
ω
=
=
ζω
−
−
, (1.31)
przy czym
τ
ζ
−
=
ζ
−
ω
=
ω
)
1
(
)
1
(
2
2
0
n
(1.32)
oznacza pulsację drgań tłumionych. Odpowiedź skokową członu (1.4) określa zależność
),
(
)
1
(
sin
cos
1
)
(
)
sin(
)
1
(
1
1
]
/
)
(
[
)
(
2
0
0
0
2
1
t
t
t
e
t
t
e
s
s
G
L
t
h
t
t
n
n
1
1
⋅
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
ζ
−
ω
ζ
+
ω
−
=
=
⋅
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
α
+
ω
ζ
−
−
=
=
ζω
−
ζω
−
−
(1.33)
gdzie
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
ζ
ζ
−
=
α
)
1
(
arctg
2
.
(1.34)
Przykładowe przebiegi odpowiedzi skokowych dla różnych wartości współczynnika tłumienia
ζ
pokazano na rys. 1.7.
0
1
2
3
4
5
6
7
8
9
10
0
0.2
0.4
0.6
0.8
1
1.2
1.4
h(
t)
t [s]
ζ
1
=0.3
ζ
2
=0.6
ζ
3
=1
Rys. 1.7. Odpowiedź skokowa członu oscylacyjnego
Widmowa charakterystyka członu drugiego rzędu (1.4) dana jest wzorem
2
2
2
arctg
2
2
2
2
2
)
(
)
2
(
)
(
)
(
)
(
ω
−
ω
ω
ζω
⋅
−
ω
φ
ω
=
ω
ζω
+
ω
−
ω
ω
=
ω
=
n
n
j
n
n
n
j
j
s
e
e
M
s
G
,
(1.35)
zaś przykładowe przebiegi funkcji
M ( )
ω oraz
)
(
ω
φ
dla różnych wartości współczynnika
tłumienia
ζ
pokazano na rys. 1.8.
0
1
2
3
4
5
0
0.5
1
1.5
ω
[rad/s]
M(
ω
)
ζ
1
=0.3
ζ
2
=0.6
ζ
3
=1
0
1
2
3
4
5
-180
-90
0
ω
[rad/s]
φ
(
ω
) [
°
]
ζ
1
=0.3
ζ
2
=0.6
ζ
3
=1
Rys. 1.8. Częstotliwościowe charakterystyki członu oscylacyjnego dla
s
1
=
τ
.
Dla odpowiedzi skokowej
h t
( )
definiuje się następujące wskaźniki (por. rys. 1.9):
a) przeregulowanie
κ
κ =
− ∞
∞
⋅
h
h
h
max
( )
( )
100%
(1.36)
gdzie
h
h t
t
max
max ( )
=
≥0
(1.37)
b) czas osiągnięcia maksimum T
κ
(czas piku, czas maksimum)
}
)
(
:
{
max
h
t
h
t
T
=
=
κ
(1.38)
c) czas ustalania T
s
∆
dla strefy kontrolnej o szerokości
2
∆
)}
(
|
)
(
)
(
|:
{
max
arg
0
∞
⋅
∆
=
∞
−
=
≥
∆
h
h
t
h
t
T
t
s
(1.39)
Odpowiednie wskaźniki definiuje się także dla amplitudowej charakterystyki
M ( )
ω
rozważanego członu dynamicznego (por. rys. 1.10)):
a)
wskaźnik oscylacyjności (szczyt rezonansowy) M
r
)
0
(
max
M
M
M
r
=
(1.40)
gdzie
M
M
max
max
( )
=
≥
ω
ω
0
(1.41)
b) pulsacja rezonansowa
ω
r
ω
ω
ω
r
M
M
=
=
{
:
( )
}
max
,
(1.42)
c) trzydecybelowe pasmo przenoszenia
ω
3dB
ω
ω
ω
3
3
0
2
dB
dB
=
=
{
:
(
)
( ) /
}
M
M
.
(1.43)
0
0
↓
↑
∆
h(
t)
t
T
s
∆
T
κ
h
max
h(
∞
)
Rys. 1.9. Definicja wskaźników dotyczących odpowiedzi skokowej
0
0
ω
M(
ω
)
ω
r
ω
3dB
M(0)
M
max
M(
ω
3dB
)
Rys. 1.10. Definicja wskaźników dotyczących charakterystyki amplitudowej
Przeregulowanie
κ
oraz czas maksimum T
κ
wiążą się z parametrami
ζ
,
0
1
< <
ζ
, oraz
τ
transmitancji operatorowej (1.4) następującymi wzorami
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
ζ
−
ζπ
−
=
κ
)
1
(
exp
2
(1.42)
)
1
(
2
ζ
−
τπ
=
κ
T
(1.43)
Czas ustalania T
s
∆
jest nieciągłą funkcją współczynnika tłumienia
ζ
, dla której można podać
następującą ciągłą funkcję majoryzującą
(
)
τ
⋅
ζ
ζ
−
∆
=
≤
∆
∆
)
1
(
ln
2
s
s
T
T
,
0
1
< <
ζ
(1.44)
Ze wzoru (1.44) wynika, iż dla
∆ = 002
. i
∆ = 005
. oraz dla dostatecznie małych wartości
ζ
obowiązują oszacowania
ζ
τ
≅
4
%
2
s
T
(1.45)
ζ
τ
≅
3
%
5
s
T
(1.46)
Wskaźniki M
r
,
ω
r
oraz
ω
3dB
opisujące amplitudową charakterystykę M ( )
ω członu (1.4)
związane są z parametrami
ζ
oraz
τ tego członu następującymi formułami:
)
1
(
2
1
2
ζ
−
ζ
=
r
M
,
0
1
2
< <
ζ
/
(1.47)
τ
ζ
−
=
ω
)
2
1
(
2
r
,
0
1
2
< <
ζ
/
(1.48)
τ
+
ζ
−
+
ζ
−
=
ω
1
)
2
1
(
2
1
2
2
2
dB
3
(1.49)
Wykresy rozważanych wskaźników pokazano na rys. 1.11a-c (wskaźniki dotyczące odpowiedzi
skokowej) oraz na rys. 1.12a-c (wskaźniki dotyczące charakterystyki amplitudowej).
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0
10
20
30
40
50
60
70
ζ
κ [%
]
(a)
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0
2
4
6
8
10
12
14
16
18
20
ζ
T
κ
/ τ
(b)
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
0
5
10
15
20
25
30
35
40
(c)
ζ
T
s
∆
/ τ
T
s5%
/
τ
T
s2%
/
τ
Rys. 1.11. Wskaźniki dotyczące odpowiedzi skokowej członu oscylacyjnego
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
1
1.5
2
2.5
3
3.5
4
4.5
5
ζ
M
r
(a)
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
ζ
ω
r
τ
(b)
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0.4
0.6
0.8
1
1.2
1.4
1.6
(c)
ζ
ω
3d
B
τ
Rys. 1.12. Wskaźniki dotyczące charakterystyki amplitudowej członu oscylacyjnego
Ze wzoru (1.47) wynika, iż
2
)
1
(
1
2
−
−
−
=
ζ
r
M
,
M
r
≥ 1
(1.50)
Zatem, korzystając ze wzoru (1.42) wyznaczyć można zależność przeregulowania
κ
odpowiedzi
skokowej (1.33) od wskaźnika oscylacyjności M
r
. Zachodzi ponadto
κ
+
π
κ
=
ζ
2
2
ln
|
ln
|
,
κ > 0
(1.51)
Na podstawie wzoru (1.44) można określić relację wskaźnika oscylacyjności M
r
w zależności
od przeregulowania , co ilustruje wykres dany na rys. 1.13.
1
1.5
2
2.5
3
3.5
4
4.5
5
0
10
20
30
40
50
60
70
80
M
r
κ
[%
]
Rys. 1.13. Związek między wskaźnikiem oscylacyjności a przeregulowaniem.
Bieguny
s
1 2
,
transmitancji operatorowej (1.4) dane są wzorem
2
2
,
1
1
ζ
−
ω
±
ζ
ω
−
=
n
n
j
s
(1.52)
zaś ich rozmieszczenie na płaszczyźnie zespolonej pokazano na rys. 1.14.
)
(
Im
ω
j
G
)
(
Re
ω
j
G
α
n
ζω
−
2
1
ζ
−
ω
n
j
Rys. 1.14. Położenie biegunów transmitancji członu oscylacyjnego
Z rysunku tego wynika, iż
ζ
ζ
−
=
α
2
1
tg
(1.53)
a zatem
ζ
=
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
ζ
ζ
−
=
α
cos
arc
)
1
(
tg
arc
2
(1.54)
Załóżmy, iż (1.4) oznacza operatorową transmitancję zamkniętego układu sterowania o
strukturalnym schemacie danym na rys. 1.15, gdzie
G s
k
s
Ts
0
1
( )
(
)
=
+
(1.55)
jest transmitancją operatorową układu otwartego, zaś
ζτ
=
2
1
k
(1.56)
ζ
τ
=
2
T
(1.57)
)
(t
r
)
(t
c
)
(
0
s
G
Rys. 1.15. Schemat strukturalny układu sterowania
Pulsacja odcięcia
ω
gc
amplitudowej charakterystyki transmitancji rozważanego układu
otwartego zdefiniowana jest wzorem
1
|
)
(
|
0
=
=
gc
j
s
s
G
ω
(1.58)
zaś jej związek z parametrami transmitancji (1.4) opisuje formuła
τ
ζ
ζ
ω
2
4
2
1
4
−
+
=
gc
(1.59)
Zapas fazy rozważanego układu sterowania wyraża się wzorem
⎟⎟
⎟
⎠
⎞
⎜⎜
⎜
⎝
⎛
ζ
−
+
ζ
ζ
=
+
π
=
∆
ω
=
2
4
0
2
1
4
2
arctg
)
(
arg
gc
j
s
p
s
G
(1.60)
Przebieg funkcji
ω ζ
gc
( ) przedstawiono na rys. 1.16, zaś funkcję
∆
p
( )
ζ ilustruje rys. 1.17.
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
0.4
0.5
0.6
0.7
0.8
0.9
1
ζ
ω
gc
τ
Rys. 1.16. Pulsacja odcięcia
ω
gc
w funkcji współczynnika tłumienia.
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
10
20
30
40
50
60
70
80
ζ
∆
p
[
°
]
Rys. 1.17. Zapas fazy
∆
p
w funkcji współczynnika tłumienia.
Niech
2
/
)
(
arg
2
/
π
−
=
ω
π
−
ω
=
ω
j
G
(1.61)
2
/
)
(
arg
2
/
|
)
(
|
π
−
=
ω
π
−
ω
=
j
G
j
G
M
(1.62)
0
,
1
|
)
(
>
ω
=
ω
ω
=
ω
j
G
gc
,
2
/
1
<
ζ
(1.63)
ϑ
ω
gc
gc
G j
= arg (
)
,
0
1
2
< <
ζ
/
(1.64)
gc
ω
ω
=
γ
dB
3
(1.65)
Jak łatwo sprawdzić zachodzą następujące relacje:
2
cos
tg
p
p
∆
∆
=
ζ
(1.66)
1
)
2
1
(
2
1
cos
2
2
2
dB
3
+
ζ
−
+
ζ
−
∆
ω
=
ω
p
gc
(1.67)
p
p
gc
∆
⋅
τ
ζ
=
τ
∆
=
ω
tg
2
)
(cos
(1.68)
2
2
4
2
1
2
1
4
ζ
−
ζ
−
+
ζ
ω
=
ω
r
gc
(1.69)
p
s
gc
T
∆
≈
ω
tg
6
%
5
(1.70)
p
s
gc
T
∆
≈
ω
tg
8
%
2
(1.71)
τ
=
ω
π
−
1
2
/
(1.72)
ζ
=
π
−
2
1
2
/
M
(1.73)
τ
ζ
−
=
ω
)
2
1
(
2
2
gc
,
0
1
2
< <
ζ
/
(1.74)
1
)
2
1
(
2
1
2
2
2
+
ζ
−
+
ζ
−
=
γ
(1.75)
Na podstawie powyższych rozważań opracowano tabelę formuł, przydatnych przy identyfikacji
modeli drugiego rzędu (Tabela 1.1).
Dane
wejściowe
Parametry modelu drugiego rzędu Ograniczenia
M
r
r
,
ω
r
r
r
M
M
M
2
1
2
−
−
=
ζ
,
r
ω
ζ
−
=
τ
2
2
1
2
1
0
<
ζ
<
2
/
2
/
,
π
−
π
−
ω
M
|
)
(
2
1
2
/
π
−
ω
=
ζ
M
,
2
/
1
π
−
ω
=
τ
ω
ω
π
gc
,
/
− 2
2
/
2
2
π
−
ω
ω
=
ζ
gc
,
2
/
1
π
−
ω
=
τ
2
1
0
<
ζ
<
ω
ω
π
3
2
dB
,
/
−
γ
+
γ
−
γ
=
ζ
2
1
2
4
2
,
2
1
π/
−
ω
=
τ
ϑ
ω
gc
gc
,
⎪
⎪
⎪
⎪
⎩
⎪
⎪
⎪
⎪
⎨
⎧
π
−
>
ϑ
ϑ
+
+
π
−
≤
ϑ
ϑ
+
−
=
ζ
,
2
2
tan
1
1
1
,
2
2
tan
1
1
1
2
2
gc
gc
gc
gc
,
gc
ω
ζ
−
=
τ
)
2
1
(
2
2
2
1
0
<
ζ
<
κ
κ
,T
κ
+
π
κ
=
ζ
2
2
ln
|
ln
|
,
κ
ζ
−
=
τ
κ
ln
T
κ > 0
Tabela 1.1. Formuły służące identyfikacji modelu drugiego rzędu
Należy podkreślić, iż metody identyfikacji oparte o pary danych pomiarowych
)
,
(
2
/
2
/
π
−
π
−
ω
M
oraz (
,
)
/
ω
ω
π
3
2
dB
−
pozwalają na identyfikację członu dynamicznego drugiego rzędu w
przypadku, gdy ζ ≥ 1
2
/
, w tym także członu dwuinercyjnego
ζ ≥ 1
.
1.2.5 Charakterystyki czasowe i częstotliwościowe członu nieminimalnofazowego
Odpowiedź impulsową układu nieminimalnofazowego (1.5) opisuje wzór
[
]
)
(
1
1
)
(
)
(
/
1
t
e
T
T
T
s
G
L
t
g
y
T
t
y
x
y
1
⋅
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
+
=
=
−
−
(1.76)
Odpowiedź skokowa tego członu dana jest wzorem
)
(
1
1
)
(
)
(
/
1
t
e
T
T
s
s
G
L
t
h
y
T
t
y
x
1
⋅
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
+
−
=
⎥⎦
⎤
⎢⎣
⎡
=
−
−
(1.77)
zatem dla
0
>
x
T
0
)
0
(
<
−
=
y
x
T
T
h
(1.78)
1
)
(
=
∞
h
(1.79)
Przykładowe znormalizowane przebiegi odpowiedzi skokowej dla różnych wartości
y
x
T
T /
pokazano na rys.1.18.
0
1
2
3
4
5
6
-6
-5
-4
-3
-2
-1
0
1
t [s]
h(
t)
T
x
/ T
y
= 2
T
x
/ T
y
= 4
T
x
/ T
y
= 6
Rys. 1.18. Odpowiedź skokowa członu dynamicznego z nieminimalnofazowym zerem.
Na podstawie odpowiedzi skokowej można dokonać identyfikacji stałych czasowych
x
T oraz
y
T .
Wprowadza się dodatkową zmienną pomocniczą
y
x
T
T
h
x
=
−
=
)
0
(
(1.80)
Z kolei na podstawie wzoru (1.77), znając
0
t (Rys. 1.19) jaki upłynął od momentu, gdy
odpowiedź
0
)
(
0
=
t
h
można wyznaczyć
)
1
ln(
0
x
t
T
y
+
=
(1.81)
Stałą czasową
x
T wyznacza się na podstawie wzoru (1.80)
h(t
0
)=0
h(0)
t
0
t
h(
t)
Rys. 1.19. Wyznaczenie chwili czasu
0
t .
Widmowa charakterystyka członu dynamicznego z nieminimalnofazowym zerem (1.5) dana jest
wzorem
)
arctan
(arctan
2
2
2
2
)
(
1
1
)
(
)
(
y
x
T
T
j
y
x
j
j
s
e
T
T
e
M
s
G
ω
+
ω
−
ω
φ
ω
=
ω
+
ω
−
=
ω
=
(1.82)
zatem
1
)
0
(
=
M
(1.83)
y
x
T
T
M
=
∞)
(
(1.84)
Przykładowe przebiegi charakterystyki amplitudowej
)
(
ω
M
oraz fazowej
)
(
ω
φ
tego członu dla
różnych wartości
y
x
T
T / pokazano na Rys. 1.20.
0
1
2
3
4
5
6
0
2
4
6
ω
[rad/s]
|M
(
ω
)|
T
x
/ T
y
= 2
T
x
/ T
y
= 4
T
x
/ T
y
= 6
0
1
2
3
4
6
0
-180
-90
0
ω
[rad/s]
φ
(
ω
) [
°
]
T
x
/ T
y
= 2
T
x
/ T
y
= 4
T
x
/ T
y
= 6
Rys. 1.20. Charakterystyki częstotliwościowe członu z nieminimalnofazowym zerem.
Aby dokonać identyfikacji stałych czasowych
x
T oraz
y
T w oparciu o charakterystyki widmowe,
w pierwszej kolejności należy wyznaczyć pulsację
2
/
π
−
ω
, dla której
2
)
(
2
/
π
−
=
ω
ϕ
π
−
(Rys. 1.21).
Wykorzystując wzór (1.82) można pokazać, że dla pulsacji
2
/
π
−
ω
charakterystyka amplitudowa
ma następującą wartość
)
(
)
(
2
/
∞
=
=
ω
π
−
M
T
T
M
y
x
(1.85)
Z powyższego faktu wynika, że pulsację
2
/
π
−
ω
można wyznaczyć także w oparciu o
charakterystykę
)
(
ω
M
. Dodatkowo ze wzoru (4.8) oraz (4.11) można wyprowadzić zależność
2
/
)
(
π
−
ω
∞
=
M
T
x
(1.86)
Stałą czasową
y
T wyznacza się w oparciu o wzór (4.10).
M(
ω
-
π/2
)
M(
∞
)
ω
-
π/2
ω
M(
ω
)
-180
-90
0
ω
-
π/2
ω
φ
(
ω
) [
°
]
Rys. 1.21. Wyznaczenie pulsacji
2
/
π
−
ω
.
1.3 Opis stanowiska laboratoryjnego
Stanowisko laboratoryjne tworzą:
1. Zestaw Analogowych Modeli Procesów Przemysłowych (ZAMPP)
Zawiera on w sobie obiekty dynamiczne opisane transmitancjami (1.1)-(1.5) oraz miernik
przesunięcia fazowego sygnału wyjściowego względem sygnału wejściowego. Widok płyty
czołowej układu ZAMPP przedstawiono na rys.1.22. Z kolei Rys. 1.23 prezentuje schematy
ideowe identyfikowanych układów.
Rys. 1.22. Widok zestawu laboratoryjnego.
Rys. 1.23. Schematy ideowe badanych układów.
2. Wielofunkcyjny Zestaw Pomiarowy typu MS-9140
Przyrząd ten zawiera między innymi:
• generator funkcji, stanowiący źródło wejściowych sygnałów periodycznych podawanych
na wejście modelu układu,
• częstościomierz, umożliwiający odczyt częstotliwości sygnałów wejściowych.
3. Oscyloskop dwukanałowy.
Umożliwia wizualizację sygnałów wejściowego i wyjściowego oraz pomiar ich amplitud i
parametrów czasowych.
Na Rys. 1.24. przedstawiono widok stanowiska pomiarowego.
Rys. 1.24. Widok stanowiska pomiarowego.
1.4 Zadania pomiarowe
a) należy pomierzyć charakterystykę amplitudową oraz fazową badanych obiektów;
b) w oparciu o charakterystyki czasowe należy oszacować wskaźniki odpowiedzi skokowej
badanych obiektów
Uwagi:
• postać sygnałów wejściowych, to znaczy ich amplitudy oraz pulsacje, należy
dobierać w ten sposób, aby spełnione były warunki umożliwiające racjonalną
identyfikację badanego układu, co sprowadza się przede wszystkim do postulatu
stosowania pobudzeń, przy których układ laboratoryjny pracuje w zakresie
liniowym;
• przy pomiarze charakterystyk amplitudowo-fazowych zanotować charakterystyczne
punkty odpowiednio dla każdego układu (np. dla układu inercyjnego pierwszego
rzędu – pulsację trzydecybelową, przesunięcie fazowe dla pulsacji trzydecybelowej
oraz wzmocnienie w stanie ustalonym);
• pomiary przeprowadzać w zakresie dolnych częstotliwości;
• przy pomiarze charakterystyk czasowych należy oszacować wskaźniki odpowiedzi
skokowej (nie mierzyć czasu ustalania
s
T
!!!).
• odpowiednie formuły służące identyfikacji parametrów poszczególnych układów
można znaleźć w poprzednich punktach owego rozdziału. Przykładowo dla układu
inercyjnego pierwszego rzędu z członem opóźniającym można posłużyć się
formułami (1.5)-(1.7).
1.5 Opracowanie wyników
W sprawozdaniu z ćwiczenia należy:
• zamieścić wykresy charakterystyk Bodego oraz Nyquista badanych obiektów;
• przedstawić zidentyfikowane modele tych układów;
• przeprowadzić dyskusję wyników (porównanie różnych metod identyfikacji);
• wyznaczyć modele optymalne w sensie odpowiednich nieliniowych zadań najmniejszych
kwadratów. W tym celu stosuje się MATLABową funkcję fmins, która rozwiązuje
odpowiednie zadanie optymalizacji.
1.6 Oprogramowanie wspomagające przetwarzanie danych pomiarowych
MATLABowa funkcja fmins pozwala na optymalizację parametrów modeli (1.1)-(1.5) na
podstawie zbioru dostępnych danych pomiarowych
n
j
j
m
j
out
U
1
)}
(
),
(
{
=
ω
Φ
ω
oraz
in
U
, gdzie
)
(
j
out
U
ω jest amplitudą wyjściowego sygnału harmonicznego, Φ
m
j
(
)
ω jest przesunięciem fazy
tego sygnału w stosunku do fazy wejściowego sygnału harmonicznego o amplitudzie
in
U
i
pulsacji
ω
j
,
n
j
,
,
1 …
=
. Rozwiązywane zadanie optymalizacji ma postać następującego
problemu z kwadratową funkcją kosztów
2
/
1
1
2
0
2
0
)
,
(
0
)]
(
sin
/
)
(
)
,
,
,
(
[Im
)]
(
cos
/
)
(
)
,
,
,
(
[Re
min
arg
)
,
(
0
⎪⎭
⎪
⎬
⎫
⎪⎩
⎪
⎨
⎧
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
ω
Φ
⋅
ω
−
ω
+
ω
Φ
⋅
ω
−
ω
=
∑
=
∗
∗
n
j
j
m
in
j
out
i
j
j
m
in
j
out
i
j
T
T
i
U
U
T
T
k
H
U
U
T
T
k
H
T
T
i
(1.87)
Funkcję fmins można wywołać z wiersza poleceń MATLABa w następujący sposób:
>> x = fmins('FunOpt',x0)
gdzie FunOpt oznacza nazwę funkcji kryterialnej według której poszukiwany jest optymalny
wektor parametrów x, zaś x0 jest punktem startowym procedury optymalizacyjnej.
Do optymalizacji parametrów identyfikowanych modeli należy użyć następujących funkcji
kryterialnych dla modelu:
• funm1
–
dla modelu inercyjnego pierwszego rzędu oraz dla układu całkującego w
pętli sprzężenia jednostkowego opisanych odpowiednio transmitancjami
(1.1) i (1.3);
• funm1t –
dla modelu pierwszego rzędu z opóźnieniem transportowym opisanego
transmitancją (1.2)
• funm2
–
dla modelu pierwszego rzędu z opóźnieniem transportowym opisanego
transmitancją (1.4)
• funmnf –
dla modelu nieminimalnofazowego opisanego transmitancją (1.4)
Ponadto należy zdefiniować w MATLABie dane pomiarowe w postaci macierzy o następującej
strukturze:
DATA_M
=[freq Uout Phase]
gdzie freq oznacza kolumnę częstoliwości - [Hz], Uout jest kolumną zawierającą amplitudą
sygnału wyjściowego - [V] (przy założeniu, że amplituda napięcia wejściowego wynosi 1V) oraz
Phase
reprezentuje kolumnę przesunięcia fazowego - [
ο
]. Następnie należy zdefiniować zmienną
DATA_M
jaką zmienną globalną poprzez polecenie:
>> global DATA_M
Dla celów graficznej wizualizacji stopienia dopasowania identyfikowanych modeli stworzono
dodatkowe m-pliki:
• gvidm1 –
dla modelu inercyjnego pierwszego rzędu oraz dla układu całkującego w
pętli sprzężenia jednostkowego opisanych odpowiednio transmitancjami
(1.1) i (1.3);
• gvidm1t –
dla modelu pierwszego rzędu z opóźnieniem transportowym opisanego
transmitancją (1.2)
• gvidm2 –
dla modelu pierwszego rzędu z opóźnieniem transportowym opisanego
transmitancją (1.4)
• gvidmnf –
dla modelu nieminimalnofazowego opisanego transmitancją (1.4)
Powyższe funkcje można wywołać z wiersza poleceń MATLABa w następujący sposób:
>> gvidm1 (x)
>> gvidm1t (x)
>> gvidm2 (x)
>> gvidmnf (x)
gdzie x oznacza znaleziony optymalny wektor parametrów. Przykładowa charakterystyka
Bodego, uzyskana za pomocą omawianych MATLABowych funkcji, przedstawiono na rys. 1.25.
10
1
10
2
10
3
10
4
10
5
-20
-15
-10
-5
0
5
ω
[rad/s]
|G
(j
ω
)|
a
1
= 0.00022199 a
2
= 5.1313e-008
10
1
10
2
10
3
10
4
10
5
-200
-150
-100
-50
0
ω
[rad/s]
Φ
(
ω
)
punkty pomiarowe
model zidentyfikowany
Rys. 1.24. Przykładowa charakterystyka częstotliwościowa modelu drugiego rzędu.