summaryrefslogtreecommitdiff
path: root/packages/sundials
diff options
context:
space:
mode:
Diffstat (limited to 'packages/sundials')
-rw-r--r--packages/sundials/LICENSE2
-rw-r--r--packages/sundials/hmatrix-sundials.cabal14
-rw-r--r--packages/sundials/src/Main.hs34
-rw-r--r--packages/sundials/src/Numeric/Sundials/Arkode/ODE.hs (renamed from packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs)31
-rw-r--r--packages/sundials/stack.yaml66
5 files changed, 37 insertions, 110 deletions
diff --git a/packages/sundials/LICENSE b/packages/sundials/LICENSE
index 18e5956..e0cc65d 100644
--- a/packages/sundials/LICENSE
+++ b/packages/sundials/LICENSE
@@ -1,4 +1,4 @@
1Copyright (c) 2018, Dominic Steinitz 1Copyright (c) 2018 Novadiscovery, Dominic Steinitz
2 2
3All rights reserved. 3All rights reserved.
4 4
diff --git a/packages/sundials/hmatrix-sundials.cabal b/packages/sundials/hmatrix-sundials.cabal
index 762537e..1ddede3 100644
--- a/packages/sundials/hmatrix-sundials.cabal
+++ b/packages/sundials/hmatrix-sundials.cabal
@@ -1,23 +1,18 @@
1-- Initial sundials.cabal generated by cabal init. For further
2-- documentation, see http://haskell.org/cabal/users-guide/
3
4name: hmatrix-sundials 1name: hmatrix-sundials
5version: 0.1.0.0 2version: 0.1.0.0
6-- synopsis: 3homepage: https://github.com/albertoruiz/hmatrix/tree/sundials
7-- description:
8homepage: https://github.com/idontgetoutmuch/hmatrix/tree/sundials
9license: BSD3 4license: BSD3
10license-file: LICENSE 5license-file: LICENSE
11author: Dominic Steinitz 6author: Dominic Steinitz
12maintainer: dominic@steinitz.org 7maintainer: dominic@steinitz.org
13-- copyright:
14category: Math 8category: Math
15build-type: Simple 9build-type: Simple
16extra-source-files: ChangeLog.md 10extra-source-files: ChangeLog.md
17cabal-version: >=1.10 11cabal-version: >=1.10
18 12
19extra-source-files: src/helpers.c, src/helpers.h 13extra-source-files:
20 14 src/helpers.c
15 src/helpers.h
21 16
22library 17library
23 build-depends: base >=4.10 && <4.11, 18 build-depends: base >=4.10 && <4.11,
@@ -41,6 +36,7 @@ executable sundials
41 template-haskell >=2.12 && <2.13, 36 template-haskell >=2.12 && <2.13,
42 containers >=0.5 && <0.6, 37 containers >=0.5 && <0.6,
43 hmatrix>=0.18, 38 hmatrix>=0.18,
39 lens,
44 plots, 40 plots,
45 diagrams-lib, 41 diagrams-lib,
46 diagrams-rasterific 42 diagrams-rasterific
diff --git a/packages/sundials/src/Main.hs b/packages/sundials/src/Main.hs
index 5e51372..d1f6755 100644
--- a/packages/sundials/src/Main.hs
+++ b/packages/sundials/src/Main.hs
@@ -1,6 +1,5 @@
1{-# OPTIONS_GHC -Wall #-} 1module Main where
2 2
3import qualified Data.Vector.Storable as V
4import Numeric.Sundials.Arkode.ODE 3import Numeric.Sundials.Arkode.ODE
5import Numeric.LinearAlgebra 4import Numeric.LinearAlgebra
6 5
@@ -9,11 +8,10 @@ import qualified Diagrams.Prelude as D
9import Diagrams.Backend.Rasterific 8import Diagrams.Backend.Rasterific
10 9
11import Control.Lens 10import Control.Lens
12import Data.List (zip4)
13 11
14 12brusselator :: a -> [Double] -> [Double]
15brusselator _t x = [ a - (w + 1) * u + v * u^2 13brusselator _t x = [ a - (w + 1) * u + v * u^(2::Int)
16 , w * u - v * u^2 14 , w * u - v * u^(2::Int)
17 , (b - w) / eps - w * u 15 , (b - w) / eps - w * u
18 ] 16 ]
19 where 17 where
@@ -24,6 +22,7 @@ brusselator _t x = [ a - (w + 1) * u + v * u^2
24 v = x !! 1 22 v = x !! 1
25 w = x !! 2 23 w = x !! 2
26 24
25stiffish :: Double -> [Double] -> [Double]
27stiffish t v = [ lamda * u + 1.0 / (1.0 + t * t) - lamda * atan t ] 26stiffish t v = [ lamda * u + 1.0 / (1.0 + t * t) - lamda * atan t ]
28 where 27 where
29 lamda = -100.0 28 lamda = -100.0
@@ -45,15 +44,14 @@ kSaxis xs = P.r2Axis &~ do
45 44
46main :: IO () 45main :: IO ()
47main = do 46main = do
48 let res = odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0]) 47 do let res = odeSolve brusselator [1.2, 3.1, 3.0] (fromList [0.0, 0.1 .. 10.0])
49 putStrLn $ show res 48 putStrLn $ show res
50 renderRasterific "diagrams/brusselator.png" 49 renderRasterific "diagrams/brusselator.png"
51 (D.dims2D 500.0 500.0) 50 (D.dims2D 500.0 500.0)
52 (renderAxis $ lSaxis $ [0.0, 0.1 .. 10.0]:(toLists $ tr res)) 51 (renderAxis $ lSaxis $ [0.0, 0.1 .. 10.0]:(toLists $ tr res))
53 52
54 let res = odeSolve stiffish [0.0] (fromList [0.0, 0.1 .. 10.0]) 53 do let res = odeSolve stiffish [0.0] (fromList [0.0, 0.1 .. 10.0])
55 putStrLn $ show res 54 putStrLn $ show res
56 renderRasterific "diagrams/stiffish.png" 55 renderRasterific "diagrams/stiffish.png"
57 (D.dims2D 500.0 500.0) 56 (D.dims2D 500.0 500.0)
58 (renderAxis $ kSaxis $ zip [0.0, 0.1 .. 10.0] (concat $ toLists res)) 57 (renderAxis $ kSaxis $ zip [0.0, 0.1 .. 10.0] (concat $ toLists res))
59
diff --git a/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs b/packages/sundials/src/Numeric/Sundials/Arkode/ODE.hs
index f432951..6d9a1b2 100644
--- a/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs
+++ b/packages/sundials/src/Numeric/Sundials/Arkode/ODE.hs
@@ -1,14 +1,14 @@
1{-# OPTIONS_GHC -Wall #-}
2
3{-# LANGUAGE QuasiQuotes #-}
4{-# LANGUAGE TemplateHaskell #-}
5{-# LANGUAGE MultiWayIf #-} 1{-# LANGUAGE MultiWayIf #-}
6{-# LANGUAGE OverloadedStrings #-} 2{-# LANGUAGE OverloadedStrings #-}
3{-# LANGUAGE QuasiQuotes #-}
7{-# LANGUAGE ScopedTypeVariables #-} 4{-# LANGUAGE ScopedTypeVariables #-}
5{-# LANGUAGE TemplateHaskell #-}
8 6
9module Numeric.Sundials.Arkode.ODE ( solveOde 7module Numeric.Sundials.Arkode.ODE
10 , odeSolve 8 ( SundialsDiagnostics(..)
11 ) where 9 , solveOde
10 , odeSolve
11 ) where
12 12
13import qualified Language.C.Inline as C 13import qualified Language.C.Inline as C
14import qualified Language.C.Inline.Unsafe as CU 14import qualified Language.C.Inline.Unsafe as CU
@@ -40,10 +40,10 @@ C.include "<stdio.h>"
40C.include "<math.h>" 40C.include "<math.h>"
41C.include "<arkode/arkode.h>" -- prototypes for ARKODE fcts., consts. 41C.include "<arkode/arkode.h>" -- prototypes for ARKODE fcts., consts.
42C.include "<nvector/nvector_serial.h>" -- serial N_Vector types, fcts., macros 42C.include "<nvector/nvector_serial.h>" -- serial N_Vector types, fcts., macros
43C.include "<sunmatrix/sunmatrix_dense.h>" -- access to dense SUNMatrix 43C.include "<sunmatrix/sunmatrix_dense.h>" -- access to dense SUNMatrix
44C.include "<sunlinsol/sunlinsol_dense.h>" -- access to dense SUNLinearSolver 44C.include "<sunlinsol/sunlinsol_dense.h>" -- access to dense SUNLinearSolver
45C.include "<arkode/arkode_direct.h>" -- access to ARKDls interface 45C.include "<arkode/arkode_direct.h>" -- access to ARKDls interface
46C.include "<sundials/sundials_types.h>" -- definition of type realtype 46C.include "<sundials/sundials_types.h>" -- definition of type realtype
47C.include "<sundials/sundials_math.h>" 47C.include "<sundials/sundials_math.h>"
48C.include "../../../helpers.h" 48C.include "../../../helpers.h"
49 49
@@ -108,7 +108,7 @@ odeSolve f y0 ts = case solveOde g (V.fromList y0) (V.fromList $ toList ts) of
108 nC = length y0 108 nC = length y0
109 g t x0 = V.fromList $ f t (V.toList x0) 109 g t x0 = V.fromList $ f t (V.toList x0)
110 110
111solveOde :: 111solveOde ::
112 (Double -> V.Vector Double -> V.Vector Double) -- ^ The RHS of the system \(\dot{y} = f(t,y)\) 112 (Double -> V.Vector Double -> V.Vector Double) -- ^ The RHS of the system \(\dot{y} = f(t,y)\)
113 -> V.Vector Double -- ^ Initial conditions 113 -> V.Vector Double -- ^ Initial conditions
114 -> V.Vector Double -- ^ Desired solution times 114 -> V.Vector Double -- ^ Desired solution times
@@ -116,7 +116,7 @@ solveOde ::
116solveOde f y0 tt = case solveOdeC (coerce f) (coerce y0) (coerce tt) of 116solveOde f y0 tt = case solveOdeC (coerce f) (coerce y0) (coerce tt) of
117 Left c -> Left $ fromIntegral c 117 Left c -> Left $ fromIntegral c
118 Right (v, d) -> Right (coerce v, d) 118 Right (v, d) -> Right (coerce v, d)
119 119
120solveOdeC :: 120solveOdeC ::
121 (CDouble -> V.Vector CDouble -> V.Vector CDouble) -- ^ The RHS of the system \(\dot{y} = f(t,y)\) 121 (CDouble -> V.Vector CDouble -> V.Vector CDouble) -- ^ The RHS of the system \(\dot{y} = f(t,y)\)
122 -> V.Vector CDouble -- ^ Initial conditions 122 -> V.Vector CDouble -- ^ Initial conditions
@@ -285,13 +285,13 @@ solveOdeC fun f0 ts = unsafePerformIO $ do
285 flag = ARKDlsGetNumRhsEvals(arkode_mem, &nfeLS); 285 flag = ARKDlsGetNumRhsEvals(arkode_mem, &nfeLS);
286 check_flag(&flag, "ARKDlsGetNumRhsEvals", 1); 286 check_flag(&flag, "ARKDlsGetNumRhsEvals", 1);
287 ($vec-ptr:(long int *diagMut))[9] = ncfn; 287 ($vec-ptr:(long int *diagMut))[9] = ncfn;
288 288
289 /* Clean up and return */ 289 /* Clean up and return */
290 N_VDestroy(y); /* Free y vector */ 290 N_VDestroy(y); /* Free y vector */
291 ARKodeFree(&arkode_mem); /* Free integrator memory */ 291 ARKodeFree(&arkode_mem); /* Free integrator memory */
292 SUNLinSolFree(LS); /* Free linear solver */ 292 SUNLinSolFree(LS); /* Free linear solver */
293 SUNMatDestroy(A); /* Free A matrix */ 293 SUNMatDestroy(A); /* Free A matrix */
294 294
295 return flag; 295 return flag;
296 } |] 296 } |]
297 if res == 0 297 if res == 0
@@ -311,4 +311,3 @@ solveOdeC fun f0 ts = unsafePerformIO $ do
311 return $ Right (m, d) 311 return $ Right (m, d)
312 else do 312 else do
313 return $ Left res 313 return $ Left res
314
diff --git a/packages/sundials/stack.yaml b/packages/sundials/stack.yaml
deleted file mode 100644
index 9c6b17c..0000000
--- a/packages/sundials/stack.yaml
+++ /dev/null
@@ -1,66 +0,0 @@
1# This file was automatically generated by 'stack init'
2#
3# Some commonly used options have been documented as comments in this file.
4# For advanced use and comprehensive documentation of the format, please see:
5# https://docs.haskellstack.org/en/stable/yaml_configuration/
6
7# Resolver to choose a 'specific' stackage snapshot or a compiler version.
8# A snapshot resolver dictates the compiler version and the set of packages
9# to be used for project dependencies. For example:
10#
11# resolver: lts-3.5
12# resolver: nightly-2015-09-21
13# resolver: ghc-7.10.2
14# resolver: ghcjs-0.1.0_ghc-7.10.2
15# resolver:
16# name: custom-snapshot
17# location: "./custom-snapshot.yaml"
18resolver: lts-10.9
19
20# User packages to be built.
21# Various formats can be used as shown in the example below.
22#
23# packages:
24# - some-directory
25# - https://example.com/foo/bar/baz-0.0.2.tar.gz
26# - location:
27# git: https://github.com/commercialhaskell/stack.git
28# commit: e7b331f14bcffb8367cd58fbfc8b40ec7642100a
29# - location: https://github.com/commercialhaskell/stack/commit/e7b331f14bcffb8367cd58fbfc8b40ec7642100a
30# extra-dep: true
31# subdirs:
32# - auto-update
33# - wai
34#
35# A package marked 'extra-dep: true' will only be built if demanded by a
36# non-dependency (i.e. a user package), and its test suites and benchmarks
37# will not be run. This is useful for tweaking upstream packages.
38packages:
39- .
40# Dependency packages to be pulled from upstream that are not in the resolver
41# (e.g., acme-missiles-0.3)
42# extra-deps: []
43
44# Override default flag values for local packages and extra-deps
45# flags: {}
46
47# Extra package databases containing global packages
48# extra-package-dbs: []
49
50# Control whether we use the GHC we find on the path
51# system-ghc: true
52#
53# Require a specific version of stack, using version ranges
54# require-stack-version: -any # Default
55# require-stack-version: ">=1.6"
56#
57# Override the architecture used by stack, especially useful on Windows
58# arch: i386
59# arch: x86_64
60#
61# Extra directories used by stack for building
62# extra-include-dirs: [/path/to/dir]
63# extra-lib-dirs: [/path/to/dir]
64#
65# Allow a newer minor version of GHC than the snapshot specifies
66# compiler-check: newer-minor \ No newline at end of file