I am not sure whether or not it addresses your needs, but have you looked at
the tensor calculus capabilities that are already in the itensor package? In
particular, look at the einhil demo, which shows how itensor can be used to
derive the Einstein field equations from the Einstein-Hilbert Lagrangian.
(The functionality is not limited to general relativity though.)
Viktor
-----Original Message-----
From: maxima-bounces at math.utexas.edu [mailto:maxima-bounces at math.utexas.edu]
On Behalf Of Mike Valenzuela
Sent: Thursday, March 22, 2012 6:19 PM
To: maxima
Subject: Development of a matrix calculus package
Hello everyone,
a while back I asked about writing rules for derivatives of matrices. There
were many, many rules. But by replacing some fundamental identities (dX/dX
--> 4th dimensional object), I think I have found a general approach which
appears to rederive about 80% of all the rules. The draw back is I need to
be able to replace some fundamental identities in Maxima's default calculus
package. If I can just rewrite the fundamental rules, then the rest of
calculus applies as normal. I know this is possible since the itensor
package does something similar.
Way too many details below, but it should help clarify what I mean by
replacing Maxima's default fundamental identities:
TASK 1: replacing fundamental identities. I was wondering how to make the
following rules, using as much of Maxima's underlying functionality as
possible.
(a) diff( X, X[i][j] ), where X is a general matrix --> J[i][j], the single
entry matrix
I may prefer the notation X[i,j], but notation is just a side note.
(b) diff( X, X[i][j] ), where X is a structured matrix --> Struct[i][j]
(c) diff( X, X[i][j] ), where X is a symmetric matrix --> Sym[i][j]
expand( Sym[i][j] ) --> delta(i,j) * J[i][i] + (1-delta(i,j) ) * ( J[i][j] +
J[j][i] )
(d) diff( f(U(X)), X[i][j] ), where X is a matrix --> trace( transpose(
diff( f(U), U ) ) . diff(U, X[i][j] ) )
(e) diff( g(X) . f(X) , X[i][j] ), where X is a matrix --> diff( g(X),
X[i][j] ) . f(X) + g(X) . diff( f(X), X[i][j] )
For instance, (a) and (e) ought to be able to process:
diff( X^3, X[i][j] ) --> diff( X.(X.X), X[i][j] ) --> J[i][j] . (X . X) + X
. diff( X . X, X[i][j]) --> J[i][j] . X . X + X . J[i][j] . X + X . X .
J[i][j]
(f) diff( v, v[i] ), where v is a general vector --> e[i]
(g) diff( transpose(f(X)), X[i][j] ) --> transpose( diff( f(X), X[i][j] ) )
(h) diff( exp( X . t ), t), where X is a matrix, t is scalar --> X . exp(X .
t) == exp(X . t) . X
The output matrices commute, and Maxima knows this!
(i) diff( trace( X ), X[i][j] ) --> I[i][j]
Which reduces to the identity matrix when indices are removed
(j) diff( trace( F(X) ), X ), where X is a general matrix --> transpose(
f(X) )
Where f(.) is the scalar derivative of F(.).
And several other fundamental rules like these. I wish to include Hadamard
and Kronecker products as well. I would appreciate help with extensions of
these rules for symmetric, hermitian, diagonal, toeplitz, etc.
---
TASK2: Attempting to remove indices
(2a) diff( transpose( v ) . X . v, v ) --> contractIndices( diff(
transpose( v ) . X . v, v[i] ), [i] ) -->
contractIndices( transpose(e[i]) . X . v + transpose( v ) . X . e[i], i )
-->
(X + transpose(X)) . v
(2b) diff( trace(X . X), X ) --> contractIndices( diff( trace(X . X),
X[i][j] ), [i,j] ) --> we use the scalar derivative
contractIndices( transpose( 1 * X[i][j] + X[i][j] * 1 ), [i,j] ) -->
2 transpose( X )
I don't know a good general procedure for index removal, but I would
appreciate anyone willing to help with the contract indices method.
Thanks in advance for any help/pointers.