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
拘束力の求め方
例:球関節の拘束.
2物体が 点pAと点pBで繋がっている
拘束の式: pA = pB
A
pA
pB
B
この式を満たすように,関節に働く力を求める
 A –p
B  0
方法1:運動方程式に拘束の式を含めて解く p
 解析法,David Baraff 89-93など
 Havok,Tokamak, Open Dynamics Engine
方法2:拘束違反に応じた力を加える
f  k (p A –p B )
 b (p A –p B )
 ペナルティー法,昔から用いられてきた
 Springhead
 接触判定時に拘束違反の量(侵入量)を調べる必要がある
Human interface Section, P&I Lab, Titech
解析法1(関節)
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
解析法2(抗力)
 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
解析法3(摩擦力)
クーロンの摩擦モデル
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)
立方体の面接触1つ:4次元
10個つむと:40次元
Human interface Section, P&I Lab, Titech
解析法の高速化
 巨大な行列Aを作らず,
2物体単位で計算する
 繰り返し計算する
 繰り返し回数が少ないと精度が落ちる
Af  b  0
f 0
f (Af  b)  0
T
巨大行列Aを解いた場合
2体単位で,繰り返し
解く場合
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 )
1の拘束の式: B1ArA  B1BrB  b1  0
Bの・・・・・・・・・・・・・・・・ ここで, rやf は6次元ベクトル
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の方法
巨大行列を作らない
葉から根に向かって1つずつ拘束力を求める
とても速いが,
構造を変えるのに多少時間がかかる
輪があると計算できない
→抗力計算には向かない
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),安定性重視
シミュレータでロボット対戦(ロボコン)をやりたい
 2台以上
 フィールドに障害物
 押し合い → 摩擦が重要
従来のものは
 更新が遅すぎる
 摩擦のシミュレーションがいい加減すぎる
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)
1
V0 : 凸形状内の任意の点
Wi :OViとOWiの内積が最小の点
Vi :三角形Wi-2 Wi-1 Wi内の点で原点に
一番近い点
2
3
4
Human interface Section, P&I Lab, Titech
凸形状2
凸形状の便利な性質2
凸形状の交差部分も凸形状
交差部分の形状が簡単に求まる
Half space
representation
D. E. Muller and F.P.Preparata:
“Finding the intersection of two convex” (1978)