{-# LANGUAGE MagicHash, UnboxedTuples #-}

{-# OPTIONS_GHC -fno-warn-missing-signatures #-}

{- |
Module      :  Numeric.GSL.Interpolation
Copyright   :  (c) Matthew Peddie 2015
License     :  GPL
Maintainer  :  Alberto Ruiz
Stability   :  provisional

Interpolation routines.

<https://www.gnu.org/software/gsl/manual/html_node/Interpolation.html#Interpolation>

The GSL routines @gsl_spline_eval@ and friends are used, but in spite
of the names, they are not restricted to spline interpolation.  The
functions in this module will work for any 'InterpolationMethod'.

-}


module Numeric.GSL.Interpolation (
  -- * Interpolation methods
  InterpolationMethod(..)
  -- * Evaluation of interpolated functions
  , evaluate
  , evaluateV
    -- * Evaluation of derivatives of interpolated functions
  , evaluateDerivative
  , evaluateDerivative2
  , evaluateDerivativeV
  , evaluateDerivative2V
    -- * Evaluation of integrals of interpolated functions
  , evaluateIntegral
  , evaluateIntegralV
) where

import Numeric.LinearAlgebra(Vector, fromList, size, Numeric)
import Foreign.C.Types
import Foreign.Marshal.Alloc(alloca)
import Foreign.Ptr(Ptr)
import Foreign.Storable(peek)
import Numeric.GSL.Internal
import System.IO.Unsafe(unsafePerformIO)

-- FIXME
import qualified Data.Vector.Storable as S
import GHC.Base (IO(..), realWorld#)

data InterpolationMethod = Linear
                         | Polynomial
                         | CSpline
                         | CSplinePeriodic
                         | Akima
                         | AkimaPeriodic
                         deriving (InterpolationMethod -> InterpolationMethod -> Bool
(InterpolationMethod -> InterpolationMethod -> Bool)
-> (InterpolationMethod -> InterpolationMethod -> Bool)
-> Eq InterpolationMethod
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: InterpolationMethod -> InterpolationMethod -> Bool
$c/= :: InterpolationMethod -> InterpolationMethod -> Bool
== :: InterpolationMethod -> InterpolationMethod -> Bool
$c== :: InterpolationMethod -> InterpolationMethod -> Bool
Eq, Int -> InterpolationMethod -> ShowS
[InterpolationMethod] -> ShowS
InterpolationMethod -> String
(Int -> InterpolationMethod -> ShowS)
-> (InterpolationMethod -> String)
-> ([InterpolationMethod] -> ShowS)
-> Show InterpolationMethod
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [InterpolationMethod] -> ShowS
$cshowList :: [InterpolationMethod] -> ShowS
show :: InterpolationMethod -> String
$cshow :: InterpolationMethod -> String
showsPrec :: Int -> InterpolationMethod -> ShowS
$cshowsPrec :: Int -> InterpolationMethod -> ShowS
Show, ReadPrec [InterpolationMethod]
ReadPrec InterpolationMethod
Int -> ReadS InterpolationMethod
ReadS [InterpolationMethod]
(Int -> ReadS InterpolationMethod)
-> ReadS [InterpolationMethod]
-> ReadPrec InterpolationMethod
-> ReadPrec [InterpolationMethod]
-> Read InterpolationMethod
forall a.
(Int -> ReadS a)
-> ReadS [a] -> ReadPrec a -> ReadPrec [a] -> Read a
readListPrec :: ReadPrec [InterpolationMethod]
$creadListPrec :: ReadPrec [InterpolationMethod]
readPrec :: ReadPrec InterpolationMethod
$creadPrec :: ReadPrec InterpolationMethod
readList :: ReadS [InterpolationMethod]
$creadList :: ReadS [InterpolationMethod]
readsPrec :: Int -> ReadS InterpolationMethod
$creadsPrec :: Int -> ReadS InterpolationMethod
Read)

methodToInt :: Integral a => InterpolationMethod -> a
methodToInt :: forall a. Integral a => InterpolationMethod -> a
methodToInt InterpolationMethod
Linear = a
0
methodToInt InterpolationMethod
Polynomial = a
1
methodToInt InterpolationMethod
CSpline = a
2
methodToInt InterpolationMethod
CSplinePeriodic = a
3
methodToInt InterpolationMethod
Akima = a
4
methodToInt InterpolationMethod
AkimaPeriodic = a
5

dim :: Numeric t => Vector t -> Int
dim :: forall t. Numeric t => Vector t -> Int
dim = Vector t -> Int
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size

-- FIXME
appVector :: (Ptr a -> a) -> Vector a -> a
appVector Ptr a -> a
f Vector a
x = IO a -> a
forall {a}. IO a -> a
unsafeInlinePerformIO (Vector a -> (Ptr a -> IO a) -> IO a
forall a b. Storable a => Vector a -> (Ptr a -> IO b) -> IO b
S.unsafeWith Vector a
x (a -> IO a
forall (m :: * -> *) a. Monad m => a -> m a
return (a -> IO a) -> (Ptr a -> a) -> Ptr a -> IO a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Ptr a -> a
f))

unsafeInlinePerformIO :: IO a -> a
unsafeInlinePerformIO (IO State# RealWorld -> (# State# RealWorld, a #)
f) = case State# RealWorld -> (# State# RealWorld, a #)
f State# RealWorld
realWorld# of
    (# State# RealWorld
_, a
x #) -> a
x

applyCFun :: String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
hsname String
cname Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt
fun InterpolationMethod
mth Vector t
xs Vector a
ys t
x
  | Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Vector a -> Int
forall t. Numeric t => Vector t -> Int
dim Vector a
ys = String -> a
forall a. HasCallStack => String -> a
error (String -> a) -> String -> a
forall a b. (a -> b) -> a -> b
$
                         String
"Error: Vectors of unequal sizes " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                         Int -> String
forall a. Show a => a -> String
show (Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" and " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show (Vector a -> Int
forall t. Numeric t => Vector t -> Int
dim Vector a
ys) String -> ShowS
forall a. [a] -> [a] -> [a]
++
                         String
" supplied to " String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
hsname
  | Bool
otherwise = IO a -> a
forall {a}. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$
      ((Ptr t -> IO a) -> Vector t -> IO a)
-> Vector t -> (Ptr t -> IO a) -> IO a
forall a b c. (a -> b -> c) -> b -> a -> c
flip (Ptr t -> IO a) -> Vector t -> IO a
forall {a} {a}. Storable a => (Ptr a -> a) -> Vector a -> a
appVector Vector t
xs ((Ptr t -> IO a) -> IO a) -> (Ptr t -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr t
xs' ->
       ((Ptr a -> IO a) -> Vector a -> IO a)
-> Vector a -> (Ptr a -> IO a) -> IO a
forall a b c. (a -> b -> c) -> b -> a -> c
flip (Ptr a -> IO a) -> Vector a -> IO a
forall {a} {a}. Storable a => (Ptr a -> a) -> Vector a -> a
appVector Vector a
ys ((Ptr a -> IO a) -> IO a) -> (Ptr a -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr a
ys' ->
        (Ptr a -> IO a) -> IO a
forall a b. Storable a => (Ptr a -> IO b) -> IO b
alloca ((Ptr a -> IO a) -> IO a) -> (Ptr a -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr a
y' -> do
          Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt
fun Ptr t
xs' Ptr a
ys'
            (Int -> t
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> t) -> Int -> t
forall a b. (a -> b) -> a -> b
$ Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs) t
x
            (InterpolationMethod -> t
forall a. Integral a => InterpolationMethod -> a
methodToInt InterpolationMethod
mth) Ptr a
y' IO CInt -> (IO CInt -> IO ()) -> IO ()
forall x y. x -> (x -> y) -> y
// String -> IO CInt -> IO ()
check String
cname
          Ptr a -> IO a
forall a. Storable a => Ptr a -> IO a
peek Ptr a
y'

foreign import ccall safe "spline_eval" c_spline_eval
  :: Ptr Double -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt

--------------------------------------------------------------------
{- | Evaluate a function by interpolating within the given dataset.  For
example:

>>> let xs = vector [1..10]
>>> let ys = vector $ map (**2) [1..10]
>>> evaluateV CSpline xs ys 2.2
4.818867924528303

To successfully @evaluateV xs ys x@, the vectors of corresponding
domain-range values @xs@ and @ys@ must have identical lengths, and
@xs@ must be monotonically increasing.  The evaluation point @x@ must
lie between the smallest and largest values in @xs@.

-}
evaluateV :: InterpolationMethod  -- ^ What method to use to interpolate
             -> Vector Double     -- ^ Data points sampling the domain of the function
             -> Vector Double     -- ^ Data points sampling the range of the function
             -> Double            -- ^ Point at which to evaluate the function
             -> Double            -- ^ Interpolated result
evaluateV :: InterpolationMethod
-> Vector Double -> Vector Double -> Double -> Double
evaluateV = String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluateV" String
"spline_eval" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval

{- | Evaluate a function by interpolating within the given dataset.  For
example:

>>> let xs = [1..10]
>>> let ys map (**2) [1..10]
>>> evaluate Akima (zip xs ys) 2.2
4.840000000000001

To successfully @evaluate points x@, the domain (@x@) values in
@points@ must be monotonically increasing.  The evaluation point @x@
must lie between the smallest and largest values in the sampled
domain.

-}
evaluate :: InterpolationMethod    -- ^ What method to use to interpolate
            -> [(Double, Double)]  -- ^ (domain, range) values sampling the function
            -> Double              -- ^ Point at which to evaluate the function
            -> Double              -- ^ Interpolated result
evaluate :: InterpolationMethod -> [(Double, Double)] -> Double -> Double
evaluate InterpolationMethod
mth [(Double, Double)]
pts =
  String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluate" String
"spline_eval" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval
  InterpolationMethod
mth ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
xs) ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
ys)
  where
    ([Double]
xs, [Double]
ys) = [(Double, Double)] -> ([Double], [Double])
forall a b. [(a, b)] -> ([a], [b])
unzip [(Double, Double)]
pts

foreign import ccall safe "spline_eval_deriv" c_spline_eval_deriv
  :: Ptr Double -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt

{- | Evaluate the derivative of a function by interpolating within the
given dataset.  For example:

>>> let xs = vector [1..10]
>>> let ys = vector $ map (**2) [1..10]
>>> evaluateDerivativeV CSpline xs ys 2.2
4.338867924528302

To successfully @evaluateDerivativeV xs ys x@, the vectors of
corresponding domain-range values @xs@ and @ys@ must have identical
lengths, and @xs@ must be monotonically increasing.  The interpolation
point @x@ must lie between the smallest and largest values in @xs@.

-}
evaluateDerivativeV :: InterpolationMethod  -- ^ What method to use to interpolate
                       -> Vector Double     -- ^ Data points @xs@ sampling the domain of the function
                       -> Vector Double     -- ^ Data points @ys@ sampling the range of the function
                       -> Double            -- ^ Point @x@ at which to evaluate the derivative
                       -> Double            -- ^ Interpolated result
evaluateDerivativeV :: InterpolationMethod
-> Vector Double -> Vector Double -> Double -> Double
evaluateDerivativeV =
  String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluateDerivativeV" String
"spline_eval_deriv" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval_deriv

{- | Evaluate the derivative of a function by interpolating within the
given dataset.  For example:

>>> let xs = [1..10]
>>> let ys map (**2) [1..10]
>>> evaluateDerivative Akima (zip xs ys) 2.2
4.4

To successfully @evaluateDerivative points x@, the domain (@x@) values
in @points@ must be monotonically increasing.  The evaluation point
@x@ must lie between the smallest and largest values in the sampled
domain.

-}
evaluateDerivative :: InterpolationMethod    -- ^ What method to use to interpolate
                      -> [(Double, Double)]  -- ^ (domain, range) points sampling the function
                      -> Double              -- ^ Point @x@ at which to evaluate the derivative
                      -> Double              -- ^ Interpolated result
evaluateDerivative :: InterpolationMethod -> [(Double, Double)] -> Double -> Double
evaluateDerivative InterpolationMethod
mth [(Double, Double)]
pts =
  String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluateDerivative" String
"spline_eval_deriv" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval_deriv
  InterpolationMethod
mth ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
xs) ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
ys)
  where
    ([Double]
xs, [Double]
ys) = [(Double, Double)] -> ([Double], [Double])
forall a b. [(a, b)] -> ([a], [b])
unzip [(Double, Double)]
pts

foreign import ccall safe "spline_eval_deriv2" c_spline_eval_deriv2
  :: Ptr Double -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt

{- | Evaluate the second derivative of a function by interpolating within the
given dataset.  For example:

>>> let xs = vector [1..10]
>>> let ys = vector $ map (**2) [1..10]
>>> evaluateDerivative2V CSpline xs ys 2.2
2.4

To successfully @evaluateDerivative2V xs ys x@, the vectors @xs@ and
@ys@ must have identical lengths, and @xs@ must be monotonically
increasing.  The evaluation point @x@ must lie between the smallest
and largest values in @xs@.

-}
evaluateDerivative2V :: InterpolationMethod  -- ^ What method to use to interpolate
                        -> Vector Double     -- ^ Data points @xs@ sampling the domain of the function
                        -> Vector Double     -- ^ Data points @ys@ sampling the range of the function
                        -> Double            -- ^ Point @x@ at which to evaluate the second derivative
                        -> Double            -- ^ Interpolated result
evaluateDerivative2V :: InterpolationMethod
-> Vector Double -> Vector Double -> Double -> Double
evaluateDerivative2V =
  String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluateDerivative2V" String
"spline_eval_deriv2" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval_deriv2

{- | Evaluate the second derivative of a function by interpolating
within the given dataset.  For example:

>>> let xs = [1..10]
>>> let ys map (**2) [1..10]
>>> evaluateDerivative2 Akima (zip xs ys) 2.2
2.0

To successfully @evaluateDerivative2 points x@, the domain (@x@)
values in @points@ must be monotonically increasing.  The evaluation
point @x@ must lie between the smallest and largest values in the
sampled domain.

-}
evaluateDerivative2 :: InterpolationMethod    -- ^ What method to use to interpolate
                       -> [(Double, Double)]  -- ^ (domain, range) points sampling the function
                       -> Double              -- ^ Point @x@ at which to evaluate the second derivative
                       -> Double              -- ^ Interpolated result
evaluateDerivative2 :: InterpolationMethod -> [(Double, Double)] -> Double -> Double
evaluateDerivative2 InterpolationMethod
mth [(Double, Double)]
pts =
  String
-> String
-> (Ptr Double
    -> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> a
applyCFun String
"evaluateDerivative2" String
"spline_eval_deriv2" Ptr Double
-> Ptr Double -> CUInt -> Double -> CInt -> Ptr Double -> IO CInt
c_spline_eval_deriv2
  InterpolationMethod
mth ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
xs) ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
ys)
  where
    ([Double]
xs, [Double]
ys) = [(Double, Double)] -> ([Double], [Double])
forall a b. [(a, b)] -> ([a], [b])
unzip [(Double, Double)]
pts

foreign import ccall safe "spline_eval_integ" c_spline_eval_integ
  :: Ptr Double -> Ptr Double -> CUInt -> Double -> Double -> CInt -> Ptr Double -> IO CInt

applyCIntFun :: String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> t
-> a
applyCIntFun String
hsname String
cname Ptr t -> Ptr a -> t -> t -> t -> t -> Ptr a -> IO CInt
fun InterpolationMethod
mth Vector t
xs Vector a
ys t
a t
b
  | Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= Vector a -> Int
forall t. Numeric t => Vector t -> Int
dim Vector a
ys = String -> a
forall a. HasCallStack => String -> a
error (String -> a) -> String -> a
forall a b. (a -> b) -> a -> b
$
                         String
"Error: Vectors of unequal sizes " String -> ShowS
forall a. [a] -> [a] -> [a]
++
                         Int -> String
forall a. Show a => a -> String
show (Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs) String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" and " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Int -> String
forall a. Show a => a -> String
show (Vector a -> Int
forall t. Numeric t => Vector t -> Int
dim Vector a
ys) String -> ShowS
forall a. [a] -> [a] -> [a]
++
                         String
" supplied to " String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
hsname
  | Bool
otherwise = IO a -> a
forall {a}. IO a -> a
unsafePerformIO (IO a -> a) -> IO a -> a
forall a b. (a -> b) -> a -> b
$
      ((Ptr t -> IO a) -> Vector t -> IO a)
-> Vector t -> (Ptr t -> IO a) -> IO a
forall a b c. (a -> b -> c) -> b -> a -> c
flip (Ptr t -> IO a) -> Vector t -> IO a
forall {a} {a}. Storable a => (Ptr a -> a) -> Vector a -> a
appVector Vector t
xs ((Ptr t -> IO a) -> IO a) -> (Ptr t -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr t
xs' ->
       ((Ptr a -> IO a) -> Vector a -> IO a)
-> Vector a -> (Ptr a -> IO a) -> IO a
forall a b c. (a -> b -> c) -> b -> a -> c
flip (Ptr a -> IO a) -> Vector a -> IO a
forall {a} {a}. Storable a => (Ptr a -> a) -> Vector a -> a
appVector Vector a
ys ((Ptr a -> IO a) -> IO a) -> (Ptr a -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr a
ys' ->
        (Ptr a -> IO a) -> IO a
forall a b. Storable a => (Ptr a -> IO b) -> IO b
alloca ((Ptr a -> IO a) -> IO a) -> (Ptr a -> IO a) -> IO a
forall a b. (a -> b) -> a -> b
$ \Ptr a
y' -> do
          Ptr t -> Ptr a -> t -> t -> t -> t -> Ptr a -> IO CInt
fun Ptr t
xs' Ptr a
ys'
            (Int -> t
forall a b. (Integral a, Num b) => a -> b
fromIntegral (Int -> t) -> Int -> t
forall a b. (a -> b) -> a -> b
$ Vector t -> Int
forall t. Numeric t => Vector t -> Int
dim Vector t
xs) t
a t
b
            (InterpolationMethod -> t
forall a. Integral a => InterpolationMethod -> a
methodToInt InterpolationMethod
mth) Ptr a
y' IO CInt -> (IO CInt -> IO ()) -> IO ()
forall x y. x -> (x -> y) -> y
// String -> IO CInt -> IO ()
check String
cname
          Ptr a -> IO a
forall a. Storable a => Ptr a -> IO a
peek Ptr a
y'

{- | Evaluate the definite integral of a function by interpolating
within the given dataset.  For example:

>>> let xs = vector [1..10]
>>> let ys = vector $ map (**2) [1..10]
>>> evaluateIntegralV CSpline xs ys 2.2 5.5
51.89853207547169

To successfully @evaluateIntegralV xs ys a b@, the vectors @xs@ and
@ys@ must have identical lengths, and @xs@ must be monotonically
increasing.  The integration bounds @a@ and @b@ must lie between the
smallest and largest values in @xs@.

-}
evaluateIntegralV :: InterpolationMethod  -- ^ What method to use to interpolate
                     -> Vector Double     -- ^ Data points @xs@ sampling the domain of the function
                     -> Vector Double     -- ^ Data points @ys@ sampling the range of the function
                     -> Double            -- ^ Lower integration bound @a@
                     -> Double            -- ^ Upper integration bound @b@
                     -> Double            -- ^ Resulting area
evaluateIntegralV :: InterpolationMethod
-> Vector Double -> Vector Double -> Double -> Double -> Double
evaluateIntegralV =
  String
-> String
-> (Ptr Double
    -> Ptr Double
    -> CUInt
    -> Double
    -> Double
    -> CInt
    -> Ptr Double
    -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> t
-> a
applyCIntFun String
"evaluateIntegralV" String
"spline_eval_integ" Ptr Double
-> Ptr Double
-> CUInt
-> Double
-> Double
-> CInt
-> Ptr Double
-> IO CInt
c_spline_eval_integ

{- | Evaluate the definite integral of a function by interpolating
within the given dataset.  For example:

>>> let xs = [1..10]
>>> let ys = map (**2) [1..10]
>>> evaluateIntegralV CSpline (zip xs ys) (2.2, 5.5)
51.909

To successfully @evaluateIntegral points (a, b)@, the domain (@x@)
values of @points@ must be monotonically increasing.  The integration
bounds @a@ and @b@ must lie between the smallest and largest values in
the sampled domain..
-}
evaluateIntegral :: InterpolationMethod    -- ^ What method to use to interpolate
                    -> [(Double, Double)]  -- ^ (domain, range) points sampling the function
                    -> (Double, Double)    -- ^ Integration bounds (@a@, @b@)
                    -> Double              -- ^ Resulting area
evaluateIntegral :: InterpolationMethod
-> [(Double, Double)] -> (Double, Double) -> Double
evaluateIntegral InterpolationMethod
mth [(Double, Double)]
pts (Double
a, Double
b) =
  String
-> String
-> (Ptr Double
    -> Ptr Double
    -> CUInt
    -> Double
    -> Double
    -> CInt
    -> Ptr Double
    -> IO CInt)
-> InterpolationMethod
-> Vector Double
-> Vector Double
-> Double
-> Double
-> Double
forall {t} {a} {a} {t} {t} {t} {t}.
(Numeric t, Numeric a, Storable a, Integral t, Num t) =>
String
-> String
-> (Ptr t -> Ptr a -> t -> t -> t -> t -> Ptr a -> IO CInt)
-> InterpolationMethod
-> Vector t
-> Vector a
-> t
-> t
-> a
applyCIntFun String
"evaluateIntegral" String
"spline_eval_integ" Ptr Double
-> Ptr Double
-> CUInt
-> Double
-> Double
-> CInt
-> Ptr Double
-> IO CInt
c_spline_eval_integ
  InterpolationMethod
mth ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
xs) ([Double] -> Vector Double
forall a. Storable a => [a] -> Vector a
fromList [Double]
ys) Double
a Double
b
  where
    ([Double]
xs, [Double]
ys) = [(Double, Double)] -> ([Double], [Double])
forall a b. [(a, b)] -> ([a], [b])
unzip [(Double, Double)]
pts