diff options
Diffstat (limited to 'lib/Numeric')
-rw-r--r-- | lib/Numeric/GSL/Differentiation.hs | 2 | ||||
-rw-r--r-- | lib/Numeric/GSL/Fourier.hs | 3 | ||||
-rw-r--r-- | lib/Numeric/GSL/Integration.hs | 4 | ||||
-rw-r--r-- | lib/Numeric/GSL/Matrix.hs | 52 | ||||
-rw-r--r-- | lib/Numeric/GSL/Minimization.hs | 34 | ||||
-rw-r--r-- | lib/Numeric/GSL/Polynomials.hs | 3 | ||||
-rw-r--r-- | lib/Numeric/GSL/Special/Internal.hs | 4 | ||||
-rw-r--r-- | lib/Numeric/GSL/Vector.hs | 12 | ||||
-rw-r--r-- | lib/Numeric/LinearAlgebra/LAPACK.hs | 34 |
9 files changed, 91 insertions, 57 deletions
diff --git a/lib/Numeric/GSL/Differentiation.hs b/lib/Numeric/GSL/Differentiation.hs index e7fea92..09236bd 100644 --- a/lib/Numeric/GSL/Differentiation.hs +++ b/lib/Numeric/GSL/Differentiation.hs | |||
@@ -35,7 +35,7 @@ derivGen c h f x = unsafePerformIO $ do | |||
35 | r <- malloc | 35 | r <- malloc |
36 | e <- malloc | 36 | e <- malloc |
37 | fp <- mkfun (\x _ -> f x) | 37 | fp <- mkfun (\x _ -> f x) |
38 | c_deriv c fp x h r e // check "deriv" [] | 38 | c_deriv c fp x h r e // check "deriv" |
39 | vr <- peek r | 39 | vr <- peek r |
40 | ve <- peek e | 40 | ve <- peek e |
41 | let result = (vr,ve) | 41 | let result = (vr,ve) |
diff --git a/lib/Numeric/GSL/Fourier.hs b/lib/Numeric/GSL/Fourier.hs index e975fbf..4b08625 100644 --- a/lib/Numeric/GSL/Fourier.hs +++ b/lib/Numeric/GSL/Fourier.hs | |||
@@ -26,7 +26,8 @@ import Foreign | |||
26 | 26 | ||
27 | genfft code v = unsafePerformIO $ do | 27 | genfft code v = unsafePerformIO $ do |
28 | r <- createVector (dim v) | 28 | r <- createVector (dim v) |
29 | c_fft code // vec v // vec r // check "fft" [v] | 29 | ww2 withVector v withVector r $ \ v r -> |
30 | c_fft code // v // r // check "fft" | ||
30 | return r | 31 | return r |
31 | 32 | ||
32 | foreign import ccall "gsl-aux.h fft" c_fft :: Int -> TCVCV | 33 | foreign import ccall "gsl-aux.h fft" c_fft :: Int -> TCVCV |
diff --git a/lib/Numeric/GSL/Integration.hs b/lib/Numeric/GSL/Integration.hs index d756417..747b34c 100644 --- a/lib/Numeric/GSL/Integration.hs +++ b/lib/Numeric/GSL/Integration.hs | |||
@@ -44,7 +44,7 @@ integrateQAGS prec n f a b = unsafePerformIO $ do | |||
44 | r <- malloc | 44 | r <- malloc |
45 | e <- malloc | 45 | e <- malloc |
46 | fp <- mkfun (\x _ -> f x) | 46 | fp <- mkfun (\x _ -> f x) |
47 | c_integrate_qags fp a b prec n r e // check "integrate_qags" [] | 47 | c_integrate_qags fp a b prec n r e // check "integrate_qags" |
48 | vr <- peek r | 48 | vr <- peek r |
49 | ve <- peek e | 49 | ve <- peek e |
50 | let result = (vr,ve) | 50 | let result = (vr,ve) |
@@ -75,7 +75,7 @@ integrateQNG prec f a b = unsafePerformIO $ do | |||
75 | r <- malloc | 75 | r <- malloc |
76 | e <- malloc | 76 | e <- malloc |
77 | fp <- mkfun (\x _ -> f x) | 77 | fp <- mkfun (\x _ -> f x) |
78 | c_integrate_qng fp a b prec r e // check "integrate_qng" [] | 78 | c_integrate_qng fp a b prec r e // check "integrate_qng" |
79 | vr <- peek r | 79 | vr <- peek r |
80 | ve <- peek e | 80 | ve <- peek e |
81 | let result = (vr,ve) | 81 | let result = (vr,ve) |
diff --git a/lib/Numeric/GSL/Matrix.hs b/lib/Numeric/GSL/Matrix.hs index e803c53..09a0be4 100644 --- a/lib/Numeric/GSL/Matrix.hs +++ b/lib/Numeric/GSL/Matrix.hs | |||
@@ -51,7 +51,8 @@ eigSg' m | |||
51 | | otherwise = unsafePerformIO $ do | 51 | | otherwise = unsafePerformIO $ do |
52 | l <- createVector r | 52 | l <- createVector r |
53 | v <- createMatrix RowMajor r r | 53 | v <- createMatrix RowMajor r r |
54 | c_eigS // matc m // vec l // matc v // check "eigSg" [cdat m] | 54 | ww3 withMatrix m withVector l withMatrix v $ \m l v -> |
55 | c_eigS // m // l // v // check "eigSg" | ||
55 | return (l,v) | 56 | return (l,v) |
56 | where r = rows m | 57 | where r = rows m |
57 | foreign import ccall "gsl-aux.h eigensystemR" c_eigS :: TMVM | 58 | foreign import ccall "gsl-aux.h eigensystemR" c_eigS :: TMVM |
@@ -84,7 +85,8 @@ eigHg' m | |||
84 | | otherwise = unsafePerformIO $ do | 85 | | otherwise = unsafePerformIO $ do |
85 | l <- createVector r | 86 | l <- createVector r |
86 | v <- createMatrix RowMajor r r | 87 | v <- createMatrix RowMajor r r |
87 | c_eigH // matc m // vec l // matc v // check "eigHg" [cdat m] | 88 | ww3 withMatrix m withVector l withMatrix v $ \m l v -> |
89 | c_eigH // m // l // v // check "eigHg" | ||
88 | return (l,v) | 90 | return (l,v) |
89 | where r = rows m | 91 | where r = rows m |
90 | foreign import ccall "gsl-aux.h eigensystemC" c_eigH :: TCMVCM | 92 | foreign import ccall "gsl-aux.h eigensystemC" c_eigH :: TCMVCM |
@@ -120,7 +122,8 @@ svd' x = unsafePerformIO $ do | |||
120 | u <- createMatrix RowMajor r c | 122 | u <- createMatrix RowMajor r c |
121 | s <- createVector c | 123 | s <- createVector c |
122 | v <- createMatrix RowMajor c c | 124 | v <- createMatrix RowMajor c c |
123 | c_svd // matc x // matc u // vec s // matc v // check "svdg" [cdat x] | 125 | ww4 withMatrix x withMatrix u withVector s withMatrix v $ \x u s v -> |
126 | c_svd // x // u // s // v // check "svdg" | ||
124 | return (u,s,v) | 127 | return (u,s,v) |
125 | where r = rows x | 128 | where r = rows x |
126 | c = cols x | 129 | c = cols x |
@@ -149,7 +152,8 @@ qr = qr' . cmat | |||
149 | qr' x = unsafePerformIO $ do | 152 | qr' x = unsafePerformIO $ do |
150 | q <- createMatrix RowMajor r r | 153 | q <- createMatrix RowMajor r r |
151 | rot <- createMatrix RowMajor r c | 154 | rot <- createMatrix RowMajor r c |
152 | c_qr // matc x // matc q // matc rot // check "qr" [cdat x] | 155 | ww3 withMatrix x withMatrix q withMatrix rot $ \x q rot -> |
156 | c_qr // x // q // rot // check "qr" | ||
153 | return (q,rot) | 157 | return (q,rot) |
154 | where r = rows x | 158 | where r = rows x |
155 | c = cols x | 159 | c = cols x |
@@ -161,7 +165,8 @@ qrPacked = qrPacked' . cmat | |||
161 | qrPacked' x = unsafePerformIO $ do | 165 | qrPacked' x = unsafePerformIO $ do |
162 | qr <- createMatrix RowMajor r c | 166 | qr <- createMatrix RowMajor r c |
163 | tau <- createVector (min r c) | 167 | tau <- createVector (min r c) |
164 | c_qrPacked // matc x // matc qr // vec tau // check "qrUnpacked" [cdat x] | 168 | ww3 withMatrix x withMatrix qr withVector tau $ \x qr tau -> |
169 | c_qrPacked // x // qr // tau // check "qrUnpacked" | ||
165 | return (qr,tau) | 170 | return (qr,tau) |
166 | where r = rows x | 171 | where r = rows x |
167 | c = cols x | 172 | c = cols x |
@@ -172,9 +177,10 @@ unpackQR (qr,tau) = unpackQR' (cmat qr, tau) | |||
172 | 177 | ||
173 | unpackQR' (qr,tau) = unsafePerformIO $ do | 178 | unpackQR' (qr,tau) = unsafePerformIO $ do |
174 | q <- createMatrix RowMajor r r | 179 | q <- createMatrix RowMajor r r |
175 | rot <- createMatrix RowMajor r c | 180 | res <- createMatrix RowMajor r c |
176 | c_qrUnpack // matc qr // vec tau // matc q // matc rot // check "qrUnpack" [cdat qr,tau] | 181 | ww4 withMatrix qr withVector tau withMatrix q withMatrix res $ \qr tau q res -> |
177 | return (q,rot) | 182 | c_qrUnpack // qr // tau // q // res // check "qrUnpack" |
183 | return (q,res) | ||
178 | where r = rows qr | 184 | where r = rows qr |
179 | c = cols qr | 185 | c = cols qr |
180 | foreign import ccall "gsl-aux.h QRunpack" c_qrUnpack :: TMVMM | 186 | foreign import ccall "gsl-aux.h QRunpack" c_qrUnpack :: TMVMM |
@@ -196,20 +202,22 @@ cholR :: Matrix Double -> Matrix Double | |||
196 | cholR = cholR' . cmat | 202 | cholR = cholR' . cmat |
197 | 203 | ||
198 | cholR' x = unsafePerformIO $ do | 204 | cholR' x = unsafePerformIO $ do |
199 | res <- createMatrix RowMajor r r | 205 | r <- createMatrix RowMajor n n |
200 | c_cholR // matc x // matc res // check "cholR" [cdat x] | 206 | ww2 withMatrix x withMatrix r $ \x r -> |
201 | return res | 207 | c_cholR // x // r // check "cholR" |
202 | where r = rows x | 208 | return r |
209 | where n = rows x | ||
203 | foreign import ccall "gsl-aux.h cholR" c_cholR :: TMM | 210 | foreign import ccall "gsl-aux.h cholR" c_cholR :: TMM |
204 | 211 | ||
205 | cholC :: Matrix (Complex Double) -> Matrix (Complex Double) | 212 | cholC :: Matrix (Complex Double) -> Matrix (Complex Double) |
206 | cholC = cholC' . cmat | 213 | cholC = cholC' . cmat |
207 | 214 | ||
208 | cholC' x = unsafePerformIO $ do | 215 | cholC' x = unsafePerformIO $ do |
209 | res <- createMatrix RowMajor r r | 216 | r <- createMatrix RowMajor n n |
210 | c_cholC // matc x // matc res // check "cholC" [cdat x] | 217 | ww2 withMatrix x withMatrix r $ \x r -> |
211 | return res | 218 | c_cholC // x // r // check "cholC" |
212 | where r = rows x | 219 | return r |
220 | where n = rows x | ||
213 | foreign import ccall "gsl-aux.h cholC" c_cholC :: TCMCM | 221 | foreign import ccall "gsl-aux.h cholC" c_cholC :: TCMCM |
214 | 222 | ||
215 | 223 | ||
@@ -223,7 +231,8 @@ luSolveR a b = luSolveR' (cmat a) (cmat b) | |||
223 | luSolveR' a b | 231 | luSolveR' a b |
224 | | n1==n2 && n1==r = unsafePerformIO $ do | 232 | | n1==n2 && n1==r = unsafePerformIO $ do |
225 | s <- createMatrix RowMajor r c | 233 | s <- createMatrix RowMajor r c |
226 | c_luSolveR // matc a // matc b // matc s // check "luSolveR" [cdat a, cdat b] | 234 | ww3 withMatrix a withMatrix b withMatrix s $ \ a b s -> |
235 | c_luSolveR // a // b // s // check "luSolveR" | ||
227 | return s | 236 | return s |
228 | | otherwise = error "luSolveR of nonsquare matrix" | 237 | | otherwise = error "luSolveR of nonsquare matrix" |
229 | where n1 = rows a | 238 | where n1 = rows a |
@@ -240,7 +249,8 @@ luSolveC a b = luSolveC' (cmat a) (cmat b) | |||
240 | luSolveC' a b | 249 | luSolveC' a b |
241 | | n1==n2 && n1==r = unsafePerformIO $ do | 250 | | n1==n2 && n1==r = unsafePerformIO $ do |
242 | s <- createMatrix RowMajor r c | 251 | s <- createMatrix RowMajor r c |
243 | c_luSolveC // matc a // matc b // matc s // check "luSolveC" [cdat a, cdat b] | 252 | ww3 withMatrix a withMatrix b withMatrix s $ \ a b s -> |
253 | c_luSolveC // a // b // s // check "luSolveC" | ||
244 | return s | 254 | return s |
245 | | otherwise = error "luSolveC of nonsquare matrix" | 255 | | otherwise = error "luSolveC of nonsquare matrix" |
246 | where n1 = rows a | 256 | where n1 = rows a |
@@ -256,7 +266,8 @@ luRaux = luRaux' . cmat | |||
256 | 266 | ||
257 | luRaux' x = unsafePerformIO $ do | 267 | luRaux' x = unsafePerformIO $ do |
258 | res <- createVector (r*r+r+1) | 268 | res <- createVector (r*r+r+1) |
259 | c_luRaux // matc x // vec res // check "luRaux" [cdat x] | 269 | ww2 withMatrix x withVector res $ \x res -> |
270 | c_luRaux // x // res // check "luRaux" | ||
260 | return res | 271 | return res |
261 | where r = rows x | 272 | where r = rows x |
262 | c = cols x | 273 | c = cols x |
@@ -269,7 +280,8 @@ luCaux = luCaux' . cmat | |||
269 | 280 | ||
270 | luCaux' x = unsafePerformIO $ do | 281 | luCaux' x = unsafePerformIO $ do |
271 | res <- createVector (r*r+r+1) | 282 | res <- createVector (r*r+r+1) |
272 | c_luCaux // matc x // vec res // check "luCaux" [cdat x] | 283 | ww2 withMatrix x withVector res $ \x res -> |
284 | c_luCaux // x // res // check "luCaux" | ||
273 | return res | 285 | return res |
274 | where r = rows x | 286 | where r = rows x |
275 | c = cols x | 287 | c = cols x |
diff --git a/lib/Numeric/GSL/Minimization.hs b/lib/Numeric/GSL/Minimization.hs index f523849..e44b2e5 100644 --- a/lib/Numeric/GSL/Minimization.hs +++ b/lib/Numeric/GSL/Minimization.hs | |||
@@ -26,7 +26,6 @@ import Data.Packed.Matrix | |||
26 | import Foreign | 26 | import Foreign |
27 | import Complex | 27 | import Complex |
28 | 28 | ||
29 | |||
30 | ------------------------------------------------------------------------- | 29 | ------------------------------------------------------------------------- |
31 | 30 | ||
32 | {- | The method of Nelder and Mead, implemented by /gsl_multimin_fminimizer_nmsimplex/. The gradient of the function is not required. This is the example in the GSL manual: | 31 | {- | The method of Nelder and Mead, implemented by /gsl_multimin_fminimizer_nmsimplex/. The gradient of the function is not required. This is the example in the GSL manual: |
@@ -85,9 +84,10 @@ minimizeNMSimplex f xi sz tol maxit = unsafePerformIO $ do | |||
85 | szv = fromList sz | 84 | szv = fromList sz |
86 | n = dim xiv | 85 | n = dim xiv |
87 | fp <- mkVecfun (iv (f.toList)) | 86 | fp <- mkVecfun (iv (f.toList)) |
88 | rawpath <- createMIO maxit (n+3) | 87 | rawpath <- ww2 withVector xiv withVector szv $ \xiv szv -> |
89 | (c_minimizeNMSimplex fp tol maxit // vec xiv // vec szv) | 88 | createMIO maxit (n+3) |
90 | "minimizeNMSimplex" [xiv,szv] | 89 | (c_minimizeNMSimplex fp tol maxit // xiv // szv) |
90 | "minimizeNMSimplex" | ||
91 | let it = round (rawpath @@> (maxit-1,0)) | 91 | let it = round (rawpath @@> (maxit-1,0)) |
92 | path = takeRows it rawpath | 92 | path = takeRows it rawpath |
93 | [sol] = toLists $ dropRows (it-1) path | 93 | [sol] = toLists $ dropRows (it-1) path |
@@ -150,9 +150,10 @@ minimizeConjugateGradient istep minimpar tol maxit f df xi = unsafePerformIO $ d | |||
150 | df' = (fromList . df . toList) | 150 | df' = (fromList . df . toList) |
151 | fp <- mkVecfun (iv f') | 151 | fp <- mkVecfun (iv f') |
152 | dfp <- mkVecVecfun (aux_vTov df') | 152 | dfp <- mkVecVecfun (aux_vTov df') |
153 | rawpath <- createMIO maxit (n+2) | 153 | rawpath <- withVector xiv $ \xiv -> |
154 | (c_minimizeConjugateGradient fp dfp istep minimpar tol maxit // vec xiv) | 154 | createMIO maxit (n+2) |
155 | "minimizeDerivV" [xiv] | 155 | (c_minimizeConjugateGradient fp dfp istep minimpar tol maxit // xiv) |
156 | "minimizeDerivV" | ||
156 | let it = round (rawpath @@> (maxit-1,0)) | 157 | let it = round (rawpath @@> (maxit-1,0)) |
157 | path = takeRows it rawpath | 158 | path = takeRows it rawpath |
158 | sol = toList $ cdat $ dropColumns 2 $ dropRows (it-1) path | 159 | sol = toList $ cdat $ dropColumns 2 $ dropRows (it-1) path |
@@ -169,7 +170,7 @@ foreign import ccall "gsl-aux.h minimizeWithDeriv" | |||
169 | 170 | ||
170 | --------------------------------------------------------------------- | 171 | --------------------------------------------------------------------- |
171 | iv :: (Vector Double -> Double) -> (Int -> Ptr Double -> Double) | 172 | iv :: (Vector Double -> Double) -> (Int -> Ptr Double -> Double) |
172 | iv f n p = f (createV n copy "iv" []) where | 173 | iv f n p = f (createV n copy "iv") where |
173 | copy n q = do | 174 | copy n q = do |
174 | copyArray q p n | 175 | copyArray q p n |
175 | return 0 | 176 | return 0 |
@@ -187,25 +188,30 @@ foreign import ccall "wrapper" | |||
187 | aux_vTov :: (Vector Double -> Vector Double) -> (Int -> Ptr Double -> Ptr Double -> IO()) | 188 | aux_vTov :: (Vector Double -> Vector Double) -> (Int -> Ptr Double -> Ptr Double -> IO()) |
188 | aux_vTov f n p r = g where | 189 | aux_vTov f n p r = g where |
189 | v@V {fptr = pr} = f x | 190 | v@V {fptr = pr} = f x |
190 | x = createV n copy "aux_vTov" [] | 191 | x = createV n copy "aux_vTov" |
191 | copy n q = do | 192 | copy n q = do |
192 | copyArray q p n | 193 | copyArray q p n |
193 | return 0 | 194 | return 0 |
194 | g = withForeignPtr pr $ \_ -> copyArray r (ptr v) n | 195 | g = withForeignPtr pr $ \p -> copyArray r p n |
195 | 196 | ||
196 | -------------------------------------------------------------------- | 197 | -------------------------------------------------------------------- |
197 | 198 | ||
198 | createV n fun msg ptrs = unsafePerformIO $ do | 199 | |
200 | createV n fun msg = unsafePerformIO $ do | ||
199 | r <- createVector n | 201 | r <- createVector n |
200 | fun // vec r // check msg ptrs | 202 | withVector r $ \ r -> |
203 | fun // r // check msg | ||
201 | return r | 204 | return r |
202 | 205 | ||
206 | {- | ||
203 | createM r c fun msg ptrs = unsafePerformIO $ do | 207 | createM r c fun msg ptrs = unsafePerformIO $ do |
204 | r <- createMatrix RowMajor r c | 208 | r <- createMatrix RowMajor r c |
205 | fun // matc r // check msg ptrs | 209 | fun // matc r // check msg ptrs |
206 | return r | 210 | return r |
211 | -} | ||
207 | 212 | ||
208 | createMIO r c fun msg ptrs = do | 213 | createMIO r c fun msg = do |
209 | r <- createMatrix RowMajor r c | 214 | r <- createMatrix RowMajor r c |
210 | fun // matc r // check msg ptrs | 215 | withMatrix r $ \ r -> |
216 | fun // r // check msg | ||
211 | return r | 217 | return r |
diff --git a/lib/Numeric/GSL/Polynomials.hs b/lib/Numeric/GSL/Polynomials.hs index 42694f0..e663711 100644 --- a/lib/Numeric/GSL/Polynomials.hs +++ b/lib/Numeric/GSL/Polynomials.hs | |||
@@ -47,7 +47,8 @@ polySolve = toList . polySolve' . fromList | |||
47 | polySolve' :: Vector Double -> Vector (Complex Double) | 47 | polySolve' :: Vector Double -> Vector (Complex Double) |
48 | polySolve' v | dim v > 1 = unsafePerformIO $ do | 48 | polySolve' v | dim v > 1 = unsafePerformIO $ do |
49 | r <- createVector (dim v-1) | 49 | r <- createVector (dim v-1) |
50 | c_polySolve // vec v // vec r // check "polySolve" [v] | 50 | ww2 withVector v withVector r $ \ v r -> |
51 | c_polySolve // v // r // check "polySolve" | ||
51 | return r | 52 | return r |
52 | | otherwise = error "polySolve on a polynomial of degree zero" | 53 | | otherwise = error "polySolve on a polynomial of degree zero" |
53 | 54 | ||
diff --git a/lib/Numeric/GSL/Special/Internal.hs b/lib/Numeric/GSL/Special/Internal.hs index a08809b..ca36009 100644 --- a/lib/Numeric/GSL/Special/Internal.hs +++ b/lib/Numeric/GSL/Special/Internal.hs | |||
@@ -45,7 +45,7 @@ type Size_t = Int | |||
45 | createSFR :: Storable a => String -> (Ptr a -> IO Int) -> (a, a) | 45 | createSFR :: Storable a => String -> (Ptr a -> IO Int) -> (a, a) |
46 | createSFR s f = unsafePerformIO $ do | 46 | createSFR s f = unsafePerformIO $ do |
47 | p <- mallocArray 2 | 47 | p <- mallocArray 2 |
48 | f p // check s [] | 48 | f p // check s |
49 | [val,err] <- peekArray 2 p | 49 | [val,err] <- peekArray 2 p |
50 | free p | 50 | free p |
51 | return (val,err) | 51 | return (val,err) |
@@ -60,7 +60,7 @@ createSFR_E10 s f = unsafePerformIO $ do | |||
60 | let sd = sizeOf (0::Double) | 60 | let sd = sizeOf (0::Double) |
61 | let si = sizeOf (0::Int) | 61 | let si = sizeOf (0::Int) |
62 | p <- mallocBytes (2*sd + si) | 62 | p <- mallocBytes (2*sd + si) |
63 | f p // check s [] | 63 | f p // check s |
64 | val <- peekByteOff p 0 | 64 | val <- peekByteOff p 0 |
65 | err <- peekByteOff p sd | 65 | err <- peekByteOff p sd |
66 | expo <- peekByteOff p (2*sd) | 66 | expo <- peekByteOff p (2*sd) |
diff --git a/lib/Numeric/GSL/Vector.hs b/lib/Numeric/GSL/Vector.hs index d94b377..65f3a2e 100644 --- a/lib/Numeric/GSL/Vector.hs +++ b/lib/Numeric/GSL/Vector.hs | |||
@@ -73,24 +73,28 @@ data FunCodeS = Norm2 | |||
73 | 73 | ||
74 | toScalarAux fun code v = unsafePerformIO $ do | 74 | toScalarAux fun code v = unsafePerformIO $ do |
75 | r <- createVector 1 | 75 | r <- createVector 1 |
76 | fun (fromEnum code) // vec v // vec r // check "toScalarAux" [v] | 76 | ww2 withVector v withVector r $ \v r -> |
77 | fun (fromEnum code) // v // r // check "toScalarAux" | ||
77 | return (r `at` 0) | 78 | return (r `at` 0) |
78 | 79 | ||
79 | vectorMapAux fun code v = unsafePerformIO $ do | 80 | vectorMapAux fun code v = unsafePerformIO $ do |
80 | r <- createVector (dim v) | 81 | r <- createVector (dim v) |
81 | fun (fromEnum code) // vec v // vec r // check "vectorMapAux" [v] | 82 | ww2 withVector v withVector r $ \v r -> |
83 | fun (fromEnum code) // v // r // check "vectorMapAux" | ||
82 | return r | 84 | return r |
83 | 85 | ||
84 | vectorMapValAux fun code val v = unsafePerformIO $ do | 86 | vectorMapValAux fun code val v = unsafePerformIO $ do |
85 | r <- createVector (dim v) | 87 | r <- createVector (dim v) |
86 | pval <- newArray [val] | 88 | pval <- newArray [val] |
87 | fun (fromEnum code) pval // vec v // vec r // check "vectorMapValAux" [v] | 89 | ww2 withVector v withVector r $ \v r -> |
90 | fun (fromEnum code) pval // v // r // check "vectorMapValAux" | ||
88 | free pval | 91 | free pval |
89 | return r | 92 | return r |
90 | 93 | ||
91 | vectorZipAux fun code u v = unsafePerformIO $ do | 94 | vectorZipAux fun code u v = unsafePerformIO $ do |
92 | r <- createVector (dim u) | 95 | r <- createVector (dim u) |
93 | fun (fromEnum code) // vec u // vec v // vec r // check "vectorZipAux" [u,v] | 96 | ww3 withVector u withVector v withVector r $ \u v r -> |
97 | fun (fromEnum code) // u // v // r // check "vectorZipAux" | ||
94 | return r | 98 | return r |
95 | 99 | ||
96 | --------------------------------------------------------------------- | 100 | --------------------------------------------------------------------- |
diff --git a/lib/Numeric/LinearAlgebra/LAPACK.hs b/lib/Numeric/LinearAlgebra/LAPACK.hs index 315be17..19516e3 100644 --- a/lib/Numeric/LinearAlgebra/LAPACK.hs +++ b/lib/Numeric/LinearAlgebra/LAPACK.hs | |||
@@ -61,7 +61,8 @@ svdAux f st x = unsafePerformIO $ do | |||
61 | u <- createMatrix ColumnMajor r r | 61 | u <- createMatrix ColumnMajor r r |
62 | s <- createVector (min r c) | 62 | s <- createVector (min r c) |
63 | v <- createMatrix ColumnMajor c c | 63 | v <- createMatrix ColumnMajor c c |
64 | f // matf x // matf u // vec s // matf v // check st [fdat x] | 64 | ww4 withMatrix x withMatrix u withVector s withMatrix v $ \x u s v -> |
65 | f // x // u // s // v // check st | ||
65 | return (u,s,trans v) | 66 | return (u,s,trans v) |
66 | where r = rows x | 67 | where r = rows x |
67 | c = cols x | 68 | c = cols x |
@@ -73,7 +74,8 @@ eigAux f st m | |||
73 | l <- createVector r | 74 | l <- createVector r |
74 | v <- createMatrix ColumnMajor r r | 75 | v <- createMatrix ColumnMajor r r |
75 | dummy <- createMatrix ColumnMajor 1 1 | 76 | dummy <- createMatrix ColumnMajor 1 1 |
76 | f // matf m // matf dummy // vec l // matf v // check st [fdat m] | 77 | ww4 withMatrix m withMatrix dummy withVector l withMatrix v $ \m dummy l v -> |
78 | f // m // dummy // l // v // check st | ||
77 | return (l,v) | 79 | return (l,v) |
78 | where r = rows m | 80 | where r = rows m |
79 | 81 | ||
@@ -115,7 +117,8 @@ eigRaux m | |||
115 | l <- createVector r | 117 | l <- createVector r |
116 | v <- createMatrix ColumnMajor r r | 118 | v <- createMatrix ColumnMajor r r |
117 | dummy <- createMatrix ColumnMajor 1 1 | 119 | dummy <- createMatrix ColumnMajor 1 1 |
118 | dgeev // matf m // matf dummy // vec l // matf v // check "eigR" [fdat m] | 120 | ww4 withMatrix m withMatrix dummy withVector l withMatrix v $ \m dummy l v -> |
121 | dgeev // m // dummy // l // v // check "eigR" | ||
119 | return (l,v) | 122 | return (l,v) |
120 | where r = rows m | 123 | where r = rows m |
121 | 124 | ||
@@ -144,7 +147,8 @@ eigS' m | |||
144 | | otherwise = unsafePerformIO $ do | 147 | | otherwise = unsafePerformIO $ do |
145 | l <- createVector r | 148 | l <- createVector r |
146 | v <- createMatrix ColumnMajor r r | 149 | v <- createMatrix ColumnMajor r r |
147 | dsyev // matf m // vec l // matf v // check "eigS" [fdat m] | 150 | ww3 withMatrix m withVector l withMatrix v $ \m l v -> |
151 | dsyev // m // l // v // check "eigS" | ||
148 | return (l,v) | 152 | return (l,v) |
149 | where r = rows m | 153 | where r = rows m |
150 | 154 | ||
@@ -166,7 +170,8 @@ eigH' m | |||
166 | | otherwise = unsafePerformIO $ do | 170 | | otherwise = unsafePerformIO $ do |
167 | l <- createVector r | 171 | l <- createVector r |
168 | v <- createMatrix ColumnMajor r r | 172 | v <- createMatrix ColumnMajor r r |
169 | zheev // matf m // vec l // matf v // check "eigH" [fdat m] | 173 | ww3 withMatrix m withVector l withMatrix v $ \m l v -> |
174 | zheev // m // l // v // check "eigH" | ||
170 | return (l,v) | 175 | return (l,v) |
171 | where r = rows m | 176 | where r = rows m |
172 | 177 | ||
@@ -177,7 +182,8 @@ foreign import ccall "LAPACK/lapack-aux.h linearSolveC_l" zgesv :: TCMCMCM | |||
177 | linearSolveSQAux f st a b | 182 | linearSolveSQAux f st a b |
178 | | n1==n2 && n1==r = unsafePerformIO $ do | 183 | | n1==n2 && n1==r = unsafePerformIO $ do |
179 | s <- createMatrix ColumnMajor r c | 184 | s <- createMatrix ColumnMajor r c |
180 | f // matf a // matf b // matf s // check st [fdat a, fdat b] | 185 | ww3 withMatrix a withMatrix b withMatrix s $ \a b s -> |
186 | f // a // b // s // check st | ||
181 | return s | 187 | return s |
182 | | otherwise = error $ st ++ " of nonsquare matrix" | 188 | | otherwise = error $ st ++ " of nonsquare matrix" |
183 | where n1 = rows a | 189 | where n1 = rows a |
@@ -201,7 +207,8 @@ foreign import ccall "LAPACK/lapack-aux.h linearSolveSVDC_l" zgelss :: Double -> | |||
201 | 207 | ||
202 | linearSolveAux f st a b = unsafePerformIO $ do | 208 | linearSolveAux f st a b = unsafePerformIO $ do |
203 | r <- createMatrix ColumnMajor (max m n) nrhs | 209 | r <- createMatrix ColumnMajor (max m n) nrhs |
204 | f // matf a // matf b // matf r // check st [fdat a, fdat b] | 210 | ww3 withMatrix a withMatrix b withMatrix r $ \a b r -> |
211 | f // a // b // r // check st | ||
205 | return r | 212 | return r |
206 | where m = rows a | 213 | where m = rows a |
207 | n = cols a | 214 | n = cols a |
@@ -251,7 +258,8 @@ cholS = cholAux dpotrf "cholS" . fmat | |||
251 | 258 | ||
252 | cholAux f st a = unsafePerformIO $ do | 259 | cholAux f st a = unsafePerformIO $ do |
253 | r <- createMatrix ColumnMajor n n | 260 | r <- createMatrix ColumnMajor n n |
254 | f // matf a // matf r // check st [fdat a] | 261 | ww2 withMatrix a withMatrix r $ \a r -> |
262 | f // a // r // check st | ||
255 | return r | 263 | return r |
256 | where n = rows a | 264 | where n = rows a |
257 | 265 | ||
@@ -270,8 +278,8 @@ qrC = qrAux zgeqr2 "qrC" . fmat | |||
270 | qrAux f st a = unsafePerformIO $ do | 278 | qrAux f st a = unsafePerformIO $ do |
271 | r <- createMatrix ColumnMajor m n | 279 | r <- createMatrix ColumnMajor m n |
272 | tau <- createVector mn | 280 | tau <- createVector mn |
273 | withForeignPtr (fptr $ fdat $ a) $ \p -> | 281 | ww3 withMatrix a withMatrix r withVector tau $ \ a r tau -> |
274 | f m n p // vec tau // matf r // check st [fdat a] | 282 | f // a // tau // r // check st |
275 | return (r,tau) | 283 | return (r,tau) |
276 | where m = rows a | 284 | where m = rows a |
277 | n = cols a | 285 | n = cols a |
@@ -292,7 +300,8 @@ hessC = hessAux zgehrd "hessC" . fmat | |||
292 | hessAux f st a = unsafePerformIO $ do | 300 | hessAux f st a = unsafePerformIO $ do |
293 | r <- createMatrix ColumnMajor m n | 301 | r <- createMatrix ColumnMajor m n |
294 | tau <- createVector (mn-1) | 302 | tau <- createVector (mn-1) |
295 | f // matf a // vec tau // matf r // check st [fdat a] | 303 | ww3 withMatrix a withMatrix r withVector tau $ \ a r tau -> |
304 | f // a // tau // r // check st | ||
296 | return (r,tau) | 305 | return (r,tau) |
297 | where m = rows a | 306 | where m = rows a |
298 | n = cols a | 307 | n = cols a |
@@ -313,7 +322,8 @@ schurC = schurAux zgees "schurC" . fmat | |||
313 | schurAux f st a = unsafePerformIO $ do | 322 | schurAux f st a = unsafePerformIO $ do |
314 | u <- createMatrix ColumnMajor n n | 323 | u <- createMatrix ColumnMajor n n |
315 | s <- createMatrix ColumnMajor n n | 324 | s <- createMatrix ColumnMajor n n |
316 | f // matf a // matf u // matf s // check st [fdat a] | 325 | ww3 withMatrix a withMatrix u withMatrix s $ \ a u s -> |
326 | f // a // u // s // check st | ||
317 | return (u,s) | 327 | return (u,s) |
318 | where n = rows a | 328 | where n = rows a |
319 | 329 | ||