SDRE法による制御器の設計法とその応用・例

SDRE法による制御器の設計法とその応用・例
SDRE法による制御器の設計法とその応用
SDRE
法による制御器の設計法とその応用
法による制御器の設計法とその応用・例
例
鈴 木
聡
東京電機大学
未来科学部・ロボットメカトロニクス学科
(21世紀
(21
世紀COE
COEプロジェクト推進室
プロジェクト推進室))
日本テクノセンター・定期セミナー
2005年
2005
年11
11月
月28
28・・29
29日
日 @ 日本テクノセンター西新宿研修室
本講義内容
1. SDRE法の概要
2 SDRE法の設計概略
2.
3.SDRE法の理論的背景
4. SDRE法の類似手法紹介
法の類似手法紹介
5. SDRE法のTopics
6. SDRE法の拡張例
法 拡張例
I SDRE法の概要
I.
SDRE法の概要
1.. 非線形システム論における位置づけ
非線形システム論における位置 け
◆ SDRE法とは
・SDRE=State Dependent Riccati Equation / 状態依存リカッチ方程式
⇔ 通常のリカッチ方程式は状態非依存
※リカッチ方程式…最適化問題の解を得るために必要な式
・非線形システム制御系設計法
非線 シ
ム制御系設計法
- 理論の箍を多少緩め, 実用性に重きをおいた手法
- 数値解(実用解)を得ることを第一目的にした手法
理論的には非線形システムに対する最適制御系設計問題での,
・理論的には非線形システムに対する最適制御系設計問題での
Hamilton-Jacobi方程式の近似解計算手法
・オンライン-ゲイン可変-フィードバック制御法
オ
イ ゲイ 可変
ドバ ク制御法
◆ フィードバック制御とは
・ 入力を使い,
入力を使 制御対象の安定性/応答性を変える手段
制御対象 安定性 応答性を変える手段
フィードバック制御
入力
誤差
目標
e
+
r
e0
u
制御器
-
出力
y
プラント
(制御対象)
利点: 安定性が強い
不確定性にロバスト
欠点: 応答がやや遅い
(出力)フィードバック
フィードフォワード制御
yr
入力
目標
r
u
制御器
出力
y
プラント
(制御対象)
利点: 応答が速い
欠点: プラントの正確な
モデルが必要
◆ 制御設計手順
・ 通常のフィードバック系設計手順
オンライン
目標
標
+
r
オフライン
入力
誤差
e
-
u
制御器
出力
プラント
(制御対象)
y
プラントモデル
パラメータ調整
ラ
タ調整
(出力)フィードバック
制御器設計
・ SDRE系設計手順
制御器
オンライン 設計
e
+
r
-
パラメ タ調整
パラメータ調整
プラントモデル
制御器
プラント
(制御対象)
オフライン
出力
y
状態依存プラント
状態依存プラ
ト
モデル
状態依存パラ
メータ
◆ 非線形システムとは
・線形システム
- 「重ね合わせの原理 」が成立するシステム
- 時不変, 連続時間
出力
入力
ui
線形
システム
yi
出力
入力
au1  bu2
線形
線
システム
ay1  by2
i  1, 2
・非線形システム=線形システム以外全て
… 不連続,
不連続 時変,
時変 非線形関数,
非線形関数 不確定要素,
不確定要素 etc,,,
- 静止/動摩擦... F    sgn( x )
(機構系)
B T 1  T0 1 
経年変化,, 温度変
温度変化… R  Ro e
((化学プラント)
学 ラ
)
- 経年変
- 三角関数, 幾何学的特性… mgl sin 
(ロボット)
- ノイズ
非線形システムのクラスは非常に広い!!!
◆ 非線形システムの制御法
・対象システムの非線形性を分析する必要あり
((a)) (非線形)モデルがあるか、
(
)
ないか? ...モデルベースド可否
(b) 動作点は狭いか、広いか? ...線形化近似可否
(a) モデル(ダイナミクス)がある場合
(a).
-近似線形化…1次近似で線形(化)モデルを導出 →線形理論適用
-厳密な線形化…モデルによる座標変換
厳密な線形化 モデルによる座標変換 →線形理論適用
線形理論適用
(b). モデルがないor立式困難な場合
- (簡易)制御器 + 実時間調整 (適応制御, オンライン最適化)
- 入出力応答設計 ((ex. ファジイ, ニューラルネット))
◆ SDRE法の制御理論としての分類
様々な制御法の亜流?
- 最適制御問題(LQR)の一種
状態フィードバック則として解を得る
( ィ
(フィードバックゲインを可変にする)
ックゲインを可変 する)
- モデルベースト制御の一種
制御対象のモデル(運動方程式)が必要
- 実時間最適化手法の一種
計算機パワ はほどほど必要
計算機パワーはほどほど必要
- リカッチ方程式の解を実時間で求めるため
- ただしサンプリング時間とトレードオフの関係
◆ SDRE法の適用対象
SDRE法適用のための条件(基本的に)
(a) 連続時間システム
x (t ), u (t )  C 
(b) プラントの数式モデルが存在
x (t )  f ( x (t ), u (t ))
“基本的に”の意味
(a)に対して…
- SDRE法は、元々、理論的困難性(大域安定の補償、解の一意性
など) を度外視して、実用解を得ることを優先した手法。
- 理論的補償なくとも制御入力を計算することは可能という意味で。
理論的補償なくとも制御入力を計算する とは 能と う意味
(b)に対して…
- プラントのモデル自体も設計要因として調整できるので、
プラントの デ 自体も設計要因として調整できるので
状態に依存させプラント風に表現できれば、SDREの適用は可能
◆ 最適制御問題(最適レギュレ
最適制御問題(最適レギュレータ問題:LQR)とは
タ問題:LQR)とは
- 評価関数(基準となるものさし)を考え、その最小化問題として
制御器を導出する方法



T
- 2次形式評価関数(基準):
(
) J   x(t )T Q
Qx
(
t
)

u
(
t
)
Ru (t ) dt
0
- 制御対象のモデル: x (t )  A  x (t )  B  u (t )
- 手段=リカッチ方程式を解く: AT P  PA  Q  PBR 1 B T P  0
P
1 T
- 解=状態フィードバック則: u (t )   R B P  x(t )
x(t ), u (t )
J
t
◆ モデルベ
モデルベースト制御とは
スト制御とは
・制御対象の動特性を記述する(数式)モデルを元に制御器を導出
⇔入出力応答ベースとは異なる
- システムのダイナミクス情報を利用した設計法
「制御対象のモデル例」 x (t )  A  x (t )  B  u (t )
- モデルの良し悪しが制御器の善し悪しに影響
- 制御対象自体の特性が変動した場合に対応困難
(パラメータ変動には弱い...)
- システム構造の特徴を活用した設計法が考案
例:逆モデル制御, モデルによる厳密線形化(ロボット)
t
◆ 実時間最適化手法
・プラントの状態を実時間で計測し、プラントモデルを常に修正、
オンラインで最適化問題を逐次解くことで入力量を計算
・計算機パワーはある程度必要(制御周期時間とトレードオフの関係)
区間外について
は感知しない
障害物
現在から有限時間先
の区間で入力最適化
考慮する区間内で最適と思
われる入力を探索し、その
初回のみ使用
…続き
評価区間
プラント
最終状態
解の形式
LQR
x(t )
LQR
t=0~∞
不変
x → 0 (t →∞)
時不変/FB則/式
RH
x ( )  0
実時間最適化(モデル予測制御)
t = 0 ~ 有限
時変
仮定不要
逐次変更/入力量
x(t )
時刻 t’’ においてT
だけ未来までの制
御を最適化
を信じて入力則
を導く
t
u  Fx
t
T
実際のプラントへの
入力はuの初回のみ
t
u  {u (t ),  u (T )}
I SDRE法の概要
I.
SDRE法の概要
2. 実装例紹介...
...倒立振子の振上げ
倒立振子の振上げ//安定化制御
◆ SDRE法による広域安定化制御
A ngles of link
2
モータ
[rad]]
1
0
-1
Hit!
-2
Hit!
-3
振子部
0
2
4
6
Input
8
10
12
0
2
4
6
[s]
8
10
12
2
[N.m]
1
0
-1
-2
・350MHz PC/AT +Windows98, RT-MaTX(http://www.matx.org/)
・制御周期2msec
Riccati方程式算出: Schur分解
・Riccati方程式算出:
◆ ‘倒立振子’=メカニカルシステム基本形
倒立振子 メカニカルシステム基本形
・ 慣性, 粘性/コリオリ力/遠心力, 重力を含む
0 
0 
 

M ( )     N ( , )     g ( )   


0
1 
1 
 : [  0 , 1 ]T …2次の最も簡単な系
 J 0  m1L20  0.5m1l12 (1  c11 )  m1l1L0c1 
M 

 m1l1L0c
J 1  m1l12 

 m1l12 s111  d 0 m1l1L0 s11 
N 

2

d1 
  0.5m1l1 s11 0
g  0
m1l1 sin 1 g  …M,N,Gが状態
T
θによって変化
◆ 倒立振子モデルの応用例
・ 歩行ロボット(受動歩行)
hip
stance-leg
m
swing-leg
h
m

st
sw
st

sw
M ( )  N ( , )  G ( ,  ) g  
 : [ st , sw ]T  R 2
 : [ st ,  sw ]T  R 2
…劣駆動
…ロボットの一般表現系
ボットの 般表現系
l1
 0 , 1
angle
J 0 , J 1 inertia
m1
mass
d 0 , d1 friction
L0 , l1 length
0
1
J1, m1, d1
J 0 ,d0
L0
◆ なぜ
なぜ‘振上げ/安定制御’か?
振上げ/安定制御 か?
・ 非線形 / 不安定システムの大域安定化問題のベンチマーク
・First
Fi t application
li ti
M i Ni hih F t (76)
Mori,Nishihara,Furuta(76)
・Pseudo State Feedback
Furuta(92)
Energy control
Chung(95), Astrom(96)
・Energy
・FB linearization
Saeki(99), Spong(acrobot:94)
・VSS/Bang-Bang
Yamakita (96),Furuta(92)
◆ 難しさの要因
1
・ 安定平衡点→不安定平衡点への状態遷移
- プラントなどのモ
プラントなどのモード移行
ド移行
・ 特異姿勢の通過
◆ なぜ,
0
この成果がすばらしいのか?
(i) `制御器’は1つ
制御器 は
⇔従来法 = 2つの制御器を切替
① 振上げ制御系 : 非線形コントロ
非線形コントローラ
ラ
② 倒立安定化制御系 : 線形コントローラ
振上げ
制御器
倒立安定
制御器
振上げモード
安定化モード
- 振上げ制御系&倒立安定制御系を別々に設計
- 切替タイミングの指針なし,
切替タイミングの指針なし 切替前後の安定補償困難
I. SDRE法の概要
SDRE法の概要
法 概要
3.SDRE法の実装手順概略
3.SDRE
法の実装手順概略
◆ SDRE法の要点
・状態依存型リカッチ方程式(SDRE)法とは
状態依存型リカッチ方程式(SDRE)法とは
状態に依存してリカッチ方程式を実時間で解き,
毎サンプリング時間で制御入力を更新する方法.
=>実時間で最適化問題の解を更新する方法(大雑把)
実時間で最適化問題の解を更新する方法(大雑把)
もう少し理論的にいうと,
もう少し理論的にいうと
SDRE法=ハミルトンヤコビ方程式(HJE)の近似解を得る代替手法
HJEってなんでしょう?
◆ 従来の
従来の一般的な制御器設計法
般的な制御器設計法…線形化モデルを使用
線形化モデルを使用
(1) x  f ( x , u ) の動作点近傍 ( x , u )  ( x e , u e ) のみ考え、
(2) 動作点周りの線形傾き(テイラー展開の1次項係数)


B
f ( x, u )
A
f ( x, u )
x  xe
x

x
を用
を用いて、
u
x
u u
u u
e
e
e
(3) 動作点まわりの簡略モデル((=線形化近似モデル)を導出し、
~
x  A  ~
x  B  u~ (~
x  x  x , u~ u u )

e
e
(4) 線形システム制御理論を
用いて、制御器を設計
u   g ( x)
f ( x, u e )
x
f ( x, u )
x  xe
u ue
u  ue
に固定
xe
x
◆ 線形化モデル設計の限界と実際
・動作点(線形化点)近傍でのみ、線形化モデルは有効
→ ( x , u )  ( x e , u e )以外では本来は安定性は保証されない
・しかし、LQR自体ある程度のロバスト性を持つので,
- 状態量の変動差(非線形性)が小さいか,
- 動作点まわりでの安定化のみしか気にしない,
ならば実用上は耐えることが多い
・動作点を広げるには?
- 複数 Ai , Bi に対して設計した制御器を切替
に対し 設計した制御器を切替
…マルチコントローラ, ゲインスケジューリング
- A, B を状態に依存させて変化
を状態 依存さ
変
→ A ( x ), B ( x )
…SDRE法
◆ 最適レギュレ
最適レギュレータ問題(Linear
タ問題(Linear
Quadratic Regulator: LQR)
・ 設計手順と操作内容
①制御対象のモデルを準備 A, B (  const )
②2次形式評価関数で重み調整: Q ( 0), R (  0)
③リカッチ方程式を解く AT P  PA  Q  PBR 1 B T P  0
④FBゲイン算出 F  R 1 B T P
⑤状態FBで閉ループ系を安定化 x  ( A  BF ) x
P
・ 逆にいうと, A  BF が安定になるような F を求める作業
・ 閉ループ系が必ず安定化されることが補償されている
閉ル プ系が必ず安定化されることが補償されている
- 安定性に関しては気にせず、性能調整に集中できる
- ただし線形時不変・可制御システムに対して.
… では時変システムに対しては?
◆ 時変システムへの展開
(SDRE法へ)
・状態依存という時変システムに拡張する
x  A  x  B  u
x  A( x )  x  B ( x )  u
- 本来の時変システム= 時間に応じて変化するシステム
{ x  A(t )  x  B (t )  u } > { x  A( x )  x  B ( x )  u }
・LQR
Q 設計は, 様々な仮定/条件の下に導かれた解
→単純に{ A, B}  { A( x ), B ( x )}と置き換えて, 大丈夫?
) B ( x ) を凍結する”
・”各時刻 t 毎の状態 x (t ) に対して、係数
に対して 係数 A( x ),
という概念を導入
◆ SDRE法の計算アルゴリズム
・毎制御周期時間でリカッチ方程式を解く
1. x  x(t) を計測 (場合によっては推定)
2. SDLRの係数行列
係数行列 A(x), B( x) を算出
3. 状態依存重み Q(x), R( x) を計算
4 SDREを解き, 正定解 P(x) を算出
4.
AT ( x)P( x)  P( x) A( x)  P(x)B(x)R1(x)B( x)P( x)  Q( x)  0
5 制御入力(時刻tでの)を算出
5.
u' (t )   R  x  B  x  P  x   x
6. 時刻 t=t~t+Δtの間 u’(t)をプラントへ出力
7. 1へ戻る.
1
T
◆ SDRE法の計算アルゴリズム
状態依存
仕様数値化
オンライン
状態量固定
+
r (t )
x (t ) サン
- プラ
t
x : r  y
x [t ]
Q(x), R(x)
計算
リカッチ解
P((x) 計算
A(x), B(x)
計算
入力
計算
u  u(x[t])
u (t )  u ( x )
t  [t , t   t ]
ホー
ルダ
プラント
(制御対象)
y (t )
状態依存線形
モデル係数計算
) B(x) の数式表現はオフラインで求めておく
・ プラント A(x),
・ 制御指標重み行列 Q( x), R( x) の数式表現も予め決めておく
・ それらに従い, 実時間で最適入力 u ( x ) を算出し, t = t’の入力値だけを使用
u ((tt )  u ( x )
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
1.設計の
1.
設計のkeypoints
keypoints
~本章の方針~
SDRE法を理解するために必要な、制御理論の基礎的概念
の説明を交えて進めます
◆ SDRE制御器の実装アルゴリズム
1. x  x(t) を計測 (場合によっては推定)
) B( x) を算出
2. SDLRの係数行列
係数行列 A(x),
A(x), B(x)の算出?
3. 状態依存重み Q(x), R( x) を計算
4. SDREを解き, 正定解 P(x) を算出
Q(x),
) R(x) の選定?
AT ( x)P( x)  P( x) A( x)  P(x)B(x)R1(x)B( x)P( x)  Q( x)  0
5. 制御入力(時刻tでの)を算出
実時間で解くには?
u' (t )   R  x  B  x  P  x   x
6. 時刻 t=t~t+Δtの間 u’(t)をプラントへ出力
7 1へ戻る.
7.
大域安定性の補償は?
1
T
◆ 制御器設計の
制御器設計の一般論
般論
1. モデリング…システムを数式で表現
シミュレータ構築の為
シミ
レ タ構築の為. 制御器設計の為(簡略版)
2. パラメータ同定…数式モデルの係数を, 実機データから推定
3. 設計仕様/指針決定…どのくらいの応答性, ロバスト性
4. 適用理論選択
5 制御器設計…制御器パラメータ調整
5.
制御器パラメ タ調整, シミュレーション
シミ レ ション
6. 実装…状態推定, 実装数値計算, 不確定性(ノイズ, 未知要素)
7. 実機検証/評価…結果に応じ1~6へ(繰返し)
SDRE法では 1 → A( x),
) B( x)
3 → Q( x), R( x)
6 → P((x)
と密接に関係
◆ SDRE制御器設計のキ
SDRE制御器設計のキーポイント
ポイント
(1) 状態依存線形表現(SDLR)の導出(ゼロ割,速度情報)
x  A( x) x  B( x)u の係数行列 A( x), B( x) の算出
(2) 状態依存重みの設定

V ( x )   ( x T Q ( x ) x  uT R ( x )u ) dt の Q( x), R( x)
0
(3) リカッチ方程式の実時間求法
AT (x)P( x)  P(x) A( x)  P( x)B( x)R1(x)B( x)P(x)  Q(x)  0
の P(x) の算出
制御系問題としての以下の基礎が必要
1.モデリング, 2.安定性, 3.フィードバック, 4.最適制御
…次以降で, 概略説明します
◆
SDRE法の為のモデリング 状態依存線形表現(SDLR)
SDRE法の為のモデリング…状態依存線形表現(SDLR)
・オフラインでの準備その1
・基本方針: ダイナミクスモデルから, 入力affine-1回微分形へ変形
x  f ' ( x,u)
x  f ( x)  g( x,u)
x  f S ( x)  x  g S ( x,u)  u
ただしf S ( x), g S ( x,u ) はsmooth
x  A( x) x  B( x)u
では、ダイナミクスモデル x  f ' ( x,u ) の求め方は?
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
2.モデリングの基礎
2.
モデリングの基礎
◆ モデリングの意義&必要性
・モデリングとは…システムを数式で表現すること
= システムの特徴/振る舞いの事前情報を得ること
・予測(シミュレーション), 特徴解析, 制御系設計に利用
・システムとは…目的とする機能を満足するように,互いに関連
する要素の集まりが作用しているもの。
LOS
LOS state
96
G1 ( s )
Y ( s)
15
24
13
ctrl.of
movement
49
42
left
34
q4
6
q6
q6
q5
mouse
q6
stop
11
q4
left
70
q5
direct
111
q2
down
q3
33
backward
73
up
11
245
9
head
mouse state
q1
right
9
q1
186
com. board 43
53
システムモデル例1
-ブロック線図-
8
8
9
q3
forward
q2
stop
27
head ctrl.
+
R(s )
68
71
q4
80
camera
+
G1 ( s )
robot
robot state
84
q1
q3
50
right
left turn
q9
17
12
11
q10
88
q7
forward
fo
wa d
9
18
17
q8
49
right turn
backward
システムモデル例2
- 状態遷移モデル -
◆
制御工学でのシステム表現
•外部から見ることができる独立した変数を入力、従属した変数を
出力と言う すべての変数は 時間的に変化し この変化は情報
出力と言う。すべての変数は、時間的に変化し、この変化は情報
をあたえるので、その情報的な側面を捕らえて、信号という。
•現在の入力の影響が将来の出力に現れるものを動的システム
と言う。通常入力は原因であり、出力は結果として現れる。この
場合は因果性が成立すると言う
場合は因果性が成立すると言う。
input: u
system
output: y
state: x
x  f S ( x)  x  g S ( x,u)  u
y  Ix
-入力、出力、状態をもつ
-信号の流れを作る
-積分/積分要素を含む
◆
物理システムのモデリング
・ニュートンの運動法則(慣性の法則, 運動の法則, 作用反作用)
u
→ 並進運動方程式
d2
m
dt
2
m
yu
y
・オイラーの運動法則
→ 回転運動方程式
d2
J 2   Fl  
dt
J

l
l
F
F  l  力 長さトルク 
- 慣性モーメント
慣性モ メント: 回転しにくさを表す量
(並進運動における質量に相当)
◆
ニュートンの運動の3法則
ニュ
トンの運動の3法則
1. 慣性の法則
物体は外力からの作用を受けない限り、現在の運動状態を
保持する。
2. 運動方程式
運動量の変化量とその方向は 外力に比例する
運動量の変化量とその方向は、外力に比例する。
3. 作用
作用・反作用の法則
反作用の法則
力の作用とそれから生じる反作用は常に逆向きで、大き
さは等しい
◆
インピーダンスモデルへの変形
インピ
ダンスモデルへの変形 / 行列表現化
・ベクトル/行列表現で簡素化する
例: 2質点系
・例:
u1
1 質量(仮想質量を含む)の点に
1.
変位を表す記号をつける。
m1
f1
k1
y1
2 質量にかかる力を考える。
2.
質量にかかる力を考える
m2
f2
k2
m1 y1  u1  f1 ( y1  y 2 )  k1 ( y1  y2 )
y2
m2 y2   f 2 y2  k2 y2
 f1 ( y1  y2 )  k1 ( y1  y2 )
3. 変数ごとに係数をまとめる
m1 y1  f1 y1  k1 y1  u1  f1 y 2  k1 y2
m2 y2  ( f1  f 2 ) y 2  (k1  k 2 ) y2  f1 y1  k1 y1
4. 変数をベクトルでまとめ, 行列を使って一つの等式にまとめる。
C
m1
 0

0   y1   f 1

m 2   y2    f 1
M
   k 1
  k1

 f 1   y 1 
f 1  f 2   y 2 
 k 1   y1   u1 

k 1  k 2   y 2   0 
K

インピーダンスモデル(2階微分方程式)
M    C    K    u
u
◆
状態空間モデル
・1階微分行列方程式…制御系設計に便利な表現
状態量:x 、入力量:u、出力量:y、の時間変数
・状態量:
 x  A  x  B  u

y  C  x  D u
y に出る.
u から入って
から入って…
D
+
B
+
x (t )

x(t )
C
+
u (t )
y (t )
+
A
状態空間表現モデルのブロック線図表記
◆
2階微分形式→1階微分形式モデルへ
・インピーダンスモデルから、状態空間モデル表現への変換
クトル表現を利用する
・ベクトル表現を利用する
・逆行列の存在性に注意する
- インピーダンスモデル(2階微分)
M  s  C  s  K  s  u
 
x :  
 
- 拡張ベクトル導入
- 状態空間モデルの状態方程式(1階微分)
I     0 
  0
  1 u

  
1
1   
    M K  M C     M 
x  A  x  B  u
◆
運動方程式の求め方
・ニュートン運動方程式法 = 力・トルクの釣り合い関係
…システムが複雑になると不便
が複雑 なると 便
・ラグランジェ法…エネルギーの釣り合い関係から導出
…システムが複雑でもエネルギ計算は容易
が複
も
ギ 算
並進の運動エネルギー
y
y
x
回転の運動エネルギー
1
m( x 2  y 2 )
2

m
l
1
1
1
m (l) 2  m l 2 2  J 2
2
2
2
x
重心位置がわかればよい。
重心位置の微分=速度 →速度の2乗 x慣性=エネルギー
ネルギ
…重心位置の微分
◆
ラグランジェ法
L  T U
d 
L R
( L) 
u

Lagrange方程式
dt y
y y
Lagrange関数
u
運動エネルギー
運動エネルギ
m
k
f
y
ポテンシャルエ
ネ ギ
ネルギー
散逸エネルギー
1
m y 2
2
1
U  k y2
2
1
R  f y 2
2
T
上記の例の場合Lagrange
g g 方程式を解くと… my  kyy  fy  u
◆
ラグランジェ法の原理
1 2 2
ml 
2
ポテンシャルエネルギ U  l (1  cos  ) mg
運動エネルギー

l
T
・外力を加えなければエネルギーの変化
はないので、総エネルギーの微分は0

mg
mg sin 
1 2 2
ml   lmg (1  cos )
2
d
V  ml 2  lmg sin   (ml 2  lmg sin  )  0
dt
l   g sin  なる微分方程式が得られる。
V  T U 
◆
例…台車型倒立振子
(A)ニュートン運動方程式による場合
(A)ニュ
トン運動方程式による場合
振子長: 2l
I  Vl sin   Hl cos  c
V
Vl sin
i 
u

( x, y )
mg

H

(1)
振子・回転運動
振
d2
d2
m 2 y  m 2 (l cos  )  V  mg (2)
dt
dt
振子・並進運動
(垂直方向)
d2
d2
m 2 x  m 2 (  l sin  )  H (3)
dt
dt
振子・並進運動
(水平方向)
M  u  H  f
(4)
一般に
般に
台車・並進運動
d2
d
m 2 (l sin  )  ml
((cos )  ml ( sin  2  cos)
d
dt
d
dt
が成立するので..
(3)式は H  m  mll ( sin  2  cos)
(3’))
(3
d2
d
( sin )
(2)式の中項は m 2 (l cos  )  ml
dt
dt
 ml ( cos 2  sin )
 V  mg (=(2)式右項)
Vでまとめると
V  ml ( cos  2  sin )  mg
(2’)
(1)に(2’)と(3’)を代入して, (中間変数である)V,Hを消去すると


I  c  ml ( cos  2  sin )  mg  l sin 

相殺

 m  ml sin  2  ml cos   l cos 
 mgl sin   ml 2( s2  c2 )  ml cos 
 mgl sin   ml 2  ml cos 
よって
( I  ml 2 )  c  mgl sin   ml cos 
(5)
一方(4)式と(3’)からHを消去して
M  u  m  ml sin  2  ml cos  f
( M  m)  f  ml cos   u  ml sin  2
(6)
・ニュートン運動方程式による手法(手順まとめ)
①個々の質点ごとに並進運動、回転運動の方程式を立て、
②
②異なる質点間の作用・反作用力を中間変数を用いて表現
表
③中間変数を代数計算で消去
・中間変数を都度設定する煩わしさがある。
・中間変数を都度設定する煩わしさがある
(B)ラグランジェ法による場合
振子の重心位置は ( x, y )  (l cos  ,   l sin  ) であるから,
・運動エネルギー
運動エネルギ
2
2
1  2 1  d
 d
  1  2
T  I  m (l cos )    (  l sin  )    M
2
2  dt
  dt
  2





2
2
1
1
1
 I 2  m  sin l 2  cos l 2  2l cos   2  M 2
2
2
2
・位置エネルギ
U  mgl cos 
・消費エネルギR  1 f  2  1 c  2
2
2
振子長: 2l
Vl sin
i 
u
V

( x, y )
mg

H

L  T U 
1 2 1
1
I  m(l 2 2  2lc   2 )  M 2  mglc
2
2
2
y : [ ,  ]T とおくと
d  L  d  L  d

     (mlc   m  M)
dt  y1  dt    dt
 ( M  m)  ml (c   s  2 )
d  L  d  L  d 

     ( I  ml 2  mlc )
dt  y 2  dt    dt
 ( I  ml 2 )  ml c   s 

L L

  ml  s  mgl  s
y2 

であるから….
Lagrange方程式に代入して,
d  L  L R


 
dt     
→ ( M  m)  f  ml cos   u  ml sin  2
i=1のとき
→ ( M  m )  ml ( c   s  2 )  f  u
i=2のとき d  L   L  R  0


dt   

(6)と同じ

相殺
2 




→ ( I  ml )  ml c   s  ml  s  mgl
g  s  c  0
(5)と同じ
→ ( I  ml 2 )  c  mgl sin   ml cos 


・幾何学的な座標計算と, 簡単な微分計算/代入操作のみで
求まるため導出手順は簡単.
◆
m自由度アーム運動方程式
m自由度ア
ム運動方程式
・ロボットやメカニカルシステムに適した表現法
- 物理的に意味のある項でまとめる
⇔インピーダンスモデルは状態量の微分次数毎でまとめたもの
・アーム運動方程式一般形
ア ム運動方程式 般形
M ( )    N ( , )    g ( )  
: 状態
  Rm
: 一般化入力
  Rm
M ( )  R mm : 一般慣性(質量, 慣性)項
N ( , )  R mm: 非線形項(粘性, 遠心力, コリオリ力)
g( )  Rm1
: 重力項
◆ 補足…行列計算
SDRE法を含む状態空間表現ベースの計算では、行列演算が
必須なので念のため復習
・ 一般に非可換: AB  BA
・ 対称行列:
・ 転置の分配:
AT  A
( AB)T  B T AT
・固有値:特性方程式の根
- | sI  A | s n  an s n 1  an 1s n  2    a2 s  a1  0
- 正方行列Aについてのみ
- 行列式で計算される
- 安定判別などに必須
安定判別など 必須
◆ 補足…行列式とは
・ 行列式(determinant)
n
| A |  (1) i  j  aij  M ij
j 1
aijは行列 Aの(i,j)成分
・ 小行列: M kl …行列のk行とl列を除いてできた (n  1)  (n  1) 次元の
行列、の行列式
行列、
行列式
例
 a11 a12 
 (1)11 a11  M 11  (1)1 2 a12  M 12
a

 21 a22   a  | [a ] | a  | [a ] | a a  a a
11
22
12
21
11 22
12 21
・ 固有値の算出や、逆行列計算仮定で必要
◆ 補足…逆行列の計算
・ AX  I を満たす正方行列 X を Aの逆行列という
・ A1と記述する
adj ( A)
(ただし | A | 0 )
| A|
・ adj
( dj i t matrix)
ti )
dj ( A):余因子行列(adjoint
i行j列の要素 Aij が Aij  (1)i  j  M ji で与えられる行列
・逆行列の定義: A1 
1 2
 a a   11
 | [a ] |  [a12 ] |  a22  a12 
・例: adj   11 12     (1) M 11 (1) M 21    22
   a a 
 a a  (1) 21 M (1) 2 2 M

|
[
a
]
|
[
a
]
|
21
11  
21
11 
12
22  
  21 22   
・ 逆行列の分配:
逆行列の分配 ( AB) 1  B 1 A1 (ただし| A | 0, | B | 0)
・ 合成分配:
( AT ) 1  ( A1 )T
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
3.状態依存線形表現
3.
状態依存線形表現
State--Dependent Linear Representation (SDLR)
State
◆
SDRE法設計でのモデル
「(a). モデル(ダイナミクス)がある場合
- 近似線形化…
近似線形化 1次近似で線形(化)モデルを導出
次 似 線
デ を導 →線形理論適用
- 厳密な線形化…モデルによる座標変換 →線形理論適用
」
の折衷的手法に相当
折
近似線形化 のようですが、近似動作点を変化させます
・ “近似線形化”…のようですが、近似動作点を変化させます
・ “厳密な線形化”…のように広範囲の領域での均一化を図ります。
f (x, ue )
f (x, ue )
f (x, u)
f (x, u)
x
xe
◆
x
xe1 xe2 xe3
x
x
メカニカルシステムの状態依存線形表現(SDLR)導出
目標: M ( )    N ( , )    g ( )   ,  ,  R m
x  A( x) x  B( x)u
・手順1: 状態ベクトルを拡張し, 2回微分形を1回微分形へ
k
m k
x : [  T , T ]T  R n 2 m   Lu, u  R , L  R
・手順2: 加速度を, 速度/角速度/入力からの計算式へ変形
…((A))
   M 1 ( ) N ( ,)    g ( )  


・手順3: さらに, 重力項を角度の線形一次形に変形
g ( )i
G : [ gij ], gij :
, i, j  1,, m
g ( ) i は g ( )のi行
j 
…(B)
( )
続き..
・手順
手順4:
4 1階微分行列形式にまとめる
(A)式→    M 1 N    M 1G    M 1 L  u
0
I
    0 

  
    M 1 ( )L  u
1
1


M
(

)

G
(

)

M
(

)

N
(

,

)
   

  
→  
SDLR
I
0


 0 
x  

x

M 1 ( )L  u
1
1
)

(

)

(

)

(

)

(

,

M
G
M
N




: A( x ) x  B ( x )u
keypoints - 状態量と入力の線形一次形に変形すること
- ゼロ割を回避するように工夫すること
- 数値的には可制御性を維持すること
数値的には可制御性を維持する と
(真の特異姿勢はサンプル点間に擬似的に納める)
◆
SDLRの具体例…回転型倒立振子
M ( )    N ( , )    g ( )   ,  ,  R 2
  [ 0, u ]T
 : [  0 , 1 ]T
l1
 J 0  m1L20  0.5m1l12 (1  c11 )  m1l1L0c1 
M 

J 1  m1l12 
 m1l1L0c1

 m1l12 s111  d 0
N 
2

  0.5m1l1 s11 0
g  0
m1l1L0 s11 

d1

m1l1 sin 1 g 
 0 , 1
angle
J 0 , J 1 inertia
m1
mass
0
d 0 , d1 friction
L0 , l1
length
T
問題点 : 平衡点 1  n では0になる
1
J1, m1, d1
J 0 ,d0
L0
I
0


 0 
x  

x


M 1L  u
1
1

M

G

M

N




x : [  0 , 1 , 0 , 1 ]T
g ( ) i
より
より、


j
0

sin
i 1 
m1 g l1 
1   
] g ij :
先の(B)式に従って変形すると
式に従 て変形すると G : [ g ij ],
0
G  0


0


g

m1 g l1  sin 1 
次の補正式を導入
1
| 1 |  …最上点近傍補正:理論値と近似値との整合化

2

00(|1| )
h(1 ) :  sin 1  0.01e 100
else …最下点で(微小)有限値化





1
0
G 
0
SDLR形式が求まる
x(t)  A( x) x  B( x)u
◆

m1 gl1  h(1 )
0
SDLRの解釈
解釈
・通常の線形化…現在の状態近傍での局所的勾配を計算
・SDLRによる線形化…現状態と収束目標状態との大域的勾配
x
原点(安定化目標状態)
x

f ( x, u )
x  xe
x
u u
x
x
A( x )  f ( x , u ) / x
e
x  f (x, u)
テイラ展開1次線形化の概念図
x  f (x, u)
SDLRによる線形化の概念図
◆ SDRE制御器の実装アルゴリズム
1. x  x(t) を計測 (場合によっては推定)
) B( x) を算出
2. SDLRの係数行列
係数行列 A(x),
3. 状態依存重み Q(x), R( x) を計算
4. SDREを解き, 正定解 P(x) を算出
AT ( x)P( x)  P( x) A( x)  P(x)B(x)R1(x)B( x)P( x)  Q( x)  0
5. 制御入力(時刻tでの)を算出
u' (t )   R  x  B  x  P  x   x
6. 時刻 t=t~t+Δtの間 u’(t)をプラントへ出力
7 1へ戻る.
7.
1
T
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
4.状態依存重み行列の決定指針
4.
状態依存重み行列の決定指針
Q ( x ) x  u T R ( x )u dt
0
における, 重み行列 Q ( x ),
) R ( x ) の選定
・2次形式風評価関数 J 
 x

T
一般的方針
般的方針
・ Q ( x )  0, R ( x )  0
・制御対象の特性に応じて選定
制御対象の特性に応じて選定
・LQRのように直感的に各状態への依存を考えて選定
・代表的な姿勢時の重みを連続的に繋げる関数を選ぶ
代表的な姿勢時 重 を連続的 繋げる関数を選ぶ
(Q,Rは連続という仮定から)
倒立振子の場合
・振上げ/安定で重要な変数=振子角度 1 に依存させて変化
・倒立安定可能なLQR重み Q , R ( const ) へ滑らかに変化
通常のLQRで準備調整
・不安定平衡点
不安定平衡点(振子角度 1  0 の時)での線形化モデル+LQR
で制御器設計
1
・最上端での重み行列を決定
Q  diag( 5000 , 3000 , 1, 1000 ),
R  10
0
全域への拡大
・減少させたい位置状態, 動きを抑制したい速度状態に対し,
大きい重み値となるように関数形状を(状態依存で)調整
①
①LQRから
から Q , R相当
Q( x)  diag(
g( q1 , q2 , 1, q4 )
1
2000
q1  1 
2.5
2
1.5
1
0.5
0
1
0.5
0
1
0.5
0
1
0.5
0
2
5000
②: 振上げ動作を重視させる為
1 e
0
3
2.5
2
1.5
1000
③b: 振子の速度抑制
4
q
100
500
0
1 e
10(1  / 2 ) 2
1 e
50(|1| / 6 ) 2
10
3
2.5
2
1.5
④特異姿勢では入力抑制
100
R(x)
r2 
3
10000
10(|1| / 6 )
R( x)  0.2  r1  r2
r1 
0
q
q4  1 
103
③a: モ
モータの移動抑制
タの移動抑制
4000
q
5000
1  e10(|1| / 9 )
q2  3 103 (2  cos 1 )
we igh t fu n c tio n s
6000
50
0
3
2.5
2
do wn < -
1.5
 1 [rad]
最下点
-> up
最上点
◆ SDRE制御器の実装アルゴリズム
1. x  x(t) を計測 (場合によっては推定)
) B( x) を算出
2. SDLRの係数行列
係数行列 A(x),
3. 状態依存重み Q(x), R(x) を計算
4. SDREを解き, 正定解 P(x) を算出
AT ( x)P( x)  P( x) A( x)  P(x)B(x)R1(x)B( x)P( x)  Q( x)  0
5. 制御入力(時刻tでの)を算出
u' (t )   R  x  B  x  P  x   x
6. 時刻 t=t~t+Δtの間 u’(t)をプラントへ出力
7 1へ戻る.
7.
1
T
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
5.実時間リカッチの求解
5.
実時間リカッチの求解
・リカッチ方程式の解を数値的に求める方法
- 数値積分法…微分リカッチ方程式の数値積分(繰り返し計算)
- 有本-ポッタ
ポッター法
法…固有値分解法.
固有値分解法 ハミルトン行列の応用
ミルトン行列の応用
- Schur分解法…固有値分解法の改良版. 実数範囲のみで計算可能
- ニュートン法…リカッチ方程式の線形化版を繰り返し計算
・SDREでは実時間でリカッチ方程式を解かねばならない.
- 各制御周期時間内で計算が終了する必要あり.
- 繰り返し計算はその点で不向き
★ リカッチ方程式については, システムの安定, 2次形式などの
概念を理解した方が良いため, 後ほど説明いたします.
II. SDRE制御器の設計
SDRE制御器の設計
制御器 設計
6.実装上の問題点の回避策
6.
実装上の問題点の回避策
◆
考えられる問題点
・SDLR(状態依存線形化モデリング)に於いて…
- 0割の回避…分母に微小値追加,
分母に微小値追加 数学理論式の置き換え
- 制御対象固有の不連続点回避…真の特異姿勢を回避
- 速度情報のsmoothing
thi …制御対象モデルの振動的変動回避
制御対象モデルの振動的変動回避
・状態依存重み行列Q(x),R(x)選定に於いて…
- 理論的必要条件…連続性, 正定性
- 補間性…代表動作点での制御則=線形制御器設計
◆
具体的回避法…倒立振子の場合
・SDLR(状態依存線形化モデリング)に於いて…
- 0割の回避…分母に微小値/数学理論値の置き換え
1
0.9
| 1 | 
 1
sin 1
 sin 
 1 (1  0) →実装上: h(1 ) : 
理論上:
1
else
1
1  
0.8
0.7
0.6
0.5
0.4
0.3
0.2
01
0.1
0
-4
-3
-2
-1
0
1
2
- 制御対象固有の不連続点回避…真の特異姿勢を回避
厳密に 1   / 2 の時は,
の時は B ( x ) : [0 0 * *]  [0 0 * 0] となってA(x),
) B(x)不可制御
1   / 2   とすることで, 制御周期点では数値的な不可制御状態を回避
((同様に最下点でも微小数値を加え,, 
1
| 1 | 

100(| | )
不連続点を回避:
h(1 ) :  sin 1  0.01e
else

)
1  

1
2
- 速度情報のsmoothing…制御対象モデルの振動的変動回避
計測直接状態量 x と, フィルタで平滑化した速度情報を含む状態量 x
を用いて、実装上はSDLRを x  A( x)  x  B ( x)  u として扱う
III..SDRE
III
SDREの理論的背景
の理論的背景
1.. 非線形システムと各種制御手法
3
4
概説 ~非線形システムへの取り組み~
分類
・偏微分方程式 or 常微分方程式 or…
・smooth (無限回微分可能) or non-smooth …
・線形近似できるクラス
………….
・affine系
・線形近似できないクラス
・双線形システム(セミアクティブサス) …
・非ホロノミックシステム(車)
x  f ( x, u )
x  f ( x)  g ( x)u
x  f ( x)  x  u
非線形制御理論のトレンド
・微分幾何学に基づく厳密な線形化(1980前後)
・ダイナミクスベースド制御
ダイナミク
ド制御((1990前後)
・出力零化制御, Chained-form(1990前後) etc, etc ……
実用志向理論のトレンド
・平衡点近傍での線形近似
x  Ax  Bu
:
A : f x x  x0 , B : f u x  x0
u u0
u u0
:
・モデルベースド非線形補償
Mq  Nq  G  
q  u
~
~
~
加速度分解制御
  Mu  Nq  G
:
x  A( ) x  B ( )u
・ゲインスケジューリング
:
・Receeding Horizon(非線形モデル予測制御)
・SDRE法
◆ 本章の(全体講義における)目的
・ SDREの理論を理解する為の基礎復習/補足説明
- 伝達関数とシステムの安定性
数
定
- ブロック線図によるシステム表現
- 状態フィードバックによる安定化の原理
- 状態空間表現と状態フィードバック
- 最適レギュレータ(LQR)…対線形システム
- リカッチ方程式
(a) 2次形式, (b) リカッチ行列方程式,
(c) リカッチ微分方程式,
リカッチ微分方程式 (d) 数値解の求め方
- 非線形最適制御問題...対非線形システム
- ハミルトン-ヤコビ方程式
- SDRE法の理論的解釈
III..SDRE
III
SDREの理論的背景
の理論的背景
2. 線形最適制御問題
線形最適制御問題..a.
..a.状態フィードバックと安定化
状態フィ ド ックと安定化
◆
伝達関数とは
・動的なシステムの入出力関係の一表現法
・時間信号をラプラス領域に変換することで、時間領域での
畳み込み積分の記述を簡単な掛け算に記述するツール
動的システム=微分方程式で記述. その代数解を得る為には積分計算必要
・インパルス応答のラプラス変換
時間領域
u (t )
システム
H
 (t )
y (t )
畳み込み積分


y (t )  h(t ) u (t   ) d
0
h (t )
:インパルス応答
インパルス応答
L[h (t )]  H ( s )
伝達関数
・伝達関数=「入力信号のLaplace変換」に対する、「出力信号の
Laplace変換」の比
U (s
( )
H (s )
Y (s
( )
Y ( s)  H ( s) U ( s)
・利点: (畳み込み)積分形式が積の形式に変換できる
→複数
複数の動的システムの連結関係が四則演算的に可能
動
結関係
則演算
能
伝達関数
H ( s) 


0
h(t ) e  st dt
◆
ラプラス変換とは
時間関数
複素関数
F (s
( )
s : 複素数
f (t )
t : 時間
に移す変換

F ( s)   f (t ) e st dt : L[ f (t )]
0
例: f (t )  e tのラプラス変換は
F ( s )  L[e
t


]   e e dt   e (  s ) t dt
t  st
0
0

 e  ( s  )t 
1
 


 s   0 s  
◆
伝達関数と安定性
・伝達関数の更なる利点: 動的システムの安定判別が容易になる
・極と零点
H (s) 
( s  z1 )( s  z 2 )
( s  p1 )( s  p2 )( s  p3 )
- H ( s )を  にする s の値: s   p1 , p2 , p3 を極という
- H ( s )を 0 にする s の値: s   z1 , z 2 を零点という
・極=時間関数に戻した時の時間変数係数
先の例の場合:
f ( t )  e  t
L[ f (t )] 
極: -α
1
s 
・(漸近)安定
h(t )  0, t  
なるとき、システムは安定という。
・線形な伝達関数 H (s ) のシステムが安定となる必要十分条件は
の全ての極が左半平面にあることである。
ある とである。
H ((s )の全ての極が左半平面
 im1 ( s  zi )
H (s)  n
 i 1 ( s  pi ) とすると
定 あ 。
Re(( pi )  0, i  1,  , n ならシステムは安定である。
d ( s )   in1 ( s  pi ) をシステムの特性方程式という。
・LQR(最適制御設計)は、得られる解(制御則)が、必ずシステム
を安定化することを保障する設計理論
◆ 準備2…ブロック線図
準備2 ブロック線図
・フィードバック安定化の概念を理解していただく為に
・動的システムを信号の流れで表現する方法
動的シ テムを信号の流れで表現する方法
直列結合
R((s )
G1 ( s )
G1 ( s )
+
並列結合
H (s )
+
R(s )
G (s )
Y (s
( )
Y ((s )
Y ( s)  H ( s)  G ( s)  R( s)
Y ( s )  G1 ( s )  G2 ( s )   R( s )
フィードバック結合
+
R(s )
E (s )
-
G (s )
一般に左記のフィード
バック系において
Y (s )
R( s ) と Y ( s ) の間の
H ((s )
伝達関数は
伝達関数
Y( s ) 
なぜなら
G( s )
R( s )
1  G( s )H ( s )
Y ( s )  G ( s ) E ( s )  G ( s ) R ( s )  H ( s )  Y ( s ) 
1  G ( s) H ( s)   Y ( s)  G ( s) R( s)
Y( s ) 
G( s )
R( s )
1  G( s )H ( s )
◆ 準備3…フィードバック安定化
・開ループシステムの極を、シフトさせる手段
R(s)
n (s )
d (s )
閉ループ系
閉ル
プ系
Y (s)
R(s)
+
開ループ系
開ル
プ系
E(s) n (s )
-
Y (s)
d (s )
 (s )
 (s )
Y ( s) 
n( s )
R( s )
d ( s)
Y ( s) 
 (s)n(s)
R(s)
d (s) (s)  n(s) (s)
特性多項式が変わる→極が移動する→安定性を変更できる
・例:1バネ-マス-ダンパ系
2. 中間変数付きで運動方程式
u
my  u  F1  F2
y
m
k
F1  ky, F2  f y
f
3 ブロック線図化
3.
+
m
+
u
1.構成要素
-
-
F1
1
ms
1
s
f
F2
k
f
k
y
4. 内側のフィードバックループを変形
+
r1
e
-
1
ms
y
f
L[ y ]  Yd ( s ) とすると
Yd ( s ) 
1
1
R1 ( s)  ffYd ( s) 
E (s) 
ms
ms
f 
1

R1 ( s )
1 
Yd 
ms
 ms 
1
Yd  ms R1 ( s )
f
1
ms
y
5. 外側のフィードバックループを変形
+
U (s
( )
-
1
ms  f
1
s
Y (s
( )
k
1
1
(ms  f ) s
U (s)
Y (s) 
U (s)  2
k
ms  fs  k
1
(ms  f ) s
★ もしもフィ
もしもフィードバックゲイン
ド ックゲインkがないとすると
がないとすると…
1
1

ms 2  fs  k
ms ( s  f )
s=0の不安定極 →非安定
III..SDRE
III
SDREの理論的背景
の理論的背景
2. 線形最適制御問題
線形最適制御問題....b.最適レギュレタ問題とリカッチ方程式
◆
伝達関数と状態空間表現
・伝達関数(Transfer Function: TF)と状態空間表現(State Space:SS)は互換
・表現可能なクラスの大きさ: SS > TF
・SS→TFは一意に決定 H ( s )  C ( sI  A) 1 B  D
・TF→SSにはいくつか手段がある.
可制御正準形
bn 1s n 1    b0
H (s)  n
s  an 1s n 1    a0
1
 0
 0 

  

 x   u
x  

1  0 

  

a



a
0
n

1

 1 
y  b0   bn 1 x
◆
状態空間表現における安定性判別
・行列Aの固有値=伝達関数の極
x  Ax が安定であるための必要十分条件
A の全ての固有値の実部が負
・不安定な系の安定化には入力を使う
x  Ax
A  Bu
B
…状態フィードバック
u  Fx
u  Fx を代入した
x  ( A  BF ) x …閉ループ系
が安定になるように F を決める
( A  BF )の全ての固有
値の実部を負にする
例題…状態空間表現に基づく安定化の確認
0 1  0 
x  
 x  1u
1
2
  

u   f1
f 2 x
なるシステムを安定にする
を求めよ
f1 , f 2 を求めよ.
閉ループ系
閉
系
  0 1  0 
x   
 x  1 f1
1
2
  

特性方程式

f 2  x

 0
x  
 f1  1
| sI  A |
 f1  1  0
f1   1
s 2  ( f 2  2) s  ( f1  1)  0
が負の実部を零(解)に持つため
が負
実部を
解 持
の必要十分条件
◆
1 
x
f 2  2
 f2  2  0
f 2  2
安定化フィードバック則の
条件が求まる
最適制御とリカッチ行列方程式
・最適制御問題…2次形式評価関数の意味で最適な解を導く
(1) 線形時不変プラント:
x (t )  A  x (t )  B  u (t ) を漸近安定化(x (  )  0)し、
(2) 2次形式評価関数:

J  0 x(t )T Qx (t )  u (t )T Ru (t ) dt を最小にする解は、
(3) リカッチ行列方程式:
AT P  PA  Q  PBR 1 B T P  0 の正定対称解 P を用いて、
(4) 状態フィードバック則:
u (t )   R 1 B T P  x(t ) として与えられる。


重要な事 ・A,Bは定数行列
・漸近安定を仮定… t  で x (t )  0
・解が`状態フィードバック’という形式
◆ リカッチ方程式の本質
・ 最適問題の直感的理解…2
2次関数の最小化
y  a  bx  cx 2 , c  0, ( 4ac  b 2 )
2
2
2 
2

b
b
b
b



  a  c x 
 c  x 
 
 a
2 

2
c
4
c
2
c
4c





y
x*  b / 2c
a  b 2 / 4c
②最小値
を最小化
①最小値を与える
位置を算出
0
x*
x
・ リカッチ方程式
- J.F.Riccati(1676-1754)の名をとった1階非線形微分方程式
x (t )  a (t )  b(t ) x (t )  c(t ) x (t ) 2
- 微分方程式論の基本 → システム制御工学では行列版頻出
◆ 準備…2次形式と正定行列
・ 2次形式(quadratic form):
- 変数 x1 , x2 ,  xn についての2次の
次のべき関数
き関数   aij xi x j
i
j
- 対称行列 A を用いて x T Ax と書ける
例: 2 x12  2 x1 x2  4 x1 x3  x22  6 x32
 2 x12  x1 x2
2
 x1 x2 x3  1
 2
 2 x1 x3  x22  x2 x1  6 x32  2 x3 x1
 1 2  x1 
1 0  x2   xT Ax
一見行列のようだが、
0 6  x3 
スカラーであることに
・正定
正定 (positive definite):
xT Ax  0
準(半)正定 (positive semi-definite): x T Ax  0
負定 (negative definite):
xT Ax  0
注意
(  x  0)
(  x  0)
(  x  0)
続き…
・ 正定行列(positive definite matrix):
- 2次形式 x T Axが正定のとき、行列 Aを正定行列といい、
A  0と書く
- 2つの対称行列について、 A  B  0 のときは A  B と書く
・ 正定行列の性質(★重要)
(a) A  0 である為の必要十分条件は A の固有値が全て正
((b)) 任意
任意の行列
行列 C  R l nに対して、
対し 、C T C  0
(c) x T C T Cx  0 と Cx  0 は等価
◆ リカッチ行列方程式はどこから?
・LQR問題の証明の過程で導かれる
- 仮定: x (t )  A  x (t )  B  u (t ) は漸近安定可能:x (  )  0
- 準備:
準備:一般的な正定対称行列
般的な正定対称行列 P に対して、


0
d T
( x Px ) dt
dt
なる式の値を計算しておく
0

d
T
T
T

x (t ) Px (t ) dt  x (t ) Px (t ) 0  x (  ) Px (  )  x (0)T Px (0)
 0 dt
d
  x (0)T Px (0) …①
x  Ax  Bu
一方
方、積分内の微分を先に計算すると
積分内の微分を先に計算すると


 x (t )

0
T
Px (t )  x (t )T Px (t )  


0
( Ax  B)
①と②から次の恒等式が導ける
 ( Ax  B )

0
T
T
Px  x T P ( Ax  Bu )dt …②
Px  x T P ( Ax  Bu )dt  x T (0) Px (0)  0 …③
③
続き
- 証明
元の2次形式評価関数:

J   {x T Qx  u T Ru} dt
0
に③式の恒等式を足す.
0

J   {xT Qx  uT Ru  ( Ax  Bu)T Px  xT P( Ax  Bu)} dt  xT (0) Px(0)
0

J   {(uT  xT PBR1)R(u  R1BT Px
P )  xT ( AT P  PA Q  PBR1BT P)x} dt
d
0
 x T (0) Px (0)
J を最小化する為には, 2次項=0であればよいので,
u   R 1B T Px
…状態フィードバック則
状態フィ ドバック則
AT P  PA  Q  PBR 1B T P  0 …リカッチ方程式
…以上がリカッチ行列方程式の出生
以上がリカ チ行列方程式の出生
◆ リカッチ行列方程式の解き方:例題
・1入力2次-可制御正準系の場合
1  0
 0
システム x  
 x  1u に対し,

a

a
1
 
 0

評価関数を J    x1
0

q
x2  1
0
0   x1  2 
 u r dt とすると,
q2   x2 

解くべきRiccati方程式は AT P  PA  Q  PBR 1B T P  0 から,
1 
 0
 0  a0 
1  0
P

P
Q
P  0 1P  0


 a  a 
1  a 
r
1
1 
1

 0
 p pb 
これならば, P   a
 とおいて, 代数的に解けそうである.
p
p
c
 b
( PT  P  0 の条件を満たすものを選ぶ)
…2次以上で解くには?
◆ 解法1…リカッチ微分方程式の数値積分
・一般化最適制御問題
(1) 線形時変プラント:
x ( t )  A(t )  x (t )  B ( t )  u (t )
(2) 一般化2次形式評価関数:
J   x (t )T Q (t1 ) x (t )  u (t )T R (t1 )u (t )  dt  x (t1 )T Q (t1 ) x (t1 )
t1
t0
積分区間が有限
(
)
終端値が0でない(拘束条件となる)
(3) リカッチ行列微分方程式:
P (t)   A(t)T P(t )  P(t ) A(t )  Q(t)  P(t )B(t)R(t )1 B(t )T P(t )


(4) 最適制御入力:終端条件P (t1 )を出発して時間軸上逆方向に
積分して、P (t ) について解いて得られる( t0  t  t1の区間)
u * (t )   R (t ) 1 B (t )T P (t )  x(t )
通常の最適制御問題に対するリカッチ方程式
= リカッチ微分方程式の定常解
…続き
(リカッチ微分方程式の数値積分)
・リカッチ行列微分方程式:
P (t )   A(t )T P(t )  P(t ) A(t )  Q(t)  P(t )B(t)R(t )1 B(t )T P(t )


の定常解とは、 P (t )  0 の時の解 P(t)
・最適レギュレータ問題の時の条件
A(t ),
) B (t ),
) Q (t ),
) R (t )  A, B, Q , R  const に置き換え,
リカッチ微分方程式を使って, 逆時間方向に数値積分する.
P(  d )  P( )  P ( )  d P(0)  O
P ( )   AT P( )  P( ) A  Q  P( )BR1BT P( )


P(   d )  P( ) と十分収束した P()をリカッチ解 P( x) とする.
★注意:このときのτは単に数値計算の為変数で、実時間 t とは無関係
★dτは数値積分ための、数値積分の刻み幅
…続き2
(リカッチ微分方程式の数値積分)
・本手法の利点
- 実装が容易…低級水準の行列計算で実現
(行列の加減乗除とfor/while文程度)
R1も慣例的に対角行列に選べば逆行列計算不要
 r1

R 



 r1 1


  R 1 



rn 


・本手法の欠点
本手法 欠点
- 収束に必要な繰り返し計算回数が変動
- 数値的に不安定になる場合もある



rn 1 
・SDRE法では制御周期毎にリカッチ方程式の解が必要
本手法では制御周期時間内で計算が完了する保障がない
…続き3
(リカッチ微分方程式の数値積分)
・別の見方…ホモトピー法(変分法)的考え方
- 変数は少しづつ変化するので, 逐次積算で代替する.
・SDRE法の場合
- 状態 x の変化は連続的
- t t   の時の A(x(t))  A(x(t  )) 等の変動量は少ないはず
- P(x(t))  P(x(t  )) の変動量も少ないはず.
・1サンプル時間前の P(x(t  )) を、リカッチ微分方程式の数値
積分開始点として使用すれば 収束までの計算回数を減らせる.
積分開始点として使用すれば、収束までの計算回数を減らせる
t
実時間
t
P( x(t  ))  P( ) | 0
繰返し計算
P(  d )  P( )  P ( )  d
P( x(t ))  P( )
P( x(t ))  P( ) | 0
◆ 補足…固有値/固有ベクトル
・リカッチ解の代数的解法=固有値法(有本-ポッター法)の理解の為
・固有値:特性方程式の根
- | sI  A | s n  an s n 1  an 1s n  2    a2 s  a1  0
- 正方行列
方行 Aについてのみ
・固有ベクトル
- A v i  i v i を満たす n  1 の1次元ベクトル v i
(i  1,  n)
- 重複固有値がない場合, 各固有値i について, ただ1つの
固有ベクトル v i が存在する
- 固有ベクトルは互いに線形独立
a1v1a2 v2 an vn  0 が成立するのは a i  0 ( i ) の時のみ
…続き
1
0
・例 A  

  2  2
 s 1 
- 特性方程式:sI  A  
 s(s  2)  2  s 2  2s  2  0

2 s  2
- 固有値: i  1  j , 1  j
- 1  1  j の固有ベクトルは (1I  A)v1  0 を満たす筈なので
とお
展開すると,
v1  [   ]T とおいて展開すると
(1 j)    0
1  j 1    0
等価
(1I  A)v1  

   0
2


(
1

j
)


0
2
1

j

   
  (1 j)
T
- 等価なので,   1 と選べば   1 j となり, v1  [ 1 1  j ]
- v2 についても同様に v2  [ 1 1  j ]T と求められる
◆ リカッチ方程式の解法2…有本-ポッターの方法
ハミルトン行列:
 A  BR 1 B T 
2 n 2 n
M 
R
T
A 
 Q
とおくと, 次式が成り立つ
| sI 2 n  M | (1) n  | sI  A |  |  sI  A |
ここで A  A  BF , F  R 1 B T P
Mの固有値は、虚軸を対称にして配置→  a  bj , a  bj
① 最適極(実部が負)を選ぶ→ 1 ,  , n (Re(i )  0, i  1~n)
② 各 i に対応する固有ベクトルwi  C 2 nを求める→ Mwi  i wi
③ 固有ベクトルを2つに分ける
u 
wi   i  (ui , vi  C n )
 vi 
④ リカッチ解Pが次式で求まる P  [v1  vn ][u1  u n ]1
◆ ハミルトン行列とリカッチ方程式
P(t )  V (t )  U (t ) 1 とおく と(以下, Pは時間関数と一般化する)
P代入
P   VU 1  VU 1UU 1 …(1) が成立する
(なぜならPU  V  P U  PU  V  P   VU 1  PUU 1 )
両辺微分
U逆行列
ここで V  ATV  QU
V (0)  P0 …(2)
U  BR1BTV  AU U (0)  I n …(3) とすると
(1)式は P   AT VU 1  Q  VU 1 BR 1 B T VU 1  VU 1 A
 AT P  Q  PBR 1 B T P  PA…リカッチ微分方程式((4))
 V (t ) 
そこで Z (t )  
ハミルトン行列
 とおくと

U
(
(t
t
)


 V   AT
 Q  V 

Z 
 MT Z





T
1

 A  U 
 U   BR B
…続き
- ハミルトン行列Mには先に示したように, その固有値は実部が
が
正と負のペアになって現れる.
- リカッチ微分方程式
リカ チ微分方程式(4)の定常解は, Z(t)
( )の収束解に相当する.
- よってZ(t) の安定モードである, 安定固有値 1 ,  , n (Re(i )  0)
とその固有ベクトル wi  C 2 n が求まれば
 v   v   v    V 
[ w1 w2  wn ]    1   2    n    

u

u
u


U

n


1
2








P(t )  V (t )  U (t ) 1 とおいたので、
 [v1  vn ][u1  u n ]1  P
◆ 有本-ポッター法改善策
・有本-ポッター法は定常解の解析には有効な手法
・実際の計算では複素数演算が必要
・次元が大きい時や, 固有値の大小差が大きい時には精度劣化
・解決法の一例: 実Shur分解による手法
- 固有ベクトルの代わりにShurベクトルを用いる
- Aの複素固有値部分を
複
有値
2x2の実数対角部分に変換し
角
変換 ,
Hの複素固有値に対
複素数演算を回避する.
*
*  応した対角部分のみ
1 *
T
2 n 2 n
に2x2の対角ブロック
M  W MW  U  R


2a * * 
それ以外の対角要
U 
- Uに対して固有値
素に
実数固有

* 2b *  素にHの実数固有
値をもつ


分解法を適用する
n 
0
◆ SDRE制御器の実装アルゴリズム
1. x  x(t) を計測 (場合によっては推定)
) B( x) を算出
2. SDLRの係数行列
係数行列 A(x),
3. 状態依存重み Q(x), R(x) を計算
4. SDREを解き, 正定解 P(x) を算出
AT ( x)P( x)  P( x) A( x)  P(x)B(x)R1(x)B( x)P( x)  Q( x)  0
5. 制御入力(時刻tでの)を算出
u' (t )   R  x  B  x  P  x   x
6. 時刻 t=t~t+Δtの間 u’(t)をプラントへ出力
7 1へ戻る.
7.
1
T
III..SDRE
III
SDREの理論的背景
の理論的背景
3 非線形最適制御問題
非線形最適制御問題…
…a.ハミルトン
ミルトン-ヤコビ方程式
ヤ ビ方程式(HJE)
◆ 非線形最適制御問題
x  f ( x, u, t)
非線形システム
非線形システム:
評価関数:
f ,  ,  :連続微分可能
T
V ( x(t ), u(), t ) :  ( x( ), u( ), ) d   ( x(T ))
t
非線形最適制御問題
V を最小化する最適入力 u* (t ), t  [t , T ] を見つける問題
=次式の最適評価値を求める問題
T
*

V ( x (t ), t )  min   ( x ( ), u( ), ) d   ( x (T ))

u(t ) 
 t
V * ( x (t ), t ) は次のHamilton-Jacobi equation(HJE)の解
T


 V * 
V *
 f ( x (t ),
  min  ( x (t ),
) u(t ),
) t )  
) u(t ),
) t )
u(t )

t
x




(1)
→解析的には解けない…
◆ 非線形最適制御問題の簡単化
簡単化その1… 非線形システムと 評価関数のクラスを限定
システム
x  f ( x, u, t )
被積分項
 ( x ( ), u( ), )
x  f ( x)  g ( x)  u , f (0)  0
x T Q ( x ) x  u T R ( x )u

V * ( x )  min  ( x T Q ( x ) x  uT R ( x )u ) dt
新評価 数
新評価関数
u(t )
0
Eq.(1)右辺の鍵括弧内
右辺 鍵括弧内
T
 V * 
   f ( x )  g ( x ) u 
 x Q ( x ) x  u R ( x ) u  

x
T
T


 V * 
1  V * 
V *
T
 g ( x ) R 1 g T ( x )
 f ( x )  
 x Q ( x ) x  
4  x 
x
 x 
(a)
* T


1 1
1 1
V 
V *
T
T
 R ( x )  u  R ( x ) g ( x )
  u  R ( x ) g ( x )
2
2
x

x



(b)
T
T



…続き
(a)=0 ⇒
T
簡単版 HJE
(b)=0 ⇒
最適入力
T
 V * 
 V *  T
1  V * 
1
T

 f ( x)  
 g( x) R ( x) g ( x)
  x Q( x) x  0
4  x 
 x 
 x 
(2)
1 1
V *
T
u   R ( x) g ( x)
2
x
(3)
・まだまだ解析解を得るのは困難…
・常套
常套(代替)手段…数値的に近似
数値的 近似
Taylor級数展開法, Galerkin近似, 非線形行列不等式の活用, GA, …
従来法: 状態xの変動幅を限定, 非線形要素の変動幅を小さ
くして近似解析解や数値近似解を求める.
SDRE法: 逆に変化する状態xに依存させ, HJEの状態依存数
値解を算出. これらを時間的に紡いで元のHJEの解とする.
III..SDRE
III
SDREの理論的背景
の理論的背景
理論的背景
3 非線形最適制御問題
非線形最適制御問題…
…b.
b HJEとSDRE
◆ 非線形最適制御問題からSDRE法へ
・ SDRE法…近年提案されたHJEの近似解を得る為の代替手法
ex.) J.R., Cloutier, ACC. (1997)
・簡単化その2… 状態依存線形表現(SDLR)で表せるシステムに限定
(SDLR=State Dependent Linear Representation)
SDLR
x  A( x ) x  B ( x )u
YH
Y.Huang
and
dW
W.Lu,
L CDC(1996)
最適指標 V *は, もし存在するとすれば, 対称行列 P ( x ) を用いて,
x T Px のようにかける. Anderson & J.B.Moore, Optimal Control(1990)
V * ( x (t ),
) t )  x T Px
P ,
f ( x)  A( x)  x, g ( x)  B ( x)
…続き
E (2)
Eq.(2)→
SDRE
AT ( x)P( x)  P( x) A( x)  P( x)B( x)R1( x)B( x)P( x)  Q( x)  0
(4)

T p
T p j
satisfying x
i
x j
xi
( x ),
x, i, j  1,, n
P( x )   p1 ( x ),, pn ( x )
T
L. Wei Min and J. Doyle, IEEE Trans. of AC(1995)
Eq.(3)→
q( )
制御則
( x)  x
u   R ( x ) 1 B ( x ) T P ( x )  x
・状態 xを固定(frozen)したら, Eq.(4)は通常のリカッチ方程式
→数値解として解ける!
(5)
IV..SDRE
IV
S
SDRE法と類似手法との比較
法と類似手法との比較
1. 非線形モデル予測制御
2. ゲインスケジューリング
3 拡張カルマンフィルタ
3.
4. ハイブリッドシステム
IV--1. 非線形モデル予測制御
IV
非線形モデル予測制御((Receeding Horizon)
・現時刻での最適入力の初回のみを使用するという点で類似
・現時刻からある有限未来時刻までの最適化問題を解く.
・何らかの最適ツールを用いて解を探索→計算量が多い.
・実時間実装には本質的に不向き
実
実装
質
(制御周期の遅いプラント
御周期
)
x((t )
t
T
t
Horizon(予測区間)
◆ LQRとモデル予測制御
・理論的前提条件の違い
理論的前提条件の違
LQR
評価区間 t = 0 ~ ∞
プラント
不変
最終状態 x → 0 (t →∞)
解の形式 時不変/FB則/式
RH
x(t )
t
実時間最適化(モデル予測制御)
t = 0 ~ 有限
時変
仮定不要
逐次変更/入力量
時刻 t’ においてT
だけ未来までの制
御を最適化
t
T
実際のプラントへの
入力はuの初回のみ
u  {u (t ),  u (T )}
IV--2. ゲインスケジューリング
IV
・制御ゲインをオンラインで変更するという点で類似
制御ゲインをオンラインで変更するという点で類似
・起源: プラントの特性変動を補償する方式として考案
・特性変動の原因
特性変動の原因(環境の変化)が測定できることが前提
→その値によってコントローラのパラメータを調整
スケジューラ
y (t )
+
r (t )
 (t )
環境条件
-
制御器
プラント
・ 長所: 特性変動
特性変動に対する即応性はよい
対する即応性はよ .
・ 短所: 原理的には開ループ方式なので, 環境条件と制御器の
パラメータの対応関係が不正確だと
ラメ タの対応関係が不正確だと, 性能/動作が保証困難
◆ LMI(線形行列不等式)によるGS設計法
・スケジューリングパラメータ   [1 ,  2 ,,  L ]T  R L を用いて,
プラントの変動を表現(上下限を仮定:   [ ,  ])
・LPV (Linear Parameter Varying)プラントモデル:
x  A( ) x  B1 ( ) w  B2 ( )u
G ( ) z  C1 ( ) x  D11( ) w  D12 ( )u
y  C2 ( ) x  D21( ) w
A( )  A0  i 1  i ( )  Ai

z (t ) : 制御量
y (t ) : 観測量
G(α)
w (t ) : 外生信号
u (t ) : 制御入力
…続き
・制御器構造

K ( ) xC  AC ( ) xC  BC ( ) y
u  CC ( ) xC  DC ( ) y
z (t )
w (t )
G(α)
u ((t )
y (t )
K(α)
・全変動に対し内部安定, 所望性能を満たす個々の制御器を設計
- L個のパラメ
個のパラメータを
タを, 各々N L分割した格子点ごとに制御器設計
  [1 ,  2 ,,  L ]Tの 1 について→ {1,1 , 1, 2 , , 1, N }
同様に他の i についても分割
→ { L ,1 ,  L , 2 , ,  L , N }
- 全ての i の組み合わせ N1  N 2    N L通りで実施
1
L
・多項式表現で表せる程度のクラスが限度(計算量大なので)
…続き
・個々の制御器の設計は H 制御問題に相当
・2次形式の行列不等式: AX  XAT  XCT CX  BBT  0
の正定対称解 X を見つける問題に帰着
Schur Complement(
p
(行列変換公式の一種)
 AX  XAT  BBT

CX

XCT 
  0
I 
線形行列不等式に変形できる
(Linear Matrix Inequarity:LMI)
・LMIはXに関して凸問題(最適化問題の特徴の
最適化問題の特徴の一つ)
- 凸問題: 局所最適解が見つかれば, それは大域解
→解X行列が求まり, 制御器が設計できる
◆ ゲインスケジューリングのまとめ
・変化するプラントを、格子点で線形システムに近似し、その代表
点で設計した個々の制御器を切り替えて適用
・状態xとパラメータ変動αとを, 理論上は分離して扱っている為,
状態に依存した非線形性が強いプラントには不向き
GS,RH
, とも理論的には堅実&正攻法
攻法. 実用解を得るのは困難
⇔SDREは理論的には発展途上だが, 実用解は得意
IV--3. 拡張カルマンフィルタ
IV
・非線形状態推定器(オブザーバ
ブザ
)の一種
・プラント非線形モデルに基づき、状態に依存した線形化モデルを

もとに計算する点がSDREと類似
f ( x)
f ( x, u )
x
- 非線形プラントモデル: x  f ( x , u )
x  xˆ
u 0
- 現在の状態量(推定値): xˆ ( t )
x  xˆ
- 現時刻の線形化モデル: x  A( xˆ )  x  B ( xˆ )  u


A( xˆ ) 
f ( x, u )
B ( xˆ ) 
f ( x, u )
x  xˆ
x  xˆ
x

u
u 0
u 0
x
- 拡張カルマンフィルタ(Extended Kalman Filter: EKF)
xˆ  A( xˆ )  xˆ  B ( xˆ )  u  K ( xˆ )  ( y  C xˆ )
推定状態量
(推定)状態依存係数 オブサーバゲイン
◆ (線形)カルマンフィルタ
・推定状態量 x̂ と真の状態量 x との二乗誤差平均を最小化
するように状態を推定 E [|| x ( t )  xˆ ( t ) ||2 ]  min
・観測モデル(状態を推定したいプラントの数式モデル)を準備
観測 デ (状態を推定 た プ
数式 デ )を準備
x  Ax  Bu   v システム雑音
y  Cx  w 観測雑音
・外乱(ノイズ)の統計的性質を仮定 共分散行列
1 t  0

(t )  
E[v(t )]  0, cov[v(t ), v( )]  V (t  )
0 else
E[w(t )]  0, cov[[w(t )), w( )]  W (t  )
cov[a, b]  E[(a  a )  (b  b )]
・カルマンフィルタは動的システムとして得られる
xˆ  A xˆ  Bu
B  K  ( y  C xˆ )
K   C T W 1 : フィルタゲイン
A    A T   V  T   C T W 1C   0
システム雑音 v(t )
観測雑音 w(t )

+

+
+
A
プラント
カルマンフィルタ
B
+
u (t )
C
y ((t )
+
+
B
x(t )
+
u (t )
x  Ax  Bu   v
y  Cx  w
K

-
xˆ (t )
C
xˆ (t )
状態推定値
A
xˆ  A xˆ  Bu  K  ( y  C xˆ )
◆…続き
・フィルタゲインKを与える  は推定誤差の共分散行列
  cov[ x (t )  xˆ (t ),
) x (t )  xˆ (t )]
・そのリカッチ方程式はLQRと双対(本質的には等価)
K l
Kalman-filter:
fil
A    A T   V  T   C T W 1C   0
{, A, CT ,W , VT , K}  {P, AT , B, R, Q,F T }
LQR:
AT P  PA  Q  PBR 1 B T P   0
◆…拡張カルマンフィルタ(EKF)とSDRE法
・EKFは,観測モデルが非線形の場合に、それらを線形化近似し、
線
線形カルマンフィルタの構造で、状態を推定する手法
ィ
構
、状 を推定す
法
・EKFの線形化モデルは推定状態近傍の局所近似モデル
⇔SDREでは、SDLRモデルを使用している点が異なる
・とはいえ,EKFとSDREは, 解くべきリカッチ方程式も双対であり,
計算アルゴリズムもほぼ同等であり, 双対的手法といえる.
・EKFは基本的に“状態推定”であり, プラントの動作に直接的に
は関与しない. →気軽に使えるため, 多用されてきた(と推測)
・SDRE法のように制御器は, 何よりも安定性の補償(システムの
安全保証)が第一なため, 最近まで実使用されなかった(と推測)
IV--4. ハイブリッドシステム
IV
・基本的にはSDRE法は連続時間システムが対象
続
が 象.
しかし時々刻々最適解を更新することで不連続システムにも
対応できる可能性有り… 倒立振子の振り上げがその例
・理論的に不連続時間システムを扱うには?
→ハイブリッドシステム
イブリッドシステム(不連続切替を扱った
不連続切替を扱った一般的枠組み
般的枠組み)
・HS = 連続時間システム + 離散時間プロセス(論理・意思決定)
前者...差動/微分方程式, 連続/離散時間システム
後者...有限オートマトン、離散事象システム
例: 車の変速器, ディスクドライブ, ステッピンモータ, 炉
•理論的に証明することに注力して研究
論的 証 する と 注力
究
論点:切替時の安定性補償、切替制御器の存在性
(‘倒立振子の振り上げ/安定化制御’がこれにあたる)
◆
HSの難しさ
・個々のシステムが安定だからといって、それらをつないだ
個 のシ テ が安定だからと
て、それらを な だHSは安定であ
るとは限らない。
x (t )  Ap ( t )  x (t )  R 2 , where p (t )  {1,2} :切替え変数
HS
  1  100
A1  
,

10
1


  1 10 
1, 2  1  j 1000
A2  

  100  1
:双方とも安定
切替則
1, if p(t  )  2 and x2 (t)  1/ k  x1(t)
p(t )  
if p(t  )  1 and x2 (t)  kx1(t)
2,
(1.2)

・ k  0.2 の場合はどんな初期値( x  0)
に対しても発散する。
◆
Fig.1 Unstable state traj. of
switched stable systems.
systems
HSの面白さ
・個々のシステムが不安定だからといって、それらをつないだ
個 のシ テ が不安定だからと
て、それらを な だHSは 不安定
であるとも限らない。
HS
x (t )  Ap ( t )  x (t )  R 2 , where p (t )  {1,2}
0 10
A1  
,
0
0


2   ( A1 )  0,  ( A2 )  0.5  j 3
1.5
A2  
 :双方とも不安定
  2  0.5
切替則
1, if p(t )  2 and x2 (t)  0.25x1(t)
p(t)  

2, if p(t )  1 and x2 (t)  0.5x1(t)
(1 3b)
(1.3b)
・HSは安定。
Fig.3 State traj. Resultinf
from switching bet.
bet A1 and A2
◆
HSのモデリング
状態方程式
x (t )  f  x(t ),
) p (t ),
) u (t ) 
(2 1a)
(2.1a)
出力方程式
y (t )  g  x(t ), p (t ), u (t ) 
(2.1b)
モード決定式
ド決定式
p (t  )    x(t ), p (t ), u (t ),  (t ) 
(2.1c)
連続時間状態量
x (t )  R n
連続時間入力 or 外因性信号(参照信号,外乱含む)
u (t )
p (t )  {1,  , M } 離散状態量。ベクトル場 f () を切替。入力となり得る。
離散事象(制御可能)
 (t )
x
x
B

 (  ,  ,  ,  ) 不連続遷移関数
1
(記述モデルand/or
記述 デ
離散事象システム含む)
u
B2
A1
   ( x)
A2
BN
x (t )  Ai  x (t )  Bi  u (t )
◆
AN
   ( x)
HSの例:オートマティック変速器
・斜度 の斜面を速度 v で走行させるための自動変速装置付車の特性
G p (t )
k
v 車速
v   v 2 sign (v)  g sin  
T
T エンジントルク(入力)
m
m
(2.5a)
G p (t ) 変速ギア比
  G p (t )  v
(2.5b)
k, m 定数
G p ( t )  {G1 , G2 , G3 , G4 }
 エンジン角速度
G1  G2  G3  G4

地面傾斜角度
離散状態推移関数
if p (t )  i  4 and v  1 / Gi   high シフトアップ
i  1

p (t )  
if p (t )  i  1  2 and v  1 / Gi 1  loqq シフトダウン (2.6)
 i
Fig.5 Automaton illustrating discrete state transistors
PIクルーズ制御 …道路の傾斜によらず速度を一定に保つ制御
k 2
TP  K p ( t ) (vref  v) :比例制御分
T  TP  TI 
v sign (v )
G p (t )
TI :積分制御分
((2.7))
vref :参照速度
閉ループ系
G p (t )
K p (t ) (vref  v)  TI   g sin 
v 
m
K p (t )
TI 
(vref  v ) TR :積分時定数
TR
(2.8a)
(2.8b)
制御条件…快適な乗り心地の為に=シフト時に車速が滑らかに繋がる様に。
(2.9a)
Gi K i  Gi 1 K i 1
G p ( t  )TI (t k )  G p ( t  )TI (t k ) (2.9b)
k
k
…速度偏差に対する補正
…累積トルク積分量に対する補正
・(2.8b)の積分器は vref が変わったとき0にリセット。
◆
HSの一般的安定性
・システムを切り替えても、
シ
ムを切り替え も Lyapunov((リアプノフ)関数が存在すればよい。
プ
)関数が存在すればよ
リアプノフ関数Vとは
・状態空間上に定義された正値スカラー関数。V(x)>0
・状態が安定するに従って、値が0に収束する関数.
機械系,電気回路のような物理システムのエネルギーの概念を,
微分方程式や差分方程式で表される動的システム(DS)に拡張したもの。
・利点: システムの安定性を, そのシステムの解を求める事なく,
間接的に調べられる。
・HSに対しては、複数のLyapunov-like関数
(
(LLF)
)を繋げたLFを用いることで安定性の
を用いる とで安定性の
証明を狙う
◆
HSに関する疑問点
a)全ての切替入力に対して安定に状態遷移を行う
全
切替入力に対し 安定に状態遷移を行うHSのクラスは?
ク
は
b)HS全体を安定に保つ切替入力はどんなものか?
(大雑把な)答え
a)システムの切替に関係なく,共通のLyapunov関数が存在すること。
関数が存在すること
b) 個々の(漸近)安定なシステムよりも入力の切替が十分遅いこと。
(直感的回答)
報告事例
・2つの安定な線形システムに対する、共通な(2次)Lyapunov関数が存在
するための必要十分条件は導出済み
→ 逆に言うと3つ以上のシステムに対しては完全には未解決
◆
安定解析の例…オートクルージングモデル
・LMI(Linear
LMI(Li
Matrix
M t i Inequality
I
lit :線形行列不等式)で安定性を確認する。
で安定性を確認する
・リアプノフ関数候補を
す
V ( x )  x T Px とすると
V ( x )   x T Qx, Q  0 を満たすような行列Qが存在すればよい
・システムは x  Axなので、上の条件は
AiT P  PAi  Q とLMI条件に変換できる
条件:平地走行(斜度α=0). v : vref  v(速度差), TI  TI(トルク指令積分値)
→ 閉ループ  v    G p ( t ) K p ( t ) / m  G p ( t ) / m   v 

 T    K / T

0
p(t )
R
 I 
  TI 
この Ai に対して, P が存在すれば安定
(A) TI がjumpしない場合
…(シフト切替前後で(2.9b)を満たす為に不連続にリセットする必要の無いとき)
・領域を分ける必要なく、単一のLLFですべての離散状態をカバーできる。
・LMI問題の条件を満たす解として
255.589 72.262 
P

 72.262 40.822 
が得られる。
・正定解Pが存在するので、状態跳躍ない場合はHSは大域的指数安定。
(B) TI がjumpする場合
・ロー→ハイへのシフトアップ:条件 TI  K p(t ) (v  vref ) を満たす軌道が、
ギア切替線 v  1/ Gp(t )  low を左から右へ横切るとき。
・ハイ→ローへのシフトダウン:条件
ハイ ロ へのシフトダウン:条件 TI  K p(t ) (v  vref ) を満たす軌道が、
を満たす軌道が
ギア切替線 v  1/ Gp(t )  high を右から左へ横切るとき。
・新しい参照速度
新し 参照速度 vref が与えられたときは
与えられた きはTI (0) は常
は常に0へリセット
リ ッ
1速から制御開始したとすると、切替条件(2.9b)から
1
 v 
  T   0
 I  i 1 
0
 v 
 v 
 : Fi  




G i 1 / G i    T I  i
  TI  i
LMI条件((LLFが切替後減少すると
が切替後減少するという条件
う条件)は F T PF  P である。
i
i
2
これを展開すると p2,2 (Gi1 / Gi )  p2,2 となる( p2,2 はPの(2,2)要素)。
なので、
、xT Px で表されるエネルギーは常に切替後も減少。
表さ る
ギ
常 切替後も減少。
Gi1 / Gi  1 な
→よって指数安定
このように理論的にHSの安定性を示せる。
の安定性を示せる
◆
ハイブリッドシステム論のまとめ
・より広い事象を表現するには適した手段
り広 事象を表現する は適 た手段
大半のシステム= 連続時間部分とそれらの不連続切替からなる
・HSの安定性に関する一般的必要十分条件は未確立
-共通リアプノフ関数, LMIによる手法が有力候補
-それ以外にも様々な定理/安定化条件が提案されている.
-それだけに本理論は奥深い
・現存するシステム設計論を大きく変える可能性も秘めている.
-例: 不安定なサブシステムの切替で全域安定化が達成できる.
・HS論はSDRE法よりも理論志向. 参考にする点は多い.
V SDRE法の研究トピック
V.
SDRE法の研究トピック
1.理論に関して
1.
理論に関して
大域安定性
・一般的には保証されない
・現時刻で固定された評価関数に対する最適化問題の
解であることのみ保証
・大域漸近安定性の議論&試みも有り.
① Y.Huang and A.Jadbabaie, IFAC(1999)
SDLRの非
の非一意性
意性
・変な線形化モデル(厳密な線形ではない)
・大域安定性を補償できない理由の つ
・大域安定性を補償できない理由の一つ
・A(x)x のとり方は任意. 自由度が残る
→逆にこの自由度を活用する研究も有り
② C.P. Mracek, J.R.Cloutier, CCA(1996)
・SDRE法の概要
・J.R. Cloutier, ``State-Dependent Riccati Equation Techniques: An Overview,'' in
Proc. of the American Control Conference, Albuquerque, New Mexico, pp.932-936,
1997.
・非線形H∞ 制御…ゲインスケジューリング+SDRE手法
・Yun
Yun Huang and Ali Jadbabaie,
Jadbabaie ``Nonlinear
Nonlinear H_inf
H inf Control: An Enhance Quasi
Quasi-LPV
LPV
Approach,‘’ in Proc. of the 14th World Congress of IFAC, Beijing, P.R., China, pp.8590, 1999….①
・J.R.
J.R. Cloutier, et.al, ``Nonlinear
Nonlinear regulation and nonlinear H
H-infinity
infinity control via the
state-dependent Ricatti equation technique:Part 1, Theory,'' in Proc. of the Int. Conf.
on Nonlinear Problems in Aviation and Aerospace, Daytona Beach, FL. 1996.
・非線形オブザーバ…拡張カルマンフィルタのSDRE版
・Curtis P. Mracek, James R. Cloutier, and Christopher A.D’Souza, ``A New Technique
f N
for
Nonlinear
li
Estimation,‘’
E i i ‘’ in
i Proc.
P
off the
h 1996 IEEE Int.
I Conf.
C f on Control
C
l
Applications , Dearborn, MI, Sept. 15-18, pp.338-343, 1996…②
・非線形最適制御問題との関連性
.
・ Y. Huang, et.al., ``Nonlinear Optimal Control: Alternatives to Hamilton-Jacobi Equation,''
in Proc. of the 35th CDC, Kobe, Japan, pp.3942--3947, 1996.
.
・ W.M.Lukes,
k ``O
``Optimal
i l regulation
l i off nonlinear
li
dynamical
d
i l systems,'''' SIAM
S A J. on Control
C
l
and Optimization, vol.7, No.1, pp.75--100, 1969.
・SDRE法での理論的安定性補償
.
・ Jeff S.Shamma, and James R.Cloutier, “Existence of SDRE stabilizing feedback,” IEEE
pp
, 2003.
Transaction on Automatic Control,, vol.48,, no.3,, March,, pp.513-517,
・ J.Willard Curtis and Randal W.Beard, “Ensuring stability of state-dependent Riccati
equation controllers via satisficing, in Proc. of the 41st IEE Conf. on Decision and Control,
LasVegas, Negada USA, pp.2645-2650.
・RH(モデル予測制御)とSDRE法との組み合わせ
.
・ Mario Sznnaier and Rodolfo Suarez, “Suboptimal control of constrained nonlinear
sysmtes via receding horizon state dependent Riccati equations,” in Proc. of the 40th IEEE
Conf. on Decision and Control, Orland, Florida USA, December, pp.3832-3836, 2001.
V SDRE法の研究トピック
V.
SDRE法の研究トピック
2.実装に関して
2.
実装に関して…
…適用事例
・飛翔体の制御
飛翔体の制御
ミサイルの誘導制御 … N.F. Palumbo, CCA(1999)
音速飛行機の航行制御
音速飛行機
航行制御 …J.R.Cloutier, CCA(1999)
高度によって変化する動特性をSDRE法で補償
Jet plane(Sonic cruiser)
・N.F. Palumbo and T.D. Jackson, ``Integrated missile guidance and control: A state
dependent Riccati differential equation approach,'' in Proc. of the IEEE CCA, Kohala
Coast-Island of Hawai
Hawai'ii, pp.243-248,
pp 243-248 1999
・James R.Cloutier and Peter H.Zipfel,”Hypersonic Guidence via the StateDependent Riccati Equation Control Method,” in Proc. of the 1999 IEEE Int. Conf.
on Control Applications, Hawaii USA, pp.219-2224, 1999.
・Praveen Shankar, Rama K.Yedavalli, and David B.Doman, “Dynamic inversion via
state-dependent Riccati equation approach: application to flight vehicles,” American
Insititute of Aeronautics and Astronautics Paper,
.
・ヘリコプターの姿勢安定制御
ハイブリッド的使用法
高度に依存する非線形特性 →SDRE
不確定性補償 →ニューラルネット
Hovering Control
・Eric A. Wan and Alexander A.Bogdanov, ``Model predictive neural control with
applicatoins to a 6 DoF Helicopter model ,' in Proc. of the American Control
Conference, Arlington, VA June 25-27, pp.488-493, 2001.
・Alexander A. Bogdanov and Eric A. Wan, `Model predictive neural control of a
high-fidelity helicopter model, ’ American Institute of Aeronautics and Astronautics
paper, 2001-4164.
.
・マニピュレータへの応用
・ M.Xin,
M Xi et.al,
l ``Robust
``R b state dependent
d
d Riccati
Ri
i equation
i based
b d robot
b manipulator
i l
control,'' in Proc. of the 2001 IEEE CCA, pp.369--374, 2001.
VI. SDRE法の拡張
SDRE法の拡張
法 拡張
1.状態依存評価関数の活用
1.
状態依存評価関数の活用
- 機構内力抑制制御…機械・制御同時設計の一例
- Projection Method…物理システムのモデリング法
◆
状態依存型評価関数と同時最適化問題
・SDRE法のポイント: 評価したい物理量や制御指標が状態依存
表
形で表現できれば制御系設計可能
SDRE法の評価関数: V ( x ) 


0
( x T Q ( x ) x  uT R ( x )u ) dt
メカトロニクス設計
・機械系
機械系(ハード
ド)設計製作→制御系(ソフト
ト)設計の傾向あり
・機構設計条件: 機構機能, モータ容量, 構造強度, 寿命/安全率
…運動設計に関係
運動設計 関係=制御系に依存
制御系 依存
・制御系と機械系の同時最適化問題
・機械系の指標をSDRE法の評価関数に考慮できたら…
→ 機械系の特性を加味した制御が実現
◆ 内力とは
・内力(internal force)…機械構造における各要素間の運動を拘束
するための内在力. 重力, コリオリ力, 関節結合力などの合力.
V
V
mg
mg
H
H
+
H
V
・内力の低減
減→機構部への力負荷低減
構部
負荷 減→穏やかな制御
◆ 内力の計算
・倒立振子のモデリング(ニュートン運動学に基づく場合)
振子長:2l
振子長
Vl sin 
u

V

( x, y )
mg

H
I  Vl sin   Hl cos   c
d2
d2
m 2 y  m 2 (l cos )  V  mg
dt
dt
2
d
d2
m 2 x  m 2 (  l sin  )  H
dt
dt
M  u  H  f
式変形を経て…
振子・回転運動
振子・並進運動
(垂直方向)
振子・並進運動
(水平方向)
台車 並進運動
台車・並進運動
H  m  ml ( sin  2  cos)
V  ml ( cos 2  sin )  mg
・ニュートン運動学に基づく運動方程式導出時には、中間
ュ トン運動学に基づく運動方程式導出時には、中間
変数である内力を手計算により式変形で導いていた。
法では 内力の一般計算式が得られる
般計算式が得られる。
・PJ法では、内力の
◆ 機械系の拘束システム表現…Projection法
…Projection Method. W.Blajer, J.Appl.Mech(1992)
非拘束運動:各リンクが別々に運動できる状態
内力(運動拘束力)
M a qa  ha  CaT  a
拘束運動: 各リンクが連結した運動状態
J1 ,C1
1
m1 g ( x, y , z )
J 0 ,C0
0
z
qa :   0 x y z 1
J0


Ma  




T
拡張状態量
o
u  C00 




 0 

, ha   0 



  m1 g 

  C11 
J 1 


m1
m1
m1
y
x
J 0 , J1
m1
C0 ,C1
u
:inertia
:mass
:friction
:torque
◆ (機構的)拘束条件の導出…回転型倒立振子を例に
重心位置の幾何学関係=運動拘束条件
x  L0 sin  0  l1 sin 1 cos 0
y  L0 cos 0  l1 sin 1 sin  0
z  l1 cos1
0
・ 方程式=等式=拘束条件
・ 微分しても等式は成立
微分し も等式は成立
0
z  l1 sin 11
1
1
0
0
1
1
( x, y , z )
L0
z
x  (L0 cos0  l1 sin1 cos0 )0  l1 cos1 cos01
y  (L sin  l sin cos )  l cos sin 
0
1
l1
0 1
o
x
y
L0 :arm length
l1 :length from pendulumpivot to c.o.g
行列形式に一般化できる…運動拘束条件
 1 0 0 l1c1c0 
 L0 c0  l1s1c0
Ca q a  0
Ca :  L0 s0  l1s1c0 0
 1 0  l1c1s0 
①

0
0
0  1  l1s1 

qa  0 x y z 1

T
si : sini ,ci : cosi
拡大状態量qa と一般化(独立)変数
拡
般 (
)変数 q    0 1
1

 L c l s c
q a  Da q
 0 0 11 0
② Da :  L0 s0  l1s1c0

0

0

①と②の条件から
Ca Da  0
T との関係
関係

l1c1c0 
 l1c1s0 

 l1s1 
1 
0
③
①から Ca は, 速度場 qa に直交することから, ラグランジェ乗数 a を
導入すると 拘束システムは次式で表現できる
導入すると、拘束システムは次式で表現できる
M a qa  ha  CaT  a
④
④から 内力の表現式が導ける
④から,
a  Ca M a1CaT  Ca qa  Ca M a1ha 
1
④ DaT を掛け, ③の条件式 Ca Da  0 を利用すると, 独立変数 q
④に
に関する(最小次数での)運動方程式が得られる
DaT M a Da  q  DaT M a D a  q  DaT ha
→Lagrange法による運動方程式と同じ結果を与える
g g 法
る運動方程式 同 結果を与 る
◆ PJ法のまとめ
拘束システム表現
M a ( qa )  qa  ha ( qa , q a , )  Ca ( qa )  
T
拘束条件
Ca ( p )  qa  0
内力
  ( Ca M a1CaT )1  ( C a q a  Ca M a1ha )
 ( Ca M a1CaT )1  ( C a Da q  Ca M a1h   Ca M a1hc )
・PJ法は, Lagrange法では得られない, 内力を計算できる
・計算が容易である(偏微分などを必要としない)
◆ SDRE法応用例:
内力制御
・内力=状態依存項と入力依存項に分離可能
a  Ca M a1CaT  Ca qa  Ca M a1ha 
1
ha  hx x  hu u  hg
  (Ca M a1CaT ) 1   C a Da Ex  Ca M a1 (hx x  hu u )  q a  Da x
  ( x)  x   ( x)  u :バイアス除去後の拘束力
 ( x) : (Ca M a1CaT ) 1  (C a Da E  Ca M a1hx )
 ( x) : (Ca M a1CaT ) 1  Ca M a1hu
x : [ 0 1 0 1 ] T , E : [ 0 I ]
・評価関数修正
x Q( x) x   S ( x)  u R( x)u dt
  x Q ( x) x  2 x N ( x)u  u R ( x)u dt
J 

0

T
T
T
0
T
T
T
Q ( x) : Q ( x)   ( x)T S ( x) ( x)
N ( x) :  ( x)T S ( x)  ( x)
R ( x) :  ( x)T S ( x)  ( x)  R ( x)
・内力抑制制御則


u   R ( x) 1 N ( x)T  B ( x)T P ( x)  x
A( x)T P( x)  P( x) A( x)  Q( x)
 P( x) B( x)  N ( x) R( x) 1 B( x)T P( x)  N ( x)T  0


・評価関数の重み行列の選定
Q,Rについては, 基本的に前回までのSDRE法での場合と同じ
Q( x )  diag ( q1 , q2 , q3 , q4 )
q1  1  5 103 1  exp(10( 1   / 9 ))
q2  3 103 ( 2  cos 1 )
q3  1
q4  1  103 1  exp(10( 1   / 6 ))
R( x )  0.2  r1  r2
r1  100 exp 10( 1   / 2 )2 
r2  10 1  exp( 50( 1   / 6 ))
・内力に関する重み行列
S ( x)  diag( s1 , s1 , s1 )
s1  50  50 exp 10 | 1 |
y
1 -
0 +
1
80
s
link-1
( x1 , y1 , z1 )
x
w e ig h t in g f u n c
100
z
l1
L0
link-0
60
z
40
20
0
y
3
2 .5
2
1 .5
1
0 .5
0
u
o
x
+内力抑制制御
通常SDRE法
大きく 激しい動き
大きく,
小さく 柔らかな動き
小さく,
・グラフ解析
angle
g of ppendulum
lower
SDRE+内力制御
通常SDRE法
振上げ
upright
i h
internal force
減少
input
0
0.5
1
1.5 [s]
・内力を評価関数に加味したSDRE制御は, 通常のSDRE
制御よりも 穏やかな振上げ制御が実現できている.
制御よりも,
穏やかな振上げ制御が実現 き
る
◆ 例:
二重倒立振子の振上げ制御の場合
Weighting functions(for discrete time control)
Q ( x)  diag
g (1, q1 , q2 ,1, q3 , q4 )
q1  10 4 ( 2  0.9 cos 1 )
q2  10 4 f (1 )
q3  103 ( 2  cos 1 )  10 4 1  exp(10( 1   / 6)) 
q4  10 2 f (1 )
1  sin | 1 | | 1 |  / 2
f (1 )  
2
| 1 |  / 2

R( x)  2(1  cos1 )
S ( x )  diag ( 0.1, 0.1, 0, 0.1, 0.1, 0 )
・シミュレーション結果
[de
egree]
200
lower
swing-up
150
100
top link
50
0
-50
angle of links
upright
p g
bottom link
0
0.2
0.4
0.6
0.8
1
1.2
1.4
[deggree]
200
1.6
2
angle of motor
150
100
50
0
0
0.2
0.4
0.6
0.8
1
1.2
1.4
1.6
200
[N
N.m]
1.8
1.8
2
input
100
0
-100
-200
0
0.2
0.4
0.6
0.8
1
time [s]
1.2
1.4
1.6
1.8
2
◆ 参考文献
[1] S.Suzuki, K.Furuta, et.al. ''Nonlinear Optimal Constraint-Force Control and Application
to Swing-up/Stabilization of Pendulum,'' Transaction of the ASME: Journal of dynamics
systems measurement and control,
systems,
control vol.126,
vol 126 2004
.
[2] S.Suzuki, et.al, ``State-dependent Sliding-sector VS-control and Application to Swingup Control of Pendulum,'' in Proc. of the 42nd CDC, Hawaii, USA
.
[8] W. Blajer, ``A Projection Method Approach to Constrained Dynamic Analysis,'' J. of
Applied Mechanics, vol.59, pp.643--649, 1992.
VI. SDRE法の拡張
SDRE法の拡張
法 拡張
2.線形システム論の非線形化
2.
線形システム論の非線形化
- セクター可変構造制御(線形システム論)
- SDRE法による非線形化
◆ スライディングセクタ可変構造理論…線形システム設計
・リカッチ方程式を使った既存線形システム設計手法を, 非線形
システムへ拡張する手段として
シ
テ
拡張する手段
SDRE法を利用する
・スライディングセクタ可変構造制御(SS-VSC)
VS i
VS-input
t
1入力時不変
x (t )  Ax(t )  Bu (t )
P ノルム
P-
x P : xT Px
PR
PR-sector
t
PR-スライディングセクタ
= PP ノルムが入力0で減少する部分空間
sliding mode
・従来のVSCのようなチャタリングがない.
・入力0で安定化が促進できる領域を活用する
が
.
→最小干渉制御(lazy control)の可能性
Furuta & Pan,
Pan Automatica(2000)
◆ スライディングセクタとリカッチ方程式
PR-スライディングセクタ(PR-SS)
S  x s( x)   ( x) , x  R n


s ( x)  S x :スライディングモード
 ( x)  0 :スカラー関数
VS-input
PR-sector
PR-SSの一設計法…リカッチ方程式を利用
S  BT P,  ( x )  x T (Q  R ) x
AT P  PA  PBBT P  Q , (Q  QT  0)
SDREにより非線形へ拡張
x  f ( x)  g ( x)u
A(x)T P(x)  P(x) A(x)  P(x) B(x) B(x)T P(x)  Q(x)
S   x B T P( x) x   ( x) , x  R n


◆ 状態依存型SS-VSC
・アルゴリズム
0. SSの幅を   R nn で決定しておく.
1. 現状態量 x(t ) に対してSDLR { A( x(t )), B ( x(t )) } を導出
2. 重み行列 Q ( x(t )) , R ( x(t )) を計算
3. SDREを解く AT ( x) P  PA( x)  PBT BP  Q ( x)( R ( x)) 1
4. SS関数 s( x )  BT Px(t ) と (x) を計算.
5. 次式でVS入力を計算, プラントに入力印加
 BT PAx
PA  k ( x)  sgn( s ( x)) x(t ) S 
u (t )  
x(t ) S 
 BT PAx  k  s ( x)

6. 1へ戻り繰り返し
6
戻り繰り返し
Suzuki & Pan, CDC(2003)
◆ シミュレーション結果
状態依存
状態依存SS-VSC
通常のSDRE法
( a ) A n g le s
(a) Angles
6
3
base-link
pendulum
p
2
4
b a s e - lin k
p e n d u lu m
2
少ない移動量
1
0
0
0
0
0.2
0.4
0.6
0.8 (b) Input
1
1.2
1.4
1.6
1.8
1
( b ) In p u t
1 .5
2
0 .5
1
[s ]
1 .5
2
100
[N.m]
100
[N.m]
0 .5
2
0
0
小さい入力
-1 0 0
-100
0
0
0.2
0.4
0.6 (c) Boundary
0.8
1 of sectors
1.2
1.4
1.6
1.8
2
[no unit]
30
ほぼセクタの内側
20
大きなバンバン状入力
|s(x)|
(x)
10
0
0
0.2
0.4
0.6
0.8
1
[s]
1.2
1.4
1.6
1.8
2
状態依存SS-VSC
SS VSCは, 倒立成功と同時にエネルギ消費低減可能
2
入力累積値  0 |  (t ) | dt 148.5  80 [ Nms]
◆ 実験結果
倒立 安定化
倒立→安定化
(a) Angles
3
base-link
pendulum
2
安定化中の外乱印加
1
0
(a) Angles
0
0.2
0.4
0.6
0.8 (b) Input
1
1.2
1.4
[N.m]
100
1.6
1.8
2
2
小さい入力
0
0
-1
-2
12
-100
100
0
0.2
0.4
0.6 (c) Boundary
0.8
1 of sectors
1.2
1.4
1.6
1.8
12.5
13
2
13.5
14
14.5
15
15.5
16
16.5
17
(b) Boundary of sectors
20
30
|s(x)|
(x)
ほぼセクターの内側
20
10
|s(x)|
(x)
SS外→復帰
15
[no unit]
[no
o unit]
base-link
pendulum
外乱
1
10
5
0
0
0.2
0.4
0.6
0.8
1
[s]
1.2
1.4
1.6
1.8
2
0
12
12.5
13
13.5
14
14.5
[s]
15
15.5
16
16.5
・振上げ~倒立の広域安定化が実現…SDREの特徴の継承
・滑らかで少ない入力で上記制御実現… SS-VSCの特徴の継承
◆ VI章のまとめ
・SDREの展開事例として以下の2つを紹介した.
((i)) 内力抑制制御
(ii) 状態依存型スライディングセクター制御法(状態依存SS-VSC)
・(i)内力抑制制御では, 機械系物理量と制御指標を同等に評価関数に
考慮する制御器設計法を示し, SDRE法の同時最適化問題への発展性
を示した.
・(ii)状態依存SS-VSCでは, SDRE法が線形システム設計論を非線形シ
ステムに拡張する手段として利用できることを示した.
17