Cusp. Finite difference library ?

Accelerated Computing CUDA CUDA Programming and Performance
#1

Hi
Is there any libray e.g. Cusp for computing finite difference as an approximation of derivative ?

I am honest, I am doing an exercise for a course and it is asked to compare my implementation to an already developed library.
I am having hard time to find anything comparable. I need something doing a 1D derivative or maybe a finite difference approximation using my stencil (cumulative sum of the product of the function times the coefficients).
Any idea ?
Cheers

#2

Have a look at OpenCurrent.

#3

Thanks for the reply.

There is a lot of staff on OpenCurrent, even a complete implementation of finite difference in 3D. However, I only have a simple problem solved with a fast coefficient computation and a sparse matrices multiplication (I was thinking at the gradient descent in the Krylov method (since I have sparce matrices with a bunch of zeros and just few coefficients).

However It is not so easy to interpret what is going on inside Cusp. How can I display (eg. with printf) the content of a cusp-array1D (the last line in my code gives an error).

Assume I have this code with Ax=b, where A is a matrix of coefficient, x is a discretization of point of a function and b is the derivative as in the picture (solving the equations I obtain the derivatives as b):

#include <cusp/csr_matrix.h>

#include <cusp/monitor.h>

#include <cusp/krylov/cg.h>

#include <cusp/gallery/poisson.h>

#include <stdlib.h>

int main(void)

   {

       // create an empty sparse matrix structure (CSR format)

       cusp::csr_matrix<int, float, cusp::device_memory> A;

// initialize matrix

       cusp::gallery::poisson5pt(A, 10, 10);

// allocate storage for solution (x) and right hand side (b)

       cusp::array1d<float, cusp::device_memory> x(A.num_rows, 0);

      cusp::array1d<float, cusp::device_memory> b(A.num_rows, 1);

// set stopping criteria:

       //  iteration_limit    = 100

       //  relative_tolerance = 1e-6

       cusp::verbose_monitor<float> monitor(b, 100, 1e-6);

// set preconditioner (identity)

       cusp::identity_operator<float, cusp::device_memory> M(A.num_rows, A.num_rows);

// solve the linear system A x = b

       cusp::krylov::cg(A, x, b, monitor, M);

for(int kk=0;kk<10;kk++)

   printf("Val of A = %i \n", x.begin() );

        return 0;

   }