This is the mail archive of the
gsl-discuss@sourceware.cygnus.com
mailing list for the GSL project.
Re: numerical differentiation?
I'm not sure which methods you are talking about. There are
different notions, depending on usage. Do you mean one
or all of the following?
/* differentiate a function at a point */
gsl_diff_function(const gsl_function * f, double x0, double *
dydx_x0);
/* differentiate equally spaced data in a table */
gsl_diff_table_unif(
const double * y_data,
size_t length,
double delta_x,
double * dydx_data
);
/* differentiate a table */
gsl_diff_table(
const double * y_data,
const double * x_data,
size_t length,
double * dydx_data
);
/* evaluate a derivative poinwise, using given data */
gsl_diff_point(
const double * y_data,
const double * x_data,
size_t length,
double x0,
double * dydx_x0
);
Making a general purpose interface for this sort of functionality
will require a little thought. The following issues arise:
o Both constant and non-constant spacing should be supported
for the table methods.
o 'gsl_diff_table' and 'gsl_diff_table_unif' are not complete
interfaces.
They leave open the question of what happens at the endpoints,
for instance whether or not the data is assumed periodic (which will
in itself be confusing when the spacing is not constant).
o 'gsl_diff_function' might have trouble with the domain of the
function.
You don't want the needed evaluations of the function to wander
outside the domain [ex. f(x) = log(x), x in [0, infinity), should
never try to evaluate f(-1.0e-56) or something]. Brian: is there any
facility hovering around gsl_function for declaring a domain?
o People who use 'gsl_diff_point' repeatedly on the same data,
but with different x0, may want an optimization based on
some kind of cached state, depending on the complexity
of the evaluation. We have to decide whether or not
that is worth supporting. There is already something like that
in the chebyshev evaluation functions, where you get much
of it for free because you have to cache the chebyshev fit
anyway. Of course, in the case of chebshev, you only get the
one evaluation method, which is only obliguely related to
difference quotients. My feeling is that it is not worth
supporting any sort of cached state for differencing, since
the implementation complexity will probably be low. But that
really depends on the set of methods to be encapsulated;
if somebody comes along with a whiz-bang but expensive method,
then we might feel bad about not having an interface flexible
enough to support it properly.
o As Brian points out, you really want to have some error estimate.
I'm not sure what the theory is for this. For 'gsl_diff_function',
presumably you do some kind of extrapolation (Richardson?), and
produce an error
estimate from that. But I don't know what you do for fixed data;
that
requires some assumption about the class of functions. Also,
error estimates will complicate the interface for something
like 'gsl_diff_table'; but maybe it could be something simple like
gsl_diff_table(
const double * y_data,
const double * x_data,
size_t length,
double * dydx_data,
double * dydx_error
);
where dydx_error==0 indicates no error estimate will be provided.
--
G. Jungman
Brian Gough wrote:
>
> Yes, that would be very useful, especially if there is an estimate of
> the error in the derivative.
>
> Don't feel constrained by omissions in the TODO file -- it is not an
> exhaustive list. As a general rule anything that's covered by
> existing proprietary libraries is desirable for GSL.
>
> Dave Morrison writes:
> > Hi all,
> >
> > Is there any interest in adding routines for numerical differentiation
> > to GSL?