diff options
author | Dominic Steinitz <dominic@steinitz.org> | 2018-04-27 14:19:59 +0100 |
---|---|---|
committer | Dominic Steinitz <dominic@steinitz.org> | 2018-04-27 14:19:59 +0100 |
commit | 149dedfc6ec8dea039a4df7ad1d31880820c52eb (patch) | |
tree | 6f5ead50e4006d3ea578d8e744849470004cdf78 /packages/sundials/src/Numeric/Sundials/ARKode | |
parent | d48892298d5e87ec12b29adc69af2fdd59f4c8bd (diff) |
More restructuring
Diffstat (limited to 'packages/sundials/src/Numeric/Sundials/ARKode')
-rw-r--r-- | packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs | 56 |
1 files changed, 40 insertions, 16 deletions
diff --git a/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs b/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs index a8d418b..13b7eb8 100644 --- a/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs +++ b/packages/sundials/src/Numeric/Sundials/ARKode/ODE.hs | |||
@@ -125,6 +125,7 @@ import Data.Maybe (isJust) | |||
125 | 125 | ||
126 | import Foreign.C.Types (CDouble, CInt, CLong) | 126 | import Foreign.C.Types (CDouble, CInt, CLong) |
127 | import Foreign.Ptr (Ptr) | 127 | import Foreign.Ptr (Ptr) |
128 | import Foreign.Storable (poke) | ||
128 | 129 | ||
129 | import qualified Data.Vector.Storable as V | 130 | import qualified Data.Vector.Storable as V |
130 | 131 | ||
@@ -139,10 +140,33 @@ import Numeric.LinearAlgebra.HMatrix (Vector, Matrix, toList, (><), | |||
139 | subMatrix, rows, cols, toLists, | 140 | subMatrix, rows, cols, toLists, |
140 | size, subVector) | 141 | size, subVector) |
141 | 142 | ||
142 | import qualified Numeric.Sundials.CLangToHaskellTypes as T | ||
143 | import Numeric.Sundials.Arkode | ||
144 | import qualified Numeric.Sundials.Arkode as B | ||
145 | import qualified Numeric.Sundials.ODEOpts as SO | 143 | import qualified Numeric.Sundials.ODEOpts as SO |
144 | import qualified Numeric.Sundials.Arkode as T | ||
145 | import Numeric.Sundials.Arkode (getDataFromContents, putDataInContents, arkSMax, | ||
146 | sDIRK_2_1_2, | ||
147 | bILLINGTON_3_3_2, | ||
148 | tRBDF2_3_3_2, | ||
149 | kVAERNO_4_2_3, | ||
150 | aRK324L2SA_DIRK_4_2_3, | ||
151 | cASH_5_2_4, | ||
152 | cASH_5_3_4, | ||
153 | sDIRK_5_3_4, | ||
154 | kVAERNO_5_3_4, | ||
155 | aRK436L2SA_DIRK_6_3_4, | ||
156 | kVAERNO_7_4_5, | ||
157 | aRK548L2SA_DIRK_8_4_5, | ||
158 | hEUN_EULER_2_1_2, | ||
159 | bOGACKI_SHAMPINE_4_2_3, | ||
160 | aRK324L2SA_ERK_4_2_3, | ||
161 | zONNEVELD_5_3_4, | ||
162 | aRK436L2SA_ERK_6_3_4, | ||
163 | sAYFY_ABURUB_6_3_4, | ||
164 | cASH_KARP_6_4_5, | ||
165 | fEHLBERG_6_4_5, | ||
166 | dORMAND_PRINCE_7_4_5, | ||
167 | aRK548L2SA_ERK_8_4_5, | ||
168 | vERNER_8_5_6, | ||
169 | fEHLBERG_13_7_8) | ||
146 | 170 | ||
147 | 171 | ||
148 | C.context (C.baseCtx <> C.vecCtx <> C.funCtx <> T.sunCtx) | 172 | C.context (C.baseCtx <> C.vecCtx <> C.funCtx <> T.sunCtx) |
@@ -451,9 +475,9 @@ solveOdeC method initStepSize jacH (absTols, relTol) fun f0 ts = unsafePerformIO | |||
451 | funIO x y f _ptr = do | 475 | funIO x y f _ptr = do |
452 | -- Convert the pointer we get from C (y) to a vector, and then | 476 | -- Convert the pointer we get from C (y) to a vector, and then |
453 | -- apply the user-supplied function. | 477 | -- apply the user-supplied function. |
454 | fImm <- fun x <$> SO.getDataFromContents dim y | 478 | fImm <- fun x <$> getDataFromContents dim y |
455 | -- Fill in the provided pointer with the resulting vector. | 479 | -- Fill in the provided pointer with the resulting vector. |
456 | SO.putDataInContents fImm dim f | 480 | putDataInContents fImm dim f |
457 | -- FIXME: I don't understand what this comment means | 481 | -- FIXME: I don't understand what this comment means |
458 | -- Unsafe since the function will be called many times. | 482 | -- Unsafe since the function will be called many times. |
459 | [CU.exp| int{ 0 } |] | 483 | [CU.exp| int{ 0 } |] |
@@ -465,8 +489,8 @@ solveOdeC method initStepSize jacH (absTols, relTol) fun f0 ts = unsafePerformIO | |||
465 | jacIO t y _fy jacS _ptr _tmp1 _tmp2 _tmp3 = do | 489 | jacIO t y _fy jacS _ptr _tmp1 _tmp2 _tmp3 = do |
466 | case jacH of | 490 | case jacH of |
467 | Nothing -> error "Numeric.Sundials.ARKode.ODE: Jacobian not defined" | 491 | Nothing -> error "Numeric.Sundials.ARKode.ODE: Jacobian not defined" |
468 | Just jacI -> do j <- jacI t <$> SO.getDataFromContents dim y | 492 | Just jacI -> do j <- jacI t <$> getDataFromContents dim y |
469 | SO.putMatrixDataFromContents j jacS | 493 | poke jacS j |
470 | -- FIXME: I don't understand what this comment means | 494 | -- FIXME: I don't understand what this comment means |
471 | -- Unsafe since the function will be called many times. | 495 | -- Unsafe since the function will be called many times. |
472 | [CU.exp| int{ 0 } |] | 496 | [CU.exp| int{ 0 } |] |
@@ -675,7 +699,7 @@ butcherTable method = | |||
675 | case getBT method of | 699 | case getBT method of |
676 | Left c -> error $ show c -- FIXME | 700 | Left c -> error $ show c -- FIXME |
677 | Right (ButcherTable' v w x y, sqp) -> | 701 | Right (ButcherTable' v w x y, sqp) -> |
678 | ButcherTable { am = subMatrix (0, 0) (s, s) $ (B.arkSMax >< B.arkSMax) (V.toList v) | 702 | ButcherTable { am = subMatrix (0, 0) (s, s) $ (arkSMax >< arkSMax) (V.toList v) |
679 | , cv = subVector 0 s w | 703 | , cv = subVector 0 s w |
680 | , bv = subVector 0 s x | 704 | , bv = subVector 0 s x |
681 | , b2v = subVector 0 s y | 705 | , b2v = subVector 0 s y |
@@ -710,25 +734,25 @@ getButcherTable method = unsafePerformIO $ do | |||
710 | 734 | ||
711 | btSQP :: V.Vector CInt <- createVector 3 | 735 | btSQP :: V.Vector CInt <- createVector 3 |
712 | btSQPMut <- V.thaw btSQP | 736 | btSQPMut <- V.thaw btSQP |
713 | btAs :: V.Vector CDouble <- createVector (B.arkSMax * B.arkSMax) | 737 | btAs :: V.Vector CDouble <- createVector (arkSMax * arkSMax) |
714 | btAsMut <- V.thaw btAs | 738 | btAsMut <- V.thaw btAs |
715 | btCs :: V.Vector CDouble <- createVector B.arkSMax | 739 | btCs :: V.Vector CDouble <- createVector arkSMax |
716 | btBs :: V.Vector CDouble <- createVector B.arkSMax | 740 | btBs :: V.Vector CDouble <- createVector arkSMax |
717 | btB2s :: V.Vector CDouble <- createVector B.arkSMax | 741 | btB2s :: V.Vector CDouble <- createVector arkSMax |
718 | btCsMut <- V.thaw btCs | 742 | btCsMut <- V.thaw btCs |
719 | btBsMut <- V.thaw btBs | 743 | btBsMut <- V.thaw btBs |
720 | btB2sMut <- V.thaw btB2s | 744 | btB2sMut <- V.thaw btB2s |
721 | let funIOI :: CDouble -> Ptr T.SunVector -> Ptr T.SunVector -> Ptr () -> IO CInt | 745 | let funIOI :: CDouble -> Ptr T.SunVector -> Ptr T.SunVector -> Ptr () -> IO CInt |
722 | funIOI x y f _ptr = do | 746 | funIOI x y f _ptr = do |
723 | fImm <- funI x <$> SO.getDataFromContents dim y | 747 | fImm <- funI x <$> getDataFromContents dim y |
724 | SO.putDataInContents fImm dim f | 748 | putDataInContents fImm dim f |
725 | -- FIXME: I don't understand what this comment means | 749 | -- FIXME: I don't understand what this comment means |
726 | -- Unsafe since the function will be called many times. | 750 | -- Unsafe since the function will be called many times. |
727 | [CU.exp| int{ 0 } |] | 751 | [CU.exp| int{ 0 } |] |
728 | let funIOE :: CDouble -> Ptr T.SunVector -> Ptr T.SunVector -> Ptr () -> IO CInt | 752 | let funIOE :: CDouble -> Ptr T.SunVector -> Ptr T.SunVector -> Ptr () -> IO CInt |
729 | funIOE x y f _ptr = do | 753 | funIOE x y f _ptr = do |
730 | fImm <- funE x <$> SO.getDataFromContents dim y | 754 | fImm <- funE x <$> getDataFromContents dim y |
731 | SO.putDataInContents fImm dim f | 755 | putDataInContents fImm dim f |
732 | -- FIXME: I don't understand what this comment means | 756 | -- FIXME: I don't understand what this comment means |
733 | -- Unsafe since the function will be called many times. | 757 | -- Unsafe since the function will be called many times. |
734 | [CU.exp| int{ 0 } |] | 758 | [CU.exp| int{ 0 } |] |