有色雑音を有する一入力一出力システムに対する雑音ゲイン関数を用いたパラメータ同定法
8
0
0
全文
(2) 30. システム制御情報学会論文誌 第 34 巻 第 2 号 (2021). あらかじめ計測できるが，フィードバックを行うための 計測装置には，未知の雑音が入るようなシステムにおい て成り立つ． 本同定アルゴリズムでは，出力雑音の分散と自己共分 散の値を求めるために，伝達関数パラメータから導出さ れる雑音ゲイン関数を用いる．雑音ゲイン関数により入 力雑音の分散と自己共分散の値から，有色入力雑音に よって出力に生じる雑音成分の分散と自己共分散の値を. Fig. 2. The target SISO system. 求めることができる． なお本同定アルゴリズムは，誤差分散行列を仮定して. えられ，出力誤差として ηt に有色ガウス雑音 Δyt が加え. システムパラメータを同定する計算と，システムパラ. られ出力 y となっている．ただし，Δut ,Δyt の平均は 0. メータを仮定して誤差分散行列を推定する計算を繰り返. とする．また，Δut と Δyt は無相関とし，Δut の分散と. す計算アルゴリズムである．以下の章において，まず対. 2 ,τΔu,k (k = 1,2,...)，Δyt の分散と自 自己共分散を σΔu. 象とする有色雑音を有するシステムにおける仮定および，. 2 2 ,τΔy,d (d = 1,2,...) とする．また，σΔu 己共分散を σΔy. 対象とするシステムについて述べる．つぎに，同定アル. 2 と σΔy を用いるより，それらを τΔu,0 と τΔy,0 と表した. ゴリズムを示し，誤差分散行列の推定値から伝達関数パ. 方が後の計算式が簡素に表現できる．したがって，次の. ラメータの推定値を求める方法，伝達関数パラメータか. ように仮定する．. ら誤差分散行列の推定値を求める方法を順に示す．さら. E[Δut ] = 0. に，本同定法の有効性を確認するために行った，有色雑 音を有するシステムの同定計算シミュレーションの内容， およびその結果について述べる．. (1). E[Δut Δut−k ] = τΔu,k (: known,k = 0,1,...) E[Δyt ] = 0. (2) (3). E[Δyt Δyt−d ] = τΔy,d (: unknown,d = 0,1,...) (4). 提案する同定法. 2. 2.1. E[Δut Δyt−d ] = 0. 有色雑音を有するシステムに対する仮定. (d = 0,1,...). (5). なお，E[∗] は t についての平均を表すものとする．た. 筆者は，白色ガウス入出力雑音を有するシステムに対. とえば十分大きな N に対して，. して，入出力関係から直接得られる誤差分散についての 方程式と，伝達関数パラメータから得られる白色雑音に. E[Δyt ] = (N + 1)−1. 対する白色ゲイン関数を用いた方法を提案し，その同定. N . Δyt. (6). t=0. アルゴリズムを構成した [4]．本研究では，そのパラメー. とする．つまり十分大きなサンプリング回数に対する時. タ同定法を，Fig. 1 のような有色ガウス雑音を有するシ. 間平均を表す．. ステムに対して拡張することを目指す．ただし，本論文. また，本論文の有色ガウス雑音 Δyt ，Δut は，入力が. では，入力雑音は有色ガウス雑音であるが，その分散と. 白色ガウス雑音である未知の線形システムの出力である. 自己共分散は既知であると仮定する．また入力雑音が白. と仮定する．したがってエルゴード性が成り立ち，確率. 色ガウス雑音の場合は自己共分散が 0 となる．. 平均と時間平均が一致すると仮定する． また，t 0 に対して，yt , ut , ηt , ξt , Δyt , Δut はす べて 0 とする．そして，τΔu,k (k = 0,1,...,n − 1) は既知 とし，τΔy,d (d = 0,1,...) は未知とする．さらに yt ，ut ，. Fig. 1. 2.2. ηt ，ξt ，Δyt ，Δut を z 変換したものをそれぞれ y(z)， u(z)，η(z)，ξ(z)，Δy(z)，Δu(z) とする． なお，伝達関数 G(z) の次数 n は既知とし，. The SISO system with colored noises. 対象とするシステム. G(z) =. 提案する同定アルゴリズムは，Fig. 2 に示されるよう. b1 z −1 + ··· + bn z −n 1 + a1 z −1 + ··· + an z −n. (7). とする．このとき，. な，有色入出力雑音を持つ SISO システムに対するもの である．出力雑音の分散と自己共分散は未知とする．シ. η(z) = G(z)ξ(z) b1 z −1 + ··· + bn z −n ξ(z) = 1 + a1 z −1 + ··· + an z −n. ステム全体の入力を ut ，出力を yt とし，t = 1,...,N にお けるそれらの値が計測されているとする．Fig. 2 の G(z) で示される安定な要素に対する入力および出力を ξt ,ηt. となる．また伝達関数パラメータベクトル a,b を，. とする． また，入力誤差として ut に有色ガウス雑音 Δut が加. –2–. (8) (9).
(3) 小松：有色雑音を有する一入力一出力システムに対する雑音ゲイン関数を用いたパラメータ同定法. a = [a1 ,···,an ]T ∈ Rn. (10). b = [b1 ,···,bn ] ∈ R. (11). T. n. 31. まず，. T z = z −1 ,...,z −n. とする．. (15). とすると，(9) 式 より，. さらに，τΔu,k が事前の計測等によって求められ，k > k0. . で τΔu,k の値は十分に小さく無視できるとする． また，G(z) は重極を持たないと仮定する．これによっ. 1 + aT z η(z) − bT zξ(z) = 0. (16). となる．さらに，. て対象とするシステムが限定されると思われるが，一般. η t = [ηt−1 ,···,ηt−n ]T ∈ Rn. に計測値を用いて同定を行う場合，A/D 変換で得られる. T. データには誤差が含まれる．またコンピュータの実数表. ξ t = [ξt−1 ,···,ξt−n ] ∈ R. n. (17) (18). 現は有限長の桁での表現であるため，それらの平均値な とすると，G(z) の入出力の関係は次式のように表される．. どを用いて求められる a1 ,...,an や b1 ,...,bn などのパラ. . メータの推定値は誤差を持つ．そのため，真の伝達関数 が重極を持っていても，ほとんどの場合，推定される伝 仮定しても問題ないと考えられる．. T T 2n+1 q 0,t = [ηt ,η T t ,−ξ t ] ∈ R. 同定アルゴリズム. 本節では，提案する同定計算のアルゴリズムについて， 計算の過程を項目に分け，それらをまとめて示す．なお， 続く節において主な計算項目について述べる．. かる．. (21). さらに，p に直交し原点を通る超平面，. p を， T p = 1,aT ,bT. qT p = 0. (12). ˆ，p ˆ ，b ˆ とし， と す る ．ま た ，a，b，p の 推 定 値 を a τΔy,d (d = 0,1,2,...) の推定値を τˆΔy,d とする．なお， τΔu,k (k = 0,1,2,...,k0 ) は既知とする． このとき，提案する同定計算の各ステップは以下のよ. τˆΔy,d の初期値を τΔu,k とする． ˆ を求める． τˆΔy,d を用いて計測データから p ˆ を用いて計測データから τˆΔy,d を求める． p ˆ の値が収束判定条件を満たすまで繰り S2,S3 を p. を考える．なお，q はこの超平面上の任意の点の位置ベ クトルとする．提案するパラメータ推定は実際に計測さ れる入出力データで構成される q t を用いて (22) 式を満 たす p を推定するものである．. yt = ηt + Δyt. (23). ξt = ut + Δut. (24). である． ここで，. 返す．. ˆ を求める． ˆ から，a, b の推定値 a ˆ, b S5 p なお S4 における収束判定条件は， 回目のループでの ˆ を用いて，ベクトル p を， ベクトル p T ˆT p = 1, p (13). y t = [yt−1 ,···,yt−n ]. T. ut = [ut−1 ,···,ut−n ]. (25). T. Δy t = [Δyt−1 ,···,Δyt−n ]. (26) T. Δut = [Δut−1 ,···,Δut−n ]. とし，ある打ち切り判定係数 μ を決め，次式で与える．. ||p − p−1 || · ||p ||−1 μ. (22). ところで，. うになる．. T. (27) (28). とすると，(23), (24) 式 より，. [yt ,y Tt ,−uTt ]. (14). = [ηt ,η Tt ,−ξ Tt ] + [Δyt ,Δy Tt ,ΔuTt ]. なお S1 において，τΔy,d の初期値として適切な値があ. (29). となる．. ればそれを用いる．ここでは，そのような値がない場合. さらに，. を想定し τΔu,k を用いている．. 2.4. (20). とすると次式が得られ，q 0,t と p は直交することがわ. qT 0,t p = 0. まず，パラメータベクトル a,b より作られるベクトル. S1 S2 S3 S4. (19). そこで，. 達関数は重極を持たなくなる．したがって，このように. 2.3. T ηt ,η T t ,−ξ t p = 0. T T q t = yt ,y T t ,−ut T T Δq t = Δyt ,Δy T t ,Δut. ˆ を求める方法 τˆΔy,d を用いて p. ˆを 同定アルゴリズムの S2 における，τˆΔy,d を用いて p 求める方法について述べる．. とすると，計測値である q t は，. –3–. (30) (31).
(4) 32. システム制御情報学会論文誌 第 34 巻 第 2 号 (2021). q t = q 0,t + Δq t. −1 T T −1 ZN QN +1 QT +1 ZN +1 = W N +1 W. (32). QN +1 QT N +1. となるが，(21) 式より q 0,t は (22) 式の超平面に含まれ. る．そして q t は，雑音 Δq t によってその超平面外にあ. = QN QT N. + q N +1 q T N +1. き，計算の高速化が可能である [3]．. ことが考えられる．. 次のようになる．. ˆ を求める方法は 以上をまとめると，τˆΔy,d を用いて p. まず雑音 Δq t の分布について考える．仮定より Δq t. ・ 計測データ ut ，yt から (30) 式を用いて q t を求める．. の分散行列 S は次式で与えられる．. S=E ただし，. . Δq t Δq T t. ⎡. . Sy 0 = 0 Su. T ZN を求める． ・ q t から (42), (43) 式を用いて ZN. T ZN の最小固有値に対する単位ベクトル n を求 ・ ZN. (33). める． ・ τˆΔy,d ，τΔu,k から (33)∼(35), (37) 式 を用いて W −1 を求める．. ⎤. ··· τΔy,n . ⎥ (n+1)×(n+1) .. . .. ⎥ ⎦∈R ··· τΔy,0 ⎡ ⎤ τΔu,0 ··· τΔu,n−1 ⎢ ⎥ .. .. .. ⎥ ∈ Rn×n Su = ⎢ . . . ⎣ ⎦ τΔu,n−1 ··· τΔu,0 ⎢ Sy = ⎢ ⎣. τΔy,0 .. . τΔy,n. (43). T となり，N に対して ZN ZN を逐次的に求めることがで. ると考えられる．そこで，その超平面から q t までの距離 ˆ とする et を考え，e2t の総和を最小とする p を，推定値 p. (42). ・ W −1 ,n から (39) 式を用いて w1 ,w2 ,···,w2n+1 を求 める．. (34). ˆ を求める． ・ w1 ,w2 ,···,w2n+1 から (40) 式を用いて p. 2.5. ˆ を用いて τˆΔy,d を求める方法 p. ˆ を用いて計測デー 同定アルゴリズムの S3 における，p. (35). タから τˆΔy,d を求める方法について述べる． まず (9) 式より，. そこで p と直交する超平面上の任意のベクトル q を用. y(z) − Δy(z) = G(z)(u(z) + Δu(z)). (44). いて，距離 et を次のように定義する．. e2t = min(q t − q)T W −2 (q t − q) q. である．そこで，. (36). φ(z) = y(z) − G(z)u(z). ただし，. W=. √. (45). とおく（なお φ(z) を φt の z 変換とする）．すると (44). S. 式は，. (37). G(z)Δu(z) + Δy(z) = φ(z). である． すると，付録 1 に示すように，すでに提案された同定 方法 [1–5] で用いられたパラメータ推定方法と同様に，. T. れる計算手順は次のようなものである．. φt = −aT φt−1:t−n + 1,aT y t:t−n − bT ut−1:t−n. (38). (48). T より ZN を作り，ZN ZN の最小固有値に対する単位固有. となる．さらに Δη(z) を，. ベクトル n を求める．つぎに，W −1 を用いて， T. [w1 ,w2 ,···,w2n+1 ] = W −1 n. Δη(z) = G(z)Δu(z). (39). Δη(z) + Δy(z) = φ(z). (40). (50). となる．ここで，Δη(z) が Δηt の z 変換であるとすると，. ˆp ˆ , b, ˆ を求める． より，a. Δηt + Δyt = φt. なお，. QN = [q 1 ,...,q N ]. (49). とすると，(46), (49) 式より，. より，w1 ,w2 ···w2n+1 が求められ，. ˆT = w−1 [w1 ,···,w2n+1 ]T ˆ = 1, a ˆ T,b p 1. (47). の表記を導入すると，(45) 式より，. まず計測データ q t と W −1 を用いて，. ZN = [q 1 ,...,q N ]T W −1 ∈ RN ×(2n+1). となる．ここで，一般に整数 k,(k < ) に対して，. φ:k = [φ ,φ−1 ,···,φk ] ∈ R−k+1. パラメータ推定を行うことが可能となる．その結果得ら. (46). (41). (51). となる． さらに，. とすると (38) 式より，. τΔη,d = E[Δηt Δηt−d ] –4–. (52).
(5) 小松：有色雑音を有する一入力一出力システムに対する雑音ゲイン関数を用いたパラメータ同定法. τφ,d = E[φt φt−d ] τ Δy = [τΔy,0 ,τΔy,1 ,···,τΔy,n ]. τ Δu = [τΔu,0 ,τΔu,1 ,...,τΔu,k0 ] τ Δη = [τΔη,0 ,τΔη,1 ,...,τΔη,n ] τ φ = [τφ,0 ,τφ,1 ,...,τφ,n ]. から (48) 式を用いて φt ，τ φ を求める．. (53) T T. T. T. (54). ˆ から (59) 式の vr ,zr を求める． ・p. (55). ・ vr ,zr から (61), (62) 式を用いて fd,k (G) を求める．. (56). ・ τ Δu ，τ φ ，F (G) より (65) 式を用いて τˆΔy,d ，τˆΔy,d を求める．. (57). ここでは (60) 式を満たす fd,k (G) が (61), (62) 式 で. となるので，(51) 式より，. 与えられることを示す． まず， = 1,2,... に対して，. (58). となる．. w =. 仮定より G(z) は重極を持たない．そこで G(z) を，. G(z) =. n vr z −1 1 − zr z −1 r=1. fd,k (G)τΔu,k. G(z) =. (60). n . vr z −1. ∞ . zr z −1. . (67) (68) (69). =1. となる．したがって (49) 式より，. と表す．ただし，次節で述べるように，. Δη(z) =. (61). ∞ . w z −. ∞ . Δut z −t. (70). t=0. =1. となる．また t 0 で Δut = 0 であるので，. であり，k = 1,2,... に対して， n n vr vs fd,k (G) = (zrk+d + zr|k−d| ) 1 − z z r s r=1 s=1. (66). r=1 n =0 ∞ vr zr−1 z − = r=1 =1 ∞ w z − =. k=0. n n vr vs d fd,0 (G) = z 1 − zr zs r r=1 s=1. vr zr−1. とすると，(59) 式より，. (59). |zr | < 1(r = 1,2,...,n) である． つぎに，τΔη,d を雑音ゲイン関数 fd,k (G) を用いて， ∞ . n r=1. とすることができる．ただし，G(z) は安定であるので. τΔη,d =. fd,k (G) の導出. 2.6. とすると，(5) 式および (49) 式より Δηt と Δyt は無相関. τ Δη + τ Δy = τ φ. 33. Δη(z) =. (62). ∞ ∞ . w Δut− z −t. (71). t=0 =1. となる．よって，. である． ところで，k > k0 に対して τΔu,k は無視できると仮定 しているので，fd,k (G)(k = 0,...,k0 ) を用いて F (G) を，. ⎤ f0,0 (G) ... f0,k0 (G) ⎥ ⎢ .. .. .. ⎥ F (G) = ⎢ . . . ⎦ ⎣ fn,0 (G) ... fn,k0 (G). Δηt =. ⎡. w Δut−. (72). =1. が得られる．これより，. (63) Δη t Δηt−d ∞ ∞ = w Δut− w Δut−d−. とすると，τ Δη の推定値 τˆ Δη を既知の τ Δu を用いて，. τˆ Δη = F (G)τ Δu. ∞ . (64). =. と表すことができる． したがって (58), (64) 式より τ Δy の推定値 τˆ Δy は，. +. =1 ∞ ∞ p=0 =1 ∞ ∞ . (73). =1. w w+p Δut− Δut−d−−p w+p w Δut−−p Δut−d−. (74). p=1 =1. τˆ Δy = τ φ − F (G)τ Δu. (65). となる．したがって，. より求められる．また，計測データ ut ，yt から (48) 式を 用いて φt ，τ φ が求められるので，(65) 式を用いて τˆ Δy ，. τˆΔy,d が求められる．. =. ˆ を用いて τˆΔy,d を求める方法は 以上をまとめると，p 次のようになる．. +. ・ 計測データ ut から τ Δu を求め，計測データ ut ，yt. τΔη,d ∞ ∞ p=0 =1 ∞ ∞ p=1. –5–. =1. w w+p τΔu,d+p w w+p τΔu,d−p. (75).
(6) 34. システム制御情報学会論文誌 第 34 巻 第 2 号 (2021). となる． つぎに，. Ar,s =. vr vs 1 − zr zs. (76). とすると，p = 0,1,2,... に対して，. . w w+p = ∞ . w w+p =. n . . vs zs−1. s=1 n n . n . vr zr+p−1. (77). r=1. Ar,s zrp. (78). r=1 s=1. =1. となる．したがって (75), (78) 式より，. = + = +. τΔη,d n n ∞ p=0 r=1 s=1 n n ∞ p=1 r=1 s=1 n n ∞ . Fig. 3. Estimated parameters of G(z). Ar,s zrp. τΔu,d+p . Ar,s zrp. τΔu,d−p. (79). Ar,s zrk−d τΔu,k. k=d r=1 ns=1 n ∞ . Ar,s zrk+d τΔu,−k. k=1−d r=1 s=1 n n d Ar,s zr τΔu,0 = r=1 s=1 n n ∞ k+d |k−d| + Ar,s zr + zr τΔu,k k=1 r=1 s=1. (80). Fig. 4. (81). Estimated variances. とした．そのときの Δu の分散と自己共分散は，シミュ レーション結果より，. となる．したがって (76), (81) 式より (61), (62) 式が得 られる．なお τΔu,k = 0(k = 1,2,...)，d = 0 とすると (81). 2 σΔu = τΔu,0 = 0.0038919. (87). 式は既発表の同定法 [4] で導出された白色雑音ゲイン関. τΔu,1 = −0.0022512. (88). 数に一致する．ただし，本論文での導出の方がより簡易. τΔu,2 = 0.0012550. (89). な方法となっている．. τΔu,3 = −0.0002976. (90). τΔu,4 = 0.0000189. (91). シミュレーション結果. 3. 3.1. となる．これより (63) 式の k0 を k0 = 3 とし，τΔu,k. 入力雑音の分散と自己共分散を既知とした 場合. (k = 0,...,3) の値を既知とした． なお，本来 τΔy,0 ,τΔu,0 の値は正であるが，計測データ. 伝達関数およびノイズのパラメータを仮定し，パラ メータ同定を行うシミュレーションを実施した．. が少ない部分では，雑音の影響のためその値が負になっ. システムは Fig. 1 において， −1. ている部分がある．そのような場合，重み行列 W に複 素数が含まれる．ただし，本シミュレーションプログラ. −2. 1.0z + 0.5z 1 − 1.5z −1 + 0.7z −2 1 − 0.2z −1 + 0.3z −2 Hu (z) = 1 + 0.4z −1 + 0.1z −2 1 + 0.3z −1 + 0.4z −2 Hy (z) = 1 + 0.1z −1 + 0.2z −2 G(z) =. (82). ムでは，転置計算を共役複素転置とすることなどによっ て計算を継続させ，繰り返し計算を収束させることがで. (83). きている． また，入力は ±1 の値をとるランダム信号とし，入出. (84). 力データをシミュレーション計算で発生させ，データ数. としたものである．また，Fig. 1 の u ， y は，白色正規. 100 ごとにその時点でのデータ全体（総数は増加する） q t (t = 1,...,N ) を用いて，パラメータ同定を行った．N に対する伝達関数パラメータ推定値の変化を Fig. 3 に示 す．また，ノイズ分散推定値の変化を Fig. 4 に示す．. 雑音から生成し，それらの分散 σ2u ，σ2y の値は，. σ2u = 0.0025. (σu = 0.05). (85). = 0.0025. (σy = 0.05). (86). σ2y. –6–.
(7) 小松：有色雑音を有する一入力一出力システムに対する雑音ゲイン関数を用いたパラメータ同定法. Table 1. Results with simulation. pG ,σ 2. T.V.. E.V.. E.V.-TV. a1 a2 b1 b2. −1.5 0.7 1.0 0.5. −1.4974 0.6976 1.0067 0.5056. 0.0026 −0.0024 0.0067 0.0056. 0.00269 0.00056 0.00036. 0.00297 0.00082 0.00059. 0.00029 0.00026 0.00022. 2 σΔy τΔy,1 τΔy,2. 35. T.V.: True values,. Fig. 5. A motor control system of a mobile robot. を計測することで，u から y への開ループシステムに対. E.V.: Estimated values. する同定を行う場合では，本同定法の雑音についての仮 定は妥当であると考えられる．. そして N = 10000 において得られた，伝達関数パラ. ところで本論文で導出された雑音ゲイン関数は，白色. メータ推定値およびノイズ分散推定値を Table 1 に示す．. Fig. 3 より各伝達関数パラメータがほぼ正しく推定さ れていることがわかる．また Fig. 4 よりノイズの分散も 正しく推定されていることがわかる．なお Table 1 より a1 ,a2 の推定誤差は 0.003 程度であり，b1 ,b2 の推定誤差 2 ,τΔy,k の は 0.007 程度であることがわかる．さらに σΔy 残差は 0.0003 程度であることがわかる．. 4.. 入力雑音を仮定した場合，つまり τu,k = 0(k = 1,2,...) の場合，ARMA モデルに対する自己共分散関数に対応 する．そしてそのような ARMA モデルに対する自己共 分散の値についての計算法が発表されている [10]．しか し，有色入力雑音に対するは示されておらず，導出され ている式は ARMA モデルの AR 部分を用いた多項式方 程式の根，つまり，本論文でのパラメータ zr に対応する パラメータと，MA 部分のパラメータ，つまり，本論文. 考察. 本同定法では，入力雑音が制御装置からの信号であり，. の b1 ,b2 ,... に対応するパラメータによって表現されてお り，本論文のようにパラメータ vr ,zr を用いた形式とは. あらかじめ計測できるが，フィードバックを行うための. 異なっている．また，本論文で導出した式の方がより簡. 出力信号の計測装置には未知の雑音が入るようなシステ. 素に表現されていると考えられる．. ムを想定し，入力雑音の分散を既知として有色ゲイン関. なお，本同定法の繰り返し計算アルゴリズムにおける. 数を用いてノイズパラメータを推定している．そこで，. 収束性については，入力雑音がないシステムの場合 [5]. そのような仮定が実際のシステムにおいて妥当であるか. と同様の検証が必要となる．. 考察する．. 5.. 例として，モータ回転角や回転速度をセンサによって. 結論. 電圧として検出され，制御装置の制御電圧がモータの駆. 今回，入力雑音の分散と自己共分散が既知であると仮. 動回路へ入力されるモータ制御装置を考える．このとき. 定した場合の，有色の入出力雑音を有する SISO システ. 入力部分がシールドされた機器内部にあるとすると，上. ムに対するシステム同定アルゴリズムを提案した．そし. 記の仮定が成り立つと考えられる．. て，一入力一出力システムに対するシミュレーションを. このようなシステムの構成の具体的な例を Fig. 5 に示. 行い，パラメータ推定値が得られることを確認した．. す．これは移動ロボットにおけるモータ制御システムで. 本同定法は有色出力雑音を考慮しており，一般的な雑. ある．MPU がモータ出力の回転速度を検出し，さらに. 音環境に適した同定法となっている．また入力雑音が白. その値に応じてモータドライバを通してモータを駆動す. 色でその分散が既知のシステムにも適用可能である．. る．MPU，D/A コンバータ，A/D コンバータ，モータ. 今後実験システムへ適用し，本同定法の有効性につい. ドライバは電磁シールド内にある．電磁シールドされた. て検証を行う必要があると考えられる．また，既発表の. ロボット内部のモータドライバ入力への電磁的雑音につ. 多入力多出力システムの場合の同定法 [3] と同様に，本. いては，あらかじめ信号の雑音成分を計測することで雑. 同定法を多入力多出力システムにも利用できると考えら. 音の統計量を求めることが可能であると考えられる．な. れる．. お，モータ本体や出力部には多くの場合減速ギアがつい. なお，フィードバックシステムに本同定法を適応する. ており，それらの軸受に働く摩擦力は通常極めて小さい．. 場合，事前にオープンループでの入出力値を計測し同定. また減速ギアの出力は十分なトルク出力を持つため，摩. 計算を行い，その結果を用いてフィードバックシステム. 擦力の影響は大きくないと考えられる．一方，センサ信. を構築する．ただし，制御を進めていく過程で，十分な. 号にはロボットが移動する環境に応じた外部の雑音が入. 入出力データが得られるならば，その時点で同定計算を. り，それらの分散などは通常未知となる．したがって，. 行い，制御対象の伝達関数を更新することにより，フィー. Fig. 5 において，MPU によって u を適当に発生させ，y –7–.
(8) 36. システム制御情報学会論文誌 第 34 巻 第 2 号 (2021). ドバックゲインの更新を行い，制御性能を向上させるシ. が得られる．したがって，そのとき，. 2 e2t = q Tt W −1 n Wp n= ||W p||. ステムを構築する閉ループ同定への応用も考えられる．. 参 考 文 献 [1] 小松：入出力誤差分散推定を用いた一入力一出力システ. T. Z = [q 1 ,q 2 ,···,q N ] W −1. N . e2t = ||Zn||2. (A6). t=1. となる．この総和を最小とする n を求めるため，ラグラ ンジュの未定数法を用いる．つまり，. ムに対するパラメータ同定法；システム制御情報学会論 文誌, Vol. 31, No. 9, pp. 336–345 (2018). J2 = (Zn)T (Zn) + λ(1 − nT n). [5] 小松：有色出力雑音を有する一入力一出力システムに対. (A7). を考え，J2 を未定定数 λ と n で偏微分し，それらを 0 と. する固有ベクトル法を用いたパラメータ同定法；シス テム制御情報学会論文誌, Vol. 33, No. 3, pp. 95–103. おくことにより，. (2020). Z T Zn = λn. [6] 和田：固有ベクトル法によるシステム同定について；シ ステム／制御／情報, Vol. 51, No. 1, pp. 34–36 (2007) [7] T. S¨ oderstr¨ om: Errors-in-varialbes models in system identification; Automatica, Vol. 43, pp. 939–958 (2007) [8] S. Beghelli, R.P. Guidorzi and U. Soverini: The Frisch scheme in dynamic system identification; Automatica, Vol. 26, pp. 171–176 (1990) [9] R. Diversi, R. Guidorzi and U. Soverini: A new criterion in EIV identification and filtering applications; 13th IFAC Symposium on System Identification (2003) [10] M. Karanasos: A new method for obtaining the autocovariance of an ARMA model: An exact form solution; Econometric Theory, Vol. 14, No. 5, pp. 622–640 (1998). (A8). を得る．したがって n は Z T Z の単位固有ベクトルとな り，J2 の最小値は Z T Z の最小固有値で与えられる．ま た，そのとき (A4) 式より，n と定数 k を用いて，. W p = kn. (A9). となる．これより，. p = kW −1 n. (A10). が得られる．ここで，. W −1 n = [w1 ,···,w2n+1 ]. T. (A11). とすると，. ˆT = w−1 [w1 ,···,w2n+1 ]T ˆ = 1, a ˆ T,b p 1. 録. (A12). ˆp ˆ , b, ˆ が求められる． より，a. 付録 1. パラメータ推定方法について ここでは，すでに提案された同定方法 [1–5] で用いら. 著 者 略 歴. れた固有ベクトル法によるパラメータ推定法 [6] につい. こ. まつ. のぶ. お. 小 松 信 雄 (正会員) 1983 年 3 月大阪府立大学大学院工学研 究科修士課程修了（航空工学）．1988 年 11 月大阪府立大学工学部航空宇宙工学科 助手，1996 年 4 月静岡理工科大学理工学 部機械工学科講師，1999 年 4 月同助教授， 2002 年 4 月大阪工業大学情報科学部情報 科学科助教授（2007 年より准教授）となり現在に至る．移動. て述べる． まず (36) 式の e2t を求めるため，(22) 式を拘束条件と して，ラグランジュの未定定数法を用いる．つまり，. (A1). を考え，J1 を未定定数 λ と q で偏微分し，それらを 0 と おくことにより (36) 式の q を求める．すると，. qT p qt − q = T t 2 W 2 p p W p. (A5). とすると，e2t の総和は，. Vol. 28, No. 7, pp. 310–319 (2015) [3] 小松：入出力誤差分散推定を用いた多入力多出力システ ムのパラメータ同定法；システム制御情報学会論文誌, Vol. 29, No. 8, pp. 362–371 (2016) [4] 小松：白色雑音ゲイン関数を用いた一入力一出力システ. J1 = (q t − q)T W −2 (q t − q) + λq T p. (A4). となる．さらに. ムのパラメータ同定アルゴリズム；システム制御情報学 会論文誌, Vol. 27, No. 5, pp. 216–229 (2014) [2] 小松：入出力誤差分散推定を用いた多入力一出力システ ムのパラメータ同定法；システム制御情報学会論文誌,. 付. (A3). 体の位置計測，誘導制御に関する研究に従事．博士（工学） 大阪府立大学．. (A2). –8–.
(9)
図
関連したドキュメント
The aqueous layer was extracted with ethyl acetate, and the combined organic layers were washed with brine, dried over anhydrous MgSO 4 , and concentrat- ed in vacuo.. The solution
One dimensional classification problem is used for simulation to show the validity of adding one randomly selected data to a pair of the boundary data.. The location of the boundary
This paper presents a new wavelet interpolation Galerkin method for the numerical simulation of MEMS devices under the effect of squeeze film damping.. Both trial and weight
３．排出水に対する規制
8, and Peng and Yao 9, 10 introduced some iterative schemes for finding a common element of the set of solutions of the mixed equilibrium problem 1.4 and the set of common fixed
高出力、高トルク、クリーン排気を追求した排ガ ス対応エンジンは、オフロード法 2014 年基準に 適合する低エミッション性能を実現。また超低騒
During stage 1, we used an adaptively preconditioned thick restarted FOM method to approximately solve the linear system and then used recycled spectral information gathered during
These recent studies have been focused on stabilization of the lowest equal-order finite element pair P 1 − P 1 or Q 1 − Q 1 , the bilinear function pair using the pressure
