summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--README62
-rw-r--r--examples/speed.hs20
-rw-r--r--examples/speed.m7
-rw-r--r--examples/tests.hs5
-rw-r--r--lib/Data/Packed/Internal/Matrix.hs12
-rw-r--r--lib/Data/Packed/Matrix.hs8
6 files changed, 77 insertions, 37 deletions
diff --git a/README b/README
index e2b7aa3..f117a26 100644
--- a/README
+++ b/README
@@ -1,3 +1,6 @@
1A simple scientific library for Haskell
2---------------------------------------
3
1REQUIREMENTS 4REQUIREMENTS
2 5
31) GNU Scientific Library (http://www.gnu.org/software/gsl) development packages 61) GNU Scientific Library (http://www.gnu.org/software/gsl) development packages
@@ -7,19 +10,20 @@ REQUIREMENTS
7 10
8INSTALLATION 11INSTALLATION
9 12
10(More detailed information is included in the "tutorial",
11available in the web page of the project.)
12
13$ runhaskell Setup.lhs configure --prefix=$HOME 13$ runhaskell Setup.lhs configure --prefix=$HOME
14$ runhaskell Setup.lhs build 14$ runhaskell Setup.lhs build
15$ runhaskell Setup.lhs haddock 15$ runhaskell Setup.lhs haddock
16$ runhaskell Setup.lhs install --user 16$ runhaskell Setup.lhs install --user
17 17
18See below for installation on Windows.
19
18USING ATLAS 20USING ATLAS
19 21
20$ ln -s /usr/lib/atlas/libblas.so.3.0 $HOME/lib/hssl-0.1/ghc-6.6.1/libcblas.so 22$ ln -s /usr/lib/atlas/libblas.so.3.0 $HOME/lib/hssl-0.1/ghc-6.6.1/libcblas.so
21$ ln -s /usr/lib/atlas/liblapack.so.3.0 $HOME/lib/hssl-0.1/ghc-6.6.1/liblapack.so 23$ ln -s /usr/lib/atlas/liblapack.so.3.0 $HOME/lib/hssl-0.1/ghc-6.6.1/liblapack.so
22 24
25(More info in the tutorial, available from the web page of the project.)
26
23TESTS 27TESTS
24 28
25$ runhaskell examples/tests.hs 29$ runhaskell examples/tests.hs
@@ -46,10 +50,12 @@ Prelude Numeric.LinearAlgebra> u <> d <> trans v
46Prelude Numeric.GSL> :q 50Prelude Numeric.GSL> :q
47Leaving GHCi. 51Leaving GHCi.
48 52
53A few illustrative programs are included in the examples folder.
54
49CHANGES 55CHANGES
50 56
51This is a new version of GSLHaskell. The package is (provisionally) 57This is a new version of GSLHaskell. The package is provisionally
52called \C{hssl} (a simple scientific library for Haskell) because only 58called "hssl" (a simple scientific library for Haskell) because only
53a small part of GSL is available and linear algebra is based on LAPACK. 59a small part of GSL is available and linear algebra is based on LAPACK.
54 60
55The code has been extensively refactored. There is a new internal representation 61The code has been extensively refactored. There is a new internal representation
@@ -57,35 +63,18 @@ which admits both C and Fortran matrices and avoids many transposes.
57 63
58There are only minor API changes: 64There are only minor API changes:
59 65
60- the matrix product operator \C{(<>)} is now overloaded only for matrix-matrix, 66- The matrix product operator (<>) is now overloaded only for matrix-matrix,
61matrix-vector and vector-matrix, with the same base type. The dot product and the scaling 67 matrix-vector and vector-matrix, with the same base type. Dot product and scaling
62of vectors or matrices is now denoted by `dot` and `scale`. Conversions from real to 68 of vectors or matrices is now denoted by `dot` or (<.>) and `scale` or (.*).
63complex objects must be explicit. 69 Conversions from real to complex objects must now be explicit.
64 70
65- Most linear algebra functions admit both real and complex objects. Utilities such as 71- Most linear algebra functions admit both real and complex objects. Utilities such as
66ident or constant are now polymorphic. 72 ident or constant are now polymorphic.
67 73
68- Runtime errors produced by GSL or LAPACK can be handled using \C{Control.Exeception.catch}. 74- Runtime errors produced by GSL or LAPACK can be handled using Control.Exeception.catch.
69 75
70Old GSLHaskell code will work with small modifications. 76Old GSLHaskell code will work with small modifications.
71 77
72ACKNOWLEDGEMENTS
73
74I thank Henning Thielemann and all the people in the Haskell mailing lists for their help.
75
76- Nico Mahlo discovered a bug in the eigendecomposition wrapper.
77
78- Frederik Eaton discovered a bug in the design of the wrappers.
79
80- Eric Kidd has created a wiki page explaining the installation on MacOS X:
81 http://www.haskell.org/haskellwiki/GSLHaskell_on_MacOS_X
82
83- Fawzi Mohamed discovered a portability bug in the lapack wrappers.
84
85- Pedro E. López de Teruel fixed the interface to lapack.
86
87- Antti Siira discovered a bug in the plotting functions.
88
89INSTALLATION ON WINDOWS 78INSTALLATION ON WINDOWS
90 79
911) Download the developer files gsl-1.8-lib.zip from 801) Download the developer files gsl-1.8-lib.zip from
@@ -109,3 +98,20 @@ INSTALLATION ON WINDOWS
109Unfortunately the lapack dll supplied by the R system does not include zgels_ and zgelss_, 98Unfortunately the lapack dll supplied by the R system does not include zgels_ and zgelss_,
110so the functions depending on them (linearSolveLS and linearSolveSVD for complex data) 99so the functions depending on them (linearSolveLS and linearSolveSVD for complex data)
111will produce a "non supported" runtime error. 100will produce a "non supported" runtime error.
101
102ACKNOWLEDGEMENTS
103
104I thank Henning Thielemann and all the people in the Haskell mailing lists for their help.
105
106- Nico Mahlo discovered a bug in the eigendecomposition wrapper.
107
108- Frederik Eaton discovered a bug in the design of the wrappers.
109
110- Eric Kidd has created a wiki page explaining the installation on MacOS X:
111 http://www.haskell.org/haskellwiki/GSLHaskell_on_MacOS_X
112
113- Fawzi Mohamed discovered a portability bug in the lapack wrappers.
114
115- Pedro E. López de Teruel fixed the interface to lapack.
116
117- Antti Siira discovered a bug in the plotting functions.
diff --git a/examples/speed.hs b/examples/speed.hs
index a937f31..865e3e1 100644
--- a/examples/speed.hs
+++ b/examples/speed.hs
@@ -1,3 +1,23 @@
1{- speed tests
2
3$ ghc --make -O speed
4In my machine:
5$ ./speed 5 100000 1
6(3><3)
7 [ -1.7877285611885504e-2, 0.0, -0.9998401885597121
8 , 0.0, 1.0, 0.0
9 , 0.9998401885597168, 0.0, -1.7877285611891697e-2 ]
100.29 CPU seconds
11
12GNU-Octave:
13./speed.m
14 -0.017877255967426 0.000000000000000 -0.999840189089781
15 0.000000000000000 1.000000000000000 0.000000000000000
16 0.999840189089763 0.000000000000000 -0.017877255967417
179.69 seconds
18
19-}
20
1import Numeric.LinearAlgebra 21import Numeric.LinearAlgebra
2import System 22import System
3import Data.List(foldl1') 23import Data.List(foldl1')
diff --git a/examples/speed.m b/examples/speed.m
index 10dbf78..2f41665 100644
--- a/examples/speed.m
+++ b/examples/speed.m
@@ -12,9 +12,10 @@ end
12t0=time(); 12t0=time();
13x = linspace(0,1,1E5); 13x = linspace(0,1,1E5);
14ac = eye(3); 14ac = eye(3);
15for a=x 15for a = x
16 ac = rot(a)*ac; 16 ac = ac*rot(a);
17end 17end
18 18
19format long
19disp(ac); 20disp(ac);
20disp(time()-t0) 21printf("%.2f seconds\n",time()-t0)
diff --git a/examples/tests.hs b/examples/tests.hs
index b44d140..9388671 100644
--- a/examples/tests.hs
+++ b/examples/tests.hs
@@ -12,6 +12,7 @@ import qualified Numeric.GSL.Matrix as GSL
12import Test.QuickCheck hiding (test) 12import Test.QuickCheck hiding (test)
13import Test.HUnit hiding ((~:),test) 13import Test.HUnit hiding ((~:),test)
14import System.Random(randomRs,mkStdGen) 14import System.Random(randomRs,mkStdGen)
15import System.Info
15 16
16type RM = Matrix Double 17type RM = Matrix Double
17type CM = Matrix (Complex Double) 18type CM = Matrix (Complex Double)
@@ -323,7 +324,9 @@ tests = do
323 quickCheck (invTest . sqm :: SqM (Complex Double) -> Bool) 324 quickCheck (invTest . sqm :: SqM (Complex Double) -> Bool)
324 putStrLn "--------- pinv ------" 325 putStrLn "--------- pinv ------"
325 quickCheck (pinvTest . sqm :: SqM Double -> Bool) 326 quickCheck (pinvTest . sqm :: SqM Double -> Bool)
326 quickCheck (pinvTest . sqm :: SqM (Complex Double) -> Bool) 327 if os == "mingw32"
328 then putStrLn "complex pinvTest skipped in this OS"
329 else quickCheck (pinvTest . sqm :: SqM (Complex Double) -> Bool)
327 putStrLn "--------- chol ------" 330 putStrLn "--------- chol ------"
328 runTestTT $ TestList 331 runTestTT $ TestList
329 [ test "cholR" cholRTest 332 [ test "cholR" cholRTest
diff --git a/lib/Data/Packed/Internal/Matrix.hs b/lib/Data/Packed/Internal/Matrix.hs
index 77906dc..605a6f3 100644
--- a/lib/Data/Packed/Internal/Matrix.hs
+++ b/lib/Data/Packed/Internal/Matrix.hs
@@ -48,12 +48,14 @@ import Data.List(transpose)
48- we could carry both the matrix and its (lazily computed) transpose. 48- we could carry both the matrix and its (lazily computed) transpose.
49 This may save some transpositions, but it is necessary to keep track of the 49 This may save some transpositions, but it is necessary to keep track of the
50 data which is actually computed to be used by functions like the matrix product 50 data which is actually computed to be used by functions like the matrix product
51 which admit both orders. Therefore, maybe it is better to have something like 51 which admit both orders.
52 viewC and viewF, which may actually perform a transpose if required.
53 52
54- but if we need the transposed data and it is not in the structure, we must make 53- but if we need the transposed data and it is not in the structure, we must make
55 sure that we touch the same foreignptr that is used in the computation. Access 54 sure that we touch the same foreignptr that is used in the computation.
56 to such pointer cannot be made by creating a new vector. 55
56- a reasonable solution is using two constructors for a matrix. Transposition just
57 "flips" the constructor. Actual data transposition is not done if followed by a
58 matrix product or another transpose.
57 59
58-} 60-}
59 61
@@ -66,7 +68,7 @@ data Matrix t = MC { rows :: Int, cols :: Int, cdat :: Vector t, fdat :: Vector
66-- MC: preferred by C, fdat may require a transposition 68-- MC: preferred by C, fdat may require a transposition
67-- MF: preferred by LAPACK, cdat may require a transposition 69-- MF: preferred by LAPACK, cdat may require a transposition
68 70
69-- | matrix transpose 71-- | Matrix transpose.
70trans :: Matrix t -> Matrix t 72trans :: Matrix t -> Matrix t
71trans MC {rows = r, cols = c, cdat = d, fdat = dt } = MF {rows = c, cols = r, fdat = d, cdat = dt } 73trans MC {rows = r, cols = c, cdat = d, fdat = dt } = MF {rows = c, cols = r, fdat = d, cdat = dt }
72trans MF {rows = r, cols = c, fdat = d, cdat = dt } = MC {rows = c, cols = r, cdat = d, fdat = dt } 74trans MF {rows = r, cols = c, fdat = d, cdat = dt } = MC {rows = c, cols = r, cdat = d, fdat = dt }
diff --git a/lib/Data/Packed/Matrix.hs b/lib/Data/Packed/Matrix.hs
index e94d038..a705975 100644
--- a/lib/Data/Packed/Matrix.hs
+++ b/lib/Data/Packed/Matrix.hs
@@ -193,6 +193,14 @@ dsp' sep as = unlines . map unwords' $ transpose mtp where
193 pad n str = replicate (n - length str) ' ' ++ str 193 pad n str = replicate (n - length str) ' ' ++ str
194 unwords' = concat . intersperse sep 194 unwords' = concat . intersperse sep
195 195
196{- | Creates a string from a matrix given a separator and a function to show each entry. Using
197this function the user can easily define any desired display function:
198
199@import Text.Printf(printf)@
200
201@disp = putStrLn . format \" \" (printf \"%.2f\")@
202
203-}
196format :: (Field t) => String -> (t -> String) -> Matrix t -> String 204format :: (Field t) => String -> (t -> String) -> Matrix t -> String
197format sep f m = dsp' sep . map (map f) . toLists $ m 205format sep f m = dsp' sep . map (map f) . toLists $ m
198 206