Библиотека CUSP, вызванная из Фортрана, не работает

Question

Я хочу повторно решить CG/BicGSTAB с помощью решателя CUSP, вызванного из Фортрана. Чтобы избежать переводов, я передаю данные Фортрана непосредственно в CUSP. Код компилируется, но прерывается во время выполнения, отмечая:

terminate called after throwing an instance of 'thrust::system::system_error'
  what():  invalid argument
terminate called recursively
Aborted (core dumped)

Не говоря уже о ядре кода, даже поток печати не происходит. Код, конечно, находится на предварительной стадии, но мне интересно, что с ним не так.

extern "C" void bicgstab_(int *device_I, int *device_J, float *device_V, float *device_x, float *device_b, int *n, int *nnz){

        int N = *n;
        int NNZ = *nnz;
        std::cout << N << " " << NNZ << " " << *device_I << std::endl;
        for(int i=0; i<N;i++)std::cout << device_I[i] << " "; std::cout << std::endl;
        for(int i=0; i<NNZ;i++)std::cout << device_J[i] << " "; std::cout << std::endl;
        for(int i=0; i<NNZ;i++)std::cout << device_V[i] << " "; std::cout << std::endl;
        for(int i=0; i<N;i++)std::cout << device_x[i] << " "; std::cout << std::endl;
        for(int i=0; i<N;i++)std::cout << device_b[i] << " "; std::cout << std::endl;

         // *NOTE* raw pointers must be wrapped with thrust::device_ptr!
    thrust::device_ptr<int> wrapped_device_I(device_I);
    thrust::device_ptr<int> wrapped_device_J(device_J);
    thrust::device_ptr<float> wrapped_device_V(device_V);
    thrust::device_ptr<float> wrapped_device_x(device_x);
    thrust::device_ptr<float> wrapped_device_b(device_b);

    // use array1d_view to wrap the individual arrays
    typedef typename cusp::array1d_view< thrust::device_ptr<int> > DeviceIndexArrayView;
    typedef typename cusp::array1d_view< thrust::device_ptr<float> > DeviceValueArrayView;
        std::cout << wrapped_device_I[3];
/*
    DeviceIndexArrayView row_indices (wrapped_device_I, wrapped_device_I + (N+1));
    DeviceIndexArrayView column_indices(wrapped_device_J, wrapped_device_J + NNZ);
    DeviceValueArrayView values (wrapped_device_V, wrapped_device_V + NNZ);
    DeviceValueArrayView x (wrapped_device_x, wrapped_device_x + N);
    DeviceValueArrayView b (wrapped_device_b, wrapped_device_b + N);

//      std::cout << device_x[0] ;
//      for(int i=0;i<NNZ;i++)std::cout << column_indices[i] << std::endl;

        // combine the three array1d_views into a csr_matrix_view
    typedef cusp::csr_matrix_view<DeviceIndexArrayView,
    DeviceIndexArrayView,
    DeviceValueArrayView> DeviceView;

    // construct a csr_matrix_view from the array1d_views
    DeviceView A(N, N, NNZ, row_indices, column_indices, values);

    // set stopping criteria:    // iteration_limit = 100    // relative_tolerance = 1e-5
    cusp::verbose_monitor<float> monitor(b, 100, 1e-5);

    // solve the linear system A * x = b with the Conjugate Gradient method
//    cusp::krylov::bicgstab(A, x, b);*/

}

Если это невозможно, я могу перейти к другому подходу, но поскольку я не уверен в правильности, я не могу принять решение. Любая помощь приветствуется.