Apply Matrix Order of Operations on Scalars

Is there a way to make variables defined by Symbolics (Or any symbolic module I suppose) to be treated with the algebraic rules and order of operations for matrices, while still displaying as scalars?

I am using Symbolics.jl to double check some math I am doing, mostly LFTs. The hand calculations use matrices of matrices, so it seemed worth using the same notation in Julia. My general use case is

using Symbolics, LinearAlgebra

@variables G W1 W2 W3 W4 K1 K2 # All variables are matrices with appropriate dimensions assumed
# Constructing generalized plant P column by column
# 0 and 1 are taking place of zero and identity matrices.  
C1 = [0, 0, 0, 0, 1, 0]
C2 = [0, W2*G, W3*G, W4, 0, G]
C3 = [0, 0, W3, 0, 0, 1]
C4 = [W1, W2*G, W3*G, W4, 0, G]
P = [C1 C2 C3 C4]

# Partitioning
P11 = P[1:4, 1:3]
P12 = P[1:4, 4]
P21 = P[5:6, 1:3]
P22 = P[5:6, 4]
K = [K2*K1 (-K2)]

# Calculating N via LFT
N = P11 + P12 * K * inv(convert.(Float64, I(2)) - P22 * K)*P21

I’m checking that what I derived by hand matches N, or various steps along the way like P22 * K.

This gives me what seems to be the correct result if everything was scalars. But for matrices, some of the order of operations may get flipped. For instance K[1] has it’s elements flipped.

julia> K[1]
# Assumes commutative. Would like to be K2*K1.

And I get that this is correct and proper for scalars, so I shouldn’t have anything to complain about.

Inversion creates dividing terms, such as

julia> inv(I(2) - P22 * K)
2×2 Matrix{Num}:
                1                  0
 (G*K1*K2) / (true + G*K2)  1 / (true + G*K2)
# Would like division to be (1 + G*K2)^-1 on left.

Further down the line, terms get collected in ways which make sense for scalars, but not for my weird use case.

julia> K2 * G * K2
# Again, assumes commutative to collect terms. Would like to be K2*G*K2.

If I re-write the script using Symbolics arrays (@variables G[1:2, 1:2], W1[1:2,1:2], etc...), everything appears correct, and the numerical results appear to match. But I’m trying to check some intermediate result which vanishes when I do this.

I imagine there is a Julia-esque way to do this which overloads the operators. I guess this could be as general as I need some scalar type to be use the matrix order of operations, not just Symbolics types.

I do not have an answer but I can offer a concise reformulation of your question: what you are looking for is capability to do noncommutative algebra in Julia.

Recognising the control theoretic motivation in your post, I recall that a package called NCAlgebra offering the needed functionality was developed by William Helton for Mathematica. I do not know if anything similar is available in Julia.

1 Like