Home My Page Projects Code Snippets Project Openings diderot
Summary Activity Tracker Tasks SCM

SCM Repository

[diderot] Annotation of /branches/vis15/src/compiler/ein/mk-operators.sml
ViewVC logotype

Annotation of /branches/vis15/src/compiler/ein/mk-operators.sml

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3496 - (view) (download)

1 : jhr 3477 (* mk-operators.sml
2 : jhr 3476 *
3 : jhr 3477 * Functions to create the various Ein operators.
4 :     *
5 : jhr 3476 * This code is part of the Diderot Project (http://diderot-language.cs.uchicago.edu)
6 :     *
7 :     * COPYRIGHT (c) 2015 The University of Chicago
8 :     * All rights reserved.
9 :     *)
10 :    
11 :     structure MkOperators : sig
12 :    
13 : jhr 3477 type dim = int
14 :     type shape = dim list
15 :     type ids = Ein.index_id list
16 : jhr 3476
17 : jhr 3477 val addRR : Ein.ein
18 :     val addTT : shape -> Ein.ein
19 :     val addTF : dim * shape -> Ein.ein
20 :     val addFF : dim * shape -> Ein.ein
21 :    
22 :     val subRR : Ein.ein
23 :     val subTT : shape -> Ein.ein
24 :     val subTF : dim * shape -> Ein.ein
25 : jhr 3478 val subFT : dim * shape -> Ein.ein
26 : jhr 3477 val subFF : dim * shape -> Ein.ein
27 :    
28 :     val mulRT : shape -> Ein.ein
29 :     val mulRR : Ein.ein
30 :     val mulRF : dim * shape -> Ein.ein
31 :     val mulST : dim * shape -> Ein.ein
32 :     val mulSS : dim -> Ein.ein
33 :     val mulSF : dim * shape -> Ein.ein
34 :    
35 :     val divTR : shape -> Ein.ein
36 :     val divRR : Ein.ein
37 :     val divFR : dim * shape -> Ein.ein
38 :     val divSS : dim -> Ein.ein
39 :     val divFS : dim * shape -> Ein.ein
40 :    
41 :     val negTT : shape -> Ein.ein
42 :     val negFF : dim * shape -> Ein.ein
43 :    
44 :     val cross2TT : Ein.ein
45 :     val cross3TT : Ein.ein
46 :     val cross2FF : Ein.ein
47 :     val cross3FF : Ein.ein
48 :    
49 : cchiw 3495 val outerTT : shape * shape -> Ein.ein
50 :     val outerFF : dim * shape * shape -> Ein.ein
51 : cchiw 3496 val outerTF : dim * shape * shape -> Ein.ein
52 :     val outerFT : dim * shape * shape -> Ein.ein
53 :    
54 : jhr 3477
55 :     val innerTT : shape * ids -> Ein.ein
56 :     val innerFF : shape * dim * ids -> Ein.ein
57 :     val innerFT : shape * dim * ids -> Ein.ein
58 :     val innerTF : shape * dim * ids -> Ein.ein
59 :    
60 :     val colonTT : shape * ids -> Ein.ein
61 :     val colonFF : dim * shape * ids -> Ein.ein
62 :     val colonFT : dim * shape * ids -> Ein.ein
63 :     val colonTF : dim * shape * ids -> Ein.ein
64 :    
65 : cchiw 3494 val normT : shape -> Ein.ein
66 :     val normF : dim * shape -> Ein.ein
67 : jhr 3477
68 :     val normalizeTT : shape -> Ein.ein
69 :     val normalizeFF : dim * shape -> Ein.ein
70 :    
71 :     val traceT : dim -> Ein.ein
72 :     val traceF : dim * shape -> Ein.ein
73 :    
74 :     val transposeT : shape -> Ein.ein
75 :     (* QUESTION: should these be index_kind? *)
76 :     val transposeF : dim * Ein.index_id * Ein.index_id -> Ein.ein
77 :    
78 : jhr 3478 val det2T : Ein.ein
79 :     val det3T : Ein.ein
80 :     val det2F : Ein.ein
81 :     val det3F_full : Ein.ein
82 :     val det3FA : Ein.ein
83 :     val det3F : Ein.ein
84 : jhr 3477
85 :     val expF : dim -> Ein.ein
86 :     val expT : Ein.ein
87 :    
88 : jhr 3478 val powF : dim * int -> Ein.ein
89 :    
90 : jhr 3477 val sqrtF : dim -> Ein.ein
91 :     val cosF : dim -> Ein.ein
92 :     val acosF : dim -> Ein.ein
93 :     val sinF : dim -> Ein.ein
94 :     val asinF : dim -> Ein.ein
95 :     val tanF : dim -> Ein.ein
96 :     val atanF : dim -> Ein.ein
97 : cchiw 3494 val sqrtR : Ein.ein
98 :     val cosR : Ein.ein
99 :     val acosR : Ein.ein
100 :     val sinR : Ein.ein
101 :     val asinR : Ein.ein
102 :     val tanR : Ein.ein
103 :     val atanR : Ein.ein
104 : jhr 3477
105 :     val modulate : dim -> Ein.ein
106 :    
107 :     val identity : dim -> Ein.ein
108 : jhr 3478 val zeros : shape -> Ein.ein
109 : jhr 3477
110 :     val slice : shape * bool list * int list * shape -> Ein.ein
111 :    
112 :     val conv : dim * shape -> Ein.ein
113 :     val probe : shape * dim -> Ein.ein
114 :    
115 :     val curl2d : Ein.ein
116 :     val curl3d : Ein.ein
117 :     val grad : shape -> Ein.ein
118 :     val dotimes : dim * shape -> Ein.ein (* ?? *)
119 :     val divergence : dim * shape -> Ein.ein
120 :    
121 : jhr 3476 end = struct
122 :    
123 :     structure E = Ein
124 :    
125 : jhr 3477 type dim = int
126 :     type shape = dim list
127 :     type ids = int list
128 :    
129 : jhr 3476 fun specialize (alpha, inc) = List.mapi (fn (i, _) => E.V(i + inc)) alpha
130 :    
131 :     fun sumIds (n, inc, alpha) = let
132 :     val vs = List.tabulate(n, fn v => E.V(v+inc))
133 :     in
134 :     ListPair.map (fn(v, i) => (v, 0, i-1)) (vs, alpha)
135 :     end
136 :    
137 :     fun sumIds2 (n, i) = List.tabulate(n, fn v => (E.V v, 0, i))
138 :    
139 :     val subst_flag = 0 (*here*) (* ????? *)
140 :    
141 :     (******************************* Addition *****************************************)
142 :    
143 :     (* Adding tensors : < X{\alpha} + Y_{\alpha}>_{\alpha} *)
144 :     fun addTT alpha = let
145 :     val expindex = specialize(alpha, 0)
146 :     in
147 :     E.EIN{
148 :     params = [E.TEN(subst_flag, alpha), E.TEN(subst_flag, alpha)],
149 :     index = alpha,
150 :     body = E.Opn(E.Add, [E.Tensor(0, expindex), E.Tensor(1, expindex)])
151 :     }
152 :     end
153 :    
154 :     val addRR = addTT []
155 :    
156 : jhr 3477 (* Tensor and Fields *)
157 :     fun addTF (dim, shape) =let
158 :     val expindex = specialize(shape, 0)
159 : jhr 3476 in
160 :     E.EIN{
161 : jhr 3477 params = [E.TEN(subst_flag, shape), E.FLD dim],
162 : jhr 3476 index = shape,
163 : jhr 3477 body = E.Opn(E.Add, [E.Lift(E.Tensor(0, expindex)), E.Field(1, expindex)])
164 : jhr 3476 }
165 :     end
166 :    
167 : jhr 3477 (* Adding Fields : < F{\alpha} + G_{\alpha}>_{\alpha} *)
168 :     fun addFF (dim, shape) =let
169 : jhr 3476 val expindex = specialize(shape, 0)
170 :     in
171 :     E.EIN{
172 : jhr 3477 params = [E.FLD dim, E.FLD dim],
173 : jhr 3476 index = shape,
174 : jhr 3477 body = E.Opn(E.Add, [E.Field(0, expindex), E.Field(1, expindex)])
175 : jhr 3476 }
176 :     end
177 :    
178 :     (********************************* Subtraction **************************************)
179 :    
180 :     fun subTT alpha = let
181 :     val expindex = specialize(alpha, 0)
182 :     in
183 :     E.EIN{
184 : jhr 3477 params = [E.TEN(subst_flag, alpha), E.TEN(subst_flag, alpha)],
185 : jhr 3476 index = alpha,
186 : jhr 3477 body = E.Op2(E.Sub, E.Tensor(0, expindex), E.Tensor(1, expindex))
187 : jhr 3476 }
188 :     end
189 :    
190 :     val subRR = subTT []
191 :    
192 :     fun subTF (dim, shape) = let
193 : jhr 3477 val expindex = specialize(shape, 0)
194 : jhr 3476 in
195 :     E.EIN{
196 : jhr 3477 params = [E.TEN(subst_flag, shape), E.FLD dim],
197 : jhr 3476 index = shape,
198 : jhr 3477 body = E.Opn(E.Add,
199 :     [E.Lift(E.Tensor(0, expindex)), E.Op1(E.Neg, E.Field(1, expindex))])
200 : jhr 3476 }
201 :     end
202 :    
203 :     fun subFT (dim, shape) = let
204 : jhr 3477 val expindex = specialize(shape, 0)
205 : jhr 3476 in
206 :     E.EIN{
207 : jhr 3477 params = [E.TEN(subst_flag, shape), E.FLD dim],
208 : jhr 3476 index = shape,
209 :     body = E.Op2(E.Sub, E.Field(1, expindex), E.Lift(E.Tensor(0, expindex)))
210 :     }
211 :     end
212 :    
213 :     fun subFF (dim, shape) = let
214 : jhr 3477 val expindex = specialize(shape, 0)
215 : jhr 3476 in
216 :     E.EIN{
217 : jhr 3477 params = [E.FLD dim, E.FLD dim],
218 : jhr 3476 index = shape,
219 : jhr 3477 body = E.Op2(E.Sub, E.Field(0, expindex), E.Field(1, expindex))
220 : jhr 3476 }
221 :     end
222 :    
223 :     (********************************** Multiplication *************************************)
224 :    
225 :     (* scalar times tensor product: <s * T_{\alpha}>_{\alpha} *)
226 :     fun mulRT alpha = let
227 : jhr 3477 val expindex = specialize(alpha, 0)
228 : jhr 3476 in
229 :     E.EIN{
230 : jhr 3477 params = [E.TEN(subst_flag, []), E.TEN(subst_flag,alpha)],
231 : jhr 3476 index = alpha,
232 : jhr 3477 body = E.Opn(E.Prod, [E.Tensor(0, []), E.Tensor(1, expindex)])
233 : jhr 3476 }
234 :     end
235 :    
236 :     val mulRR = mulRT []
237 :    
238 : jhr 3477 fun mulRF (dim, shape) =let
239 : jhr 3476 val expindex = specialize(shape, 0)
240 :     in
241 :     E.EIN{
242 : jhr 3477 params = [E.TEN(subst_flag, []), E.FLD dim],
243 : jhr 3476 index = shape,
244 : jhr 3477 body = E.Opn(E.Prod, [E.Lift( E.Tensor(0, [])), E.Field(1,expindex)])
245 : jhr 3476 }
246 :     end
247 :    
248 :     fun mulST (dim, shape) =let
249 :     val expindex = specialize(shape, 0)
250 :     in
251 :     E.EIN{
252 : jhr 3477 params = [E.TEN(subst_flag,shape), E.FLD dim],
253 : jhr 3476 index = shape,
254 : jhr 3477 body = E.Opn(E.Prod, [E.Lift( E.Tensor(0,expindex)), E.Field(1, [])])
255 : jhr 3476 }
256 :     end
257 :    
258 :     fun mulSS dim = E.EIN{
259 : jhr 3477 params = [E.FLD dim, E.FLD dim],
260 : jhr 3476 index = [],
261 : jhr 3477 body = E.Opn(E.Prod, [E.Field(0, []), E.Field(1, [])])
262 : jhr 3476 }
263 :    
264 : jhr 3477 fun mulSF(dim, shape) =let
265 :     val expindex= specialize(shape, 0)
266 : jhr 3476 in E.EIN{
267 : jhr 3477 params = [E.FLD dim, E.FLD dim],
268 : jhr 3476 index = shape,
269 : jhr 3477 body = E.Opn(E.Prod, [E.Field(0, []), E.Field(1, expindex)])
270 : jhr 3476 }
271 :     end
272 :    
273 :     (************************************ Division ************************************)
274 :    
275 : jhr 3477 fun divTR alpha = let
276 :     val expindex = specialize(alpha, 0)
277 :     in
278 : jhr 3476 E.EIN{
279 : jhr 3477 params = [E.TEN(subst_flag, alpha), E.TEN(subst_flag, [])],
280 : jhr 3476 index = alpha,
281 : jhr 3477 body = E.Op2(E.Div, E.Tensor(0, expindex), E.Tensor(1, []))
282 :     }
283 :     end
284 : jhr 3476
285 : jhr 3477 val divRR = divTR []
286 : jhr 3476
287 : jhr 3477 fun divFR (dim, shape) = let
288 :     val expindex = specialize(shape, 0)
289 :     in
290 :     E.EIN{
291 :     params = [E.FLD dim, E.TEN(subst_flag, [])],
292 :     index = shape,
293 :     body = E.Op2(E.Div, E.Field(0, expindex), E.Lift(E.Tensor(1, [])))
294 :     }
295 :     end
296 :    
297 : jhr 3476 fun divSS dim = E.EIN{
298 : jhr 3477 params = [E.FLD dim, E.FLD dim],
299 :     index = [],
300 :     body = E.Op2(E.Div, E.Field(0, []), E.Field(1, []))
301 :     }
302 : jhr 3476
303 : jhr 3477 fun divFS(dim, shape) = let
304 :     val expindex = specialize(shape, 0)
305 :     in
306 :     E.EIN{
307 :     params = [E.FLD dim, E.FLD dim],
308 :     index = shape,
309 :     body = E.Opn(E.Prod, [E.Field(0, expindex), E.Op2(E.Div, E.B(E.Const 1), E.Field(1, []))])
310 :     }
311 :     end
312 : jhr 3476
313 : jhr 3477 (************************************* Negation **********************************)
314 :    
315 : jhr 3476 fun negTT alpha = let
316 : jhr 3477 val expindex= specialize(alpha, 0)
317 : jhr 3476 (*changed tensor lift variable here *)
318 : jhr 3477 in
319 :     E.EIN {
320 :     params = [E.TEN(subst_flag, alpha)], index = alpha,
321 :     body = E.Op1(E.Neg, E.Tensor(0, expindex))
322 :     }
323 :     end
324 : jhr 3476
325 : jhr 3477 fun negFF (dim, shape) = let
326 :     val expindex = specialize(shape, 0)
327 :     in
328 :     E.EIN{
329 :     params = [E.FLD dim], index = shape,
330 :     body = E.Op1(E.Neg, E.Field(0, expindex))
331 :     }
332 :     end
333 : jhr 3476
334 : jhr 3477 (****************************** cross product ***********************************)
335 : jhr 3476
336 : jhr 3477 (* 2-d cross product Eps_{ij}U_i V_j *)
337 : jhr 3476 val cross2TT = E.EIN{
338 : jhr 3477 params = [E.TEN(subst_flag, [2]), E.TEN(subst_flag, [2])],
339 :     index = [],
340 :     body = E.Sum([(E. V 0, 0, 1), (E.V 1, 0, 1)],
341 :     E.Opn(E.Prod, [E.G(E.Eps2(0, 1)), E.Tensor(0, [E.V 0]), E.Tensor(1, [E.V 1])]))
342 :     }
343 : jhr 3476
344 : jhr 3477 (* crossProduct is on 3D vectors ..vec3 t8=t0 × t1; *)
345 :     val cross3TT = E.EIN{
346 :     params = [E.TEN(subst_flag, [3]), E.TEN(subst_flag, [3])],
347 :     index = [3],
348 :     body = E.Sum([(E.V 1, 0, 2), (E.V 2, 0, 2)],
349 :     E.Opn(E.Prod, [
350 :     E.G(E.Epsilon(0, 1, 2)),
351 :     E.Tensor(0, [E.V 1]),
352 :     E.Tensor(1, [E.V 2])
353 :     ]))
354 :     }
355 :    
356 :     (* Field Cross Product *)
357 : jhr 3476 val cross2FF = E.EIN{
358 : jhr 3477 params = [E.FLD(2), E.FLD(2)], index = [],
359 :     body = E.Sum([(E.V 0, 0, 1), (E.V 1, 0, 1)],
360 :     E.Opn(E.Prod, [E.G(E.Eps2(0, 1)), E.Field(0, [E.V 0]), E.Field(1, [E.V 1])]))
361 :     }
362 : jhr 3476
363 : jhr 3477 (* Field Cross Product *)
364 : jhr 3476 val cross3FF = E.EIN{
365 :     params = [E.FLD(3), E.FLD(3)],index= [3],
366 : jhr 3477 body = E.Sum([(E.V 1, 0, 2), (E.V 2, 0, 2)],
367 :     E.Opn(E.Prod, [E.G(E.Epsilon(0, 1, 2)),
368 :     E.Field(0, [E.V 1]), E.Field(1, [E.V 2])]))
369 :     }
370 : jhr 3476
371 : jhr 3477 (******************** outer product ********************************)
372 : cchiw 3495 fun outerTT(alpha,beta) =let
373 :     val expindexA= specialize(alpha,0)
374 :     val expindexB= specialize(beta,length(alpha))
375 :     in
376 :     E.EIN{
377 :     params = [E.TEN(subst_flag,alpha), E.TEN(subst_flag,beta)],
378 :     index= alpha@beta,
379 :     body= E.Opn(E.Prod,[E.Tensor(0, expindexA), E.Tensor(1, expindexB)])
380 :     }
381 :     end
382 : jhr 3476
383 : cchiw 3495 (*Assumes same dimension vector field *)
384 :     fun outerFF(dim,alpha,beta) =let
385 :     val expindexA= specialize(alpha,0)
386 :     val expindexB= specialize(beta,length(alpha))
387 :     in E.EIN{
388 :     params = [E.FLD(dim),E.FLD(dim)],
389 :     index = alpha@beta,
390 :     body = E.Opn(E.Prod,[E.Field(0, expindexA),E.Field(1, expindexB)])
391 :     }
392 :     end
393 : jhr 3477
394 : cchiw 3496 fun outerTF(dim,alpha,beta) =let
395 :     val expindexA= specialize(alpha,0)
396 :     val expindexB= specialize(beta,length(alpha))
397 :     in
398 :     E.EIN{
399 :     params = [E.TEN(subst_flag,alpha),E.FLD(dim)],
400 :     index = alpha@beta,
401 :     body = E.Opn(E.Prod,[E.Lift(E.Tensor(0, expindexA)),E.Field(1, expindexB)])
402 :     }
403 :     end
404 :    
405 :     fun outerFT(dim,alpha,beta) =let
406 :     val expindexA= specialize(alpha,0)
407 :     val expindexB= specialize(beta,length(alpha))
408 :     in
409 :     E.EIN{
410 :     params = [E.FLD(dim),E.TEN(subst_flag,alpha)],
411 :     index = alpha@beta,
412 :     body = E.Opn(E.Prod,[E.Field(0, expindexA),E.Lift(E.Tensor(1, expindexB))])
413 :     }
414 :     end
415 : jhr 3477
416 : cchiw 3496
417 : jhr 3477 (*************************** inner product **********************************)
418 : jhr 3476 (* generic inner product: <T_{\alpha i} * T_{i \beta}>_{\alpha \beta} *)
419 : jhr 3477 fun innerTT (shape1, i::beta) = let
420 :     val alpha = List.take(shape1, length shape1 - 1)
421 :     val expindexA = specialize(alpha, 0)
422 :     val expindexB = specialize(beta, length alpha)
423 :     val s' = E.V(length alpha + length beta)
424 :     val s'' = [(s', 0, i-1)]
425 :     in
426 : jhr 3476 E.EIN{
427 : jhr 3477 params = [E.TEN(subst_flag, shape1), E.TEN(subst_flag, i::beta)],
428 : jhr 3476 index = alpha@beta,
429 : jhr 3477 body = E.Sum(s'', E.Opn(E.Prod, [
430 :     E.Tensor(0, expindexA@[s']), (* T_{\alpha i} *)
431 :     E.Tensor(1, [s']@expindexB ) (* T'_{i \beta} *)
432 :     ]))
433 :     }
434 :     end
435 : jhr 3476 | innerTT _ = raise Fail "Wrong shape for inner product"
436 :    
437 : jhr 3477 (* generic inner product: <T_{\alpha i} * T_{i \beta}>_{\alpha \beta} *)
438 :     fun innerFF (shape1, dim, i::beta) = let
439 :     val alpha = List.take(shape1, length shape1 - 1)
440 :     val expindexA = specialize(alpha, 0)
441 :     val expindexB = specialize(beta, length alpha)
442 :     val sid = E.V(length alpha + length beta)
443 :     in
444 :     E.EIN{
445 :     params = [E.FLD dim, E.FLD dim],
446 :     index = alpha@beta,
447 :     body = E.Sum([(sid, 0, i-1)],
448 :     E.Opn(E.Prod, [
449 :     E.Field(0, expindexA@[sid]), (* F_{\alpha i} *)
450 :     E.Field(1, [sid]@expindexB) (* F'_{i \beta} *)
451 :     ]))
452 :     }
453 :     end
454 :     | innerFF _ = raise Fail "Wrong shape for innerProductField"
455 : jhr 3476
456 : jhr 3477 fun innerFT (shape1, dim, i::beta) = let
457 :     val alpha = List.take(shape1, length shape1-1)
458 :     val expindexA = specialize(alpha, 0)
459 :     val expindexB = specialize(beta, length alpha)
460 :     val sid = E.V(length alpha + length beta)
461 :     in
462 :     E.EIN{
463 :     params = [E.FLD dim, E.TEN(subst_flag, i::beta)],
464 :     index = alpha@beta,
465 :     body = E.Sum([(sid, 0, i-1)],
466 :     E.Opn(E.Prod, [
467 :     E.Field(0, expindexA@[sid]), (* F_{\alpha i} *)
468 :     E.Lift(E.Tensor(1, [sid]@expindexB )) (* F'_{i \beta} *)
469 :     ]))
470 :     }
471 :     end
472 :     | innerFT _ = raise Fail "Wrong shape for innerProductFieldTensor"
473 : jhr 3476
474 : jhr 3477 fun innerTF (shape1, dim, i::beta) = let
475 :     val alpha = List.take(shape1, length shape1 - 1)
476 :     val expindexA = specialize(alpha, 0)
477 :     val expindexB = specialize(beta, length alpha)
478 :     val sid = E.V(length(alpha) + length beta)
479 :     in
480 :     E.EIN{
481 :     params = [E.TEN(subst_flag,shape1), E.FLD dim],index = alpha@beta,
482 :     body = E.Sum([(sid, 0, i-1)],
483 :     E.Opn(E.Prod, [
484 :     E.Lift(E.Tensor(0, expindexA@[sid])), (* F_{\alpha i} *)
485 :     E.Field(1, [sid]@expindexB) (* F'_{i \beta} *)
486 :     ]))
487 :     }
488 :     end
489 :     | innerTF _ = raise Fail "Wrong shape for innerProductTensorField"
490 : jhr 3476
491 : jhr 3477 (*************************** colon product **********************************)
492 : jhr 3476
493 : jhr 3477 (* <T_{\alpha i j} * B{i j \beta }>_\alpha \beta *)
494 :     fun colonTT (shape1, i::j::beta) = let
495 :     val lenAlpha = length shape1 - 2
496 :     val alpha = List.take(shape1, lenAlpha)
497 :     val expindexA = specialize(alpha, 0)
498 :     val expindexB = specialize(beta, lenAlpha)
499 :     val sumi = lenAlpha + length beta
500 :     val s' = [E.V sumi, E.V(sumi+1)]
501 :     val sx = [(E.V sumi, 0, i-1), (E.V(sumi+1), 0, j-1)]
502 :     in
503 : jhr 3476 E.EIN{
504 : jhr 3477 params = [E.TEN(subst_flag, shape1), E.TEN(subst_flag, i::j::beta)],
505 : jhr 3476 index = alpha@beta,
506 : jhr 3477 body = E.Sum(sx,
507 :     E.Opn(E.Prod, [E.Tensor(0, expindexA@s'), E.Tensor(1, s'@expindexB)]))
508 :     }
509 :     end
510 : jhr 3476
511 : jhr 3477 (* <F_{\alpha i j} * G_{i j \beta }>_\alpha \beta *)
512 :     fun colonFF (dim, shape1, i::j::beta) = let
513 :     val lenAlpha = length shape1 - 2
514 :     val alpha = List.take(shape1, lenAlpha)
515 :     val expindexA = specialize(alpha, 0)
516 :     val expindexB = specialize(beta, lenAlpha)
517 :     val sumi = lenAlpha + length beta
518 :     val s' = [E.V sumi, E.V(sumi+1)]
519 :     val sx = [(E.V sumi, 0, i-1), (E.V(sumi+1), 0, j-1)]
520 :     in
521 : jhr 3476 E.EIN{
522 : jhr 3477 params = [E.FLD dim, E.FLD dim],
523 : jhr 3476 index = alpha@beta,
524 : jhr 3477 body = E.Sum(sx,
525 :     E.Opn(E.Prod, [E.Field(0, expindexA@s'), E.Field(1, s'@expindexB)]))
526 :     }
527 :     end
528 : jhr 3476
529 :    
530 : jhr 3477 (* <F_{\alpha i j} * T_{i j \beta }>_\alpha \beta *)
531 :     fun colonFT (dim, shape1, i::j::beta) = let
532 :     val lenAlpha = length shape1 - 2
533 :     val alpha = List.take(shape1, lenAlpha)
534 :     val expindexA = specialize(alpha, 0)
535 :     val expindexB = specialize(beta, lenAlpha)
536 :     val sumi = lenAlpha + length beta
537 :     val s' = [E.V sumi, E.V(sumi+1)]
538 :     val sx = [(E.V sumi, 0, i-1), (E.V(sumi+1), 0, j-1)]
539 :     in
540 : jhr 3476 E.EIN{
541 : jhr 3477 params = [E.FLD dim, E.TEN(subst_flag,shape1)],
542 :     index = alpha@beta,
543 :     body = E.Sum(sx,
544 :     E.Opn(E.Prod, [E.Field(0, expindexA@s'), E.Lift(E.Tensor(1, s'@expindexB))]))
545 :     }
546 :     end
547 : jhr 3476
548 : jhr 3477 (* <T_{\alpha i j} * G{i j \beta }>_\alpha \beta *)
549 :     fun colonTF (dim, shape1, i::j::beta) = let
550 :     val lenAlpha = length shape1 - 2
551 :     val alpha = List.take(shape1, lenAlpha)
552 :     val expindexA = specialize(alpha, 0)
553 :     val expindexB = specialize(beta, lenAlpha)
554 :     val sumi = lenAlpha + length beta
555 :     val s' = [E.V sumi, E.V(sumi+1)]
556 :     val sx = [(E.V sumi, 0, i-1), (E.V(sumi+1), 0, j-1)]
557 :     in
558 : jhr 3476 E.EIN{
559 : jhr 3477 params = [E.TEN(subst_flag, i::j::beta), E.FLD dim],
560 :     index = alpha@beta,
561 :     body = E.Sum(sx,
562 :     E.Opn(E.Prod, [E.Lift(E.Tensor(0, expindexA@s')), E.Field(1, s'@expindexB)]))
563 :     }
564 :     end
565 : jhr 3476
566 : jhr 3477 (******************** Norm ********************************)
567 : jhr 3476
568 :     (*n.t.s tensor norms/normalize use high-il ops
569 :     *we can norm any size and use einapp
570 :     * normalize only implemented for vectors and use mid il-op
571 :     *)
572 : jhr 3477 (* get norm, but without the sqrt. Implemented as a summation over a modulate *)
573 : cchiw 3495 fun normT alpha = let
574 : cchiw 3494 val expindex= specialize(alpha,0)
575 :     val sx= sumIds(length(alpha),0,alpha)
576 :     in
577 : jhr 3476 E.EIN{
578 : cchiw 3494 params = [E.TEN(subst_flag,alpha)],
579 :     index = [],
580 :     body = E.Op1(E.Sqrt,E.Sum(sx, E.Opn(E.Prod,[E.Tensor(0, expindex), E.Tensor(0, expindex)])))
581 :     }
582 :     end
583 :    
584 : cchiw 3495 fun normF(dim,alpha) =let
585 : cchiw 3494 val expindex= specialize(alpha,0)
586 :     val sx= sumIds(length(alpha),0,alpha)
587 :     in
588 :     E.EIN{
589 :     params = [E.FLD(dim)],
590 :     index = [],
591 :     body = E.Op1(E.Sqrt,E.Sum(sx, E.Opn(E.Prod,[E.Field(0, expindex),E.Field(0,expindex)])))
592 :     }
593 :     end
594 : jhr 3476
595 : cchiw 3494 fun normalizeTT alpha = let
596 : jhr 3477 val expindex = specialize(alpha, 0)
597 :     val len = length alpha
598 :     val expindexDot = specialize(alpha, len)
599 :     val sx = sumIds(len, len, alpha)
600 :     val f = E.Tensor(0, expindex)
601 :     val g = E.Tensor(1, expindexDot)
602 :     in
603 : jhr 3476 E.EIN{
604 : jhr 3477 params = [E.TEN(subst_flag, alpha), E.TEN(subst_flag,alpha)],
605 :     index = alpha,
606 :     body = E.Opn(E.Prod, [
607 :     f, E.Op2(E.Div, E.B(E.Const 1), E.Op1(E.Sqrt, E.Sum(sx, E.Opn(E.Prod, [g, g]))))
608 :     ])
609 :     }
610 :     end
611 : jhr 3476
612 : cchiw 3494 fun normalizeFF (dim, alpha as i::_) = let
613 : jhr 3477 val expindex = specialize(alpha, 0)
614 :     val len = length alpha
615 :     val expindexDot = specialize(alpha, len)
616 :     val sx = sumIds(len, len, alpha)
617 :     val f = E.Field(0, expindex)
618 :     val g = E.Field(1, expindexDot)
619 :     in
620 : jhr 3476 E.EIN{
621 : jhr 3477 params = [E.FLD dim, E.FLD dim],
622 : jhr 3476 index = [dim],
623 : jhr 3477 body = E.Opn(E.Prod, [
624 :     f, E.Op2(E.Div, E.B(E.Const 1), E.Op1(E.Sqrt, E.Sum(sx, E.Opn(E.Prod, [g, g]))))
625 :     ])
626 :     }
627 :     end
628 : jhr 3476
629 : jhr 3477 (************************* trace *************************)
630 :    
631 :     (* Trace: <M_{i, i}> This one Sx represents both i's*)
632 : jhr 3476 fun traceT dim = E.EIN{
633 : jhr 3477 params = [E.TEN(1, [dim, dim])], index = [],
634 :     body = E.Sum([(E.V 0, 0, dim-1)], E.Tensor(0, [E.V 0, E.V 0]))
635 :     }
636 : jhr 3476
637 : jhr 3477 (* Trace: <Sigma_i F_{\alpha i, i}> This one Sx represents both i's *)
638 :     fun traceF (dim, alpha) = let
639 :     val expindex = specialize(alpha, 0)
640 :     val s = E.V(length alpha)
641 :     in
642 : jhr 3476 E.EIN{
643 : jhr 3477 params = [E.FLD dim],
644 : jhr 3476 index = alpha,
645 : jhr 3477 body = E.Sum([(s, 0, dim-1)], E.Field(0, expindex@[s, s]))
646 :     }
647 :     end
648 : jhr 3476
649 : jhr 3477 (************************* tranpose *************************)
650 : jhr 3476
651 : jhr 3477 fun transposeT alpha = E.EIN{
652 :     params = [E.TEN(subst_flag, alpha)],
653 :     index = List.rev alpha,
654 :     body = E.Tensor(0, [E.V 1, E.V 0])
655 :     }
656 : jhr 3476
657 : jhr 3477 (* Transpose Field F_{ji} *)
658 :     fun transposeF (dim, i, j) = E.EIN{
659 :     params = [E.FLD dim],
660 :     index = [i, j],
661 :     body = E.Field(0, [E.V 1, E.V 0])
662 :     }
663 : jhr 3476
664 : jhr 3477 (************************* determinant *************************)
665 : jhr 3476
666 : jhr 3478 val det2T = E.EIN{
667 : jhr 3477 params = [E.TEN(0, [2, 2])],
668 :     index = [],
669 :     body = E.Op2(E.Sub,
670 :     E.Opn(E.Prod, [E.Tensor(0, [E.C 0, E.C 0]), E.Tensor(0, [E.C 1, E.C 1])]),
671 :     E.Opn(E.Prod, [E.Tensor(0, [E.C 0, E.C 1]), E.Tensor(0, [E.C 1, E.C 0])]))
672 :     }
673 : jhr 3476
674 : jhr 3478 val det3T = let
675 : jhr 3477 val a = E.Tensor(0, [E.C 0, E.C 0])
676 :     val b = E.Tensor(0, [E.C 0, E.C 1])
677 :     val c = E.Tensor(0, [E.C 0, E.C 2])
678 :     val d = E.Tensor(0, [E.C 1, E.C 0])
679 :     val e = E.Tensor(0, [E.C 1, E.C 1])
680 :     val f = E.Tensor(0, [E.C 1, E.C 2])
681 :     val g = E.Tensor(0, [E.C 2, E.C 0])
682 :     val h = E.Tensor(0, [E.C 2, E.C 1])
683 :     val i = E.Tensor(0, [E.C 2, E.C 2])
684 :     in
685 :     E.EIN{
686 :     params = [E.TEN(0, [3, 3])],
687 :     index = [],
688 :     body = E.Op2(E.Sub,
689 :     E.Opn(E.Add, [
690 :     E.Opn(E.Prod, [a, e, i]), E.Opn(E.Prod, [b, f, g]), E.Opn(E.Prod, [c, d, h])
691 :     ]),
692 :     E.Opn(E.Add, [
693 :     E.Opn(E.Prod, [c, e, g]), E.Opn(E.Prod, [b, d, i]), E.Opn(E.Prod, [a, f, h])
694 :     ]))
695 :     }
696 :     end
697 :    
698 : jhr 3478 val det2F = E.EIN{
699 : jhr 3477 params = [E.FLD 2],
700 :     index = [],
701 :     body = E.Op2(E.Sub,
702 :     E.Opn(E.Prod, [E.Field(0, [E.C 0, E.C 0]), E.Field(0, [E.C 1, E.C 1])]),
703 :     E.Opn(E.Prod, [E.Field(0, [E.C 0, E.C 1]), E.Field(0, [E.C 1, E.C 0])]))
704 :     }
705 :    
706 : jhr 3478 val det3F_full = let
707 : jhr 3477 val a = E.Field(0, [E.C 0, E.C 0])
708 :     val b = E.Field(0, [E.C 0, E.C 1])
709 :     val c = E.Field(0, [E.C 0, E.C 2])
710 :     val d = E.Field(0, [E.C 1, E.C 0])
711 :     val e = E.Field(0, [E.C 1, E.C 1])
712 :     val f = E.Field(0, [E.C 1, E.C 2])
713 :     val g = E.Field(0, [E.C 2, E.C 0])
714 :     val h = E.Field(0, [E.C 2, E.C 1])
715 :     val i = E.Field(0, [E.C 2, E.C 2])
716 :     in
717 :     E.EIN{
718 :     params = [E.FLD 3],
719 :     index = [],
720 :     body = E.Op2(E.Sub,
721 :     E.Opn(E.Add, [
722 :     E.Opn(E.Prod, [a, e, i]), E.Opn(E.Prod, [b, f, g]), E.Opn(E.Prod, [c, d, h])
723 :     ]),
724 :     E.Opn(E.Add, [
725 :     E.Opn(E.Prod, [c, e, g]), E.Opn(E.Prod, [b, d, i]), E.Opn(E.Prod, [a, f, h])
726 :     ]))
727 :     }
728 :     end
729 : jhr 3476
730 : jhr 3478 val det3FA = E.EIN{
731 : jhr 3477 params = [E.FLD 3],
732 :     index = [],
733 :     body = E.Sum([(E.V 0, 0, 2), (E.V 1, 0, 2), (E.V 2, 0, 2)],
734 :     E.Opn(E.Prod, [
735 :     E.G(E.Epsilon(0, 1, 2)),
736 :     E.Field(0, [E.C 0, E. V 0]),
737 :     E.Field(0, [E.C 1, E. V 1]),
738 :     E.Field(0, [E.C 2, E. V 2])
739 :     ]))
740 :     }
741 : jhr 3476
742 : jhr 3478 val det3F = E.EIN{
743 : jhr 3477 params = [E.FLD 3],
744 :     index = [],
745 :     body = E.Sum([(E.V 0, 0, 2)],
746 :     E.Opn(E.Prod, [
747 :     E.Field(0, [E.C 0, E. V 0]),
748 :     E.Sum([(E.V 1, 0, 2)],
749 :     E.Opn(E.Prod, [
750 :     E.Field(0, [E.C 1, E. V 1]),
751 :     E.Sum([(E.V 2, 0, 2)],
752 :     E.Opn(E.Prod, [E.G(E.Epsilon(0, 1, 2)), E.Field(0, [E.C 2, E. V 2])]))
753 :     ]))
754 :     ]))
755 :     }
756 : jhr 3476
757 : jhr 3477 (************************* Exponential **************************)
758 :     fun expF dim = E.EIN{params = [E.FLD dim], index = [], body = E.Op1(E.Exp, E.Field(0, []))}
759 :     val expT = E.EIN{params = [E.TEN(0, [])], index = [], body = E.Op1(E.Exp, E.Tensor(0, []))}
760 : jhr 3476
761 : jhr 3477 (************************* Lifted single-argument math functions *************************)
762 :     local
763 : cchiw 3494 fun liftFn op1 dim = E.EIN{
764 : jhr 3477 params = [E.FLD dim],
765 :     index = [],
766 : cchiw 3494 body = E.Op1(op1, E.Field(0, []))
767 : jhr 3477 }
768 : cchiw 3494 fun op1Fn op0= E.EIN{
769 :     params = [E.TEN(0,[])],
770 :     index = [],
771 :     body = E.Op1(op0, E.Field(0, []))
772 :     }
773 : jhr 3477 in
774 : jhr 3478 fun powF (dim, n) = E.EIN{
775 :     params = [E.FLD dim],
776 :     index = [], body = E.Op1(E.PowInt n, E.Field(0, []))
777 :     }
778 : jhr 3477 val sqrtF = liftFn E.Sqrt
779 :     val cosF = liftFn E.Cosine
780 :     val acosF = liftFn E.ArcCosine
781 :     val sinF = liftFn E.Sine
782 :     val asinF = liftFn E.ArcSine
783 :     val tanF = liftFn E.Tangent
784 :     val atanF = liftFn E.ArcTangent
785 : cchiw 3494 val sqrtR =op1Fn E.Sqrt
786 :     val acosR = op1Fn E.ArcCosine
787 :     val sinR = op1Fn E.Sine
788 :     val asinR = op1Fn E.ArcSine
789 :     val tanR = op1Fn E.Tangent
790 :     val atanR = op1Fn E.ArcTangent
791 : jhr 3477 end (* local *)
792 : jhr 3476
793 : jhr 3477 (************************* other tensor ops *************************)
794 :     fun modulate dim = E.EIN{
795 :     params = [E.TEN(subst_flag, [dim]), E.TEN(subst_flag, [dim])],
796 :     index = [dim],
797 :     body = E.Opn(E.Prod, [E.Tensor(0, [E.V 0]), E.Tensor(1, [E.V 0])])
798 :     }
799 : jhr 3476
800 : jhr 3477 fun identity dim = E.EIN{
801 :     params = [], index = [dim, dim], body = E.G(E.Delta(E.V(0), E.V(1)))
802 :     }
803 : jhr 3476
804 : jhr 3478 fun zeros shape = E.EIN{
805 :     params = [], index = shape, body = E.B(E.Const 0)
806 :     }
807 :    
808 : jhr 3477 (* QUESTION: why do we need the const list? The indices are implicit in the position of
809 :     * of the mask element! Likewise, the result type can be determined from the argTy and
810 :     * mask.
811 :     *)
812 :     fun slice (argTy, mask, const, rstTy) = let
813 :     fun iter ([], _, cnt) = []
814 :     | iter (true::es, c::cs, cnt) = E.C c :: iter(es, cs, cnt)
815 :     | iter (false::es, cs, cnt) = E.V cnt :: iter(es, cs, cnt+1)
816 :     val ix = iter(mask, const, 0)
817 :     in
818 :     E.EIN{
819 :     params = [E.TEN(subst_flag, argTy)],
820 :     index = rstTy,
821 :     body = E.Tensor(0, ix)
822 :     }
823 :     end
824 : jhr 3476
825 : jhr 3477 (******************** other field ops ********************************)
826 : jhr 3476
827 : jhr 3477 (* FLD here is bounded to image field, and dimension of h *)
828 :     fun conv (dim, shape) =let
829 :     val expindex = specialize(shape, 0)
830 :     in
831 :     E.EIN{
832 :     params = [E.IMG(dim,shape), E.KRN],
833 :     index = shape,
834 :     body = E.Conv(0,expindex, 1, [])
835 :     }
836 :     end
837 : jhr 3476
838 : jhr 3477 (* Probe: <F(x)>_{\alpha} *)
839 :     fun probe (alpha, dim) = let
840 :     val expindex = specialize(alpha, 0)
841 :     in
842 :     E.EIN{
843 :     params = [E.FLD dim, E.TEN(0, [])], index = alpha,
844 :     body = E.Probe(E.Field(0, expindex), E.Tensor(1, []))
845 :     }
846 :     end
847 : jhr 3476
848 : jhr 3477 (***************************** derivative ****************************)
849 : jhr 3476
850 : jhr 3477 (* \EinExp{\sum_{ij}\mathcal{E}_{ij} \frac{F_j}{\partial x_i} *)
851 :     val curl2d = E.EIN{
852 :     params = [E.FLD 2],
853 :     index = [],
854 :     body = E.Sum([(E.V 0, 0, 1), (E.V 1, 0, 1)],
855 :     E.Opn(E.Prod, [
856 :     E.G(E.Eps2(0, 1)),
857 :     E.Apply(E.Partial[E.V 0], E.Field(0, [E.V 1]))
858 :     ]))
859 :     }
860 : jhr 3476
861 : jhr 3477 val curl3d = E.EIN{
862 :     params = [E.TEN(1, [3])],
863 :     index = [3],
864 :     body = E.Sum([(E.V 1, 0, 2), (E.V 2, 0, 2)],
865 :     E.Opn(E.Prod, [
866 :     E.G(E.Epsilon(0, 1, 2)),
867 :     E.Apply(E.Partial[E.V 1], E.Field(0, [E.V 2]))
868 :     ]))
869 :     }
870 : jhr 3476
871 : jhr 3477 (*< d F / d_i>_i *)
872 :     fun grad (alpha as a::_) = let
873 :     val expindex = specialize(alpha, 0)
874 :     in
875 :     E.EIN{
876 :     params = [E.FLD a],
877 :     index = alpha,
878 :     body = E.Apply(E.Partial expindex, E.Field(0, []))
879 :     }
880 :     end
881 :    
882 :     (*< Sigma d F_alpha / d x_i>ALpha i CHANGE HERE *)
883 :     fun dotimes (dim, alpha) = let
884 :     val n = length alpha
885 :     (*
886 :     fun expIndex(n, inc) = List.tabulate(n, (fn(x)=>E.V (x+inc)))
887 :     val i' = expIndex(n, 0)
888 :     *)
889 :     val i' = List.tabulate (n, fn x => E.V x)
890 :     in
891 :     E.EIN{
892 :     params = [E.FLD dim], index = alpha@[dim],
893 :     body = E.Apply(E.Partial[E.V n], E.Field(0, i'))
894 :     }
895 :     end
896 :    
897 :     (* <d F_i /d_i> *)
898 :     fun divergence (dim, alpha) = let
899 :     val expindex = specialize(alpha, 0)
900 :     val sumI = length alpha
901 :     val sumIndex = E.V sumI
902 :     val sumIndexL = [sumIndex]
903 :     val S = expindex@sumIndexL
904 :     in
905 :     E.EIN{
906 :     params = [E.FLD dim],
907 :     index = alpha,
908 :     body = E.Sum([(sumIndex, 0, dim-1)], E.Apply(E.Partial sumIndexL, E.Field(0, S)))
909 :     }
910 :     end
911 :    
912 : jhr 3476 end (* mkOperators *)

root@smlnj-gforge.cs.uchicago.edu
ViewVC Help
Powered by ViewVC 1.0.0