38 人が閲覧(直近 30 日) シミュレーション

PBDとは|Position Based DynamicsとXPBDの仕組み・実装コード・Unity/Blenderでの使い方

PBDとは|Position Based DynamicsとXPBDの仕組み・実装コード・Unity/Blenderでの使い方

PBD(Position Based Dynamics)は、力から加速度を求めて積分する代わりに、質点の位置を直接補正して拘束条件を満たす物理シミュレーション手法です。Matthias Müllerらが2006年のワークショップVRIPHYSで発表し、布・ロープ・髪・ソフトボディをゲームで安定して動かす方法として広まりました。後継のXPBD(Extended Position Based Dynamics)は、PBDの弱点だった「硬さが反復回数と時間刻みで変わる」問題を解いた拡張で、HoudiniのVellum、Blender 5.2のXPBD Solverノード、NVIDIAらのNewtonなど現行のツールが採用しています。

この記事では、PBDの処理手順と距離拘束の式をMüllerらのPBD原論文とXPBD論文に沿って説明し、Node.jsで実行した2本のコードの実測値で、PBDとXPBDの差、ガウス・ザイデル型とヤコビ型の差、反復とサブステップの配分、質量比の影響を示します。そのうえでUnity・Unreal Engine・Blender・Pythonでの使い方と、XPBDの後に出たVBD・AVBDとの違いを整理します。なお「PbD」はプライバシー・バイ・デザイン、「.pbd」は別形式のファイル拡張子、「PBD」は企業の略称としても使われますが、この記事が扱うのは物理シミュレーションのPBDだけです。

まとめ:PBDとXPBDの要点

  • PBDの考え方:速度で位置を予測し、拘束(距離・曲げ・衝突など)を満たすように予測位置を反復で補正し、補正後の位置から速度を逆算する。時間刻みを大きくしても発散しにくく、実装が短い。
  • PBDの弱点:拘束の硬さ k が反復回数と時間刻みに依存する。後述の実測では、同じ k=0.5 の鎖が反復5回で30.56%、80回で1.76%伸びた。
  • XPBD:拘束ごとにコンプライアンス α(剛性の逆数)とラグランジュ乗数 λ を持たせ、材料の硬さを α で指定できるようにした。実測では反復20回と80回の伸びが19.64%と19.62%で、理論値19.6%に一致した。
  • 計算量の配り方:反復回数よりサブステップ数を増やす。伸びない鎖では、拘束を処理する回数が合計20回でも「1ステップ×反復20回」の2.31%に対し「20サブステップ×反復1回」は0.12%だった。並列化しやすいヤコビ型は、この鎖ではガウス・ザイデル型の4倍の反復で同じ精度だった。
  • 苦手な条件:質量比の大きい系。末端に2,000倍の重りを吊ると、伸びないはずの鎖が反復20回で502%伸びた。質量比が大きい・剛性差が大きい系はVBDやAVBDを検討する。
  • 使える環境:UnityはObi、Unreal EngineはChaos ClothのXPBD拘束(利用経路ごとに状態を確認)、BlenderはXPBD Solverノード、HoudiniはVellum、PythonはNewtonの SolverXPBD やGenesisの PBDSolver。
  • 使わない場面:応力やひずみの数値を設計値として使う解析。材料モデルとメッシュを持つ有限要素法と実験値との照合が要る。

Position Based Dynamicsの定義と力ベース手法との違い

一般的な物理シミュレーションは力ベースです。ばねの力や重力を合計し、ニュートンの運動方程式で加速度を求め、速度、位置の順に積分します。この方式で硬いばねを陽解法で解くと、時間刻みが大きいときに補正が行き過ぎて振動が増え、最後は発散します。

PBDは速度の層を飛ばして位置を直接操作します。原論文の要旨は、剛体シミュレータの多くが速度を直接操作する(インパルスベース)のに対し、PBDは「速度の層も省いて位置に直接作用する」と説明し、主な利点に制御しやすさと、陽的積分の行き過ぎ(overshooting)を避けられる点を挙げています。著者4人は当時PhysXを開発していたAGEIAの所属で、論文ではこの手法でゲーム向け物理ライブラリのリアルタイム布シミュレータを作ったと述べています。

論文は積分方式の位置づけも書いています。反復1回なら陽解法に近く、反復を増やすと拘束系をいくらでも硬くでき、陰解法に近い振る舞いになります。ボトルネックは衝突判定からソルバーへ移ります。

PBDのアルゴリズム:1ステップの処理手順

原論文のアルゴリズムを、1フレーム分の処理として並べると次のとおりです。w は質量の逆数で、w=0 の点は固定点になります。

  1. 外力(重力など)で速度を更新する:v ← v + Δt·w·fext
  2. 速度を減衰させる(dampVelocities)
  3. 予測位置を求める:p ← x + Δt·v(陽的オイラーの1ステップ)
  4. 現在位置 x から予測位置 p への移動を調べ、衝突拘束を生成する
  5. すべての拘束について、予測位置 p を補正する処理を決められた回数だけ繰り返す(solverIterations)
  6. 補正後の位置から速度を逆算する:v ← (p − x) / Δt
  7. 位置を確定する:x ← p
  8. 摩擦や反発を速度に反映する(velocityUpdate)

中心は手順3・5・6・7の4つです。論文によると、6と7は陽解法のように未来へ外挿するのではなく、ソルバーが求めた「物理的に妥当な位置」へ点を移すため、不安定になりうるのはソルバー自身だけで、その安定性は時間刻みではなく拘束関数の形で決まります。

距離拘束と衝突拘束の作り方

距離拘束の補正式

布やロープの基本になる距離拘束を例に、手順5の中身を見ます。2点 p1、p2 の距離を静止長 d に保つ拘束は、次の式で表します。

C(p1, p2) = |p1 − p2| − d = 0

C が0でなければ、2点を結ぶ方向 n = (p1 − p2) / |p1 − p2| に沿って、質量の逆数の比で補正量を配分します。

  • Δp1 = −w1 / (w1 + w2) · C · n
  • Δp2 = +w2 / (w1 + w2) · C · n

重い点(w が小さい点)ほど動かず、固定点(w=0)は動きません。どちらの点も有限の質量を持つ場合、質量で重み付けした補正量の和 m1Δp1 + m2Δp2 は0になるので、この内部拘束だけで2点の重心が動くことはありません(原論文も運動量が保存されることを示しています)。曲げ・体積の各拘束も、拘束関数 C とその勾配を定義すれば同じ枠組みで扱えます。

衝突拘束と空間ハッシュによる近傍探索

衝突も拘束として解きます。手順4で「床より下に入ったら床の高さへ戻す」「2粒子の距離が半径の和より小さければ離す」といった不等式拘束を毎ステップ作り直し、手順5で他の拘束と一緒に解きます。めり込んだ点を正しい位置へ射影するだけで貫通を解消できる点は、原論文の要旨でもPBDの利点として挙げられています。

問題は、どの粒子の組を調べるかです。全組み合わせを調べると粒子数 n に対して n² に比例して重くなります。粒子の半径がそろっている場合は、空間を一定サイズのセルに区切り、セル座標をハッシュ値にして粒子を登録し、自分と隣接するセルの粒子だけを調べる空間ハッシュが定番です。Müllerの解説シリーズTen Minute Physicsは、第11回でこの空間ハッシュを、第15回でそれを使った布の自己衝突を扱っています。後述のFleXの論文も、粒子半径をシーンごとに固定して一様グリッドで衝突を検出する設計です。

JavaScriptによる鎖の実装例:PBDとXPBDの比較

次のコードは、長さ5cmの区間19本でつないだ質点20個の鎖を上端で固定して吊るし、600フレーム(60fpsで10秒)後の全長が静止長から何%伸びたかを出力します。mode=’pbd’ は硬さ k で補正量を縮める原論文の方式、mode=’xpbd’ はコンプライアンス α を使うXPBDの方式です。Node.js v26.5.0でそのまま実行できます。

// 鎖(質点N個)を上端固定で吊るし、静止後の全長の伸びを測る
function simulate({ n = 20, iters = 10, substeps = 1, mode = 'pbd', k = 0.5, alpha = 1e-4, steps = 600 }) {
  const dt = 1 / 60, h = dt / substeps, L = 0.05, g = -9.8, damp = Math.pow(0.99, 1 / substeps);
  const x = [], y = [], px = [], py = [], vx = [], vy = [], w = [], lambda = [];
  for (let i = 0; i < n; i++) { x.push(i * L); y.push(0); px.push(0); py.push(0); vx.push(0); vy.push(0); w.push(i === 0 ? 0 : 1); }
  for (let s = 0; s < steps * substeps; s++) {
    for (let i = 0; i < n; i++) { if (w[i] > 0) vy[i] += g * h; px[i] = x[i] + vx[i] * h; py[i] = y[i] + vy[i] * h; }
    for (let j = 0; j < n - 1; j++) lambda[j] = 0;
    for (let it = 0; it < iters; it++) {
      for (let j = 0; j < n - 1; j++) {
        const a = j, b = j + 1, wSum = w[a] + w[b];
        if (wSum === 0) continue;
        const dx = px[a] - px[b], dy = py[a] - py[b], d = Math.hypot(dx, dy);
        if (d === 0) throw new Error(`particles ${a} and ${b} overlap`);
        const C = d - L, nx = dx / d, ny = dy / d;
        let dl;
        if (mode === 'pbd') dl = -k * C / wSum;
        else { const at = alpha / (h * h); dl = (-C - at * lambda[j]) / (wSum + at); lambda[j] += dl; }
        px[a] += w[a] * dl * nx; py[a] += w[a] * dl * ny;
        px[b] -= w[b] * dl * nx; py[b] -= w[b] * dl * ny;
      }
    }
    for (let i = 0; i < n; i++) { vx[i] = (px[i] - x[i]) / h * damp; vy[i] = (py[i] - y[i]) / h * damp; x[i] = px[i]; y[i] = py[i]; }
  }
  let len = 0;
  for (let j = 0; j < n - 1; j++) len += Math.hypot(x[j] - x[j + 1], y[j] - y[j + 1]);
  return ((len / ((n - 1) * L) - 1) * 100).toFixed(2) + '%';
}
for (const iters of [5, 20, 80]) {
  console.log(`iters=${iters}  PBD(k=0.5)=${simulate({ iters, mode: 'pbd' })}  XPBD(alpha=1e-4)=${simulate({ iters, mode: 'xpbd' })}`);
}
console.log(`substeps=20,iters=1  XPBD=${simulate({ iters: 1, substeps: 20, mode: 'xpbd' })}  PBD=${simulate({ iters: 1, substeps: 20, mode: 'pbd' })}`);
console.log(`substeps=1,iters=20  XPBD=${simulate({ iters: 20, substeps: 1, mode: 'xpbd' })}`);

PBDとXPBDの差は、内側ループの dl を求める1行だけです。質点の質量は1kg(w=1)、長さの単位はm、コンプライアンス α の単位はm/N です。速度に掛けている damp は振動を止めるための簡易な減衰で、原論文の手順2(dampVelocities)の代わりです。サブステップ数を変えても1フレームあたりの減衰が0.99倍でそろうよう、0.99 のサブステップ数乗根にしています。2点が完全に重なると方向を定義できないため、その場合は例外を投げます。実行結果は次のとおりでした。

iters=5  PBD(k=0.5)=30.56%  XPBD(alpha=1e-4)=24.04%
iters=20  PBD(k=0.5)=7.49%  XPBD(alpha=1e-4)=19.64%
iters=80  PBD(k=0.5)=1.76%  XPBD(alpha=1e-4)=19.62%
substeps=20,iters=1  XPBD=19.62%  PBD=0.38%
substeps=1,iters=20  XPBD=19.64%

PBDの弱点:硬さが反復回数と時間刻みで変わる問題

上の結果のPBD列が、この手法の最大の弱点です。k=0.5 という同じ設定の鎖が、反復5回では30.56%伸び、20回で7.49%、80回で1.76%まで硬くなりました。1フレームを20サブステップに割って反復1回ずつにすると0.38%で、ほぼ伸びないロープになります。反復回数や時間刻みを変えるたびに、アーティストが決めた「布の硬さ」が別物になるということです。

原論文もこの問題を認めています。k を補正量に掛けるだけだと、1本の距離拘束の残差は ns 回の反復後に Δp(1 − k)ns となり、効き方が反復回数に対して非線形になります。そこで論文は k の代わりに k′ = 1 − (1 − k)1/ns を掛け、残差を Δp(1 − k) にそろえる補正を示しました。ただし同じ段落で、この補正後も材料の硬さは時間刻みに依存したままであり、固定時間刻みのリアルタイム環境なら問題にならない、と限定しています。

ゲームのように時間刻みも反復回数も固定で、見た目が破綻しなければよい用途なら、PBDのままでも困りません。硬さをパラメータとして資産に保存し、品質設定で反復回数を変えるような作り方をするなら、次のXPBDを使うべきです。

XPBDとは:コンプライアンスで硬さを反復回数から切り離す拡張

XPBDは、Miles Macklin、Matthias Müller、Nuttapong Chentanez(いずれもNVIDIA)が2016年のMotion in Games(MIG ’16)で発表した「XPBD: Position-Based Simulation of Compliant Constrained Dynamics」で提案されました。論文は要旨で、PBDの長年の問題である「反復回数と時間刻みに依存する拘束の硬さ」に取り組むと書いています。

XPBDの更新式とPBDからの移行

論文のAlgorithm 1は、元のPBDに3行(λ の初期化、Δλ の計算、λ の更新)を足しただけの形をしています。拘束ごとにコンプライアンス α(剛性の逆数)を持たせて α̃ = α / Δt² とし、拘束ごとの累積ラグランジュ乗数 λ を各ステップの初めに0へ戻してから積み上げます。1本の拘束を解くたびに求める増分は次の式です。

Δλj = (−Cj − α̃jλj) / (∇Cj M−1 ∇CjT + α̃j)

位置の補正は Δx = M−1∇CTΔλ です。論文が指摘するとおり、α=0(無限に硬い拘束)のときこの式は、硬さ係数 k=1 のPBDの補正と一致します。既存のPBD実装に数行足すだけで移行できるのはこのためです。上のコードの else 節がこの式そのもので、距離拘束では ∇C M−1 ∇CT が w1 + w2 になります。

鎖の伸びの実測値と理論値の照合

実測では、α=1e-4 の鎖の伸びは反復20回で19.64%、80回で19.62%、20サブステップ×反復1回で19.62%でした。区間 j にかかる張力はその下にぶら下がる質点の重さなので、α×張力の総和から静止時の伸びを手計算すると 1e-4 × 9.8 × (1+2+…+19) = 0.186m、静止長0.95mに対して19.6%です。反復20回以上とサブステップ20回の結果はこの理論値とよく一致しており、この鎖の静止時の伸びについては、反復回数や時間刻みを変えても α が決める値に落ち着くことが確認できます。反復5回の24.04%だけが外れているのは、反復が少なすぎて収束しきっていないためです。α は材料の硬さを反復回数から切り離しますが、収束しきらない分の数値誤差まで消すわけではありません。

拘束減衰(damping)の式

上のコードの damp は全速度を一律に縮める簡易版で、回転や剛体的な移動まで減速させます。XPBD論文の第5節は、拘束方向の速度だけを減衰させる方法を示しています。減衰の強さを β(コンプライアンスと違い、逆数ではなく減衰の強さそのもの)として β̃ = Δt²β と置き、γj = α̃jβ̃j / Δt を使って増分を次の式(論文の式26)に置き換えます。

Δλj = (−Cj − α̃jλj − γj∇Cj(xi − xn)) / ((1 + γj)∇CjM−1∇CjT + α̃j)

xi − xn はステップ開始時からの移動量で、∇Cj を掛けると拘束を伸び縮みさせる向きの動きだけが取り出されます。布の変形に伴う振動を抑えたい場合は、拘束方向の減衰を検討します。ただし、この式ではα=0にするとγも0になり、有限のβによる追加の減衰項は消えるため、伸びない拘束には別途減衰の設計が必要です。

ソルバーの組み方の実測:ヤコビ型・サブステップ・質量比

拘束の解き方には、拘束を1本ずつ解いてすぐ位置に反映するガウス・ザイデル型(原論文の方式)と、全拘束の補正量をいったん粒子ごとに集め、平均してから反映するヤコビ型があります。ヤコビ型はGPUで全拘束を並列に処理できる一方、収束が遅くなります。どれだけ遅いのかを、先ほどの鎖を伸びない設定(α=0)にして測りました。以下のヤコビ型はα=0での比較用です。αを正にする場合は、位置補正の平均化と乗数の更新が整合するよう実装を見直す必要があります。理想は伸び0%で、値が小さいほど拘束を正確に満たしています。ratio は末端の質点だけを何倍重くするかの指定です。

// 硬い鎖(XPBDのalpha=0)で、ガウス・ザイデル型とヤコビ型、反復とサブステップの配分、末端の重りの質量比を比べる
function simulate({ n = 20, iters = 20, substeps = 1, solver = 'gs', alpha = 0, omega = 1, ratio = 1, steps = 600 }) {
  const dt = 1 / 60, h = dt / substeps, L = 0.05, g = -9.8, damp = Math.pow(0.99, 1 / substeps), at = alpha / (h * h);
  const x = [], y = [], px = [], py = [], vx = [], vy = [], w = [], lambda = [];
  for (let i = 0; i < n; i++) { x.push(i * L); y.push(0); px.push(0); py.push(0); vx.push(0); vy.push(0); w.push(i === 0 ? 0 : i === n - 1 ? 1 / ratio : 1); }
  for (let s = 0; s < steps * substeps; s++) {
    for (let i = 0; i < n; i++) { if (w[i] > 0) vy[i] += g * h; px[i] = x[i] + vx[i] * h; py[i] = y[i] + vy[i] * h; }
    for (let j = 0; j < n - 1; j++) lambda[j] = 0;
    for (let it = 0; it < iters; it++) {
      const ax = new Array(n).fill(0), ay = new Array(n).fill(0), cnt = new Array(n).fill(0);
      for (let j = 0; j < n - 1; j++) {
        const a = j, b = j + 1, wSum = w[a] + w[b];
        const dx = px[a] - px[b], dy = py[a] - py[b], d = Math.hypot(dx, dy);
        const C = d - L, nx = dx / d, ny = dy / d;
        const dl = (-C - at * lambda[j]) / (wSum + at);
        lambda[j] += dl;
        if (solver === 'gs') {
          // ガウス・ザイデル型:解いた補正をすぐ位置に反映する
          px[a] += w[a] * dl * nx; py[a] += w[a] * dl * ny;
          px[b] -= w[b] * dl * nx; py[b] -= w[b] * dl * ny;
        } else {
          // ヤコビ型:補正を粒子ごとに集めるだけで、位置はまだ動かさない
          ax[a] += w[a] * dl * nx; ay[a] += w[a] * dl * ny; cnt[a]++;
          ax[b] -= w[b] * dl * nx; ay[b] -= w[b] * dl * ny; cnt[b]++;
        }
      }
      // ヤコビ型:集めた補正を拘束数で割って平均し、過緩和係数omegaを掛けて反映する
      if (solver === 'jacobi') for (let i = 0; i < n; i++) if (cnt[i] > 0) { px[i] += omega * ax[i] / cnt[i]; py[i] += omega * ay[i] / cnt[i]; }
    }
    for (let i = 0; i < n; i++) { vx[i] = (px[i] - x[i]) / h * damp; vy[i] = (py[i] - y[i]) / h * damp; x[i] = px[i]; y[i] = py[i]; }
  }
  let len = 0;
  for (let j = 0; j < n - 1; j++) len += Math.hypot(x[j] - x[j + 1], y[j] - y[j + 1]);
  return ((len / ((n - 1) * L) - 1) * 100).toFixed(2) + '%';
}
for (const [substeps, iters] of [[1, 20], [4, 5], [20, 1]]) {
  console.log(`substeps=${substeps} iters=${iters}  GS=${simulate({ substeps, iters, solver: 'gs' })}  Jacobi=${simulate({ substeps, iters, solver: 'jacobi' })}`);
}
for (const iters of [80, 200, 320]) console.log(`substeps=1 iters=${iters}  GS=${simulate({ iters, solver: 'gs' })}  Jacobi=${simulate({ iters, solver: 'jacobi' })}`);
for (const omega of [1.5, 2]) console.log(`substeps=1 iters=20  Jacobi(omega=${omega})=${simulate({ iters: 20, solver: 'jacobi', omega })}`);
for (const ratio of [200, 2000]) console.log(`ratio=1:${ratio}  1x20=${simulate({ ratio, iters: 20 })}  1x200=${simulate({ ratio, iters: 200 })}  20x1=${simulate({ ratio, iters: 1, substeps: 20 })}`);
substeps=1 iters=20  GS=2.31%  Jacobi=9.67%
substeps=4 iters=5  GS=0.61%  Jacobi=2.44%
substeps=20 iters=1  GS=0.12%  Jacobi=0.49%
substeps=1 iters=80  GS=0.47%  Jacobi=2.31%
substeps=1 iters=200  GS=0.12%  Jacobi=0.84%
substeps=1 iters=320  GS=0.04%  Jacobi=0.47%
substeps=1 iters=20  Jacobi(omega=1.5)=6.39%
substeps=1 iters=20  Jacobi(omega=2)=5.00%
ratio=1:200  1x20=53.90%  1x200=5.24%  20x1=2.69%
ratio=1:2000  1x20=502.18%  1x200=51.33%  20x1=26.57%

ガウス・ザイデル型とヤコビ型の収束差

同じ反復20回で、ガウス・ザイデル型は2.31%、ヤコビ型は9.67%でした。ヤコビ型で2.31%に届くのは反復80回、ガウス・ザイデル型の80回(0.47%)に届くのは320回で、この鎖ではヤコビ型にガウス・ザイデル型の4倍の反復が必要でした。中間の粒子は2本の拘束に挟まれており、補正を2で割って平均するため1回の反復で進む量が小さくなります。

この遅さを補うのが過緩和係数 ω です。GPU向けのPBDを提案した「Unified Particle Physics for Real-Time Applications」(SIGGRAPH 2014、FleXの元になった論文)は、補正量を拘束数で割る平均化に ω を掛ける形を示し、論文中のシミュレーションはすべて 1 ≤ ω ≤ 2 で動かしたとしています。実測では反復20回のヤコビ型が ω=1.5 で6.39%、ω=2 で5.00%まで改善しましたが、ガウス・ザイデル型の2.31%には届きません。GPUでもガウス・ザイデル型の収束を保ちたい場合は、点を共有しない拘束どうしをグループ(色)に分け、グループ内だけを並列に解くグラフ彩色を使います。

反復回数よりサブステップ数を増やす理由

同じ20回分の拘束計算を、1ステップ×反復20回に使うと2.31%、4サブステップ×反復5回で0.61%、20サブステップ×反復1回で0.12%でした。この0.12%は、1ステップで反復を200回回した場合と小数第2位まで同じです。拘束を処理する回数は10分の1ですが、位置予測や速度更新などの回数も変わるため、総計算量や実行時間が10分の1になるとは限りません。2019年のSCAで発表された「Small Steps in Physics Simulation」(Macklinほか、NVIDIA)の主張、すなわち1回の時間刻みで n 回反復するより時間刻みを n 分割して各1回反復するほうが拘束誤差と数値減衰が小さい、という結果と一致します。

サブステップごとに衝突判定をやり直すと費用が膨らむため、同論文のAlgorithm 1は衝突判定をフレームの先頭で1回だけ行い、サブステップ内ではその結果の拘束を解く構成です。論文は限界も書いています。Δt² の項があるため、サブステップを重ねると32ビット浮動小数点の精度の限界に達することがあり、大きな座標値に小さな位置の変化を足しても値が変わらない場合がある、というものです。原点から遠い座標で細かいステップを刻むときは、浮動小数点数の仕組みと誤差を踏まえて座標系の置き方を決めてください。

質量比が大きい系での破綻

末端の質点を200倍に重くすると、伸びないはずの鎖が反復20回で53.90%、2,000倍では502.18%伸びました。重い点はほとんど動かないため、軽い点の側だけで誤差を吸収しようとして収束が極端に遅くなるためです。反復を200回に増やしても2,000倍で51.33%、20サブステップでも26.57%残ります。サブステップのほうが効くことは、Unity向けアセットObiの公式マニュアルも「大きな質量比にはサブステップのほうがはるかに効果的」と書いています。重い荷物を細いロープで吊る、重い剛体を布で包むといった場面がこれに当たり、サブステップでも足りない場合は後述のVBD・AVBDが候補になります。

PBDの応用:布・髪・ソフトボディ・流体・剛体

布・ロープ・髪の距離拘束と曲げ拘束

距離拘束で布の伸びを、隣り合う三角形の二面角を保つ曲げ拘束で折れにくさを表します。Müllerらの2006年の原論文も、この2種類の拘束と自己衝突で布を実装しています。髪やロープは質点の列に距離拘束と曲げ拘束を並べたもので、上の鎖のコードがそのまま出発点になります。

流体の密度拘束とPosition Based Fluids

MacklinとMüller(NVIDIA)がSIGGRAPH 2013で発表し、ACM Transactions on Graphics 32(4)に掲載された「Position Based Fluids」(PBF)は、PBDの枠組みに密度の拘束を加えた流体の手法です。粒子 i ごとに Ci = ρi / ρ0 − 1(ρ0 は静止密度、ρi は近傍粒子からSPHのカーネルで推定した密度)という拘束を置き、密度を一定に保つことで非圧縮性を表します。論文の要旨は、現代のSPHソルバーに近い非圧縮性と収束性を持ちながら、PBDの安定性を引き継いで大きな時間刻みを使え、リアルタイム用途に向くと述べています。

剛体・関節のXPBDによる処理

剛体は2020年のSCAで発表された「Detailed Rigid Body Simulation with Extended Position Based Dynamics」(Müllerほか、NVIDIA)でXPBDに取り込まれました。論文は、接触と複数種類の関節、ソフトボディとの相互作用を扱う剛体ソルバーを実装するための詳細をすべて示したとしています。速度レベルで拘束を線形化する従来の剛体シミュレータと違い、常に最新の拘束方向で解くため、曲面に高速でぶつかる物体も追跡しやすくなります。前述のSmall Stepsのサブステップ方式と組み合わせる前提の設計です。

Unity・Unreal Engine・Blender・PythonでPBDとXPBDを使う方法

環境 機能・ライブラリ 方式 状態(2026年10月)
Unity 標準Cloth NvClothの位置拘束 方式の切り替え不可
Unity Obi(アセット) XPBD 有償アセット
Unreal Engine Chaos Cloth PBD/XPBD拘束 利用経路による(本文参照)
Blender XPBD Solverノード XPBD 5.2で搭載
Houdini Vellum XPBD Houdini 22.0
Python Newton XPBD・VBDほか 1.7.0・Apache-2.0
Python Genesis PBD 1.4.4
C++ PhysX 5 PBD(粒子) BSD-3-Clause
C++ PositionBasedDynamics PBD・XPBD MIT・タグは2.2.0

Unity:標準ClothとObi

Unity標準のClothコンポーネントは、Unity 2019.2のリリースノートのとおり、PhysX 4.x への移行に合わせてNVIDIAのNvClothソルバーに切り替わりました。Unity 6.6のマニュアルにある設定項目は1秒あたりの反復回数(Solver Frequency)などで、PBDかXPBDかを選ぶ項目はありません。硬さを反復回数から切り離したいなら、標準Clothではなく別の手段を選ぶことになります。

「xpbd unity」で探す人が実際に使えるのは、アセットのObiです。開発元Virtual Methodの公式フォーラムで開発者は「Obi is a full XPBD implementation」と回答しており、同じ回答で、接触拘束はコンプライアンスを持たないため通常のPBDに帰着するとも説明しています。マニュアルは拘束の反復回数を1にしてサブステップ数で品質を調整する設定を勧めており、前章の実測と同じ考え方です。自作する場合は、計算部分をC#へ書き換え、配列の初期化と状態の保持を更新処理から分離します。FixedUpdate ごとに1ステップを進め、時間刻みには Time.fixedDeltaTime を使い、計算した位置を描画へ反映します。Unityでゲームの物理表現から作り込む開発は、ゲームアプリ開発として外部に任せる選択肢もあります。

Unreal Engine:Chaos ClothのXPBD拘束

Unreal EngineのChaos Clothは、使う機能によってXPBD拘束の扱いが違います。Unreal Engine 5.3で加わったPanel Cloth Editorについて、公式ドキュメントは従来のクロスエディタと同じChaosソルバーを使いつつ、XPBDの拘束を任意で選べる拡張機能があると説明しています。一方、Unreal Engine 5.8のAPIドキュメントでは、UChaosClothSharedSimConfig の bUseXPBDConstraints が「反復回数と無関係に剛性を解くXPBD拘束を有効にする」と説明されています。ただし同じ欄に「Experimental, this feature might be removed without warning, not for production use」とあります。製品に組み込むならこの設定には頼らず、Panel Cloth EditorのXPBD拘束を候補にして、使うエンジンのバージョンでの扱いを公式ドキュメントで確かめてから採用します。

Blender・Houdini:XPBD SolverノードとVellum

Blenderは5.2でGeometry NodesにXPBD Solverノードを加え、5.2のリリースノートは、Cloth DynamicsとHair Dynamicsのアセットが内部でこのノードを使うと説明しています。ノードの入力はサブステップ数(Substeps)と拘束の反復回数(Constraint Iterations)で、マニュアルは出力される残差が1より大きければ反復回数を増やすよう勧め、サブステップは衝突判定を繰り返すため費用が増えると注意しています。開発者向けドキュメントは、布の機能は概念実証の段階で、高度な布シミュレーションに必要な機能がまだ欠けていると書いています。

映像制作で実績があるのはHoudiniのVellumです。Houdini 22.0のVellumの概要は「拡張Position Based Dynamicsのアプローチを使うシミュレーションの枠組み」と説明し、布・髪・ソフトボディ・風船・粒を作れるとしています。トラブルシューティングのページには、XPBDは布を硬くしすぎない一方で、反復が足りないと収束しないという注意があります。製品の選び方はHoudiniの特徴とBlenderとの違いで整理しています。

Python・GPU:Newton・Genesis・PhysX

PythonからGPUで動かすなら、Disney Research、Google DeepMind、NVIDIAが始めたLinux FoundationのプロジェクトNewtonが候補です。2026年10月14日公開のv1.7.0(Apache-2.0)は、ソルバーとして SolverXPBD と後述のVBDを実装した SolverVBD を同じAPIで提供しており、同じモデルで両者を比べられます。SolverXPBD の反復回数の既定値は2で、ソースのdocstringは、接触の補正を接触数で割る既定の重み付けにより、接触では運動量が保存されないと明記しています。

ロボット学習向けのGenesis(genesis-world 1.4.4)は、PBDSolverで布・弾性体・液体・粒子の4種類の材料を扱います。公式ドキュメントはXPBDではなくPBDと書いているため、反復回数を変えると硬さが変わる前提で設定を固定してください。C++ではPhysX 5に粒子をPBDで解く PxPBDParticleSystem があり、NVIDIAのOmniverse Physicsは布を旧来のPBDの粒子布から、XPBDで有限要素法の弾性モデルを解く面の変形体へ移行しています。PhysXを使うIsaac Simで布を扱うときは、使っている版がどちらのモデルかを先に確認してください。実装を読んで学ぶなら、Jan BenderらのC++ライブラリPositionBasedDynamics(MIT)が剛体・変形体・流体の例を一通り含みます。タグ付きリリースは2022年12月の2.2.0が最後で、リポジトリへのpushは2026年9月も続いています。FleXのリポジトリは2021年4月から更新が止まっているため、新規採用より手法の参考資料として読むのが妥当です。

XPBDの後に出た手法:VBD・AVBDとの違い

XPBDは拘束1本ずつに乗数 λ を持たせて解く方式です。これに対し、SIGGRAPH 2024の「Vertex Block Descent」(VBD、Chenほか)は、陰的オイラー法を最小化問題として書き、頂点1つずつ位置を更新するガウス・ザイデル反復で解きます。論文の要旨は、無条件に安定で、反復回数を制限すれば決まった計算予算に収められるとしています。著者のプロジェクトページは、XPBDが扱えなかった1:2000の質量比をVBDは問題なく扱えたと比較しています。前章の鎖で2,000倍の重りが502%伸びたのは、この弱点を再現したものです。

VBDにも弱点があり、無限に硬い拘束(伸びない関節など)をそのままでは扱えません。SIGGRAPH 2025の「Augmented Vertex Block Descent」(AVBD、Gilesほか)は、拡張ラグランジュ法を入れて硬い拘束を数値的に安定して扱えるようにし、剛体の積み重ねと摩擦、自由度を制限した関節、硬い物体と柔らかい物体の相互作用を、GPUで数百万個の物体までリアルタイムに解けると報告しています。AVBD論文は、VBDとXPBDはどちらも反復法だが、XPBDは拘束1本ずつ双対空間で解く点が違うと整理しています。研究側では、PBDを線形のコンプライアンス拘束で書けない非線形の力にも使えるよう一般化する研究(Towards Generalized Position-Based Dynamics、2025年11月)も出ています。

選び方ははっきりしています。Vellum・Blender・Obi・Chaos Clothのように使うツールが決まっているならXPBDで不足はなく、質量比による伸びが問題なら、まずサブステップを増やし、許容する拘束誤差と処理時間に収まるかを確認します。重い剛体と細いロープ、剛体の高い積み重ねのように質量比や剛性比が数百倍を超える系を自前で組むなら、VBD・AVBDを実装したライブラリ(Newtonの SolverVBD など)で同じシーンを比べてから決めてください。

有限要素法・SPH・質点ばねモデルとの比較と採用しない場面

有限要素法やSPHは物体をどう離散化するかの方法で、PBDやXPBDは拘束をどう解くかの方法です。同じ層の選択肢ではないため、前述のOmniverse PhysicsのようにXPBDで有限要素法の弾性モデルを解く組み合わせや、SPHの密度推定を使うPBFのような組み合わせが成り立ちます。

手法 解くもの 時間刻みへの強さ 精度を左右する条件 主な用途
PBD 位置(拘束) 強い 拘束モデル・反復数・時間刻み ゲームの布・ロープ
XPBD 位置+乗数λ 強い 拘束モデル・収束誤差・質量比 ゲーム・VFX・医療訓練
VBD・AVBD 頂点ごとの位置 強い 反復数・並列化の色数 大質量比の剛体・弾性体
質点ばね(陽解法) 力→加速度 弱い(硬いと発散) ばねモデル・積分法・時間刻み 教育・単純な弾性体
有限要素法(FEM) 応力・ひずみ 陰解法なら強い 材料モデル・メッシュ・解法 構造解析・CAE
SPH 粒子の密度・圧力 中(CFL条件) カーネル・粒子解像度・圧力解法 流体

質点ばねモデルは、ばね定数を上げるほど安定に必要な時間刻みが小さくなります。PBDはこの制約から自由になった代わりに、硬さの物理的な意味を失いました。XPBDはその意味をコンプライアンスとして取り戻しましたが、論文自身が「陰的な運動方程式を厳密に解いているとは言えない」と書いているとおり、近似解法です。

判断を分けるのは、結果の数値を保証する必要があるかどうかです。部材の応力や変形量を設計値として使う検証に、本記事のような距離拘束だけのモデルを使うべきではありません。材料モデルとメッシュ、収束判定、実験値との照合までそろえたCAEで使われる有限要素法の解析が必要です。XPBDを使う場合でも、拘束を有限要素法の弾性エネルギーから作り、収束するまで反復して誤差を評価しない限り、設計値の根拠にはなりません。反対に、60fpsなら1フレーム約16.7ms以内に結果を返す必要があり、多少の誤差より破綻しないことが優先されるゲーム・VR・インタラクティブな演出では、XPBDが第一候補になります。流体でPBFを選ぶ理由は、大きな時間刻みで安定させたい場合と、布や剛体と同じ拘束ソルバーで流体を連成させたい場合です。

よくある質問

PBDとXPBDの違いは何ですか?

XPBDは拘束ごとにコンプライアンス α と累積ラグランジュ乗数 λ を追加したPBDの拡張です。PBDでは設定した硬さの効き方が反復回数と時間刻みで変わりますが、XPBDでは材料の硬さを α で指定でき、反復が収束していれば結果は α で決まる値に近づきます。α=0 のXPBDは、硬さ係数 k=1 のPBDと同じ補正になります。

PBDは何の略ですか?

物理シミュレーションの文脈ではPosition Based Dynamics(位置ベース動力学)の略です。同じ綴りでも、個人情報保護のPrivacy by Design(PbD)や企業名の略称、「.pbd」という別形式のファイル拡張子を指す場合があり、これらはPosition Based Dynamicsとは無関係です。

UnityでXPBDを使うにはどうすればよいですか?

Unity標準のClothはNvClothのソルバーで動き、XPBDへの切り替え項目はありません。XPBDを使うなら、開発者が「full XPBD implementation」と説明しているアセットのObiを導入するか、本記事の鎖のコードのような数十行のソルバーをC#で自作します。

反復回数とサブステップ数はどちらを増やすべきですか?

サブステップ数です。本記事の伸びない鎖では、拘束を処理する回数を合計20回にそろえると、1ステップ×反復20回が2.31%、20サブステップ×反復1回が0.12%の伸びでした。「Small Steps in Physics Simulation」(2019年)も同じ結論で、Obiはサブステップを増やす調整を勧めています。一方、BlenderのXPBD Solverは残差が大きい場合に拘束の反復回数を増やすよう案内しているため、残差と処理時間を見て配分を決めます。

PBFとは何ですか?

Position Based Fluidsの略で、PBDの枠組みに粒子の密度を一定に保つ拘束を加えて非圧縮性の流体を表す手法です。MacklinとMüllerが2013年に発表しました。GenesisのPBDSolverの液体材料のように、PBD系のライブラリが流体を扱うときの基本になっています。

関連記事

お気に入りに入れた記事の一覧

一創の公式 SNS

X では毎朝の新着記事と IT ニュースの速報を、YouTube・TikTok・Instagram・LinkedIn・Facebook では記事の解説動画を配信しています。

この記事は以下の記事からリンクされています

資料請求

RELATED POSTS 関連記事

目次