1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
|
{- speed tests
$ ghc --make -O speed
In my machine:
$ ./speed 5 100000 1
(3><3)
[ -1.7877285611885504e-2, 0.0, -0.9998401885597121
, 0.0, 1.0, 0.0
, 0.9998401885597168, 0.0, -1.7877285611891697e-2 ]
0.29 CPU seconds
GNU-Octave:
./speed.m
-0.017877255967426 0.000000000000000 -0.999840189089781
0.000000000000000 1.000000000000000 0.000000000000000
0.999840189089763 0.000000000000000 -0.017877255967417
9.69 seconds
-}
import Numeric.LinearAlgebra
import System
import Data.List(foldl1')
import System.CPUTime
import Text.Printf
import Debug.Trace
debug x = trace (show x) x
timing act = do
t0 <- getCPUTime
act
t1 <- getCPUTime
printf "%.2f CPU seconds\n" $ (fromIntegral ((t1 - t0) `div` (10^10)) / 100 :: Double)
op a b = trans $ (trans a) <> (trans b)
op2 a b = trans $ (trans a) + (trans b)
rot' :: Double -> Matrix Double
rot' a = matrix [[ c,0,s],
[ 0,1,0],
[-s,0,c]]
where c = cos a
s = sin a
matrix = fromLists
rot :: Double -> Matrix Double
rot a = (3><3) [ c,0,s
, 0,1,0
,-s,0,c ]
where c = cos a
s = sin a
fun n r = foldl1' (<>) (map r [0,delta..1]) where delta = 1 /(fromIntegral n-1)
main = do
args <- getArgs
let [p,n,d] = map read args
let ms = replicate n ((ident d :: Matrix Double))
let mz = replicate n (diag (constant (0::Double) d))
timing $ case p of
0 -> print $ foldl1' (<>) ms
1 -> print $ foldl1' (<>) (map trans ms)
2 -> print $ foldl1' op ms
3 -> print $ foldl1' op2 mz
4 -> print $ fun n rot'
5 -> print $ fun n rot
|
