Curl Noise with OpenCL

どうも再び@jirohclです。そしてGPGPU Advent Calender@12日目です！
人数が少なく、下手をすればもう一度まわってくるらしいのですが完全にネタも尽きてはたしてどうなる事やら。
まぁ愚痴はアンドロメダ星雲にでもすっ飛ばしておくとして、今回はこんなものを実装してみましたネタです。

Curl Noiseと聞いてピンと来る方も多いと思いますが、今をときめくスクウェアエニックス様が作ったAgni's Philosophyでも使われていたパーティクルの技法です。
主な特徴としては、パーティクル同士の判定を行わないために超高速、剛体との衝突がなんかそれっぽく見える、といった所でしょうか。

これをどうOpenCLに持っていくかというと、パーティクルの更新が都合よく並列化に向いていたためupdate関数をカーネル化してやりました。

__kernel void update(__global float2*        in,
		     __global float2*        rx,
		     __global float *         t,
		     __global float2*        rt,
		     __global unsigned char* rti){
	int const id = get_global_id(0);
	float2 mid = in[id];
	mid = mid+0.5f*curl_vel(in[id],rx,t[0],rt,rti)*0.04f;
	in[id] += curl_vel(mid,rx,t[0],rt,rti)*0.04f;
}

ここで特筆する点は乱数テーブル(rt)をホストから渡している所でしょうか。
この乱数テーブルは円状に配置された点をswapしてテーブル化してあります。

次に

float2 curl_vel(float2 x,__global float2* rx,float t,__global float2* rt, __global unsigned char* rti)
{
	float const delta_x = (1e-4);
	float2 v;
	v.x =  -(potential((float2)(x.x,x.y+delta_x),rx,t,rt,rti)
	       - potential((float2)(x.x,x.y-delta_x),rx,t,rt,rti))/(2*delta_x);
	v.y = (potential((float2)(x.x+delta_x,x.y),rx,t,rt,rti)
	     - potential((float2)(x.x-delta_x, x.y),rx,t,rt,rti))/(2*delta_x);
	return v;
}

速度の計算です、パーティクルの現在位置から微小距離ずらし、ポテンシャル関数を用いて速度を求めます。

float potential(float2 px,__global float2* rx,float t,__global float2* rt, __global unsigned char* rti)
{
	float2 wind_velocity = (float2)(0.1f,0.5f);
	float const inv_env_sqr = 1/sqrt(0.17f);
	float numer= inv_env_sqr*cross2(px,wind_velocity);
	float denom= inv_env_sqr;
        //---------------hard coded rigid X( ---
	float2 dx=px-rx[0];//centrer
	float psi=cross2(px,rx[1]);
	float phi=((sqrt(dx.x*dx.x)+sqrt(dx.y*dx.y))-0.06f);//radius;
	float m=1/(sqrt(phi)+1e-18);
	numer+=m*psi;
	denom+=m;
	//--------------------------------------
	float base_psi=numer/denom;

	float d=phi/0.1f;
	float g=smooth_step(1.1-px.x);
	float turb=g*d*0.0026f*noise((px-t*wind_velocity)/0.1f,rt,rti);
	return base_psi+turb;
}

そしてこのポテンシャル関数の中で、回転するノイズ関数、剛体や風等の周りからの力の付与を計算することで流体を表現するらしいです。

実行結果


実装に疲れたので説明はおざなりですが中々面白い技術なので皆様試してみてはいかがでしょうか？

Windows用の実行ファイルになり、かつ他の環境で動くかどうか怪しいですが一応デモも置いておきます。
※PCが壊れた等の責任は負いかねますので悪しからず