Calcul formel : Mode d'emploi - Exemples en Maple

Claude Gomez, Bruno Salvy, Paul Zimmermann

Masson, 1995

Chapitre VIII, section 1.4, exercice 1, page 202.

Philippe.Dumas@inria.fr
http://algo.inria.fr/dumas/Maple/

Page du Projet Algorithmes | Page de Philippe Dumas | Page Maple de Philippe Dumas

Retour en page principale
Table des matières
Index
Maple V.4 worksheet
Maple V.5 worksheet


Nous employons la seconde des deux formules. De plus nous considérons que le calcul du premier terme de l'expression est effectué par evalf. C'est ce que semble signifier l'énoncé.

Le terme général de la série a l'expression suivante.

> u:=1/(2*n+1)^2/binomial(2*n,n);

[Maple Math]

Nous calculons le quotient de deux termes consécutifs de façon à minimiser le nombre d'opérations.

> r:=convert(u/subs(n=n-1,u),GAMMA);

[Maple Math]

> simplify(r,GAMMA);

[Maple Math]

Nous voyons que le passage du terme de rang n-1 au terme de rang n demande sept multiplications ou divisions. Le calcul des n premiers termes et de la somme partielle d'indice n a donc un coût qui est un O(n). Si nous appliquions brutalement la formule, le coût serait un O(n^2).

L'énoncé demande de passer en paramètre le nombre de chiffres désiré. Il est plus dans l'esprit du logiciel d'employer la variable d'environnement Digits pour contrôler la précision du calcul. D'autre part, nous n'effectuons pas de calcul si le nombre de chiffres demandé est faible. Dans ce cas le résultat est tout simplement stocké en mémoire. C'est ce qui est fait pour les constantes connues du système. On pourra s'en convaincre par l'expérience suivante.

> interface(verboseproc=2):
readlib(`evalf/constant/Pi`);
readlib(`evalf/constant/bigPi`);

Si le nombre de chiffres demandé dépasse 50, nous augmentons la précision des calculs pour éviter les erreurs d'arrondis, puis nous sommons les termes de la série, jusqu'à ce qu'ils deviennent négligeables.

> catalan:=proc()
local t,s,n,oldDigits;
if Digits<=50 then
evalf(.915965594177219015054603514932384110774149374281672)
else
oldDigits:=Digits;
Digits:=floor(1.1*Digits);
t:=1;
s:=1;
for n while t>10^(-oldDigits-5) do
t:=evalf(t/2*n*(2*n-1)/((2*n+1)^2));
s:=s+t
od;
evalf(evalf(Pi/8*ln(2+3^(1/2))+3/8*s),oldDigits)
fi
end:

Comparons le résultat fourni avec celui que donne la procédure du système.

> Digits:=200:

> catalan();

[Maple Math]
[Maple Math]
[Maple Math]

> evalf(Catalan);

[Maple Math]
[Maple Math]
[Maple Math]

> Digits:=10:

Tout semble correct.

Le test d'arrêt de la boucle n'a pas été justifié. Il ne suffit pas que le terme d'indice N d'une série soit plus petit que 10^(-D) pour que le reste d'indice N vérifie la même inégalité. Par exemple le terme 1/N^2 est plus petit que 10^(-4) pour N > 100, mais le reste d'indice 101 vaut environ 0.009. Ici les termes décroissent tellement vite vers 0, que le reste d'indice N est plus petit que le terme d'indice N. Pour nous en convaincre nous regardons le comportement du terme général de la série.

> asympt(u,n,3);

[Maple Math]

> combine(%,power,symbolic);

[Maple Math]

Nous avons obtenu un équivalent du terme général de la série. Comme il s'agit de séries convergentes, les restes des deux séries sont donc équivalents. Pour cette nouvelle série, la suite des termes décroît à partir d'un certain rang. Le reste d'indice N est donc plus petit que l'intégrale de N à l'infini de l'équivalent obtenu (et plus grand que l'intégrale de N+1 à l'infini).

> int(sqrt(Pi)/4/4^n/n^(3/2),n=N..infinity);

[Maple Math]

> asympt(%,N,2);

[Maple Math]

Nous voyons que cette intégrale vaut environ 1/2/ln(2) fois le terme d'indice N, ce qui justifie la méthode employée.

> evalf(1/2/ln(2));

[Maple Math]

Retour en page principale