diff options
Diffstat (limited to 'packages/gsl/src/Numeric')
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Fitting.hs | 16 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Fourier.hs | 5 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/IO.hs | 2 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Internal.hs | 7 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/LinearAlgebra.hs | 2 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Minimization.hs | 17 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/ODE.hs | 16 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Polynomials.hs | 7 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Random.hs | 16 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Root.hs | 18 | ||||
-rw-r--r-- | packages/gsl/src/Numeric/GSL/Vector.hs | 3 |
11 files changed, 59 insertions, 50 deletions
diff --git a/packages/gsl/src/Numeric/GSL/Fitting.hs b/packages/gsl/src/Numeric/GSL/Fitting.hs index 0a92373..db9d82f 100644 --- a/packages/gsl/src/Numeric/GSL/Fitting.hs +++ b/packages/gsl/src/Numeric/GSL/Fitting.hs | |||
@@ -1,3 +1,5 @@ | |||
1 | {-# LANGUAGE FlexibleContexts #-} | ||
2 | |||
1 | {- | | 3 | {- | |
2 | Module : Numeric.GSL.Fitting | 4 | Module : Numeric.GSL.Fitting |
3 | Copyright : (c) Alberto Ruiz 2010 | 5 | Copyright : (c) Alberto Ruiz 2010 |
@@ -50,7 +52,7 @@ module Numeric.GSL.Fitting ( | |||
50 | fitModelScaled, fitModel | 52 | fitModelScaled, fitModel |
51 | ) where | 53 | ) where |
52 | 54 | ||
53 | import Numeric.LinearAlgebra | 55 | import Numeric.LinearAlgebra.HMatrix |
54 | import Numeric.GSL.Internal | 56 | import Numeric.GSL.Internal |
55 | 57 | ||
56 | import Foreign.Ptr(FunPtr, freeHaskellFunPtr) | 58 | import Foreign.Ptr(FunPtr, freeHaskellFunPtr) |
@@ -80,13 +82,13 @@ nlFitting :: FittingMethod | |||
80 | nlFitting method epsabs epsrel maxit fun jac xinit = nlFitGen (fi (fromEnum method)) fun jac xinit epsabs epsrel maxit | 82 | nlFitting method epsabs epsrel maxit fun jac xinit = nlFitGen (fi (fromEnum method)) fun jac xinit epsabs epsrel maxit |
81 | 83 | ||
82 | nlFitGen m f jac xiv epsabs epsrel maxit = unsafePerformIO $ do | 84 | nlFitGen m f jac xiv epsabs epsrel maxit = unsafePerformIO $ do |
83 | let p = dim xiv | 85 | let p = size xiv |
84 | n = dim (f xiv) | 86 | n = size (f xiv) |
85 | fp <- mkVecVecfun (aux_vTov (checkdim1 n p . f)) | 87 | fp <- mkVecVecfun (aux_vTov (checkdim1 n p . f)) |
86 | jp <- mkVecMatfun (aux_vTom (checkdim2 n p . jac)) | 88 | jp <- mkVecMatfun (aux_vTom (checkdim2 n p . jac)) |
87 | rawpath <- createMatrix RowMajor maxit (2+p) | 89 | rawpath <- createMatrix RowMajor maxit (2+p) |
88 | app2 (c_nlfit m fp jp epsabs epsrel (fi maxit) (fi n)) vec xiv mat rawpath "c_nlfit" | 90 | app2 (c_nlfit m fp jp epsabs epsrel (fi maxit) (fi n)) vec xiv mat rawpath "c_nlfit" |
89 | let it = round (rawpath @@> (maxit-1,0)) | 91 | let it = round (rawpath `atIndex` (maxit-1,0)) |
90 | path = takeRows it rawpath | 92 | path = takeRows it rawpath |
91 | [sol] = toRows $ dropRows (it-1) path | 93 | [sol] = toRows $ dropRows (it-1) path |
92 | freeHaskellFunPtr fp | 94 | freeHaskellFunPtr fp |
@@ -99,7 +101,7 @@ foreign import ccall safe "nlfit" | |||
99 | ------------------------------------------------------- | 101 | ------------------------------------------------------- |
100 | 102 | ||
101 | checkdim1 n _p v | 103 | checkdim1 n _p v |
102 | | dim v == n = v | 104 | | size v == n = v |
103 | | otherwise = error $ "Error: "++ show n | 105 | | otherwise = error $ "Error: "++ show n |
104 | ++ " components expected in the result of the function supplied to nlFitting" | 106 | ++ " components expected in the result of the function supplied to nlFitting" |
105 | 107 | ||
@@ -114,9 +116,9 @@ err (model,deriv) dat vsol = zip sol errs where | |||
114 | sol = toList vsol | 116 | sol = toList vsol |
115 | c = max 1 (chi/sqrt (fromIntegral dof)) | 117 | c = max 1 (chi/sqrt (fromIntegral dof)) |
116 | dof = length dat - (rows cov) | 118 | dof = length dat - (rows cov) |
117 | chi = norm2 (fromList $ cost (resMs model) dat sol) | 119 | chi = norm_2 (fromList $ cost (resMs model) dat sol) |
118 | js = fromLists $ jacobian (resDs deriv) dat sol | 120 | js = fromLists $ jacobian (resDs deriv) dat sol |
119 | cov = inv $ trans js <> js | 121 | cov = inv $ tr js <> js |
120 | errs = toList $ scalar c * sqrt (takeDiag cov) | 122 | errs = toList $ scalar c * sqrt (takeDiag cov) |
121 | 123 | ||
122 | 124 | ||
diff --git a/packages/gsl/src/Numeric/GSL/Fourier.hs b/packages/gsl/src/Numeric/GSL/Fourier.hs index 734325b..d824b4f 100644 --- a/packages/gsl/src/Numeric/GSL/Fourier.hs +++ b/packages/gsl/src/Numeric/GSL/Fourier.hs | |||
@@ -16,14 +16,13 @@ module Numeric.GSL.Fourier ( | |||
16 | ifft | 16 | ifft |
17 | ) where | 17 | ) where |
18 | 18 | ||
19 | import Data.Packed | 19 | import Numeric.LinearAlgebra.HMatrix |
20 | import Numeric.GSL.Internal | 20 | import Numeric.GSL.Internal |
21 | import Data.Complex | ||
22 | import Foreign.C.Types | 21 | import Foreign.C.Types |
23 | import System.IO.Unsafe (unsafePerformIO) | 22 | import System.IO.Unsafe (unsafePerformIO) |
24 | 23 | ||
25 | genfft code v = unsafePerformIO $ do | 24 | genfft code v = unsafePerformIO $ do |
26 | r <- createVector (dim v) | 25 | r <- createVector (size v) |
27 | app2 (c_fft code) vec v vec r "fft" | 26 | app2 (c_fft code) vec v vec r "fft" |
28 | return r | 27 | return r |
29 | 28 | ||
diff --git a/packages/gsl/src/Numeric/GSL/IO.hs b/packages/gsl/src/Numeric/GSL/IO.hs index 0d6031a..936f6bf 100644 --- a/packages/gsl/src/Numeric/GSL/IO.hs +++ b/packages/gsl/src/Numeric/GSL/IO.hs | |||
@@ -14,7 +14,7 @@ module Numeric.GSL.IO ( | |||
14 | fileDimensions, loadMatrix, fromFile | 14 | fileDimensions, loadMatrix, fromFile |
15 | ) where | 15 | ) where |
16 | 16 | ||
17 | import Data.Packed | 17 | import Numeric.LinearAlgebra.HMatrix hiding(saveMatrix, loadMatrix) |
18 | import Numeric.GSL.Vector | 18 | import Numeric.GSL.Vector |
19 | import System.Process(readProcess) | 19 | import System.Process(readProcess) |
20 | 20 | ||
diff --git a/packages/gsl/src/Numeric/GSL/Internal.hs b/packages/gsl/src/Numeric/GSL/Internal.hs index a1c4e0c..a269224 100644 --- a/packages/gsl/src/Numeric/GSL/Internal.hs +++ b/packages/gsl/src/Numeric/GSL/Internal.hs | |||
@@ -22,14 +22,13 @@ module Numeric.GSL.Internal( | |||
22 | aux_vTom, | 22 | aux_vTom, |
23 | createV, | 23 | createV, |
24 | createMIO, | 24 | createMIO, |
25 | module Data.Packed.Development, | 25 | module Numeric.LinearAlgebra.Devel, |
26 | check, | 26 | check, |
27 | Res,TV,TM,TCV,TCM | 27 | Res,TV,TM,TCV,TCM |
28 | ) where | 28 | ) where |
29 | 29 | ||
30 | import Data.Packed | 30 | import Numeric.LinearAlgebra.HMatrix |
31 | import Data.Packed.Development hiding (check) | 31 | import Numeric.LinearAlgebra.Devel hiding (check) |
32 | import Data.Complex | ||
33 | 32 | ||
34 | import Foreign.Marshal.Array(copyArray) | 33 | import Foreign.Marshal.Array(copyArray) |
35 | import Foreign.Ptr(Ptr, FunPtr) | 34 | import Foreign.Ptr(Ptr, FunPtr) |
diff --git a/packages/gsl/src/Numeric/GSL/LinearAlgebra.hs b/packages/gsl/src/Numeric/GSL/LinearAlgebra.hs index 17e2258..cb78bf4 100644 --- a/packages/gsl/src/Numeric/GSL/LinearAlgebra.hs +++ b/packages/gsl/src/Numeric/GSL/LinearAlgebra.hs | |||
@@ -15,7 +15,7 @@ module Numeric.GSL.LinearAlgebra ( | |||
15 | fileDimensions, loadMatrix, fromFile | 15 | fileDimensions, loadMatrix, fromFile |
16 | ) where | 16 | ) where |
17 | 17 | ||
18 | import Data.Packed | 18 | import Numeric.LinearAlgebra.HMatrix hiding (RandDist,randomVector,saveMatrix,loadMatrix) |
19 | import Numeric.GSL.Internal hiding (TV,TM,TCV,TCM) | 19 | import Numeric.GSL.Internal hiding (TV,TM,TCV,TCM) |
20 | 20 | ||
21 | import Foreign.Marshal.Alloc(free) | 21 | import Foreign.Marshal.Alloc(free) |
diff --git a/packages/gsl/src/Numeric/GSL/Minimization.hs b/packages/gsl/src/Numeric/GSL/Minimization.hs index 056d463..00e0619 100644 --- a/packages/gsl/src/Numeric/GSL/Minimization.hs +++ b/packages/gsl/src/Numeric/GSL/Minimization.hs | |||
@@ -1,3 +1,6 @@ | |||
1 | {-# LANGUAGE FlexibleContexts #-} | ||
2 | |||
3 | |||
1 | {- | | 4 | {- | |
2 | Module : Numeric.GSL.Minimization | 5 | Module : Numeric.GSL.Minimization |
3 | Copyright : (c) Alberto Ruiz 2006-9 | 6 | Copyright : (c) Alberto Ruiz 2006-9 |
@@ -56,7 +59,7 @@ module Numeric.GSL.Minimization ( | |||
56 | ) where | 59 | ) where |
57 | 60 | ||
58 | 61 | ||
59 | import Data.Packed | 62 | import Numeric.LinearAlgebra.HMatrix hiding(step) |
60 | import Numeric.GSL.Internal | 63 | import Numeric.GSL.Internal |
61 | 64 | ||
62 | import Foreign.Ptr(Ptr, FunPtr, freeHaskellFunPtr) | 65 | import Foreign.Ptr(Ptr, FunPtr, freeHaskellFunPtr) |
@@ -99,7 +102,7 @@ uniMinimizeGen m f xmin xl xu epsrel maxit = unsafePerformIO $ do | |||
99 | rawpath <- createMIO maxit 4 | 102 | rawpath <- createMIO maxit 4 |
100 | (c_uniMinize m fp epsrel (fi maxit) xmin xl xu) | 103 | (c_uniMinize m fp epsrel (fi maxit) xmin xl xu) |
101 | "uniMinimize" | 104 | "uniMinimize" |
102 | let it = round (rawpath @@> (maxit-1,0)) | 105 | let it = round (rawpath `atIndex` (maxit-1,0)) |
103 | path = takeRows it rawpath | 106 | path = takeRows it rawpath |
104 | [sol] = toLists $ dropRows (it-1) path | 107 | [sol] = toLists $ dropRows (it-1) path |
105 | freeHaskellFunPtr fp | 108 | freeHaskellFunPtr fp |
@@ -137,13 +140,13 @@ minimize method eps maxit sz f xi = v2l $ minimizeV method eps maxit (fromList s | |||
137 | ww2 w1 o1 w2 o2 f = w1 o1 $ \a1 -> w2 o2 $ \a2 -> f a1 a2 | 140 | ww2 w1 o1 w2 o2 f = w1 o1 $ \a1 -> w2 o2 $ \a2 -> f a1 a2 |
138 | 141 | ||
139 | minimizeV method eps maxit szv f xiv = unsafePerformIO $ do | 142 | minimizeV method eps maxit szv f xiv = unsafePerformIO $ do |
140 | let n = dim xiv | 143 | let n = size xiv |
141 | fp <- mkVecfun (iv f) | 144 | fp <- mkVecfun (iv f) |
142 | rawpath <- ww2 vec xiv vec szv $ \xiv' szv' -> | 145 | rawpath <- ww2 vec xiv vec szv $ \xiv' szv' -> |
143 | createMIO maxit (n+3) | 146 | createMIO maxit (n+3) |
144 | (c_minimize (fi (fromEnum method)) fp eps (fi maxit) // xiv' // szv') | 147 | (c_minimize (fi (fromEnum method)) fp eps (fi maxit) // xiv' // szv') |
145 | "minimize" | 148 | "minimize" |
146 | let it = round (rawpath @@> (maxit-1,0)) | 149 | let it = round (rawpath `atIndex` (maxit-1,0)) |
147 | path = takeRows it rawpath | 150 | path = takeRows it rawpath |
148 | sol = flatten $ dropColumns 3 $ dropRows (it-1) path | 151 | sol = flatten $ dropColumns 3 $ dropRows (it-1) path |
149 | freeHaskellFunPtr fp | 152 | freeHaskellFunPtr fp |
@@ -191,7 +194,7 @@ minimizeD method eps maxit istep tol f df xi = v2l $ minimizeVD | |||
191 | 194 | ||
192 | 195 | ||
193 | minimizeVD method eps maxit istep tol f df xiv = unsafePerformIO $ do | 196 | minimizeVD method eps maxit istep tol f df xiv = unsafePerformIO $ do |
194 | let n = dim xiv | 197 | let n = size xiv |
195 | f' = f | 198 | f' = f |
196 | df' = (checkdim1 n . df) | 199 | df' = (checkdim1 n . df) |
197 | fp <- mkVecfun (iv f') | 200 | fp <- mkVecfun (iv f') |
@@ -200,7 +203,7 @@ minimizeVD method eps maxit istep tol f df xiv = unsafePerformIO $ do | |||
200 | createMIO maxit (n+2) | 203 | createMIO maxit (n+2) |
201 | (c_minimizeD (fi (fromEnum method)) fp dfp istep tol eps (fi maxit) // xiv') | 204 | (c_minimizeD (fi (fromEnum method)) fp dfp istep tol eps (fi maxit) // xiv') |
202 | "minimizeD" | 205 | "minimizeD" |
203 | let it = round (rawpath @@> (maxit-1,0)) | 206 | let it = round (rawpath `atIndex` (maxit-1,0)) |
204 | path = takeRows it rawpath | 207 | path = takeRows it rawpath |
205 | sol = flatten $ dropColumns 2 $ dropRows (it-1) path | 208 | sol = flatten $ dropColumns 2 $ dropRows (it-1) path |
206 | freeHaskellFunPtr fp | 209 | freeHaskellFunPtr fp |
@@ -217,6 +220,6 @@ foreign import ccall safe "gsl-aux.h minimizeD" | |||
217 | --------------------------------------------------------------------- | 220 | --------------------------------------------------------------------- |
218 | 221 | ||
219 | checkdim1 n v | 222 | checkdim1 n v |
220 | | dim v == n = v | 223 | | size v == n = v |
221 | | otherwise = error $ "Error: "++ show n | 224 | | otherwise = error $ "Error: "++ show n |
222 | ++ " components expected in the result of the gradient supplied to minimizeD" | 225 | ++ " components expected in the result of the gradient supplied to minimizeD" |
diff --git a/packages/gsl/src/Numeric/GSL/ODE.hs b/packages/gsl/src/Numeric/GSL/ODE.hs index 7549a65..3258b83 100644 --- a/packages/gsl/src/Numeric/GSL/ODE.hs +++ b/packages/gsl/src/Numeric/GSL/ODE.hs | |||
@@ -1,3 +1,6 @@ | |||
1 | {-# LANGUAGE FlexibleContexts #-} | ||
2 | |||
3 | |||
1 | {- | | 4 | {- | |
2 | Module : Numeric.GSL.ODE | 5 | Module : Numeric.GSL.ODE |
3 | Copyright : (c) Alberto Ruiz 2010 | 6 | Copyright : (c) Alberto Ruiz 2010 |
@@ -32,7 +35,7 @@ module Numeric.GSL.ODE ( | |||
32 | odeSolve, odeSolveV, ODEMethod(..), Jacobian | 35 | odeSolve, odeSolveV, ODEMethod(..), Jacobian |
33 | ) where | 36 | ) where |
34 | 37 | ||
35 | import Data.Packed | 38 | import Numeric.LinearAlgebra.HMatrix |
36 | import Numeric.GSL.Internal | 39 | import Numeric.GSL.Internal |
37 | 40 | ||
38 | import Foreign.Ptr(FunPtr, nullFunPtr, freeHaskellFunPtr) | 41 | import Foreign.Ptr(FunPtr, nullFunPtr, freeHaskellFunPtr) |
@@ -68,7 +71,7 @@ odeSolve | |||
68 | -> Vector Double -- ^ desired solution times | 71 | -> Vector Double -- ^ desired solution times |
69 | -> Matrix Double -- ^ solution | 72 | -> Matrix Double -- ^ solution |
70 | odeSolve xdot xi ts = odeSolveV RKf45 hi epsAbs epsRel (l2v xdot) (fromList xi) ts | 73 | odeSolve xdot xi ts = odeSolveV RKf45 hi epsAbs epsRel (l2v xdot) (fromList xi) ts |
71 | where hi = (ts@>1 - ts@>0)/100 | 74 | where hi = (ts!1 - ts!0)/100 |
72 | epsAbs = 1.49012e-08 | 75 | epsAbs = 1.49012e-08 |
73 | epsRel = 1.49012e-08 | 76 | epsRel = 1.49012e-08 |
74 | l2v f = \t -> fromList . f t . toList | 77 | l2v f = \t -> fromList . f t . toList |
@@ -107,14 +110,14 @@ odeSolveV' | |||
107 | -> Vector Double -- ^ desired solution times | 110 | -> Vector Double -- ^ desired solution times |
108 | -> Matrix Double -- ^ solution | 111 | -> Matrix Double -- ^ solution |
109 | odeSolveV' method mbjac h epsAbs epsRel f xiv ts = unsafePerformIO $ do | 112 | odeSolveV' method mbjac h epsAbs epsRel f xiv ts = unsafePerformIO $ do |
110 | let n = dim xiv | 113 | let n = size xiv |
111 | fp <- mkDoubleVecVecfun (\t -> aux_vTov (checkdim1 n . f t)) | 114 | fp <- mkDoubleVecVecfun (\t -> aux_vTov (checkdim1 n . f t)) |
112 | jp <- case mbjac of | 115 | jp <- case mbjac of |
113 | Just jac -> mkDoubleVecMatfun (\t -> aux_vTom (checkdim2 n . jac t)) | 116 | Just jac -> mkDoubleVecMatfun (\t -> aux_vTom (checkdim2 n . jac t)) |
114 | Nothing -> return nullFunPtr | 117 | Nothing -> return nullFunPtr |
115 | sol <- vec xiv $ \xiv' -> | 118 | sol <- vec xiv $ \xiv' -> |
116 | vec (checkTimes ts) $ \ts' -> | 119 | vec (checkTimes ts) $ \ts' -> |
117 | createMIO (dim ts) n | 120 | createMIO (size ts) n |
118 | (ode_c (method) h epsAbs epsRel fp jp // xiv' // ts' ) | 121 | (ode_c (method) h epsAbs epsRel fp jp // xiv' // ts' ) |
119 | "ode" | 122 | "ode" |
120 | freeHaskellFunPtr fp | 123 | freeHaskellFunPtr fp |
@@ -126,7 +129,7 @@ foreign import ccall safe "ode" | |||
126 | ------------------------------------------------------- | 129 | ------------------------------------------------------- |
127 | 130 | ||
128 | checkdim1 n v | 131 | checkdim1 n v |
129 | | dim v == n = v | 132 | | size v == n = v |
130 | | otherwise = error $ "Error: "++ show n | 133 | | otherwise = error $ "Error: "++ show n |
131 | ++ " components expected in the result of the function supplied to odeSolve" | 134 | ++ " components expected in the result of the function supplied to odeSolve" |
132 | 135 | ||
@@ -135,6 +138,7 @@ checkdim2 n m | |||
135 | | otherwise = error $ "Error: "++ show n ++ "x" ++ show n | 138 | | otherwise = error $ "Error: "++ show n ++ "x" ++ show n |
136 | ++ " Jacobian expected in odeSolve" | 139 | ++ " Jacobian expected in odeSolve" |
137 | 140 | ||
138 | checkTimes ts | dim ts > 1 && all (>0) (zipWith subtract ts' (tail ts')) = ts | 141 | checkTimes ts | size ts > 1 && all (>0) (zipWith subtract ts' (tail ts')) = ts |
139 | | otherwise = error "odeSolve requires increasing times" | 142 | | otherwise = error "odeSolve requires increasing times" |
140 | where ts' = toList ts | 143 | where ts' = toList ts |
144 | |||
diff --git a/packages/gsl/src/Numeric/GSL/Polynomials.hs b/packages/gsl/src/Numeric/GSL/Polynomials.hs index b1be85d..246e301 100644 --- a/packages/gsl/src/Numeric/GSL/Polynomials.hs +++ b/packages/gsl/src/Numeric/GSL/Polynomials.hs | |||
@@ -16,9 +16,8 @@ module Numeric.GSL.Polynomials ( | |||
16 | polySolve | 16 | polySolve |
17 | ) where | 17 | ) where |
18 | 18 | ||
19 | import Data.Packed | 19 | import Numeric.LinearAlgebra.HMatrix |
20 | import Numeric.GSL.Internal | 20 | import Numeric.GSL.Internal |
21 | import Data.Complex | ||
22 | import System.IO.Unsafe (unsafePerformIO) | 21 | import System.IO.Unsafe (unsafePerformIO) |
23 | 22 | ||
24 | #if __GLASGOW_HASKELL__ >= 704 | 23 | #if __GLASGOW_HASKELL__ >= 704 |
@@ -47,8 +46,8 @@ polySolve :: [Double] -> [Complex Double] | |||
47 | polySolve = toList . polySolve' . fromList | 46 | polySolve = toList . polySolve' . fromList |
48 | 47 | ||
49 | polySolve' :: Vector Double -> Vector (Complex Double) | 48 | polySolve' :: Vector Double -> Vector (Complex Double) |
50 | polySolve' v | dim v > 1 = unsafePerformIO $ do | 49 | polySolve' v | size v > 1 = unsafePerformIO $ do |
51 | r <- createVector (dim v-1) | 50 | r <- createVector (size v-1) |
52 | app2 c_polySolve vec v vec r "polySolve" | 51 | app2 c_polySolve vec v vec r "polySolve" |
53 | return r | 52 | return r |
54 | | otherwise = error "polySolve on a polynomial of degree zero" | 53 | | otherwise = error "polySolve on a polynomial of degree zero" |
diff --git a/packages/gsl/src/Numeric/GSL/Random.hs b/packages/gsl/src/Numeric/GSL/Random.hs index f1f49e5..139c921 100644 --- a/packages/gsl/src/Numeric/GSL/Random.hs +++ b/packages/gsl/src/Numeric/GSL/Random.hs | |||
@@ -21,11 +21,13 @@ module Numeric.GSL.Random ( | |||
21 | ) where | 21 | ) where |
22 | 22 | ||
23 | import Numeric.GSL.Vector | 23 | import Numeric.GSL.Vector |
24 | import Numeric.LinearAlgebra(cholSH) | 24 | import Numeric.LinearAlgebra.HMatrix hiding ( |
25 | import Numeric.Container hiding ( | ||
26 | randomVector, | 25 | randomVector, |
27 | gaussianSample, | 26 | gaussianSample, |
28 | uniformSample | 27 | uniformSample, |
28 | Seed, | ||
29 | rand, | ||
30 | randn | ||
29 | ) | 31 | ) |
30 | import System.Random(randomIO) | 32 | import System.Random(randomIO) |
31 | 33 | ||
@@ -40,10 +42,10 @@ gaussianSample :: Seed | |||
40 | -> Matrix Double -- ^ covariance matrix | 42 | -> Matrix Double -- ^ covariance matrix |
41 | -> Matrix Double -- ^ result | 43 | -> Matrix Double -- ^ result |
42 | gaussianSample seed n med cov = m where | 44 | gaussianSample seed n med cov = m where |
43 | c = dim med | 45 | c = size med |
44 | meds = konst 1 n `outer` med | 46 | meds = konst 1 n `outer` med |
45 | rs = reshape c $ randomVector seed Gaussian (c * n) | 47 | rs = reshape c $ randomVector seed Gaussian (c * n) |
46 | m = rs `mXm` cholSH cov `add` meds | 48 | m = rs <> cholSH cov + meds |
47 | 49 | ||
48 | -- | Obtains a matrix whose rows are pseudorandom samples from a multivariate | 50 | -- | Obtains a matrix whose rows are pseudorandom samples from a multivariate |
49 | -- uniform distribution. | 51 | -- uniform distribution. |
@@ -55,10 +57,10 @@ uniformSample seed n rgs = m where | |||
55 | (as,bs) = unzip rgs | 57 | (as,bs) = unzip rgs |
56 | a = fromList as | 58 | a = fromList as |
57 | cs = zipWith subtract as bs | 59 | cs = zipWith subtract as bs |
58 | d = dim a | 60 | d = size a |
59 | dat = toRows $ reshape n $ randomVector seed Uniform (n*d) | 61 | dat = toRows $ reshape n $ randomVector seed Uniform (n*d) |
60 | am = konst 1 n `outer` a | 62 | am = konst 1 n `outer` a |
61 | m = fromColumns (zipWith scale cs dat) `add` am | 63 | m = fromColumns (zipWith scale cs dat) + am |
62 | 64 | ||
63 | -- | pseudorandom matrix with uniform elements between 0 and 1 | 65 | -- | pseudorandom matrix with uniform elements between 0 and 1 |
64 | randm :: RandDist | 66 | randm :: RandDist |
diff --git a/packages/gsl/src/Numeric/GSL/Root.hs b/packages/gsl/src/Numeric/GSL/Root.hs index b9f3b94..724f32f 100644 --- a/packages/gsl/src/Numeric/GSL/Root.hs +++ b/packages/gsl/src/Numeric/GSL/Root.hs | |||
@@ -1,3 +1,5 @@ | |||
1 | {-# LANGUAGE FlexibleContexts #-} | ||
2 | |||
1 | {- | | 3 | {- | |
2 | Module : Numeric.GSL.Root | 4 | Module : Numeric.GSL.Root |
3 | Copyright : (c) Alberto Ruiz 2009 | 5 | Copyright : (c) Alberto Ruiz 2009 |
@@ -39,7 +41,7 @@ module Numeric.GSL.Root ( | |||
39 | rootJ, RootMethodJ(..), | 41 | rootJ, RootMethodJ(..), |
40 | ) where | 42 | ) where |
41 | 43 | ||
42 | import Data.Packed | 44 | import Numeric.LinearAlgebra.HMatrix |
43 | import Numeric.GSL.Internal | 45 | import Numeric.GSL.Internal |
44 | import Foreign.Ptr(FunPtr, freeHaskellFunPtr) | 46 | import Foreign.Ptr(FunPtr, freeHaskellFunPtr) |
45 | import Foreign.C.Types | 47 | import Foreign.C.Types |
@@ -69,7 +71,7 @@ uniRootGen m f xl xu epsrel maxit = unsafePerformIO $ do | |||
69 | rawpath <- createMIO maxit 4 | 71 | rawpath <- createMIO maxit 4 |
70 | (c_root m fp epsrel (fi maxit) xl xu) | 72 | (c_root m fp epsrel (fi maxit) xl xu) |
71 | "root" | 73 | "root" |
72 | let it = round (rawpath @@> (maxit-1,0)) | 74 | let it = round (rawpath `atIndex` (maxit-1,0)) |
73 | path = takeRows it rawpath | 75 | path = takeRows it rawpath |
74 | [sol] = toLists $ dropRows (it-1) path | 76 | [sol] = toLists $ dropRows (it-1) path |
75 | freeHaskellFunPtr fp | 77 | freeHaskellFunPtr fp |
@@ -100,7 +102,7 @@ uniRootJGen m f df x epsrel maxit = unsafePerformIO $ do | |||
100 | rawpath <- createMIO maxit 2 | 102 | rawpath <- createMIO maxit 2 |
101 | (c_rootj m fp dfp epsrel (fi maxit) x) | 103 | (c_rootj m fp dfp epsrel (fi maxit) x) |
102 | "rootj" | 104 | "rootj" |
103 | let it = round (rawpath @@> (maxit-1,0)) | 105 | let it = round (rawpath `atIndex` (maxit-1,0)) |
104 | path = takeRows it rawpath | 106 | path = takeRows it rawpath |
105 | [sol] = toLists $ dropRows (it-1) path | 107 | [sol] = toLists $ dropRows (it-1) path |
106 | freeHaskellFunPtr fp | 108 | freeHaskellFunPtr fp |
@@ -132,13 +134,13 @@ root method epsabs maxit fun xinit = rootGen (fi (fromEnum method)) fun xinit ep | |||
132 | 134 | ||
133 | rootGen m f xi epsabs maxit = unsafePerformIO $ do | 135 | rootGen m f xi epsabs maxit = unsafePerformIO $ do |
134 | let xiv = fromList xi | 136 | let xiv = fromList xi |
135 | n = dim xiv | 137 | n = size xiv |
136 | fp <- mkVecVecfun (aux_vTov (checkdim1 n . fromList . f . toList)) | 138 | fp <- mkVecVecfun (aux_vTov (checkdim1 n . fromList . f . toList)) |
137 | rawpath <- vec xiv $ \xiv' -> | 139 | rawpath <- vec xiv $ \xiv' -> |
138 | createMIO maxit (2*n+1) | 140 | createMIO maxit (2*n+1) |
139 | (c_multiroot m fp epsabs (fi maxit) // xiv') | 141 | (c_multiroot m fp epsabs (fi maxit) // xiv') |
140 | "multiroot" | 142 | "multiroot" |
141 | let it = round (rawpath @@> (maxit-1,0)) | 143 | let it = round (rawpath `atIndex` (maxit-1,0)) |
142 | path = takeRows it rawpath | 144 | path = takeRows it rawpath |
143 | [sol] = toLists $ dropRows (it-1) path | 145 | [sol] = toLists $ dropRows (it-1) path |
144 | freeHaskellFunPtr fp | 146 | freeHaskellFunPtr fp |
@@ -169,14 +171,14 @@ rootJ method epsabs maxit fun jac xinit = rootJGen (fi (fromEnum method)) fun ja | |||
169 | 171 | ||
170 | rootJGen m f jac xi epsabs maxit = unsafePerformIO $ do | 172 | rootJGen m f jac xi epsabs maxit = unsafePerformIO $ do |
171 | let xiv = fromList xi | 173 | let xiv = fromList xi |
172 | n = dim xiv | 174 | n = size xiv |
173 | fp <- mkVecVecfun (aux_vTov (checkdim1 n . fromList . f . toList)) | 175 | fp <- mkVecVecfun (aux_vTov (checkdim1 n . fromList . f . toList)) |
174 | jp <- mkVecMatfun (aux_vTom (checkdim2 n . fromLists . jac . toList)) | 176 | jp <- mkVecMatfun (aux_vTom (checkdim2 n . fromLists . jac . toList)) |
175 | rawpath <- vec xiv $ \xiv' -> | 177 | rawpath <- vec xiv $ \xiv' -> |
176 | createMIO maxit (2*n+1) | 178 | createMIO maxit (2*n+1) |
177 | (c_multirootj m fp jp epsabs (fi maxit) // xiv') | 179 | (c_multirootj m fp jp epsabs (fi maxit) // xiv') |
178 | "multiroot" | 180 | "multiroot" |
179 | let it = round (rawpath @@> (maxit-1,0)) | 181 | let it = round (rawpath `atIndex` (maxit-1,0)) |
180 | path = takeRows it rawpath | 182 | path = takeRows it rawpath |
181 | [sol] = toLists $ dropRows (it-1) path | 183 | [sol] = toLists $ dropRows (it-1) path |
182 | freeHaskellFunPtr fp | 184 | freeHaskellFunPtr fp |
@@ -189,7 +191,7 @@ foreign import ccall safe "multirootj" | |||
189 | ------------------------------------------------------- | 191 | ------------------------------------------------------- |
190 | 192 | ||
191 | checkdim1 n v | 193 | checkdim1 n v |
192 | | dim v == n = v | 194 | | size v == n = v |
193 | | otherwise = error $ "Error: "++ show n | 195 | | otherwise = error $ "Error: "++ show n |
194 | ++ " components expected in the result of the function supplied to root" | 196 | ++ " components expected in the result of the function supplied to root" |
195 | 197 | ||
diff --git a/packages/gsl/src/Numeric/GSL/Vector.hs b/packages/gsl/src/Numeric/GSL/Vector.hs index af79f32..0cd99eb 100644 --- a/packages/gsl/src/Numeric/GSL/Vector.hs +++ b/packages/gsl/src/Numeric/GSL/Vector.hs | |||
@@ -14,8 +14,7 @@ module Numeric.GSL.Vector ( | |||
14 | fwriteVector, freadVector, fprintfVector, fscanfVector | 14 | fwriteVector, freadVector, fprintfVector, fscanfVector |
15 | ) where | 15 | ) where |
16 | 16 | ||
17 | import Data.Packed | 17 | import Numeric.LinearAlgebra.HMatrix hiding(randomVector, saveMatrix) |
18 | import Numeric.LinearAlgebra(RandDist(..)) | ||
19 | import Numeric.GSL.Internal hiding (TV,TM,TCV,TCM) | 18 | import Numeric.GSL.Internal hiding (TV,TM,TCV,TCM) |
20 | 19 | ||
21 | import Foreign.Marshal.Alloc(free) | 20 | import Foreign.Marshal.Alloc(free) |