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]
K1*K2
# 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
G*(K2^2)
# 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.