The implementation of an algorithm in Macsyma: computing the formal solutions of differential systems in the neighborhood of regular singular point

We discuss in this paper the problems arised in the implementation in Macsyma of a direct algorithm for computing the formal solutions of differential systems in the neighborhood of regular singular point.

The differential system to be considered is of the form x^h dy/dx = A(x) y (1) with A(x) = A₀ + A₁ + … is an n by n matrices of formal series.

This paper deals with system of the form (1) with h = 1. Another paper [3] has considered the systems with h > 1. So these two papers give complete consideration of systems in the form (1) in a neighborhood of both regular and irregular singular point without reducing the systems to super-irreducible forms, i.e. without knowing (for h > 1) a singular point is regular or irregular.

The first step in the algorithm is to transform the leading matrix A₀ to its Jordan form, for which we can use the algorithm in [1][11] that computes at first the Frobenius form and then reduces it to its Jordan form. The reduction from the Frobenius matrix to its Jordan matrix is based on a theorem in Wilkinson [13] which reduces each companion matrix to its Jordan matrix. One of the important aspects is that we don't compute the eigenvalues of the matrix but just keep them as the roots of some corresponding polynomials, i.e. they are looked as algebraic numbers.

One simple case in the computations is that when the eigenvalues of the leading matrix A₀ have no integer differences, which is called the generic case.

The formal solutions of the differential system can be obtained by solving a linear system of the form A X - X B = C with A and B in Jordan form. The implementation of this resolution is immediate. Hence this gives the fundamental system of formal solutions wanted.

For the general case, we should determine at first the integer differences of the eigenvalues of the matrix A₀. A classical method is applied and implemented in Macsyma for our situation. Then we search for some transformations to reduce the differential system to the generic case.

Then an interesting question arised: how to compute the Jordan normal form of the new leading matrix whose entries are rational functions of the eigenvalues which are algebraic numbers. We have studied these problems and proposed two algorithms for two different cases. These algorithms are implemented in Macsyma.