現代制御論

制御

概要 †

• 最適レギュレータをなるべく速く実装するためのページ
• 以下が可能
  • octaveによる極配置、最適レギュレータ
  • 実機とのパラメータ合わせ
  • 実機のモータを擬似力制御する方法（まとも度：ひずみゲージ＞電流制御・疑似電流制御＝電圧＋速度制御、高剛性バネ＋変位センサ）
• ゲイン多様体による制御器の構造の把握

参考 †

• octaveスクリプト
• 倒立振子モデル

疑問 †

• 独立制御可能性を示すためにはどうすればいいか？（直進と回転など）
• 制御不要の時にフィードバックゲインの一部を0にすることは最適レギュレータに対してどう影響するか。

基本式 †

• \dot{x}=Ax(t)+Bu(t));
• u(t)=-F(x(t)-r(t)));（ただしrは平衡点。平衡点であれば、以下の記述はr=0として一般性を失わない）

状態フィードバック †

• \dot{x}=Ax(t)+Bu(t)=(A-BF)x(t));に対して、A-BFの固有値実部が全て負で漸近安定である。
  • 固有値は、時定数の逆数を表す。（証明は実際にxを解くによる。よく固有値を-3[/t]とか-4[/t]とかになるように制御する）

フィードバックゲインの修正 †

制御しなくても何とかなる変数 †

• Aのi列が0ベクトルかつBのi行が0ベクトルならば、変数x_i);は制御不要である。（最適制御のゲイン多様体を見ることでも分かる）
• 制御不要な変数を制御しない方法
  • 真面目にやるには、不要な変数を削除して最適レギュレータを再度組む
    • パラメータが増えるデメリットがある。
• 魔改造的には、制御の必要ないフィードバックゲインを0にすればよい。
  • 例えば、倒立振子ロボットでの速度制御をしたい場合は、フィードバックゲインのうちxのフィードバックを支配する要素を0にする。
  • パラメータを増やす必要がないので、実装も楽。

フィードバックゲインの決定法 †

• 直接法、可制御正順系を介した直接法、アッカーマン法、最適レギュレータ、折り返し法がある。
• アッカーマン法と、最適レギュレータが現実的。これらは関数が用意されている。

直接法 †

1. 可制御性のチェック
2. フィードバック係数ベクトルfを要素で表現して、A-Bfを計算
3. 特性方程式|sI-A+bf|を計算 (1)
4. 目的とする固有値e_i);を恣意的に選び、\Pi (s-e_i));を計算 (2)
5. (1)と(2)の係数比較により、fを算出

可制御正準型を介した直接法 †

1. 可制御性のチェック
2. 可制御正準型に変換行列Tで変換
3. 可制御正準型なので、Aの最下行に特性方程式が現れる (1)
4. 目的とする固有値e_i);を恣意的に選び、\Pi (s-e_i));を計算 (2)
5. (2)-(1)により、可制御正準系でのフィードバック係数ベクトル\tilde{f'});を算出
6. f=\tilde{f}T^{-1});により、fを算出

アッカーマン法 †

1. 可制御性のチェック
2. 目的とする固有値e_i);を恣意的に選び、P(s)=\Pi (s-e_i));を計算 (1)
3. 可制御性行列U_c);を用いて、f=(0 ... 0 1)U_c^{-1}P(A));

• 実装
  • 最適レギュレータで倒立振子を制御するプログラムのparam.mをAchermannで検索（モデルはここを参照、ただしxの次元の順序が2,3で逆になってるので注意）
  • octave run.mでシミュレータが走る。

最適レギュレータ †

1. 可制御性のチェック
2. 状態フィードバックのうち、恣意的に与えたQ, Rに対してJ=\frac{1}{2} \int_0^\infty x^t Q x + u^t R u dt);を最小化するフィードバックゲインは、A^tP+PA+Q-PBR^{-1}B^tP=0);を満たす正定行列Pを用いて、F=R^{-1}B^tP);である。

• 実装
  • 最適レギュレータで倒立振子を制御するプログラムのparam.mをLinear Quadratic Controllerで検索。（モデルはここを参照、ただしxの次元の順序が2,3で逆になってるので注意）
  • octave run.mでシミュレータが走る。
  • モータのモデル化についてはモータへ

折り返し法 †

1. 可制御性のチェック
2. 固有値を折り返すための直線を\mimetex(Re \lambda = - \alpha);で指定する（Aの固有値にかぶらないように）。
3. フィードバックゲインは、恣意的に与えたRに対して(A+\alpha I)^tP+P(A+\alpha I)-PBR^{-1}B^tP=0);を満たす半正定な最大解P_+);（他の解との差が半正定となる解）を用いて、F=R^{-1}B^tP_+);である。

• 備考
  • Q=2 \alpha P_+);とした最適レギュレータに一致。

ゲイン多様体 †

最適レギュレータのゲイン多様体 †

• 最適レギュレータのパラメータのうち1つを固定して他をランダムに振ることで得られた、制御ゲインを多様体とみなして解析することで色々と分かることがある。

• https://github.com/hamko/sample/blob/master/opt_reg/manifold.m
  • 倒立振子の例

• gnuplotで

p "data" u 5:8 #など

• 
  • xが0に張り付く平面を許容するので、xは制御不要である。
  • 一方、xの強い制御をしたければ、thetaゲインもあげなければならなくなる
  • thetaには最低限必要なゲインが存在する。
• 
  • xは制御不要
  • 一方、thetaと同様に、xの強い制御をしたければ、xdotゲインもあげなければならない
  • xdotの制御には最低限必要なゲインが存在する
• 図はないが、x-thetadotも同様。
• 
  • xdotとthetaは強く相関し、パラメータ選定において拘束されている
  • thetaとxdotの魔改造には注意が必要

推定誤差の制御偏差への伝播 †

• 標準偏差の敏感度[cm/deg]=状態推定が角度[deg]ずれしたときに、どれくらい位置[cm]に偏差を残すか
  • =k_theta[N/deg]/k_x[N/cm]

plot "<awk '$7 > 2.2015*$5-220 && $5 < -10 {print}' mani" u 0:(($7/60)/($5/100))

安定ゲイン多様体との比較 †

• https://github.com/hamko/sample/blob/master/opt_reg/stability.m
• 安定領域は，A-bf（fはフィードバックゲイン）の特性多項式の解の実部が全て0未満に相当する．
  • 数学的にはラウス・フルビッツの安定判別法によって計算可能だが面倒なので、フィードバックゲインをランダムで振って根を求めちゃったほうが早い。
  • 当然，安定ゲインは最適レギュレータゲインを内包する．
• x-theta
  • 
• x-xdot
  • 
• x-thetadot
  • 
• x-thetaの制御のトレードオフは，ダイナミクスの構造によるもの
  • ここを魔改造すると，ロバスト性が失われる
• 速度成分は，最適レギュレータによる収束性の制約によるもの
  • ここを魔改造すると，収束性能が失われる．
• 総じて，最適レギュレータの解でゲインが強すぎになってきたら，thetaのゲインを安定限界方向に減らすのが良い．

オブザーバ †

• 完全なオブザーバ付き最適制御
  • 自分で勝手に極をもってきてしまうと、重みとの兼ね合いがよくわからなくなるので、実際のmotorに合わせるのがよい
  • オブザーバが遅くなってしまう場合は、その時点で立ち戻ればよい。ただしその場合は、どれくらいロバストネスを上げるべきか？などといった上位の視点が必要になる
  • あまりにオブザーバが早すぎると過学習してビクンビクンしてしまうので、状態推定の時定数はそんなに早くないほうがよい。

教科書と実機とのギャップ †

操作量uがトルク・力の時のモータでの扱い †

• 基本的には、Duty出力直前まで、マイコン内部でもSI単位系の力[N]で計算して、最終段に力-duty変換をかける。
• モータの実出力トルクは
  1. 電流に比例する→電流センサをつける
  2. 電圧に比例する無回転トルク-回転数に比例する粘性抵抗→これらを分解し、電圧に比例する無回転トルクを操作量、粘性抵抗を損失としてモデル化する

理論と実機の単位の違い †

• 状態・操作量の単位補正行列U_x, U_u);を計算し、フィードバックゲインK);をロボット用のフィードバックゲインK_{r} = U_u^{-1} K U_x);に変換する
  • x=(x, \theta));で表されるとする。x_1);のロボットの内部単位である1mmは、SI単位系では0.001である。x_2);のロボット内部単位である1radは、SI単位系では1である。このような数値を集めた単位補正行列U_x=diag([0.001\,1]));を用いて、x = U_x x_{r});
  • 同様に、トルクのロボットの内部単位である1 dutyは、SI単位系でk[Nm]だとする。単位補正行列U_u=diag([k]));を用いて、u = U_u u_{r});
  • フィードバックゲインの式u=-Kx);は、U_u u_{r} = -K U_x x_{r});なので、u_{r} = - U_u^{-1} K U_x x_{r});
  • したがって、ロボットのプログラミング時のフィードバックゲインは、K_{r} = U_u^{-1} K U_x);と表され、u_{r} = -K_{r} x_{r});である。
  • 状態・操作量はSI単位系や、それに近い（U_x, U_uなどで容易に変換可能）になるように作るべき。U_x, U_uに全ての単位変換を押し付けてはいけない
• K_robot=U_x*Kだが、マイコンに書き込むプログラムでは別の変数として保持し、マイコンのプリプロセッサによって掛け算させるべき
  • ハードウェアアブストラクションと制御器設計を分離するため（U_xも容易に変わる、センサの換装や推定器の変更などによる）。

実機のduty-トルク特性が非線形だったら †

• 上の単位の違いの吸収の式で、U_u);が非線形になるだけ。

x=0以外の目標値 †

• u(t)=-F(x(t)-r(t)));とすればよいだけ。（ただしrは平衡点）

出力オフセット †

• 静止摩擦力はオフセットdutyを使って保証する
• 出力0付近で出力方向が振動して発振の危険があるので，単位時間あたりの出力制限を入れてそれを防ぐ

モデルの切り替え †

• 例えば、倒立振子が空中に浮いたら？
  • どの制御が一番効いているか？xっぽい。
  • 複数モデルについてモデル化するのではなく切り替えなくても、両方でいけるパラメータを探す方法でも可能

LQRの使えない時 †

• 硬い力学系でサンプリング時間も十分短い場合には有効
  • 本質的に振動的
  • サンプリングが粗い
  • 系が高速に切り替わる

Tips †

• 倒立振子の数式とoctaveスクリプト
• 軌道の高速化=タイムステップの変更とみなして、ダイナミクスとLQRパラメータをいじる方法がある [Jur van den Berg+, 2010]