Как правильно построить CUSP Coo матрицу из пропущенных массивов

Я пытаюсь интегрировать CUSP в существующий код Fortran. Прямо сейчас я просто пытаюсь настроить базовую настройку Thrust / CUSP для подачи массивов из Fortran и использовать их для построения матрицы CUSP (сейчас это круто).
Я смог получить оболочку, подобную подпрограмме C, скомпилировать в библиотеку и связать ее с кодом на Фортране благодаря этой теме: неразрешенные ссылки, используя-ifort-с-NVCC-и-параболическая

И я могу убедиться, что Fortran правильно вводит указатели на массивы благодаря помощи из предыдущего потока: Генерация CUSP coo_matrix из переданных массивов FORTRAN

К сожалению, я до сих пор не могу заставить CUSP использовать их для генерации и печати матрицы.
Код и вывод показаны ниже:

выход

$ ./fort_cusp_test
testing 1 2 3
n: 3, nnz: 9
i,  row_i,  col_j,        val_v
0,      1,      1,   1.0000e+00
1,      1,      2,   2.0000e+00
2,      1,      3,   3.0000e+00
3,      2,      1,   4.0000e+00
4,      2,      2,   5.0000e+00
5,      2,      3,   6.0000e+00
6,      3,      1,   7.0000e+00
7,      3,      2,   8.0000e+00
8,      3,      3,   9.0000e+00
initialized row_i into thrust
initialized col_j into thrust
initialized val_v into thrust
defined CUSP integer array view for row_i and col_j
defined CUSP float array view for val_v
loaded row_i into a CUSP integer array view
loaded col_j into a CUSP integer array view
loaded val_v into a CUSP float array view
defined CUSP coo_matrix view
Built matrix A from CUSP device views
sparse matrix <3, 3> with 9 entries
libc++abi.dylib: terminating with uncaught exception of type thrust::system::system_error: invalid argument

Program received signal SIGABRT: Process abort signal.

Backtrace for this error:
#0  0x10d0fdff6
#1  0x10d0fd593
#2  0x7fff8593ff19
Abort trap: 6

fort_cusp_test.f90

program fort_cuda_test

implicit none

! interface
!    subroutine test_coo_mat_print_(row_i,col_j,val_v,n,nnz) bind(C)
!       use, intrinsic :: ISO_C_BINDING, ONLY: C_INT,C_FLOAT
!       implicit none
!       integer(C_INT),value :: n, nnz
!       integer(C_INT) :: row_i(:), col_j(:)
!       real(C_FLOAT) :: val_v(:)
!    end subroutine test_coo_mat_print_
! end interface

integer*4   n
integer*4   nnz

integer*4, target :: rowI(9),colJ(9)
real*4, target :: valV(9)

integer*4, pointer ::   row_i(:)
integer*4, pointer ::   col_j(:)
real*4, pointer ::   val_v(:)

n     =  3
nnz   =  9
rowI =  (/ 1, 1, 1, 2, 2, 2, 3, 3, 3/)
colJ =  (/ 1, 2, 3, 1, 2, 3, 1, 2, 3/)
valV =  (/ 1, 2, 3, 4, 5, 6, 7, 8, 9/)

row_i => rowI
col_j => colJ
val_v => valV

write(*,*) "testing 1 2 3"
call test_coo_mat_print (rowI,colJ,valV,n,nnz)

end program fort_cuda_test

cusp_runner.cu

#include <stdio.h>
#include <cusp/coo_matrix.h>
#include <iostream>
// #include <cusp/krylov/cg.h>
#include <cusp/print.h>

#if defined(__cplusplus)
extern "C" {
#endif

void test_coo_mat_print_(int * row_i, int * col_j, float * val_v, int * N, int * NNZ ) {

int n, nnz;

n = *N;
nnz = *NNZ;

printf("n: %d, nnz: %d\n",n,nnz);

printf("%6s, %6s, %6s, %12s \n","i","row_i","col_j","val_v");
for(int i=0;i<n;i++) {
printf("%6d, %6d, %6d, %12.4e\n",i,row_i[i],col_j[i],val_v[i]);
}
//if ( false ) {
//wrap raw input pointers with thrust::device_ptr
thrust::device_ptr<int> wrapped_device_I(row_i);
printf("initialized row_i into thrust\n");
thrust::device_ptr<int> wrapped_device_J(col_j);
printf("initialized col_j into thrust\n");
thrust::device_ptr<float> wrapped_device_V(val_v);
printf("initialized val_v into thrust\n");

//use array1d_view to wrap individual arrays
typedef typename cusp::array1d_view< thrust::device_ptr<int> > DeviceIndexArrayView;
printf("defined CUSP integer array view for row_i and col_j\n");
typedef typename cusp::array1d_view< thrust::device_ptr<float> > DeviceValueArrayView;
printf("defined CUSP float array view for val_v\n");

DeviceIndexArrayView row_indices(wrapped_device_I, wrapped_device_I + nnz);
printf("loaded row_i into a CUSP integer array view\n");
DeviceIndexArrayView column_indices(wrapped_device_J, wrapped_device_J + nnz);
printf("loaded col_j into a CUSP integer array view\n");
DeviceValueArrayView values(wrapped_device_V, wrapped_device_V + nnz);
printf("loaded val_v into a CUSP float array view\n");

//combine array1d_views into coo_matrix_view
typedef cusp::coo_matrix_view<DeviceIndexArrayView,DeviceIndexArrayView,DeviceValueArrayView> DeviceView;
printf("defined CUSP coo_matrix view\n");

//construct coo_matrix_view from array1d_views
DeviceView A(n,n,nnz,row_indices,column_indices,values);
printf("Built matrix A from CUSP device views\n");

cusp::print(A);
printf("Printed matrix A\n");
//}
}
#if defined(__cplusplus)
}
#endif

Makefile

Test:
nvcc -Xcompiler="-fPIC" -shared cusp_runner.cu -o cusp_runner.so -I/Developer/NVIDIA/CUDA-6.5/include/cusp
gfortran -c fort_cusp_test.f90
gfortran fort_cusp_test.o cusp_runner.so -L/Developer/NVIDIA/CUDA-6.5/lib -lcudart -o fort_cusp_test

clean:
rm *.o *.so fort_cusp_test

Функциональная версия cusp_runner.cu:

#include <stdio.h>
#include <cusp/coo_matrix.h>
#include <iostream>
// #include <cusp/krylov/cg.h>
#include <cusp/print.h>

#if defined(__cplusplus)
extern "C" {
#endif

void test_coo_mat_print_(int * row_i, int * col_j, float * val_v, int * N, int * NNZ ) {

int n, nnz;

n = *N;
nnz = *NNZ;

printf("n: %d, nnz: %d\n",n,nnz);

printf("printing input (row_i, col_j, val_v)\n");
printf("%6s, %6s, %6s, %12s \n","i","row_i","col_j","val_v");
for(int i=0;i<nnz;i++) {
printf("%6d, %6d, %6d, %12.4e\n",i,row_i[i],col_j[i],val_v[i]);
}

printf("initializing thrust device vectors\n");
thrust::device_vector<int> device_I(row_i,row_i+nnz);
printf("device_I initialized\n");
thrust::device_vector<int> device_J(col_j,col_j+nnz);
printf("device_J initialized\n");
thrust::device_vector<float> device_V(val_v,val_v+nnz);
printf("device_V initialized\n");

cusp::coo_matrix<int, float, cusp::device_memory> A(n,n,nnz);
printf("initialized empty CUSP coo_matrix on device\n");

A.row_indices = device_I;
printf("loaded device_I into A.row_indices\n");
A.column_indices = device_J;
printf("loaded device_J into A.column_indices\n");
A.values = device_V;
printf("loaded device_V into A.values\n");

cusp::print(A);
printf("Printed matrix A\n");
//}
}
#if defined(__cplusplus)
}
#endif

2

Решение

Ваш код толчка / стороны CUSP для обработки указателей абсолютно неверен. Это:

thrust::device_ptr<int> wrapped_device_I(row_i);

не делает то, что вы думаете, что делает. То, что вы фактически сделали, это преобразование адреса хоста в адрес устройства. Это незаконно, если вы не работаете с управляемой памятью CUDA, и я не вижу доказательств этого в этом коде. Что вы хотите сделать, это выделить память и скопировать массивы Fortran в GPU перед запуском. Сделать что-то вроде:

thrust::device_ptr<int> wrapped_device_I = thrust::device_malloc<int>(nnz);
thrust::copy(row_i, row_i + nnz, wrapped_device_I);
[отказ от ответственности: полностью не проверено, используйте на свой страх и риск]

Для каждого из векторов COO. Однако я бы предложил заменить большую часть кода в разделе настройки графического процессора test_coo_mat_print_ с thrust::vector экземпляры вместо. Помимо простоты в использовании, вы получаете бесплатное освобождение памяти, когда они выходят из области видимости, так что вероятность утечек памяти намного меньше. Так что-то вроде:

thrust::device_vector<int> device_I(row_i, row_i + nnz);

заботится обо всем в один звонок.

В заключение, если вы разрабатываете многоязыковые кодовые базы, спроектируйте их так, чтобы код на каждом языке был полностью независимым и имел собственный собственный тестовый код. Если бы вы сделали это в этом случае, вы бы обнаружили, что часть C ++ не работает независимо от каких-либо проблем с Fortran. Это сделало бы отладку намного проще.

3

Другие решения