実時間物理シミュレーション技術 東京工業大学精密工学研究所 長谷川晶一.

Slides:



Advertisements
Similar presentations
1 高速フーリエ変換 (fast Fourier transform). 2 高速フーリエ変換とは? – 簡単に言うとフーリエ変換を効率よく計算 する方法 – アルゴリズムの設計技法は分割統治法に基 づいている 今回の目的は? – 多項式の積を求める問題を取り上げ、高速 フーリエ変換のアルゴリズムを用いた解法.
Advertisements

Absolute Orientation. Absolute Orientation の問題 二つの座標系の間における剛体 (rigid body) 変換を復元す る問題である。 例えば: 2 台のステレオカメラから得られた3次元情報の間の関 係を推定する問題。 2 台のステレオカメラから得られた3次元情報の間の関.
1 運動方程式の例2:重力. 2 x 軸、 y 軸、 z 軸方向の単位ベクトル(長さ1)。 x y z O 基本ベクトルの復習 もし軸が動かない場合は、座標で書くと、 参考:動く電車の中で基本ベクトルを考える場合は、 基本ベクトルは時間の関数になるので、 時間で微分して0にならない場合がある。
Determining Optical Flow. はじめに オプティカルフローとは画像内の明る さのパターンの動きの見かけの速さの 分布 オプティカルフローは物体の動きの よって変化するため、オプティカルフ ローより速度に関する情報を得ること ができる.
物理シミュレーションの現状と未来 Physics simulation, state of the art and future
コンピュータビジョン特論 第8回対象追跡 2006年11月22日 加藤丈和.
HOG特徴に基づく 単眼画像からの人体3次元姿勢推定
CGアニメーションの原理 基本技術 対象物体の動きや変形の設定方法 レンダリング技術
電磁気学C Electromagnetics C 7/27講義分 点電荷による電磁波の放射 山田 博仁.
自己重力多体系の 1次元シミュレーション 物理学科4年 宇宙物理学研究室  丸山典宏.
・力のモーメント ・角運動量 ・力のモーメントと角運動量の関係
伝達事項 皆さんに数学と物理の全国統一テストを受けても らいましたが、この時の試験をまた受けていただ きます。
コリオリ力の復習資料 見延 庄士郎(海洋気候物理学研究室)
実時間物理シミュレーション技術 東京工業大学精密工学研究所 長谷川晶一.
定在波型熱音響エンジンにおける 臨界温度比推定のための適応制御系の 安定性に関する実験と理論の比較 長岡技術科学大学
エージェントモデル シミュレーション.
剛体の物理シミュレーション は難しい? 佐藤研助手 長谷川晶一.
スペクトル法による数値計算の原理 -一次元線形・非線形移流問題の場合-
重力3体問題の数値積分Integration of 3-body encounter.
上坂吉則 尾関和彦 文一総合出版 宮崎大輔2003年6月28日(土)
大阪工業大学 情報科学部 情報システム学科 宇宙物理研究室 B 木村悠哉
プロセス制御工学 6.PID制御 京都大学  加納 学.
次に 円筒座標系で、 速度ベクトルと加速度ベクトルを 求める.
工業力学 補足・復習スライド 第13回:偏心衝突,仕事 Industrial Mechanics.
1.Atwoodの器械による重力加速度測定 2.速度の2乗に比例する抵抗がある場合の終端速度 3.減衰振動、強制振動の電気回路モデル
衝突判定法をソート・検索アルゴリズムの観点から眺める
周期境界条件下に配置されたブラックホールの変形
エンタテインメントのための 力覚提示 東京工業大学 精密工学研究所 長谷川晶一.
最短路問題のための LMS(Levelwise Mesh Sparsification)
高山建志 五十嵐健夫 テクスチャ合成の新たな応用と展開 k 情報処理 vol.53 No.6 June 2012 pp
ストークスの定理と、 渦度・循環の関係を 直感で理解する方法
(ラプラス変換の復習) 教科書には相当する章はない
3次元剛体運動の理論と シミュレーション技法
動力学(Dynamics) 運動方程式のまとめ 2008.6.17
7-3.高度な木 (平衡木) AVL木 平衡2分木。回転操作に基づくバランス回復機構により平衡を保つ。 B木
物理学セミナー 2004 May20 林田 清 ・ 常深 博.
合成伝達関数の求め方(1) 「直列結合 = 伝達関数の掛け算」, 「並列結合 = 伝達関数の足し算」であった。
ロボティクス・メカトロニクスの基礎 東京大学 生産技術研究所 倉爪 亮.
Bottom-UpとTop-Down アプローチの統合による 単眼画像からの人体3次元姿勢推定
スペクトル法の一部の基礎の初歩への はじめの一歩
力覚計算と更新周期 電気通信大学知能機械工学科 長谷川晶一.
カオス水車のシミュレーションと その現象解析
6. ラプラス変換.
有効座席(出席と認められる座席) 左 列 中列 右列.
電磁気学C Electromagnetics C 7/17講義分 点電荷による電磁波の放射 山田 博仁.
物理学Ⅰ - 第 11 回 - 前回のまとめ 回転軸の方向が変化しない運動 回転運動のエネルギーとその応用 剛体の回転運動の方程式
システム制御基礎論 システム工学科2年後期.
物理学Ⅰ - 第 9 回 -.
ルンゲクッタ法 となる微分方程式の解を数値的に解く方法.
決定木 Decision Tree DT 明治大学 理工学部 応用化学科 データ化学工学研究室 金子 弘昌.
動力学(Dynamics) 力と運動方程式 2008.6.10
変換されても変換されない頑固ベクトル どうしたら頑固になれるか 頑固なベクトルは何に使える?
材質感提示のための振動を用いた力覚インタラクション環境の提案
第 6 章 :フィードバック制御系の安定性 6.1 フィードバック系の内部安定性 6.2 ナイキストの安定定理
実時間物理シミュレーション技術 東京工業大学精密工学研究所 長谷川晶一.
DPFのマスモジュールにおける残留ガス雑音の研究II
物理学Ⅰ - 第 7 回 - アナウンス 中間試験 第8回講義(6/16)終了前30分間 第7回講義(本日)(運動量)までの内容 期末試験
抗力への振動付加による 高剛性とすべり感提示
4. システムの安定性.
Bottom-UpとTop-Down アプローチの組み合わせによる 単眼画像からの人体3次元姿勢推定
ニュートン力学(高校レベル) バージョン.2 担当教員:綴木 馴.
宿題を提出し,宿題用解答用紙を 1人2枚まで必要に応じてとってください 配布物:ノート 2枚 (p.85~89), 小テスト用解答用紙 1枚
ガウス分布における ベーテ近似の理論解析 東京工業大学総合理工学研究科 知能システム科学専攻 渡辺研究室    西山 悠, 渡辺澄夫.
バネモデルの シミュレータ作成 精密工学科プログラミング基礎 資料.
力覚インタラクションのための 物理ベースモデリング
宿題を提出してください. 配布物:ノート 3枚 (p.49~60), 中間アンケート, 解答用紙 3枚 (1枚は小テスト,2枚は宿題用)
Cプログラミング演習 ニュートン法による方程式の求解.
* Ehime University, Japan
市松模様を使用した カメラキャリブレーション
逆運動学(Inverse Kinematics) 2007.5.15
Presentation transcript:

実時間物理シミュレーション技術 東京工業大学精密工学研究所 長谷川晶一

内容 物理シミュレーションの意義 シミュレータの仕組みと選定 制御 デモ なぜ,今物理シミュレーションなのか? 物理シミュレーションでできること シミュレータの仕組みと選定 接触力計算法,高速化法 制御 シミュレータ内の物体の動かし方 デモ 我々が開発している物理シミュレータ付きVR開発環境Springheadと力覚インタフェースSPIDARのデモ

なぜ物理シミュレーションが必要か? 入力に対する反応の多様さ →ゲーム世界の多様さ・楽しさ 入力の進化 アナログパッド,画像,力覚 入力に対する反応の多様さ →ゲーム世界の多様さ・楽しさ 入力の進化 アナログパッド,画像,力覚 ユーザの選択肢が増加 たくさんの反応を用意しなければならない 従来の延長 たくさんの手間と時間 フレーム問題 シミュレーション 自動的に多様でリアルな反応

物理シミュレーションが役立つ例 動きの生成 例:歩行 当たり判定 ダメージ計算 効果音 動きの生成 例:歩行 シミュレーションなし: 足の動きと体の動きを別々に計算 シミュレーションあり: 足を動かすと自然と体が動く 当たり判定 なし: 当たり判定専用計算 あり: ポリゴンモデルの接触位置が毎ステップ求まる ダメージ計算 なし: 当たったもの,当たり方毎に値を用意 あり: 加わった力の大きさからダメージを計算 効果音 なし: イベントごとに音を用意 あり: 加わった力,材質から自動的に音を発生

物理シミュレーションの役割 ゲーム世界内の物体が,現実世界と同様に自律的に動くようになる. 現実の物体と同様に,計測や制御をすることができる. 従来のゲームでも物理は使われていた? 従来は,特別に,物理法則を作りこんでいた. 物理シミュレータでは,物体は,すべて物理法則にしたがって動作する.

物理シミュレータの選定 いくつもの物理シミュレータ いくつかの方式と特徴 物理シミュレータの中身を知る必要がある. Havok, MathEngine, Tokamak (商用) Open Dynamics Engine, Springhead (オープンソース) いくつかの方式と特徴 スピードと精度のトレードオフ 何の精度を重視するか 物理シミュレータの中身を知る必要がある.

物理シミュレーションとは … 物理法則(現実世界)は微分方程式で記述できる シミュレーションは微分方程式の数値解の計算 m たとえば 質点の運動 剛体の運動 流体の運動 運動方程式: m f 差分方程式にすると: x 順番に求めて行く: … x(0) x(Dt) x(2Dt)

剛体運動のシミュレーション 剛体の運動の話だけにします. 剛体 硬いもの,変形しないもの 剛体だと考えてよいものが多い 剛体でないもの 積み木,ボール,ロボット,人体・・・ 剛体でないもの スポンジ,粘土,水・・・

剛体の運動 運動方程式 ならば,速度一定・角運動量一定 m: 質量 I: 慣性テンソル f: 外力 r: 外力の作用点の位置 v: 速度 ω:角速度 (すべて絶対座標系) 運動方程式 ならば,速度一定・角運動量一定

剛体に働く力 kx mg fn ft 重力→ f=mg… 定数 バネ→ f=kx… 位置に比例 拘束力 拘束力の計算が難しい 力の大きさは不明 剛体同士の位置・速度関係が決まっている 蝶番:2物体の相対位置が一定 抗力:2物体が互いに侵入しない 静止摩擦力:物体が滑らない 拘束力の計算が難しい kx mg fn ft

拘束力の求め方 例:球関節の拘束. 2物体が 点pAと点pBで繋がっている B 拘束の式: pA = pB この式を満たすように,関節に働く力を求める 解析法:運動方程式に拘束の式を含めて解く David Baraff 89-93など Havok,Tokamak, Open Dynamics Engine ペナルティ法:拘束違反に応じた力を加える 昔からいろいろな用途に使われてきた Springhead 接触判定時に拘束違反の量(侵入量)を調べる必要がある B pB pA A

解析法1(関節) f A pA B pB 拘束: 運動方程式: からfを求められる

解析法2(抗力) fn 拘束: B 運動方程式: A 抗力は,反発だけ: pB pA 利点: 物体同士が侵入しない. 欠点: 遅い,跳ね返り係数を考慮していない.

解析法3(摩擦力) クーロンの摩擦モデル 拘束: 運動方程式: B 抗力は,反発力: 摩擦力の条件: A 場合分けを無くすため としてしまっている

ペナルティ法(抗力) . 侵入量 r,相対速度 r 拘束を解かない. 拘束を侵した分だけ罰(ペナルティ)として力を加える. 繰り返すうちに拘束が満たされる...はず. 力が直接決まるので,計算量はo(n) 接触部にバネとダンパを入れたと考えられる 侵入量 r,相対速度 r . バネ ダンパ 利点: 高速,跳ね返り係数を考慮できる. 欠点: 物体同士が多少侵入してしまう.

ペナルティ法(摩擦力) 静止摩擦:ずれに比例した力を加える 接触の履歴を利用 クーロンの摩擦モデルをそのまま利用できる 静止摩擦 動摩擦

ペナルティ法(面接触) 面で接触する場合 正確な接触力の計算 力は接触部分全体から発生 抗力・摩擦力の分布モデルを考える 静止摩擦←→動摩擦は一度に切り替わるとする 発生する力を接触部全体について積分 意外と簡単な式になる 解析法に比べ,とくに摩擦力が正確になる

計算量と高速化 ここまでで,物理の話は終わりです. ここからは,リアルタイム動作に必要な, 計算量と高速化について話します. 関節を持つ剛体 剛体同士の接触 のシミュレーションができるようになりました. ここからは,リアルタイム動作に必要な, 計算量と高速化について話します.

行列を解くので,行列の次数nに対して,計算量はo(n3) 解析法の計算量と接触数 接触している 関節で繋がれる 多物体が 場合多くの拘束が働く. f1n f2n 行列Aの次元=抗力 fi の数 行列を解くので,行列の次数nに対して,計算量はo(n3) 立方体の面接触1つ:4次元 10個つむと:40次元

解析法の高速化 巨大な行列Aを作らず, 2物体単位で計算する 繰り返し計算する 繰り返し回数が少ないと精度が落ちる 巨大行列Aを解いた場合 2体単位で,繰り返し 解く場合 Open Dynamics Engineのドキュメントより

ペナルティ法は高速 60 ペナルティ法(Springhead) 50 解析法(Open Dynamics Engine) 40 30 計算時間[ms] 20 10 5 10 15 ブロック数

構造を利用した,解析法の高速化 Articulated Body 構造に特徴 多数の剛体を関節でつないだもの 動物や人間の体,ロボットなど 繋がっていない剛体の組が多い 剛体が輪になっていない→木構造 輪になっている例 全部繋がっている例

解析法で求めると... D C B A 3 2 1 並べてみると 0だらけ,すかすか行列,sparse行列 1の拘束の式: Bの・・・・・・・・・・・・・・・・ 並べてみると 0だらけ,すかすか行列,sparse行列 普通に解くより早い方法があるのでは?

解析法の高速化 Featherstoneの方法 とても速いが, 巨大行列を作らない 葉から根に向かって1つずつ拘束力を求める 構造を変えるのに多少時間がかかる 輪があると計算できない →抗力計算には向かない

Featherstoneの方法 リンクiとより葉側のリンク全体の運動方程式が, のように,書ける. Mi: リンクiとより葉側のリンクの 質量・慣性による項 Zi: 外力,重力,コリオリ力による項 Mi,Ziは加速度 によらない 葉のMiは既知なので, 葉から根に向かってMiを 計算する. 根ではfiは0なので,運動方程式から加速度 を求める. 根から葉に向かって, を求める 葉のリンク リンクi 葉のリンク 根側のリンク リンクi とより葉側のリンク 根のリンク

Featherstoneの方法 D C B A 先ほどのsparse行列で考えると, 3 2 1 1. 葉から根に向かってMiを計算 2.加速度 を求める.

Featherstoneの方法 3.根から葉に向かって, を求める

解析法の高速化 ループを含めた高速化方法 山根克 98 「構造可変なリンク系の動力学計算法とヒューマンフィギュアの運動計算」 Open HRP (Humanoid Robotics Platform)に 使われているようです (なぜかOpen HRPは非常に重い)

Springhead 我々が開発している物理シミュレータ 開発の動機 力覚インタフェースに使いたい 超高速更新(>300Hz),安定性重視 シミュレータでロボット対戦(ロボコン)をやりたい 2台以上 フィールドに障害物 押し合い → 摩擦が重要 従来のものは 更新が遅すぎる 摩擦のシミュレーションがいい加減すぎる

Springheadの選択 接触力計算はペナルティ法 ロボットなど関節を持つ物はFeatherstone法 高速性,安定性 接触領域に分布モデルを考える→摩擦の精度 多少剛体同士が侵入するが,あきらめる. ロボットなど関節を持つ物はFeatherstone法 構造が変化しないので非常に高速. ループのある機構は少ない. ペナルティ法とは簡単に組み合わせられる. シミュレーション法・高速化法の特性を考慮して, シミュレータを選んでください.

ゲームの作成 シミュレータの限界 限界を見極めてモデリングする必要がある 安定性の限界 計算速度の限界 極端に重いものと軽いもの 衝突時に速度が大きくなりすぎる 質量と慣性モーメントの比率がおかしいもの 回転速度が極端に大きくなることがある 極端に大きさが異なるもの 衝突判定の精度が問題となる 時間刻み 衝突を見落としたり,バネダンパが不安定になる 計算速度の限界 物体数,ポリゴン数,時間刻み デモしながら

ゲームの作成 物理以外に必要なもの 例: 物理シミュレーションを活かして作ると効果的 ゲームのルール,ダメージの計算,... 接触力や関節に働く力の計測が効果的 人物,動物,車などの動き 自然の動きではなく,意図を持った動き → コントロール=制御が必要

物体の制御 シミュレータ内の物体は,速度・慣性を持つ. 自動的に運動する 強制的に位置を指定すると・・・・ やわらかく,物体を動かすには? 不自然な動き,非常に大きな力が発生 シミュレーションの意味がない やわらかく,物体を動かすには? 物理を無視せず,実世界と同じく力を加えて動かす. =制御をする

PD制御 f m xg x0 目標位置に質点を持っていくには? 誤差に応じた力を加えてやる 例えば f=k(xg – x)では? 振動し続ける   減衰し,振動が止まる f m x0 xg 初期位置 目標位置 誤差に比例:バネ Proportional制御 微分に比例:ダンパ Differential制御

PD制御の性質と調節 k b バネ係数kとダンパ係数bでPD制御の性質が決まる バネ・ダンパと考えられるので, m バネ係数 k 大きいほど早く目標に近づこうとする 小さいほど柔らかい =力が加わったときに動きやすい ダンパ係数b 負だと振動がどんどん大きくなる 0だと振動が止まらない 大きいほど振動しにくい 大きいほど撃力が伝わる kが大きいほど早く動く. b>0ならば,収束する(振動が収まる).   連続系ではこうなるが,シミュレータは離散系 そううまくはいかない. m 目標位置 k b

シミュレータ上でのPD制御 シミュレータ上でのPD制御 時間が離散 m x v 位置 速度 目標位置

シミュレータ上でのPD制御 安定性の確認 行列A An →0 ならば,このPD制御は安定になる

シミュレータ上でのPD制御 安定性の確認(つづき) |Aの固有値| < 1 ならば, An →0 なので,Aの固有値を調べる

シミュレータ上でのPD制御 安定な k,b の範囲, 早く静止する k,b b k 2 b=k 実際に良く使うのはこの辺

PD制御の実験 b=k いろいろな k, b での挙動をご覧ください b 2 k 2

ペナルティ法とPD制御 . ペナルティ法もバネ・ダンパモデル バネ・ダンパ係数は? 侵入量 r,相対速度 r では多物体が接触したときに安定しない. Springheadでは, を使用 侵入量 r,相対速度 r . バネ ダンパ

衝突 実世界の物体は 互いに侵入しない. 跳ね返る. 再現するためには 衝突検知 剛体間に働く力の決定

衝突検知 衝突の検出 衝突しているかどうかすべての剛体について調べる. 剛体の形状は,多面体で表現されている. 衝突検知を簡単にするため,凸形状に分割しておく

凸形状 凸形状の便利な性質 凸形状 非凸形状 凸形状 距離が極小となる点が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) 凸形状 非凸形状

凸形状(GJK) 凸形状A上の点から,凸形状B上の点へのベクトルを 原点を始点に並べると ベクトルの終点の集合も凸形状になる

凸形状(GJK2) 1 2 3 4 V0 : 凸形状内の任意の点 Wi :OViとOWiの内積が最小の点 Vi :三角形Wi-2 Wi-1 Wi内の点で原点に 一番近い点 2 3 4

凸形状2 凸形状の便利な性質2 凸形状の交差部分も凸形状 交差部分の形状が簡単に求まる D. E. Muller and F.P.Preparata: “Finding the intersection of two convex” (1978) Half space representation