Calcul formel : Mode d'emploi - Exemples en Maple
Claude Gomez, Bruno Salvy, Paul Zimmermann
Masson, 1995
Chapitre V, section 2.5, exercice 3, page 137.
Philippe.Dumas@inria.fr
http://algo.inria.fr/dumas/Maple/
|
|
Il est crucial de remarquer que la matrice donnée est un polynôme en une matrice nilpotente N que nous définissons tout d'abord. Comme la taille 3 n'a rien de particulier, nous généralisons un peu la question en introduisant un paramètre d.
> d:=3:
> N:=array(1..d,1..d,sparse):
> for i to d-1 do N[i,i+1]:=1 od:
> print(N);
> M:=evalm(1+2*N+3*N^2);
L'expression 1+2*N+3*N^2 est le début d'une série entière. C'est le caractère nilpotent de N qui fait que cette série est tronquée. Évaluons cette série entière (qu'il ne serait pas mauvais de connaître par coeur).
> sum((n+1)*x^n,n=0..infinity);
Si nous substituons N à x dans la fraction rationnelle obtenue, nous retrouvons bien M. Notez qu'il n'y a pas de problème avec le calcul de l'inverse ; puisque N est nilpotente, sa seule valeur propre est 0 et N-I (noté ici N-1) est inversible.
> evalm((N-1)^(-2));
Du coup les puissances de M s'exprime très simplement par la formule et ceci est valable que k soit entier naturel ou entier négatif, puisque M est inversible. Par exemple le cube de M peut s'obtenir comme suit (et bien sûr par evalm(M^3)).
> evalm((N-1)^(-2*3));
Grâce au facteur 2 qui figure dans l'exposant, la formule a même un sens pour les demi-entiers. Ceci nous fournit une racine carrée de M.
> R[0]:=evalm((N-1)^(-1));
> evalm(R[0]^2);
On pourrait aussi dire que tout ceci fonctionne parce que M est une exponentielle. Nous le vérifions en employant encore le fait qu'on peut substituer une matrice nilpotente dans une série entiére.
> l:=convert(series(ln(1+x),x,d),polynom);
> L:=evalm(subs(x=M-1,l));
> linalg[exponential](L);
Nous avons ainsi obtenu la formule . Nous pouvons ainsi donner un sens à la quantité M^t avec t réel ou même complexe.
Revenons sur le calcul de racine carrée. Le texte demande d'évaluer la racine carrée de M. Il s'agit bien sûr d'un lapsus digiti, car il y a au moins deux racines carrées. Nous avons déjà obtenu une racine carrée. Elle va servir de référence pour les déterminer toutes. Nous voulons résoudre l'équation . Nous procédons au changement d'inconnue , ce qui signifie que nous définissons . Ceci ne pose pas de problème car les matrices considérées sont nécessairement inversibles comme M. L'équation se récrit ou encore (car nous pouvons simplifier par une matrice inversible) . Nous devons discuter cette équation en S. En toute généralité cela serait un peu délicat, mais l'outil de calcul dont nous disposons permet de régler rapidement la question dans le cas proposé en procédant à une résolution explicite.
> S:=matrix(d,d,s):
> print(S);
> A:=evalm(S&*R[0]&*S);
> sys:={seq(seq(A[i,j]=R[0][i,j],j=1..d),i=1..d)}:
> inc:={seq(seq(s(i,j),j=1..d),i=1..d)}:
> sol:=[solve(sys,inc)];
> SOL[1]:=subs(sol[1],eval(S));
> SOL[2]:=subs(sol[2],eval(S));
Il n'y a pas d'autres solutions que les deux solutions évidentes, donc les racines carrées de la matrice proposée sont données par les formules et