summaryrefslogtreecommitdiff
path: root/examples/pru.hs
diff options
context:
space:
mode:
Diffstat (limited to 'examples/pru.hs')
-rw-r--r--examples/pru.hs45
1 files changed, 25 insertions, 20 deletions
diff --git a/examples/pru.hs b/examples/pru.hs
index fb33962..f855bce 100644
--- a/examples/pru.hs
+++ b/examples/pru.hs
@@ -52,12 +52,14 @@ main = do
52 print $ foldl part t [("p",1),("r",2),("q",0)] 52 print $ foldl part t [("p",1),("r",2),("q",0)]
53 print $ foldl part t $ reverse [("p",1),("r",2),("q",0)] 53 print $ foldl part t $ reverse [("p",1),("r",2),("q",0)]
54 54
55t = T [(4,(Covariant,"p")),(2,(Covariant,"q")),(3,(Contravariant,"r"))] $ fromList [1..24::Double] 55t = T [IdxDesc 4 Covariant "p",IdxDesc 2 Covariant "q" ,IdxDesc 3 Contravariant "r"]
56 $ fromList [1..24::Double]
56 57
57 58
58 59t1 = T [IdxDesc 4 Covariant "p",IdxDesc 4 Contravariant "q" ,IdxDesc 2 Covariant "r"]
59t1 = T [(4,(Covariant,"p")),(4,(Contravariant,"q")),(2,(Covariant,"r"))] $ fromList [1..32::Double] 60 $ fromList [1..32::Double]
60t2 = T [(4,(Covariant,"p")),(4,(Contravariant,"q"))] $ fromList [1..16::Double] 61t2 = T [IdxDesc 4 Covariant "p",IdxDesc 4 Contravariant "q"]
62 $ fromList [1..16::Double]
61 63
62 64
63 65
@@ -72,15 +74,18 @@ e i n = fromList [ delta k i | k <- [1..n]]
72diagl = diag.fromList 74diagl = diag.fromList
73 75
74scalar x = T [] (fromList [x]) 76scalar x = T [] (fromList [x])
75tensorFromVector idx v = T {dims = [(dim v,idx)], ten = v} 77tensorFromVector (tp,nm) v = T {dims = [IdxDesc (dim v) tp nm]
76tensorFromMatrix idxr idxc m = T {dims = [(rows m,idxr),(cols m,idxc)], ten = cdat m} 78 , ten = v}
79tensorFromMatrix (tpr,nmr) (tpc,nmc) m = T {dims = [IdxDesc (rows m) tpr nmr,IdxDesc (cols m) tpc nmc]
80 , ten = cdat m}
77 81
78td = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ diagl [1..4] :: Tensor Double 82td = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ diagl [1..4] :: Tensor Double
79 83
80tn = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ (2><3) [1..6] :: Tensor Double 84tn = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ (2><3) [1..6] :: Tensor Double
81tt = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ (2><3) [1..6] :: Tensor Double 85tt = tensorFromMatrix (Contravariant,"i") (Covariant,"j") $ (2><3) [1..6] :: Tensor Double
82 86
83tq = T [(3,(Covariant,"p")),(2,(Covariant,"q")),(2,(Covariant,"r"))] $ fromList [11 .. 22] :: Tensor Double 87tq = T [IdxDesc 3 Covariant "p",IdxDesc 2 Covariant "q" ,IdxDesc 2 Covariant "r"]
88 $ fromList [11 .. 22] :: Tensor Double
84 89
85r1 = contraction tt "j" tq "p" 90r1 = contraction tt "j" tq "p"
86r1' = contraction' tt "j" tq "p" 91r1' = contraction' tt "j" tq "p"
@@ -101,7 +106,7 @@ pru = do
101 print $ foldl part t2 [("j'",0),("p",1),("r",1)] 106 print $ foldl part t2 [("j'",0),("p",1),("r",1)]
102 107
103scsig t = scalar (signature (nms t)) `prod` t 108scsig t = scalar (signature (nms t)) `prod` t
104 where nms = map (snd.snd) . dims 109 where nms = map idxName . dims
105 110
106antisym' t = addT $ map (scsig . flip tridx t) (perms (names t)) 111antisym' t = addT $ map (scsig . flip tridx t) (perms (names t))
107 112
@@ -115,9 +120,9 @@ antisym' t = addT $ map (scsig . flip tridx t) (perms (names t))
115-} 120-}
116 121
117auxrename (T d v) = T d' v 122auxrename (T d v) = T d' v
118 where d' = [(n,(c,show (pos q))) | (n,(c,q)) <- d] 123 where d' = [IdxDesc n c (show (pos q)) | IdxDesc n c q <- d]
119 pos n = i where Just i = elemIndex n nms 124 pos n = i where Just i = elemIndex n nms
120 nms = map (snd.snd) d 125 nms = map idxName d
121 126
122antisym t = T (dims t) (ten (antisym' (auxrename t))) 127antisym t = T (dims t) (ten (antisym' (auxrename t)))
123 128
@@ -156,16 +161,16 @@ l1 = vector "p" [0,0,0,1]
156l2 = vector "q" [1,0,0,1] 161l2 = vector "q" [1,0,0,1]
157l3 = vector "r" [0,1,0,1] 162l3 = vector "r" [0,1,0,1]
158 163
159leviCivita n = antisym $ foldl1 prod $ zipWith tcovector (map show [1..]) (toRows (ident n)) 164leviCivita n = antisym $ foldl1 prod $ zipWith tcovector (map show [1,2..]) (toRows (ident n))
160 165
161contractionF t1 t2 = contraction t1 n1 t2 n2 166contractionF t1 t2 = contraction t1 n1 t2 n2
162 where n1 = fn t1 167 where n1 = fn t1
163 n2 = fn t2 168 n2 = fn t2
164 fn = snd . snd . head . dims 169 fn = idxName . head . dims
165 170
166 171
167dualV vs = foldl' contractionF (leviCivita n) vs 172dualV vs = foldl' contractionF (leviCivita n) vs
168 where n = fst . head . dims . head $ vs 173 where n = idxDim . head . dims . head $ vs
169 174
170 175
171dual1 = foldl' contractionF (leviCivita 3) [u,v] 176dual1 = foldl' contractionF (leviCivita 3) [u,v]
@@ -184,16 +189,16 @@ x3 = vector "r" [-3,-1,-1]
184x4 = vector "s" [12,0,3] 189x4 = vector "s" [12,0,3]
185 190
186raise (T d v) = T (map raise' d) v 191raise (T d v) = T (map raise' d) v
187 where raise' (n,(Covariant,s)) = (n,(Contravariant,s)) 192 where raise' idx@IdxDesc {idxType = Covariant } = idx {idxType = Contravariant}
188 raise' (n,(Contravariant,s)) = (n,(Covariant,s)) 193 raise' idx@IdxDesc {idxType = Contravariant } = idx {idxType = Covariant}
189 194
190dualMV t = prod (foldl' contract1b (lc <*> t) ds) (scalar (recip $ fromIntegral $ fact (length ds))) 195dualMV t = prod (foldl' contract1b (lc <*> t) ds) (scalar (recip $ fromIntegral $ fact (length ds)))
191 where 196 where
192 lc = leviCivita n 197 lc = leviCivita n
193 nms1 = map (snd.snd) (dims lc) 198 nms1 = map idxName (dims lc)
194 nms2 = map ((++"'").snd.snd) (dims t) 199 nms2 = map ((++"'").idxName) (dims t)
195 ds = zip nms1 nms2 200 ds = zip nms1 nms2
196 n = fst . head . dims $ t 201 n = idxDim . head . dims $ t
197 202
198-- intersection of two lines :-) 203-- intersection of two lines :-)
199-- > raise $ dualMV $ raise $ dualMV (x1/\x2) /\ dualV [x3,x4] 204-- > raise $ dualMV $ raise $ dualMV (x1/\x2) /\ dualV [x3,x4]
@@ -211,7 +216,7 @@ y4 = vector "s" [12,0,0,3]
211 216
212asBase r n = filter (\x-> (x==nub x && x==sort x)) $ sequence $ replicate r [1..n] 217asBase r n = filter (\x-> (x==nub x && x==sort x)) $ sequence $ replicate r [1..n]
213 218
214partF t i = part t (name,i) where name = snd . snd . head . dims $ t 219partF t i = part t (name,i) where name = idxName . head . dims $ t
215 220
216--partL = foldl' partF 221--partL = foldl' partF
217 222
@@ -219,7 +224,7 @@ niceAS t = filter ((/=0.0).fst) $ zip vals base
219 where vals = map ((`at` 0).ten.foldl' partF t) (map (map pred) base) 224 where vals = map ((`at` 0).ten.foldl' partF t) (map (map pred) base)
220 base = asBase r n 225 base = asBase r n
221 r = length (dims t) 226 r = length (dims t)
222 n = fst . head . dims $ t 227 n = idxDim . head . dims $ t
223 228
224z1 = vector "p" [0,0,0,1] 229z1 = vector "p" [0,0,0,1]
225z2 = vector "q" [1,0,0,1] 230z2 = vector "q" [1,0,0,1]