Real-Time Fluid Dynamics for Gamesを読んでみる(2)

Jos StamさんのStable Fluids関連の論文を読んできましたが、これで最後になります。
今回は、GDC’03で発表された、Real-Time Fluid Dynamics for Gamesのvelocity fieldの計算部分になります。

Evolving Velocities

velocityソルバーの説明にをします。もう一度、Fig.1の式を見てください。
densityソルバーの右辺の項については、既に我々は知っているので、それをvelocityの式にある3つの項にに応用します。外力、粘性による拡散、self-advection(訳注:自己移流?)に関する項です。
self-advectionに関しては、なじみが無いかもしれませんが、単に、velocity filedが自身の流れに沿って動くことと解釈できます。
もっと重要なことは、densityソルバーのために開発したルーチンをvelocity filedの更新に応用することが出来ると言うことです。
以下にコードを示します。force field(訳注:sourceと外力)が、配列u0, v0に格納されているとします。

void vel_step() (訳注:オリジナルの論文を参照して下さい)

densityの更新ルーチンと似ていることに注目して下さい。ほとんどの計算で、我々はvelocityの各コンポーネントのために、単に関数の呼び出しを複製しているだけです。
しかしながら、densityソルバーの説明で登場しなかった、新しい、project()というルーチンがあります。
このルーチンは、velocity filedを質量保存にさせます。これは、この法則を守ったリアルなfluidsでは重要な特性です。
ビジュアル的には、これは、多数の渦に流れを沿わせ、リアリスティックなグニャグニャとした流れを生成します。従って、これはこのソルバーの重要な項目です。

Fig.7 (訳注:オリジナルの論文を参照して下さい)

project()ルーチンの呼び出し直前では、velocity fieldはめったに質力保存の状態になっていません。要点は、この最後のステップで、velocity fieldを質量保存の状態にするということです。
これを実現するために、Hodge decompositionと呼ばれる、純粋な数学的成果「どんなvelocity fieldも、mass conserving fieldとgradient fieldの和であるということ」を利用します。
これを説明したものが、Fig.7(上)になります。mass conserving fieldが、我々の求めたいと思っている、ナイスなグニャグニャした渦巻き状になっていることに注目して下さい。対して、Fig.7の右上に示される、gradient fieldがワーストケースになっています。幾つかの点で、流れが、内向き、あるいは外向きになっています。
gradient fieldは、実は、とある高さの関数の、最も急な、下降勾配の向きを示すものになっています。
丘や谷がある地形の全ての場所に、最も急な下向きの勾配が示された矢印がある状態を想像してください。gradient fieldを計算するということは、height fieldを計算するということと同義です。
ひとたび、height filedを計算すれば、我々は、velocity filedから、そのgradient fieldを減算することができ、Fig.7(下)に示すように、質量保存のものを手に入れることが出来ます。
ここでは、数学的な詳細にはあまり触れませんが、height fieldの計算には、Poisson方程式と呼ばれる、線形システムが用いられています。このシステムも疎なので、ここでも我々は、densityの拡散を解く際に用いた、Gauss-Seidel法を用いています。
以下がprojection処理のコードになります。

void project(){} (訳注:オリジナルの論文を参照して下さい)

project()ルーチンが2度呼び出されていることに注目してください。これは、velocity fieldが質量保存であるときに、advect()ルーチンの振る舞いを、より正確にするために行っています。
境界の取り扱いについての説明が残っています。具体的には、コードの諸所で見られるset_bnd()ルーチンについてです。
我々は、fluidが、強固な壁の箱で覆われていて、fluidが外に流れ出ないと仮定しています。これは、単に、垂直な壁では、velocityの水平成分がゼロ、水平な壁では、velocityの垂直成分がゼロになるべきだということを意味します。
densityや他のfieldについては、単に連続的であると仮定しています。
以下のコードがこれ等の条件の実装になります。

void set_bnd() {} (訳注:オリジナルの論文を参照して下さい)

もちろん他の境界条件も適用可能です。たとえば、境界の一方から出た流れが、反対側から再入するような、循環している境界条件を定義することが出来ます。このケースを制御するのは、とても簡単です。読者のエクササイズとして残しておきます。このケースでは、advect()ルーチンも変更されるべきです。
他の可能性としては、風洞などで見られる流れ込みをシミュレートするために、幾つかの境界で、fixed velocityを適用したケースです。我々は、読者に他の境界条件について探求することをお勧めします。この章のまとめとして、我々のプロトタイプで、どのように我々のコードが使用されているかを示します。

while ( simulating) {} (訳注:オリジナルの論文を参照して下さい)

[ここまでで。感想]

先に、Siggraph’99の氏の論文に目を通しておいた方は、すんなり読めてしまうと思います。
Projectionの処理で使用されている、Poisson方程式は、やはり疎線形システムになるので、density fieldの拡散項の計算で使用した、Gauss-Seidel法を用いて解いています。projectionの処理を説明している部分で、たびたびmass conserving:質量保存という言葉が出てきますが、これはSig’99の論文中で示されていた、連続の式の事に他なりません。それ以外では特筆する部分がありません。

[全体の感想]

都合4回の記事に分けて、氏のStable Fluidsに関する発表を読んでみました。特にGDC’03の発表は、勤めて平易に書かれているので、読むのは簡単です。前提となる背景知識について先に触れておくと、さらに簡単に読むことが出来ます。
これ等は、現在のGPUならば、かなり高速に計算できることは、容易に想像がつきます。従って、ある程度限定された空間にのみ適用するならば、非常に高速にvelocity fieldの計算を行うことが可能だと思われます。
問題は、シミュレーションが適用できるのは、あくまで設定したvoxelグリッドの範囲で、また、グリッドの粒度が、流速が低下した際のシミュレーションの結果の綺麗さにつながると思う一方で、メモリ使用量と計算コストを考えるとそれほど細かいグリッドを採用することは難しいのかもしれないということです。
また、実際に計算していないので分かりませんが、細かいグリッド内で、速い流速を扱い、時間間隔dtを大きくとると、シミュレーション自体が意味消失(blow upはしないが、求めている結果にはならない)してしまうかも知れません。しかし、流速の影響を最も受けるのは、semi-Lagrange法のパーティクルトレーサーの部分なので、この部分をカスタマイズすればうまく行くのかもしれません。
ゲームに対しての組み込みという面では、外力をどのように適用するかと、干渉する物体をfield内でどのように取り扱うかという問題も残っていると思います。これに関しては、set_bnd()関数のように、逐次的に境界条件を適用してやる方法で解決できるかもしれません。

最後になりますが、本発表のHistorical Notesの章に触れられている通り、この手法は、
U. S.patent # 6,266,071 B1
で保護されているそうです。
数学的背景知識は、ずっと昔からあるものばかりで、加えて先行論文も結構あるので、いったいどの部分が特許となっているのかはよく分かりません。このあたりは専門家に相談した方が良いと思います。

コメントを残す

以下に詳細を記入するか、アイコンをクリックしてログインしてください。

WordPress.com ロゴ

WordPress.com アカウントを使ってコメントしています。 ログアウト / 変更 )

Twitter 画像

Twitter アカウントを使ってコメントしています。 ログアウト / 変更 )

Facebook の写真

Facebook アカウントを使ってコメントしています。 ログアウト / 変更 )

Google+ フォト

Google+ アカウントを使ってコメントしています。 ログアウト / 変更 )

%s と連携中