An algorithm for computing Jordan chains and inverting analytic matrix functions (Q2383013)

The main motivation for this paper, according to the author, is solving linear systems of ODE's with linear differential operator with constant coefficients of the form \[ A_l{d^lx\over{dt^l}}+\cdots+ A_0x=0.\tag{1} \] A solution of the equation (1) can be expressed as a vector function \[ x(t)=\sum_{j=0}^k {{t^{k-j}}\over{(k-j)!}}x_j, \] provided that the following identity holds: \[ \sum_{m=0}^k{{t^{k-m}e^{\lambda_0t}} \over{(k-m)!}}\left[\sum_{j=0}^m{1\over{j!}}A^{(j)}(\lambda_0)x_{m-j}\right]=0, \] with \(A(\lambda)=A_0+\cdots+{\lambda}^lA_l\). The sequence of vectors \(\{x_j\}\) forms, generally speaking, the so-called Jordan chain of the matrix \(A(\lambda)\) at \(\lambda_0\). Hence the proposed algorithm is looking for convenient Jordan chains (see Algorithm 3.1). This procedure starts from the kernel of \(A_0\) and then deals with certain extensions of matrices created at consecutive steps, thus extensions of system of columns to full basis are needed. As a byproduct, algorithms for finding the Laurent series for \(A(\lambda)^{-1}\) are given (see the Algorithms 4.1, 4.6). The floating-point and the `symbolic' versions of the algorithms are discussed and compared. Reviewer's remark. I agree with the author that the proposed algorithm (rather complicated and fastidious, mostly in the case of large dimensions) can be useful when a discussion of the solution, or when the `symbolic' computing is needed. However, in the situations most frequently met in computational practice such as, for example, when evaluation of the solution of the initial value problem for (1) is desired, very simple finite difference algorithms are much more effective. Nonetheless, the problem presented in this paper seems to be interesting from the theoretical point of view.