(57)【特許請求の範囲】
【請求項１】
複数の入力層素子及び複数の中間層素子を有し、全ての入力層素子に中間層素子が結合されてなる全結合部分と、複数の入力層素子のうち関連の強い一部の入力層素子に中間層素子が結合されてなる疎結合部分とを備える解析可能な階層型構造のニューラルネットワークを対象として、気象条件、曜日や平日・休日の区別、プラントの操業状態、プラント構成機器の運転パターンを少なくとも含む入力因子の実績値と各種負荷の実績値とを用いて学習させ、学習済みのニューラルネットワークに、予測対象日の気象条件を予報値により他の前記入力因子とともに入力して対象プラントの電力負荷、熱負荷、空気負荷等のプラント負荷を予測し、
この負荷予測値に対する負荷種別毎の需給バランスを満足し、かつ、各プラント構成機器の入出力限界、起動・停止時間などの特性上の制約を考慮しながら、少なくとも所定期間のプラントの運用費用の最小化を目的関数とし、離散値である各制御時間毎の各プラント構成機器の起動・停止状態及び連続値である各制御時間毎の各プラント構成機器への燃料注入量を状態変数として、プラントの最適化問題を非線形混合整数計画問題として定式化すると共に、前記状態変数の初期値を計画値として生成し、
探索過程における各プラント構成機器の起動・停止状態及び燃料注入量を定常プラントシミュレータに入力して各プラント構成機器の入出力状態を前記目的関数の計算用に出力させ、前記計画値をモダンヒューリスティック最適化手法により逐次修正しながら最終的に大域最適解に近い各プラント構成機器の起動・停止状態及び燃料注入量を決定するプラントの最適運用方法であって、
前記定常プラントシミュレータは、各制御時間毎のプラント構成機器の起動・停止状態及び燃料注入量が計画値として与えられたときに、各プラント構成機器相互の接続状態、各プラント構成機器の線形または非線形な入出力特性及び各プラント構成機器の運転パターンを用いて、各制御時間毎の各プラント構成機器の定常的な入出力状態を燃料注入側の機器から末端に向かい逐次計算して出力することを特徴とするプラントの最適運用方法。
【請求項２】
請求項１に記載したプラントの最適運用方法において、前記最適化問題をParticle Swarm Optimizationを利用して解くことを特徴とするプラントの最適運用方法。
【請求項３】
請求項１または２に記載したプラントの最適運用方法を用いて過去に作成された当日の運用計画を対象として、
当日の現時点におけるプラント構成機器の運転状態等をプラントの固定的な初期状態とし、この初期状態を各機器の特性上の制約として考慮すると共に当日の残りの時間の負荷予測値に基づく需給バランスを制約条件として考慮しながら当日の運用計画を補正することを特徴とするプラントの最適運用方法。
【請求項４】
プラントの設計段階において、各プラント構成機器の最適容量の決定問題を、設定可能な離散値としての機器容量値の中から最適値を決定する組合せ最適化問題として定式化すると共に、所定期間内の対象プラントの電力負荷、熱負荷、空気負荷等の負荷パターンと各プラント構成機器の容量の設計値とを設定し、
請求項１または２に記載したプラントの最適運用方法により決定される各プラント構成機器の起動・停止状態及び燃料注入量を用いて、少なくともプラント運用費用に関する評価値が最良になるように組合せ最適化手法を用いて前記設計値を逐次変更しながら各プラント構成機器の最適容量を決定していくことを特徴とするプラントの最適設計方法。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明は、プラントの効率的な操業を可能にするために、プラントにおける電力負荷、熱負荷、空気負荷等の各種負荷をニューラルネットワークにより予測し、定常プラントシミュレータによりプラント構成機器の定常的な入出力状態を計算してプラント全体の挙動を模擬すると共に、プラントの運用費用や排出ガスを最小化するような最適運用を非線形混合整数計画問題として定式化し、この最適化問題を前記定常プラントシミュレータとの間の情報の受け渡しにより解くようにしたプラントの最適運用方法、並びに、この最適運用方法を利用して一定期間内のプラント運用費用や排出ガス等に関する評価値を最良化するようにプラント構成機器の最適容量を組み合わせ最適化問題として定式化し、これを求解するようにしたプラントの最適設計方法に関する。
【０００２】
【従来の技術及び発明が解決しようとする課題】
（１）プラントの各種負荷予測に関して
ここでは、プラントにおける各種負荷の中で、熱負荷を例として述べる。
例えば、熱源機器及び蓄熱プラントを有する蓄熱システムにおいては、翌日に必要となるであろう熱負荷の一部または全部を、夜間に冷水、氷、温水などの形で蓄熱している。ここで、翌日、当日の熱負荷を正確に予測できれば、必要な熱量を過不足なく蓄熱でき、熱源機器や蓄熱プラントの効率的な運転が可能となる。
【０００３】
従来の熱負荷予測方法としては、基準熱負荷を利用する方法、熱負荷の時系列的な推移を利用する方法、気温などの外的要因と熱負荷との相関モデルを利用する方法が知られている。
【０００４】
上記従来技術のうち、基準熱負荷を利用する方法は、月別などで基準となる熱負荷を設定し、この基準熱負荷を気温や実績熱負荷を用いて補正することで予測を行っている。しかし、同一月の中でも平日と土曜日、休日、操業状態、気象条件等によって実際の熱負荷は大きく変わるため、高精度な予測が行えないという問題がある。
また、熱負荷の時系列的な推移を利用する方法としては、例えば、過去３点の熱負荷の自己回帰モデルなどが知られている。しかし、この方法は、時間毎の変動が大きい場合や突発的な変動がある場合には適用できない。
【０００５】
また、気温などの外的要因と熱負荷との相関モデルを利用する方法には、重回帰式を用いて線形モデルを構築する方法や、近年では、ニューラルネットワークを用いて非線形モデルを構築する方法がある。
しかしながら、線形重回帰式を用いる方法では、対象となる熱負荷の非線形性を考慮することができない。
【０００６】
一方、一般的な従来型のニューラルネットワークを用いれば非線形モデル化が可能であるが、種々の入力因子や中間層ニューロン数を試行錯誤的に設定しなければならず、必ずしも最適なニューラルネットワークが得られるとは限らない。また、ニューラルネットワークの内部がブラックボックスであるため、ある入力が与えられた時にどのような理由でその予測値が得られたかという予測理由の説明ができず、運用者に不安感を与えてしまう場合がある。
このような問題は、熱負荷ばかりでなく、電力負荷や空気負荷など、他のプラント負荷にも共通する問題である。
【０００７】
このように従来の負荷予測方法は、何れも予測精度や予測理由の説明の点で問題がある。すなわち、重回帰式などの線形モデルを利用する方法では予測対象の非線形性を考慮できず、非線形モデルを利用する方法では、最適な非線形モデルが得られていない場合があるため予測精度を保証できない。また、予測理由の説明も的確に行えないので、運用者に不安感を与えてしまう等の問題があった。
【０００８】
上述したように、電力負荷、熱負荷、空気負荷等の各種プラント負荷の予測モデルを構築するにあたっては、予測精度の観点から非線形関係のモデル化を可能とし、かつ、運用者に与える安心感、信頼性の観点から予測理由の説明が可能であることが望ましい。
しかし、従来では、予測精度を高めたい場合にはニューラルネットワークを適用し、予測理由の説明という要求を満足させたい場合には重回帰式を適用することが多く、予測精度向上と理由説明という二つの要求を同時に満足させる予測方法は実現されていない現状である。
【０００９】
そこで、請求項１記載の発明においては、ニューラルネットワークによる非線形の予測モデルにより高精度な負荷予測を可能とし、更に、予測モデルから出力される負荷予測値と入力因子との因果関係、相関関係等による予測理由の説明を容易にしたプラント負荷の予測方法を用いるものである。
【００１０】
（２）定常プラントシミュレータ及びプラントの最適運用方法に関して
従来、各種プラント構成機器の起動・停止まで考慮したうえで最適化手法を用いてプラントの最適運用を求める方法は種々提案されているが、プラント構成機器の入出力特性や運転パターンを線形な一般的な方程式として扱い、混合整数計画問題として定式化したうえで問題を解くか、あるいは、プラント構成機器の入出力特性の非線形性（非連続、微分不可能等も含む）や運転パターンを考慮したうえで、エキスパートシステムやファジイ推論等の非線形も扱えるヒューリスティックな手法によりプラントの最適運用を求める方法が一般的であった。
【００１１】
しかしながら、各種プラント構成機器の入出力特性は非線形性を有するほか、時間帯によって異なった運転パターンをとることから、最適化手法で直接扱えるように一般的な方程式として定式化することは通常、困難である。
従って、プラントの最適運用は最適化問題の中でも最も求解が困難な非線形混合整数計画問題として定式化することができる。しかし、プラント構成機器の起動・停止などの離散値と構成機器の燃料注入量などの連続値とを総合的に考慮したうえで、各プラント構成機器の特性上の制約や負荷種別ごとの需給バランス等の制約条件、及び、プラント運用費用や排出ガスの最小化等の目的関数と共に、非線形性を有する問題に対して数学的に最適化できる方法はこれまで存在しなかった。
【００１２】
そこで請求項１に記載した発明においては、プラント構成機器の非線形な入出力特性や複雑な運転パターンを考慮したうえで、プラント全体の挙動を模擬し、各構成機器の定常的な入出力状態を出力するようにした定常プラントシミュレータを使用し、特に、請求項１〜３に記載した発明においては、プラントの最適運用を非線形混合整数計画問題として定式化し、前記定常プラントシミュレータの出力を用いて、プラントの非線形特性及び運転パターンの変更等を考慮した上で、大域最適解に近い解をParticle Swarm Optimization等の最適化手法を利用して求めることによりプラント構成機器の起動・停止状態及び燃料注入量を決定するようにしたプラントの最適運用方法を提供しようとするものである。
【００１３】
（３）プラントの最適設計方法に関して
従来、プラントの設計方法としては、上記（２）で述べた混合整数計画問題による定式化やヒューリスティックな手法により構築された最適運用方法を用いて、プラントの負荷予測値を様々なシナリオにより変更しつつ最適化を何回も実行しながら、各プラント構成機器の適切な容量を決定する方法が知られている。
【００１４】
このような従来の最適設計方法では、混合整数計画問題として定式化した最適運用方法を用いた場合は、プラントの非線形特性や運転パターンの変更を考慮できないため、実際のプラントに適合した現実的な解を生成することが困難である。
また、ヒューリスティックな手法を用いた場合は、プラントの非線形特性や運転パターンの変更を考慮できるものの、制約条件を満たす解を生成できるという程度であり、より大域最適解に近い解を生成できる保証がまったくない。
更に、負荷予測自体の精度も満足できるものとは言えないため、この負荷予測値を基準とした最適設計にも自ずから限界があった。
【００１５】
そこで請求項４に記載した発明は、プラントの非線形特性や運転パターンの変更等を考慮した上で、前記最適運用方法を利用して一定期間内のプラント運用費用や排出ガス等に関する評価値を最良化するようにプラント構成機器の最適容量を組み合わせ最適化問題として定式化し、大域最適解に近い解を求めるようにしたプラントの最適設計方法を提供しようとするものである。
【００１６】
【課題を解決するための手段】
上記課題を解決するため、請求項１に記載したプラントの最適運用方法は、複数の入力層素子及び複数の中間層素子を有し、全ての入力層素子に中間層素子が結合されてなる全結合部分と、複数の入力層素子のうち関連の強い一部の入力層素子に中間層素子が結合されてなる疎結合部分とを備える解析可能な階層型構造のニューラルネットワークを対象として、気象条件、曜日や平日・休日の区別、プラントの操業状態、プラント構成機器の運転パターンを少なくとも含む入力因子の実績値と各種負荷の実績値とを用いて学習させ、学習済みのニューラルネットワークに、予測対象日の気象条件を予報値により他の前記入力因子とともに入力して対象プラントの電力負荷、熱負荷、空気負荷等のプラント負荷を予測し、
この負荷予測値に対する負荷種別毎の需給バランスを満足し、かつ、各プラント構成機器の入出力限界、起動・停止時間などの特性上の制約を考慮しながら、少なくとも所定期間のプラントの運用費用の最小化を目的関数とし、離散値である各制御時間毎の各プラント構成機器の起動・停止状態及び連続値である各制御時間毎の各プラント構成機器への燃料注入量を状態変数として、プラントの最適化問題を非線形混合整数計画問題として定式化すると共に、前記状態変数の初期値を計画値として生成し、
探索過程における各プラント構成機器の起動・停止状態及び燃料注入量を定常プラントシミュレータに入力して各プラント構成機器の入出力状態を前記目的関数の計算用に出力させ、前記計画値をモダンヒューリスティック最適化手法により逐次修正しながら最終的に大域最適解に近い各プラント構成機器の起動・停止状態及び燃料注入量を決定するプラントの最適運用方法であって、
前記定常プラントシミュレータは、各制御時間毎のプラント構成機器の起動・停止状態及び燃料注入量が計画値として与えられたときに、各プラント構成機器相互の接続状態、各プラント構成機器の線形または非線形な入出力特性及び各プラント構成機器の運転パターンを用いて、各制御時間毎の各プラント構成機器の定常的な入出力状態を燃料注入側の機器から末端に向かい逐次計算して出力するものである。
【００１９】
請求項２に記載したプラントの最適運用方法は、請求項１に記載したプラントの最適運用方法において、前記最適化問題をParticle Swarm Optimizationを利用して解くものである。
【００２０】
請求項３に記載したプラントの最適運用方法は、請求項１または２に記載したプラントの最適運用方法を用いて過去（例えば前日）に作成された当日の運用計画を対象として、当日の現時点におけるプラント構成機器の運転状態等をプラントの固定的な初期状態とし、この初期状態を各機器の特性上の制約として考慮すると共に当日の残りの時間の負荷予測値に基づく需給バランスを制約条件として考慮しながら当日の運用計画を補正するものである。
【００２１】
請求項４に記載したプラントの最適設計方法は、プラントの設計段階において、各プラント構成機器の最適容量の決定問題を、設定可能な離散値としての機器容量値の中から最適値を決定する組合せ最適化問題として定式化すると共に、所定期間内の対象プラントの電力負荷、熱負荷、空気負荷等の負荷パターンと各プラント構成機器の容量の設計値とを設定し、請求項１または２に記載したプラントの最適運用方法により決定される各プラント構成機器の起動・停止状態及び燃料注入量を用いて、少なくともプラント運用費用に関する評価値が最良になるように組合せ最適化手法を用いて前記設計値を逐次変更しながら各プラント構成機器の最適容量を決定していくものである。
【００２２】
【発明の実施の形態】
以下、図に沿って各発明の実施形態を説明する。
１．まず、請求項１に記載したプラント負荷の予測方法に関する実施形態を説明する。
本実施形態では、対象となるプラントの電力負荷、熱負荷、空気負荷等の各種負荷を予測する当たり、従来から用いられている重回帰式を用いた線形モデルや一般的なニューラルネットワークによる非線形モデルを用いるのではなく、例えば、本出願人の先願である特願２０００−１６６５２８「ニューラルネットワーク及びその学習方法、解析方法並びに異常判定方法」、同２０００−２３０６６５「ニューラルネットワークの学習方法」、同２０００−２２７０５７「ニューラルネットワークの構造」等に記載された解析（任意の入力因子と出力との因果関係、相関関係の解析）が可能なニューラルネットワークを用いることとした。
【００２３】
以下、本実施形態における予測モデルとして用いる階層型ニューラルネットワークの一般的な構造及びその学習方法を、特願２０００−１６６５２８、同２０００−２３０６６５に基づいて説明する。
【００２４】
まず、ニューラルネットワーク構造を説明する。
通常の階層型ニューラルネットワークは、入力層素子と中間層素子が全て結合している（全結合部分という）が、本実施形態において使用する階層型ニューラルネットワークは、図１に示すように任意の入力層素子と任意の中間層素子間だけが結合している。
【００２５】
すなわち、このニューラルネットワークは、全ての入力層素子と結合している中間層素子からなる全結合部分１１と、一部の入力層素子と結合している中間層素子からなる疎結合部分１２とからなっている。
このように、一部の入力層素子との間の重み（重み係数または結合係数）の値を０にした疎結合部分１２を設けることで、従来の階層型ニューラルネットワークと完全な互換性を持つことができる。
つまり、全結合部分のみからなる従来のニューラルネットワークの一部の重みの値を０とすることで図１に示すようなニューラルネットワークを実現できるため、プログラム上は従来のニューラルネットワークをそのまま使用することができる。
【００２６】
例えば、
ｙ＝ａｘ１＋ｂｘ２＋ｃｘ１ｘ２＋ｄ
（ｘ１，ｘ２は入力変数、ｙは出力変数、ａ，ｂ，ｃ，ｄは係数）
として表現されるような従来の回帰式では、
▲１▼出力の要因が明確である（ｘ１，ｘ２の要因によりｙが変化することが明確である）、
▲２▼各要因の働きが明確である（ｘ１，ｘ２の独立した成分（回帰式の右辺第１項、第２項）とその相互成分（同第３項）から成り立つ）、
▲３▼各要因の入出力に対する影響度合いが明確である（ｘ１，ｘ２，ｘ１ｘ２の影響度合いはそれぞれａ，ｂ，ｃであり、固定分としてｄがある）、
といった特徴があり、これらの特徴によって内部解析が容易になっている。
【００２７】
従って、図１に示したように解析可能な疎結合部分１２と精度を保証する全結合部分１１とからなる階層型ニューラルネットワークを用いれば、回帰式の解析の容易さを既存のニューラルネットワークにより実現することができ、前述したような従来の階層型ニューラルネットワークとの互換性によって既存のシステムへの適用も容易に行うことができる。
【００２８】
ここで、図１に示したニューラルネットワークの学習方法には種々の方法があるが、一例として、特願２０００−１６６５２８に記載された学習方法を図２のフローチャートにて沿って説明する。
まず、第１ステップＡ１は、通常のニューラルネットワークの重み初期化処理である。具体的には、通常のニューラルネットワークの各層の素子間の全ての重みに対して初期値を小さい数の乱数で与える。ここで、プログラム的に、入力層と中間層との間の任意の結合がないニューラルネットワーク構造を定義してもよい。この場合は、以下の第２，第４ステップの処理は不要である。
【００２９】
第２ステップＡ２は、重みを初期化したニューラルネットワークに対し、図１に示したニューラルネットワーク構造に変更するための処理である。すなわち、任意の入力層素子と中間層素子との間の結合を削除する。ここで、結合を削除するための最も簡単な方法は、任意の重みの値を０に置き換える方法である。
【００３０】
第３ステップＡ３は、通常のニューラルネットワークの重み修正量の計算である。学習誤差を評価する評価関数が小さくなるように、入力層素子と中間層素子との間の重みの修正量を計算する。この場合の評価関数の一例を、以下の数式１に示す。
【００３１】
【数１】
Ｊ＝１／２・（ｏ−ｔ）２
【００３２】
なお、数式１において、Ｊ：評価関数、ｏ：ニューロ出力、ｔ：教師値（学習目標値）である。
【００３３】
第４ステップＡ４は、図１に示したニューラルネットワーク構造のための重み修正量の計算ステップである。第３ステップＡ３の計算により、本来的に重みが０である結合に対して修正量が計算され、元の重み（ｗ＝０）に修正量△ｗが加えられる結果、重みが再構築されることがある。それを防止するために、第４ステップＡ４ではこのような重みの修正量を強制的に０にする。
【００３４】
第５ステップＡ５は、重みの修正処理である。第３ステップ、第４ステップを経て計算された最終的な修正量に従って、入力層素子と中間層素子との間の重みを修正する。
重みの修正量を△ｗｉｊ、重みをｗｉｊ、学習係数をαとすると、数式２により重みを修正することができる。
【００３５】
【数２】
ｗｉｊ＝ｗｉｊ＋α△ｗｉｊ
【００３６】
なお、第３ステップＡ３以降の処理は、学習誤差が規定値以下になって学習終了が確認されるまで繰り返し行われる（ステップＡ６）。ここで、学習終了の判断は、評価関数や全学習データに対する誤差が規定値以下になったか否か、もしくは、学習回数が所定回数に達したか否かにより判断することができる。
【００３７】
上述した学習方法では、図１における疎結合部分１２と全結合部分１１とを区別せずに一括して同時に学習している。
その結果、本来ならば疎結合部分１２の中間層素子につながるはずの部分（入力層素子）が全結合部分１１に結合されてしまったり、疎結合部分１２の結合が小さく、全結合部分１１の結合が大きくなりすぎたりして、ニューラルネットワークの構造が解析困難になる場合がある。
【００３８】
上記の点に鑑み提案された特願２０００−２３０６６５記載の学習方法を、以下に説明する。
すなわち、図３に示すように、第１ステップＢ１として、一の疎結合部分に単独でさせる学習を全ての疎結合部分について行う。次に、第２ステップＢ２として、第１ステップＢ１により学習済みの全ての疎結合部分を結合して第１ネットワークを生成する。
【００３９】
次いで、第３ステップＢ３として、第２ステップＢ２により生成された第１ネットワークに学習を行って第２ネットワークを生成する。
そして、第４ステップＢ４として、第３ステップＢ３により生成された第２ネットワークに全結合部分を結合して第３ネットワークを生成し、第５ステップＢ５として、第４ステップＢ４により生成された第３ネットワークに学習を行って学習誤差を低減し、最終的なニューラルネットワークを生成する。
【００４０】
図４は第１ステップＢ１を経たニューラルネットワークを示しており、１２Ａ，１２Ｂ，１２Ｃはそれぞれ疎結合部分である。これらの疎結合部分１２Ａ，１２Ｂ，１２Ｃは関連の強い入力因子同士が一つの疎結合部分に含まれるように定義されており、図示例では関連の強い入力１，２が一つの疎結合部分１２Ａに含まれている。
各疎結合部分１２Ａ，１２Ｂ，１２Ｃは所定の学習データを与えて各々が単独で学習を行っており、各入力因子と出力との関係が強く学習されている。ここで、各疎結合部分１２Ａ，１２Ｂ，１２Ｃは比較的小規模のネットワークであるため、学習速度も速く、局所解へ陥る危険性も少ない。
【００４１】
図５は、第２ステップＢ２により、図４の疎結合部分１２Ａ，１２Ｂ，１２Ｃを一つに結合して形成された第１ネットワーク１２Ｄを示している。
なお、このネットワークについて後述の第３ステップＢ３により学習を行えば、ネットワークの構造上は、図５は第２ネットワークを示したものでもあると言える。
【００４２】
この場合、複数の疎結合部分を単純に結合するだけではニューラルネットワークは正常に動作しない。
例えば、疎結合部分が図示するように三つあり、学習データの出力値の範囲が０．１〜０．９であったとする。
【００４３】
第１ステップＢ１において各疎結合部分が良好に学習できたとすると、入出力特性を示す図６のシグモイド関数（ニューラルネットワークを構成する各素子（ニューロン）の入出力特性としては、通常、シグモイド関数の入出力特性が使用される）によれば、各疎結合部分の入力値の範囲は−２．２〜＋２．２となる。この場合、三つの疎結合部分を単純に結合すると、入力値の範囲は−６．６〜＋６．６となってしまい、その結果、出力値の範囲も０〜１に変わってしまう。
これを回避するために、入力値の範囲は複数の疎結合部分を結合した後も−２．２〜＋２．２に保つ必要がある。
【００４４】
そこで、結合の方法は各種考えられるが、ここでは以下の方法を用いる。
（１）中間層と出力層との間の結合係数の値を疎結合部分の数ｎで割る。
（２）各疎結合部分を結合する。
（３）すべてのパターンを想起し、出力値が所定の範囲に入っていなければ、その範囲内に収まるように結合係数を修正する。例えば、出力値の範囲が０．１〜０．９である場合には入力値が−２．２〜＋２．２になるように、出力値の範囲が０．４〜０．６である場合には入力値が−０．４〜＋０．４になるようにする。
これらの対応は、図６のシグモイド関数の入出力特性によって決定される。但し、この誤差が小さい場合には、次の第３ステップＢ３の学習だけでも修正は可能である。
【００４５】
図３の第２ステップＢ２により生成された第１ネットワーク（図５におけるネットワーク１２Ｄ）は、複数の疎結合部分を機械的に結合したものに相当するので、学習誤差が比較的大きい。
そこで、第３ステップＢ３では、各疎結合部分のバランス調整を行うために、第１ネットワークに再度、学習させて第２ネットワークを生成する。
また、第５ステップＢ５では、学習誤差を更に低減させるために再度、学習を行う。
【００４６】
第３ステップＢ３により生成された第２ネットワークは、疎結合部分１２Ａ，１２Ｂ，１２Ｃのみから構成されている。つまり、各入力同士の相互作用が欠落しているため、誤差の低下にも限界がある。
そこで、第４ステップＢ４では、各入力同士の相互作用を反映させた全結合部分１１を追加して第３ネットワーク（図１のニューラルネットワーク）を形成し、更に第５ステップＢ５ではこの第３ネットワークを対象として再度、学習を行う。これにより、最終的に学習誤差が低減され、かつ解析容易なニューラルネットワークが生成される。
【００４７】
次に、プラントの熱負荷を予測する場合を例にとって、本発明の実施形態を説明する。
（１）予測モデル構築
まず、非線形性を有する熱負荷を予測するために、図７に示すような予測モデルを構築する。この予測モデルは図１の予測モデルの変形例であり、図７に示すような全結合部分１１及び疎結合部分１２からなる階層型ニューラルネットワークである。
【００４８】
学習時に予測モデルに入力される入力因子としては、出力である熱負荷と関連がある因子、例えば、予測対象日の毎時気温、最高気温、最低気温、最小湿度、天候、日射量などの気象条件に関する因子、及び、曜日や平日・休日（土曜日は平日・休日の何れかに含まれる）の区別、イベントなどの有無及び種別、プラントの操業状態、プラント構成機器（空調機器等）の運転パターンのごとく熱負荷パターンを大まかに決定付ける因子の全部または一部を用いる。これらの入力データには実績値が用いられる。
また、出力には熱負荷の実績値を用いることとし、上記入力因子及び出力をニューラルネットワークに与えて前述した学習方法により学習させ、予測モデルを構築する。
【００４９】
図７の例では、入力層に与える入力因子を最高気温、最小湿度、曜日とし、翌日の熱負荷を予測して出力する構成となっている。
また、中間層は、最高気温が入力される入力層素子のみに結合して気温成分を出力する素子と、最小湿度が入力される入力層素子のみに結合して湿度成分を出力する素子と、曜日が入力される入力層素子のみに結合して曜日成分を出力する素子と、すべての入力層素子に結合して相互作用成分を出力する素子とからなっている。ここで、出力層素子は単一である。
【００５０】
なお、翌日の毎時の熱負荷を予測する場合には、ニューラルネットワークの出力層素子をそれぞれ各時間の熱負荷に割り当てて２４個の出力層素子を有するニューラルネットワークを用いれば良い。この場合、学習時には、入力層素子に例えば一定時間ごとの気温や湿度等を入力し、１番目の出力層素子に０時の熱負荷実績値を与えて学習させ、実際の予測時には、前記入力因子の予報値を入力層素子に与えて１番目の出力層素子から０時の熱負荷予測値を出力させればよい。
このように出力層素子を２４個にすれば２４点全体で１日の熱負荷パターンをモデル化できるため、時間ごとにニューラルネットワークを構築する（予測モデルを２４個用いる）方法に比べて予測モデル数を削減することができ、計算機負荷の削減が可能になる。
【００５１】
（２）予測実行
上記（１）で構築された予測モデルを用いて、翌日の熱負荷の予測計算を実行する。このとき、入力因子としての最高気温、最小湿度等の気象条件に関する入力データには、翌日の予報値を用いる。また、プラントの操業状態、プラント構成機器の運転パターン等を入力する場合には、それらの計画値を用いれば良い。
【００５２】
（３）予測理由説明
予測結果が得られた理由を説明する。予測理由の説明には、重回帰モデルでは、偏回帰係数、ニューラルネットワークでは入出力の感度を用いるのが一般的であるが、本実施形態では全結合部分１１及び疎結合部分１２を有するニューラルネットワークを用いることにより、どのような理由で当該予測値が得られたかという予測理由の説明が可能である。
図２の例では、熱負荷予測値を最高気温に関連した成分、最小湿度に関連した成分、及び曜日に関連した成分と、これらの最高気温、最小湿度、曜日が相互に関連する相互作用成分とに分割して提示することができる。
【００５３】
以下、予測理由の説明に当たり、任意のデータを入力したときの、中間層素子から出力層素子に伝達される情報量から入力因子の出力への影響を解析する方法について述べる。
この解析方法は、例えば特願２０００−１６６５２８の請求項１１，１２等に記載されている。
【００５４】
始めに、特願２０００−１６６５２８の請求項１１に記載された解析方法について説明する。
図８に示すニューラルネットワークにおいて、中間層から出力層に伝達される情報は、中間層出力Ｏと重みｖとの積和である。つまり、この情報は数式３によって表される。
【００５５】
【数３】
0004059014.tif 000002
【００５６】
ここで、｜ｖｉＯｉ｜が最も大きい値を出力する中間層素子は、出力への影響が最も強く、更にその中間層素子に結合されている入力因子の影響も強いと言える。
例えば、図８において、中間層素子１による｜ｖ１Ｏ１｜が最も大きい場合には入力１の影響が強く、中間層素子２による｜ｖ２Ｏ２｜が最も大きい場合には入力１，２の相互作用分が強いと言える。
従って、各中間層素子から出力層素子に伝達される情報量の大きさを検出することで、中間層素子、出力層素子及び入力層素子の間の結合の強さを知ることができ、ニューラルネットワークにおける入力因子の出力への影響を解析することができる。
【００５７】
図７の例について言えば、例えば中間層の気温成分の出力と、当該成分と出力層との間の重みとの積が、他の成分に比べて大きい場合には、予測値（翌日の熱負荷）に占める気温成分（つまり最高気温）の影響が大きいと考えられるから、その予測値が得られた理由は最高気温が大きく影響しているという説明が可能になる。
また、中間層の相互作用成分の出力と、当該成分と出力層との間の重みとの積が、他の成分に比べて大きい場合には、予測値に占める相互作用成分（つまり最高気温、最小湿度、曜日）の影響が大きいと考えられるから、その予測値が得られた理由を一義的に断定することはできず、最高気温、最小湿度、曜日の相互作用によるものであるという説明が可能になる。
【００５８】
次に、特願２０００−１６６５２８の請求項１２に記載された解析方法について説明する。
この解析方法は、任意のデータを入力したときの、入力データと中間層素子から出力層素子へ伝達される情報との相関関係から、入力因子の出力への影響を解析する方法である。
【００５９】
すなわち、図８のニューラルネットワークにおいて、入力１，入力２の値を常に等しくして（入力１，入力２）＝（０，０），（０．２，０．２），……，（０．８，０．８），（１，１）の０．２刻みの複数のデータを入力し、そのときに中間層素子１〜３から出力層素子へ伝達される情報量が図９のようになったとする。つまり、中間層素子１については出力層素子へ伝達される情報量が次第に増加し、中間層素子２についてはほぼ一定であり、中間層素子３については次第に減少していくとする。これから、以下のことが判る。
【００６０】
▲１▼入力１のみに結合されている中間層素子１は正の相関、つまり、入力１は出力に対し正の相関がある。
▲２▼中間層素子２は出力にほとんど影響を与えない。つまり、中間層素子２に結合している入力１，２の相互作用はほとんどない。
▲３▼入力２のみに結合されている中間層素子３は負の相関、つまり、入力２は出力に対し負の相関がある。
また、以上の▲１▼〜▲３▼により、次の▲４▼が導かれる。
▲４▼出力を大きくするためには、入力１に大きな値を入力し、入力２に小さな値を入力すればよい。
【００６１】
上記のように、入力値に対する中間層素子の挙動（それぞれの中間層素子から出力層へ伝達される情報量）を検出することにより、未知データ（学習していない入力ｘ１，ｘ２のパターン）に対しても、各入力因子や中間層素子の働きが判っているために出力値を容易に推定することができる。
【００６２】
図７の例に則して言えば、最高気温や最小湿度等の入力因子について、中間層の気温成分や湿度成分から出力層へ伝達される情報に関する図９のような相関関係が予め判っていれば、過去に学習していない最高気温や最小湿度等のデータが与えられた場合にも、翌日の熱負荷がどのような値になるかをある程度予測することができ、また、当該予測値が得られる理由についても説明することができる。
【００６３】
以上のように本実施形態によれば、予測モデルの非線形化と予測理由の説明とを両立させながら、熱負荷等の非線形性を有する負荷を高精度に予測することが可能である。なお、本発明は、熱負荷以外のプラントの電力負荷、空気負荷等の予測にも適用可能であることは言うまでもない。
【００６４】
２．次に、請求項１における定常プラントシミュレータに関する実施形態を説明する。ここで、定常プラントシミュレータは、プラント構成機器相互の接続状態、線形または非線形な入出力特性、時間帯に応じて変化する運転パターン（運転パターンの変化により入出力特性が変化する）等のプラントの状態の計算に必要な特性を保存したデータベースを用いて、一定の時間間隔をおいた各制御時間毎のプラント構成機器に対する燃料注入量、プラント構成機器の起動・停止状態が計画値として後述する最適化部から与えられた時に、各制御時間毎の各構成機器毎の定常的な入出力状態（ここでは、過渡的な入出力状態については考えず、機器の出力が安定した定常時のみを考慮する）を計算する。そして、ある構成機器の出力を次段の構成機器の入力として逐次計算していくことにより、最終的にプラント全体の定常的な入出力状態を模擬するものである。このようにして模擬されたプラント構成機器の入出力状態は、プラントの最適運用を決定するための目的関数の計算用に、後述する最適化部に送られることになる。
【００６５】
例えば、図１０に示すように、ガスタービン２１、排ガスボイラ２２、スチームタービン２３、コンプレッサ２４、ボイラ２５、吸収式冷凍機２６、ガス焚冷凍機２７、蓄熱層２８から構成され、空気負荷、電気負荷、蒸気負荷、熱（空調）負荷を得るためのエネルギープラントを想定したとする。定常プラントシミュレータでは、上記各機器２１〜２８の図示するような構成状態、入出力特性、運転パターンが与えられ、更に、制御時間毎のガスタービン２１、ボイラ２５、ガス焚冷凍機２７に対する燃料注入量及び各機器２１〜２８の起動・停止状態が与えられた時に、例えばガスタービン２１から出力される排ガスの量やボイラ２５から出力されるスチームの量を計算する。そして、これらを次段の排ガスボイラ２２や吸収式冷凍機２６等の入力として、各機器２２，２６，……の出力を逐次計算していくことにより、このエネルギープラント全体の定常的な入出力状態を模擬する。
【００６６】
この際、各プラント構成機器の非線形な入出力特性は、以下のような方法を用いて表現する。
▲１▼最小二乗法により、その係数が調整された２次までの多項式
▲２▼応答曲面法により、その係数が調整された一般的な非線形方程式
▲３▼ニューラルネットワーク
【００６７】
なお、▲２▼の応答曲面法とは、例えば、Andre I. Khuri, and John A. Cornell, Response Surfaces Designs and Analyses, Second Edition, Marcel Dekker Inc., 1996等に説明されている。
また、ニューラルネットワークは、D. E. Rumelhart and J. L. McCelland(Eds.), Parallel Distributed Processing, MIT Press, 1986.等に説明されているところの、多層パーセプトロンを用いたモデルのように、ニューラルネットワークに対する入力と出力との写像関数を自動的に生成する全てのニューラルネットワークを包含する。場合によっては、請求項１の発明の実施形態として説明した構造化ニューラルネットワークを用いても良く、その場合には、各プラント構成機器の入出力関係についての理由説明が容易になる利点が得られる。
【００６８】
ここで、前記の応答曲面法は、タグチメソッドと同様に、製品のばらつきなどを考慮して製品開発を行う品質工学分野において、特にアメリカで実用化されている手法である。
応答曲面法は、適切な実測値の組合せに対する限られた数の応答を求め、予測変数から予測される応答の近似関数を作成する方法であり、この応答曲面法により作成した近似関数を用いて各機器の入出力特性を表現することができる。ここで、応答曲面を求める方法は、ソフトウェア・シミュレーションを用いても、実際に試作したハードウェアによる実験に基づいても、あるいはプラント実測値を用いても良い。
【００６９】
近似関数の関数形に制限はないが、近似関数としては以下のような形態が考えられる。
▲１▼１次関数
数式４の線形関数を利用する場合であり、統計で利用されている重回帰に相当する。
【００７０】
【数４】
ｙ＝β０＋β１ｘ１＋β２ｘ２＋……
【００７１】
ここで、
βi：係数
ｘi：予測変数（設計変数）
である。
【００７２】
▲２▼２次多項式
数式５の２次多項式を利用する場合である。
【００７３】
【数５】
ｙ＝β０＋β１ｘ１＋β２ｘ２＋……＋β１１ｘ１２＋β２２ｘ２２＋……＋β１２ｘ１ｘ２＋β１３ｘ１ｘ３＋……
【００７４】
同様にして、３次以上の多項式を用いることも可能である。また、２次以上の多項式には、以下に示すように直交実験とするために直交多項式を用いる場合もある。
【００７５】
なお、近似関数は、実験点を選定し、その実験点を用いて作成する。
【００７６】
実験点の選定には、以下のような方法がある。
▲１▼直交計画
タグチメソッドで用いる直交表を用いる場合であり、通常は１次多項式（線形重回帰）による近似に用いられる。しかし、２次以上の多項式回帰モデルでは、L27などの多水準の直交表を用いる。また、この際に、近似関数として直交関数を用いることにより、近似関数の作成において必要な係数の選定を独立に行えるようになるので、有効である。
【００７７】
▲２▼中央複合計画(Central Composite Design: CCD)
ＣＣＤは、例えば２変数の場合に図１１に示す点を実験点とする方法であり、２次多項式による近似で有効であることが知られている。
【００７８】
▲３▼Ｄ最適水準
計算機支援の実験計画として、以下の数式６で示されるＤ最適性指標を最大化するｐ個の実験点を用いる方法である。
【００７９】
【数６】
0004059014.tif 000003
【００８０】
ここで、Ｘ：実験点の行列(ｎ×ｋ)
ｎ：実験点の組の総数
ｋ： 設計変数の数
ｐ：近似に必要な実験点数
である。
【００８１】
なお、以下に示すように、実験点の選定と用いる近似関数とには密接な関係がある。
・直交計画による実験点の選定と直交関数による近似は、計算が簡単であり、更に高次の多項式にも対応できるという利点を有している。しかし、例えば、Ｌ２７直交表＋２次直交関数のＤｅｆｆは0.462であり、余りよい近似とはいえない。
・ＣＣＤは２次多項式だけに有効である。ＣＣＤ＋２次多項式のＤｅｆｆは0.996であり、かなり良い近似となる。
・上記以外の場合には、Ｄ最適規準が必要である。
【００８２】
２次以上の多項式も変数変換することにより１次多項式にすることができるため、１次多項式、つまり統計の分野で線形重回帰式に利用されている方法を用いて近似関数を作成することができる。
（１）モデルの係数の求め方
数式４，５のモデルの各係数βｉは、重回帰で用いられるように最小２乗法で求めることができる。
（２）モデルの適切性の評価
モデルが適切かどうかは、数式７に示す決定係数（回帰モデルの適合性を示す指標）Ｒ２により評価することができる。
数式７において、ＳＳＲは回帰二乗和、Ｓｙｙは応答ｙの平均値まわりの変動である。回帰式が応答に完全に一致する場合にはＳＳＲ，Ｓｙｙが完全に同一になり、決定係数Ｒ２は１となる。この決定係数Ｒ２は０〜１の値をとり、一般には１に近いほど適切な回帰モデルということができる。
なお、モデルの各係数はｔ検定により有効性を確認することができる。
【００８３】
【数７】
0004059014.tif 000004
【００８４】
上記のように、この実施形態にかかる定常プラントシミュレータは、図１０に示したような各機器の線形または非線形な入出力特性を最小二乗法や応答曲面法、ニューラルネットワーク等により表現し、これらの入出力特性と、時間帯に応じた各機器の運転パターン（始動後に定常状態に達するまでの時間等の特性も含む）等をデータベースとして予め備えている。そして、各制御時間毎の燃料注入量や各機器の起動・停止状態を入力情報として、燃料の注入側からプラントの末端方向に向かって各機器の入出力状態（入出力されるガスやスチーム、水等の量やこれらの入出力時間）を逐次計算することにより、プラント全体の定常的な入出力状態を模擬する。
【００８５】
３．次に、請求項１，２，３に記載したプラントの最適運用方法の実施形態を説明する。
本実施形態では、最適化部が、計画案として各制御時間毎の各プラント構成機器の起動・停止状態、燃料注入量を定常プラントシミュレータに渡し、定常プラントシミュレータ側では、これらの入力情報を用いて前述した動作により求めた各機器の入出力状態を最適化部に返すと共に、最適化部では、各機器の入出力状態を用いて所定の目的関数（プラント運転費用の最小化やガス排出量の最小化等）を満足するような各機器の起動・停止状態及び燃料注入量を各種の最適化手法により探索していくものである。
【００８６】
図１２は、定常プラントシミュレータを用いてプラントの最適運用方法を決定する際の概念図を示している。なお、図１２における最適化部３０及び定常プラントシミュレータ４０は、何れも計算機のハードウェア及びソフトウェアにより実現されるものである。
まず、請求項１記載の発明により各種プラント負荷の予測値が得られたら、各負荷種別毎の負荷予測値に対するエネルギー供給（需給バランス）を満足し、各機器の特性上の制約を考慮しながら（これらを制約条件とする）、所定期間、例えば１日のプラントの運用費用（燃料費を含む）及びガス排出量（例えばＣＯ２排出量）を最小化することを目的関数として（ペナルティとして需給バランス不均衡量及び機器特性制約逸脱量を設定する）、プラントの最適化問題を非線形混合整数計画問題として定式化する。
【００８７】
そして、各制御時間毎の各プラント構成機器の起動・停止状態（起動か停止かという離散値）及び燃料注入量（連続値）の初期値を生成し、前記定常プラントシミュレータに対して、各探索過程におけるプラント構成機器の起動・停止状態、及び、燃料注入量を逐次与えることにより、定常プラントシミュレータから各プラント構成機器の入出力状態を出力として得ると共に、起動・停止状態及び燃料注入量の計画値を最適化手法により逐次修正しながら、最終的に大域最適解に近い各機器の起動・停止状態及び燃料注入量を決定することにより、プラントの最適運用方法を決定する。
【００８８】
以下に、最適化の定式化を示す。
(１)状態変数
状態変数は、以下のプラント量とする。
▲１▼各制御時間毎の各プラント構成機器（例えばガスタービン，ボイラ，冷凍機，コジェネ等）の起動・停止状態
▲２▼各制御時間毎の燃料注入量（上記各機器への燃料注入量）
【００８９】
(２)目的関数
目的関数は、以下のような項からなるとする。
▲１▼プラント運用費用の最小化(ｆ１とする)
▲２▼ＣＯ２排出量の最小化(ｆ２とする)
▲３▼ペナルティ(需給バランス不均衡量，機器特性制約逸脱量) (ｆ３とする)実際には、上記各項からなる以下の目的関数を用いる。
ｆ＝ｗ１ｆ１＋ｗ２ｆ２＋ｗ３ｆ３
ここで、ｗｉ：各項に対する重み
であり、任意に設定される。
【００９０】
(３)制約条件
▲１▼各負荷種別毎の需給バランス
電力系、熱系、空気系など、それぞれの負荷種別に対する需要と供給とのバランスに関する制約であり、この需給バランスを決定するに当たって負荷予測値が考慮される。
▲２▼各機器の特性上の制約
各プラント構成機器の入出力限界、起動・停止時間などの特性上の制約である。
【００９１】
(４)求解アルゴリズム
定常プラントシミュレータによる計算結果（各制御時間毎の各プラント構成機器の入出力状態）を利用しながら、起動・停止状態及び燃料注入量の最適解を探索していく最適化手法としては、モダンヒューリスティック最適化手法を用いる。具体的には、Particle Swarm Optimization（以下、ＰＳＯという）とその改良手法、遺伝的アルゴリズム（以下、ＧＡという）とその改良手法、タブーサーチ（以下，ＴＳという）とその改良手法、Ant Colony Optimization（以下、ＡＣＯという）とその改良手法などを用いる。
ここで、ＰＳＯ及びその改良手法とは、J. Kennedy and R. Eberhart, Swarm Intelligence, Morgan Kaufmann Publishers, 2001で解説されているように、Eberhart氏らにより開発された群れ理論を基にしたＰＳＯのGbestモデル、Lbestモデル、あるいは、Angeline氏により開発されたHybrid Particle Swarm OptimizationなどのＰＳＯの様々なバリエーションを含む手法を指す。
ＧＡ及びその改良手法とは、D. E. Goldberg, Genetic Algorithms in Search, Optimization, and Machine Learning, Addison-Wesley, 1989の中で述べられているSimple Genetic Algorithm（以下、ＳＧＡという）とその改良手法を意味し、ＴＳ及びその改良手法とは、F. Glover, "Tabu Search Part I", ORSA Journal of Computing, Vol. 1, NO. 3, Summer 1989で述べられているタブーサーチとその改良手法を意味し、ＡＣＯ及びその改良手法とは、A. Colorni, M. Dorigo, and V. Maniezzo, "Distributed Optimization by Ant Colonies", Proc. of First European Conference on Artificial Life, pp.134-142, Cambridge, MA: MIT Press 1991.で述べられているAnt Colony Optimization (ACO)とその改良手法を意味する。
【００９２】
以下では、例として、ＰＳＯ，ＳＧＡ及びＴＳを用いた場合のアルゴリズムを示す。
（１）ＰＳＯ
▲１▼状態表現
ＰＳＯでは、各制御時間毎の各プラント構成機器への燃料注入量、及び、起動・停止状態を変更できる各構成機器の起動・停止状態からなるデータ列によって状態を表現する。
▲２▼求解アルゴリズム
ＰＳＯを非線形混合整数計画問題に適用できるようにする方法は、福山他による「電圧信頼度を考慮したParticle Swarm Optimizationによる電圧無効電力制御方式の検討」電気学会論文誌Ｂ １１９巻１２号(１９９９年１２月)等で提案されている方法を用いる。
【００９３】
ＰＳＯは、簡単化した社会モデルのシミュレーションを通して開発されたモダンヒューリスティック手法の１つであり、鳥の群れの動きを連続変数の２次元空間で表現することを通して開発された。ＰＳＯでは、各エージェント（上述した鳥）の位置（状態量）をｘ，ｙ座標で表現し、前記位置（状態量）の変化分に相当する速度（ベクトル）をｖｘ（ｘ方向の速度）、ｖｙ（ｙ方向の速度）で表現する。
これらの位置と速度情報から、次時点の各エージェントの位置を更新することができる。この概念に基づき、鳥の群れ全体が何らかの目的関数を最適化するような行動をとると考えると、以下のような最適化が考えられる。
【００９４】
つまり、各エージェントは各々の探索における目的関数のそれまでの個々の最良値（ｐbest）と、その位置（状態量）を示すｘ，ｙ座標とを覚えている。また、各エージェントは、ｐbestのうち集団の中で最も最良のもの、すなわち集団のそれまでの目的関数の最良値（ｇbest）情報を共有している。
そして、各エージェントは、現在の自己のｘ，ｙ座標と速度ｖｘ，ｖｙ、及び、ｐbestとｇbestとの距離に応じて、ｐbest，ｇbestの存在する位置に方向を変更しようとする。この変更しようとする行動は速度を修正することで表現される。現在の速度とｐbest及びｇbestを用いて、各エージェントの速度は数式８により修正される。
【００９５】
【数８】
0004059014.tif 000005
【００９６】
数式８において、ｖｉ：エージェントｉの速度、rand()：０〜１までの一様乱数、ｓｉｋ：エージェントｉの探索ｋ回目の位置（探索点）、ｐbestｉ：エージェントｉのｐbest、ｗ：エージェントの速度に対する重み関数、cｉ：各項に対する重み係数である。
【００９７】
上記数式８を用いることにより、各エージェントのこれまでの最良解及び集団の最良解に確率的に近付くような速度が求められ、これにより各エージェントの現在の位置（探索点）を数式９により修正する。
ここで、本発明に即していえば、各エージェントがプラント構成機器に相当し、数式８，数式９における各エージェントの位置ｓｉが各プラント構成機器の起動・停止状態（離散値）や燃料注入量（連続値）に相当し、また、各エージェントの速度ｖｉがそれらの変化分に相当する。
【００９８】
【数９】
ｓｉｋ＋１＝ｓｉｋ＋ｖｉｋ＋１
【００９９】
ＰＳＯは遺伝的アルゴリズム（Genetic Algorithm：以下、ＧＡ）等と同様に複数の探索点を持った多点探索で、各探索点のｐbest及び集団のｇbestを用いて各探索点を確率的に変更していくことにより、大域最適解（最良解）を得る方法である。また、これまでの速度を維持しようとする大域探索（数式８の右辺第１項）と、ｐbest，ｇbestを用いてそれらに近付こうとする局所探索（それぞれ数式８の右辺第２，３項）とをバランスよく行う機構をもった探索手法である。
更に、ＰＳＯでは探索の各ステップで目的関数値を評価する必要があるが、評価の回数は問題の規模によらずエージェント数のみで良いというメリットがあるため、プラントの最適運用といった大規模問題への適用が可能である。
【０１００】
ｇbestモデルのＰＳＯを利用した場合のプラント最適運用方法の全体的なアルゴリズムを、以下に示す。なお、図１３はそのフローチャートである。
【０１０１】
Step.1  データの入力（図１３のＣ１）
▲１▼負荷予測値の入力
・各負荷種別毎の各制御時間毎の負荷予測値を入力する。この負荷予測値は、請求項１記載の発明によって得たプラントの電力負荷、熱負荷、空気負荷等の予測値である。
▲２▼プラント情報の入力（これらはデータベースとして予め格納されている）
・各プラント構成機器毎の入出力特性（線形、非線型）と各制御時間毎の運転パターン
・各構成機器相互の接続関係とそれに対応する特性式
・各燃料種別毎の燃料費
これらの情報をプラント情報として入力する。
▲３▼Particle Swarm Optimizationに関する情報の入力
・エージェント数、各最適化パラメータ値、最大探索回数を入力する。
【０１０２】
Step.2  初期値の生成（図１３のＣ２）
▲１▼各エージェント毎に、各制御時間毎の以下の値をランダムに計算する。
・各プラント構成機器毎の起動・停止状態（離散値）
・上記プラント構成機器のうち、燃料注入機器の燃料注入量（連続値）
▲２▼各エージェント毎のプラント状態の計算及び評価値の計算
・▲１▼で計算した各プラント構成機器毎の起動・停止状態及び燃料注入量を用いて、請求項２に係る定常プラントシミュレータにより、各制御時間毎の各プラント構成機器の入出力状態を求める。
・更に、図１２に示した目的関数及び制約条件を用いて、各エージェント毎に評価値を計算する。
ここで、評価値とは、プラント運用費用やＣＯ２排出量の総和である。
▲３▼ｐbest及びｇbestの初期設定
・▲２▼で計算した各エージェント毎の評価値を、現在の各エージェント毎のｐbest値とする。
・上記のｐbestのうちの最良値（最小値）をｇbestとする。
これにより、ｐbest及びｇbestの初期値が設定されることになる。
【０１０３】
Step.3  各エージェントの位置（探索点）の修正（図１３のＣ３）
▲１▼前述した数式８、数式９により、離散値及び連続値について各エージェントの位置（探索点）を修正する。
このときの修正の概念図を図１４に示す。図１４において、ｖｋ＋１は数式９の右辺第２項に相当している。
【０１０４】
Step.4  各エージェントの評価（図１３のＣ４）
▲１▼各エージェント毎のプラント状態の計算及び評価値の計算
・Step.2で計算した各プラント構成機器毎の起動・停止状態及び燃料注入量を用いて、定常プラントシミュレータにより、各制御時間毎の各プラント構成機器の入出力状態を求める。
・図１２の目的関数及び制約条件を用いて、各エージェント毎に評価値を計算する。
▲２▼ｐbest及びｇbestの修正
・▲１▼で計算した各エージェント毎の評価値が現在の各エージェント毎のｐbest値より良かったら、前記評価値を各エージェント毎のｐbestとする。
・上記のｐbestのうちの最良値が現在のｇbest値より良かったら、その最良値をｇbestとする。
これらの処理を通じて、各プラント構成機器の起動・停止状態や燃料注入量の最適解を求めていく。
【０１０５】
Step.5  終了条件のチェック（図１３のＣ５）
▲１▼探索回数が入力した最大探索回数に達したら終了する。そうでない場合は、Step.3へ戻る。
【０１０６】
（２）ＳＧＡ
▲１▼状態表現
ＳＧＡでは、全て離散値としなければならないことから、本来、連続値である各燃料注入量は離散化する最小燃料費刻みで離散化した値を利用し、離散化した数値の最小値から最大値までを整数と対応させる。内部で整数で表現し、連続値である実際の燃料注入量へは整数との対応表から変換する。
一方、もともと離散値である各制御時間毎の機器の起動・停止状態は、起動を１、停止を０として整数により表現する。つまり、燃料注入量及び起動・停止状態をすべて整数で表現して対応表を作成すると共に、燃料注入量を検討する機器（燃料注入機器）の数及び起動・停止状態を考慮する機器の数の合計の長さを持つ遺伝子で状態を表現する。従って、各遺伝子座は、各燃料注入量あるいは機器の起動・停止状態である。
【０１０７】
▲２▼求解アルゴリズム
Step.1  前提条件と初期条件設定（図１５のＥ１）
・請求項１の負荷予測方法で求めた負荷予測値を入力する。
・ストリング数、交差確率、突然変異確率、最大世代数を設定する。
【０１０８】
Step.2  初期値の生成（図１５のＥ２）
・各ストリングの遺伝子座について、ランダムに整数値を選択して初期値とする。
・現在の世代数を１とする。
【０１０９】
Step.3  各ストリングの評価と選択（図１５のＥ３）
・各ストリングの遺伝子座の整数値と対応表（前述した▲１▼の状態表現における対応表）を用いて、燃料注入量あるいは機器の起動・停止状態を決定する。
・プラントモデル及び負荷予測値を利用して各遺伝子に対する評価値を計算する。ここで、評価値とは、プラント運用費用やＣＯ２排出量の総和である。
・各ストリングの評価値を用いて、ルーレット・ホイール・セレクションによりストリングの選択を行う。
【０１１０】
Step.4  ストリング操作（図１５のＥ４）
・ストリング集合に対して、交差確率、突然変異確率を用いて、交差及び突然変異を実行することにより、上記評価値を最良とするようなストリングを決定する。そして、このストリングの遺伝子座の整数値と対応表を用いて、燃料注入量や機器の起動・停止状態の最適解を決定していく。
【０１１１】
Step.5  終了判定（図１５のＥ５，Ｅ６）
・世代が事前に決定した最大世代数の設定値に達したら、終了する。
・設定値に達しない場合には、世代数に１を足してStep.3へ戻る。
【０１１２】
(２)ＴＳ
▲１▼状態表現
ＴＳにおいても、全て離散値としなければならないことから、本来、連続量である各燃料注入量は離散化する最小燃料費刻みで離散化した値を利用し、離散化した数値の最小値から最大値までを整数と対応させる。内部で整数で表現し、連続量である実際の燃料注入量へは整数との対応表から変換する。
また、もともと離散値である各制御時間毎の機器の起動停止状態は、起動を１、停止を０として整数により表現する。つまり、ＳＧＡと同様に燃料注入量及び起動・停止状態をすべて整数で表現して対応表を作成すると共に、燃料注入機器の数及び起動・停止状態を考慮する機器の数の合計の長さを持つ整数配列によって表現する。従って、配列の各要素は、各燃料注入量あるいは機器の起動・停止状態である。
【０１１３】
▲２▼求解アルゴリズム
Step.1  前提条件と初期条件設定（図１６のＦ１）
・請求項１の負荷予測方法により求めた負荷予測値を入力する。
・タブー長及び最大探索回数を設定する。
【０１１４】
Step.2  初期値の生成（図１６のＦ２）
・状態表現の配列の各要素について、ランダムに整数値を選択して初期値とする。
・現在の状態をタブーリストに入れる。
・現在の探索回数を１とする。
【０１１５】
Step.3  隣接状態の生成と次状態の決定（図１６のＦ３）
・現在の状態の配列の各要素に対して、整数値の＋１及び−１の値（上下限値の場合は、値を生成しない）を隣接状態として生成する。
例えば配列されている整数値が（２，３，４）であった場合、各要素について＋１及び−１を加算した（１，３，４），（３，３，４），（２，２，４），（２，４，４），（２，３，３），（２，３，５）が隣接状態となる。
・各隣接状態の整数値の配列に対して、対応表を参照することにより燃料注入量あるいはプラント構成機器の起動・停止状態を決定する。
・プラントモデル及び負荷予測値を利用して、各遺伝子に対する評価値を計算する。ここで、評価値とは、プラント運用費用やＣＯ２排出量の総和である。
・各隣接状態の中で、タブーでない最も評価がよいものを次状態とする処理により、上記評価値を最良とするような隣接状態を決定する。そして、この隣接状態の整数値と対応表を用いて、燃料注入量や機器の起動・停止状態の最適解を決定していく。
【０１１６】
Step.4  終了判定（図１６のＦ４，Ｆ５）
・探索回数が事前に設定した最大探索回数に達したら、終了する。
・設定値に達しない場合には、現在の探索回数に１を足してStep.3へ戻る。
【０１１７】
なお、上記の最適運用方法を利用して当日の運用計画を補正（再計画）することが望ましい。
すなわち、前述したように、所定の状態変数（▲１▼各制御時間毎の各プラント構成機器の起動・停止状態、▲２▼各制御時間毎の燃料注入量）を用い、目的関数（ｆ＝ｗ１ｆ１＋ｗ２ｆ２＋ｗ３ｆ３）を最小化するように各種制約条件（▲１▼各負荷種別毎の需給バランス、▲２▼各機器の特性上の制約）を考慮しながら最適化手法によりプラント構成機器の起動・停止状態、燃料注入量を決定した結果、図１７（ａ）に示すような当日の運用計画が前日に得られたとする（前日に計画されたという意味で前日計画と称する）。
【０１１８】
この前日計画において考慮された制約条件のうち各負荷種別毎の需給バランスは、前日における負荷予測値を用いたものであり、当日の負荷値は予測値からずれて推移している場合が多い。このため、当日には、新たな負荷予測値に基づく需給バランスを制約条件として考慮する必要があり、言い換えれば、前日計画を補正して当日補正計画を作成する必要がある。
【０１１９】
このため、図１７（ｂ）に示すように現時点が当日のｔ時（ｔ時間経過後）であるとすると、当日の残りの時間（２４−ｔ時間）を考察時間（補正計画の対象時間）とし、現時点におけるプラント構成機器の運転状態や既に運転指令済みで変更不可能な機器の運転状態をプラントの固定的な初期状態とし、この状態を各機器の特性上の制約として考慮した上で、当日の残りの時間（２４−ｔ時間）の負荷予測値に基づく需給バランスを制約条件として考慮しつつ、最適運用方法により当日補正計画を作成すればよい。
以上の着想が、請求項３に記載した発明に相当する。
【０１２０】
４．次いで、請求項４に記載したプラントの最適設計方法の実施形態を説明する。
この最適設計方法とは、プラントの設計段階において、各プラント構成機器の最適容量の決定問題を、設定可能な機器容量値（離散値）の中から最適値を決定する組合せ最適化問題として定式化するとともに、年間の各月の代表日の負荷パターンを実績値等を用いて経験的に設定し、請求項１または２に記載した最適運用方法を利用して、年間のプラント運用費用やＣＯ２排出量に関する評価値を最良化するような各プラント構成機器の容量を求めることをいう。
この最適設計方法は、以下の手順により実現する。
【０１２１】
Step.1  前提条件及び初期条件設定（図１８のＤ１）
・プラントのシステム構成(プラント構成機器と相互の接続状態)を決定する。
・所定期間、例えば年間の各月の代表日の各制御時間毎の負荷値を負荷パターンとして設定する。対象プラントが既存プラントであれば、これまでの実績値を用いて負荷パターンを設定する。また、新規プラントであれば、既存プラントの実績負荷値などを参考にして設定する。
・各プラント構成機器の容量の初期値を計画値として設定する。
【０１２２】
Step.2  Step.1で設定した値を用いて、請求項３または４のプラント最適運用方法により、各代表日に対する各機器の起動・停止状態及び燃料注入量を決定して最適運用方法を求める（図１８のＤ２）。
【０１２３】
Step.3  １年を通したプラント運用費用やＣＯ２排出量に関する評価値（各代表日の評価値の和）が最良となるまで、各プラントの構成機器の容量を組合せ最適化手法を用いて逐次変更する（図１８のＤ３）。
【０１２４】
ここで、組合せ最適化手法としては、遺伝的アルゴリズムとその改良手法、及びタブーサーチとその改良手法、Ant Colony Optimizationとその改良手法などを利用する。
上記遺伝的アルゴリズムとは、E. Goldberg, Genetic Algorithms in Search, Optimization, and Machine Learning, Addison-Wesley, 1989の中で述べられているSimple Genetic Algorithm（以下、ＳＧＡという）とその改良手法、タブーサーチとはF. Glover, "Tabu Search Part I", ORSA Journal of Computing, Vol. 1, NO. 3, Summer 1989で述べられているタブーサーチ（以下、ＴＳという）とその改良手法、Ant Colony OptimizationとはA. Colorni, M. Dorigo, and V. Maniezzo, "Distributed Optimization by Ant Colonies", Proc. of First European Conference on Artificial Life, pp.134-142, Cambridge, MA: MIT Press 1991.で述べられているAnt Colony Optimization （以下、ＡＣＯという）及びその改良手法を意味する。
以下では、例としてＳＧＡ及びＴＳを用いた場合のアルゴリズムを前述の図１５，図１６を参照しつつ説明する。
【０１２５】
(１)ＳＧＡ
▲１▼状態表現
ＳＧＡでは、最適容量を求めるプラント構成機器の数だけの遺伝子長とし、各遺伝子座は各構成機器の選択可能な容量を小さい方から整数で表現した整数値とする。従って、実際の容量は、各整数値と実際の容量との対応表から求める。例えば、１が１００ｋＶＡ、２が１５０ｋＶＡ、３が２００ｋＶＡなどの対応となる。
【０１２６】
▲２▼求解アルゴリズム
Step.1  前提条件及び初期条件設定（図１５のＥ１）
・プラントのシステム構成(プラントの構成機器と相互の接続状態)を決定する。
・年間の各月の代表日の各制御時間毎の負荷値を設定する。既存プラントであれば、これまでの実績値を用いて設定する。また、新規プラントであれば、既存プラントの実績負荷値などを参考にして設定する。
・各プラント構成機器の容量の初期値を計画値として設定する。
・ストリング数、交差確率、突然変異確率、最大世代数を設定する。
【０１２７】
Step.2  初期値の生成（図１５のＥ２）
・各ストリングの遺伝子座について、ランダムに整数値を選択して初期値とする。
・現在の世代数を１とする。
【０１２８】
Step.3  各ストリングの評価と選択（図１５のＥ３）
・各ストリングの遺伝子座の整数値を対応表を用いて設定し、容量を決定する。
・各構成機器容量が決定したプラントモデル、及び、最初に設定した各月の代表日の各制御時間毎の負荷値を利用して各機器の起動・停止状態及び燃料注入量を決定することにより最適運用方法を求め、１年を通した各月の代表日に対するプラント運用費用やＣＯ２排出量に関する評価値の和を各ストリングの評価値とする。
・各ストリングの評価値を用いて、ルーレット・ホイール・セレクションによりストリングの選択を行う。
【０１２９】
Step.4  ストリング操作（図１５のＥ４）
・ストリング集合に対して、交差確率、突然変異確率を用いて、交差及び突然変異を実行することにより、１年を通した評価値の和が最良となるような各構成機器の容量を決定する。
【０１３０】
Step.5  終了判定（図１５のＥ５，Ｅ６）
・世代数が事前に決定した最大世代数の設定値に達したら、終了する。
・設定値に達しない場合には、世代数に１を足してStep.3へ戻る。
【０１３１】
(２)ＴＳ
▲１▼状態表現
ＴＳでは、最適容量を求めるプラント構成機器の数だけの配列長とし、配列の各要素は各構成機器の選択可能な容量を小さい方から整数で表現した整数値とする。従って、実際の容量は、各整数値と容量との対応表から求める。例えば、１が１００ｋＶＡ、２が１５０ｋＶＡ、３が２００ｋＶＡなどの対応となる。
【０１３２】
▲２▼求解アルゴリズム
Step.1  前提条件及び初期条件設定（図１６のＦ１）
・プラントのシステム構成(プラントの構成機器と相互の接続状態)を決定する。
・年間の各月の代表日の各制御時間毎の負荷値を設定する。既存プラントであれば、これまでの実績値を用いて設定する。また、新規プラントであれば、既存プラントの実績負荷などを参考にして設定する。
・各プラント構成機器の容量の初期値を計画値として設定する。
・タブー長及び最大探索回数を設定する。
【０１３３】
Step.2  初期値の生成（図１６のＦ２）
・状態表現の配列の各要素について、ランダムに整数値を選択して初期値とする。
・現在の状態をタブーリストに入れる。
・現在の探索回数を１とする。
【０１３４】
Step.3 隣接状態の生成と次状態の決定（図１６のＦ３）
・現在の状態の配列の各要素に対して、整数値の＋１及び−１の値（上下限値の場合は、値を生成しない）を隣接状態として生成する。例えば配列されている整数値が（２，３，４）であった場合、各要素について＋１及び−１を加算した（１，３，４），（３，３，４），（２，２，４），（２，４，４），（２，３，３），（２，３，５）が隣接状態となる。
・各隣接状態の整数値の配列に対して、対応表を参照することにより各構成機器の容量を決定する。
・各構成機器容量が決定したプラントモデル、及び、最初に設定した各月の代表日の各制御時間毎の負荷値を利用して各機器の起動・停止状態及び燃料注入量を決定することにより最適運用方法を求め、１年を通した各月の代表日に対するプラント運用費用やＣＯ２排出量に関する評価値の和を各隣接状態の評価値とする。
・各隣接状態の中でタブーでない最も評価がよいものを次状態とする処理により、１年を通した評価値の和が最良となるような各構成機器の最適容量を決定する。
【０１３５】
Step.4  終了判定（図１６のＦ４，Ｆ５）
・探索回数が事前に設定した最大探索回数に達したら、終了する。
・設定値に達しない場合には、現在の探索回数に１を足してStep.3へ戻る。
【０１３６】
【発明の効果】
以上述べたように、請求項１記載の発明によれば、非線形性を有するプラントの各種負荷を高精度に予測可能として、プラント構成機器の効率的な運転、省エネルギー化が可能になると共に、予測理由の説明を容易にして運用者の不安感を払拭することができる。
【０１３７】
また、請求項１に記載した発明によれば、プラント構成機器の非線形な入出力特性や複雑な運転パターンを考慮したうえで、プラント全体の挙動を正確に模擬する定常プラントシミュレータを実現することができる。
【０１３８】
更に、請求項１，２に記載した発明によれば、プラントの最適運用を非線形混合整数計画問題として定式化し、前記定常プラントシミュレータの出力を用いることで、大域最適解に近いプラント構成機器の起動・停止状態及び燃料注入量を決定することができる。
請求項３記載の発明によれば、負荷予測値と当日の実際の負荷値とのずれを考慮して運用計画を補正することが可能であり、プラントの運用を一層最適化することが可能である。
【０１３９】
請求項４に記載した発明によれば、前記最適運用方法を利用して一定期間内のプラント運用費用や排出ガスに関する評価値を最良化するプラント構成機器の最適容量を、組み合わせ最適化問題を解くことにより容易に求めることができる。
【図面の簡単な説明】
【図１】請求項１の発明の実施形態が適用される階層型ニューラルネットワークの一般的な構造を示す図である。
【図２】ニューラルネットワークの学習方法を示すフローチャートである。
【図３】ニューラルネットワークの学習方法を示すフローチャートである。
【図４】図３における第１ステップのネットワーク構造を示す図である。
【図５】図３における第１ネットワークの説明図である。
【図６】シグモイド関数の入出力関係を示す図である。
【図７】請求項１の発明の実施形態が適用される階層型ニューラルネットワークの具体的な構造を示す図である。
【図８】入力因子の出力への影響を解析する方法を説明するためのニューラルネットワークの構造を示す図である。
【図９】入力値と出力層へ伝達する情報との相関関係を示す図である。
【図１０】エネルギープラントの一例を示す構成図である。
【図１１】応答曲面法における中央複合計画の実験点を示す概念図である。
【図１２】プラントシミュレータを用いてプラントの最適運用方法を決定する際の概念図である。
【図１３】請求項１に記載した発明の実施形態を示すフローチャートである。
【図１４】ＰＳＯにおける修正の概念図である。
【図１５】ＳＧＡにおける求解アルゴリズムを示すフローチャートである。
【図１６】ＴＳにおける求解アルゴリズムを示すフローチャートである。
【図１７】運用計画の補正を示す概念図である。
【図１８】請求項４に記載した発明の実施形態を示すフローチャートである。
【符号の説明】
１１，１１Ａ 全結合部分
１２，１２Ａ，１２Ｂ，１２Ｃ 疎結合部分
１２Ｄ 第１ネットワーク
２１ ガスタービン
２２ 排ガスボイラ
２３ スチームタービン
２４ コンプレッサ
２５ ボイラ
２６ 吸収式冷凍機
２７ ガス焚冷凍機
２８ 蓄熱層
３０ 最適化部
４０ 定常プラントシミュレータ
【図１】
0004059014.tif 000006
【図２】
0004059014.tif 000007
【図３】
0004059014.tif 000008
【図４】
0004059014.tif 000009
【図５】
0004059014.tif 000010
【図６】
0004059014.tif 000011
【図７】
0004059014.tif 000012
【図８】
0004059014.tif 000013
【図９】
0004059014.tif 000014
【図１０】
0004059014.tif 000015
【図１１】
0004059014.tif 000016
【図１２】
0004059014.tif 000017
【図１３】
0004059014.tif 000018
【図１４】
0004059014.tif 000019
【図１５】
0004059014.tif 000020
【図１６】
0004059014.tif 000021
【図１７】
0004059014.tif 000022
【図１８】
0004059014.tif 000023