JPH04266101A - 多変数モデル予測制御装置 - Google Patents

Abstract

(57)【要約】本公報は電子出願前の出願データであるた
め要約のデータは記録されません。

Description

【発明の詳細な説明】

【０００１】

【産業上の利用分野】本発明は、制御対象の動特性モデ
ルに基づいて制御応答の未来の動きを予測しそれを考慮
しながら操作量を算出するモデル予測制御装置に関する
。

【０００２】

【従来の技術】近年、プロセス制御の分野で、モデル予
測制御装置がしばしば用いられる。モデル予測制御は、
■むだ時間の長いプロセスに対し安定した制御応答を実
現できる。■未来目標値を用いたフィードフォワード制
御で追従性を改善できる。■制御対象の正確な動特性モ
デルを必要とせず、例えばステップ応答から、制御系を
容易に構成できる。■予測モデルにプラントの物理的法
則や非線形特性を含めることにより、きめの細かい制御
が期待できる等の特徴がある。

【０００３】このため、種々の予測制御方式が提案され
ている。これ等の例が、例えば、 （１）　西谷：モデル予測制御の応用、計測と制御ｖｏ
ｌ．２８，Ｎｏ．１１，ｐｐ．９９６−１００４（１９
８９）（２）　Ｄ．Ｗ．　Ｃｌａｒｋｅ　＆　Ｃ．Ｍｏ
ｈｔａｄｉ：　Ｐｒｏｐｅｒｔｉｅｓ　ｏｆ　Ｇｅｎｅ
ｒａｌｉｚｅｄ　Ｐｒｅｄｉｃｔｉｖｅ　Ｃｏｎｔｒｏ
ｌ，　Ａｕｔｏｍａｔｉｃａ　２５−６　ｐｐ．８５９
　（１９８９） に解説されている。特に、（２）　においては、多種の
モデル予測制御方式を包含した一般化予測制御（Ｇｅｎ
ｅｒａｌｉｚｅｄ　Ｐｒｅｄｉｃｔｉｖｅ　Ｃｏｎｔｒ
ｏｌ：　ＧＰＣ）が提案されている。これは、未来目標
値ｙ＊　が与えられたとき、プロセス（制御対象）のモ
デルに基づいて制御応答未来値を予測し、制御要求を表
す評価関数 　　　　　　　　　　Ｎｐ 　　Ｊ　　＝　　Σ｛Ｄ（ｚ−１）（ｙ（ｋ＋ｉ）−ｙ
＊　（ｋ＋ｉ））｝２　　　　　　　　　　　　　ｉ＝
Ｌ　　　 　　　　　　　　　　　　　　　　　　　　　　　　　
　　　Ｎｕ　　　　　　　　　　　　　　　　　　　　
　　　　　　＋λΣ（Δｕ（ｋ＋ｉ−１））２　　　　
　　　　　（１）　　　　　　　　　　　　　　　　　
　　　　　　　　　　　　　ｉ＝１　　　を最小化する
操作量増分Δｕ（ｋ）を求める方式である。

【０００４】ここで、（１）式中のパラメータのうち、
Ｌは予測開始時間、Ｎｐは予測長、Ｎｕは制御長、λは
重み係数、Ｄ（ｚ−１）は極配置多項式を表している。

【０００５】

【発明が解決しようとする課題】上述した従来のモデル
予測制御方式は、例えば、図５に示されるように１入出
力プロセス６に対し、モデル予測制御演算部９と制御定
数算出手段８とが設けられた構成としている。すなわち
、操作量と制御量が共に一つである１入出力系のプロセ
スを対象としている。しかし、化学、鉄鋼、セメント、
抄紙、食品などのプロセスの多くは、温度、圧力、流量
、液位など複数の制御変数が互いに干渉している多変数
系（多入出力系）であり、これらの対象に有効な制御方
式が望まれている。モデル予測制御をこれらの対象に適
用する場合、次の問題点が生じる。

【０００６】第１に、これらの多入出力系では、操作量
と制御量の数が異なる場合がある。従来より提案されて
いる入手の容易な伝達関数を用いた多変数制御系設計法
では、操作量と制御量の数が等しいことを前提としてい
る。例えば、計測制御学会学術図書、線形制御系の設計
理論、第６章（ｐ．１８６〜２２１）にこのような一例
が記載されている。従って、これをモデル予測制御系の
設計にそのまま用いることはできない。よって、操作量
と制御量の数が異なる場合にも適用できるモデル予測制
御方式が必要とされている。

【０００７】第２に、モデル予測制御系では、上述の評
価関数（１）式に含まれる評価関数パラメータＬ、Ｎｐ
、Ｎｕ、λ、Ｄ（ｚ−１）の選び方により制御系の特性
、特に安定性や特性変動に対するロバスト性、いわゆる
ナイキスト判定法における安定性余裕が大きく変わるた
め、制御装置の始動時には上記パラメータを適切に調節
（チューニング）する必要がある。

【０００８】しかし、従来のモデル予測制御装置では評
価関数パラメータと制御特性との関係が明確にされてお
らず、これらの評価関数パラメータは操作員が試行錯誤
を繰り返しながら、経験的に定めていた。このため、制
御系が十分安定になるように調整するのに手間を要し、
制御装置の起動に時間を要している。

【０００９】そこで、本発明のモデル予測制御装置は、
評価関数パラメータが自動的に最適値に設定されるモデ
ル予測制御装置を提供することを目的とする。

【００１０】

【課題を解決するための手段】上記目的を達成するため
本発明のモデル予測制御装置は、図１に示されるように
ｑ本の操作量とｐ本の制御量（ｐ≠ｑ，すなわち操作量
と制御量の数が異なる場合もありうる）を持つ多入出力
プロセス１と、複数の操作量を算出するモデル予測制御
演算部５と、制御定数を算出する制御定数算出手段４と
に加えて、評価関数パラメータを適切な値に調整するた
めの評価関数パラメータ調整手段３を備えている。

【００１１】制御定数算出手段４には、操作量と制御量
の数が異なる多入出力プロセスの動的モデルが、例えば
パルス伝達関数行列の形で与えられる。制御定数算出手
段４は、これを変形して未来の制御量を予測するための
予測式を導出する。また、複数の制御偏差の２乗和と複
数の操作量増分の２乗和の項からなる２次形式評価関数
を与え、評価関数の値が最小となる最適操作量を算出す
るモデル予測制御演算式を導出する。この演算式の制御
パラメータをモデル予測制御演算部５に設定する。

【００１２】評価関数パラメータ調節手段３は、上述の
評価関数パラメータのうち、（１）予測長Ｎｐ，予測開
始時間Ｌは制御対象のステップ応答、またはインパルス
応答に基づき決定し、（２）制御長Ｎｕは制御演算の数
値的安定性を確保するために、制御演算式の導出過程で
得られる行列のランクに基づき決定し、（３）重み係数
λと閉ループ極配置多項式Ｄ（ｚ−１）は、制御系の安
定性、特性変動に対するロバスト性を確保するために、
制御系の一巡伝達関数行列の周波数応答、感度関数、相
補感度関数を算出し、操作員が予め指定した安定余裕、
ロバスト安定余裕に関する条件を調べながら、それを満
たすような値に決定する。この原理を図２及び図３を用
いて以下に説明する。

【００１３】図１に示されるモデル予測制御系は、多入
出力プロセス１とモデル予測制御演算部５からなる閉ル
ープ系を有する。この閉ループ系の特性が制御応答特性
や安定性を支配する。多入出力プロセスの周波数応答行
列をＧ（ｊω）、モデル予測制御演算部の周波数応答行
列をＣ（ｊω）、単位行列Ｉ、夫々の不安定極の数をＮ
Ｇ　個、ＮＣ　個とすると、一般化ナイキスト定理によ
れば、 ｄｅｔ（Ｉ＋Ｇ（ｊω）Ｃ（ｊω）） または、 ｄｅｔ（Ｉ＋Ｃ（ｊω）Ｇ（ｊω）） の周波数ω＝−∞〜∞に対する複素平面上のベクトル軌
跡が図２のように原点を通らずに、原点に対し反時計回
りにＮＧ　＋ＮＣ　回だけ回転すれば閉ループ系は安定
である。この性質を利用すれば、一度安定になるように
設計したモデル予測制御系で制御定数を連続的に変化さ
せる場合、上述のベクトル軌跡が原点からある距離εｍ
ｉｎ　を保つ限り、すなわち ｄｅｔ（Ｉ＋Ｃ（ｊω）Ｇ（ｊω））＞０である限り、
閉ループ系の安定性は保たれる。この条件は、行列式　
ｄｅｔ（…）の代わりに最小および最大特異値に関する
条件 σｍｉｎ　（Ｉ＋Ｃ（ｊω）Ｇ（ｊω））＞０または、 σｍａｘ　（Ｉ＋Ｃ（ｊω）Ｇ（ｊω））−１＜∞に置
き換えられる。ここに、σｍｉｎ　及びσｍａｘ　は夫
々最小及び最大特異値を抽出すること表す。

【００１４】そこで、安定余裕の尺度として上記εｍｉ
ｎ　の逆数に相当するＳｍａｘ　＜∞を与えておき、σ
ｍａｘ　（Ｉ＋Ｃ（ｊω）Ｇ（ｊω））−１≦Ｓｍａｘ
　が成り立つように制御定数を調整すれば、制御系の安
定性が保たれる。これは、結局、図３（Ａ）に示される
感度関数（Ｉ＋Ｃ（ｊω）Ｇ（ｊω））−１の最大特異
値の周波数応答曲線のピーク値を閾値以内に抑えること
に相当する。

【００１５】次に、制御対象の特性変動に対するロバス
ト安定性を保つための条件を示す。制御対象Ｇ（ｊω）
がＧ（ｊω）（１＋Δ（ｊω））に変化した場合、相補
感度関数 　　　　　　Ｔ（ｊω）＝（Ｉ＋Ｃ（ｊω）Ｇ（ｊω）
）−１Ｃ（ｊω）Ｇ（ｊω）に対し、 σｍａｘ　｛Ｔ（ｊω）｝　　≦　　Ｔｍａｘ　（ω）
ただし、 Ｔｍａｘ　（ω）　　＝　　（σｍａｘ　｛Δ（ｊω）
｝）−１が成り立つならば、制御系の安定性は確保され
る。そこで、あらかじめ予想される制御対象の特性変動
に対応した関数Ｔｍａｘ　（ω）を与えておき、図３（
Ｂ）に示すように上述の条件が成立すれば、制御系のロ
バスト安定性が確保できる。

【００１６】従って、図３に示すような感度関数、相補
感度関数に関する条件をチェックしながら、制御定数、
あるいは、評価関数パラメータを調整する限り、モデル
予測制御系の安定性、ロバスト安定性は確保される。以
上が本発明の提案する評価関数パラメータ調整手段３に
おけるパラメータ調整方法の原理である。

【００１７】こうして、決定された評価関数パラメータ
に基づいてモデル予測制御装置の制御演算式に含まれる
制御定数を決定する。

【００１８】

【作用】制御演算式の導出に用いる未来の制御量を予測
するための予測式は、操作量と制御量の数が異なる場合
も適用できる形式なので、導出される制御演算式も、操
作量と制御量の数が異なる多入出力プロセスに適用する
ことが出来る。

【００１９】また、上記制御演算式を導出するに際して
、参照される評価関数の評価関数パラメータは、評価関
数パラメータ調整手段によって制御系の安定性、対象の
特性変動に対しロバスト安定性を確保するように設定さ
れるので、操作員の試行錯誤に頼らずにパラメータ調整
が出来、モデル予測制御装置を迅速に立ち上げることが
出来る。また、制御系の制御応答特性を良好に保つこと
が出来る。

【００２０】

【実施例】本発明のモデル予測制御装置の実施例を図１
を参照して説明する。図１において制御対象プロセス１
は、ｑ本の操作量とｐ本の制御量を持つ多入出力プロセ
スである。この多入出力プロセス１を制御するモデル予
測制御装置２は、プロセスの制御量ｙ１　〜ｙｐ　を目
標値ｙ＊１　〜ｙ＊ｐ　に追従させるための操作量ｕ１
　〜ｕｑ　を出力する。

【００２１】このモデル予測制御装置２は、与えられた
安定余裕パラメータ等に基づいて適切な評価関数パラメ
ータを決定する評価関数パラメータ調整手段３と、この
評価関数パラメータと、動特性モデルとして与えられる
制御対象のパルス伝達関数行列から制御定数を決定する
制御定数算出手段４と、決定された制御定数に基づき最
適な操作量を逐次算出するモデル予測制御演算部５とか
ら構成される。以下、各機能を説明する。

【００２２】多入出力プロセスに対する動特性モデルと
して次式のパルス伝何時関数行列をモデル予測制御装置
２に与える。

【００２３】

【数２】 ただし、 Ａｉｊ（ｚ−１）＝１＋ａ１　ｉｊＺ−１＋…＋ａｎ　
ｉｊｚ−ｎＢｉｊ（ｚ−１）＝ｂ０　ｉｊ＋ｂ１　ｉｊ
ｚ−１＋…＋ｂｍ　ｉｊｚ−ｍである。ここで、ｐは制
御量の数、ｑは操作量の数で、ｐ≠ｑでもよい。

【００２４】また制御目的は、次式の評価関数が最小に
なるように制御量ｙ１　を目標値ｙ１　＊に追従させる
ことである。

【００２５】 　　　　　　　　Ｎｐ　　ｐ　 　　Ｊ＝　　Σ　　Σ｛Ｄ（ｚ−１）（ｙｉ　（ｋ＋ｊ
）−ｙｉ　＊（ｋ＋ｊ））｝２　　　　　　　　　　ｊ
＝Ｌ　　ｉ＝１　　　　　　　　　　　　　　　　　　
　　　　　Ｎｕ　　ｑ　　　　　　　　　　　　　　　
　　　　＋λΣ　　Σ（Δｕｉ　（ｋ＋ｊ−１））２　
　　　　　　　　　　　　…（３）　　　　　　　　　
　　　　　　　　　　ｊ＝１　　ｉ＝１　　　ただし、
Δｕｉ　（ｋ）＝ｕｉ　（ｋ）−ｕｉ　（ｋ−１）は操
作量増分である。また、（３）式中のパラメータのうち
、ＬおよびＮｐは予測長、Ｎｕは制御長、λは重み係数
、Ｄ（ｚ−１）は極配置多項式で、これ等をまとめて評
価関数パラメータと称する。なお、パラメータＬ、Ｎｐ
、Ｎｕ、λ、Ｄ（ｚ−１）を各操作量ｕ１　〜ｕｑ　、
各制御量ｙ１　〜ｙｐ　毎に変えることが出来るが、こ
の実施例では区別しない場合を示す。

【００２６】まず、制御定数算出手段４における制御定
数の算出手順を示す。制御定数算出手段４は、与えられ
た制御対象モデルのパルス伝達関数行列（２）式と、評
価関数パラメータに基づき、次の計算を行う。

【００２７】（２）式右辺行列中の第１〜ｐ行について
、行の中で分母を通分する。例えば、第ｉ行については
、 　　　　　　Ｂｉ１（ｚ−１）／Ａｉ１（ｚ−１）…Ｂ
ｉｑ（ｚ−１）／Ａｉｑ（ｚ−１）を、 　　　　　　Ｂ′ｉ１（ｚ−１）／Ａｉ　（ｚ−１）…
Ｂｉｑ′（ｚ−１）／Ａｉ　（ｚ−１）のように通分す
る。この結果、（２）式は、

【００２８】

【数３】 の形に変形できる。次に、各行（ｉ＝１〜ｐ）について
、およびｊ＝１〜Ｎｐについて、次の方程式（これをＤ
ｉｏｐｈａｎｔｉｎｅ　方程式という）　　　　Ｄ（ｚ
−１）＝Ｅｊｉ（ｚ−１）（１−ｚ−１）Ａｉ　（ｚ−
１）＋ｚ−ｊＦｊｉ（ｚ−１）　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　…（５）を解き、多項式 　　　　Ｅｊｉ（ｚ−１）＝１＋ｅ１　ｊｉｚ−１＋…
＋ｅｊ−１　ｊｉｚ−ｊ＋１（ｊ−１次モニック）　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　…（６ａ）　　　　Ｆｊｉ（ｚ−
１）＝ｆ０　ｊｉ＋ｆ１　ｊｉｚ−１＋…＋ｆｎ　ｊｉ
ｚ−ｎ（ｎ次）　　…（６ｂ）を求める。ここで、閉ル
ープ多項式は、例えば、　　　　Ｄ（ｚ−１）＝１−ρ
ｚ−１　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　　　　　　　…（７）と与える。 ρはフィルタ係数であり、応答速度を表す。 ρが０に近い値であると応答が急であり、１に近い値で
あると応答は緩やかになる。

【００２９】また、この結果、ｊステップ先の制御量ｙ
１　〜ｙｐの予測値は、

【００３０】

【数４】 で与えられる。ただし、Ｇｊ−１　、Ｇ０　はｐ×ｑの
定数行列、Ｈｊ　（ｚ−１）はｚ−１に関するｐ×ｑの
多項式行列で、次の関係式から求められる。

【００３１】

【数５】 なお、Ｇ０　…Ｇｊ−１　は制御対象（２）式のステッ
プ応答行列に相当する。

【００３２】次に（８）式を次式のように簡略に表現す
ることにする。 　　　　ｙ（ｋ＋ｊ）＝Ｇｊ−１　Δｕ（ｋ＋ｊ−１）
＋…＋Ｇ０　Δｕ（ｋ）　　　　　　　　　　　　＋Ｆ
ｊ　（ｚ−１）ｙ（ｋ）＋Ｈｊ　（ｚ−１）Δｕ（ｋ−
１）　　…（１０）これを、ｊ＝１〜Ｎｐについてまと
めて表現すると、

【００３３】

【数６】 となる。そこで、さらに（１１）をまとめて、　　　　
ｙ＝ＧΔｕ＋Ｆ（ｚ−１）ｙ（ｋ）＋Ｈ（ｚ−１）Δｕ
（ｋ−１）　　　　…（１２）と表現する。このとき、
（３）式の評価関数は、　　Ｊ＝（Ｄ（ｚ−１）ｙ−Ｄ
（ｚ−１）ｙ＊）Ｔ　（Ｄ（ｚ−１）ｙ−Ｄ（ｚ−１）
ｙ＊）　　　　　　　　　　　　　　＋λΔｕＴ　Δｕ
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　　　　…（１３）と表すことができる
。ただし、Ｔ　は転値行列を意味し、

【００３４】

【数７】 で、ｙ＊１　（ｋ＋ｊ）は、第ｉ番目の制御量ｙ１　に
対するｊステップ先の目標値である。

【００３５】このとき、評価関数（３）式（または（１
３）式）を最小化する最適操作量は、 　　　　Δｕ＝（ＧＴ　Ｇ＋λＩ）−１ＧＴ　｛Ｄ（ｚ
−１）ｙ＊−Ｆ（ｚ−１）ｙ（ｋ）　　　　　　　　　
　　　　　−Ｈ（ｚ−１）Δｕ（ｋ−１）｝　　　　　
　　　　　　　　　　　　　　　　　…（１４）で与え
られる。実際の操作量は、ベクトルΔｕのはじめのｑ要
素 　　　　Δｕ（ｋ＋１）＝Ｇｑ＊｛Ｄ（ｚ−１）ｙ＊−
Ｆ（ｚ−１）ｙ（ｋ）　　　　　　　　　　　　　　−
Ｈ（ｚ−１）Δｕ（ｋ−１）｝　　　　　　　　　　　
　　　　　　　　　　　…（１５）（ただし、Ｇｑ＊は
（ＧＴ　Ｇ＋λＩ）−１ＧＴ　の上からｑ行を取りだし
たｑ×（ｐ×Ｎｐ）のマトリクスである。）を取り出し
、 　　　　　　ｕ（ｋ）＝ｕ（ｋ−１）＋Δｕ（ｋ）　　
　　　　　　　　　　　　　　　　　　　　　　…（１
６）により、計算される。

【００３６】従って、制御定数算出手段４は導出した（
１５）式から制御パラメータ（制御定数）としてＧｑ＊
、Ｆ（ｚ−１）、Ｈ（ｚ−１）を抽出し、これ等パラメ
ータをモデル予測制御演算部５に提供する。

【００３７】次に、モデル予測制御演算部５では、制御
定数算出手段４から制御定数Ｇｑ＊、Ｆ（ｚ−１）、Ｈ
（ｚ−１）を受け、プロセスの制御量ｙ１　〜ｙｐ　と
目標値ｙ＊１〜ｙ＊ｐ　を用いて、（１５）、（１６）
式に基づいて制御演算を行う。

【００３８】次に、評価関数パラメータ調整手段３にお
ける、評価関数パラメータ（予測長ＬおよびＮｐ、制御
長Ｎｕ、重み係数λ、極配置多項式Ｄ（ｚ−１））の調
整手段について説明する。

【００３９】（ｉ）予測長Ｎｐは以下の様に決定する。

【００４０】 　　Ｎｐ＝ｍａｘ　（パルス伝達関数行列（２）式ｉｊ
要素のステップ応答最終値　　　　　　　　ｉｊ　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　　　の９０％に達する時間）　　または、 　　Ｎｐ＝ｍａｘ　（パルス伝達関数行列（２）式ｉｊ
要素のインパルス応答　　　　　　　　ｉｊ　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　ピーク値に達する時間）　　（ｉｉ）予測
開始長Ｌは以下の様に決定する。

【００４１】・対象プロセス（２）式のむだ時間がわか
っている場合 　　Ｌ＝ｍａｘ　（パルス伝達関数行列（２）式ｉｊ要
素のむだ時間）　　　　　　ｉｊ 　　　　・対象プロセス（２）式のむだ時間が明かでな
い場合 Ｌ＝１ （ｉｉｉ）制御長Ｎｕは、

【００４２】

【数８】 　　に対し、Ｎｕ＝　ｉｎｔ（ｒａｎｋＧ／ｑ）＋１　
　　　　　　　　　　　　　　　　　　　　　　　　　
…（１８）ただし、Ｎｕ＞Ｎｐの場合はＮｕ＝Ｎｐとす
る。ｉｎｔ　（…）は整数部を取り出す演算である。

【００４３】ｒａｎｋＧの計算は、Ｇを特異値分解し、
　　　　　　Ｇ　　＝　　Ｕ　　ｄｉａｇ｛σ１　…σ
ｋ　…σＮｐ｝ＶＴ　　　　　　　　　　　　　　　　
　　　…（１９）その特異値σから、 　　　　　　　　　　　　　　　　　　　　　　　　　
Ｋ　　　　　　　　Ｎｐ×ｐ　　　　　ｒａｎｋＧ＝　
ｍｉｎ｛ｋ；　　Σσｉ　／　　Σσｉ　≧１−ε｝　
　　　　　　　　　　　　　…（２０）　　　　　　　
　　　　　　　　　　　　　　　　　ｉ＝１　　　　　
　　ｉ＝１　　　　　　　　　　　　　閾値εは、例え
ば０．１とする。 により決定する。上記Ｕ及びＶは特異ベクトルである。

【００４４】（ｉｖ）重み係数λ及び閉ループ極配置多
項式Ｄ（ｚ−１）は図４のフローチャートに従って調整
を行う。以下に手順を示す。

【００４５】ステップ１：まず、λ、Ｄ（ｚ−１）の初
期値を例えば次式のように与える。 　　　　　　λ＝ε　　　　　　　　　　　　　　　　
（　　ε≧０　　、例えばε＝１０４　）　　　　　　
Ｄ（ｚ−１）＝１−ρｚ−１（０≦ρ＜１、例えばρ＝
０．９９９９）ステップ２：この評価関数パラメータを
制御定数算出手段４に渡し、（２）〜（１５）式の手順
で制御定数を決定する。

【００４６】ステップ３：制御対象（２）式と制御則（
１５）式から得られる閉ループ系の特性多項式　　　　
　ｄｅｔ｛Δ（ｚ−１）Ａ（ｚ−１）　　　　　　＋Ｂ
（ｚ−１）（Ｉ＋Ｇｑ＊Ｈ（ｚ−１））−１Ｇｑ＊Ｆ（
ｚ−１）｝＝０　　…（２１）の根を計算し、その実部
が全て負であるか否かを判定する。全てが負であれば、
閉ループ系は安定であるので、次のステップに進む。そ
うでない場合は、ステップ１に戻り、評価関数パラメー
タ初期値を変更する。ただし、

【００４７】

【数９】 である。

【００４８】ステップ４：制御則（１５）式から得られ
るコントローラの特性多項式 　　　　　ｄｅｔ（Ｉ＋Ｇｑ＊Ｈ（ｚ−１））＝０　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　…（２２）の根を計算し、その実部が全て負であ
るか否かを判定する。全てが負であれば、コントローラ
は安定であるので、次のステップに進む。そうでない場
合は、ステップ１に戻り、評価関数パラメータ初期値を
変更する。

【００４９】ステップ５：制御系の一巡伝達関数行列周
波数応答Ｌ（ｊω）を以下の計算で求める。 　　Ｌ（ｊω）＝Δ（ｚ−１）−１（Ｉ＋Ｇｑ＊Ｈ（ｚ
−１））−１　　　　　　　　Ｇｑ＊Ｆ（ｚ−１）Ａ（
ｚ−１）−１Ｂ（ｚ−１）｜ｚ＝ｅｊωτ　　　　　　
…（２３）ただし、τは制御周期である。続いて、感度
関数　　　　　　Ｓ（ｊω）＝（Ｉ＋Ｌ（ｊω））−１
　　　　　　　　　　　　　　　　　　　　　　　　　
　…（２４）を求める。

【００５０】ステップ６：感度関数の最大特異値が以下
の条件を満たすか否かを判定する。

【００５１】 　　　　　　Ｓｍａｘ　−εＳ　≦σｍａｘ　｛Ｓ（ｊ
ω）｝Ｓｍａｘ　　　　　　　　　　　　　　　　　　
　…（２５）ただし、Ｓｍａｘ　、εＳ　は、あらかじ
め与えた、安定余裕に関する閾値である。この条件が成
立すれば、次のステップへ、不成立ならばステップ９へ
進む。

【００５２】ステップ７：相補感度関数　　　　　　Ｔ
（ｊω）＝（Ｉ＋Ｌ（ｊω））−１Ｌ（ｊω）　　　　
　　　　　　　　　　　　…（２６）を求める。

【００５３】ステップ８：相補感度関数の最大特異値が
以下の条件を満たす否かを判定する 。　　　　　　σｍａｘ　｛Ｔ（ｊω）｝　　≦　　Ｔ
ｍａｘ　（ω）　　　　　　　　　　　　　　　　　　
　　…（２７）ただし、Ｔｍａｘ　（ω）は、あらかじ
め与えた、ロバスト安定余裕に関する閾値関数で、たと
えば、　　　　Ｔｍａｘ　（ω）＝｜Ｋ／（１＋Ｔｊω
）ｎ　｜　　　　　　　　　　　　　　　　　　　　　
　…（２８）　　　　　　　　　　　　　　　　　　　
　　　　　　　　　（Ｋ＞０、Ｔ＞０、整数ｎ≧１）　
　等を用いる。この条件が成立すれば、評価関数パラメ
ータλ、Ｄ（ｚ−１）の調整を終了する。不成立ならば
ステップ９へ進む。

【００５４】ステップ９：評価関数パラメータλ、Ｄ（
ｚ−１）＝１−ρｚ−１を以下の要領で調整する。 　　　　　　　　λ←λ×０．９９　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　…（２９）または、 　　　　　　　　ρ←０．９９×ρ　　　　　　　　　
　　　　　　　　　　　　　　　　　　　　　　　　　
　　　　　　…（３０）なお、感度関数及び相補感度関
数に関する安定条件を共に満足する必要は必ずしもなく
、（２９）式、（３０）式のいずれか一方のみ実行する
場合も有り得る。

【００５５】ステップ１０：調整し直した評価関数パラ
メータを制御定数算出手段４に渡し、ステップ２と同様
に　（２）〜（１５）式の手順で制御定数を決定する。

【００５６】以上のようにして重み係数λ及び閉ループ
極配置多項式Ｄ（ｚ−１）は制御系が指定した安定余裕
、ロバスト安定余裕を満たすまで繰り返し調整される。

【００５７】こうして、評価関数パラメータ調整手段３
において、予測開始時間Ｌ、予測長Ｎｐ、制御長Ｎｕ、
重み係数λ、極配置多項式Ｄ（ｚ−１）等の評価関数パ
ラメータが適切な値に調整される。

【００５８】なお、（２９）、（３０）式のパラメータ
調整方法は、一例であり、（２５）、（２７）式の条件
が満たされるようなパラメータを、他の非線形最適化計
算法を用いて探索することも可能である。

【００５９】本発明のモデル予測制御装置を、石油化学
プロセスの蒸留塔温度制御系に適用した例について説明
する。

【００６０】図６に示される蒸留塔温度系の特性として
次のような伝達関数が測定された。ここで、操作量ｕ１
　，ｕ２　はそれぞれ還流量調節弁開度、スチーム流量
調節弁開度で、制御量ｙ１　，ｙ２　はそれぞれ、塔頂
温度差、塔底温度である。

【００６１】

【数１０】 この伝達関数をサンプリング周期２秒で離散化したパル
ス伝達関数に対し、本モデル予測制御装置により、評価
関数パラメータを調節して、温度設定値変更に対する制
御応答を測定した。調整前後の制御系の応答を夫々図７
及び図８に示す。調整前は、評価関数パラメータＬ＝１
、Ｎｐ＝１０、Ｎｕ＝５ λ＝１、Ｄ（ｚ−１）＝１−０．１ｚ−１であり、制御
量ｙ１（蒸留塔の塔頂温度差）の目標値ステップ変化に
対する応答は、オフセットがしばらく残り、整定が遅い
。これに対し、評価関数パラメータ調整手段３を動作さ
せてパラメータを調整した結果、評価関数パラメータは
、 Ｌ＝１、Ｎｐ＝１０、Ｎｕ＝５ λ＝０．０５、Ｄ（ｚ−１）＝１ となり、調整後の制御量ｙ１　（蒸留塔の塔頂温度差）
の目標値ステップ変化に対する応答は、整定時間が速く
、また、制御量ｙ１　とｙ２　の間の干渉の影響も十分
に低減されている。この結果、充分な安定性を保ちなが
ら、速応性、非干渉性に優れた多変数モデル予測制御系
が実現できた。

【００６２】

【発明の効果】以上説明したように本発明のモデル予測
制御装置は、従来、経験的に設定していた評価関数パラ
メータを、評価関数パラメータ調整手段を設けてロバス
ト性を考慮した一定のアルゴリズムを実行して適切な評
価関数パラメータを見出すことにより、評価関数パラメ
ータを適切な値に調整するので、制御装置を迅速に起動
出来る。また、予め指定した安定余裕条件を満たすよう
に調整するので、常に十分な安定性を持ち、プラントの
特性変動に影響されにくいロバストなモデル予測制御系
が得られる。また、操作量と制御量の数が共に複数であ
る多入出力のプラント等にも適用可能なモデル予測制御
装置が得られる。

【図面の簡単な説明】

【図１】本発明の多変数モデル予測制御装置の構成を示
すブロック図。

【図２】多変数制御系の安定性を判別するための一般化
ナイキスト線図。

【図３】制御系の安定性とロバスト性を判別するための
感度特性及び相補感度特性を示す特性図。

【図４】評価関数パラメータ調整手段の処理手順を表す
フローチャート。

【図５】従来のモデル予測制御装置の構成を示すブロッ
ク図。

【図６】蒸留塔温度制御系の構成を示すブロック図。

【図７】従来のモデル予測制御による蒸留塔の温度制御
応答特性を示す図。

【図８】本発明によるモデル予測制御による蒸留塔の温
度制御応答特性を示す図。

【符号の説明】

１　　多入出力プロセス（制御対象） ２　　多変数モデル予測制御装置 ３　　評価関数パラメータ調整手段 ４　　制御定数算出手段 ５　　モデル予測制御演算部 ６　　プロセス（制御対象） ７　　従来形モデル予測制御装置 ８　　制御定数算出手段 ９　　モデル予測制御演算部

Claims (7)

【特許請求の範囲】
1. 【請求項１】操作量と制御量が共に複数である多入出力
   プラントの動特性モデルに基づき複数の制御量未来値ｙ
   ｉ　を予測し、未来目標値ｙｉ　＊との偏差信号および
   操作量増分Δｕｉ　に関する以下の２次形式評価関数を
   最小化するような複数の最適操作量を制御量と未来目標
   値から算出するモデル予測制御装置であって、予測開始
   時間Ｌ　、予測長Ｎｐ、制御長Ｎｕ、重み係数λ、閉ル
   ープの極配置多項式Ｄ（ｚ−１）を評価関数パラメータ
   とする評価関数、 　　　　　　　　Ｎｐ　　ｐ　 　　Ｊ＝　　Σ　　Σ｛Ｄ（ｚ−１）（ｙｉ　（ｋ＋ｊ
   ）−ｙｉ　＊（ｋ＋ｊ））｝２　　　　　　　　　　ｊ
   ＝Ｌ　　ｉ＝１　　　　　　　　　　　　　　　　　　
   　　　　　Ｎｕ　　ｑ　　　　　　　　　　　　　　　
   　　　　＋λΣ　　Σ（Δｕｉ　（ｋ＋ｊ−１））２　
   　　　　　　　　　　　　　　　　　　　　　ｊ＝１　
   　ｉ＝１　　　　　　　　　　　（ただし、Δｕｉ　（
   ｋ）＝ｕｉ　（ｋ）−ｕｉ　（ｋ−１）である）　　に
   おける予測長Ｎｐを制御対象のステップ応答形状または
   インパルス応答形状に基づいて決定することを特徴とす
   る多変数モデル予測制御装置。
2. 【請求項２】前記予測長Ｎｐの決定において、制御対象
   のステップ応答が定常値の９０％に達する時間または制
   御対象のインパルス応答がピークに達する時間に基づい
   て前記予測長Ｎｐとすることを特徴とする請求項１記載
   の多変数モデル予測制御装置。
3. 【請求項３】操作量と制御量が共に複数である多入出力
   プラントの動特性モデルに基づき複数の制御量未来値ｙ
   ｉ　を予測し、未来目標値ｙｉ　＊との偏差信号および
   操作量増分Δｕｉ　に関する以下の２次形式評価関数を
   最小化するような複数の最適操作量を制御量と未来目標
   値から算出するモデル予測制御装置であって、予測開始
   時間Ｌ　、予測長Ｎｐ、制御長Ｎｕ、重み係数λ、閉ル
   ープの極配置多項式Ｄ（ｚ−１）を評価関数パラメータ
   とする評価関数、 　　　　　　　　Ｎｐ　　ｐ　 　　Ｊ＝　　Σ　　Σ｛Ｄ（ｚ−１）（ｙｉ　（ｋ＋ｊ
   ）−ｙｉ　＊（ｋ＋ｊ））｝２　　　　　　　　　　ｊ
   ＝Ｌ　　ｉ＝１　　　　　　　　　　　　　　　　　　
   　　　　　Ｎｕ　　ｑ　　　　　　　　　　　　　　　
   　　　　＋λΣ　　Σ（Δｕｉ　（ｋ＋ｊ−１））２　
   　　　　　　　　　　　　　　　　　　　　　ｊ＝１　
   　ｉ＝１　　　　　　　　　　　（ただし、Δｕｉ　（
   ｋ）＝ｕｉ　（ｋ）−ｕｉ　（ｋ−１）である）　　に
   おける前記予測開始時間Ｌ　は、制御対象の各要素のむ
   だ時間のうち最大のものに選定されることを特徴とする
   多変数モデル予測制御装置。
4. 【請求項４】操作量と制御量が共に複数である多入出力
   プラントの動特性モデルに基づき複数の制御量未来値ｙ
   ｉ　を予測し、未来目標値ｙｉ　＊との偏差信号および
   操作量増分Δｕｉ　に関する以下の２次形式評価関数を
   最小化するような複数の最適操作量を制御量と未来目標
   値から算出するモデル予測制御装置であって、予測開始
   時間Ｌ　、予測長Ｎｐ、制御長Ｎｕ、重み係数λ、閉ル
   ープの極配置多項式Ｄ（ｚ−１）を評価関数パラメータ
   とする評価関数、 　　　　　　　　Ｎｐ　　ｐ　 　　Ｊ＝　　Σ　　Σ｛Ｄ（ｚ−１）（ｙｉ　（ｋ＋ｊ
   ）−ｙｉ　＊（ｋ＋ｊ））｝２　　　　　　　　　　ｊ
   ＝Ｌ　　ｉ＝１　　　　　　　　　　　　　　　　　　
   　　　　　Ｎｕ　　ｑ　　　　　　　　　　　　　　　
   　　　　＋λΣ　　Σ（Δｕｉ　（ｋ＋ｊ−１））２　
   　　　　　　　　　　　　　　　　　　　　　ｊ＝１　
   　ｉ＝１　　　　　　　　　　　（ただし、Δｕｉ　（
   ｋ）＝ｕｉ　（ｋ）−ｕｉ　（ｋ−１）である）　　に
   おける前記制御長Ｎｕは、制御対象のステップ応答行列
   Ｇ０　，Ｇ１　，Ｇ２　，…をブロック要素とする行列
   【数１】 の階数（ｒａｎｋ）に基づいて決定されることを特徴と
   する多変数モデル予測制御装置。
5. 【請求項５】前記制御長Ｎｕの決定法として、制御対象
   のステップ応答行列をブロック要素とする行列Ｇの特異
   値分解 Ｇ　　＝Ｕ　　ｄｉａｇ｛σ１　…σｋ　…σＮｐ｝Ｖ
   Ｔ　（ただし、ｄｉａｇ｛σ１　…σＮｐ｝は対角行列
   、Ｕは特異ベクトルの行列、ＶＴ　は特異ベクトルの転
   値行列である）により求めた特異値σを用いて、しきい
   値ε＜＜１に対し、 により行列Ｇの階数（ｒａｎｋ）を定め、Ｎｕ＝ｉｎｔ
   　（ｒａｎｋＧ／ｑ）＋１（ただし、Ｎｕ＞Ｎｐの場合
   はＮｕ＝Ｎｐ，　ｉｎｔ（…）は整数部を取り出す演算
   を表す）により制御長Ｎｕを決定することを特徴とする
   請求項４記載の多変数モデル予測制御装置。
6. 【請求項６】操作量と制御量が共に複数である多入出力
   プラントの動特性モデルに基づき複数の制御量未来値ｙ
   ｉ　を予測し、未来目標値ｙｉ　＊との偏差信号および
   操作量増分Δｕｉ　に関する以下の２次形式評価関数を
   最小化するような複数の最適操作量を制御量と未来目標
   値から算出するモデル予測制御装置であって、予測開始
   時間Ｌ　、予測長Ｎｐ、制御長Ｎｕ、重み係数λ、閉ル
   ープの極配置多項式Ｄ（ｚ−１）を評価関数パラメータ
   とする評価関数、 　　　　　　　　Ｎｐ　　ｐ　 　　Ｊ＝　　Σ　　Σ｛Ｄ（ｚ−１）（ｙｉ　（ｋ＋ｊ
   ）−ｙｉ　＊（ｋ＋ｊ））｝２　　　　　　　　　　ｊ
   ＝Ｌ　　ｉ＝１　　　　　　　　　　　　　　　　　　
   　　　　　Ｎｕ　　ｑ　　　　　　　　　　　　　　　
   　　　　＋λΣ　　Σ（Δｕｉ　（ｋ＋ｊ−１））２　
   　　　　　　　　　　　　　　　　　　　　　ｊ＝１　
   　ｉ＝１　　　　　　　　　　　（ただし、Δｕｉ　（
   ｋ）＝ｕｉ　（ｋ）−ｕｉ　（ｋ−１）である）　　に
   おける前記重み係数λまたは前記閉ループの極配置多項
   式Ｄ（ｚ−１）は、感度特性における感度関数Ｓ（ｊω
   ）の許容最大値Ｓｍａｘ　、安定余裕パラメータεＳ　
   に対し感度関数Ｓ（ｊω）の最大特異値σｍａｘ　｛Ｓ
   （ｊω）｝が、 Ｓｍａｘ　−εＳ　≦σｍａｘ　｛Ｓ（ｊω）｝≦Ｓｍ
   ａｘ　を満足するように定められることを特徴とする多
   変数モデル予測制御装置。
7. 【請求項７】操作量と制御量が共に複数である多入出力
   プラントの動特性モデルに基づき複数の制御量未来値ｙ
   ｉ　を予測し、未来目標値ｙｉ　＊との偏差信号および
   操作量増分Δｕｉ　に関する以下の２次形式評価関数を
   最小化するような複数の最適操作量を制御量と未来目標
   値から算出するモデル予測制御装置であって、予測開始
   時間Ｌ　、予測長Ｎｐ、制御長Ｎｕ、重み係数λ、閉ル
   ープの極配置多項式Ｄ（ｚ−１）を評価関数パラメータ
   とする評価関数、 における前記重み係数λまたは前記閉ループの極配置多
   項式Ｄ（ｚ−１）は、相補感度特性における相補感度関
   数Ｔ（ｊω）の許容最大値Ｔｍａｘ　（ω）に対し、相
   補感度関数Ｔ（ｊω）の最大特異値σｍａｘ　｛Ｔ（ｊ
   ω）｝が、σｍａｘ　｛Ｔ（ｊω）｝≦Ｔｍａｘ　（ω
   ）を満足するように定められることを特徴とする多変数
   モデル予測制御装置。