Ancilliary files for the start of CVODE support

3 files changed, 24 insertions, 10 deletions

diff --git a/packages/sundials/hmatrix-sundials.cabal b/packages/sundials/hmatrix-sundials.cabal
index 388f1db..4cc02c6 100644
--- a/packages/sundials/hmatrix-sundials.cabal
+++ b/packages/sundials/hmatrix-sundials.cabal
@@ -25,10 +25,12 @@ library
                        template-haskell >=2.12 && <2.13,
                        containers >=0.5 && <0.6,
                        hmatrix>=0.18
-  extra-libraries:     sundials_arkode
+  extra-libraries:     sundials_arkode,
+                       sundials_cvode
   other-extensions:    QuasiQuotes
   hs-source-dirs:      src
-  exposed-modules:     Numeric.Sundials.ARKode.ODE
+  exposed-modules:     Numeric.Sundials.ARKode.ODE,
+                       Numeric.Sundials.CVode.ODE
   other-modules:       Types,
                        Arkode
   c-sources:           src/helpers.c src/helpers.h
@@ -39,6 +41,7 @@ test-suite hmatrix-sundials-testsuite
   main-is:             Main.hs
   other-modules:       Types,
                        Numeric.Sundials.ARKode.ODE,
+                       Numeric.Sundials.CVode.ODE,
                        Arkode
   build-depends:       base >=4.10 && <4.11,
                        inline-c >=0.6 && <0.7,
@@ -52,6 +55,7 @@ test-suite hmatrix-sundials-testsuite
                        lens,
                        hspec
   hs-source-dirs:      src
-  extra-libraries:     sundials_arkode
+  extra-libraries:     sundials_arkode,
+                       sundials_cvode
   c-sources:           src/helpers.c src/helpers.h
   default-language:    Haskell2010
diff --git a/packages/sundials/src/Arkode.hsc b/packages/sundials/src/Arkode.hsc
index 9db37b5..558ce9e 100644
--- a/packages/sundials/src/Arkode.hsc
+++ b/packages/sundials/src/Arkode.hsc
@@ -10,6 +10,7 @@ import Foreign.C.Types
 #include <nvector/nvector_serial.h>
 #include <sunmatrix/sunmatrix_dense.h>
 #include <arkode/arkode.h>
+#include <cvode/cvode.h>
 
 
 #def typedef struct _generic_N_Vector SunVector;
@@ -40,6 +41,11 @@ getContentPtr ptr = (#peek SunVector, content) ptr
 getData :: Storable a => Ptr b -> IO a
 getData ptr = (#peek SunContent, data) ptr
 
+cV_ADAMS :: Int
+cV_ADAMS = #const CV_ADAMS
+cV_BDF :: Int
+cV_BDF = #const CV_BDF
+
 arkSMax :: Int
 arkSMax = #const ARK_S_MAX
 
diff --git a/packages/sundials/src/Main.hs b/packages/sundials/src/Main.hs
index 729d35a..3904b09 100644
--- a/packages/sundials/src/Main.hs
+++ b/packages/sundials/src/Main.hs
@@ -1,6 +1,7 @@
 {-# OPTIONS_GHC -Wall #-}
 
-import           Numeric.Sundials.ARKode.ODE
+import qualified Numeric.Sundials.ARKode.ODE as ARK
+import qualified Numeric.Sundials.CVode.ODE  as CV
 import           Numeric.LinearAlgebra
 
 import           Plots as P
@@ -97,24 +98,24 @@ kSaxis xs = P.r2Axis &~ do
 main :: IO ()
 main = do
 
-  let res1 = odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0])
+  let res1 = ARK.odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0])
   renderRasterific "diagrams/brusselator.png"
                    (D.dims2D 500.0 500.0)
                    (renderAxis $ lSaxis $ [0.0, 0.1 .. 10.0]:(toLists $ tr res1))
 
-  let res1a = odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0])
+  let res1a = ARK.odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0])
   renderRasterific "diagrams/brusselatorA.png"
                    (D.dims2D 500.0 500.0)
                    (renderAxis $ lSaxis $ [0.0, 0.1 .. 10.0]:(toLists $ tr res1a))
 
-  let res2 = odeSolve stiffish [0.0] (fromList [0.0, 0.1 .. 10.0])
+  let res2 = ARK.odeSolve stiffish [0.0] (fromList [0.0, 0.1 .. 10.0])
   renderRasterific "diagrams/stiffish.png"
                    (D.dims2D 500.0 500.0)
                    (renderAxis $ kSaxis $ zip [0.0, 0.1 .. 10.0] (concat $ toLists res2))
 
-  let res2a = odeSolveV (SDIRK_5_3_4') Nothing 1e-6 1e-10 stiffishV (fromList [0.0]) (fromList [0.0, 0.1 .. 10.0])
+  let res2a = ARK.odeSolveV (ARK.SDIRK_5_3_4') Nothing 1e-6 1e-10 stiffishV (fromList [0.0]) (fromList [0.0, 0.1 .. 10.0])
 
-  let res2b = odeSolveV (TRBDF2_3_3_2') Nothing 1e-6 1e-10 stiffishV (fromList [0.0]) (fromList [0.0, 0.1 .. 10.0])
+  let res2b = ARK.odeSolveV (ARK.TRBDF2_3_3_2') Nothing 1e-6 1e-10 stiffishV (fromList [0.0]) (fromList [0.0, 0.1 .. 10.0])
 
   let maxDiff = maximum $ map abs $
                 zipWith (-) ((toLists $ tr res2a)!!0) ((toLists $ tr res2b)!!0)
@@ -123,7 +124,10 @@ main = do
     it "for two different RK methods" $
       maxDiff < 1.0e-6
 
-  let res3 = odeSolve lorenz [-5.0, -5.0, 1.0] (fromList [0.0, 0.01 .. 10.0])
+  let res2c = CV.odeSolveV (CV.BDF) Nothing 1e-6 1e-10 stiffishV (fromList [0.0]) (fromList [0.0, 0.1 .. 10.0])
+  putStrLn $ show res2c
+
+  let res3 = ARK.odeSolve lorenz [-5.0, -5.0, 1.0] (fromList [0.0, 0.01 .. 10.0])
 
   renderRasterific "diagrams/lorenz.png"
                    (D.dims2D 500.0 500.0)