I have an implementation for implicitly stored Kronecker products laying on my disc from 3 years ago. You can check it out here https://gitlab.com/mfeuerle/mastersthesis (ignore the tamg module, that is a train wreck, the kron module is the topic).
I spent quite some time back then implementing, documenting, and writing pytests for it, so it is already in quite a good shape, although not touched for three years. I recently figured that this might be a good addition to scipy it self as i think more people could find that useful. I could get it ready for a pull-request within a few days of work. But before doing so, i wanted to ask if there is actually any interest in adding these features.
So whats this all about:
When taking the the Kronecker of an n×m matrix A and an p×q matrix B, the resulting Kronecker product A⊗B is of dimension np×mq. In the full matrix we basically store the matrix B not one but nm times, just scaled differently. That is a lot of redundancy and a heavy burden on memory consumption, becoming infeasible quite fast. So the idea of the implicit Kronecker module is to store A⊗B not as a full matrix but implicitly by just storing the matrices A and B. Assuming n=m=p=q and dense matrices, the memory is reduced from O(n^4) to O(n^2).
But this i not only beneficial memory wise, as most operations on A⊗B can be reduced to separate operations on A and B, e.g.
(A⊗B)x = BXA^T, whereXis a reshape of the vectorxinto a matrix(A⊗B)^-1 = A^-1⊗B^-1det(A⊗B) = det(A)^n det(b)^p(ifn=mandp=q)- …
So for many operations on Kronecker products, the full matrix is not need and the operations can be performed much faster using A and B separately, e.g. the matrix-vector multiplication is reduced from O(n^4) to O(n^3).
In the implicit Kronecker module i implemented this approach of representing Kronecker products of arbitrary order A1⊗A2⊗....
I also implemented an implicit storage of sums A1+A2+... where A1,A2,... are each Kronecker products and for block-wise matrices kron.block([[A11,A12],[A21,A22]]) where Aij are each Kronecker products or (sums of Kronecker products).
Although for sums and blocks of Kronecker products, many of the neat properties get lost, the efficient storage and the fast matrix-vector multiplication stays. So all these matrices can be used very efficiently for iterative solvers like cg, gmres, etc. using the scipy.sparse.linalg.LinearOperator class.
In my implementation i focused on the matrix-vector implementation, as this works for all of these matrices. So that is what i could provide. But having the format at hand, i think it would be quite easy to just add a typ check inside the scipy.linalg / scipy.sparse.linalg functions solve,inv,det,eig,... to also provide efficient implementations for the other properties of pure Kronecker products (no sums or blocks).
My Question:
Do you think that this is something useful to add to scipy? In that case i could try to get it ready for a pull-request.