On Tue, Jan 27, 2009 at 10:46 PM, Liam Healy lhealy@common-lisp.net wrote:
GSL's gsl_spline_ functions do exactly what you are looking for: "remember" xarr and yarr so you don't have to keep sending them. That's why they're missing from the arglist of integral-evaluate-spline in GSLL. The gsl_interp_ functions require the arrays, and that's why the GSLL equivalents take them as arguments.
Please see and compare: http://www.gnu.org/software/gsl/manual/html_node/Evaluation-of-Interpolating... double gsl_interp_eval_integ (const gsl_interp * interp, const double xa[], const double ya[], double a, double b, gsl_interp_accel * acc) which becomes #'intergral-evaluate-interpolation in GSLL and http://www.gnu.org/software/gsl/manual/html_node/Higher_002dlevel-Interface.... Function: double gsl_spline_eval_integ (const gsl_spline * spline, double a, double b, gsl_interp_accel * acc) which becomes #'integral-evaluate-spline in GSLL. They are two different functions with two different sets of arguments, though they do the same thing.
OK, but is there an error in evaluation.lisp integral-evaluate-spline?
It has as arguments x, low, high, acceleration. What is the `x' doing there? The call for gsl_spline_eval_integ is: gsl_spline_eval_integ (const gsl_spline * spline, double a, double b, gsl_interp_accel * acc). If so, should the definition be:
(defmfun integral-evaluate-spline (spline low high acceleration) "gsl_spline_eval_integ" (((mpointer spline) :pointer) (low :double) (high :double) ((mpointer acceleration) :pointer)) :c-return :double :documentation ; FDL "Find the numerical integral of an interpolated function over the range [low, high], using the spline object spline, data arrays xa and ya and the accelerator acceleration.")
(That may have been the original stumbling block that lead me to wonder down the incorrect path)
Liam
On Tue, Jan 27, 2009 at 11:29 AM, Mirko Vukovic mirko.vukovic@gmail.com wrote:
Hello,
In order to use integral-evaluate-spline, I had to modify its definition. The one below works, but I am not sure if there is a better way. Here is the definition: (defmfun integral-evaluate-spline (spline xa ya low high acceleration) "gsl_spline_eval_integ" (((mpointer spline) :pointer) ((c-pointer xa) :pointer) ((c-pointer ya) :pointer) (low :double) (high :double) ((mpointer acceleration) :pointer)) :c-return :double :documentation ; FDL "Find the numerical integral of an interpolated function over the range [low, high], using the spline object spline, data arrays xa and ya and the accelerator acceleration.")
I removed a reference to a double-float `x' and added the x&y array pointers into the arguments, as they are necessary for gsl's `gsl_interp_eval_integ'.
They are needed for this function and they are there, in integral-evaluate-interpolation. You copied and pasted the definition of integral-evaluate-spline, then put in arguments that don't belong there.
You are right. I messed up there.
(I am really not sure what is the purpose of the ..._e functions).
Here is a bit of test code that exercises the defintion by integrating pi*sin(pi x) between 0 & 0.5 which should equal 1.0:
(let* ((acc (make-acceleration)) (xarr (make-marray 'double-float :initial-contents (loop for i from 0.0d0 below 1d0 by 0.2d0 collect i))) (yarr (make-marray 'double-float :initial-contents (loop for x across (cl-array xarr) collect (sin (* x pi))))) (spline (make-spline *cubic-spline-interpolation* xarr yarr))) (* pi (integral-evaluate-spline spline xarr yarr 0d0 0.5d0 acc)))
I wonder if there is a way for `make-spline' to "remember" xarr and yarr so they don't have to be passed to `integral-evaluate-spline'.
As always, thank for the heavy lifting :-)
Mirko
Gsll-devel mailing list Gsll-devel@common-lisp.net http://common-lisp.net/cgi-bin/mailman/listinfo/gsll-devel