Position Based Fluid
をテンプレートにして作成
[
トップ
|
新規
|
一覧
|
検索
|
最終更新
|
ヘルプ
]
開始行:
*位置ベース流体(Position Based Fluids) [#k29c9bfc]
粒子法を用いた流体シミュレーションにおいて重要な要素とな...
非圧縮性の強制である．
非圧縮性とは流体の密度が流れによって空間的にも時間的にも...
水のような液体はもちろん，
流速が音速より十分小さいならば空気も非圧縮性流体と考える...
(ちなみに音速に近くなると空気は圧縮する，
音速に達した飛行機が発する衝撃波いわゆるソニックブームは...
非圧縮性のための様々な方法(MPS,WCSPH,PCISPH,IISPHなど)が...
それらに比べて非常に高速かつ安定した計算が可能なのが位置...
である．
PBFは位置ベース法(Position-based method)¬e{Muller2007:...
の一種である．
位置ベース法では従来の力学ベースの方法と異なり，位置(物体...
力学ベースでは支配方程式に従い加速度->速度->位置と変更し...
それぞれで積分に伴う誤差の蓄積が発生するが位置に直接制約...
ただし，頂点や粒子の位置関係に基づく制約をかけるためどち...
あまり物理的ではないのではないかという議論もある．
**位置ベース法の基本 [#v82f65a4]
PBFの前に基本となる位置ベース法の考え方を説明する．
まず，移動する計算点&ref(pbf.eq1.gif,nolink,70%);に対して...
&ref(pbf.eq3.gif,nolink,70%);の条件が常に満たされていると...
粒子が微少量&ref(pbf.eq4.gif,nolink,70%);移動した後もこの...
&ref(pbf.eq5.gif,nolink,70%);となる．
これをテイラー展開してやれば，
#ref(pbf.eq6.gif,nolink,70%)
となる．2次以降の項を無視すれば&ref(pbf.eq4.gif,nolink,70...
#ref(pbf.eq7.gif,nolink,70%)
&ref(pbf.eq8.gif,nolink,70%);は制約として定義しているので...
いま，&ref(pbf.eq1.gif,nolink,70%);を粒子の中心座標とする...
粒子数を&ref(pbf.eq9.gif,nolink,70%);としたら3次元空間の...
一方で&ref(pbf.eq11.gif,nolink,70%);が各粒子ごとに計算さ...
上記の方程式の数は&ref(pbf.eq9.gif,nolink,70%);なので未知...
(疑似逆行列などを使ってできないことはない)．
そこで&ref(pbf.eq4.gif,nolink,70%);に対して別の制約をかけ...
計算点&ref(pbf.eq1.gif,nolink,70%);が運動方程式に従って移...
ここでの&ref(pbf.eq4.gif,nolink,70%);はその運動とは別の制...
そうすると元の運動(平行移動と回転)を阻害しないような移動...
言い換えると平行移動と回転運動量を保存するような移動方向...
制約&ref(pbf.eq11.gif,nolink,70%);がこれらの運動の状態に...
その傾き&ref(pbf.eq13.gif,nolink,70%);はこれらの運動方向...
つまり，&ref(pbf.eq4.gif,nolink,70%);の方向を&ref(pbf.eq1...
これはちょうど制約付きの最適化問題になるのでラグランジュ...
#ref(pbf.eq15.gif,nolink,70%)
とする．
これを上の&ref(pbf.eq4.gif,nolink,70%);に関する方程式に代...
#ref(pbf.eq16.gif,nolink,70%)
となる．求めた&ref(pbf.eq14.gif,nolink,70%);を&ref(pbf.eq...
それで計算点位置を修正してやればよい．
複数の計算点が絡む場合は1回の修正では制約を満たすことは難...
上記の式による制約を条件を満たすまで反復処理する．
**位置ベース法の粒子法への適用 [#b2340bb0]
さて，話を粒子法の場合に戻そう．
まず，非圧縮性条件が満たすべき制約(Constraint)を考える．
流体が圧縮しないということはその密度&ref(pbf.eq17.gif,nol...
つまり，条件としては&ref(pbf.eq19.gif,nolink,70%);である...
&ref(pbf.eq20.gif,nolink,70%); ・・・(1)&br;
となる．
ここでは粒子法の一種であるSPH法を使ってある粒子&ref(pbf.e...
#ref(pbf.eq23.gif,nolink,70%)
ここで&ref(pbf.eq24.gif,nolink,70%);は粒子&ref(pbf.eq21.g...
近傍粒子数を&ref(pbf.eq9.gif,nolink,70%);とすると粒子&ref...
関数となるので，&ref(pbf.eq14.gif,nolink,70%);の計算式は...
&ref(pbf.eq28.gif,nolink,70%); ・・・(2)&br;
&ref(pbf.eq29.gif,nolink,70%);は密度&ref(pbf.eq22.gif,nol...
&ref(pbf.eq30.gif,nolink,70%);も同様で，SPH法における勾配...
#ref(pbf.eq32.gif,nolink,70%)
となる．
この式で注意しないといけないのは&ref(pbf.eq31.gif,nolink,...
&ref(pbf.eq35.gif,nolink,70%);で微分したときに0にならない...
そうすると上記の式はさらに書き下すことができて，&br;
&ref(pbf.eq38.gif,nolink,70%); ・・・(3)&br;
となる．&ref(pbf.eq36.gif,nolink,70%);の時は&ref(pbf.eq39...
&ref(pbf.eq41.gif,nolink,70%);の場合は&ref(pbf.eq42.gif,n...
(マイナスの符号がついているのは&ref(pbf.eq43.gif,nolink,7...
また，粒子質量&ref(pbf.eq25.gif,nolink,70%);がすべての粒...
最終的に&ref(pbf.eq4.gif,nolink,70%);の計算ではすべて打ち...
式(2)と式(3)により各粒子のラグランジュ乗数&ref(pbf.eq45.g...
これらから各粒子の位置修正量&ref(pbf.eq46.gif,nolink,70%)...
このとき近傍粒子への影響と近傍粒子からの影響が非対称にな...
&ref(pbf.eq47.gif,nolink,70%); ・・・(4)&br;
先ほども説明したように質量一定ならばこちらの&ref(pbf.eq48...
&ref(pbf.eq45.gif,nolink,70%);と&ref(pbf.eq49.gif,nolink,...
SOR法の加速緩和係数と同じ役割を果たすようで，ない方が収束...
**張力の安定化 [#za0cd3e4]
SPH法では近傍粒子数が少ない場合に負の圧力が発生して粒子同...
(particle clusteringやparticle stackingと呼ばれている)現...
PBFではこれを解決するためにMonaghan
¬e{Monaghan2000:J. Monaghan, "SPH without a Tensile In...
が提案した以下のような項を式(4)に追加している．
#ref(pbf.eq51.gif,nolink,70%)
式中のパラメータの値としては&ref(pbf.eq52.gif,nolink,70%)...
式(4)の&ref(pbf.eq55.gif,nolink,70%);の部分を&ref(pbf.eq5...
Monaghanの元の論文では&ref(pbf.eq57.gif,nolink,70%);は加...
下記の私の実装では&ref(pbf.eq58.gif,nolink,70%);を掛けて...
**境界粒子 [#fb3d5633]
前述のparticle stackingの問題は固体境界付近でも発生する．
粒子と固体境界の単純な衝突処理だけだと固体内に粒子が存在...
境界付近に粒子が集まってしまう．
それを解決する一つの方法が固体境界内に位置固定の粒子(境界...
通常，境界粒子は有効半径に合わせて複数層配置するが，
Akinciらの方法¬e{Akinci2012:N. Akinci, M. Ihmsen, G. A...
を使えば，仮想的な体積を事前に計算しておくことで1層だけで...
境界に沿って粒子を配置した後，境界粒子だけを使って以下の...
#ref(pbf.eq59.gif,nolink,70%)
そして，PBFでも計算時に境界粒子の質量を&ref(pbf.eq60.gif,...
(詳しい式などは[[こちらの論文:http://slis.tsukuba.ac.jp/~...
*実装例 [#fc7280bb]
Visual Studio 2012 + CUDA7.5の環境で作成したコードは以下．
#ref(rx_pbf.zip)
-ビルドするのに必要なライブラリ
freeglut, GLEW, CUDA, (OpenVDB)
各ライブラリについては[[ライブラリのインストール]]を参照．
-簡単な説明＆注意事項
--GUIツールキットとしてFLTKを用いている．
--rx_fltk_glcanvas.hの55行目
#define RX_USE_GPU
をコメントアウトするとCPUで計算する．
--OpenVDBはポリゴンモデルを読み込んでその内部に境界粒子を...
--モデルの読み込み自体は[[3Dモデルファイルの入出力]]のrx_...
--Visual Studioから実行する場合は，作業ディレクトリを"../...
--"s"キー or Start/Stopボタンでシミュレーションスタート/...
--F1〜F9キー or Sceneメニューからシーン切替可能．binフォ...
これがシーンを記述するファイル．デフォルトでは4つのシーン...
--CUDAはv7.5で動作を確認している．
--ポリゴンモデルを固体境界とした場合に境界粒子配置のため...
終了行:
*位置ベース流体(Position Based Fluids) [#k29c9bfc]
粒子法を用いた流体シミュレーションにおいて重要な要素とな...
非圧縮性の強制である．
非圧縮性とは流体の密度が流れによって空間的にも時間的にも...
水のような液体はもちろん，
流速が音速より十分小さいならば空気も非圧縮性流体と考える...
(ちなみに音速に近くなると空気は圧縮する，
音速に達した飛行機が発する衝撃波いわゆるソニックブームは...
非圧縮性のための様々な方法(MPS,WCSPH,PCISPH,IISPHなど)が...
それらに比べて非常に高速かつ安定した計算が可能なのが位置...
である．
PBFは位置ベース法(Position-based method)¬e{Muller2007:...
の一種である．
位置ベース法では従来の力学ベースの方法と異なり，位置(物体...
力学ベースでは支配方程式に従い加速度->速度->位置と変更し...
それぞれで積分に伴う誤差の蓄積が発生するが位置に直接制約...
ただし，頂点や粒子の位置関係に基づく制約をかけるためどち...
あまり物理的ではないのではないかという議論もある．
**位置ベース法の基本 [#v82f65a4]
PBFの前に基本となる位置ベース法の考え方を説明する．
まず，移動する計算点&ref(pbf.eq1.gif,nolink,70%);に対して...
&ref(pbf.eq3.gif,nolink,70%);の条件が常に満たされていると...
粒子が微少量&ref(pbf.eq4.gif,nolink,70%);移動した後もこの...
&ref(pbf.eq5.gif,nolink,70%);となる．
これをテイラー展開してやれば，
#ref(pbf.eq6.gif,nolink,70%)
となる．2次以降の項を無視すれば&ref(pbf.eq4.gif,nolink,70...
#ref(pbf.eq7.gif,nolink,70%)
&ref(pbf.eq8.gif,nolink,70%);は制約として定義しているので...
いま，&ref(pbf.eq1.gif,nolink,70%);を粒子の中心座標とする...
粒子数を&ref(pbf.eq9.gif,nolink,70%);としたら3次元空間の...
一方で&ref(pbf.eq11.gif,nolink,70%);が各粒子ごとに計算さ...
上記の方程式の数は&ref(pbf.eq9.gif,nolink,70%);なので未知...
(疑似逆行列などを使ってできないことはない)．
そこで&ref(pbf.eq4.gif,nolink,70%);に対して別の制約をかけ...
計算点&ref(pbf.eq1.gif,nolink,70%);が運動方程式に従って移...
ここでの&ref(pbf.eq4.gif,nolink,70%);はその運動とは別の制...
そうすると元の運動(平行移動と回転)を阻害しないような移動...
言い換えると平行移動と回転運動量を保存するような移動方向...
制約&ref(pbf.eq11.gif,nolink,70%);がこれらの運動の状態に...
その傾き&ref(pbf.eq13.gif,nolink,70%);はこれらの運動方向...
つまり，&ref(pbf.eq4.gif,nolink,70%);の方向を&ref(pbf.eq1...
これはちょうど制約付きの最適化問題になるのでラグランジュ...
#ref(pbf.eq15.gif,nolink,70%)
とする．
これを上の&ref(pbf.eq4.gif,nolink,70%);に関する方程式に代...
#ref(pbf.eq16.gif,nolink,70%)
となる．求めた&ref(pbf.eq14.gif,nolink,70%);を&ref(pbf.eq...
それで計算点位置を修正してやればよい．
複数の計算点が絡む場合は1回の修正では制約を満たすことは難...
上記の式による制約を条件を満たすまで反復処理する．
**位置ベース法の粒子法への適用 [#b2340bb0]
さて，話を粒子法の場合に戻そう．
まず，非圧縮性条件が満たすべき制約(Constraint)を考える．
流体が圧縮しないということはその密度&ref(pbf.eq17.gif,nol...
つまり，条件としては&ref(pbf.eq19.gif,nolink,70%);である...
&ref(pbf.eq20.gif,nolink,70%); ・・・(1)&br;
となる．
ここでは粒子法の一種であるSPH法を使ってある粒子&ref(pbf.e...
#ref(pbf.eq23.gif,nolink,70%)
ここで&ref(pbf.eq24.gif,nolink,70%);は粒子&ref(pbf.eq21.g...
近傍粒子数を&ref(pbf.eq9.gif,nolink,70%);とすると粒子&ref...
関数となるので，&ref(pbf.eq14.gif,nolink,70%);の計算式は...
&ref(pbf.eq28.gif,nolink,70%); ・・・(2)&br;
&ref(pbf.eq29.gif,nolink,70%);は密度&ref(pbf.eq22.gif,nol...
&ref(pbf.eq30.gif,nolink,70%);も同様で，SPH法における勾配...
#ref(pbf.eq32.gif,nolink,70%)
となる．
この式で注意しないといけないのは&ref(pbf.eq31.gif,nolink,...
&ref(pbf.eq35.gif,nolink,70%);で微分したときに0にならない...
そうすると上記の式はさらに書き下すことができて，&br;
&ref(pbf.eq38.gif,nolink,70%); ・・・(3)&br;
となる．&ref(pbf.eq36.gif,nolink,70%);の時は&ref(pbf.eq39...
&ref(pbf.eq41.gif,nolink,70%);の場合は&ref(pbf.eq42.gif,n...
(マイナスの符号がついているのは&ref(pbf.eq43.gif,nolink,7...
また，粒子質量&ref(pbf.eq25.gif,nolink,70%);がすべての粒...
最終的に&ref(pbf.eq4.gif,nolink,70%);の計算ではすべて打ち...
式(2)と式(3)により各粒子のラグランジュ乗数&ref(pbf.eq45.g...
これらから各粒子の位置修正量&ref(pbf.eq46.gif,nolink,70%)...
このとき近傍粒子への影響と近傍粒子からの影響が非対称にな...
&ref(pbf.eq47.gif,nolink,70%); ・・・(4)&br;
先ほども説明したように質量一定ならばこちらの&ref(pbf.eq48...
&ref(pbf.eq45.gif,nolink,70%);と&ref(pbf.eq49.gif,nolink,...
SOR法の加速緩和係数と同じ役割を果たすようで，ない方が収束...
**張力の安定化 [#za0cd3e4]
SPH法では近傍粒子数が少ない場合に負の圧力が発生して粒子同...
(particle clusteringやparticle stackingと呼ばれている)現...
PBFではこれを解決するためにMonaghan
¬e{Monaghan2000:J. Monaghan, "SPH without a Tensile In...
が提案した以下のような項を式(4)に追加している．
#ref(pbf.eq51.gif,nolink,70%)
式中のパラメータの値としては&ref(pbf.eq52.gif,nolink,70%)...
式(4)の&ref(pbf.eq55.gif,nolink,70%);の部分を&ref(pbf.eq5...
Monaghanの元の論文では&ref(pbf.eq57.gif,nolink,70%);は加...
下記の私の実装では&ref(pbf.eq58.gif,nolink,70%);を掛けて...
**境界粒子 [#fb3d5633]
前述のparticle stackingの問題は固体境界付近でも発生する．
粒子と固体境界の単純な衝突処理だけだと固体内に粒子が存在...
境界付近に粒子が集まってしまう．
それを解決する一つの方法が固体境界内に位置固定の粒子(境界...
通常，境界粒子は有効半径に合わせて複数層配置するが，
Akinciらの方法¬e{Akinci2012:N. Akinci, M. Ihmsen, G. A...
を使えば，仮想的な体積を事前に計算しておくことで1層だけで...
境界に沿って粒子を配置した後，境界粒子だけを使って以下の...
#ref(pbf.eq59.gif,nolink,70%)
そして，PBFでも計算時に境界粒子の質量を&ref(pbf.eq60.gif,...
(詳しい式などは[[こちらの論文:http://slis.tsukuba.ac.jp/~...
*実装例 [#fc7280bb]
Visual Studio 2012 + CUDA7.5の環境で作成したコードは以下．
#ref(rx_pbf.zip)
-ビルドするのに必要なライブラリ
freeglut, GLEW, CUDA, (OpenVDB)
各ライブラリについては[[ライブラリのインストール]]を参照．
-簡単な説明＆注意事項
--GUIツールキットとしてFLTKを用いている．
--rx_fltk_glcanvas.hの55行目
#define RX_USE_GPU
をコメントアウトするとCPUで計算する．
--OpenVDBはポリゴンモデルを読み込んでその内部に境界粒子を...
--モデルの読み込み自体は[[3Dモデルファイルの入出力]]のrx_...
--Visual Studioから実行する場合は，作業ディレクトリを"../...
--"s"キー or Start/Stopボタンでシミュレーションスタート/...
--F1〜F9キー or Sceneメニューからシーン切替可能．binフォ...
これがシーンを記述するファイル．デフォルトでは4つのシーン...
--CUDAはv7.5で動作を確認している．
--ポリゴンモデルを固体境界とした場合に境界粒子配置のため...
ページ名:
![[PukiWiki]](image/pukiwiki.png "[PukiWiki]")
