Accessing MatrixSymbol

[baggepinnen]

I can't seem to figure out how to make use of the sympy MatrixSymbol type, described at http://docs.sympy.org/latest/modules/matrices/expressions.html I would like to manipulate matrix expressions on a block matrix level, not have the expressions be manipulated on the level of scalar matrix entries.

To be more specific, how do I in SymPy.jl create a matrix with specified dimensions, without having it be an array of scalars, like (in python)

>>> X = MatrixSymbol('X', 3, 3)
>>> Matrix(X)
Matrix([
[X[0, 0], X[0, 1], X[0, 2]],
[X[1, 0], X[1, 1], X[1, 2]],
[X[2, 0], X[2, 1], X[2, 2]]])

[mzaffalon]

Interesting question. What is your use case?

[baggepinnen]

I want to manipulate the expression for the KL divergence between to multivariate gaussian distributions that each are the joint distribution of two groups of variables. I explicitly want to see how the expression simplifies if one of the subgroups of variables are the same in the two distributions. If I do this with scalars, sympy assumes that multiplication is commutative and the results won't hold for matrices.

Sent from mobile device

On Dec 3, 2017 22:31, "Michele Zaffalon" notifications@github.com wrote:

Interesting question. What is your use case?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/JuliaPy/SymPy.jl/issues/200#issuecomment-348816321, or mute the thread https://github.com/notifications/unsubscribe-auth/ADnx85cDmVvXu-Yr__K_wIyr2zdpqWqgks5s8xNEgaJpZM4QzrMK .

[mzaffalon]

I do not know if I understand it, but why can you not define an array of Any where each entry is one SymPy matrix and define some appropriate operations on the Array{Any}, specifically matrix-matrix multiplications? Unary operations can be carried out with the dotted version.

[baggepinnen]

I'm not sure if defining a Base method on a Base type is a good idea, it can alter the behavior in other packages? In either case, if I define the block matrix multiplication, wouldn't SymPy go through and carry out the matrix multiplications between the individual blocks on a scalar level? What I'm looking for is essentially a type that has a size but individual elements are not "defined". It would behave according to

A = SomeMatrixType(3,3, invertible=true)
B = SomeMatrixType(3,3, invertible=true)
C = SomeMatrixType(2,2)
A*inv(A) == I  # If A invertible 
A*B*A != A^2*B # These expressions are equal if A and B are scalar, but not in general equal if A,B <: Matrix
A*C == Error   # Size does not match

The MatrixSymbol object in sympy seems to work like this

>>> from sympy import MatrixSymbol, Identity
>>> A = MatrixSymbol('A', 3, 4) # A 3 by 4 Matrix
>>> B = MatrixSymbol('B', 4, 3) # A 4 by 3 Matrix
>>> A.shape
(3, 4)
>>> 2*A*B + Identity(3)
I + 2*A*B

[jverzani]

Here is a painful way to use this class:

X = sympy"MatrixSymbol"("X",3,3)
Y = sympy"MatrixSymbol"("Y",3,3)

So far so good. But to use them is cumbersome:

PyCall.py"$(X[:T]()) * $Y * $X"

The latter should be simply X'*Y*X but to do that requires some changes in SymPy:

• create a type in Julia to hold these (say SymMatrixSymbol).
• create a constructor for this type (to avoid sympy"MatrixSymbol"
• Set up a pytype_mapping, as in

pytype_mapping(sympy[:MatrixSymbol], SymMatrixSymbol)

• define the operations for this class

With this, matrices of SymMatrixSymbol values should be doable as in SymPy

I have a branch (https://github.com/jverzani/SymPy.jl/blob/v0.6-7/src/permutations.jl) where I do this for permutations. I can integrate this in at the same time.

[jverzani]

Well, I tried, but didn't get this done. Following the above comment leads to a situation where simple things work fine, but ultimately fail as the two types (SymMatrixSymbol and Sym) don't play nicely together with the generic linear algebra functions.

For your use case, you might just create your variables as scalars and use commutative=false when defining them.