diff options
Diffstat (limited to 'lib/GSLHaskell.hs')
-rw-r--r-- | lib/GSLHaskell.hs | 327 |
1 files changed, 0 insertions, 327 deletions
diff --git a/lib/GSLHaskell.hs b/lib/GSLHaskell.hs deleted file mode 100644 index 254a957..0000000 --- a/lib/GSLHaskell.hs +++ /dev/null | |||
@@ -1,327 +0,0 @@ | |||
1 | {-# OPTIONS_GHC -fglasgow-exts -fallow-undecidable-instances #-} | ||
2 | ----------------------------------------------------------------------------- | ||
3 | {- | | ||
4 | Module : GSLHaskell | ||
5 | Copyright : (c) Alberto Ruiz 2006 | ||
6 | License : GPL-style | ||
7 | |||
8 | Maintainer : Alberto Ruiz (aruiz at um dot es) | ||
9 | Stability : provisional | ||
10 | Portability : uses -fffi and -fglasgow-exts | ||
11 | |||
12 | Old GSLHaskell interface. | ||
13 | |||
14 | -} | ||
15 | ----------------------------------------------------------------------------- | ||
16 | |||
17 | module GSLHaskell( | ||
18 | module Data.Packed.Vector, | ||
19 | module Data.Packed.Matrix, | ||
20 | module LinearAlgebra.Algorithms, | ||
21 | module LinearAlgebra.Linear, | ||
22 | module LAPACK, | ||
23 | module GSL.Integration, | ||
24 | module GSL.Differentiation, | ||
25 | module GSL.Fourier, | ||
26 | module GSL.Polynomials, | ||
27 | module GSL.Minimization, | ||
28 | module GSL.Matrix, | ||
29 | module GSL.Special, | ||
30 | module Graphics.Plot, | ||
31 | module Complex, | ||
32 | Mul,(<>), readMatrix, size, dispR, dispC, format, gmap, Joinable, (<|>),(<->), | ||
33 | fromArray2D, GSLHaskell.pnorm, | ||
34 | ) where | ||
35 | |||
36 | import GSL.Integration | ||
37 | import GSL.Differentiation | ||
38 | import GSL.Fourier | ||
39 | import GSL.Polynomials | ||
40 | import GSL.Minimization | ||
41 | import Graphics.Plot | ||
42 | import Complex | ||
43 | import GSL.Special(setErrorHandlerOff, | ||
44 | erf, | ||
45 | erf_Z, | ||
46 | bessel_J0_e, | ||
47 | exp_e10_e, | ||
48 | gamma) | ||
49 | import Data.Packed.Vector | ||
50 | import Data.Packed.Matrix | ||
51 | import Data.Packed.Matrix hiding ((><)) | ||
52 | import GSL.Vector | ||
53 | import qualified LinearAlgebra.Algorithms | ||
54 | import LAPACK | ||
55 | import GSL.Matrix | ||
56 | import LinearAlgebra.Algorithms hiding (pnorm) | ||
57 | import LinearAlgebra.Linear hiding (Mul,(<>)) | ||
58 | import Data.Packed.Internal.Matrix(multiply) | ||
59 | import Complex | ||
60 | import Numeric(showGFloat) | ||
61 | import Data.List(transpose,intersperse) | ||
62 | import Foreign(Storable) | ||
63 | import Data.Array | ||
64 | import LinearAlgebra.Instances | ||
65 | |||
66 | |||
67 | class Mul a b c | a b -> c where | ||
68 | infixl 7 <> | ||
69 | {- | An overloaded operator for matrix products, matrix-vector and vector-matrix products, dot products and scaling of vectors and matrices. Type consistency is statically checked. Alternatively, you can use the specific functions described below, but using this operator you can automatically combine real and complex objects. | ||
70 | |||
71 | @v = 'fromList' [1,2,3] :: Vector Double | ||
72 | cv = 'fromList' [1+'i',2] | ||
73 | m = 'fromLists' [[1,2,3], | ||
74 | [4,5,7]] :: Matrix Double | ||
75 | cm = 'fromLists' [[ 1, 2], | ||
76 | [3+'i',7*'i'], | ||
77 | [ 'i', 1]] | ||
78 | \ | ||
79 | \> m \<\> v | ||
80 | 14. 35. | ||
81 | \ | ||
82 | \> cv \<\> m | ||
83 | 9.+1.i 12.+2.i 17.+3.i | ||
84 | \ | ||
85 | \> m \<\> cm | ||
86 | 7.+5.i 5.+14.i | ||
87 | 19.+12.i 15.+35.i | ||
88 | \ | ||
89 | \> v \<\> 'i' | ||
90 | 1.i 2.i 3.i | ||
91 | \ | ||
92 | \> v \<\> v | ||
93 | 14.0 | ||
94 | \ | ||
95 | \> cv \<\> cv | ||
96 | 4.0 :+ 2.0@ | ||
97 | |||
98 | -} | ||
99 | (<>) :: a -> b -> c | ||
100 | |||
101 | |||
102 | instance Mul Double Double Double where | ||
103 | (<>) = (*) | ||
104 | |||
105 | instance Mul Double (Complex Double) (Complex Double) where | ||
106 | a <> b = (a:+0) * b | ||
107 | |||
108 | instance Mul (Complex Double) Double (Complex Double) where | ||
109 | a <> b = a * (b:+0) | ||
110 | |||
111 | instance Mul (Complex Double) (Complex Double) (Complex Double) where | ||
112 | (<>) = (*) | ||
113 | |||
114 | --------------------------------- matrix matrix | ||
115 | |||
116 | instance Mul (Matrix Double) (Matrix Double) (Matrix Double) where | ||
117 | (<>) = multiply | ||
118 | |||
119 | instance Mul (Matrix (Complex Double)) (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
120 | (<>) = multiply | ||
121 | |||
122 | instance Mul (Matrix (Complex Double)) (Matrix Double) (Matrix (Complex Double)) where | ||
123 | c <> r = c <> liftMatrix comp r | ||
124 | |||
125 | instance Mul (Matrix Double) (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
126 | r <> c = liftMatrix comp r <> c | ||
127 | |||
128 | --------------------------------- (Matrix Double) (Vector Double) | ||
129 | |||
130 | instance Mul (Matrix Double) (Vector Double) (Vector Double) where | ||
131 | (<>) = mXv | ||
132 | |||
133 | instance Mul (Matrix (Complex Double)) (Vector (Complex Double)) (Vector (Complex Double)) where | ||
134 | (<>) = mXv | ||
135 | |||
136 | instance Mul (Matrix (Complex Double)) (Vector Double) (Vector (Complex Double)) where | ||
137 | m <> v = m <> comp v | ||
138 | |||
139 | instance Mul (Matrix Double) (Vector (Complex Double)) (Vector (Complex Double)) where | ||
140 | m <> v = liftMatrix comp m <> v | ||
141 | |||
142 | --------------------------------- (Vector Double) (Matrix Double) | ||
143 | |||
144 | instance Mul (Vector Double) (Matrix Double) (Vector Double) where | ||
145 | (<>) = vXm | ||
146 | |||
147 | instance Mul (Vector (Complex Double)) (Matrix (Complex Double)) (Vector (Complex Double)) where | ||
148 | (<>) = vXm | ||
149 | |||
150 | instance Mul (Vector (Complex Double)) (Matrix Double) (Vector (Complex Double)) where | ||
151 | v <> m = v <> liftMatrix comp m | ||
152 | |||
153 | instance Mul (Vector Double) (Matrix (Complex Double)) (Vector (Complex Double)) where | ||
154 | v <> m = comp v <> m | ||
155 | |||
156 | --------------------------------- dot product | ||
157 | |||
158 | instance Mul (Vector Double) (Vector Double) Double where | ||
159 | (<>) = dot | ||
160 | |||
161 | instance Mul (Vector (Complex Double)) (Vector (Complex Double)) (Complex Double) where | ||
162 | (<>) = dot | ||
163 | |||
164 | instance Mul (Vector Double) (Vector (Complex Double)) (Complex Double) where | ||
165 | a <> b = comp a <> b | ||
166 | |||
167 | instance Mul (Vector (Complex Double)) (Vector Double) (Complex Double) where | ||
168 | (<>) = flip (<>) | ||
169 | |||
170 | --------------------------------- scaling vectors | ||
171 | |||
172 | instance Mul Double (Vector Double) (Vector Double) where | ||
173 | (<>) = scale | ||
174 | |||
175 | instance Mul (Vector Double) Double (Vector Double) where | ||
176 | (<>) = flip (<>) | ||
177 | |||
178 | instance Mul (Complex Double) (Vector (Complex Double)) (Vector (Complex Double)) where | ||
179 | (<>) = scale | ||
180 | |||
181 | instance Mul (Vector (Complex Double)) (Complex Double) (Vector (Complex Double)) where | ||
182 | (<>) = flip (<>) | ||
183 | |||
184 | instance Mul Double (Vector (Complex Double)) (Vector (Complex Double)) where | ||
185 | a <> v = (a:+0) <> v | ||
186 | |||
187 | instance Mul (Vector (Complex Double)) Double (Vector (Complex Double)) where | ||
188 | (<>) = flip (<>) | ||
189 | |||
190 | instance Mul (Complex Double) (Vector Double) (Vector (Complex Double)) where | ||
191 | a <> v = a <> comp v | ||
192 | |||
193 | instance Mul (Vector Double) (Complex Double) (Vector (Complex Double)) where | ||
194 | (<>) = flip (<>) | ||
195 | |||
196 | --------------------------------- scaling matrices | ||
197 | |||
198 | instance Mul Double (Matrix Double) (Matrix Double) where | ||
199 | (<>) a = liftMatrix (a <>) | ||
200 | |||
201 | instance Mul (Matrix Double) Double (Matrix Double) where | ||
202 | (<>) = flip (<>) | ||
203 | |||
204 | instance Mul (Complex Double) (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
205 | (<>) a = liftMatrix (a <>) | ||
206 | |||
207 | instance Mul (Matrix (Complex Double)) (Complex Double) (Matrix (Complex Double)) where | ||
208 | (<>) = flip (<>) | ||
209 | |||
210 | instance Mul Double (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
211 | a <> m = (a:+0) <> m | ||
212 | |||
213 | instance Mul (Matrix (Complex Double)) Double (Matrix (Complex Double)) where | ||
214 | (<>) = flip (<>) | ||
215 | |||
216 | instance Mul (Complex Double) (Matrix Double) (Matrix (Complex Double)) where | ||
217 | a <> m = a <> liftMatrix comp m | ||
218 | |||
219 | instance Mul (Matrix Double) (Complex Double) (Matrix (Complex Double)) where | ||
220 | (<>) = flip (<>) | ||
221 | |||
222 | ----------------------------------------------------------------------------------- | ||
223 | |||
224 | size :: Vector a -> Int | ||
225 | size = dim | ||
226 | |||
227 | gmap :: (Storable a, Storable b) => (a->b) -> Vector a -> Vector b | ||
228 | gmap f v = liftVector f v | ||
229 | |||
230 | -- shows a Double with n digits after the decimal point | ||
231 | shf :: (RealFloat a) => Int -> a -> String | ||
232 | shf dec n | abs n < 1e-10 = "0." | ||
233 | | abs (n - (fromIntegral.round $ n)) < 1e-10 = show (round n) ++"." | ||
234 | | otherwise = showGFloat (Just dec) n "" | ||
235 | -- shows a Complex Double as a pair, with n digits after the decimal point | ||
236 | shfc n z@ (a:+b) | ||
237 | | magnitude z <1e-10 = "0." | ||
238 | | abs b < 1e-10 = shf n a | ||
239 | | abs a < 1e-10 = shf n b ++"i" | ||
240 | | b > 0 = shf n a ++"+"++shf n b ++"i" | ||
241 | | otherwise = shf n a ++shf n b ++"i" | ||
242 | |||
243 | dsp :: String -> [[String]] -> String | ||
244 | dsp sep as = unlines . map unwords' $ transpose mtp where | ||
245 | mt = transpose as | ||
246 | longs = map (maximum . map length) mt | ||
247 | mtp = zipWith (\a b -> map (pad a) b) longs mt | ||
248 | pad n str = replicate (n - length str) ' ' ++ str | ||
249 | unwords' = concat . intersperse sep | ||
250 | |||
251 | format :: (Field t) => String -> (t -> String) -> Matrix t -> String | ||
252 | format sep f m = dsp sep . map (map f) . toLists $ m | ||
253 | |||
254 | disp m f = putStrLn $ "matrix ("++show (rows m) ++"x"++ show (cols m) ++")\n"++format " | " f m | ||
255 | |||
256 | dispR :: Int -> Matrix Double -> IO () | ||
257 | dispR d m = disp m (shf d) | ||
258 | |||
259 | dispC :: Int -> Matrix (Complex Double) -> IO () | ||
260 | dispC d m = disp m (shfc d) | ||
261 | |||
262 | -- | creates a matrix from a table of numbers. | ||
263 | readMatrix :: String -> Matrix Double | ||
264 | readMatrix = fromLists . map (map read). map words . filter (not.null) . lines | ||
265 | |||
266 | ------------------------------------------------------------- | ||
267 | |||
268 | class Joinable a b c | a b -> c where | ||
269 | joinH :: a -> b -> c | ||
270 | joinV :: a -> b -> c | ||
271 | |||
272 | instance Joinable (Matrix Double) (Vector Double) (Matrix Double) where | ||
273 | joinH m v = fromBlocks [[m,reshape 1 v]] | ||
274 | joinV m v = fromBlocks [[m],[reshape (size v) v]] | ||
275 | |||
276 | instance Joinable (Vector Double) (Matrix Double) (Matrix Double) where | ||
277 | joinH v m = fromBlocks [[reshape 1 v,m]] | ||
278 | joinV v m = fromBlocks [[reshape (size v) v],[m]] | ||
279 | |||
280 | instance Joinable (Matrix Double) (Matrix Double) (Matrix Double) where | ||
281 | joinH m1 m2 = fromBlocks [[m1,m2]] | ||
282 | joinV m1 m2 = fromBlocks [[m1],[m2]] | ||
283 | |||
284 | instance Joinable (Matrix (Complex Double)) (Vector (Complex Double)) (Matrix (Complex Double)) where | ||
285 | joinH m v = fromBlocks [[m,reshape 1 v]] | ||
286 | joinV m v = fromBlocks [[m],[reshape (size v) v]] | ||
287 | |||
288 | instance Joinable (Vector (Complex Double)) (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
289 | joinH v m = fromBlocks [[reshape 1 v,m]] | ||
290 | joinV v m = fromBlocks [[reshape (size v) v],[m]] | ||
291 | |||
292 | instance Joinable (Matrix (Complex Double)) (Matrix (Complex Double)) (Matrix (Complex Double)) where | ||
293 | joinH m1 m2 = fromBlocks [[m1,m2]] | ||
294 | joinV m1 m2 = fromBlocks [[m1],[m2]] | ||
295 | |||
296 | infixl 3 <|>, <-> | ||
297 | |||
298 | {- | Horizontal concatenation of matrices and vectors: | ||
299 | |||
300 | @\> 'ident' 3 \<-\> i\<\>'ident' 3 \<|\> 'fromList' [1..6] | ||
301 | 1. 0. 0. 1. | ||
302 | 0. 1. 0. 2. | ||
303 | 0. 0. 1. 3. | ||
304 | 1.i 0. 0. 4. | ||
305 | 0. 1.i 0. 5. | ||
306 | 0. 0. 1.i 6.@ | ||
307 | -} | ||
308 | (<|>) :: (Joinable a b c) => a -> b -> c | ||
309 | a <|> b = joinH a b | ||
310 | |||
311 | -- | Vertical concatenation of matrices and vectors. | ||
312 | (<->) :: (Joinable a b c) => a -> b -> c | ||
313 | a <-> b = joinV a b | ||
314 | |||
315 | ---------------------------------------------------------- | ||
316 | |||
317 | fromArray2D :: (Field e) => Array (Int, Int) e -> Matrix e | ||
318 | fromArray2D m = (r><c) (elems m) | ||
319 | where ((r0,c0),(r1,c1)) = bounds m | ||
320 | r = r1-r0+1 | ||
321 | c = c1-c0+1 | ||
322 | |||
323 | |||
324 | pnorm :: (Normed t1, Num t) => t -> t1 -> Double | ||
325 | pnorm 0 = LinearAlgebra.Algorithms.pnorm Infinity | ||
326 | pnorm 1 = LinearAlgebra.Algorithms.pnorm PNorm1 | ||
327 | pnorm 2 = LinearAlgebra.Algorithms.pnorm PNorm2 \ No newline at end of file | ||