適用例

GPUはSIMDマシンなので,CUDAの潜在力を利用するにはSIMDの方法で問題を提示しなければならない.各スレッドが1つの要素を独立して計算できる方法で分割できる計算はGPUに理想的である.

アルゴリズムによっては並列で書けなかったり,CUDA上で使えなかったり(アーキテクチャの制約のため)する.そのような場合に,それらの計算を行うためにGPUを使用する別の方法を導入するための研究が行われている.

このセクションでは,Mathematica 内におけるCUDAプログラミングの使い方を示す.以下のすべての例において,CUDAソース,バイナリまたはライブラリを Mathematica にロードすることができるようになるCUDAFunctionLoadが使われる.

CUDAFunctionLoadCUDA関数を Mathematica にロードする

CUDAFunctionLoadを使うと,CUDAのソース,バイナリ,ライブラリを Mathematica にロードすることが可能になる

画像処理

このセクションでは画像処理操作を行うCUDAの適用例を示す.CUDALink にはCUDAImageConvolveCUDABoxFilterCUDAErosionCUDADilationCUDAOpeningCUDAClosingCUDAImageAdd等,画像処理を行う関数がいくつか組み込まれている.

画像の二値化

Binarizeは入力画像を取り,画素が閾値より上の場合は白,それ以外は黒である二値画像を出力する.

CUDALink アプリケーションをまだインポートしていない場合はインポートする.

In[1]:=
Click for copyable input

入力画像を定義する.GPU上の必要なメモリ量を減らすために,画像を表すのにを使う.

In[2]:=
Click for copyable input

CUDAFunctionをロードする.型はパラメータリスト中でで表される.

In[3]:=
Click for copyable input
Out[3]=

入力画像を定義し,出力にメモリを割り当てる.

In[4]:=
Click for copyable input

閾値をとして二値化関数を呼び出す.

In[8]:=
Click for copyable input
Out[8]=

出力画像を表示する.

In[9]:=
Click for copyable input
Out[9]=

結果は Mathematica のものと一致する.

In[10]:=
Click for copyable input
Out[10]=

割り当てられたメモリをアンロードする.

In[11]:=
Click for copyable input

再帰的ガウシアン

再帰的ガウシアンは入力データにガウスフィルタを適用して近似する.これは速度を取ると精度を犠牲にする.以下はDericheスタイルの再帰的ガウスフィルタを実装する.

In[12]:=
Click for copyable input
Out[12]=

2つの関数,関数と関数をロードする.関数はRGBA画像を取り,データはまたはとして渡される.

In[13]:=
Click for copyable input

の値を指定する.ガウシアンはである.

In[15]:=
Click for copyable input

以下のパラメータはの値に基づいて計算される.

In[16]:=
Click for copyable input

次は入力画像を定義し,出力のためのメモリを割り当てる.

In[17]:=
Click for copyable input

関数を実行する.垂直ガウシアン,転置,もう一度垂直ガウシアン,もう一度転置を行って結果を得る.

In[21]:=
Click for copyable input

出力画像を得る.

In[25]:=
Click for copyable input
Out[25]=

割り当てられたメモリは解放しなければならない.

In[26]:=
Click for copyable input

ボックスフィルタ

ぼっくするフィルタはカーネルがBoxMatrixの時に最適化された畳込みである.ここでこれを実装する.

In[27]:=
Click for copyable input
Out[27]=

ボックスフィルタを適用するのに必要なCUDA関数をロードする.

In[28]:=
Click for copyable input

入力パラメータを設定する.

In[30]:=
Click for copyable input

半径はに設定される.

In[34]:=
Click for copyable input

関数を呼び出す.

In[35]:=
Click for copyable input

画像を得る.

In[37]:=
Click for copyable input
Out[37]=

割り当てられたメモリを解放する.

In[38]:=
Click for copyable input

画像の調整

これはCUDAにおけるImageAdjustの実装である.

In[39]:=
Click for copyable input

ソース文字列からCUDAFunctionをロードし,定数値にFloatを使用する.

In[40]:=
Click for copyable input
Out[40]=

CUDAFunctionをラップし,ImageAdjustと似たシンタックスで作る.

In[41]:=
Click for copyable input

関数はImageAdjustと同じように使える.

In[43]:=
Click for copyable input
Out[43]=

キャニー(Canny)エッジ検出

キャニーエッジ検出は10種類程度のフィルタを組み合せ,画像のエッジを見付ける.MathematicaEdgeDetectは同様の機能を提供する.以下がその実装である.

In[44]:=
Click for copyable input
Out[44]=

エッジ検出を行うのに必要なCUDA関数をロードする.

In[45]:=
Click for copyable input

入力画像とパラメータを設定する.

In[57]:=
Click for copyable input

ここでCUDAマネージャにホストとデバイスを加える.これらは入力と出力の両方を保持する.

In[60]:=
Click for copyable input

次に計算に使うためのデバイスのみの一時メモリを定義する.キャニーエッジ検出には多くのフィルタを使うので,多くのメモリが必要である.

In[62]:=
Click for copyable input

キャニーエッジ検出関数を呼び出す.

In[68]:=
Click for copyable input

出力を画像として表示する.

In[80]:=
Click for copyable input
Out[80]=

割り当てられたメモリをアンロードする.

In[81]:=
Click for copyable input

線形代数とリスト処理

このセクションでは線形代数操作を行うCUDA適用例を示す.ここで取り上げる関数の多くはCUDATransposeCUDADotを使って実行できる.

行列の転置

行列の転置は多くのアルゴリズムにおいて不可欠である.CUDALinkCUDATransposeという形ですぐに使える実装を提供する.しかし,ユーザは独自のものを実装することもできる.

行列の転置のためにCUDAFunctionをロードし,実数行列を取ってその転置を出力する新しい関数を定義する.

In[82]:=
Click for copyable input

行列を設定する.

In[83]:=
Click for copyable input
Out[83]//MatrixForm=

行列を転置する.

In[84]:=
Click for copyable input
Out[84]//MatrixForm=

Mathematica の結果と一致する.

In[85]:=
Click for copyable input
Out[85]//MatrixForm=

行列とベクトルの乗算

行列とベクトルの乗算は線形代数,有限要素解析等において一般的な操作である.CUDAFunctionをロードし,行列とベクトルの乗算を実装する.

In[86]:=
Click for copyable input

入力行列,ベクトルを設定する.

In[87]:=
Click for copyable input

上で定義した関数を呼び出し,MatrixFormを使って結果を表示する.

In[89]:=
Click for copyable input
Out[89]//MatrixForm=

Mathematica の結果と一致する.

In[90]:=
Click for copyable input
Out[90]//MatrixForm=

行列と行列の乗算

行列と行列の乗算は多くのアルゴリズムにおいて重要な関数である.ソースファイルからCUDAFunctionをロードし,ブロック次元をに設定する.

In[91]:=
Click for copyable input

入力値を設定し,出力のためのメモリを割り当てる.

In[93]:=
Click for copyable input

計算を行う.

In[96]:=
Click for copyable input

MatrixFormを使って結果を表示する.

In[97]:=
Click for copyable input
Out[97]//MatrixForm=

Mathematica の結果と一致する.

In[98]:=
Click for copyable input
Out[98]//MatrixForm=

出力バッファをアンロードする.

In[99]:=
Click for copyable input

ドット積

2つのベクトルのドット積は線形代数において一般的な操作である.ベクトルの集合を取り,それぞれのドット積を返す関数を実装する.

In[100]:=
Click for copyable input
Out[100]=

関数を Mathematica にロードする.

In[101]:=
Click for copyable input
Out[101]=

のランダムベクトルを生成する.

In[102]:=
Click for copyable input

関数を呼び出し,結果を返す.

In[105]:=
Click for copyable input
Out[105]=

結果は Mathematica の出力と一致する.

In[106]:=
Click for copyable input
Out[106]=

凸包

凸包はGPU上で並列化が難しいため,この例ではGPUプログラミングへのハイブリッドアプローチを取る.つまり,GPUに適した計算をGPU上で行い,残りをCPUで行う.

この凸包の実装はAndrewのアルゴリズムの教科書通りの実装で,簡単になるように設計されており,効率的ではない.

Mathematica 内で述語を定義する.

In[107]:=
Click for copyable input

端点を結ぶ与えられた線の上か下かどちらかに点集合を分割する関数を定義する.

In[108]:=
Click for copyable input

上記関数をCUDA関数としてロードする.

In[109]:=
Click for copyable input

2つの分割点集合を取り,それらの凸包を求める.

In[110]:=
Click for copyable input

上記関数を呼び出す,リストはが処理する前にCUDASortでソートされる.

In[111]:=
Click for copyable input

検証するために,一様に分布したランダムな20000の点を作成する.

In[112]:=
Click for copyable input

凸包を計算する.

In[114]:=
Click for copyable input

結果を可視化する.凸包の点の間に線が引かれる.

In[115]:=
Click for copyable input
Out[115]=

上記アルゴリズムはすべて,またはほとんどの点が凸包に含まれる極端な場合も扱う.単位円板じょうに一様に分布した点を生成する.

In[116]:=
Click for copyable input

凸包の点を計算する.

In[118]:=
Click for copyable input

凸包を可視化する.

In[119]:=
Click for copyable input
Out[119]=

これは直列コードでしか書けないようなアルゴリズムを,Mathematica とCUDAプログラミングを組み合せて部分的に並列にする重要な例である.

乱数生成

レイトレーシングから偏微分方程式の解放に至るまでの多くのアルゴリズムで,入力に乱数が必要となる.しかしこれは各コアが他のコアと同じ状態である複数のコアのシステムの多くでは難しい問題である.それを避けるために,通常は並列乱数生成器は一日の中の時間等のエントロピー値を乱数生成器のシードに利用するが,これらの呼出しはCUDAでは利用できない.

次のセクションでは乱数生成の3クラスのアルゴリズムを取り上げる.最初は一様乱数生成器である(擬似乱数生成器と準乱数生成器を例示する).2つ目は乱数生成にハッシング関数の一様分布を利用する乱数生成器,そして最後は正規乱数生成器である.

擬似乱数生成器

擬似乱数生成器はランダムであるように見える数を生成する決定論アルゴリズムである.これはシード値を利用して乱数を生成する.

このセクションでは,簡単な線形合同乱数生成器(Park-Millerのアルゴリズム)と,より複雑なメルセンヌツイスタ(Mersenne Twister)を取り上げる.

Park-Miller

Park-Miller乱数生成器は次の再帰方程式で定義される.

これは Mathematica で簡単に実装できる.

In[120]:=
Click for copyable input

に共通の値を使う.NestListを使って1000個の数のリストを生成し,プロットすることができる.

In[121]:=
Click for copyable input
Out[121]=

百万個の数の生成にかかる時間である.

In[122]:=
Click for copyable input
Out[122]=

代りにSeedRandomMethodオプションを使うこともできる.これも前回と同じように使える.

In[123]:=
Click for copyable input
Out[123]=

Mathematica での実装と比較すると,これは300倍も速い.

In[124]:=
Click for copyable input
Out[124]=

CUDAでの実装は Mathematica で書かれたものと似ている.実装は CUDALink で分配され,以下に場所を示す.

In[125]:=
Click for copyable input
Out[125]=

CUDA関数を Mathematica にロードする.

In[126]:=
Click for copyable input
Out[126]=

出力メモリを割り当てる.

In[127]:=
Click for copyable input
Out[127]=

CUDAFunctionを呼び出す.

In[128]:=
Click for copyable input
Out[128]=

結果はランダムである.

In[129]:=
Click for copyable input
Out[129]=

時間を測定すると, Mathematica の組込みメソッドより2倍,純粋な Mathematica 実装より600倍速いことが分かる.

In[130]:=
Click for copyable input
Out[130]=
In[131]:=
Click for copyable input
Out[131]=

Compileに対する計算時間もCUDAのものと同程度である.Compile文から整数オーバーフロー検出のないCコードを生成する.

In[132]:=
Click for copyable input

計算時間を求める.少し違いがある.

In[133]:=
Click for copyable input
Out[133]=

メモリの割当てに必要な時間を無視すると(生成された乱数はGPU上で再利用されるので,これは正当である),速度は10倍改善されたことが分かる.

In[134]:=
Click for copyable input
Out[136]=

メルセンヌツイスタ

メルセンヌツイスタは乱数の生成にシフトレジスタを利用する.実装は簡単なので,GPUにうまくマップされる.以下のファイルに実装が含まれている.

In[137]:=
Click for copyable input
Out[137]=

ファイルから関数をロードする.

In[138]:=
Click for copyable input

メルセンヌツイスタ入力パラメータを定義する.つい巣他派シード値を必要とする.これらの値はオフラインで計算してファイルに保管したり,Mathematica で生成したりすることができる.ここでは後者を示す.

In[139]:=
Click for copyable input

出力メモリを割り当てる.出力は上書きされるので,Mathematica からGPU上にメモリをロードする必要はない.

In[146]:=
Click for copyable input
Out[146]=

CUDAFunctionをパラメータを与えて呼び出す.

In[147]:=
Click for copyable input
Out[147]=

出力をプロットするとランダムであることが示される.

In[148]:=
Click for copyable input
Out[148]=

生成する乱数を入力として取り,パラメータを設定して必要な割当てを行い,ランダムな出力メモリを返す Mathematica 関数を書くことができる.

In[149]:=
Click for copyable input

最初の10000要素のプロットを生成する.

In[150]:=
Click for copyable input
Out[150]=

一億個の数を精製するのに乱数生成器が必要とする時間を測定する.

In[151]:=
Click for copyable input
Out[151]=

これは Mathematica の乱数生成器が必要とする時間と同程度である.

In[152]:=
Click for copyable input
Out[152]=

乱数が他の問題のシードであると考えると,Mathematica の所要時間がCUDA実装のものより優れていても,アルゴリズム全体としてパフォーマンスが改善される.

準乱数生成器

このセクションでは準乱数生成器について述べる.擬似乱数生成器と異なり,これらの列は一様ではなく,数値メソッドにおいて便利であることがある潜在的な構造を持つ.例えば,これらの列は通常多次元モンテカルロ(Monte Carlo)積分でより素早く収束する.

ハルトン(Halton)列

ハルトン列は単位区間で一様な準乱数を精製する.コードは任意の次元で動作するが,ここでは一次元空間でのみ動作するvan der Corput列について述べる.これは比較の意味で適切である.

ハルトン(またはvan der Corput)列は決定論的であるが,単位区間において相違が少ない.これは空間を一様に埋めるので,モンテカルロ積分など,適用分野によっては擬似乱数生成器より好ましいことがある.

指定の と基底

基底 における一次元ハルトン(またはvan der Corput)の値.

長さ の列は次のように書ける.

van der Corput列は,基底の表現に数が与えられると,数を小数点の反対側にミラーし,列の値はとなる.

Mathematica では,IntegerDigitsを使って列が求められる.

In[153]:=
Click for copyable input

基底をに設定すると,列の最初の1000要素が計算できる.

In[154]:=
Click for copyable input

結果をプロットすると,空間が一様に埋められることが分かる.

In[155]:=
Click for copyable input
Out[155]=

相違の少ない列の特性は,列の中の隣の要素は,前の要素がどこに位置するのかを知っているというものである.これはManipulateで見ることができる.

In[156]:=
Click for copyable input
Out[156]=

CUDA実装については,IntegerDigitsに相当する独自のものを実装しなければならないが,これは難しくはない.まず実装のソースコードをロードする.

In[157]:=
Click for copyable input
Out[157]=

CUDAFunctionをロードする.

In[158]:=
Click for copyable input
Out[158]=

出力用にメモリを割り当てる.乱数は1024個だけ生成する.

In[159]:=
Click for copyable input
Out[159]=

次元1で関数を実行する.

In[160]:=
Click for copyable input
Out[160]=

結果をプロットする.

In[161]:=
Click for copyable input
Out[161]=

Sobol列

Sobol列も相違の少ない列である.これは次のCUDAファイルに実装されている.

In[162]:=
Click for copyable input
Out[162]=

ブロック次元をとして関数をロードする.

In[163]:=
Click for copyable input
Out[163]=

入力パラメータをロードする.Sobol列が必要とする方向ベクトルは前もって計算され,ファイルに保管される.

In[164]:=
Click for copyable input

Sobol関数にパラメータを渡して実行する.

In[169]:=
Click for copyable input

列の中の最初の10000個の値をプロットする.空間が一様に埋められることが分かる(準乱数生成器).

In[170]:=
Click for copyable input
Out[170]=

完了したら,メモリはアンロードしなければならない.

In[171]:=
Click for copyable input

Niederreiter列

Niederreiter列もよく使われる,相違の少ない列である.これは以下のCUDAファイルに実装されている.

In[172]:=
Click for copyable input
Out[172]=

乱数ジェネレータで使うシードデータをロードする.

In[173]:=
Click for copyable input

CUDAFunctionをロードする.

In[174]:=
Click for copyable input
Out[174]=

計算に必要なメモリを割り当て,適切なパラメータを設定する.次元3の列を生成し,解像度はとする.

In[175]:=
Click for copyable input

関数を実行する.

In[181]:=
Click for copyable input

GPUから出力を取得する.

In[182]:=
Click for copyable input

列の最初の10000個の値をプロットする.Sobol列のときのように空間が一様に埋められているのが分かる.

In[183]:=
Click for copyable input
Out[183]=

完了したら,メモリはアンロードしなければならない.

In[184]:=
Click for copyable input

ハッシング乱数生成器

ハッシングを利用する乱数生成器は質のより低い乱数を生成するが,高速に生成する.これは多くの適用分野において十分である.

小型暗号化アルゴリズムハッシング

小型暗号化アルゴリズム(TEA)は非常に簡単なハッシングアルゴリズムで,以下のファイルに実装されている.

In[185]:=
Click for copyable input
Out[185]=

CUDAFunctionをロードする.

In[186]:=
Click for copyable input
Out[186]=

出力用のメモリを割り当てる.

In[187]:=
Click for copyable input
Out[187]=

CUDAFunction関数を呼び出す.

In[188]:=
Click for copyable input
Out[188]=

結果をプロットする.

In[189]:=
Click for copyable input
Out[189]=

割り当てられたメモリを削除する.

In[190]:=
Click for copyable input

MD5ハッシング

乱数生成器にはその他の一般的なハッシング法が使える.以下はよく知られたハッシングアルゴリズムである,MD5アルゴリズムの実装である.

In[191]:=
Click for copyable input
Out[191]=

ソースからCUDAFunctionをロードする.

In[192]:=
Click for copyable input
Out[192]=

出力メモリをロードする.

In[193]:=
Click for copyable input
Out[193]=

CUDAFunctionを呼び出す.

In[194]:=
Click for copyable input
Out[194]=

結果をプロットする.

In[195]:=
Click for copyable input
Out[195]=

割り当てられたメモリを削除する.

In[196]:=
Click for copyable input

正規乱数

以下のアルゴリズムは正規分布の乱数を生成する.

逆累積正規分布

以下は正規分布の乱数を生成する方法を実装する.

In[197]:=
Click for copyable input
Out[197]=

CUDAFunctionをロードする.

In[198]:=
Click for copyable input
Out[198]=

100000個の乱数のためのメモリを割り当てる.

In[199]:=
Click for copyable input

CUDAFunctionを呼び出す.

In[201]:=
Click for copyable input
Out[201]=

メモリを Mathematica 内に取得する.

In[202]:=
Click for copyable input

Histogramを使って結果をプロットする.

In[203]:=
Click for copyable input
Out[203]=

メモリをアンロードする.

In[204]:=
Click for copyable input

Box-Muller

Box-Mullerは一式の一様分布の乱数を与えられると,正規分布の数を生成する方法である.CUDAの実装は以下のファイルにある.

In[205]:=
Click for copyable input
Out[205]=

CUDAFunctionをロードする.

In[206]:=
Click for copyable input
Out[206]=

入力引数を設定する.

In[207]:=
Click for copyable input

メルセンヌツイスタ(2セクション前に定義したもの)を使って一様に分布した乱数を生成する.

In[212]:=
Click for copyable input
Out[212]=

一様な乱数のリストを正規分布乱数に変換する.

In[213]:=
Click for copyable input
Out[213]=

Histogramを使うとベル状の曲線が示される.

In[214]:=
Click for copyable input
Out[214]=

T割り当てられたメモリを削除する.

In[226]:=
Click for copyable input

乱数ジェネレータの適用例

乱数は多くの分野に適用できる.ここでは,モンテカルロ積分( と任意の関数を近似して)とブラウン(Brownian)運動のシミュレーションという2つの重要な適用例を示す.

の近似

の値はモンテカルロ積分を使って近似できる.まず単位正方形で一様乱数を生成する.次に単位円の最初の四分円内の点の数を数える.それから結果を点の数で割る.これがとなる.

単位円内の点の数を数え,リダクションする実装である.

In[216]:=
Click for copyable input
Out[216]=

CUDAFunctionをロードする.

In[217]:=
Click for copyable input
Out[217]=

1000000個の点を使う.

In[218]:=
Click for copyable input

前述のメルセンヌツイスタアルゴリズムを使って乱数を生成する.

In[221]:=
Click for copyable input
Out[221]=

出力メモリを割り当てる.

In[222]:=
Click for copyable input
Out[222]=

計算を行う.

In[223]:=
Click for copyable input
Out[223]=

出力メモリを取得する.

In[224]:=
Click for copyable input
Out[224]=

結果は Mathematica のものと一致する.

In[225]:=
Click for copyable input
Out[225]=

計算速度は格段に速い.

In[226]:=
Click for copyable input
Out[226]=

Mathematica と比較する.

In[227]:=
Click for copyable input
Out[227]=

モンテカルロ積分

モンテカルロ積分は多くの分野で使われる.ここではSqrt[x]を0から1まで積分する.

In[228]:=
Click for copyable input

関数をロードする.

In[229]:=
Click for copyable input
Out[229]=

乱数にSobel準乱数生成器を使う.

In[230]:=
Click for copyable input

長さに乱数生成器の数を使う.

In[238]:=
Click for copyable input
Out[238]=

出力用にメモリを割り当てる.

In[239]:=
Click for copyable input
Out[239]=

関数を呼び出す.

In[240]:=
Click for copyable input
Out[240]=

最初のいくつかの要素が適切であるかどうか確認する.

In[241]:=
Click for copyable input
Out[241]=

出力を足し合わせる.これはCUDAFoldを使って行う.

In[242]:=
Click for copyable input
Out[242]=

結果はNIntegrateのものと一致する.

In[243]:=
Click for copyable input
Out[243]=

割り当てられたメモリをアンロードする.

In[244]:=
Click for copyable input

ブラウン運動

シミュレーションのためのメモリを割り当てる.

In[245]:=
Click for copyable input
Out[245]=
In[246]:=
Click for copyable input
In[247]:=
Click for copyable input

擬似乱数列の値は正規分布している.

Out[248]=
In[249]:=
Click for copyable input
In[250]:=
Click for copyable input
Out[250]=

コードの生成

CUDALinkMathematica に統合されているので,CUDAカーネルコードの生成にSymbolicC等の Mathematica 機能を利用することができる.CUDALink をまだインポートしていない場合はインポートする.

In[251]:=
Click for copyable input

個の例ではSymbolicCパッケージが必要である.

In[252]:=
Click for copyable input

いくつかの一般的な Mathematica 構文をSymbolicC表現に変換して定義する.

In[253]:=
Click for copyable input

検証するために Mathematica 文を渡すと,SymbolicC出力が得られる.

In[271]:=
Click for copyable input
Out[271]=

C文字列に変換するにはToCCodeStringメソッドを使う.

In[272]:=
Click for copyable input
Out[272]=

以上により,Mathematica 関数(純粋であってもそうでなくても)を取り,適切なCUDAカーネルソースを生成する関数が書けるようになる.

In[273]:=
Click for copyable input

純粋関数をに渡すと,カーネルコードが返される.

In[276]:=
Click for copyable input
Out[276]=

Mathematica 関数と入力リストを与えられると,CUDAカーネルコードを生成し,コードをCUDAFunctionとしてロードし,CUDAFunctionを実行し,結果を返す関数を定義する.

In[277]:=
Click for copyable input

入力リストの各要素に2を足す純粋関数でを検証する.

In[279]:=
Click for copyable input
Out[279]=

が変換したあらゆる構文がによってサポートされる.以下では各要素を2乗する.

In[280]:=
Click for copyable input
Out[280]=

モンテカルロ積分を行う.

In[281]:=
Click for copyable input
Out[281]=

関数は定義してに渡すことができる.次は色反転関数を定義する.

In[282]:=
Click for copyable input

色反転関数を呼び出す.

In[283]:=
Click for copyable input
Out[283]=

以上は簡単な例であるが,より複雑なものの基盤として使える.

New to Mathematica? Find your learning path »
Have a question? Ask support »