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
Have a look at OpenCurrent.
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;
}