Calcul formel : Mode d'emploi - Exemples en Maple
Claude Gomez, Bruno Salvy, Paul Zimmermann
Masson, 1995
Chapitre VIII, section 3.8, exercice 5, page 225.
Philippe.Dumas@inria.fr
http://algo.inria.fr/dumas/Maple/
|
|
Idée et amorce de preuve | Méthode de Kummer
Idée et amorce de preuve. Une idée simple, voire simpliste, vient à l'esprit. Nous développons le terme général de la série entière et nous sommons chaque terme. Pour tester cela, calculons un petit développement.
> u:=(1+1/n)^(n^2);
> dau:=asympt((1+1/n)^(n^2),n,4);
> daterm:=expand(dau*x^n/n!);
Ensuite sommons.
> map(sum,daterm,n=1..infinity);
La réussite est limitée, mais elle nous fait espérer un équivalent en .
Pour simplifier posons . Le développement du terme général de la série nous pousse à considérer les séries.
> Sum(1/n^k*z^n/n!,n=1..infinity);
Nous voulons majorer les séries de ce type pour prouver que les termes du développement ont des contributions de plus en plus faibles. Dans une telle série les termes vont d'abord en croissant puis en décroissant. Le changement a lieu quand le rapport de deux termes consécutifs vaut environ l'unité. Or ce rapport a pour équivalent . Nous allons couper la série en deux en tenant compte de ceci. Notons la partie entière de . La somme se récrit
en coupant au rang .
Dans la première somme, le terme le plus grand est obtenu pour et il y a . termes. Pour la seconde somme le terme le plus grand est obtenu pour et les termes suivants se majorent par le produit de ce premier terme et d'une suite géométrique de raison . On a ainsi un majorant simple que l'on peut évaluer asymptotiquement.
> asympt(2*z*z^z/z!+z^(2*z)/(2*z)!*2/2^k,z,3)/z^k;
On note que la seconde somme est négligeable devant la première. En tout cas nous avons un majorant qui se comporte en . Évidemment pour ce n'est pas une majoration de grande qualité, mais pour elle devient intéressante et suffisante pour nos besoins.
Méthode de Kummer. Nous avons vu que Maple n'arrivait pas à sommer les termes du développement. Ceci n'est pas choquant en soi, mais nous devons le pousser dans le sens qui nous convient. Pour cela nous employons la méthode de Kummer qui consiste à remplacer le terme à sommer par un équivalent que l'on sait sommer. Par exemple, ici on a envie de remplacer dans le second terme du développement le par pour que le système puisse sommer (et fasse apparaître une exponentielle). De même nous voulons que soit récrit ou que devienne . Bien sûr il faut à chaque fois introduire des termes correctifs. On aura par exemple
Plutôt que de bricoler chaque terme, nous écrivons une procédure. Cette procédure prend en entrée une expression rationnelle en une indéterminée à coefficients rationnels, éventuellement augmentée d'un terme d'erreur de la forme où est un entier relatif. Il est possible de renforcer le code pour vérifier que cette spécification est respectée. Nous commençons par repérer l'indéterminée (#1). Ensuite nous tenons compte de l'eventuel terme d'erreur. S'il est présent, il va déterminer la précision du développement. Sinon la précison dépend de la valeur de la variable d'environnement Order. Dans cette phase, nous traitons le terme d'erreur et nous rangeons dans une liste otherterms les autres termes de l'expression. Dans la phase de calcul proprement dite (#3) nous traitons successivement chaque terme et nous engrangeons les termes qui apparaissent dans la variable result. En sortie (#4) nous renvoyons les nouveaux termes, qui sont de la forme demandée, et un terme d'erreur, qui est lui-même approprié à la sommation. Finalement, cette procédure fournit un développement asymptotique dans une échelle différente de l'échelle des puissances, mais qui s'avère tout aussi utile.
> kummer:=proc(expr)
local n,errorterm,otherterms,accuracy,
Errorterm,result,i,j,dj,c,newterm,k;
i:=indets(expr,name); #1
if nops(i)>1 then
ERROR(`kummer expects its argument to be a function
of one variable, but received`,expr)
elif nops(i)=0 then
RETURN(expr)
else
n:=op(i)
fi;
if has(expr,O) then #2
if type(expr,function) and op(0,expr)=O then
errorterm:=expr
else
errorterm:=select(has,expr,O)
fi;
otherterms:=normal(expr-errorterm);
accuracy:=degree(op(1,errorterm),n);
if accuracy>0 then
Errorterm:=mul(n-k,k=0..accuracy-1)
elif accuracy=0 then
Errorterm:=1
else
Errorterm:=1/mul(n+k,k=1..-accuracy)
fi
else
otherterms:=expr;
accuracy:=-Order;
Errorterm:=1/mul(n+k,k=1..-accuracy)
fi;
if type(otherterms,`+`) then
otherterms:=[op(otherterms)]
else
otherterms:=[otherterms]
fi;
result:=0; #3
for i in otherterms do
j:=normal(i,expanded);
dj:=degree(numer(j),n)-degree(denom(j),n);
c:=limit(j/n^dj,n=infinity);
while (j<>0) and (dj>accuracy) do
dj:=degree(numer(j),n)-degree(denom(j),n);
c:=limit(j/n^dj,n=infinity);
if dj>0 then
newterm:=c*mul(n-k,k=0..dj-1)
elif dj=0 then
newterm:=c;
else
newterm:=c/mul(n+k,k=1..-dj);
fi;
j:=normal(j-newterm,expanded);
dj:=degree(numer(j),n)-degree(denom(j),n);
result:=result+newterm;
od
od;
result+O(1)*Errorterm #4
end:
Appliquons la procédure à l'exemple traité. Nous sortons d'abord la partie rationnelle utile. Puis nous lui appliquons la procédure.
> da:=asympt(normal(daterm/exp(1)^n/x^n/exp(-1/2)*n!),n);
> kummer(da);
Nous reconstituons les termes du développement, puis nous sommons.
> map(proc(expr) exp(n)*x^n*exp(-1/2)/n!*expr end,%);
> map(sum,%,n=0..infinity);
Nous simplifions l'écriture du résultat.
> map(combine,%,power,symbolic);
> map(expand,%);
Les termes qui ne comportent pas l'exponentielle sont négligeables. Nous les évacuons et nous obtenons ainsi le développement asymptotique demandé. Avec ce qui précède, il n'est pas difficile de justifier sa validité.
> select(has,%,exp(exp(1)*x));
La procédure kummer permet de traiter les sommes
Si est un polynôme, elle permet d'obtenir la somme sous forme explicite. Si est une fraction rationnelle, elle permet d'obtenir un développement asymptotique de la somme. Sur la méthode de Kummer, on pourra consulter le livre de Konrad Knopp [Knopp90]. La procédure kummer ci-dessus est inspirée d'une procédure de même nom (mais plus simple et plus limitée) qui figure dans [DuGo97, page 405].