Computing the rank of a symbolic matrix in Julia

No worries: I had suspected that \mathcal{C}(s) and \mathcal{O}(s) for the PBH test matrices were misnomers, if you will, but the point was clear and of use. Grateful for the SNF survey and suggestion, but I cannot deviate from symbolic forms, because the scope is to advance eigenvalue conditions for the VAR representation of x_{t} in y_{t} under all possible cases. Can I please ask you to try this code in order to check how long the machine may employ to deliver the eigenvalues? Thankful.

using SymPy, LinearAlgebra
A=[Sym[“a$i$j” for i in 1:3, j in 1:3] zeros(3, 1); zeros(1, 4)]; B=Sym[“b$i” for i in 1:4];
M=[1 0 0 0]; C=MA; D=MB;
Obs=[C; CA; CA^2; CA^3]; ro=Obs.rank();
T=[Obs[1:3, 1:4]; 0 0 0 1];
Ad=inv(T)AT; Bd = inv(T)B; Cd = CT;
Am=Ad[1:3, 1:3]; Bm=Bd[1:3]; Cm=Cd[1:1, 1:3];
Fm=Am-Bm
inv(D)*Cm; eigvals(Fm)

To j_verzani: did you manage to run the code, by any chance? Cheers.

Can you reformat? There is some missing operator with

C=M A; D=M B;

But it likely doesn’t matter. SymPy seems to be fine with eigenvalues of arbitrary 3x3 matrices, but 4x4 are a problem.

Sorry about that, it was supposed to show the multiplication signs. In pasting, the quotation marks were moreover transmuted into smart quotes; B was also defined inefficiently. The correct version be should today’s, which I presented to the other user; at any rate, my machine managed to return the eigenvalues after an odd hour, in line with your technical insight. Thanks for the symbolic solution, anyway.

Eigenvalues (asymmetric)|690x201

How can I get SymPy to compute controllability matrix \mathcal{C}=[B \; \cdots \; A^{n_{x}-1}B] through an element loop? This command returns a 1\times3 row vector, as opposed to a 3\times3 matrix.

using SymPy, LinearAlgebra, Statistics, Compat
A=[Sym[“a$i$j” for i in 1:2, j in 1:2] zeros(2, 1); zeros(1, 3)]; B=Sym[“b$i” for i in 1:3];
M=[1 0 0]; C=MA; D=MB; nx=size(A, 1)
Con=[A^(i)*B for i in 0:nx-1]

The element wise (dot) alternative by contrast returns the wrong matrix. Thanks.

As used here, SymPy is just using Julia’s AbstractArray type. So this is a linear algebra question. Maybe, reduce(hcat, ([A^i for i in 0:(n-1)]..., B)) is what you want. (BTW, the code you had earlier seems to have an issue with your “D” matrix, so I couldn’t test.)