(57)【要約】【課題】  多変数システムを制御するための最適化アルゴリズムを提供する。【解決手段】  動的システムの実時間コントロールが、可能な限り所望の応答を生じることに近づくコントロール変数を決定することにより提供される。物理的制限の追加的約因(additional consideration)は、実時間で解くことが要求される不等式制約を有する凸二次計画をもたらす。新規なアクティブセットアルゴリズムは、実時間の要求を満足する凸二次計画効率を効率的に解くために述べられる。システムの物理的制限が、実際の最適なコントロールが変化するかもしれないとしても）時間が経っても比較的変わらずに維持される最適なアクティブセットに形を変えるという重要な知見に基づき、前の時間周期における最後の反復から得られたアクティブセットについてのスターティング予測が反復を大幅に減少させ、従って、二次計画が実時間内に収束して解かれることを可能にする。【選択図】  図１
2004288161.tif 000002
【特許請求の範囲】
【請求項１】
ａ）前記システムの現在の状態を示す複数のセンサ信号であってサンプリングされた信号を受信するステップと、
ｂ）複数のコマンドを受信するステップと、
ｃ）前記コマンドおよび前記センサ信号に基づき前記システムの所望の動的応答を決定するステップと、
ｄ）複数の時間ステップのそれぞれにおいて、前記システムを制御することの問題を定式化して、二次計画問題に対する解として前記所望の動的応答を達成するステップと、
ｅ）最適なアクティブセットを探索する反復アルゴリズムを用いて各時間ステップにおいて前記二次計画問題を解くステップと、
ｆ）前記複数の時間ステップの後に続く各時間ステップにおいて、前記複数の時間ステップの前の時間ステップの最後のアクティブセットで最善のアクティブセットに対する探索を初期化するステップと、
を具備する多変数システムを制御するための方法。
【請求項２】
前記アクティブセットが、最適化された解で拘束力を有する制約のセットを具備することを特徴とする請求項１に記載された方法。
【請求項３】
前記ステップｅ）が、負のラグランジェ乗数を備える少なくとも一つの制約を落とす手順を更に具備することを特徴とする請求項２に記載された方法。
【請求項４】
前記ステップｅ）が、負のラグランジェ乗数を備える全ての制約を落とす手順を更に具備することを特徴とする請求項３に記載された方法。
【請求項５】
ＱＲ分解を用いて縮退のイベントにおける堅牢性を改善することを特徴とする請求項１に記載された方法。
【請求項６】
前記反復アルゴリズムの反復数が、時間予測可能とする目的で固定されたことを特徴とする請求項１に記載された方法。
【請求項７】
次の反復が前の反復に関する改善であるような方法で全ての反復で二次計画問題を定式化するステップを更に具備することを特徴とする請求項１に記載された方法。
【請求項８】
Ｈのコレスキー分解を前計算するステップと、
前記ステップｅ）の反復アルゴリズムの複数の反復において、前記Ｈの前計算されたコレスキー分解を用いるステップと、
を更に具備することを特徴とする請求項１に記載された方法。
【請求項９】
Ｈの一般化された二乗根マトリックスのＱＲ分解を経て前記コレスキー分解を決定するステップを更に具備することを特徴とする請求項１に記載された方法。
【請求項１０】
多変数のコントロールシステムであって、
前記システムの現在の状態を示す複数のセンサと、
エフェクタコマンドに応答して前記システムのダイナミクスを変えるための複数のエフェクタと、
所望の動的応答を決定する動的フィードフォワードモジュールと、
所望のシステム応答と実際のシステム応答を比較して必要に応じて出力変更コマンドを調整する動的モデルフォロアと、
前記システムを制御することの問題を定式化することにより前記動的フィードフォワードおよび動的モデルフォロアモジュールによる決定に基づき前記システムの出力を変えるために前記エフェクタコマンドを決定して、複数の時間ステップのそれぞれにおいて二次計画問題に対する解として前記所望の動的応答を達成する動的インバージョンモジュールと、
最適なアクティブセットを探索する反復アルゴリズムを用いて各時間ステップにおいて二次計画問題に対して最適化された解を見つけ出す二次計画法モジュールと、を具備し、
前記二次計画法モジュールが、前の時間ステップの最後のアクティブセットで後に続く各時間ステップにおける最善のアクティブセットの探索を初期化する多変数コントロールシステム。
【請求項１１】
前記二次計画法モジュールが、前記制約のサブセットのラグランジェ乗数に基づき前記アクティブセットを予測することを特徴とする請求項１０に記載されたシステム。
【請求項１２】
前記二次計画法モジュールが、負のラグランジェ乗数を備える全ての制約を落とすことを特徴とする請求項１１に記載されたシステム。
【請求項１３】
前記アクティブセットが、最適化された解で拘束力を有する制約のセットを具備することを特徴とする請求項１０に記載されたシステム。
【請求項１４】
ａ）システムの現在の状態を示す複数のセンサ信号であってサンプリングされた信号を受信するステップと、
ｂ）複数のコマンドを受信するステップと、
ｃ）前記コマンドおよび前記センサ信号に基づき前記システムの所望の動的応答を決定するステップと、
ｄ）複数の時間ステップのそれぞれにおいて、前記システムを制御することの問題を定式化して、二次計画問題に対する解として前記所望の動的応答を達成するステップと、
ｅ）最適なアクティブセットを探索する反復アルゴリズムを用いて各時間ステップにおいて前記二次計画問題を解くステップと、
ｆ）前記複数の時間ステップの後に続く各時間ステップにおいて、前記複数の時間ステップの前の時間ステップの最後のアクティブセットで最善のアクティブセットに対する探索を初期化するステップと、
をコンピュータに実行させるコンピュータプログラムを格納するコンピュータ読み取り可能な記録媒体。
【請求項１５】
前記アクティブセットが、最適化された解で拘束力を有する制約のセットを具備することを特徴とする請求項１４に記載されたコンピュータ読み取り可能な記録媒体。
【請求項１６】
前記ステップｆ）が、前記制約のサブセットのラグランジェ乗数に基づき前記アクティブセットを予測するステップを更に具備することを特徴とする請求項１５に記載されたコンピュータ読み取り可能な記録媒体。
【請求項１７】
前記ステップｆ）が、負のラグランジェ乗数を備える全ての制約を落とすステップを更に具備することを特徴とする請求項１６に記載されたコンピュータ読み取り可能な記録媒体。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、概してコントロールシステムに関し、更に詳しくは二次計画法(quadratic programming)を用いた動的インバージョンコントロールシステム(dynamic inversion control system)に関する。
本発明は、米国政府契約Ｎ００４２１−０１−２−０１３１の下での作業を遂行する中で着想された。
【背景技術】
【０００２】
時間的に変化するシステムは、時間瞬間ごとに出力および性能指数(performance index)を入力に関連付ける線形代数方程式のセットによりモデル化される（非線形動的システムでは、これらの数式は時刻ごとに変化する）。これらの入力は、形式的に、コントロール変数として知られている。なぜなら、これら変数の値は、オペレータのコントロールの下にあるからである。その目的は、システムの性能指数を最大化するコントロール変数の最適なシーケンス(sequence)を時間内に見つけ出すことである。加えて、重大なシステム故障を防止するためには、システムの物理的制約が尊重されなければならない。コントロール変数の最適な数列を決定する問題は、最適化問題、即ち時間における瞬間ごとの新たな最適化問題のシーケンスを時間内に解くことを必要とする。このような各問題は、ほとんど瞬時に解かれなければならない。高速なダイナミクスを備えたシステムにおける動的インバージョンシステムに対しては、これはほんの数ミリセカンドしか許されない。従って、最適化技術が解を迅速に提示する要件がこれらのアプリケーションにおけるキーポイントになる。
【０００３】
このような問題のシーケンスは、過去においては動的システムを制御するために直接的または間接的に提起された。しかしながら、それらを解くための技術は、アプリケーションが実時間クリティカルではなかったために実時間では動作しないか又は実時間で動作していたが、しかし、下流の石油精製作業のコントロールおよびポートフォリオ最適化におけるように、実時間インターバル(real-time interval)にははるかに余裕があった（恐らく、ミリセカンドではなく時間のオーダーで）。他の技術は実時間で運用されるが、しかしその技術は、解の最適性については厳しい要求がなされず、且つシステムの物理的制約に違反することには注意が払われないようなものであった。この技術の例としては、ある種のフライトコントロールアルゴリズムがある。
【発明の開示】
【課題を解決するための手段】
【０００４】
本発明は、最適化アルゴリズムを提供し、それは、重大な故障(failure)または機能停止(breakdown)を回避するためにシステムの物理的制約を尊重しつつ、動的システムの制御用のオンボードコンピュータに導入される。それは、特に、比較的速いダイナミクスを備えるシステムのコントロールに重要な意味を持ち、それはコントロール変数に対する迅速な解を必要とする。そのアルゴリズムは、以下のように、実時間コントロールアプリケーションに対してそれを特に有用なものとする特別な特徴を備えている。
【０００５】
創意に富んだアルゴリズムは、唯一の最適解が存在し、極小値または極大値が存在せず、且つ、その解が相応な計算で見出し得ることを保証するために、通常は解が凸形式(convex form)に制限され、ときおり二次計画法に制限される。
【０００６】
上記アルゴリズムは、また、実時間凸最適化問題(real time convex optimization problems)を解くためのアクティブセットメソッド(active set method)を使用し、そこでは、前の解(prior solution)が“フィードバック”されて現在の問題において最適なアクティブセットの探索を初期化する。これは、或るアップデートから次のアップデートへのごく小さな変化しか存在しない程の十分に高い周波数でコントロール解が定期的にアップデートされるコントロールシステムにおけるアプリケーションを予見させる。
【０００７】
本発明は、また、どの反復(iteration)も前のものを越えた改善であるような最適なアクティブセットの反復探索が数式化されたアクティブセットメソッドを提供する。これは、必要であればアルゴリズムを中断することができ、または、与えられた反復が予期しない数値の問題のために失敗したとしても、良い解を前の反復から供給することができることを意味する。
【０００８】
上記アルゴリズムは、また、アクティブセットフィードバックと共に、アクティブセットメソッドを使用する。そこでは、最適なアクティブセットの探索における反復の数が、実コントロールソフトウェアにより要求されたときに、予測可能な計算時間を提供するために固定される。数式の線形セットは、制約における縮退(degeneracy)があり得るにもかかわらず、良好な数値特性を提供するために、コレスキー(Cholesky)分解およびＱＲ分解の組み合わせを用いて解かれる。
【０００９】
主要な挑戦は、二次計画法問題を、システムが準静的(quasi-static)と見なされる所定の実時間インターバル(real-time interval)内で解くことである。いくつかの高速なシステムアプリケーションにおいては、この時間インターバルは、動的インバーションにおいてＱＰを解くためにほんの数ミリ秒が割り当てられるとして、たったの１秒の４０分の１である。ＱＰを解くための全ての既存の（非実時間の）アルゴリズムは反復(iterative)のものであり、そこでは、すべての一連の反復が（常にではないが）実際の解についてのより良い推測をもたらす。アクセラレーテッドＱＰ(Accelerated QP)のフレームワークは、明確になるという理由で、内点メソッド(interior point method)とは対照的に、最適性で拘束力のある制約のセット（アクティブセット）を推測することに基づく二次計画法技術によって提供される。このアクセラレーテッドＱＰアルゴリズムにおいて、すべての反復は目的関数(objective function)を単調に改善する。これらの反復のどれも、線形代数計算を含み、そして一定の浮動小数点演算命令実行回数(flops)を獲得し、それは、代わりにＣＰＵを消耗する。常にＱＰを一定の（実）時間インターバル内で解くために、従って一定数の反復内でＱＰを解く必要がある。解がアクティブセットメソッドを用いて見つけられる反復の理論上の最大数が存在する一方、既存のコントロールハードウェアを用いて所定の時間インターバル内でそれらの反復を常に実行することは可能ではない。
【００１０】
従って、アクセラレーテッドＱＰは、可能な限り少ない反復で各ＱＰを解き、反復または計算の数を一定とし、そして各反復に要するＣＰＵ時間が最小となるように反復を計算上効率的とすることに主眼を置く。
【００１１】
最小の反復(minimal iteration)を達成することの背後にあるキーとなる考え：
各反復では、アクティブセットＱＰアルゴリズムは、どの制約が解に拘束力があるのかを推測する。この推測は、反復で拘束力(binding)があると評価された制約（制約のワーキングセットと呼ばれる）のラグランジェ乗数(Lagrange multipliers)の符号(sign)に基づいている。負の符号のラグランジェ乗数は、対応する制約の内部で移動可能であり、そしてコスト関数(cost function)の値を改善することが可能であることを示している。従来のＱＰアルゴリズムは、ワーキングセットから最も負の乗数を備える制約(the most negative Lagrange multiplier)のみを落とす一方、アクセラレーテッドＱＰは、ワーキングセットから負の乗数を備える全ての制約を落とす。これは、従来のメソッドが結局は一つずつ制約を落とすことになり、さらなる反復を招くためである。負の乗数を備える全ての制約を落とすことは、次のＱＰ反復を自由にして大きな部分空間において移動することを可能にし、従ってより少ない反復で解への収束を一層速める。このことは重要な手がかりである。なぜならば、演算を一定時間インターバル内で完了しなければならない実時間セッティング(real-time setting)においては、ＱＰは一定数の反復の後に停止されるからである。
【００１２】
収束に必要とされるＱＰ反復数を最小化することにおける他の重要な考えは、‘ホットスタート(hot start)’の考えである。二次計画法のシーケンスは、発展する動的システムによって要求されるように、時間内で解かれる。ＱＰにおける制約は、システムの物理的状態と関連づけられるので、解に直接的に拘束力のある制約は、システムを制限する物理的要素に対応する。これらの制限要素は、或る時間インターバルから次のインターバルにかけて、少ししか変化しないか又は全く変化しないので、現在の時間ステップにおけるＱＰからの拘束力のある制約のセットは、優秀な、そして大抵は的確な、次の時間ステップにおけるＱＰのための最適性で拘束力のある制約のセットの予測材料(predictor)である（たとえＱＰを定義するパラメータが異なっているとしても）。この重要な知見は、我々が少ない反復で最適に全てのＱＰを解くことを可能とし、それは、動的インバージョンに対するMatlabシミュレーションについてワーストケースで３−５反復である。
【００１３】
もし、短い実時間インターバルが、最適性に至るＱＰの解を可能としなければ、それはこのようにして得られた最良の解とは依然として程遠い状態で終わる。さて、最適に解かれないＱＰは次の時間ステップにおけるＱＰとほとんど類似し、非最適な解についてのワーキングセットは、次のＱＰについての最適な拘束力のあるセットに極めて近い。従って時間ステップでのＱＰに対する解は最適ではないかもしれないにしても、結局のところ、スターティング推測(starting guess)はとても良好であるから、時間ステップは、ＱＰが最適に解かれるときに到達される。この方法では、‘ホットスタート(hot start)’スキームは、アルゴリズムが、およそ少ない時間ステップの間に実時間制限によって引き起こされる‘非最適性から回復’することを可能にする。動的インバージョンについてのコンピュータシミュレーションにおいては、最適性は、たった一つの余剰の時間ステップにおいて達成されることが分かる。
【００１４】
アクティブセットメソッドは、当然に、‘ホットスタート’情報を有効に利用することができ、それは、内点メソッド(Interior Point Method)とは対照的に、それがアクセラレーテッドＱＰの土台を形成する理由である。
【００１５】
各ＱＰ反復をスピードアップすることの背後にあるキーとなる考え：
アクティブセット探索の各反復をスピードアップすることの背後にあるキーとなる考えは、来たるセクションにおいて詳細に述べられるように、有効な線形代数および性能プログラミングメソッド(performance programming method)を基にしている。
【００１６】
本発明の他の効果は、本明細書における添付の図面と関連づけて検討すれば、同様のことが次の詳細な説明を参照することにより一層良く理解されるにつれ、容易に理解されるであろう。
ここで、図１は、本発明の多変数コントロールシステムを例示する。図２は、制限に到達しない場合の図１のシステムの応答のグラフを示す。図３は、制限に到達した場合の本システムの応答のグラフを示す。図４は、図１のＱＰソルバー(QP solver)の計算負荷を示す。
【発明を実施するための最良の形態】
【００１７】
図１は、本発明の実時間二次計画アルゴリズム(real-time quadratic programming algorithm)によって動的インバージョンを実施するコントロールシステム２０を例示する。本コントロールシステム２０は、ここでは、簡潔に要約されており、本特許出願と同一の発明者により同日に出願された“制約付き動的インバージョンコントロールアルゴリズム(Constrained Dynamic Inversion Control Algorithm)”と題された現在継続中の米国シリアル番号＿において更に詳細に述べられており、そしてそれは、本実施形態においてそっくりそのまま参照することにより組み込まれる。制約付き動的インバージョンアルゴリズムは、動的システムを制御するために離散時間型サンプル値コンピュータ(discrete time sampled-data computer)２２上で実施されるように意図されている。このコンピュータ２２は、ＣＰＵ、メモリ２４および他の公知のハードウェアを具備し、且つ本実施形態で説明されるように機能してアクチュエータコマンド(actuator command)ｕｃを生成するようにプログラムされている。このアクチュエータコマンドｕｃは、複数のアクチュエータ２６を制御し、そしてそれは電子的コマンドｕｃを受け取って、それらを、複数のエフェクタ２８を調節する物理的力に変え、複数のアクチュエータ２６のそれぞれは、複数のエフェクタ２８の一つに関連づけられている。エフェクタ２８は、プラント３０のダイナミクスを変える。複数のセンサ３２は、コンピュータ２２にフィードバックｙを供給する。コンピュータ２２は、センサ３２からの情報に基づき、且つシステム出力コマンド３４に基づき、アクチュエータコマンドの次のセットを決定し、それは人間のオペレータまたは高次のコントロールシステムから到来するものであってもよい。
【００１８】
コンピュータ２２は、アクチュエータコマンドｕｃを作成するために、コマンドｙｃおよびフィードバックｙを定期的にサンプリングする。ｙ（ｎ）およびｙｃ（ｎ）は、時間ｎで生じるｎ番目のサンプリングでのｙおよびｙｃの値を指す。アクチュエータ２６は、ｕｃ（ｎ）を、時間ｎから時間ｎ＋１にかけて保持される一定のコマンドとして解釈する。ｙ、ｙｃ、およびｕは、全て、複数の自由度を表すベクトルであってもよい。本実施形態で用いられるように、現在の時間は時間ｎであり、且つ次のコントロールサンプルの時間は時間ｎ＋１である。コントロールアルゴリズムは、コントロールソフトウェアによって実施される。それは、命令された値ｙｃにシステム出力ｙを一致させるアクチュエータコマンドｕを決定する。このコマンドは、人間のオペレータまたは高次のコントロールシステムから到来するものであってもよい。これは、コントロールシステムに基づくデジタルコンピュータの良く知られた形式である。本発明の要旨は、メモリ２４に格納され且つコンピュータ２２により実行される新規な制御にある。
【００１９】
コントロールロジック（ソフトウェア）アーキテクチャは、コントロールロジックデザインおよび実施を、動的フィードフォワード（ＤＦＦ; Dynamic Feed-Forward）４０、動的モデルフォロア（ＤＭＦ; Dynamic Model Follower）４２、および動的インバージョン（ＤＩ; Dynamic Inversion）４４に分割する。ＤＦＦ４０モジュールは、次のコントロールサンプル期間までにシステム出力がどの程度変化すべきかのコマンドのコンポーネントを生成する。そのデザインは、予定通り(nominal)のコマンド応答ダイナミクス(command response dynamics)を確立する。ＤＩ４４モジュールは、要求どおりにシステム出力を最適に変化させるエフェクタ２８のコマンドを決定する。ＤＩ４４は、従って、それらを、簡単で分離したデジタル積分器として現すことにより、システムダイナミクスを“インバート(invert)”させる。システムダイナミクスが構成または動作点によって変化する場合でさえ、ＤＩ４４を備えるダイナミクスは、システムの動作エンベロープ(operating envelope)の全体にわたって同一であり、且つ一つのシステムから他のシステムにかけて同一であるように見える。ＤＭＦ４２モジュールは、所望のシステム応答および実際のシステム応答を比較して、モデリングエラーおよび障害の影響を低減する必要に応じて、出力変調コマンド(output change command)を調整する。ＤＭＦ４２は、モデルエラーおよび障害の影響がどのようにして迅速且つ完全に低減されるかを定める。
【００２０】
ＤＩ４４モジュールは、図１に示すように、更に四つのモジュールに再分割され、次のように全部で六つのモジュールになる。
状態変数モデル(state variable model)４８は、アクチュエータコマンドｕｃ（ｎ）に応じて、ひとつ先（別名、次のステップ）のコントロールサンプルであるシステム変数ｙ（ｎ＋１）の値を予測するために線形化したモデルを生成する。この予測は、また、現在の出力ｙ（ｎ）のフィードバック量(feedback measure)を使用する。
性能指数ブロック５０は、このモデルを、線形制約(linear constraints)に関する二次最適化問題のように、アクチュエータコマンドｕｃ（ｎ）の決定を定式化するために使用する。このブロックに対する一つの入力は、コマンドｖ（ｎ）であり、それは、システム出力ｙが次の時間ｎ＋１までに変化すべき量である。性能指数は、どの程度良好に達成されるかの二次量(quadratic measure)を含む。
制約ブロック５２は、エフェクタに関する制限(limits)および制約(constraints)、およびシステム応答の種々の変数に関する制限および制約を、独立変数としてｕｃ（ｎ）を有する線形不等式の形式に変換する。
二次計画法ブロック５４は、性能指数および制約式によって釣り合いがとられた最適化問題の特殊なタイプを解く。二次計画法ブロック５４は、本発明の中心である。
【００２１】
最適化ベース(optimization-based)のコントローラは多くの利点を提供する。第一に、それらは、多変数システムにおいて、制約を実施するための一般的な理にかなった方法を提供する。制限付き最適化(constrained optimization)の一般的な計算(mathematics)は、制約が満たされるようにシステムの振る舞いを変化させる。第二に、それらは、柔軟性がある。モデル、制約、および性能指数は、コントローラサンプル時間(controller sample times)の間に変化することができる。コントロールアクチュエータ故障および他のシステム故障は、モデルに対する変更(changes to the model)および課された制約(imposed constraints)としてコントローラに反映されてもよい。これらが一旦変わると、最適化ベースのコントローラは、機能的制御(functional control)を用いてできる最良の制御を行おうと試みる。最後に、最適化ベース(optimization-based)のコントローラは、より高レベルの自律性で作動する。それらは、様々なそして変わりやすい順序(priorities)で、多数の、ことによると相反するゴール(goal)の形式でコマンドを受信することができる。最適化はゴール間のトレードオフを可能にする。これらの長所を提供することは本来正確な最適解を必要としないということに注意を払うことは重要である。それどころか、最適化は、より高レベルのアブストラクト(abstract)でコントロールソフトウェアをデザインするための数学的フレームワークとして使用されているので、より柔軟性があり且つ再使用可能である。
【００２２】
コントローラは、ハード(hard)な実時間で自動的に作動する。このため、関連する最適化ソフトウェアは、高度に信頼できるものでなければならず、且つ一定の時間インターバル内にその計算を完了しなければならない。制御理論上の理由から、このインターバルは、一般に、システム動作に関する最速の時間定数(fastest time constant)の２０分の１から１０分の１である。例えば、高性能の航空機は、実質的に、０．３秒以内の突然の突風(sudden gust)の影響を低減することを必要とする。このことは、コントローラが、およそ２５ｍｓでその計算を完了する必要があることを意味する。従って、最適化ソフトウェアは、計算効率が良く、実際的に実行可能であることを必要とする。
【００２３】
実時間が制約された最適化に対する高レベルの戦略は、最適化を凸二次計画法（ＱＰ）問題に見せかけるのに必要な近似(approximation)を行うことである。どの凸ＱＰ問題も、ユニークな世界的な最適条件(global optimum)を有している。ＱＰの確かな探索手順を解くためのアルゴリズムは予測可能な数値特性を含んでいる。ＱＰ問題は、二次指標(quadratic indices)を用いて性能を測定すること、制限を部分的に線形不等制約として現すこと、および部分的に、もしかするとパラメータが変化する、線形差分方程式を用いてシステムをモデル化することにより形成され得る。これらは、ほとんどのアプリケーションで緩やかな近似である。
【００２４】
以下の説明は、浮動少数点演算について一定の割当量(budget)が与えられた、時間における各点でのこのようなＱＰの事例の数式および解に焦点を当てる。明らかに、（浮動少数点演算の数が、ワーキングセット(working set)に対する組み合わせの可能性(combinatorial possibilities)の数に等しい反復数を可能とする程度に大きくなければ）アルゴリズムが一定数の浮動少数点演算(floating point operation)内で解に到達することを保証することは困難である。従って、許容時間において、または浮動少数点演算（ｆｌｏｐ）の割当量内で極めて良好な解に到達することが特別に重要となる。
【００２５】
問題の記述(Problem Description)
コントロール変数ｕ（ｔ）、即ち
【００２６】
【数１】
2004288161.tif 000003
と、状態変数ξ（ｔ）、即ち
【００２７】
【数２】
2004288161.tif 000004
とを有する動的システムを考える。システムのダイナミクスについて線形モデルを考える。
【００２８】
【数３】
2004288161.tif 000005
ここで、ｕｃはコントローラに発行されたコントロールコマンドを表す。システムｙ（ｔ）、即ち
【００２９】
【数４】
2004288161.tif 000006
の応答は、状態変数およびコントロール変数の一次結合(linear combination)として定義される。
【００３０】
【数５】
2004288161.tif 000007
【００３１】
目的は、可能な限り最も良好にｄｙ／ｄｔが所望のベクトルｖと一致するようなｕｃ（ｔ）を見つけ出すことである。このことは、次のように提示される。
【００３２】
【数６】
2004288161.tif 000008
ここで、‖・‖ｗは、重み付けされたｌ２ノルム(norm)を表し、しばしばエリプチィックノルム(elliptic norm)として参照される。ここで、Ｗは、正の重みの対角マトリックスである。上記線形システムは、次のように（Δｔ＝１として）、ディスクリタイズ(discretized)できる。
【００３３】
【数７】
2004288161.tif 000009
【００３４】
コントローラダイナミクスが、ｕｃ＝ｕｔとなるように、“速い(fast)”、即ちＢ２＝−Ｉ且つＭ＝Ｉと仮定する。これは、主として計算を簡単化する図示目的のためになされ、且つ取り扱うのにもはや全く困難でない。２に１を代入すると、次式を生じる。
【００３５】
【数８】
2004288161.tif 000010
時間点ｔで、目的は、それから次のように最小化することである。
【００３６】
【数９】
2004288161.tif 000011
または、３を代入して次式を得る。
【００３７】
【数１０】
2004288161.tif 000012
【００３８】
時間ｔでの上述の目的関数(objective function)において、時間ｔ−１からの量が既に知られており、且つｕｔは、目的が最適化された変数であることに注意されたい。これは、一つの時間ユニットよりも長さが大きい時間ウィンドウを使用するモデル予測制御(Model predictive Control)とは対照的である。
目的関数４は、ｕｔにおいて凸二次的(convex quadratic)である。加えて、結び付けられた制約と、おそらくはＱＰの解が結果としてシステムの故障原因となるアクチュエータコマンドに終わることを妨げるｙｔ，Δｙｔ，Δｕｔを含む他の線形制約が存在する。これらは、ｕｔにおける線形不等制約(linear inequality constraint)として書き直すことができ、そして数理計画法の数式に加えられる。
【００３９】
たまたま、検討中のシステムは、しばしば極度の状態で動作されるものである。従って、これらの不等性制約を考慮に入れることは重要である。これは、実際に、不等性制約のないＱＰに対する類似解は、極度の状態下でのそれらの制約に違反し、そしてシステムを危機的故障に導くからである。
システムデザインは、応答(response)よりも多くのアクチュエータ、従ってコントロール変数を有し、従ってマトリックスＤは過小決定(under-determined)される。
変数の自由度は、ｕｔについて所望の値ｕｄを満足させることにおいて第二の目的を達成するために使用される。この第二の目的は次のように表される。
【００４０】
【数１１】
2004288161.tif 000013
【００４１】
４における目的(objective)は、結合した目的(combined objective)を得るため、上述の第二の目的でペナルティを課される(penalized)。微分可能性(differentiability)のために二つのパートを二乗して加算すると、結合した目的関数は、次のようになる。
【００４２】
【数１２】
2004288161.tif 000014
ここで、ρは第二の目的の相対的重要性を表す小さなスカラーであり、代表的には１０−３である。
ペナルティ項としての第二の目的(secondary objective)の加算は、第一項について規則化する効果を有する。それは、ｕｔにおいて二次の項における対角要素の大きさを増加し、且つアイゲン値(eigenvalue)を的確に正にする。従ってＱＰは的確に凸ＱＰ(convex QP)になる。
【００４３】
第二目的の役割
冗長なコントローラは、次のような第二のゴールを達成するために利用できる付加的な自由を提供する。
アクチュエータ故障のイベントにおいてシステムの堅牢性を増加する。
より低い温度でエンジンを作動させ、従ってエンジンの寿命を延ばす。
コントロール変数レベルを予め決定された最良の設定の近くに維持する。
エンジンの燃料好率を増加する。
素早い操作(maneuver)についての余裕を増加する。例えば、素早く応答するアクチュエータは、遅いアクチュエータよりも、そのバウンド(bound)から離してセットできる。そして中心のより速いアクチュエータは、要求の多い敏速な操作にすばやく応答するために使用できる。
面倒な表記法を避けるため、最適化変数をｕｔに代えてｘと名前を付け変え、縮小された問題は、単に不等制約を用いて次の厳格な凸ＱＰの形式にすることができる。
【００４４】
【数１３】
2004288161.tif 000015
チャレンジは、堅牢な方法で許可された実時間インターバル内に良好な解に到達することに存在する。
更に、Ｈ，ｃ，Ａおよびｂの全ては、一つの時間ステップから次の時間ステップで変わるが、添え字ｔは、表記法を簡単にするために落とされる（ここでは特別な時間点でのＱＰのみを考慮する）。
【００４５】
凸ＱＰの解決
凸ＱＰを解くための技術の二つの幅広いカテゴリーは、どのように不等制約がアルゴリズムにおいて取り扱われるかという点において異なる。
アクティブセットメソッド(Active Set Method)
二つのうちの古いものは、アクティブセットメソッドであり、そこでは、不等式(inequality)のうちのいくつかは、どの反復でも等式(equality)として取り扱われ、そして次のステップは、これらの等式が違反されないようにとられる。等式として保持され、または拘束力(binding)があり、またはアクティブ(active)である不等式のセットは、アクティブセット(active set)として参照される。全ての不等式制約に関する実行可能性は、どの反復でも維持される。線形計画法についてのシンプレックスメソッド(Simplex Method)は、アクティブセットメソッドの良く知られた例である。理論上、アクティブセットアルゴリズムは、結局は有界凸多面体(polytope)の全ての頂点を訪れること、即ち、制約の全ての可能性のある組み合わせをくまなく探索することになる（線形計画法についての病理学のクレーミンティ(Klee Minty)の例により立証されたように）。しかしながら、実際には収束は極めて速い。
【００４６】
内点(interior point)アプローチ
他の選択肢は、内点メソッドであり、それは、ＱＰの最適性についてのかき乱された一次の必要条件(perturbed first-order necessary condition)にニュートンのメソッドを適用することと考えることができる。反復(iterates)は、どの反復(iteration)でも完全に実行可能のままである（即ち、どの不等式も拘束力があるようになることが許されていない）。凸ＱＰについて、反復(iterate)は、収束の最中にセントラルパス(central path)として知られているものに極めて近い状態であり得た。多項式時間収束(polynomial time convergence)は、内点メソッドについて立証することができる。即ち、それらは、標準シンプレックスメソッドがクレーミンティの例に関してなすように、最悪条件の指数関数時間挙動(the worst case exponential time behavior)を示さない。内点メソッドは、多くの不等式制約に関する大きな問題を取り扱うのに大変都合よくできている。
【００４７】
アクティブセット法がより適切である理由
この場合、ある洞察が物理的システムによって提供される。このシステムは、極めて頻繁に極度の状態下で動作し、そしてこのシステムが極度の状態の或るセットから極度の状態の他のセットに速やかにジャンプするとは予期されない。極度の状態が、代わりに、どの不等式が解で拘束力を有するのかを支配するので、特定の時間ステップでＱＰの解に結び付けられた不等式のセットは、次の時間ステップにおいてＱＰについて最適に、とても良好な、時には完璧なアクティブセットの評価(estimate)を提供する。次の時間ステップにおけるＱＰを初期化するために現在の時間ステップで解を使用できるので、これは、問題に対する内点メソッドを越えてアクティブセットメソッドを選択する良い理由である。
【００４８】
問題についてのアクティブセットメソッド
初期時間点ｔ０で、ＱＰは、アクティブセットメソッドを用いて最適に解かれる。任意の時間点ｔ＞ｔ０では、浮動少数点演算命令実行回数の限られた割当量(budget)が存在し、その割当量は、結果として、多くても、実行され得る反復の最大数になる。このようなＱＰの解決のための手順を以下に略述する。
【００４９】
初期実行可能なスタート点の発見
アクティブセットメソッドが制約に関して実行可能な反復を生じる事実は、開始すべき実行可能なスタート点を有することに基づいている。一般のＱＰについては、実行可能な点を見つけるための“フェーズ１”の線形問題を解く必要がある。幸いなことに、検討中の物理的システムは、実行可能なセット内でｘ＝０という仮定を許容する。ｘ＝０は、現在の時間点において新たなコントロールコマンドを何ら発行しないことに対応し、従って前の時間点の状態を維持することに対応し、そしてそれは最適ではないが、少なくとも制約を満足する。
【００５０】
後に続く反復(subsequent iteration)
ワーキングセットを形成する。アクティブセットは、現在の反復で等式のようにその右手側を満足する不等式制約のセットである一方、ワーキングセットは、アクティブ制約のサブセットである。このサブセットは、厳格なサブセットである必要はない。アルゴリズムについて、最初のｎのアクティブ制約はワーキングセットに含まれる。ここで、ｎは変数の数である。実際には、アクティブ制約の数は、ｎよりも十分に小さく、故に全てのアクティブ制約はワーキングセットの中にある。
ｘｋ、現在のｋ番目の反復が与えられると、等式制約付きのＱＰを解くことにより次のステップｓｋを決定する。
【００５１】
【数１４】
2004288161.tif 000016
ここで、Ｅｋは、その行(row)がワーキングセットに制約を含むＡの行に対応するマトリックスであり、一方、ｒｋのコンポーネントは、ｂから右手側のエントリーの対応するサブセットを含む。プレーンワーク(plain work)において、この上述の問題は、次の点ｘｋ＋ｓｋにとって、二次目的(quadratic objective)に関して最適であるように見える一方、ワーキングセットにおける制約を等式として持続させる。しかしながら、ｘｋ＋ｓｋは、可能領域の外にあってもよい。
ｘｋは、実行可能(feasible)であり、且つＥｋは、前の反復においてアクティブであった行を含むので、Ｅｋｘｋ＝ｒｋである。これを上式に代入し、そして目的に改めて目を向ければ、次のステップを決定するための次の等式制約付きＱＰが導かれる。
【００５２】
【数１５】
2004288161.tif 000017
【００５３】
Ｈの正の確定性(definiteness)は、我々がクローズド形式(closed form)で（８）の唯一の全体的ミニマイザ(unique global minimizer)を表現することを可能にするが、Ｅｋがランク不十分(rank-deficient)である場合の事例に特別な注意を払う必要がある。これは、次のセクションにおいて詳しく述べられるであろう。Ｈが反転可能(invertible)であると仮定すると（それは、厳格に正の重さを有する第二の目標の付加によって保証される）、（８）に対するクローズド形式の解は次のようになる。
【００５４】
【数１６】
2004288161.tif 000018
ここで、λｋは、制約Ｅｋｓｋ＝０のラグランジェ乗数を表す。
従ってこの形式の解の実施は、我々が、反復(iteration)の始まりで、Ｈ（その反転の代用）のコレスキー(Cholesky)分解を前計算(pre-compute)し、そしてどの反復(iteration)でもこの既に計算された分解を再使用することを可能にし、従ってあらゆる反復で元の線形システムを解くことに比較して、反復ごとの全体の計算コストを低減させる。従って、この手順は、Ｈがあらゆる反復で極めて大きな線形システムの解を避けることにおいて反転可能(invertible)であるという事実を利用している。
【００５５】
ワーキングセットから乗数λｋの負のコンポーネントに対応する制約を落とす。負の値の乗数は、不等式制約の厳格な内部で移動し且つ目的の値を改善することが可能であるということを指す。故に、それらの不等式制約は最適条件で拘束力があり得ない。
ステップｓｋが与えられると、ｘｋ＋αｓｋが実行可能となるように、α∈[0.1]の最小値を決定する。換言すると、そのステップは、アップデートされた点が実行可能であるようにスケールバック(scale back)される。これは、アクティブ制約の新たなセットをもたらす。
現在の反復をアップデートする。
【００５６】
【数１７】
2004288161.tif 000019
もし、ワーキングセットから制約が落とされなかったならば、即ち全ての乗数が負でなければ、且つもし何ら新たな制約がアクティブにならなければ、即ちステップが可能領域内にあり且つα＝１であるならば、現在の点ｚｋ＋１は最適である。
もしそうでなければ、ワーキングセットを選択するために初期ステップに戻る。
【００５７】
上述の反復が有限量の時間で全体的な最適値に収束することは証明できる一方、解は、許された時間フレーム内に供給されなければならない。従って、上述したプリセットされた反復の最大数のみの実行が可能となり、それは、プロセッサが、許された実時間インターバル内でどの位多くの浮動小数点演算命令実行回数（ｆｌｏｐｓ）を処理することができるかに依存する。また浮動小数点演算の割当量(budget)の一部は、可能領域内に置くためにスケールバック(scale back)された制約されない目的の全体的最適値を計算するためにとっておかれる。この点ｘ＊は次のように計算される。
【００５８】
【数１８】
2004288161.tif 000020
ここで、αは、ｘ＊が可能領域内に存在するような[0,1]において最も大きな値である。このように、ｘ＝０という先の仮定が与えられれば、αの値が常に存在する。
【００５９】
さて、スケール(scale)され制約されない解として参照される上述の解、およびＱＰ反復によって得られた最良点が与えられると、これまでは最良解はリターンされる。ＱＰアルゴリズムは、あらゆる反復で目的(objective)を改善する実行可能な反復を生じる。凸性に基づくこの正式な証明は次のサブセクションにおいて提供される。
【００６０】
ＱＰが時間インターバルのいくつかにおいて適切に解かれないかもしれないということはあり得る。しかしながら、現在の時間インターバルにおける最後の反復が、次の時間インターバルにおけるＱＰに対する良い導きまたはホットスタート(hot start)を提供することが期待され、従って収束は、続く時間インターバルにおいて結局は達成されるであろう。それは、許された時間インターバルにおいてどのＱＰについても収束が達成されるので、必要ではないようである。
ＱＰ反復は、どの反復でも単調に目的関数を改善し、それは厳密に証明できる。
【００６１】
代わりに他の制約の一次結合として表現されていた冗長制約を導入することにより、縮退(degeneracy)が引き起こされる。しばしば、このような依存性の制約は、同時にアクティブになり得、そして最後のハワーキングセットで終わる。これは、結果としてＥｋにおける行ランクの損失になり、そして計算において処理されることを必要とする。ただ一つのアクティブ制約のみをワーキングセットに加え且つ与えられた反復(given iteration)においてただ一つのみを落とすところの線形計画法に対するシンプレックスメソッドにおいて、縮退は、アクティブセットのいくつかの選択の中から‘サイクリング(cycling)’をもたらし、そしてブランド(Bland)のルールのようなシンプルなルールで調整され得る。しかしながら、アクティブセットは、このＱＰアルゴリズムについて別に構成され、そしてサイクリングは、ここで取り扱う課題ではない。それにもかかわらず、Ｅｋがランクを失う場合に（９）において乗数の計算を取り扱うことは必要ではなく、これについて、以下に概略を述べる。
【００６２】
根本的な考え方は、ＥｋＨ−１ＥｋＴのＱＲ分解を使用することである。非線形計画法においては、ＥＴのＱＲ分解を実施し、そしてＥのゼロ空間(null space)の基盤を作り出すための分解を使用し、且つそれをλおよびｓｋの計算において使用することが、より一般的であるので、これは型破りのように見える。しかしながら、このアプローチは、この特別なアルゴリズムに対し、わずかしか浮動小数点命令実行回数を要求しない。
【００６３】
【数１９】
2004288161.tif 000021
【００６４】
マトリックスＢｋの列について解くために先に計算されたＨのコレスキー(Cholesky)分解（Ｌ）を用いてバックソルブ(back solve)およびフォワードソルブ(forward solve)の手順を使用する。
【００６５】
【数２０】
2004288161.tif 000022
ｍ×ｍマトリックスＭｋのＱＲ分解を実施する。ここで、ｍはワーキングセットにおける制約の数である。そうすると、次のようになる。
【００６６】
【数２１】
2004288161.tif 000023
ＱＴｒｋを形成する。Ｑは直交的(orthogonal)であるから、乗数は次式を解くことにより得られる。
【００６７】
【数２２】
2004288161.tif 000024
【００６８】
これは、Ｒは上三角(upper triangular)であるから、Ｏ（ｍ２）オペレーションを含むフォワードソルブのみである。これは、縮退の影響が作用し始める場所である。縮退ワーキングセットは、ゼロになるため依存性制約に対応するＲにおける対角上のエントリーを生じる。そこで、フォワードソルブの期間にゼロピボット(zero pivot)に遭遇したときには、対応するλｋのコンポーネントはその上ゼロに設定される。数学上、乗数の値は、任意の数に設定できる。しかしながら、それらはゼロに設定されて、これらの乗数が現在の下位の問題(subproblem)において役割を果たす冗長制約に相当することを表す。
一旦、λｋが決定されると、ｘにおけるステップは次のように得られる。
【００６９】
【数２３】
2004288161.tif 000025
【００７０】
計算効率についての性能プログラミング(performance programming)の利用
計算効率の最大の利点は、オブジェクト指向のＣ＋＋の実施において、上述の数式における線形代数計算を実施することにおいて得られる。オブジェクトのコピーを作り出すような一般の落とし穴を避けることと共に、アクセスオペレータのようにたびたび呼び出されるファンクションが自動的にインライン化(inlined)され、そしてＣ＋＋の他の特長が使用されて、結果が以下にリストされる効率になる。以下の一式の効率最大化スキームを採用する商業的または別に利用可能な二次計画法についてはコードが存在しないと見られている。
【００７１】
マトリックスのエレメントをアクセスする間、同時にメモリの隣接するブロックをアクセスすることはより効率的である。ＱＰアルゴリズムに含まれる線形代数は、慎重に分析され、そして、いくつかのマトリックスを行方向(row-wise)にアクセスすることが最善であることが決定され、一方、他のものについては、エレメントは列方向(column-wise)にアクセスされなければならない。非プログラム言語で予め定義された二次元アレー(two dimensional array)は、二つのタイプのアクセス（例えば、Ｃにおいては２−Ｄアレーが行方向に保存され、フォートランでは列方向に保存される）を可能にする。アクセラレーテッドＱＰアルゴリズム(accelerated QP algorithm)が、エレメントを別に保存するマトリックスを利用できることは、Ｃ＋＋においてオブジェクトを定義することのみを用いている。エレメントがどのように保存されるかに基づき、カスタムルーチン(custom routine)がまた、ストレージスキーム(storage scheme)を完全に利用するマトリックスベクトル乗算(matrix-vector multiplication)について実施される。また上及び下三角および対称マトリックス(Upper and Lower Triangular and Symmetric Matrix)のような他のマトリックスクラスが実施され、ストレージスキーム(storage scheme)を完全に利用するためアクセレレートＱＰアルゴリズムを可能にする。
【００７２】
ＱＰアルゴリズムはマトリックスを使用し、そのマトリックスの行は、元の問題に対する制約マトリックスの行のサブセットである。Ｃ＋＋でのアクセラレーテッドＱＰの実施は、‘パラサイトマトリックス(parasite matrix)’を可能にし、それはそれ自身にメモリを含まないが、ペアレントマトリックス(parent matrix)の行を簡単に指す。ホストマトリックスに関する全ての線形代数演算の可能なことは、パラサイトマトリックスに関して可能なことである。これは多くのコピーを回避することを助ける。
【００７３】
計算コストの大部分を形成するワーキング制約について乗数を決定するために、すべてのＱＰ反復で線形システムを解くことが必要である。安定性の目的のため、これは、ＱＲ分解を用いて行われる。ここでの重要な効率増進は、Ｑマトリックスがシステムを解くために形成される必要がないことを認識することにある。必要とされることの全ては、Ｑの転置にベクトルを乗算することである。これは、Ｑを全く形成することなしに、ＱＲ分解についてハウスホルダー変換ベクトル(householder transformation vector)を用いることによりなし得る。従って非常に多くの浮動小数点命令実行回数を節約する。ここで気づくべき他の重要なことは、とても頻繁に解かれるシステムの次元は、ほんの１×１または２×２であったと言うことである。なぜなら、システムの次元は、ワーキングセットにおける制約の数に等しいからである。線形システムが、１×１のケースについて平凡に解かれ、そしてさらに、２×２のケースについての行列式に基づきクローズド形式の数式(closed form formula)を用いて直接的に解かれ、これらの二つのケースは、別々に取り扱われ、従ってＱＰ分解を用いるよりも一層少ない命令実行回数で足りる。
【００７４】
上述のアルゴリズムを開発することにおいて、それを安定させ、信頼できるものにし、そして堅牢なものとすることに主要な注意が向けられたことを強調することは重要である。従ってアクセラレーテッドＱＰアルゴリズムは、迅速であるのみならず、とても堅牢である。従来のものとは反対に、ＱＲ分解アプローチは、縮退(degenerate)であるかもしれなケース、即ち、同時に拘束力(biding)を有するようになる複数の制約が線形依存性を有するようなケースを取り扱うために、ラグランジェ乗数を計算するために使用される。この状況は、同じ物理的制限を示す冗長制約によって引き起こされてもよい。本発明が適用された本システムは、堅牢性のために組み込まれたこのような冗長制約を備えてもよい。
【００７５】
更に、アルゴリズムにおいて、ワーキングセットにおける制約に対応する行(rows)についてのいくつかのマトリックスベクトル乗算が回避される。なぜなら、それらがワーキングセットにあるという事実のおかげで、結果が知られているからである。これは、さもなければ徐々に増える小さな（しかし重大な）浮動少数点エラーを回避する一方、追加のボーナスとして浮動小数点命令実行回数を下げる。
【図面の簡単な説明】
【００７６】
【図１】本発明の多変数コントロールシステムを例示する図である。
【図２】制限に到達しない場合の図１のシステムの応答のグラフを示す図である。
【図３】制限に到達した場合の本システムの応答グラフを示す図である。
【図４】図１のＱＰソルバーの計算負荷を示す図である。
【符号の説明】
【００７７】
２０  コントロールシステム
２２  コンピュータ
２４  メモリ
２６  アクチュエータ
２８  エフェクタ
３０  プラント
３２  センサ
３４  出力コマンド
４０  動的フィードフォワード（ＤＦＦ）
４２  動的モデルフォロア（ＤＭＦ）
４４  動的インバージョン（ＤＩ）
４８  状態変数モデル
５０  性能指数ブロック
５２  制約ブロック
５４  二次計画法ブロック
【図１】
2004288161.tif 000026
【図２】
2004288161.tif 000027
【図３】
2004288161.tif 000028
【図４】
2004288161.tif 000029
【外国語明細書】
2004288161.tif 000030
2004288161.tif 000031
2004288161.tif 000032
2004288161.tif 000033
2004288161.tif 000034
2004288161.tif 000035
2004288161.tif 000036
2004288161.tif 000037
2004288161.tif 000038
2004288161.tif 000039
2004288161.tif 000040
2004288161.tif 000041
2004288161.tif 000042
2004288161.tif 000043
2004288161.tif 000044
2004288161.tif 000045
2004288161.tif 000046
2004288161.tif 000047
2004288161.tif 000048
2004288161.tif 000049
2004288161.tif 000050
2004288161.tif 000051
2004288161.tif 000052
2004288161.tif 000053
2004288161.tif 000054
2004288161.tif 000055