\documentclass[XUPS,XML,SOM,Unicode,francais, NoFloatCountersInSection, NoEqCountersInSection, screen]{cedram}
\usepackage{xups91-03}
\setcounter{tocdepth}{2}
\begin{document}
\frontmatter

\title[Procédés formels et numériques de sommation de séries]{Procédés formels et numériques\\ de sommation de séries\\ solutions d'équations différentielles}

\author[\initial{J.} \lastname{Thomann}]{\firstname{Jean} \lastname{Thomann}}

\address{Centre de Calcul du CNRS,
BP20 Cro, 67037 Strasbourg Cedex, France}

\thanks{Journées X-UPS 1991. Séries divergentes et procédés de resommation. Prépublication du Centre de mathématique de l'École polytechnique, 1991}

\maketitle
\vspace*{-\baselineskip}
\tableofcontents
\mainmatter

\section{Préliminaires}
Une série formelle $\widehat f(x)=a_0+a_1x+a_2x^2\cdots +a_nx^n+\cdots$
est définie si tous les coefficients $a_j \quad (j\in \NN)$ peuvent
être définis par un ou plusieurs algorithmes \emph{numériques},
\emph{formels} ou \emph{formels-numériques}:
\begin{itemize}
\item
\emph{numériques} au sens où chaque $a_j$ peut être évalué
par un calcul numérique classique (par exemple une quadrature),
avec nécessairement une erreur de représentation.

\item
\emph{formels} au sens où, par exemple:
$$ (7+2j)(5+2j)a_j+2(7+2j)a_{j+1}-(4+2j)a_{j+2}=0 \quad (j=0,1,2,
\ldots) $$
et où nous connaissons les conditions initiales: $a_0=1$ et $a_1=5$.
\end{itemize}

Dans ce cas, les $a_j$ sont solutions d'une équation de récurrence
avec des coefficients polynomiaux. Ces polynômes ont des coefficients
qui appartiennent à $\ZZ$, $\QQ$ ou
qui peuvent être des nombres algébriques
(par exemple $i$, $\sqrt 2$, $\ldots$, ou beaucoup plus
compliqués).

Si les valeurs initiales appartiennent à ces ensembles,
nous pouvons calculer \emph{exactement} tous les $a_j$ (même
si leur représentation est symbolique pour les nombres algébriques).

Ces définitions sont d'une importance primordiale parce que les
valeurs
de $a_j$, si elles sont calculées par des méthodes numériques
faisant
appel à l'arithmétique de la virgule flottante, peuvent dévier
très
rapidement de leurs vraies valeurs, pour des raisons d'instabilité
bien connues notamment dans les récurrences à 3 termes.
\begin{itemize}
\item
\emph{formels-numériques}, par exemple, au cas où les conditions
initiales sont exprimées par des nombres en virgule flottante et
où l'équation de récurrence est \emph{formelle} au sens ci-dessus.
\end{itemize}

Cette distinction ou synthèse \emph{formelle-numérique} sera utile
dans
tout calcul concernant les équations différentielles, les solutions
sous formes de séries, d'approximants rationnels, etc.

La règle sera: utiliser les algorithmes formels, donc exacts, le plus
longtemps possible, puis passer à un calcul numérique {\it
contrôlé.}

\section{Somme d'une série convergente}
Connaissant une série formelle $\widehat f(x)$, le but est de lui
associer une
fonction \emph{unique} $f(x)$, somme de la série formelle en $x$.

Si tout va bien, nous pouvons utiliser directement les critères de
convergence (Cauchy, d'Alembert, etc.) et constater la convergence
d'une série dans un disque de rayon non nul, comme dans le cas
trivial:
$$ \widehat f(x)=1-x+x^2-x^3+\cdots +(-1)^nx^n+\cdots $$
convergente dans un disque de rayon $< 1$.

Mais aussitôt, on aimerait connaître la somme $f(x)$ dans le
domaine
le plus grand possible du plan complexe; d'où les méthodes de
prolongement analytique bien connues (Weierstrass, Borel, Mittag-Leffler, approximants rationnels, etc.).

Dans notre cas simple, la fonction rationnelle $\spfrac{1}{1+x}$ nous
donne immédiatement une fonction holomorphe dans tout le plan
complexe privé du point singulier $-1$.

D'où l'idée d'effectuer un prolongement analytique par des fonctions
rationnelles avec l'aide d'algorithmes qui sont plus
faciles à mettre en
oeuvre que les autres méthodes mentionnées et qui peuvent être
entièrement formels.

La série formelle
$$\widehat f(x)=1-{3\over 4}x+{39\over 32}x^2-{267\over 128}x^3
+{7563\over 2048}x^4-{54789\over 8192}x^5+\cdots $$
semble plus difficile à analyser, à sommer et à prolonger que la
précédente.

Par contre, si nous savons que cette série est solution formelle de
l'équation fonctionnelle:\enlargethispage{\baselineskip}
$$ (1+2x)f^2-1-{1\over 2}x=0,$$
nous pouvons constater qu'une singularité se situe en $-{1/ 2}$,
et
que $\widehat f$ est convergent dans un disque de rayon $< {1/ 2}$

D'autre part, l'approximant de Padé $\lbrack 1/1 \rbrack:
{1+{7/ 8}x\over 1+{13/ 8}x}$ (approximant rationnel admettant
une série de Taylor coïncidant avec $\widehat f$ jusqu'à l'ordre
$1+1=2$
inclus et calculé directement à partir de 3 termes de la série
$\widehat f$)
représente la somme $f(x)=\sqrt { 1+{1/ 2}x \over 1+2x}$
avec une précision déjà remarquable, de sorte que $f(1)=0.707...$
et $\lbrack 1/1\rbrack(1)=0.714...$, ainsi que
$f(\infty)=0.5$
et $\lbrack 1/1\rbrack (\infty)=0.54$. (Voir \cite{BG}.)

Dans notre cas simple, il est facile de voir que:
$$ \widehat f(x)=1-x+x^2\cdots +(-1)^nx^n+\cdots $$
est une série
formelle, solution
formelle de l'équation différentielle:
$$ (x+1)y'+y=0.$$
Cette série formelle a été obtenue en reportant la série
$$ y(x)=a_0+a_1x+a_2x^2\cdots +a_nx^n+\cdots $$
ainsi que sa dérivée formelle dans l'équation différentielle et
en
identifiant les coefficients $a_j$ des mêmes termes en $x^j$.

D'où une équation de récurrence:
$$ a_{j+1}+a_j=0\qquad (j=0,1,2...)$$
associée à l'équation différentielle génératrice et
permettant le
calcul exact de tous les $a_j$ à partir de $a_0=1$.

Cette dualité \emph{équation différentielle-équation de
récurrence}
(utilisée par Frobenius et
Pincherle) a été développée et est un élément de
base des logiciels de calcul formel des solutions d'équations
différentielles linéaires (par exemple le code DESIR développé
au
LMC de Grenoble \cite{2,10}).\enlargethispage{\baselineskip}

Il est important de remarquer que nous pouvons lire la singularité
$-1$ sur l'équation différentielle.
En effet, l'équation caractéristique $x+1=0$ a comme racine $-1$.

\section{Somme d'une série divergente}
Nous pouvons utiliser le même procédé sur l'équation
différentielle:
$$ x^2y'+y=x$$
pour trouver la série formelle solution:
$$x\widehat f(x)=x\lbrack 1-1!x+2!x^2+\cdots +(-1)^nn!x^n+\cdots\rbrack .$$
Cette fois-ci, la série $\widehat f(x)$ est divergente; ce qui pouvait
être prévu en regardant l'équation différentielle qui
présente une
singularité irré\-gu\-lière à l'origine.

La théorie des développements asymptotiques Gevrey, solutions
d'équations différentielles linéaires \cite{4},
nous permet de dire qu'à
l'intérieur d'un secteur \emph{assez petit} $V$ (d'ouverture $<
\pi /k$), la sommation au \emph{plus petit terme} de $\widehat f(x)$
donne une approximation d'une solution $f(x)$ holomorphe dans $V$
telle que:
$$ \vert f(x)-\sum_{p< L_V\vert x\vert^{-k}}a_px^p\vert
\leq C_V\exp(-K_V\vert x\vert^{-k}),$$
où $L_V$, $C_V$ et $K_V$ sont des constantes liées à $V$.
$k$ vaut 1 dans notre cas particulier.

Cette somme tronquée donne une solution à une fonction
exponentiellement décroissante près.

Dans les cas usuels, notamment les séries solutions d'équations
de fonctions spéciales, ce procédé revient à tronquer la
série juste
avant le plus petit terme $\vert a_nx^n\vert$.

C'est ainsi qu'on trouve la \emph{somme} de la série, c'est à dire
une solution de l'équation différentielle: par exemple en $0.1$:
$f(0.1)=0.91528$ (la valeur exacte de cette fonction connue étant
$0.91563$).

Au sens de l'analyse numérique classique, ce procédé permet
d'avoir
une très bonne approximation dans un petit voisinage de~$O$,
puisque l'erreur décroît exponentiellement si on tend vers~$O$
sur n'importe quel rayon issu de l'origine.

Si on entreprend une analyse plus fine, on peut dans certains cas
représenter la solution d'une équation différentielle linéaire,
dont le développement de Taylor est une série formelle $\widehat f(x)$,
sous la forme:
$$ f(x)={k\over x^k}\int_d\varphi (t)\exp(-{t^k/ x^k}),
t^{k-1}dt$$
où
$k$ est un entier positif.

$\varphi(t)$ est le prolongement analytique sur une droite $d$
issue de $0$ de la somme de la série formelle
\[
\widehat\varphi(t)=
\sum_{n\ge 0}{a_n\over\Gamma(1+{n/ k})}t^n
\]
convergente
dans un disque de rayon non nul.

Cette représentation est vraie si $\widehat f(x)$ est \emph{$k$-sommable}
d'après la théorie développée par J.-P.\,Ramis \cite{4}.

Si $k$ est entier (sinon on peut s'y ramener par une ramification,
si l'équation est à coefficients polynomiaux), cette intégrale
existe dans toute direction sauf un nombre fini de \emph{directions
singulières} et $f(x)$ est l'unique solution de l'équation
différentielle dans un secteur d'ouverture $>{\pi/ k}$, dont le
développement de Taylor est $\widehat f(x)$.

Toutes les séries solutions d'équations différentielles du
premier et
du second
ordre (fonctions spéciales, etc.) homogènes sont
\emph{$k$-som\-ma\-bles} et bien d'autres encore...

La valeur de $k$ peut être lue directement à partir des coefficients de
l'équation
différentielle avec l'aide du polygone de Newton.

\begin{definition*}[polygone N-R-M (Newton, Ramis, Malgrange) en $x=0$]
Considérons un opérateur différentiel
\[
{\sum_{i=0}^{n}\sum_{j=0}^{\infty}p_{i,j}^{}x^j \left(\frac{d}{dx}\right)^i}.
\]
Si $Q^+(u,v)=\lbrace(x,y)\in \RR^2,\; x\leq u,\,y\geq v\rbrace$, est
le second quadrant de~$\RR^2$ translaté en $(u,v)$, on définit:
$$
M^+(L)=\bigcup_{p_{i,j}\not=0}Q^+(i,j-i). $$
Le \emph{polygone} {N-R-M} de $L$
est l'enveloppe convexe inférieure de
$M^+(L)$.
\end{definition*}

Toutes ses pentes finies sont rationnelles. Par exemple, le
polygone N-R-M de l'équation:
$$ 4x^5y''+2x^2y'+y=0 $$
est représenté sur la figure 1.\enlargethispage{\baselineskip}

\begin{figure}[htb]
\begin{center}
\begin{picture}(2.5,4.5)(0,0)

%axe horizontal

\put(0,0){\vector(1,0){2.5}}
\put(2.4,-.3){${}_u$}

%axe vertical

\put(0,0){\vector(0,1){4.5}}
\put(-.3,4.4){${}_v$}

%coordonnees

\put(0,0){\circle*{.15}}
\put(0,-.3){${}_0$}

\put(1,0){\circle*{.03}}
\put(1,-.3){${}_1$}

\put(2,0){\circle*{.03}}
\put(2,-.3){${}_2$}

\put(-.3,1){${}_1$}

\put(-.3,2){${}_2$}

\put(-.3,3){${}_3$}

\put(1,1){\circle*{.15}}

\put(2,3){\circle*{.15}}


\multiput(0,1)(.2,-.2){3}{\circle*{.03}}

\multiput(0,1.5)(.2,-.2){4}{\circle*{.03}}
\multiput(0,2)(.2,-.2){5}{\circle*{.03}}
\multiput(0,2.5)(.2,-.2){6}{\circle*{.03}}
\multiput(0,3)(.2,-.2){7}{\circle*{.03}}
\multiput(0,3.5)(.2,-.2){8}{\circle*{.03}}
\multiput(0,4)(.2,-.2){8}{\circle*{.03}}
\multiput(.2,4.3)(.2,-.2){8}{\circle*{.03}}


\multiput(.6,4.4)(.2,-.2){7}{\circle*{.03}}
\multiput(1,4.5)(.2,-.2){5}{\circle*{.03}}

\thicklines
\put(0,0){\line(1,1){1}}
\put(1,1){\line(1,2){1}}
\put(2,3){\line(0,1){1.5}}

\end{picture}
\end{center}
\caption{}
\end{figure}

Pour l'équation d'Euler $x^2y'+y=x$, la série $\widehat f(x)$ est
\emph{$1$-sommable} parce que le polygone de Newton n'a qu'une seule
pente égale à $1$.

Le fait que la pente soit différente de $0$ nous indique aussi que
la singularité à l'origine est \emph{irrégulière} et que par
conséquent
la série est en général divergente.

La série $\widehat f(x)$ a une somme unique définie par
$$ f(x)={1\over x}\int_d\varphi(t)\exp(-{t/ x})dt$$
où $d$ est toute direction
issue de $O$ sauf la direction $\RR^-$.

Pourquoi la direction singulière est-elle $\RR^-$?
Si nous calculons formellement $\widehat\varphi(t)$, définie
précédemment
et
appelée transformée formelle de Borel, nous obtenons:
$$ 1-t+t^2-\cdots +(-1)^nt^n\cdots ={1\over{1+t}} $$
et nous connaissons son prolongement sur toute droite $d$, sauf
celle qui correspond à $\RR^-$ puisque $\spfrac{1}{1+t}$ a un pôle en
$-1$.

Si $d$ prend toutes les directions non singulières, les fonctions
définies par l'intégrale et holomorphes dans un demi-plan
bissecté par
$d$, se recollent pour donner une seule fonction $f(x)$ définie
dans un secteur d'ouverture $3\pi$.

Dans le demi-plan bissecté par $\RR^-$, nous voyons que nous obtenons
deux déterminations de $f(x)$ en faisant tourner $d$ depuis la position
où $d$ forme un angle
$-\pi+\varepsilon$ avec $\RR^+$ jusqu'à la position où $d$ forme un
angle $+\pi-\varepsilon$ avec $\RR^+$; et par un calcul de
résidu on trouve que la différence des deux déterminations est
$2i\pi\exp(1/x)$.

Nous voyons donc que la notion de disque de convergence est remplacée
par la notion de secteurs de convergence.

La transformée de Borel formelle permet de \emph{transposer}
l'étude de $f(x)$ dans un nouveau plan complexe où les séries
$\widehat\varphi(t)$ sont cette fois-ci convergentes comme précédemment
avec les mêmes problèmes de prolongement analytique.

Dans cette transposition il y a correspondance entre \emph{directions
singulières} de $\widehat f(x)$ et \emph{directions sur
lesquelles se trouvent
les singularités} de $\varphi(t)$.

\section{Un exemple particulier}
Nous allons illustrer cette démarche sur un exemple un peu plus
complexe.

La série formelle
$$\widehat g(x)=1+5x+105/4x^2+525/4x^3
\cdots$$ est l'une des 2 solutions de l'équation différentielle:
$$ 4x^3y''+2(14x^2+2x-1)y'+5(7x+2)y=0.$$
Le polygone de Newton de cette équation, représenté sur la figure 2,
est composé d'un segment de pente nulle et d'un segment de pente 2.
Nous pouvons en déduire que la série formelle $\widehat g(x)$,
correspondant
à la pente nulle, est \emph{$2$-sommable}.

\begin{figure}[htb]
\begin{center}
\begin{picture}(2.5,4)(0,0)

%axe horizontal

\put(0,1){\vector(1,0){2.5}}
\put(2.4,.7){${}_u$}

%axe vertical

\put(0,0){\vector(0,1){4}}
\put(-.3,3.9){${}_v$}

%coordonnees

\put(0,1){\circle*{.03}}
\put(0,.7){${}_0$}

\put(1,1){\circle*{.03}}
\put(1,.7){${}_1$}

\put(2,1){\circle*{.03}}
\put(2,.7){${}_2$}

\put(-.3,2){${}_1$}
\put(0,2){\circle*{.03}}

\put(-.5,0){${}_{-1}$}
\put(0,0){\circle*{.03}}

\put(1,0){\circle*{.15}}

\put(2,2){\circle*{.15}}


\thicklines
\put(0,0){\line(1,0){1}}
\put(1,0){\line(1,2){1}}
\put(2,2){\line(0,1){2}}

\end{picture}
\end{center}
\caption{}
\end{figure}

Les coefficients de la série formelle sont obtenus par l'équation
de récurrence associée à l'équation différentielle
génératrice:
$$ (7+2j)(5+2j)a_j+2(7+2j)a_{j+1}-(4+2j)a_{j+2}=0\quad (j=0,1,2\ldots)
, $$
où $a_0=1$ et $a_1=5$.

Si nous scindons $\widehat g(x)$ en 2 sous-séries: $\widehat\psi_1(u)+x\widehat\psi_2(u)$, où $u=x^2$, $\psi_1(u)$ et
$\psi_2(u)$ sont chacune \emph{$1$-sommable}.

Les
coefficients des transformées de Borel formelles: $\widehat\varphi_1(t)$
et $\widehat\varphi_2(t)$ seront solutions respectivement des équations
de récurrence:
\begin{multline*}
(5+4j)(7+4j)(9+4j)(11+4j)a_j\\
-2(7+4j)(9+4j)(11+4j)(j+1)a_{j+1}\\
+(6+4j)(8+4j)(j+1)(j+2)a_{j+2}=0
\end{multline*}
et
\begin{multline*}
(7+4j)(9+4j)(11+4j)(13+4j)a_j\\
-2(9+4j)(11+4j)(13+4j)(j+1)a_{j+1}\\
+(8+4j)(10+4j)(j+1)(j+2)a_{j+2}=0.
\end{multline*}

Si nous calculons formellement
l'équation différentielle génératrice de la première
équation
nous trouvons:
\begin{multline*}
16(16t^2-8t+1)t^2{d^4U\over dt^4
+(3584t^2-1248t+72)t{d^3U\over dt^3} }\\
+(13920t^2-2904t+48){d^2U\over dt^2}
+(15840t-1386){dU\over dt}
+3465U=0.
\end{multline*}
Nous savons que $\widehat \varphi_1(t)$ est solution de cette équation.

L'équation caractéristique: $16t^2-8t+1=(4t-1)^2=0$ a une racine
double
en $1/4$, qui est donc la seule singularité à distance finie de
l'équation différentielle, donc la seule possible à distance
finie pour
$\varphi_1(t)$.
On peut situer de la même manière la seule singularité
possible
à distance
finie de $\varphi_2(t)$,
en $1/4$.

En analysant formellement cette équation différentielle,
par exem\-ple par le code DESIR, nous
constatons que nous nous trouvons en face d'une singularité
essentielle.

Nous savons donc dans quelles directions nous pouvons effectuer le
prolongement analytique de $\varphi_1(t)$ afin de pouvoir calculer:
$$ \psi_1(u)={1\over u}\int_d \varphi_1(t)\exp(-{t/ u})dt. $$
Ainsi, la direction $d$ pourra être toute direction sauf $\RR^+$ car elle contient~$1/4$.

\section{Calculs effectifs.}
Nous pouvons donc passer au calcul effectif de ces transformées de
Laplace:
\[
f(x)={1/ x}\int_d\varphi(t)\exp(-{t/ x})dt,
\]
où $\varphi(t)$ est la somme d'une série formelle $\widehat\varphi(t)$
convergente dans un disque de rayon non nul.

Si $\widehat f(x)$ est \emph{$k$-sommable}, nous pouvons toujours, comme
dans l'exemple précédent, scinder $\widehat f(x)$ en sous-séries
\emph{$1$-sommables} et nous ramener à ce cas.

Rappelons que les coefficients de la série $\widehat f(x)$, donc de la
série~$\widehat\varphi(t)$ sont des nombres algébriques si $\widehat f(x)$
est solution d'une équation différentielle linéaire à
coefficients
polynomiaux, ces polynômes ayant eux-mêmes des coefficients qui
sont des nombres algébriques.

Nous connaissons l'équation de récurrence exacte que vérifient
les coefficients de $\widehat \varphi(t)$,
et par conséquent l'équation différentielle que
vérifie $\widehat\varphi(t)$.
Donc ses singularités sont connues ainsi que les directions
singulières correspondantes.

De plus la théorie de la \emph{$k$-sommabilité} nous assure de la
croissance au plus exponentielle de $\varphi(t)$ quand $t$
tend vers l'infini dans toute direction non singulière.

\subsection{Décomposition spectrale formelle}
Une première méthode consiste à prolonger $\varphi(t)$ par un
approximant rationnel de type Padé dans les directions non
singulières.

Dans le cas où il n'y a qu'une seule direction singulière
(équations différentielles homogènes du second ordre, fonctions
spéciales, etc.), nous pouvons approcher
$\varphi (t)$ par ${P_M(t)/ Q_N(t)}$ où
$P_M(t)$ et $Q_N(t)$ sont des polynômes respectivement de degré
$M=N+j$, où $j\geq -1$, et de degré $N$.
On choisit $Q_N(t)$ de telle sorte que les racines de $Q_N(t)$ soient
situées sur la coupure $T$ correspondant à la demi-droite issue de
la singularité la plus proche de $O$, d'affixe $a$, et de direction,
la direction singulière.

À la différence des approximants de Padé, pour les approximants de
\emph{type Padé}, on choisit les dénominateurs $Q_N(t)$, puis
on détermine $P_M(t)$ de sorte que le développement de Taylor
de ${P_M(t)/ Q_N(t)}$ coïncide avec la série formelle
jusqu'à
l'ordre $M$ inclus.

On prendra $Q_N(t)=t^Nv_N({1/ t})$ où $v_N$ sera un polynôme
orthogonal (Legendre, Tchebycheff) ayant ses racines situées sur
le segment inverse de $T$.

Dans ce cas, il est démontré \cite{3} que:
$$ {1\over x}\int_d {P_M(t)\over {Q_N(t)}}\exp({-t/ x})dt\quad
\hbox {\rm tend vers}\quad {1\over x}\int_d\varphi (t)\exp({-t/
x})dt,$$
quand $N$ tends vers l'infini et $M=N+j$ où $j\geq -1$.

En décomposant ${P_{M}(t)/ {Q_N(t)}}$ en
éléments simples, nous obtiendrons une expression de la forme
\begin{align*}
f(x)&=\Pol(x)
+\sum_{i=1}^N{A_i\over x}\int_d {1\over {t-t_i}}\,e^{-{t/
x}}dt \\
&=\Pol(x)
+\sum_{i=1}^N-{A_i\over{t_i}}\EXPI(-{x/{t_i}}).
\end{align*}
où $\Pol(x)$ est un polynôme en $x$ et
où $\EXPI$ désigne:
$$ \EXPI(x)={1\over x}\int_d {1\over {1+t}}e^{-{t/ x}}dt,$$
qui est la solution de l'équation d'Euler.

D'ailleurs dans ce cas simple, nous avons bien exprimé directement
$\varphi(t)$ par la fonction rationnelle $\spfrac{1} {1+t}$, qui peut
être considérée comme un approximant de type Padé.

Bien plus, comme nous le verrons en 5.4, ces fonctions $\EXPI(x)$
peuvent être calculées par de vrais approximants de Padé, obtenus
par un algorithme purement formel. Les coefficients ${A_i/ t_i}$
étant des nombres algébriques, cette méthode peut donc être
menée
au bout par des algorithmes formels.

L'évaluation de $\psi_1(u)=\psi_1(x^2)$ pour $x=0.15$ donne le
résultat $0.622338$ en utilisant cette méthode.

\Subsection{Prolongement numérique de la transformée de Borel}
Une autre manière de concrétiser le prolongement analytique de~$\varphi(t)$ est mi-formelle, mi-numérique.
En effet, le calcul formel nous fournit l'équation différentielle
vérifiée par
\[
\varphi(t)=\sum_{n=\mu}^\infty{a_n\over{n!}}t^n.
\]
Les conditions initiales de $\varphi(t)$
en $0$ sont facilement lisibles sur la série elle-même.
Cette équation différentielle est transformée en un
système différentiel d'ordre 1 par un algorithme formel.
À partir de ce système formel, nous pouvons générer un programme
numérique classique, par exemple en Fortran ou en Pascal, qui va
résoudre le système différentiel par un sous-programme de
bibliothèque de résolution par la méthode de Runge-Kutta, puis
l'intégration de
\[
\frac{1}{x}\int_d\varphi(t)\exp(-{t/ x})dt
\]
par un sous-programme de bibliothèque de quadrature optimale par la
méthode de Gauss-Laguerre.

Ces méthodes numériques sont stables
dans le plan complexe, sauf à proximité des singularités de
$\varphi(t)$, qui sont connues.

C'est ainsi que $\psi_1(x^2)$ vaut $0.622150$ au même endroit que
précédemment.

\subsection{Séries de factorielles}
Une troisième méthode numérique
particulièrement efficace consiste à transformer
$\widehat f(x)$ en série de factorielles généralisées (ou séries
de
\emph{facultés} généralisées).
Les \it nombres de Stirling généralisés de $1^{\rm ere}$ espèce
$S_t(n,k)$ \/\rm sont définis par la formule de récurrence:
$$
S_t(n+1,k)=S_t(n,k-1)-t_nS_t(n,k) \qquad (n,k\geq 1)
$$
$S_t(0,0)=1,\quad S_t(n,0)=0$ si $n\not= 0$ et $S_t(n,k)=0$ si
$k\geq n+1$ \hfill \break
où $n,k\in N$ et $\lbrace t_1,t_2,\ldots t_n,\ldots\rbrace$ est une
suite quelconque de nombres complexes, que nous supposerons situés sur
$1/d$.

Toute série formelle:
$$ \widehat f(x)=a_0+a_1x+a_2x^2+\ldots +a_nx^n+\ldots $$
peut être transformée formellement, telle que:
$$
{1\over z}\widehat f({1/ z})=
\sum_{n\geq 0}{b_{n+1}\over{z(z+t_1)\ldots
(z+t_n)}}, $$
où $$
b_{n+1}=\sum_{k=0}^n
(-1)^{k+n}S_t(n,k)a_k, $$ et
où $z=1/x$ (\cite{4}).

On appellera SFG (comme série de factorielle généralisée)
formelle, la série: $$\sum_{n\geq0}{b_{n+1}\over
z(z+t_1)(z+t_2)\ldots (z+t_n)}.$$

\begin{remarque*}
Par récurrence, on peut construire les formules de passage entre SFG
formelles, correspondant à des suites $\lbrace t_1,t_2,\ldots
t_n,\ldots
\rbrace$ différentes.
\end{remarque*}


On a vu
que, si une série $\widehat f(x)$ était 1-sommable dans
une direction~$d$, d'angle $\theta$ avec le demi-axe positif, sa somme
s'exprimait par:
$$
f(x)={1\over x}\int_d\varphi (t)e^{-{t/ x}}dt={1\over
x^*}\int_0^\infty \varphi (te^{i\theta})e^{-{t/ x^*}}dt\quad
\hbox{où}\ x^*=xe^{-i\theta}
$$
ou par:
$$
{1\over z}f({1/ z})=e^{i\theta}\int _0^\infty \varphi(te^{i\theta})
e^{-tz^*}dt \quad x={1/ z}\quand x^*={1/ z^*}
$$ \par
En situant $\lbrace t_1,t_2,\ldots t_n\rbrace$
sur la demi-droite de direction $-\theta$,
on a l'équivalence formelle:
\begin{align*}
{1\over z}
\widehat f({1/ z})&=\sum_{n\geq 0}{b_{n+1}\over {z(z+t_1)
\ldots (z+t_n)}}\\
&=\sum_{n\geq 0}{b^*_{n+1}\over {ze^{i\theta}(ze^{i\theta}+\tau_1)
\ldots (ze^{i\theta}+\tau_n)}},
\end{align*}
où
\begin{align*}
b_{n+1}&=\sum_{k=0}^n (-1)^{k+n} S_t(n,k)a_k,\\
b^*_{n+1}&=
\sum_{k=0}^n (-1)^{k+n} S_{\tau}(n,k)a_k
e^{i(k+1)\theta},
\end{align*}
et $\tau_l=t_l\exp(i\theta)$ $(l=1,2,\ldots, n)$.

Les $\tau_l$ sont alors situés sur la demi-droite $\RR^+$.

En posant:
$$ u_1={a_0\over z},\ u_2={a_1\over z^2},\ldots, u_m={a_{m-1}\over z^m}
\qquad {\rm et} \qquad y=ze^{i\theta},
$$
il vient
$$
{1\over z}\widehat f({1/ z})=\sum_{n\geq 0}
v_{n+1}=\sum_{n\geq 0}{\sum_{k=0}^n
u_{k+1}y^{k+1}\vert S_\tau(n,k)\vert \over {y(y+\tau_1)
\ldots (y+\tau_n)}}.
$$

L'évaluation directe des termes $v_{n+1}$ est impossible à cause
de la
croissance explosive des $S_\tau(n,k)$.
En utilisant la récurrence définissant les nombres de
Stirling généralisés
de 5.3, on peut construire un algorithme récursif où les termes:
$$
v_{n+1}^{(j)}={\tau_{n-1}v_n^{(j)}+yv_n^{(j+1)}\over {y+\tau_n}}
$$
sont calculés par colonnes dans le tableau:
$$\arraycolsep4.5pt
\begin{matrix}
v_1^{(1)}=u_1\cr &v_2^{(1)}={yu_2\over {y+\tau_1}}\cr
v_1^{(2)}=u_2 & &v_3^{(1)}=
{\tau_1v_2^{(1)}+yv_2^{(2)}\over {y+\tau_2}}\cr
&v_2^{(2)}={yu_3\over {y+\tau_1}}
& &v_4^{(1)}={\tau_2v_3^{(1)}+yv_3^{(2)}
\over {y+\tau_3}} \cr
v_1^{(3)}=u_3 &&v_3^{(2)}=
{\tau_1v_2^{(2)}+yv_2^{(3)}\over {y+\tau_2}}&&
\ddots\cr
& v_2^{(3)}={yu_4\over {y+\tau_1}} \cr
v_1^{(4)}=u_4 \cr
\;\vdots
\end{matrix}
$$

Les termes diagonaux: $v_1^{(1)},v_2^{(1)},\ldots, v_n^{(1)},v_{n+1}^{(1)
}\ldots $ sont les termes de la SFG cherchée:
$ v_1,v_2,\ldots, v_n,v_{n+1}$. D'après la théorie de la $k$-sommabilité développée par
J.-P.\,Ramis~\cite{4}, reprenant un résultat de Watson,
nous
avons correspondance entre les résultats suivants: \par
\begin{itemize}
\item
d'une part: $$
f(x)={1\over x}\int _d \varphi
(t)e^{-{t/ x}}dt$$ \par
\item
d'autre part: $$
{1\over z}f({1/ z})=\sum_{n\geq 0}{b
_{n+1}\over{z(z+t_1)\cdots
(z+t_n)}} , $$
où les coefficients $b_{n+1}$
sont définis précédemment, et où la SFG est
uniformément convergente dans le demi-plan $\mathrm{Re}(ze^{i\theta})\geq
\lambda$ ($\lambda> 0$), à \emph{condition que}:
$$ t_\ell=\ell\omega \exp(-i\theta) \quad (\ell=1,\ldots, n), $$
où $\omega \in \RR$, $\omega > 0$, est tel que $\omega >
\omega _0 > 0$. Ici, $\omega_0$ est la valeur critique liée au type de $\widehat f$ et aux
singularités
de $\varphi$.

La valeur de $\omega_0$, qui dépend de $\theta$, peut s'évaluer
simplement
à partir des situations des singularités de $\varphi$
dans le plan complexe.
\end{itemize}

La démonstration, d'après Watson, Nörlund, Nevanlinna, de la convergence uniforme de
$$
\sum_{n\geq0}{b^*_{n+1}\over z^*(z^*+\omega)\ldots (z^*+n\omega)}
\quad \hbox{où } z^*=z\exp(i\theta)$$
utilise la transformation conforme
$s=\exp(-\omega t)$.

L'image réciproque du cercle $\lbrace\vert 1-s\vert =1\rbrace$ est la
courbe C constituée des $t=\omega (u+iv)$, où: $u=-\log(2\cos(v))$,
$v\in{}\rbrack -\pi/2,\pi/2\lbrack$. L'étude de
$\varphi(te^{i\theta})$ dans
l'image réciproque $U$ du disque ouvert \hbox{$\lbrace\vert 1-s\vert< 1
\rbrace $} est ainsi remplacée par l'étude du développement de
Taylor
de $\varphi(-\log(s)e^{i\theta}/\omega)$
au centre de ce disque. Ce changement de
variable permet d'ailleurs d'écrire la transformée de Laplace sous
forme
de transformée de Mellin, puis de série de factorielle convergente.

Tout ceci dépend de la condition d'holomorphie de
$\varphi (te^{i\theta})$ dans $U$.
La valeur maximale $\omega_0$ de $\omega$ est telle que C passe par la
première singularité rencontrée de $\varphi(te^{i\theta})$.

Comme nous pouvons localiser les singularités de $\varphi(t)$ dans
le plan complexe, nous pouvons donc calculer facilement
$\omega_0$ par un algorithme numérique.

\begin{remarque*}
Des expériences numériques montrent que la
convergence
de la SFG est d'autant plus rapide que $\omega$ est proche de $\omega
_0$.
\end{remarque*}


\subsection{Approximation rationnelle directe}
Si nous examinons les séries de factorielles généralisées, nous
constatons que l'approximant de type Padé $P_n(x)/ Q_n(x)$
de la série formelle $\widehat f(x)$, où:
$$ Q_n(x)=\left(x+{1\over t_1}\right)
\left(x+{1\over t_2}\right)\cdots \left(x+{1\over t_n}\right),$$
coïncide formellement avec la somme des $n+1$ premiers termes de la
SFG divisée par $x$ (ou multipliée par $z=1/x$).

Si $\widehat f(x)$ n'a qu'une direction singulière, nous pouvons
démontrer
que, d'une manière plus générale \cite{9},
${1\over x}\int_d \varphi (t)\exp(-t/x)dt$ où $d$ forme avec
$\RR^+$ un angle $\theta$ et avec
$\sum (D)$ un angle $\theta_s$ supérieur à~$\pi /2$,
peut être approché indéfiniment quand $n$ tend vers l'infini,
pour~$x$ fixé, par l'approximant
de type Padé ${P_m(x)/ Q_n(x)}$ de $\widehat f(x)$, où $m=n+j \quad
(j\ge 0)$ et où nous avons choisi le dénominateur:
$$ Q_n(x)=\Bigl (x+{1\over t_1^{(n)}}\Bigr)
\Bigl (x+{1\over t_2^{(n)}}\Bigr)\cdots
\Bigl (x+{1\over t_n^{(n)}}\Bigr). $$
Les $t_j^{(n)}\quad (j=1,2,\ldots n)$
sont cette fois dépendants de $n$:
\[
t_j^{(n)}=A\tau_j^{(n)}\exp(-i\theta);
\]
les $\tau_j^{(n)}$ ($j=1,2
\ldots n$) sont les racines
du polynôme de Laguerre $L_n(x)$ de degré $n$ et $A$ est un scalaire
$>{1/ {2\vert a\cos (\theta_s)\vert }}$, $a$ étant la
singularité
de $\varphi (t)$ la plus proche de $O$.

En pratique, cette approximation rationnelle revient à accumuler,
de plus en plus près de l'origine, les racines $-{1/ t_j}$ de
$Q_n$ sur $-d$, d'une part si $n$ augmente, mais aussi si $d$ tend
vers les positions limites, où $d$ forme un angle de valeur absolue
${\pi/ 2}$ avec la direction singulière.
Entre ces positions $d$ peut prendre toutes les directions dans le
demi-plan complémentaire au demi-plan bissecté par cette direction
singulière.

Autrement dit, nous pouvons encore construire une approximation
rationnelle de $f(x)$, qui réalise une généralisation du
prolongement
analytique, en \emph{simulant} la direction singulière par une
accumulation de pôles dans cette direction.

L'évaluation de $\psi_1(x^2)$ en $x=0.15$ donne cette fois-ci:
$0.622166$.

Soit la fonctionnelle linéaire $c$ telle que:
$$ c(1)=a_0,\ c(x)=a_1,\ldots c(x^n)=a_n\ldots $$
Un choix de $v(x)$, s'il existe, tel que:
$$ c(x^\ell v(x))=0 \qquad {\rm pour} \ \ell=0,1\ldots, n-1, $$
nous conduit à un choix des $t_\ell^{(n)}$ ($\ell=1,2,\ldots n$)
du dénominateur $Q_n(x)=x^nv(1/x)$ tels que l'approximant de type
Padé ${P_n(x)/ Q_n(x)}$ devienne l'approximant de Padé
tout court.

Malheureusement, dans ce cas, on ne peut plus rien conclure sur la
convergence de ces approximants de Padé vers $f(x)$, sauf si la
série $\widehat f(x)$ est du type Stieltjes.

Dans le cas de la série d'Euler $1-1!x+2!x^2\cdots$, l'approximant
de Padé, construit à partir du dénominateur où $v_n(x)$ est un
polynôme de Laguerre, converge bien vers la somme $f(x)$ sauf si
$x\in \RR^-$.

Pour $x=0.1$, comme dans le paragraphe 3, nous obtenons $f(0.1)=
0.91563333$ avec 8 chiffres significatifs exacts.

Dans les cas où on peut établir la convergence des approximants de
Padé vers $f(x)$, on dispose ainsi d'un algorithme formel permettant
de
calculer cet approximant, alors qu'en 5.3 et 5.4, on dispose d'un
algorithme essentiellement numérique.

\subsection{Conclusion}
Nous constatons que l'évaluation de la fonction
$\psi_1(x^2)$ pour $x=0.15$, par les 3 méthodes décrites, donne
un résultat stable sur 3 chiffres significatifs, ce qui peut être
considéré comme acceptable en tenant compte de
la divergence très forte de $\widehat g(x)$,
donc de $\widehat\psi_1(x^2)$.

D'ailleurs la seule façon de vérifier la validité des
résultats est
de confronter, en diverses régions
du plan complexe, les valeurs obtenues
par différentes méthodes.

Un logiciel complet de sommation de séries divergentes, solutions
d'équations différentielles, doit être interactif, et utiliser
toutes
les ressources données par ces équations, pour analyser les
régions
de stabilité des différentes méthodes de calcul utilisées, et
confronter
ces méthodes par recoupement des résultats.

L'outil graphique interactif est alors indispensable et le seul
moyen d'avoir une vision d'ensemble (comme dans le logiciel
développé
par F.\,Richard-Jung \cite{8}).

Cette analyse est caractérisée par le \emph{contrôle formel}, qui,
suivant des choix prévus ou donnés interactivement, engendre
les différentes parties numériques.

\backmatter
\nocite{*}
\bibliographystyle{jepplain+eid}
\bibliography{xups91-03}
\end{document}