仮眠プログラマーのつぶやき

自分がプログラムやっていて、思いついたことをつぶやいていきます。

自作ゲームやツール、ソースなどを公開しております。
①ポンコツ自動車シュライシュラー
DOWNLOAD
②流体力学ソース付き
汚いほうDOWNLOAD
綺麗なほうDOWNLOAD
③ミニスタヲズ
DOWNLOAD
④地下鉄でGO
DOWNLOAD
⑤ババドン
DOWNLOAD
⑥圧縮拳(ツール)
DOWNLOAD
⑦複写拳
DOWNLOAD
⑧布シミュレーション
DOWNLOAD
⑨minecraft巨大電卓地形データ
DOWNLOAD
⑩フリュードランダー
デジゲー博頒布α版
DOWNLOAD
⑪パズドラルート解析GPGPU版
DOWNLOAD
⑫ゲーム「流体de月面着陸」
DOWNLOAD

2011年07月

HSPで高速並列演算

(※並列にできるのは整数型、float型の四則演算のみです。命令や関数の並列処理は難or不可です。)


 前にブログに書いたようにhgimg3やHSPDXfixなどのプラグインは標準命令と比べて描画がとても早い。
 標準命令ではスペックに限界があり、例えば自分の作りたいシューティングゲームなどがとても60FPS出せない場合はこれらのプラグインを使うと解決されることが多い。

 たいていこういうプラグインは、ハードウェアに処理をさせているため高速に描画できます、とか書いてある。

o0380008511359222298


 その「ハードウェア」ってのはつまるところGPU(広義のグラフィックボードやビデオカード)のことで、端的に言えばこいつもCPUみたく四則演算ができる演算装置だ。ただCPUみたく高度な処理はきなくて、全画面のピクセルの書き換えとか、同時に何万もの同じ処理をしたいときなどに向いているとされている。

 GPUは本来グラフィック関連の処理に特化されたチップなので、その用途以外で使われることは考えられていなかったようなのだが、GPGPUと言って最近このGPUがCPUにかわって様々な処理をするようになってきた。 動画変換や計算科学(流体力学やブラックホールシミュレーションなど)といった、グラフィックに関係ないところに使われるようになってきて、それもまたCPUで処理するより何十、何百倍も早いそうなのだ。

 http://dic.nicovideo.jp/a/gpgpu (ニコニコ大百科がまとまってて分かりやすい!)
 後にも書くが、流体力学の計算をCPUとGPUでやらせたときの処理時間比が私のパソコンでは最大180倍になるという結果になった!
これは2^n×2^nの計算領域でCPUとGPUで流体計算させたときの1ループの処理時間(秒)

o0703046311359223918



 あまりにもGPUの計算時間が早くて、いまいちすごさがパッとしないので処理時間の比をとってみた


o0725041611359223917


 横軸が若干省略されてしまっているが、
512^2と2048^2の間にある1024^2の計算格子領域で一番GPU/CPU比が高く181倍も高速化ができた。
大岡山に存在するとされている某スパコンは、主にこのGPUで計算させているのでCPUのみのスパコンと比べすごく演算が早いとのこと。(もしCPUだけで同じ性能にしようとすると莫大なお金がかかる)
周波数が高いわけではなく、何万もの計算を並列に行なえるから結果的に早いのだそうだ。


・・・さて長い前置きは置いておいて、 私はHSPコンテストで分子シミュレーター
o0640048011359227067


とか
重力シミュレーター
o0431032111359227068

とか


まぁ普通にくだらないクソソフトを毎年量産しては応募しているのだが・・・
何がクソかというと、なんといっても処理が重いのだ!

これはプログラムが非効率なのもあるが、多数の粒子や相互作用などを計算すること自体がそもそも物理的に膨大な計算量を必要とする。

非力な私のノーパソのCPUでは、一億桁×一億桁の計算に何日もかかるというのに、分子のシミュレーションをリアルタイムでやろうとするなんぞ分不相応というか背伸びのしすぎであった・・・

今よりもCPUの性能が上がればいいのだが、生憎CPUのクロック数はもう頭打ちと言われていてせいぜい3.6Ghzが一般CPUの限界だから、今の1.5倍にしかならない。
HSPはマルチスレッド対応ではないため、私のやりたいシミュレーション系のプログラムはもうCPUだけでは限界なのだ。 多重起動の擬似マルチスレッド化とかして工夫はしてきたものやはりCPUに頼っているだけでは根本解決にはならない
そこで、本格的にGPGPUプログラミングに手を出そうと考えていた。

調べていくうちにGPGPUプログラミングのためにはnVIDIA製のGPU搭載デスクトップPCを買い、CUDAを勉強するのが一番手っ取り早いということがわかったのだが・・・

1 買う金がない
2 デスクトップ置く場所がない
3 CUDAを勉強する時間がない

の理由から諦めかけていた。
HSPでGPGPUプログラムができれば全て解決んだけど・・と思って何気なくHSPのサンプルを眺めていたら・・ これはもしかしたら、 HSPでGPGPUぽいことができるんじゃないのかと、そんな確信に近いある予感がひらめいた!!

きっかけはeasy3Dのサンプルのe3dhsp3_CbCr.hspを見た時だった
実行してみると、左には普通の3Dモデルのサンプルがリアルタイムで動いていて、右にはそのモデルの色だけがセピア色に変換されて同じく動いている

o0230017311356902500


 そうセピア変換とは、昔私がブログに無駄に長いソースを貼っつけたアレである。

まず画面の全ドットのrgb成分をyuvに逐一変換してuvに一定値を代入してまたrgbに変換・・・というやつである。
これをeasy3Dの中で実現しているということは、GPUにrgbをyuvに変換する並列処理をさせているということに他ならない!(CPUでリアルタイムなんかまず無理であろうし)

ということはソースかdllのプログラムのどこかにこのrgb成分をyuv成分に変換する、掛け、割り、足し、引きの計算式が書かれているはずだ!
そこのソースを改変すれば、自分の好きな手順、順番で四則演算ができるということではないのか!? この仮説が正しければeasy3Dのプラグインを使って高速並列計算が可能ということになる!
高速シミュレーションの道を開くために、執念でdllを解析しようかとか考えていたら、思いのほかもっと近くにその四則演算部分のプログラムを発見した
e3dhsp3.dllのと同じフォルダにある「E3D_HLSL」だ!
o0480024011361295228


そのフォルダにposteffect.fxが入っていて、こいつがセピア変換の四則演算を記述してあるソースだった。
調べたらこの.fxという第二のソースはHLSLといってHigh Level Shader Languageというもうひとつの言語らしい
おや、fxファイルをメモ帳で開くとプログラムソースが見れる!
ということは編集もできるではないか!?

o0730063411361295227



このファイルに自分で式を打って正常に動けば、自分の思い通りの並列処理が可能ってこと!?
試しに計算式を消したり追加したり実験してみたところ、コンパイルが正常に行なわれた

驚くべきことに、HSPで並列処理ができることがこれでわかった!!(正確にはHLSLでだけど・・・)


さらに驚くべきことに、easy3Dのソースを調べていくと律儀にも、私みたいなユーザーがposteffect.fxを改変できるようにと、なんとすでにカスタム領域が備わっていたのだ!
あとからおちゃっこさんのブログを見てわかったのだが「カスタムHLSL機能の需要は?」という名目で2年以上も前にこの自由に定義できるカスタムHLSLの機能を作っていたらしい。
つまりHSP側で並列計算させたいデータをeasy3Dの命令でGPUに送り、GPU側の処理の操作はHLSL(posteffect.fx)のカスタム領域で自分の好きな処理をさせられ、処理結果をまたHSP側で取り出すといった一連の操作が、我々ユーザーのためにすでに可能にされていたということだ。

HLSLのカスタムを可能にしてくれたおちゃっこさんには本当に感謝しないといけない!!
おちゃっこさんは需要なさそうと気にしていらっしゃったようだけど、これぞまさに私の求めていたもの!
この機能に出会えて本当に良かった!




HSPで並列演算ができるとわかった今、シミュレーションプログラムの幅が大きく広がること間違いなしである!!
これでもう数値解析やシミュレーションのデバッグに何日もかけないで済む!
HSPコンテストに糞ソフトと呼ばれるような生半可なシミュレーションソフトは出さないぞ。(呼んでいるのは主に自分)
さらにさらに驚くべきことに自分の非力なノーパにもそれなりのGPU・・オンボードだが・・がそなわっていたらしい!
o0399042211361313932

今の今まで気づかなかった・・そもそもGPUなんて付いてないものかと・・
最大単精度浮動小数点バッファまで対応してるみたいだし、シェーダモデル4.0で、なんだか考えられる中で最高の環境が気づかぬ間に手元にそろっていたようだ。
で、 具体的な並列処理の方法だが、自分の手順をさかのぼってみて、例をおりまぜてまとめてみた。





まずHSP3.2は入っている前提で

・easy3dのサイトからEasy3DForHSP3_ver5233をダウンロードする
・解凍してreadme読む ・readmeの言われたとおりにする
・次に解凍した中のE3D_HLSLのフォルダもhsp3フォルダにコピー
・下のリンクからダウンロードしてきたposteffect.fxを、hsp3フォルダのE3D_HLSLの中にあるposteffect.fxに上書き

いきなりHLSLを自分でカスタムするのは難しいので、他人のカスタムしたものを参考に少しづつ変えて勉強していくのが一番いいだろう。またHLSLのサイトにいけば命令とか演算子が全部見られる。
サンプル 上のサンプルは、左の2枚の画像を加算、減算、割り算、掛け算しその結果を右に出力するというものだ。1~8のキー操作を受け付けて8種類の処理が見れるようになっている。

o0777040011361295229

話を戻して、例えば足し算の並列演算をしたい場合なら、

・#include "e3dhsp3.as"して、E3DCreateRenderTargetTextureで、足す方を格納する2枚のテクスチャと、足した結果を出力するためのテクスチャをひとつ、合計3つ、同じ大きさで作る。
 ようはc=a+bという計算をさせたい場合、aというテクスチャとbというテクスチャとcというテクスチャを作ってください、ということだ。
 テクスチャの大きさは、並列に処理したい要素数によって変わってくる。65536要素の計算を並列に処理したいときは256*256の大きさのテクスチャをつくればいい。
・E3DShaderConstUserTexで足す2枚のテクスチャをGPUに送る
・E3DCallUserShaderでshadernoを0、passnoを1にして、足し算実行。これでE3DCallUserShaderの描画先のスワップチェインID(テクスチャ)に足し算結果が出力される。(ここらへんはサンプルの中にソースあり)





といった流れが、足し算の並列処理の方法だ。
passnoを変えれば割り算などを指定できる。
当然posteffect.fxを書き換えれば3つ4つ連続の足し算や、割った余りを求めるものなど、無限の処理が可能となる。
出力された計算結果(テクスチャ)はE3DShaderConstUserTexで指定しE3DCallUserShaderでつぎの処理に使いまわしていく事ができる。

上のスクリーンショットにも書いたが、E3DCreateRenderTargetTextureで作成するテクスチャのbit数を16bitや32bitにすれば、浮動小数点同士の足し算掛け算を並列処理できることとなる。

問題は、並列処理するときのデータが全てテクスチャ形式で行なわれているという点である。

はっきりいってHSP側で作ったDOUBLE型の配列変数を浮動小数点テクスチャに変換するというのはかなり至難の業である。(整数テクスチャならば話は簡単で、bmp形式でデータのやりとりを行なえばよい。)
というのもE3DCreateRenderTargetTextureで作ったテクスチャに1ドットずつ特定の値を書き込む命令がない。(というかeasy3Dはおそらくこの使いかたは想定していなかった)。だから何らかの方法で工夫しないといけないのだ。
それに、作れるテクスチャの精度が最大32bitなので標準命令のdoble型(64bit)と比べ精度が落ちる。グラフィックチップによっては16bitが限界だったりするから、CPUでのdoble型の計算をGPUで精度を落とさずやるには、16bitを4回計算みたいな感じでやるなどさらなる工夫が必要となる。

そして浮動小数点テクスチャから配列変数に戻すのも同様に大変である。
そこらへんを考慮した設計でないと、予想と結果が大きく異なってしまう可能性がある。


ということは、HSP側で配列変数を作る必要は本当にあるのかという疑問がわく。
並列演算出力結果を例えば視覚的に見るだけで良い場合(流体解析や分子動力学など)なら、posteffect.fx内で可視化プログラムを作ってしまえば、easy3Dのテクスチャレンダリング命令でスクリーンにいっきに可視化できる。わざわざHSP側で作った配列変数に値を戻さなくても、だ。
理想はそのような形で、最初から最後までテクスチャとして数値を保存しておく方が楽で良いのだが、厳密な数値を文字列として出力しなければいけないなど、どうしてもHSP側の変数に数値を戻す作業が欠かせない場合がでてくるだろう。 浮動小数点バッファを小数のまま保存できたりHSPのメモリに転送できる命令があればいいが、ない以上、1枚の浮動小数点バッファを何枚ものbmpに情報を分解して保存するアルゴリズムを構築するしかないのかもしれない。少なくとも私のお思いつく方法は今のとこそれしかない・・・
それだけがHSPでの並列処理においての大きな問題である。逆にもし単精度浮動小数点テクスチャと標準命令のdouble変数が両方向に対して変換可能になれば、まぁ64→32bitで精度が悪くなるのは仕方ないとして、タイトル通りの「HSPで高速並列演算」が可能になるというわけだ。これからは浮動小数点の構造の勉強をしてなんとか数値の取り込み、取り出しの問題を解決していきたいと思う・・・。








最後に、標準命令だけで作った流体力学プログラムを、HLSLに移植して並列処理化してみた。 計算は32bit浮動小数点バッファを使用、格子数は256*256で、CPUの80倍以上高速に処理できた。

o0512025611361324083

(一応解説すると、上から流体が流れてきてて、画面上の方にあるみにくいけど黒い動かないのが障害物。画面右は圧力を可視化した物。緑=高圧、赤=低圧)

これは、出力がスクリーンのみなので、posteffect.fxのほうに可視化プログラムが記述されており、可視化したテクスチャをeasy3DのE3DRenderSpriteでスクリーンにコピーしている、といったところ。
また入力のほうは、障害物の場所情報とか、0(なし) 1(あり) の2進数で表せるような情報しかないため、HSP側でbmpの生成→easy3Dでbmpの読み込み、bmpをスプライト(整数)に貼りつけ→32bit浮動小数点テクスチャにスプライト貼りつけ、・・と、数値の劣化なくGPUに送ることが出来ている。
並列処理を利用した数値解析シミュレーション特に流体力学はHSPコンテスト2011に向けて温めているもののひとつなのでネタバレはこのくらいにしておくとしよう。


次回:流体力学のプログラムを作りたい!その1













追記
HSP→浮動小数点バッファ への32bit浮動小数のデータ転送は
userFL4_0 ~ userFL4_9 のユーザー定義定数を使って正確に 行なえることが分かった。

また32bit浮動小数点を、4つの8bit*3(r,g,b)バッファに分解して転送して、HLSL側で数値を復元する方法も試みたが、
仮数部の最小1bit部分でたまに誤差が生まれる模様・・・


また浮動小数データの分割&整数値化による GPUバッファ→HSP へのデータ転送は全くうまくいかず・・・orz

パーリンノイズアルゴリズム 後編

以下の説明は引き続きバリューノイズの説明です。パーリンノイズは近日中にコード載せます

完成系
02c458d3.png

を作るには今まで説明に使ってきたような、テクスチャをいくつか使用して作る方法と、関数のみで作る方法の2種類がある。
関数のみでやる場合、作業用のテクスチャは1バイトも使わなくてすみ、並列処理化が容易であるためGPUでの処理もできる。
関数のみで作る方法は前の記事の外部リンク(リンク)で見ることができた。
自分でもHSPで関数のみのほうのアルゴリズムを作ってみたのだがどっかにやってしまった・・・
結局テクスチャを使うほうが圧倒的に早かったのだ

ちなみに前の前の記事で載せたgifたちは全て完成系を3D展開して、sin やらabsやら1/fの変換をほどこしたものである。
そして3Dのスライス面をgifのフレームに対応させてある

後編というわりに内容がほとんどないけどパーリンノイズの話は終わり
アルゴリズムのアの字も書けなかった・・ソース発見したらはるかも







追記、ソース発見!!

HSPスクリプトエディタに貼りつければ単独実行
(※処理に時間がかかります・・3分とか)

#module

#defcfunc Interpolate double a,double b,double c;aとbのあいだをcos補間で出す(内挿処理)
f = (1.0 - cos(c * 3.1415927)) * 0.5
return a*(1.0-f) + b*f

#defcfunc Noise int x,int y;ただのランダム数値を出す
z=y*9799+(y*y+435547)*(x-904120)*319-8129471
r=(z-32459)*(z+235152)+z+9022093
repeat (abs(z))\19+13+((r\11)*82)\12
r=(r+1114122123)*(r-108612193)+1182655124+r\12392874
loop
return 1.0*(r*(r*r*15731+789221)+1376312589)/2147483647.0

#defcfunc SmoothedNoise int x,int y;ガウシアンフィルタ3*3
return Noise(x, y) / 4+( Noise(x-1, y) +Noise(x+1, y) +Noise(x, y-1) +Noise(x, y+1) ) / 8+( Noise(x-1, y-1)+Noise(x+1, y-1)+Noise(x-1, y+1)+Noise(x+1, y+1) ) / 16

#defcfunc InterpolatedNoise double a,double b;加算合成する前の画像のx,yのドット色を出す
integer_X = int(a)
fractional_X = a - integer_X
integer_Y = int(b)
fractional_Y = b - integer_Y
v1 = SmoothedNoise(integer_X, integer_Y)
v2 = SmoothedNoise(integer_X + 1, integer_Y)
v3 = SmoothedNoise(integer_X, integer_Y + 1)
v4 = SmoothedNoise(integer_X + 1, integer_Y + 1)
i1 = Interpolate(v1 , v2 , fractional_X)
i2 = Interpolate(v3 , v4 , fractional_X)
return Interpolate(i1 , i2 , fractional_Y)

#defcfunc pnoise double a,double b;n枚の画像をブレンド率変えて加算合成
persistence=0.6
total = 0.0
p = persistence
n = Number_Of_Octaves@
frequency=1.0
amplitude=0.65
repeat n
total=total+(InterpolatedNoise(a*frequency,b*frequency)) *amplitude
frequency*=2.0
amplitude*=p
loop
return total

#global

;ここまで関数とかの定義

;スクリーンにレンダリング
Number_Of_Octaves=6
ddim k,640,480
screen 0,640,480,0

repeat 640
d=cnt
repeat 480
k.d.cnt=pnoise(0.01*d,0.01*cnt)
a=limit(218.0*(0.5+k.d.cnt),0,255)
color a,a,a
pset d,cnt
loop
loop
;bmpsave "a.bmp"



結果
f61c76fa.jpg 

プロフィール

toropippi

記事検索
アクセスカウンター

    QRコード
    QRコード
    • ライブドアブログ