特開 2024-126120 制御装置、制御方法、及びプログラム

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】公開特許公報(A)

(11)【公開番号】P2024126120

(43)【公開日】2024-09-20

(54)【発明の名称】制御装置、制御方法、及びプログラム

(51)【国際特許分類】

   G05B 13/02 20060101AFI20240912BHJP

【ＦＩ】

G05B13/02 Z

【審査請求】有

【請求項の数】9

【出願形態】ＯＬ

(21)【出願番号】P 2023034303

(22)【出願日】2023-03-07

(11)【特許番号】

(45)【特許公報発行日】2023-05-30

(71)【出願人】

【識別番号】000005234

【氏名又は名称】富士電機株式会社

(74)【代理人】

【識別番号】100107766

【弁理士】

【氏名又は名称】伊東 忠重

(74)【代理人】

【識別番号】100070150

【弁理士】

【氏名又は名称】伊東 忠彦

(72)【発明者】

【氏名】丹下 吉雄

【テーマコード（参考）】

5H004

【Ｆターム（参考）】

5H004GA05

5H004GA07

5H004GB01

5H004HA07

5H004HB01

5H004KA33

5H004KA35

(57)【要約】

【課題】制御性能の低下を抑制できる技術を提供すること。
【解決手段】本開示の一態様による制御装置は、制御対象に対する操作量を出力し、前記制御対象の制御量を目標値に追従させる制御装置であって、第１の制御周期Ｔ_ｃ毎に、前記操作量と前記制御量とに基づいて、モデル予測制御により前記制御対象に対する次の操作量ｕ_ａを出力する主制御部と、前記第１の制御周期Ｔ_ｃよりも短い第２の制御周期Ｔ_ｆ毎に、前記制御対象に対する次の操作量ｕ_ｂを出力する高速補完制御部と、を有し、前記高速補完制御部は、前記主制御部から前記操作量ｕ_ａが出力された直後の時刻ｔ_０における前記制御量と、所定の先読み長Ｔ_ｐと、前記目標値とに基づいて、前記先読み長Ｔ_ｐ後における目標偏差の予測値である予測目標偏差ｅ_ｆを計算し、前記予測目標偏差ｅ_ｆと所定のゲインｋ_Ｉとに基づいて計算した瞬時的な操作変化量の平均と前記操作量ｕ_ａとの和を前記操作量ｕ_ｂとして計算する。
【選択図】図１

【特許請求の範囲】

【請求項1】

制御対象に対する操作量を出力し、前記制御対象の制御量を目標値に追従させる制御装置であって、
第１の制御周期Ｔ_ｃ毎に、前記操作量と前記制御量とに基づいて、モデル予測制御により前記制御対象に対する次の操作量ｕ_ａを出力する主制御部と、
前記第１の制御周期Ｔ_ｃよりも短い第２の制御周期Ｔ_ｆ毎に、前記制御対象に対する次の操作量ｕ_ｂを出力する高速補完制御部と、を有し、
前記高速補完制御部は、
前記主制御部から前記操作量ｕ_ａが出力された直後の時刻ｔ_０における前記制御量と、所定の先読み長Ｔ_ｐと、前記目標値とに基づいて、前記先読み長Ｔ_ｐ後における目標偏差の予測値である予測目標偏差ｅ_ｆを計算し、
前記予測目標偏差ｅ_ｆと所定のゲインｋ_Ｉとに基づいて計算した瞬時的な操作変化量の平均と前記操作量ｕ_ａとの和を前記操作量ｕ_ｂとして計算する、制御装置。

【請求項2】

前記高速補完制御部は、
前記予測目標偏差ｅ_ｆに対して前記ゲインｋ_Ｉを乗じることで、前記瞬時的な操作変化量を計算する、請求項１に記載の制御装置。

【請求項3】

前記先読み長Ｔ_ｐと前記ゲインｋ_Ｉは、
前記モデル予測制御に用いられる予測モデルに基づいて計算される、請求項２に記載の制御装置。

【請求項4】

前記高速補完制御部は、
前記時刻ｔ_０における前記制御量から前記先読み長Ｔ_ｐ後に前記目標値に到達するような線形な制御軌道を用いて、前記目標値と、現在時刻ｔにおける前記制御量から前記制御軌道により予測される時刻ｔ_０＋Ｔ_ｐの前記制御量との差を前記予測目標偏差ｅ_ｆとして計算する、請求項１乃至３の何れか一項に記載の制御装置。

【請求項5】

前記高速補完制御部は、
前記時刻ｔ_０における前記制御量から現在時刻ｔにおける前記制御量までの線形な制御軌道により予測される時刻ｔ_０＋Ｔ_ｐの前記制御量と、前記目標値と差を前記予測目標偏差ｅ_ｆとして計算する、請求項１乃至３の何れか一項に記載の制御装置。

【請求項6】

前記高速補完制御部は、
前記時刻ｔ_０における初期値を０として、前記第２の制御周期Ｔ_ｆ毎にインクリメントされるカウンタｋを用いて、前記カウンタｋに関する前記予測目標偏差ｅ_ｆ（ｋ）と前記瞬時的な操作変化量の平均とを計算する、請求項１に記載の制御装置。

【請求項7】

前記高速補完制御部は、
前記主制御部から出力される前記操作量ｕ_ａに適用される上下限範囲を用いて、前記上下限範囲を満たす前記操作量ｕ_ｂを計算する、請求項１に記載の制御装置。

【請求項8】

制御対象に対する操作量を出力し、前記制御対象の制御量を目標値に追従させる制御装置が、
第１の制御周期Ｔ_ｃ毎に、前記操作量と前記制御量とに基づいて、モデル予測制御により前記制御対象に対する次の操作量ｕ_ａを出力する主制御手順と、
前記第１の制御周期Ｔ_ｃよりも短い第２の制御周期Ｔ_ｆ毎に、前記制御対象に対する次の操作量ｕ_ｂを出力する高速補完制御手順と、を実行し、
前記高速補完制御手順は、
前記主制御手順で前記操作量ｕ_ａが出力された直後の時刻ｔ_０における前記制御量と、所定の先読み長Ｔ_ｐと、前記目標値とに基づいて、前記先読み長Ｔ_ｐ後における目標偏差の予測値である予測目標偏差ｅ_ｆを計算し、
前記予測目標偏差ｅ_ｆと所定のゲインｋ_Ｉとに基づいて計算した瞬時的な操作変化量の平均と前記操作量ｕ_ａとの和を前記操作量ｕ_ｂとして計算する、制御方法。

【請求項9】

制御対象に対する操作量を出力し、前記制御対象の制御量を目標値に追従させる制御装置に、
第１の制御周期Ｔ_ｃ毎に、前記操作量と前記制御量とに基づいて、モデル予測制御により前記制御対象に対する次の操作量ｕ_ａを出力する主制御手順と、
前記第１の制御周期Ｔ_ｃよりも短い第２の制御周期Ｔ_ｆ毎に、前記制御対象に対する次の操作量ｕ_ｂを出力する高速補完制御手順と、を実行させ、
前記高速補完制御手順は、
前記主制御手順で前記操作量ｕ_ａが出力された直後の時刻ｔ_０における前記制御量と、所定の先読み長Ｔ_ｐと、前記目標値とに基づいて、前記先読み長Ｔ_ｐ後における目標偏差の予測値である予測目標偏差ｅ_ｆを計算し、
前記予測目標偏差ｅ_ｆと所定のゲインｋ_Ｉとに基づいて計算した瞬時的な操作変化量の平均と前記操作量ｕ_ａとの和を前記操作量ｕ_ｂとして計算する、プログラム。

【発明の詳細な説明】

【技術分野】

【0001】

本開示は、制御装置、制御方法、及びプログラムに関する。

【背景技術】

【0002】

制御対象の制御量を目標値に追従させることを目的とする制御方式の１つとして、モデル予測制御が知られている。例えば、ＰＬＣ（Programmable Logic Controller）やＤＣＳ（Distributed Control System）等のエッジデバイス上でモデル予測制御を実現しつつ、モデル予測制御に必要な制御パラメータも自動調整する技術が提案されている（特許文献１）。

【先行技術文献】

【特許文献】

【0003】

【特許文献1】特許第７０１４３３０号公報

【発明の概要】

【発明が解決しようとする課題】

【0004】

しかしながら、モデル予測制御の制御周期は、制御対象の時定数に応じて適切に定める必要がある。予測の効果を十分に発揮するためには時定数に応じて或る程度長い制御周期を定める必要があるが、長過ぎる制御周期は制御性能が低下することがあるためである。

【0005】

本開示は、上記の点に鑑みてなされたもので、制御性能の低下を抑制できる技術を提供することを目的とする。

【課題を解決するための手段】

【0006】

本開示の一態様による制御装置は、制御対象に対する操作量を出力し、前記制御対象の制御量を目標値に追従させる制御装置であって、第１の制御周期Ｔ_ｃ毎に、前記操作量と前記制御量とに基づいて、モデル予測制御により前記制御対象に対する次の操作量ｕ_ａを出力する主制御部と、前記第１の制御周期Ｔ_ｃよりも短い第２の制御周期Ｔ_ｆ毎に、前記制御対象に対する次の操作量ｕ_ｂを出力する高速補完制御部と、を有し、前記高速補完制御部は、前記主制御部から前記操作量ｕ_ａが出力された直後の時刻ｔ_０における前記制御量と、所定の先読み長Ｔ_ｐと、前記目標値とに基づいて、前記先読み長Ｔ_ｐ後における目標偏差の予測値である予測目標偏差ｅ_ｆを計算し、前記予測目標偏差ｅ_ｆと所定のゲインｋ_Ｉとに基づいて計算した瞬時的な操作変化量の平均と前記操作量ｕ_ａとの和を前記操作量ｕ_ｂとして計算する。

【発明の効果】

【0007】

制御性能の低下を抑制できる技術が提供される。

【図面の簡単な説明】

【0008】

【図1】第一の実施形態に係る制御装置の全体構成の一例を示す図である。

【図2】プラント応答関数の動作の一例を説明するための図である。

【図3】プラント応答関数の計算処理の一例を説明するためのフローチャートである。

【図4】モデルパラメータの推定処理の一例を説明するためのフローチャートである。

【図5】先読み応答補正部の動作の一例を説明するための図である。

【図6】操作変化量計算部の動作の一例を説明するための図である。

【図7】制御パラメータの計算処理の一例を説明するためのフローチャートである。

【図8】先読み長の計算概要の一例を説明するための図である。

【図9】先読み長の計算処理の一例を説明するためのフローチャートである。

【図10】高速補完制御部の動作の一例を説明するための図である。

【図11】補完操作量の計算処理の一例を説明するためのフローチャートである。

【図12】予測目標偏差の一例を説明するための図（その１）である。

【図13】予測目標偏差の一例を説明するための図（その２）である。

【図14】高速補完制御部を用いた場合の応答の一例を説明するための図である。

【図15】第一の実施形態に係る制御装置のハードウェア構成の一例を示す図である。

【図16】第二の実施形態に係る制御装置の全体構成の一例を示す図である。

【図17】実施例におけるプラント応答を示す図である。

【図18】実施例において高速補完制御が無効である場合の制御量、目標値及び操作量を示す図である。

【図19】実施例において高速補完制御が有効である場合の制御量、目標値及び操作量を示す図（その１）である。

【図20】実施例において高速補完制御が有効である場合の制御量、目標値及び操作量を示す図（その２）である。

【発明を実施するための形態】

【0009】

以下、本発明の一実施形態について、図面を参照しながら詳細に説明する。以下では、モデル予測制御の制御周期が長過ぎることによる制御性能の低下を抑制できる制御装置１０について説明する。なお、以下で説明する制御装置１０は、例えば、ＰＣ（パーソナルコンピュータ）等と比較して計算資源が乏しいエッジデバイス（ＰＬＣやＤＣＳ等）であることを想定する。

【0010】

ここで、一般に、時定数が長い制御対象に対してあまりに短い制御周期を定めた場合、予測の効果が十分に発揮できないことがある。また、多くのメモリを消費してしまったり、モデルが正しく同定できなかったりする傾向がある。他方で、長過ぎる制御周期を定めた場合、制御タイミングの間に生じた目標値（制御量に対する目標値）の変更や外乱等の影響を制御に反映できず、制御性能が低下してしまうことがある。このため、制御対象の時定数に応じて適切に長い制御周期を定める必要があるものの、制御周期が長過ぎることによって制御性能が低下してしまうことがあるというジレンマが存在する。

【0011】

そこで、以下で説明する制御装置１０では、モデル予測制御の制御器（後述する主制御部１０４）に加えて、モデル予測制御の制御周期Ｔ_ｃ（以下、「主制御周期Ｔ_ｃ」ともいう。）よりも短い制御周期Ｔ_ｆ（以下、「高速制御周期Ｔ_ｆ」ともいう。）で制御対象を高速に制御する制御器（後述する高速補完制御部１０５）を設けた構成としている。これにより、モデル予測制御の主制御周期Ｔ_ｃが長い場合であっても、高速制御周期Ｔ_ｆ毎に制御対象が高速に制御されるため、例えば、制御量に対する目標値の変更や変動の早い外乱の影響等にも素早く追従が可能となり、制御性能の低下を抑制することができる。また、制御装置１０の計算資源が乏しく、主制御周期Ｔ_ｃを長くせざるを得ない場合であっても、制御性能の低下を抑制することができる（言い換えれば、制御装置１０の計算資源が乏しい場合であっても、制御性能の低下を抑えつつ、長い主制御周期Ｔ_ｃを定めたモデル予測制御を行うことができる。更に、モデル予測制御では抑制が不十分な制御偏差（例えば、オーバーシュート等による制御偏差）の改善も期待できる。なお、制御量の例としては制御対象の温度、目標値の例としては設定温度等が挙げられる。ただし、これらはいずれも一例であって、制御量及び目標値は、温度及び設定温度に限られるものではなく、制御対象における任意の制御量及び当該制御量の目標となる目標値を用いることが可能である。

【0012】

以下では、制御対象としてプラントを想定し、制御装置１０が制御対象プラント２０を制御する場合について説明する。ただし、制御対象はプラントに限られるものではなく、任意の機器、装置、設備等を制御対象とすることが可能である。

【0013】

［第一の実施形態］
以下、第一の実施形態について説明する。

【0014】

＜制御装置１０の全体構成例＞
まず、第一の実施形態に係る制御装置１０の全体構成例について、図１を参照しながら説明する。図１は、第一の実施形態に係る制御装置１０の全体構成の一例を示す図である。

【0015】

図１に示すように、第一の実施形態に係る制御装置１０は、モデルパラメータ推定部１０１と、計測部１０２と、差分器１０３と、主制御部１０４と、高速補完制御部１０５と、第１のタイマ１０６と、第２のタイマ１０７と、切替器１０８とを有する。これら各部は、例えば、制御装置１０にインストールされた１以上のプログラムが、ＣＰＵ（Central Processing Unit）やＭＰＵ（Micro-Processing Unit）等のプロセッサに実行させる処理により実現される。

【0016】

第１のタイマ１０６は、主制御周期Ｔ_ｃ毎に、モデルパラメータ推定部１０１及び主制御部１０４の動作トリガーとして動作する。なお、主制御周期Ｔ_ｃとは、モデル予測制御により制御対象プラント２０を制御する周期のことであり、その値は予め設定される。

【0017】

第２のタイマ１０７は、高速制御周期Ｔ_ｆ毎に、計測部１０２及び高速補完制御部１０５の動作トリガーとして動作する。なお、高速制御周期Ｔ_ｆとは、ｎを予め決められた任意の自然数としてｎＴ_ｆ＝Ｔ_ｃと表される周期であり、高速補完制御部１０５により制御対象プラント２０を制御する周期のことである。高速制御周期Ｔ_ｆの値も主制御周期Ｔ_ｃと同様に予め設定される。

【0018】

計測部１０２は、高速制御周期Ｔ_ｆ毎に、制御対象プラント２０の制御量ｙを計測（観測）し、計測した制御量ｙを出力する。なお、制御対象プラント２０の制御量ｙは、操作量ｕと外乱ｖとに応じて決定される。外乱ｖの例としては、制御量ｙが温度である場合には外気温の低下又は上昇等が挙げられる。

【0019】

また、計測部１０２は、高速制御周期Ｔ_ｆ毎に、切替器１０８から出力された操作量ｕを取得（観測）し、取得した操作量ｕを出力する。なお、操作量ｕの例としては、制御量ｙが温度である場合に熱媒体の流量を制御するバルブの開閉量等が挙げられる。

【0020】

以下、制御量ｙが計測された時刻ｔを明示する場合、ｙ（ｔ）と表す。同様に、操作量ｕが取得された時刻ｔを明示する場合、ｕ（ｔ）と表す。その他の変数についても同様に、時刻ｔを明示する場合はその変数の後に「（ｔ）」を表記するものとする。

【0021】

差分器１０３は、目標値ｒと、制御量ｙとの差（偏差）を目標偏差ｅ_０として出力する。時刻ｔにおける目標偏差ｅ_０（ｔ）は、ｅ_０（ｔ）＝ｒ（ｔ）－ｙ（ｔ）で計算される。なお、以下では、一例として、目標値は一定、すなわちｒ（ｔ）＝定数であるものとする。ただし、これは一例であって、目標値ｒ（ｔ）は時刻ｔに応じて異なり得るものであってもよい。

【0022】

主制御部１０４は、主制御周期Ｔ_ｃ毎に、モデル予測制御により制御対象プラント２０を制御する。すなわち、主制御部１０４は、主制御周期Ｔ_ｃ毎に、制御対象プラント２０のプラント応答モデルを表すプラント応答関数｛Ｓ_θ（ｔ）｝を用いて、当該制御対象プラント２０に対する操作量ｕ_ａを計算し、出力する。以下では、一例として、主制御部更新カウンタと呼ばれるカウンタがメモリ内に記憶されており、操作量ｕ_ａが計算及び出力される毎に、当該主制御部更新カウンタが更新（例えば、１が加算）されるものとする。なお、プラント応答関数｛Ｓ_θ（ｔ）｝には、モデルパラメータ推定部１０１によって推定対象となるモデルパラメータθが含まれる。

【0023】

ここで、主制御部１０４には、制御パラメータ計算部１１１と、先読み応答補正部１１２と、操作変化量計算部１１３と、加算器１１４とが含まれる。

【0024】

制御パラメータ計算部１１１は、プラント応答関数｛Ｓ_θ（ｔ）｝と、調整係数α，βとを入力し、制御パラメータとして制御ゲインｋ_Ｉと先読み長Ｔ_ｐとを計算し、操作変化量計算部１１３と先読み応答補正部１１２とにそれぞれ出力する。ただし、制御パラメータ計算部１１１は、主制御部１０４に含まれていなくてもよい。

【0025】

先読み応答補正部１１２は、プラント応答関数｛Ｓ_θ（ｔ）｝と、目標偏差ｅ_０（ｔ）と、過去の操作量ｕ_ａの変化量ｄｕの時系列データである操作変化量時系列｛ｄｕ（ｔ）｝と、先読み長Ｔ_ｐとに基づいて、目標偏差ｅ_０（ｔ）を補正した補正目標偏差ｅ^＊（ｔ）を計算する。

【0026】

操作変化量計算部１１３は、補正目標偏差ｅ^＊（ｔ）と制御ゲインｋ_Ｉとに基づいて、操作変化量ｄｕ（ｔ）を計算する。操作変化量計算部１１３は、例えば、ｄｕ（ｔ－３Ｔ_ｃ）、ｄｕ（ｔ－２Ｔ_ｃ）、ｄｕ（ｔ－Ｔ_ｃ）という順序で操作変化量ｄｕ（ｔ）を算出し、出力する。なお、操作変化量ｄｕは、主制御周期Ｔ_ｃ毎に操作量ｕ_ａが変化した量である。

【0027】

加算器１１４は、計測部１０２から出力された操作量ｕと、操作変化量計算部１１３から取得された操作変化量ｄｕとを加算して、新たな操作量ｕ_ａを算出する。すなわち、加算器１１４は、任意の整数ｎ_１に対してｔ＝ｎ_１Ｔ_ｃを満たす時刻ｔ毎に、新たな操作量ｕ_ａ（ｔ）＝ｕ（ｔ－Ｔ_ｃ）＋ｄｕ（ｔ）を算出する。そして、加算器１１４は、この新たな操作量ｕ_ａを切替器１０８に出力する。

【0028】

モデルパラメータ推定部１０１は、主制御周期Ｔ_ｃ毎に、制御量ｙと操作量ｕとを入力し、プラント応答関数｛Ｓ_θ（ｔ）｝のモデルパラメータの推定値θ_ｅｓｔを算出し、出力する。以下では、プラント応答関数｛Ｓ_θ（ｔ）｝に設定されているモデルパラメータθを「モデルパラメータ設定値θ」ともいう。なお、モデルパラメータの推定値θ_ｅｓｔを算出する際に、モデルパラメータの初期値θ_０が入力されてもよい。

【0029】

高速補完制御部１０５は、高速制御周期Ｔ_ｆ毎に、制御対象プラント２０を制御する。すなわち、高速補完制御部１０５は、高速制御周期Ｔ_ｆ毎に、制御対象プラント２０に対する操作量ｕ_ｂを計算し、出力する。ただし、或る整数ｎ_２，ｎ_３に対してｔ＝ｎ_２Ｔ_ｆかつｔ＝ｎ_３Ｔ_ｃとなる時刻ｔに関しては、ｕ_ｂ（ｔ）＝ｕ_ａ（ｔ）となる。

【0030】

切替器１０８は、操作量ｕ_ｂ又はｕ_ａとｕ_ｂの両方を入力し、操作量ｕを制御対象プラント２０に出力する。すなわち、切替器１０８は、任意の整数をｎ_４，ｎ_５として、ｔ≠ｎ_４Ｔ_ｃかつｔ＝ｎ_５Ｔ_ｆを満たす時刻ｔ毎に操作量ｕ_ｂを入力して操作量ｕ＝ｕ_ｂを出力し、ｔ＝ｎ_４Ｔ_ｃかつｔ＝ｎ_５Ｔ_ｆを満たす時刻ｔ毎に操作量ｕ_ａ及びｕ_ｂを入力して操作量ｕ＝ｕ_ａ＝ｕ_ｂを出力する。

【0031】

なお、第一の実施形態に係る制御装置１０は、例えば、主制御周期Ｔ_ｃ毎に計算されたモデルパラメータや制御パラメータ（制御ゲイン、先読み長）をディスプレイ等の表示装置に逐次表示する表示部を有していてもよいし、これらのモデルパラメータや制御パラメータを通信ネットワークにより接続されるＰＣ等といった端末の表示装置上に逐次表示させる表示制御部を有していてもよい。これにより、制御装置１０のユーザや管理者（以降、「ユーザ等」とも表す。）は、モデルパラメータや制御パラメータのトレンドを知ることが可能となり、プラント応答モデルの同定具合や制御パラメータの調整具合を確認することができるようになる。また、プラント応答モデルに急激な変動が起きる場合にその兆候を発見することも可能となる。

【0032】

＜プラント応答関数｛Ｓ_θ（ｔ）｝の動作＞
次に、プラント応答関数｛Ｓ_θ（ｔ）｝の動作について、図２を参照しながら説明する。図２は、プラント応答関数｛Ｓ_θ（ｔ）｝の動作の一例を説明するための図である。

【0033】

図２に示すように、プラント応答関数｛Ｓ_θ（ｔ）｝は、モデルパラメータ設定値θと時刻ｔとが入力されると、初期時刻０から時間ｔ経過後の時刻ｔにおける単位ステップ応答Ｓ_θ（ｔ）を出力する。なお、単位ステップ応答とは、操作量を単位ステップ入力とした場合における応答（つまり、制御対象プラント２０のプラント応答モデルの出力）のことである。

【0034】

＜プラント応答関数｛Ｓ_θ（ｔ）｝の計算＞
次に、プラント応答関数｛Ｓ_θ（ｔ）｝の単位ステップ応答Ｓ_θ（ｔ）を計算する処理について、図３を参照しながら説明する。図３は、プラント応答関数｛Ｓ_θ（ｔ）｝の計算処理の一例を説明するためのフローチャートである。なお、図３では、モデルパラメータ設定値θと時刻ｔとが入力されたものとする。

【0035】

ここで、本実施形態では、プラント応答関数Ｓ_θ（ｔ）の計算モデルとして、制御量ｙについては過去Ｎ点の自己回帰、操作量ｕについては現在値と過去Ｍ点の移動平均を用いたＡＲＭＡ（autoregressive moving average）モデルを採用した場合について説明する。なお、Ｎ及びＭは、例えば、ユーザ等によって予め設定される。ただし、これは一例であって、例えば、ＡＲＭＡＸモデル等の他のモデルを採用してもよい。

【0036】

このとき、モデルパラメータ設定値θを

【数1】

とする。モデルパラメータ設定値θの各要素は、以下の数２に示すＡＲＭＡモデルの係数である。

【数2】

ここで、ｋはインデックスである。ただし、これは一例であって、例えば、過去Ｌ点（Ｌは予め設定された１以上の整数）の外乱ｖを考慮したモデルが用いられてもよい。

【0037】

ステップＳ１１：主制御部１０４は、時刻を表す変数をτ、インデックスｋにおける状態ベクトルをφ（ｋ）として、時刻τと状態ベクトルφ（０）とを初期化する。また、主制御部１０４は、インデックスｋをｋ＝０に初期化する。

【0038】

ここで、状態ベクトルφ（ｋ）は、

【数3】

と表される。

【0039】

なお、例えば、過去Ｌ点の外乱ｖを考慮する場合には、状態ベクトルφ（ｋ）を、ｙ（ｋ）及びｕ（ｋ）に加えてＬ個の外乱ｖ（ｋ）を要素として持つ状態ベクトルに拡大すればよい。

【0040】

主制御部１０４は、例えば、τ＝０と初期化すると共に、

【数4】

と初期化する。すなわち、φ（０）はｕ（ｋ）に相当する要素のみ１、他の要素は０に初期化される。これは、単位ステップ応答を模擬する際の初期値として設定していることを意味する。

【0041】

ステップＳ１２：次に、主制御部１０４は、モデルパラメータ設定値θと状態ベクトルφ（ｋ）とに基づいて制御量予測値ｙ（ｋ）を計算する。主制御部１０４は、例えば、ｙ（ｋ）＝φ（ｋ）^Τθにより制御量予測値ｙ（ｋ）を計算する。ここで、Τは転置を表す。

【0042】

ステップＳ１３：次に、主制御部１０４は、制御量予測値ｙ（ｋ）と状態ベクトルφ（ｋ）とに基づいて、次のインデックス（つまり、ｋ＋１）における状態ベクトルφ（ｋ＋１）を更新する。主制御部１０４は、例えば、以下により状態ベクトルφ（ｋ＋１）を更新する。

【数5】

ここで、ｙ（ｋ）には上記のステップＳ１２で計算した制御量予測値を設定し、ｕ（ｋ＋１）には１を設定する。また、ｙ（ｋ－１），・・・，ｙ（ｋ－Ｎ＋１）及びｕ（ｋ），・・・，ｕ（ｋ－Ｍ＋１）には、状態ベクトルφ（ｋ）と同じ値を設定する（すなわち、状態ベクトルφ（ｋ）に含まれるｙ（ｋ－１），・・・，ｙ（ｋ－Ｎ＋１）及びｕ（ｋ），・・・，ｕ（ｋ－Ｍ＋１）をそれぞれ設定する。）。

【0043】

ステップＳ１４：次に、主制御部１０４は、時刻τをτ＋ΔＴに更新すると共に、インデックスｋをｋ＋１に更新する。ここで、ΔＴはＡＲＭＡモデルの１ステップの時間幅を表す。プラント応答関数の時間刻みΔＴは、制御対象プラント２０の時定数に応じて任意に設定することができるが、例えば、主制御周期Ｔ_ｃと同一にすることが好適である。

【0044】

ステップＳ１５：次に、主制御部１０４は、τ≧ｔであるか否かを判定する。そして、τ≧ｔであると判定されなかった場合（ステップＳ１５でＮＯ）、主制御部１０４は、ステップＳ１２に戻る。これにより、τ≧ｔとなるまで、ステップＳ１２～ステップＳ１４が繰り返し実行される。

【0045】

一方で、τ≧ｔであると判定された場合（ステップＳ１５でＹＥＳ）、主制御部１０４は、処理を終了する。これにより、最終的に計算されたｙ（ｋ）が単位ステップ応答Ｓ_θ（ｔ）として得られる（つまり、Ｓ_θ（ｔ）＝ｙ（ｋ）が、ステップ応答関数｛Ｓ_θ（ｔ）｝の時刻ｔにおける単位ステップ応答として計算される。）。

【0046】

＜モデルパラメータθの推定＞
次に、プラント応答関数｛Ｓ_θ（ｔ）｝のモデルパラメータθを推定する処理（つまり、モデルパラメータ推定値θ_ｅｓｔを算出する処理）について、図４を参照しながら説明する。図４は、モデルパラメータθの推定処理の一例を説明するためのフローチャートである。図４では、現在時刻ｔにおける制御量ｙ（ｔ）と操作量ｕ（ｔ）とが入力されたものとして、或るインデックスｋにおけるモデルパラメータ推定値θ_ｅｓｔを算出する場合について説明する。なお、このインデックスｋは、図３でプラント応答関数の計算に用いたインデックスｋとは独立した値であり、モデルパラメータ推定処理の実行時インデックスを表す。

【0047】

ステップＳ２１：モデルパラメータ推定部１０１は、モデルパラメータ推定値θ_ｅｓｔ（ｋ）と共分散行列Ｐ（ｋ）とを初期化するか否かを判定する。ここで、初期化すると判定される場合としては、例えば、モデルパラメータ推定値θ_ｅｓｔを初めて算出する場合（つまり、モデルパラメータ推定値θ_ｅｓｔの初回算出時）、ユーザ等により初期化指示が行われた場合等が挙げられる。

【0048】

モデルパラメータ推定値θ_ｅｓｔ（ｋ）と共分散行列Ｐ（ｋ）とを初期化すると判定した場合（ステップＳ２１でＹＥＳ）、モデルパラメータ推定部１０１は、ステップＳ２２に進む。一方で、モデルパラメータ推定値θ_ｅｓｔ（ｋ）と共分散行列Ｐ（ｋ）とを初期化すると判定しなかった場合（ステップＳ２１でＮＯ）、モデルパラメータ推定部１０１は、ステップＳ２３に進む。

【0049】

ステップＳ２２：上記のステップＳ２１でモデルパラメータ推定値θ_ｅｓｔ（ｋ）と共分散行列Ｐ（ｋ）とを初期化すると判定された場合、モデルパラメータ推定部１０１は、θ_ｅｓｔ（０）＝θ_０及びＰ（０）＝Ｉと初期化し、更に制御量前回値ｙ_－１←ｙ（ｔ）と初期化すると共に、ｋ←０と初期化する。ここで、θ_０としては、既にプラント応答関数｛Ｓ_θ（ｔ）｝に設定されているモデルパラメータ設定値θを用いてもよいし、ユーザ等が事前に想定した初期値を用いてもよいし、全ての要素が０のベクトル等の固定の初期値を用いてもよい。また、Ｉとしては、単位行列としてもよいし、予め決められた任意の行列としてもよい。

【0050】

ステップＳ２３：モデルパラメータ推定部１０１は、状態ベクトルφ（ｋ－１）と操作量ｕ（ｔ）と制御量前回値ｙ_－１とに基づいて、状態ベクトルφ（ｋ）を更新する。すなわち、モデルパラメータ推定部１０１は、以下により状態ベクトルφ（ｋ）を更新する。

【数6】

ここで、ｙ（ｋ－１）にはｙ_－１を設定し、ｕ（ｋ）にはｕ（ｔ）を設定する。また、ｙ（ｋ－２），・・・，ｙ（ｋ－Ｎ）及びｕ（ｋ－１），・・・，ｕ（ｋ－Ｍ）には、状態ベクトルφ（ｋ－１）と同じ値を設定する（すなわち、状態ベクトルφ（ｋ－１）に含まれるｙ（ｋ－２），・・・，ｙ（ｋ－Ｎ）及びｕ（ｋ－１），・・・，ｕ（ｋ－Ｍ）をそれぞれ設定する。）。更に、モデルパラメータ推定部１０１は、制御量前回値ｙ_－１←ｙ（ｔ）と更新する。

【0051】

なお、例えば、過去Ｌ点の外乱ｖを考慮する場合には、上述したように、状態ベクトルφ（ｋ）を、ｙ（ｋ）及びｕ（ｋ）に加えてＬ個の外乱ｖ（ｋ）を要素として持つ状態ベクトルに拡大し、外乱に対する共分散行列とモデルパラメータの推定とを行うようにすればよい。

【0052】

ステップＳ２４：次に、モデルパラメータ推定部１０１は、モデルパラメータ推定値θ_ｅｓｔ（ｋ－１）と状態ベクトルφ（ｋ）と制御量ｙ（ｔ）とに基づいて、予測誤差ε（ｋ）を計算する。なお、モデルパラメータ推定値θ_ｅｓｔ（ｋ－１）は前回（つまり、ｋ－１のとき）推定したモデルパラメータθの推定値である。

【0053】

モデルパラメータ推定部１０１は、例えば、以下により予測誤差ε（ｋ）を計算する。

【0054】

ｙ（ｋ）＝ｙ（ｔ）
ε（ｋ）＝ｙ（ｋ）－φ（ｋ）^Τθ_ｅｓｔ（ｋ－１）

【0055】

ステップＳ２５：次に、モデルパラメータ推定部１０１は、共分散行列Ｐ（ｋ）を更新する。モデルパラメータ推定部１０１は、例えば、以下により共分散行列Ｐ（ｋ）を更新する。

【数7】

ここで、Ｐ（ｋ－１）は前回（つまり、ｋ－１のとき）得られた共分散行列である。また、λは０＜λ≦１を取る忘却係数であり、予め設定された値である。忘却係数λは過去データを忘却するための係数であり、０に近いほど急激に過去データの影響が減少し、１に近いほど過去データの影響が保持される。

【0056】

ステップＳ２６：次に、モデルパラメータ推定部１０１は、モデルパラメータ推定値θ_ｅｓｔを更新するか否かを判定する。ここで、モデルパラメータ推定値θ_ｅｓｔを更新すると判定される場合としては、例えば、上記のステップＳ２２での初期化時から所定の期間が経過するまでモデルパラメータ推定値θ_ｅｓｔが更新されなかった場合、ユーザ等により更新指示が行われた場合等が挙げられる。

【0057】

モデルパラメータ推定値θ_ｅｓｔを更新すると判定した場合（ステップＳ２６でＹＥＳ）、モデルパラメータ推定部１０１は、ステップＳ２７に進む。一方で、モデルパラメータ推定値θ_ｅｓｔを更新すると判定しなかった場合（ステップＳ２６でＮＯ）、モデルパラメータ推定部１０１は、ステップＳ２７を実行せずに終了する。

【0058】

ステップＳ２７：上記のステップＳ２６でモデルパラメータ推定値θ_ｅｓｔを更新すると判定された場合、モデルパラメータ推定部１０１は、モデルパラメータ推定値θ_ｅｓｔ（ｋ）を更新する。モデルパラメータ推定部１０１は、例えば、以下によりモデルパラメータ推定値θ_ｅｓｔ（ｋ）を更新する。

【数8】

ステップＳ２８：そして、モデルパラメータ推定部１０１は、モデルパラメータ推定処理の実行時インデックスｋをｋ＋１に更新する。

【0059】

＜先読み応答補正部１１２の動作＞
次に、先読み応答補正部１１２の動作について、図５を参照しながら説明する。図５は、先読み応答補正部１１２の動作の一例を説明するための図である。

【0060】

図５に示すように、先読み応答補正部１１２は、プラント応答関数｛Ｓ_θ（ｔ）｝と、目標偏差ｅ_０（ｔ）と、操作変化量時系列｛ｄｕ（ｔ）｝と、先読み長Ｔ_ｐとが入力されると、制御量が過去の操作変化量によって現在時刻ｔからＴ_ｐ経過後に変化すると予測される値を先読み応答補正値ｙ_ｎ（ｔ）として算出する。なお、先読み長Ｔ_ｐは制御パラメータ計算部１１１により計算される。

【0061】

そして、先読み応答補正部１１２は、先読み応答補正値ｙ_ｎ（ｔ）により目標偏差ｅ_０（ｔ）を補正した補正目標偏差ｅ^＊（ｔ）を出力する。ここで、補正目標偏差ｅ^＊（ｔ）は、ｅ^＊（ｔ）＝ｒ（ｔ）－（ｙ（ｔ）＋ｙ_ｎ（ｔ））＝ｅ_０（ｔ）－ｙ_ｎ（ｔ）で計算される。なお、先読み応答補正値ｙ_ｎ（ｔ）は、例えば、国際公開第２０１６／０９２８７２号や特開２０２０－２１４１１号公報等に記載されている手法により算出することができる。

【0062】

＜操作変化量計算部１１３の動作＞
次に、操作変化量計算部１１３の動作について、図６を参照しながら説明する。図６は、操作変化量計算部１１３の動作の一例を説明するための図である。

【0063】

図６に示すように、操作変化量計算部１１３は、補正目標偏差ｅ^＊（ｔ）と制御ゲインｋ_Ｉとが入力されると、補正目標偏差ｅ^＊（ｔ）に対して制御ゲインｋ_Ｉを乗じて操作変化量ｄｕ（ｔ）を出力する。すなわち、操作変化量ｄｕ（ｔ）は、ｄｕ（ｔ）＝ｋ_Ｉ×ｅ^＊（ｔ）で算出される。

【0064】

ただし、補正目標偏差ｅ^＊（ｔ）に対して制御ゲインｋ_Ｉを乗じた結果が上限値ｄｕ_ｍａｘを超える場合、操作変化量計算部１１３は、ｄｕ_ｍａｘを操作変化量ｄｕ（ｔ）とする。同様に、補正目標偏差ｅ^＊（ｔ）に対して制御ゲインｋ_Ｉを乗じた結果が下限値ｄｕ_ｍｉｎを下回る場合、操作変化量計算部１１３は、ｄｕ_ｍｉｎを操作変化量ｄｕ（ｔ）とする。これにより、上下限範囲に対するリミッターを設けることができる。なお、操作量ｕ（ｔ）と操作変化量ｄｕとが加算器１１４によって加算された後の操作量ｕ_ａが所定の上下限範囲を逸脱しないように、ｄｕ_ｍａｘ及びｄｕ_ｍｉｎを都度設定してもよい。

【0065】

＜制御パラメータの計算＞
次に、制御パラメータとして制御ゲインｋ_Ｉと先読み長Ｔ_ｐを計算する処理について、図７を参照しながら説明する。図７は、制御パラメータの計算処理の一例を説明するためのフローチャートである。なお、図７では、プラント応答関数｛Ｓ_θ（ｔ）｝と調整係数α，βとが入力されたものとする。

【0066】

ステップＳ３１：制御パラメータ計算部１１１は、プラント応答関数｛Ｓ_θ（ｔ）｝と調整係数βとに基づいて、先読み長Ｔ_ｐを計算する。制御パラメータ計算部１１１は、プラント応答関数｛Ｓ_θ（ｔ）｝を用いて、十分長い値として予め設定された最終時刻Ｔ_ｍａｘのプラント応答関数値Ｓ_θ（Ｔ_ｍａｘ）に対して調整係数βを乗じた値にプラント応答関数値が等しくなる時点を先読み長Ｔ_ｐ時点とする。すなわち、制御パラメータ計算部１１１は、以下により先読み長Ｔ_ｐを計算する。

【0067】

Ｆｉｎｄ Ｔ_ｐ，ｗｈｅｒｅ Ｓ_θ（Ｔ_ｐ）＝β×Ｓ_θ（Ｔ_ｍａｘ）
ここで、βは、０＜β≦１であることが好ましい。なお、閉ループの応答が遅い場合には、調整係数βをより小さい値とすることで応答の早さを調整することができる。

【0068】

また、制御パラメータ計算部１１１は、先読み長Ｔ_ｐ時点でのゲインｇ_ｐを以下により計算する。

【0069】

ｇ_ｐ＝Ｓ_θ（Ｔ_ｐ）
このように、先読み長Ｔ_ｐに基づいて制御ゲインｇ_ｐを計算する。これにより、例えば、むだ時間と時定数が長い場合であっても、より早く応答させることが可能となる。

【0070】

なお、Ｔ_ｍａｘは、プラント応答関数が収束するのに十分な値としてもよいし、メモリ上の制約等から定めた値であってもよい。又は、先読み応答補正部１１２で予測可能な期間の最大値をＴ_ｍａｘとしてもよい。

【0071】

ステップＳ３２：そして、制御パラメータ計算部１１１は、先読み長Ｔ_ｐ時点でのゲインｇ_ｐと調整係数αとに基づいて、制御ゲインｋ_Ｉを計算する。制御パラメータ計算部１１１は、以下により制御ゲインｋ_Ｉを計算する。

【数9】

すなわち、先読み長Ｔ_ｐ時点でのゲインｇ_ｐの逆数に対して調整係数αを乗じた値を制御ゲインｋ_Ｉとする。ここで、αは、０＜α≦１であることが好ましい。また、εはゲインが大きくなりすぎるのを回避するための微小な非負の値であり、ε≧０を満たす。

【0072】

＜先読み長Ｔ_ｐの計算概要＞
次に、上記のステップＳ３１における先読み長Ｔ_ｐの計算概要について、図８を参照しながら説明する。図８は、先読み長Ｔ_ｐの計算概要の一例を説明するための図である。

【0073】

図８に示すように、プラント応答関数Ｓ_θ（ｔ）の値が時刻０からＴ_ｍａｘまで徐々に変化していくが、或る時刻Ｔ_ｐにおいてＳ_θ（Ｔ_ｐ）＝β×Ｓ_θ（Ｔ_ｍａｘ）が成り立つ。そこで、このような時刻Ｔ_ｐを先読み長とする。これにより、先読み長Ｔ_ｐ時点におけるプラントの挙動が最終的なプラントの挙動と連動するため、無駄な操作量の変化を抑制することが可能となり、オーバーシュートやアンダーシュートの発生を抑制することができる。また、それに加えて、むだ時間や逆応答が生じる応答初期の期間を回避した先読み長を決定することが可能となる。

【0074】

なお、従来技術として公知のＰＩＤ（Proportional-Integral-Differential）制御においては、むだ時間の大きさに応じて微分ゲインを大きくする必要があり、ノイズに対して弱いという欠点がある。例えば、Ziegler-Nicholsのステップ応答法による制御則では、ｋ_ｐ＝１．２／（ＲＬ）、Ｔ_Ｉ＝２Ｌ、Ｔ_Ｄ＝０．５Ｌ、Ｒ＝ｇ_ｐ／Ｔ_ｐが用いられており、むだ時間Ｌに比例して微分時間Ｔ_Ｄが長く、微分強度が大きくなるため、結果的にノイズやインパルス状の外乱に脆弱な制御となってしまう。これに対して、本実施形態では、観測値の微分信号を用いないため、ノイズやむだ時間に対してロバストな制御を実現することができる。

【0075】

＜先読み長Ｔ_ｐの計算詳細＞
次に、上記のステップＳ３１における先読み長Ｔ_ｐの計算詳細について、図９を参照しながら説明する。図９は、先読み長Ｔ_ｐの計算詳細の一例を説明するためのフローチャートである。

【0076】

ステップＳ４１：まず、制御パラメータ計算部１１１は、先読み長Ｔ_ｐを探索するための変数をＴ_１として、Ｔ_１←０とする。

【0077】

ステップＳ４２：次に、制御パラメータ計算部１１１は、Ｔ_１≧Ｔ_ｍａｘであるか否かを判定する。

【0078】

Ｔ_１≧Ｔ_ｍａｘであると判定した場合（ステップＳ４２でＹＥＳ）、制御パラメータ計算部１１１は、ステップＳ４６に進む。一方で、Ｔ_１≧Ｔ_ｍａｘであると判定しなかった場合（ステップＳ４２でＮＯ）、制御パラメータ計算部１１１は、ステップＳ４３に進む。

【0079】

ステップＳ４３：制御パラメータ計算部１１１は、Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすか否かを判定する。

【0080】

Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすと判定した場合（ステップＳ４３でＹＥＳ）、制御パラメータ計算部１１１は、ステップＳ４５に進む。一方で、Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすと判定しなかった場合（ステップＳ４３でＮＯ）、制御パラメータ計算部１１１は、ステップＳ４４に進む。

【0081】

ステップＳ４４：制御パラメータ計算部１１１は、Ｔ_１←Ｔ_１＋ΔＴと更新し、ステップＳ４２に戻る。これにより、Ｔ_１≧Ｔ_ｍａｘ又はＳ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）のいずれかを満たすまで、ステップＳ４４が繰り返し実行される。

【0082】

ステップＳ４５：Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすと判定した場合、制御パラメータ計算部１１１は、Ｔ_ｐ←Ｔ_１とする。すなわち、制御パラメータ計算部１１１は、Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすＴ_１を先読み長Ｔ_ｐとする。

【0083】

ステップＳ４６：Ｔ_１≧Ｔ_ｍａｘであると判定しなかった場合、制御パラメータ計算部１１１は、Ｔ_ｐ←γ×Ｔ_ｍａｘとする。ここで、０＜γ＜１である。すなわち、Ｓ_θ（Ｔ_１）≧β×Ｓ_θ（Ｔ_ｍａｘ）を満たすＴ_１が探索できなかった場合は、Ｔ_ｍａｘに対してγを乗じた値を先読み長Ｔ_ｐとする。

【0084】

＜高速補完制御部１０５の動作＞
次に、高速補完制御部１０５の動作について、図１０を参照しながら説明する。図１０は、高速補完制御部１０５の動作の一例を説明するための図である。

【0085】

図１０に示すように、主制御周期Ｔ_ｃ毎に主制御部１０４から操作量ｕ_ａが出力されるのに対して、高速補完制御部１０５は、高速制御周期Ｔ_ｆ毎に、操作量ｕ_ｂを出力する。ただし、或る整数ｎ_２，ｎ_３に対してｔ＝ｎ_２Ｔ_ｆかつｔ＝ｎ_３Ｔ_ｃとなる時刻ｔに関しては、ｕ_ｂ（ｔ）＝ｕ_ａ（ｔ）となる。これは、主制御部１０４及び高速補完制御部１０５はいずれも現在時刻ｔにおける制御量ｙ（ｔ）を基準に計算されるためである。以下では、主制御部１０４から出力される操作量ｕ_ａのことを「主制御操作量ｕ_ａ」ともいい、高速補完制御部１０５から出力される操作量ｕ_ｂのことを「補完操作量ｕ_ｂ」ともいう。

【0086】

＜補完操作量ｕ_ｂの計算＞
次に、補完操作量ｕ_ｂを計算する処理について、図１１を参照しながら説明する。図１１は、補完操作量ｕ_ｂの計算処理の一例を説明するためのフローチャートである。なお、図１１では、目標値ｒと、現在時刻ｔにおける制御量ｙ（ｔ）とが入力されたものとして、或るインデックスｋにおける補完操作量ｕ_ｂを計算する場合について説明する。なお、このインデックスｋは、図３でプラント応答関数の計算に用いたインデックスｋや図４でモデルパラメータθの推定に用いたインデックスｋとは独立した値であり、補完操作量ｕ_ｂの計算処理の実行時インデックスを表す。

【0087】

ステップＳ５１：高速補完制御部１０５は、主制御部更新カウンタの値が変化（更新）されたか否かを判定する。言い換えれば、高速補完制御部１０５は、主制御部更新カウンタの値の変化により、主制御部１０４から主制御操作量ｕ_ａが出力されたか否かを判定する。なお、主制御部更新カウンタの値の変化により主制御操作量ｕ_ａが出力されたか否かを判定することは一例であって、主制御部更新カウンタを用いずに、主制御操作量ｕ_ａが出力されたか否かを判定してもよい。例えば、高速補完制御部１０５は、前回の主制御操作量ｕ_ａが出力されてから主制御周期Ｔ_ｃが経過したか否かにより、主制御操作量ｕ_ａが出力されたか否かを判定してもよい。

【0088】

主制御部更新カウンタの値が変化したと判定した場合（ステップＳ５１でＹＥＳ）、高速補完制御部１０５は、ステップＳ５２に進む。一方で、主制御部更新カウンタの値が変化したと判定しなかった場合（ステップＳ５１でＮＯ）、高速補完制御部１０５は、ステップＳ５３に進む。

【0089】

ステップＳ５２：高速補完制御部１０５は、ｙ_０＝ｙ（ｔ）、ｍ_ｄｕ（０）＝０、ｋ＝０と初期化する。ここで、ｍ_ｄｕ（ｋ）は、後述するように、インデックスｋにおける平均操作変化量である。

【0090】

ステップＳ５３：高速補完制御部１０５は、インデックスｋをｋ←ｋ＋１と更新する。なお、インデックスｋはカウンタの役割を持つため、例えば、「補完カウンタ」等と呼ばれてもよい。

【0091】

ステップＳ５４：高速補完制御部１０５は、予測目標偏差ｅ_ｆ（ｋ）を計算する。予測目標偏差ｅ_ｆ（ｋ）とは、現在時刻ｔで予測される先読み長Ｔ_ｐ後の暫定的な目標偏差のことである。高速補完制御部１０５は、例えば、以下の計算方法１又は計算方法２により予測目標偏差ｅ_ｆ（ｋ）を計算する。なお、予測目標偏差ｅ_ｆ（ｋ）の例については後述する。

【0092】

・予測目標偏差ｅ_ｆ（ｋ）の計算方法１
高速補完制御部１０５は、以下により予測目標偏差ｅ_ｆ（ｋ）を計算する。

【数10】

ここで、ｒ（ｔ）は現在時刻ｔの目標値、ｙ（ｔ）は現在時刻ｔの制御量である。

【0093】

・予測目標偏差ｅ_ｆ（ｋ）の計算方法２
高速補完制御部１０５は、以下により予測目標偏差ｅ_ｆ（ｋ）を計算する。

【数11】

ここで、ｒ（ｔ）は現在時刻ｔの目標値、ｙ（ｔ）は現在時刻ｔの制御量である。

【0094】

ステップＳ５５：高速補完制御部１０５は、平均操作変化量ｍ_ｄｕ（ｋ）を計算する。平均操作変化量ｍ_ｄｕ（ｋ）は以下で計算される。

【数12】

ここで、ｄｕ_ｆ（ｋ）は、瞬時的な操作変化量を表す瞬時操作変化量である。このように、平均操作変化量ｍ_ｄｕ（ｋ）は、ｋ＝０以降の瞬時操作変化量ｄｕ_ｆ（ｋ）の平均となっている。

【0095】

ステップＳ５６：高速補完制御部１０５は、補完操作量ｕ_ｂを計算する。インデックスｋにおける補完操作量ｕ_ｂ（ｋ）は、主制御操作量ｕ_ａと、平均操作変化量ｍ_ｄｕ（ｋ）とを用いて、以下で計算される。

【数13】

ここで、ｕ_ｍｉｎ及びｕ_ｍａｘは操作量に対して課された上下限制約である。このｕ_ｍｉｎ及びｕ_ｍａｘは、例えば、主制御部１０４から出力される操作量ｕ_ａに対して課されている上下限制約をそのまま用いてもよい。

【0096】

以上により、現在時刻ｔにおける補完操作量ｕ_ｂ（ｔ）＝ｕ_ｂ（ｋ）が得られる。

【0097】

＜計算方法１による予測目標偏差ｅ_ｆ（ｋ）＞
次に、上記のＳ５４で計算方法１によって計算される予測目標偏差ｅ_ｆ（ｋ）について、図１２を参照しながら説明する。図１２は、予測目標偏差ｅ_ｆ（ｋ）の一例を説明するための図（その１）である。

【0098】

図１２に示すように、計算方法１による予測目標偏差ｅ_ｆ（ｋ）は、主制御部更新カウンタが更新された直後の時刻ｔ＝ｔ_０における制御量ｙ_０＝ｙ（ｔ_０）から先読み長Ｔ_ｐに目標値ｒ（ｔ_０＋Ｔ_ｐ）に到達するような線形な軌道（制御軌道）を考えたときに、目標値ｒ（ｔ_０＋Ｔ_ｐ）と、現在時刻ｔにおける制御量ｙ（ｔ）（つまり、高速補完制御部１０５の動作時の制御量ｙ（ｔ））から当該軌道により時刻ｔ_０＋Ｔ_ｐの制御量を予測したときの値（予測制御量）との差として計算される。言い換えれば、計算方法１による予測目標偏差ｅ_ｆ（ｋ）は、目標値ｒ（ｔ_０＋Ｔ_ｐ）と、現在時刻ｔにおける制御量ｙ（ｔ）に当該軌道を平行移動させたときの当該軌道上の時刻ｔ_０＋Ｔ_ｐの制御量（予測制御量）との差として計算される。

【0099】

＜計算方法２による予測目標偏差ｅ_ｆ（ｋ）＞
次に、上記のＳ５４で計算方法２によって計算される予測目標偏差ｅ_ｆ（ｋ）について、図１３を参照しながら説明する。図１３は、予測目標偏差ｅ_ｆ（ｋ）の一例を説明するための図（その２）である。

【0100】

図１３に示すように、計算方法２による予測目標偏差ｅ_ｆ（ｋ）は、主制御部更新カウンタが更新された直後の時刻ｔ＝ｔ_０における制御量ｙ_０＝ｙ（ｔ_０）から現在時刻ｔにおける制御量ｙ（ｔ）までの線形な軌道により時刻ｔ_０＋Ｔ_ｐの制御量を予測したときの値（予測制御量）と、目標値ｒ（ｔ_０＋Ｔ_ｐ）との差として計算される。言い換えれば、計算方法２による予測目標偏差ｅ_ｆ（ｋ）は、当該軌道を時刻ｔ_０＋Ｔ_ｐまで延長したときの制御量（予測制御量）と、目標値ｒ（ｔ_０＋Ｔ_ｐ）との差として計算される。

【0101】

＜高速補完制御部１０５を用いた場合の応答＞
次に、高速補完制御部１０５を用いた場合の応答について、図１４を参照しながら説明する。図１４は、高速補完制御部１０５を用いた場合の応答の一例を説明するための図である。

【0102】

図１４に示すように、主制御部１０４により先読み長Ｔ_ｐ先の時点で目標偏差が０となるような制御が行われる場合、高速補完制御部１０５を用いた場合の応答でも、先読み長Ｔ_ｐ先の時点で目標偏差が０となるような軌道に沿う制御が行われる。ただし、高速補完制御部１０５を用いた場合の応答は、主制御部１０４のみを用いた軌道とは必ずしも一致せず、例えば、よりオーバーシュートが抑制された軌道になることが考えられる。

【0103】

＜制御装置１０のハードウェア構成＞
次に、本実施形態に係る制御装置１０のハードウェア構成について、図１５を参照しながら説明する。図１５は、第二の実施形態に係る制御装置１０のハードウェア構成の一例を示す図である。

【0104】

図１５に示すように、第二の実施形態に係る制御装置１０は、入力装置２０１と、表示装置２０２と、外部Ｉ／Ｆ２０３と、通信Ｉ／Ｆ２０４と、プロセッサ２０５と、メモリ装置２０６とを有する。これらの各ハードウェアは、それぞれがバス２０７により通信可能に接続される。

【0105】

入力装置２０１は、例えば、タッチパネルや各種ボタン等である。表示装置２０２は、例えば、表示パネル等である。なお、制御装置１０は、入力装置２０１及び表示装置２０２のうちの少なくとも一方を有していなくてもよい。

【0106】

外部Ｉ／Ｆ２０３は、記録媒体２０３ａ等の外部装置とのインタフェースである。記録媒体２０３ａとしては、例えば、ＳＤメモリカード（Secure Digital memory card）やＵＳＢ（Universal Serial Bus）メモリカード等が挙げられる。

【0107】

通信Ｉ／Ｆ２０４は、制御装置１０を通信ネットワークに接続するためのインタフェースである。プロセッサ２０５は、例えば、ＣＰＵやＭＰＵ等の各種演算装置である。メモリ装置２０６は、例えば、ＳＳＤ（Solid State Drive）、ＲＡＭ（Random Access Memory）、ＲＯＭ（Read Only Memory）、フラッシュメモリ等の各種記憶装置である。

【0108】

なお、図１５に示すハードウェア構成は一例であって、制御装置１０は、他のハードウェア構成を有していてもよい。例えば、制御装置１０は、図示したハードウェア以外にも種々のハードウェアを有していてもよい。

【0109】

［第二の実施形態］
以下、第二の実施形態について説明する。なお、第二の実施形態では、第一の実施形態との相違点について説明し、第一の実施形態と同様としてよい構成要素についてはその説明を省略する。

【0110】

＜制御装置１０の全体構成例＞
第二の実施形態に係る制御装置１０の全体構成例について、図１６を参照しながら説明する。図１６は、第二の実施形態に係る制御装置１０の全体構成の一例を示す図である。

【0111】

図１６に示すように、第二の実施形態に係る制御装置１０は、事前に推定したモデルパラメータθをプラント応答関数｛Ｓ_θ（ｔ）｝に設定しておくことで、第一の実施形態で説明した全体構成例からモデルパラメータ推定部１０１を除いたものである。

【0112】

［実施例］
次に、実施例について説明する。本実施例では、一例として、第二の実施形態に係る制御装置１０によって制御対象プラント２０を制御する場合について説明する。

【0113】

本実施例では、主制御周期Ｔ_ｃ＝５［ｓｅｃ］、高速制御周期Ｔ_ｆ＝０．２［ｓｅｃ］と設定した。すなわち、高速補完制御部１０５は、主制御部１０４の２５倍の速さで動作するものとした。また、操作量の上限ｕ_ｍａｘ＝９０、下限ｕ_ｍｉｎ＝０とし、調整係数はα＝１．０、β＝０．５とした。

【0114】

本実施例における制御対象プラント２０のステップ応答を図１７に示す。図１７に示すように、本実施例では、ステップ応答が収束するまでにおよそ１５０［ｓｅｃ］掛かっている。このため、例えば、主制御周期をＴ_ｃ＝０．５［ｓｅｃ］とした場合、ステップ応答の収束までの予測にはおよそ３００点の予測バッファが必要になる。一方で、本実施例では主制御周期をＴ_ｃ＝５［ｓｅｃ］としているため、ステップ応答の収束までの予測には３０点の予測バッファで済む。したがって、主制御周期をＴ_ｃ＝５［ｓｅｃ］と長くしておくことにより、エッジデバイス等である制御装置１０の予測バッファを約１／１０に節約し、メモリ削減と共に計算時間も短縮することができる。

【0115】

本実施例における制御装置１０との比較例として、第二の実施形態に係る制御装置１０の高速補完制御部１０５を無効とした場合における制御量、目標値及び操作量を図１８に示す。図１８では、上段に制御量（ＰＶ）、目標値（ＳＰ）を、下段に操作量（ＭＶ）を示している。また、右図は、左図の一部の時間を切り出して拡大したものである。図１８に示すように、高速補完制御を無効とし、主制御部１０４のみを用いた制御では、主制御周期Ｔ_ｃ＝５［ｓｅｃ］毎に階段上に操作量が変化している。

【0116】

予測目標偏差ｅ_ｆ（ｋ）を計算方法１により計算するときの本実施例における制御装置１０による制御量、目標値及び操作量（つまり、高速補完制御部１０５による高速補完制御が有効である場合における制御量、目標値及び操作量）を図１９に示す。図１９では、上段に制御量（ＰＶ）、目標値（ＳＰ）を、下段に操作量（ＭＶ）を示している。また、右図は、左図の一部の時間を切り出して拡大したものである。図１９に示すように、図１８と比較すると、制御量の目標値への追従性が向上しており、オーバーシュートが低減している。また、操作量の拡大図を見ると、操作量は主制御周期Ｔ_ｃでの変化に加え、高速補完制御部１０５による高速制御周期Ｔ_ｆでの変化が生じていることが確認できる。したがって、高速補完制御部１０５による高速補完制御を有効にすることにより、主制御部１０４のみを用いる場合よりも制御性能を向上させることができる。

【0117】

また、予測目標偏差ｅ_ｆ（ｋ）を計算方法２により計算するときの本実施例における制御装置１０による制御量、目標値及び操作量を図２０に示す。図２０では、上段に制御量（ＰＶ）、目標値（ＳＰ）を、下段に操作量（ＭＶ）を示している。また、右図は、左図の一部の時間を切り出して拡大したものである。図２０に示すように、図１８と比較すると、制御量の目標値への追従性が向上しており、オーバーシュートがほぼゼロになっている。また、操作量の拡大図を見ると、操作量は主制御周期Ｔ_ｃでの変化に加え、高速補完制御部１０５による高速制御周期Ｔ_ｆでの変化が生じていることが確認できる。したがって、高速補完制御部１０５による高速補完制御を有効にすることにより、主制御部１０４のみを用いる場合よりも制御性能を向上させることができる。

【0118】

なお、図１９と図２０を比較すると、図２０では操作量の変動が急速となっているが、これは予測目標偏差ｅ_ｆ（ｋ）の計算方法２では制御量の時間方向の傾きに大きく影響されるためである。

【0119】

予測目標偏差ｅ_ｆ（ｋ）の計算方法１と計算方法２ではオーバーシュート量と操作量の変動幅が異なっている。このため、操作量変動が大きくてもオーバーシュート量を抑えたい用途では計算方法２を、操作量変動をできるだけ抑制しつつオーバーシュート量もバランスよく抑制したい用途では計算方法１を用いることが好適である。例えば、ＰＷＭ（Pulse Width Modulation）制御におけるＯＮ時間を操作量とするような電気的制御では操作量の変動が大きくとも問題が生じるケースは少ないため、計算方法２を用いることができる。一方で、流量の電動弁制御の開度を操作量とするような機械的制御では操作量の短期的な変動が機械的ストレスになり得るため、計算方法１を用いるのが好適である。

【0120】

［まとめ］
以上のように、第一及び第二の実施形態に係る制御装置１０は、主制御部１０４によりも短い制御周期Ｔ_ｆで制御対象プラント２０を高速する制御する高速補完制御部１０５を有している。これにより、主制御部１０４の制御周期Ｔ_ｃが長い場合であっても、制御周期Ｔ_ｆ毎に制御対象プラント２０が高速に制御されるため、例えば、制御量に対する目標値の変更や変動の早い外乱の影響等にも素早く追従が可能となり、制御性能の低下を抑制することができる。また、制御装置１０の計算資源が乏しく、主制御周期Ｔ_ｃを長くせざるを得ない場合であっても、制御性能の低下を抑制することができる。更に、モデル予測制御では抑制が不十分な制御偏差（例えば、オーバーシュート等による制御偏差）の改善も期待できる。

【0121】

なお、第一の実施形態のようにモデルパラメータ推定部１０１によってモデルパラメータθの推定を行う場合、高速補完制御部１０５が有効だと、高速補完制御を含めた閉ループの同定となってしまうため、高速補完制御部１０５が無効のときのみモデルパラメータ推定部１０１によるモデルパラメータθの推定を行うものとしてもよい。

【0122】

また、主制御部１０４の動作は第一の実施形態で説明したものに限られず、制御方式としてモデル予測制御を採用したものであれば任意の動作を採用することが可能であるが、高速補完制御部１０５が制御量ｙ（ｔ）を起点として、時刻ｔ_０＋Ｔ_ｐに目標値ｒに到達する軌道（制御軌道）を暫定的に設定して、当該軌道に沿うように制御ゲインｋ_Ｉによって補完制御するため、主制御部１０４も同様な動作をするものであることが好ましい。

【0123】

本発明は、具体的に開示された上記の実施形態に限定されるものではなく、特許請求の範囲の記載から逸脱することなく、種々の変形や変更、既知の技術との組み合わせ等が可能である。

【符号の説明】

【0124】

１０ 制御装置
２０ 制御対象プラント
１０１ モデルパラメータ推定部
１０２ 計測部
１０３ 差分器
１０４ 主制御部
１０５ 高速補完制御部
１０６ 第１のタイマ
１０７ 第２のタイマ
１０８ 切替器
１１１ 制御パラメータ計算部
１１２ 先読み応答補正部
１１３ 操作変化量計算部
１１４ 加算器
２０１ 入力装置
２０２ 表示装置
２０３ 外部Ｉ／Ｆ
２０３ａ 記録媒体
２０４ 通信Ｉ／Ｆ
２０５ プロセッサ
２０６ メモリ装置
２０７ バス

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】

【図17】

【図18】

【図19】

【図20】