0. はじめに
- もともと所有するFortranコードをGPGPU移植するにあたって、Cuda Fortranが最も移植コストが少ないと考えた。
- Monte Carlo法を有するコードなのでcuRANDを使いたいと思って調べたところ、Cuda Fortranのライブラリから直接読み込めるようであった。
- 乱数長が事前にわからないので、curand_init/curand_uniformを使ってカーネル内から直接cuRANDを使用したかった。
- しかしその解説記事が存在しておらず、ChatGPT (2024年5月15日時点のGPT-4o)を用いても、情報古いのか、CUDA Cで書いたコードをCUDA Fortranカーネルから呼び出すように言われた。
- 最終的に実装できたのでその方法とTipsを示す。
環境
OS： Ubuntu 22.4
コンパイラ： nvfortran 2024 / 24.3
GPU: RTX 1660 Super
1. 実装
はじめにchatGPTに聞いて出てきたコードが以下
質問：
CUDA Fortranで各スレッドが独自の乱数を生成するためには、cuRANDライブラリを使用します。以下に、10個の乱数を別々のスレッドで生成するためのテストコードを示します。
module curand_module
use cudafor
implicit none
include 'curand_kernel.h'
contains
attributes(global) subroutine setup_kernel(state, seed)
type(curandState), device :: state(*)
integer, value :: seed
integer :: id
id = threadIdx%x + (blockIdx%x - 1) * blockDim%x
call curand_init(seed, id, 0, state(id))
end subroutine setup_kernel
attributes(global) subroutine generate_kernel(state, result)
type(curandState), device :: state(*)
real(kind=4), device :: result(*)
integer :: id
id = threadIdx%x + (blockIdx%x - 1) * blockDim%x
result(id) = curand_uniform(state(id))
end subroutine generate_kernel
end module curand_module
program main
use cudafor
use curand_module
implicit none
integer, parameter :: n = 10
integer, parameter :: seed = 1234
type(curandState), device, allocatable :: dev_states(:)
real(kind=4), device, allocatable :: dev_random_numbers(:)
real(kind=4), allocatable :: random_numbers(:)
! デバイスメモリのアロケート
allocate(dev_states(n))
allocate(dev_random_numbers(n))
allocate(random_numbers(n))
! 乱数状態の初期化
call setup_kernel<<<1, n>>>(dev_states, seed)
! 乱数生成
call generate_kernel<<<1, n>>>(dev_states, dev_random_numbers)
! 結果をホストにコピー
random_numbers = dev_random_numbers
! 結果を表示
print *, 'Generated random numbers:'
print *, random_numbers
! メモリの解放
deallocate(dev_states)
deallocate(dev_random_numbers)
deallocate(random_numbers)
end program main
include 'curand_kernel.h'
となっているし、適切に出力できていない。
これでできないですと返したら、
1.CUDA Cコードで乱数生成関数を実装
2.FortranモジュールでCUDA C関数のインターフェースを定義
3.Fortranメインプログラムで乱数生成を呼び出す
を提案してきた。いくつか他の聞き方をしてもイマイチ理解していないので、自分で調べることにした。
コード改修
上記のコードをベースに動かすために重要だったポイントを記す。
だいたいここ (https://docs.nvidia.com/hpc-sdk/compilers/fortran-cuda-interfaces/index.html#cfrand-runtime) のサンプルコードを見たらわかった。
type(curandState)は使えない
代わりに、type(curandStateXORWOW) が使用できる。
curand_deviceをインクルードする必要がある
use curand_deviceを最初に入れる。
コンパイルオプションに-gpu=nonvvmを入れる必要がある
旧情報だと-Mcuda=nonvvmだが、最新版コンパイラなら-gpu=nonvvmと記載する。
2. 最終版コード
上記以外にもいくつか修正を行いながら最終的に問題なく動くと思われる状態になったコードを以下に示す。
モジュール
curand.cuf
module curand_module
use cudafor
! use curand
use curand_device
implicit none
contains
attributes(global) subroutine setup_kernel(n, state, seed)
integer, value :: n
type(curandStateXORWOW), device :: state(*)
integer, value :: seed
integer :: id
id = threadIdx%x + (blockIdx%x - 1)*blockDim%x
if (id <= n) then
call curand_init(seed, id, 0, state(id))
end if
end subroutine setup_kernel
attributes(global) subroutine generate_kernel(n, state, result)
integer, value :: n
type(curandStateXORWOW), device :: state(*)
double precision, device :: result(*)
integer :: id
id = threadIdx%x + (blockIdx%x - 1)*blockDim%x
if (id <= n) then
result(id) = curand_uniform(state(id))
end if
end subroutine generate_kernel
end module curand_module
メインプログラム
main.f90
program main
use curand_module
implicit none
integer, parameter :: n = 1000000
integer :: seed
type(curandStateXORWOW), device, allocatable :: dev_states(:)
double precision, device, allocatable :: dev_random_numbers(:)
double precision, allocatable :: random_numbers(:)
integer :: values(8) ! to hold date_and_time values
double precision :: start_time, end_time, elapsed_time
double precision :: mean, stddev, sumrand, sumsq
integer :: i, n_thread, n_block
integer(kind=c_size_t) :: mem_size
! 現在の時間を取得してシード値を生成
call date_and_time(values=values)
seed = values(5)*10000 + values(6)*100 + values(7) ! 時・分・秒を組み合わせてシード値を生成
print *, 'Seed value: ', seed
! デバイスメモリのアロケート
allocate (dev_states(n))
allocate (dev_random_numbers(n))
allocate (random_numbers(n))
mem_size = size(dev_states) * c_sizeof(dev_states(1))
print *, 'Memory size of dev_states: ', mem_size, ' bytes'
mem_size = size(dev_random_numbers) * c_sizeof(dev_random_numbers(1))
print *, 'Memory size of dev_random_numbers: ', mem_size, ' bytes'
! 乱数生成開始時刻の記録
call cpu_time(start_time)
n_thread = 256
n_block = (n + n_thread - 1) / n_thread
! 乱数状態の初期化
call setup_kernel <<<n_block, n_thread>>> (n, dev_states, seed)
! 乱数生成
call generate_kernel <<<n_block, n_thread>>> (n, dev_states, dev_random_numbers)
! 乱数生成終了時刻の記録
call cpu_time(end_time)
! 結果をホストにコピー
random_numbers = dev_random_numbers
! 経過時間を計算
elapsed_time = end_time - start_time
! 乱数の平均値と標準偏差を計算
sumrand = 0.0d0
sumsq = 0.0d0
do i = 1, n
sumrand = sumrand + random_numbers(i)
sumsq = sumsq + random_numbers(i)**2
end do
mean = sumrand/n
stddev = sqrt((sumsq/n) - mean**2)
! 結果を表示
print *, 'Generated random numbers (first 10):'
print *, random_numbers(1:10)
print *, 'Elapsed time (seconds): ', elapsed_time
print *, 'Mean: ', mean
print *, 'Standard Deviation: ', stddev
! メモリの解放
deallocate (dev_states)
deallocate (dev_random_numbers)
deallocate (random_numbers)
end program main
なお上記で配列のメモリを確認したり平均・標準偏差をとっているのは、乱数が適切に生成されているかたしかめるため。
コンパイルコマンド
nvfortran -cuda -gpu=nonvvm curand.cuf main.cuf -o main
3. 実行結果
Seed value: 142021
Memory size of dev_states: 48000000 bytes
Memory size of dev_random_numbers: 8000000 bytes
Generated random numbers (first 10):
0.1877072751522064 0.2312117666006088 0.6245507597923279
0.6655896306037903 0.7520658969879150 0.4479118287563324
0.8056853413581848 0.8860529661178589 0.4329888820648193
0.4152245223522186
Elapsed time (seconds): 4.4918060302734375E-004
Mean: 0.4998809951801293
Standard Deviation: 0.2885744465300747