I've seen a function in numpy (http://docs.scipy.org/doc/numpy/reference/generated/numpy.gradient.html) for calculating the gradient of a 2-dimensional array. Numpy also provides 2 examples on that website (I'm more interested in the second one).

I now want to use such a function in C++. Does anyone know if there is such a function in GSL? Or can someone explain to me what the numpy function does, so that I can write my own in C++?

Thanks for any help.

looking at the 1-D example(s) on the page you referenced:

grad[0] = (vals[1] - vals[0]) / dx;
grad[i] = (vals[i+1] - vals[i-1]) / (2*dx);  // for i in [1,N-2]
grad[N-1] = (vals[N-1] - vals[N-2]) / dx;

for the 2-D case, you have to figure out which direction in the input is X vs. Y (should be able to tell by specifying different values for dx and dy and seeing which output array changes in what way), but essentially the input is a 2-D array of scalar values and the output is a 3-D array, or two 2-D arrays matching the size/shape of the input. Each of the resulting 2-D arrays is the gradient in one direction: in this case, it looks like the Y gradients are output first, and the X gradients second, but again, it depends on the definition of X and Y.

Hope this helps!