[ANN] PeriodicMatrices.jl - Handling of periodic time-varying matrices

Package Announcements
, , , ,
1

PeriodicMatrices.jl provides the basic tools to handle periodic time-varying matrices. This package is the basis for developing the PeriodicMatrixEquations.jl and PeriodicSystems.jl packages.

For a real periodic matrix A(t) with period T, the dependence of the time variable t can be either continuous or discrete.

A continuous-time periodic matrix can be specified in one of the following forms:

  • periodic matrix function, with A(t) a matrix function of the real variable t āˆˆ [0, T);

  • periodic symbolic matrix, with A(t) a symbolic matrix as defined in the Symbolics.jl package depending on the (symbolic) real variable t āˆˆ [0, T);

  • harmonic matrix series, with A(t) defined as

                       p
     A(t) = A_0 +  āˆ‘ ( Ac_i*cos(iĻ‰t)+As_i*sin(iĻ‰t) ) ,
                  i=1

    where Ļ‰ = 2Ļ€/T and A_0, Ac_i, As_i for i = 1,..., p are real matrices;

  • periodic matrix time series with constant dimensions on a uniform time grid;

  • periodic matrix time series with constant dimensions on a non-uniform time grid;

  • Fourier matrix series approximation, with A(t) a Fourier series representation (similar to the harmonic matrix series representation) as defined in the ApproxFun.jl package.

A discrete-time periodic matrix can be specified in the following forms:

  • periodic matrix time series with time-varying dimensions on a uniform time grid;

  • periodic matrix time series with time-varying dimensions on a non-uniform time grid;

  • periodic matrix time series with constant dimensions on a uniform time grid;

  • periodic matrix time series with constant dimensions on a non-uniform time grid.

All possible conversions between the above representations are supported.
The usage of the FFTW, Symbolics.jl, ApproxFun.jl and Interpolations packages to implement some conversions is highly acknowledged.

Several operations on periodic matrices are implemented, such as, inversion, transposing, norms, derivative/time-shifting, trace. All operations with two periodic matrices such as addition/substraction, multiplication, horizontal/vertical concatenation, block-diagonal appending, allow different, but commensurate, periods.

Several advanced computational functions are provided to compute the characteristic multipliers and characteristic exponents of periodic matrices, using methods based on the periodic Schur decomposition of matrix products (provided in the SLICOT library or in the PeriodicSchurDecompositions.jl package) or structure exploitung fast algorithms requiring no external supporting packages. These functions are instrumental to apply Floquet theory to study the properties of solutions of various classes of differential equations (e.g., Mathieu, Hill, Meissner) and the stability of linear periodic systems (see PeriodicSystems.jl package). The implementations of several functions rely on the high performance ODE solvers available in the OrdinaryDiffEq and IRKGaussLegendre packages.

In developing this package I benefitted from the generous support of the Julia developers community, with special mention of @kellertuer, Ralph A. Smith, @ChrisRackauckas, mikelehu, @mkitti, @dlfivefifty.

I am open to constructive critics and suggestions. I am aware of many possible improvements for this package, addressing issues such as, the reduction of the number of allocations, improvements of the definitions of periodic matrix objects, extension with new periodic matrix representations, better support of large periodic matrices, to mention a few of them.

11 Likes
2

Congrats for making your work available to the community! Much appreciated.

Please apologize for being only poorly literate on periodic matrices.

I do, however, wonder how PeriodicMatrices.jl implements operations on truncated harmonic series. This curiosity stems from the fact that we wish to implement the same functionality for the project described at GitHub - ziolai/sediment-transport-rivers .

In particular, I imagine that the multiplication of periodic matrices results in higher order harmonics. A trial example is sin^2 (x) = 0.5 - 0.5 cos(2x). How does PeriodicMatrices.jl handle this case? Does the package truncate the new series or does it instead allow to store the higher harmonic as well?

Many thanks in advance for sharing your thoughts on this matter.

3

Periodic matrices can be represented simply as function matrices or symbolic function matrices, in which cases there is no information stored about the used number of harmonics. However, the HarmonicArray representation stores all harmonic terms (packed in a complex 3-dimensional array) and yes, operations such as a product of two such matrices may change the number of harmonics.
For your example, the result for squaring sin^2(t) is converted to a symbolic representation for better visualization:

using PeriodicMatrices
using Symbolics
A = PeriodicFunctionMatrix(t-> sin(t)^2,2pi) 
Ah = convert(HarmonicArray,A);
convert(PeriodicSymbolicMatrix,Ah*Ah)

and the result is:

PeriodicSymbolicMatrix{:c, Num, Matrix{Num}}(Num[0.375 - 0.5cos(2.0t) + 0.125cos(4.0t);;], 6.283185307179586, 1)

which is 0.375 - 0.5cos(2t) + 0.125cos(4t), as expected.
So, no truncation is performed. Truncation can be however explicitly performed using the function hrtrunc.

4

Note that ApproxFun.jl also implements operations on truncated Fourier series representations of periodic functions (along with several other spectral representations).