[Gsll-devel] Fwd: integral-evaluate-spline crashes for xarr containing containing negative numbers

[Mirko Vukovic]

A correction regarding when the crash occurs:  The x-array lower bound
has to be exactly 0d0.  Cannot be positive or negative.

Mirko


---------- Forwarded message ----------
From: Mirko Vukovic <mirko.vukovic at gmail.com>
Date: Tue, Jan 27, 2009 at 4:16 PM
Subject: integral-evaluate-spline crashes for xarr containing
containing  negative numbers
To: gsll-devel at common-lisp.net


Hello.

To illustrate the problem consider the following routine

(defun spline-integration-example (&optional (lower-limit 0d0))
 "Integral of  pi * sin (pi x) between 0 and 0.5.  Should come close to 1.0"
 (let* ((acc (make-acceleration))
        (xarr
         (make-marray
          'double-float
          :initial-contents
          (loop for i from lower-limit 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))))

Doing (spline-integration-example) gives 1.0000108...
Doing (spline-integration-example -0.1d0) causes an error:

Unhandled memory fault at #x0.
  [Condition of type SB-SYS:MEMORY-FAULT-ERROR]

Restarts:
 0: [RETRY] Retry SLIME REPL evaluation request.
 1: [ABORT] Return to SLIME's top level.
 2: [ABORT] Exit debugger, returning to top level.


What is interesting are is the environment prior to the call to
integral-evaluate-spline (#7 below):
Backtrace:
 0: (SB-SYS:MEMORY-FAULT-ERROR)
 1: ("foreign function: #x41A9A2")
 2: ("foreign function: #x41AA80")
 3: ("foreign function: #x2AE6C85E1D1B")
 4: ("foreign function: #x2AE6C85E4683")
 5: ("foreign function: #x2AE6C85E541A")
 6: ("foreign function: #x2AE6C85E5F93")
 7: (INTEGRAL-EVALUATE-SPLINE ..)
     Locals:
       SB-DEBUG::ARG-0 = #<SPLINE {10039177D1}>
       SB-DEBUG::ARG-1 = #<VECTOR-DOUBLE-FLOAT #(-0.1d0 0.1d0
0.30000000000000004d0 0.5d0 0.7d0 ..)>
       SB-DEBUG::ARG-2 = #<VECTOR-DOUBLE-FLOAT
#(-0.3090169943749474d0 0.3090169943749474d0 ..)>
       SB-DEBUG::ARG-3 = 0.5d0
       SB-DEBUG::ARG-4 = 3.395932351202108d-313
       SB-DEBUG::ARG-5 = #<ACCELERATION {1003910A71}>
 8: (SPLINE-INTEGRATION-EXAMPLE -0.1d0)
 9: (SB-INT:SIMPLE-EVAL-IN-LEXENV (SPLINE-INTEGRATION-EXAMPLE -0.1d0)
#<NULL-LEXENV>)
 --more--

Arguments 3 & 4 are incorrect !.  They should be 0d0 and 0.5d0.  I
don't have a clue how that happens.

btw, I am using my definition of integral-evaluate-spline, which I
posted earlier today:

(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.")

Thanks to the properties of integrals, this problem is not a
show-stopper.  Just shifting the integration axis will do the job.
But it is annoying :-)

Hope Liam can help finding the fix.

Mirko