Document

Human interface Section, P&I Lab, Titech
実時間物理シミュレーション技術
東京工業大学精密工学研究所
長谷川晶一
Human interface Section, P&I Lab, Titech
内容
物理シミュレーションの意義
なぜ，今物理シミュレーションなのか？
物理シミュレーションでできること
シミュレータの仕組みと選定
接触力計算法，高速化法
制御
シミュレータ内の物体の動かし方
デモ
我々が開発している物理シミュレータ付きVR開発環境
Springheadを用いたデモ
なぜ物理シミュレーションが必要
か？
Human interface Section, P&I Lab, Titech
入力に対する反応の多様さ
→ゲーム世界の多様さ・楽しさ
入力の進化
アナログパッド，画像，力覚
ユーザの選択肢が増加
たくさんの反応を用意しなければならない
 従来の延長
 シミュレーション
たくさんの手間と時間
フレーム問題
自動的に多様でリアルな反応
Human interface Section, P&I Lab, Titech
物理シミュレーションが役立つ例
動きの生成
例：歩行
シミュレーションなし：足の動きと体の動きを別々に計算
シミュレーションあり：足を動かすと自然と体が動く
当たり判定
なし：当たり判定専用計算
あり：ポリゴンモデルの接触位置が毎ステップ求まる
ダメージ計算
なし：当たったもの，当たり方毎に値を用意
あり：加わった力の大きさからダメージを計算
効果音
なし：イベントごとに音を用意
あり：加わった力，材質から自動的に音を発生
Human interface Section, P&I Lab, Titech
物理シミュレーションの役割
ゲーム世界内の物体が，現実世界と同様に
自律的に動くようになる．
現実の物体と同様に，計測や制御をすること
ができる．
従来のゲームでも物理は使われていた？
従来は，特別に，物理法則を作りこんでいた．
物理シミュレータでは，物体は，すべて物理法則
にしたがって動作する．
Human interface Section, P&I Lab, Titech
物理シミュレータの選定
いくつもの物理シミュレータ
Havok， MathEngine, Tokamak （商用）
Open Dynamics Engine, Springhead （オープンソース）
いくつかの方式と特徴
スピードと精度のトレードオフ
何の精度を重視するか
物理シミュレータの中身を知る必要がある．
Human interface Section, P&I Lab, Titech
物理シミュレーションとは
物理法則(現実世界)は微分方程式で記述される
たとえば
質点の運動
剛体の運動
流体の運動
f  mx
 
mv  f , Iω
シミュレーションは微分方程式の数値解の計算
f  mx
運動方程式：
f
m
f t
差分方程式にすると： x (t  t )  x (t ) 
m
x(t  t )  x(t )  x (t )t
順番に求めて行く：
x(t )  x (t )t  x(0)
x(2t )  ...
x(3t )  ......
x
…
x(0) x(t) x(2t)
Human interface Section, P&I Lab, Titech
剛体運動のシミュレーション
剛体の運動の話だけにします．
剛体
硬いもの，変形しないもの
剛体だと考えてよいものが多い
積み木，ボール，ロボット，人体・・・
剛体でないもの
スポンジ，粘土，水・・・
Human interface Section, P&I Lab, Titech
剛体の運動
v: 速度
ω:角速度
m: 質量
I: 慣性テンソル
f: 外力
r: 外力の作用点の位置
(すべて絶対座標系)
運動方程式
dmv
 f
dt
dIω

dt
mv (t  t )  mv (t )  ft
I (t  t )ω (t  t )  I (t )ω (t )  t
f  0,  0 ならば，速度一定・角運動量一定
Human interface Section, P&I Lab, Titech
剛体に働く力
重力→ f=mg… 定数
バネ→ f=kx… 位置に比例
拘束力
mg
kx
力の大きさは不明
剛体同士の位置・速度関係が決まっている
 蝶番：2物体の相対位置が一定
 抗力：2物体が互いに侵入しない
 静止摩擦力：物体が滑らない
拘束力の計算が難しい
fn
ft
Human interface Section, P&I Lab, Titech
拘束力の求め方
例：球関節の拘束．
２物体が点pAと点pBで繋がっている
拘束の式： pA = pB
A
pA
pB
B
この式を満たすように，関節に働く力を求める
 A –p
B  0
方法１：運動方程式に拘束の式を含めて解く p
 解析法，David Baraff 89-93など
 Havok，Tokamak, Open Dynamics Engine
方法２：拘束違反に応じた力を加える
f  k (p A –p B )
 b (p A –p B )
 ペナルティー法，昔から用いられてきた
 Springhead
 接触判定時に拘束違反の量（侵入量）を調べる必要がある
Human interface Section, P&I Lab, Titech
解析法１（関節）
f
 B  p
 A  0
拘束： p
B
(rA,  θ A  (p A  rA )  θ A  (p A  rA ))
 (r  θ  (p  r )  θ  (p  r ))  0
pA pB
B,
A
IA
mB
M
B
B
B
B
B
Br  b  0
運動方程式：
 mA






B
1
 rA  

  



 θ A    (p A  rA )  
 r   
f
1
 B  







I B  θ B   (p B  rB )  
r 
C
f
Af  b  0
からfを求められる
Human interface Section, P&I Lab, Titech
解析法２（抗力）
 B  p
 A )  n  0
拘束： (p
B
pB
fn
運動方程式: Mr  Cfn
af  b  0
pA
A
(Br  b)  n  0
抗力は，反発だけ：
これを満たす f を見つける→
利点：
欠点：
f 0
線形計画法,2次計画法・・・
物体同士が侵入しない．
遅い，跳ね返り係数を考慮していない．
Human interface Section, P&I Lab, Titech
解析法３（摩擦力）
クーロンの摩擦モデル
f S  0 f N (静止摩擦)
f d   f N (動摩擦)
拘束：
B
運動方程式:
 B  p
 A )  n  0
(p
Mr  Cf
抗力は，反発力:
f n  0
摩擦力の条件： f n   / 1   2
A
場合分けを無くすため 0   としてしまっている
Human interface Section, P&I Lab, Titech
ペナルティ法（抗力）
 拘束を解かない．
Af  b  0
 拘束を侵した分だけ罰(ペナルティ)として力を加える．
.
侵入量 r，相対速度 r
f  Kr  Br  c
バネダンパ
 繰り返すうちに拘束が満たされる．．．はず．
力が直接決まるので，計算量はo(n)
接触部にバネとダンパを入れたと考えられる
利点：
欠点：
高速，跳ね返り係数を考慮できる．
物体同士が多少侵入してしまう．
Human interface Section, P&I Lab, Titech
ペナルティ法（摩擦力）
静止摩擦：ずれに比例した力を加える
接触の履歴を利用
クーロンの摩擦モデルをそのまま利用できる
f d   f N (動摩擦)
f S  0 f N (静止摩擦)
静止摩擦
動摩擦静止摩擦
Human interface Section, P&I Lab, Titech
ペナルティ法（面接触）
面で接触する場合
力は接触部分全体から発生
正確な接触力の計算
抗力・摩擦力の分布モデルを考える
静止摩擦←→動摩擦は一度に切り替わるとする
発生する力を接触部全体について積分
意外と簡単な式になる
解析法に比べ，とくに摩擦力が正確になる
Human interface Section, P&I Lab, Titech
計算量と高速化
ここまでで，物理の話は終わりです．
関節を持つ剛体
剛体同士の接触
のシミュレーションができるようになりました．
ここからは，リアルタイム動作に必要な，
計算量と高速化について話します．
Human interface Section, P&I Lab, Titech
解析法の計算量と接触数
接触している
多物体が
場合多くの拘束が働く．
関節で繋がれる
f1n
f2n
Af  b  0
f 0
f (Af  b)  0
T
行列Aの次元＝抗力 fi の数
行列を解くので，行列の次数nに対して，計算量はo(n3)
立方体の面接触１つ：４次元
１０個つむと：４０次元
Human interface Section, P&I Lab, Titech
解析法の高速化
 巨大な行列Aを作らず，
2物体単位で計算する
 繰り返し計算する
 繰り返し回数が少ないと精度が落ちる
Af  b  0
f 0
f (Af  b)  0
T
巨大行列Aを解いた場合
２体単位で，繰り返し
解く場合
Open Dynamics Engineのドキュメントより
Human interface Section, P&I Lab, Titech
ペナルティ法は高速
60
ペナルティ法（Springhead）
解析法（Open Dynamics Engine）
計算時間[ms]
50
40
30
20
10
0
0
5
10
ブロック数
15
Human interface Section, P&I Lab, Titech
構造を利用した，解析法の高速化
Articulated Body
多数の剛体を関節でつないだもの
動物や人間の体，ロボットなど
構造に特徴
 繋がっていない剛体の組が多い
剛体が輪になっていない→木構造
輪になっている例全部繋がっている例
Human interface Section, P&I Lab, Titech
解析法で求めると．．．
3
2
B1
A
D
Aの運動方程式： MArA  CA1 (f1 )
１の拘束の式： B1ArA  B1BrB  b1  0
Bの・・・・・・・・・・・・・・・・ここで， rやｆは６次元ベクトル
r  vx
v y
v z
 x
 y
 z T
並べてみると
MA





B
 1A



MB
MC
MD
B1B
B 2B
B 2C
B 3C
B 3D
C A1
 rA   0 
  

CB1 CB2
 rB   0 
CC2 CC3  rC   0 
  

C D3  rD    0 
 f    b 
 1   1 
 f 2    b 2 
  

f

b
 3   3 
0だらけ，すかすか行列，sparse行列
普通に解くより早い方法があるのでは？
Human interface Section, P&I Lab, Titech
解析法の高速化
Featherstoneの方法
巨大行列を作らない
葉から根に向かって１つずつ拘束力を求める
とても速いが，
構造を変えるのに多少時間がかかる
輪があると計算できない
→抗力計算には向かない
Human interface Section, P&I Lab, Titech
Featherstoneの方法
 リンクiとより葉側のリンク全体
の運動方程式が，
fi  Miri  Zi
のように，書ける．
Mi: リンクiとより葉側のリンクの
質量・慣性による項
Zi: 外力，重力，コリオリ力による項
Mi,Ziは加速度 r によらない
1. 葉のMiは既知なので，
葉から根に向かってMiを
計算する．
fi
2. 根ではfiは0なので，運動方程
式から加速度r を求める．
3. 根から葉に向かって，fi , r
を求める
リンクi
ri
葉のリンク
葉のリンク
リンクi とより葉側のリンク
根のリンク
Human interface Section, P&I Lab, Titech
Featherstoneの方法
先ほどのsparse行列で考えると，
3
2
B1
A
D
1. 葉から根に向かってMiを計算
MA





B
 1A



MB
MC
MD
B1B
B 2B
B 2C
B 3C
B 3D
C A1
 rA   0 
  



CB1 CB2
r
0
 B  

CC2 CC3  rC   0 
  

C D3  rD    0 
 f    b  f 3を消去
 1   1 
 f 2    b 2 
  

f

b
 3   3  r を消去
2.加速度rA を求める．
D
Human interface Section, P&I Lab, Titech
Featherstoneの方法
3.根から葉に向かって，fi , r を求める
MA





B
 1A



MB
MC
MD
B1B
B 2B
B 2C
B 3C
B 3D
C A1
 rA   0 
  

CB1 CB2
 rB   0 
CC2 CC3  rC   0 
  

C D3  rD    0 
 f    b 
 1   1 
 f 2    b 2 
  

f

b
 3   3 
rAが求まったので，
f1を求める
Human interface Section, P&I Lab, Titech
解析法の高速化
ループを含めた高速化方法
山根克 98 「構造可変なリンク系の動力学計算
法とヒューマンフィギュアの運動計算」
Open HRP (Humanoid Robotics Platform)に
使われているようです
（なぜかOpen HRPは非常に重い）
Human interface Section, P&I Lab, Titech
Springhead
我々が開発している物理シミュレータ
開発の動機
力覚インタフェースに使いたい
 超高速更新(>300Hz)，安定性重視
シミュレータでロボット対戦（ロボコン）をやりたい
 ２台以上
 フィールドに障害物
 押し合い → 摩擦が重要
従来のものは
 更新が遅すぎる
 摩擦のシミュレーションがいい加減すぎる
Human interface Section, P&I Lab, Titech
Springheadの選択
接触力計算はペナルティ法
高速性，安定性
接触領域に分布モデルを考える→摩擦の精度
多少剛体同士が侵入するが，あきらめる．
ロボットなど関節を持つ物はFeatherstone法
構造が変化しないので非常に高速．
ループのある機構は少ない．
ペナルティ法とは簡単に組み合わせられる．
シミュレーション法・高速化法の特性を考慮して，
シミュレータを選んでください．
Human interface Section, P&I Lab, Titech
ゲームの作成
シミュレータの限界
限界を見極めてモデリングする必要がある
安定性の限界
 極端に重いものと軽いもの
• 衝突時に速度が大きくなりすぎる
 質量と慣性モーメントの比率がおかしいもの
• 回転速度が極端に大きくなることがある
 極端に大きさが異なるもの
• 衝突判定の精度が問題となる
 時間刻み
• 衝突を見落としたり，バネダンパが不安定になる
計算速度の限界
 物体数，ポリゴン数，時間刻み
Human interface Section, P&I Lab, Titech
ゲームの作成
物理以外に必要なもの
物理シミュレーションを活かして作ると効果的
例：
ゲームのルール，ダメージの計算，．．．
接触力や関節に働く力の計測が効果的
人物，動物，車などの動き
自然の動きではなく，意図を持った動き
→ コントロール＝制御が必要
Human interface Section, P&I Lab, Titech
物体の制御
シミュレータ内の物体は，速度・慣性を持つ．
自動的に運動する
強制的に位置を指定すると・・・・
不自然な動き，非常に大きな力が発生
シミュレーションの意味がない
やわらかく，物体を動かすには？
物理を無視せず，実世界と同じく力を加えて動かす．
＝制御をする
Human interface Section, P&I Lab, Titech
PD制御
 目標位置に質点を持っていくには？
誤差に応じた力を加えてやる
例えば f=k(xg – x)では？
mx  k ( xg  x) 
m
x0
f
xg
初期位置目標位置
x  ( x0  xg ) cos( k / m t )  xg
運動方程式：f  mx
振動し続ける
誤差に比例:バネ
Proportional制御
 f  k ( xg  x)  b( x g  x ) では？
x  ( x0  xg ) e (b / 2)t cos( k / m t )  xg
減衰し，振動が止まる
微分に比例：ダンパ
Differential制御
Human interface Section, P&I Lab, Titech
PD制御の性質と調節
 バネ係数kとダンパ係数bでPD制御の性質が決まる
 バネ・ダンパと考えられるので，
k
 バネ係数 k
 大きいほど早く目標に近づこうとする
 小さいほど柔らかい
＝力が加わったときに動きやすい
目標位置
m
 ダンパ係数b
 負だと振動がどんどん大きくなる
 0だと振動が止まらない
 大きいほど振動しにくい
 大きいほど撃力が伝わる
 kが大きいほど早く動く．
 b>0ならば，収束する（振動が収まる）．
連続系ではこうなるが，シミュレータは離散系
そううまくはいかない．
b
Human interface Section, P&I Lab, Titech
シミュレータ上でのPD制御
シミュレータ上でのPD制御
時間が離散
m
位置 x
速度 v
0
目標位置
Human interface Section, P&I Lab, Titech
シミュレータ上でのPD制御
安定性の確認
行列A
An →0 ならば，このPD制御は安定になる
Human interface Section, P&I Lab, Titech
シミュレータ上でのPD制御
安定性の確認（つづき）
|Aの固有値| < 1 ならば， An →0 なので，Aの固有値を調べる
Human interface Section, P&I Lab, Titech
シミュレータ上でのPD制御
安定な k,b の範囲，早く静止する k,b
b
b=k
2
実際に良く使うのはこの辺
k
Human interface Section, P&I Lab, Titech
PD制御の実験
いろいろな k, b での挙動をご覧ください
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
Human interface Section, P&I Lab, Titech
衝突
実世界の物体は
互いに侵入しない．
跳ね返る．
再現するためには
衝突検知
剛体間に働く力の決定
Human interface Section, P&I Lab, Titech
衝突検知
衝突の検出
衝突しているかどうかすべての剛体について調
べる．
剛体の形状は，多面体で表現されている．
衝突検知を簡単にするため，凸形状に分割しておく
Human interface Section, P&I Lab, Titech
凸形状
凸形状の便利な性質
凸形状距離が極小となる点が1点
凸形状
最近傍点が簡単に求まる
非凸形状
GJK algorithm
E. G. Gilbert, D. W. Johnson and S. S. Keerthi
A Fast Procedure for Computing the Distance between
Complex Objects in Three-Dimensional Space (1988)
Human interface Section, P&I Lab, Titech
凸形状（GJK）
凸形状A上の点から，
凸形状B上の点へのベ
クトルを
原点を始点に並べると
ベクトルの終点の集合
も凸形状になる
Human interface Section, P&I Lab, Titech
凸形状（GJK2）
１
V0 : 凸形状内の任意の点
Wi :OViとOWiの内積が最小の点
Vi :三角形Wi-2 Wi-1 Wi内の点で原点に
一番近い点
２
３
４
Human interface Section, P&I Lab, Titech
凸形状２
凸形状の便利な性質２
凸形状の交差部分も凸形状
交差部分の形状が簡単に求まる
Half space
representation
D. E. Muller and F.P.Preparata:
“Finding the intersection of two convex” (1978)

Download Report