cuFFTパラメータ編 #C++

概要

cuFFTで主に使用するパラメータの紹介

はじめに

最初に言います。
「cuFFTまじでむずい！！」
少し扱う機会があったので、勉強をしてみたのですが最初使い方が本当にわかりませんでした。
今もわからない部分はありますが、忘れたくないので勉強したことを記事にしていこうと思います。
間違っている部分があれば、遠慮なく教えていただけると幸いです。

cufftComplexに画素値を代入してみる

cuFFTとは

まずcuFFTとは何かと言いますと
「GPU上で離散フーリエ変換を行うCUDAのライブラリ」
です。
CUDAのドキュメントにもIntroductionのところに説明が書いてあります。

英語：
The cuFFT library provides a simple interface for computing FFTs on an NVIDIA GPU, which allows users to quickly leverage the floating-point power and parallelism of the GPU in a highly optimized and tested FFT library.

日本語(DeepL)：
cuFFTライブラリは、NVIDIA GPU上でFFTを計算するためのシンプルなインターフェースを提供し、高度に最適化されテストされたFFTライブラリでGPUの浮動小数点演算能力と並列性を迅速に活用することを可能にします。

cuFFTドキュメント

cuFFTで主に使うであろうもの

cuFFTの大まかな流れとしては

cuFFT planの保存とアクセスに使用するハンドルを定義
planの詳細なパラメータを定義
FFTを実行

みたいな感じです。途中のホストからデバイスにデータをコピーするといった処理は省略していますが、ざっくり書くと上記のようになります。
少し補足をすると、「plan」とは「CUFFTプランの保存とアクセスに使用されるハンドル型」です。わかりやすく言い換えると、フーリエ変換をするときにこのplanを介して行うみたいな感じです。

plan関数

cufftPlan1D

1次元信号をフーリエ変換するときに定義する「plan」で、パラメータとしてはこんな感じ。

.cpp

cufftPlan1d(cufftHandle *plan, int nx, cufftType type, int batch);

Input
plan	ハンドル
nx	データのサイズ（信号の長さ）
type	フーリエ変換の方法を選択（下記参照）
batch	Inputのフーリエ変換計算回数：基本はnxと同じ値

cufftPlan2D

2次元のデータをフーリエ変換するときに定義する「plan」で、パラメータとしてはこんな感じ。
画像のフーリエ変換などに使用するはず。

.cpp

cufftPlan2d(cufftHandle *plan, int nx, int ny, cufftType type);

Input
plan	ハンドル
nx	データのサイズ（データの大きさ：画像だとwidth）
ny	データのサイズ（データの大きさ：画像だとheight）
type	フーリエ変換の方法を選択（下記参照）

※cufftPlan2dにはbatchのオプションはないので、もしbatchを使いたければcufftPlanManyを使用すること

cufftPlan3d

3次元データをフーリエ変換するときに定義する「plan」で、パラメータとしてはこんな感じ。

.cpp

cufftPlan3d(cufftHandle *plan, int nx, int ny, int nz, cufftType type);

Input
plan	ハンドル
nx	xのデータのサイズ
ny	yのデータのサイズ
nz	zのデータのサイズ
type	フーリエ変換の方法を選択（下記参照）

cufftPlanMany

これはわかりにくいので詳しく説明します。まずパラメータがこんな感じです。

.cpp

cufftPlanMany(cufftHandle *plan, int rank, int *n, int *inembed,
        int istride, int idist, int *onembed, int ostride,
        int odist, cufftType type, int batch);

Input
plan	ハンドル
rank	入力データの次元（1次元～3次元）
n	入力データの大きさ・フーリエ変換のカーネルサイズ（おそらくこの解釈で合ってるはず）
inembed	入力するpixelの幅（画像であればwidth） or NULL
istride	入力するブロックの間隔 or NULL
idist	入力するブロックの大きさ or NULL
onembed	出力するpixelの幅（画像であればwidth） or NULL
ostride	出力するブロックの間隔 or NULL
odist	出力するブロックの大きさ or NULL
type	フーリエ変換の方法を選択（下記参照）
batch	Inputのフーリエ変換計算回数

とパラメータがたくさんあります。cufftPlanManyはAdvanced Data Layoutというもので、ユーザーが自由にデータの入出力を調整することができます。例えば、640×300の画像データを7×7に分割して入力したり、640×300の画像データを7×7の分割データで出力したりすることができます。
簡潔に言えば、cufftPlan1dとかの設定をより細かくできるみたいな感じです。
それでは、それぞれのパラメータの説明に入ります。

rank

これは単純に、入力データが何次元かを指定するものです。1次元～3次元で選択します。

n

このnはわかりにくいです。まずは公式ドキュメントの内容を見てみます。

英語：
Array of size rank, describing the size of each dimension, n[0] being the size of the outermost and n[rank-1] innermost (contiguous) dimension of a transform.

日本語（DeepL）：
大きさ rank の配列で，各次元の大きさを表す．n[0] は変換の一番外側の次元， n[rank-1] は一番内側の（連続した）次元の大きさである．

「うーん。。。つまりどういうことだ？」
最初見たときに↑みたいな感じになりました。
このnが意味するのはフーリエ変換に使用するkernel サイズを意味しているみたいです。
例えば、2次元データである画像に対してフーリエ変換を行うときに、

.cpp

int rank = 2;
int n[rank] = {7, 7}; // nの配列サイズはrank

と設定してあげれば、7×7のカーネルでフーリエ変換ができるということです。

inembed

これはフーリエ変換を行うブロックのPixel幅を設定します。
公式ドキュメントでは以下のように説明されています。

英語：
The inembed and onembed parameters define the number of elements in each dimension in the input array and the output array respectively. The inembed[rank-1] contains the number of elements in the least significant (innermost) dimension of the input data excluding the istride elements; the number of total elements in the least significant dimension of the input array is then istride*inembed[rank-1]. The inembed[0] or onembed[0] corresponds to the most significant (that is, the outermost) dimension and is effectively ignored since the idist or odist parameter provides this information instead. Note that the size of each dimension of the transform should be less than or equal to the inembed and onembed values for the corresponding dimension, that is n[i] ≤ inembed[i], n[i] ≤ onembed[i]，where i∈{0,…,rank−1}

日本語（DeepL）：
inembedとonembedパラメータは，それぞれ入力配列と出力配列の各次元の要素数を定義します．inembed[rank-1] には，入力データの最下位（最内周）次元の要素数（istride 要素を除く）が入ります．したがって，入力配列の最下位次元の総要素数は istride*inembed[rank-1] となります．inembed[0] または onembed[0] は，最も大きい（つまり一番外側の）次元に対応し， idist または odist パラメータが代わりにこの情報を提供するので，実質的には無視されます．変換の各次元のサイズは，対応する次元の inembed と onembed の値以下でなければならないことに注意してください．

ドキュメントではinembedは入力配列の要素数を定義すると記載されています。
またまた画像を例えに用いますと、inembedは入力する画像の幅を設定するということになります。
ここで、公式ドキュメントに記載してある式を以下に示します。

1D：
input[ b * idist + x * istride]
output[ b * odist + x * ostride]

2D：
input[b * idist + (x * inembed[1] + y) * istride]
output[b * odist + (x * onembed[1] + y) * ostride]

3D：
input[b * idist + ((x * inembed[1] + y) * inembed[2] + z) * istride]
output[b * odist + ((x * onembed[1] + y) * onembed[2] + z) * ostride]

上記の式はInput、Outputされる値を表したものになります。
この式を見るとお気づきの方もいるかもしれませんが、inembed[0]の値は無視されています。なので、inembedで意味のある値は**inembed[1]**以上の値になります。
ただ、that is n[i] ≤ inembed[i], n[i] ≤ onembed[i]という記載があるのが少し気にかかりはします。inembed[0]は意味を持たないのに、なぜ条件に含まれているのかが不思議です。

inembedは下記のように定義します。

上記の式の意味合いなどは別の記事を記載いたしますので少々お待ちを。。。（長くなってしまいそうなので）

.cpp

int inembed[rank] = {1, width};

nとinembedを2Dで使う際のイメージ図を作ってみました。

idist

これは、入力する要素の数を定義します。1個ずつでよければ1を設定し、7×7で入力したければ49を設定します。ただ、公式ドキュメントでは

idist または odist パラメータが代わりにこの情報を提供するので，実質的には無視されます

と書かれているため、重要なのはodistの定義の方になります。

istride

このistrideは深層学習とかで畳込み層にて使用するstrideとは意味合いが少し違います。
これは図を見てもらえばわかりやすいかなと思うので先に図を載せます。1の場合と2の場合を載せています。
図の補足として、四角内の数字は入力順を表しています。

言葉で表すと、入力する形を一つのブロック(kernel×kernel)とし、そのブロックの入力間隔をistrideで調整するということになるのかなと思います。
図では

istride=1であれば、入力はブロック順の通り
istride=2であれば、入力は一個飛ばし

のようになります。
それでは次に出力に使用するパラメータの説明に入っていきましょう。

onembed

onembedはinembedとそんなに違いはありません。出力のwidthを定義する感じです。inembedと同様に0番目の値は無視されるので、2Dで使用すると

.cpp

int onembed[rank] = {1, width}

のようになります。

odist

odistは出力する要素の数を定義します。

.cpp

int odist = 1; // 要素1個を1ブロックとする場合

int odist = 49; // 要素49個を1ブロックとする場合（7×7）

ostride

これもistrideとほぼ同じ役割です。出力するブロックの間隔を定義します。

batch

フーリエ変換の計算回数を定義します。正確には並列処理で実行する個数の定義を行います。batchはGPUのスペックが許す限りは大きくしていいのかなと思います。

英語：
where batch denotes the number of transforms that will be executed in parallel, rank is the number of dimensions of the input data

日本語（DeepL）：
ここで、batch は並列に実行される変換の数を表し、rank は入力データの次元数を表す

typeについて

planを作成する際にはフーリエ変換の方法を選択する必要があります。

type
CUFFT_R2C	Float型の実数からFloat型の複素数に変換
CUFFT_C2R	Float型の複素数からFloat型の実数に変換
CUFFT_C2C	Float型の複素数からFloat型の複素数に変換
CUFFT_D2Z	Double型の実数からDouble型の複素数に変換
CUFFT_Z2D	Double型の複素数からDouble型の実数に変換
CUFFT_Z2Z	Double型の複素数からDouble型の複素数に変換

上記の中から実装環境に合わせて選択します。

Plan作成のエラーチェック

上記で紹介してきたplanには戻り値があり、それでplanの作成に成功したか失敗したか知ることができます。
まず戻り値の内容としては公式ドキュメントより以下のようになってます。

typedef enum cufftResult_t {
    CUFFT_SUCCESS        = 0,  //  The cuFFT operation was successful
    CUFFT_INVALID_PLAN   = 1,  //  cuFFT was passed an invalid plan handle
    CUFFT_ALLOC_FAILED   = 2,  //  cuFFT failed to allocate GPU or CPU memory
    CUFFT_INVALID_TYPE   = 3,  //  No longer used
    CUFFT_INVALID_VALUE  = 4,  //  User specified an invalid pointer or parameter
    CUFFT_INTERNAL_ERROR = 5,  //  Driver or internal cuFFT library error
    CUFFT_EXEC_FAILED    = 6,  //  Failed to execute an FFT on the GPU
    CUFFT_SETUP_FAILED   = 7,  //  The cuFFT library failed to initialize
    CUFFT_INVALID_SIZE   = 8,  //  User specified an invalid transform size
    CUFFT_UNALIGNED_DATA = 9,  //  No longer used
    CUFFT_INCOMPLETE_PARAMETER_LIST = 10, //  Missing parameters in call
    CUFFT_INVALID_DEVICE = 11, //  Execution of a plan was on different GPU than plan creation
    CUFFT_PARSE_ERROR    = 12, //  Internal plan database error 
    CUFFT_NO_WORKSPACE   = 13  //  No workspace has been provided prior to plan execution
    CUFFT_NOT_IMPLEMENTED = 14, // Function does not implement functionality for parameters given.
    CUFFT_LICENSE_ERROR  = 15, // Used in previous versions.
    CUFFT_NOT_SUPPORTED  = 16  // Operation is not supported for parameters given.
} cufftResult;

肝心のどうやって値が返ってくるかわかりにくい。。笑
取得できる値としては、0,1,2,...の数字が取得できます。この数字なんですが、cufftResultという型なので

.cpp

cufftResult cufftR = cufftPlanMany(&plan, rank, n, inembed, istride, idist, onembed, ostride, odist, CUFFT_C2C, batch);

のようする必要があります。あとは、Swich文やif文でリザルトの内容を表示するものを作ってあげればわかりやすくなるかなと思います。
もちろんCUFFT_SUCCESS であればplanの作成は成功してます。

フーリエ変換関数

最後にフーリエ変換を実際に行う関数の説明を行っていきます。この関数もplanと同じ方法で戻り値の確認をすることができます。

cufftExecR2C

実数から虚数に変換するフーリエ変換です。

.cpp

cufftExecR2C(cufftHandle plan, cufftReal *idata, cufftComplex *odata);

Input
plan	ハンドル
idata	入力データ
odata	出力データ

出力データのところにidataを指定するとインプレイス変換され、idataに出力データが保存されます。

cufftExecC2R

虚数から実数に変換するフーリエ変換です。R2Cでフーリエ変換した後に、C2Rで逆フーリエ変換するパターンが多いのかなと思います。

.cpp

cufftExecC2R(cufftHandle plan, cufftComplex *idata, cufftReal *odata);

Input
plan	ハンドル
idata	入力データ
odata	出力データ

出力データのところにidataを指定するとインプレイス変換され、idataに出力データが保存されます。

cufftExecC2C

複素数から複素数に変換するフーリエ変換です。

.cpp

cufftExecC2C(cufftHandle plan, cufftComplex *idata, cufftComplex *odata, int direction);

Input
plan	ハンドル
idata	入力データ
odata	出力データ
direction	CUFFT_FORWARD or CUFFT_INVERSEのどちらかを指定する（FORWARDがフーリエ変換、INVERSEが逆フーリエ変換）

出力データのところにidataを指定するとインプレイス変換され、idataに出力データが保存されます。
cufftComplexの型は、xに実数部、yに虚数部が設定されています。↓非常に簡略化してるけどこんな感じで値は参照可能。

.cpp

cufftComplex *idata;
x = idata[10].x;
y = idata[10].y;

あとはダブル型の関数もありますが、全く一緒なのでリンクだけ貼っておきます。

cuFFT Execution

終わりに

コードの実装も書くと記事が長くなるので一旦ここでこの記事は締めます。
これからもこれ必要だなと感じたら説明を加えて行きます。

cuFFTパラメータ編

概要

はじめに

関連記事

cuFFTとは

cuFFTで主に使うであろうもの

plan関数

cufftPlan1D

cufftPlan2D

cufftPlan3d

cufftPlanMany

rank

n

inembed

idist

istride

onembed

odist

ostride

batch

typeについて

Plan作成のエラーチェック

フーリエ変換関数

cufftExecR2C

cufftExecC2R

cufftExecC2C

終わりに