/*
  Computes the derivative of the Function and returns it as a new Function.
  The newly created function shares the X vector with the previous one.

  WARNING: this is a very naive 3-points algorithm.
*/
static VALUE function_derivative(VALUE self)
{
  long size = function_sanity_check(self);
  const double *x = Dvector_Data_for_Read(get_x_vector(self),NULL);
  const double *y = Dvector_Data_for_Read(get_y_vector(self),NULL);
  VALUE derivative = Dvector_Create();
  long i = 0;
  /* First value */
  Dvector_Push_Double(derivative, (y[i+1] - y[i]) /(x[i+1] - x[i]));
  i++;
  while(i < (size - 1))
    {
      Dvector_Push_Double(derivative, 
                          .5 * (
                                (y[i+1] - y[i]) /(x[i+1] - x[i]) + 
                                (y[i] - y[i-1]) /(x[i] - x[i-1])
                                ));
      i++;
    }
  Dvector_Push_Double(derivative, (y[i] - y[i-1]) /(x[i] - x[i-1]));
  return Function_Create(get_x_vector(self), derivative);
}