Formal Solutions of Completely Integrable Pfaffian Systems With Normal Crossings
aa r X i v : . [ c s . S C ] O c t Formal Solutions of Completely IntegrablePfaffian Systems With Normal Crossings
Barkatou, Moulay A.
XLIM UMR 7252 , DMI, University of Limoges; CNRS123, Avenue Albert Thomas, 87060 Limoges, France [email protected]
Jaroschek, Maximilian
Max Planck Institute for Informatics,Saarbruecken, Germany [email protected]
Maddah, Suzy S.
Fields Institute222 College St, Toronto, ON M5T 3J1 Ontario, Canada [email protected]
Abstract
In this paper, we present an algorithm for computing a fundamental matrix of formal solutions ofcompletely integrable Pfaffian systems with normal crossings in several variables. This algorithmis a generalization of a method developed for the bivariate case based on a combination of severalreduction techniques and is partially implemented in the computer algebra system Maple . Key words:
Linear systems of partial differential equations, Pfaffian systems, Formalsolutions, Rank reduction, Hukuhara-Turrittin’s normal form, Normal crossings.
1. Introduction
Pfaffian systems arise in many applications [20], including the studies of aerospace,celestial mechanics [16], and statistics [24]. So far, the most important systems for appli-cations are those with so-called normal crossings [31]. Our Maple package PfaffInt can be downloaded at:
It con-tains functionalities illustrated by examples for the splitting, column reduction, rank reduction, and thecomputation of exponential parts of multivariate completely integrable systems with normal crossings
Preprint submitted to Journal of Symbolic Computation 18 October 2018 univariate completely integrable Pfaffian system with normal crossings reduces toa singular linear system of ordinary differential equations (ODS, in short), which havebeen studied extensively (see [5, 34] and references therein). Moreover, unlike the generalcase of several variables considered herein, algorithms to related problems leading to thecomputation of formal solutions have been developed by various authors (see [6, 11, 12, 13]and references therein). The
Maple package
Isolde [14] and
Mathemagix package
Lindalg [28] are dedicated to the symbolic resolution of such systems.More recently, bivariate systems were treated by the first and third author of thispaper in [1]. This paper refines the results of the bivariate case and generalizes them totreat the more general multivariate case.To get an intuition of the kind of systems we consider, we informally study the followingsimple bivariate completely integrable Pfaffian system with normal crossings. A formaldefinition of these systems will be given in Section 2.
Example 1. [1, Example 2] Given the following bivariate system over the ring of formalpower series in ( x , x ) with complex coefficients: x ∂∂x F = A F = x + x + x x − x + x − x ! F,x ∂∂x F = A F = x − x − x − x − x − x − ! F, we are interested in constructing the formal objects F that satisfy the system. Theexistence of a fundamental matrix of solutions and its general form follows from wellknown theoretical results (see Corollary 6). The proof, however, is not constructive. Forsimplicity, we assume we already know that a fundamental matrix of formal solutions inour particular case is of the formΦ( x , x ) x C x C e q ( x − /s ) e q ( x − /s ) , (1)where Φ( x , x ) is a matrix with formal power series entries, C and C are matriceswith entries in C , and q , q are polynomials in C [ z ] , C [ z ] respectively. We now wantto determine Φ , C , C , q , q , s , s . For this purpose, we use the algorithm presentedin [1]. The main idea is to compute one part of the solution by considering an associatedODS in only one variable and then use this information to compute the other parts of thesolution by transforming and decoupling the system into smaller and simpler systems: • First, we construct two associated systems whose equations are derived by settingeither x = 0 or x = 0: x ∂∂x F = A ( x , F = x + x − x + x ! F,x ∂∂x F = A (0 , x ) F = x − x − x − x − x − x − ! F. We show in Section 4 that the formal invariants q , q , s and s can be computedfrom these associated systems. Via Isolde or Lindalg we compute s = s = 1 and q (1 /x ) = − x , q (1 /x ) = ( x + x ) , and (1) becomesΦ( x , x ) x C x C e − x e x + x . Next, we apply the so-called eigenvalue shifting F = e − x e x + x G (for a new unknownvector G ), to facilitate the next step. The shifting yields: x ∂∂x G = x + x x − x − x ! G,x ∂∂x G = x x − − x ! G. • After the eigenvalue-shifting we apply another transformation that reduces the ordersof the singularities in x and x to their minimal integer values. By setting G = T H where T = x x − x , we get: x ∂∂x H = − − x ! H,x ∂∂x H = − − x − ! H. • Finally, via some linear algebra (see [25, Chapter 3] for general cases) we compute thetransformation T = x + 2 x − , and setting H = T U results in the system x ∂∂x U = C U = − ! U,x ∂∂x U = C U = − − ! U. We can now read off C and C . We collect the applied transformations and get afundamental matrix of solutions: T T | {z } =:Φ x C x C e − x e x + x , where C = − and C = − − .Unlike this simple example, the steps of computation can be far more involved anddemand multiple levels of recursion. In order to generalize this algorithm to more thantwo variables, the following nontrivial questions have to be addressed: • Can the information on the formal invariants still be obtained from the associatedODS systems? 3
Can a rank reduction algorithm be developed without relying on properties of principalideal domains as in the univariate and bivariate case?The results of [1] which have an immediate generalization to the multivariate setting arerefined herein and supported by fully transparent proofs and illustrating figures (Theo-rem 13 which answers positively and unconditionally the first question, the structure ofthe main algorithm described in Section 3, and Theorem 18). However, the answer tothe second question is more elaborate (see Section 5) and requires the discussion of twoproblems which are not discussed in [1]: • The major obstacle to a generalization of the results on the bivariate case lies in theprocess of finding integral relations among generators of certain modules over the ringof multivariate power series. In Section 5.4, we propose a solution that relaxes thecondition of working over a principal ideal domain to working over local rings andshow how to utilize Nakayama’s Lemma in the formal reduction process if the modulesunder investigation are free. • We discuss an algorithmic difficulty which arises as not all formal power series undermanipulation admit a finite representation, even if the input Pfaffian system is givenin a finite form (see Example 23). Although this problem arises in the bivariate caseas well, it has not been addressed before (neither in [1] nor in [15, 25]). We provide areasoning to check the correctness of our algorithm.We thus present the first comprehensive description of the state of the art algorithmicapproach for solving completely integrable Pfaffian systems with normal crossings inthe multivariate setting. Our investigation also involves the multivariate versions of thetransformations used classically in the well-studied univariate case (e.g. shearing transfor-mations in Section 5.3, column reduction transformations in Section 5.4, and propertiesof transformations in Proposition 9). Not only does this discussion serve the manipula-tion of such transformations within our proposed formal reduction, but it also plays arole in future generalizations of many other algorithms available for univariate systems(e.g. the alternative rank reduction algorithm of Section 6.2 and the notion of simplesystems as suggested in the conclusion).This paper is divided as follows: In Section 2, we recall the basic definitions and thenecessary theory for our algorithm. This includes the general form of the solutions, thenotion of equivalence between systems, the classification of singularities, and a descrip-tion of the necessary transformations whose generalization to the multivariate case isstraightforward. In Section 3, we give the general structure of our proposed algorithmwhich relies on two major components: The first is associating to our system a set ofODS’s from which its formal invariants can be efficiently derived. This is detailed inSection 4. The second component is the rank reduction which we give in Section 5. Themain algorithm is then given in Section 6 before concluding in Section 7.
2. Preliminaries
The systems considered in this paper are those whose associated differential form isa 1-form. More explicitly, let R := C [[ x , . . . , x n ]] be the ring of formal power series4n x , x , . . . , x n over the field of complex numbers C . A Pfaffian system with normalcrossings is a system of linear partial differential equations of the form x p i +1 i ∂∂x i F = A i F, ≤ i ≤ n, (2)where the A i ’s are d × d matrices with entries in R. The system (2) is completely deter-mined by the A i ’s and p i ’s and we conveniently denote it by [ A ]. Each of the p i ’s is aninteger and the number p ′ i := max(0 , p i ) is called the Poincar´e rank of the i th component A i . The n -tuple p := ( p ′ , . . . , p ′ n )is called the Poincar´e rank of the system [ A ]. If p i ≤ − i ∈ { , . . . , n } thenthe origin is an ordinary (non-singular) point of the system and the system is said tobe regular. In this paper, we tackle the rather more interesting singular systems. Thesingular locus of a system with normal crossings is a union of hyperplanes of coordi-nates x x . . . x n = 0. A Pfaffian system is called completely integrable, if the followingcommutation rule holds for all i, j ∈ { , . . . , n } : A i A j − A j A i = x p i +1 i ∂∂x i A j − x p j +1 j ∂∂x j A i . (3)Subsequently, whenever we refer to a Pfaffian system, we assume it is a completelyintegrable system with normal crossings. For the remainder of this paper we once andfor all fix a Pfaffian system [ A ] for which all the A i are non-zero and there is at least onestrictly positive p i . All subsequent definitions and theorems are stated in this setting,disregarding systems for which the origin is an ordinary point. Our notations follow a set of guidelines in order to help the reader remember themultitude of different objects involved in our work. Single letter identifiers are usuallychosen to be the initial letter of the mathematical term attached to the referenced object,like d for dimension and R for a ring. For a vector v , its i th component is given by v i andfor a univariate power series s the i th coefficient is denoted by s i . We do not distinguishbetween row and column vectors. Upper case letters are used for algebraic structures,matrices and the unknown in a Pfaffian system. A family of matrices is given with lowerindices, e.g. ( M i,j ) i,j ≥ , and for a matrix M i,j , blocks are given with upper indices, e.g. M i,j = M i,j M i,j M i,j M i,j , where the size of the different blocks are clear from the context. By x we denote the collec-tion of variables x , . . . , x n and we use ¯ x i to refer to the variables x , . . . , x i − , x i +1 , . . . , x n .One can expand the A i in system [ A ] as a formal power series with respect to x i : x p i +1 i ∂∂x i F = ( A i, + A i, x i + A i, x i + . . . ) F, where the A i,j are elements of C [[¯ x i ]]. We denote this ring by R ¯ x i . The first coefficient A i, = A ( x i = 0) in such an expansion can be regarded as non-zero without any loss of5enerality, otherwise p i can be readjusted. We call A i, the leading matrix coefficient ofthe i th component.Aside from the rings R and R ¯ x i , we will frequently have to work in other algebraicstructures. We denote by K := Frac(R) (respectively K ¯ x i ) the fraction field of R (respec-tively R ¯ x i ). Let L be the set of monomials given by, L = { x β = x β x β . . . x β n n for β = ( β , β , . . . , β n ) ∈ N n } . Clearly, L is closed under multiplication and contains the unit element. Then, one candefine R L := L − R, the localization of R at L , i.e. the ring of series with only finitelymany terms having monomials of strictly negative exponents. Unlike in the univariatecase, K and R L do not refer to the same algebraic structure, e.g. ( x + x ) − is anelement of K which is not an element of R L . In fact, there exists no β ∈ N such that x β ( x + x ) − ∈ R: ( x + x ) − = /x − ( − x /x ) whose formal expansion with respect to x is P i ≥ ( − i x − i − x i , which has infinitely many poles in x . For a further characterizationof K, one may refer to [3].Finally, in the sequel it will be necessary to introduce ramifications of the form x i = t α i i for new variables t i and positive integers α i . We will therefore write R t for C [[ t , . . . , t n ]]and allow analogous notations for all structures introduced so far. The identity andzero matrices of dimension d are denoted by I d and O d , and we set: GL d (R) = { M ∈ Mat d × d (R) | det(M) is invertible in R } . As we have seen in the introductory example, we will make use of transformations tobring a system into particular forms. Such a transformation acts on a Pfaffian systemas follows: A linear transformation (also called gauge transformation) F = T G , where T ∈ GL d (K), applied to (2) results in the system x ˜ p i +1 i ∂∂x i G = ˜ A i G, ≤ i ≤ n, (4)where ˜ A i x ˜ p i +1 i = T − ( A i x p i +1 i T − ∂∂x i T ) , ≤ i ≤ n. (5)We say that system (4) is equivalent to system (2) and we write T [ A ] := [ ˜ A ]. It can beeasily verified that complete integrability is inherited by an equivalent system. Subse-quently, to stay in the same class of systems under study, special care will be taken sothat the transformations used in our considerations do not alter the normal crossings. Infact, a major difficulty within the symbolic manipulation of system (2) arises from (5).It is evident that any transformation alters all the components simultaneously. In par-ticular, the equivalent system does not necessarily inherit the normal crossings even forvery simple examples. 6 xample 2. [15, Section 4] Consider the following completely integrable Pfaffian systemwith normal crossings of Poincar´e rank (3 , x ∂∂x F = A ( x , x ) F = x + x x − − x + x ! F,x ∂∂x F = A ( x , x ) F = x x − − x ! F. This system appears within the reduction of the system of Example 1 in the introduction.As we have seen, there exists a transformation which drops p to zero. This can also beattained by the transformation F = x − x x G, which is computed by the univariate-case Moser-based rank reduction algorithm, uponregarding the first component as an ODS in x and x as a transcendental constant. Thisresults in the equivalent system x x ∂∂x G = ˜ A ( x , x ) G = − x − x ! G,x ∂∂x G = ˜ A ( x , x ) G = − x − x − x ! G. We can see that such a transformation achieves the goal of reducing the Poincar´e rank ofthe first component. However, it alters the normal crossings as it introduces the factor x on the left hand side of the first component. Moreover, it elevates the Poincar´e rank ofthe second component.In order to preserve the normal crossings, we restrict the class of transformations thatwe use in our algorithm: Definition 3.
Let T ∈ GL d (K). We say that the transformation F = T G (respec-tively T ) is weakly compatible with system [ A ] if T [ A ] := ˜ A is again a completelyintegrable Pfaffian system with normal crossings. In particular, ˜ A i ∈ R d × d for every i ∈ { , . . . , n } . Clearly, any constant or unimodular invertible matrix is an example of such transfor-mations.In the sequel, we will also need to resort to transformations with stronger properties:
Definition 4.
Let T ∈ GL d (K). We say that the transformation F = T G (respec-tively T ) is compatible with system [ A ] if it is weakly compatible with [ A ] and thePoincar´e rank of each individual component of T [ A ] does not exceed that of the respec-tive component of [ A ]. 7 .4. Fundamental Matrix of Formal Solutions Before studying how to construct formal solutions to a given system, the questionarises if and how many solutions exist. The language of stable modules over the ringof power series is used in [18, Theorem 1] and [22, Main Theorem] independently toestablish the following theorem which gives an answer to this question.
Theorem 5.
There exist strictly positive integers α i , 1 ≤ i ≤ n , and an invertible matrix T ∈ R d × dt such that, upon setting x i = t α i i , the transformation T ( t ) yields the followingequivalent system: t α i ˜ p i +1 i ∂∂t i G = ˜ A i ( t i ) G, ≤ i ≤ n, where ˜ A i ( t i ) = Diag( ˜ A i ( t i ) , ˜ A i ( t i ) , . . . , ˜ A jji ( t i )) , and for every ℓ ∈ { , . . . , j } we have that ˜ A ℓℓi ( t i ) is a square matrix of dimension d ℓ ofthe form ˜ A ℓℓi ( t i ) = w ℓℓi ( t i ) I d ℓ + t α i ˜ p i i ( c ℓℓi I d ℓ + N ℓℓi ) , where • d + d + · · · + d j = d ; • w ℓℓi ( t i ) = P α i ˜ p i − m =0 λ imℓ t mi is a polynomial in t i , with coefficients in C ; • c ℓℓi ∈ C and N ℓℓi is a constant (with respect to all derivations ∂/∂t i ) d ℓ -square matrixhaving nilpotent upper triangular form; • for any fixed ℓ ∈ { , . . . , j } , the matrices { N ℓℓi } i =1 ,...,n are permutables; • for all ℓ ∈ { , . . . , j − } , there exists i ∈ { , . . . , n } such that w ℓℓi ( t i ) = w ( ℓ +1)( ℓ +1) i ( t i ) or c ℓℓi − c ( ℓ +1)( ℓ +1) i Z . This theorem guarantees the existence of a transformation which takes system (2)to the so-called Hukuhara-Turrittin’s normal form from which the construction of afundamental matrix of formal solutions (6) is straightforward. In fact, we have:
Corollary 6.
Given system (2), a fundamental matrix of formal solutions exists and isof the form Φ( x /s , . . . , x /s n n ) n Y i =1 x C i i n Y i =1 exp( Q i ( x − /s i i )) , (6)where Φ is an invertible matrix with entries in R t and for each i ∈ { , . . . , n } we have: • s i is a positive integer; • the diagonal matrix Q i ( x − /s i i ) = Diag (cid:16) q i, ( x − /s i i ) , q i, ( x − /s i i ) , . . . , q i,d ( x − /s i i ) (cid:17) contains polynomials in x − /s i i over C without constant terms. We refer to Q i ( x − /s i i )as the x i -exponential part. Under the notations of Theorem 5, it is obtained by formallyintegrating w ℓℓi t αi ˜ pi +1 i ; • C i is a constant matrix which commutes with Q i ( x − /s i i ).8 singular system [ A ] is said to be regular singular whenever, for every i ∈ { , . . . , d } , Q i ( x − /s i i ) is a zero matrix. Otherwise, system (2) is said to be irregular singular and theentries of Q i ( x − /s i i ), 1 ≤ i ≤ n , determine the main asymptotic behavior of the actualsolutions as x i → Definition 7.
Let i ∈ { , . . . , n } . If Q i ( x − /s i i ) is a nonzero matrix then we set m i,j tobe the minimum order in x i within the terms of q i,j ( x − /s i i ) for 1 ≤ j ≤ d . The x i - formalexponential growth order ( x i -exponential order, in short) of A i is the rational number ω ( A i ) = − min ≤ j ≤ d m i,j . The n -tuple of rational numbers ω ( A ) = ( ω ( A ) , . . . , ω ( A n )) then defines the exponentialorder of system [ A ]. Otherwise, we set ω ( A i ) = 0.If two systems are equivalent then they have the same x i -exponential parts, and con-sequently the same x i exponential orders, for all 1 ≤ i ≤ n , under any transformation T ∈ GL d (K). Example 8 (Example 1 cont.) . From our investigations in the example of Section 1,we see that for the given fundamental system of formal solutions, we have non-zeroexponential parts with ω ( A ) = 1 and ω ( A ) = 2 and so the system is irregular singular(although s = s = 1).The above theoretical results on existence do not establish the formal reduction itself,that is the algorithmic procedure which computes explicitly the α i ’s and a transformationwhich takes the system to a normal form that allows the construction of such solutions.This will be our interest in the following sections.The computation of the formal invariants is a difficult task in the univariate case [6,9]. However, we will prove in Section 4 that in the multivariate case, these invariantscan be computed from associated univariate systems. Unlike the univariate case, themain difficulties of the algorithm lie in rank reduction. Before proceeding to describe thealgorithms we propose, we give a property of the transformations which can be deployed: Proposition 9.
Consider a completely integrable Pfaffian system [ A ] with normal cross-ings. Let T ∈ GL d (K) and set T [ A ] = [ ˜ A ]. If T is a transformation which is weaklycompatible with system [ A ] then T ∈ GL d ( R L ). Proof.
It follows from (5) that ∂∂x i T = A i x p i +1 i T − T ˜ A i x ˜ p i +1 i , ≤ i ≤ n. Thus, we have (see, e.g. [5, Proposition 1, proof, pp 6]): ∂∂x i det( T ) = ( tr( A i ) x p i +1 − tr( ˜ A i ) x ˜ p i +1 ) det( T ) , ≤ i ≤ n. (7)Therefore, det( T ) itself is a solution of a completely integrable Pfaffian system withnormal crossings. By Corollary 6, det( T ) has the form (6). Since T ∈ K d × d then det( T ) isfree of logarithmic and exponential terms. Hence, det( T ) corresponds to a log-free regular9 nputsystem FirstcomponentLastcomponent System w.lower dim.System w.lower dim. ≥ ≥ Fig. 1. Computing a fundamental matrix of formal solutions by working with one of the com-ponents, e.g. the first component. The other components would follow the chosen component inthe uncoupling. solution of (7). Thus, det( T ) ∈ R L . The same argument serves to prove that det( T − )is an element of R L as well, upon remarking that T − [ ˜ A ] = [ A ]. Hence, det( T ) − ∈ R L and consequently T ∈ R d × dL . ✷ However, the converse of Proposition (9) is not true, which complicates the task ofconstructing adequate transformations in the reduction process (see, e.g., Example 2 orthe shearing transformations of Section 5.3).
3. Structure of the Main Algorithm
If one is only interested in the asymptotic behavior of the solutions of system [ A ],then one can compute the formal invariants from associated univariate systems as weprove in Section 4. If the singularity is regular or one is interested in computing a fullfundamental matrix of formal solutions, as given by (6), then, besides computing theseinvariants, further involved steps are required, as illustrated in Example 1. The recursivealgorithm we propose generalizes that of the univariate case given by the first authorin [6]. At each level of recursion with input [ A ], we consider the leading matrix coefficients A i, = A i ( x i = 0) (we use both notation interchangeably) and distinguish between threemain cases:(1) There exists at least one index i ∈ { , . . . , n } such that A i, has at least two distincteigenvalues.(2) All of the leading matrix coefficients have exactly one eigenvalue and there existsat least one index i ∈ { , . . . , n } such that A i, has a nonzero eigenvalue.(3) For all i ∈ { , . . . , n } , A i, is nilpotent.10n order to identify the properties of the eigenvalues of A i ( x i = 0), it suffices toconsider the constant matrix A i ( x = 0) due to the following well-known proposition (see,e.g., [22, Proposition 1, pp 8] or [8, Proposition 2.2] for a proof within the context ofeigenrings): Proposition 10.
The eigenvalues of A i, , 1 ≤ i ≤ n , belong to C .Then, based on the above classification, a linear or an exponential transformation willbe computed as described in the following subsections. Whenever there exists an index i ∈ { , . . . , n } such that A i, has at least two distincteigenvalues, the system can be uncoupled into subsystems of lower dimensions as shownin Theorem 11. For a constructive proof, one may refer to [23, Section 5.2, pp 233]. Theorem 11.
Suppose that for some i ∈ { , . . . , n } , the leading matrix coefficient A i, has at least two distinct eigenvalues. Then there exists a unique transformation T ∈ GL d (R) of the form T ( x ) = T T T T = I d ′ T ( x ) T ( x ) I d − d ′ , where 0 < d ′ < d , such that the transformation F = T G yields the equivalent system x p i +1 i ∂∂x i G = ˜ A i ( x ) OO ˜ A i ( x ) G, ≤ i ≤ n. and ˜ A i ( x ) , ˜ A i ( x ), i ∈ { , . . . , n } are of dimensions d ′ and d − d ′ respectively.The theorem can be restated by saying that if one of the components of the systemhas a leading matrix coefficient with at least two distinct eigenvalues, then it can beuncoupled. All of the other components are uncoupled simultaneously. In the sequel, weaim to determine changes of the independent variables x i (ramifications), and constructtransformations, which will allow the reduction of any input system to a system forwhich the leading matrix coefficient of at least one of its components has at least twodistinct eigenvalues. This allows us to either arrive at a system with lower Poincar´e rankor uncouple it into several subsystems of lower dimensions. The recursion stops wheneverwe arrive at regular singular ( p = (0 , . . . , d = 1) subsystems. The formerhave been already investigated in [25, Chapter 3] and the resolution of the latter isstraightforward.We remark that, by Proposition 10, it suffices that there exists i ∈ { , . . . , n } suchthat the constant matrix A i ( x = 0 , . . . , x i = 0 , . . . , x n = 0) has at least two distincteigenvalues. 11 .2. Unique Eigenvalue: Shifting For any i ∈ { , . . . , n } such that A i, has a unique nonzero eigenvalue γ i ∈ C , applyingthe so-called eigenvalue shifting F = exp (cid:18)Z x i γ i z − p i − i dz i (cid:19) G, yields a system [ ˜ A ] whose i th component has a nilpotent leading matrix coefficient: x ˜ p i +1 i ∂∂x i G = ˜ A i ( x ) G, where ˜ A i ( x ) = A i ( x ) − γ i I d . The other components of the system are not modified by this transformation which isclearly compatible with system [ A ].Hence, due to the uncoupling and shifting, we can assume without loss of generalitythat for all i ∈ { , . . . , n } , the leading matrix coefficients A i, are nilpotent. In the univariate case, n = 1, the nilpotency of A , suggests at least one of the fol-lowing two steps, as proposed by the first author in [6]: Rank reduction and computationof the exponential order ω ( A ). The former reduces p to its minimal integer value. Itis possible that p drops to zero, i.e. we arrive at a regular singular system, or that theleading matrix coefficient of the resulting system has at least two distinct eigenvalues,in which case we can again uncouple the system. Otherwise, ω ( A ) = ℓ/m is to be com-puted, where ℓ and m are coprime. Then, by setting x = t m and applying rank reductionagain, it is proven that we arrive at a system whose leading matrix coefficient has twodistinct eigenvalues. Therefore, the system can be uncoupled (see Figure 1).The bivariate case, n = 2, is studied by the first and third authors of this paper in [1].For rank reduction, the properties of principal ideal domains were used. To determinethe formal exponential order ω ( A ), associated univariate systems were defined. In thispaper, we show that on the one hand, this approach to determine the formal exponentialorder remains valid in the multivariate setting, as we will see in the next section. On theother hand, the generalization of the rank reduction algorithm to the multivariate caseis nontrivial and is discussed in Section 5. The multivariate formal reduction algorithmis then summed up in Section 6.
4. Computing the Formal Invariants
In the univariate case, where the system is given by a single matrix A , ω ( A ) canbe computed from the characteristic polynomial of A , i.e. det( λI d − A ), based on theanalysis of a Newton polygon associated with the system [6, Theorem 1]. In this sectionwe show that one need not search for a generalization of this algorithm to the multivariatecase as the formal invariants of [ A ], i.e. the exponential parts and ω ( A ), can be obtainedfrom an associated univariate system. We do not only give a method to retrieve theseinvariants but we also reduce computations to computations with univariate rather thanmultivariate formal series. 12 nputsystem FirstcomponentLastcomponent FirstassociatedODSLastassociatedODS Exp. partin first var.Exp. partin last var. Exp. part Fig. 2. Computing the exponential part from associated ODS’s
Definition 12.
Given a Pfaffian system [ A ], we call the following the associated ODS of [ A ]: x p i +1 i ddx i F i = A i ( x i ) F i , ≤ i ≤ n, where A i ( x i ) := A ( x = 0 , . . . , x i − = 0 , x i , x i +1 = 0 , . . . , x n = 0) . Theorem 13.
For every i ∈ { , . . . , n } , the x i -exponential part of a Pfaffian system isequal to the exponential part of the i th component of its associated ODS.To establish this result, we rely on a triangular form weaker than the Hukuhara-Turrittin’s normal form given in Theorem 5. This weaker form suffices to give insightinto the computation of (6).The following theorem is an reformulation of a theorem which was first given in [17,Proposition 3, pp 654] for the bivariate case, and then generalized in [18, Theorem 2.3]to the general multivariate case. Theorem 14.
Consider the Pfaffian system [ A ]. There exists a positive integer α , anda transformation T ∈ GL d (K t ) (where x = t α and x i = t i , ≤ i ≤ n ), such that thetransformation F = T G yields the equivalent system: t α ˆ p +11 ∂∂t G = ˆ A ( t , x , . . . , x n ) G,x ˆ p i +1 i ∂∂x i G = ˆ A i ( x , . . . , x n ) G, ≤ i ≤ n, (8)where ˜ A ( t , x , . . . , x n ) = Diag( ˆ A , ˆ A , . . . , ˆ A jj ) , ˆ A i ( x , . . . , x n ) = Diag( ˆ A i , ˆ A i , . . . , ˆ A jji ) , ≤ i ≤ n, and for all ℓ ∈ { , . . . , j } and i ∈ { , . . . , n } the entries of ˆ A ℓℓi lie in R ¯ x . The ˆ A ℓℓ ’s are ofthe form ˆ A ℓℓ = w ℓℓ ( t ) I d ℓ + t α ˆ p ( ˆ N ℓℓ ( x , . . . , x n ) + c ℓℓ I d ℓ ) , where • d + d + · · · + d j = d ; • w ℓℓ ( t ) and c ℓℓ are as in Theorem 5; • If ℓ, ℓ ′ ∈ { , . . . , j − } and ℓ = ℓ ′ , then w ℓℓ ( t ) = w ℓ ′ ℓ ′ ( t ) or c ℓℓ − c ℓ ′ ℓ ′ Z ;13 ˆ N ℓℓ ( x , . . . , x n ) is a nilpotent d ℓ -square matrix whose entries lie in R ¯ x .Moreover, T can be chosen as a product of transformations in GL d (R t ) and transforma-tions of the form Diag( t β , . . . , t β d ), where β , . . . , β d are non-negative integers. ✷ Proof of Theorem 13.
Upon the change of independent variable x = t α , the transfor-mation F = T G yields system (8) for which the first component is given by t α ˆ p +11 ∂∂t G = ˜ A ( t , x , . . . , x n ) G. with the notations and properties as in Theorem 14. It then follows from (5) that t α ˆ p +11 ∂∂t T = α A ( x = t α ) T − T ˆ A . (9)Due to the particular choice of T in Theorem 14, we can set x i = 0, 2 ≤ i ≤ n in (9). Inparticular, the relation between the leading terms A ( x = t α ) := A ( x = t α , x = 0 , . . . , x n = 0) , ˆ A ( t ) := ˆ A ( t , x = 0 , . . . , x n = 0) , T ( t ) := T ( t , x = 0 , . . . , x n = 0) , is given by t α ˆ p +11 ∂∂t T = α A ( x = t α ) T − T ˆ A . Hence, the systems given by α A ( x = t α ) (respectively A ) and ˆ A are equivalent.It follows that they have the same formal invariants. Clearly, the same result can beobtained for any of the other components via permutation with the first component.Noting that the x i -exponential part is independent of ¯ x i completes the proof. ✷ For univariate systems, the true Poincar´e rank p true ( A ) is defined as the smallestinteger greater or equal than the exponential order ω ( A ) of the system [ A ]. It is knownthat this integer coincides with the minimal value for p which can be obtained uponapplying any non-ramified linear transformation to A . With the help of Theorem 13 wecan establish the analogous result for multivariate systems. We first give the followingdefinition: Definition 15.
Let [ A ] be a Pfaffian system and for any i ∈ { , . . . , n } , let p true ( A i ) bethe minimal integer value which bounds the exponential order in x i , i.e. p true ( A i ) − < ω ( A i ) ≤ p true ( A i ) . Then p true ( A ) = ( p true ( A ) , . . . , p true ( A n )) is called the true Poincar´e rank of A .It is shown by Deligne and van den Essen separately in [19, 21], in the multivariatesetting that a necessary and sufficient condition for system [ A ] to be regular singularis that each individual component A i , considered as a system of ordinary differentialequations in x i , with the remaining variables held as transcendental constants, is regularsingular. As a consequence, system [ A ] is regular singular if and only if its true Poincar´erank is (0 , , . . . , n = 1 (e.g. [13, 26]) can be applied separately to each of the individual components. Thefollowing corollary follows directly from Theorem 13, showing that the i th component14f the true Poincar´e rank of system [ A ] is equal to the true Poincar´e rank of the i th associated (univariate) ODS. Corollary 16.
For all 1 ≤ i ≤ n we have p true ( A i ) = p true ( A i ) . Proof.
For the proof, it suffices to remark that p true ( A i ) − < ω ( A i ) = ω ( A i ) ≤ p true ( A i ) . ✷ From this Corollary it does not yet follow for the multivariate case, as in the univariatecase, that it is possible to apply a compatible transformation to system [ A ] such thatall the p i simultaneously equal the p true ( A i ). We investigate this possibility in the nextsection.In summary, the formal exponential order, the true Poincar´e rank, and most impor-tantly the Q i ’s in (6), can be obtained efficiently by computations with univariate ratherthan multivariate series, making use of existing algorithms and packages. As mentionedin the introduction, this exponential part is of central importance in applications since itdetermines the asymptotic behavior of the solution in the neighborhood of an irregularsingularity. To compute a full fundamental matrix of formal solutions, we still have todetermine suitable rank reduction transformations. Transformations which reduce therank of the associated systems do not suffice, since they are not necessarily compatible.We therefore proceed to develop a multivariate rank reduction algorithm.
5. Rank Reduction
In this section, we are interested in the rank reduction of Pfaffian systems, morespecifically, the explicit computation of a transformation which, given system [ A ], yieldsan equivalent system whose Poincar´e rank is the true Poincar´e rank. We show that, undercertain conditions, the true Poincar´e ranks of the subsystems of [ A ] can be attainedsimultaneously via a transformation compatible with [ A ].We first generalize Moser’s reduction criterion [30] to multivariate systems. We thenestablish an extension of the algorithm we gave in [1] for the bivariate case to the mul-tivariate setting. The main problem in the treatment of multivariate systems is that theentries of the A i do not necessarily lie in a principal ideal domain. This is a commonproblem within the study of systems of functional equations. The same obstacle arisesin [15] and in the analogous theory of formal decomposition of commuting partial lineardifference operators established in [32]. For univariate systems, Moser’s criterion characterizes systems for which rank reduc-tion is possible. To adapt this criterion to our setting, we follow [7, 30] and define thegeneralized Moser rank and the Moser invariant of a system [ A ] as the following n -tuplesof rational numbers: Stronger bounds are given in [6, Remark 3] ( A ) = ( m ( A ) , . . . , m ( A n )) , where m ( A i ) = max (cid:18) , p i + rank( A i, ) d (cid:19) ,µ ( A ) = ( µ ( A ) , . . . , µ ( A n )) , where µ ( A i ) = min( { m ( T ( A i )) | T ∈ GL d (R L ) } ) . We remark that µ ( A ) is well-defined due to Corollary 16. Definition 17.
Consider the partial order ≤ on Q n for which ℓ < k holds if and only if ℓ j ≤ k j for all 1 ≤ j ≤ n and there is at least one index for which the inequality is strict.The system [ A ] is called reducible if µ ( A ) < m ( A ). Otherwise it is said to be irreducible.In other words, system [ A ] is irreducible whenever each of its components is. In par-ticular, it is easy to see from this definition that a system [ A ] is regular singular if andonly if µ ( A i ) ≤ i ∈ { , . . . , n } , i.e. the true Poincar´e rank is a zero n -tuple, whichcoincides with Deligne’s and van den Essen’s criterion for regular singular systems. With the help of compatible transformations and the criterion established in Sec-tion 5.1, we study the rank reduction of some component A i of [ A ] which is givenby (2). We will see that rank reduction can be carried out for each component inde-pendently without affecting the individual Poincar´e ranks of the other components. Wefix i ∈ { , . . . , n } and we recall that one can expand the components of [ A ] with respectto x i . In particular, we have, x p i +1 i ∂∂x i F = A i ( x ) F = ( A i, (¯ x i ) + A i, (¯ x i ) x i + A i, (¯ x i ) x i + . . . ) F. We set r := rank( A i, ) . For all i we can assume without loss of generality that A i, is not the zero matrix andthus the reducibility of system [ A ] coincides with the existence of an equivalent systemsuch that for some i the rank of the leading matrix coefficient ˜ A i, is less than r . Weestablish the following theorem: Theorem 18.
Consider a Pfaffian system [ A ] and suppose that m ( A i ) > i ∈ { , . . . , n } . If µ ( A i ) < m ( A i ) then the polynomial θ i ( λ ) := x ir det( λI + A i, x i + A i, ) | x i =0 (10)vanishes identically in λ . Proof.
Suppose that there exists a transformation T ( x ) ∈ R d × dL which reduces m ( A i ) forsome i ∈ { , . . . , n } . That is, setting T [ A ] := [ ˜ A ], we have:m( ˜ A i ) < m( A i ) . (11)The i th component of system [ A ] can also be viewed as a system of ordinary differen-tial equations (ODS) in x i upon considering the x j ’s for j ∈ { , . . . , n } , j = i , to betranscendental constants. Hence, by (11) and [30, Theorem 1], θ i ( λ ) = 0. ✷ A i /x i is used in (10) to detect the truePoincar´e rank of the i th component via the valuation of x i . It turns out that the valuationis only influenced by A i, and A i, . Even though the true Poincar´e rank can be determinedfrom the associated ODS, the criterion is essential as it leads to the construction of thetransformation T .The converse of Theorem 18 also hold true under certain conditions (see Theorem 26):If θ i ( λ ) vanishes for some index i ∈ { , . . . , n } then we can construct a compatibletransformation T ∈ GL d (R L ) which reduces m ( A i ). We will establish this result anddescribe the steps of the algorithm after establishing a series of intermediate ones. Wewill need two kinds of transformations, shearing transformations and column reductions,which we explain in the next two subsections. Consider the expansion A i = P ∞ k =0 A i,k x ki of A i with respect to x i for a fixed i ∈ { , . . . , n } . Shearing transformations are polynomial transformations that, roughlyspeaking, are used to exchange blocks between the A i,k ’s. The ones we consider here areof the form S = Diag( x β i , . . . , x β d i ) , with β j ∈ { , } for all j ∈ { , . . . , d } . We illustrate the shearing effect of such a trans-formation in an easy example. Example 19.
We apply a shearing transformation to a univariate system [ A ] given by x p +11 ∂∂x F = A F where A = P ∞ k =0 A ,k x k , with A , = − and A , = −
58 9 0 08 6 2 45 6 3 3 . The transformation S = Diag( x , x , ,
1) yields the equivalent system [ ˜ A ] given by x p +11 ∂∂x F = ˜ A F where ˜ A = S − A S − x p Diag(1 , , , . As we are interested in the effect of S on the first few terms of A , we look into S − A S which exchanges the upper right and lower left 2 × A as exhibited in thefollowing diagram: x + x + x + x + . . .x + x + x + x + . . . x + x + x + x + . . . ∗ Consequently, A , and A , become˜ A , = −
50 0 0 00 0 0 00 0 0 0 , ˜ A , = ∗ ∗ ∗ ∗− . Note that the lower left zero entries in ˜ A , come from A , − , which is a zero matrix.In return, the upper right block of A , is sent to A , − . Since it is a zero block, thistransformation does not introduce denominators. The upper right entries in ˜ A , comefrom A , . With this transformation, we reduced the rank of A , from 2 to 1.More generally, let [ A ] be a multivariate Pfaffian system. The transformation F = SG ,where S is a shearing transformation in x i , yields an equivalent system with: ( ˜ A i = S − A i S − x p i i Diag( β , . . . , β d ) , ˜ A j = S − A j S, ≤ j = i ≤ n, where S − A j S = A j, A j, x β − β i . . . A j, d x β d − β i A j, x β − β i A j, . . . A j, d x β d − β i ... ... . . . ... A j,d x β − β d i A j,d x β − β d i . . . A j,dd for all 1 ≤ j ≤ n. The shearing in Example 19 reduced the rank of the leading matrix coefficient andwas compatible with the system, i.e. it did not introduce undesired denominators of x i ,because of the column reduced form of A i, . The input system is not always given insuch a form for A i, , and so we investigate in the following subsection how to achieve it. To enable rank reduction, we alternate between the shearing transformation and trans-formations which reduce some columns of a leading matrix coefficient to zero. For thiswe discuss in this section the following problem.(P) Given a square matrix A = [ v , . . . , v d ] ∈ Mat d × d (R) (where v i denotes the i thcolumn) of rank r < d when considered as an element of Mat d × d (K), find a transfor-mation T ∈ GL d (R) such that the last d − r columns of T AT − are zero.Before considering the algorithmic aspects, we first discuss the existence of such atransformation. As the next example shows, the desired transformation does not neces-sarily exist for any matrix A . 18 xample 20. The matrix x x is obviously of rank 1 over K. There is, however, no transformation T ∈ GL (R) suchthat T AT − contains only one non-zero column.Consider the finitely generated R d -submodule M := h v , . . . , v d i . We call it the col-umn module of A . In order to construct a suitable transformation for bivariate Pfaffiansystems (for which the leading matrix coefficients are univariate) the authors of [1] usethe fact that C [[ x ]] and C [[ x ]] are principal ideal domains and hence that every finitelygenerated submodule of a free module over this ring is free. We generalize this for themultivariate case by showing in Corollary 22 that the freeness of the column module M is a necessary and sufficient condition for the existence of a transformation that meetsour requirements. This is a direct consequence of Nakayama’s Lemma for local rings. Theorem 21. [29, Theorem 2.3, pp 8] Let R be a local ring, M its maximal idealand let M be a finitely generated R-module. Then v , . . . , v r ∈ M form a minimal set ofgenerators for M if and only if their images ¯ v , . . . , ¯ v r under the canonical homomorphism M → M/ M M form a basis of the vector space M/ M M over the field R / M .The central consequence of Theorem 21 for us is that if M is free, a module basis of M can be chosen among the columns of A . We adapt the theorem to our situation toshow that we can bring A into a column-reduced form if and only if its column moduleis free. Corollary 22.
Let A ∈ Mat d × d (R) be of rank r over K and let M be the modulegenerated by the columns of A . If M is free, then there exists a subset B of the columnsin A with r elements such that B is a module basis of M . Furthermore, B is also aK-vector space basis of the column space of A . Proof.
By Theorem 21 we can find a basis B = { b , . . . , b k } of M among the columns of A .By definition, the b i are linearly independent over R, so they are also linearly independentover K (otherwise, multiplying a linear relation in K with a common denominator yieldsa relation in R). Since B is a basis of the column module, it also contains a generatingset of the K-vector space generated by the columns of A . ✷ In theory, Corollary 22 would allow the computation of a unimodular column reductiontransformation as required in ( P ) simply via Gaussian elimination. Assume we are givena matrix A and already know a subset B = ( b , . . . , b r ) of the columns of A which formsa basis of the column module. Let v be a column vector of A which is not in B . Then,since B is a vector space basis, there exist c , . . . , c r ∈ K such that c b + · · · + c r b r = v. By assumption, B is also a module basis, so there also exist d , . . . , d r ∈ R with d b + · · · + d r b r = v. b i are linearly independent, and therefore the cofactors of v with respect to B areunique. It follows that c i = d i for all 1 ≤ i ≤ r .The main algorithmic difficulty stems from the fact that not all formal power seriesadmit a finite representation, and even if the initial system is given in a finite form, thesplitting transformation as in Theorem 11 does not preserve finiteness. In particular, weface two main problems when working with truncated power series:( P
1) Detecting the correct rank and the linear independent columns of A ( P
2) If we know the independent columns, a column reduction transformation computedafter truncation is not uniquely determined.These computational problems arise for general multivariate and for bivariate systems,but were not addressed in previous algorithmic works on this topic [1, 15, 25]. Before wepropose our resolution, we illustrate both problems in the following example:
Example 23.
Consider the matrix x x x + x x x x . Here, the first three columns v , v , v are linearly independent and generate the columnmodule. A linear combination of the fourth column v is given by1 · v + 0 · v + 1 · v = v . When truncating at order 1, the system is given as x x x x x The original rank cannot be determined from the truncated matrix. Furthermore, evenif we know that v , v , v are linearly independent, there are several linear combinationsof the fourth column after truncation:1 · v + 0 · v + 1 · v = v . · v + 1 · v + 0 · v = v . The cofactors of the second linear combination are not the truncated cofactors of thefirst. It can not be extended with higher order terms to a suitable linear combinationover the formal power series ring without truncation.We can solve both ( P
1) and ( P
2) with the help of minors of the original system. Let r be the rank of A . Then there exists a nonzero r × r submatrix B of A whose determinantis nonzero. Let k be the order of the determinant. If we take the truncated system˜ A = A rem x k +1 , the same submatrix ˜ B in ˜ A will have a non-zero determinant modulo x k +1 and we can therefore identify in ˜ A which columns in A are linearly independent.This resolves ( P
1) as long as the truncation order k is chosen big enough.Next assume that for instance the first r columns of A are linearly independent, i.e.we can choose B such that its columns correspond to v , . . . , v r . Let k be as above, ℓ be20 positive integer and let v be a column vector that is linearly dependent on the columnsof B . Then there exist c , . . . , c r ∈ C [[ x , . . . , x n ]] such that B · ( c , . . . , c r ) = v. By Cramer’s rule, we know that the c i are given by c i = det( B i )det( B ) , (12)where B i is the matrix obtained by replacing the i th column of B by v . RewritingEquation (12) gives det( B ) c i − det( B i ) = 0 , (13)and this equation allows the computation of c i by coefficient comparison. In particular,we are guaranteed to obtain the correct c i up to order ℓ if in (13) we replace B by ˜ B , itstruncation at order ℓ + k + 1, and B i by ˜ B i , the truncation of B i at order ℓ + k + 1. Thisresolves ( P k such that we canfind a submatrix of maximal dimension with non-zero determinant. We have to remark,however, that by the nature of formal power series, it is in general not possible to tella priori if a given truncation is high enough. Furthermore, we emphasize that it is ingeneral not possible to draw a conclusion about the freeness of the column module fromthe integral relations among the truncated column vectors, since any linear combinationof the form c v + · · · + c d − v d − − v d = 0 mod x k can require a non-unit cofactor c d forhigher truncation orders. However, if no integral relations can be found with the abovemethod, also the column module without truncation cannot be free. Both observationslead to the following practical approach. The full algorithm is carried out with a giventruncation order. If the output is correct (compared to the invariant exponential partwhich can be obtained by Theorem 13), we are done. If not, we increase the truncationorder until we get a correct output or arrive at a point where no integral relations canbe found anymore. This procedure necessarily terminates, since there exists a suitabletruncation order.One should note that not every K-vector space basis of the column space of A is also amodule basis. So, in the worst case, (cid:0) dr (cid:1) submatrices have to be tested to obtain a modulebasis. We consider again a multivariate system [ A ] as in (2). We fix i ∈ { , . . . , n } andinvestigate the rank reduction of its i th component given by x p i +1 i ∂∂x i F = A i F = ( A i, + A i, x i + A i, x i + A i, x i + . . . ) F, (14)where the matrices A i,j have their entries in R ¯ x i and the algebraic rank of A i, is denotedby r . We recall that we defined ¯ x i := ( x , . . . , x i − , x i +1 , . . . , x n ) and R ¯ x i := C [[¯ x i ]].The establishment of the converse of Theorem 18 for the reduction in x i follows es-sentially the steps of that of the bivariate case which was given in [1]. The constructionrequires successive application of transformations in GL d ( R ¯ x i ) and shearing transforma-tions in x i . We remark that when applying a transformation T ∈ GL d (R ¯ x i ) on the i th component, (5) reduces to ˜ A i = T − A i T. A i, is free, then one can compute a transformation U ∈ GL d (R ¯ x i ) such that U − A i, U = B OB O has rank r , entries in R ¯ x i and with diagonal blocks of sizes r × r and ( d − r ) × ( d − r )respectively. Let v be the rank of B . If also the column module of B is free, then onecan compute a transformation U ∈ GL r (R ¯ x i ) such that U − B U = E OE O has rank v , entries in R ¯ x i and with diagonal blocks of sizes v × v and ( r − v ) × ( r − v )respectively. We set U := Diag( U , I d − r ) · U . Then the leading coefficient ˜ A i, of theequivalent system U [ A i ] has the following form:˜ A i, = ˜ A i, O O ˜ A i, O r − v O ˜ A i, ˜ A i, O d − r (15)with diagonal blocks of sizes v × v , ( r − v ) × ( r − v ) and ( d − r ) × ( d − r ) respectively forsome 0 ≤ v < r and where ˜ A i, ˜ A i, and ˜ A i, O ˜ A i, O ˜ A i, ˜ A i, are r × v and d × r matrices of full column ranks v and r respectively. Clearly, U iscompatible with system A since it is unimodular.From now on, we assume that the leading coefficient A i, of (14) is in form (15). Inparticular, we require the column module of A i, and the column module of B as givenabove to be free. We then partition A i, in accordance with A i, and set G A i ( λ ) := A i, O A i, A i, O A i, A i, A i, A i, + λI d − r . (16)The polynomial det( G A i ( λ )) vanishes identically in λ if and only if θ i ( λ ) given by (10)does. In fact, let D ( x i ) = Diag( x i I r , I d − r ). Then we can write x i − A i ( x ) = N ( x ) D − ( x i )where N ( x ) ∈ R d × d , and set D = D ( x i = 0), N = N ( x i = 0). Then we have22et( G A i ( λ )) = det( N + λD ) = det( N + λD ) | x i =0 = (det( Ax i + λI d ) det( D )) | x i =0 = (det( A i, x i + A i, + λI d ) x ir ) | x i =0 = θ i ( λ ) . Moreover, G A i ( λ ) has an additional important application within the construction of adesired transformation as we show in the following proposition: Proposition 24.
Suppose that m ( A i ) > G A i ( λ )) is identical to zero. If therow module of G A i ( λ = 0) is free then there exists a transformation Q (¯ x i ) in GL d (R ¯ x i )with det( Q ) = ±
1, compatible with system A , such that the matrix G ˜ A i ( λ )) has theform G ˜ A i ( λ ) = A i, O U U A i, O U U V V W + λI d − r − ̺ W M O M W + λI ̺ , (17)where 0 ≤ ̺ ≤ d − r , and rank A i, U A i, U M M = rank A i, U A i, U , (18)rank A i, U A i, U < r. (19) Proof.
If the row module of G A i ( λ = 0) is free then the transformation Q (¯ x i ) can beconstructed as in the proof of [1, Proposition 3] for the bivariate case. ✷ We remark that in the particular case of v = 0, (15) is given by˜ A i, (¯ x i ) = O r O ˜ A i, O d − r , with rank( ˜ A i, ) = r. Consequently, it can be easily verified that (17) is given by G ˜ A i ( λ ) = O r U V W + λI d − r , and ̺ = 0 . roposition 25. If m ( A i ) > G A i ( λ )) ≡ A i of A in (2) is reducible and reduction can becarried out with the shearing F = S ( x i ) G where S ( x i ) = Diag( x i I r , I d − r − ̺ , x i I ̺ ) if ̺ = 0 S ( x i ) = Diag( x i I r , I d − r ) otherwise.Furthermore, this shearing is compatible with system A . Proof.
Given system (2). For any j ∈ { , . . . , n } we partition A j according to (17) A j = A j A j A j A j A j A j A j A j A j A j A j A j A j A j A j A j , ≤ j ≤ n, where A j , A j , A j , A j are square matrices of dimensions v, r − v, d − r − ̺, and ̺ respectively. It is easy to verify that the equivalent system S [ A ] ≡ ˜ A given by (4) admitsthe form ˜ A i = A i A i x − i A i A i A i A i x − i A i A i x i A i x i A i A i x i A i A i A i x − i A i A i − x p i i Diag( I r , O d − r − ̺ , I ̺ )˜ A j = A j A j x − j A j A j A j A j x − i A j A j x j A j x j A j A j x j A j A j A j x − j A j A j , ≤ j = i ≤ n. Hence, the leading matrix coefficient of the equivalent i th -component is given by˜ A i, (¯ x i ) = A i, O U OA i, O U OO O O OM O M O A i, ) < r since (18) and (19) are satisfied.It remains to prove the compatibility of S ( x i ) with the system (2), in particular, that thenormal crossings are preserved. It suffices to prove that the entries of A j , ≤ j = i ≤ n, which will be multiplied by x − i upon applying S ( x i ), namely, the entries of A j , A j , and A j are zero matrices modulo x i otherwise poles in x i will be introduced. This can berestated as requiring A j ( x i = 0) , A j ( x i = 0) , and A j ( x i = 0) to be zero submatrices.This requirement is always satisfied due to the integrability condition and the resultingequality, obtained by setting x i = 0, which we restate here x p j +1 j ∂∂x j A i, = A j ( x i = 0) A i, − A i, A j ( x i = 0) , ≤ j = i ≤ n. (20)This equality induces a structure of A j ( x i = 0) which depends on that of A i, . Since G A i ( λ ) is as in (17), then, before applying the shearing transformation, A i, (¯ x i ) hasthe following form (21) and A j ( x i = 0) can be partitioned accordingly. So we have for1 ≤ j = i ≤ nA i, (¯ x i ) = A i, O O OA i, O ( r − v )( r − v ) O OV V O ( d − r − ̺ )( d − r − ̺ ) OM O O O ̺̺ , (21) A j ( x i = 0) = A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ) ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) A j ( x i = 0) . (22)Inserting (21) and (22) in (20), one can obtain the desired results by equating the entriesof (20). More explicitly, upon investigating the entries in (Column 3), (Rows 1 and 2,Column 2), and (Row 4, Column 2), we observe the following respectively: • We have that A i, OA i, OV V M O · A j ( x i = 0) A j ( x i = 0) = O n,n − r − ̺ . The former matrix is of full rank r by construction thus A j ( x i = 0) A j ( x i = 0) is a zeromatrix. 25 We also get A i, A i, · A j ( x i = 0) = O r,r − v . The former is of full rank v by construction thus A j ( x i = 0) is a zero matrix. • Finally, A j ( x i = 0) · V − M · A j ( x i = 0) = O ̺, ( r − v ) . Since A j ( x i = 0) is null and V is of full column rank r − v by construction then A j ( x i = 0) is a zero matrix aswell.This completes the proof. ✷ We can thus establish the following theorem:
Theorem 26.
Consider a Pfaffian system [ A ] and suppose that m ( A i ) > i ∈ { , . . . , n } . If θ i ( λ ) given by (10) vanishes for some index i ∈ { , . . . , n } ,then under the conditions required to attain (15) and Proposition 24, we can construct acompatible transformation T ∈ GL d (R L ) which reduces m ( A i ) (and consequently m ( A )).In this case, T can be chosen to be a product of transformations in GL d (R ¯ x i ) andpolynomial transformations of the form Diag( x β i , . . . , x β d i ) where β , . . . , β d are non-negative integers. Proof.
Under the required conditions, we can assume that A i, has the form (15). Let G A i ( λ ) be given by (16). Then det( G A i ( λ )) vanishes identically in λ if and only if θ i ( λ )does. Then the system S [ Q [ A ]] where S, Q are as in Propositions 24 and 25 respectively,has the desired property. ✷ For a given index i ∈ { , . . . , n } , the algebraic rank of the leading matrix coefficientcan be decreased as long as θ i ( λ ) vanishes identically in λ . In case the leading matrixcoefficient eventually reduces to a zero matrix, the Poincar´e rank drops at least by one.This process can be repeated until the Moser rank of system [ A ] equals to its Moserinvariant. Due to the compatibility of T in Theorem 18, rank reduction can be appliedto any of the components of [ A ] without altering the Moser rank of the others. Hence,by Corollary 16, the true Poincar´e rank of system [ A ] can be attained by a successiveapplication of the rank reduction to each of its components.Finally, we remark that the conditions of Theorem 26 are always satisfied in thebivariate case n = 2 of arbitrary dimension. Example 27.
Consider the completely integrable Pfaffian system with normal crossingsgiven by x ∂∂x F = A F = ( x x x + 1)( x − x ( x − x x (1 − x + x x x − x x x ) x x x (1 − x ) ! F,x ∂∂x F = A F = (2 + 3 x )( x x x + 1) x (2 + 3 x ) − x x (3 x x x + 2 x x x + x + 3 x + 2) − x x x (2 + 3 x ) ! F,x ∂∂x F = A F = − x x ! F.
26 fundamental matrix of formal solutions is given by:Φ( x , x , x ) x C x C x C e q ( x − /s ) e q ( x − /s ) e q ( x − /s ) . (23)If we only seek to compute the exponential parts q , q , q in (23), then from theassociated ODS, we compute: s = 1 and q ( x ) = x ,s = 1 and q ( x ) = − − x x ,s = 1 and q ( x ) = 0 . Furthermore, if we wish to compute a fundamental matrix of formal solutions, thenfollowing the steps of our formal reduction algorithm, we look at the leading coefficientsof each of the three components of the given system. If one of these coefficients has twodistinct eigenvalues, the we can apply the splitting lemma (Theorem 11). Indeed, since A , ( x , x ) has this property, we can compute such a transformation F = T G up to anyorder. In particular, up to order 10, we compute: T = x x x − x x x + x x x − x − x x , which yields the following diagonalized system (up to order 10): x ∂∂x G = A G = x − x x x f ( x ) ! G,x ∂∂x G = A G = (2 + 3 x ) 00 x x x f ( x ) ! G,x ∂∂x G = A G = x x x f ( x ) ! G, with f ( x ) := ( x x x − x x x + 1). Hence, the system can be uncoupled into twosubsystems of linear scalar equations and integrated to construct G , and consequently, F .An example of the reduction process for a system in two variables was given in Exam-ple 1. However, examples of dimensions two and three do not cover the richness of thetechniques presented. So, to illustrate the full process, we treat an example of dimensionsix. Due to the size of the system and the number of necessary computation steps, weare not able to include it directly in this paper. It is available in several formats at
6. Formal Reduction Algorithm
We now give the full algorithm in pseudo-code and we refer to more detailed descrip-tions within the article whenever necessary.27 emark 28.
Throughout the article, we adopted the field of complex numbers C as thebase field for the simplicity of the presentation. However, any computable commutativefield K with Q ⊆ K ⊆ C can be considered instead. In this case, the restrictions on theextensions of the base field discussed in [6] apply as well and are taken into considerationwithin our Maple implementation.Given system [ A ], we discuss the eigenvalues of the leading matrix coefficients A i, , i ∈{ , . . . , n } , of its n components. If for all of these components uncoupling is unattainable,then we fix i ∈ { , . . . , n } and proceed to compute the exponential order ω ( A i ) from theassociated ODS. Suppose that ω ( A i ) = ℓm with ℓ, m coprime positive integers. One canthen set t = x i /m (re-adjustment of the independent variable), and perform again rankreduction to get an equivalent system whose i th component has Poincar´e rank equalto ℓ and leading matrix coefficient with at least m distinct eigenvalues. Consequently,block-diagonalization can be re-applied to uncouple the i th -component. By Section 3.1,this uncoupling results in an uncoupling for the whole system. As mentioned before, thisprocedure can be repeated until we attain either a scalar system, i.e. a system whose n components are scalar equations, or a system whose Poincar´e rank is given by (0 , . . . , lgorithm 1: fmfs pfaff( p, A )Input: p = ( p , . . . , p n ) , A ( x ) = ( A , . . . , A n ) of (2). Output:
A fundamental matrix of formal solutions (6). { C i } ≤ i ≤ n ← O n ; { Q i } ≤ i ≤ n ← O d ; Φ ← I d WHILE d = 1 or p i > i ∈ { , . . . , n } DOIF A i, has at least two distinct eigenvaluesSplit system [ A ] as in Section 3.1; Update Φ Fmfs pfaff ( p, ˜ A ); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n Fmfs pfaff ( p, ˜ A ); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n ELSE IF A i, has one non-zero eigenvalueUpdate Q i from the eigenvalues of A i, A ( x ) ← Follow Section 3.2 ( A i, is now nilpotent) fmfs pfaff ( p, ˜ A ( x )); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n ELSE
Apply rank reduction of Section 5; Update Φ; Update p ; Update A i, IF p i > A i, has at least two distinct eigenvaluesSplit system as in Section 3.1; Update Φ fmfs pfaff ( p, ˜ A ( x )); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n fmfs pfaff ( p, ˜ A ( x )); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n ELSE IF A i, has one non-zero eigenvalueUpdate Q i from the eigenvalues of A i, A ( x ) ← Follow Section 3.2; ( A i, is now nilpotent) Fmfs pfaff ( p, ˜ A ( x )); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n ELSE
Follow Section 4 to compute ω ( A i ) = ℓm x i ← x im Apply rank reduction of Section 5Update Φ; Update p ( p i ← ℓ ); Update A i, Update Q i from the eigenvalues of A i, A ( x ) ← Follow Section 3.2; ( A i, is now nilpotent) fmfs pfaff ( p, A ( x )); Update Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n END IFEND IFEND WHILERETURN p , A , Φ, { C i } ≤ i ≤ n , { Q i } ≤ i ≤ n .29 lgorithm 2: rankReduce( p, A )Input: p . . . , p n , A . . . . , A n of (2). Output: T ( x ) ∈ GL d (R L ) and an irreducible equivalent system { T ( A ) } whose Poincar´e rank is its true Poincar´e rankand the rank of its leading coefficient matrices is minimal( µ ( T ( A )) = m ( T ( A ))) T ← I d FOR every i from 1 to n DO T i ← I d ; p i ← Poincar´e rank of A i U ( ¯ x i ) ← yields the form (15) IF U ( ¯ x i ) cannot be determined RETURN ERROR “Column module not free.”
END IF A i ← U − A i U ; T i ← T i U WHILE det( G A i ( λ )) = 0 and p i > DO Q ( ¯ x i ) , ̺ ← Proposition 24 IF Q ( ¯ x i ) cannot be determined RETURN ERROR “Row module not free.”
END IF S ( x i ) ← Proposition 25 P ← QS ; T i ← T i PA i ← P − A i P − x p i i S − ∂S∂x i p i ← Poincar´e rank of A i U ( ¯ x i ) ← yields the form (15) A i ← U − A i U ; T i ← T i U END WHILEFOR every j = i from 1 to n DO A j ← T − i A j T i − x p j +1 j T − i ∂∂x j T i END FOR T ← T T i END FORRETURN ( T, p , . . . , p n , A , . . . , A n ). In the case of univariate systems, Levelt’s investigations of the existence of stationarysequences of free lattices lead to an algorithm which reduces the Poincar´e rank to itsminimal integer value [26]. This algorithm was then generalized to the bivariate caseby the first author of this paper et al. in [15]. The theoretical basis of this algorithmdiffers substantially from the algorithm given herein based on Moser’s criterion. Thefinal result of both approaches however, i.e. the algorithm itself, is based on applyingcolumn reductions and shearing transformations in both algorithms, though in a differentmanner. In fact, the algorithms coincide for the particular case of ̺ = 0. The limitationin [15] within the generalization to the multivariate case is in guaranteeing the freenessconditions for the leading matrix coefficient A i, as stated in Section 5.5. The additionalcondition of the freeness of the row module as in Proposition 24 is not required. Sincethe linear algebra problem is resolved in Section 5.4, this results in Algorithm 3.Although both algorithms have an identical cost [25, pp 108], experimental results forthe univariate case and certain bivariate systems (singularly-perturbed linear differential30ystems) suggest that the lattice-based algorithm complicates dramatically the coeffi-cients of the system under reduction, even if Moser’s criterion is adjoined to avoid someunnecessary computations [2, Section 4.3]. Hence, Algorithm 2 can be used as long asthe required freeness conditions hold. Nevertheless, if the freeness of the row module of G A i ( λ = 0) is not satisfied, then Algorithm 3 can be used as long as the column modulesof A i, and B is free. There remains however, the question on the equivalence of theseconditions. Algorithm 3: rankReduce alt( p, A )Input: p . . . , p n , A . . . . , A n of (2). Output: T ( x ) ∈ GL d (R L ) and an irreducible equivalent system { T ( A ) } whose Poincar´e rank is its true Poincar´e rankand the rank of its leading coefficient matrices is minimal( µ ( T ( A )) = m ( T ( A ))) T ← I d FOR every i from 1 to n DO T i ← I d ; p i ← Poincar´e rank of A i WHILE j < d − p i > DO U ( ¯ x i ) ← yields the form (15) IF U ( ¯ x i ) cannot be determined RETURN ERROR “Column module not free.”
END IF r = rank( A i, ) S ( x i ) ← Proposition 25 with ̺ = 0 (i.e. S ( x i ) ← Diag( x i I r , I d − r )) P ← U SA i ← P − A i P − x p i i S − ∂S∂x i ˜ p i ← Poincar´e rank of A i IF ˜ p i < p i THEN j ← ELSE j ← j + 1 END IF p i ← ˜ p i T i ← T i P END WHILEFOR every j = i from 1 to n DO A j ← T − i A j T i − x p j +1 j T − i ∂T i ∂x j END FOR T ← T T i END FORRETURN ( T, A i { ≤ i ≤ n } ).
7. Conclusion
In this article, we studied completely integrable Pfaffian systems with normal crossingsin several variables. We showed that one can associate a set of univariate linear singulardifferential systems from which the formal invariants of the former can be retrieved. This31educes computations to computations over a univariate field via
Isolde or Lindalg , andlimits the numbers of coefficients necessary for the computations. We then complementedour work with a rank reduction algorithm based on generalizing Moser’s criterion andthe algorithm given by Barkatou in [7]. The former is applicable to any bivariate system.However, for multivariate systems, it demands that several explicitly described conditionsare met.One field of investigation is the possibility of weakening the conditions required inthe multivariate setting for the rank reduction. Another one is the adaptation of thetechniques developed herein for rank reduction to generalize the notion of simple systems(see [11] for m = 1). This notion, in the univariate case, gives another approach toconstruct a basis of the space of regular solutions [12], and the obstacles encountered aresimilar to those in rank reduction. An additional field is to study closed form solutions [10]in the light of associated ODS introduced in Section 4.In future work we aim to thoroughly study the theoretical complexity of the proposedalgorithm as well as give detailed information of the effects on the truncation of the inputsystem during the computation. This work is not straightforward and even in the case ofregular systems, it is not studied in the existing work (i.e. [25, Chapter 3] and referencestherein).Systems arising from applications do not necessarily or directly fall into the classof completely integrable Pfaffian systems with normal crossings. Investigations in moregeneral classes can be found in [31, 24, 33] and references therein.An additional field of investigation is the formal reduction in the difference case usingthe approaches proposed herein. Praagman established in [32] a formal decomposition of m commuting partial linear difference operators. This study was intended as an analogto that established by Levelt, van den Essen, G´erard, Charri`ere, Deligne, and others [17,18, 19, 21]. Acknowledgments
We are grateful to the anonymous referees whose comments and suggestions improvedthe readability of this paper. The second author would also like to thank the TechnischeUniversit¨at Wien, where he is employed by the time of submission and supported by theERC Starting Grant 2014 SYMCAR 639270. He furthermore would like to thank MatteoGallet for valuable discussions. The third author would like to thank the University ofLimoges since a part of this work was done during the period of her doctoral studies atDMI.
References [1] H. Abbas, M. Barkatou, and S.S. Maddah. Formal Solutions of a Class of PfaffianSystems in Two Variables. In
Proceedings of the International Symp. on Symb. andAlgeb. Comp. , ACM Press, pp 312-319, 2014.[2] H. Abbas, M. A. Barkatou, and S.S. Maddah. On the Reduction of Singularly-Perturbed Linear Differential Systems. In
Proceedings of the 39th International Sym-posium on Symbolic and Algebraic Computation , ACM Press, pp 320-327, 2014.[3] A. Aparicio Monforte and M. Kauers, Formal Laurent Series in Several Variables,
Expositiones Mathematicae , 31(4):350–367, December 2013,324] D. G. Babbitt and V. S. Varadarajan. Formal Reduction of Meromorphic DifferentialEquations: A Group Theoretic View. In
Pacific Journal of Mathematica , 109 (1),pp 1- 80, 1983.[5] W. Balser.
Formal Power Series and Linear Systems of Meromorphic Ordinary Dif-ferential Equations . Springer, 2000.[6] M. Barkatou. An algorithm to compute the exponential part of a formal fundamentalmatrix solution of a linear differential system.
Journal of App. Alg. in Eng. Comm.and Comp. , 8(1):1-23, 1997.[7] M. Barkatou. A Rational Version of Moser’s Algorithm. In
Proceedings of the Inter-national Symposium on Symbolic and Algebraic Computation , pages 297-302. ACMPress, July 1995.[8] M. Barkatou. Factoring systems of linear functional equations using eigenrings. InI. S. Kotsireas and E. V. Zima, editors,
Latest Advances in Symbolic Algorithms , pp22–42. World Scientific, 2007.[9] M. Barkatou. Symbolic methods for solving systems of linear ordinary differentialequations. Tutorial at the
International Symposium on Symbolic and Algebraic Com-putation , pages 7-8. ACM Press, July 2010.[10] M. A. Barkatou, T. Cluzeau, C. El Bacha, and J. A. Weil. Computing closed formsolutions of integrable connections. In
Proceedings of the 37th International Sympo-sium on Symbolic and Algebraic Computation , pp 43-50, ACM press, 2012.[11] M. Barkatou and C. El Bacha. On k-simple forms of first-order linear differentialsystems and their computation.
Journal of Symb. Comp. ,54:36-58, 2013.[12] M. Barkatou and E. Pfluegel. An algorithm computing the regular formal solutionsof a system of linear differential equations.
Journal of Symbolic Computation , 28(4),569-587, 1999.[13] M. Barkatou and E. Pfluegel. Computing super-irreducible forms of systems of lin-ear differential equations via moser-reduction: a new approach. In
Proceedings ofthe 2007 international symposium on Symbolic and algebraic computation , pp. 1-8.ACM, 2007.[14] M. Barkatou and E. Pfluegel. ISOLDE, Integration of Systems of Ordinary LinearDifferential Equations. Available at: http://isolde.sourceforge.net/[15] M. Barkatou and N. LeRoux. Rank Reduction of a class of Pfaffian Systems in TwoVariables. In
Proceedings of the International Symposium on Symbolic and AlgebraicComputation , pages 204-211. ACM Press, 2006.[16] R. Broucke. On Pfaff’s equations of motion in dynamics; Applications to SatelliteTheory.
Journal of Celestial Mechanics , 18(3): 207-222, 1978.[17] H. Charri`ere. Triangulation Formelle de certains Syst`emes de Pfaff Compl`etementInt´egrables et Application `a l’etude C ∞ des Syst`emes non Lin´eaires. Annali dellaScuola Normale Superiore di Pisa-Classe di Scienze , 7(4): 625 - 714, 1980.[18] H. Charri`ere and R. G´erard. Formal Reduction of Integrable Linear Connexion hav-ing a certain kind of Irregular Singularities.
Analysis , 1:85-115, 1981.[19] P. Deligne. Equations Diff´erentielles `a Points Singuliers R´eguliers.
Lecture Notes InMathematics , Volume 163, Springer-Verlag,1970.[20] D. Delphenich. The role of integrability in a large class of physical systems. arXivpreprint avilable at: arXiv:1210.4976 (2012).[21] A. van den Essen. Regular Singularities along Normal Crossings. In Gerard Ramis,editor,
Syst`emes de Pfaff et Equations Diff´erentielles dans le Champ Complexe ,volume 712 of
Lecture Notes in Mathematics , pages 88-130. Springer-Verlag, 1979.3322] A. van den Essen and A.H.M. Levelt. Irregular Singularities in Several Variables.
Memoirs of AMS , 40(270), 1982.[23] R. G´erard and Y. Sibuya. Etude de certains syst`emes de Pfaff avec singularit´es. InGerard Ramis, editor,
Systemes de Pfaff et Equations Diff´erentielles dans le ChampComplexe , volume 712 of
Lecture Notes in Mathematics , pages 131-288. Springer-Verlag, 1979.[24] H. Hashiguchi, Y. Numata, N. Takayama, and A. Takemura. The Holonomic Gradi-ent Method for the Distribution Function of the Largest Root of a Wishart Matrix.
Journal of Multivariate Analysis , volume 117, pp 296 - 312, 2013.[25] N. LeRoux. Solutions formelles d’´equations aux d´eriv´ees partielles. Ph.D. Thesis.University of Limoges. 2006.[26] A.H.M. Levelt. Stabilizing Differential Operators: a method for Computing Invari-ants at Irregular Singularities.
Differential Equations and Computer Algebra, M.Singer (ed.) , pages 181-228, 1991.[27] D. A. Lutz and R. Schaefke. On the Identification and Stability of Formal Invariantsfor Singular Differential Equations. In
Linear Algebra and its Applications , 72, pp1- 46, 1985.[28] S. S. Maddah. Mathemagix Lindalg Package for Symbolic Resolution of Linear Sys-tems of Differential Equations. Available at: http : [29] H. Matsumura. Commutative Ring Theory . Cambridge Studies in Advanced Math-ematics, Cambridge University Press, ISBN 9780521367646, 1989. Available at: http : //books.google.de/books ? id = yJwN rABugDEC [30] J. Moser. The Order of a Singularity in Fuchs’ Theory. Mathematische Zeitschrift ,72:379- 398, 1960.[31] D. Novikov, S. Yakovenko. Lectures on meromorphic flat connections.
Normal forms,bifurcations and finiteness problems in diff-l eqn. , pages 387-430, 2004.[32] C. Praagman. Formal Decomposition of n Commuting Partial Linear DifferenceOperators.
Rijksuniversiteit Groningen, Mathematisch Instituut , 1983.[33] K. Takano. Local Equivalence of Linear Pfaffian Systems with Simple Poles.
Funk-cialaj Ekvacioj , 24, 373-382, 1981.[34] W. Wasow.