(57)【要約】【課題】少ない計算量で、係数行列の変化に対する固有値の変化を精度良く推定することができるようにする。【解決手段】計算された係数行列Ａの固有ベクトルと、算出された行列ΔＡとから得られる行列ΔΛ（＝Ｐ−１ΔＡＰ）から、行列ΔΛの近似行列ΔΛ’を生成する（１７０）。生成された近似行列ΔΛ’と行列Λ（＝Ｐ−１ＡＰ）との和の行列（Λ＋ΔΛ’）について、非零の非対角要素を左上に集めるように行と列とを並び替える（１７２）。並び替えられた行列のうちの、非零の非対角要素が集められた部分からなる部分行列の固有値及び固有ベクトルを算出し、算出した部分行列の固有値及び固有ベクトルに基づいて、行列（Λ＋ΔΛ’）を並び替えた行列の固有値及び固有ベクトルを算出する。そして、元の行列（Λ＋ΔΛ’）の固有値及び固有ベクトルを、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルの推定値とする（１７４）。【選択図】図９
2010074995.tif 000002
【特許請求の範囲】
【請求項１】
状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、
前記係数行列Ａに含まれるパラメータを変化させたときの係数行列Ａの変化分を算出する変化分算出手段と、
前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段と、
を含む固有値推定装置。
【数１】
2010074995.tif 000037
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【請求項２】
前記固有値推定手段は、前記行列ΔΛの近似行列と前記固有値行列Δとの和の行列に対して、０でない非対角要素を集めるように行及び列の並び替えを行い、０でない非対角要素が集められた部分行列の固有値及び固有ベクトルを算出し、前記算出した固有値及び固有ベクトルに基づいて、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する請求項１記載の固有値推定装置。
【請求項３】
該着目モードと異なる対象モードの各々について、前記固有値行列Λと、前記算出された前記行列ΔΛとの和の行列における、該着目モード及び前記対象モードの各々に対応する対角要素、並びに該着目モード及び前記対象モードに対応する非対角要素で構成される行列の固有値を算出し、前記着目モードに対応する対角要素と、前記算出された行列の固有値のうちの前記着目モードに対応する固有値との差の絶対値を算出する差算出手段と、
前記差算出手段によって算出された前記差の絶対値が大きくなる順に、前記対象モードを前記強相関モードとして所定個選択するモード選択手段とを更に含む請求項１又は２記載の固有値推定装置。
【請求項４】
該着目モードと異なる対象モードの各々について、前記算出された前記行列ΔΛにおける、該着目モード及び前記対象モードに対応する非対角要素の積の絶対値を算出する積算出手段と、
前記積算出手段によって算出された前記積の絶対値が大きくなる順に、前記対象モードを前記強相関モードとして所定個選択するモード選択手段とを更に含む請求項１又は２記載の固有値推定装置。
【請求項５】
複数種類の選択個数の各々について、前記差算出手段によって算出された前記差の絶対値が大きくなる順に、前記対象モードを前記強相関モードとして前記選択個数だけ選択したときに、選択されない前記対象モードに対する前記差の絶対値の総和を、誤差として推定する誤差推定手段と、
前記誤差推定手段によって推定された前記誤差が所定値以下となるように、選択個数を決定する個数決定手段とを更に含み、
前記モード選択手段は、前記対象モードを前記強相関モードとして前記決定された選択個数だけ選択する請求項３記載の固有値推定装置。
【請求項６】
電力系統の安定化を図るための安定化制御を行う電力系統安定化装置に設定する制御パラメータと、入力される前記電力系統の系統情報とに基づいて、前記電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、
前記制御パラメータを変化させたときの前記係数行列Ａの変化分を算出する変化分算出手段と、
前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段と、
前記固定値推定手段によって推定された固有値に基づいて、前記電力系統の安定度を評価する安定度評価手段と、
前記安定度評価手段によって評価された前記安定度の評価値を最大化又は最小化するように、前記制御パラメータを変更するパラメータ変更手段と、
を含む電力系統制御装置。
【数２】
2010074995.tif 000038
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【請求項７】
電力系統を運転操作するために設定される操作パラメータと、入力される前記電力系統の系統情報と、予め設定された想定事故とに基づいて、前記電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、
前記操作パラメータを変化させたときの前記係数行列Ａの変化分を算出する変化分算出手段と、
前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段と、
前記固定値推定手段によって推定された固有値に基づいて、前記電力系統の定態安定度余裕を算出する安定度余裕算出手段と、
前記安定度余裕算出手段によって算出された前記定態安定度余裕及び前記操作パラメータをパラメータとする予め定められた制約条件を満たす範囲内で、前記電力系統が最適な状態となるように、前記操作パラメータを変更するパラメータ変更手段と、
を含む電力系統予防制御装置。
【数３】
2010074995.tif 000039
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【請求項８】
プラントを運転操作するために設定される操作パラメータと、入力される前記プラントのプラント情報とに基づいて、前記プラントのプラント状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、
前記操作パラメータを変化させたときの行列Ａの変化分を算出する変化分算出手段と、
前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段と、
前記固定値推定手段によって推定された固有値に基づいて、前記プラントが望ましい運転状態であるか評価するための予め定められた評価関数の値を算出する評価関数算出手段と、
前記評価関数算出手段によって算出された前記評価関数の値が最適な値になるように、前記操作パラメータを変更するパラメータ変更手段と、
を含むプラント制御装置。
【数４】
2010074995.tif 000040
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【請求項９】
コンピュータを、
状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段、
前記係数行列Ａに含まれるパラメータを変化させたときの係数行列Ａの変化分を算出する変化分算出手段、
前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段
として機能させるためのプログラム。
【数５】
2010074995.tif 000041
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、固有値推定装置、電力系統制御装置、電力系統予防制御装置、プラント制御装置、及びプログラムに関する。
【背景技術】
【０００２】
電力供給コストの低減を進めるためには、送変電線設備を最大限に利用することが必要であるが、電力系統の安定限界をより正確に把握しておくことが極めて重要となっている。
【０００３】
電力系統の定態安定度は、線形状態微分方程式を表す係数行列Ａの固有値を求めることにより定量的に評価できる。
【０００４】
また、係数行列の変化に対する固有値の変化を推定する方法として固有値感度解析法を用いた、大規模電力系統の定態安定度解析手法が知られている（非特許文献１）。また、固有値感度解析法を用いた定態安定度の監視・予防制御方式が知られている（非特許文献２）。
【非特許文献１】内田、長尾、「電力系統の定態安定度向上対策（その１）  ‐ＰＳＳの設置箇所選定と定数最適化論理‐」、電力中央研究所報告、１９８３、１８３０４０
【非特許文献２】内田、高崎、吉村、「電力系統の安定運用のための予防制御論理（その２）  ‐固有値感度による予防制御の基本論理‐」、電力中央研究所報告、１９８６、Ｔ８６０３９
【発明の開示】
【発明が解決しようとする課題】
【０００５】
しかしながら、上記の非特許文献１、２に記載の技術で用いられる固有値感度解析法は線形推定であるのに対し、固有値の変化が係数行列の変化に対し強い非線形性を示すことがあるため、パラメータの変化に対しその都度、固有値を再計算する必要があり、多大な量の計算が必要となってしまう、という問題がある。
【０００６】
本発明は、上記の問題点を解決するためになされたもので、少ない計算量で、係数行列の変化に対する固有値の変化を精度良く推定することができる固有値推定装置、電力系統制御装置、電力系統予防制御装置、プラント制御装置、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記の目的を達成するために第１の発明に係る固有値推定装置は、状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、前記係数行列Ａに含まれるパラメータを変化させたときの係数行列Ａの変化分を算出する変化分算出手段と、前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段とを含んで構成されている。
【０００８】
【数１】
2010074995.tif 000003
【０００９】
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【００１０】
第２の発明に係るプログラムは、コンピュータを、状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段、前記係数行列Ａに含まれるパラメータを変化させたときの係数行列Ａの変化分を算出する変化分算出手段、前記変化分算出手段によって算出された前記係数行列Ａの変化分、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分だけ変化したときの係数行列Ａの固有値及び固有ベクトルを推定する固有値推定手段として機能させるためのプログラムである。
【００１１】
【数２】
2010074995.tif 000004
【００１２】
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【００１３】
第１の発明及び第２の発明によれば、固有値算出手段によって、状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出し、変化分算出手段によって、係数行列Ａに含まれるパラメータを変化させたときの係数行列Ａの変化分ΔＡを算出する。
【００１４】
そして、固有値推定手段によって、変化分算出手段によって算出された係数行列Ａの変化分ΔＡ、及び固有値算出手段によって算出された係数行列Ａの固有ベクトルに基づいて、行列ΔΛを算出し、行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした行列ΔΛの近似行列と、係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。
【００１５】
このように、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができる。
【００１６】
第１の発明に係る固有値推定手段は、行列ΔΛの近似行列と固有値行列Δとの和の行列に対して、０でない非対角要素を集めるように行及び列の並び替えを行い、０でない非対角要素が集められた部分行列の固有値及び固有ベクトルを算出し、算出した固有値及び固有ベクトルに基づいて、変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定することができる。これによって、０でない非対角要素が集められた部分行列の固有値及び固有ベクトルを算出することにより、計算量を減らすことができる。
【００１７】
上記の固有値推定装置は、着目モードと異なる対象モードの各々について、固有値行列Λと、算出された行列ΔΛとの和の行列における、着目モード及び対象モードの各々に対応する対角要素、並びに着目モード及び対象モードに対応する非対角要素で構成される行列の固有値を算出し、着目モードに対応する対角要素と、算出された行列の固有値のうちの着目モードに対応する固有値との差の絶対値を算出する差算出手段と、差算出手段によって算出された差の絶対値が大きくなる順に、対象モードを強相関モードとして所定個選択するモード選択手段とを更に含むことができる。これによって、着目モードと強い相互作用を持つモードを強相関モードとして選択することができる。
【００１８】
上記の固有値推定装置は、着目モードと異なる対象モードの各々について、算出された行列ΔΛにおける、着目モード及び対象モードに対応する非対角要素の積の絶対値を算出する積算出手段と、積算出手段によって算出された積の絶対値が大きくなる順に、対象モードを強相関モードとして所定個選択するモード選択手段とを更に含むことができる。これによって、着目モードと強い相互作用を持つモードを強相関モードとして選択することができる。
【００１９】
上記の差算出手段を含む固有値推定装置は、複数種類の選択個数の各々について、差算出手段によって算出された差の絶対値が大きくなる順に、対象モードを強相関モードとして選択個数だけ選択したときに、選択されない対象モードに対する差の絶対値の総和を、誤差として推定する誤差推定手段と、誤差推定手段によって推定された誤差が所定値以下となるように、選択個数を決定する個数決定手段とを更に含み、モード選択手段は、対象モードを強相関モードとして、決定された選択個数だけ選択することができる。これによって、誤差が所定値以下となるように、強相関モードを選択することができる。
【００２０】
第３の発明に係る電力系統制御装置は、電力系統の安定化を図るための安定化制御を行う電力系統安定化装置に設定する制御パラメータと、入力される前記電力系統の系統情報とに基づいて、前記電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、前記制御パラメータを変化させたときの前記係数行列Ａの変化分ΔＡを算出する変化分算出手段と、前記変化分算出手段によって算出された前記係数行列Ａの変化分ΔＡ、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する固有値推定手段と、前記固定値推定手段によって推定された固有値に基づいて、前記電力系統の安定度を評価する安定度評価手段と、前記安定度評価手段によって評価された前記安定度の評価値を最大化又は最小化するように、前記制御パラメータを変更するパラメータ変更手段とを含んで構成されている。
【００２１】
【数３】
2010074995.tif 000005
【００２２】
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【００２３】
第３の発明に係る電力系統制御装置によれば、固有値算出手段によって、電力系統の安定化を図るための安定化制御を行う電力系統安定化装置に設定する制御パラメータと、入力される電力系統の系統情報とに基づいて、電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する。変化分算出手段によって、制御パラメータを変化させたときの係数行列Ａの変化分ΔＡを算出する。
【００２４】
そして、固有値推定手段によって、変化分算出手段によって算出された係数行列Ａの変化分ΔＡ、及び固有値算出手段によって算出された係数行列Ａの固有ベクトルに基づいて、行列ΔΛを算出し、行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした行列ΔΛの近似行列と、係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。
【００２５】
そして、安定度評価手段によって、固定値推定手段によって推定された固有値に基づいて、電力系統の安定度を評価し、パラメータ変更手段によって、安定度評価手段によって評価された安定度の評価値を最大化又は最小化するように、制御パラメータを変更する。
【００２６】
このように、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができ、電力系統の安定度を精度良く評価することができる。
【００２７】
第４の発明に係る電力系統予防制御装置は、電力系統を運転操作するために設定される操作パラメータと、入力される前記電力系統の系統情報と、予め設定された想定事故とに基づいて、前記電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、前記操作パラメータを変化させたときの前記係数行列Ａの変化分ΔＡを算出する変化分算出手段と、前記変化分算出手段によって算出された前記係数行列Ａの変化分ΔＡ、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する固有値推定手段と、前記固定値推定手段によって推定された固有値に基づいて、前記電力系統の定態安定度余裕を算出する安定度余裕算出手段と、前記安定度余裕算出手段によって算出された前記定態安定度余裕及び過負荷に関する予め定められた制約条件を満たす範囲内で、前記電力系統が最適な状態となるように、前記操作パラメータを変更するパラメータ変更手段とを含んで構成されている。
【００２８】
【数４】
2010074995.tif 000006
【００２９】
ただし、ΔＡは、係数行列Ａの変化分であり、Ｐは、係数行列Ａの固有ベクトルを並べた行列である。
【００３０】
第４の発明に係る電力系統予防制御装置によれば、固有値算出手段によって、電力系統を運転操作するために設定される操作パラメータと、入力される電力系統の系統情報と、予め設定された想定事故とに基づいて、電力系統の系統状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する。変化分算出手段によって、操作パラメータを変化させたときの係数行列Ａの変化分ΔＡを算出する。
【００３１】
そして、固有値推定手段によって、変化分算出手段によって算出された係数行列Ａの変化分ΔＡ、及び固有値算出手段によって算出された係数行列Ａの固有ベクトルに基づいて、行列ΔΛを算出し、行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした行列ΔΛの近似行列と、係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。
【００３２】
そして、安定度余裕算出手段によって、固定値推定手段によって推定された固有値に基づいて、電力系統の定態安定度余裕を算出し、パラメータ変更手段によって、安定度余裕算出手段によって算出された定態安定度余裕及び過負荷に関する予め定められた制約条件を満たす範囲内で、電力系統が最適な状態となるように、操作パラメータを変更する。
【００３３】
このように、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができ、電力系統の定態安定度余裕を精度良く算出することができる。
【００３４】
第５の発明に係るプラント制御装置は、プラントを運転操作するために設定される操作パラメータと、入力される前記プラントのプラント情報とに基づいて、前記プラントのプラント状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する固有値算出手段と、前記操作パラメータを変化させたときの行列Ａの変化分ΔＡを算出する変化分算出手段と、前記変化分算出手段によって算出された前記係数行列Ａの変化分ΔＡ、及び前記固有値算出手段によって算出された前記係数行列Ａの固有ベクトルに基づいて、以下の式で表される行列ΔΛを算出し、前記行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした前記行列ΔΛの近似行列と、以下の式で表される前記係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、前記変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する固有値推定手段と、前記固定値推定手段によって推定された固有値に基づいて、前記プラントが望ましい運転状態であるか評価するための予め定められた評価関数の値を算出する評価関数算出手段と、前記評価関数算出手段によって算出された前記評価関数の値が最適な値になるように、前記操作パラメータを変更するパラメータ変更手段とを含んで構成されている。
【００３５】
【数５】
2010074995.tif 000007
【００３６】
ただし、ΔＡは、前記係数行列Ａの変化分であり、Ｐは、前記係数行列Ａの固有ベクトルを並べた行列である。
【００３７】
第５の発明に係るプラント制御装置によれば、固有値算出手段によって、プラントを運転操作するために設定される操作パラメータと、入力されるプラントのプラント情報とに基づいて、プラントのプラント状態を示す状態変数ベクトルを変数とした線形微分方程式の係数行列Ａの固有値及び固有ベクトルを算出する。変化分算出手段によって、操作パラメータを変化させたときの行列Ａの変化分ΔＡを算出する。
【００３８】
そして、固有値推定手段によって、変化分算出手段によって算出された係数行列Ａの変化分ΔＡ、及び固有値算出手段によって算出された係数行列Ａの固有ベクトルに基づいて、行列ΔΛを算出し、行列ΔΛのうち、予め設定された着目モード及び他のモードより該着目モードに強い相互作用を持つ予め選択された強相関モードに対応する非対角要素以外の非対角要素を０とした行列ΔΛの近似行列と、係数行列Ａの固有値行列Λとの和の行列について、固有値及び固有ベクトルを算出することにより、変化分ΔＡだけ変化したときの係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。
【００３９】
そして、評価関数算出手段によって、固定値推定手段によって推定された固有値に基づいて、プラントが望ましい運転状態であるか評価するための予め定められた評価関数の値を算出し、パラメータ変更手段によって、評価関数算出手段によって算出された評価関数の値が最適な値になるように、操作パラメータを変更する。
【００４０】
このように、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができ、プラントの運転状態を精度良く評価することができる。
【発明の効果】
【００４１】
以上説明したように、本発明の固有値推定装置、電力系統制御装置、電力系統予防制御装置、プラント制御装置、及びプログラムによれば、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができる、という効果が得られる。
【発明を実施するための最良の形態】
【００４２】
以下、図面を参照して本発明の実施の形態を詳細に説明する。なお、第１の実施の形態では、電力系統安定化装置によって電力系統の安定化を図る電力系統システムに本発明を適用した場合を例に説明する。
【００４３】
図１に示すように、第１の実施の形態に係る電力系統システム１０は、多数の発電機を含む電力系統１２と、電力系統１２の安定化を図るために励磁を調整する電力系統安定化装置（ＰＳＳ）１４と、入力される電力系統１２の系統構成や運転状態を示す系統情報、及び発電機の潮流計算結果を読み込んで、電力系統安定化装置１４の制御パラメータとしてのＰＳＳ定数を最適化して、電力系統安定化装置１４に設定するコンピュータ１６とを備えている。
【００４４】
コンピュータ１６は、ＣＰＵ、後述する制御パラメータ最適化処理ルーチンを実現するためのプログラムを記憶したＲＯＭ、データを一時的に記憶するＲＡＭ、及びＨＤＤ等の記憶装置で構成されている。
【００４５】
ここで、電力系統１２の安定度について説明する。まず、図２に示すように、無限大母線（電圧の位相と大きさが不変）に同期機が接続されている場合の系統を一機無限大母線系統と称する。この一機無限大母線系統は、安定度を論じる場合の最も基本となる系統である。タービン発電機のような非凸極機の場合は、直軸リアクタンスxdと横軸リアクタンスxeは等しく、同期リアクタンスxdの背後電圧が一定に保たれる。上記図２では、無限大母線の電圧E∞を基準とし、発電機の内部電圧の大きさをEq、内部相差角をδ、端子電圧の大きさをVt、端子相差角をδtとしている。
【００４６】
また、電力系統の安定度は、発電機のローターに磁束を介して加わる電磁トルクの作用によるローターの回転速度や位相の変動としてとらえることができる。従来から考えられてきた電力系統の安定度は、運動方程式が支配的な主として有効電力とローターの位相角の動揺として現れる、いわゆる電力モードに主眼をおいている。しかし、磁束変動の微分方程式が支配的な、電圧の動揺として大きく現れるモード、いわゆる電圧モードが不安定の要因となっている。数学的にはこの二つのモード以外にも多種の振動モードが存在するが、一機無限大系統ではこの両モード以外は系統の安定性にはほとんど影響がない。
【００４７】
次に、電力系統１２の基本方程式について説明する。
【００４８】
電力系統の基本方程式は、一般に以下の（１）式に示す非線形微分方程式で表すことができる。
【００４９】
【数６】
2010074995.tif 000008
【００５０】
ここで、ｘは、電力系統の系統状態を示す状態変数ベクトルであり、αは任意のパラメータである。電力系統の状態変数としては、例えば、発電機のローターの回転速度、位相、磁束などが該当する。
【００５１】
定態安定度解析は、上記（１）式の平衡点ｘ０（ｆ（ｘ，α）＝０となる点ｘ＝ｘ０）のまわりでの微少変化に対する電力系統の安定性を論じるものである。すなわち、以下の（２）式で微小変化量Δｘに変数の変換を行い、上記（１）式を（３）式で線形近似し、この線形微分方程式の安定性を問題とする。
【００５２】
【数７】
2010074995.tif 000009
【００５３】
ここで係数行列Ａは、以下の（４）式で表されるｎ×ｎの行列である。
【００５４】
【数８】
2010074995.tif 000010
【００５５】
本実施の形態では、いずれも状態変数の微少変動量に対する線形微分方程式を基本方程式とするので、煩雑さを避けるため、Δｘをあらためてｘと書き直し、以下の（５）式を基本方程式とする。
【００５６】
【数９】
2010074995.tif 000011
【００５７】
上記（５）式は、定常運転状態で線形化した微分方程式を表わし、状態変数ベクトルを変数とした線形微分方程式である。電力系統１２の定態安定度として、上記の（５）式の係数行列Ａの固有値λｉ（ｉ＝１，・・・，ｎ）の実部が全て負のとき安定と判断され、そうでないとき不安定であると判断される。
【００５８】
また、本実施の形態では、電力系統１２の定態安定度を適切に評価できる評価関数として、以下の（６）式を用いる。
【００５９】
【数１０】
2010074995.tif 000012
【００６０】
なお、上記（６）式の評価関数に限定するものではなく、電力系統１２の定態安定度を適切に評価できる評価関数であれば、他の評価関数を用いてもよい。
【００６１】
次に、係数行列Ａの固有値について説明する。
【００６２】
まず、係数行列Ａの右固有ベクトルを並べた行列Ｐ（以下、モード行列と呼ぶ）を用いると、以下の（７）式が得られる。
【００６３】
【数１１】
2010074995.tif 000013
【００６４】
上記（７）式を用いて、上記（５）式を、以下の（８）式のように変形する。
【００６５】
【数１２】
2010074995.tif 000014
【００６６】
ここで、Ｐ−１ＡＰは、以下の（９）式に示すように、係数行列Ａの固有値（λ１、λ２、・・・、λｎ）が対角線上に並ぶ対角行列Λとなる。
【００６７】
【数１３】
2010074995.tif 000015
【００６８】
ただし、行列Ｐ−１は、行列Ａの左固有ベクトルを並べた行列である。
【００６９】
なお、固有値λｉが実数の場合は、この一組の固有値λｉと固有ベクトルＰｉのセットをモードｉと呼ぶ。また、固有値λｉが複素数の場合は、この固有値λｉと固有ベクトルＰｉ、及びこれらと各々共役な固有値λｉ＋１と固有ベクトルＰｉ＋１のセットを、モード（ｉ，ｉ＋１）と呼ぶこととする。ただし、特に混乱を生じない場合には、後者の複素数のモード（ｉ，ｉ＋１）も単にモードｉと呼ぶことにする。
【００７０】
また、以下の記述では煩雑さを避け、分かりやすさを優先するため、２つのモードを関連づけて説明するとき、まず両モードとも実数のモードの場合についてのみ記述を行うことにする。どちらかあるいは両者が複素数の場合にも容易に拡張できるため、その場合の説明は省略する。
【００７１】
次に、固有値の変化の推定方法について説明する。
【００７２】
まず、係数行列Ａに含まれるパラメータα＝α０がΔαだけ変化すると、以下の（１０）式に示すように、係数行列Ａが係数行列Ａのチルダに変化する。
【００７３】
【数１４】
2010074995.tif 000016
【００７４】
係数行列Ａが上記（１０）式のように変化したとき、変化後の係数行列Ａのチルダに対して、係数行列Ａの固有ベクトルからなるモード行列Ｐによる相似変換を行うと、以下の（１１）式が得られる。
【００７５】
【数１５】
2010074995.tif 000017
【００７６】
行列（Λ＋ΔΛ）は相似変換で得られる行列であるため、行列（Λ＋ΔΛ）の固有値は元の係数行列Ａのチルダ（＝Ａ＋ΔＡ）の固有値と厳密に一致する。
【００７７】
ここで、行列（Λ＋ΔΛ）の近似を考える。行列ΔΛのすべての要素は一般に非零であり、非対角要素はモード間の相互作用を表すが、例えば、着目するモードをモード（２，３）とし、これらのモードと強い相互作用を持つモードをモード（６，７）として、着目するモードの各々と、着目するモードと強い相互作用を持つ２個のモードの各々とに対応する非対角要素Δλｉｊを残し、それ以外の非対角要素を０とすれば、以下の（１２）式に示すような、行列ΔΛの近似行列ΔΛ’が得られる。
【００７８】
【数１６】
2010074995.tif 000018
【００７９】
また、上記（１２）式で表される近似行列ΔΛ’において、非零の非対角要素が左上に集まるように行と列を並べ替えると、以下の（１３）式に示すような行列が得られる。
【００８０】
【数１７】
2010074995.tif 000019
【００８１】
着目するモードの新しい固有値の推定値を得るために、上記（１３）式のように並び替えられた近似行列ΔΛ’と同じように並び替えた行列Λとの和の行列のうち、４×４の部分行列の固有値を求めれば良いことになる。
【００８２】
このとき、以下のように、並び替えられた行列（Λ＋ΔΛ’）の部分行列の固有値及び固有ベクトルを求める。
【００８３】
はじめに、行列（Ａ＋ΔＡ）の固有値及び固有ベクトルの算出方法について説明する。上記（７）式、（１１）式から、以下の（１４）式が得られる。
【００８４】
【数１８】
2010074995.tif 000020
【００８５】
また、行列Λ＋ΔΛの固有ベクトルを並べてできる行列Ｑを用いて、以下の（１５）式の変換を行う。
【００８６】
【数１９】
2010074995.tif 000021
【００８７】
上記（１４）式に上記（１５）式を代入して、左からＱ−１を乗算すると、以下の（１６）式が得られる。
【００８８】
【数２０】
2010074995.tif 000022
【００８９】
従って、Ｑ−１（Λ＋ΔΛ）Ｑから、行列（Λ＋ΔΛ）の固有値が得られる。
【００９０】
また、上記（７）式に、上記（１５）式を代入すると、以下の（１７）式が得られる。
【００９１】
【数２１】
2010074995.tif 000023
【００９２】
上記（１７）式の行列ＰＱから、行列（Ａ＋ΔＡ）の固有ベクトルが得られる。
【００９３】
また、行列（Λ＋ΔΛ’）について、行列（Λ＋ΔΛ’）の固有ベクトルを並べてできる行列Ｑ’を用いると、以下の（１８）式が得られる。
【００９４】
【数２２】
2010074995.tif 000024
【００９５】
また、上記（１６）式と同様に、以下の（１９）式が得られる。
【００９６】
【数２３】
2010074995.tif 000025
【００９７】
従って、Ｑ’ −１（Λ＋ΔΛ’）Ｑ’から、行列（Λ＋ΔΛ’）の固有値が得られる。
【００９８】
また、上記（１７）式と同様に、以下の（２０）式が得られる。
【００９９】
【数２４】
2010074995.tif 000026
【０１００】
上記（１７）式の行列ＰＱ’から、行列（Ａ＋ΔＡ’）の固有ベクトルが得られる。
【０１０１】
上記（１８）式〜（２０）式において、Ｐ、Ｑ、Ｑ−１について、上記（１３）式の並び替えと同様に並び替え、かつ、同様の左上の部分行列を用いることにより、少ない計算量で、上記（１９）式、（２０）式と同様の式から、並び替えられた行列（Λ＋ΔΛ’）の部分行列の固有値及び固有ベクトルが得られる。
【０１０２】
また、並び替えられた行列（Λ＋ΔΛ’）の部分行列の固有値及び固有ベクトルから、元の行列（Λ＋ΔΛ’）の固有値及び固有ベクトルが得られ、その固有値及び固有ベクトルが、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルの推定値となる。
【０１０３】
一般的に記述すれば以下のとおりである。正確に求めたいモードを着目モードと呼ぶことにし、着目モードをモードｉとする。着目モードが複数個あっても容易に拡張することができるので、ここでは、説明を複雑にしないため、着目モードを１個としておく。この着目モードｉに強い相互作用を持つ、すなわち強い影響を与えるモードを選択し結合する。これらの強い影響を与えるモードを強相関モードと呼ぶこととする。この強相関モードがｊ１，ｊ２，・・・，ｊｍのｍ個選択されたとする。着目モードとｍ個の強相関モードの各々とに対応する非対角要素Δλｉｊだけを残し、それ以外の非対角要素を０として近似する。そして、近似された行列に対して、０でない非対角要素を集めるように行及び列の並び替えを行い、０でない非対角要素が集められた部分行列の固有値を算出するようにすることにより、少ない計算量で精度良く近似固有値を求めることができる。
【０１０４】
次に、強相関モードを選択する方法について説明する。
【０１０５】
本実施の形態では、モードｉ，ｊ間の影響度ｄｉｊ（ｊ＝１，・・・，ｎ，ｊ≠ｉ）を求め、ｄｉｊが大きくなる順にｍ個のモードｊ１，ｊ２，・・・，ｊｍを選ぶ。
【０１０６】
モード間の影響度は、以下のように算出する。なお、着目モードがモードｉであり、結合対象となるモードがモードｊである場合について説明する。
【０１０７】
まず、モードｉとｊの両モードを結合して固有値を求める。すなわち、以下の（２１）式で表される行列の固有値を求める。
【０１０８】
【数２５】
2010074995.tif 000027
【０１０９】
なお、上記（２１）式で表される行列Ｍｉｊは、行列（Λ＋ΔΛ）における、モードｉ，ｊの各々に対応する対角要素、及びモードｉ，ｊに対応する非対角要素で構成される行列である。
【０１１０】
また、行列Ｍｉｊの固有値のうち、モードｉに対応する要素（λｉ＋Δλｉｉ）に対応する固有値をμｉｊとして、以下の（２２）式に従って、モードｉとｊとを結合したことにより生じる推定固有値の差の絶対値ｄｉｊを、モードｉ，ｊ間のモード間影響度として算出する。
【０１１１】
【数２６】
2010074995.tif 000028
【０１１２】
着目モード以外の各モードについて、上記（２２）式に従ってモード間影響度を算出して、モード間影響度が大きくなる順にｍ個のモードが、強相関モードとして選択される。
【０１１３】
なお、煩雑さを避けるため両モードとも実数の固有値の場合について説明したが、複素数の場合にも共役となる要素を考慮して容易に拡張できるので省略する。その場合には、組み合わせによって行列Ｍｉｊが３×３あるいは４×４の行列となる。
【０１１４】
次に、結合するモードの選択個数の決定方法について説明する。
【０１１５】
変化後の係数行列（Ａ＋ΔＡ）の真の固有値が不明な状態で、上述した方法によって推定される固有値について発生する誤差の概略値をあらかじめ知ることができれば、結合すべきモードを何個選べば十分であるかの判断を事前に行うことができる。そこで、結合すべきモードの選択個数を以下のように決定する。
【０１１６】
強相関モードを選択する際に算出される影響度ｄｉｊの値は、モードｊを結合しなかったときに発生する誤差の大きさを示したものであると考えることができる。全体でＬ個のモードがあり、着目モードを除くＬ‐１個のモードのうち、モード間影響度ｄｉｊの大きい方からｍ個のモードｊ１，ｊ２，・・・，ｊｍを選ぶとする。このとき発生する誤差は、Ｌ‐ｍ‐１個の結合されなかったモードによる誤差と考えることが出来るので、以下の（２３）式に従って、結合されなかったモードに対するモード間影響度ｄｉｊの総和を、誤差の概略推定値Ｅｓ（ｍ）として使用することができる。
【０１１７】
【数２７】
2010074995.tif 000029
【０１１８】
上記（２３）式に従って、複数種類の選択個数ｍの各々について、誤差の概略推定値Ｅｓ（ｍ）を算出し、誤差の概略推定値Ｅｓ（ｍ）が所定値以下となるときの選択個数ｍを、結合すべきモードの個数として決定する。
【０１１９】
次に、比較例としての従来の固有値感度解析法を利用した固有値の推定方法について説明する。
【０１２０】
上記（１）式において、パラメータα＝α０がｄαだけ無限小の変化をしたときの固有値λｉの変化率、すなわち固有値感度は、以下の（２４）式で求めることができる。
【０１２１】
【数２８】
2010074995.tif 000030
【０１２２】
ここで、ｑｉは左固有ベクトルであり、その大きさは、以下の（２５）式に示すように規格化しておくこととする。
【０１２３】
【数２９】
2010074995.tif 000031
【０１２４】
上記（２４）式で表される固有値感度を利用すれば、α＝α０が、αのチルダ（＝α０＋Δα）に変化したことによって係数行列Ａが係数行列Ａのチルダ（＝Ａ＋ΔＡ）に変化したときの、固有値の変化Δλｉを、以下の（２６）式に従って推定することができる。
【０１２５】
【数３０】
2010074995.tif 000032
【０１２６】
上記（２６）式による推定は、ΔλｉとΔＡの間に線形の関係が成り立っていることを前提としているので、線形推定法と呼ぶこともできる。
【０１２７】
上記（２６）式に従って、以下の（２７）式に示すような、Ｐ−１ΔＡＰの非対角要素を全て０とした近似行列が得られる。
【０１２８】
【数３１】
2010074995.tif 000033
【０１２９】
上記（２７）式により、行列（Ｐ−１ＡＰ＋Ｐ−１ΔＡＰ）の非対角要素はすべて０となるため、計算量が大幅に削減される。ただし、ΔＡが大きくなると誤差が大きくなる欠点がある。
【０１３０】
次に、本実施の形態における固有値推定方法を用いたシミュレーションの結果について説明する。
【０１３１】
まず、モデル系統と解析条件について説明する。
【０１３２】
電気学会ＥＡＳＴ３０（ｐｅａｋ）標準モデル系統（非特許文献：「電力系統の標準モデル」、電気学会技術報告第８９８号、１９９９年）を用いて、潮流の変化に伴う固有値の変化を、本実施の形態における固有値推定方法を用いて推定し、その精度について検討を行った。潮流の変化については、ノード番号２００１に接続されている発電機１番（Ｇ１）の出力ＰＧ１をベースケースの３．８ｐ．ｕ．から最大４．８ｐ．ｕ．まで０．２ｐ．ｕ．ずつ増加させた。また、各ケースについて、本実施の形態における固有値推定方法による推定固有値を求め、真の固有値との比較を行った。
【０１３３】
スィング発電機は発電機２５番（ノード番号３０１０）とし、ここで変化分を吸収させた。界磁巻線及びｄ，ｑ軸に制動巻線を各１個とした３次のＰａｒｋ式に、ロータの運動方程式を含めて、発電機モデルを５次の微分方程式で表わし、ＡＶＲは直流励磁機方式の３次の標準定数モデル（Ｌａｔ＝１）とした。ガバナタービンモデルは省略し、機械入力は一定とした。従って、発電機１機あたり８次の微分方程式で表わし、電力系統全体で、２４０次の微分方程式を用いて表わした。
【０１３４】
以下の表1に、各ケースでのこの電力系統の固有値を示す。なお、ＰＧ１＝４．８ｐ．ｕ．では電力系統の定態安定度が不安定となった。これは系統の最も周期が長く不安定となりやすい振動モードの固有値であり、主動揺モードと呼ぶべきモードであるが、この主動揺モードを、本解析における着目モードとした。
【０１３５】
【表１】
2010074995.tif 000034
【０１３６】
また、本実施の形態における固有値推定方法による固有値の推定は、ベースケースで求めた固有ベクトルを用いて行った。また、潮流計算および係数行列Ａの計算については、線形近似などにより簡略計算することはせず、すべてのケースで厳密計算を行った。
【０１３７】
図３に、上記で説明したモード間影響度の算出方法を用いて計算した着目モードｉに対するモードｊのモード間影響度の大きさｄｉｊを示す。図３の横軸はモード番号であり、モード間影響度の大きい順番に並べ直してある。また、３本の曲線は、ＰＧ１＝４．０， ４．４， ４．８ ［ｐ．ｕ．］の各ケースに対応している。
【０１３８】
次に、本実施の形態における固有値推定方法による推定固有値の精度について説明する。
【０１３９】
図４（Ａ）、（Ｂ）に、着目モードの推定固有値の実部と虚部を示す。横軸ＰＧ１はＧ１の有効電力出力であり、パラメータｎは推定に用いたモードの総次数である。すなわち着目モードと強相関モードで、複素数モードの次数を２、実数のモードの次数を１としたときの次数の総和である。
【０１４０】
ｎ＝２のときは、着目モードのみを用いて推定した場合、すなわち、比較例として上記で説明した固有値感度による線形推定法を用いた場合であり、モード間の相互作用を全く無視した場合と同等である。また、ｎ＝２４０のときは固有値の厳密解と完全に一致する。ＰＧ１＝４．０ｐ．ｕ．では、線形推定法でもそれほど大きな誤差は生じないが、ＰＧ１＝４．６ｐ．ｕ．以上では、非線形性が強く線形推定では大きな誤差となる。本実施の形態の固有値推定方法ではｎ＝６でも大幅な改善がなされており、ｎ＝３０（元の次数の１／８）で実用上支障が無い程度に精度が改善されている。
【０１４１】
また、図５に推定固有値の誤差を示す。横軸は推定に用いたモードの総次数ｎである。ＰＧ１が大きいほど誤差は大きくなるが、ＰＧ１＝４．８の場合でもｎ＝３０程度とすれば誤差が１％程度になることが分かった。
【０１４２】
図６に、推定固有値の誤差について上記（２３）式による概略推定値と実際に発生した誤差の相関を示した。両者には十分な相関があり、上記（２３）式による誤差の推定法が有用かつ有効であることが分かった。
【０１４３】
次に、本実施の形態に係る固有値推定方法を用いた場合の計算時間について説明する。以下の表２に示すように、行列のサイズを１６、１８、２４０としたときの計算時間を計測した。本実施の形態に係る固有値推定方法を用いることにより、短い計算時間で固有値を推定することができることが分かった。
【０１４４】
【表２】
2010074995.tif 000035
【０１４５】
次に、第１の実施の形態に係る電力系統システム１０の作用について説明する。
【０１４６】
まず、自動あるいはオペレータの指示によって、系統監視制御システムから、コンピュータ１６に、電力系統１２の系統情報として、系統構成や運転状態を示すデータが入力され、コンピュータ１６は、系統構成や運転状態を示すデータに基づいて潮流計算を行い、計算結果と系統情報とを記憶しておく。
【０１４７】
そして、コンピュータ１６において、図７に示す制御パラメータ最適化処理ルーチンが実行される。
【０１４８】
まず、ステップ１００において、系統情報及び潮流計算結果を読み込み、ステップ１０２において、ＰＳＳ定数αに初期値を設定する。
【０１４９】
そして、ステップ１０４において、ＰＳＳ定数αに基づいて、線形微分方程式の係数行列Ａを計算し、ステップ１０６で、係数行列Ａの固有値及び固有ベクトルを計算する。
【０１５０】
次のステップ１０８では、変数Ｌに初期値の１を設定し、ＰＳＳ定数αを修正する。そして、ステップ１１０において、上記ステップ１０８又は後述するステップ１２４で修正されたαの変化に応じて、係数行列Ａの変化分ΔＡを算出する。
【０１５１】
ステップ１１２において、変数Ｌが１であるか否かを判定し、変数Ｌが初期値１である場合には、ステップ１１４において、強相関モードを選択するが、一方、変数Ｌが２以上である場合には、ステップ１１６へ移行する。
【０１５２】
上記ステップ１１４は、図８に示す強相関モード選択処理ルーチンによって実現される。以下、強相関モード選択処理ルーチンについて説明する。
【０１５３】
まず、ステップ１４０において、上記ステップ１０４、１０６で計算された係数行列Ａ、係数行列Ａの固有値、固有ベクトルと、上記ステップ１１０で算出された行列ΔＡとから得られる行列（Ｐ−１ＡＰ＋Ｐ−１ΔＡＰ）に基づいて、予め設定された着目モードと結合対象モードとを結合した行列Ｍを作成し、行列Ｍの固有値を計算する。
【０１５４】
そして、ステップ１４２において、上記ステップ１４０で計算された行列Ｍの固有値を用いて、着目モードと対象モードとのモード間影響度を算出する。ステップ１４４では、着目モード以外の全てのモードを、結合対象モードとして、モード間影響度を算出したか否かを判定し、モード間影響度を算出していないモードがある場合には、結合対象モードとして、算出していないモードを設定して、ステップ１４０へ戻る。一方、着目モード以外の全てのモードを結合対象モードとして、モード間影響度を算出した場合には、ステップ１４６へ進む。
【０１５５】
ステップ１４６では、上記ステップ１４２で算出されたモード間影響度が大きい順に、着目モード以外のモードを並べ、ステップ１４８において、結合対象モードの選択個数を示す変数ｋを初期値の１に設定する。なお、着目モード以外のモードの個数をＬ−１個とする。
【０１５６】
そして、ステップ１５０において、モード間影響度が大きい順にｋ個のモードを選択した場合に、残りの（Ｌ−ｋ−１）個のモードによる誤差を推定する。次のステップ１５２では、ｋがＬ−１未満であるか否かを判定し、ｋがＬ−１未満である場合には、ステップ１５４において、変数ｋをインクリメントして、ステップ１５０へ戻る。一方、ｋがＬ−１に達した場合には、ステップ１５６へ進む。
【０１５７】
ステップ１５６では、上記ステップ１５０で推定された誤差に基づいて、誤差が所定値以下となるときの選択個数ｍを決定する。そして、ステップ１５８において、モード間影響度が大きい方からｍ個のモードを強相関モードとして選択して、強相関モード選択処理ルーチンを終了する。
【０１５８】
そして、上記図７のステップ１１６において、行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。ステップ１１６は、図９に示す固有値推定処理ルーチンによって実現される。以下に、固有値推定処理ルーチンについて説明する。
【０１５９】
まず、ステップ１７０において、上記ステップ１０４、１０６で計算された係数行列Ａの固有ベクトルと、上記ステップ１１０で算出された行列ΔＡとから得られる行列ΔΛ（＝Ｐ−１ΔＡＰ）から、行列ΔΛの近似行列ΔΛ’を生成する。
【０１６０】
そして、ステップ１７２において、上記ステップ１７０で生成された近似行列ΔΛ’と行列Λ（＝Ｐ−１ＡＰ）との和の行列（Λ＋ΔΛ’）について、非零の非対角要素を左上に集めるように行と列とを並び替える。
【０１６１】
次のステップ１７４では、上記ステップ１７２で並び替えられた行列のうちの、非零の非対角要素が集められた部分からなる部分行列の固有値及び固有ベクトルを算出し、算出した部分行列の固有値及び固有ベクトルに基づいて、行列（Λ＋ΔΛ’）を並び替えた行列の固有値及び固有ベクトルを算出する。そして、算出した固有値及び固有ベクトルから得られる元の行列（Λ＋ΔΛ’）の固有値及び固有ベクトルを、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルの推定値とする。
【０１６２】
そして、図７のステップ１１８において、上記ステップ１１６で推定された変化後の係数行列（Ａ＋ΔＡ）の固有値を用いて、上記（６）式に示す評価関数ｆに従って、評価値を算出する。次のステップ１２０では、上記ステップ１１６で推定された変化後の係数行列（Ａ＋ΔＡ）の固有値λｉ及び固有ベクトルｑｉ，ｐｉを用いて、以下の（２０）式、（２１）式に従って、Δλｉ、Δｆを算出し、Δｆを用いて、Δｆ／Δαを算出する。
Δλｉ＝（ｑｉΔＡｐｉ）／ｑｉｐｉ      ・・・（２０）
Δｆ＝ΣΔｆ（Δλｉ）      ・・・（２１）
【０１６３】
そして、ステップ１２２において、上記ステップ１１８で算出された評価関数ｆの値に基づいて、ＰＳＳ定数αは準最適であるか否かを判定する。ＰＳＳ定数αが準最適ではないと判定された場合には、ステップ１２４において、上記ステップ１２０で算出されたΔｆ／Δαに基づいて、最大傾斜の方向にαを修正し、ステップ１２６で、変数Ｌをインクリメントして、上記ステップ１１０へ戻る。
【０１６４】
上記ステップ１２２でＰＳＳ定数αが準最適であると判定されるまで、上記ステップ１１０〜ステップ１２０を繰り返すことにより、図１０に示すように、ＰＳＳ定数αを最適値の近くまで修正することができる。
【０１６５】
上記ステップ１２２において、ＰＳＳ定数αが準最適であると判定された場合には、ステップ１２８において、上記ステップ１１８で算出された評価関数ｆの値に基づいて、ＰＳＳ定数αが最適化されたか否かを判定する。評価関数ｆの値が最小化されていない場合には、ＰＳＳ定数αが最適化されていないと判定し、上記ステップ１０４へ戻るが、一方、評価関数ｆの値が最小化された場合には、ＰＳＳ定数αが最適化されたと判定し、ステップ１３０において、最適化されたＰＳＳ定数を、電力系統安定化装置１４に設定して、制御パラメータ最適化処理ルーチンを終了する。
【０１６６】
以上説明したように、第１の実施の形態に係る電力系統システムによれば、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列（Λ＋ΔΛ）を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができる。また、電力系統の安定度を精度良く評価することができる。
【０１６７】
また、近似された行列（Λ＋ΔΛ’）の０でない非対角要素が集められた部分行列の固有値及び固有ベクトルを算出することにより、計算量を減らすことができる。
【０１６８】
また、一般に電力系統の各種パラメータが変化すると定態安定度解析で用いる係数行列が変化し、さらに係数行列の変化に対して固有値が変化する。係数行列の変化も固有値の変化も非線形な変化となるが、係数行列の変化に伴う固有値の変化は極めて強い非線形性を示すことがある。本実施の形態の固有値推定方法では、着目モードと強相関モードとに対応する非対角要素を残し、それ以外の非対角要素を０とすることにより、この非線形な特性を持つ固有値の変化を精度良く推定することが可能となる。
【０１６９】
また、係数行列Ａの全ての固有値の計算には行列のサイズｎの３乗に比例する計算時間を必要とするため、電力系統のような大規模システムでは、固有値の計算時間が長大となるが、本実施の形態の固有値推定方法に用いる行列のサイズは、本計算例では、元の行列のおよそ１／８程度で十分であり、固有値推定に要する計算時間を、元の行列の計算時間の１％以下にすることができる。
【０１７０】
なお、上記の実施の形態では、入力された電力系統の系統情報に基づいて、制御パラメータを最適化する場合を例に説明したが、これに限定されるものではなく、電力系統の系統情報をリアルタイムに取得するようにして、電力系統のリアルタイムな状態変化に応じて、最適な状態に制御パラメータをリアルタイムに修正するようにしてもよい。これによって、適応制御、すなわち時々刻々の系統の状態に応じて常に最も安定な状態に電力系統を維持する制御が可能となる。
【０１７１】
次に第２の実施の形態について説明する。第２の実施の形態では、想定された事故等に対し電力系統が安定になるように予防制御を行う電力系統システムに、本発明を適用した場合を例に説明する。なお、第１の実施の形態と同様の構成となる部分については、同一符号を付して説明を省略する。
【０１７２】
図１１に示すように、第２の実施の形態に係る電力系統システム２１０は、電力系統１２と、電力系統１２の運転操作を行なうための端末である系統監視制御システム２１４と、電力系統１２の系統構成や運転状態を示す系統情報、及び発電機の潮流データを読み込んで、想定される事故に対する定態安定度余裕を算出し、全ての発電機が安定となるような運転操作パラメータを計算するコンピュータ２１６とを備えている。
【０１７３】
コンピュータ２１６は、ＣＰＵ、後述する予防制御処理ルーチンを実現するためのプログラムを記憶したＲＯＭ、データを一時的に記憶するＲＡＭ、及びＨＤＤ等の記憶装置で構成されている。
【０１７４】
次に、予防制御の方法について説明する。
【０１７５】
まず、電力系統１２の複数の発電機に対する通常の潮流データを取得した後、想定される系統事故あるいは補修作業による送電線停止時の潮流計算を行い、各々の状態での固有値解析を行う。ここで、運転操作パラメータをわずかに変化させて各発電機の出力をわずかに変化させたときの固有値の変化量も求めておく。この固有値を用いて定態安定度余裕を算定すると共に、予防制御のための最適な発電機出力配分が得られる運転操作パラメータを算定する。
【０１７６】
ここで、各発電機の出力指定値（運転操作パラメータ）を増減させたときに電力系統が不安定（固有値の実部が正）あるいは警戒状態（固有値の実部が予め定められた負の値よりも大）となる出力指定値の変化量Δα＋、Δα‐を、上述した固有値推定法で求め、この値を発電機定格出力で割った値を、各発電機の定態安定度余裕とする。
【０１７７】
また、予防制御では、想定される事故等に対し電力系統が不安定になると予想される場合には、安定性が維持させるように、安定化に最も効果的な発電機出力配分が得られる操作パラメータを求め、求められた操作パラメータに予め修正する。また、想定される事故等に対し電力系統が安定であると予想される場合には、安定性を維持できる範囲で最も経済的な発電機出力配分となる操作パラメータを求め、求められた操作パラメータに予め修正する。
【０１７８】
次に、安定度余裕の計算方法について説明する。ここで、説明を簡単にするために、電力系統のある安定な状態における発電機ｊの出力指定値αＧｊをαとし、不安定となる固有値をλとする。αが、α＋Δαｋ（ｋ＝１、２、３、・・・）と変化したときの固有値λｋを、上記第１の実施の形態で説明した固有値推定方法によって求める。このときの固有値の実部及びαの変化を図１２に示す。上記図１２では、縦軸を固有値λの実部ρ、横軸をα＋Δαｋとしている。
【０１７９】
上記図１２において、不安定となる境界（ρ＝０）とαの変化に応じた固有値の実部の変化を表わす曲線との交点を求めることにより、不安定状態までの安定度余裕を求めることができる。警戒状態までの余裕も同様に求めることができる。
【０１８０】
複数の固有値が不安定となる場合は、より小さな変化に対して不安定となる固有値によって安定度余裕を求める。また、Δαを負の値とするときにも不安定となることがあるので、正の場合をΔα＋、負の場合をΔα−と表わす。
【０１８１】
次に、第２の実施の形態に係る電力系統システム２１０の作用について説明する。
【０１８２】
まず、オペレータによって、コンピュータ２１６に、解析基本データと想定される全ての事故を示す条件とが入力され、コンピュータ２１６に記憶される。例えば、「３０分後に補修作業のためにある送電線Ａの停止が予定されていたところ、別の送電線Ｂを雷雲が通過するとの予報が入った」という状況において、送電線Ｂ事故停止、及び送電線Ａの作業停止の各々を示す条件が入力される。
【０１８３】
そして、コンピュータ２１６において、図１２に示す予防制御処理ルーチンが実行される。なお、第１の実施の形態と同様の処理については、同一符号を付して詳細な説明を省略する。
【０１８４】
まず、ステップ２５０において、基本解析データや想定事故条件を読み込み、ステップ２５２で、操作パラメータである運転指令値αに初期値を設定する。
【０１８５】
そして、ステップ２５３で、電力系統１２から、現在及び一定時間後に関する系統構成や運転状態を示す系統情報と、現在及び一定時間後に関する潮流データとを取得する。
【０１８６】
次のステップ１０４では、運転指令値αに基づいて、線形微分方程式の係数行列Ａを計算し、ステップ１０６で、係数行列Ａの固有値及び固有ベクトルを計算する。
【０１８７】
次のステップ２５４では、上記ステップ２５０で読み込まれた想定事故条件が示す各想定事故に対する安定度余裕を算出する。上記ステップ２５４は、図１３に示す処理ルーチンによって実現される。
【０１８８】
まず、ステップ２８０において、上記ステップ２５０で読み込まれた想定事故条件が示す想定事故のうちの何れか１つの想定事故Ｊを選択する。そして、ステップ２８２において、運転指令値の中から１つの運転指令値αＩを選択する。
【０１８９】
次のステップ２８４では、想定される運転指令値αＩの最大変化幅Δαだけ運転指令値αＩを変化させたときの係数行列Ａの変化分ΔＡを算出し、ステップ１１４において、強相関モードを選択する。
【０１９０】
そして、ステップ２８６において、上記で説明したように、運転指令値αＩを徐々に増大させながら、増大した運転指令値αＩ＋Δα＋の各々について、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定することにより、αＩを増大させたときの想定事故Ｊにおける安定度余裕を算出する。また、運転指令値αＩを徐々に減少させながら、減少した運転指令値αＩ＋Δα−の各々について、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定することにより、αＩを減少させたときの想定事故Ｊにおける安定度余裕を算出する。
【０１９１】
次のステップ２８８では、全ての運転指令値αＩについて、想定事故Ｊにおける安定度余裕を算出したか否かを判定し、運転指令値αのうち、想定事故Ｊにおける安定度余裕を算出していない運転指令値αＩがある場合には、上記ステップ２８２へ戻り、他の運転指令値αＩについて想定事故Ｊにおける安定度余裕を算出する。一方、上記ステップ２８８において、全ての運転指令値αＩについて、想定事故Ｊにおける安定度余裕を算出したと判定された場合には、ステップ２９０へ移行する。
【０１９２】
ステップ２９０では、全ての想定事故における安定度余裕を算出したか否かを判定し、読み込んだ想定事故条件が示す想定事故のうち、安定度余裕を算出していない想定事故がある場合には、上記ステップ２８０へ戻り、他の想定事故における安定度余裕を算出する。一方、上記ステップ２９０において、全ての想定事故における安定度余裕を算出したと判定された場合には、処理ルーチンを終了する。
【０１９３】
そして、予防制御処理ルーチンのステップ２５６において、運転指令値αに基づいて、例えば運転コストを評価するための評価関数を算出すると共に、評価関数の偏微分係数を算出する。次のステップ２５８では、上記ステップ２５４で算出された安定度余裕に基づいて、全ての想定事故において安定であるか否かを判定する。
【０１９４】
全ての想定事故に対して、安定であると判定されなかった場合には、ステップ２６０において、最大傾斜法により運転指令値αを修正して、ステップ２５４へ戻る。
【０１９５】
上記ステップ２５８において、全ての想定事故に対して、安定であると判定された場合には、ステップ２６２において、上記ステップ２５６で算出した評価関数の値に基づいて、運転指令値αが準最適であるか否かを判定し、運転指令値αが準最適でないと判定された場合には、ステップ２６０へ移行する。一方、運転指令値αが準最適であると判定された場合には、ステップ２６４において、運転指令値αの修正量Δαは小さいか否かを判定し、修正量Δαが大きい場合には、ステップ２６６で、運転指令値αをα＋Δαに修正して、ステップ１０４へ戻る。
【０１９６】
上記ステップ２６４で、運転指令値αの修正量Δαが小さいと判定された場合には、ステップ２６８において、運転指令値αをオペレータに対して表示して、予防制御処理ルーチンを終了する。
【０１９７】
オペレータは、コンピュータ２１６で表示された運転指令値αを、系統監視制御システム２１４に設定して、運転指令値を最適な値に変更する。
【０１９８】
このように、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列（Λ＋ΔΛ）を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列の変化に対する固有値の変化を精度良く推定することができる。また、電力系統の定態安定度余裕を精度良く算出することができる。
【０１９９】
また、本実施の形態の固有値推定方法を用いることによって、高速かつ高精度な固有値の推定が可能となり、高速性を要求される電力系統の定態安定度の信頼度監視及び予防制御が可能となる。
【０２００】
次に、第３の実施の形態について説明する。第３の実施の形態では、プラントの運転状態を目的とする状態にするように制御するプラントシステムに本発明を適用した場合を例に説明する。なお、第１の実施の形態と同様の構成となる部分については、同一符号を付して説明を省略する。
【０２０１】
図１４に示すように、第３の実施の形態に係るプラントシステム３１０は、プラント３１２と、プラント３１２の運転操作を行なうための端末であるプラント監視制御装置３１４と、プラント３１２の構成や運転状態を示すプラント情報を読み込んで、最適な状態となるような運転操作パラメータを計算し、プラント監視制御装置３１４に対して、計算した運転操作パラメータに変更するように指令するコンピュータ３１６とを備えている。
【０２０２】
コンピュータ３１６は、ＣＰＵ、後述する制御パラメータ最適化処理ルーチンを実現するためのプログラムを記憶したＲＯＭ、データを一時的に記憶するＲＡＭ、及びＨＤＤ等の記憶装置で構成されている。
【０２０３】
次に、第３の実施の形態に係るプラントシステム３１０の作用について説明する。
【０２０４】
まず、オペレータによって、コンピュータ３１６に、プラント３１２に関する基本データと、プラント３１２の運転に関する制約条件と、プラント３１２の運転状態が目的とする望ましい状態であるかを評価するための評価関数とが入力され、コンピュータ３１６に記憶される。
【０２０５】
そして、コンピュータ３１６において、図１５に示す制御パラメータ最適化処理ルーチンが実行される。なお、第１の実施の形態と同様の処理については、同一符号を付して詳細な説明を省略する。
【０２０６】
まず、ステップ３５０において、基本データや、制約条件、評価関数を読み込み、ステップ２５２で、操作パラメータαに初期値を設定する。
【０２０７】
そして、ステップ３５４で、プラント３１２から構成や運転状態を示すプラント情報を取得する。次のステップ１０４では、操作パラメータαに基づいて、線形微分方程式の係数行列Ａを計算し、ステップ１０６で、係数行列Ａの固有値及び固有ベクトルを計算する。
【０２０８】
次のステップ１０８では、変数Ｌに初期値の１を設定し、運転指令値αを修正する。そして、ステップ１１０において、上記ステップ１０８又は後述するステップ１２４で修正されたαの変化に応じて係数行列Ａの変化分ΔＡを算出する。
【０２０９】
ステップ１１２において、変数Ｌが１であるか否かを判定し、変数Ｌが初期値１である場合には、ステップ１１４において、強相関モードを選択するが、一方、変数Ｌが２以上である場合には、ステップ１１６へ移行する。
【０２１０】
そして、ステップ１１６において、変化後の係数行列（Ａ＋ΔＡ）の固有値及び固有ベクトルを推定する。
【０２１１】
次のステップ３５６では、上記ステップ１１６で推定された変化後の係数行列（Ａ＋ΔＡ）の固有値に基づいて、上記ステップ３５０で読み込んだ評価関数ｆの値を算出し、ステップ１２０において、上記ステップ１１６で推定された変化後の係数行列（Ａ＋ΔＡ）の固有値λｉ及び固有ベクトルｑｉ，ｐｉを用いて、Δλｉ、Δｆを算出し、Δｆ／Δαを算出する。
【０２１２】
そして、ステップ３５８において、上記ステップ３５６で算出された評価関数ｆの値に基づいて、操作パラメータαは最適であるか否かを判定する。評価関数ｆの値が、最大（又は最小）ではなく、または、操作パラメータαが上記ステップ３５０で読み込んだ制約条件を満たさない場合には、操作パラメータαは最適でないと判定されて、ステップ１２４へ移行する。ステップ１２４では、上記ステップ１２０で算出されたΔｆ／Δαに基づいて、最大傾斜の方向にαを修正し、ステップ１２６で、変数Ｌをインクリメントして、上記ステップ１１０へ戻る。
【０２１３】
上記ステップ３５８において、評価関数ｆの値が、最大（又は最小）であり、かつ、操作パラメータαが上記ステップ３５０で読み込んだ制約条件を満たす場合には、操作パラメータαは最適であると判定して、ステップ３６０へ移行する。
【０２１４】
ステップ３６０では、修正された操作パラメータαに変更するように、プラント監視制御装置３１４に指令を出す。そして、ステップ３６２において、処理を終了するか否かを判定し、終了しない場合には、ステップ３５４へ戻るが、処理を終了する場合には、制御パラメータ最適化処理ルーチンを終了する。
【０２１５】
以上説明したように、第３の実施の形態に係るプラントシステムによれば、着目モード及び強相関モードに対応する非対角要素以外の非対角要素を０として、行列（Λ＋ΔΛ）を近似することにより、少ない計算量で、固有値及び固有ベクトルを算出することができると共に、係数行列Ａの変化に対する固有値の変化を精度良く推定することができる。また、プラントの運転状態を精度良く評価することができる。また、時々刻々変化するプラントの状態に応じて最適な制御を行うことができる。
【０２１６】
次に、第４の実施の形態について説明する。第４の実施の形態では、選択個数の決定方法が、第１の実施の形態と異なっている。
【０２１７】
第４の実施の形態では、以下に説明するように、強相関モードを選択する。
【０２１８】
本実施の形態では、着目モードをｉとし、結合対象モードをｊとして、モードｉ，ｊ間の影響度ｄｉｊ（ｊ＝１，・・・，ｎ，ｊ≠ｉ）を、以下の（２８）式に従って算出する。
【０２１９】
【数３２】
2010074995.tif 000036
【０２２０】
上記（２８）式は、行列（Λ＋ΔΛ）における、モードｉ，ｊに対応する非対角要素の積の絶対値を表している。
【０２２１】
そして、ｄｉｊが大きくなる順にｍ個のモードｊ１，ｊ２，・・・，ｊｍを強相関モードとして選択する。
【０２２２】
なお、他の構成及び作用については、上記の実施の形態と同様であるため、説明を省略する。
【０２２３】
なお、上記の実施の形態では、係数行列Ａの固有値及び固有ベクトルを厳密計算する毎に、強相関モードを選択する場合を例に説明したが、これに限定されるものではなく、αを修正する毎に、強相関モードを選択するようにしてもよい。この場合には、αを修正する毎に、各モードについてモード間影響度を算出して、モード間影響度が大きい方から、強相関モードを選択すればよい。
【０２２４】
また、本発明に係るプログラムを、記録媒体に格納して提供することも可能である。
【図面の簡単な説明】
【０２２５】
【図１】本発明の第１の実施の形態に係る電力系統システムの構成を示す概略図である。
【図２】一機無限大母線系統を示す図である。
【図３】各モードに対する着目モードとのモード間影響度を示すグラフである。
【図４】（Ａ）推定に用いたモードの総次数の各々における、有効電力出力と着目モードの推定固有値の実部との関係を示すグラフ、及び（Ｂ）推定に用いたモードの総次数の各々における、有効電力出力と着目モードの推定固有値の虚部との関係を示すグラフである。
【図５】推定に用いたモードの総次数と推定固有値の誤差との関係を示すグラフである。
【図６】概略推定値と実際に発生した誤差とを示すグラフである。
【図７】第１の実施の形態に係る電力系統システムのコンピュータにおける制御パラメータ最適化処理ルーチンの内容を示すフローチャートである。
【図８】第１の実施の形態に係る電力系統システムのコンピュータにおける強相関モード選択処理ルーチンの内容を示すフローチャートである。
【図９】第１の実施の形態に係る電力系統システムのコンピュータにおける固有値推定処理ルーチンの内容を示すフローチャートである。
【図１０】ＰＳＳ定数αを最適値の近くまで修正する様子を示すイメージ図である。
【図１１】本発明の第２の実施の形態に係る電力系統システムの構成を示す概略図である。
【図１２】第２の実施の形態に係る電力系統システムのコンピュータにおける予防制御処理ルーチンの内容を示すフローチャートである。
【図１３】第２の実施の形態に係る電力系統システムのコンピュータにおける安定度余裕を算出する処理ルーチンの内容を示すフローチャートである。
【図１４】本発明の第３の実施の形態に係るプラントシステムの構成を示す概略図である。
【図１５】第３の実施の形態に係るプラントシステムのコンピュータにおける制御パラメータ最適化処理ルーチンの内容を示すフローチャートである。
【符号の説明】
【０２２６】
１０、２１０  電力系統システム
１２   電力系統
１４   電力系統安定化装置
１６、２１６、３１６ コンピュータ
２１４  系統監視制御システム
３１０  プラントシステム
３１２  プラント
３１４  プラント監視制御装置
【図１】
2010074995.tif 000042
【図２】
2010074995.tif 000043
【図３】
2010074995.tif 000044
【図４】
2010074995.tif 000045
【図５】
2010074995.tif 000046
【図６】
2010074995.tif 000047
【図７】
2010074995.tif 000048
【図８】
2010074995.tif 000049
【図９】
2010074995.tif 000050
【図１０】
2010074995.tif 000051
【図１１】
2010074995.tif 000052
【図１２】
2010074995.tif 000053
【図１３】
2010074995.tif 000054
【図１４】
2010074995.tif 000055
【図１５】
2010074995.tif 000056