モデリング装置および該方法
【課題】本発明は、先験的知識を利用することによって不充分な入出力データであっても、より精度の高いモデルを構築し得るモデリング装置および該方法を提供する。
【解決手段】本発明のモデリング装置Ｓは、システムのモデルをそのパラメータの値を決定することで構築する装置であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られているシステムに関する先験的情報を制約条件として設定する制約条件入力設定部２１および制約条件記憶部４１と、所定の入力データを実績入力データとしてシステムに入力した場合に出力された実績出力データと、前記実績入力データをモデルに入力した場合に出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部１３と、前記制約条件を満たし、かつ、前記評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル更新決定部１４とを備える。
【解決手段】本発明のモデリング装置Ｓは、システムのモデルをそのパラメータの値を決定することで構築する装置であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られているシステムに関する先験的情報を制約条件として設定する制約条件入力設定部２１および制約条件記憶部４１と、所定の入力データを実績入力データとしてシステムに入力した場合に出力された実績出力データと、前記実績入力データをモデルに入力した場合に出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部１３と、前記制約条件を満たし、かつ、前記評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル更新決定部１４とを備える。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、システムのモデルを作成するモデリング装置およびモデリング方法に関し、特に、先験的知識による先験的情報を利用することができるモデリング装置およびモデリング方法に関する。
【背景技術】
【０００２】
近年の計算機（コンピュータ）の発展に伴い、システムの制御系の設計には、コンピュータを用いたシミュレーションが活用されている。実世界のシステムをコンピュータ上の仮想世界に取り入れるためには、実世界のシステムをコンピュータで取り扱い可能な形式で表現した、実世界と仮想世界とを結び付けるモデルで実世界のシステムを記述する必要がある。このモデルは、通常、微分方程式、差分方程式、伝達関数、状態方程式および論理式等の数学モデルである。シミュレーション結果が実世界のシステムを正確に表すものであるか否かは、このモデルの良し悪し（モデルの精度）によって決定され、このモデルを構築する作業は、モデリングと呼ばれる。
【０００３】
モデリング方法は、ホワイトボックスモデリングとブラックボックスモデリングとに大別され、さらに、これらホワイトボックスモデリングとブラックボックスモデリングとを折衷したグレーボックスモデリングがある。
【０００４】
このホワイトボックスモデリングは、いわゆる第１原理に基づくモデリングであり、モデリングの対象であるシステムを支配する物理法則や化学法則等の自然法則（第１原理）と物理化学的な諸量とに基づいてモデリングを行う方法であり、対象が物理システムの場合には物理モデリングと呼ばれる。この方法は、システム（対象）が比較的小規模でこのシステムに作用している物理法則等の自然法則がよく分かっている場合に有効であるが、一般にシステムの大規模化あるいは複雑化に従ってその適用が難しくなる。なお、大規模システムの場合に、この大規模システムを小規模なサブシステムに分割し、各サブシステムに本方法を適用してモデリングを行い、それらを統合することによって大規模システムに本方法を適用することもある。
【０００５】
ブラックボックスモデリングは、物理法則や化学法則等の自然法則や先験的知識を用いることなく、システムをブラックボックスと見なして、システム（対象）で観測された有限個の入出力データから、所定の目的（仕様）に応じたモデルを構築する方法であり、制御工学の分野では、システム同定と呼ばれる。システムが比較的大規模あるいは複雑であり、物理法則等の自然法則のみによってシステムの動特性が記述することができない場合に有効であるが、同じ入出力データが与えられても、モデリングの目的（仕様）やモデリングを行う技術者によって得られるモデルが異なる場合があり、主観的要素が入り込む余地がある。
【０００６】
グレーボックスモデリングは、これらホワイトボックスモデリングとブラックボックスモデリングとの中間に位置づけられる方法であり、例えば、モデルの構築に、観測された入出力データに加えて、物理法則等の自然法則が部分的に利用される。
【０００７】
モデリングの対象となるシステムには、例えば、多関節ロボット（マニピュレータ）等のようにシステムの基本的な挙動が運動方程式によって記述することができ、ホワイトボックスモデリングを採用可能なシステムもあるが、一般的には、例えば熱処理や化学処理等の伴う燃焼炉や化学プラント等のように、その動特性を自然法則から導出し難いケースが多く、システムのモデリング（システムのモデル化）には、ブラックボックスモデルイングやグレーボックスモデルイングを用いることが多い。
【０００８】
このブラックボックスモデリングやグレーボックするモデリングでは、入出力データを用いることによってモデルが構築されるため、モデルの良し悪しがその入出力データに依存し、高精度なモデリングを可能とする質および量ともに充分な入出力データが事前に得られることが望ましい。しかしながら、現実には、実際のシステムに入力可能な波形や大きさ（レベル）等が限られていたり、これら入出力データを採取する実験機会が限られていたり等するために、必ずしも前記質および量ともに充分な入出力データが得られない場合が多い。さらに、実際のシステムから採取された入出力データには、ノイズも存在し、このノイズが場合によっては、モデリングに無視できないほどの大きな悪影響を及ぼすこともあり得る。また、システムによっては、例えば、入力信号に対する持続励振性（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)
【発明の概要】
【発明が解決しようとする課題】
【００１３】
本発明は、上述の事情に鑑みて為された発明であり、その目的は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができるモデリング装置およびモデリング方法を提供することである。
【課題を解決するための手段】
【００１４】
本発明者は、種々検討した結果、上記目的は、以下の本発明により達成されることを見出した。すなわち、本発明の一態様にかかるモデリング装置は、システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング装置であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備えることを特徴とする。
【００１５】
そして、本発明の他の一態様にかかるモデリング方法は、システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング方法であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定ステップと、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算ステップと、前記制約条件設定ステップで設定された制約条件を満たし、かつ、前記評価演算ステップで演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算ステップとを備えることを特徴とする。
【００１６】
このような構成のモデリング装置およびモデリング方法では、実績出力データとモデル出力データとの誤差が評価関数で評価され、この評価結果が最良となるようにモデルのパラメータの値が求められ、モデルが構築される。そして、先験的知識を表す先験的情報が制約条件として設定され、このモデルを構築する際に、制約条件も満たすように、モデルのパラメータの値が求められる。このため、このような構成のモデリング装置およびモデリング方法は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【００１７】
また、他の一態様では、上述のモデリング装置において、前記先験的情報は、前記システムにおける、定常ゲインに関する上限値および下限値のうちの少なくとも一方であることを特徴とする。
【００１８】
この構成によれば、定常ゲインに関する先験的情報も制約条件として加味され、定常ゲインに関する先験的情報を満たすシステムのモデルを構築することができる。
【００１９】
また、他の一態様では、これら上述のモデリング装置において、前記評価関数は、前記実績出力データと前記モデル出力データとの誤差に重み付けを行う重みを有し、前記重みは、その値が大きいほど前記評価結果に大きく影響を与えるものである場合に、前記重みの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されることを特徴とする。
【００２０】
この構成によれば、所定の閾値以上のノイズを含む実績出力データの、評価結果に与える影響が低減され、システムのモデルをより高精度に構築することができる。
【００２１】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、前記モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うことを特徴とする。
【００２２】
一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。この構成によれば、この初期値自体も未知パラメータとして取り扱われるので、初期値に左右されることなく、システムのモデルをより高精度に構築することができる。
【００２３】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄ（ＰＳＯ）によって前記モデルのパラメータの値を求めることを特徴とする。
【００２４】
この構成によれば、モデルのパラメータの値を求める最適化方法としてＰＳＯを用いたモデリング装置を提供することができる。
【００２５】
また、他の一態様では、上述のモデリング装置において、前記モデル演算部は、所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値を前記Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄの初期値に加えることを特徴とする。
【００２６】
この構成によれば、ＰＳＯに用いる初期値として、予め決められた所定の最適化方法によって求められたモデルのパラメータの値も用いるので、より真に近い値を加えた初期値からＰＳＯ法の演算処理を開始することができるから、不必要な探索を行うことなく、システムのモデルをより高精度に構築することができる。
【００２７】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによる更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いることを特徴とする。
【００２８】
この構成によれば、これら定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いてＰＳＯの更新を行うので、不必要な探索を行うことなく、これらのうちの少なくとも１つを制約条件として満たすシステムのモデルを構築することができる。
【発明の効果】
【００２９】
本発明にかかるモデリング装置およびモデリング方法は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【図面の簡単な説明】
【００３０】
【図１】実施形態におけるモデリング装置の構成を示すブロック図である。
【図２】実施形態におけるモデリング装置の概略動作を示すフローチャートである。
【図３】実施形態のモデリング装置における制約条件の設定動作を示すフローチャートである。
【図４】実施形態のモデリング装置におけるモデルの演算動作を示すフローチャートである。
【図５】ゴミ焼却炉において、給塵機の速度をステップ状に変化させた場合における蒸気発生量の変化を示す図である。
【図６】制約条件としてのステップ応答の上限および下限を示す図である。
【図７】実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【図８】比較例としての従来の手法によりモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【図９】定常ゲインの制約条件を考慮した場合であって、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【発明を実施するための形態】
【００３１】
以下、本発明にかかる実施の一形態を図面に基づいて説明する。なお、各図において同一の符号を付した構成は、同一の構成であることを示し、適宜、その説明を省略する。
【００３２】
図１は、実施形態におけるモデリング装置の構成を示すブロック図である。実施形態のモデリング装置Ｓは、モデリングすべきシステムのモデルを前記モデルのパラメータの値を決定することによって構築する装置であって、前記モデルは、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備える装置である。
【００３３】
このようなモデリング装置Ｓは、例えば、図１に示すように、演算部１と、入力設定部２と、出力部３と、記憶部４と、これら演算部１、入力設定部２、出力部３および記憶部４の間で相互にデータを交換することができるようにこれらを互いに接続するバス５とを備えて構成される。
【００３４】
入力設定部２は、モデリング装置Ｓの演算開始指示等の各種コマンドや、例えば実績入出力データ等のシステムのモデル化に必要な各種データをモデリング装置Ｓに入力する機器であり、例えば、キーボードやマウス等である。そして、本実施形態では、入力設定部２は、制約条件を入力設定する制約条件入力設定部２１と、モデル次数および同定パラメータを入力設定するモデル次数同定パラメータ入力設定部２２と、最適化計算パラメータを入力設定する最適化計算パラメータ入力設定部２３と、評価関数を入力設定する評価関数入力設定部２４とを備えている。これら制約条件、モデル次数、同定パラメータおよび最適化計算パラメータは、システムのモデル化の際に利用されるように、モデリング装置Ｓの記憶部４に格納されることで、モデリング装置Ｓに設定される。
【００３５】
制約条件入力設定部２１は、モデルの制約条件（先験的情報）を入力してモデリング装置Ｓに設定するためのものである。制約条件は、システムのモデル化の際に、このモデル化によって得られたモデルが満たすべき条件であり、本実施形態では、モデリング対象のシステムに関する先験的知識によって与えられる。先験的知識は、任意の手段によってシステムから得られる情報であり、例えば、システムの設計値や理論式等から計算によって推定される情報であって良く、また例えば、システムを設計した場合に過去の事例に照らして経験的に推測される情報であって良く、また例えば、システムを実際に組み、この実際のシステムを稼動させてこれを観察することによって得られる情報であって良い。この得られた先験的知識は、モデリング装置Ｓが取り扱うことができる形式で表現（表記、記述）され、先験的情報とされる。
【００３６】
本実施形態では、予め決められた規範（基準）となるデータであって時間に依存する規範入力を入力データとしてシステムに与えた場合に、システムから出力される出力データ（計算によって得られる出力データ、および、システムで観測される実績出力データを含む）が制約条件（時間応答制約条件）とされる。規範入力は、例えば、ステップ状の入力（ステップ入力）、ランプ状の入力（ランプ入力）、三角関数波（例えばｓｉｎ波）や三角波等の周期入力およびそれらの周波数掃引波（スイープ波、ｓｗｅｅｐ波）等の、時間に関する所定の関数である。時間応答制約条件は、例えば前記出力データのそのものであって良く、また例えば前記出力データの上限値および／またはその下限値であって良く、さらにこれらに加えてその時間範囲が設定されて良い。すなわち、規範入力の入力時点を時刻０として、どの時刻ｔＳからどの時刻ｔＥまでの前記出力データが制約条件とされるのか、あるいは、どの時刻ｔＳからどの時刻ｔＥまでの前記出力データにおける上限値および／または下限値が制約条件とされるのか、が設定される（０≦ｔＳ＜ｔＥ）。この時間応答制約条件が入力設定される場合、制約条件入力設定部２１は、時間応答制約条件入力設定部として機能する。
【００３７】
なお、Ａおよび／またはＢは、ＡおよびＢのうちの少なくとも一方の意味であり、ＡあるいはＢ、または、ＡとＢとの両方の意味である。
【００３８】
そして、本実施形態では、システムの定常ゲインの符号や、その範囲（定常ゲインの上限値および下限値）が制約条件（定常ゲイン制約条件）とされる。この定常ゲイン制約条件が入力設定される場合、制約条件入力設定部２１は、定常ゲイン制約条件入力設定部として機能する。
【００３９】
モデル次数未知パラメータ入力設定部２２は、モデルの次数および前記モデルの同定パラメータを入力してモデリング装置Ｓに設定するためのものである。モデルの同定パラメータは、前記モデルにおけるパラメータのうち未知のパラメータとして同定すべきパラメータである。モデルの同定パラメータは、後述のシステムのモデルを作成する演算によってその値が決定され、これによってシステムのモデルが構築される。これらモデルの次数および同定パラメータは、例えば制御系設計に使用するためやシミュレーションに使用するため等の、モデリングの目的によって適宜に選定される。また、モデルの精度、モデル化に使用される入出力データの質や量、および、モデル化の情報処理量（計算量）等も勘案されて良い。
【００４０】
最適化計算パラメータ入力設定部２３は、最適化計算パラメータを入力してモデリング装置Ｓに設定するためのものである。最適化計算パラメータは、最適化の計算を実行する際に必要とされる調整パラメータである。最適化計算パラメータは、本実施形態では、最適化計算手法として、後述するようにＰＳＯ（Particle Swarm Optimization）法が用いられていることから、このＰＳＯ法の演算に必要なパラメータであり、粒子ｐの個数Ｐ、解の更新回数（繰り返し回数）Ｋ、初期値の値域、後述の式におけるｗ、ｃ１、ｃ２である。最適化計算パラメータの入力設定の際に、モデルの次数や同定パラメータの個数が勘案されて良い。例えば、モデルの次数が大きくなるに従って、あるいは、同定パラメータの個数が多くなるに従って、探索すべき解空間が広くなるので、粒子ｐの個数Ｐや解の更新回数Ｋが大きく設定されることが好ましい。
【００４１】
評価関数入力設定部２４は、評価関数を入力してモデリング装置Ｓに設定するためのものである。評価関数は、同定パラメータのパラメータ値を決めた場合に、このパラメータ値によって規定されるモデルがシステムをどの程度よく表現しているか否かの度合い（モデルの精度）を示す値（評価値）を求める関数である。評価関数は、予め１つだけ設定され、プログラムに予め組み込まれていても良いが、本実施形態では、より適切にモデルを評価するために、評価関数入力設定部２４から入力可能とされている。さらに、オペレータの便宜を図るために、本実施形態では、予め複数の評価関数が用意されており、これら複数の評価関数の中から評価関数入力設定部２４を介して選定することで、評価関数が入力設定される。
【００４２】
出力部３は、入力部２から入力されたコマンドやデータおよびモデリング装置Ｓの演算結果（システムのモデル）等を出力する機器であり、例えばＣＲＴディスプレイ、ＬＣＤ、有機ＥＬディスプレイ又はプラズマディスプレイ等の表示装置やプリンタ等の印字装置等である。
【００４３】
記憶部４は、データやプログラムを記憶する例えばＲＯＭ（Read Only Memory）およびＥＥＰＲＯＭ（Electrically Erasable Programmable Read Only Memory）等の不揮発性の記憶素子やハードディスク等を備えるとともに、演算部１が実行するモデリングプログラムや制御プログラム等の実行中における各データを一時的に記憶するいわゆるワーキングメモリであり、例えば揮発性の記憶素子であるＲＡＭ（Random Access Memory）を備えて構成される。記憶部４は、機能的に、制約条件記憶部４１と、モデル次数同定パラメータ記憶部４２と、最適化計算パラメータ記憶部４３と、評価関数記憶部４４と、モデル出力等記憶部４５と、実績データ記憶部４６とを備える。
【００４４】
制約条件記憶部４１は、制約条件を記憶するものであり、制約条件入力設定部２１から入力された制約条件が制約条件記憶部４１に格納され、記憶される。制約条件記憶部４１は、本実施形態では、時間応答制約条件を記憶する時間応答制約条件記憶部４１１と、定常ゲイン制約条件を記憶する定常ゲイン制約条件記憶部４１２とを備えている。
【００４５】
モデル次数同定パラメータ記憶部４２は、モデルの次数および前記モデルの同定パラメータを記憶するものであり、モデル次数同定パラメータ入力設定部２２から入力されたモデルの次数および同定パラメータがモデル次数同定パラメータ記憶部４２に格納され、記憶される。
【００４６】
最適化計算パラメータ記憶部４３は、最適化計算パラメータを記憶するものであり、最適化計算パラメータ入力設定部２３から入力された最適化計算パラメータが最適化計算パラメータ記憶部４３に格納され、記憶される。
【００４７】
評価関数記憶部４４は、評価関数を記憶するものであり、評価関数入力設定部２４から入力された評価関数が評価関数記憶部４４に格納され、記憶される。本実施形態では、上述したように、複数の所定の評価関数が予め記憶されており、これら複数の評価関数の中から１つの評価関数が評価関数入力設定部２４から選定され、この選定されたことを示す情報も評価関数記憶部４４に記憶される。
【００４８】
モデル出力等記憶部４５は、モデルを規定する所定のパラメータ値に対し、演算部１で演算されたモデル出力、評価結果および制約条件を満たすか否かの情報を互いに対応付けて記憶するものである。前記所定のパラメータ値は、本実施形態では、ＰＳＯ法によって与えられる値である。そして、本実施形態では、制約条件を満たすか否かの情報は、評価結果を示す評価値に変えられ、評価結果としてモデル出力等記憶部４５に記憶される。このように制約条件を満たすか否かの情報が評価結果を示す評価値に変換されて評価結果としてモデル出力等記憶部４５に記憶されることで、モデルのパラメータの値を決定する際に、評価結果を参酌するだけで済み、このモデルのパラメータの値を決定する処理が簡便化される。
【００４９】
実績データ記憶部４６は、実績データを記憶するものである。実績データは、モデリング対象のシステムに入力（印加）される所定の入力データ（実績入力データ）、および、この入力データの入力によってシステムから出力される出力データ（実績出力データ）であり、これら実績入力データおよび実績出力データ（実績入出力データ）は、互いに対応付けて実績データ記憶部４３に記憶される。前記所定の入力データは、任意のデータであって良く、例えば、ステップ状の入力データ、ランダム入力データおよび擬似白色系列（例えばＭ系列等）データ等の時系列データである。出力データもこれに対応して時系列データとなる。前記所定の入力データは、制御工学の分野では、同定入力と呼ばれる。時系列入力データおよび時系列出力データは、１組であっても良く、また複数の組であっても良い。複数の組である場合には、時系列入出力データの各組は、同日に得られたデータであっても良く、また、異なる日に得られたデータであって良い。そして、これら複数の組の時系列入出力データは、モデリングの際に、これらの中から適宜に選択されるようにモデリング装置Ｓが構成されてもよい。
【００５０】
演算部１は、例えば、マイクロプロセッサおよびその周辺回路等を備えて構成され、プログラムに従い入力設定部２、出力部３および記憶部４を当該機能に応じてそれぞれ制御することによってモデリング装置Ｓ全体の制御を司るものであって、モデリング対象のシステムのモデルを前記モデルのパラメータの値を決定することによって構築するものである。演算部１は、機能的に、モデル出力演算部１１と、制約条件評価部１２と、評価演算部１３と、モデル更新決定部１４とを備える。
【００５１】
モデル出力演算部１１は、モデル次数、同定パラメータおよび最適化計算パラメータを設定することによって設定された前記システムのモデルに、入力データを入力（印加）することで前記モデルの出力データ（モデル出力データ）を求めるものである。この求めたモデル出力データは、モデル出力等記憶部４５に記憶され、格納される。
【００５２】
制約条件評価部１２は、モデル出力演算部１１によって得られたモデル出力データが制約条件を満たしているか否かを評価するものである。この評価したモデル出力データが制約条件を満たしているか否かの情報は、モデル出力等記憶部４５に記憶され、格納される。
【００５３】
評価演算部１３は、予め決められた所定の入力データを実績入力データとしてシステムに入力した場合にこのシステムから出力された実績出力データと、前記実績入力データを前記システムのモデルに入力した場合にこのモデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求めるものである。この求めた評価結果は、モデル出力等記憶部４５に記憶され、格納される。
【００５４】
モデル更新決定部１４は、前記システムのモデルを予め決められた所定の更新則に従って更新するとともに、制約条件入力設定部２１で設定された制約条件を満たし、かつ、評価演算部１３で演算された評価結果が最良となるように、前記システムのモデルにおけるそのパラメータの値を求め、前記システムのモデルを決定するものである。
【００５５】
このようなモデリング装置Ｓは、モデリング方法のプログラム（ソフトウェア）を実装することで、例えば、デスクトップ型パソコンやノート型パソコン等のコンピュータによって構成可能である。
【００５６】
なお、必要に応じてモデリング装置Ｓは、図１に破線で示す外部記憶部６や通信インターフェース部７をさらに備えてもよい。
【００５７】
外部記憶部６は、例えば、フレキシブルディスク、ＣＤ−ＲＯＭ（Compact Disc Read Only Memory）、ＣＤ−Ｒ（Compact Disc Recordable）およびＤＶＤ−Ｒ（Digital Versatile Disc Recordable）等の記録媒体との間でデータを読み込みおよび／または書き込みを行う装置であり、例えば、フレキシブルディスクドライブ、ＣＤ−ＲＯＭドライブ、ＣＤ−Ｒドライブ及びＤＶＤ−Ｒドライブ等である。通信インターフェース部７は、通信網（ネットワーク）に接続され、この通信網を介して他の装置との間で通信信号を送受信するための機器である。
【００５８】
各プログラムが格納されていない場合には、これらを記録した記録媒体から外部記憶部６を介して記憶部４にインストールされるように構成してもよく、また、これらプログラムを管理するサーバコンピュータ（不図示）から通信網および通信インターフェース部７を介して各プログラムがダウンロードされるように構成してもよい。また、システムのモデル化の演算に当たってモデリング装置Ｓに入力すべき例えば実績入出力データ等のデータは、このデータを記憶した記録媒体によって外部記憶部６を介してモデリング装置Ｓに入力されるように構成してもよく、また、このデータを管理するサーバコンピュータから通信網および通信インターフェース部７を介してモデリング装置Ｓに入力されるように構成してもよい。
【００５９】
次に、本実施形態におけるモデリング装置の動作について説明する。図２は、実施形態におけるモデリング装置の概略動作を示すフローチャートである。図３は、実施形態のモデリング装置における制約条件の設定動作を示すフローチャートである。図４は、実施形態のモデリング装置におけるモデルの演算動作を示すフローチャートである。
【００６０】
本実施形態のモデリング装置Ｓでは、概略、図２に示すように、まず、ステップＳ１１では、制約条件が設定され、続いて、ステップＳ１２では、モデルを構築すべく、制約条件を満たし、かつ、評価結果が最良となるように、モデルのパラメータの値を決定することによって、モデルが演算される。続いて、ステップＳ１３では、この演算されたモデルの精度が予め決められた許容範囲以内であるか否かが判定され、この判定の結果、モデルの精度が許容範囲以内ではない場合（ＮＯ）には処理がステップＳ１１に戻され、一方、モデルの精度が許容範囲内である場合（ＹＥＳ）にはステップＳ１４が実行され、この演算されたモデルがシステムのモデルとして出力部３に出力され、出力部３で出力される。
【００６１】
次に、前記ステップＳ１１についてより具体的に説明する。図３において、まず、ステップＳ２１では、入力設定部２の制約条件入力設定部２１から、時間応答波形における制約条件（前記時間応答制約条件）が入力され、記憶部４の制約条件記憶部４１における時間応答制約条件記憶部４１１に格納され、記憶される。これによって、後述するように、演算部１の制約条件評価部１２が時間応答制約条件記憶部４１１に記憶されている時間応答制約条件を用いることで、この時間応答制約条件を満たすように、システムのモデルが構築される。
【００６２】
以下に、説明の簡単化のため、１入力１出力系のシステムの場合について、より具体的に説明する。もちろん、以下の説明は、一般性を損なうことなく、多入力多出力系のシステムに拡張可能である。
【００６３】
まず、例えば燃焼炉プラント等の実際のシステムから入出力データ（入力データＵｄａｔａおよび出力データＹｄａｔａ）が実績入出力データ（実績入力データＵｄａｔａおよび実績出力データＹｄａｔａ）として取得され、実績データ記憶部４６に格納され、記憶されているものとする。システムのモデルが連続した時間と見なされる時間変数ｔを含む連続時間状態空間モデル式または連続時間伝達関数モデル式であることから、実績入出力データＵｄａｔａ、Ｙｄａｔａは、時間的に連続しているデータに見えるように、より短いサンプリング間隔で取得されることが好ましい。ここでは、表記の簡単化のために、１秒間隔（１秒周期）でサンプリングされているものとする。サンプリング開始の時刻０から時刻ｔｍａｘまでのｔｍａｘ秒間における実績入出力データＵｄａｔａ、Ｙｄａｔａは、時刻ｔ秒の入力値をｕ（ｔ）とし、時刻ｔ秒の出力値をｙ（ｔ）とする場合には、それぞれ、次の式１−１、式１−２となる。
Ｕｄａｔａ＝［ｕ（０） ｕ（１） ・・・ ｕ（ｔｍａｘ）］Ｔ ・・・（１−１）
Ｙｄａｔａ＝［ｙ（０） ｙ（１） ・・・ ｙ（ｔｍａｘ）］Ｔ ・・・（１−２）
ここで、上付き添え字Ｔは、いわゆる転置を表す記号である。
【００６４】
時刻０から時刻ｔＥまでのｔＥ秒間における式２の規範入力Ｕｒをシステムに加えた場合に得られる出力データに関し、時刻ｔＳ１から時刻ｔＥ１までの間における上限Ｙｍａｘおよび時刻ｔＳ２から時刻ｔＥ２までの間における下限Ｙｍｉｎを式３−１および式３−２とすると、これら式３−１によって記述される上限Ｙｍａｘおよび式３−２によって記述される下限Ｙｍｉｎが制約条件入力設定部２１からモデリング装置Ｓに与えられる。
Ｕｒ＝［ｕｒ（０） ｕｒ（１） ・・・ ｕｒ（ｔE）］Ｔ ・・・（２）
Ｙｍａｘ＝［ｙｍａｘ（ｔｓ１） ｙｍａｘ（ｔｓ１＋１） ・・・ ｙｍａｘ（ｔＥ１）］Ｔ ・・・（３−１）
Ｙｍｉｎ＝［ｙｍｉｎ（ｔｓ２） ｙｍｉｎ（ｔｓ２＋１） ・・・ ｙｍｉｎ（ｔＥ２）］Ｔ ・・・（３−２）
ここで、ｕｒ（ｔ）は、時刻ｔ秒における入力値であり、ｙｍａｘ（ｔ）は、時刻ｔにおける出力値の上限値であり、そして、ｙｍｉｎ（ｔ）は、時刻ｔにおける出力値の下限値である。また、０≦ｔＳ１＜ｔＥ１≦ｔＥであり、０≦ｔＳ２＜ｔＥ２≦ｔＥである。
【００６５】
なお、規範入力の一具体例として、単位ステップ入力は、Ｕｒ＝［０ １ １ ・・・ １］Ｔであり、単位インパルス入力は、Ｕｒ＝［０ １ ０ ・・・ ０］Ｔであり、またランプ入力は、Ｕｒ＝［０ １ ２ ３ ・・・］Ｔである。また、規範入力は、１種類であって良く、また複数の種類であっても良い。
【００６６】
続いて、ステップＳ２２では、入力設定部２の制約条件入力設定部２１から、定常ゲインにおける制約条件（前記定常ゲイン制約条件）が入力され、記憶部４の制約条件記憶部４１における定常ゲイン制約条件記憶部４１２に格納され、記憶される。例えば、システムの定常ゲインＧがＧｍｉｎ≦Ｇ≦Ｇｍａｘとすると、これら下限値Ｇｍｉｎおよび上限値Ｇｍａｘが制約条件入力設定部２１からモデリング装置Ｓに与えられる。これによって、後述するように、演算部１の制約条件評価部１２が定常ゲイン制約条件記憶部４１２に記憶されている定常ゲイン制約条件を用いることで、この定常ゲイン制約条件を満たすように、システムのモデルが構築される。
【００６７】
次に、前記ステップＳ１２についてより具体的に説明する。図４において、まず、ステップＳ３１では、入力設定部２のモデル次数同定パラメータ入力設定部２２から、モデルの次数および前記モデルの同定パラメータが入力され、記憶部４のモデル次数同定パラメータ記憶部４２に格納され、記憶される。
【００６８】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、例えば、モデリング対象のシステムに対する連続時間状態空間モデル式が式４（式４−１、式４−２）であるとする。
ｄｘ（ｔ）／ｄｔ＝Ａｘ（ｔ）＋Ｂｕ（ｔ） ・・・（４−１）
ｙｍ（ｔ）＝Ｃｘ（ｔ）＋Ｄｕ（ｔ） ・・・（４−２）
ここで、ｘ（ｔ）は、ｎ次の状態ベクトルであり、ｙｍ（ｔ）は、モデルの出力（モデル出力、モデルによる計算値）であり、Ａ、Ｂ、ＣおよびＤは、それぞれｎ行ｎ列、ｎ行１列、１行ｎ列および１行１列の定数行列である。ここでは、説明の簡単化のため、モデリング対象のシステムが厳密にプロパーである、すなわち、Ｄ＝０とする。
【００６９】
モデルの次数は、このｎの値であり、この値ｎがモデルの次数としてモデル次数同定パラメータ入力設定部２２から入力される。この例では、ｎ＝２、すなわち、モデルの次数が２である２次のモデルを求めるものとする。また、定数行列Ａ、ＢおよびＣは、例えば可制御正準形式等の所定の形式でよいが、冗長なパラメータが無い方が好ましいことから、この例では、可観測正準形式であるとし、それぞれ式５−１、式５−２および式５−３となる。
【００７０】
【数１】
【００７１】
ここで、ａ０、ａ１、ｂ０およびｂ１は、モデルを規定する未知パラメータである。これらａ０、ａ１、ｂ０およびｂ１の未知パラメータが、モデル次数同定パラメータ入力設定部２２から入力され、同定パラメータとして指定される。すなわち、モデリング対象のシステムのモデルを構築することは、これら同定パラメータとして指定された未知パラメータａ０、ａ１、ｂ０およびｂ１の各値を求めることである。
【００７２】
なお、式４（式４−１、式４−２）および式５（式５−１〜式５−３）の求解は、式６の連続時間伝達関数モデル式Ｇ（ｓ）を求めること、すなわち、式６の連続時間伝達関数モデル式Ｇ（ｓ）における未知パラメータａ０、ａ１、ｂ０およびｂ１を求めることと等価である。
Ｇ（ｓ）＝（ｂ１ｓ＋ｂ０）／（ｓ２＋ａ１ｓ＋ａ０） ・・・（６）
【００７３】
続いて、ステップＳ３２では、入力設定部２の最適化計算パラメータ入力設定部２３から、最適化計算パラメータが入力され、記憶部４の最適化計算パラメータ記憶部４３に格納され、記憶される。
【００７４】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、未知パラメータ（同定パラメータ）ａ０、ａ１、ｂ０およびｂ１の同定は、本実施形態では、例えば、ＰＳＯ（Particle Swarm Optimization）法が用いられる。このＰＳＯ法は、解空間に複数Ｐ個の粒子ｐ（粒子群）をランダム（無作為）にばらまき、これら解空間にばらまかれた粒子群が所定の更新則に従って解空間内を移動しながら、最適解を探索する手法である。このＰＳＯ法を本実施形態に適用する場合では、解空間は、未知パラメータ（同定パラメータ）ａ０、ａ１、ｂ０およびｂ１を各座標軸とする座標空間（この例では４次元空間）であり、粒子ｐは、４次元ベクトルで表される前記座標空間上の点となり、更新回数ｋ番目におけるＰ個のうちのｉ番目の粒子ｐ（ｉ，ｋ）は、式７によって表現されることになる。
ｐ（ｉ，ｋ）＝［ａ０ ａ１ ｂ０ ｂ１］Ｔ ・・・（７）
【００７５】
以下、適宜、このｐ（ｉ，ｋ）の［ａ０ ａ１ ｂ０ ｂ１］Ｔをパラメータセットと呼称することとする。
【００７６】
このような最適化手法では、粒子ｐの個数Ｐ（０＜ｉ≦Ｐ）、更新回数Ｋ（０＜ｋ≦Ｋ）および所定の更新則におけるパラメータ（後述のｗ、ｃ１、ｃ２）が最適化計算パラメータであり、最適化計算パラメータ入力設定部２３から入力される。
【００７７】
続いて、ステップＳ３３では、入力設定部２から評価関数が入力され、設定される。
【００７８】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、評価関数は、例えば、本実施形態では、実績入力データＵｄａｔａに対応する実績出力データＹｄａｔａ（＝［ｙ（０） ｙ（１） ・・・ ｙ（ｔｍａｘ）］Ｔ）と、前記実績入力データＵｄａｔａをモデル（上述の例では、式４（式４−１、式４−２）および式５（式５−１〜式５−３）または式６）に入力した場合にこのモデルから出力されるモデル出力データＹｉ（＝［ｙｉ（０） ｙｉ（２） ・・・ ｙｉ（ｔｍａｘ）］Ｔ）との誤差を評価する関数である。このような評価関数Ｊｉは、例えば、式８−１によって与えられる。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜ ・・・（８−１）
ここで、｜｜ｖ｜｜は、いわゆる、ベクトルｖのユークリッドノルムを表し、ｖが縦ベクトルである場合では、｜｜ｖ｜｜＝√（ｖＴｖ）である。また、Ｗは、出力の予測誤差に対していずれの部分を重要視して評価するかを定める重み行列であり、例えば対角行列とする。
【００７９】
また例えば、評価関数Ｊｉは、例えば、式８−２によって与えられる無限ノルムである。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜∞ ・・・（８−２）
ここで、｜｜ｖ｜｜∞は、いわゆる、ベクトルｖの無限ノルム（最大値ノルム）を表し、ｖ＝［ｖ１ ｖ２ ・・・ ｖｎ］Ｔである場合では、｜｜ｖ｜｜∞＝ｍａｘ（｜ｖ１｜、｜ｖ２｜、・・・、｜ｖｎ｜）である。
【００８０】
また例えば、評価関数Ｊｉは、例えば、式８−３によって与えられる。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜１ ・・・（８−３）
ここで、｜｜ｖ｜｜１は、いわゆる、ベクトルｖの１ノルムを表し、ｖ＝［ｖ１ ｖ２ ・・・ ｖｎ］Ｔである場合では、｜｜ｖ｜｜１＝Σ｜ｖｉ｜である。ただし、Σは、ｉ＝１からｉ＝ｎまでの和である。
【００８１】
評価関数記憶部４４には、例えば、このような複数の評価関数Ｊｉ（式８−１〜式８−３）が記憶されており、これら複数の評価関数Ｊｉの中から、評価関数入力設定部２４を介して１つの評価関数Ｊｉが選定され、入力設定される。
【００８２】
なお、このような評価関数Ｊが、実績出力データとモデル出力データとの誤差に重み付けを行う重みＷを有し、この重みＷが、その値が大きいほど評価関数Ｊによる評価結果に大きく影響を与えるものである場合に、この重みＷの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されてよい。このような前記区間およびその重みＷならびに前記他の区間およびその重みＷは、例えば、入力設定部２の最適化計算パラメータ入力設定部２３から入力設定される。このように構成することによって、所定の閾値以上のノイズを含む実績出力データの、評価結果に与える影響が低減され、システムのモデルをより高精度に構築することができる。
【００８３】
続いて、ステップＳ３４では、演算部１のモデル出力演算部１１によって、モデル出力が演算され、記憶部４のモデル出力等記憶部４５に記憶される。
【００８４】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、まず、Ｐ個の粒子ｐ（ｉ，０）（更新回数ｋ＝０、ｉ＝０〜Ｐ）、例えば１００個の粒子ｐ（ｉ，０）（ｉ＝０〜１００）の初期値が例えば乱数を生成することによってランダムに設定される。これによってＰ個の各パラメータセットｐ（ｉ，０）の値が設定される。そして、ｉ番目の粒子ｐ（ｉ，０）の移動速度ｖ（ｉ，０）の初期値が例えばｖ（ｉ，０）＝［０ ０ ０ ０］Ｔに設定される。
【００８５】
ここで、制約条件として、モデリング対象のシステムが安定な系であることが、先験的知識として判明している場合には、ａ０＞０かつａ１＞０の範囲で、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８６】
また、制約条件として、モデリング対象のシステムにおける定常ゲインが正の値であることが、先験的知識として判明している場合には、ｂ０の符号とａ０の符号が等しくなるように、すなわち、ｂ０×ａ０＞０となるように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８７】
また、制約条件として、モデリング対象のシステムに不安定なゼロ点がないことが、先験的知識として判明している場合には、ｂ１の符号とｂ０の符号が等しくなるように、すなわち、ｂ１×ｂ０＞０となるように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８８】
また、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下であることが、先験的知識として判明している場合には、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８９】
このようにｉ番目の粒子ｐ（ｉ，０）の初期値を設定する際に、制約条件を加味することによって、解が存在しない領域を探索する無駄を省くことが可能となる。
【００９０】
このようにＰ個の粒子ｐの値をそれぞれ設定することで、１個の粒子ｐに対応して１個のモデルが規定される。すなわち、ｉ番目の粒子ｐ（ｉ，０）に対応して１個のモデルが規定される。上述の例では、１００個の各粒子ｐのそれぞれに対応して１００個のモデルが規定されることになる。
【００９１】
このｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデル（この例では式４および式５によって表現されるモデル、あるいは、式６によって表現されるモデル）に対し、実績データ記憶部４６に記憶されている実績入力データＵｄａｔａが代入され、モデル出力データＹｉが得られ、この得られたモデル出力データＹｉが記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このモデル出力データＹｉの計算には、公知の制御系シミュレーションソフトウェア、例えば、Ｔｈｅ Ｍａｔｈ Ｗｏｒｋ社のＭＡＴＬＡＢ／Ｓｉｍｕｌｉｎｋ等が用いられる。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、モデル出力が演算され、記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９２】
続いて、ステップＳ３５では、制約条件評価部１２によって、制約条件が評価され、この評価結果が記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９３】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ステップＳ３４と同様に、このｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデルに対し、所定の規範入力Ｕｒが代入され、前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉ（＝［ｙｒｉ（０） ｙｒｉ（１） ・・・ ｙｒｉ（ｔｍａｘ）］Ｔ）が得られる。そして、この得られた前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが規範入力に対する制約条件を満たしているか否かが判断される。例えば、制約条件として、時刻ｔＳ１から時刻ｔＥ１までの間における上限値Ｙｍａｘおよび時刻ｔＳ２から時刻ｔＥ２までの間における下限値Ｙｍｉｎが設定されている場合には、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが時刻ｔＳ１から時刻ｔＥ１までの間で上限Ｙｍａｘ以下で、かつ、時刻ｔＳ２から時刻ｔＥ２までの間で下限Ｙｍｉｎ以上であるか否かが判断される。この判断の結果、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが時刻ｔＳ１から時刻ｔＥ１までの間で上限Ｙｍａｘ以下で、かつ、時刻ｔＳ２から時刻ｔＥ２までの間で下限Ｙｍｉｎ以上ではない場合には、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが制約条件を満たさないものとされ、この制約条件の評価結果が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、制約条件が評価され、この評価結果が記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９４】
なお、上述では、粒子ｐ（ｉ，０）の初期値の設定の際に、定常ゲインＧの制約条件を考慮したが、定常ゲインＧの下限値Ｇｍｉｎおよび上限値Ｇｍａｘが制約条件として設定されている場合に、実績入力データＵｄａｔａに対応するモデル出力データＹｉが下限値Ｇｍｉｎ以上であって、かつ、上限値Ｇｍａｘ以下であるか否かが判断されることで、定常ゲインＧの制約条件が評価されても良い。
【００９５】
続いて、ステップＳ３６では、評価演算部１３によって、評価関数による評価値が演算され、記憶部４のモデル出力等記憶部４５に記憶される。
【００９６】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデルに対し、前記評価関数Ｊｉによる評価値が演算され、この評価関数Ｊｉの評価値が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、評価関数による評価値が演算され、記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９７】
ここで、上記制約条件の評価結果を記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶する際に、上述のように、この制約条件の評価結果そのものを記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶してもよいが、本実施形態では、この制約条件の評価結果が評価関数Ｊｉの評価値で記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。例えば、制約条件を満たさない場合に、評価が悪いことを示す評価値が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。より具体的には、例えば、評価関数Ｊｉが、評価が高いほどその評価値が小さくなるように構成されている場合では、制約条件を満たさない場合に、例えばモデル出力データＹｉがその値として取り得ない大きな評価値が、記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。また例えば、評価関数Ｊｉが、その取り得る値の範囲が負であるように構成されている場合では、制約条件を満たさない場合に、正の値（例えば固定値や制約条件の外れ度合いに従った値等）が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このような評価関数Ｊｉは、例えば、式８−４によって与えられる。
Ｊｉ＝ｔａｎ−１（｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜）−π／２ ・・・（８−４）
【００９８】
続いて、ステップＳ３７では、モデル更新決定部１４によって、最適化方法の繰り返し計算における更新が終了しているか否かが判断され、この判断の結果、更新が終了していない場合（ＮＯ）では、ステップＳ３８の処理が実行され、一方、更新が終了している場合（ＹＥＳ）では、ステップＳ３９の処理が実行される。更新回数Ｋは、任意の回数で良いが、ＰＳＯ法では、解の精度や演算処理量等を考慮することによって適宜な回数に設定され、例えば、Ｋ＝１００回とされる。
【００９９】
ステップＳ３８では、モデル更新決定部１４によって、所定の更新則に従って更新され、処理がステップＳ３４に戻される。
【０１００】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ｋ回目（ｋ＝０，１，２，・・・、Ｋ）の更新時点において、Ｐ個の各粒子ｐの中からその評価関数Ｊの値が最も小さい粒子ｐが探索され、第１最良評価粒子ｐ（ｂｅｓｔ，ｋ）とされる。例えば、ｋ回目の更新時点において、ｑ番目の粒子ｐ（ｑ、ｋ）の評価関数Ｊｑの値が最も小さい場合には、ｐ（ｂｅｓｔ，ｋ）＝ｐ（ｑ，ｋ）である。
【０１０１】
また、ｉ番目の粒子ｐにおいて、０回目（初期値）からｋ−１回目までの中で評価関数Ｊの値が最も小さい粒子ｐが探索され、第２最良評価粒子ｐ（ｉ，ｂｅｓｔ）とされる。例えば、ｒ回目の更新時点において、ｉ番目の粒子ｐ（ｉ、ｒ）の評価関数Ｊｉの値が最も小さい場合には、ｐ（ｉ，ｂｅｓｔ）＝ｐ（ｉ，ｒ）である。なお、ｋ＝０の場合では、それぞれの初期値そのものが各粒子ｐのｐ（ｉ，ｂｅｓｔ）である。
【０１０２】
そして、例えば、式９（式９−１、式９−２）によって与えられる更新則に従って、各粒子ｐの速度および位置がそれぞれ決定され、更新される。
ｖ（ｉ，ｋ＋１）＝ｗｖ（ｉ，ｋ）
＋ｃ１Ｒ１ｉｋ（ｐ（ｉ，ｂｅｓｔ）−ｐ（ｉ，ｋ））
＋ｃ２Ｒ２ｉｋ（ｐ（ｂｅｓｔ，ｋ）−ｐ（ｉ，ｋ）） ・・・（９−１）
ｐ（ｉ，ｋ＋１）＝ｐ（ｉ，ｋ）＋ｖ（ｉ，ｋ＋１） ・・・（９−２）
ここで、ｗ、ｃ１、ｃ２は、それぞれ、調整パラメータであり、最適化計算パラメータとして上述したように最適化計算パラメータ入力設定部２３から入力され、設定される。例えば、ｗ＝０．９、ｃ１＝０．８およびｃ２＝０．８等である。また、Ｒ１ｉｋおよびＲ２ｉｋは、それぞれ対角行列であり、対角要素のそれぞれは、１〜０の範囲で一様に分布する乱数である。この乱数値は、固定値ではなく、各粒子ｐごとに、かつ、各更新ｋごとに、それぞれランダムな値が使用される。
【０１０３】
ここで、制約条件として、モデリング対象のシステムが安定な系であることが、先験的知識として判明している場合には、ａ０＞０かつａ１＞０の範囲に入るように更新することが好ましい。より具体的には、更新後の粒子ｐがａ０＜０となった場合には、式１０−１に示すように、更新前のａ０と０とを内分する点に更新される。すなわち、０〜１の乱数αにおいて、
ａ０（ｉ，ｋ＋１）＝（１−α）・ａ０（ｉ，ｋ）＋α・０ ・・・（１０−１）
ここで、ｉ番目の粒子ｐについて、ａ０（ｉ，ｋ＋１）は、更新後のａ０であり、ａ０（ｉ，ｋ）は、更新前のａ０である。また、粒子ｐの速度ｖは、ａ０に対応する部分ｖ（ｉ，ｋ＋１）｛ａ０｝について、式１０−２となる。
ｖ（ｉ，ｋ＋１）｛ａ０｝＝ａ０（ｋ＋１）−ａ０（ｋ） ・・・（１０−２）
【０１０４】
また、制約条件として、モデリング対象のシステムにおける定常ゲインが正の値であることが、先験的知識として判明している場合には、ｂ０の符号とａ０の符号が等しくなるように、すなわち、ｂ０×ａ０＞０となるように、更新することが好ましい。
【０１０５】
また、制約条件として、モデリング対象のシステムに不安定なゼロ点がないことが、先験的知識として判明している場合には、ｂ１の符号とｂ０の符号が等しくなるように、すなわち、ｂ１×ｂ０＞０となるように、更新することが好ましい。
【０１０６】
また、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下の範囲であることが、先験的知識として判明している場合には、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、更新することが好ましい。
【０１０７】
このように粒子ｐを更新する際に、制約条件を加味することによって、解が存在しない領域を探索する無駄を省くことが可能となり、よりよい解をより早く探索することができる。
【０１０８】
さらに、制約条件として、同定パラメータの符号ａ０、ａ１、ｂ０、ｂ０１が先験的知識として判明している場合には、これら同定パラメータが指数形式で表現されても良い。例えば、式６の伝達関数が安定かつ不安定零点なしかつゲインが負であることが先験的知識として判明している場合では、ａ０の符号およびａ１の符号が正であって、かつ、ｂ０の符号およびｂ１の符号が負でなければならないことから、ａ０＝ｅｘｐ（λ０）、ａ１＝ｅｘｐ（λ１）、ｂ０＝−ｅｘｐ（μ０）、ｂ１＝ｅｘｐ（μ１）とおき、式７のパラメータセットがｐ（ｉ，ｋ）＝［λ０ λ１ μ０ μ１］Ｔとなる。したがって、同定パラメータ［ａ０ ａ１ ｂ０ ｂ１］Ｔの座標空間ではなく、このパラメータ［λ０ λ１ μ０ μ１］Ｔの座標空間で解が探索される。このように構成することによって、−∞＜λ０＜＋∞、−∞＜λ１＜＋∞、−∞＜μ０＜＋∞、−∞＜μ１＜＋∞の実数全体での演算処理が可能となり、符号や範囲をチェックすることなく演算処理が可能となる。
【０１０９】
このようにＰＳＯ法による更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つが用いられてよい。このようにモデリング装置Ｓを構成することによって、不必要な探索を行うことなく、これらのうちの少なくとも１つを制約条件として満たすシステムのモデルを構築することができる。
【０１１０】
こうしてＰ個の各パラメータセットｐ（ｉ，ｋ）が所定の更新則に従って更新され、処理がステップＳ３４に戻される。
【０１１１】
一方、ステップＳ３９では、モデル更新決定部１４によって、制約条件を満たし、かつ最良の評価値を与えるように、同定パラメータセットの値が決定され、このシステムのモデルが構築される。例えば、上述のＰＳＯ法では、Ｐ×Ｋ個の各パラメータセットｐ（ｉ，ｋ）の中から、制約条件を満たし、かつ最良の評価値を与えるパラメータセットｐ（ｉ，ｋ）が選択され、この選択されたパラメータセットｐ（ｉ，ｋ）の値が、モデリング対象のシステムのモデルを規定する同定パラメータの値として決定され、このシステムのモデルが決定される。
【０１１２】
このように本実施形態におけるモデリング装置Ｓおよびモデリング方法では、実績出力データとモデル出力データとの誤差が評価関数Ｊで評価され、この評価結果が最良となるようにモデルのパラメータの値が求められ、モデルが構築される。そして、先験的知識を表す先験的情報が制約条件として設定され、このモデルを構築する際に、制約条件も満たすように、モデルのパラメータの値が求められる。このため、本実施形態におけるモデリング装置Ｓおよびモデリング方法は、上述の非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【０１１３】
また、本実施形態におけるモデリング装置Ｓおよびモデリング方法では、システムにおける、定常ゲインＧに関する上限値Ｇｍａｘおよび下限値Ｇｍｉｎのうちの少なくとも一方が、本実施形態では定常ゲインＧに関する上限値Ｇｍａｘおよび下限値Ｇｍｉｎが制約条件とされている。このため，本実施形態におけるモデリング装置Ｓおよびモデリング方法は、定常ゲインＧに関する先験的情報も制約条件として加味され、定常ゲインＧに関する先験的情報を満たすシステムのモデルを構築することができる。
【０１１４】
次に、一実施例として、実施形態におけるモデリング装置Ｓを用いて例えばごみ焼却炉（ガス化溶融炉）のモデルを作成した場合について、説明する。
【０１１５】
図５は、ゴミ焼却炉において、給塵機の速度をステップ状に変化させた場合における蒸気発生量の変化を示す図である。図５（Ａ）は、給塵機の速度を示し、その横軸は、秒単位で表す経過時間を示し、その縦軸は、速度を示す。図５（Ｂ）は、蒸気発生量を示し、その横軸は、秒単位で表す経過時間を示し、その縦軸は、蒸気発生量を示す。図６は、制約条件としてのステップ応答の上限および下限を示す図である。図６の横軸は、秒単位で表す経過時間を示し、その縦軸は、蒸気発生量を示す。図７は、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図７（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図７（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。図８は、比較例としての従来の手法によりモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図８（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図８（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。図９は、定常ゲインの制約条件を考慮した場合であって、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図９（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図９（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。
【０１１６】
ゴミ焼却炉では、ゴミを焼却炉に供給する給塵機を入力とし、焼却炉で発生する蒸気発生量を出力とすると、図５（Ａ）に示すように、給塵機の速度がステップ状に変化させた場合に、蒸気発生量は、図５（Ｂ）に示すように、基本的に、給塵機の稼動開始に伴って徐々に増大するとともに給塵機の稼動停止に伴って徐々に減少するプロファイルであるが、このプロファイルに外乱の影響によりノイズが重畳しており、その波形が大きく乱れている。この原因として、例えば、給塵機の速度が一定であっても給塵機によって実際に焼却炉に投入されるゴミの量は、必ずしも一定ではなく、変動し、蒸気発生量の外乱となると考えられる。しかしながら、実際に稼働中のゴミ焼却炉（生産設備やプラント等の制御対象）では、ゴミ焼却炉のモデル化のために実績入出力データを取得する実験機会や、ゴミ焼却炉に投入可能なゴミ量等が限定され、このようなノイズを含むデータしか取得できない場合もある。
【０１１７】
このようなゴミ焼却炉において、先験的情報として、時刻０において給塵機の速度を１単位だけステップ状に変化させる単位ステップ入力で変化させた場合に、蒸気発生量（ステップ応答）は、図６に示すように、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線、ステップ応答上限）と、時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線、ステップ応答下限）との間に、存在することが分かっているものとする。このような先験的知識は、通常、ゴミ焼却炉を実際に稼動させるオペレータの経験や、ゴミ焼却炉の過去実績データ等が理解されるものであるが、ゴミ焼却炉の設計値等からも理解される場合もある。例えば、ゴミのカロリー（熱量）における変動範囲と、給塵機を１単位だけ変化させた場合におけるゴミの供給量変化に基づき、蒸気発生量の変化は、概略計算が可能である。
【０１１８】
なお、このような推定計算は、このゴミ焼却炉のようなプラントだけでなく、種々の装置で可能である。例えば、リンク機構として、マニピュレータでは、入力トルクを１単位だけ変化させた場合に、アームの回転速度がどのくらいの時間でどのくらいの範囲で変化するか、その慣性や粘性等のパラメータから、時間応答の上限および下限が概略計算することができる。
【０１１９】
図５に示す実績データおよび図６に示す制約条件の下に、ゴミ燃焼炉を本実施形態のモデリング装置Ｓによってモデル化することによって作成されたゴミ燃焼炉のモデルにおけるステップ応答を図７に示す。この本実施形態のモデリング装置Ｓによるゴミ燃焼炉のモデルでは、そのステップ応答（蒸気発生量）は、図７（Ａ）に示すように、図５（Ｂ）に示す実際の蒸気発生量の変化によく一致（フィッティング）しており、そして、図７（Ｂ）に示すように、ステップ応答の上限および下限の制約条件を満たしている。
【０１２０】
一方、図５に示す実績データおよび図６に示す制約条件の下に、ゴミ燃焼炉を従来の手法によってモデル化することによって作成されたゴミ燃焼炉のモデルにおけるステップ応答を図８に示す。この従来の手法によるゴミ燃焼炉のモデルでは、図７と比較すると、図７に示すステップ応答（蒸気発生量）ほど、実際の蒸気発生量の変化に一致していない。特に、図８（Ａ）に破線の○で囲んだ外乱の大きな部分の実績出力データにモデルが合わせ込もうとされており、そのため、図８に示すステップ応答では、図７に示すステップ応答に較べて緩やかな立ち上がりとなってしまってしい、図７に示すステップ応答に較べてフィッティングの精度が落ちている。さらに、図８（Ｂ）に示すように、ステップ応答の上限および下限の制約条件を満たしていない部分がある。この結果、従来の手法では、再度、モデル化に使用するデータ区間を変更したり、フィルタリング等のデータの前処理を行ったり、あるいは、実績入出力データを取り直ししたり等、先験的情報の制約条件を満たすモデルが作成されるまで、試行錯誤的にモデル化作業を繰り返す必要がある。
【０１２１】
しかしながら、本実施形態のモデリング装置Ｓは、図７と図８とを比較すると分かるように、このようなモデル化作業が軽減される。
【０１２２】
ここで、上述のゴミ燃焼炉のモデル化では、図６（Ｂ）に示すように、ステップ応答が１２００秒で略整定（飽和）していることから、制約条件の時間範囲をこの１２００秒で区切っている。このため、１２００秒を越える部分で、モデルの実用上問題の無いレベルであるが、ステップ応答がステップ応答の上限を越えてしまっている。そこで、厳密に、ステップ応答がステップ応答の上限および下限を逸脱しないようにするために、定常ゲインＧの上下限値が制約条件とされてもよい。すなわち、上述したように、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下であることが、先験的知識として判明している場合に、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定され、そして、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、更新される。このような制約条件を課すことにより、ステップ応答（蒸気発生量）は、図９（Ａ）に示すように、図５（Ｂ）に示す実際の蒸気発生量の変化によく一致（フィッティング）しているだけでなく、図９（Ｂ）に示すように、１２００秒を越えた範囲も含めて、ステップ応答の上限および下限の制約条件を満たしている。
【０１２３】
なお、上述では、給塵機の速度を入力とし、蒸気発生量を出力としたが、入出力の関係は、これに限定されるものではなく、種々の関係であって良い。例えば、入力としては、押込空気量や、空気供給量や、バーナ油量や、蒸気弁開度や、触媒投入量等が挙げられ、出力としては、砂層温度や、炉内温度や、ボイラドラム圧力や、炉内圧力や、酸素濃度や、一酸化炭素濃度や、ＮＯｘ濃度等が挙げられる。
【０１２４】
本実施形態のモデリング装置Ｓの適用例としてゴミ焼却炉のモデル化について説明したが、本実施形態のモデリング装置Ｓは、操作量を入力とするとともに制御量を出力とすることで、様々な制御対象のモデル化に適用可能である。本実施形態のモデリング装置Ｓは、例えば、圧延機のモデル化に適用可能であり、例えば、圧下位置やロール速度等を入力とし、スタンド間張力や板幅・板厚等を出力とすることができる。また例えば、本実施形態のモデリング装置Ｓは、マニピュレータのモデル化に適用可能であり、トルク等を入力とし、アーム速度や位置等を出力とすることができる。
【０１２５】
なお、上述の実施形態において、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うように、モデリング装置Ｓが構成されてもよい。一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。このため、このように構成されることによって、この初期値自体も未知パラメータとして取り扱われるので、初期値に左右されることなく、システムのモデルをより高精度に構築することが可能となる。
【０１２６】
例えば、上述の１入力１出力系のシステムの場合において、初期状態が０（ゼロベクトル）であることを前提にモデルが構築されている。すなわち、モデリング対象のシステムが整定している状態（平衡状態）において入力データが印加されたものと見なされており、モデルの構築に用いた例えば実績入出力データ等の入出力データは、初期状態が０であるとされている。
【０１２７】
式４において、状態ｘ（ｔ）の挙動は、一般に、式１１となることが知られている。
ｘ（ｔ）＝Φ（ｔ）ｘ（０）＋∫Φ（ｔ−τ）Ｂｕ（τ）ｄτ ・・・（１１）
ここで、積分区間は、０からｔまでであり、Φ（ｔ）は、遷移マトリックスと呼ばれ、式１２によって表される。
Φ（ｔ）＝ｅｘｐ（Ａｔ）
＝Ｉ＋Ａｔ＋Ａ２ｔ２／２！＋Ａ３ｔ３／３！＋・・・ ・・・（１２）
【０１２８】
この式１１で表される状態ｘ（ｔ）を式４−２に代入すると、式１３が得られ、初期状態ｘ（０）＝０でなければ、出力ｙｍ（ｔ）は、初期状態ｘ（０）の影響を受けることが分かる。
ｙｍ（ｔ）＝ＣΦ（ｔ）ｘ（０）＋Ｃ∫Φ（ｔ−τ）Ｂｕ（τ）ｄτ＋Ｄｕ（ｔ） ・・・（１３）
【０１２９】
モデルを構築するための入出力データを得る実際の実験において、例えばマニピュレータ等の機械系のシステムでは、システムが静止した状態から実験を開始することが可能であり、この場合では初期状態を０と見なせるが、例えばゴミ燃焼炉や化学プラント等のプラント系のシステムでは、稼動中に実験を行わなければならず、状態が絶えず変動しており、通常、初期状態が０とはならない。
【０１３０】
このため、上述のように、初期状態自体もパラメータとして扱い、モデルを構築することが好ましい。
【０１３１】
より具体的には、例えば、上述の１入力１出力系のシステムの場合において、式４および式５、あるいは、式６では、モデリング対象のシステムのモデル次数が２次であるから、初期状態ｘ（０）＝［ｘ１ ｘ２］Ｔとおくことができる。このため、Ｐ個のうちのｉ番目の粒子ｐ（ｉ，ｋ）は、式７に代えて、式１４によって表現されることになる。
ｐ（ｉ，ｋ）＝［ａ０ ａ１ ｂ０ ｂ１ ｘ１ ｘ２］Ｔ ・・・（１４）
【０１３２】
また、この場合におけるｉ番目の粒子ｐの移動速度ｖ（ｉ，０）の初期値は、ｖ（ｉ，０）＝［０ ０ ０ ０ ０ ０］Ｔに設定される。
【０１３３】
このような初期値の設定を行って、以下、上述と同様の処理によってシステムのモデルが構築される。
【０１３４】
また、上述の実施形態において、予め決められた所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値をＰＳＯ法の初期値に加えて用いるように、モデリング装置Ｓが構成されてもよい。前記所定の最適化方法として、ＰＳＯ法、最小二乗法および部分空間法等の各種の最適化手法が挙げられる。なお、最小二乗法や部分空間法が用いられる場合には、離散時間で一旦状態空間モデル式や伝達関数モデル式を求め、これを連続時間のモデルに変換される。また、ＰＳＯ法を用いる場合では、評価関数のみを単に最小化するモデルが求められても良いし、その際に、Ｊｉ＝｜｜Ｗ（Ｙｉ−Ｙｄａｔａ）｜｜＜ε（εは所定の閾値である）となる粒子ｐを保存しておいても良い。このＰＳＯ法を用いる場合では、初期値の決定およびモデルの構築の２回、ＰＳＯ法が使われることになる。一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。このため、ＰＳＯ法による場合でもランダムに振る範囲が試行錯誤によって決定される。このため、この試行錯誤の手間を軽減し、最適なモデルを求めるためには、最初のモデルが入出力データにある程度フィッティングしていることが望ましい。この構成では、ＰＳＯに用いる初期値として、予め決められた所定の最適化方法によって求められたモデルのパラメータの値も用いるので、より真に近い値を含む初期値からＰＳＯ法の演算処理を開始することができるから、不必要な探索を行うことなく、システムのモデルをより高精度に構築することができる。
【０１３５】
本発明を表現するために、上述において図面を参照しながら実施形態を通して本発明を適切且つ十分に説明したが、当業者であれば上述の実施形態を変更および／または改良することは容易に為し得ることであると認識すべきである。したがって、当業者が実施する変更形態または改良形態が、請求の範囲に記載された請求項の権利範囲を離脱するレベルのものでない限り、当該変更形態または当該改良形態は、当該請求項の権利範囲に包括されると解釈される。
【符号の説明】
【０１３６】
Ｓ モデリング装置
１ 演算部
２ 入力設定部
４ 記憶部
１１ モデル出力演算部
１２ 制約条件評価部
１３ 評価演算部
１４ モデル更新決定部
２１ 制約条件入力設定部
２２ モデル次数同定パラメータ入力設定部
２３ 最適化計算パラメータ入力設定部
２４ 評価関数入力設定部
４１ 制約条件記憶部
４２ モデル次数同定パラメータ記憶部
４３ 最適化計算パラメータ記憶部
４４ 評価関数記憶部
４５ モデル出力等記憶部
４６ 実績データ記憶部
【技術分野】
【０００１】
本発明は、システムのモデルを作成するモデリング装置およびモデリング方法に関し、特に、先験的知識による先験的情報を利用することができるモデリング装置およびモデリング方法に関する。
【背景技術】
【０００２】
近年の計算機（コンピュータ）の発展に伴い、システムの制御系の設計には、コンピュータを用いたシミュレーションが活用されている。実世界のシステムをコンピュータ上の仮想世界に取り入れるためには、実世界のシステムをコンピュータで取り扱い可能な形式で表現した、実世界と仮想世界とを結び付けるモデルで実世界のシステムを記述する必要がある。このモデルは、通常、微分方程式、差分方程式、伝達関数、状態方程式および論理式等の数学モデルである。シミュレーション結果が実世界のシステムを正確に表すものであるか否かは、このモデルの良し悪し（モデルの精度）によって決定され、このモデルを構築する作業は、モデリングと呼ばれる。
【０００３】
モデリング方法は、ホワイトボックスモデリングとブラックボックスモデリングとに大別され、さらに、これらホワイトボックスモデリングとブラックボックスモデリングとを折衷したグレーボックスモデリングがある。
【０００４】
このホワイトボックスモデリングは、いわゆる第１原理に基づくモデリングであり、モデリングの対象であるシステムを支配する物理法則や化学法則等の自然法則（第１原理）と物理化学的な諸量とに基づいてモデリングを行う方法であり、対象が物理システムの場合には物理モデリングと呼ばれる。この方法は、システム（対象）が比較的小規模でこのシステムに作用している物理法則等の自然法則がよく分かっている場合に有効であるが、一般にシステムの大規模化あるいは複雑化に従ってその適用が難しくなる。なお、大規模システムの場合に、この大規模システムを小規模なサブシステムに分割し、各サブシステムに本方法を適用してモデリングを行い、それらを統合することによって大規模システムに本方法を適用することもある。
【０００５】
ブラックボックスモデリングは、物理法則や化学法則等の自然法則や先験的知識を用いることなく、システムをブラックボックスと見なして、システム（対象）で観測された有限個の入出力データから、所定の目的（仕様）に応じたモデルを構築する方法であり、制御工学の分野では、システム同定と呼ばれる。システムが比較的大規模あるいは複雑であり、物理法則等の自然法則のみによってシステムの動特性が記述することができない場合に有効であるが、同じ入出力データが与えられても、モデリングの目的（仕様）やモデリングを行う技術者によって得られるモデルが異なる場合があり、主観的要素が入り込む余地がある。
【０００６】
グレーボックスモデリングは、これらホワイトボックスモデリングとブラックボックスモデリングとの中間に位置づけられる方法であり、例えば、モデルの構築に、観測された入出力データに加えて、物理法則等の自然法則が部分的に利用される。
【０００７】
モデリングの対象となるシステムには、例えば、多関節ロボット（マニピュレータ）等のようにシステムの基本的な挙動が運動方程式によって記述することができ、ホワイトボックスモデリングを採用可能なシステムもあるが、一般的には、例えば熱処理や化学処理等の伴う燃焼炉や化学プラント等のように、その動特性を自然法則から導出し難いケースが多く、システムのモデリング（システムのモデル化）には、ブラックボックスモデルイングやグレーボックスモデルイングを用いることが多い。
【０００８】
このブラックボックスモデリングやグレーボックするモデリングでは、入出力データを用いることによってモデルが構築されるため、モデルの良し悪しがその入出力データに依存し、高精度なモデリングを可能とする質および量ともに充分な入出力データが事前に得られることが望ましい。しかしながら、現実には、実際のシステムに入力可能な波形や大きさ（レベル）等が限られていたり、これら入出力データを採取する実験機会が限られていたり等するために、必ずしも前記質および量ともに充分な入出力データが得られない場合が多い。さらに、実際のシステムから採取された入出力データには、ノイズも存在し、このノイズが場合によっては、モデリングに無視できないほどの大きな悪影響を及ぼすこともあり得る。また、システムによっては、例えば、入力信号に対する持続励振性（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)
【発明の概要】
【発明が解決しようとする課題】
【００１３】
本発明は、上述の事情に鑑みて為された発明であり、その目的は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができるモデリング装置およびモデリング方法を提供することである。
【課題を解決するための手段】
【００１４】
本発明者は、種々検討した結果、上記目的は、以下の本発明により達成されることを見出した。すなわち、本発明の一態様にかかるモデリング装置は、システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング装置であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備えることを特徴とする。
【００１５】
そして、本発明の他の一態様にかかるモデリング方法は、システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング方法であって、前記モデルが、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定ステップと、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算ステップと、前記制約条件設定ステップで設定された制約条件を満たし、かつ、前記評価演算ステップで演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算ステップとを備えることを特徴とする。
【００１６】
このような構成のモデリング装置およびモデリング方法では、実績出力データとモデル出力データとの誤差が評価関数で評価され、この評価結果が最良となるようにモデルのパラメータの値が求められ、モデルが構築される。そして、先験的知識を表す先験的情報が制約条件として設定され、このモデルを構築する際に、制約条件も満たすように、モデルのパラメータの値が求められる。このため、このような構成のモデリング装置およびモデリング方法は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【００１７】
また、他の一態様では、上述のモデリング装置において、前記先験的情報は、前記システムにおける、定常ゲインに関する上限値および下限値のうちの少なくとも一方であることを特徴とする。
【００１８】
この構成によれば、定常ゲインに関する先験的情報も制約条件として加味され、定常ゲインに関する先験的情報を満たすシステムのモデルを構築することができる。
【００１９】
また、他の一態様では、これら上述のモデリング装置において、前記評価関数は、前記実績出力データと前記モデル出力データとの誤差に重み付けを行う重みを有し、前記重みは、その値が大きいほど前記評価結果に大きく影響を与えるものである場合に、前記重みの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されることを特徴とする。
【００２０】
この構成によれば、所定の閾値以上のノイズを含む実績出力データの、評価結果に与える影響が低減され、システムのモデルをより高精度に構築することができる。
【００２１】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、前記モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うことを特徴とする。
【００２２】
一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。この構成によれば、この初期値自体も未知パラメータとして取り扱われるので、初期値に左右されることなく、システムのモデルをより高精度に構築することができる。
【００２３】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄ（ＰＳＯ）によって前記モデルのパラメータの値を求めることを特徴とする。
【００２４】
この構成によれば、モデルのパラメータの値を求める最適化方法としてＰＳＯを用いたモデリング装置を提供することができる。
【００２５】
また、他の一態様では、上述のモデリング装置において、前記モデル演算部は、所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値を前記Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄの初期値に加えることを特徴とする。
【００２６】
この構成によれば、ＰＳＯに用いる初期値として、予め決められた所定の最適化方法によって求められたモデルのパラメータの値も用いるので、より真に近い値を加えた初期値からＰＳＯ法の演算処理を開始することができるから、不必要な探索を行うことなく、システムのモデルをより高精度に構築することができる。
【００２７】
また、他の一態様では、これら上述のモデリング装置において、前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによる更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いることを特徴とする。
【００２８】
この構成によれば、これら定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いてＰＳＯの更新を行うので、不必要な探索を行うことなく、これらのうちの少なくとも１つを制約条件として満たすシステムのモデルを構築することができる。
【発明の効果】
【００２９】
本発明にかかるモデリング装置およびモデリング方法は、前記非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【図面の簡単な説明】
【００３０】
【図１】実施形態におけるモデリング装置の構成を示すブロック図である。
【図２】実施形態におけるモデリング装置の概略動作を示すフローチャートである。
【図３】実施形態のモデリング装置における制約条件の設定動作を示すフローチャートである。
【図４】実施形態のモデリング装置におけるモデルの演算動作を示すフローチャートである。
【図５】ゴミ焼却炉において、給塵機の速度をステップ状に変化させた場合における蒸気発生量の変化を示す図である。
【図６】制約条件としてのステップ応答の上限および下限を示す図である。
【図７】実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【図８】比較例としての従来の手法によりモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【図９】定常ゲインの制約条件を考慮した場合であって、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。
【発明を実施するための形態】
【００３１】
以下、本発明にかかる実施の一形態を図面に基づいて説明する。なお、各図において同一の符号を付した構成は、同一の構成であることを示し、適宜、その説明を省略する。
【００３２】
図１は、実施形態におけるモデリング装置の構成を示すブロック図である。実施形態のモデリング装置Ｓは、モデリングすべきシステムのモデルを前記モデルのパラメータの値を決定することによって構築する装置であって、前記モデルは、連続した時間と見なされる時間変数を含む数式であり、予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備える装置である。
【００３３】
このようなモデリング装置Ｓは、例えば、図１に示すように、演算部１と、入力設定部２と、出力部３と、記憶部４と、これら演算部１、入力設定部２、出力部３および記憶部４の間で相互にデータを交換することができるようにこれらを互いに接続するバス５とを備えて構成される。
【００３４】
入力設定部２は、モデリング装置Ｓの演算開始指示等の各種コマンドや、例えば実績入出力データ等のシステムのモデル化に必要な各種データをモデリング装置Ｓに入力する機器であり、例えば、キーボードやマウス等である。そして、本実施形態では、入力設定部２は、制約条件を入力設定する制約条件入力設定部２１と、モデル次数および同定パラメータを入力設定するモデル次数同定パラメータ入力設定部２２と、最適化計算パラメータを入力設定する最適化計算パラメータ入力設定部２３と、評価関数を入力設定する評価関数入力設定部２４とを備えている。これら制約条件、モデル次数、同定パラメータおよび最適化計算パラメータは、システムのモデル化の際に利用されるように、モデリング装置Ｓの記憶部４に格納されることで、モデリング装置Ｓに設定される。
【００３５】
制約条件入力設定部２１は、モデルの制約条件（先験的情報）を入力してモデリング装置Ｓに設定するためのものである。制約条件は、システムのモデル化の際に、このモデル化によって得られたモデルが満たすべき条件であり、本実施形態では、モデリング対象のシステムに関する先験的知識によって与えられる。先験的知識は、任意の手段によってシステムから得られる情報であり、例えば、システムの設計値や理論式等から計算によって推定される情報であって良く、また例えば、システムを設計した場合に過去の事例に照らして経験的に推測される情報であって良く、また例えば、システムを実際に組み、この実際のシステムを稼動させてこれを観察することによって得られる情報であって良い。この得られた先験的知識は、モデリング装置Ｓが取り扱うことができる形式で表現（表記、記述）され、先験的情報とされる。
【００３６】
本実施形態では、予め決められた規範（基準）となるデータであって時間に依存する規範入力を入力データとしてシステムに与えた場合に、システムから出力される出力データ（計算によって得られる出力データ、および、システムで観測される実績出力データを含む）が制約条件（時間応答制約条件）とされる。規範入力は、例えば、ステップ状の入力（ステップ入力）、ランプ状の入力（ランプ入力）、三角関数波（例えばｓｉｎ波）や三角波等の周期入力およびそれらの周波数掃引波（スイープ波、ｓｗｅｅｐ波）等の、時間に関する所定の関数である。時間応答制約条件は、例えば前記出力データのそのものであって良く、また例えば前記出力データの上限値および／またはその下限値であって良く、さらにこれらに加えてその時間範囲が設定されて良い。すなわち、規範入力の入力時点を時刻０として、どの時刻ｔＳからどの時刻ｔＥまでの前記出力データが制約条件とされるのか、あるいは、どの時刻ｔＳからどの時刻ｔＥまでの前記出力データにおける上限値および／または下限値が制約条件とされるのか、が設定される（０≦ｔＳ＜ｔＥ）。この時間応答制約条件が入力設定される場合、制約条件入力設定部２１は、時間応答制約条件入力設定部として機能する。
【００３７】
なお、Ａおよび／またはＢは、ＡおよびＢのうちの少なくとも一方の意味であり、ＡあるいはＢ、または、ＡとＢとの両方の意味である。
【００３８】
そして、本実施形態では、システムの定常ゲインの符号や、その範囲（定常ゲインの上限値および下限値）が制約条件（定常ゲイン制約条件）とされる。この定常ゲイン制約条件が入力設定される場合、制約条件入力設定部２１は、定常ゲイン制約条件入力設定部として機能する。
【００３９】
モデル次数未知パラメータ入力設定部２２は、モデルの次数および前記モデルの同定パラメータを入力してモデリング装置Ｓに設定するためのものである。モデルの同定パラメータは、前記モデルにおけるパラメータのうち未知のパラメータとして同定すべきパラメータである。モデルの同定パラメータは、後述のシステムのモデルを作成する演算によってその値が決定され、これによってシステムのモデルが構築される。これらモデルの次数および同定パラメータは、例えば制御系設計に使用するためやシミュレーションに使用するため等の、モデリングの目的によって適宜に選定される。また、モデルの精度、モデル化に使用される入出力データの質や量、および、モデル化の情報処理量（計算量）等も勘案されて良い。
【００４０】
最適化計算パラメータ入力設定部２３は、最適化計算パラメータを入力してモデリング装置Ｓに設定するためのものである。最適化計算パラメータは、最適化の計算を実行する際に必要とされる調整パラメータである。最適化計算パラメータは、本実施形態では、最適化計算手法として、後述するようにＰＳＯ（Particle Swarm Optimization）法が用いられていることから、このＰＳＯ法の演算に必要なパラメータであり、粒子ｐの個数Ｐ、解の更新回数（繰り返し回数）Ｋ、初期値の値域、後述の式におけるｗ、ｃ１、ｃ２である。最適化計算パラメータの入力設定の際に、モデルの次数や同定パラメータの個数が勘案されて良い。例えば、モデルの次数が大きくなるに従って、あるいは、同定パラメータの個数が多くなるに従って、探索すべき解空間が広くなるので、粒子ｐの個数Ｐや解の更新回数Ｋが大きく設定されることが好ましい。
【００４１】
評価関数入力設定部２４は、評価関数を入力してモデリング装置Ｓに設定するためのものである。評価関数は、同定パラメータのパラメータ値を決めた場合に、このパラメータ値によって規定されるモデルがシステムをどの程度よく表現しているか否かの度合い（モデルの精度）を示す値（評価値）を求める関数である。評価関数は、予め１つだけ設定され、プログラムに予め組み込まれていても良いが、本実施形態では、より適切にモデルを評価するために、評価関数入力設定部２４から入力可能とされている。さらに、オペレータの便宜を図るために、本実施形態では、予め複数の評価関数が用意されており、これら複数の評価関数の中から評価関数入力設定部２４を介して選定することで、評価関数が入力設定される。
【００４２】
出力部３は、入力部２から入力されたコマンドやデータおよびモデリング装置Ｓの演算結果（システムのモデル）等を出力する機器であり、例えばＣＲＴディスプレイ、ＬＣＤ、有機ＥＬディスプレイ又はプラズマディスプレイ等の表示装置やプリンタ等の印字装置等である。
【００４３】
記憶部４は、データやプログラムを記憶する例えばＲＯＭ（Read Only Memory）およびＥＥＰＲＯＭ（Electrically Erasable Programmable Read Only Memory）等の不揮発性の記憶素子やハードディスク等を備えるとともに、演算部１が実行するモデリングプログラムや制御プログラム等の実行中における各データを一時的に記憶するいわゆるワーキングメモリであり、例えば揮発性の記憶素子であるＲＡＭ（Random Access Memory）を備えて構成される。記憶部４は、機能的に、制約条件記憶部４１と、モデル次数同定パラメータ記憶部４２と、最適化計算パラメータ記憶部４３と、評価関数記憶部４４と、モデル出力等記憶部４５と、実績データ記憶部４６とを備える。
【００４４】
制約条件記憶部４１は、制約条件を記憶するものであり、制約条件入力設定部２１から入力された制約条件が制約条件記憶部４１に格納され、記憶される。制約条件記憶部４１は、本実施形態では、時間応答制約条件を記憶する時間応答制約条件記憶部４１１と、定常ゲイン制約条件を記憶する定常ゲイン制約条件記憶部４１２とを備えている。
【００４５】
モデル次数同定パラメータ記憶部４２は、モデルの次数および前記モデルの同定パラメータを記憶するものであり、モデル次数同定パラメータ入力設定部２２から入力されたモデルの次数および同定パラメータがモデル次数同定パラメータ記憶部４２に格納され、記憶される。
【００４６】
最適化計算パラメータ記憶部４３は、最適化計算パラメータを記憶するものであり、最適化計算パラメータ入力設定部２３から入力された最適化計算パラメータが最適化計算パラメータ記憶部４３に格納され、記憶される。
【００４７】
評価関数記憶部４４は、評価関数を記憶するものであり、評価関数入力設定部２４から入力された評価関数が評価関数記憶部４４に格納され、記憶される。本実施形態では、上述したように、複数の所定の評価関数が予め記憶されており、これら複数の評価関数の中から１つの評価関数が評価関数入力設定部２４から選定され、この選定されたことを示す情報も評価関数記憶部４４に記憶される。
【００４８】
モデル出力等記憶部４５は、モデルを規定する所定のパラメータ値に対し、演算部１で演算されたモデル出力、評価結果および制約条件を満たすか否かの情報を互いに対応付けて記憶するものである。前記所定のパラメータ値は、本実施形態では、ＰＳＯ法によって与えられる値である。そして、本実施形態では、制約条件を満たすか否かの情報は、評価結果を示す評価値に変えられ、評価結果としてモデル出力等記憶部４５に記憶される。このように制約条件を満たすか否かの情報が評価結果を示す評価値に変換されて評価結果としてモデル出力等記憶部４５に記憶されることで、モデルのパラメータの値を決定する際に、評価結果を参酌するだけで済み、このモデルのパラメータの値を決定する処理が簡便化される。
【００４９】
実績データ記憶部４６は、実績データを記憶するものである。実績データは、モデリング対象のシステムに入力（印加）される所定の入力データ（実績入力データ）、および、この入力データの入力によってシステムから出力される出力データ（実績出力データ）であり、これら実績入力データおよび実績出力データ（実績入出力データ）は、互いに対応付けて実績データ記憶部４３に記憶される。前記所定の入力データは、任意のデータであって良く、例えば、ステップ状の入力データ、ランダム入力データおよび擬似白色系列（例えばＭ系列等）データ等の時系列データである。出力データもこれに対応して時系列データとなる。前記所定の入力データは、制御工学の分野では、同定入力と呼ばれる。時系列入力データおよび時系列出力データは、１組であっても良く、また複数の組であっても良い。複数の組である場合には、時系列入出力データの各組は、同日に得られたデータであっても良く、また、異なる日に得られたデータであって良い。そして、これら複数の組の時系列入出力データは、モデリングの際に、これらの中から適宜に選択されるようにモデリング装置Ｓが構成されてもよい。
【００５０】
演算部１は、例えば、マイクロプロセッサおよびその周辺回路等を備えて構成され、プログラムに従い入力設定部２、出力部３および記憶部４を当該機能に応じてそれぞれ制御することによってモデリング装置Ｓ全体の制御を司るものであって、モデリング対象のシステムのモデルを前記モデルのパラメータの値を決定することによって構築するものである。演算部１は、機能的に、モデル出力演算部１１と、制約条件評価部１２と、評価演算部１３と、モデル更新決定部１４とを備える。
【００５１】
モデル出力演算部１１は、モデル次数、同定パラメータおよび最適化計算パラメータを設定することによって設定された前記システムのモデルに、入力データを入力（印加）することで前記モデルの出力データ（モデル出力データ）を求めるものである。この求めたモデル出力データは、モデル出力等記憶部４５に記憶され、格納される。
【００５２】
制約条件評価部１２は、モデル出力演算部１１によって得られたモデル出力データが制約条件を満たしているか否かを評価するものである。この評価したモデル出力データが制約条件を満たしているか否かの情報は、モデル出力等記憶部４５に記憶され、格納される。
【００５３】
評価演算部１３は、予め決められた所定の入力データを実績入力データとしてシステムに入力した場合にこのシステムから出力された実績出力データと、前記実績入力データを前記システムのモデルに入力した場合にこのモデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求めるものである。この求めた評価結果は、モデル出力等記憶部４５に記憶され、格納される。
【００５４】
モデル更新決定部１４は、前記システムのモデルを予め決められた所定の更新則に従って更新するとともに、制約条件入力設定部２１で設定された制約条件を満たし、かつ、評価演算部１３で演算された評価結果が最良となるように、前記システムのモデルにおけるそのパラメータの値を求め、前記システムのモデルを決定するものである。
【００５５】
このようなモデリング装置Ｓは、モデリング方法のプログラム（ソフトウェア）を実装することで、例えば、デスクトップ型パソコンやノート型パソコン等のコンピュータによって構成可能である。
【００５６】
なお、必要に応じてモデリング装置Ｓは、図１に破線で示す外部記憶部６や通信インターフェース部７をさらに備えてもよい。
【００５７】
外部記憶部６は、例えば、フレキシブルディスク、ＣＤ−ＲＯＭ（Compact Disc Read Only Memory）、ＣＤ−Ｒ（Compact Disc Recordable）およびＤＶＤ−Ｒ（Digital Versatile Disc Recordable）等の記録媒体との間でデータを読み込みおよび／または書き込みを行う装置であり、例えば、フレキシブルディスクドライブ、ＣＤ−ＲＯＭドライブ、ＣＤ−Ｒドライブ及びＤＶＤ−Ｒドライブ等である。通信インターフェース部７は、通信網（ネットワーク）に接続され、この通信網を介して他の装置との間で通信信号を送受信するための機器である。
【００５８】
各プログラムが格納されていない場合には、これらを記録した記録媒体から外部記憶部６を介して記憶部４にインストールされるように構成してもよく、また、これらプログラムを管理するサーバコンピュータ（不図示）から通信網および通信インターフェース部７を介して各プログラムがダウンロードされるように構成してもよい。また、システムのモデル化の演算に当たってモデリング装置Ｓに入力すべき例えば実績入出力データ等のデータは、このデータを記憶した記録媒体によって外部記憶部６を介してモデリング装置Ｓに入力されるように構成してもよく、また、このデータを管理するサーバコンピュータから通信網および通信インターフェース部７を介してモデリング装置Ｓに入力されるように構成してもよい。
【００５９】
次に、本実施形態におけるモデリング装置の動作について説明する。図２は、実施形態におけるモデリング装置の概略動作を示すフローチャートである。図３は、実施形態のモデリング装置における制約条件の設定動作を示すフローチャートである。図４は、実施形態のモデリング装置におけるモデルの演算動作を示すフローチャートである。
【００６０】
本実施形態のモデリング装置Ｓでは、概略、図２に示すように、まず、ステップＳ１１では、制約条件が設定され、続いて、ステップＳ１２では、モデルを構築すべく、制約条件を満たし、かつ、評価結果が最良となるように、モデルのパラメータの値を決定することによって、モデルが演算される。続いて、ステップＳ１３では、この演算されたモデルの精度が予め決められた許容範囲以内であるか否かが判定され、この判定の結果、モデルの精度が許容範囲以内ではない場合（ＮＯ）には処理がステップＳ１１に戻され、一方、モデルの精度が許容範囲内である場合（ＹＥＳ）にはステップＳ１４が実行され、この演算されたモデルがシステムのモデルとして出力部３に出力され、出力部３で出力される。
【００６１】
次に、前記ステップＳ１１についてより具体的に説明する。図３において、まず、ステップＳ２１では、入力設定部２の制約条件入力設定部２１から、時間応答波形における制約条件（前記時間応答制約条件）が入力され、記憶部４の制約条件記憶部４１における時間応答制約条件記憶部４１１に格納され、記憶される。これによって、後述するように、演算部１の制約条件評価部１２が時間応答制約条件記憶部４１１に記憶されている時間応答制約条件を用いることで、この時間応答制約条件を満たすように、システムのモデルが構築される。
【００６２】
以下に、説明の簡単化のため、１入力１出力系のシステムの場合について、より具体的に説明する。もちろん、以下の説明は、一般性を損なうことなく、多入力多出力系のシステムに拡張可能である。
【００６３】
まず、例えば燃焼炉プラント等の実際のシステムから入出力データ（入力データＵｄａｔａおよび出力データＹｄａｔａ）が実績入出力データ（実績入力データＵｄａｔａおよび実績出力データＹｄａｔａ）として取得され、実績データ記憶部４６に格納され、記憶されているものとする。システムのモデルが連続した時間と見なされる時間変数ｔを含む連続時間状態空間モデル式または連続時間伝達関数モデル式であることから、実績入出力データＵｄａｔａ、Ｙｄａｔａは、時間的に連続しているデータに見えるように、より短いサンプリング間隔で取得されることが好ましい。ここでは、表記の簡単化のために、１秒間隔（１秒周期）でサンプリングされているものとする。サンプリング開始の時刻０から時刻ｔｍａｘまでのｔｍａｘ秒間における実績入出力データＵｄａｔａ、Ｙｄａｔａは、時刻ｔ秒の入力値をｕ（ｔ）とし、時刻ｔ秒の出力値をｙ（ｔ）とする場合には、それぞれ、次の式１−１、式１−２となる。
Ｕｄａｔａ＝［ｕ（０） ｕ（１） ・・・ ｕ（ｔｍａｘ）］Ｔ ・・・（１−１）
Ｙｄａｔａ＝［ｙ（０） ｙ（１） ・・・ ｙ（ｔｍａｘ）］Ｔ ・・・（１−２）
ここで、上付き添え字Ｔは、いわゆる転置を表す記号である。
【００６４】
時刻０から時刻ｔＥまでのｔＥ秒間における式２の規範入力Ｕｒをシステムに加えた場合に得られる出力データに関し、時刻ｔＳ１から時刻ｔＥ１までの間における上限Ｙｍａｘおよび時刻ｔＳ２から時刻ｔＥ２までの間における下限Ｙｍｉｎを式３−１および式３−２とすると、これら式３−１によって記述される上限Ｙｍａｘおよび式３−２によって記述される下限Ｙｍｉｎが制約条件入力設定部２１からモデリング装置Ｓに与えられる。
Ｕｒ＝［ｕｒ（０） ｕｒ（１） ・・・ ｕｒ（ｔE）］Ｔ ・・・（２）
Ｙｍａｘ＝［ｙｍａｘ（ｔｓ１） ｙｍａｘ（ｔｓ１＋１） ・・・ ｙｍａｘ（ｔＥ１）］Ｔ ・・・（３−１）
Ｙｍｉｎ＝［ｙｍｉｎ（ｔｓ２） ｙｍｉｎ（ｔｓ２＋１） ・・・ ｙｍｉｎ（ｔＥ２）］Ｔ ・・・（３−２）
ここで、ｕｒ（ｔ）は、時刻ｔ秒における入力値であり、ｙｍａｘ（ｔ）は、時刻ｔにおける出力値の上限値であり、そして、ｙｍｉｎ（ｔ）は、時刻ｔにおける出力値の下限値である。また、０≦ｔＳ１＜ｔＥ１≦ｔＥであり、０≦ｔＳ２＜ｔＥ２≦ｔＥである。
【００６５】
なお、規範入力の一具体例として、単位ステップ入力は、Ｕｒ＝［０ １ １ ・・・ １］Ｔであり、単位インパルス入力は、Ｕｒ＝［０ １ ０ ・・・ ０］Ｔであり、またランプ入力は、Ｕｒ＝［０ １ ２ ３ ・・・］Ｔである。また、規範入力は、１種類であって良く、また複数の種類であっても良い。
【００６６】
続いて、ステップＳ２２では、入力設定部２の制約条件入力設定部２１から、定常ゲインにおける制約条件（前記定常ゲイン制約条件）が入力され、記憶部４の制約条件記憶部４１における定常ゲイン制約条件記憶部４１２に格納され、記憶される。例えば、システムの定常ゲインＧがＧｍｉｎ≦Ｇ≦Ｇｍａｘとすると、これら下限値Ｇｍｉｎおよび上限値Ｇｍａｘが制約条件入力設定部２１からモデリング装置Ｓに与えられる。これによって、後述するように、演算部１の制約条件評価部１２が定常ゲイン制約条件記憶部４１２に記憶されている定常ゲイン制約条件を用いることで、この定常ゲイン制約条件を満たすように、システムのモデルが構築される。
【００６７】
次に、前記ステップＳ１２についてより具体的に説明する。図４において、まず、ステップＳ３１では、入力設定部２のモデル次数同定パラメータ入力設定部２２から、モデルの次数および前記モデルの同定パラメータが入力され、記憶部４のモデル次数同定パラメータ記憶部４２に格納され、記憶される。
【００６８】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、例えば、モデリング対象のシステムに対する連続時間状態空間モデル式が式４（式４−１、式４−２）であるとする。
ｄｘ（ｔ）／ｄｔ＝Ａｘ（ｔ）＋Ｂｕ（ｔ） ・・・（４−１）
ｙｍ（ｔ）＝Ｃｘ（ｔ）＋Ｄｕ（ｔ） ・・・（４−２）
ここで、ｘ（ｔ）は、ｎ次の状態ベクトルであり、ｙｍ（ｔ）は、モデルの出力（モデル出力、モデルによる計算値）であり、Ａ、Ｂ、ＣおよびＤは、それぞれｎ行ｎ列、ｎ行１列、１行ｎ列および１行１列の定数行列である。ここでは、説明の簡単化のため、モデリング対象のシステムが厳密にプロパーである、すなわち、Ｄ＝０とする。
【００６９】
モデルの次数は、このｎの値であり、この値ｎがモデルの次数としてモデル次数同定パラメータ入力設定部２２から入力される。この例では、ｎ＝２、すなわち、モデルの次数が２である２次のモデルを求めるものとする。また、定数行列Ａ、ＢおよびＣは、例えば可制御正準形式等の所定の形式でよいが、冗長なパラメータが無い方が好ましいことから、この例では、可観測正準形式であるとし、それぞれ式５−１、式５−２および式５−３となる。
【００７０】
【数１】
【００７１】
ここで、ａ０、ａ１、ｂ０およびｂ１は、モデルを規定する未知パラメータである。これらａ０、ａ１、ｂ０およびｂ１の未知パラメータが、モデル次数同定パラメータ入力設定部２２から入力され、同定パラメータとして指定される。すなわち、モデリング対象のシステムのモデルを構築することは、これら同定パラメータとして指定された未知パラメータａ０、ａ１、ｂ０およびｂ１の各値を求めることである。
【００７２】
なお、式４（式４−１、式４−２）および式５（式５−１〜式５−３）の求解は、式６の連続時間伝達関数モデル式Ｇ（ｓ）を求めること、すなわち、式６の連続時間伝達関数モデル式Ｇ（ｓ）における未知パラメータａ０、ａ１、ｂ０およびｂ１を求めることと等価である。
Ｇ（ｓ）＝（ｂ１ｓ＋ｂ０）／（ｓ２＋ａ１ｓ＋ａ０） ・・・（６）
【００７３】
続いて、ステップＳ３２では、入力設定部２の最適化計算パラメータ入力設定部２３から、最適化計算パラメータが入力され、記憶部４の最適化計算パラメータ記憶部４３に格納され、記憶される。
【００７４】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、未知パラメータ（同定パラメータ）ａ０、ａ１、ｂ０およびｂ１の同定は、本実施形態では、例えば、ＰＳＯ（Particle Swarm Optimization）法が用いられる。このＰＳＯ法は、解空間に複数Ｐ個の粒子ｐ（粒子群）をランダム（無作為）にばらまき、これら解空間にばらまかれた粒子群が所定の更新則に従って解空間内を移動しながら、最適解を探索する手法である。このＰＳＯ法を本実施形態に適用する場合では、解空間は、未知パラメータ（同定パラメータ）ａ０、ａ１、ｂ０およびｂ１を各座標軸とする座標空間（この例では４次元空間）であり、粒子ｐは、４次元ベクトルで表される前記座標空間上の点となり、更新回数ｋ番目におけるＰ個のうちのｉ番目の粒子ｐ（ｉ，ｋ）は、式７によって表現されることになる。
ｐ（ｉ，ｋ）＝［ａ０ ａ１ ｂ０ ｂ１］Ｔ ・・・（７）
【００７５】
以下、適宜、このｐ（ｉ，ｋ）の［ａ０ ａ１ ｂ０ ｂ１］Ｔをパラメータセットと呼称することとする。
【００７６】
このような最適化手法では、粒子ｐの個数Ｐ（０＜ｉ≦Ｐ）、更新回数Ｋ（０＜ｋ≦Ｋ）および所定の更新則におけるパラメータ（後述のｗ、ｃ１、ｃ２）が最適化計算パラメータであり、最適化計算パラメータ入力設定部２３から入力される。
【００７７】
続いて、ステップＳ３３では、入力設定部２から評価関数が入力され、設定される。
【００７８】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、評価関数は、例えば、本実施形態では、実績入力データＵｄａｔａに対応する実績出力データＹｄａｔａ（＝［ｙ（０） ｙ（１） ・・・ ｙ（ｔｍａｘ）］Ｔ）と、前記実績入力データＵｄａｔａをモデル（上述の例では、式４（式４−１、式４−２）および式５（式５−１〜式５−３）または式６）に入力した場合にこのモデルから出力されるモデル出力データＹｉ（＝［ｙｉ（０） ｙｉ（２） ・・・ ｙｉ（ｔｍａｘ）］Ｔ）との誤差を評価する関数である。このような評価関数Ｊｉは、例えば、式８−１によって与えられる。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜ ・・・（８−１）
ここで、｜｜ｖ｜｜は、いわゆる、ベクトルｖのユークリッドノルムを表し、ｖが縦ベクトルである場合では、｜｜ｖ｜｜＝√（ｖＴｖ）である。また、Ｗは、出力の予測誤差に対していずれの部分を重要視して評価するかを定める重み行列であり、例えば対角行列とする。
【００７９】
また例えば、評価関数Ｊｉは、例えば、式８−２によって与えられる無限ノルムである。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜∞ ・・・（８−２）
ここで、｜｜ｖ｜｜∞は、いわゆる、ベクトルｖの無限ノルム（最大値ノルム）を表し、ｖ＝［ｖ１ ｖ２ ・・・ ｖｎ］Ｔである場合では、｜｜ｖ｜｜∞＝ｍａｘ（｜ｖ１｜、｜ｖ２｜、・・・、｜ｖｎ｜）である。
【００８０】
また例えば、評価関数Ｊｉは、例えば、式８−３によって与えられる。
Ｊｉ＝｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜１ ・・・（８−３）
ここで、｜｜ｖ｜｜１は、いわゆる、ベクトルｖの１ノルムを表し、ｖ＝［ｖ１ ｖ２ ・・・ ｖｎ］Ｔである場合では、｜｜ｖ｜｜１＝Σ｜ｖｉ｜である。ただし、Σは、ｉ＝１からｉ＝ｎまでの和である。
【００８１】
評価関数記憶部４４には、例えば、このような複数の評価関数Ｊｉ（式８−１〜式８−３）が記憶されており、これら複数の評価関数Ｊｉの中から、評価関数入力設定部２４を介して１つの評価関数Ｊｉが選定され、入力設定される。
【００８２】
なお、このような評価関数Ｊが、実績出力データとモデル出力データとの誤差に重み付けを行う重みＷを有し、この重みＷが、その値が大きいほど評価関数Ｊによる評価結果に大きく影響を与えるものである場合に、この重みＷの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されてよい。このような前記区間およびその重みＷならびに前記他の区間およびその重みＷは、例えば、入力設定部２の最適化計算パラメータ入力設定部２３から入力設定される。このように構成することによって、所定の閾値以上のノイズを含む実績出力データの、評価結果に与える影響が低減され、システムのモデルをより高精度に構築することができる。
【００８３】
続いて、ステップＳ３４では、演算部１のモデル出力演算部１１によって、モデル出力が演算され、記憶部４のモデル出力等記憶部４５に記憶される。
【００８４】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、まず、Ｐ個の粒子ｐ（ｉ，０）（更新回数ｋ＝０、ｉ＝０〜Ｐ）、例えば１００個の粒子ｐ（ｉ，０）（ｉ＝０〜１００）の初期値が例えば乱数を生成することによってランダムに設定される。これによってＰ個の各パラメータセットｐ（ｉ，０）の値が設定される。そして、ｉ番目の粒子ｐ（ｉ，０）の移動速度ｖ（ｉ，０）の初期値が例えばｖ（ｉ，０）＝［０ ０ ０ ０］Ｔに設定される。
【００８５】
ここで、制約条件として、モデリング対象のシステムが安定な系であることが、先験的知識として判明している場合には、ａ０＞０かつａ１＞０の範囲で、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８６】
また、制約条件として、モデリング対象のシステムにおける定常ゲインが正の値であることが、先験的知識として判明している場合には、ｂ０の符号とａ０の符号が等しくなるように、すなわち、ｂ０×ａ０＞０となるように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８７】
また、制約条件として、モデリング対象のシステムに不安定なゼロ点がないことが、先験的知識として判明している場合には、ｂ１の符号とｂ０の符号が等しくなるように、すなわち、ｂ１×ｂ０＞０となるように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８８】
また、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下であることが、先験的知識として判明している場合には、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定されることが好ましい。
【００８９】
このようにｉ番目の粒子ｐ（ｉ，０）の初期値を設定する際に、制約条件を加味することによって、解が存在しない領域を探索する無駄を省くことが可能となる。
【００９０】
このようにＰ個の粒子ｐの値をそれぞれ設定することで、１個の粒子ｐに対応して１個のモデルが規定される。すなわち、ｉ番目の粒子ｐ（ｉ，０）に対応して１個のモデルが規定される。上述の例では、１００個の各粒子ｐのそれぞれに対応して１００個のモデルが規定されることになる。
【００９１】
このｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデル（この例では式４および式５によって表現されるモデル、あるいは、式６によって表現されるモデル）に対し、実績データ記憶部４６に記憶されている実績入力データＵｄａｔａが代入され、モデル出力データＹｉが得られ、この得られたモデル出力データＹｉが記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このモデル出力データＹｉの計算には、公知の制御系シミュレーションソフトウェア、例えば、Ｔｈｅ Ｍａｔｈ Ｗｏｒｋ社のＭＡＴＬＡＢ／Ｓｉｍｕｌｉｎｋ等が用いられる。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、モデル出力が演算され、記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９２】
続いて、ステップＳ３５では、制約条件評価部１２によって、制約条件が評価され、この評価結果が記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９３】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ステップＳ３４と同様に、このｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデルに対し、所定の規範入力Ｕｒが代入され、前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉ（＝［ｙｒｉ（０） ｙｒｉ（１） ・・・ ｙｒｉ（ｔｍａｘ）］Ｔ）が得られる。そして、この得られた前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが規範入力に対する制約条件を満たしているか否かが判断される。例えば、制約条件として、時刻ｔＳ１から時刻ｔＥ１までの間における上限値Ｙｍａｘおよび時刻ｔＳ２から時刻ｔＥ２までの間における下限値Ｙｍｉｎが設定されている場合には、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが時刻ｔＳ１から時刻ｔＥ１までの間で上限Ｙｍａｘ以下で、かつ、時刻ｔＳ２から時刻ｔＥ２までの間で下限Ｙｍｉｎ以上であるか否かが判断される。この判断の結果、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが時刻ｔＳ１から時刻ｔＥ１までの間で上限Ｙｍａｘ以下で、かつ、時刻ｔＳ２から時刻ｔＥ２までの間で下限Ｙｍｉｎ以上ではない場合には、この前記所定の規範入力Ｕｒに対応するモデル出力データＹｒｉが制約条件を満たさないものとされ、この制約条件の評価結果が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、制約条件が評価され、この評価結果が記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９４】
なお、上述では、粒子ｐ（ｉ，０）の初期値の設定の際に、定常ゲインＧの制約条件を考慮したが、定常ゲインＧの下限値Ｇｍｉｎおよび上限値Ｇｍａｘが制約条件として設定されている場合に、実績入力データＵｄａｔａに対応するモデル出力データＹｉが下限値Ｇｍｉｎ以上であって、かつ、上限値Ｇｍａｘ以下であるか否かが判断されることで、定常ゲインＧの制約条件が評価されても良い。
【００９５】
続いて、ステップＳ３６では、評価演算部１３によって、評価関数による評価値が演算され、記憶部４のモデル出力等記憶部４５に記憶される。
【００９６】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ｉ番目の粒子ｐ（ｉ，０）に対応して規定されるモデルに対し、前記評価関数Ｊｉによる評価値が演算され、この評価関数Ｊｉの評価値が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このようにＰ個の各パラメータセットｐ（ｉ，ｋ）のそれぞれに対し、評価関数による評価値が演算され、記憶部４のモデル出力等記憶部４５にパラメータセットｐ（ｉ，ｋ）に対応付けて記憶される。
【００９７】
ここで、上記制約条件の評価結果を記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶する際に、上述のように、この制約条件の評価結果そのものを記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶してもよいが、本実施形態では、この制約条件の評価結果が評価関数Ｊｉの評価値で記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。例えば、制約条件を満たさない場合に、評価が悪いことを示す評価値が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。より具体的には、例えば、評価関数Ｊｉが、評価が高いほどその評価値が小さくなるように構成されている場合では、制約条件を満たさない場合に、例えばモデル出力データＹｉがその値として取り得ない大きな評価値が、記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。また例えば、評価関数Ｊｉが、その取り得る値の範囲が負であるように構成されている場合では、制約条件を満たさない場合に、正の値（例えば固定値や制約条件の外れ度合いに従った値等）が記憶部４のモデル出力等記憶部４５にｉ番目の粒子ｐ（ｉ，０）に対応付けて記憶される。このような評価関数Ｊｉは、例えば、式８−４によって与えられる。
Ｊｉ＝ｔａｎ−１（｜｜Ｗ・（Ｙｉ−Ｙｄａｔａ）｜｜）−π／２ ・・・（８−４）
【００９８】
続いて、ステップＳ３７では、モデル更新決定部１４によって、最適化方法の繰り返し計算における更新が終了しているか否かが判断され、この判断の結果、更新が終了していない場合（ＮＯ）では、ステップＳ３８の処理が実行され、一方、更新が終了している場合（ＹＥＳ）では、ステップＳ３９の処理が実行される。更新回数Ｋは、任意の回数で良いが、ＰＳＯ法では、解の精度や演算処理量等を考慮することによって適宜な回数に設定され、例えば、Ｋ＝１００回とされる。
【００９９】
ステップＳ３８では、モデル更新決定部１４によって、所定の更新則に従って更新され、処理がステップＳ３４に戻される。
【０１００】
上述の１入力１出力系のシステムの場合において、より具体的に説明すると、ｋ回目（ｋ＝０，１，２，・・・、Ｋ）の更新時点において、Ｐ個の各粒子ｐの中からその評価関数Ｊの値が最も小さい粒子ｐが探索され、第１最良評価粒子ｐ（ｂｅｓｔ，ｋ）とされる。例えば、ｋ回目の更新時点において、ｑ番目の粒子ｐ（ｑ、ｋ）の評価関数Ｊｑの値が最も小さい場合には、ｐ（ｂｅｓｔ，ｋ）＝ｐ（ｑ，ｋ）である。
【０１０１】
また、ｉ番目の粒子ｐにおいて、０回目（初期値）からｋ−１回目までの中で評価関数Ｊの値が最も小さい粒子ｐが探索され、第２最良評価粒子ｐ（ｉ，ｂｅｓｔ）とされる。例えば、ｒ回目の更新時点において、ｉ番目の粒子ｐ（ｉ、ｒ）の評価関数Ｊｉの値が最も小さい場合には、ｐ（ｉ，ｂｅｓｔ）＝ｐ（ｉ，ｒ）である。なお、ｋ＝０の場合では、それぞれの初期値そのものが各粒子ｐのｐ（ｉ，ｂｅｓｔ）である。
【０１０２】
そして、例えば、式９（式９−１、式９−２）によって与えられる更新則に従って、各粒子ｐの速度および位置がそれぞれ決定され、更新される。
ｖ（ｉ，ｋ＋１）＝ｗｖ（ｉ，ｋ）
＋ｃ１Ｒ１ｉｋ（ｐ（ｉ，ｂｅｓｔ）−ｐ（ｉ，ｋ））
＋ｃ２Ｒ２ｉｋ（ｐ（ｂｅｓｔ，ｋ）−ｐ（ｉ，ｋ）） ・・・（９−１）
ｐ（ｉ，ｋ＋１）＝ｐ（ｉ，ｋ）＋ｖ（ｉ，ｋ＋１） ・・・（９−２）
ここで、ｗ、ｃ１、ｃ２は、それぞれ、調整パラメータであり、最適化計算パラメータとして上述したように最適化計算パラメータ入力設定部２３から入力され、設定される。例えば、ｗ＝０．９、ｃ１＝０．８およびｃ２＝０．８等である。また、Ｒ１ｉｋおよびＲ２ｉｋは、それぞれ対角行列であり、対角要素のそれぞれは、１〜０の範囲で一様に分布する乱数である。この乱数値は、固定値ではなく、各粒子ｐごとに、かつ、各更新ｋごとに、それぞれランダムな値が使用される。
【０１０３】
ここで、制約条件として、モデリング対象のシステムが安定な系であることが、先験的知識として判明している場合には、ａ０＞０かつａ１＞０の範囲に入るように更新することが好ましい。より具体的には、更新後の粒子ｐがａ０＜０となった場合には、式１０−１に示すように、更新前のａ０と０とを内分する点に更新される。すなわち、０〜１の乱数αにおいて、
ａ０（ｉ，ｋ＋１）＝（１−α）・ａ０（ｉ，ｋ）＋α・０ ・・・（１０−１）
ここで、ｉ番目の粒子ｐについて、ａ０（ｉ，ｋ＋１）は、更新後のａ０であり、ａ０（ｉ，ｋ）は、更新前のａ０である。また、粒子ｐの速度ｖは、ａ０に対応する部分ｖ（ｉ，ｋ＋１）｛ａ０｝について、式１０−２となる。
ｖ（ｉ，ｋ＋１）｛ａ０｝＝ａ０（ｋ＋１）−ａ０（ｋ） ・・・（１０−２）
【０１０４】
また、制約条件として、モデリング対象のシステムにおける定常ゲインが正の値であることが、先験的知識として判明している場合には、ｂ０の符号とａ０の符号が等しくなるように、すなわち、ｂ０×ａ０＞０となるように、更新することが好ましい。
【０１０５】
また、制約条件として、モデリング対象のシステムに不安定なゼロ点がないことが、先験的知識として判明している場合には、ｂ１の符号とｂ０の符号が等しくなるように、すなわち、ｂ１×ｂ０＞０となるように、更新することが好ましい。
【０１０６】
また、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下の範囲であることが、先験的知識として判明している場合には、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、更新することが好ましい。
【０１０７】
このように粒子ｐを更新する際に、制約条件を加味することによって、解が存在しない領域を探索する無駄を省くことが可能となり、よりよい解をより早く探索することができる。
【０１０８】
さらに、制約条件として、同定パラメータの符号ａ０、ａ１、ｂ０、ｂ０１が先験的知識として判明している場合には、これら同定パラメータが指数形式で表現されても良い。例えば、式６の伝達関数が安定かつ不安定零点なしかつゲインが負であることが先験的知識として判明している場合では、ａ０の符号およびａ１の符号が正であって、かつ、ｂ０の符号およびｂ１の符号が負でなければならないことから、ａ０＝ｅｘｐ（λ０）、ａ１＝ｅｘｐ（λ１）、ｂ０＝−ｅｘｐ（μ０）、ｂ１＝ｅｘｐ（μ１）とおき、式７のパラメータセットがｐ（ｉ，ｋ）＝［λ０ λ１ μ０ μ１］Ｔとなる。したがって、同定パラメータ［ａ０ ａ１ ｂ０ ｂ１］Ｔの座標空間ではなく、このパラメータ［λ０ λ１ μ０ μ１］Ｔの座標空間で解が探索される。このように構成することによって、−∞＜λ０＜＋∞、−∞＜λ１＜＋∞、−∞＜μ０＜＋∞、−∞＜μ１＜＋∞の実数全体での演算処理が可能となり、符号や範囲をチェックすることなく演算処理が可能となる。
【０１０９】
このようにＰＳＯ法による更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つが用いられてよい。このようにモデリング装置Ｓを構成することによって、不必要な探索を行うことなく、これらのうちの少なくとも１つを制約条件として満たすシステムのモデルを構築することができる。
【０１１０】
こうしてＰ個の各パラメータセットｐ（ｉ，ｋ）が所定の更新則に従って更新され、処理がステップＳ３４に戻される。
【０１１１】
一方、ステップＳ３９では、モデル更新決定部１４によって、制約条件を満たし、かつ最良の評価値を与えるように、同定パラメータセットの値が決定され、このシステムのモデルが構築される。例えば、上述のＰＳＯ法では、Ｐ×Ｋ個の各パラメータセットｐ（ｉ，ｋ）の中から、制約条件を満たし、かつ最良の評価値を与えるパラメータセットｐ（ｉ，ｋ）が選択され、この選択されたパラメータセットｐ（ｉ，ｋ）の値が、モデリング対象のシステムのモデルを規定する同定パラメータの値として決定され、このシステムのモデルが決定される。
【０１１２】
このように本実施形態におけるモデリング装置Ｓおよびモデリング方法では、実績出力データとモデル出力データとの誤差が評価関数Ｊで評価され、この評価結果が最良となるようにモデルのパラメータの値が求められ、モデルが構築される。そして、先験的知識を表す先験的情報が制約条件として設定され、このモデルを構築する際に、制約条件も満たすように、モデルのパラメータの値が求められる。このため、本実施形態におけるモデリング装置Ｓおよびモデリング方法は、上述の非特許文献１や特許文献１とは異なる手法で、先験的知識を利用することによって、前記質および量の不充分な入出力データであっても、より精度の高いモデルを構築することができる。
【０１１３】
また、本実施形態におけるモデリング装置Ｓおよびモデリング方法では、システムにおける、定常ゲインＧに関する上限値Ｇｍａｘおよび下限値Ｇｍｉｎのうちの少なくとも一方が、本実施形態では定常ゲインＧに関する上限値Ｇｍａｘおよび下限値Ｇｍｉｎが制約条件とされている。このため，本実施形態におけるモデリング装置Ｓおよびモデリング方法は、定常ゲインＧに関する先験的情報も制約条件として加味され、定常ゲインＧに関する先験的情報を満たすシステムのモデルを構築することができる。
【０１１４】
次に、一実施例として、実施形態におけるモデリング装置Ｓを用いて例えばごみ焼却炉（ガス化溶融炉）のモデルを作成した場合について、説明する。
【０１１５】
図５は、ゴミ焼却炉において、給塵機の速度をステップ状に変化させた場合における蒸気発生量の変化を示す図である。図５（Ａ）は、給塵機の速度を示し、その横軸は、秒単位で表す経過時間を示し、その縦軸は、速度を示す。図５（Ｂ）は、蒸気発生量を示し、その横軸は、秒単位で表す経過時間を示し、その縦軸は、蒸気発生量を示す。図６は、制約条件としてのステップ応答の上限および下限を示す図である。図６の横軸は、秒単位で表す経過時間を示し、その縦軸は、蒸気発生量を示す。図７は、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図７（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図７（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。図８は、比較例としての従来の手法によりモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図８（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図８（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。図９は、定常ゲインの制約条件を考慮した場合であって、実施形態のモデリング装置でモデリングされたゴミ焼却炉のモデルにおけるステップ応答を説明するための図である。図９（Ａ）は、給塵機の速度を示し、その太い実線は、モデル出力としての蒸気発生量を示し、その細い実線は、図５（Ｂ）に示す実績データとしての蒸気発生量を示す。図９（Ｂ）は、モデルのステップ応答（実線）、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線）および時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線）が示されている。
【０１１６】
ゴミ焼却炉では、ゴミを焼却炉に供給する給塵機を入力とし、焼却炉で発生する蒸気発生量を出力とすると、図５（Ａ）に示すように、給塵機の速度がステップ状に変化させた場合に、蒸気発生量は、図５（Ｂ）に示すように、基本的に、給塵機の稼動開始に伴って徐々に増大するとともに給塵機の稼動停止に伴って徐々に減少するプロファイルであるが、このプロファイルに外乱の影響によりノイズが重畳しており、その波形が大きく乱れている。この原因として、例えば、給塵機の速度が一定であっても給塵機によって実際に焼却炉に投入されるゴミの量は、必ずしも一定ではなく、変動し、蒸気発生量の外乱となると考えられる。しかしながら、実際に稼働中のゴミ焼却炉（生産設備やプラント等の制御対象）では、ゴミ焼却炉のモデル化のために実績入出力データを取得する実験機会や、ゴミ焼却炉に投入可能なゴミ量等が限定され、このようなノイズを含むデータしか取得できない場合もある。
【０１１７】
このようなゴミ焼却炉において、先験的情報として、時刻０において給塵機の速度を１単位だけステップ状に変化させる単位ステップ入力で変化させた場合に、蒸気発生量（ステップ応答）は、図６に示すように、時刻ｔＳ１から時刻ｔＥ１までの間における上限のステップ応答Ｙｍａｘ（破線、ステップ応答上限）と、時刻ｔＳ２から時刻ｔＥ２までの間における下限のステップ応答Ｙｍｉｎ（破線、ステップ応答下限）との間に、存在することが分かっているものとする。このような先験的知識は、通常、ゴミ焼却炉を実際に稼動させるオペレータの経験や、ゴミ焼却炉の過去実績データ等が理解されるものであるが、ゴミ焼却炉の設計値等からも理解される場合もある。例えば、ゴミのカロリー（熱量）における変動範囲と、給塵機を１単位だけ変化させた場合におけるゴミの供給量変化に基づき、蒸気発生量の変化は、概略計算が可能である。
【０１１８】
なお、このような推定計算は、このゴミ焼却炉のようなプラントだけでなく、種々の装置で可能である。例えば、リンク機構として、マニピュレータでは、入力トルクを１単位だけ変化させた場合に、アームの回転速度がどのくらいの時間でどのくらいの範囲で変化するか、その慣性や粘性等のパラメータから、時間応答の上限および下限が概略計算することができる。
【０１１９】
図５に示す実績データおよび図６に示す制約条件の下に、ゴミ燃焼炉を本実施形態のモデリング装置Ｓによってモデル化することによって作成されたゴミ燃焼炉のモデルにおけるステップ応答を図７に示す。この本実施形態のモデリング装置Ｓによるゴミ燃焼炉のモデルでは、そのステップ応答（蒸気発生量）は、図７（Ａ）に示すように、図５（Ｂ）に示す実際の蒸気発生量の変化によく一致（フィッティング）しており、そして、図７（Ｂ）に示すように、ステップ応答の上限および下限の制約条件を満たしている。
【０１２０】
一方、図５に示す実績データおよび図６に示す制約条件の下に、ゴミ燃焼炉を従来の手法によってモデル化することによって作成されたゴミ燃焼炉のモデルにおけるステップ応答を図８に示す。この従来の手法によるゴミ燃焼炉のモデルでは、図７と比較すると、図７に示すステップ応答（蒸気発生量）ほど、実際の蒸気発生量の変化に一致していない。特に、図８（Ａ）に破線の○で囲んだ外乱の大きな部分の実績出力データにモデルが合わせ込もうとされており、そのため、図８に示すステップ応答では、図７に示すステップ応答に較べて緩やかな立ち上がりとなってしまってしい、図７に示すステップ応答に較べてフィッティングの精度が落ちている。さらに、図８（Ｂ）に示すように、ステップ応答の上限および下限の制約条件を満たしていない部分がある。この結果、従来の手法では、再度、モデル化に使用するデータ区間を変更したり、フィルタリング等のデータの前処理を行ったり、あるいは、実績入出力データを取り直ししたり等、先験的情報の制約条件を満たすモデルが作成されるまで、試行錯誤的にモデル化作業を繰り返す必要がある。
【０１２１】
しかしながら、本実施形態のモデリング装置Ｓは、図７と図８とを比較すると分かるように、このようなモデル化作業が軽減される。
【０１２２】
ここで、上述のゴミ燃焼炉のモデル化では、図６（Ｂ）に示すように、ステップ応答が１２００秒で略整定（飽和）していることから、制約条件の時間範囲をこの１２００秒で区切っている。このため、１２００秒を越える部分で、モデルの実用上問題の無いレベルであるが、ステップ応答がステップ応答の上限を越えてしまっている。そこで、厳密に、ステップ応答がステップ応答の上限および下限を逸脱しないようにするために、定常ゲインＧの上下限値が制約条件とされてもよい。すなわち、上述したように、制約条件として、モデリング対象のシステムにおける定常ゲインがＧｍｉｎ以上かつＧｍａｘ以下であることが、先験的知識として判明している場合に、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、ｉ番目の粒子ｐ（ｉ，０）の初期値がランダムに設定され、そして、ａ０×Ｇｍｉｎ≦ｂ０≦ａ０×Ｇｍａｘを満たすように、更新される。このような制約条件を課すことにより、ステップ応答（蒸気発生量）は、図９（Ａ）に示すように、図５（Ｂ）に示す実際の蒸気発生量の変化によく一致（フィッティング）しているだけでなく、図９（Ｂ）に示すように、１２００秒を越えた範囲も含めて、ステップ応答の上限および下限の制約条件を満たしている。
【０１２３】
なお、上述では、給塵機の速度を入力とし、蒸気発生量を出力としたが、入出力の関係は、これに限定されるものではなく、種々の関係であって良い。例えば、入力としては、押込空気量や、空気供給量や、バーナ油量や、蒸気弁開度や、触媒投入量等が挙げられ、出力としては、砂層温度や、炉内温度や、ボイラドラム圧力や、炉内圧力や、酸素濃度や、一酸化炭素濃度や、ＮＯｘ濃度等が挙げられる。
【０１２４】
本実施形態のモデリング装置Ｓの適用例としてゴミ焼却炉のモデル化について説明したが、本実施形態のモデリング装置Ｓは、操作量を入力とするとともに制御量を出力とすることで、様々な制御対象のモデル化に適用可能である。本実施形態のモデリング装置Ｓは、例えば、圧延機のモデル化に適用可能であり、例えば、圧下位置やロール速度等を入力とし、スタンド間張力や板幅・板厚等を出力とすることができる。また例えば、本実施形態のモデリング装置Ｓは、マニピュレータのモデル化に適用可能であり、トルク等を入力とし、アーム速度や位置等を出力とすることができる。
【０１２５】
なお、上述の実施形態において、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うように、モデリング装置Ｓが構成されてもよい。一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。このため、このように構成されることによって、この初期値自体も未知パラメータとして取り扱われるので、初期値に左右されることなく、システムのモデルをより高精度に構築することが可能となる。
【０１２６】
例えば、上述の１入力１出力系のシステムの場合において、初期状態が０（ゼロベクトル）であることを前提にモデルが構築されている。すなわち、モデリング対象のシステムが整定している状態（平衡状態）において入力データが印加されたものと見なされており、モデルの構築に用いた例えば実績入出力データ等の入出力データは、初期状態が０であるとされている。
【０１２７】
式４において、状態ｘ（ｔ）の挙動は、一般に、式１１となることが知られている。
ｘ（ｔ）＝Φ（ｔ）ｘ（０）＋∫Φ（ｔ−τ）Ｂｕ（τ）ｄτ ・・・（１１）
ここで、積分区間は、０からｔまでであり、Φ（ｔ）は、遷移マトリックスと呼ばれ、式１２によって表される。
Φ（ｔ）＝ｅｘｐ（Ａｔ）
＝Ｉ＋Ａｔ＋Ａ２ｔ２／２！＋Ａ３ｔ３／３！＋・・・ ・・・（１２）
【０１２８】
この式１１で表される状態ｘ（ｔ）を式４−２に代入すると、式１３が得られ、初期状態ｘ（０）＝０でなければ、出力ｙｍ（ｔ）は、初期状態ｘ（０）の影響を受けることが分かる。
ｙｍ（ｔ）＝ＣΦ（ｔ）ｘ（０）＋Ｃ∫Φ（ｔ−τ）Ｂｕ（τ）ｄτ＋Ｄｕ（ｔ） ・・・（１３）
【０１２９】
モデルを構築するための入出力データを得る実際の実験において、例えばマニピュレータ等の機械系のシステムでは、システムが静止した状態から実験を開始することが可能であり、この場合では初期状態を０と見なせるが、例えばゴミ燃焼炉や化学プラント等のプラント系のシステムでは、稼動中に実験を行わなければならず、状態が絶えず変動しており、通常、初期状態が０とはならない。
【０１３０】
このため、上述のように、初期状態自体もパラメータとして扱い、モデルを構築することが好ましい。
【０１３１】
より具体的には、例えば、上述の１入力１出力系のシステムの場合において、式４および式５、あるいは、式６では、モデリング対象のシステムのモデル次数が２次であるから、初期状態ｘ（０）＝［ｘ１ ｘ２］Ｔとおくことができる。このため、Ｐ個のうちのｉ番目の粒子ｐ（ｉ，ｋ）は、式７に代えて、式１４によって表現されることになる。
ｐ（ｉ，ｋ）＝［ａ０ ａ１ ｂ０ ｂ１ ｘ１ ｘ２］Ｔ ・・・（１４）
【０１３２】
また、この場合におけるｉ番目の粒子ｐの移動速度ｖ（ｉ，０）の初期値は、ｖ（ｉ，０）＝［０ ０ ０ ０ ０ ０］Ｔに設定される。
【０１３３】
このような初期値の設定を行って、以下、上述と同様の処理によってシステムのモデルが構築される。
【０１３４】
また、上述の実施形態において、予め決められた所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値をＰＳＯ法の初期値に加えて用いるように、モデリング装置Ｓが構成されてもよい。前記所定の最適化方法として、ＰＳＯ法、最小二乗法および部分空間法等の各種の最適化手法が挙げられる。なお、最小二乗法や部分空間法が用いられる場合には、離散時間で一旦状態空間モデル式や伝達関数モデル式を求め、これを連続時間のモデルに変換される。また、ＰＳＯ法を用いる場合では、評価関数のみを単に最小化するモデルが求められても良いし、その際に、Ｊｉ＝｜｜Ｗ（Ｙｉ−Ｙｄａｔａ）｜｜＜ε（εは所定の閾値である）となる粒子ｐを保存しておいても良い。このＰＳＯ法を用いる場合では、初期値の決定およびモデルの構築の２回、ＰＳＯ法が使われることになる。一般に、モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法では、この最適化手法によって構築されるシステムのモデルの精度は、初期値に左右されることが多い。このため、ＰＳＯ法による場合でもランダムに振る範囲が試行錯誤によって決定される。このため、この試行錯誤の手間を軽減し、最適なモデルを求めるためには、最初のモデルが入出力データにある程度フィッティングしていることが望ましい。この構成では、ＰＳＯに用いる初期値として、予め決められた所定の最適化方法によって求められたモデルのパラメータの値も用いるので、より真に近い値を含む初期値からＰＳＯ法の演算処理を開始することができるから、不必要な探索を行うことなく、システムのモデルをより高精度に構築することができる。
【０１３５】
本発明を表現するために、上述において図面を参照しながら実施形態を通して本発明を適切且つ十分に説明したが、当業者であれば上述の実施形態を変更および／または改良することは容易に為し得ることであると認識すべきである。したがって、当業者が実施する変更形態または改良形態が、請求の範囲に記載された請求項の権利範囲を離脱するレベルのものでない限り、当該変更形態または当該改良形態は、当該請求項の権利範囲に包括されると解釈される。
【符号の説明】
【０１３６】
Ｓ モデリング装置
１ 演算部
２ 入力設定部
４ 記憶部
１１ モデル出力演算部
１２ 制約条件評価部
１３ 評価演算部
１４ モデル更新決定部
２１ 制約条件入力設定部
２２ モデル次数同定パラメータ入力設定部
２３ 最適化計算パラメータ入力設定部
２４ 評価関数入力設定部
４１ 制約条件記憶部
４２ モデル次数同定パラメータ記憶部
４３ 最適化計算パラメータ記憶部
４４ 評価関数記憶部
４５ モデル出力等記憶部
４６ 実績データ記憶部
【特許請求の範囲】
【請求項１】
システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング装置において、
前記モデルは、連続した時間と見なされる時間変数を含む数式であり、
予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、
所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、
前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備えること
を特徴とするモデリング装置。
【請求項２】
前記先験的情報は、前記システムにおける、定常ゲインに関する上限値および下限値のうちの少なくとも一方であること
を特徴とする請求項１に記載のモデリング装置。
【請求項３】
前記評価関数は、前記実績出力データと前記モデル出力データとの誤差に重み付けを行う重みを有し、
前記重みは、その値が大きいほど前記評価結果に大きく影響を与えるものである場合に、前記重みの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されること
を特徴とする請求項１または請求項２に記載のモデリング装置。
【請求項４】
前記モデル演算部は、前記モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うこと
を特徴とする請求項１ないし請求項３の何れか１項に記載のモデリング装置。
【請求項５】
前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによって前記モデルのパラメータの値を求めること
を特徴とする請求項１ないし請求項４のいずれか１項に記載のモデリング装置。
【請求項６】
前記モデル演算部は、所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値を前記Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄの初期値に加えること
を特徴とする請求項５に記載のモデリング装置。
【請求項７】
前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによる更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いること
を特徴とする請求項５または請求項６に記載のモデリング装置。
【請求項８】
システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング方法において、
前記モデルは、連続した時間と見なされる時間変数を含む数式であり、
予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定ステップと、
所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算ステップと、
前記制約条件設定ステップで設定された制約条件を満たし、かつ、前記評価演算ステップで演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算ステップとを備えること
を特徴とするモデリング方法。
【請求項１】
システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング装置において、
前記モデルは、連続した時間と見なされる時間変数を含む数式であり、
予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定部と、
所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算部と、
前記制約条件設定部で設定された制約条件を満たし、かつ、前記評価演算部で演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算部とを備えること
を特徴とするモデリング装置。
【請求項２】
前記先験的情報は、前記システムにおける、定常ゲインに関する上限値および下限値のうちの少なくとも一方であること
を特徴とする請求項１に記載のモデリング装置。
【請求項３】
前記評価関数は、前記実績出力データと前記モデル出力データとの誤差に重み付けを行う重みを有し、
前記重みは、その値が大きいほど前記評価結果に大きく影響を与えるものである場合に、前記重みの値は、所定の閾値以上のノイズを含む実績出力データに対応する区間に対し、前記区間を除く他の区間よりも相対的に小さくなるように設定されること
を特徴とする請求項１または請求項２に記載のモデリング装置。
【請求項４】
前記モデル演算部は、前記モデルのパラメータにおける初期値を設定してから前記モデルのパラメータの値を求める最適化方法によって前記モデルのパラメータの値を求める場合に、前記初期値自体を前記モデルのパラメータとして扱うこと
を特徴とする請求項１ないし請求項３の何れか１項に記載のモデリング装置。
【請求項５】
前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによって前記モデルのパラメータの値を求めること
を特徴とする請求項１ないし請求項４のいずれか１項に記載のモデリング装置。
【請求項６】
前記モデル演算部は、所定の最適化方法によって前記モデルのパラメータの値を求め、この求めた前記モデルのパラメータの値を前記Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄの初期値に加えること
を特徴とする請求項５に記載のモデリング装置。
【請求項７】
前記モデル演算部は、Ｐａｒｔｉｃｌｅ Ｓｗａｒｍ Ｏｐｔｉｍｉｚａｔｉｏｎ Ｍｅｔｈｏｄによる更新を行う場合に、定常ゲインの符号、定常ゲインの範囲、前記システムの安定性に関する情報および前記システムの不安定なゼロ点に関する情報のうちの少なくとも１つを用いること
を特徴とする請求項５または請求項６に記載のモデリング装置。
【請求項８】
システムのモデルを前記モデルのパラメータの値を決定することによって構築するモデリング方法において、
前記モデルは、連続した時間と見なされる時間変数を含む数式であり、
予め得られている前記システムに関する先験的情報を制約条件として設定する制約条件設定ステップと、
所定の入力データを実績入力データとして前記システムに入力した場合に前記システムから出力された実績出力データと、前記実績入力データを前記モデルに入力した場合に前記モデルから出力されたモデル出力データとの誤差を評価する評価関数の評価結果を求める評価演算ステップと、
前記制約条件設定ステップで設定された制約条件を満たし、かつ、前記評価演算ステップで演算された評価結果が最良となるように、前記モデルのパラメータの値を求めるモデル演算ステップとを備えること
を特徴とするモデリング方法。
【図１】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【図２】
【図３】
【図４】
【図５】
【図６】
【図７】
【図８】
【図９】
【公開番号】特開２０１２−１９０２７４（Ｐ２０１２−１９０２７４Ａ）
【公開日】平成２４年１０月４日（２０１２．１０．４）
【国際特許分類】
【出願番号】特願２０１１−５３４３６（Ｐ２０１１−５３４３６）
【出願日】平成２３年３月１０日（２０１１．３．１０）
【出願人】（０００００１１９９）株式会社神戸製鋼所 (5,860)
【Ｆターム（参考）】
【公開日】平成２４年１０月４日（２０１２．１０．４）
【国際特許分類】
【出願日】平成２３年３月１０日（２０１１．３．１０）
【出願人】（０００００１１９９）株式会社神戸製鋼所 (5,860)
【Ｆターム（参考）】
[ Back to top ]