fc2ブログ

制御工学その18

5.3 垂直駆動アームのPID制御



例2.4の垂直駆動アームを対象とし、PID制御により角度制御を実際に行ってみる。
今回の対象は

伝達関数

P(s) = b0/(s^2+a1s+a0)

で表すことが出来、かつ

a0 = Mlg/J
a1 = c/J
b0 = 1/J

である。

5.3.1 P制御

Pコントローラ式

C(s) = kp

を用いると、目標値、入力外乱から制御量への伝達関数はそれぞれ

G_yr = P(s)C(s)/1+P(s)C(s) = b0/(s^2+a1s+a0)*kp / 1 + b0/(s^2+a1s+a0)*kp = b0kp / s^2+a1s+a0+b0kp

G_yd = P(s)/1+P(s)C(s) = b0 / s^2+a1s+a0+b0kp

これは2次遅れ系なので、2次遅れ標準形に治すと

ωn = sqrt(a0+b0*kp)
ζ = a1/2sqrt(a0+b0*kp)
K1 = b0*kp/(a0+b0*kp)
K2 = b0/(a0+b0*kp)

つまり

G_yr = K1*ωn^2 / (s^2+2ζωn+ωn^2) = b0*kp/(a0+b0*kp)*sqrt(a0+b0*kp)^2 / (s^2 + 2*a1/2sqrt(a0+b0*kp)*sqrt(a0+b0*kp) + sqrt(a0+b0*kp)^2) = b0kp / s^2+a1s+a0+b0kp

G_yd = K2*ωn^2 / (s^2+2ζωn+ωn^2) = b0/(a0+b0*kp)*sqrt(a0+b0*kp)^2 / (s^2 + 2*a1/2sqrt(a0+b0*kp)*sqrt(a0+b0*kp) + sqrt(a0+b0*kp)^2) = b0 / s^2+a1s+a0+b0kp0

となる。ここで、kpに注目すれば、2次遅れ要素の固有角周波数ωnと減衰係数ζに同時に現れているため、どちらもkpにより操作することができるが、
逆にどちらか一つしか指定することができない。
つまり、速応性か減衰性のどちらか一つしか指定できない。

また、定常位置偏差

e_∞ = 1/(1+lim[s→0]{P(s)C(s)}) = 1/(1+b0kp/a0) = a0/a0+b0kp

が必ず発生し、
ステップ状の入力外乱に対する入力外乱応答

G_yd(0) = K2 = b0/(a0+b0*kp)

がkpがいかなる値となっても0にはならない。
しかし、定常位置偏差、入力外乱応答はkpが大きくなればなるほど0に近づいてはいく。

シミュレートしてみると、kpを大きく、つまりωnを大きくすればするほど
入力外乱応答、定常位置偏差は0に近づいていくが、振動的な応答になる。

問題5.1

ζを指定するときは

ζM = a1/2sqrt(a0+b0*kp)

より

kp = a1^2/4b0ζM^2 - a0/b0

とすればいい。


5.2.2 P-D制御

P制御では速応性か減衰性のいずれか一方しか指定できなかったので、今回はP-D制御で解決させたい。

P(I)-Dコントローラの式

u(s) = kp*e(s) + kI/s*e(s) - kDs*y(s)

を変形して

u(s) = (Kp+kI/s)*(r(s)-y(s)) - kDs*y(s) = (kp+kI/s)*r(s) - (kp+kI/s+kDs)*y(s)

  = (kp+kI/s+kDs)*(((kp+KI/s)/(kp+kI/s+kDs))*r(s) - y(s))

となるので、ここで

C1(s) = kp+kI/s+kDs
C2(s) = (kps+KI)/(kDs^2+kps+kI)

とすると,目標値、入力外乱から制御量への伝達関数は

G_yr = P(s)C1(s)C2(s)/(1+P(s)C1(s))
G_yd = P(s)/(1+P(s)C1(s))

となる。今回はP-DなのでkI=0として計算すると

G_yr = b0/(s^2+a1s+a0)*(kp+kDs)*(kps)/(kDs^2+kps) / 1+b0/(s^2+a1s+a0)*(kp+kDs)

= b0/(s^2+a1s+a0)*(kp) / (s^2+a1s+a0+b0(kp+kDs))/(s^2+a1s+a0)

= kpb0/(s^2+(a1+kDb0)s+a0+kpb0)

G_yd = b0/(s^2+(a1+kDb0)s+a0+kpb0)

となるので、これを二次遅れ標準形に直すと

ωn = sqrt(a0+b0*kp)
ζ = (a1+b0*kD)/2sqrt(a0+b0*kp)
K1 = b0*kp/(a0+b0*kp)
K2 = b0/(a0+b0*kp)

つまり

G_yr = K1*ωn^2 / (s^2+2ζωn+ωn^2) = b0*kp/(a0+b0*kp)*sqrt(a0+b0*kp)^2 / (s^2 + 2*(a1+b0*kp)/2sqrt(a0+b0*kp)*sqrt(a0+b0*kp) + sqrt(a0+b0*kp)^2) = b0kp / s^2+(a1+kDb0)s+a0+b0kp

G_yd = K2*ωn^2 / (s^2+2ζωn+ωn^2) = b0/(a0+b0*kp)*sqrt(a0+b0*kp)^2 / (s^2 + 2*(a1+b0*kp)/2sqrt(a0+b0*kp)*sqrt(a0+b0*kp) + sqrt(a0+b0*kp)^2) = b0 / s^2+(a1+kDb0)s+a0+b0kp0

となるので、

kp = (ωM^2-a0)/b0
kD = (2*ζM*ωM-a1)/b0

と選べば、任意の固有角周波数、減衰係数を設定することが出来る。
つまり、過渡特性の改善が行われる。
しかし、K1、K2は同じなので、定常特性の改善はできない。

問題5.2

本文中にやってるやんけ



5.3.3 PI-D制御

P制御やP-D制御では、コントローラに積分器を含んでいないため、定常位置偏差や入力外乱定常値が0にならないという問題があった。
今回はPI-D制御で解決させたい。また、その時の部分モデルマッチング法により、パラメータを決定する。

PI-Dのコントローラは前節で示したとおり、

u(s) = (kp+kI/s+kDs)*(((kp+KI/s)/(kp+kI/s+kDs))*r(s) - y(s))

であり、

C1(s) = kp+kI/s+kDs
C2(s) = (kps+KI)/(kDs^2+kps+kI)

とすると,目標値、入力外乱から制御量への伝達関数は

G_yr = P(s)C1(s)C2(s)/(1+P(s)C1(s))
G_yd = P(s)/(1+P(s)C1(s))

である。今回はkI=0ではないのでそのまま計算すると

G_yr = b0/(s^2+a1s+a0)*(kp+kI/s+kDs)*(kps+KI)/(kDs^2+kps+kI) / 1+b0/(s^2+a1s+a0)*(kp+kI/s+kDs)

= b0/(s^2+a1s+a0)*(kps+kI)/s / (s^2+a1s+a0)+b0((kp+kI/s+kDs))/(s^2+a1s+a0)

= b0(kps+kI)/s / s^2+a1s+a0+b0kp+b0kI/s+b0kDs = = b0(kps+kI)/s / (s^3+a1s^2+a0s+b0kps+b0kI+b0kDs^2)/s

=b0(kps+kI)/(s^3+(a1+b0kD)s^2+(a0+b0kp)s+b0kI)

G_yd = b0s/(s^3+(a1+b0kD)s^2+(a0+b0kp)s+b0kI)

である。ここから、パラメータを何とかして決定する。

(a)目標値応答に注目した部分的モデルマッチング法

目標値応答に注目するので、G_yrを規範モデルGMと近似的に一致させる方法である。

そのままでは難しいので、伝達関数の逆数

1/G_yr をs=0でテイラー級数展開したもの

1+a0/b0kI*s+((a1+b0kD)/b0kI - a0kP/b0kI^2)*s^2 + (1/b0kI - kp(a1+b0kD)/b0kI^2 + a0kp^2/b0kI^3)*s^3 ...

計算略!(一次導関数くらいでめんどくさくなった

と、規範モデルGMの逆数1/GMを近似的に一致させればいいので

例えば規範モデルを2次遅れ系

GM(s) = ωM^2/s^2+2ζM*ωM*s+ωM^2

とすると、1/GMと1/G_yrのテイラー展開のsに関する係数が一致すればいいので


a0/b0kI = 2ζM/ωM

(a1+b0kD)/b0kI - a0kP/b0kI^2 = 1

1/b0kI - kp(a1+b0kD)/b0kI^2 + a0kp^2/b0kI^3 = 0


となるように、kI、kD、kpを定めればよい。

kI = ωM*a0/2*ζM*b0

kp = ωM^2/b0

kD = 2*ζM*ωM/b0 + a0/2*ζM*ωM*b0 + a1/b0

と選べばよい。

(b)入力外乱応答に注目した部分的モデルマッチング法

同様に、入力外乱応答に注目した部分的モデルマッチング法は、
G_yd(s)を規範モデルs*GM/kIと近似させる方法である。?

前回同様、

1/GMとs/kIG_yd(s)を一致させたいので、

s/kIG_yd(s)をテイラー級数展開し、(今回はテイラー展開しなくてもなんとかなる)

s/kIG_yd(s) = (s^3+(a1+b0kD)s^2+(a0+b0kp)s+b0kI)/b0kI



突然現れた3次の規範モデル

GM(s) = ωM^3/(s3+α2ωMs^2+α1ωM^2*s+ωM^3)

の逆数と係数が一致すればいいので

kI = ωM^3/b0

kp = (α1ωM^2-a0)/b0

kD = (α2ωM-a1)/b0

とすればよい。

この値で実際にシミュレートしてみると、目標値応答に注目した場合、目標値応答は適度な速応性と減衰性を兼ね備えており、
定常位置偏差も0であるが、入力外乱の影響を完全に抑制するまで時間がかかる。
また、入力外乱応答に注目した場合、入力外乱の影響は速やかに抑制できるが、目標値応答におけるオーバーシュートが大きい。

問題5.3

G_yr式にパラメータをぶち込んでみる。

G_yr = b0(kps+kI)/(s^3+(a1+b0kD)s^2+(a0+b0kp)s+b0kI)

kI = ωM*a0/2*ζM*b0

kp = ωM^2/b0

kD = 2*ζM*ωM/b0 + a0/2*ζM*ωM*b0 + a1/b0

とすると

G_yr = (b0*(ωM^2/b0)*s + b0*ωM*a0/2*ζM*b0) / s^3 + a1*s^2+b0*(2*ζM*ωM/b0 + a0/2*ζM*ωM*b0 + a1/b0)*s^2 + a0s + b0*ωM^2/b0*s + b0*ωM*a0/2*ζM*b0

=ωM^2*s + ωM*a0/2*ζM / s^3 + (a1 + 2*ζM*ωM + a0/2*ζM*ωM + a1)*s^2 + (a0 + ωM^2)*s + ωM*a0/2*ζM

=スキップ


5.3.4 I-PD制御

PI-D制御で、入力外乱に注目すれば、入力外乱の影響は速やかに抑制できたが、オーバーシュートが大きいという問題があった。

そこで、今回は入力外乱の影響を素早く抑制でき、しかも目標値応答を改善するI-PD制御を設計する。

I-PDコントローラの式を用いると
3回目なので略して

G_yr = b0kI/(s^3 + (a1+b0kD)*s^2 + (a0+b0kp)*s + b0kI)
G_yd = b0s/(s^3 + (a1+b0kD)*s^2 + (a0+b0kp)*s + b0kI)

となる。このとき、G_yrはPI-D制御で入力外乱応答に注目した時のGMの式

GM = [PI-D]G_yd*kI/s

と等しいため、G_yrを3次の規範モデルと等しくさせるためのパラメータは
PI-Dで入力外乱応答に注目した時のパラメータに等しい。

実際に同じパラメータでシミュレートしてみると、I-PDは優れた入力外乱抑制特性を保持しつつ、目標値応答を改善している。

問題5.4

対象となる台車系の制御対象は

P(s) = 1/(Ms^2+cs)

I-PDコントローラは

C1(s) = kp+kI/s+kDs
C2(s) = kI/(kDs^2+kps+kI)

よって、

G_yr = P(s)C1(s)C2(s)/(1+P(s)C1(s)) = 1/(Ms^2+cs)*(kp+kI/s+kDs)*(kI/(kDs^2+kps+kI)) / 1 + 1/(Ms^2+cs)*kp+kI/s+kDs

= 1/(Ms^2+cs)*kI/s / (Ms^2+cs+kp+kI/s+kDs)/(Ms^2+cs)

= kI/s / Ms^2+cs+kp+kI/s+kDs = kI/s / Ms^3+cs^2+kps+kI+kDs^2 /s

= kI/Ms^3+(c+kD)*s^2+kps+kI


*別方向からのアプローチ
P(s)を二次遅れ標準形に直すと

P(s) = (1/M) / s^2 + c/Ms)

a1 = c/M
a0 = 0
b0 = 1/M

前節で求めた標準式は

G_yr = b0kI/(s^3 + (a1+b0kD)*s^2 + (b0kp)*s + b0kI)

分母をb0でくくりだすと

G_yr = b0kI/b0(s^3/b0 + (a1+kD)*s^2 + (kp)*s + kI)

   = kI/s^3/(s^3/b0 + (a1+kD)*s^2 + (kp)*s + kI)

となる。

問題に戻って、3次の規範モデルは

GM(s) = ωM^3/(s3+α2ωMs^2+α1ωM^2*s+ωM^3)

G_yrを最も簡単な形の状態で、それぞれの逆数を考える。

1/G_yr = (s^3 + (a1+b0kD)*s^2 + (b0kp)*s + b0kI)/b0kI

1/GM = (s3+α2ωMs^2+α1ωM^2*s+ωM^3)/ωM^3

係数比較すると

1/ωM^3 = 1/b0kI

α2/ωM^2 = (a1+b0kD)/b0kI

α1/ωM = kp/kI

より

kI = ωM^3/b0
kp = α1/ωM * ωM^3/b0 = α1ωM^2/b0
kD = (b0kI*α2/ωM^2 - a1)/b0 = (ωM^3*α2/ωM^2 - a0) / b0 = (α2ωM - a0)/b0


つまり

kI = ωM^3*M
kp = α1*ωM^2*M
kD = (α2*ωM^2 - c/M) * M = α2*ωM*M - c

となる。

まとめ

垂直駆動アームを対象に実際に様々な制御タイプを考えた。
P制御、P-D制御、PI-D制御、I-PD制御と進んでいくにつれ

パラメータの決定はめんどくさくなるが、
I-PD制御では、入力外乱応答の抑制を速やかに行ないながら、目標値応答におけるオーバーシュートが小さい、すなわち減衰性を抑えることが出来た。


わからないところ

PI-D制御のパラメータ決定に関して、入力外乱応答に注目した際に、どうして通常の伝達関数G_ydにkI/sとかいう係数がつくのだろうか。
確かにあとの方から見るとすごく綺麗にまとまるんだけど、どこから出てきたんだろうか。積分動作だけ外れているとかそういうことなのか?

問題5.3はとちゅうで計算力が尽きた。頑張ればなんとかなる気がするけどとてもめんどくさい。
スポンサーサイト



コメントの投稿

非公開コメント

No title

大変参考になりました。このエントリの問題番号などは、何かの書籍の内容と一致しているのでしょうか。もしよろしければ教えて頂けると幸いです。

承認待ちコメント

このコメントは管理者の承認待ちです
検索フォーム
RSSリンクの表示
リンク
ブロとも申請フォーム

この人とブロともになる

QRコード
QR