0. はじめに
もともと所有するFortranコードをGPGPU移植するにあたって、Cuda Fortranが最も移植コストが少ないと考えた。キー配列を用いた配列の並べ替えやReduction(総和)の計算が必要になる計算コードで、CUDA CのThrustで実装できることはわかっていたが、Cuda FortranでThrustはまだサポートされていない。
CUDA Cで書いたコードをCUDA Fortranカーネルから呼び出すことで実装を行った。
環境
OS: Ubuntu 22.4
コンパイラ: nvfortran 2024 / 24.3
GPU: RTX 1660 Super
1. 実装
はじめにchatGPT(2024年5月15日時点のGPT-4o)に聞いて出てきたコードが以下
質問:
thrust::stable_sort_by_keyをCudaFortranで使うための方法を教えてください。なお、コンパイラはNVIDIA SDKを使用します。
#include <thrust/sort.h>
#include <thrust/device_vector.h>
extern "C" void sort_by_key(int* keys, int* values, int size) {
thrust::device_vector<int> d_keys(keys, keys + size);
thrust::device_vector<int> d_values(values, values + size);
thrust::stable_sort_by_key(d_keys.begin(), d_keys.end(), d_values.begin());
thrust::copy(d_keys.begin(), d_keys.end(), keys);
thrust::copy(d_values.begin(), d_values.end(), values);
}
module thrust_sort
use iso_c_binding
implicit none
interface
subroutine sort_by_key(keys, values, size) bind(C, name="sort_by_key")
use iso_c_binding
integer(c_int), dimension(*), intent(inout) :: keys
integer(c_int), dimension(*), intent(inout) :: values
integer(c_int), intent(in) :: size
end subroutine sort_by_key
end interface
contains
subroutine sort_keys_and_values(keys, values, size)
use iso_c_binding
integer(c_int), dimension(*), intent(inout) :: keys
integer(c_int), dimension(*), intent(inout) :: values
integer(c_int), intent(in) :: size
call sort_by_key(keys, values, size)
end subroutine sort_keys_and_values
end module thrust_sort
program test_sort
use thrust_sort
implicit none
integer(c_int), dimension(5) :: keys = [3, 1, 4, 1, 5]
integer(c_int), dimension(5) :: values = [9, 2, 8, 7, 6]
call sort_keys_and_values(keys, values, 5)
print *, 'Sorted keys: ', keys
print *, 'Sorted values: ', values
end program test_sort
nvcc -c thrust_sort.cu -o thrust_sort.o
nvfortran thrust_sort.F90 test_sort.F90 thrust_sort.o -o test_sort_app
かなり優秀にコードを出力できたとは思うが、コンパイルは通らなかった。
なお、そもそもの連携法について解説はこちら (https://www.slideshare.net/ssuserf87701/pgi-cuda-fortrangpu) などが詳しい。
コンパイルに関して修正
最初に出てきたエラーは「ライブラリへのリンクができていない」というものだった。試行錯誤するうちにnvfortranのオプションにc++ライブラリへのリンクをつけること、PTX(Parallel Thread Execution)のバージョンが実行時とコンパイル時でことなることが問題であるということがわかってきた
(参考サイトを忘れてしまいました。。。)
2. 最終コード
ソートや総和計算は型が違う場合それぞれ定義する必要がるので、bool,int,float,doubleの4パターンで用意した。
CUDA Cコード
#include <thrust/device_ptr.h>
#include <thrust/device_vector.h>
#include <thrust/sort.h>
#include <thrust/reduce.h>
#include <thrust/functional.h>
// #include <thrust/execution_policy.h>
extern "C"
{
// ######################################### //
// ############### Reduction ############### //
// ######################################### //
// Reduce for bool array (true =1, false = 0)
int thrust_reduce_bool(bool *d_array, int N)
{
thrust::device_ptr<bool> d_ptr(d_array);
int result = thrust::reduce(d_ptr, d_ptr + N, 0, thrust::plus<int>());
return result;
}
// Reduce for int array
int thrust_reduce_int(int *d_array, int N)
{
thrust::device_ptr<int> d_ptr(d_array);
int result = thrust::reduce(d_ptr, d_ptr + N, 0, thrust::plus<int>());
return result;
}
// Reduce for float array
float thrust_reduce_float(float *d_array, int N)
{
thrust::device_ptr<float> d_ptr(d_array);
float result = thrust::reduce(d_ptr, d_ptr + N, 0.0f, thrust::plus<float>());
return result;
}
// Reduce for double array
double thrust_reduce_double(double *d_array, int N)
{
thrust::device_ptr<double> d_ptr(d_array);
double result = thrust::reduce(d_ptr, d_ptr + N, 0.0, thrust::plus<double>());
return result;
}
// ######################################### //
// ############## SORT_BY_KEY ############## //
// ######################################### //
// Sort int array by int key array
void thrust_sort_int_by_intkey(int *keys, int *values, int n)
{
thrust::device_ptr<int> dev_keys(keys);
thrust::device_ptr<int> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
// Sort float array by int key array
void thrust_sort_float_by_intkey(int *keys, float *values, int n)
{
thrust::device_ptr<int> dev_keys(keys);
thrust::device_ptr<float> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
// Sort double array by bool key array
void thrust_sort_double_by_intkey(int *keys, double *values, int n)
{
thrust::device_ptr<int> dev_keys(keys);
thrust::device_ptr<double> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
// Sort int array by bool key array
void thrust_sort_int_by_boolkey(bool *keys, int *values, int n)
{
thrust::device_ptr<bool> dev_keys(keys);
thrust::device_ptr<int> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
// Sort float array by bool key array
void thrust_sort_float_by_boolkey(bool *keys, float *values, int n)
{
thrust::device_ptr<bool> dev_keys(keys);
thrust::device_ptr<float> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
// Sort double array by bool key array
void thrust_sort_double_by_boolkey(bool *keys, double *values, int n)
{
thrust::device_ptr<bool> dev_keys(keys);
thrust::device_ptr<double> dev_values(values);
thrust::stable_sort_by_key(dev_keys, dev_keys + n, dev_values);
}
}
Fortranで呼び出すためのモジュール
module thrust_module
use iso_c_binding
implicit none
interface
! ################ REDUCTION ################# !
function dev_reduce_bool(values, n) bind(C, name="thrust_reduce_bool") result(result)
use iso_c_binding
implicit none
logical(c_bool), device :: values(*)
integer(c_int), value :: n
integer(c_int) :: result
end function dev_reduce_bool
function dev_reduce_int(values, n) bind(C, name="thrust_reduce_int") result(result)
use iso_c_binding
implicit none
integer(c_int), value :: n
integer(c_int), device :: values(*)
integer(c_int) :: result
end function dev_reduce_int
function dev_reduce_float(values, n) bind(C, name="thrust_reduce_float") result(result)
use iso_c_binding
integer(c_int), value :: n
real(c_float), device :: values(*)
real(c_float) :: result
end function dev_reduce_float
function dev_reduce_double(values, n) bind(C, name="thrust_reduce_double") result(result)
use iso_c_binding
integer(c_int), value :: n
double precision(c_double), device :: values(*)
double precision(c_double) :: result
end function dev_reduce_double
! ################ REDUCTION ################# !
! ################ SORT BY KEY ################### !
subroutine dev_sort_int_by_intkey(keys, values, n) bind(C, name="thrust_sort_int_by_intkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
integer(c_int), device :: keys(*)
integer(c_int), device :: values(*)
end subroutine dev_sort_int_by_intkey
subroutine dev_sort_float_by_intkey(keys, values, n) bind(C, name="thrust_sort_float_by_intkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
integer(c_int), device :: keys(*)
real(c_float), device :: values(*)
end subroutine dev_sort_float_by_intkey
subroutine dev_sort_double_by_intkey(keys, values, n) bind(C, name="thrust_sort_double_by_intkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
integer(c_int), device :: keys(*)
double precision(c_double), device :: values(*)
end subroutine dev_sort_double_by_intkey
! --------------------------------------------------
subroutine dev_sort_int_by_boolkey(keys, values, n) bind(C, name="thrust_sort_int_by_boolkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
logical(c_bool), device :: keys(*)
integer(c_int), device :: values(*)
end subroutine dev_sort_int_by_boolkey
subroutine dev_sort_float_by_boolkey(keys, values, n) bind(C, name="thrust_sort_float_by_boolkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
logical(c_bool), device :: keys(*)
real(c_float), device :: values(*)
end subroutine dev_sort_float_by_boolkey
subroutine dev_sort_double_by_boolkey(keys, values, n) bind(C, name="thrust_sort_double_by_boolkey")
use iso_c_binding
implicit none
integer(c_int), value :: n
logical(c_bool), device :: keys(*)
double precision(c_double), device :: values(*)
end subroutine dev_sort_double_by_boolkey
! ################ SORT BY KEY ################### !
end interface
end module thrust_module
Mainプログラム
program main
use cudafor
use thrust_module
implicit none
integer, parameter :: n = 10
logical(c_bool), device :: d_keys_bool(n)
integer, device :: d_keys(n)
integer, device :: d_int_values(n)
real, device :: d_float_values(n)
double precision, device :: d_double_values(n)
logical(c_bool) :: keys_bool(n)
integer :: keys(n)
integer :: int_values(n)
real :: float_values(n)
double precision :: double_values(n)
integer :: i
integer :: sum_result_bool
integer :: sum_result_int
real :: sum_result_float
double precision :: sum_result_double
! -----------bool 型キーを用いて整列する-------------
! 初期データを設定
keys_bool = (/1, 1, 0, 1, 0, 1, 0, 1, 0, 1/)
int_values = (/3, 1, 4, 1, 5, 9, 2, 6, 5, 3/)
float_values = (/3.0, 1.0, 4.0, 1.1, 5.0, 9.0, 2.0, 6.0, 5.1, 3.1/)
double_values = (/3.0d0, 1.0d0, 4.0d0, 1.1d0, 5.0d0, 9.0d0, 2.0d0, 6.0d0, 5.1d0, 3.1d0/)
! デバイスメモリにコピー
print *, '-------------------------------'
print *, 'Sort by Bool keys: ', keys_bool
d_keys_bool = keys_bool
d_int_values = int_values
d_float_values = float_values
d_double_values = double_values
call dev_sort_int_by_boolkey(d_keys_bool, d_int_values, n)
d_keys_bool = keys_bool
call dev_sort_float_by_boolkey(d_keys_bool, d_float_values, n)
d_keys_bool = keys_bool
call dev_sort_double_by_boolkey(d_keys_bool, d_double_values, n)
! 結果をホストにコピー
keys_bool = d_keys_bool
keys = d_keys
int_values = d_int_values
float_values = d_float_values
double_values = d_double_values
! 結果を表示
print *, 'Sorted keys: ', keys_bool
print *, 'Sorted int values: ', int_values
print *, 'Sorted float values: ', float_values
print *, 'Sorted double values: ', double_values
print *, '-------------------------------'
! -----------int 型キーを用いて整列する-------------
! 初期データを設定
keys = (/3, 1, 4, 1, 5, 9, 2, 6, 5, 3/)
int_values = (/3, 1, 4, 1, 5, 9, 2, 6, 5, 3/)
float_values = (/3.0, 1.0, 4.0, 1.1, 5.0, 9.0, 2.0, 6.0, 5.1, 3.1/)
double_values = (/3.0d0, 1.0d0, 4.0d0, 1.1d0, 5.0d0, 9.0d0, 2.0d0, 6.0d0, 5.1d0, 3.1d0/)
! デバイスメモリにコピー
d_keys = keys
d_int_values = int_values
d_float_values = float_values
d_double_values = double_values
! Thrustを使ってキーによる安定ソートを実行
print *, '-------------------------------'
print *, 'Sort by Integer keys: ', keys
call dev_sort_int_by_intkey(d_keys, d_int_values, n)
d_keys = keys
call dev_sort_float_by_intkey(d_keys, d_float_values, n)
d_keys = keys
call dev_sort_double_by_intkey(d_keys, d_double_values, n)
! 結果をホストにコピー
keys_bool = d_keys_bool
keys = d_keys
int_values = d_int_values
float_values = d_float_values
double_values = d_double_values
! 結果を表示
print *, 'Sorted keys: ', keys
print *, 'Sorted int values: ', int_values
print *, 'Sorted float values: ', float_values
print *, 'Sorted double values: ', double_values
print *, '-------------------------------'
!-----------------------------------
! -------- Reduceで総和を求める----------
sum_result_bool = dev_reduce_bool(d_keys_bool, n)
sum_result_int = dev_reduce_int(d_int_values, n)
sum_result_float = dev_reduce_float(d_float_values, n)
sum_result_double = dev_reduce_double(d_double_values, n)
print *, '-------------------------------'
print *, "Sum(bol): ", sum_result_bool
print *, "Sum(int): ", sum_result_int
print *, "Sum(flt): ", sum_result_float
print *, "Sum(dbl): ", sum_result_double
print *, '-------------------------------'
end program main
boolでソート、intでソート、総和の計算が入っている。
ソートするたびにキー配列もソートされることに注意。
コンパイル・実行コマンド
nvcc -c -arch=sm_70 thrust_sort.cu -o thrust_sort.o
nvfortran -cuda -gpu=cc70 thrust_sort.o thrust_module.cuf main.cuf -o main -lstdc++
./main
3. 実行結果
-------------------------------
Sort by Bool keys: T T F T F T F T F T
Sorted keys: F F F F T T T T T T
Sorted int values: 4 5 2 5
3 1 1 9 6 3
Sorted float values: 4.000000 5.000000 2.000000
5.100000 3.000000 1.000000 1.100000
9.000000 6.000000 3.100000
Sorted double values: 4.000000000000000 5.000000000000000
2.000000000000000 5.100000000000000 3.000000000000000
1.000000000000000 1.100000000000000 9.000000000000000
6.000000000000000 3.100000000000000
-------------------------------
-------------------------------
Sort by Integer keys: 3 1 4 1
5 9 2 6 5 3
Sorted keys: 1 1 2 3 3
4 5 5 6 9
Sorted int values: 1 1 2 3
3 4 5 5 6 9
Sorted float values: 1.000000 1.100000 2.000000
3.000000 3.100000 4.000000 5.000000
5.100000 6.000000 9.000000
Sorted double values: 1.000000000000000 1.100000000000000
2.000000000000000 3.000000000000000 3.100000000000000
4.000000000000000 5.000000000000000 5.100000000000000
6.000000000000000 9.000000000000000
-------------------------------
-------------------------------
Sum(bol): 6
Sum(int): 39
Sum(flt): 39.30000
Sum(dbl): 39.30000000000000
-------------------------------