[ANN] DescriptorSystems.jl: A package for manipulation of generalized state-space (descriptor) system representations

Package Announcements
,
18

I thank you for your comments on MatrixEquations.jl, which motivated me to try to further enhance the performance of Lyapunov solvers.
I registered today a new version of MatrixEquations, where the updated solvers for Lyapunov-type matrix equations exhibit a substantial increase of execution speed (up to two times),
while the memory allocation behaviour has been drastically improved (up to 50 times reduction of the allocated memory). Further enhancements are certainly possible (e.g., by storing
intermediary results in working arrays), but the expected increase in speed (probably less than 10%) does not justify, in my opinion, the use of additional storage.
Old Fortran programmers, like me, would never waste unnecessarily memory for insignificant increase in speed!
Also, the usage of solvers for small order equations based on static arrays led in my experiments to a marginal increase of speed (and involves certainly extra storage allocations).
An alternative to achieve the limiting speed would be to integrate the solvers available in the SLICOT Fortran library (now freely available, under a BSD license).

Here are some performance figures for the improved solvers, run on 500x500 examples:

julia> n = 500; a = rand(n,n); e = rand(n,n); q = rand(n,n); q = q'*q;

julia> @btime ControlMatrixEquations.lyapd(a,q);
  225.434 ms (26 allocations: 13.51 MiB)

julia> @btime MatrixEquations.lyapd(a,q);
  241.155 ms (27 allocations: 11.60 MiB)

julia> @btime ControlMatrixEquations.lyapd(a,q,e);
  810.124 ms (40 allocations: 24.98 MiB)

julia> @btime MatrixEquations.lyapd(a,e,q);
  820.314 ms (37 allocations: 19.26 MiB)

I would like to add that the performance of the Matlab solver dlyap, which is based on the Fortran codes from SLICOT, is probably the limit of achievable speed
(a nice challenge for further speedup of your code).
Here are, for comparison purposes, the time measures computed with Matlab 2020b:


>> a = rand(500); e = rand(500); q = rand(500); q = q'*q;
>> tic; x = dlyap(a,q); toc
Elapsed time is 0.115986 seconds.
>> tic; x = dlyap(a,e,q); toc
Elapsed time is 0.215429 seconds.
5 Likes