2009年6月20日土曜日

粒子法のプログラム第6回(境界条件と粒子位置の更新)


今回は、粒子の境界条件と粒子位置の更新について説明します。

境界条件と粒子位置の更新のコードは以下のようになります。

void calculate_position( Particles* ps )
{
vec accel(3), g(3), norm(3);
double speed, diff, adj;

g(0) = 0.0; g(1) = -9.8; g(2) = 0.0;
FOR_EACH_PARTICLES( ps, p )
{
accel = p->f * SPH_PMASS;

speed = norm2( accel );
if ( speed > SPH_LIMIT*SPH_LIMIT ) {
accel *= SPH_LIMIT / sqrt(speed);
}

// Z-axis walls
diff = 2.0 * RADIUS - ( p->r(2) - MIN(2) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 0.0; norm(2) = 1.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(2) - p->r(2) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 0.0; norm(2) = -1.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

// X-axis walls
diff = 2.0 * RADIUS - ( p->r(0) - MIN(0) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 1.0; norm(1) = 0.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(0) - p->r(0) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = -1.0; norm(1) = 0.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

// Y-axis walls
diff = 2.0 * RADIUS - ( p->r(1) - MIN(1) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = 1.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}
diff = 2.0 * RADIUS - ( MAX(1) - p->r(1) ) * SPH_SIMSCALE;
if ( diff > EPSILON )
{
norm(0) = 0.0; norm(1) = -1.0; norm(2) = 0.0;
adj = SPH_EXTSTIFF * diff - SPH_EXTDAMP * inner_prod( norm, p->v );
accel += adj * norm;
}

accel += g;
p->v += accel * DT;
p->r += p->v * DT / SPH_SIMSCALE;

// stick on x-y plane
p->r(2) = 0.0;

}
}
それぞれの粒子について、粒子にかかる力から加速度を計算したあと、各壁面について境界条件を適用し、重力加速度を加えます。そうやって求まった加速度から、粒子の速度を求め、粒子の位置を更新します。

境界条件には、ペナルティ法を用いています。ペナルティ法とは、粒子が壁面へめりこんだ場合に、そのめり込みの程度によって壁面との衝突による反発力を与える方法です。

このようにして、境界条件を処理し、粒子の位置を更新しています。

次回は計算結果の出力です。最終回です。

粒子法のプログラム
粒子法のプログラム第1回(概要)
粒子法のプログラム第2回(プログラムの大枠)
粒子法のプログラム第3回(データ構造)
粒子法のプログラム第4回(密度と圧力の計算)
粒子法のプログラム第5回(力の計算)
粒子法のプログラム第6回(境界条件と粒子位置の更新)
粒子法のプログラム最終回(粒子の出力)

その他の解説エントリ
SPHによる巻き波のシミュレーション1
SPHによる巻き波のシミュレーション2
SPHによる巻き波のシミュレーション3
このあとやりたいこと

固液連成シミュレーションに関するエントリ
粒子法による固液連成シミュレーション

流体シミュレーションに関するエントリ
【粒子法】粒子を流体としてレンダリング
3次元の粒子法シミュレーション
粒子法のシーンを2倍のサイズにしてみたが…
粒子法のシーンを2倍のサイズにしてみた
Haskell、OCamlでSPH法
カメラ位置を変えて流体をレンダリング
Bunny-shaped fluid simulation

剛体シミュレーションに関するエントリ
粒子ベース剛体シミュレーション(プレビュー)
粒子ベース多体衝突シミュレーション
引き続き、粒子ベース剛体シミュレーション
Falling Rigid Bunnies

動画
シミュレーションの結果をレンダリングして作った動画です。流体シミュレーションや剛体シミュレーションの動画を見ることができます。

動画の一覧


--

2009年6月14日日曜日

Vivienne Westwoodの長財布


財布って、どんなのを使ってます?

私は父親からアメリカ土産にもらった革の三つ折り財布を10年以上使っているのですが、最近、長財布もいいなぁと思っていました。

それで今日、表参道に行って、Vivienne Westwoodの長財布を買ってきました。

買ったのは、写真にあるADVANのネイビー。表面に味のあるアドバン加工を施したレザーをベースに、全面に大胆なORB型の型押しを施してあります。ORBはVivienneのいろんなアイテムに施されていて、そこには「伝統をもって未来を創る」という彼女のコンセプトが込められています。

あと、外側のネイビーと内側のイエローのコンビネーションがとってもクール!この絶妙なコンビネーションが素晴らしい!!

ブティックで買い物をしたのが初めてだったので、ドキドキしたー

大事に長く使おう!