Development of a matrix calculus package



Thank you, I'll take a look. I've used the itensor package in the past, but
it seemed quite limited as far multiplicative inverses, transposes, and it
appears to assume everything is a general tensor. But maybe that demo will
show more than I've come to expect of the itensor package.


On Thu, Mar 22, 2012 at 8:00 PM, Viktor T. Toth <vttoth at vttoth.com> wrote:

> 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.
>
>