(57)【要約】【課題】１入力複数出力系の一次遅れシステムにおいて、質や量の不十分な入出力波形データからでも、精度と信頼性の高いモデルパラメータを推定演算してモデルを構築する。【解決手段】入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算装置１０において、当該モデルパラメータ推定演算装置１０は、上記入出力データに基づいて、上記１入力複数出力系の一次遅れシステムの時定数と定常利得との少なくとも１つを含む先験的情報を、制約条件と評価関数のうちの少なくとも一方として考慮して上記モデルパラメータを推定演算する。ここで、上記先験的情報は、時定数と定常利得の上限値と下限値の少なくとも一方である。【選択図】図３
2008287343.tif 000002
【特許請求の範囲】
【請求項１】
入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算装置において、
上記入出力データに基づいて、上記１入力複数出力系の一次遅れシステムの時定数と定常利得との少なくとも１つを含む先験的情報を、制約条件と評価関数のうちの少なくとも一方として考慮して上記モデルパラメータを推定演算する演算手段を備えたことを特徴とするモデルパラメータ推定演算装置。
【請求項２】
上記先験的情報は、時定数と定常利得の上限値と下限値の少なくとも一方であることを特徴とする請求項１記載のモデルパラメータ推定演算装置。
【請求項３】
上記演算手段は、上記入出力データに基づいて、上記先験的情報を線形の制約情報として用いて、かつ上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの誤差を示す評価関数が最小となるように、モデルパラメータを推定演算することを特徴とする請求項１又は２記載のモデルパラメータ推定演算装置。
【請求項４】
上記評価関数は、上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの二乗誤差であり、
上記演算手段は、二次計画問題法を用いて、上記線形の制約条件を満たしかつ上記評価関数を最小化するモデルパラメータを推定演算することを特徴とする請求項１乃至３のうちのいずれか１つに記載のモデルパラメータ推定演算装置。
【請求項５】
入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算方法において、
上記入出力データに基づいて、上記１入力複数出力系の一次遅れシステムの時定数と定常利得との少なくとも１つを含む先験的情報を、制約条件と評価関数のうちの少なくとも一方として考慮して上記モデルパラメータを推定演算する演算ステップを含むことを特徴とするモデルパラメータ推定演算方法。
【請求項６】
上記先験的情報は、時定数と定常利得の上限値と下限値の少なくとも一方であることを特徴とする請求項５記載のモデルパラメータ推定演算方法。
【請求項７】
上記演算ステップは、上記入出力データに基づいて、上記先験的情報を線形の制約情報として用いて、かつ上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの誤差を示す評価関数が最小となるように、モデルパラメータを推定演算することを特徴とする請求項５又は６記載のモデルパラメータ推定演算方法。
【請求項８】
上記評価関数は、上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの二乗誤差であり、
上記演算ステップは、二次計画問題法を用いて、上記線形の制約条件を満たしかつ上記評価関数を最小化するモデルパラメータを推定演算することを特徴とする請求項５乃至７のうちのいずれか１つに記載のモデルパラメータ推定演算方法。
【請求項９】
請求項５乃至８のうちのいずれか１つに記載のモデルパラメータ推定演算方法の演算ステップを含むことを特徴とするモデルパラメータ推定演算処理プログラム。
【請求項１０】
請求項９記載のモデルパラメータ推定演算処理プログラムを記録したことを特徴とするコンピュータにより読取可能な記録媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、動的モデルを有するシステムに対する理論的・経験的知識を活用するシステム同定法を用いて、当該動的モデルのモデルパラメータを推定演算するモデルパラメータの推定演算装置及び方法、モデルパラメータ推定演算プログラム並びに当該プログラムを記録したコンピュータにより読取可能な記録媒体に関する。
【背景技術】
【０００２】
あるシステムの制御系を設計したり、その挙動をシミュレートしたりするような場合には、システムの挙動を模擬できるモデルが必要となる。例えば、モデルベースのコントローラを設計する場合にはモデルが不可欠であるし、ＰＩＤ（Proportional Integral Derivative）補償器を用いる場合でも事前に制御系の挙動をシミュレートする際にモデルが必要となる。
【０００３】
一般に、モデルを構築する方法は、以下の２つ（ないしは３つ）に大別される。
【０００４】
１つは、いわゆる第一原理に基づいたモデリング（ホワイトボックスモデリングとか、物理モデリングと呼ばれる）であり、システムを支配する物理化学的な自然法則（運動方程式、電磁界方程式、物質収支、熱収支、化学反応方程式など）に基づいてモデリングを行う方法である。大規模なシステムの場合は、サブシステムに分割してそれぞれモデリングを行い、それらを統合してシステム全体のモデルを構築するアプローチである。
【０００５】
もう１つは、ブラックボックスモデリングと呼ばれるもので、システムに関する物理化学的な法則や先験的知識を利用しないアプローチである。システムをブラックボックスと見なして、観測された（有限個の）入出力波形データからモデルのパラメータを推定する方法であり、制御工学の分野ではシステム同定という。システムが比較的大規模あるいは複雑であり、物理化学法則のみによってその動特性が記述できないような場合に有効である。
【０００６】
さらには、（モデリング手法を３つに大別する場合には）３つ目の分類として、グレーボックスモデリングが存在する。これは、ホワイトボックスモデリングとブラックボックスモデリングの中間に位置づけられるもので、例えば、入出力波形データに加えて、物理化学的な自然法則を部分的に利用するようなモデリング手法である。
【０００７】
さて、モデリングの対象となるシステムには、多リンク構造のロボット（マニピュレータ）のように基本的な挙動が運動方程式で表現できるものもあるが、熱や化学反応を伴うプラント（例えば、化学プラントや燃焼炉など）のように、その動特性が自然法則から導き難いものもある。この場合、入出力波形データに基づいてモデルを構築するブラックボックスモデリングを採用することになるが、システムがプロパーかつ安定で、そのステップ応答に無視できないような振動や逆応答がない場合には、モデルの構造を「むだ時間＋一次遅れ系」で近似して差し支えない（制御やシミュレーションの目的に照らして考えた場合に問題とならない）ケースも多い。
【０００８】
例えば、非特許文献２によれば、モデルの未知パラメータは、入出力波形データから何らかの方法で推定しなければならないが、その場合、モデルを実システムに近づけようとして、モデルの次数をいたずらに増加させても、モデルの精度は向上するとは限らない。なぜなら、次数が増加することによって、推定すべきパラメータ数も多くなるため、雑音や無視した動特性・非線形特性のためにパラメータの推定精度が劣化するからである。
【０００９】
モデルの構造を「むだ時間＋一次遅れ系」とした場合、むだ時間の値は、ステップ応答などの入力波形に対する出力波形の遅れを観察したり、アクチュエーターの動作遅れやセンサの検出遅れ時間を足し合わせたり、入出力波形データを解析（例えば、相互相関関数など）したりすることによって見積もることができる。また、むだ時間のない、あるいは、むだ時間を無視できるシステムも多い。従って、入出力波形データ（むだ時間がある場合には、むだ時間の分だけ入力あるいは出力波形データの時間をシフトしたデータ）から一次遅れ系のモデルパラメータを高精度に推定することが重要である。
【００１０】
ここで、一次遅れ系について説明する。一次遅れ系の伝達関数は次式で表される。
【００１１】
【数１】
2008287343.tif 000003
【００１２】
ここで、Ｋは定常利得であり、Ｔは時定数と呼ばれる。図２は従来技術に係る一次遅れ系のステップ応答波形を示す入出力波形データを示す波形図である。ここで、Ｙｃ（ｓ）、Ｕｃ（ｓ）は、それぞれシステムの出力波形データｙｃ（ｔ）、入力波形データｕｃ（ｔ）の時系列をラプラス変換によって周波数領域で表したものである。ステップ応答波形の性質として、Ｔ秒後に出力値が定常利得の約６３．２％になり、応答の立ち上がりの傾きがＫ／Ｔになることが知られている。
【００１３】
さて、（１）式の連続時間システムをサンプリングタイムΔｔ（≪Ｔ）の離散時間システムで近似すると、以下のようになる。（１）式の右辺分母を左辺に持ってくると、次式を得る。
【００１４】
【数２】
2008287343.tif 000004
【００１５】
次いで、上記（２）式の微分要素ｓを前進差分に置き換えて近似すると次式を得る。
【００１６】
【数３】
2008287343.tif 000005
【００１７】
さらに、上記（３）式を整理して、次式を得る。
【００１８】
【数４】
2008287343.tif 000006
【００１９】
ここで、次式のように表現する。
【数５】
2008287343.tif 000007
【００２０】
ここで、ｎは整数である。そして、
【数６】
2008287343.tif 000008
【数７】
2008287343.tif 000009
【００２１】
とおくと、上記（４）式は次式で表される。
【００２２】
【数８】
2008287343.tif 000010
【００２３】
なお、実際の入出力波形データは外乱やノイズ等の影響を受けており、モデル化誤差も存在するので、（６）式には何らかの誤差が存在する。従って、厳密には（６）式に誤差項ｅ（ｋ）を導入して次式を得る。
【００２４】
【数９】
2008287343.tif 000011
【００２５】
これは、ＡＲＸ（Auto-Regressive model with eXogenious input）モデルとして知られるモデル形式である。なお、むだ時間がＬ秒（Ｌ＝τΔｔ）ある場合には、（１）式は次式で表される。
【００２６】
【数１０】
2008287343.tif 000012
【００２７】
（７）式は次式で表される。
【００２８】
【数１１】
2008287343.tif 000013
【００２９】
入力波形データをむだ時間の分だけシフトすれば、すなわち
【００３０】
【数１２】
2008287343.tif 000014
とすれば、一般性を失うことなく（７）式で表現できるので、以降、（７）式のモデルに基づいて説明する。
【００３１】
入出力波形データから（７）式のモデルパラメータａ及びｂを求める場合（ブラックボックスモデリング）、従来技術（非特許文献１のｐｐ．５５〜５８、参考文献２のｐ。７１〜ｐ。７７）では、最小二乗法を用いる方法が知られている。求めるべきモデルパラメータをベクトルθ
【００３２】
【数１３】
2008287343.tif 000015
とし、
【００３３】
【数１４】
2008287343.tif 000016
と定義すると、（７）式は
【００３４】
【数１５】
2008287343.tif 000017
と変形できる。いま、誤差項ｅ（ｋ）の時系列が平均値０、有限な分散σ２を有する白色雑音であると仮定すると、（１０）式の１段先予測値は、
【００３５】
【数１６】
2008287343.tif 000018
となることが知られており、このときの予測誤差は次式で表される。
【００３６】
【数１７】
2008287343.tif 000019
【００３７】
入出力波形データにできるだけフィッティングするようなモデルパラメータθを求めるために、評価規範として予測誤差の二乗和を考える。すなわち、時刻０からＮまでの入出力波形データが得られた場合、
【数１８】
2008287343.tif 000020
を最小化するようなモデルパラメータθを求める。同定入力が２次以上のＰＥ性（Persistently Exciting：持続励振性）であり、システムが安定かつ可観測であれば、（１３）式を最小化するモデルパラメータθは最小二乗計算によって求めることができることが知られている。
【００３８】
さて、上記のようなブラックボックスモデリングでは、モデル化誤差が白色雑音であることを仮定していたり、入力信号に対するＰＥ性などのように満たすべき条件が存在したりする。しかし、現実のプロセスにおいては、システム同定用入出力波形データを採取する実験機会が限られていたり、印加できる同定入力の波形や大きさが限られている場合も多い。また、入出力波形データには、少なからず外乱やノイズも存在し、場合によっては、無視できないほど大きな外乱が存在する場合もある。例えば、ごみ焼却炉のように、不特定で多種多様のものを燃焼させるようなプロセスでモデルを構築する場合、ごみ供給機の速度入力に対して、実際に炉内に投入されるごみの種類や量はばらついてしまう。従って、ごみ供給機の速度を入力とし、炉内の温度を出力とするようなモデルを構築する場合、入出力波形データには大きな外乱が存在することになる。しかも、データ採取の機会も限られているため、外乱の影響を打ち消すのに十分な量のデータを取得することもできない。
【００３９】
このように、入出力波形データに無視できない外乱やノイズが存在し、印加できる同定入力の波形や大きさが制約を受け、入出力波形データを採取する実験機会も限られているようなプロセスでは、高精度にモデリングするのに質的・量的に十分なデータを揃えることができないため、モデリング作業は試行錯誤の繰り返しに頼らざるを得ず、モデルを構築するまでに莫大な時間を要していた。すなわち、データの取得から、データ区間の選定、サンプリングタイムの調整、データの前処理（異常値の除去、フィルタリング等を含む。）、モデル次数の選定などを何度も何度も繰り返す必要があった。また、このように入出力波形データの質と量が限られている中で、高精度で信頼性のあるモデルを得ることは非常に困難であった。
【００４０】
【特許文献１】特開２００４−２７２９１６号公報。
【非特許文献１】足立修一著，「ＭＡＴＬＡＢによる制御のためのシステム同定」，東京電機大学出版局，ｐｐ．１〜５，１９９６年１１月３０日発行。
【非特許文献２】片山徹著，「システム同定入門（システム制御情報ライブラリー９）」，朝倉書店，ｐｐ．４〜７，１９９４年５月２５日発行。
【非特許文献３】Herbert J. A. F. Tulleken, "Grey-box Modelling and Identification Using Physical Knowledge and Bayesian Techniques", Automatica, Vol. 29, No. 2, pp.285-308, 1993.
【非特許文献４】久保幹雄ほか編集，「応用数理計画ハンドブック」，朝倉書店，ｐｐ．６５４〜６６３，２００２年４月３０日。
【非特許文献５】玉置久，「ＥＥＴｅｘｔ  システム最適化」，オーム社，ｐｐ．２〜５，ｐｐ．７３〜８０，平成１７年１１月２５日。
【発明の開示】
【発明が解決しようとする課題】
【００４１】
以上の問題点を解決するために、先験的知識を補うことによって、モデルを高精度化しようとしたものに以下の先行技術がある。
【００４２】
第１の先行技術として、非特許文献３において、モデルの安定性と定常利得の符号に関する先験的知識を制約条件として考慮したモデリング技術が開示されている。このモデリング技術の問題点は以下の通りである。プロセスに対して考慮すべき先見的情報は、定常利得の符号や安定性だけではない。一次遅れ系の場合、定常利得や時定数のおよその値（取り得る値の範囲）は、理論的・経験的に既知の場合も多いが、本手法では、これらの先験的知識を考慮に入れることができない。従って、本手法では、定常利得や時定数に関する先験的知識に合致するモデルが得られるとは限らず、モデルの精度や信頼性が不十分であるという問題点があった。
【００４３】
また、第２の先行技術として、特許文献１においては以下の手法が開示されている。まず、既存のモデリング手法を使って、入出力波形データから初期モデルを求める。この初期モデルのモデルパラメータとの誤差（重みつき誤差）が小さくなるように、かつ、定常利得や安定性、システムの応答に関する制約条件を満たすようにモデルパラメータを変更し、この変更されたモデルパラメータを新たなモデルとする。しかしながら、この手法では、変更後のモデルパラメータが、初期モデルのモデルパラメータと値が近いからといって、モデルの動特性が近いとは限らない。例えば、次式の１入力１出力で次数１の状態空間モデルを用いる。
【００４４】
【数１９】
2008287343.tif 000021
【００４５】
ここで、出力波形データｙ（ｋ）、入力波形データｕ（ｋ）、状態パラメータｘ（ｋ）、モデルパラメータａ０，ｂ０，ｃ０，ｄ０（それぞれスカラー）に対して、初期モデルが（ａ０，ｂ０，ｃ０，ｄ０）＝（０．１，１．０，１．０，０．０）だったとする。この初期モデルに対し、ｂ０，ｃ０のパラメータが異なる２つのモデル、モデル１（ａ０，ｂ０，ｃ０，ｄ０）＝（０．１，０．９，０．９，０．０）、モデル２（ａ０，ｂ０，ｃ０，ｄ０）＝（０．１，０．２，５．０，０．０）を定義する。初期モデルとパラメータ空間上で値が近いのは、モデル１の方であるが、初期モデルと動特性が近いのはモデル２の方である（この場合、初期モデルとモデル２は全く同じモデルを表している）。
【００４６】
特許文献１においては、そもそも、状態空間表現の場合、モデルのパラメータは一意に決まらない（無数の値をとり得る）ので、モデルのパラメータの近さを評価関数とすることは意味がない。
【００４７】
また、初期モデルを求めるときのみ入出力波形データを用いており、先見的情報に基づく制約条件をできるだけ満たすようにモデルパラメータを最適化計算する際には、初期モデルのみを規範としており、入出力波形データへの適合性を全く考慮していないという問題点があった。入出力波形データへの適合性は初期モデルを通して間接的に考慮されているだけであり、入出力波形データがモデリングに反映されているか否かは、初期モデルの精度に依存している。
【００４８】
本発明の目的は以上の問題点を解決し、１入力複数出力系の一次遅れシステムにおいて、質や量の不十分な入出力波形データからでも、精度と信頼性の高いモデルパラメータを推定演算してモデルを構築することができるモデルパラメータ推定演算装置及び方法、モデルパラメータ推定演算プログラム並びに当該プログラムを記録したコンピュータにより読取可能な記録媒体を提供することにある。
【課題を解決するための手段】
【００４９】
第１の発明に係るモデルパラメータ推定演算装置は、入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算装置において、
上記入出力データに基づいて、上記１入力複数出力系の一次遅れシステムの時定数と定常利得との少なくとも１つを含む先験的情報を、制約条件と評価関数のうちの少なくとも一方として考慮して上記モデルパラメータを推定演算する演算手段を備えたことを特徴とする。
【００５０】
上記モデルパラメータ推定演算装置において、上記先験的情報は、時定数と定常利得の上限値と下限値の少なくとも一方であることを特徴とする。
【００５１】
また、上記モデルパラメータ推定演算装置において、上記演算手段は、上記入出力データに基づいて、上記先験的情報を線形の制約情報として用いて、かつ上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの誤差を示す評価関数が最小となるように、モデルパラメータを推定演算することを特徴とする。
【００５２】
さらに、上記モデルパラメータ推定演算装置において、上記評価関数は、上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの二乗誤差であり、上記演算手段は、二次計画問題法を用いて、上記線形の制約条件を満たしかつ上記評価関数を最小化するモデルパラメータを推定演算することを特徴とする。
【００５３】
第２の発明に係るモデルパラメータ推定演算方法は、入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算方法において、
上記入出力データに基づいて、上記１入力複数出力系の一次遅れシステムの時定数と定常利得との少なくとも１つを含む先験的情報を、制約条件と評価関数のうちの少なくとも一方として考慮して上記モデルパラメータを推定演算する演算ステップを含むことを特徴とする。
【００５４】
上記モデルパラメータ推定演算方法において、上記先験的情報は、時定数と定常利得の上限値と下限値の少なくとも一方であることを特徴とする。
【００５５】
また、上記モデルパラメータ推定演算方法において、上記演算ステップは、上記入出力データに基づいて、上記先験的情報を線形の制約情報として用いて、かつ上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの誤差を示す評価関数が最小となるように、モデルパラメータを推定演算することを特徴とする。
【００５６】
さらに、上記モデルパラメータ推定演算方法において、上記評価関数は、上記入出力データに含まれる出力データと、モデルパラメータのモデルで推定演算したときの出力データとの二乗誤差であり、上記演算ステップは、二次計画問題法を用いて、上記線形の制約条件を満たしかつ上記評価関数を最小化するモデルパラメータを推定演算することを特徴とする。
【００５７】
第３の発明に係るモデルパラメータ推定演算処理プログラムは、上記モデルパラメータ推定演算方法の演算ステップを含むことを特徴とする。
【００５８】
第４の発明に係るコンピュータにより読み取り可能な記録媒体は、上記モデルパラメータ推定演算処理プログラムを記録したことを特徴とする。
【発明の効果】
【００５９】
従って、本発明に係るモデルパラメータ推定演算装置及び方法によれば、以下の特有の効果を有する。
（１）時定数に関する先験的情報（時定数の上下限）定常利得の少なくとも一方に関する先験的情報（定常利得の上下限）を考慮してモデリングすることにより、質が悪く量の少ない入出力波形データからでも、精度と信頼性の高いモデルを求めることができる。
（２）精度と信頼性の高いモデルを使って制御系を設計したり、シミュレーションを実施したりできるので、制御系やシミュレーションの精度・信頼性も向上する。
（３）時定数や定常利得に関する先験的情報を、制約条件と評価関数のうちの少なくとも一方として活用した場合、得られたモデルは必然的にこれらの先験的情報に合致したものとなる。従来法では、入出力波形データからモデルを求めてから、そのモデルが先験的情報に合致するか検証し、合致しない場合には、再度、モデリング工程を（最初から、あるいは、途中から）やり直す必要があった。特に、外乱やノイズの大きなプロセス、印加できる入力波形が限られたプロセス、採取できるデータ数が少ないプロセスにおいては、先験的情報に合致したモデルを構築することが困難であり、入出力波形データの採取・データ区間の選定・サンプリングタイムの調整・データの前処理などを試行錯誤しながら繰り返す必要があり、モデリング工程に多大な時間と労力を要した。先験的情報を予め考慮に入れることによって、モデリング工程における試行錯誤・やり直しの回数を大幅に減らすことが可能となり、短期間で精度と信頼性の高いモデルを構築できるようになる。その結果、制御系構築やシミュレーションによる検討のリードタイムが短縮化され、早い段階からそれらの効果を享受することができる。
【発明を実施するための最良の形態】
【００６０】
以下、本発明に係る実施形態について図面を参照して説明する。なお、以下の各実施形態において、同様の構成要素については同一の符号を付している。
【００６１】
図１は本発明の一実施形態に係るモデルパラメータ推定演算装置１０の構成を示すブロック図である。本実施形態のモデルパラメータ推定演算装置１０は、図１に示すように、ディジタル計算機を含むように構成され、例えば図２のモデルパラメータ推定演算処理を実行することにより、モデルパラメータ推定演算結果を表示して出力する。モデルパラメータ推定演算装置１０は、入出力データに基づいて１入力複数出力系の一次遅れシステムのモデルパラメータを推定して演算するモデルパラメータ推定演算装置であり、
（ａ）当該モデルパラメータ推定演算装置１０の動作及び処理を演算及び制御するコンピュータのＣＰＵ（中央演算処理装置）２０と、
（ｂ）オペレーションプログラムなどの基本プログラム及びそれを実行するために必要なデータを格納するＲＯＭ（読み出し専用メモリ）２１と、
（ｃ）ＣＰＵ２０のワーキングメモリとして動作し、図３のモデルパラメータ推定演算処理において必要なパラメータやデータを一時的に格納するＲＡＭ（ランダムアクセスメモリ）２２と、
（ｄ）例えばハードディスクメモリで構成され、入力パラメータのデータやシミュレーション結果のデータなどのデータを格納するデータメモリ２３と、
（ｅ）例えばハードディスクメモリで構成され、ＣＤ−ＲＯＭドライブ装置４５を用いて読みこんだ図３のモデルパラメータ推定演算処理プログラムを格納するプログラムメモリ２４と、
（ｆ）所定のデータや指示コマンドを入力するためのキーボード４１に接続され、キーボード４１から入力されたデータや指示コマンドを受信して所定の信号変換などのインターフェース処理を行ってＣＰＵ２０に伝送するキーボードインターフェース３１と、
（ｇ）ＣＲＴディスプレイ４３上で指示コマンドを入力するためのマウス４２に接続され、マウス４２から入力されたデータや指示コマンドを受信して所定の信号変換などのインターフェース処理を行ってＣＰＵ２０に伝送するマウスインターフェース３２と、
（ｈ）ＣＰＵ２０によって処理されたデータや設定指示画面などを表示するＣＲＴディスプレイ４３に接続され、表示すべき画像データをＣＲＴディスプレイ４３用の画像信号に変換してＣＲＴディスプレイ４３に出力して表示するディスプレイインターフェース３３と、
（ｉ）ＣＰＵ２０によって処理されたデータ及び所定のモデルパラメータ推定演算結果などを印字するプリンタ４４に接続され、印字すべき印字データの所定の信号変換などを行ってプリンタ４４に出力して印字するプリンタインターフェース３４と、
（ｊ）図３のモデルパラメータ推定演算処理プログラムが記憶されたＣＤ−ＲＯＭ４５ａからモデルパラメータ推定演算処理プログラムのプログラムデータを読み出すＣＤ−ＲＯＭドライブ装置４５に接続され、読み出されたモデルパラメータ推定演算処理プログラムのプログラムデータを所定の信号変換などを行ってプログラムメモリ２４に転送するドライブ装置インターフェース３５と、
（ｋ）入力信号発生器５１からのアナログ入力信号をＡ／Ｄ変換するＡ／Ｄ変換器２５ａと、上記アナログ入力信号に応答してモデリング対象システム５０から出力されるアナログ出力信号を検出して出力する出力信号検出器５２からのアナログ出力信号をＡ／Ｄ変換するＡ／Ｄ変換器２５ｂとを含む入力インターフェース２５とを備える。
ここで、これらの回路２０−２５及び３１−３４はバス３０を介して接続される。なお、Ａ／Ｄ変換器２５ａ，２５ｂによりＡ／Ｄ変換された入力データ及び出力データはバス３０を介してデータメモリ２３に格納される。
【００６２】
次いで、本実施形態に係るモデルパラメータ推定演算装置１０において用いるモデルパラメータ推定演算処理の推定演算方法について以下に詳述する。
【００６３】
ある入力からｉ番目の出力（ｉ＝１，２，…，Ｍ；ここで、自然数Ｍは出力の個数であり、以降各記号の添字ｉは、ｉ番目の出力に関するものであることを示す。）までの
（Ａ）モデルパラメータ：
【数２０】
2008287343.tif 000022
【数２１】
2008287343.tif 000023
【数２２】
2008287343.tif 000024
（Ｂ）入出力波形データ：
【数２３】
2008287343.tif 000025
に対して、モデリングの際に最小化すべき評価関数を次式のように表す。
【００６４】
【数２４】
2008287343.tif 000026
【００６５】
いま、先験的情報として、時定数Ｔｉ（＞０）の上下限
【数２５】
2008287343.tif 000027
と、定常利得Ｋｉの上下限
【数２６】
2008287343.tif 000028
が与えられるとする。また、（６）式の離散システムの安定性の条件は、よく知られているように、次式で表される。
【００６６】
【数２７】
2008287343.tif 000029
【００６７】
（１５）、（１６）、（１７）式の制約条件の下で、（１３）式の評価関数を最小化することを考える。このままでは、モデルパラメータ
【数２８】
2008287343.tif 000030
に関して、制約条件が線形でないので、次式のようにおく。
【００６８】
【数２９】
2008287343.tif 000031
【数３０】
2008287343.tif 000032
【００６９】
このとき、（１３）式は、次式で表される。
【００７０】
【数３１】
2008287343.tif 000033
【００７１】
ここで、
【数３２】
2008287343.tif 000034
【００７２】
【数３３】
2008287343.tif 000035
【００７３】
ここで、
【数３４】
2008287343.tif 000036
（時定数は正のためである。）、サンプリングタイム
【数３５】
2008287343.tif 000037
であるから、（１８）式を用いると（１５）式の制約条件は次式で表される。
【００７４】
【数３６】
2008287343.tif 000038
【００７５】
また、（１６）式の各辺にＸｉ（＞０）をかけて、（１９）式を代入すると、（１６）式の制約条件は次式で表される。
【００７６】
【数３７】
2008287343.tif 000039
【００７７】
さらに、
【数３８】
2008287343.tif 000040
を（１７）式の制約条件に代入して整理すると、次式を得る。
【００７８】
【数３９】
2008287343.tif 000041
【００７９】
さらに、（１８）式を代入して整理すると、次式を得る。
【００８０】
【数４０】
2008287343.tif 000042
【００８１】
ここで、（２１）式、（２２）式、（２４）式をまとめると、次式を得る。
【００８２】
【数４１】
2008287343.tif 000043
【００８３】
なお、（２４）式の左側の制約条件
【数４２】
2008287343.tif 000044
は、（２１）式の左側の制約条件
【数４３】
2008287343.tif 000045
を満たせば自然に成り立つので、（２５）式では考慮していない。パラメータ
【数４４】
2008287343.tif 000046
に対して、（２５）式の線形制約条件を満たし、（２０）式の評価関数を最小化するものを求める問題は、例えば非特許文献４−５にあるように、二次計画問題として広く知られ、当該二次計画問題法に係る既知の数値計算アルゴリズムを用いて、すなわち、市販の数値計算ツールを用いて、パラメータ
【数４５】
2008287343.tif 000047
を求めることができる。そして、（１８）式及び（１９）式から、
【数４６】
2008287343.tif 000048
【数４７】
2008287343.tif 000049
によって、モデルパラメータである時定数Ｔｉと定常利得Ｋｉの値を求めることができる。
【００８４】
図３は図１のモデルパラメータ推定演算装置１０によって実行されるモデルパラメータ推定演算処理を示すフローチャートである。図３のフローチャートを使用して、本実施形態に係るモデルパラメータ推定演算処理プログラムのアルゴリズムについて以下に説明する。
【００８５】
［ステップＳ１］
モデリング対象システムに対して入力信号発生器５１からアナログ入力信号を印加し、これに応答してモデリング対象システム５０から出力されるアナログ出力信号を出力信号検出器５２により検出する。アナログ入力信号は入力インターフェース２５内のＡ／Ｄ変換器２５ａに入力されてＡ／Ｄ変換された後、データメモリ２３に格納される。また、アナログ出力信号は入力インターフェース２５内のＡ／Ｄ変換器２５ｂに入力されてＡ／Ｄ変換された後、データメモリ２３に格納される。このとき、制御やシミュレーションの目的に従って、何を入力波形データとして、何を出力波形データとするか決めておく。例えば、制御系を設計する際には、何を制御するために（制御量：出力波形データ）、何を操作するか（操作入力：入力波形データ）決めておく必要がある。入出力波形データとして、複数の候補がある場合には、それらすべてのデータを採取しておいてもよい。
【００８６】
具体的には、例えばごみ処理プラントの制御においては、給塵装置のごみ供給速度をモデルの入力とし、炉内温度やボイラ蒸気発生量をモデルの出力とすることや、供給空気量をモデルの入力とし、砂層温度をモデルの出力とすることが考えられる。転炉の炉圧制御系においては、ＲＳＥ（リング・スリット・エレメント）と呼ばれる排気系のダンパー開度をモデルの入力とし、転炉炉口部の圧力をモデルの出力とすることが考えられる。鋼板の加熱炉の温度予測シミュレーションにおいては、加熱バーナーの燃料供給量をモデル入力とし、鋼板の温度をモデル出力とすることが考えられる。
【００８７】
モデリングの際、プロセスに印加する入力波形として、ステップ入力、正弦波の重ね合わせ、疑似白色２値信号（Ｍ系列など）、ランダム入力信号などが使われる。なお、印加する入力データには、理論上満たすべき条件（ＰＥ性など）が存在するので、理想的にはこれらの条件を満たしていることが望ましい。しかしながら、モデリング対象となるプロセスには、装置上や操業上、印加できる入力の大きさや波形に制約がある場合も多く、必ずしも、理想的な入力データを印加できるとは限らないのが現状である。
【００８８】
［ステップＳ２］
データメモリ２３に格納された入力データ及び出力データを含む入出力波形データの中から、モデリングに使用するデータの区間を選定する。入出力波形データには、外乱やノイズ、異常値（センサの検出誤差やデータ記録時の異常）などが含まれている場合があり、モデリングの際には、これら外乱やノイズ、異常値の影響がより小さい部分のデータを抜き出して使用した方がよい。一方で、モデリングに使用するデータの長さが短すぎても、外乱やノイズの影響を大きく受けて、モデルパラメータの推定精度が劣化することになる。モデリングに使用するデータ区間は、これらのバランスを考えて選定して、選定結果を例えばキーボード４１又はマウス４２を用いて入力する。
【００８９】
［ステップＳ３］
モデルのサンプリングタイムΔｔを設定して例えばキーボード４１又はマウス４２を用いて入力する。サンプリングタイムΔｔは、モデリング対象となるシステムの応答速度やモデリングしたい周波数帯域などに応じて設定する。
【００９０】
［ステップＳ４］
ステップＳ３で設定したサンプリングタイムΔｔに基づいてデシメーション（データのフィルタリングと間引き処理をいう。）したり、異常値を除去したり、外乱やノイズの影響を低減するためにフィルタリングをしたりするなど、入出力波形データの前処理を行う。
【００９１】
［ステップＳ５］
先験的情報（理論的・経験的検知）に基づいて、定常利得の上下限値と時定数の上下限値を設定し、例えばキーボード４１又はマウス４２を用いて入力する。この際、必ずしも上限値と下限値の両方を設定する必要は必ずしもなく、上限値だけ、下限値だけでもかまわない。また、定常利得のみに制約を設けても、時定数のみに制約を設けてもよい。ただし、定常利得や時定数のおよその値が分かっている場合には、それぞれに上下限値を設定した方が望ましい。これはモデルの精度・信頼性向上のためである。
【００９２】
［ステップＳ６］
（２５）式の制約条件の下で、（２０）式の評価関数を最小化する二次計画問題を解く。二次計画問題の解法は広く知られているので、これらの問題を解く際には、市販の最適化計算ツールを使用してもよいし、自分で問題を解くプログラムを構築してもよい。なお、時定数の下限値を制約として考慮しない場合には、（２５）式の不等式の１行目を取り除けばよい。時定数の上限値を制約として考慮しない場合には、（２５）式の２行目の右辺の
【数４８】
2008287343.tif 000050
を０と置き換えればよい。これは、（２４）式の安定性の制約条件を考慮するためである。定常利得の下限値を制約として考慮しない場合には、（２５）式の不等式の３行目を取り除けばよい。同様に、定常利得の上限値を制約として考慮しない場合には、（２５）式の不等式の４行目を取り除けばよい。そして、
【００９３】
［ステップＳ７］
前ステップで求めたパラメータ
【数４９】
2008287343.tif 000051
から、（２６）式と（２７）式によって、時定数Ｔと定常利得Ｋの値を計算する。
【００９４】
［ステップＳ８］
計算されたモデルパラメータを有するモデルが妥当であるか検証する。例えば、モデリングに用いていない部分のデータを使って、モデルを用いてシミュレートした出力波形と、実システムから採取した出力波形データを比較して検証する（クロスバリデーション）。
【００９５】
［ステップＳ９］
前ステップＳ８においてモデルを検証した結果、制御やシミュレーションの目的に照らし合わせてモデリングの精度が所定のしきい値以上であって十分であるならば、そのモデルを採用する。もし、精度的に不十分であると判断されたなら、ステップＳ２からステップＳ４のいずれかに戻って、以降のステップの処理をやり直す。それでも、満足のいくモデルが求まらない場合には、ステップＳ１に戻って、入出力波形データを採取するところからやり直す。ここで、処理の戻り先はケースバイケースで決定してもよい。
【００９６】
［ステップＳ１０］
得られたモデルパラメータを採用してデータメモリ２３に格納した後、得られたモデルのボード線図（利得線図と位相線図）を自動的にプロットして例えばＣＲＴディスプレイ４３に表示出力する。ここで、周波数領域でみた場合にも妥当なモデルが得られているか、すなわち制御上やシミュレーション上で重要な周波数帯域で十分な特性が得られているかについて検討したり、得られたモデルと実際の入力波形データを使って時間応答シミュレーションを行い、実際の出力波形と比較したりする。例えば、実際の出力波形とシミュレーション波形の平均二乗誤差（Root Mean Square）が許容値より大きければステップＳ１〜Ｓ４のいずれかに戻ってやり直してもよいし、複数のデータでそれぞれモデルを求めておいて最も精度が良いものを最終的に採用してもよい。あるいは、モデルの残差を解析して、すなわち残差の自己相関関数を求めたり、残差と入力との相互相関関数を求めたりして、残差が白色で、かつ、入力と独立かどうか検証する。
【００９７】
以上のステップに基づいてモデリングすることによって、時定数や定常利得に関する先験的知識を活用した高精度で信頼性の高いモデリングを実現することができる。
【実施例１】
【００９８】
本実施形態に係るモデルパラメータ推定演算装置１０のモデルパラメータ推定演算方法の有効性を示すために以下に実施例を示す。
【００９９】
図４は実施例１に係るシステムにおける入力波形データと外乱加算後の入力波形データとを示す波形図である。また、図５は実施例１に係るシステムにおける出力波形データを示す波形図である。さらに、図６は実施例１に係るシステムにおける、真のモデルの出力波形データと、本実施形態のモデルパラメータ推定演算装置１０による推定演算結果である出力波形データと、従来技術に係るモデルパラメータ推定演算結果である出力波形データとを示す波形図である。
【０１００】
実施例１において、真のモデリング対象は、（１）式において定常利得Ｋ＝１、時定数Ｔ＝１であるとする。すなわち、モデリング対象の真の伝達関数が次式で表されるとする。
【０１０１】
【数５０】
2008287343.tif 000052
【０１０２】
ここで、サンプリングタイムΔｔ＝０．１［秒］とすると、（６）式のように離散化した真のモデルは次式で表される。
【０１０３】
【数５１】
2008287343.tif 000053
【０１０４】
いま、先験的情報として、時定数Ｔが
【数５２】
2008287343.tif 000054
の範囲にあり、定常利得Ｋが
【数５３】
2008287343.tif 000055
の範囲にあることが分かっているものとする。この場合、（２５）式の制約条件は次式で表される。
【０１０５】
【数５４】
2008287343.tif 000056
【０１０６】
（２９）式の特性を持つモデリング対象に対して、図４の●印で表される入力（周期１秒、振幅１の正弦波、ｕ（０）〜ｕ（９）、ｕ（０）＝ｕ（５）＝０）を入力波形データとして印加したとする。しかし、入力に外乱が加わった結果、実際のモデリング対象には、図４の◆印で表される入力が加わったとする。実際には、どのような外乱が入力に加わったかは測定できないものとする。なお、ここで加えた外乱は、平均値０、標準偏差１．２の正規分布に従うランダム系列を市販の演算ソフトで発生させ、そのうちの最初の１０点を採用したものである。図４の◆印の入力を加えると、（２９）式の特性を持つシステムは、図５のような応答波形（ｙ（１）〜ｙ（１０））を出力する（ただし、初期値ｙ（０）＝０とする）。
【０１０７】
さて、実際のモデリングの際には、入力に加わった外乱は未知であるので、図４の●印で表される正弦波入力波形データと、図５の◆印の出力波形データを用いて、モデリングをすることになる。（３２）式の制約条件の下で、（２０）式の評価関数（Ｎ＝１０）を最小化するモデルパラメータ
【数５５】
2008287343.tif 000057
を二次計画問題の解法を利用して求め、（２６）式と（２７）式を利用して時定数と定常利得を求めると、時定数Ｔ＝１．１００４、定常利得Ｋ＝１．２のモデルを求めることができる。これらの値は、（３０）式と（３１）式の制約条件を満たしていることに注意する。このモデルに対し、図４の●印の正弦波入力を加えると、シミュレーションの結果、出力波形は図６の□印となる。また、（２９）式の特性を持つ真のモデル（時定数Ｔ＝１、定常利得Ｋ＝１）に、図４の●印の正弦波入力を加えると、その出力波形は、図６の●印となる。これらふたつの波形（図６の□印と●印）は、ほとんど一致しており、図４の◆印のような大きな入力外乱が加わったような場合でも、時定数と定常利得に関する先験的知識を活用することによって、高精度にモデリングができることが分かる。
【０１０８】
一方、従来技術、すなわち、（３０）式と（３１）式のような制約条件をつけず、（１３）式の評価関数のみを最小化するようにモデルパラメータを求めた場合、図４の●印の正弦波入力と図５の◆印の出力波形データからは、時定数Ｔ＝−０．４８２３９、定常利得Ｋ＝−０．７７３８８のモデルが求まってしまい、真のモデルとは、時定数と定常利得が全く異なるものになってしまう。なお、従来技術で求めたモデルに対し、図４の●印の正弦波入力を加えると、出力波形は図６の△印となり、真モデルの波形（図６の●印）とは、全く異なるものになる。なお、モデルが（２３）式の安定性の条件を満たすためには、サンプリングタイムΔｔ＝０．１なので、時定数がＴ≧０．０５でなければならない。時定数がＴ＝−０．４８２３９であるということは、不安定なモデルが求まっていることを意味する。
【０１０９】
さらに、実施例１において外乱の大きさを変化させて、比較した結果を表１に示す。
【０１１０】
【表１】
2008287343.tif 000058
【０１１１】
ここで加えた外乱は、平均値０、標準偏差σの正規分布に従うランダム系列を市販の演算ソフトで発生させ、そのうちの最初の１０点を採用したものであり、標準偏差σの値を０．０（すなわち外乱なし）〜２．０まで、０．１刻みで増加させたものである。外乱がない場合（σ＝０．０）には、従来法でも真のモデルパラメータが求まるが、入力外乱の標準偏差σが０．２より大きくなると、従来法では、先験的知識に合致しない（すなわち（３０）式と（３１）式の制約条件を満たさない）モデルが求まってしまう。中には、時定数Ｔ＜０．０５であるような、不安定なモデルが求まっているケースもあり、定常利得の符号が逆になってしまっているものもある。このように、本実施形態に係るモデルパラメータ推定演算方法は、時定数と定常利得に関する先験的知識を活用することによって、高精度で信頼性の高いモデルが求まることが確認できる。
【０１１２】
なお、表１の本実施形態に係るモデルパラメータ推定の結果、定常利得Ｋが０．８か１．２であり、時定数Ｔが０．８＜Ｔ＜１．２であるものは、先験的知識として定常利得のみを活用しても（定常利得の制約条件のみを考慮しても）同じ結果が得られたものである。裏を返せば、定常利得のみに制約条件をつけてもモデリング精度と信頼性を高める効果があった事例となっている。反対に、時定数Ｔが０．８か１．２であり、定常利得Ｋが０．８＜Ｋ＜１．２であるものは、先験的知識として時定数のみを活用しても（時定数の制約条件のみを考慮しても）同じ結果が得られたものである。裏を返せば、時定数のみに制約条件をつけてもモデリング精度と信頼性を高める効果があった事例となっている。
【実施例２】
【０１１３】
図７は実施例２に係るシステムにおける入力波形データと外乱加算後の入力波形データとを示す波形図である。また、図８は実施例２に係るシステムにおける出力波形データを示す波形図である。さらに、図９は実施例２に係るシステムにおける、真のモデルの出力波形データと、本実施形態のモデルパラメータ推定演算装置１０による推定演算結果である出力波形データと、従来技術に係るモデルパラメータ推定演算結果である出力波形データとを示す波形図である。
【０１１４】
実施例２においては、（２９）式の特性を持つモデリング対象に対して、図７の●印で表される入力（周期１秒、振幅１の正弦波、ｕ（０）〜ｕ（９）、ｕ（０）＝ｕ（５）＝０）を入力波形データとして印加したとする。しかし、入力に外乱が加わった結果、実際のモデリング対象には、図７の◆印で表される入力が加わったとする（実際には、どのような外乱が入力に加わったかは測定できないものとする）。なお、ここで加えた外乱は、平均値０、標準偏差１．７の正規分布に従うランダム系列を市販の演算ソフトで発生させ、そのうちの最初の１０点を採用したものである。図７の◆印の入力を加えると、（２９）式の特性を持つシステムは、図８のような応答波形（ｙ（１）〜ｙ（１０））を出力する（ただし、初期値ｙ（０）＝０とする）。
【０１１５】
さて、実際のモデリングの際には、入力に加わった外乱は未知であるので、図７の●印で表される正弦波入力波形データと、図８の◆印の出力波形データを用いて、モデリングをすることになる。（３２）式の制約条件の下で、（２０）式の評価関数（Ｎ＝１０）を最小化するモデルパラメータ
【数５６】
2008287343.tif 000059
を二次計画問題の解法を利用して求め、（２６）式と（２７）式を利用して時定数と定常利得を求めると、時定数Ｔ＝１．２、定常利得Ｋ＝０．８のモデルが求まめることができる。このモデルに対し、図７の●印の正弦波入力を加えると、シミュレーションの結果、出力波形は図９の□印となる。また、（２９）式の特性を持つ真のモデル（時定数Ｔ＝１、定常利得Ｋ＝１）に、図７の●印の正弦波入力を加えると、その出力波形は、図９の●印となる。これらふたつの波形（図９の□印と●印）は、およそ一致しており、図７の◆印のように非常に大きな入力外乱が加わったような場合でも、時定数と定常利得に関する先験的知識を活用することによって、精度良くモデリングができることが分かる。
【０１１６】
一方、従来技術に係る方法を用いて、すなわち、（３０）式と（３１）式のような制約条件をつけず、（１３）式の評価関数のみを最小化するようにモデルパラメータを求めた場合、図７の●印の正弦波入力と図８の◆印の出力波形データからは、時定数Ｔ＝−０．６９４１２、定常利得Ｋ＝０．１５５７８のモデルが求まってしまい、真のモデルとは、時定数と定常利得が全く異なるものになってしまう（これを従来例という。）。
【０１１７】
なお、従来例で求めたモデルに対し、図７の●印の正弦波入力を加えると、出力波形は図９の△印となり、真モデルの波形（図９の●印）とは、全く異なるものになる。なお、時定数がＴ＝−０．６９４１２であるということは、安定性の条件、すなわち、時定数Ｔ≧０．０５を満たしておらず、不安定なモデルが求まっていることを意味する。
【０１１８】
さらに、非特許文献３に開示された方法のように、モデルの安定性と定常利得の符号に制約を付けた場合を考える。ここでは、定常利得の符号だけでなく、（３１）式のように定常利得の上下限も分かっている場合を考える（これを、変形例とする。）。定常利得の符号だけでなく、上下限値を制約として与えるということは、それだけ、モデルの精度と信頼性の向上が期待できることに注意する。変形例は、（３２）式の制約条件の不等式のうち、１行目と２行目を除いたケース（時定数に関する制約を考慮しないケース）に相当する。変形例では、図７の●印の正弦波入力と図８の◆印の出力波形データからは、時定数Ｔ＝３５．４２７１、定常利得Ｋ＝１．２のモデルが求まってしまい、真のモデルとは、時定数が全く異なるものになってしまう。なお、変形例で求めたモデルに対し、図７の●印の正弦波入力を加えると、出力波形は図９の◇印となり、真モデルの波形（図９の●印）とは、大きくずれたものになってしまう。このように、定常利得と安定性に関する制約条件では、モデルの精度と信頼性が失われる場合があり、時定数に関する制約条件を考慮することも重要であることが分かる。
【０１１９】
さらに、実施例２において外乱の大きさを変化させて、比較した結果を表２に示す。
【０１２０】
【表２】
2008287343.tif 000060
【０１２１】
ここで加えた外乱は、平均値０、標準偏差σの正規分布に従うランダム系列を市販の演算ソフトで発生させ、そのうちの最初の１０点を採用したものであり、標準偏差σの値を０．０（すなわち外乱なし）〜２．０まで、０．１刻みで増加させたものである。外乱がない場合（標準偏差σ＝０．０）には、従来法でも真のモデルパラメータが求まるが、入力外乱の標準偏差σが大きくなると、従来技術では、先験的知識に合致しない（すなわち（３０）式や（３１）式の制約条件を満たさない）モデルが求まってしまう。従来例の中には、時定数Ｔ＜０．０５であるような、不安定なモデルが求まっているケースもあり、定常利得の符号が逆になってしまっているものもある。さらに、変形例の中には、（３０）式の時定数の制約条件を満たしていないモデルが多々存在する。このように、本技術は、時定数と定常利得に関する先験的知識を活用することによって、高精度で信頼性の高いモデルが求まることが確認できる。
【０１２２】
なお、表２の本実施形態に係るモデルパラメータ推定の結果、定常利得Ｋが０．８か１．２であり、時定数Ｔが０．８＜Ｔ＜１．２であるものは、先験的知識として定常利得のみを活用しても（定常利得の制約条件のみを考慮しても）同じ結果が得られたものである。裏を返せば、定常利得のみに制約条件をつけてもモデリング精度と信頼性を高める効果があった事例となっている（この場合には、変形例でも同等の結果が得られる）。反対に、時定数Ｔが０．８か１．２であり、定常利得Ｋが０．８＜Ｋ＜１．２であるものは、先験的知識として時定数のみを活用しても（時定数の制約条件のみを考慮しても）同じ結果が得られたものである。裏を返せば、時定数のみに制約条件をつけてもモデリング精度と信頼性を高める効果があった事例となっている（時定数に関する先験的情報を活用することによって、変形例に比べて、精度と信頼性の高いモデルを求めることができることが分かる）。
【０１２３】
なお、時定数は定常利得に関する先験的情報を制約条件として考慮したが、評価関数として考慮しても同様の効果があると考えられる。例えば、時定数の上下限の平均値からの誤差を評価関数に入れてもよいし、定常利得の上下限の平均値からの誤差を評価関数に入れてもよい。評価関数として考慮する場合には，入出力データへのフィッティング度合いを評価する式（例えば（２０）式）とのバランス（どの評価項目を重視するか）を調整するための重み係数を付け加えた方がよい。以下に具体例を示す。
【０１２４】
【数５７】
2008287343.tif 000061
【０１２５】
ここで、
【数５８】
2008287343.tif 000062
又は
【数５９】
2008287343.tif 000063
又は
【数６０】
2008287343.tif 000064
【０１２６】
【数６１】
2008287343.tif 000065
【数６２】
2008287343.tif 000066
【数６３】
2008287343.tif 000067
【数６４】
2008287343.tif 000068
【０１２７】
ここで、ｗｉ１，ｗｉ２，ｗｉ３は重み係数である。
【０１２８】
さらには、先験的情報を評価関数と制約条件の両方で考慮してもよい。制約条件を付けることによって、最悪の場合でも、モデルパラメータをその範囲に納めることができるし、モデルパラメータが制約条件の境界に張り付くことを防止する効果もある。実際に、上記評価関数を用いて、表１の実施例と同じ条件（同じデータ）で計算した結果を表３に示す。
【０１２９】
【表３】
2008287343.tif 000069
【０１３０】
なお、上記の式において、
【数６５】
2008287343.tif 000070
を用いて、ｗｉ１＝ｗｉ２＝ｗｉ３＝１とした。表１と表３を比較して分かるように、表３の事例では、表１の事例に比べて、時定数Ｔや定常利得Kが制約条件の境界（上下限値）に張り付く頻度が少なくなっていることが確認できる。モデルパラメータの真値が制約条件の境界には存在せず、真値が上下限値の中ほどに存在するような場合には、このような評価関数と制約条件を考慮することによって、より真値に近いパラメータが得られる可能性を高くすることができる。
【０１３１】
＜変形例＞
以上の実施形態又は実施例においては、図３のモデルパラメータ推定演算処理プログラムデータをＣＤ−ＲＯＭ４５ａに格納して実行するときにプログラムメモリ２４にロードして実行しているが、本発明はこれに限らず、ＣＤ−Ｒ、ＣＤ−ＲＷ、ＤＶＤ、ＭＯなどの光ディスク又は光磁気ディスクの記録媒体、もしくは、フロッピーディスクなどの磁気ディスクの記録媒体など種々の記録媒体に格納してもよい。これらの記録媒体は，コンピュータで読み取り可能な記録媒体である。また、図３のモデルパラメータ推定演算処理のデータを予めプログラムメモリ２４に格納して当該処理を実行してもよい。
【産業上の利用可能性】
【０１３２】
以上詳述したように、本発明に係るモデルパラメータ推定演算装置及び方法によれば、以下の特有の効果を有する。
（１）時定数に関する先験的情報（時定数の上下限）定常利得の少なくとも一方に関する先験的情報（定常利得の上下限）を考慮してモデリングすることにより、質が悪く量の少ない入出力波形データからでも、精度と信頼性の高いモデルを求めることができる。
（２）精度と信頼性の高いモデルを使って制御系を設計したり、シミュレーションを実施したりできるので、制御系やシミュレーションの精度・信頼性も向上する。
（３）時定数や定常利得に関する先験的情報を、制約条件と評価関数のうちの少なくとも一方として活用した場合、得られたモデルは必然的にこれらの先験的情報に合致したものとなる。従来法では、入出力波形データからモデルを求めてから、そのモデルが先験的情報に合致するか検証し、合致しない場合には、再度、モデリング工程を（最初から、あるいは、途中から）やり直す必要があった。特に、外乱やノイズの大きなプロセス、印加できる入力波形が限られたプロセス、採取できるデータ数が少ないプロセスにおいては、先験的情報に合致したモデルを構築することが困難であり、入出力波形データの採取・データ区間の選定・サンプリングタイムの調整・データの前処理などを試行錯誤しながら繰り返す必要があり、モデリング工程に多大な時間と労力を要した。先験的情報を予め考慮に入れることによって、モデリング工程における試行錯誤・やり直しの回数を大幅に減らすことが可能となり、短期間で精度と信頼性の高いモデルを構築できるようになる。その結果、制御系構築やシミュレーションによる検討のリードタイムが短縮化され、早い段階からそれらの効果を享受することができる。
【図面の簡単な説明】
【０１３３】
【図１】本発明の一実施形態に係るモデルパラメータ推定演算装置１０の構成を示すブロック図である。
【図２】従来技術に係る一次遅れ系のステップ応答波形を示す入出力波形データを示す波形図である。
【図３】図１のモデルパラメータ推定演算装置１０によって実行されるモデルパラメータ推定演算処理を示すフローチャートである。
【図４】実施例１に係るシステムにおける入力波形データと外乱加算後の入力波形データとを示す波形図である。
【図５】実施例１に係るシステムにおける出力波形データを示す波形図である。
【図６】実施例１に係るシステムにおける、真のモデルの出力波形データと、本実施形態のモデルパラメータ推定演算装置１０による推定演算結果である出力波形データと、従来技術に係るモデルパラメータ推定演算結果である出力波形データとを示す波形図である。
【図７】実施例２に係るシステムにおける入力波形データと外乱加算後の入力波形データとを示す波形図である。
【図８】実施例２に係るシステムにおける出力波形データを示す波形図である。
【図９】実施例２に係るシステムにおける、真のモデルの出力波形データと、本実施形態のモデルパラメータ推定演算装置１０による推定演算結果である出力波形データと、従来技術に係るモデルパラメータ推定演算結果である出力波形データとを示す波形図である。
【符号の説明】
【０１３４】
１０…モデルパラメータ推定演算装置、
２０…ＣＰＵ、
２１…ＲＯＭ、
２２…ＲＡＭ、
２３…データメモリ、
２４…プログラムメモリ、
２５…入力インターフェース、
２５ａ，２５ｂ…Ａ／Ｄ変換器、
３０…バス、
３１…キーボードインターフェース、
３２…マウスインターフェース、
３３…ディスプレイインターフェース、
３４…プリンタインターフェース、
３５…ドライブ装置インターフェース、
４１…キーボード、
４２…マウス、
４３…ＣＲＴディスプレイ、
４４…プリンタ、
４５…ＣＤ−ＲＯＭドライブ装置、
４５ａ…ＣＤ−ＲＯＭ、
５０…モデリング対象システム、
５１…入力信号発生器、
５２…出力信号検出器。
【図１】
2008287343.tif 000071
【図２】
2008287343.tif 000072
【図３】
2008287343.tif 000073
【図４】
2008287343.tif 000074
【図５】
2008287343.tif 000075
【図６】
2008287343.tif 000076
【図７】
2008287343.tif 000077
【図８】
2008287343.tif 000078
【図９】
2008287343.tif 000079