1
1
;;; -*- Mode:LISP; Package:MACSYMA -*-
2
;; (based on itensor.116 ,117)
2
;; ** (c) Copyright 1981 Massachusetts Institute of Technology **
4
;; This program is free software; you can redistribute it and/or
5
;; modify it under the terms of the GNU General Public License as
6
;; published by the Free Software Foundation; either version 2 of
7
;; the License, or (at your option) any later version.
9
;; This program is distributed in the hope that it will be
10
;; useful, but WITHOUT ANY WARRANTY; without even the implied
11
;; warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12
;; PURPOSE. See the GNU General Public License for more details.
16
;; The Itensor package was downcased, cleaned up, and moving frames
17
;; functionality was added by Viktor Toth (http://www.vttoth.com/).
19
;; As of November, 2004, the naming conventions in this package now
20
;; correspond with the naming conventions in commercial MACSYMA.
5
25
(macsyma-module itensor) ;; added 9/24/82 at UCB
6
; ** (c) Copyright 1981 Massachusetts Institute of Technology **
8
; Various functions in Itensr have been parceled out to separate files. A
27
(cond (($get '$itensor '$version) (merror "ITENSOR already loaded"))
28
(t ($put '$itensor '$v20041126 '$version))
31
; Various functions in Itensor have been parceled out to separate files. A
9
32
; function in one of these files will only be loaded in (automatically) if
10
; explicitly used in the Macsyma. (It is necessary to have first loaded in
11
; ITENSR FASL for this autoloading to take place.) The current status of
33
; explicitly used in the Maxima. (It is necessary to have first loaded in
34
; ITENSOR FASL for this autoloading to take place.) The current status of
12
35
; these separate files are:
14
37
; Filename Macsyma Functions
15
38
; -------- -----------------
16
39
; CANTEN FASL CANTEN, CONCAN, IRPMON
17
; GENER FASL GENERATE, MAKEBOX, AVERAGE, CONMETDERIV, FLUSH1DERIV,
40
; GENER FASL IC_CONVERT, MAKEBOX, AVERAGE, CONMETDERIV, FLUSH1DERIV,
19
42
; SYMTRY FASL CANFORM, DECSYM, DISPSYM, REMSYM
22
(putprop '$GENERATE '((dsk tensor) gener fasl) 'autoload)
23
(putprop '$DECSYM '((dsk tensor) symtry fasl) 'autoload)
24
(putprop '$CANFORM '((dsk tensor) symtry fasl) 'autoload)
25
(putprop '$CANTEN '((dsk tensor) canten fasl) 'autoload)
26
(putprop '$MAKEBOX '((dsk tensor) gener fasl) 'autoload)
27
(putprop '$GEODESIC '((dsk tensor) gener fasl) 'autoload)
28
(putprop '$CONMETDERIV '((dsk tensor) gener fasl) 'autoload))
45
(putprop '$ic_convert '((dsk tensor) gener fasl) 'autoload)
46
(putprop '$decsym '((dsk tensor) symtry fasl) 'autoload)
47
(putprop '$canform '((dsk tensor) symtry fasl) 'autoload)
48
(putprop '$canten '((dsk tensor) canten fasl) 'autoload)
49
(putprop '$makebox '((dsk tensor) gener fasl) 'autoload)
50
(putprop '$igeodesic_coords '((dsk tensor) gener fasl) 'autoload)
51
(putprop '$conmetderiv '((dsk tensor) gener fasl) 'autoload))
30
(putprop '$GENERATE (concat vaxima-main-dir '|//tensor//gener|) 'autoload)
31
(putprop '$DECSYM (concat vaxima-main-dir '|//tensor//symtry| )'autoload)
32
(putprop '$CANFORM (concat vaxima-main-dir '|//tensor//symtry| )'autoload)
33
(putprop '$CANTEN (concat vaxima-main-dir '|//tensor//canten| )'autoload)
34
(putprop '$MAKEBOX (concat vaxima-main-dir '|//tensor//gener| )'autoload)
35
(putprop '$GEODESIC (concat vaxima-main-dir '|//tensor//gener| )'autoload)
36
(putprop '$CONMETDERIV (concat vaxima-main-dir '|//tensor//gener| )'autoload))
53
(putprop '$ic_convert (concat vaxima-main-dir '|//tensor//gener|) 'autoload)
54
(putprop '$decsym (concat vaxima-main-dir '|//tensor//symtry| )'autoload)
55
(putprop '$canform (concat vaxima-main-dir '|//tensor//symtry| )'autoload)
56
(putprop '$canten (concat vaxima-main-dir '|//tensor//canten| )'autoload)
57
(putprop '$makebox (concat vaxima-main-dir '|//tensor//gener| )'autoload)
58
(putprop '$igeodesic_coords (concat vaxima-main-dir '|//tensor//gener| )'autoload)
59
(putprop '$conmetderiv (concat vaxima-main-dir '|//tensor//gener| )'autoload))
40
(autof '$GENERATE '|gener|)
41
(autof '$DECSYM '|symtry|)
42
(autof '$CANFORM '|symtry|)
43
(autof '$CANTEN '|canten|)
44
(autof '$MAKEBOX '|gener|)
45
(autof '$GEODESIC '|gener|)
46
(autof '$CONMETDERIV '|gener|)
47
(autof '$NAME '|canten|)
48
(autof '$CONTI '|canten|)
49
(autof '$COVI '|canten|)
50
(autof '$DERI '|canten|)
63
(autof '$ic_convert '|gener|)
64
(autof '$decsym '|symtry|)
65
(autof '$canform '|symtry|)
66
(autof '$canten '|canten|)
67
(autof '$makebox '|gener|)
68
(autof '$igeodesic_coords '|gener|)
69
(autof '$conmetderiv '|gener|)
70
(autof '$name '|canten|)
53
74
(eval-when (eval compile)
54
(defmacro fixp (x) `(typep ,x 'fixnum)))
58
#+maclisp ($UUO) ;Restore calls to SDIFF so it can be redefined
60
(DECLARE-TOP (SPECIAL SMLIST $DUMMYX $COORDINATES $METRIC $COUNTER $DIM
61
$CONTRACTIONS $COORD $ALLSYM $METRICCONVERT)
62
(*LEXPR $RENAME $DIFF $COORD $REMCOORD $LORENTZ))
64
(SETQ $DUMMYX '$% ;Prefix for dummy indices
65
$COUNTER 0. ;Dummy variable numeric indexs
66
SMLIST '(MLIST SIMP) ;Simplified MLIST header
67
$COORDINATES NIL ;Used when differentiating w.r.t. a number
68
$COORD '((MLIST SIMP)) ;Objects treated liked coordinates in DIFF
69
$ALLSYM T ;If T then all indexed objects symmetric
70
$METRICCONVERT T) ;Flag used by $GENERATE
72
;(DEFUN IFNOT MACRO (CLAUSE) (CONS 'OR (CDR CLAUSE)))
73
(DEFmacro IFNOT (&rest CLAUSE) `(or ,@ clause))
75
;(DEFUN M+OR*OR^P MACRO (CL)
76
(defmacro M+OR*OR^P (&whole cl &rest ign) ign
79
'(MEMQ (CAAR X) '(MTIMES MPLUS MEXPT))))
81
(DEFMFUN $DUMMY nil ;Sets arguments to dummy indices
82
(progn (setq $COUNTER (1+ $COUNTER))
83
(concat $DUMMYX $COUNTER)))
85
(DEFPROP $KDELTA ((/ . / )) CONTRACTIONS)
75
(defmacro fixp (x) `(typep ,x 'fixnum))
78
#+maclisp ($UUO) ;Restore calls to SDIFF so it can be redefined
81
(special smlist $idummyx $vect_coords $imetric $icounter $dim
82
$contractions $coord $allsym $metricconvert $iframe_flag
83
$itorsion_flag $inonmet_flag)
84
(*lexpr $rename $diff $idiff $coord $remcoord $lorentz_gauge)
87
(setq $idummyx '$% ;Prefix for dummy indices
88
$icounter 0. ;Dummy variable numeric index
89
smlist '(mlist simp) ;Simplified mlist header
90
$vect_coords nil ;Used when differentiating w.r.t. a number
91
$coord '((mlist simp)) ;Objects treated liked coordinates in diff
92
$allsym nil ;If T then all indexed objects symmetric
93
$metricconvert t ;Flag used by $ic_convert
98
;(defun ifnot macro (clause) (cons 'or (cdr clause)))
99
(defmacro ifnot (&rest clause) `(or ,@ clause))
101
;(defun m+or*or^p macro (cl)
102
(defmacro m+or*or^p (&whole cl &rest ign) ign
105
'(memq (caar x) '(mtimes mplus mexpt))))
107
(defmfun $idummy nil ;Sets arguments to dummy indices
108
(progn (setq $icounter (1+ $icounter))
109
(concat $idummyx $icounter)))
111
(defprop $kdelta ((/ . / )) contractions)
113
(defun isprod (x) (or (equal x '(mtimes)) (equal x '(mtimes simp))
114
(equal x '(mtimes simp ratsimp)))
117
;; Remove occurrences of ratsimp from elements of x
122
((eq (car x) 'ratsimp) (derat (cdr x)))
123
(t (cons (derat (car x)) (derat (cdr x))))
130
((and (numberp (car l)) (< (car l) 0)) (plusi (cdr l)))
131
((atom (car l)) (cons (car l) (plusi (cdr l))))
132
((and (isprod (caar l)) (eq (cadar l) -1)) (plusi (cdr l)))
133
(t (cons (car l) (plusi (cdr l))))
140
((and (numberp (car l)) (< (car l) 0)) (cons (neg (car l)) (plusi (cdr l))))
141
((atom (car l)) (minusi (cdr l)))
143
(and (isprod (caar l)) (eq (cadar l) -1))
144
(cons (caddar l) (minusi (cdr l)))
151
(defun covi (rp) (plusi (cdadr rp)))
152
(defun conti (rp) (append (minusi (cdadr rp)) (cdaddr rp)))
153
(defun deri (rp) (cdddr rp))
154
(defun name (rp) (caar rp))
155
(defmfun $covi (rp) (cond ((rpobj rp) (cons smlist (covi rp)))
156
(t (merror "Not an RPOBJ"))
159
(defmfun $conti (rp) (cond ((rpobj rp) (cons smlist (conti rp)))
160
(t (merror "Not an RPOBJ"))
163
(defmfun $deri (rp) (cond ((rpobj rp) (cons smlist (deri rp)))
164
(t (merror "Not an RPOBJ"))
167
(defmfun $name (rp) (cond ((rpobj rp) (caar rp)) (t (merror "Not an RPOBJ"))))
87
169
;KDELTA has special contraction property because it contracts with any indexed
90
(meval '(($DECLARE) %KDELTA $CONSTANT)) ;So derivative will be zero
172
(meval '(($declare) %kdelta $constant)) ;So derivative will be zero
173
(meval '(($declare) $kdelta $constant)) ;So derivative will be zero
174
(meval '(($declare) %levi_civita $constant))
175
(meval '(($declare) $levi_civita $constant))
92
(SETQ $DIM 4. $CONTRACTIONS '((MLIST SIMP)))
177
(setq $dim 4. $contractions '((mlist simp)))
94
(DEFMFUN $DEFCON N ;Defines contractions: A contracts with B to form C
96
(ADD2LNC A $CONTRACTIONS)
99
(CONS (COND ((= N 1.) '(/ . / ))
100
((= N 3.) (CONS (ARG 2.) (ARG 3.)))
101
(T (merror "DEFCON takes 1 or 3 arguments")))
102
(ZL-GET A 'CONTRACTIONS))
179
(defmfun $defcon n ;Defines contractions: A contracts with B to form C
181
(add2lnc a $contractions)
184
(cons (cond ((= n 1.) '(/ . / ))
185
((= n 3.) (cons (arg 2.) (arg 3.)))
186
(t (merror "DEFCON takes 1 or 3 arguments")))
187
(zl-get a 'contractions))
107
(DEFMSPEC $DISPCON (A) (SETQ A (CDR A))
192
(defmspec $dispcon (a) (setq a (cdr a))
108
193
;;Displays contraction definitions
110
(AND (EQ (CAR A) '$ALL) (SETQ A (CDR $CONTRACTIONS)))
115
(COND ((SETQ TMP (ZL-GET E 'CONTRACTIONS))
117
(MAPCAR #'(LAMBDA (Z)
195
(and (eq (car a) '$all) (setq a (cdr $contractions)))
200
(cond ((setq tmp (zl-get e 'contractions))
202
(mapcar #'(lambda (z)
126
(T '((MLIST SIMP)))))
211
(t '((mlist simp)))))
130
(DEFMSPEC $REMCON (A) (SETQ A (CDR A))
215
(defmspec $remcon (a) (setq a (cdr a))
131
216
;;Removes contraction definitions
132
(AND (EQ (CAR A) '$ALL) (SETQ A (CDR $CONTRACTIONS)))
133
(CONS SMLIST (MAPC (FUNCTION (LAMBDA (E) (ZL-REMPROP E 'CONTRACTIONS)
134
(DELQ E $CONTRACTIONS)))
217
(and (eq (car a) '$all) (setq a (cdr $contractions)))
218
(cons smlist (mapc (function (lambda (e) (zl-remprop e 'contractions)
219
(delq e $contractions)))
138
223
;; Helper to obtain contractions on both the noun and verb form of E
139
(COND ((AND (SYMBOLP E) (EQ (GETCHAR E 1) '%)) (ZL-GET ($VERBIFY E) 'CONTRACTIONS))
140
(T (ZL-GET E 'CONTRACTIONS))
224
(cond ((and (symbolp e) (eq (getchar e 1) '%)) (zl-get ($verbify e) 'contractions))
225
(t (zl-get e 'contractions))
144
(DEFUN RPOBJ (E) ;"True" if an indexed object and not a matrix
145
(COND ((AND (NOT (ATOM E)) (EQ (CAAR E) 'MQAPPLY)) (RPOBJ (CDR E)))
148
(NOT (EQ (CAAR E) '$MATRIX))
150
(COND ((CDDR E) ($LISTP (CADDR E)))
151
(T (NCONC E '(((MLIST SIMP))))))))))
229
(defun rpobj (e) ;"True" if an indexed object and not a matrix
230
(cond ((and (not (atom e)) (eq (caar e) 'mqapply)) (rpobj (cdr e)))
233
(not (eq (caar e) '$matrix))
235
(cond ((cddr e) ($listp (caddr e)))
236
(t (nconc e '(((mlist simp)))) t ))))))
152
237
;Transforms F([...]) into F([...],[])
154
239
;RPOBJ is the predicate for indexed objects. In the case of no contravariant
157
242
(deff $tenpr #'rpobj)
159
(DEFMFUN $METRIC (V) (SETQ $METRIC V) ($DEFCON V) ($DEFCON V V '$KDELTA))
161
(DEFUN MYSUBST0 (NEW OLD) ;To reuse subparts of old expression
162
(COND ((ALIKE1 NEW OLD) OLD) (T NEW)))
164
(DEFUN COV (A B) ;COV gives covariant form of metric
165
(COND ((BOUNDP '$METRIC)
166
(MEVAL (LIST (NCONS $METRIC)
169
(T (merror "Name of metric must be specified"))))
171
(DEFUN CONTR (A B) ;CONTR gives contraviant form of metric
172
(COND ((BOUNDP '$METRIC)
173
(MEVAL (LIST (NCONS $METRIC)
176
(T (merror "Name of metric must be specified"))))
178
(DEFUN DIFFCOV (A B D)
179
(COND ((BOUNDP '$METRIC)
180
(MEVAL (LIST (NCONS $METRIC)
244
(defmfun $imetric (v) (setq $imetric v) ($defcon v) ($defcon v v '$kdelta))
246
(defun mysubst0 (new old) ;To reuse subparts of old expression
247
(cond ((alike1 new old) old) (t new)))
249
(defun cov (a b) ;COV gives covariant form of metric
250
(cond ((boundp '$imetric)
251
(meval (list (ncons $imetric)
254
(t (merror "Name of metric must be specified"))))
256
(defun contr (a b) ;contr gives contraviant form of metric
257
(cond ((boundp '$imetric)
258
(meval (list (ncons $imetric)
261
(t (merror "Name of metric must be specified"))))
263
(defun diffcov (a b d)
264
(cond ((boundp '$imetric)
265
(meval (list (ncons $imetric)
187
(T (merror "Name of metric must be specified"))))
189
;(DEFMFUN $CHR1 (L1) ;Christoffel symbol of first kind
191
; (SETQ A (CADDDR L1) B (CADR L1) C (CADDR L1))
192
; (RETURN (LIST '(MTIMES)
193
; '((RAT SIMP) 1. 2.)
195
; (SDIFF (COV B A) C)
196
; (SDIFF (COV C A) B)
199
; (SDIFF (COV B C) A)))))))
203
; ((> NARGS 2) (RETURN (MEVAL (CONS '$COVDIFF (CONS ($CHR1 (ARG 1) (ARG 2)) (CDDR (LISTIFY NARGS)))))))
204
((> NARGS 2) (RETURN (MEVAL (CONS '$DIFF (CONS ($CHR1 (ARG 1) (ARG 2)) (APPLY #'APPEND (MAPCAR #'(LAMBDA (E) (LIST E 1)) (CDDR (LISTIFY NARGS)))))))))
205
((> NARGS 1) (AND (EQ 1 (LENGTH (ARG 2))) (RETURN ($CHR1 (ARG 1))))
206
(merror "CHR1 cannot have contravariant indices"))
208
;; ((NULL (ARG 2)) (merror "CHR1 has no contravariant indices"))
209
;; (T (RETURN ($CHR1 (ARG 1))))
213
(SETQ A (CADDDR (ARG 1)) B (CADR (ARG 1)) C (CADDR (ARG 1)))
214
(RETURN (LIST '(MTIMES)
217
;; (SDIFF (COV B A) C)
218
;; (SDIFF (COV C A) B)
222
;; (SDIFF (COV B C) A))))))
228
;(DEFMFUN $CHR2 (L1 L2) ;Christoffel symbol of second kind
230
; (SETQ A (CADR L1) B (CADDR L1) C (CADR L2))
231
; (return (do ((flag) (l (append (cdr l1) (cdr l2))))
232
; (flag (LIST '(MTIMES)
234
; ($CHR1 (LIST SMLIST A B D))))
236
; (and (not (memq d l)) (setq flag t))))))
240
((> NARGS 2) (RETURN (MEVAL (CONS '$DIFF (CONS ($CHR2 (ARG 1) (ARG 2)) (APPLY #'APPEND (MAPCAR #'(LAMBDA (E) (LIST E 1)) (CDDR (LISTIFY NARGS)))))))))
242
(SETQ A (CADR (ARG 1)) B (CADDR (ARG 1)) C (CADR (ARG 2)))
243
(return (do ((flag) (l (append (cdr (ARG 1)) (cdr (ARG 2)))))
244
(flag (LIST '(MTIMES)
246
($CHR1 (LIST SMLIST A B D))))
248
(and (not (memq d l)) (setq flag t))))))
252
(DEFMFUN $curvature (L1 L2)
255
(SETQ I (CADR L1) K (CADDR L1) H (CADDDR L1) J (CADR L2))
256
(RETURN (LIST '(MPLUS)
257
(SDIFF (LIST '($CHR2 SIMP)
263
(SDIFF (LIST '($CHR2 SIMP)
283
(DEFUN COVSUBST (X Y RP) ;Substitutes X for Y in the covariant part of RP
284
(CONS (CAR RP) (CONS (SUBST X Y (CADR RP)) (CDDR RP))))
286
(DEFUN CONSUBST (X Y RP) ;Substitutes X for Y in the contravariant part of RP
289
(CONS (SUBST X Y (CADDR RP)) (CDDDR RP)))))
291
(DEFUN DERSUBST (X Y RP) ;Substitutes X for Y in the derivative indices of RP
292
(NCONC (LIST (CAR RP) (CADR RP) (CADDR RP))
293
(SUBST X Y (CDDDR RP))))
295
(DECLARE-TOP (SPECIAL X TEMP D))
297
(DEFMFUN $COVDIFF NARGS
299
(AND (< NARGS 2) (merror "COVDIFF must have at least 2 args"))
301
AGAIN (SETQ X (ARG I) E (COVDIFF E) I (1+ I))
302
(AND (> I NARGS) (RETURN E))
308
((OR (ATOM E) (EQ (CAAR E) 'RAT)) (SDIFF E X))
310
(SETQ TEMP (MAPCAR #'(LAMBDA (V) (LIST '(MTIMES)
322
((OR (CDADR E) (CDDDR E))
328
(NCONC (MAPCAR #'(LAMBDA (V)
338
(MAPCAR #'(LAMBDA (V)
352
((EQ (CAAR E) 'MTIMES) (SIMPLUS (COVDIFFTIMES (CDR E) X) 1 T))
353
((EQ (CAAR E) 'MPLUS)
354
(SIMPLIFYA (CONS '(MPLUS)
355
(MAPCAR 'COVDIFF (CDR E)))
357
((EQ (CAAR E) 'MEXPT)
358
(SIMPTIMES (LIST '(MTIMES)
362
(LIST '(MPLUS) -1. (CADDR E)))
363
($COVDIFF (CADR E) X))
366
((EQ (CAAR E) 'MEQUAL)
367
(LIST (CAR E) (COVDIFF (CADR E)) (COVDIFF (CADDR E))))
369
(MERROR "Not acceptable to COVDIFF: ~M"
372
(DEFUN COVDIFFTIMES (L X)
374
(SETQ OUT (NCONS '(MPLUS)))
375
LOOP (SETQ SP (CAR L) L (CDR L))
377
(LIST (SIMPTIMES (CONS '(MTIMES)
378
(CONS ($COVDIFF SP X)
382
(COND ((NULL L) (RETURN OUT)))
383
(SETQ LEFT (NCONC LEFT (NCONS SP)))
386
(DECLARE-TOP (UNSPECIAL R TEMP D))
389
(cond ((equal n 0) (merror "LORENTZ requires at least one argument"))
272
(t (merror "Name of metric must be specified"))))
274
(defmfun $ichr1 nargs ; Christoffel-symbol of the first kind
278
(> nargs 2) ; Derivative indices present; use idiff() to resolve
284
($ichr1 (arg 1) (arg 2))
287
(mapcar #'(lambda (e) (list e 1)) (cddr (listify nargs)))
296
(and (eq 1 (length (arg 2))) (return ($ichr1 (arg 1))))
297
(merror "ichr1 cannot have contravariant indices")
299
(t ; G_abc = 1/2*(g_ba,c+g_ca,b-g_bc,a)
300
(setq a (cadddr (arg 1)) b (cadr (arg 1)) c (caddr (arg 1)))
309
(list '(mtimes) -1. (diffcov b c a))
318
(defmfun $ichr2 nargs ; Christoffel-symbol of the second kind
322
(> nargs 2) ; Derivative indices present; use idiff() to resolve
328
($ichr2 (arg 1) (arg 2))
331
(mapcar #'(lambda (e) (list e 1)) (cddr (listify nargs)))
338
(t ; G_ab^c=g^cd*G_abd
339
(setq a (cadr (arg 1)) b (caddr (arg 1)) c (cadr (arg 2)))
342
((flag) (l (append (cdr (arg 1)) (cdr (arg 2)))))
344
(list '(mtimes) (contr c d) ($ichr1 (list smlist a b d)))
347
(and (not (memq d l)) (setq flag t))
355
(defmfun $icurvature (l1 l2)
357
(setq r ($idummy) i (cadr l1) k (caddr l1) h (cadddr l1) j (cadr l2))
361
(idiff (list (diffop) (list smlist i k) l2) h)
364
(idiff (list (diffop) (list smlist i h) (list smlist j)) k)
368
(list (diffop) (list smlist i k) (list smlist r))
369
(list (diffop) (list smlist r h) l2)
374
(list (diffop) (list smlist i h) (list smlist r))
375
(list (diffop) (list smlist r k) l2)
382
(list '($ifb) (list smlist k h) (list smlist r))
383
(list '($icc2) (list smlist r i) (list smlist j))
393
(defun covsubst (x y rp) ;Substitutes X for Y in the covariant part of RP
394
(cons (car rp) (cons (subst x y (cadr rp)) (cddr rp))))
396
(defun consubst (x y rp) ;Substitutes X for Y in the contravariant part of RP
399
(cons (subst x y (caddr rp)) (cdddr rp)))))
401
(defun dersubst (x y rp) ;Substitutes X for Y in the derivative indices of RP
402
(nconc (list (car rp) (cadr rp) (caddr rp))
403
(subst x y (cdddr rp))))
405
;; COVARIANT DIFFERENTIATION
406
;; As of November, 2004, COVDIFF now takes into account the value of
407
;; iframe_flag. If true, COVDIFF uses the coefficients icc2 in place
408
;; of the Christoffel-symbols ichr2.
410
(defun diffop () ; ichr2 or icc2 depending on iframe_flag
413
(or $iframe_flag $itorsion_flag $inonmet_flag)
420
(declare-top (special x temp d))
422
(defmfun $covdiff nargs
425
(and (< nargs 2) (merror "COVDIFF must have at least 2 args"))
427
again (setq x (arg i) e (covdiff e) i (1+ i))
428
(and (> i nargs) (return e))
433
(defun covdiff (e) ; The covariant derivative...
436
( ; is the partial derivative for scalars (*** torsion?)
437
(or (atom e) (eq (caar e) 'rat))
446
(list (diffop) (list smlist d x) (list smlist v))
460
(or (cdadr e) (cdddr e))
461
(cons (list '(mtimes) -1. (cons '(mplus)
504
(eq (caar e) 'mtimes) ; (a*b)'
506
(covdifftimes (cdr e) x)
511
(eq (caar e) 'mplus) ; (a+b)'=a'+b'
515
(mapcar 'covdiff (cdr e))
521
(eq (caar e) 'mexpt) ; (a^b)'=b*a^(b-1)*a'
529
(list '(mplus) -1. (caddr e))
531
($covdiff (cadr e) x)
537
(eq (caar e) 'mequal)
538
(list (car e) (covdiff (cadr e)) (covdiff (caddr e)))
540
(t (merror "Not acceptable to COVDIFF: ~M" (ishow e)))
544
(defun covdifftimes (l x)
546
(setq out (ncons '(mplus)))
547
loop (setq sp (car l) l (cdr l))
551
(cons '(mtimes) (cons ($covdiff sp x) (append left l)))
556
(cond ((null l) (return out)))
557
(setq left (nconc left (ncons sp)))
562
(declare-top (unspecial r temp d))
564
(defun vecdiff (v i j d) ;Add frame bracket contribution when iframe_flag:true
571
(list (list v) '((mlist)) (list '(mlist) i) j)
574
(list (list v) '((mlist)) (list '(mlist) d))
578
(list '(%ifb) (list '(mlist) d j) (list '(mlist) i))
585
(list (list v) '((mlist)) (list '(mlist) i) j)
590
(defun liediff (v e n)
592
((not (symbolp v)) (merror "~M is not a symbol" v))
594
(or (atom e) (eq (caar e) 'rat)) ; Scalar field
595
; v([],[%1])*idiff(e,%1)
597
((dummy (implode (nconc (exploden $idummyx) (exploden n)))))
599
'(mtimes) (list (list v) '((mlist)) (list '(mlist) dummy))
605
(rpobj e) ; Tensor field
607
; Dummy implementation for logic tests
608
; (list '(%liediff) v e)
610
; Shall the dummy index be in ICOUNTER sequence? Probably yes.
611
; (let ((dummy (implode (nconc (exploden $idummyx) (exploden n)))))
617
($iframe_flag ($idummy))
627
'(mtimes) ; e([...],[...],%1)*v([],[%1])
628
(list (list v) '((mlist)) (list '(mlist) dummy))
633
#'(lambda (s) ; e([..%1..],[...])*v([],[%1],k)
636
(cond ((atom (car s)) 1) (t -1))
643
(subseq (cdadr e) 0 (- (length (cdadr e)) (length s)))
645
(cond ((atom (car s)) dummy)
646
(t (list '(mtimes simp) -1 dummy))
658
(cond ((atom (car s)) dummy) (t (caddr (car s))))
659
(cond ((atom (car s)) (car s)) (t dummy))
667
#'(lambda (s) ; +e([...],[...],..%1..)*v([],[%1],k)
671
(list (car e) (cadr e) (caddr e))
672
(subseq (cdddr e) 0 (- (length (cdddr e)) (length s)))
675
(vecdiff v dummy (car s) dummy2)
681
#'(lambda (s) ; -e([...],[..%1..])*v([],[k],%1)
685
(list (car e) (cadr e)
689
(subseq (cdaddr e) 0 (- (length (cdaddr e)) (length s)))
696
(vecdiff v (car s) dummy dummy2)
705
(eq (caar e) 'mtimes) ; Leibnitz rule
706
; Lv(cadr e)*(cddr e)+(cadr e)*Lv(cddr e)
709
(cons '(mtimes) (cons (liediff v (cadr e) n) (cddr e)))
716
(cond ((cdddr e) (cons '(mtimes) (cddr e))) (t (caddr e)))
724
(eq (caar e) 'mplus) ; Linearity
725
; We prefer mapcar to iteration, but the recursive code also works
728
; (liediff v (cadr e) n)
729
; (liediff v (cond ((cdddr e) (cons '(mplus) (cddr e))) (t (caddr e))) n)
731
(cons '(mplus) (mapcar #'(lambda (u) (liediff v u n)) (cdr e)))
733
(t (merror "~M is not a tensorial expression liediff can handle" e))
737
(defmfun $liediff (v e) (liediff v e 1))
739
(defmfun $rediff (x) (meval '(($ev) x $idiff)))
741
;;(defmfun $evundiff (x) ($rediff ($undiff x)))
742
(defmfun $evundiff (x) (meval (list '($ev) ($undiff x) '$nouns)))
753
(list '(%idiff) (list (car x) (cadr x) (caddr x)))
754
(putinones (cdddr x))
762
(simplifya (cons (ncons (caar x)) (mapcar '$undiff (cdr x))) t)
771
((cdr e) (cons (car e) (cons 1. (putinones (cdr e)))))
772
(t (list (car e) 1.))
778
(defmfun $lorentz_gauge n
779
(cond ((equal n 0) (merror "LORENTZ_GAUGE requires at least one argument"))
390
780
((equal n 1) (lorentz (arg 1) nil))
391
781
(t (lorentz (arg 1)
392
782
((lambda (l) (cond ((sloop for v in l
393
783
always (symbolp v)) l)
395
"Invalid tensor name(s) in argument to LORENTZ"))))
785
"Invalid tensor name(s) in argument to LORENTZ_GAUGE"))))
396
786
(listify (f- 1 n)))))))
398
788
;Lorentz contraction of E: indexed objects with a derivative index matching a
399
789
;contravariant index become 0. If L is NIL then do this for all indexed objects
400
790
;otherwise do this only for those indexed objects whose names are members of L.
403
793
(cond ((atom e) e)
405
795
(cond ((and (or (null l) (memq (caar e) l))
548
1044
;CONTRACT3 here if CL is known not to be null.
549
1045
;If L3 is empty then there is nothing left to
551
LOOP2(AND (NULL CL) (RETURN (NCONC L1 L3)))
553
(AND (NULL L3) (RETURN (NCONC L1 CL)))
554
(SETQ F (CAR CL) CL (CDR CL))
555
(COND ((SETQ SF (CONTRACT3 F L3)) (SETQ L3 SF))
556
(T (SETQ L3 (CONS F L3))))
559
(DEFUN CONTRACT1 (F G) ;This does the actual contraction of F with G.
560
(PROG (A B C D E CF) ;If f has any derivative indices then it can't
561
;contract G. If F is Kronecker delta then see
562
;which of the covariant, contravariant, or
563
;derivative indices matches those in G.
564
(WHEN (CDDDR F) (RETURN NIL))
571
((OR (EQ (CAAR F) '%KDELTA) (EQ (CAAR F) '$KDELTA))
572
(AND (> (LENGTH A) 1) (RETURN NIL))
573
(SETQ A (CAR A) B (CAR B))
575
(SIMPLIFYA (COND ((AND (CDR C) (AND (NOT (NUMBERP B)) (MEMQ B (CDR C))))
576
(SETQ C (SUBST A B (CDR C)))
577
(AND (NOT (MEMQ (CAAR G) CHRISTOFFELS))
579
(SETQ A (CONTRACT2 C (CDR D)))
581
D (CONS SMLIST (CDR A))))
586
((AND E (AND (NOT (NUMBERP B)) (MEMQ B E)))
587
(NCONC (LIST (CAR G) C D)
588
(itensor-SORT (SUBST A B E))))
589
((AND (CDR D) (AND (NOT (NUMBERP A)) (MEMQ A (CDR D))))
590
(SETQ D (SUBST B A (CDR D)))
592
(SETQ A (CONTRACT2 (CDR C) D))
594
C (CONS SMLIST (CAR A))))
601
;; VTT: No tensor should be able to contract LC or KDELTA.
602
(AND (OR (EQ (CAAR G) '$KDELTA) (EQ (CAAR G) '%KDELTA) (EQ (CAAR G) '$LC) (EQ (CAAR G) '%LC)) (RETURN NIL))
604
;If G has derivative indices then F must be
605
(and e ;constant in order to contract it.
606
(NOT (MGET (CAAR F) '$CONSTANT))
608
;Contraction property of F is a list of (A.B)'S
609
;; (COND ((SETQ CF (ZL-GET (CAAR F) 'CONTRACTIONS)))
610
(COND ((SETQ CF (GETCON (CAAR F))))
612
;If G matches an A then use the B for name of result.
613
;If an A is a space use name of G for result.
614
MORE (COND ((EQ (CAAR CF) '/ ) (SETQ CF (CAR G)))
615
((EQ (CAAR CF) (CAAR G))
616
(SETQ CF (NCONS (CDAR CF))))
617
(T (OR (SETQ CF (CDR CF)) (RETURN NIL)) (GO MORE)))
618
(SETQ C (CDR C) D (CDR D))
619
;If CONTRACT2 of F's contravariant and G's
620
;covariant or F's covariant and G's
621
;contravariant indicies is NIL then return NIL.
622
(COND ((AND B C (SETQ F (CONTRACT2 B C)))
623
(SETQ B (CAR F) C (CDR F)))
624
((AND A D (SETQ F (CONTRACT2 A D)))
625
(SETQ A (CAR F) D (CDR F)))
627
;Form combined indices of result
628
(AND D (SETQ B (APPEND B D)))
629
(AND C (SETQ A (APPEND C A)))
630
;Zl-remove repeated indices
631
(AND (SETQ F (CONTRACT2 A B)) (SETQ A (CAR F) B (CDR F)))
632
;; VTT: Special handling of Christoffel symbols. We can only contract them
633
;; when we turn CHR1 into CHR2 or vice versa; other index combinations are
634
;; illegal. This code checks if the index pattern is a valid one and replaces
635
;; CHR1 with CHR2 or vice versa as appropriate.
636
(COND ((OR (EQ (CAR CF) '$CHR1) (EQ (CAR CF) '%CHR1))
637
(COND ((AND (EQ (LENGTH A) 2) (EQ (LENGTH B) 1))
638
(SETQ CF (CONS (COND ((EQ (CAR CF) '$CHR1) '$CHR2) (T '%CHR2)) (CDR CF))))
639
((NOT (AND (EQ (LENGTH A) 3) (EQ (LENGTH B) 0))) (RETURN NIL)))
641
((OR (EQ (CAR CF) '$CHR2) (EQ (CAR CF) '%CHR2))
642
(COND ((AND (EQ (LENGTH A) 3) (EQ (LENGTH B) 0))
643
(SETQ CF (CONS (COND ((EQ (CAR CF) '$CHR2) '$CHR1) (T '%CHR1)) (CDR CF)))
645
((NOT (AND (EQ (LENGTH A) 2) (EQ (LENGTH B) 1))) (RETURN NIL)))
649
(SETQ F (MEVAL (LIST CF (CONS SMLIST A) (CONS SMLIST B))))
653
; (SETQ F (SDIFF F (CAR E))))
656
(SETQ F (SDIFF F (CAR E))))
665
(NCONC (LIST '(%DERIVATIVE)
666
(LIST (CAR X) (CADR X) (CADDR X)))
667
(PUTINONES (CDDDR X))))
670
(MYSUBST0 (SIMPLIFYA (CONS (NCONS (CAAR X))
671
(MAPCAR '$UNDIFF (CDR X)))
675
(COND ((CDR E) (CONS (CAR E) (CONS 1. (PUTINONES (CDR E)))))
676
(T (LIST (CAR E) 1.))))
678
(DEFMFUN $KDELTA (L1 L2)
679
(COND ((NULL (AND ($LISTP L1)
681
(= (LENGTH L1) (LENGTH L2))))
682
(merror "Improper arg to DELTA: ~M"
683
(LIST '(%KDELTA) L1 L2)
685
(T (DELTA (CDR L1) (CDR L2)))))
1047
loop2(and (null cl) (return (nconc l1 l3)))
1049
(and (null l3) (return (nconc l1 cl)))
1050
(setq f (car cl) cl (cdr cl))
1051
(cond ((setq sf (contract3 f l3))
1052
;;*** contract3 may also return a negative result
1053
(setq sf (mapcar #'(lambda (x)
1054
(cond ((atom x) x) (
1055
(and (or (equal (car x) '(mtimes)) (equal (car x) '(mtimes simp))) (eq (cadr x) -1))
1056
(setq l1 (cons -1 l1)) (caddr x)) (t x))
1060
(t (setq l3 (cons f l3))))
1063
;; Create a 'normalized' (i.e., old-style) rpobj
1064
(defmfun $renorm (e &optional (force nil))
1066
(and (not (rpobj e)) (merror "Not an RPOBJ: ~M" e))
1067
(and $allsym (setq force t))
1068
(setq c (cdaddr e) v nil)
1070
((i (reverse (cdadr e)) (cdr i)))
1072
(or (null i) (and (atom (car i)) (not force))) ; Terminating condition
1073
(setq v (append (reverse i) v)) ; Remaining covariant indices
1076
((atom (car i)) (setq v (cons (car i) v)))
1077
(t (setq c (cons (caddar i) c)))
1081
(cons (car e) (append (list (cons smlist v) (cons smlist c)) (cdddr e)))
1086
;; As above, but unconditionally. Not needed.
1087
;(defun renorm (e) (append (list (car e) ($covi e) ($conti e)) (cdddr e)))
1089
;; Add a minus sign to all elements in a list
1091
(cond ((null l) nil)
1092
(t (cons (list '(mtimes simp) -1 (car l)) (neglist (cdr l))))
1096
;; Create an 'abnormal' (i.e., new-style) rpobj
1098
(append (list (car e)
1099
(append ($covi e) (neglist (conti e)))
1105
;; Test for membership using EQUAL, to catch member lists
1106
(defun memlist (e l)
1107
(cond ((null l) nil)
1108
((equal e (car l)) l)
1109
(t (memlist e (cdr l)))
1113
;; Substitute using EQUAL, to catch member lists
1114
(defun substlist (b a l)
1116
((equal a (car l)) (cons b (cdr l)))
1117
(t (cons (car l) (substlist b a (cdr l))))
1121
;; And delete an element from a list, again using EQUAL
1122
(defun dellist (e l)
1124
((equal e (car l)) (dellist e (cdr l)))
1125
(t (cons (car l) (dellist e (cdr l))))
1129
;; Removes items not in i from l.
1130
(defun removenotin (i l)
1132
((memq (car l) i) (cons (car l) (removenotin i (cdr l))))
1133
(t (removenotin i (cdr l)))
1137
;; Removes items not in i from l. But the ones in l have a minus sign!
1138
(defun removenotinm (i l)
1140
((atom (car l)) (cons (car l) (removenotinm i (cdr l))))
1141
((and (isprod (caar l)) (eq (cadar l) -1)
1142
(not (memq (caddar l) i))) (removenotinm i (cdr l)))
1143
(t (cons (car l) (removenotinm i (cdr l))))
1147
;; Removes indices duplicated once with and once without a minus sign
1148
(defun contractinside (c)
1150
((i (minusi c) (cdr i)))
1152
(and (memlist (car i) c) (memlist (list '(mtimes simp) -1 (car i)) c)
1153
(setq c (delete (car i) (dellist (list '(mtimes simp) -1 (car i)) c)))
1159
;; This does the actual contraction of f with g. If f has any derivative
1160
;; indices then it can't contract g. If f is Kronecker delta then see which of
1161
;; the covariant, contravariant, or derivative indices matches those in g.
1162
(defun contract1 (f g)
1163
(prog (a b c d e cf sgn)
1164
(when (cdddr f) (return nil))
1165
(setq a (derat (cdadr f)) b (cdaddr f)
1166
c (derat (cadr g)) d (caddr g) e (cdddr g)
1168
(cond ; This section is all Kronecker-delta code
1170
(or (eq (caar f) '%kdelta) (eq (caar f) '$kdelta))
1172
; We normalize the indices first
1173
(setq b (append (minusi a) b) a (plusi a))
1175
;We cannot contract with higher-order or malformed Kronecker deltas
1176
(and (or (/= (length a) 1) (/= (length b) 1 )) (return nil))
1178
(setq a (car a) b (car b))
1183
(and (cdr c) (not (numberp b)) (memq b (cdr c)))
1184
(setq c (subst a b (cdr c)))
1186
(not (memq (caar g) christoffels))
1188
(setq a (contract2 c (cdr d)))
1189
(setq c (car a) d (cons smlist (cdr a)))
1191
(setq c (contractinside c))
1192
(nconc (list (car g) (cons smlist c) d) e)
1195
(and e (not (numberp b)) (memq b e))
1196
(nconc (list (car g) c d)
1198
($iframe_flag (subst a b e))
1199
(t (itensor-sort (subst a b e)))
1204
(and (cdr d) (not (numberp a)) (memq a (cdr d)))
1205
(setq d (subst b a (cdr d)))
1208
(setq a (contract2 (cdr c) d))
1209
(setq d (cdr a) c (cons smlist (car a)))
1211
(nconc (list (car g) c (cons smlist d)) e)
1214
(and (cdr c) (not (numberp a))
1215
(memlist (list '(mtimes simp) -1 a) (cdr c))
1217
(setq c (substlist (list '(mtimes simp) -1 b)
1218
(list '(mtimes simp) -1 a)
1222
(setq c (contractinside c))
1223
(nconc (list (car g) (cons smlist c) d) e)
1233
;No tensor can contract Kronecker-deltas or Levi-Civita symbols.
1235
(or (eq (caar g) '$kdelta) (eq (caar g) '%kdelta)
1236
(eq (caar g) '$levi_civita) (eq (caar g) '%levi_civita)
1237
(eq (caar g) '$icurvature) (eq (caar g) '%icurvature)
1242
;If g has derivative indices then F must be constant in order to contract it
1243
(and e (not (mget (caar f) '$constant)) (return nil))
1245
;Contraction property of f is a list of (a.b)'s
1247
((setq cf (getcon (caar f))))
1251
; Determine the sign of the result based on the expression's symmetry
1252
; properties. We use CANFORM to sort indices in the canonical order
1253
; and then extract the resulting expression's sign.
1255
(cond ((eq -1 (cadr ($canform (list '(mtimes simp) f g)))) -1) (t 1))
1258
;If g matches an a then use the b for name of result. If an a is a space
1259
;use name of G for result.
1267
(eq (caar cf) (caar g))
1268
(setq cf (ncons (cdar cf)))
1271
(or (setq cf (cdr cf)) (return nil))
1275
(setq c (cdr c) d (cdr d))
1277
;If CONTRACT2 of f's contravariant and g's covariant or f's covariant and
1278
;g's contravariant indices is nil then return nil
1281
(and b c (setq f (contract2 b c)))
1282
(setq b (car f) c (cdr f))
1285
(and a d (setq f (contract2 a d)))
1286
(setq a (car f) d (cdr f))
1289
(and a (minusi c) (setq f (contract2 a (minusi c))))
1290
; (cdr f) now contains the free indices in (minusi c).
1291
; what we need to do is find the corresponding items in c, and remove
1292
; all other negative indices (i.e., those that were dropped by
1294
; What we need to do is remove items from c one by one, and substitute
1295
; an item from (car f), which we should remove from (car f):
1296
; for i thru length(c)
1297
; if c[i] not in (cdr f)
1298
; if (car f) is nil, remove c[i]
1299
; otherwise subst c[i]
1301
; Now set c to what we made of c, a to whatever is left of (cdr f)
1309
((null i) (setq a (removenotin j a) c (reverse k)))
1312
(or (atom (car i)) (member (caddar i) (cdr f)))
1313
(setq k (cons (car i) k))
1317
(setq k (cons (car j) k) j (cdr j))
1323
(and (minusi a) c (setq f (contract2 (minusi a) c)))
1330
;; ((null i) (setq c (reverse k) a (append (plusi a) j)))
1336
(mapcar #'(lambda (x) (list '(mtimes simp) -1 x)) j)
1341
((member (car i) (cdr f)) (setq k (cons (car i) k)))
1344
(setq k (cons (list '(mtimes simp) -1 (car j)) k) j (cdr j))
1351
;Form combined indices of result
1352
(and d (setq b (append b d)))
1353
(and c (setq a (append c a)))
1354
;Zl-remove repeated indices
1355
;; (and (setq f (contract2 a b)) (setq a (car f) b (cdr f)))
1356
;; (setq a (contractinside a))
1358
;VTT: Special handling of Christoffel symbols. We can only contract them
1359
;when we turn ICHR1 into ICHR2 or vice versa; other index combinations are
1360
;illegal. This code checks if the index pattern is a valid one and replaces
1361
;ICHR1 with ICHR2 or vice versa as appropriate.
1364
(member (car cf) christoffels1)
1367
(and (eq (length a) 2) (eq (length b) 1))
1370
(elt christoffels2 (position (car cf) christoffels1))
1376
(not (and (eq (length a) 3) (eq (length b) 0)))
1382
(member (car cf) christoffels2)
1385
(and (eq (length a) 3) (eq (length b) 0))
1388
(elt christoffels1 (position (car cf) christoffels2))
1394
(not (and (eq (length a) 2) (eq (length b) 1)))
1399
((member (car cf) christoffels) (return nil))
1402
(setq f (meval (list cf (cons smlist a) (cons smlist b))))
1407
(setq f (idiff f (car e)))
1410
(return (cond ((eq sgn -1) (list '(mtimes) sgn f)) (t f)))
1414
;; In what amounts to quite an abuse of the Kronecker delta concept, we
1415
;; permit an exceptional index combination of two contravariant indices.
1416
;; This helps lc2kdt convert Levi-Civita symbols in a manner that does
1417
;; not require resorting to numeric indices, causing all sorts of problems
1418
;; with RENAME and CONTRACT.
1419
(defmfun $kdelta (l1 l2)
1420
(setq l2 (append l2 (minusi l1)) l1 (plusi l1))
1423
(and ($listp l1) ($listp l2) (= ($length l1) 0) (= ($length l2) 2))
1425
((eq (cadr l2) (caddr l2)) 1)
1427
(and (numberp (cadr l2)) (numberp (caddr l2)))
1429
((= (cadr l2) (caddr l2)) t)
1433
(t (list '(%kdelta) l1 l2))
1437
(and ($listp l1) ($listp l2) (= ($length l1) 2) (= ($length l2) 0))
1439
((eq (cadr l1) (caddr l1)) 1)
1441
(and (numberp (cadr l1)) (numberp (caddr l1)))
1443
((= (cadr l1) (caddr l1)) t)
1447
(t (list '(%kdelta) l1 l2))
1451
(null (and ($listp l1) ($listp l2) (= (length l1) (length l2))))
1452
(merror "Improper arg to DELTA: ~M" (list '(%kdelta) l1 l2))
1454
(t (delta (cdr l1) (cdr l2)))
687
1458
;kdels defines the symmetric combination of the Kronecker symbols
689
(DEFMFUN $KDELS (L1 L2)
690
(COND ((NULL (AND ($LISTP L1)
692
(= (LENGTH L1) (LENGTH L2))))
1460
(defmfun $kdels (l1 l2)
1461
(cond ((null (and ($listp l1)
1463
(= (length l1) (length l2))))
693
1464
(merror "Improper arg to DELTA: ~M"
694
(LIST '(%KDELS) L1 L2)
1465
(list '(%kdels) l1 l2)
696
(T (DELTA (CDR L1) (CDR L2) 1))))
698
;;(DECLARE-TOP (FIXNUM I))
700
;;(DEFUN DELTA (LOWER UPPER &optional (eps -1))
701
;; (COND ((NULL LOWER) $DIM)
702
;; ((NULL (CDR LOWER))
703
;; (COND ((EQUAL (CAR UPPER) (CAR LOWER))
704
;; (COND ((NUMBERP (CAR UPPER)) 1.) (T $DIM)))
705
;; ((AND (NUMBERP (CAR UPPER)) (NUMBERP (CAR LOWER))) 0.)
706
;; (T (LIST '(%KDELTA)
707
;; (CONS SMLIST LOWER)
708
;; (CONS SMLIST UPPER)))))
709
;; (T (DO ((I (LENGTH LOWER) (1- I))
713
;; (F (NCONS (CAR UPPER)))
715
;; (SIGN (ODDP (LENGTH LOWER))))
717
;; (SIMPLUS (CONS '(MPLUS) RESULT) 1. T))
718
;; (SETQ TERM (LIST (DELTA (NCONS (CAR SL)) F eps)
719
;; (DELTA (CDR SL) R eps)))
720
;; (SETQ SL (CDR (APPEND SL (NCONS (CAR SL)))))
722
;; (CONS (SIMPTIMES (CONS '(MTIMES)
1467
(t (delta (cdr l1) (cdr l2) 1))))
1469
;;(declare-top (fixnum i))
1471
;;(defun delta (lower upper &optional (eps -1))
1472
;; (cond ((null lower) $dim)
1473
;; ((null (cdr lower))
1474
;; (cond ((equal (car upper) (car lower))
1475
;; (cond ((numberp (car upper)) 1.) (t $dim)))
1476
;; ((and (numberp (car upper)) (numberp (car lower))) 0.)
1477
;; (t (list '(%kdelta)
1478
;; (cons smlist lower)
1479
;; (cons smlist upper)))))
1480
;; (t (do ((i (length lower) (1- i))
1484
;; (f (ncons (car upper)))
1486
;; (sign (oddp (length lower))))
1488
;; (simplus (cons '(mplus) result) 1. t))
1489
;; (setq term (list (delta (ncons (car sl)) f eps)
1490
;; (delta (cdr sl) r eps)))
1491
;; (setq sl (cdr (append sl (ncons (car sl)))))
1493
;; (cons (simptimes (cons '(mtimes)
731
(DEFUN DELTA (LOWER UPPER &optional (eps -1))
732
(COND ((NULL LOWER) $DIM)
734
(COND ((EQUAL (CAR UPPER) (CAR LOWER))
735
(COND ((NUMBERP (CAR UPPER)) 1.) (T $DIM)))
736
((AND (NUMBERP (CAR UPPER)) (NUMBERP (CAR LOWER))) 0.)
737
(T (LIST '(%KDELTA) (CONS SMLIST LOWER) (CONS SMLIST UPPER)))))
738
(T (DO ((LEFT NIL (APPEND LEFT (NCONS (CAR RIGHT))))
739
(RIGHT LOWER (CDR RIGHT))
741
((NULL RIGHT) (SIMPLUS (CONS '(MPLUS) RESULT) 1. T))
742
(SETQ RESULT (CONS (SIMPTIMES
743
(LIST '(MTIMES) (DELTA (NCONS (CAR RIGHT)) (NCONS (CAR UPPER)) eps)
744
(DELTA (APPEND LEFT (CDR RIGHT)) (CDR UPPER) eps)
745
(COND ((ODDP (LENGTH LEFT)) eps) (T 1))
1502
(defun delta (lower upper &optional (eps -1))
1503
(cond ((null lower) $dim)
1505
(cond ((equal (car upper) (car lower))
1506
(cond ((numberp (car upper)) 1.) (t $dim)))
1507
((and (numberp (car upper)) (numberp (car lower))) 0.)
1508
(t (list '(%kdelta) (cons smlist lower) (cons smlist upper)))))
1509
(t (do ((left nil (append left (ncons (car right))))
1510
(right lower (cdr right))
1512
((null right) (simplus (cons '(mplus) result) 1. t))
1513
(setq result (cons (simptimes
1514
(list '(mtimes) (delta (ncons (car right)) (ncons (car upper)) eps)
1515
(delta (append left (cdr right)) (cdr upper) eps)
1516
(cond ((oddp (length left)) eps) (t 1))
750
(DECLARE-TOP (NOTYPE I))
1521
(declare-top (notype i))
752
(DECLARE-TOP (SPECIAL $OUTCHAR $DISPFLAG LINELABLE FOOBAR DERIVLIST))
1523
(declare-top (special $outchar $dispflag linelable foobar derivlist))
754
1526
;Displays P([L1],[L2],I1,I2,...) by making the elements of L2 into a single
755
1527
;atom which serves as the exponent and the elements of L1 and I1,I2,... into a
756
1528
;single atom with a comma in between which serves as the subscript.
759
(progn (makelabel $LINECHAR)
761
(displa (list '(MLABLE) LINELABLE (ishow (specrepcheck f))))
762
; (setq $DISPFLAG nil)
1531
(progn (makelabel $linechar)
1533
(displa (list '(mlable) linelable (ishow (specrepcheck (derat f)))))
1534
; (setq $dispflag nil)
767
((LAMBDA (FOOBAR) ;FOOBAR intialized to NIL
769
((RPOBJ F) ;If an indexed object ...
771
(COND ((OR (CDADR F) (CDDDR F)) ;If covariant or
772
(CONS (LIST (CAAR F) ;derivative indices
774
(NCONS (MAKNAM (CONS '$ (SPLICE (CDADR F)
777
(COND ((CDADDR F) ;If contravariant indices
780
(CONS '(MTIMES SIMP) ;Make indices appear
781
(CDADDR F)))) ;as exponents for
782
(T FOOBAR))) ;proper display
784
(CONS (CAR F) (MAPCAR 'ISHOW (CDR F))))))
785
NIL)) ;Map onto subparts of F
787
(DEFUN SPLICE (L1 L2)
788
(COND (L2 (SETQ L2 (CONS '|,| (SPLICE1 L2)))
789
(AND L1 (SETQ L2 (NCONC (SPLICE1 L1) L2)))
794
(COND ((NULL (CDR L))(SPLICE2 (CAR L)))
795
(T (NCONC (SPLICE2 (CAR L))(CONS '| | (SPLICE1 (CDR L)))))))
798
(COND ((FIXP X)(EXPLODE X))
799
(T (CDR (EXPLODEc X)))))
1539
((lambda (foobar) ;FOOBAR intialized to NIL
1541
((rpobj f) ;If an indexed object ...
1543
(cond ((or (covi f) (cdddr f)) ;If covariant or
1544
(cons (list (caar f) ;derivative indices
1546
(ncons (maknam (cons '$ (splice (covi f)
1549
(cond ((conti f) ;If contravariant indices
1552
(cons '(mtimes simp) ;Make indices appear
1553
(conti f)))) ;as exponents for
1554
(t foobar))) ;proper display
1556
(cons (car f) (mapcar 'ishow (cdr f))))))
1557
nil)) ;Map onto subparts of F
1559
(defun splice (l1 l2)
1560
(cond (l2 (setq l2 (cons '|,| (splice1 l2)))
1561
(and l1 (setq l2 (nconc (splice1 l1) l2)))
1566
(cond ((null (cdr l))(splice2 (car l)))
1567
(t (nconc (splice2 (car l))(cons '| | (splice1 (cdr l)))))))
1570
(cond ((fixp x)(explode x))
1571
(t (cdr (explodec x)))))
1572
; (t (cdr (explodec (print-invert-case x))))))
802
(PROG (EXP Z COUNT V)
803
(COND ((NULL (CDR E)) (RETURN (STOTALDIFF (CAR E))))
804
((NULL (CDDR E)) (NCONC E '(1.))))
805
(SETQ EXP (CAR E) Z (SETQ E (APPEND E NIL)))
806
LOOP (COND ((OR (NULL DERIVLIST) (ZL-MEMBER (CADR Z) DERIVLIST))
1575
(prog (exp z count v)
1576
(cond ((null (cdr e)) (return (stotaldiff (car e))))
1577
((null (cddr e)) (nconc e '(1.))))
1578
(setq exp (car e) z (setq e (append e nil)))
1579
loop (cond ((or (null derivlist) (zl-member (cadr z) derivlist))
808
1581
;DERIVLIST is set by $EV
810
LOOP2(COND ((CDR Z) (GO LOOP))
811
((NULL (CDR E)) (RETURN EXP))
813
DOIT (COND ((NULL (CDDR Z))
1583
loop2(cond ((cdr z) (go loop))
1584
((null (cdr e)) (return exp))
1586
doit (cond ((null (cddr z))
814
1587
(merror "Wrong number of args to DERIVATIVE"))
815
((NOT (FIXP (SETQ COUNT (CADDR Z)))) (GO NOUN))
1588
((not (fixp (setq count (caddr z)))) (go noun))
817
1590
(merror "Improper count to DIFF: ~M"
819
LOOP1(SETQ V (CADR Z))
825
(COND ((ATOM $COORDINATES)
826
(MEVAL1 (LIST (LIST $COORDINATES 'SIMP 'ARRAY)
828
((EQ (CAAR $COORDINATES) 'MLIST)
830
(LENGTH $COORDINATES)))
1592
loop1(setq v (cadr z))
1598
(cond ((atom $vect_coords)
1599
(meval1 (list (list $vect_coords 'simp 'array)
1601
((eq (caar $vect_coords) 'mlist)
1603
(length $vect_coords)))
832
1605
"Coordinate list too short for derivative index"))
833
(T (NTH V $COORDINATES))))
835
(COND ((ZEROP COUNT) (RPLACD Z (CDDDR Z)) (GO LOOP2))
836
((ZEROP1 (SETQ EXP (SDIFF EXP V))) (RETURN 0.)))
837
(SETQ COUNT (1- COUNT))
839
NOUN (RETURN (DIFF%DERIV (CONS EXP (CDR E))))))
841
(DEFUN CHAINRULE1 (E X) ; --YS 15.02.02
843
(COND ((AND (ATOM E) (EQ (SETQ Y (CAR (MGET E 'DEPENDS)))
844
(CADR $COORD))) (RETURN (SUBST X Y (CHAINRULE E Y))))
845
(T (RETURN (CHAINRULE E X))))))
1606
(t (nth v $vect_coords))))
1608
(cond ((zerop count) (rplacd z (cdddr z)) (go loop2))
1609
((zerop1 (setq exp (sdiff exp v))) (return 0.)))
1610
(setq count (1- count))
1612
noun (return (diff%deriv (cons exp (cdr e))))))
1614
(defun chainrule1 (e x) ; --ys 15.02.02
1616
(cond ((and (atom e) (eq (setq y (car (mget e 'depends)))
1617
(cadr $coord))) (return (subst x y (chainrule e y))))
1618
(t (return (chainrule e x))))))
1620
(defun diffexpt1 (e x)
1621
;; RETURN: n*v^n*rename(v'/v) where e=v^n
1622
(list '(mtimes) (caddr e) e
1624
(list '(mtimes) (list '(mexpt) (cadr e) -1)
847
1631
;Redefined so that the derivative of any indexed object appends on the
848
1632
;coordinate index in sorted order unless the indexed object was declared
849
1633
;constant in which case 0 is returned.
850
1634
#+Franz (sstatus translink nil) ; make sdiff take hold
851
1635
#+Franz (sstatus translink t)
855
((OR (ATOM E) (MEMQ 'ARRAY (CDAR E)))
857
((MGET (CAAR E) '$CONSTANT) 0.) ;New line added
858
((EQ (CAAR E) 'MRAT) (RATDX E X))
859
((EQ (CAAR E) 'MPLUS)
860
(SIMPLUS (CONS '(MPLUS) (SDIFFMAP (CDR E) X))
1637
(cond ((mnump e) 0.)
1639
((or (atom e) (memq 'array (cdar e)))
1641
((mget (caar e) '$constant) 0.) ;New line added
1642
((eq (caar e) 'mrat) (ratdx e x))
1643
((eq (caar e) 'mplus)
1644
(simplus (cons '(mplus) (sdiffmap (cdr e) x))
863
((EQ (CAAR E) 'MEQUAL)
864
(LIST (CAR E) (SDIFF (CADR E) X) (SDIFF (CADDR E) X)))
865
((EQ (CAAR E) '$MATRIX)
868
(FUNCTION (LAMBDA (Y)
870
(SDIFFMAP (CDR Y) X))))
872
((EQ (CAAR E) 'MTIMES)
873
(ADDN (SDIFFTIMES (CDR E) X) T))
874
((EQ (CAAR E) 'MEXPT) (DIFFEXPT E X))
875
((RPOBJ E) (DIFFRPOBJ E X)) ;New line added
876
((AND (BOUNDP '$METRIC) (EQ (CAAR E) '%DETERMINANT);New line added
877
(EQ (CADR E) $METRIC))
879
(setq dummy ($dummy))
880
(COND ((EQ DUMMY X) (setq dummy ($dummy))))
881
(LIST '(MTIMES SIMP) 2. E
882
(LIST '($CHR2 SIMP) (CONS SMLIST (LIST DUMMY X))
883
(CONS SMLIST (NCONS DUMMY)))))
886
(COND ((FIXP X) (LIST '(%DERIVATIVE) E X))
888
(T (LIST '(%DERIVATIVE E X)))))
1647
((eq (caar e) 'mequal)
1648
(list (car e) (sdiff (cadr e) x) (sdiff (caddr e) x)))
1649
((eq (caar e) '$matrix)
1652
(function (lambda (y)
1654
(sdiffmap (cdr y) x))))
1656
((eq (caar e) 'mtimes)
1657
(addn (sdifftimes (cdr e) x) t))
1658
((eq (caar e) 'mexpt) (diffexpt1 e x))
1659
;; ((rpobj e) (diffrpobj e x)) ;New line added
1660
;; ((and (boundp '$imetric) (eq (caar e) '%determinant);New line added
1661
;; (eq (cadr e) $imetric))
1663
;; (setq dummy ($idummy))
1664
;; (cond ((eq dummy x) (setq dummy ($idummy))))
1665
;; (list '(mtimes simp) 2. e
1666
;; (list '($ichr2 simp) (cons smlist (list dummy x))
1667
;; (cons smlist (ncons dummy)))))
1669
((not (depends e x))
1670
(cond ((fixp x) (list '(%derivative) e x))
1672
(t (list '(%derivative) e x))))
889
1673
;This line moved down
890
((EQ (CAAR E) 'MNCTIMES)
891
(SIMPLUS (LIST '(MPLUS)
897
(SDIFF (CADDR E) X)))
900
((EQ (CAAR E) 'MNCEXPT) (DIFFNCEXPT E X))
901
((EQ (CAAR E) '%INTEGRATE) (DIFFINT E X))
902
((EQ (CAAR E) '%DERIVATIVE)
903
(COND ((OR (ATOM (CADR E))
904
(MEMQ 'ARRAY (CDAADR E)))
906
((FREEL (CDR E) X) 0.)
907
(T (DIFF%DERIV (LIST E X 1.)))))
908
((MEMQ (CAAR E) '(%SUM %PRODUCT)) (DIFFSUMPROD E X))
909
(T (SDIFFGRAD E X))))
911
(defun DIFFRPOBJ (e x) ;Derivative of an indexed object
912
(cond ((and (memq (caar e) $COORD) (null (cdadr e))
913
(equal (length (cdaddr e)) 1) (null (cdddr e)))
914
(delta (ncons x) (cdaddr e)))
915
(t (NCONC (LIST (CAR E) (CADR E) (CADDR E))
916
(COND ((NULL (CDDDR E)) (NCONS X))
917
(T (itensor-SORT (APPEND (CDDDR E) (NCONS X)))))))))
924
(IFNOT (AND A (CDR A)) (RETURN (LIST '(%LC) L1)))
926
LOOP1(IFNOT (FIXP (CAR A)) (RETURN (LIST '(%LC) L1)))
927
(AND (SETQ A (CDR A)) (GO LOOP1))
928
LOOP3(SETQ A (CAR B) B (CDR B) C B)
929
LOOP2(COND ((= (CAR C) A) (RETURN 0.))
930
((< (CAR C) A) (SETQ SIGN (NOT SIGN))))
931
(AND (SETQ C (CDR C)) (GO LOOP2))
932
(AND (CDR B) (GO LOOP3))
933
(RETURN (COND (SIGN -1.) (T 1.)))))
934
(DEFMFUN $LC (L1 &optional (L2 nil))
936
((EQ L2 nil) ($LC0 L1))
937
((LIKE L1 '((MLIST)))
938
(PROG (l) (SETQ l nil)
939
(DO ((I ($LENGTH L2) (1- I))) ((< I 1)) (SETQ l (CONS I l)))
940
(RETURN (LIST '($KDELTA SIMP) (CONS SMLIST l) L2))
1674
((eq (caar e) 'mnctimes)
1675
(simplus (list '(mplus)
1681
(sdiff (caddr e) x)))
1684
((eq (caar e) 'mncexpt) (diffncexpt e x))
1685
((eq (caar e) '%integrate) (diffint e x))
1686
((eq (caar e) '%derivative)
1687
(cond ((or (atom (cadr e))
1688
(memq 'array (cdaadr e)))
1690
((freel (cdr e) x) 0.)
1691
(t (diff%deriv (list e x 1.)))))
1692
((memq (caar e) '(%sum %product)) (diffsumprod e x))
1693
(t (sdiffgrad e x))))
1695
; VTT: several of these functions have been copied verbatim from comm.lisp and
1696
; comm2.lisp, in order to implement indicial differentiation as distinct from
1697
; differentiation with respect to an external variable.
1699
(defun idiffmap (e x) (mapcar #'(lambda (term) (idiff term x)) e))
1701
(defun idifftimes (l x)
1702
(prog (term left out)
1703
loop (setq term (car l) l (cdr l))
1704
(setq out (cons (muln (cons (idiff term x) (append left l)) t) out))
1705
(if (null l) (return out))
1706
(setq left (cons term left))
1709
(defun idiffexpt1 (e x)
1710
;; RETURN: n*v^n*rename(v'/v) where e=v^n
1711
(list '(mtimes) (caddr e) e
1713
(list '(mtimes) (list '(mexpt) (cadr e) -1)
1720
(defun idiffexpt (e x)
1721
(if (mnump (caddr e))
1722
(mul3 (caddr e) (power (cadr e) (addk (caddr e) -1)) (idiff (cadr e) x))
1723
(mul2 e (add2 (mul3 (power (cadr e) -1) (caddr e) (idiff (cadr e) x))
1724
(mul2 (simplifya (list '(%log) (cadr e)) t)
1725
(idiff (caddr e) x))))))
1727
(defmfun idiffint (e x)
1729
(cond ((null (cdddr e))
1730
(cond ((alike1 x (caddr e)) (cadr e))
1731
((and (not (atom (caddr e))) (atom x) (not (free (caddr e) x)))
1732
(mul2 (cadr e) (idiff (caddr e) x)))
1733
((or ($constantp (setq a (idiff (cadr e) x)))
1734
(and (atom (caddr e)) (free a (caddr e))))
1736
(t (simplifya (list '(%integrate) a (caddr e)) t))))
1737
((alike1 x (caddr e)) (addn (idiffint1 (cdr e) x x) t))
1738
(t (addn (cons (if (equal (setq a (idiff (cadr e) x)) 0)
1740
(simplifya (list '(%integrate) a (caddr e)
1741
(cadddr e) (car (cddddr e)))
1743
(idiffint1 (cdr e) x (caddr e)))
1746
(defun idiffint1 (e x y)
1747
(let ((u (idiff (cadddr e) x)) (v (idiff (caddr e) x)))
1748
(list (if (pzerop u) 0 (mul2 u (maxima-substitute (cadddr e) y (car e))))
1749
(if (pzerop v) 0 (mul3 v (maxima-substitute (caddr e) y (car e)) -1)))))
1751
(defun idiff%deriv (e) (let (derivflag) (simplifya (cons '(%idiff) e) t)))
1755
(cond ((null e) (wna-err '$idiff))
1756
((null (cdr e)) (return (stotaldiff (car e))))
1757
((null (cddr e)) (nconc e '(1))))
1758
(setq exp (car e) z (setq e (copy-top-level e)))
1759
loop (if (or (null derivlist) (zl-member (cadr z) derivlist)) (go doit))
1760
; DERIVLIST is set by $EV
1762
loop2(cond ((cdr z) (go loop))
1763
((null (cdr e)) (return exp))
1765
doit (cond ((nonvarcheck (cadr z) '$idiff))
1766
((null (cddr z)) (wna-err '$idiff))
1767
((not (eq (ml-typep (caddr z)) 'fixnum)) (go noun))
1768
((minusp (setq count (caddr z)))
1769
(merror "Improper count to IDIFF:~%~M" count)))
1770
loop1(cond ((zerop count) (rplacd z (cdddr z)) (go loop2))
1771
((equal (setq exp (idiff exp (cadr z))) 0) (return 0)))
1772
(setq count (f1- count))
1774
noun (return (idiff%deriv (cons exp (cdr e))))))
1777
(defmfun idiffncexpt (e x)
1778
((lambda (base* pow)
1779
(cond ((and (mnump pow) (or (not (eq (ml-typep pow) 'fixnum)) (< pow 0))) ; POW cannot be 0
1780
(idiff%deriv (list e x 1)))
1781
((and (atom base*) (eq base* x) (free pow base*))
1782
(mul2* pow (list '(mncexpt) base* (add2 pow -1))))
1783
((ml-typep pow 'fixnum)
1784
((lambda (deriv ans)
1785
(do ((i 0 (f1+ i))) ((= i pow))
1786
(setq ans (cons (list '(mnctimes) (list '(mncexpt) base* i)
1787
(list '(mnctimes) deriv
1788
(list '(mncexpt) base* (f- pow 1 i))))
1791
(idiff base* x) nil))
1792
((and (not (depends pow x)) (or (atom pow) (and (atom base*) (free pow base*))))
1793
((lambda (deriv index)
1796
(list '(mnctimes) (list '(mncexpt) base* index)
1797
(list '(mnctimes) deriv
1798
(list '(mncexpt) base*
1799
(list '(mplus) pow -1 (list '(mtimes) -1 index)))))
1800
index 0 (list '(mplus) pow -1)) nil))
1801
(idiff base* x) (gensumindex)))
1802
(t (idiff%deriv (list e x 1)))))
1803
(cadr e) (caddr e)))
1805
(defmfun idiffsumprod (e x)
1806
(cond ((or (not (atom x)) (not (free (cadddr e) x)) (not (free (car (cddddr e)) x)))
1807
(idiff%deriv (list e x 1)))
1808
((eq (caddr e) x) 0)
1809
(t (let ((u (idiff (cadr e) x)))
1810
(setq u (simplifya (list '(%sum)
1811
(if (eq (caar e) '%sum) u (div u (cadr e)))
1812
(caddr e) (cadddr e) (car (cddddr e)))
1814
(if (eq (caar e) '%sum) u (mul2 e u))))))
1816
(defun idiffgrad (e x)
1817
(let ((fun (caar e)) grad args)
1818
(cond ((and (eq fun 'mqapply) (oldget (caaadr e) 'grad))
1819
(idiffgrad (cons (cons (caaadr e) nil) (append (cdadr e) (cddr e)))
1821
((or (eq fun 'mqapply) (null (setq grad (oldget fun 'grad))))
1822
(if (not (depends e x)) 0 (idiff%deriv (list e x 1))))
1823
((not (= (length (cdr e)) (length (car grad))))
1824
(merror "Wrong number of arguments for ~:M" fun))
1825
(t (setq args (idiffmap (cdr e) x))
1830
(do ((l1 (cdr grad) (cdr l1))
1831
(args args (cdr args)) (l2))
1832
((null l1) (cons '(mlist) (nreverse l2)))
1833
(setq l2 (cons (cond ((equal (car args) 0) 0)
1839
(defmfun $idiff n (let (derivlist) (ideriv (listify n))))
1841
(defmfun idiff (e x)
1845
((or (atom e) (memq 'array (cdar e)))
1846
;; (ichainrule e x))
1847
;; (idiff%deriv (list e x 1)))
1849
((mget (caar e) '$constant) 0.) ;New line added
1850
((eq (caar e) 'mrat) (ratdx e x))
1851
((eq (caar e) 'mplus)
1852
(simplus (cons '(mplus) (idiffmap (cdr e) x))
1855
((eq (caar e) 'mequal)
1856
(list (car e) ($idiff (cadr e) x) ($idiff (caddr e) x)))
1857
((eq (caar e) '$matrix)
1860
(function (lambda (y)
1862
(idiffmap (cdr y) x))))
1864
((eq (caar e) 'mtimes)
1865
(addn (idifftimes (cdr e) x) t))
1866
((eq (caar e) 'mexpt) (idiffexpt1 e x))
1867
((rpobj e) (diffrpobj e x))
1868
((and (boundp '$imetric) (eq (caar e) '%determinant)
1869
(eq (cadr e) $imetric))
1871
(setq dummy ($idummy))
1872
(cond ((eq dummy x) (setq dummy ($idummy))))
1873
(list '(mtimes simp) 2. e
1874
(list '($ichr2 simp) (cons smlist (list dummy x))
1875
(cons smlist (ncons dummy)))))
1877
((eq (caar e) 'mnctimes)
1878
(simplus (list '(mplus)
1884
($idiff (caddr e) x)))
1887
((eq (caar e) 'mncexpt) (idiffncexpt e x))
1888
((eq (caar e) '%integrate) (idiffint e x))
1889
((eq (caar e) '%derivative)
1890
(cond ((or (atom (cadr e))
1891
(memq 'array (cdaadr e)))
1892
;; (ichainrule e x))
1893
;; (idiff%deriv (list e x 1)))
1895
;; ((freel (cdr e) x) 0.)
1896
(t (idiff%deriv (list e x 1.)))))
1897
((memq (caar e) '(%sum %product)) (idiffsumprod e x))
1902
(defun diffrpobj (e x) ;Derivative of an indexed object
1904
( ; Special case: functions declared with coord()
1906
(memq (caar e) $coord) (null (cdadr e))
1907
(equal (length (cdaddr e)) 1) (null (cdddr e))
1909
(delta (ncons x) (cdaddr e))
1911
(t ; Everything else
1913
(list (car e) (cadr e) (caddr e))
1919
( ; Derivative indices do not commute when frames are used
1920
(or $iframe_flag $itorsion_flag)
1921
(append (cdddr e) (ncons x))
1924
(itensor-sort (append (cdddr e) (ncons x)))
1936
(ifnot (and a (cdr a)) (return (list '(%levi_civita) l1)))
1938
loop1(ifnot (fixp (car a)) (return (list '(%levi_civita) l1)))
1939
(and (setq a (cdr a)) (go loop1))
1940
loop3(setq a (car b) b (cdr b) c b)
1941
loop2(cond ((= (car c) a) (return 0.))
1942
((< (car c) a) (setq sign (not sign))))
1943
(and (setq c (cdr c)) (go loop2))
1944
(and (cdr b) (go loop3))
1945
(return (cond (sign -1.) (t 1.)))))
1946
(defmfun $levi_civita (l1 &optional (l2 nil))
1948
((eq l2 nil) ($lc0 l1))
1949
((like l1 '((mlist)))
1950
(prog (l) (setq l nil)
1951
(do ((i ($length l2) (1- i))) ((< i 1)) (setq l (cons i l)))
1952
(return (list '($kdelta simp) (cons smlist l) l2))
942
((LIKE L2 '((MLIST)))
943
(PROG (l) (SETQ l nil)
944
(DO ((I ($LENGTH L1) (1- I))) ((< I 1)) (SETQ l (CONS I l)))
945
(RETURN (LIST '($KDELTA SIMP) L1 (CONS SMLIST l)))
1954
((like l2 '((mlist)))
1955
(prog (l) (setq l nil)
1956
(do ((i ($length l1) (1- i))) ((< i 1)) (setq l (cons i l)))
1957
(return (list '($kdelta simp) l1 (cons smlist l)))
947
(T (MERROR "Mixed-index Levi-Civita symbols not supported"))
1959
(t (merror "Mixed-index Levi-Civita symbols not supported"))
951
1963
;; simplification rules for the totally antisymmetric LC symbol
955
(COND ((ATOM E) (MATCHERR)))
956
(COND ((ATOM (CAR E)) (MATCHERR)))
957
(COND ((NOT (OR (EQ (CAAR E) '$LC) (EQ (CAAR E) '%LC))) (MATCHERR)))
958
(COND ((NOT ($LISTP (SETQ L1 (CADR E)))) (MATCHERR)))
959
(COND ((NOT (ALIKE1 '((MLIST SIMP)) (SETQ L2 (CADDR E)))) (MATCHERR)))
960
(COND ((CDDDR E) (MATCHERR)))
961
(SETQ N ($LENGTH L1))
963
(DO ((I N (1- I))) ((< I 1)) (SETQ l (CONS ($DUMMY) L)))
964
(RETURN (LIST '(MTIMES SIMP) -1 ($KDELTA L1 (CONS SMLIST L))
965
(LIST (CONS (CAAR E) '(SIMP)) (CONS SMLIST L) (NCONS SMLIST))
966
(LIST '(MEXPT SIMP) (MEVAL (LIST 'MFACTORIAL N)) -1))
975
(COND ((ATOM E) (MATCHERR)))
976
(COND ((ATOM (CAR E)) (MATCHERR)))
977
(COND ((NOT (OR (EQ (CAAR E) '$LC) (EQ (CAAR E) '%LC))) (MATCHERR)))
978
(COND ((NOT (ALIKE1 '((MLIST SIMP)) (SETQ L1 (CADR E)))) (MATCHERR)))
979
(COND ((NOT ($LISTP (SETQ L2 (CADDR E)))) (MATCHERR)))
980
(COND ((CDDDR E) (MATCHERR)))
981
(SETQ N ($LENGTH L2))
983
(DO ((I N (1- I))) ((< I 1)) (SETQ l (CONS ($DUMMY) L)))
984
(RETURN (LIST '(MTIMES SIMP) -1 ($KDELTA (CONS SMLIST L) L2)
985
(LIST (CONS (CAAR E) '(SIMP)) (NCONS SMLIST) (CONS SMLIST L))
986
(LIST '(MEXPT SIMP) (MEVAL (LIST 'MFACTORIAL N)) -1))
992
(ADD2LNC '$LC_L $RULES)
993
(ADD2LNC '$LC_U $RULES)
1967
(cond ((atom e) (matcherr)))
1968
(cond ((atom (car e)) (matcherr)))
1969
(cond ((not (or (eq (caar e) '$levi_civita) (eq (caar e) '%levi_civita))) (matcherr)))
1970
(cond ((not ($listp (setq l1 ($covi e)))) (matcherr)))
1971
(cond ((not (alike1 '((mlist simp)) (setq l2 ($conti e)))) (matcherr)))
1972
(cond ((cdddr e) (matcherr)))
1973
(setq nn ($length l1))
1975
(do ((i nn (1- i))) ((< i 1)) (setq l (cons ($idummy) l) n $icounter))
1976
(return (values (list '(mtimes simp) ($kdelta l1 (cons smlist l))
1977
(list (cons (caar e) '(simp)) (cons smlist l) (ncons smlist))
1978
(list '(mexpt simp) (meval (list 'mfactorial nn)) -1)) t)
1987
(cond ((atom e) (matcherr)))
1988
(cond ((atom (car e)) (matcherr)))
1989
(cond ((not (or (eq (caar e) '$levi_civita) (eq (caar e) '%levi_civita))) (matcherr)))
1990
(cond ((not (alike1 '((mlist simp)) (setq l1 ($covi e)))) (matcherr)))
1991
(cond ((not ($listp (setq l2 ($conti e)))) (matcherr)))
1992
(cond ((cdddr e) (matcherr)))
1993
(setq nn ($length l2))
1995
(do ((i nn (1- i))) ((< i 1)) (setq l (cons ($idummy) l) n $icounter))
1996
(return (values (list '(mtimes simp) ($kdelta (cons smlist l) l2)
1997
(list (cons (caar e) '(simp)) (ncons smlist) (cons smlist l))
1998
(list '(mexpt simp) (meval (list 'mfactorial nn)) -1)) t)
2004
(add2lnc '$lc_l $rules)
2005
(add2lnc '$lc_u $rules)
995
(DECLARE-TOP (SPECIAL E EMPTY $FLIPFLAG))
997
(SETQ $FLIPFLAG NIL EMPTY '((MLIST SIMP) ((MLIST SIMP)) ((MLIST SIMP))))
1001
((NUMBERP (CAR L)) (NONUMBER (CDR L)))
1003
(T (CONS (CAR L) (NONUMBER (CDR L))))
2007
(declare-top (special e empty $flipflag))
2009
(setq $flipflag nil empty '((mlist simp) ((mlist simp)) ((mlist simp))))
2013
((numberp (car l)) (nonumber (cdr l)))
2015
(t (cons (car l) (nonumber (cdr l))))
1007
(DEFUN REMOVEINDEX (E L)
1008
(COND ((NULL L) NIL)
1010
(COND ((EQ E (CAR L)) (CDR L))
1011
(T (CONS (CAR L) (REMOVEINDEX E (CDR L))))
2019
(defun removeindex (e l)
2020
(cond ((null l) nil)
2022
(cond ((eq e (car l)) (cdr l))
2023
(t (cons (car l) (removeindex e (cdr l))))
1013
(T (REMOVEINDEX (CDR E) (REMOVEINDEX (CAR E) L)))
2025
(t (removeindex (cdr e) (removeindex (car e) l)))
1018
(PROG (TOP BOTTOM X Y P Q R)
1019
(SETQ TOP NIL BOTTOM NIL)
1020
(COND ((RPOBJ E) (SETQ TOP (NONUMBER (CDADDR E)) BOTTOM (NONUMBER (APPEND (CDADR E) (CDDDR E)))))
1022
((MEMQ (CAAR E) '(MTIMES MNCTIMES MNCEXPT))
1024
(SETQ X (INDICES V) BOTTOM (APPEND BOTTOM (CADR X)) TOP (APPEND TOP (CAR X)))
1027
((MEMQ (CAAR E) '(MPLUS MEQUAL))
1028
(SETQ TOP (INDICES (CADR E)) BOTTOM (CADR TOP) TOP (CAR TOP))
1029
(SETQ P (INTERSECT TOP BOTTOM) Q (REMOVEINDEX P BOTTOM) P (REMOVEINDEX P TOP))
1030
(DOLIST (V (CDDR E))
1031
(SETQ X (INDICES V) Y (CADR X) X (CAR X))
1032
(SETQ R (INTERSECT X Y) X (REMOVEINDEX R X) Y (REMOVEINDEX R Y))
1033
(WHEN (NOT (AND (SAMELISTS X P) (SAMELISTS Y Q))) (MERROR "Improper indices in ~M" V))
1034
(SETQ TOP (UNION TOP R) BOTTOM (UNION BOTTOM R))
1037
((MEMQ (CAAR E) '($SUM %SUM))
1038
(SETQ TOP (LIST (CADDR E)) BOTTOM (LIST (CADDR E)))
1040
((MEMQ (CAAR E) '(%DERIVATIVE $DIFF))
1041
(DO ((I 1 (1+ I))) ((> I (COND ((CADDDR E) (CADDDR E)) (T 1))))
1042
(SETQ BOTTOM (CONS (CADDR E) BOTTOM)))
1044
;; (T (MERROR "Improper argument to INDICES: ~M" E))
2030
(prog (top bottom x y p q r)
2031
(setq top nil bottom nil)
2035
(setq top (nonumber (conti e))
2036
bottom (nonumber (append (covi e) (cdddr e))))
2040
(and (eq (caar e) 'mexpt) (eq (caddr e) -1))
2041
(setq x (indices (cadr e)) bottom (append bottom (car x))
2042
top (append top (cadr x)))
2045
(and (memq (caar e) '(%derivative $diff))
2046
(or (eq (length e) 3) (eq (cadddr e) 1)))
2047
(setq x (indices (cadr e)) bottom (append bottom (cadr x))
2048
top (append top (car x)))
2049
(setq x (indices (caddr e)) bottom (append bottom (car x))
2050
top (append top (cadr x)))
2053
(memq (caar e) '(mtimes mnctimes mncexpt))
2055
(setq x (indices v) bottom (append bottom (cadr x))
2056
top (append top (car x)))
2060
(memq (caar e) '(mplus mequal))
2061
(setq top (indices (cadr e)) bottom (cadr top) top (car top))
2062
(setq p (intersect top bottom) q (removeindex p bottom)
2063
p (removeindex p top))
2064
(dolist (v (cddr e))
2065
(setq x (indices v) y (cadr x) x (car x))
2066
(setq r (intersect x y) x (removeindex r x) y (removeindex r y))
2068
(not (and (samelists x p) (samelists y q)))
2069
(merror "Improper indices in ~M" v)
2071
(setq top (union top r) bottom (union bottom r))
2075
(memq (caar e) '($sum %sum))
2076
(setq top (list (caddr e)) bottom (list (caddr e)))
2079
(memq (caar e) '(%idiff $idiff))
2080
;;; This code would count derivative indices as covariant. However, it is
2081
;;; not needed. If the user wants to count derivative indices, those should
2082
;;; be part of the tensor expression; if the expression is undiff'd, there
2083
;;; must be a reason!
2085
;; ((f (cddr e) (cddr f)))
2089
;; ((> i (cond ((cadr f) (cadr f)) (t 1))))
2090
;; (setq bottom (cons (car f) bottom))
2093
(setq x (indices (cadr e)) bottom (append bottom (cadr x))
2094
top (append top (car x)))
2097
; (memq (caar e) '(%derivative $diff))
2098
; (setq x (indices (cadr e)) bottom (append bottom (cadr x))
2099
; top (append top (car x)))
2102
(return (list top bottom))
1046
(RETURN (LIST TOP BOTTOM))
1050
(DEFMFUN $INDICES (E)
1051
(PROG (TOP BOTTOM X)
1052
(SETQ TOP (INDICES E) BOTTOM (CADR TOP) TOP (CAR TOP) X (INTERSECT TOP BOTTOM))
1053
(SETQ TOP (REMOVEINDEX X TOP) BOTTOM (REMOVEINDEX X BOTTOM))
1054
(RETURN (CONS SMLIST (LIST (CONS SMLIST (APPEND TOP BOTTOM)) (CONS SMLIST X))))
1058
(DEFUN SAMELISTS (A B) ;"True" if A and B have the same distinct elements
1059
(AND (= (LENGTH A) (LENGTH B))
1064
(COND ((NULL L) (RETURN T))
1066
(T (RETURN NIL))))))
2106
(defmfun $indices (e)
2107
(prog (top bottom x)
2108
;; (setq top (indices e) bottom (cadr top) top (car top) x (intersect top bottom))
2109
(setq top (indices e) bottom (cadr top) top (car top) x (cond ($flipflag (intersect bottom top)) (t (intersect top bottom))))
2110
(setq top (removeindex x top) bottom (removeindex x bottom))
2111
(return (cons smlist (list (cons smlist (append top bottom)) (cons smlist x))))
2115
(defun samelists (a b) ;"True" if A and B have the same distinct elements
2116
(and (= (length a) (length b))
2121
(cond ((null l) (return t))
2123
(t (return nil))))))
1068
(DEFMFUN $FLUSH n ;Replaces the given (as arguments to FLUSH) indexed
2125
(defmfun $flush n ;Replaces the given (as arguments to FLUSH) indexed
1069
2126
(prog (l) ;objects by zero if they have no derivative indices.
1070
2127
(cond ((< n 2) (merror "FLUSH takes at least 2 arguments"))
1112
2169
(mapcar (function
1113
2170
(lambda (q) ($flushnd q name n)))
1114
2171
(cdr e))) e))))
1116
(DECLARE-TOP (FIXNUM INDEX N) (SPECIAL INDEX N DUMX))
1118
;(DEFMFUN $RENAME NARGS ((LAMBDA (INDEX) (RENAME (ARG 1)))
1119
; (COND ((= NARGS 1) 1) (T (ARG 2))))) ;Sets INDEX to 1 or 2nd argument of
1121
(DEFMFUN $RENAME NARGS
1122
(cond ((= NARGS 1) (setq INDEX 1)) (t (setq INDEX (arg 2)))) (rename (arg 1)))
1124
(DEFUN RENAME (E) ;Renames dummy indices consistently
1127
((OR (RPOBJ E) (EQ (CAAR E) 'MTIMES));If an indexed object or a product
1129
(SIMPTIMES (REORDER (COND (L (SUBLIS (itensor-CLEANUP L (SETQ N INDEX)) E))(T E))) 1 T))
1130
(CDADDR ($INDICES E)) ;Gets list of dummy indices
2173
(declare-top (fixnum index n) (special index n dumx))
2175
(defmfun $rename nargs
2176
(cond ((= nargs 1) (setq index 1)) (t (setq index (arg 2)))) (rename (arg 1)))
2178
(defun rename (e) ;Renames dummy indices consistently
2181
((or (rpobj e) (eq (caar e) 'mtimes););If an indexed object or a product
2182
(and (memq (caar e) '(%derivative $diff)) ; or a derivative expression
2183
(or (eq (length e) 3) (eq (cadddr e) 1)))
2186
(simptimes (reorder (cond (l (sublis (itensor-cleanup l (setq n index)) e))(t e))) 1 t))
2187
(cdaddr ($indices e)) ;Gets list of dummy indices
1132
(T ;Otherwise map $RENAME on each of the subparts e.g. a sum
1133
(MYSUBST0 (SIMPLIFYA (CONS (NCONS (CAAR E))
1134
(MAPCAR 'RENAME (CDR E)))
2189
(t ;Otherwise map $RENAME on each of the subparts e.g. a sum
2190
(mysubst0 (simplifya (cons (ncons (caar e))
2191
(mapcar 'rename (cdr e)))
1139
(DEFUN REORDER (E) ;Reorders contravariant, covariant, derivative indices
1140
(MYSUBST0 ;Example: F([A,B],[C,D],E,F)
1146
(NCONC (LIST (CAR X) ;($F SIMP)
1148
(COND ($ALLSYM (itensor-SORT (COPY (CDADR X))))
1149
(T (CDADR X)))) ;($A $B)
1152
(itensor-SORT (COPY (CDADDR X))))
1153
(T (CDADDR X))))) ;($C $D)
1154
(itensor-SORT (COPY (CDDDR X))))) ;($E $F)
1156
(COND ((EQ (CAAR E) 'MTIMES) (CDR E))
2196
(defun reorder (e) ;Reorders contravariant, covariant, derivative indices
2197
(mysubst0 ;Example: F([A,B],[C,D],E,F)
2203
(setq x ($renorm x))
2204
(nconc (list (car x) ;($f simp)
2206
(cond ($allsym (itensor-sort (copy (cdadr x))))
2207
(t (cdadr x)))) ;($a $b)
2210
(itensor-sort (copy (cdaddr x))))
2211
(t (cdaddr x))))) ;($c $d)
2212
(cond ($iframe_flag (cdddr x))
2213
(t (itensor-sort (copy (cdddr x))))))) ;($e $f)
2215
(cond ((eq (caar e) 'mtimes) (cdr e))
1160
(DEFUN itensor-CLEANUP (A N)((LAMBDA (DUMX)(CLEANUP1 A)) NIL)) ;Sets DUMX to NIL
2219
;;(defun itensor-cleanup (a n)((lambda (dumx)(cleanup1 a)) nil)) ;Sets DUMX to NIL
2220
(defun itensor-cleanup (a nn) (setq n nn dumx nil) (cleanup1 a))
1163
(AND A (SETQ DUMX (IMPLODE (NCONC (EXPLODEN $DUMMYX) ;Keep proper order of
1164
(EXPLODEN N))) N (1+ N)) ;indices
1165
(COND ((EQ DUMX (CAR A)) (CLEANUP1 (CDR A)))
1166
(T (CONS (CONS (CAR A) DUMX) (CLEANUP1 (CDR A)))))))
1167
;Make list of dotted pairs indicating substitutions i.e. ((A . #1) (B . #2))
1169
(DECLARE-TOP (NOTYPE N INDEX)(UNSPECIAL N DUMX INDEX))
1171
(DEFUN itensor-SORT (L) (COND ((CDR L) (SORT L 'LESS)) (T L)))
2223
(and a (setq dumx (implode (nconc (exploden $idummyx) ;Keep proper order of
2224
(exploden n))) n (1+ n)) ;indices
2225
(cond ((eq dumx (car a)) (cleanup1 (cdr a)))
2226
(t (cons (cons (car a) dumx) (cleanup1 (cdr a)))))))
2227
;Make list of dotted pairs indicating substitutions i.e. ((a . #1) (b . #2))
2229
(declare-top (notype n index)(unspecial n dumx index))
2231
(defun itensor-sort (l) (cond ((cdr l) (sort l 'less)) (t l)))
1172
2232
;Sort into ascending order
1174
(DEFMFUN $REMCOMPS (TENSOR)
1175
(ZL-REMPROP TENSOR 'EXPR) (ZL-REMPROP TENSOR 'CARRAYS)
1176
(ZL-REMPROP TENSOR 'TEXPRS) (ZL-REMPROP TENSOR 'INDEXED)
1177
(ZL-REMPROP TENSOR 'TSUBR) '$DONE)
2234
(defmfun $remcomps (tensor)
2235
(zl-remprop tensor 'expr) (zl-remprop tensor 'carrays)
2236
(zl-remprop tensor 'texprs) (zl-remprop tensor 'indexed)
2237
(zl-remprop tensor 'indexed) (zl-remprop tensor 'tsubr)
2238
(and (functionp tensor) (fmakunbound tensor))
1179
(DEFMFUN $INDEXED (TENSOR)
1181
(AND (ZL-GET TENSOR 'EXPR)
2241
(defmfun $indexed_tensor (tensor)
2243
(and (zl-get tensor 'expr)
1182
2244
(merror "~M has expr" tensor))
1184
(AND (SETQ FP (ZL-GET TENSOR 'SUBR))
1185
(PROGN (SETQ NEW (GENSYM))(PUTPROP NEW FP 'SUBR)
1186
(ZL-REMPROP TENSOR 'SUBR)(PUTPROP TENSOR NEW 'TSUBR)))
1187
(PUTPROP TENSOR T 'INDEXED)
1188
(PUTPROP TENSOR (SUBST TENSOR 'G '(LAMBDA N (TENSOREVAL (QUOTE G)(LISTIFY N)))) 'EXPR)
1189
(eval (subst tensor 'g (quote (defmfun g n (tensoreval 'g (listify n))))))
1194
(AND L (FIXP (CAR L)) (OR (NULL (CDR L)) (ALLFIXED (CDR L)))))
1196
;;(DEFUN TENSOREVAL (TENSOR INDXS)
1197
;; ((LAMBDA (DER CON)
1198
;; (AND (CDR INDXS) (SETQ CON (CDADR INDXS) DER (CDDR INDXS)))
1199
;; (SETQ TENSOR (SELECT TENSOR (CDAR INDXS) CON))
1200
;; (COND (DER (APPLY '$DIFF (CONS TENSOR (PUTINONES DER))))
1201
;; (T TENSOR))) NIL NIL))
1202
(DEFUN TENSOREVAL (TENSOR INDXS)
1204
(AND (CDR INDXS) (SETQ CON (CDADR INDXS) DER (CDDR INDXS)))
1205
(SETQ TENSOR (SELECT TENSOR (CDAR INDXS) CON DER))
1208
;;(DEFMFUN $COMPONENTS (TENSOR COMP)
1209
;; ((LAMBDA (LEN1 LEN2 NAME PROP)
1210
;; (COND ((OR (NOT (RPOBJ TENSOR))(CDDDR TENSOR))
1211
;; (merror "Improper 1st arg to COMPONENTS: ~M"
1214
;; (SETQ LEN1 (LENGTH (CDADR TENSOR)) LEN2 (LENGTH (CDADDR TENSOR)))
1215
;; (AND (NOT (ATOM COMP))(EQ (CAAR COMP) '$MATRIX)
1216
;; (COND ((= (f+ LEN1 LEN2) 2)(SETQ NAME (GENSYM))
1217
;; (SET NAME COMP)(SETQ COMP NAME))
1219
;; (merror "Needs two indices for COMPONENTS from matrix:~%~M"
1221
;; (COND ((AND (EQ (ML-TYPEP COMP) 'SYMBOL) (> (f+ LEN1 LEN2) 0))
1222
;; (SETQ PROP 'CARRAYS))
1223
;; ((SAMELISTS (SETQ NAME (APPEND (CDADR TENSOR) (CDADDR TENSOR)))
1224
;; (CDADR ($INDICES COMP)))
1225
;; (SETQ PROP 'TEXPRS COMP (CONS COMP NAME)))
1226
;; (T (merror "Args to COMPONENTS do not have the same free indices")))
1227
;; (SETQ TENSOR (CAAR TENSOR) LEN1 (CONS LEN1 LEN2))
1228
;; (COND ((AND (SETQ NAME (ZL-GET TENSOR PROP))
1229
;; (SETQ LEN2 (ZL-ASSOC LEN1 NAME))) (RPLACD LEN2 COMP))
1230
;; (T (PUTPROP TENSOR (CONS (CONS LEN1 COMP) NAME) PROP)))
1231
;; (OR (ZL-GET TENSOR 'INDEXED) ($INDEXED TENSOR))
1232
;; '$DONE) NIL NIL NIL NIL))
1233
(DEFMFUN $COMPONENTS (TENSOR COMP)
1234
((LAMBDA (LEN1 LEN2 LEN3 NAME PROP)
1235
(COND ((NOT (RPOBJ TENSOR))
2246
(and (setq fp (zl-get tensor 'subr))
2247
(progn (setq new (gensym))(putprop new fp 'subr)
2248
(zl-remprop tensor 'subr)(putprop tensor new 'tsubr)))
2249
(putprop tensor t 'indexed)
2250
(putprop tensor (subst tensor 'g '(lambda nn (tensoreval (quote g)(listify nn)))) 'expr)
2251
(eval (subst tensor 'g (quote (defmfun g nn (tensoreval 'g (listify nn))))))
2256
(and l (fixp (car l)) (or (null (cdr l)) (allfixed (cdr l)))))
2258
(defun tensoreval (tensor indxs)
2260
(and (cdr indxs) (setq con (cdadr indxs) der (cddr indxs)))
2261
(setq tensor (select tensor (cdar indxs) con der))
2264
(defmfun $components (tensor comp)
2265
((lambda (len1 len2 len3 name prop)
2266
(cond ((not (rpobj tensor))
1236
2267
(merror "Improper 1st arg to COMPONENTS: ~M"
1239
(SETQ LEN1 (LENGTH (CDADR TENSOR)) LEN2 (LENGTH (CDADDR TENSOR)) LEN3 (LENGTH (CDDDR TENSOR)))
1240
(AND (NOT (ATOM COMP))(EQ (CAAR COMP) '$MATRIX)
1241
(COND ((= (f+ (f+ LEN1 LEN2) LEN3) 2)(SETQ NAME (GENSYM))
1242
(SET NAME COMP)(SETQ COMP NAME))
2270
; (setq len1 (length (cdadr tensor)) len2 (length (cdaddr tensor)) len3 (length (cdddr tensor)))
2271
(setq len1 (length (covi tensor)) len2 (length (conti tensor)) len3 (length (deri tensor)))
2272
(and (not (atom comp))(eq (caar comp) '$matrix)
2273
(cond ((= (f+ (f+ len1 len2) len3) 2)(setq name (gensym))
2274
(set name comp)(setq comp name))
1244
2276
(merror "Needs two indices for COMPONENTS from matrix:~%~M"
1246
(COND ((AND (EQ (ML-TYPEP COMP) 'SYMBOL) (> (f+ (f+ LEN1 LEN2) LEN3) 0))
1247
(SETQ PROP 'CARRAYS))
1248
((SAMELISTS (SETQ NAME (APPEND (CDADR TENSOR) (CDADDR TENSOR) (CDDDR TENSOR)))
1249
(CDADR ($INDICES COMP)))
1250
(SETQ PROP 'TEXPRS COMP (CONS COMP NAME)))
1251
(T (merror "Args to COMPONENTS do not have the same free indices")))
1252
(SETQ TENSOR (CAAR TENSOR) LEN1 (LIST LEN1 LEN2 LEN3))
1253
(COND ((AND (SETQ NAME (ZL-GET TENSOR PROP))
1254
(SETQ LEN2 (ZL-ASSOC LEN1 NAME))) (RPLACD LEN2 COMP))
1255
(T (PUTPROP TENSOR (CONS (CONS LEN1 COMP) NAME) PROP)))
1256
(OR (ZL-GET TENSOR 'INDEXED) ($INDEXED TENSOR))
1257
'$DONE) NIL NIL NIL NIL NIL))
1259
;;(DEFUN SELECT (TENSOR L1 L2)
1260
;; ((LAMBDA (PROP SUBS INDEX)
1261
;; (COND ((AND (ALLFIXED SUBS) (SETQ PROP (ZL-GET TENSOR 'CARRAYS))
1262
;; (SETQ PROP (ZL-ASSOC INDEX PROP)))
1263
;; (COND ((ALIKE1 (SETQ PROP (CONS (LIST (CDR PROP) 'ARRAY) SUBS))
1264
;; (SETQ SUBS (MEVAL PROP))) 0)
1266
;; ((SETQ PROP (ZL-ASSOC INDEX (ZL-GET TENSOR 'TEXPRS)))
1267
;; (SUBLIS (MAPCAR (FUNCTION CONS)(CDDR PROP) SUBS) (CADR PROP)))
1268
;; ((SETQ PROP (ZL-GET TENSOR 'TSUBR))
1269
;; (APPLY PROP (LIST (CONS SMLIST L1)(CONS SMLIST L2))))
1270
;; (T (LIST (LIST TENSOR 'SIMP)(CONS SMLIST L1)(CONS SMLIST L2)))))
1271
;; NIL (APPEND L1 L2)(CONS (LENGTH L1)(LENGTH L2))))
1273
;;vtt: inconstant was an attempt to remove constant indices, but it really doesn't work out.
1274
;;(DEFUN INCONSTANT (L)
1277
;; (($CONSTANTP (CAR L)) (AND (NOT (EQ NIL (CDR L))) (INCONSTANT (CDR L))))
1278
;; (T (CONS (CAR L) (AND (NOT (EQ NIL (CDR L))) (INCONSTANT (CDR L)))))
1282
(DEFUN SELECT (TENSOR L1 L2 L3)
1283
((LAMBDA (PROP SUBS INDEX)
1284
(COND ((AND (ALLFIXED SUBS) (SETQ PROP (ZL-GET TENSOR 'CARRAYS))
1285
(SETQ PROP (ZL-ASSOC INDEX PROP)))
1286
(COND ((ALIKE1 (SETQ PROP (CONS (LIST (CDR PROP) 'ARRAY) SUBS))
1287
(SETQ SUBS (MEVAL PROP))) 0)
1289
((SETQ PROP (ZL-ASSOC INDEX (ZL-GET TENSOR 'TEXPRS)))
1290
;;;VTT (SUBLIS (MAPCAR (FUNCTION CONS)(CDDR PROP) SUBS) (CADR PROP)))
1291
;;; (SUBLIS (MAPCAR (FUNCTION CONS)(CDDR PROP) SUBS) ($RENAME (CADR PROP) (COND ((BOUNDP 'N) N) (T 1)))))
1292
(SUBLIS (MAPCAR (FUNCTION CONS)(CDDR PROP) SUBS) (CAR (CONS ($RENAME (CADR PROP) (1+ $COUNTER)) (SETQ $COUNTER (1- (COND ((BOUNDP 'N) N) (T 1))))))))
1293
((SETQ PROP (ZL-GET TENSOR 'TSUBR))
1294
;; (APPLY PROP (LIST (CONS SMLIST (INCONSTANT L1))(CONS SMLIST (INCONSTANT L2))(CONS SMLIST L3))))
1295
;; ((NOT (EQ L3 NIL)) (APPLY '$DIFF (SELECT TENSOR (INCONSTANT L1) (INCONSTANT L2) (CDR L3)) (LIST (CAR L3))))
1296
;; (T (APPEND (LIST (LIST TENSOR 'SIMP)(CONS SMLIST (INCONSTANT L1))(CONS SMLIST (INCONSTANT L2))) L3))))
1297
;; NIL (APPEND (INCONSTANT L1) (INCONSTANT L2) L3)(LIST (LENGTH (INCONSTANT L1))(LENGTH (INCONSTANT L2))(LENGTH L3))))
1298
(APPLY PROP (LIST (CONS SMLIST L1)(CONS SMLIST L2)(CONS SMLIST L3))))
1299
((NOT (EQ L3 NIL)) (APPLY '$DIFF (SELECT TENSOR L1 L2 (CDR L3)) (LIST (CAR L3))))
1300
(T (APPEND (LIST (LIST TENSOR 'SIMP)(CONS SMLIST L1)(CONS SMLIST L2)) L3))))
1301
NIL (APPEND L1 L2 L3)(LIST (LENGTH L1)(LENGTH L2)(LENGTH L3))))
1304
(DEFMFUN $ENTERTENSOR nargs
2278
(cond ((and (eq (ml-typep comp) 'symbol) (> (f+ (f+ len1 len2) len3) 0))
2279
(setq prop 'carrays))
2280
; ((samelists (setq name (append (cdadr tensor) (cdaddr tensor) (cdddr tensor)))
2281
((samelists (setq name (append (covi tensor) (conti tensor) (deri tensor)))
2282
(cdadr ($indices comp)))
2283
(setq prop 'texprs comp (cons comp name)))
2284
(t (merror "Args to COMPONENTS do not have the same free indices")))
2285
(setq tensor (caar tensor) len1 (list len1 len2 len3))
2286
(cond ((and (setq name (zl-get tensor prop))
2287
(setq len2 (zl-assoc len1 name))) (rplacd len2 comp))
2288
(t (putprop tensor (cons (cons len1 comp) name) prop)))
2289
(or (zl-get tensor 'indexed) ($indexed_tensor tensor))
2290
'$done) nil nil nil nil nil))
2292
(defun select (tensor l1 l2 l3)
2295
(setq l2 (append (minusi l1) l2) l1 (plusi l1))
2304
(setq prop (zl-get tensor 'carrays))
2305
(setq prop (zl-assoc idx prop))
2310
(setq prop (cons (list (cdr prop) 'array) subs))
2311
(setq subs (meval prop))
2319
(setq prop (zl-assoc idx (zl-get tensor 'texprs)))
2321
(mapcar (function cons)(cddr prop) subs)
2322
($rename (cadr prop) (cond ((boundp 'n) n) (t 1)))
2326
(setq prop (zl-get tensor 'tsubr))
2329
(list (cons smlist l1) (cons smlist l2) (cons smlist l3))
2334
(apply '$idiff (select tensor l1 l2 (cdr l3)) (list (car l3)))
2339
(list (list tensor 'simp) (cons smlist l1) (cons smlist l2))
2345
nil (append l1 l2 l3) (list (length l1)(length l2)(length l3))
2352
(defmfun $entertensor nargs
1305
2353
(prog (fun contr cov deriv)
1307
(merror "ENTERTENSOR takes 0 or 1 arguments only"))
1309
(mtell "Enter tensor name: ")
1310
(setq fun (meval (retrieve nil nil))))
1311
((setq fun (arg 1))))
2357
(merror "ENTERTENSOR takes 0 or 1 arguments only")
2361
(mtell "Enter tensor name: ")
2362
(setq fun (meval (retrieve nil nil)))
2364
((setq fun (arg 1)))
1312
2366
(mtell "Enter a list of the covariant indices: ")
1313
2367
(setq cov (checkindex (meval (retrieve nil nil)) fun))
1314
2368
(cond ((atom cov) (setq cov (cons smlist (ncons cov)))))
added
removed
share/tensor/itensor.lisp
