(57)【特許請求の範囲】
【請求項１】
モデル予測制御方法において、
前記モデル予測制御方法における予測区間のうち入力を変化させる区間である制御区間において、前記モデル予測制御方法における予測区間で予測値を算出する周期であるサンプリング周期を１番目のサンプリング間隔から線形に増加させる工程と、
前記サンプリング周期における１番目のサンプリング間隔と制御周期とを等しくする工程と、
前記サンプリング周期におけるサンプリング間隔の増加に伴って、モデル予測制御の評価関数で用いる重みを線形に増加させる工程と、
を備えるモデル予測制御方法。
【請求項２】
モデル予測制御プログラムにおいて、
コンピュータに、
前記モデル予測制御方法における予測区間のうち入力を変化させる区間である制御区間において、前記モデル予測制御方法における予測区間で予測値を算出する周期であるサンプリング周期を１番目のサンプリング間隔から線形に増加させる処理と、
前記サンプリング周期における１番目のサンプリング間隔と制御周期とを等しくする処理と、
前記サンプリング周期におけるサンプリング間隔の増加に伴って、モデル予測制御の評価関数で用いる重みを線形に増加させる処理と、
を実行させるモデル予測制御プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、モデル予測制御方法及びモデル予測制御プログラムに関し、特に所定のサンプリング周期で制約条件（重み）を満たすように制御対象の未来の入力値を算出して、１番目の入力値を現在の入力値とする工程を、所定の制御周期で繰り返すことによって断続的に現在の入力値を算出するモデル予測制御方法及びモデル予測制御プログラムに関する。
【背景技術】
【０００２】
一般的なモデル予測制御を以下に説明する。但し、本明細書においては、各式において小文字の太字はベクトルを表し、大文字の太字は行例を表す。但し、太字ではない通常の文字はスカラー量を表す。
【０００３】
制御対象の連続システムが式１、２で表され、離散システムが式３、４で表される。詳細には、連続システムである式１、２をΔｔ（スカラー量）で離散化し、状態変数ｘｋ（ベクトル）と、入力値ｕｋ（スカラー量）及び出力値ｗｋ（スカラー量）を用いて表した結果が、離散システムである式３、４である。
【数１】
0005765222.tif 000002
【０００４】
ここで、ｘ（スカラー量）は状態変数、Ａ（行列）は係数行列、ｂ（ベクトル）は入力ベクトル、ｃ（ベクトル）は出力ベクトル、ｋ（スカラー量）は時刻ステップ、ｔ（スカラー量）は連続時間、Δｔ（スカラー量）はサンプリング間隔を示す。
【０００５】
また、ハードウェアの制限や制御対象の運転条件等から、入力値や出力値に対して制御条件が課される場合の制約条件は、
【０００６】
【数２】
0005765222.tif 000003
式５、６に表すように、線形不等式の形に定式化される。
【０００７】
ここで、Ｓ１（行列）は入力制約の係数行列、ｄ１（ベクトル）は入力制約の係数ベクトル、Ｓ２（行列）は出力制約の係数行列、ｄ２（ベクトル）は出力制約の係数ベクトルを示す。
【０００８】
但し、ｕｋ（ベクトル）は現在からＮΔｔ（スカラー量）未来までの入力値を並べたベクトルであり、以下の式７で表される。ここで、Ｎは自然数を示す。
【０００９】
【数３】
0005765222.tif 000004
【００１０】
また、ｗｋ＋１（ベクトル）はＮΔｔ（スカラー量）未来までの出力値を並べて、以下の式８に示すように、ベクトルにしたものである。
【００１１】
【数４】
0005765222.tif 000005
【００１２】
モデル予測制御は現在から有限時間（ＮΔｔ（スカラー量））までの間で、式９で表される２次形式の評価関数を最小化するように入力値を決定する。
【００１３】
【数５】
0005765222.tif 000006
【００１４】
但し、上添え字のｒｅｆは目標値を表す。即ち、ｗｒｅｆｋ＋１は、式１０で示すように目標出力のＮΔｔ（スカラー量）未来までの時系列を意味する。
【００１５】
【数６】
0005765222.tif 000007
【００１６】
また、Ｑ（行列）、Ｒ（行列）は重み行列として、式１１、１２で示すように各予測点における重み値を対角に並べた対角行列である。
【００１７】
【数７】
0005765222.tif 000008
【００１８】
モデル予測制御は、式５及び式６の制約条件の範囲内で式９の評価関数を最小化する、入力時系列ｕｋ（ベクトル）を求め、ｕｋ(ベクトル)の１番目の値ｕｋ（スカラー量）を現在の入力値として用いる制御方法である。そして、Δｔ後に同様に最適解を求め、入力値時系列の１番目の値を入力値として用いる、という処理を繰り返す。
【００１９】
つまり、図１２に示すように、制御周期毎に未来の予測期間ＮΔｔ（スカラー量）を後退させながら、制御条件付き最適化問題（式５、６及び式９）を解いてその時点の入力値を決定していく制御方法である。ちなみに、図１２では太線部分が１番目の予測点である。
【００２０】
以下、上記の定式化に従って最適化問題を導く。
式３、４より、
【数８】
0005765222.tif 000009
【００２１】
式１３が導き出されるので、出力値時系列は、
【数９】
0005765222.tif 000010
式１４で表される。
【００２２】
従って、式９は式１５のように定式化される。ここで、Ｈｋ（行列）は半正定行列を示す。
【数１０】
0005765222.tif 000011
【００２３】
同様にして、式５、６についても、式１６に示すように、ｕｋ（ベクトル）に関する式で表すことができる。
【数１１】
0005765222.tif 000012
【００２４】
最終的に最適解問題は、式１５、１６より、
【数１２】
0005765222.tif 000013
式１７に表すように、凸２次計画問題に定式化され、この凸２次計画問題の求解結果ｕｋ（ベクトル）の第１番目要素が現在の入力値として実際に用いられる。
【００２５】
以上が一般的なモデル予測制御の制御則である。
なお、上記の制御則では予測区間の間、常に入力値を変化させ続けるように構成したが、モデル予測制御では、図１３に示すように、予測区間内で入力値を変化させる区間（制御区間）と入力値を変化させない区間が導入される場合もある。ここで、Ｍは自然数を示し、Ｎ以下であることを満たす。
【００２６】
予測区間内で制御区間外にある予測点については、制御区間の最終値に設定されることが一般的である。即ち、式１８であり、
【数１３】
0005765222.tif 000014
ｕｋ（ベクトル）を
【数１４】
0005765222.tif 000015
式１９と設定し直して改めて定式化を行うと、式３、４より、
【数１５】
0005765222.tif 000016
式２０となる。
【００２７】
このため、出力値時系列は、
【数１６】
0005765222.tif 000017
式２１と表される。後は制御区間と予測区間が等しい場合のモデル予測制御と同様に式１５及び式１６のように定式化され、最終的に式１７のように最適化問題が導かれる。
【００２８】
このように一般的なモデル予測制御は、モデルに基づく未来の予測・制御条件の設定・最適解の求解という３つの特徴を含んでいるため、非常に適用可能な範囲が広く有用な制御法である。
【００２９】
しかし、一番のネックは計算コストが高いことになる。これは、モデル予測制御が最適化問題（凸２次計画問題）の求解ステップを含んでいることに起因している。それ故に、多くは時定数の遅いプロセス系の制御に用いられることが多かった。
【００３０】
そのため、モデル予測制御の計算コストを低減するために、例えば特許文献１の技術が提案されている。
特許文献１の技術は、未来の制御量（入力値）を演算する区間内のサンプリング周期のサンプリング間隔を変更している。つまり、上述したモデル予測制御では、予測区間内のサンプリング間隔は一定であるが、特許文献１の技術は、制御量を変動させる区間ではサンプリング間隔を短くし、制御量を変動させない区間ではサンプリング間隔を長くしても、制御の精度がそれ程落ちないという着想から、未来の制御量を演算する区間内でサンプリング間隔を切り替えている。
【先行技術文献】
【特許文献】
【００３１】
【特許文献１】特開２００６−７２７９１号公報
【発明の概要】
【発明が解決しようとする課題】
【００３２】
例えばロボットの制御のように時定数が短く、かつ、速応性が要求されるような制御では、モデル予測制御の制御区間を長く取る必要がある。一般的に制御区間が短い方が、閉ループの安定性は高いので、制御区間を長くするということは、安定性が低下する替わりに速応性を向上させることを意味する。
【００３３】
ロボットの制御では、時定数が短いので、制御周期も早くする必要があり、モデル予測制御の計算時間も短くする必要がある。そこで、特許文献１の技術のように、計算コストを下げるために、未来の制御量を演算する区間内のサンプリング間隔を変化させる場合、サンプリング間隔の変化のさせ方を適当に設定しないと、閉ループの安定性が低く、ロボットが振動を起こしたり、不安定になって発散してしまったり、する可能性がある。
【００３４】
しかし、特許文献１には、サンプリング間隔の変化のさせ方、即ち、どのようなサンプリング間隔にするのか、という内容に関しては開示されていない。
【００３５】
本発明の目的は、このような問題を解決するためになされたものであり、閉ループの安定性が高く、計算コストを低減できるモデル予測制御方法及びモデル予測制御プログラムを提供することである。
【課題を解決するための手段】
【００３６】
本発明の一形態に係るモデル予測制御方法は、制御区間において、サンプリング周期を１番目のサンプリング間隔から線形に増加させる工程と、前記サンプリング周期における１番目のサンプリング間隔と制御周期とを等しくする工程と、前記サンプリング周期におけるサンプリング間隔の増加に伴って、重みを線形に増加させる工程と、を備える。
【００３７】
本発明の一形態に係るモデル予測制御プログラムは、コンピュータに、制御区間において、サンプリング周期を１番目のサンプリング間隔から線形に増加させる処理と、前記サンプリング周期における１番目のサンプリング間隔と制御周期とを等しくする処理と、前記サンプリング周期におけるサンプリング間隔の増加に伴って、重みを線形に増加させる処理と、を実行させる。
【発明の効果】
【００３８】
以上、説明したように、本発明によると、閉ループの安定性が高く、計算コストを低減できるモデル予測制御方法及びモデル予測制御プログラムを提供することができる。
【図面の簡単な説明】
【００３９】
【図１】本発明の実施の形態１に係るモデル予測制御方法を説明するための図である。
【図２】本発明の実施の形態１に係るモデル予測制御方法において、サンプリング周期のサンプリング間隔を線形に増加させた際の、予測点とサンプリング間隔との関係を示す図である。
【図３】本発明の実施の形態２に係るモデル予測制御方法を説明するための図である。
【図４】本発明の実施例に係るモデル予測制御方法での、２足歩行ロボットの重心モデルを模式的に示す図である。
【図５】本発明の実施例に係るモデル予測制御方法での、目標ＺＭＰとＺＭＰ制約との関係を示す図である。
【図６】本発明の実施例に係るモデル予測制御方法での、ＺＭＰ軌道と重心軌道との関係を示す図である。
【図７】本発明の実施例に係るモデル予測制御方法での、重心速度を示す図である。
【図８】本発明の実施例に係るモデル予測制御方法での、重心加速度を示す図である。
【図９】本発明の実施例に係るモデル予測制御方法において、２足歩行ロボットに外乱が生じた際のＺＭＰ軌道と重心軌道との関係を示す図である。
【図１０】本発明の実施例に係るモデル予測制御方法において、２足歩行ロボットに外乱が生じた際の重心速度を示す図である。
【図１１】本発明の実施例に係るモデル予測制御方法において、２足歩行ロボットに外乱が生じた際の重心加速度を示す図である。
【図１２】一般的なモデル予測制御を説明するための図である。
【図１３】制御区間と予測区間とが異なる、一般的なモデル予測制御を説明するための図である。
【発明を実施するための形態】
【００４０】
以下、本発明を実施するための最良の形態について、添付図面を参照しながら説明する。但し、本発明が以下の実施の形態に限定される訳ではない。また、説明を明確にするため、以下の記載及び図面は、適宜、簡略化されている。
【００４１】
＜実施の形態１＞
本実施の形態のモデル予測制御方法及びモデル予測制御プログラムは、例えばロボット等の時定数が短い制御対象を制御するために好適に実行される。本実施の形態のモデル予測制御方法は、当該制御対象の制御装置がモデル予測制御プログラムを実行することで実現される。但し、当該モデル予測制御方法は、ハードウェア資源を用いて実現しても良い。ちなみに、本実施の形態のモデル予測制御方法は、制御区間と予測区間とを等しくしている。
【００４２】
先ず、本実施の形態のモデル予測制御方法の概略を説明する。当該モデル予測制御方法は、図１及び図２に示すように、制御区間において、サンプリング周期を１番目のサンプリング間隔から線形に増加させて、各々の予測点での入力値を算出し、１番目の予測点での入力値を現在の入力値として用いる。ちなみに、図１では太線部分が1番目の予測点である。また、図２の縦軸Δｔはサンプリング間隔を示す、横軸ｊは何番目の予測点かを示し、ｊ番目のサンプリング間隔ΔｔｊをΔｔ１＋α（ｊ−１）として算出している。ここで、Δｔ１は1番目のサンプリング間隔を示し、αは係数を示す。
【００４３】
このとき、１番目の予測点での制約条件と、制御周期後の現在状態が対応するように、１番目のサンプリング間隔と制御周期とを等しくする。また、上記のようにサンプリング間隔の増加に伴って、制約条件（重み）を線形に増加させる。
【００４４】
これにより、詳細な作用効果は後述するが、最適化問題のサイズを削減することができ、制御区間を大きく取った場合でも計算コストを大幅に低減できる。
【００４５】
次に、本実施の形態のモデル予測制御方法を詳細に説明する。一般的なモデル予測制御方法では、重み行列Ｑ、Ｒの要素を個別に大きく変化させることはせず、各要素の値を概ね同じ値にするのが一般的である。すなわち、以下の式２２、２３に示すように設定する。
【数１７】
0005765222.tif 000018
【００４６】
しかし、上記のように予測時間内のサンプリング間隔を変化させた場合、各要素の重みを一定にすると、各サンプリング点（予測点）の影響度が異なるため、適切に制御することができない。これは、予測点毎に式１、２の連続システムが異なる（Ａ（行列）、ｂ（ベクトル）がΔｔ（スラカー量）の関数になっている）ためである。サンプリング間隔が大きい程、同じ入力値ｕｋ（スカラー量）でもｘｋ（ベクトル）→ｘｋ＋１（ベクトル）の変化量は大きくなるので、本実施の形態では重み行列の各要素を、
【数１８】
0005765222.tif 000019
式２４、２５のように線形に増加させる。ここで、β及びγは定数を示す。
【００４７】
即ち、重み行列を
【数１９】
0005765222.tif 000020
式２６、２７とし、対角要素の値が線形に増加するように設定する。
【００４８】
以下、上記の構成で定式化を行う。
各サンプリング間隔Δｔｊ（スカラー量）に対応した離散システム行列をＡｊ（行列）、ｂｊ（スカラー量）とすると、式３、４は
【数２０】
0005765222.tif 000021
式２８、２９となる。従って、
【数２１】
0005765222.tif 000022
式３０となるので、出力値時系列は、
【数２２】
0005765222.tif 000023
式３１と表される。後は一般的なモデル予測制御と同様にして、式９は、
【数２３】
0005765222.tif 000024
式３２のように定式化される。
【００４９】
制約条件も同様にして、式５、６についても、
【数２４】
0005765222.tif 000025
として、ｕｋ（ベクトル）に関する式３３で表すことができる。
【００５０】
最終的に最適化問題は、式３２、３３より、
【数２５】
0005765222.tif 000026
式３４に表すように、凸２次計画問題に定式化され、この凸２次計画問題の求解結果ｕｋ（ベクトル）の第１番目要素を現在の入力値として実際に用いる。
【００５１】
以上が本実施の形態のモデル予測制御の制御則である。本実施の形態のモデル予測制御は、上述のように制御区間においてサンプリング間隔を線形に増加させている。その作用効果について、本出願人は以下のように推測している。
【００５２】
目標の軌道（ｒｅｆ）の変化速度に対して予測区間（制御区間）のサンプリング間隔が長い場合、予測点間に目標軌道の山や谷が埋もれてしまう。その結果、時間が進むにつれて、前の時間の予測では入らなかった目標軌道の山や谷が、後の時間になってから現れてきたり、また逆に消えたりということが起る。この時、制御周期毎の入力値は振動的になり、状態量も目標軌道の周りで振動的な挙動を示す。
【００５３】
一方、モデル予測制御では予測区間のうち現在に近い時点の方が現在の入力値に対する影響が大きい。即ち、遠い未来の目標値が変わっても、現在の入力値はそれ程大きく変化しないが、近い未来の目標値の変化に対しては、現在の入力値を大きく変化させて対応することになる。
【００５４】
また、モデル予測制御では制御区間の予測点で制約条件を満たすように求解が行われるが、予測点間では制約条件を満たすことは保証されていない。制御区間のサンプリング間隔を変化させた場合、今回の予測点がしばしば以前の求解時の予測点間に設定されるが、予測点間では制約条件を満たす保証がないので、前回の軌道と今回の求解結果とが大きく異なる場合があり（前回の求解時には制約が守られていなかった点が、今回の求解では急に満たすように求解されるようになるため）、この時も振動の要因となる。制御区間のサンプリング間隔を急激に変化させるようなパターンを用いた場合、予測点間で制約を破る解が出力されることが特に多く、振動が発生し易くなる。
【００５５】
こうした技術的知見から、現在から予測区間（制御区間）の未来に向かって、線形にサンプリング間隔を増加させることにより、振動を発生させることなく、即ち閉ループを安定させて、予測区間の点の数を効果的に減らして計算コストを削減することが可能になる。
【００５６】
また、本実施の形態のモデル予測制御は、重みを線形に増加させているが、その作用効果について、本出願人は以下のように推測している。
【００５７】
モデル予測制御では、制御区間の予測点間では入力値が一定という仮定が前提となっている。つまり、予測点Ａから次の予測点Ｂまでの間では、離散システムにずっと予測点Ａにおける入力が作用しているとして定式化が行われる。
【００５８】
従って、サンプリング間隔を変化させた場合には、同じ入力値が作用し続ける時間が入力毎に異なってしまい、適切に最適化を行うことができない。そこで、サンプリング間隔を線形に増加させたのに合わせて重み行列を線形に増加させることにより、サンプリング間隔の大きな点の入力値ほど評価関数上の影響度が大きくなるように設定することができ、適切に最適化を行うことが可能となる。
【００５９】
＜実施の形態２＞
実施の形態１のモデル予測制御方法は、制御区間と予測区間とを等しくしたが、図３に示すように、一般的なモデル予測制御方法と同様に制御区間と予測区間とが異なるようにしても良い。なお、図３でも太線部分が1番目の予測点である。
【００６０】
この場合、ｕｋ（ベクトル）を式１９と設定し直して改めて定式化を行うと、
【数２６】
0005765222.tif 000027
式３５となるので、出力値時系列は、
【数２７】
0005765222.tif 000028
式３６と表される。その後の導出は実施の形態１のモデル予測制御方法と同様であるので説明を省略する。これにより、より長い区間の予測区間を確保することができる。
【００６１】
以上、本発明に係るモデル予測制御方法及びモデル予測制御プログラムの実施の形態を説明したが、上記に限らず、本発明の技術的思想を逸脱しない範囲で、変更することが可能である。
【００６２】
＜実施例＞
本発明のモデル予測制御方法及びモデル予測制御プログラムを２足歩行ロボットの重心軌道生成問題に適用した実施例を説明する。
図４に示すように、重心高さを一定とすると、ＺＭＰ（Zero Moment Point）方程式は良く知られているように、
【数２８】
0005765222.tif 000029
式３７となる。ここで、ｘｚｍｐは２足歩行ロボットのＺＭＰにおけるｘ座標、ｘｇは当該重心位置のｘ座標、ｈは当該重心位置の高さ、ｇは重力加速度を示す。
【００６３】
入力値をｕ（ｔ）（ベクトル）＝ｘ（・・・）ｇ（スカラー量）、出力値をｗ（ｔ）(スカラー量)＝ｘｚｍｐ（スカラー量）とすると連続システムは、
【数２９】
0005765222.tif 000030
式３８、３９となる。これをΔｔｊ（スカラー量）で離散化すると、
【数３０】
0005765222.tif 000031
式４０、４１となる。
【００６４】
例として、目標のＺＭＰ軌道とＺＭＰの制約条件（出力値の制約条件）を図５のように設定した。
【００６５】
本実施例では、予測区間と制御区間は同一にして１．６［ｓ］、制御周期は１［ｍｓ］として、サンプリング周期は図１と同様の設定を用いた（α＝０．００５）。
【００６６】
また、重みに関しては、Ｑ（スカラー量）＝１、Ｒ（スカラー量）＝１×１０−６、β(スカラー量)＝０．００５、γ（スカラー量）＝５×１０−９とした。
【００６７】
本実施例の結果を図６〜図８に示した。図６は、図５にＺＭＰ軌道と重心軌道とを重ね合わせてプロットした図である。本図より、目標ＺＭＰ軌道に追従しながら滑らかな重心軌道を生成できていることが分かる。また、図７及び図８は、重心速度・重心加速度の結果をプロットしたグラフだが、これらのグラフから、速度・加速度の次元でも振動等を発生することなく軌道生成が行えることが分かる。
【００６８】
次に制約条件が守られていることを確認するために、上記と同様の条件で、２．６［ｓ］の時に仮想的に突発外力を入力し、速度をステップ上に変化させた。この時のＺＭＰ軌道と重心軌道とを図９に、重心速度と加速度とを図１０及び図１１に示した。図９から分かるように、ＺＭＰは制約条件を満たしながら、目標ＺＭＰ軌道に戻っていることが確認できる。
【図１】
0005765222.tif 000032
【図２】
0005765222.tif 000033
【図３】
0005765222.tif 000034
【図４】
0005765222.tif 000035
【図５】
0005765222.tif 000036
【図６】
0005765222.tif 000037
【図７】
0005765222.tif 000038
【図８】
0005765222.tif 000039
【図９】
0005765222.tif 000040
【図１０】
0005765222.tif 000041
【図１１】
0005765222.tif 000042
【図１２】
0005765222.tif 000043
【図１３】
0005765222.tif 000044