![[PukiWiki]](image/pukiwiki.png "[PukiWiki]")
ODE(Open Dynamics Engine)とは†ODEは剛体ダイナミクスをシミュレーションするためのオープンソースで高性能なC/C++ライブラリです．その特徴は次のとおり
インストール†
シミュレーション空間†ODEのシミュレーション空間はWorldやSpace,Body,Geometryなどで構成される．
重力や何も設定しないときのBodyの空気抵抗などはWorldに対して設定する． また，Jointに関するパラメータとしてERPやCFMなどを設定する．これらは，
ODEのシミュレーション空間の構築例は以下 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 また，空間を破棄するときは， 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();
}
}
剛体の追加†BodyをODEのシミュレーション空間に登録する手順は以下．
立方体形状を登録するサンプルコードを以下に示す． 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つを指定できる．
タイムステップ†ODEのシミュレーションステップを進めるには，dWorldQuickStep関数などを用いる． 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関数の例は以下． //! 衝突コールバック関数用のデータ
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での描画†ODEは独自の描画フレームワークdrawstuffを含んでいるが， より詳細な描画コントロールを行い時はOpenGLなどを使うこともできる． 基本的には，
Box,Sphere,Capsule,Cylinder形状に対してOpenGLで描画を行うコード例は以下である． 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の変換行列について†OpenGLでプリミティブを描画するときは，dGeomGetPositionで現在の位置，dBodyGetRotationで回転行列などを 取得し，glTranslatefやglMultiMatrix()を使って描画する． このとき，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参照)．また，配列内の並びは row major である(これもODE WikiのFAQ参照)． 一方，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 であるので，変換が必要となる． つまり，以下のような関数を用いればよい． 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;
}
プリミティブ†球体(Sphere)†球体ジオメトリを作成する． 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)†dGeomID dCreatePlane(dSpaceID space, dReal a, dReal b, dReal c, dReal d); 平面の方程式 ax+by+cz+d=0 の各係数を指定．n = (a, b, c)が平面の法線，d は原点からの符号付距離． 円筒(Cylinder)†両端面がフラットな円筒形状． 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)†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)†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);
で行う． 関数など†位置と回転の取得†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． サンプル†立方体落下(OpenGL描画)†
リンク† |
