6 files changed, 169 insertions, 0 deletions

diff --git a/examples/ButcherTableau.hs b/examples/ButcherTableau.hs
new file mode 100644
index 0000000..e667ef1
--- /dev/null
+++ b/examples/ButcherTableau.hs
@@ -0,0 +1,47 @@
+{-# OPTIONS_GHC -Wall #-}
+
+import           Numeric.Sundials.ARKode.ODE
+import           Numeric.LinearAlgebra
+
+import           Data.List (intercalate)
+
+import           Text.PrettyPrint.HughesPJClass
+
+
+butcherTableauTex :: ButcherTable -> String
+butcherTableauTex (ButcherTable m c b b2) =
+  render $
+  vcat [ text ("\n\\begin{array}{c|" ++ (concat $ replicate n "c") ++ "}")
+       , us
+       , text "\\hline"
+       , text bs <+> text "\\\\"
+       , text b2s <+> text "\\\\"
+       , text "\\end{array}"
+       ]
+  where
+    n = rows m
+    rs = toLists m
+    ss = map (\r -> intercalate " & " $ map show r) rs
+    ts = zipWith (\i r -> show i ++ " & " ++ r) (toList c) ss
+    us = vcat $ map (\r -> text r <+> text "\\\\") ts
+    bs  = " & " ++ (intercalate " & " $ map show $ toList b)
+    b2s = " & " ++ (intercalate " & " $ map show $ toList b2)
+
+main :: IO ()
+main = do
+
+  let res = butcherTable (SDIRK_2_1_2 undefined)
+  putStrLn $ show res
+  putStrLn $ butcherTableauTex res
+
+  let resA = butcherTable (KVAERNO_4_2_3 undefined)
+  putStrLn $ show resA
+  putStrLn $ butcherTableauTex resA
+
+  let resB = butcherTable (SDIRK_5_3_4 undefined)
+  putStrLn $ show resB
+  putStrLn $ butcherTableauTex resB
+
+  let resC = butcherTable (FEHLBERG_6_4_5 undefined)
+  putStrLn $ show resC
+  putStrLn $ butcherTableauTex resC
diff --git a/examples/ChangeLog.md b/examples/ChangeLog.md
new file mode 100644
index 0000000..a7e4c4c
--- /dev/null
+++ b/examples/ChangeLog.md
@@ -0,0 +1,5 @@
+# Revision history for hmatrix-examples
+
+## 0.1.0.0  -- YYYY-mm-dd
+
+* First version. Released on an unsuspecting world.
diff --git a/examples/LICENSE b/examples/LICENSE
new file mode 100644
index 0000000..18e5956
--- /dev/null
+++ b/examples/LICENSE
@@ -0,0 +1,30 @@
+Copyright (c) 2018, Dominic Steinitz
+
+All rights reserved.
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are met:
+
+    * Redistributions of source code must retain the above copyright
+      notice, this list of conditions and the following disclaimer.
+
+    * Redistributions in binary form must reproduce the above
+      copyright notice, this list of conditions and the following
+      disclaimer in the documentation and/or other materials provided
+      with the distribution.
+
+    * Neither the name of Dominic Steinitz nor the names of other
+      contributors may be used to endorse or promote products derived
+      from this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
+"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
+LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
+A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
+OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
+SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
+LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
+DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
+THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
+(INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
diff --git a/examples/Setup.hs b/examples/Setup.hs
new file mode 100644
index 0000000..9a994af
--- /dev/null
+++ b/examples/Setup.hs
@@ -0,0 +1,2 @@
+import Distribution.Simple
+main = defaultMain
diff --git a/examples/examples.cabal b/examples/examples.cabal
new file mode 100644
index 0000000..5e85a3d
--- /dev/null
+++ b/examples/examples.cabal
@@ -0,0 +1,29 @@
+name:                examples
+version:             0.1.0.0
+synopsis:            Example usage of the various hmatrix packages
+homepage:            https://github.com/albertoruiz/hmatrix
+license:             BSD3
+license-file:        LICENSE
+author:              Dominic Steinitz
+maintainer:          dominic@steinitz.org
+copyright:           Dominic Steinitz, Novadiscovery
+category:            Math
+build-type:          Simple
+extra-source-files:  ChangeLog.md
+cabal-version:       >=1.10
+
+executable sundials
+  main-is:             sundials.hs
+  build-depends:       base >=4.10 && <4.11,
+                       hmatrix,
+                       hmatrix-sundials,
+                       hmatrix-gsl
+  default-language:    Haskell2010
+
+executable butcherTableau
+  main-is:             ButcherTableau.hs
+  build-depends:       base >=4.10 && <4.11,
+                       hmatrix,
+                       hmatrix-sundials,
+                       pretty
+  default-language:    Haskell2010
diff --git a/examples/sundials.hs b/examples/sundials.hs
new file mode 100644
index 0000000..99f662d
--- /dev/null
+++ b/examples/sundials.hs
@@ -0,0 +1,56 @@
+{-# OPTIONS_GHC -fno-warn-missing-signatures #-}
+
+{-# LANGUAGE ViewPatterns #-}
+
+import           Numeric.Sundials.ARKode.ODE
+import           Numeric.LinearAlgebra
+import           Graphics.Plot
+
+vanderpol mu = do
+    let xdot nu _t [x,v] = [v, -x + nu * v * (1-x*x)]
+        xdot _ _ _ = error "vanderpol RHS not defined"
+        ts = linspace 1000 (0,50)
+        sol = toColumns $ odeSolve (xdot mu) [1,0] ts
+    mplot (ts : sol)
+    mplot sol
+
+
+harmonic w d = do
+    let xdot u dd _t [x,v] = [v, a*x + b*v] where a = -u*u; b = -2*dd*u
+        xdot _ _ _ _ = error "harmonic RHS not defined"
+        ts = linspace 100 (0,20)
+        sol = odeSolve (xdot w d) [1,0] ts
+    mplot (ts : toColumns sol)
+
+
+kepler v a = mplot (take 2 $ toColumns sol) where
+    xdot _t [x,y,vx,vy] = [vx,vy,x*k,y*k]
+        where g=1
+              k=(-g)*(x*x+y*y)**(-1.5)
+    xdot _ _ = error "kepler RHS not defined"
+    ts = linspace 100 (0,30)
+    sol = odeSolve xdot [4, 0, v * cos (a*degree), v * sin (a*degree)] ts
+    degree = pi/180
+
+
+main = do
+    vanderpol 2
+    harmonic 1 0
+    harmonic 1 0.1
+    kepler 0.3 60
+    kepler 0.4 70
+    vanderpol' 2
+
+-- example of odeSolveV with jacobian
+vanderpol' mu = do
+    let xdot nu _t (toList->[x,v]) = fromList [v, -x + nu * v * (1-x*x)]
+        xdot _ _ _ = error "vanderpol' RHS not defined"
+        jac _ (toList->[x,v]) = (2><2) [ 0          ,          1
+                                       , -1-2*x*v*mu, mu*(1-x**2) ]
+        jac _ _ = error "vanderpol' Jacobian not defined"
+        ts = linspace 1000 (0,50)
+        hi = pure $ (ts!1 - ts!0) / 100.0
+        sol = toColumns $ odeSolveV (SDIRK_5_3_4 jac) hi 1E-8 1E-8 (xdot mu) (fromList [1,0]) ts
+    mplot sol
+
+