選炭ジグのシステム同定 尾 崎 正 憲 （Masanori Ozaki） (指導 : 太屋岡篤憲) 本稿では、現在，オペレータの勘や経験により行っている選炭ジグ装置の制御を、計算機によるによる自動制御 で代行させるための選炭ジグ装置のモデリングについて検討した選炭ジグ装置の伝達関数は一次遅れ＋無駄時間で 表されるが、本研究では、一次遅れ系を離散時間系に置き換え、システム同定を行なった。まず同定入力となる M 系列の生成し、生成した M 系列が白色性を持っているか確認のために自己相関関数とパワースペクトル密度を求 めるプログラムを作成した。また、一次遅れ系のパラメータ推定を IV 法、LS 法の二つの方法で同定するプログラ ムを作成し、さらに、雑音に対して IV 法が有効性であることを確認した。 キーワード：選炭ジグ装置、Ｍ系列、最小二乗推定、補助変数法、モデルの評価 最小二乗推定法 (LS 法) や補助変数法 (IV 法) を用い 1. は じ め に てパラメータの推定を行なう。 STEP5 モデルの評価 古くから石炭の選別装置は、石炭と不純物が混ざった原 得られたモデルの妥当性の評価を行なう。 炭を液体中に入れ、水流を用いて上下運動させることで、 石炭と不純物の比重の違いを利用して選別する”選炭ジグ装 置”が用いられてきた。しかし、装置の特性は不明な点が多 く、その制御は熟練オペレータの経験に頼って行なってい るため、製品品位を一定に保つのが困難である。 本研究は、熟練オペレータが担っていた経験による制御 を、機械に代行させる制御系を構築するため、選炭ジグ装 置のシステムを同定することを目的とする。 2. システム同定 〈2・1〉 システム同定とは システム同定とは、対象とする動的システムの入出力デー タから、ある目的の下で、対象と同一であることを証明で 図 1 選炭ジグ装置 きるような数学モデルを作成することである。図 1 に、本 研究の同定対象である選炭ジグ装置を示す。この装置にお ける操作入力は装置下部に付けられた硬抜き装置への指令 信号で、制御量は水槽中に入れられた石炭と不純物の間の 比重を有するフロートの高さである。 〈2・2〉 システム同定の手順 図 3 に，システム同定 のフローチャートを示す。以下，各ステップの内容を説明 する。 STEP1 実験の計画 ハードウェア、同定入力（本研究の場合のＭ系列）や サンプリング周期などの設定を行なう。 STEP2 データの収集 同定実験によって入出力データの採取を行い、その際、 前処理として異常値の除去を行なう。 STEP3 モデル構造の決定 モデルの形や、パラメトリックモデルを選んだ場合は 図2 次数の決定を行なう。 STEP4 システム同定 平成 18 年度 卒業研究論文集 49 システム同定フローチャート 〈2・3〉 同定入力の選定 〈2・4〉 可同定性 〈2・3・1〉 同定入力の周波数特性 可同定性とは、パラメータの推定値が求められるための システム同定を行うためには、同定入力は対象のもつ全 条件を満たしていることである。 てのモードを励起できるような広い範囲の周波数成分を含 ・可同定性の条件 （ 1 ） 入力に関する条件 む必要がある。 〈2・3・2〉 M 系列について （ a ） 入力 u(k) が有界であること。発散してはな システム同定を行なう際、システムの性質を調べるため らない。 に入力としては，すべての周波数成分を持つ白色雑音が望 （ b ） 入力 u(k) が持続的に励振していること。 ましいが、白色雑音は，物理的に実現が困難である。本研 （ 2 ） 出力に関する条件 究では、同定入力として，ほぼ白色性を有する M 系列を用 （ a ） 出力 y(k) が有界であること。すなわち、シ いる。M 系列信号は、図 3 に示すような、2 値（1/0 また ステムが安定であって、伝達関数の極が単位 は 1/-1）の周期信号である。M 系列は、人為的にある規則 円内に存在していること。 に基づいた不規則信号を生成することで、様々の乱数を生 （ b ） 出力が全て 0 になることはない。すなわ 成でき、同定の入力に適している。 ち、伝達関数の零点がすべて 0 になってはい 〈2・3・3〉 M 系列の生成 けない。 周期 T = 2p − 1 の M 系列は式 (1) で生成される。 Xn = Xn−p ⊕ Xn−q (p > q)· · · · · · · · · · · · · (1) 〈2・5〉 一次遅れ系のシステム同定 ただし、1 を真、0 を偽とする排他的論理和である。また、 〈2・5・1〉 実際のシステム 初期値として xn (n = 0, 1, 2,・ ・ ・, k) にはすべて０でなけれ 実際の入出力データを観測するシステム、すなわち一次 ば何を与えてもよい。 遅れ系の離散時間モデルを状態変数を用いて表すと式 (2) 〈2・3・4〉 M 系列の性質 となり、これを使って実際の入出力データである観測デー 以下に，M 系列の性質を述べる． タを作成した。 と「1」の数の割合が一定で、n 段のシフトを使った x(k + 1) = −ax(k) + bu(k) + w(k) とすれば、1 の数が２ n-1 個、0 の数が (２ n-1) − y(k) = x(k) + v(k) （ 1 ） どんな周期の M 系列でも 1 周期中に含まれる「0」 · · · · · · (2) 1 になる。すなわち 0 の数がひとつだけ少ない。こ 式 (2) より，状態変数 x を消去すると式 (3) となる。 れはレジスタの状態としてオールゼロを除いたこと y(k) = −ay(k − 1) + bu(k − 1) + e(k)· · · · · · · · (3) から容易に推測できる。0 と 1 の数のバランスがと れているので 0 を-1 に置き換える PN 系列では直流 ただし、 分がきわめて少なくなって好都合になる。 e(k) = v(k) + av(k − 1) + w(k)· · · · · · · · · · · · · · (4) （ 2 ） 系列中に 0 または 1 が連続して現れるとき、その 連続の長さを run と呼ぶ。系列 1 周期中の長さ m x:状態変数 u(k):入力 (選炭ジグ装置のズリの排出量) y(k):出力 (フロートの高さ) v(k):出力側の観測雑音 w(k):入力側の観測雑音 v(k)v 、w(k) はお互いに無相関な平均値 0、分散 r,q の の run の発生頻度は、長さ m ＋ 1 の run の頻度の ちょうど 2 倍になっている。これは（擬似乱数では ない）真の乱数がもつ性質と同じである。 （ 3 ） 発生されたある M 系列と、レジスタの初期値だけ が違う同じ M 系列のデータを加算すると、もとの M 系列を時間的にずらせただけの同一の系列ができ 雑音 る。ただしここでいう加算とは桁上げのない特別な a,b:推定するパラメータ 加算で、排他的論理和と同じと考えておけばよい。 また、出力結果を得るための入力データは、白色性をも つ M 系列信号 (1 or -1) を用いる。 〈2・5・2〉 一次遅れ系の離散モデル 一次遅れ系の伝達関数 (G(s) = K 1+T s ) の離散モデルは， 以下の式 (5) で表される。 図3 y(kTs ) = −ay(kTs − Ts ) + bu(kTs − Ts )· · · · · (5) M 系列の例 ただし、入力:u(kTs ), 出力:y(kTs ) , 推定するパラメータ:a,b 50 選炭ジグのシステム同定 〈2・5・4〉 一次遅れ系の補助変数 (IV 法) また、一次遅れ系の離散時間における伝達関数は式 (6) 予測誤差と過去のデータに相関があると LS 法では真値 で示せる。 G(z) = に一致せずバイアスを持ってしまう。IV 法はこの LS 法の −1 Y (z) bz b = = · · · · · · · · · (6) −1 U (z) 1 + az z+a 問題点を改良してもので、補助変数 ζ(k) を導入すること で、回帰ベクトルとは相関を持ち、雑音とは無相関にして 式 (5) より離散型の入出力データからパラメータ a,b を 最小二乗推定を行なうことができる。 決定する。 補助変数 ζ(k) を式 (11)、パラメータベクトル θ を式 (12) とすると、 〈2・5・3〉 一次遅れ系の最小二乗法 (LS 法) LS 法とは測定で得られた数値の組を、適当なモデルか ζ(k) = (−u(k − 2) u(k − 1))T · · · · · · · · · · · · (11) ら想定される特定の関数を用いて近似するときに、想定す る関数が測定値に対してよい近似となるように、残差の二 乗和を最小とするような係数を決定する方法である。本研 θ = (a1 · · · an b1 · · · bm )T · · · · · · · · · · · · · · · · · (12) 究の場合、この LS 法により，実際の出力とパラメータ推 定した値での出力結果において、両者の２乗誤差を最小に これを適用すると、 するようにパラメータを推定する同定することが可能とな る。式 (3) の一次遅れの離散時間系において、パラメータ Ã ベクトル θ 、回帰ベクトル φ(t) とすると、 θ̂ = θ = (a1 · · · an b1 · · · bm ) · · · · · · · · · · · · · · · · · · ·(7) T !−1 N 1 X ζ(k)(−y(k − 1)u(k − 1)) N k=1 × φ(k) = (−y(k − 1) u(k − 1))T · · · · · · · · · · · · · (8) k=1 ˆ ) は式 (9) のように と定義する。また、LS 法の推定値 θ(N よって、 定義される。 Ã θ̂ = !−1 N 1 X φ(t)(−y(k − 1)u(k − 1)) N Ã k=1 × N 1 X φ(t)y(k)· · · · (9) N â b̂ ! よって式 (9) に式 (8) を代入すると式 (10) に一次遅れ系の 最小二乗推定値式が得られる。 â b̂ ! 1 X 2 y (k − 1) N k=1 N 1 X − u(k − 1)y(k − 1) N k=1 −1 N 1 X y(k − 1)u(k − 1) − N k=1 N X 1 2 u (k − 1) N k=1 N 1 X y(k − 1)u(k) − N k=1 · · · · (10) × N 1 X u(k − 1)u(k) N N 1 X u(k − 2)y(k − 1) N k=1 N 1 X − u(k − 1)y(k − 1) N k=1 −1 N 1 X − u(k − 2)u(k − 1) N k=1 N X 1 u2 (k − 1) N k=1 N 1 X u(k − 2)y(k) − N k=1 · · · · (14) × N 1 X u(k − 1)y(k) N = k=1 Ã N 1 X ζ(k)y(k)· · (13) N N = k=1 〈2・6〉 モデルの妥当性の評価 システム同定では、同定結果であるとパラメータによる 出力が、システムの出力とどの程度一致しているか調べる 必要がある。 方法 1 同定値によるステップ応答またはインパルス応答 での比較。 k=1 方法 2 適合率 (%) を用いて数値で比較。適合率は式 (15) に示す。 方法 3 周波数領域での評価 (ボード線図、ナイキスト線図) 平成 18 年度 卒業研究論文集 51 v uN uX t (ŷ(k) − y(k))2 k=1 × 100 (15) v Fit(%) = 1 − uN X u 2 t (y(k) − ȳ(k)) て示す。 k=1 ただし、ȳ = N 1 X y(k) N k=1 3. 結 果 〈3・1〉 同定入力である M 系列の生成 本研究では、M 系列の生成と M 系列に様々な乱数が含 まれているかをチェックするため、図 4 のような M 系列、 図 5 プログラム実行画面 (IV 法) 自己相関関数、パワースペクトル密度関数を求めるプログ ラムを VisualBasic.NET で作成した。図 4 の (a) は M 系 列、(b) は自己相関関数、(c) はパワースペクトル密度関数 である。図 4(b) より作成した M 系列は、周期間では相関 がないことがわかる。図 4(c) のパワースペクトル密度は一 定にはなってないが、M 系列は完璧な白色性を持つもので 表1 はないため、ばらつきはあるが全体的にみるとほぼ白色性 真値 をもっていることがわかる。 同定結果 適合率 [%] a(200) b(200) -0.9 2 − LS 法 雑音小 r = 0 -0.881621380 1.989705531 91.004455 IV 法 雑音小 r = 0 雑音中 r = 2 雑音大 r = 4 -0.895318604 1.993268131 97.455349 -0.942523507 2.025839958 73.074441 -0.992557430 2.066773434 28.928566 〈3・2・1〉 最小二乗推定法による推定結果 まず、本研究では、式 (10) で表される LS 法を用いて、 一次遅れ系のシステム同定を行なった。図 6 に、出力は示 す。表 1 の推定結果より、a,b ともに真値に対して正確な 結果が得られていることがわかる。また、適合率を見ても 出力に対してパラメータ推定結果を用いた波形は 91% の適 合を示している。 図 4 研究結果 (M 系列) 〈3・2〉 一次遅れ系におけるパラメータ推定 本研究では、VisualBasic.NET を用いて生成した M 系 列を同定入力としてシステムの出力を抽出し、これらの入 出力データをもとに、LS 法と IV 法の２種類のシステム同 定法を用いてパラメータ推定を行なった。図 5 に、プログ ラムの実行画面は示す。また、推定したパラメータを用い た出力とシステムの出力を同時に示すことで、２つの出力 図6 波形を比較した。表 1 に、各同定法での同定結果はまとめ 52 LS 法 (雑音小 r = 0) 選炭ジグのシステム同定 〈3・2・2〉 補助変数法での同定結果 次に、本研究では、式 (14) で表される IV 法を用いて一 次遅れ系のシステム同定を行なった。雑音小 r = 0 の場合、 システムの出力とパラメータ推定結果を用いた波形は図 7 に示し、雑音中 r = 2 の場合を、図 9 に示す。さらに、雑 音小 r = 4 の場合は図 10 に示した。表 1 の推定結果より、 雑音小では、a,b ともに真値に対して極めて正確な結果が 得られていることがわかる。また、適合率を見ても出力に 対してパラメータ推定結果を用いた波形は 97% の適合を示 している。これらの結果より、本研究の場合、同定法とし ては LS 法より IV 法が適していることがわかる。また雑音 図 9 IV 法 (雑音大 r = 4) r = 4 の場合は、推定したパラメータや適合率は共に真値 から離れているが、極めて大きな雑音を考慮すれば、今回 のプログラムはある程度の雑音に対しても同定することが 4. 同定結果の評価 できる。 次、推定したパラメータで決定した出力波形とシステム の出力波形を比較検討する。また最小二乗推定法と補助変 数法の同定結果の比較検討も行なう。 〈4・1〉 最小二乗推定法と補助変数法の比較 ここではまず、最小二乗推定法と補助変数法の比較を行 なう。雑音が小さい場合、表 1 における適合率を見ると最 小二乗推定法における適合率は 91%、これに対して、補助 変数法の適合率は 97%を示していることがわかる。適合率 とは推定パラメータを用いて出力した波形がどの程度シス テムの出力と適合しているものか判断する要素のため、こ の結果から見ても補助変数法がより優れていることがわか る。最小二乗推定法は、残差の白色化を目的とたシステム 同定法に対して、補助変数法は残差の無相関化を目的とし 図7 IV 法 (雑音小 r = 0) ている。そのため、予測誤差と過去のデータに相関がある と最小二乗推定法では真値に一致せずバイアスを持ってし まう。そこに、補助変数 ζ(k) を加えることで回帰ベクトル とは相関を持ち、雑音とは無相関になるようにパラメータ 推定を行なうため、同じ条件下で同定を行なうと当研究の ような結果が出たと考えられる。 〈4・2〉 各同定結果のステップ応答比較 図 10 に、雑音が小さい (r = 0) という条件下において、 最小二乗推定法と補助変数法の二つの同定法で推定したパ ラメータを用いて、同定値によるステップ応答を求めた。 また、図 11 では雑音を変化させた場合の IV 法での同定値 によるステップ応答を求めた。ただし、雑音大 r=4 の場合 は、真値を大きく外れるため対象外とした。図 10 より、真 値のステップ応答に対して補助変数法を用いた同定値を利 用したほうが近い値で収束していることがわかる。この結 図8 果や表 1 に示す推定パラメータと適合率比較からも、本研 IV 法 (雑音中 r = 2) 究については、補助変数法で推定したパラメータを用いた モデルが理想的なモデリングの構築が可能であると考えら れる。また、図 11 より、大きな雑音を含むほど同定が困難 になることも示せた。 平成 18 年度 卒業研究論文集 53 申しあげます。 文 献 （ 1 ） システム同定：http://www1.odn.ne.jp/t-takig/study/controlref/sys id.pdf （ 2 ） 片山徹：システム同定入門、朝倉書店 （ 3 ） 相良節夫、秋月影雄、中沢高好 片山徹：システム同定、コロナ社 （ 4 ） 小森大輔：これからはじめる VisualBasic.NET 秀和システム （ 5 ） 河野勉：VisualBasic.NET コントロール関数編 技術評論社 （ 6 ） 若山芳三朗：学生のための VisualBasic.NET 東京電機大学出 版局 図 10 図 11 同定値によるステップ応答比較 同定値によるステップ応答比較 (IV 法) 5. 終 わ り に 本研究では、IV 法を用いて大きな雑音を含む一次遅れ系 の同定を行えた。しかし、IV 法は、N 個のデータが追加観 測されると改めて全体の計算を行なわなければならない。 今後は、効率よい同定を行なうため新しいデータが追加さ れるごとに前の推定値を修正していく、逐次補助変数法を 用いて、推定するパラメータの様子をステップごとに確認 できる同定を検討する。また、周波数応答 (ボード線図や ナイキスト線図) を用いたモデルの評価や、無駄時間を含 めた同定を行なうことも課題である。 6. 謝 辞 本研究を進めるにあたり、様々な方から有益なご助言や 励ましのお言葉を数多く頂きました。ここに深く謝意を表 します。特に５年次から太屋岡研究室の指導教官として暖 かく見守り、多くの貴重な御指導や御助言を賜りました太 屋岡先生には厚く感謝の気持ちを申し上げます。また同研 究室で、研究や発表練習において様々な御助言や御指導を いただいた専攻科生や同研究室の５年生の皆様方にも感謝 54
© Copyright 2024 ExpyDoc
