- 追加された行はこの色です。
- 削除された行はこの色です。
- ODE へ行く。
----
#contents
----
*[[ODE(Open Dynamics Engine):http://ode.org/]]とは [#s5d0d643]
ODEは剛体ダイナミクスをシミュレーションするためのオープンソースで高性能なC/C++ライブラリです．その特徴は次のとおり
-安定
-プラットフォーム非依存
-先進的なジョイントタイプの取り扱い
-摩擦まで考慮した衝突検出・応答
*インストール [#e04e30b2]
#include(build_ode,notitle)
*シミュレーション空間 [#u10e7ebc]
ODEのシミュレーション空間はWorldやSpace,Body,Geometryなどで構成される．
-World : ODE内の動力学ワールド．剛体(Body)を登録し，動力学に従い動作させる．
-Space : 形状データ(Geometry)を登録する空間．GeometryはBodyと関連づけられ，衝突判定に用いられる．
-Joint and Constraint : 2つのBody間の相対位置，方向のことで，単に接触している2物体の他に現実のジョイント(ドアのヒンジや車軸など)も含む．
重力や何も設定しないときのBodyの空気抵抗などはWorldに対して設定する．
また，Jointに関するパラメータとしてERPやCFMなどを設定する．これらは，
-ERP(error reduction parmeter) : 接続されているジョイントがユーザ操作やシミュレーションにより離れてしまった場合，
2接続点が近づくような力を負荷する．このときの力をコントロールするためのパラメータ．[0,1]．
-CFM(constraint force mixing) : 剛体間のめり込みを許す "ソフトな" ジョイントを許すためのパラメータ．
0なら完全に "ハードな" ジョイントとなり，正の0でない値なら "ソフトな" ジョイントとなる．
ODEのシミュレーション空間の構築例は以下
#code(C){{
dWorldID g_dWorld; //!< ODEワールド
dSpaceID g_dSpace; //!< ODE空間
dJointGroupID g_dContactGroupID; //!< ジョイントグループ
bool g_bWorldCreated = false;
void InitODE(void)
{
dInitODE2(0);
// ODE空間の破棄(既に存在していれば)
if(g_bWorldCreated){
dJointGroupDestroy(g_dContactGroupID);
dSpaceDestroy(g_dSpace);
dWorldDestroy(g_dWorld);
}
// ODE空間の生成
g_dWorld = dWorldCreate();
g_dSpace = dHashSpaceCreate(0);
g_dContactGroupID = dJointGroupCreate(0);
// ODE空間パラメータの設定
dWorldSetGravity(g_dWorld, 0, -9.8, 0); // 重力
dWorldSetCFM(g_dWorld, 1e-5); // グローバルCFM(constraint force mixing) : Typical Value [10^-9, 1.0]
dWorldSetContactMaxCorrectingVel(g_dWorld, 0.1);
dWorldSetERP(g_dWorld, 0.8); // グローバルERP
dWorldSetLinearDamping(g_dWorld, 0.00001); // 空気抵抗(平行移動速度)
dWorldSetAngularDamping(g_dWorld, 0.005); // 空気抵抗(角速度)
dWorldSetMaxAngularSpeed(g_dWorld, 200); // 最大角速度
dWorldSetContactSurfaceLayer(g_dWorld, 0.05); // 接触層の厚さ(この層内にいる物体は接触しているとする)
dWorldSetAutoDisableFlag(g_dWorld, 1); // しばらく動いていないbodyを自動でdisableする
g_bWorldCreated = true;
}
}}
ちなみに，デフォルトのパラメータ値はそれぞれ，
Gravity = (0, 0, 0)
CFM = 10^-5 (Typical values 0.1〜0.8)
ERP = 0.2 (Typical values 10^-9〜1)
LinearDamping = 0
AngularDamping = 0
ContactSurfaceLayer = 0
AutoDisableFlag = 0
また，空間を破棄するときは，
#code(C){{
void CleanODE(void)
{
if(g_bWorldCreated){
for(int i = 0; i < (int)g_vRObjs.size(); ++i){
dBodyDestroy(g_vRObjs[i]);
}
dGeomDestroy(g_dGround);
dJointGroupDestroy(g_dContactGroupID);
dSpaceDestroy(g_dSpace);
dWorldDestroy(g_dWorld);
dCloseODE();
}
}
}}
*剛体の追加 [#f6e02f66]
BodyをODEのシミュレーション空間に登録する手順は以下．
+dMassSet* でBodyの質量をその大きさから算出(計算はODE側で実装されているので密度と大きさを与える)
+dCreate* で形状(Geometry)を作成
+dBodyCreate で剛体オブジェクトを作成し，質量を設定(dBodySetMass)，形状を登録(dGeomSetBody)
+dBodySetPosition，dBodySetRotation，dBodySetLinearVelなどで位置，姿勢，速度などを設定
立方体形状を登録するサンプルコードを以下に示す．
#code(C){{
dBodyID SetODECube(Vec3 pos, double len)
{
// 立方体の重さを算出
dMass mass;
dMassSetBox(&mass, 1.0, (dReal)len, (dReal)len, (dReal)len);
dMassAdjust(&mass, 1.0);
// 形状の作成
dGeomID geom = dCreateBox(g_dSpace, (dReal)len, (dReal)len, (dReal)len);
// ボディの作成
dBodyID body = dBodyCreate(g_dWorld);
dBodySetMass(body, &mass);
dGeomSetBody(geom, body);
// 位置を設定
dBodySetPosition(body, pos[0], pos[1], pos[2]);
return body;
}
}}
ODEでは剛体の特性として以下の3つを指定できる．
-動力学(動作)剛体(Dynamics (moving) rigidbodies)
--BodyとGeometry
--毎シミュレーションステップで力学に基づきその位置と姿勢が更新される
-静的剛体(Static rigidbodies)
--Geometryのみ
--衝突が起こっても動かない
-運動学剛体(Kinematic rigidbodies)
--BodyとGeometry，dBodySetKinematicで設定．BodyがKinematicかどうかはdBodyIsKinematicで確認できる．
--ユーザ操作に従って動くだけで，上記二つとの衝突では影響を与えるが，影響は受けない．
*タイムステップ [#hd08003e]
ODEのシミュレーションステップを進めるには，dWorldQuickStep関数などを用いる．
#code(C){{
void StepODE(void)
{
NearCallData data;
data.world = &g_dWorld;
data.contact = &g_dContactGroupID;
dSpaceCollide(g_dSpace, &data, &ODENearCallback);
dWorldQuickStep(g_dWorld, g_fDt);
dJointGroupEmpty(g_dContactGroupID);
}
}}
ODEのステップを進める関数としては，dWorldStepの方がより精度が高いが，
大きな線形システム(サイズmとする)を解くため，最悪mの3乗のオーダの計算時間がかかる．
これに対して，dWorldQuickStepは反復的な手法を用いているので，m*Nのオーダで済む(Nは反復回数)．
ただし，精度は低くなることに注意が必要である．反復数Nは
dWorldSetQuickStepNumIterationsで設定，dWorldGetQuickStepNumIterationsで取得できる(デフォルトは20)．
さらに，dWorldQuickStepでは緩和法の一種であるSORを用いており，
その緩和加速係数を，dWorldSetQuickStepW/dWorldGetQuickStepWで設定/取得できる(デフォルトは1.3)．
ステップ関数では第二引数でタイムステップ幅を指定できるが，
ODEはタイムステップ幅が固定されていることを想定して設計されているため，
可変タイムステップは使用しない方がよい．
また，上記コード例では，NearCallback(剛体同士が近づいた時に呼び出される関数)も設定している．
NearCallback関数の例は以下．
#code(C){{
//! 衝突コールバック関数用のデータ
struct NearCallData
{
dWorldID* world;
dJointGroupID* contact;
};
/*!
* ODEの衝突処理
* 2つの物体の距離が近いときにODEから呼ばれるコールバック関数
* @param[in] data
* @param[in] o1 物体1のジオメトリID
* @param[in] o2 物体2のジオメトリID
*/
static void ODENearCallback(void *data, dGeomID o1, dGeomID o2)
{
// 衝突判定する物体のボディIDを取得
dBodyID b1 = dGeomGetBody(o1);
dBodyID b2 = dGeomGetBody(o2);
if(b1 && b2 && dAreConnected(b1, b2)){
return;
}
// 接触判定
dContact contact[RX_MAX_CONTACTS];
int n = dCollide(o1, o2, RX_MAX_CONTACTS, &contact[0].geom, sizeof(dContact));
if(n > 0){ // 1点以上での接触があった場合
for(int i = 0; i < n; ++i){
// 接触ごとにパラメータを設定
contact[i].surface.mode = dContactBounce | dContactSlip1 | dContactSlip2/* | dContactSoftERP | dContactSoftCFM | dContactApprox1*/;
contact[i].surface.mu = 1.0; // クーロン摩擦係数
contact[i].surface.bounce = 0.1; // 弾力性[0,1]
contact[i].surface.bounce_vel = 0.001;
contact[i].surface.slip1 = 0.0;
contact[i].surface.slip2 = 0.0;
// contact[i].surface.soft_erp = 0.1;
// contact[i].surface.soft_cfm = 0.01;
dJointID c = dJointCreateContact(*(((NearCallData*)data)->world), *(((NearCallData*)data)->contact), contact+i);
dJointAttach(c, dGeomGetBody(o1), dGeomGetBody(o2));
}
}
}
}}
*OpenGLでの描画 [#se80147c]
ODEは独自の描画フレームワークdrawstuffを含んでいるが，
より詳細な描画コントロールを行い時はOpenGLなどを使うこともできる．
基本的には，
+dBodyGetPosition,dBodyGetRotation でBody の位置，姿勢を取得し，glTranslatef, glMultMatrixf でOpenGL に設定
+Geometry をdBodyGetFirstGeom 関数でBody から取得
+dGeomGetClass 関数で形状の種類をGeometry から取得
+それぞれの形状に合わせてパラメータ(例えばBox ならdGeomBoxGetLengths で辺の長さ) を取得し，OpenGL で描画
の手順で行う．
Box,Sphere,Capsule,Cylinder形状に対してOpenGLで描画を行うコード例は以下である．
#code(C){{
void DrawODEBody(dBodyID body)
{
glPushMatrix();
Vec3 pos;
GLfloat m[16];
pos = ODE_GetPosition(body);
GetRotationMatrix(body, m);
glTranslatef(pos[0], pos[1], pos[2]);
glMultMatrixf(m);
// 形状の取得
dGeomID geom = dBodyGetFirstGeom(body);
// 形状タイプの取得
int type = dGeomGetClass(geom);
// 形状ごとの描画
if(type == dBoxClass){
dReal sl[4];
dGeomBoxGetLengths(geom, sl);
glScalef(sl[0], sl[1], sl[2]);
glutSolidCube(1.0);
}
else if(type == dSphereClass){
dReal rad = dGeomSphereGetRadius(geom);
glutSolidSphere(rad, 32, 16);
}
else if(type == dCapsuleClass){
dReal rad, len;
dGeomCapsuleGetParams(geom, &rad, &len);
DrawCapsule(rad, len, 20);
}
else if(type == dCylinderClass){
dReal rad, len;
dGeomCylinderGetParams(geom, &rad, &len);
DrawCylinder(rad, len, 20);
}
glPopMatrix();
}
}}
ここで，GetRotationMatrixは下で説明しているODEからOpenGLの変換行列を得る関数である．
また，DrawCapsule,DrawCylinderはそれぞれ，Capsule形状，Cylinder形状をOpenGLで描画する関数である
([[OpenGL - サンプルコード]]を参照)．
***ODEとOpenGLの変換行列について [#n88357bd]
OpenGLでプリミティブを描画するときは，dGeomGetPositionで現在の位置，dBodyGetRotationで回転行列などを
取得し，glTranslatefや[[glMultiMatrix():http://opengl.org/documentation/specs/man_pages/hardcopy/GL/html/gl/multmatrix.html]]を使って描画する．
このとき，ODEの回転行列は4x3であり，以下のような順序で配列に格納されている．
| m[0] m[1] m[2] m[3] |
| m[4] m[5] m[6] m[7] |
| m[8] m[9] m[10] m[11] |
回転行列が4x3になっているのは，SIMD friendlyにするためらしい([[ODE WikiのFAQ:http://opende.sourceforge.net/mediawiki-1.6.10/index.php/FAQ#Why_is_dVector3_the_same_size_as_dVector4.3F_Why_is_a_dMatrix_not_3x3.3F]]参照)．また，配列内の並びは row major である(これも[[ODE WikiのFAQ:http://opende.sourceforge.net/mediawiki-1.6.10/index.php/FAQ#I_have_a_matrix_does_dBodygetRotation_return_ABCDEFGHI_or_ADGBEHCFI_.4F]]参照)．
一方，OpenGLのモデルビュー変換行列は4x4で，
| m[0] m[4] m[8] m[12] |
| m[1] m[5] m[9] m[13] |
| m[2] m[6] m[10] m[14] |
| m[3] m[7] m[11] m[15] |
のように column major であるので，変換が必要となる．
つまり，以下のような関数を用いればよい．
#code(C){{
void GetRotationMatrix(dBodyID body, GLfloat m[16])
{
m[0] = dBodyGetRotation(body)[0];
m[1] = dBodyGetRotation(body)[4];
m[2] = dBodyGetRotation(body)[8];
m[3] = 0.0f;
m[4] = dBodyGetRotation(body)[1];
m[5] = dBodyGetRotation(body)[5];
m[6] = dBodyGetRotation(body)[9];
m[7] = 0.0f;
m[8] = dBodyGetRotation(body)[2];
m[9] = dBodyGetRotation(body)[6];
m[10] = dBodyGetRotation(body)[10];
m[11] = 0.0f;
m[12] = dBodyGetRotation(body)[3];
m[13] = dBodyGetRotation(body)[7];
m[14] = dBodyGetRotation(body)[11];
m[15] = 1.0f;
}
}}
*プリミティブ [#h7985a05]
***球体(Sphere) [#ra1ad75d]
球体ジオメトリを作成する．
dGeomID dCreateSphere(dSpaceID space, dReal radius);
球体の半径を指定する．これ以降の形状でも同様だが，これらの命令はジオメトリを作成するものなので，その位置などはジオメトリが関連付けされたボディ(dBodyID)で
void dBodySetPosition(dBodyID, dReal x, dReal y, dReal z);
void dBodySetRotation(dBodyID, const dMatrix3 R);
void dBodySetQuaternion dBodyID, const dQuaternion q);
void dBodySetLinearVel(dBodyID, dReal x, dReal y, dReal z);
void dBodySetAngularVel(dBodyID, dReal x, dReal y, dReal z);
などを用いて設定する．
球体のパラメータの設定と取得は，
void dGeomSphereSetRadius(dGeomID sphere, dReal radius);
dReal dGeomSphereGetRadius(dGeomID sphere);
で行う．また，ある点(x,y,z)が球体内部にあるかどうかなどの判定に用いることができるDepth値の取得には，
dReal dGeomSpherePointDepth(dGeomID sphere, dReal x, dReal y, dReal z);
を用いることができる．ここで，Depth値は球体内部で正，外部で負，球体表面でゼロとなる．
!ボックス(Box)
dGeomID dCreateBox(dSpaceID space, dReal lx, dReal ly, dReal lz);
ボックスのx,y,z軸長さを指定
***平面(Plane) [#c83ed287]
dGeomID dCreatePlane(dSpaceID space, dReal a, dReal b, dReal c, dReal d);
平面の方程式 ax+by+cz+d=0 の各係数を指定．n = (a, b, c)が平面の法線，d は原点からの符号付距離．
***円筒(Cylinder) [#v4eff95c]
両端面がフラットな円筒形状．
dGeomID dCreateCylinder(dSpaceID space, dReal radius, dReal length);
パラメータの設定と取得は
void dGeomCylinderSetParams(dGeomID cylinder, dReal radius, dReal length);
void dGeomCylinderGetParams(dGeomID cylinder, dReal *radius, dReal *length);
***円筒(Capped Cylinder, Capsule) [#c7fec56c]
dGeomID dCreateCapsule(dSpaceID space, dReal radius, dReal length);
球体を一方向にのばしたような形状の円筒(カプセル形状)で半径と長さを指定．dCreateCCylinder でも可
#define dCreateCCylinder dCreateCapsule
初期状態ではジオメトリのローカル座標のz軸方向に平行．質量設定の計算は，
void dMassSetCappedCylinder(dMass *, dReal density, int direction, dReal radius, dReal length);
void dMassSetCappedCylinderTotal(dMass *, dReal total_mass, int direction, dReal radius, dReal length);
の二つで行う．上が密度指定，下が全体質量指定．direction は円筒の length 方向がどの軸に平行かを設定(x=1, y=2, z=3)．
!光線(Ray)
dGeomID dCreateRay(dSpaceID space, dReal length);
他のジオメトリと異なり体積を持たない線を記述する．衝突時の接触点は1ヶ所になる．
***三角形メッシュ(Triangle Mesh) [#p4162b24]
dGeomID dCreateTriMesh(dSpaceID space, dTriMeshDataID Data,
dTriCallback *Callback,
dTriArrayCallback * ArrayCallback,
dTriRayCallback* RayCallback);
三角形メッシュは，メッシュデータ(dTriMeshDataID)を生成してから登録する．メッシュデータの生成は，
dTriMeshDataID dGeomTriMeshDataCreate();
void dGeomTriMeshDataBuild(dTriMeshDataID g, const void* Vertices,
int VertexStride, int VertexCount,
const void* Indices, int IndexCount,
int TriStride, const void* Normals);
で行う．
*関数など [#p83fe163]
***位置と回転の取得 [#laafdd3c]
const dReal* dBodyGetPosition(dBodyID);
const dReal* dBodyGetRotation(dBodyID);
const dReal* dBodyGetQuaternion(dBodyID);
の三つ．
dBodyGetPositionは位置の取得で，返値は dReal[3] ．(x,y,z)座標値．
dBodyGetRotationは回転行列(4x3)の取得で，返値は dReal[4*3] ．配置はソース中では R[(i)*4+(j)] (i=0,1,2, j=0,1,2,3)と定義されている．
dBodyGetQuaternionは四元数の取得で，返値は dReal[4] ．q = (cos(q0/2), sin(q0/2)*q1, sin(q0/2)*q2, sin(q0/2)*q3) = q0+q1*i+q2*j+q3*k．
*サンプル [#wfa7f9d2]
**立方体落下(OpenGL描画) [#bbf36d8a]
-&ref(box.zip);
--Visual Studio 2005用プロジェクトファイル同梱
--ODEに登録したBodyからShapeを参照して，OpenGLで描画する．
*リンク [#z3c20809]
-[[Open Dynamics Engine(本家):http://www.ode.org/]]
-[[Open Dynamics Engine v0.5 User Guide:http://ode.org/ode-latest-userguide.html]]
-[[ODE Wiki:http://opende.sf.net/wiki]]
![[PukiWiki]](image/pukiwiki.png "[PukiWiki]")
