Master 1 - Calcul Formel

Année Universitaire 2016-2017

Notes de Cours

Thomas Cluzeau
Université de Limoges - CNRS ; XLIM UMR 7252
123 avenue Albert Thomas, 87060 Limoges CEDEX
[Link]@[Link]
[Link]
2
Table des matières

1 Introduction 5
1.1 Qu’est-ce que le calcul formel ? . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.2 Structures algébriques effectives . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.3 Analyse d’un algorithme . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.4 Principes algorithmiques généraux . . . . . . . . . . . . . . . . . . . . . . . . 8
1.4.1 Évaluation-Interpolation . . . . . . . . . . . . . . . . . . . . . . . . . 8
1.4.2 Méthodes modulaires . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.4.3 Diviser pour régner . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

2 Algorithmique rapide pour les opérations de base en calcul formel 11

2.1 Multiplication de polynômes . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.1.1 Algorithme naïf . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.1.2 Algorithme de Karatsuba . . . . . . . . . . . . . . . . . . . . . . . . 12
2.1.3 Transformée de Fourier rapide (FFT) . . . . . . . . . . . . . . . . . . 13
2.2 Inversion de séries formelles et division euclidienne . . . . . . . . . . . . . . . 15
2.2.1 L’anneau des séries formelles . . . . . . . . . . . . . . . . . . . . . . . 15
2.2.2 Calcul formel dans l’anneau des séries formelles . . . . . . . . . . . . 16
2.2.3 Application à la division euclidienne . . . . . . . . . . . . . . . . . . 19
2.3 Multiplication de matrices . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
2.3.1 Algorithme naïf . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
2.3.2 Algorithme de Strassen . . . . . . . . . . . . . . . . . . . . . . . . . . 21

3 Pgcd et résultant 25
3.1 Rappels : anneaux factoriels et anneaux euclidiens . . . . . . . . . . . . . . . 25
3.1.1 Divisibilité, pgcd et ppcm . . . . . . . . . . . . . . . . . . . . . . . . 25
3.1.2 Division euclidienne . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
3.2 Calcul du pgcd et des coefficients de Bézout . . . . . . . . . . . . . . . . . . 29
3.2.1 Dans un anneau euclidien : le cas des entiers . . . . . . . . . . . . . . 29
3.2.2 Le théorème des restes chinois . . . . . . . . . . . . . . . . . . . . . . 31
3.2.3 Dans un anneau de polynômes . . . . . . . . . . . . . . . . . . . . . . 35
3.3 Résultant de deux polynômes . . . . . . . . . . . . . . . . . . . . . . . . . . 40
3.3.1 Matrice de Sylvester et forme échelonnée . . . . . . . . . . . . . . . . 40
3.3.2 Propriétés du résultant . . . . . . . . . . . . . . . . . . . . . . . . . . 43

3
 3.3.3 Calcul du résultant . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44
3.3.4 Résultant en plusieurs variables et élimination . . . . . . . . . . . . . 44

4 Factorisation de polynômes à une variable 47

4.1 Factorisation sans carré . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.1.1 Définition et existence . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.1.2 Calcul de la factorisation sans carré . . . . . . . . . . . . . . . . . . . 48
4.2 Factorisation sur Z . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.2.1 Rappels : factorisation sur un corps fini . . . . . . . . . . . . . . . . . 52
4.2.2 Algorithme utilisant un grand nombre premier . . . . . . . . . . . . . 52
4.2.3 Algorithme utilisant la remontée de Hensel . . . . . . . . . . . . . . . 53

5 Intégration symbolique 57
5.1 Algèbre différentielle et fonctions élémentaires . . . . . . . . . . . . . . . . . 57
5.2 Le cas particulier des fractions rationnelles . . . . . . . . . . . . . . . . . . . 59
5.2.1 Calcul de la partie polynomiale . . . . . . . . . . . . . . . . . . . . . 59
5.2.2 Calcul de la partie rationnelle . . . . . . . . . . . . . . . . . . . . . . 60
5.2.3 Calcul de la partie logarithmique . . . . . . . . . . . . . . . . . . . . 64
5.3 Le cas général des fonctions élémentaires . . . . . . . . . . . . . . . . . . . . 66

6 Équations différentielles linéaires et séries formelles 69

6.1 Résolution de L(y) = 0 par des séries entières . . . . . . . . . . . . . . . . . 69
6.1.1 Structure des solutions . . . . . . . . . . . . . . . . . . . . . . . . . . 69
6.1.2 Existence de séries solutions au voisinage d’un point ordinaire . . . . 71
6.1.3 Calcul des séries solutions au voisinage d’un point ordinaire . . . . . 72
6.2 Résolution de L(y) = 0 par des séries quasi-entières . . . . . . . . . . . . . . 75
6.2.1 Séries quasi-entières . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
6.2.2 Équation indicielle et calcul des séries quasi-entières solutions . . . . 76

4
Chapitre 1

Introduction

1.1 Qu’est-ce que le calcul formel ?

En mathématiques, lorsqu’un problème consiste à étudier l’existence d’un objet répondant
à certains critères, il existe deux façons différentes d’y apporter une réponse. La première
consiste à prouver l’existence de l’objet en question sans se préoccuper de sa construction
alors que la seconde produit un procédé de construction de l’objet (ce qui prouve donc son
existence). C’est ce qu’on appelle une preuve constructive. On distingue alors deux types
de procédés constructifs : les procédés numériques qui fournissent une approximation du
résultat et les procédés formels qui fournissent un résultat exact. Ceux sont ces derniers qui
seront étudiés dans ce cours.
Le terme calcul formel (en anglais computer algebra ou symbolic computation) désigne
l’étude de la représentation et de la manipulation effective d’objets mathématiques dans
un ordinateur. En d’autres termes, le calcul formel consiste à étudier dans quelle mesure
et à quelle vitesse un ordinateur peut faire des mathématiques. Pour un problème donné,
nous serons donc amené à produire un algorithme (i.e., une suite de règles opératoires per-
mettant, en un nombre finie d’étapes, d’arriver avec certitude au résultat) et en évaluer les
performances.
Les systèmes de calcul formel manipulent principalement trois types d’objets à savoir
des nombres, e.g., des entiers relatifs (arithmétique), des polynômes et des matrices (algèbre
linéaire). Un premier aspect du calcul formel consiste donc à proposer des algorithmes de
calcul efficaces pour la manipulation de ces structures de bases. Un deuxième aspect se
propose de ramener d’autres problèmes à ces derniers comme par exemple montrer que
la plupart des opérations d’algèbre linéaire peuvent se réduire à calculer des produits de
matrices (voir le théorème 2.3.1 dans la section 2.3). Les procédés que nous construirons
reposeront toujours sur de l’algèbre donc avant de construire un algorithme il faudra toujours
se placer dans un certain modèle/cadre algébrique.
Le calcul formel est une discipline récente qui a été surtout développée, à la fois par
des mathématiciens, des informaticiens et des physiciens, à partir de la seconde moitié du
xxe siècle et l’apparition des ordinateurs. Aujourd’hui les grands systèmes de calcul formel

5
comme Maple ou Mathematica sont utilisés dans l’enseignement et la recherche mais
aussi dans l’industrie.

1.2 Structures algébriques effectives

Étant donnée une structure algébrique (e.g., groupe, anneau, corps, espace vectoriel) la
première question qui se pose avant de construire des algorithmes manipulant des objects de
cette structure est le problème de la calculabilité : suis-je capable de faire des calculs dans
cette structure et d’apprendre à un ordinateur à les faire ? Le principal problème est le test
d’égalité c’est-à-dire le problème de décider si une expression donnée dans cette structure
vaut 0 ou non. En effet pour effectuer certains calculs, comme par exemple une division, il
est souvent crucial de déterminer si une expression représente 0 ou non. Notons que ce test à
zéro n’est pas décidable dans l’ensemble des nombres réels R. En calcul formel, on se ramène
donc toujours à faire des calculs dans une structure effective où le test à zéro est décidable.
Définition 1.2.1. Une structure algébrique est dite effective si l’on dispose :
1. d’une structure de données pour en représenter les éléments ;

2. d’algorithmes pour en effectuer les opérations et pour y tester l’égalité à zéro.

Le lemme suivant regroupe les résultats d’effectivité sur lesquels nous nous appuierons
dans ce cours.
Lemme 1.2.1. On a les résultats suivants :
1. L’anneau Z des entiers relatifs est effectif ;

2. Si A est un anneau effectif, alors son corps des fractions est effectif, A[X] est effectif
et pour N ∈ N, A[X]/(X N +1 ) est effectif ;

3. Si K est un corps effectif, alors sa cloture algébrique K, l’espace vectoriel Kn et l’anneau

Mn (K) sont effectifs.

1.3 Analyse d’un algorithme

Le calcul formel a pour but de produire des algorithmes qui calculent des objects mathéma-
tiques. Rappelons tout d’abord que la correction d’un algorithme est définie comme suit :
Définition 1.3.1. Un algorithme est dit correct si :
1. chaque étape est bien définie ;

2. l’algorithme se termine en un nombre fini d’étapes ;

3. le résultat renvoyé est toujours le bon, i.e., celui que l’on attend.

6
 Le problème de la calculabilité abordé ci-dessus indique uniquement la faisabilité de
certaines opérations. Pour un problème mathématique donné, il peut exister plusieurs algo-
rithmes menant au même résultat. Le calcul formel s’intéresse alors à la comparaison entre
ces différents algorithmes pour pouvoir exprimer le fait que pour un problème donné certains
algorithmes sont « meilleurs » que d’autres. Pour ce faire, on évalue l’efficacité d’un algo-
rithme en estimant sa complexité. Estimer la complexité d’un algorithme permet d’évaluer
sa qualité indépendamment de la manière dont il est programmé.
En pratique, une fois l’algorithme programmé sur un ordinateur, ses performances sont
essentiellement mesurées par le temps d’exécution du programme (et la quantité de mémoire
nécessaire à son exécution). Le temps d’exécution va dépendre à la fois du nombre d’opé-
rations effectuées et du coût de chaque opération. Le programmeur contrôle essentiellement
le nombre d’opérations effectuées donc c’est avant tout sur ce point qu’il doit agir pour
optimiser les performances de son programme. Une estimation a priori du temps d’exécu-
tion peut donc être obtenue en analysant le programme c’est-à-dire en comptant le nombre
d’opérations effectuées, le nombre d’appels récursifs, . . . .
Dans le cadre de ce cours, c’est justement le nombre d’opérations arithmétiques effectuées
par un algorithme qui nous servira de mesure de complexité. Cette mesure est évaluée en fonc-
tion d’une (ou plusieurs) grandeurs apparaissant dans l’entrée de l’algorithme (e.g., le degré
d’un polynôme, la taille d’une matrice, . . . ). De plus nous n’avons pas toujours une mesure
exacte du nombre d’opérations mais souvent un ordre de grandeur (ou une borne supérieure).
Nous utilisons alors la notation O(.) pour exprimer la complexité de nos algorithmes.

Définition 1.3.2. Soient f et g deux fonctions. On note f (n) = O(g(n)) s’il existe K > 0
et N > 0 tels que pour tout n > N , |f (n)| ≤ K |g(n)|.

Un algorithme sera dit linéaire (resp. quadratique, cubique, logarithmique, polynomial)

s’il utilise O(n) (resp. O(n2 ), O(n3 ), O(log(n)), O(nk ) avec k ∈ N) opérations arithmétiques.
Notons que ceci donne une estimation asymptotique de la complexité (l’inégalité est vérifiée
pour n > N ) et que la valeur de la constante K peut jouer un rôle important dans les
performances de l’algorithme.

Exemple 1.3.1. P Soit K un corps effectif et considérons le problème de l’évaluation d’un

polynôme P = ni=0 pi X i ∈ K[X] en un point a ∈ K.
En faisant les calculs naïvement, nous devons calculer a2 , a3 , . . . , an puis les n produits pi ai
et enfin faire la somme. Cela nous coutera donc (n − 1) + n = 2 n − 1 multiplications dans
K et n additions dans K soit 3 n − 1 opérations dans K.
Le schéma de Horner produira le même résultat en utilisant la formule suivante :

P (a) = p0 + a (p1 + a (p2 + a (· · · + a (pn−1 + a pn ) · · · ))).

En nous basant sur cette formule, nous pouvons donc calculer P (a) en posant u = pn puis
en calculant, pour i = 1, . . . , n, u = a u + pn−i . Cela nous coutera donc n additions et n
multiplications soit 2 n opérations dans K.
Ainsi, nous voyons que l’algorithme basé sur le schéma de Horner semble meilleur que celui

7
basé sur le calcul naïf. Dans cet exemple particulier, nous obtenons une complexité linéaire
en O(n) pour les deux algorithmes mais parfois le nombre d’opérations peut changer signi-
ficativement en passant d’une méthode à l’autre (Cf le calcul du déterminant d’une matrice
carrée à coefficients dans un corps effectif ). Notons aussi qu’ici l’algorithme naïf nécessite
l’utilisation de n cases mémoire alors que celui de Horner n’en a besoin que d’une seule.
Dans la suite du cours, nous nous restreindrons à évaluer la complexité arithmétique de
nos algorithmes c’est-à-dire le nombre total d’opérations arithmétiques dans la structure
algébrique effective sous-jacente A. En pratique, cette mesure reflète le coût réel de l’algo-
rithme lorsque le coup des opérations dans A est prépondérant et que chaque opération dans
A a un coût constant. Ceci est le cas par exemple lorsque A = Fp est un corps fini mais pas
lorsque A = Z ou Q où le coût d’une opération dépend de la taille des entiers qui entrent
en jeu. Dans le cas d’algorithmes manipulant des objects sur Z, la complexité arithmétique
n’est pas toujours pertinente car il faut prendre en compte la taille des entiers que l’on ma-
nipule qui peut grossir considérablement au cours des calculs. Un meilleur moyen d’évaluer
la performance a priori d’un algorithme est alors d’évaluer sa complexité binaire à savoir le
nombre d’opérations binaires effectuées par l’algorithme.

1.4 Principes algorithmiques généraux

Dans cette section, nous allons donner brièvement trois principes algorithmiques qui seront
utilisés par la suite pour développer des algorithmes efficaces/rapides en calcul formel. Cette
liste n’est évidemment pas exhaustive et il existe bien d’autres principes très utiles pour
certains problèmes comme par exemple la récursivité, la technique pas de bébés/pas de
géants, . . . .

1.4.1 Évaluation-Interpolation
Le principe d’évaluation-interpolation est utilisé pour développer des algorithmes manipulant
des polynômes. Soit une opération op à effectuer sur deux polynômes P et Q à coefficients
dans un anneau effectif A. Le principe de l’évaluation-interpolation pour calculer R = P op Q
est alors le suivant :
1. On choisit un certain nombre de points d’évaluation a1 , . . . , an ∈ A ;
2. On évalue les polynômes P et Q en a1 , . . . , an (évaluation) ;
3. On en déduit l’évaluation de R en a1 , . . . , an ;
4. À partir des valeurs R(a1 ), . . . , R(an ), on reconstruit R (interpolation).
L’efficacité de la stratégie repose alors sur le choix des points a1 , . . . , an dans l’anneau A. Le
nombre n de points d’évaluation doit être suffisamment grand pour garantir la faisabilité de
l’étape 4. De plus les points eux mêmes doivent être choisis pour optimiser le coût des étapes
2 et 4. Par exemple, choisir comme points d’évaluation des racines de l’unité, où des points
en progression arithmétique peut s’avérer particulièrement intéressant.

8
1.4.2 Méthodes modulaires
En calcul formel, beaucoup d’algorithmes manipulant des polynômes, des fractions ration-
nelles ou des matrices à coefficients entiers souffrent de la croissance des expressions inter-
médiaires. En d’autres termes, les entiers apparaissants au cours des calculs sont de très
grande taille par rapport à ceux qui figurent à la fois dans l’entrée et dans la sortie de l’algo-
rithme. Pour illustrer ce propos, considérons le problème du calcul du pgcd par l’algorithme
d’Euclide (voir l’algorithme 7 au chapitre 3).

Exemple 1.4.1. Considérons les polynômes P = 7 X 5 − 22 X 4 + 55 X 3 + 94 X 2 − 87 X + 56

et Q = 62 X 4 −97 X 3 +73 X 2 +4 X +83. Partant de R0 = P , R1 = Q, l’algorithme d’Euclide
produit alors la suite de restes suivante :
113293 409605 2 183855 272119
R2 = rem(R0 , R1 ) = 3844
X3 + 3844
X − 1922
X+ 3844
,

18423282923092 15239170790368 10966361258256

R3 = rem(R1 , R2 ) = 12835303849
X2 − 12835303849
X+ 12835303849
,

R4 = rem(R2 , R3 ) = − 216132274653792395448637
44148979404824831944178
X− 631179956389122192280133
88297958809649663888356
,

20556791167692068695002336923491296504125
R5 = rem(R3 , R4 ) = 3639427682941980248860941972667354081
.

On voit donc que les coefficients dans les calculs intermédiaires font intervenir des entiers
qui croissent de manière exponentielle alors que le résultat est pgcd(P, Q) = 1.

Les méthodes modulaires permettent de remédier à ce problème. Au lieu d’effectuer les

calculs sur Z, nous effectuons les calculs dans Fp = Z/(p Z) pour un (ou plusieurs) nombre(s)
premier(s) p. Deux cas sont alors à différencier :

1. Si notre algorithme consiste à décider si une propriété est vraie ou fausse, ou à calculer
un degré ou une dimension, alors l’algorithme appliqué aux entrées réduites modulo p
donne un algorithme probabiliste répondant à la question. De plus, si on peut maitriser
les « mauvais » nombres premiers pour lesquels la réponse est fausse, alors on obtient
un algorithme déterministe ;

2. Lorsqu’une borne sur la taille des entiers du résultat est connue, on peut alors utiliser
une stratégie basée sur le théorème des restes chinois (voir Chapitre 3) qui implique
qu’un entier inférieur à un produit de nombres premiers p1 · · · pk peut être reconstruit
à partir de ses réductions modulo p1 , . . . , pk . On effectue alors le calcul modulo suffi-
samment de nombres premiers et on reconstruit le résultat. Cette stratégie est alors du
même type que l’évaluation-interpolation.

1.4.3 Diviser pour régner

Un autre principe important utilisé (peut-être même le plus important) pour la construction
d’algorithmes efficaces en calcul formel est le paradigme « diviser pour régner ». Il consiste

9
à résoudre un problème en le réduisant à m instances du même problème sur des entrées de
taille divisée par p (en général p = 2) puis à reconstruire le résultat. Pour estimer le coût
total C(n) d’un tel algorithme, on doit connaitre le coût T (n) de la recombinaison et du
découpage initial ainsi que le coût κ de l’algorithme utilisé lorsque la taille s ≥ p de l’entrée
est suffisamment petite. Le coût total obéit alors souvent à une relation de la forme suivante :

T (n) + m C(dn/pe), si n ≥ s,
C(n) ≤ (1.1)
κ sinon.

L’efficacité de la technique dépend alors fortement de la fonction T et on a le résultat (admis)

suivant qui sera suffisant dans le cadre ce cours.

Théorème 1.4.1. Soit C(n) satisfaisant l’inégalité (1.1) avec m > 0, κ > 0 et où T est une
fonction telle qu’il existe q et r avec 1 < q ≤ r vérifiant q T (n) ≤ T (p n) ≤ r T (n) pour n
assez grand1 . Alors lorsque n → ∞, on a :

 O(T (n)), si q > m,


C(n) = O(T (n) log(n)), si q = m, (1.2)
  
 O nlogp (m) T (n) , si q < m.

nlogp (q)

Exemple 1.4.2. Au chapitre 2, nous verrons que le coût K(n) de l’algorithme de Karatsuba
pour le calcul du produit de deux polynômes de degré au plus n − 1 vérifie :

K(n) ≤ 4 n + 3 K(dn/2e), K(1) = 1.

Avec les notations précédentes, on a donc m = 3, p = s = 2, κ = 1 et T (n) = 4 n vérifie

q T (n) ≤ T (p n) ≤ r T (n) avec q = r = p = 2. On a donc q < m de sorte que (1.2) implique
 
log2 (3) 4n
K(n) = O n = O(nlog2 (3) ).
nlog2 (2)

1
Notons que cette condition est vérifiée lorsque T (n) = nα logβ (n) avec α > 0.

10
Chapitre 2

Algorithmique rapide pour les opérations

de base en calcul formel

Dans ce chapitre, nous allons nous intéresser à trois opérations de base en calcul formel à
savoir la multiplication de polynômes, la division des séries formelles (et la division eucli-
dienne des polynômes) et la multiplication de matrices. Pour chaque problème, nous allons
donner des algorithmes rapides permettant d’effectuer l’opération en question. Notons qu’il
est d’autant plus important d’avoir des algorithmes rapides pour ces opérations que la com-
plexité de beaucoup d’autres algorithmes manipulant les mêmes objects peut s’exprimer en
fonction de la complexité de l’algorithme utilisé pour effectuer ces trois opérations.

2.1 Multiplication de polynômes

Soient P et Q deux polynômes univariés de degré au plus n et à coefficients dans un anneau
effectif A. Le but de cette section est de donner des algorithmes rapides pour calculer le
produit P × Q = P Q. En d’autres termes, nous allons étudier, en fonction n, en combien
d’opérations dans A nous pouvons effectuer le produit P Q. Notons que l’on obtient des
résultats très semblables pour la complexité binaire du produit d’entiers. En effet, intuiti-
vement nous pouvons décomposer un entier sur une base B (e.g., B = 10, 2, . . .) et obtenir
« en gros » les mêmes algorithmes modulo la gestion des retenues.

2.1.1 Algorithme naïf

Soient à multiplier les polynômes
Xn n
X
P = pi X i = p0 + p1 X + · · · + pn X n , Q= q i X i = q0 + q 1 X + · · · + qn X n ,
i=0 i=0
avec p0 , . . . , pn , q0 , . . . , qn ∈ A. Le produit P Q s’écrit alors
2n
X X
PQ= hi X i , ∀ i = 0, . . . , 2 n, hi = pj q k . (2.1)
i=0 j+k=i

11
Proposition 2.1.1. En utilisant la formule (2.1), le produit P Q de deux polynômes P et Q
de degré au plus n peut s’effectuer en utilisant n2 additions et (n + 1)2 multiplications dans
A. La complexité arithmétique de l’algorithme naïf est donc O(n2 ) opérations dans A : c’est
un algorithme quadratique.

Démonstration. Pour i = 0, . . . , n, calculer hi nécessite i additions et (i + 1) multiplications.

De même, pour i = 0, . . . , n − 1, calculer h2 n−i nécessite iP
additions et (i + 1) multiplications.
n−1
En conséquence le nombre totalPd’additions est égal à 2 i=0 i + n = n2 et le nombre total
n−1
de multiplications est égal à 2 i=0 (i + 1) + (n + 1) = (n + 1)2 .

2.1.2 Algorithme de Karatsuba

L’algorithme pour la multiplication de polynômes dû à Karatsuba repose sur la remarque
suivante dans le cas de polynômes P = p0 + p1 X et Q = q0 + q1 X de degré 1. Le produit
P Q s’écrit alors
P Q = p0 q0 + (p0 q1 + p1 q0 ) X + p1 q1 X 2 ,
de sorte que l’algorithme naïf effectuera les 4 produits p0 q0 , p0 q1 , p1 q0 et p1 q1 et 1 addition.
En revanche, en remarquant que

(p0 q1 + p1 q0 ) = (p0 + p1 ) (q0 + q1 ) − p0 q0 − p1 q1 ,

on voit que l’on peut calculer P Q en effectuant seulement des 3 produits p0 q0 , p1 q1 et

(p0 + p1 ) (q0 + q1 ) et 4 additions. Nous perdons certes 3 additions par rapport à l’algorithme
naïf mais le gain d’une multiplication va au final entrainer une amélioration de la complexité.

Pour deux polynômes quelconques, le principe de l’algorithme de Karatsuba consiste à

utiliser l’observation précédente de manière récursive en utilisant une stratégie diviser pour
régner (voir la sous-section 1.4.3). Soient à multiplier deux polynômes P et Q de degré au
plus n − 1 avec n = 2k . On écrit alors
n n
P = P0 + P 1 X 2 , Q = Q0 + Q1 X 2 , (2.2)

où P0 , P1 , Q0 et Q1 sont des polynômes de degré au plus n/2−1 = 2k−1 −1. On applique alors
la stratégie expliquée ci-dessus pour les polynômes de degré 1 et on aura 3 appels récursifs
sur des polynômes de degré au plus n/2−1 pour calculer P0 Q0 , P1 Q1 et (P0 +P1 ) (Q0 +Q1 ).
On obtient l’algorithme 1.

Proposition 2.1.2. L’algorithme 1 de Karatsuba calcule le produit P Q de deux polynômes

de degré au plus n − 1 avec n = 2k en O(n1.59 ) opérations dans A.

Démonstration. Soit K(n) le nombre d’opérations arithmétiques effectuées par l’algorithme

de Karatsuba. Un appel sur des polynômes de degré au plus n − 1 fait 3 appels récursifs sur
des polynômes de degré au plus n/2 − 1 ainsi que 2 additions de polynômes de degré au plus
n/2−1 et 2 additions de polynômes de degré au plus n−1. Finalement, la dernière étape peut

12
Algorithme 1 Algorithme de Karatsuba
Entrées : P, Q ∈ A[X] de degré au plus n − 1 avec n = 2k .
Sorties : Le produit P Q.
Si n = 1 Alors
- Retourner P Q ;
Sinon
- Décomposer P et Q sous la forme (2.2) ;
- Calculer H0 = P0 Q0 et H2 = P1 Q1 en utilisant deux appels récursifs ;
- Calculer P0 + P1 et Q0 + Q1 ;
- Calculer G = (P0 + P1 ) (Q0 + Q1 ) en utilisant un appel récursif ;
- Calculer H1 = G − H0 − H2 ;
n
- Retourner H0 + H1 X 2 + H2 X n .
Finsi

se faire à l’aide d’une seule addition de polynômes de degré au plus n − 1 (car H0 et H2 X n

sont à supports monomiaux disjoints et peuvent donc être additionnés gratuitement). On a
donc K(n) ≤ 3 K(n/2) + 4 n avec K(1) = 1 (coût de la multiplication de deux polynômes
de degré 0) d’où le résultat d’après l’exemple 1.4.2 (log2 (3) ≈ 1.584962501).

Remarque 2.1.1. Lorsque n − 1 n’est pas de la forme 2k , on complète les polynômes avec
des coefficients nuls pour se ramener au cas n − 1 = 2k . Puisque le premier entier plus grand
que n−1 de la forme 2k est au plus égal à 2 (n−2), la complexité reste en O(n1.59 ) opérations
dans A.

2.1.3 Transformée de Fourier rapide (FFT)

L’algorithme basé sur la transformée de Fourier rapide est à ce jour le meilleur algorithme
connu pour la multiplication de polynômes univariés. Cet algorithme est basé sur une stra-
tégie évaluation-interpolation (voir la sous-section 1.4.1).

Étant donné un polynôme P = p0 + p1 X + · · · + pn−1 X n−1 de degré au plus n − 1 avec

n = 2k , on pose
n n
Ppair = p0 + p2 X + p4 X 2 + · · · + pn−2 X 2 −1 , Pimpair = p1 + p3 X + p5 X 2 + · · · + pn−1 X 2 −1 ,

de sorte que
P (X) = Ppair (X 2 ) + X Pimpair (X 2 ). (2.3)
Évaluer P en n points a0 , . . . , an−1 revient donc à évaluer Ppair et Pimpair en a20 , . . . , a2n−1 .

L’idée consiste maintenant à utiliser comme points d’évaluation des racines de l’unité.

Définition 2.1.1. Un élément ω ∈ A est une racine ne de l’unité si ω n = 1. C’est une racine
primitive ne de l’unité si de plus ω m 6= 1 pour tout m < n.

13
 Plaçons nous tout d’abord sur C où nous disposons de la racine primitive ne de l’unité
donnée par ω = exp 2 ni π . La suite des carrés 1, ω 2 , (ω 2 )2 , . . . , (ω 2 )n−1 correspond alors

à 1, ω 2 , ω 4 , . . . , ω n−2 , 1, ω 2 , ω 4 , . . . , ω n−2 . Évaluer Ppair (respectivement Pimpair ) en les n

points 1, ω 2 , . . . , ω 2 n−2 peut donc se faire en évaluant Ppair (resp. Pimpair ) en les n/2 points
1, ω 2 , . . . , ω n−2 .
On obtient alors l’algorithme 2 ci-dessous utilisant une stratégie diviser pour régner pour
évaluer P en 1, ω, . . . , ω n−1 ce qu’on appelle aussi calculer la transformée de Fourier discrete
(DFT) de P .

Algorithme 2 Transformée de Fourier Rapide (FFT)

Entrées : P = n−1 2iπ
P i k


i=0 pi X ∈ C[X] avec n = 2 et ω = exp n
.
n−1
Sorties : P (1), . . . , P (ω ).
Si n = 1 Alors
- Retourner p0 ;
Sinon
- Décomposer P sous la forme (2.3) ;
- Évaluer Ppair et Pimpair en 1, ω 2 , . . . , ω n−2 à l’aide de deux appels récursifs ;
- Retourner les valeurs de P (1), . . . , P (ω n−1 ) obtenus à partir de (2.3).
Finsi

Proposition 2.1.3. Soit ω = exp 2 ni π une racine primitive ne de l’unité. L’algorithme 2

ci-dessus, appelé, transformée de Fourier rapide (FFT) calcule la transformée de Fourier

discrete (DFT) de P ∈ C[X] de degré au plus n − 1 avec n = 2k en O(n log(n)) opérations
dans C.

Démonstration. Soit F (n) le nombre d’opérations arithmétiques effectuées par l’algorithme

FFT. L’algorithme utilise deux appels récursifs sur des polynômes de degré au plus n/2 − 1
à évaluer en n/2 points plus un nombre d’opérations proportionnel à n à la dernière étape.
On a donc F (n) ≤ 2 F (n/2) + C n, où C est une constante et F (1) = 0. Le résultat est donc
une application directe du théorème 1.4.1 avec p = m = q = r = 2 et κ = 0.

On peut montrer (voir TD) que l’opération inverse, i.e., l’interpolation qui consiste, étant
données les valeurs de P en 1, ω, . . . , ω n−1 à retrouver les coefficients du polynôme P se
ramène à calculer une transformée de Fourier discrete sur les points 1, ω −1 , . . . , (ω −1 )n−1 .
On peut donc donner l’algorithme de multiplication de polynômes basé sur la FFT.

Proposition 2.1.4. L’algorithme 3 de multiplication par FFT ci-dessus calcule le produit

de deux polynômes de degré au plus n − 1 avec n = 2k en O(n log(n)) opérations dans C.

Démonstration. Les étapes d’évaluation et d’interpolation se font à l’aide de 3 FFT. Les

autres étapes ont un coût en O(n) opérations dans C. Le résultat découle donc de la propo-
sition 2.1.3.

14
Algorithme 3 Multiplication de deux polynômes par FFT
Entrées : P, Q ∈ C[X] de degré au plus n2 − 1 avec n = 2k et ω = exp 2iπ


n
.
Sorties : Le produit P Q.
- Pré-calculer les puissances ω 2 , . . . , ω n−1 ;
- Évaluation : utiliser la FFT pour évaluer P et Q en 1, ω, . . . , ω n−1 ;
- Calculer les produits P (1) Q(1), P (ω) Q(ω), . . . , P (ω n−1 ) Q(ω n−1 ) ;
- Interpolation : utiliser la FFT pour calculer P Q.

Les deux algorithmes donnés dans cette sous-section reposent sur la notion de
racine de
l’unité. Lorsque nous sommes dans C, nous disposons toujours de ω = exp 2 ni π comme
racine primitive ne de l’unité. Si l’anneau A des coefficients de nos polynômes n’est pas
C, alors deux cas sont à distinguer : si A vérifie qu’il existe des racines primitives ne de
l’unité pour tout n = 2k , alors on peut encore multiplier deux polynômes de degré au plus
n − 1 en O(n log(n)). Par contre, si ce n’est pas le cas, nous devons construire ces racines
ce qui entraine un surcoût et la multiplication utilise alors O(n log(n) log(log(n))) opérations.

Résumé : nous avons donc vu dans cette section que le produit de deux polynômes de
degré au plus n et à coefficients dans A peut se calculer en :

1. O(n2 ) opérations dans A par l’algorithme naïf ;

2. O(n1.59 ) opérations dans A par l’algorithme de Karatsuba ;

3. O(n log(n) log(log(n))) opérations dans A en utilisant la FFT.

En pratique, il faut faire attention aux constantes qui se cachent dans la notation O(.) car
elles peuvent être déterminantes pour l’efficacité pratique de l’algorithme. Dans les meilleures
implémentations actuelles, l’algorithme de Karatsuba l’emporte sur l’algorithme naïf pour
des degrés de l’ordre de 20 et l’algorithme basé sur la FFT gagne pour des degrés de l’ordre
de 100. Notons que pour certains problèmes (e.g., en cryptographie) nous devons manipuler
des polynômes de degrés de l’ordre de 100000.

2.2 Inversion de séries formelles et division euclidienne

2.2.1 L’anneau des séries formelles
Définition 2.2.1. Soit A un anneau. On note A[[X]] l’ensemble des séries formelles de la
forme X
F = f0 + f1 X + f2 X 2 + · · · = fi X i , fi ∈ A.
i≥0

Le terme f0 s’appelle le terme constant de la série F .

15
 Une série formelle peut s’identifier avec la suite (fi )i∈N de ses coefficients qui est donc
une suite d’éléments de A. L’ensemble A[[X]] est muni d’une addition et d’une multiplication
données par : X X X
fi X i + gi X i = (fi + gi ) X i ,
i≥0 i≥0 i≥0
! !
X X X X
fi X i gi X i = hi X i , hi = fj gk .
i≥0 i≥0 i≥0 j+k=i

L’élément neutre pour le produit est alors la série formelle 1 + i>0 0 X i où 1 P

P
= 1A est
multiplication. Si on calcule le produit de FP= 1+X + i>1 0 X i
l’élément neutre de A pour la P
par G = 1 − X + X 2 − X 3 + i>3 (−1)i X i , alors on trouve F G = 1 + i>0 0 X i . La série G
est alors l’inverse de F pour la multiplication et on note G = F −1 . D’une manière générale,
nous avons le résultat suivant :

Proposition 2.2.1. L’ensemble A[[X]] des séries formelles est un anneau qui est commutatif
si A est commutatif. Les éléments inversibles de A[[X]] sont les séries formelles de terme
constant inversible dans A.

Démonstration. La première partie de la proposition est claire. Si F une série formelle de la

forme F = 1 + X G, où G est une série formelle, alors on a F −1 = 1 − X G + X 2 G2 − · · · .
Maintenant si le terme constant de F est un élément inversible a ∈ A, on écrit F = a (1+X G)
avec G = a−1 (F − a)/X et on a donc F −1 = (1 + X G)−1 a−1 . On a donc montré que toutes
les séries formelles dont le terme constant est inversible dans A sont inversibles dans A[[X]].
Inversement si F est inversible dans A[[X]], si G désigne son inverse et si f0 (resp. g0 ) désigne
le terme constant de F (resp. G), alors on a f0 g0 = 1 dans A de sorte que f0 est inversible
dans A.

2.2.2 Calcul formel dans l’anneau des séries formelles

Du point de vue du calcul formel se pose alors la question suivante : comment faire des cal-
culs dans A[[X]] lorsque A est un anneau effectif ? En pratique, on ne peut bien évidemment
pas manipuler des objets « infini » mais on doit considérer des troncatures de ces objets
c’est-à-dire, dans le cas des séries formelles, un certain nombre des premiers termes des sé-
ries. Dans la suite de cette section, se donner une série formelle (resp. calculer une série
formelle) F signifiera donc se donner (resp. calculer) les N premiers termes de la série. On
notera F = P + O(X N ) où P est donc un polynôme de degré au plus N − 1. Les séries
formelles deviennent alors des polynômes de sorte que les algorithmes vus dans la section
précédente peuvent être utilisés pour calculer efficacement le produit de deux séries formelles.

On va maintenant s’intéresser au problème suivant : étant donnée F ∈ A[[X]] de terme

constant inversible dans un anneau effectif A, comment calculer efficacement son inverse ? La
méthode standard pour faire cette opération consiste à utiliser la division selon les puissances
croissantes.

16
Définition 2.2.2. Soit P, Q ∈ A[X] tels que Q(0) 6= 0. Alors pour tout entier n ∈ N, il
existe un unique couple de polynômes (Qn , Rn ) ∈ A[X]2 avec deg(Qn ) ≤ n tel que

P = Qn Q + X n+1 Rn .

Le polynôme Qn (resp. Rn ) s’appelle quotient (resp. le reste) à l’ordre n de la division selon

les puissances croissantes de P par Q.
Notons que comme pour une division (euclidienne) classique de polynômes, les quotients
et restes successifs dans la division selon les puissances croissantes peuvent s’obtenir en
« posant » la division mais, contrairement à la division classique, on ordonne les polynômes
selon leurs puissances croissantes et non décroissantes.
Exemple 2.2.1. Soient P = X + 1 et Q = X 2 + 1. On a alors

P = (1 + X − X 2 ) Q + X 3 (−1 + X),

de sorte que les quotients et restes respectifs d’ordre 2 dans la division selon les puissances
croissantes de P par Q sont donnés par 1 + X − X 2 et −1 + X.
Étant donnée une série formelle F = P + O(X N ) ∈ A[[X]] de terme constant égal à 1,
on peut effectuer la division selon les puissances croissantes de 1 par P (en ne conservant à
chaque étape n que les termes de degré au plus N − (n + 1) dans Rn ). L’inverse de F est
alors donné par QN −1 + O(X N ) où QN −1 est le quotient à l’ordre N − 1 de cette division.
Exemple 2.2.2. Considérons la série formelle F = 1 + 2 X − 3 X 2 + X 3 + O(X 4 ). En posant
la division selon les puissances croissantes de 1 par P = 1 + 2 X − 3 X 2 + X 3 , on trouve que
le quotient à l’ordre 4 − 1 = 3 est Q3 = 1 − 2 X + 7 X 2 − 21 X 3 . On peut alors vérifier que
Q3 + O(X 4 ) est bien l’inverse cherché car (1 + 2 X − 3 X 2 + X 3 ) (1 − 2 X + 7 X 2 − 21 X 3 ) =
1 + O(X 4 ).
Proposition 2.2.2. Soit F = P +O(X N ) ∈ A[[X]] de terme constant égal à 1. Calculer l’in-
verse G = Q + O(X N ) ∈ A[[X]] de F en utilisant la division selon les puissances croissantes
nécessite O(N 2 ) opérations dans A.
Démonstration. Exercice.
Nous allons maintenant voir un algorithme plus efficace pour effectuer cette opération en
utilisant l’itération de Newton. La méthode de Newton (très utilisée notamment en analyse
numérique) permet de calculer les zéros de fonctions réelles φ : R → R de classe C 1 en
choisissant un point x0 ∈ R puis en utilisant l’itération de Newton
φ(xk )
xk+1 = N (xk ) = xk − , (2.4)
φ0 (xk )
pour calculer les points suivants. Sous de bonnes hypothèses, cette suite va converger vers
un zéro de φ et la convergence sera quadratique, i.e., le nombre de décimales exactes double
à chaque itération.

17
 Soit F = P + O(X N ) ∈ A[[X]] de terme constant égal à 1. L’idée de l’algorithme consiste
à utiliser la méthode de Newton pour résoudre l’équation φ(G) = 0 où φ : A[[X]] → A[[X]]
est donnée par φ(G) = 1/G − F . En effet l’unique zéro de φ est alors l’inverse de la série
formelle F que l’on cherche. L’itération (2.4) devient alors :
1
Gk
−P
Gk+1 = N (Gk ) = Gk − = 2 Gk − G2k P = Gk (2 − P Gk ). (2.5)
− G12
k

Proposition 2.2.3. Soit F = P +O(X N ) ∈ A[[X]] de terme constant égal à 1. Alors la suite
définie par G0 = 1 et (2.5) vérifie qu’à chaque itération k, il existe un polynôme Rk ∈ A[X]
k k
tel que P Gk = 1 + X 2 Rk . Autrement dit, 1/F = Gk + O(X 2 ).
Démonstration. Raisonnons par récurrence sur k. La propriété est vraie pour k = 0 car
k
F (0) = P (0) = 1 et G0 = 1. Supposons maintenant P Gk = 1 + X 2 Rk . On a
k k k+1
P Gk+1 = 2 P Gk − G2k P 2 = 2 (1 + X 2 Rk ) − (1 + X 2 Rk )2 = 1 + X 2 (−Rk )2 ,

d’où le résultat.
Exemple 2.2.3. Reprenons l’exemple 2.2.2. L’itération de Newton appliquée au polynôme
P = 1 + 2 X − 3 X 2 + X 3 donne

G0 = 1, G1 = −X 3 + 3 X 2 − 2 X + 1, G2 = −21 X 3 + 7 X 2 − 2 X + 1 + O(X 4 ),

de sorte que l’on retrouve le résultat obtenu à l’exemple 2.2.2 à l’aide de division selon les
puissances croissantes.
On obtient alors l’algorithme suivant pour calculer l’inverse d’une série formelle :

Algorithme 4 Inversion d’une série formelle par l’itération de Newton

Entrées : N = 2k ∈ N∗ et F = P + O(X N ) ∈ A[[X]] de terme constant égal à 1.
Sorties : G = Q + O(X N ) ∈ A[[X]] tel que F G = 1 + O(X N ).
Si N = 1 Alors
- Retourner 1 ;
Sinon
- Poser G0 = 1 ;
- Calculer récursivement G1 , . . . , Gk (les N premiers termes) à l’aide de (2.5) ;
- Retourner Gk .
Finsi

Proposition 2.2.4. Soit M une fonction telle que le produit de deux polynômes de degré
au plus n puisse se calculer en O(M(n)) avec M(n) ≥ 2 M(n/2)1 . L’algorithme 4 calcule
l’inverse d’une série formelle à l’aide de l’opérateur de Newton en O(M(N )) opérations dans
A.
1
Notons que cette propriété est vérifiée par les algorithmes rapides donnés dans la section 2.1.

18
Démonstration. Soit I(N ) le nombre d’opérations arithmétiques effectuées par l’algorithme.
Connaissant Gk , le calcul de Gk+1 nécessite 2 M(N ) opérations dans A de sorte que nous
avons I(N ) ≤ I(N/2) + 2 M(N ). En itérant et en utilisant M(n) ≥ 2 M(n/2), on obtient
alors
  k−1   k−1
N X N X 1
I(N ) ≤ I + 2 M(N ) ≤ 2 M i
≤ 2 M(N ) i
≤ 4 M(N ),
2 i=0
2 i=0
2

d’où le résultat.
Lorsque N n’est pas une puissance de 2 on peut encore compléter avec des 0 (voir la
remarque 2.1.1) et obtenir la même estimation de complexité.

2.2.3 Application à la division euclidienne

Intéressons nous maintenant au problème de la division euclidienne de deux polynômes :
étant donnés deux polynômes A, B ∈ A[X] avec deg(A) = n, deg(B) = m et B unitaire,
nous voulons calculer (Q, R) ∈ A[X]2 tel que

A = Q B + R, deg(R) < deg(B).

La méthode classique consiste à poser la division. Notons

A = a X n + A0 , B = X m + B0 , deg(A0 ) < n, deg(B0 ) < m.

Si n < m, nous n’avons rien à faire sinon nous effectuons la soustraction

A − a X n−m B = a X n + A0 − a X n−m (X m + B0 ) = A0 − a X n−m B0 (2.6)

pour éliminer le terme de plus haut degré et obtenir un polynôme de degré strictement
inférieur à n. On itère ensuite ce procédé jusqu’à ce que le polynôme obtenu soit de degré
strictement inférieur à m.
Proposition 2.2.5. Soient A, B ∈ A[X] de degrés respectifs n et m. Effectuer la division
euclidienne de A par B en posant la division nécessite O(m (n − m)) opérations dans A.
Démonstration. À chaque itération (2.6), l’algorithme effectue O(m) opérations. Le nombre
maximal d’itérations étant égal à n − m + 1, on obtient le résultat annoncé.
Notons que le pire cas est obtenu pour n = 2 m où l’on a une complexité quadratique en
O(m2 ) opérations dans A.

Nous allons maintenant donner un algorithme rapide pour la division euclidienne basé
sur l’inversion de séries formelles (voir la sous-section précédente 2.2.2). L’idée consiste à
ré-écrire l’égalité A = Q B + R sous la forme
A R
=Q+ ,
B B
19
de sorte que Q peut s’obtenir en calculant le développement asymptotique de A/B à l’infini.
Notons que la contrainte deg(R) < deg(B) implique que R/B tend vers 0 à l’infini. Pour
réduire le problème à un calcul de division dans l’anneau des séries formelles, on ramène
l’infini en 0 en faisant le changement de variable X = 1/T . Si l’on suppose que A est de
degré n et B unitaire de degré m ≤ n, on obtient alors

T n A(1/T ) n−m T n R(1/T )

=T Q(1/T ) + m ,
T m B(1/T ) T B(1/T )

où on a multiplié par des puissances de T pour avoir des polynômes. Notons que le terme
constant de T m B(1/T ) est alors 1 (car B est unitaire) de sorte que nous allons pouvoir
l’inverser dans l’anneau A[[X]]. On obtient l’algorithme suivant :

Algorithme 5 Division euclidienne rapide

Entrées : A, B ∈ A[X] avec deg(A) = n, deg(B) = m ≤ n et B unitaire.
Sorties : Q, R ∈ A[X] avec A = Q B + R et deg(R) < deg(B).
- Calculer T n A(1/T ) et T m B(1/T ) ;
- Calculer l’inverse de T m B(1/T ) (les n − m + 1 premiers termes) ;
- Calculer le produit de T n A(1/T ) par 1/(T m B(1/T )) (les n − m + 1 premiers termes) ;
- En déduire Q ;
- Calculer R = A − Q B (les m premiers termes) ;
- Retourner (Q, R).

Proposition 2.2.6. Soit A, B ∈ A[X] avec deg(A) = n, deg(B) = m ≤ n et B unitaire.

Soit M une fonction telle que le produit de deux polynômes de degré au plus n puisse se
calculer en O(M(n)) opérations arithmétiques avec M(n) ≥ 2 M(n/2)2 . L’algorithme 5 ci-
dessus effectue la division euclidienne de A par B en O(M(n)) opérations dans A.

Démonstration. Regardons le nombre d’opérations nécessaire à chaque étape de l’algorithme.

Les étapes 1 et 4 ne nécessitent aucune opération car elles consistent simplement à inverser
l’ordre des coefficients des polynômes. D’après la proposition 2.2.4 la deuxième étape peut
s’effectuer en O(M(n − m)) opérations. L’étape 3 nécessite O(M(n − m)) opérations et
l’étape 5 requiert O(M(m)) opérations d’où le résultat.

2.3 Multiplication de matrices

Les algorithmes traitant de problèmes d’algèbre linéaire sont omniprésents en calcul formel.
La complexité arithmétique de ces algorithmes reposent très souvent sur celle de la multipli-
cation de matrices. Voir le théorème 2.3.1 à la fin de cette sous-section.
2
Notons que cette propriété est vérifiée par les algorithmes rapides donnés dans la section 2.1.

20
2.3.1 Algorithme naïf
Soient A = (ai,j )1≤i,j≤n et B = (bi,j )1≤i,j≤n deux matrices carrées de taille n à coefficients
dans un corps effectif K. On note A, B ∈ Mn (K). Le produit C = (ci,j )1≤i,j≤n des matrices
A et B est alors défini par
n
X
∀(i, j) ∈ {1, . . . , n}, ci,j = ai,k bk,j .
k=1

La multiplication naïve des matrices consiste à utiliser les formules ci-dessus pour calculer
chaque ci,j .

Proposition 2.3.1. Soient A, B ∈ Mn (K). L’algorithme de multiplication naïve calcule

C = A B en O(n3 ) opérations dans K.

Démonstration. Le calcul d’un coefficient ci,j nécessite n multiplications et n − 1 additions

dans K. Comme nous avons n2 coefficients ci,j à calculer, nous obtenons le résultat annoncé.

2.3.2 Algorithme de Strassen

Nous allons maintenant donner un algorithme plus efficace dû à Strassen qui utilise une tech-
nique analogue à celle de l’algorithme de Karatsuba pour la multiplication de polynômes.
L’idée consiste donc à voir comment gagner un produit dans K sur le produit de matrices
de taille 2 (quitte à effectuer plus d’additions) pour pouvoir tirer partie de ce gain dans une
stratégie diviser pour régner.

Considérons les deux matrices carrées de taille 2 suivantes :

   
a b x y
A= , B= , a, b, c, d, x, y, z, t ∈ K.
c d z t

Alors les entrées de la matrice C = A B = (ci,j )1≤i,j≤2 peuvent s’obtenir comme suit


 c1,1 = a x + b z = a (x + z) + (b − a) z,
c1,2 = a y + b t = d (y + t) + (d − a) (z − y) + (b − d) (z + t) − (b − a) z,


 c2,1 = c x + d z = a (x + z) + (d − a) (z − y) + (c − a) (x + y) − (c − d) y,
c2,2 = c y + d t = d (y + t) + (c − d) y,


ce qui ne requiert que 7 multiplications dans K. En utilisant cette observation, on obtient

l’algorithme suivant :

Proposition 2.3.2. L’algorithme 6 de Strassen calcule le produit de deux matrices carrées

de taille n en O(n2.81 ) opérations arithmétiques.

21
Algorithme 6 Algorithme de Strassen pour la multiplication de matrices
Entrées : A, B ∈ Mn (K) avec n = 2k .
Sorties : C ∈ Mn (K) telle que C = A B.
Si n = 1 Alors
- Retourner A B ;
Sinon
- Décomposer A et B par blocs de taille n2 :
   
a b x y
A= , B= , a, b, c, d, x, y, z, t ∈ M n2 (K);
c d z t

- À l’aide de 7 appels récursifs, calculer



 p1 = a (x + z),
p2 = d (y + t),




 p3 = (d − a) (z − y),


p4 = (b − d) (z + t),
p5 = (b − a) z,




p = (c − a) (x + y),

 6



p7 = (c − d) y;

- Calculer 

 c1,1 = p1 + p5 ,
c1,2 = p2 + p3 + p4 − p5 ,


 c2,1 = p1 + p3 + p6 − p7 ,
c2,2 = p2 + p7 ;


- Retourner C = (ci,j )1≤i,j≤2 .

Finsi

22
Démonstration. Soit S(n) le nombre d’opérations arithmétiques effectuées par l’algorithme
de Strassen. Un appel sur des matrices de taille n fait 7 appels récursifs sur des matrices de
taille n/2 ainsi qu’un nombre constant c = 18 d’additions de matrices de taille n/2. On a
donc S(n) ≤ 7 S(n/2) + c (n/2)2 avec S(1) = 1. On applique alors le théorème 1.4.1 avec
m = 7, p = 2, q = r = 4 et κ = 1 ce qui implique que
 n 2 
log2 (7) c ( 2 )
S(n) = O n = O(nlog2 (7) ),
nlog2 (4)

d’où le résultat (log2 (7) ≈ 2.807354922).

Remarque 2.3.1. Lorsque la taille n des matrices à multiplier n’est pas une puissance de
2, alors on rajoute artificiellement des lignes et des colonnes aux matrices pour s’y ramener
ce qui n’entraine pas de surcoût.

Nous concluons ce chapitre par le théorème suivant (admis dans le cadre de ce cours)
qui montre que la complexité du produit matriciel contrôle celui de toutes les opérations
d’algèbre linéaire :

Théorème 2.3.1. Soit K un corps et supposons que deux matrices carrées de taille n à
coefficients dans K peuvent être multipliées en O(MM(n)) opérations dans K. Alors, pour
une matrice A ∈ Mn (K) donnée, il est possible de calculer

• le déterminant det(A) de A,

• l’inverse A−1 de A (si A est inversible),

• la solution x ∈ Kn du système linéaire A x = b, pour tout b ∈ Kn (si A est inversible),

• le polynôme caractéristique de A,

• une décomposition LU de A (avec la matrice de permutation P associée si nécessaire),

• le rang de A et une forme échelonnée réduite de A,

• une base du noyau de A,

en Õ(MM(n)) opérations dans K, où la notation Õ(.) cache des facteurs logarithmiques.

23
24
Chapitre 3

Pgcd et résultant

3.1 Rappels : anneaux factoriels et anneaux euclidiens

Nous allons rappeler ici quelques résultats autour de l’existence et l’unicité (à unité près) du
pgcd (resp. ppcm) et d’une division euclidienne dans un anneau A avec en point de mire
les cas particuliers A = Z et A = K[X] qui nous intéresserons ensuite. Certains résultats
sont rappelés sans preuve. Il suffit de consulter un bon ouvrage d’algèbre pour retrouver les
démonstrations.

Soit A un anneau commutatif (i.e., ∀a, b ∈ A, a b = b a), unitaire (i.e., qui admet un
élément neutre 1 pour la multiplication) et intègre (i.e., A ne possède pas de diviseur de zéro
c’est-à-dire ∀a, b ∈ A, a b = 0 ⇒ (a = 0 ou b = 0)). On note A× l’ensemble des éléments
inversibles (aussi appelés unités) de A :

A× = {a ∈ A | ∃b ∈ A, a b = 1}.

Les exemples classiques pour ce cours seront A = Z, l’anneau des entiers relatifs qui vérifie
Z× = {−1, 1} et A = K[X], l’anneau des polynômes univariés à coefficients dans un corps K
avec K[X]× = K \ {0}.

3.1.1 Divisibilité, pgcd et ppcm

Définition 3.1.1. Soient a, b ∈ A. On dit que b divise a ou que a est divisible par b et on
note b | a s’il existe c ∈ A avec a = b c. On dit que a et b sont associés et on note a ∼ b si on
a à la fois b | a et a | b.

On peut montrer que ∼ est une relation d’équivalence et que a ∼ b ⇔ ∃u ∈ A× , a = u b.

Si on note (a) = {a x | x ∈ A} l’idéal principal engendré par a ∈ A, alors on a aussi a ∼ b ⇔
(a) = (b).

/ A× et si a = b c avec
Définition 3.1.2. Soit a ∈ A. On dit que a est irréductible si a ∈
b, c ∈ A implique b ∈ A× ou c ∈ A× .

25
 Si l’idéal (a) est un idéal premier (i.e., si A/(a) √
est intègre), alors a est irréductible
√
√ (Exercice : dans A = Z[i 5], l’idéal principal (1 + i 5) n’est
mais la réciproque est fausse
pas premier mais 1 + i 5 est irréductible). Dans A = Z, les éléments irréductibles sont
les entiers de la forme ±p, où p est un nombre premier, dans A = K[X] l’ensemble des
éléments irréductibles dépend du corps K. Si K = C, alors les éléments irréductibles sont
les polynômes de degré 1 mais si K = R alors il faut ajouter les polynômes de degré 2 de
discriminant strictement négatif.

Définition 3.1.3. Soient a, b, d ∈ A. On dit que d est un plus grand commun diviseur
(pgcd) de a et b et on note d = pgcd(a, b) si

• d | a et d | b,

• si c ∈ A vérifie aussi c | a et c | b, alors c | d.

Un pgcd n’existe
√ pas toujours dans un anneau commutatif
√ unitaire et intègre (Exercice :
dans A = Z[i 5], les éléments a = 6 et b = 2 (1 + i 5) n’ont pas de pgcd) mais s’il existe
alors il est unique modulo la relation d’équivalence ∼, i.e., si d1 et d2 sont deux pgcd de a
et b, alors il existe u ∈ A× tel que d1 = u d2 . Deux éléments a et b d’un anneau A sont dits
premiers entre eux (ou étrangers) si pgcd(a, b) = 1.

Définition 3.1.4. Soient a, b, m ∈ A. On dit que m est un plus petit commun multiple
(ppcm) de a et b et on note m = ppcm(a, b) si

• a | m et b | m,

• si c ∈ A vérifie aussi a | c et b | c, alors m | c.

Tout comme le pgcd, si le ppcm existe, alors il est unique.

Pour avoir l’existence du pgcd, nous devons nous placer dans le cas d’un anneau factoriel.

Définition 3.1.5. Soit A un anneau commutatif unitaire et intègre. L’anneau A est dit
factoriel si :

• tout élément non-nul est soit inversible soit produit d’éléments irréductibles,

• Si p1 · · · pn = q1 · · · qm avec les pi et les qi irréductibles dans A, alors n = m et il

existe une permutation σ de {1, . . . , n} telle que, pour tout i ∈ {1, . . . , n}, pi ∼ qσ(i) ,
i.e., il existe ui ∈ A× avec pi = ui qσ(i) .
√
Les anneaux Z et K[X] sont factoriels mais pas Z[i 5].

Soit A un anneau factoriel, un système représentatif d’éléments irréductibles de A est une

partie P de A formée d’éléments irréductibles tels que :

• pour tout q ∈ A irréductible, il existe p ∈ P tel que q ∼ p,

26
 • les éléments de P sont deux à deux non associés.

Par exemple, pour A = Z, on peut peut prendre pour P l’ensemble des nombres premiers et
pour A = K[X], on peut prendre l’ensemble des polynômes irréductibles unitaires. Une fois
un système représentatif d’éléments irréductibles P de A fixé, tout élément de a ∈ A s’écrit
de manière unique
n
Y
a=u pei i , p1 , . . . , pn ∈ P avec ∀i 6= j, pi 6= pj , e1 , . . . , en ∈ N, u ∈ A× . (3.1)
i=1

Par conséquent, dans un anneau factoriel, si un élément irréductible p divise un produit a b,

alors p | a ou p | b de sorte que p irréductible implique que l’idéal (p) est premier. Dans un
anneau factoriel, nous avons donc l’équivalence : a irréductible ⇔ (a) premier.

Théorème 3.1.1. Soit A un anneau factoriel. Alors deux éléments non nuls a et b de A
admettent un pgcd unique à multiplication par unité près.

Avec les notations de (3.1), écrivons a = u ni=1 pei i et b = v ni=1 pfi i . Mon-
Q Q
Démonstration. Q
trons que d = ni=1 pgi i avec, pour i = 1, . . . , n, gi = min(ei , fi ) est un pgcd de a et b.
Q nous avons clairement d | a et d | b. De plus si c | a et c | b, alors nécessairement
En effet,
c = w ni=1 phi i avec w ∈ A× et pour i = 1, . . . , n, hi ∈ N (hi ≤ ei et hi ≤ fi ) de sorte que
c | d car pour i = 1, . . . , n, gi = min(ei , fi ).

On obtient un résultat analogue pour l’existence du ppcm dans un anneau factoriel

(prendre gi = max(ei , fi ) au lieu de min(ei , fi ) dans la preuve).

Le théorème suivant sera admis pour ce cours :

Théorème 3.1.2. Si un anneau A est factoriel, alors A[X] est factoriel.

Par récurrence, on voit donc que si A est factoriel, alors l’anneau des polynômes en n
variables A[X1 , . . . , Xn ] est factoriel.

3.1.2 Division euclidienne

Définition 3.1.6. Soit A un anneau commutatif unitaire et intègre. Alors l’anneau A est
dit euclidien s’il existe une application φ : A \ {0} → N appelée stathme euclidien telle que :

• pour tout a, b ∈ A \ {0}, φ(a b) ≥ φ(a) ;

• pour tout a, b ∈ A\{0}, il existe q, r ∈ A tels que a = q b+r avec r = 0 ou φ(r) < φ(b).

Les éléments q et r de A sont alors appelés respectivement quotient et reste de la division

euclidienne de a par b.

27
 L’anneau A = Z est euclidien en prenant pour stathme la valeur absolue, i.e., φ : Z\{0} →
N, r 7→ |r|. Notons que dans ce cas, le couple (q, r) n’est pas unique (8 = 2 × 3 + 2 =
3 × 3 − 1). L’anneau A = K[X] est euclidien en prenant pour stathme la fonction degré,
i.e., φ : K[X] \ {0} → N, P 7→ degX (P ). La fonction degré vérifie de plus l’inégalité
degX (P1 − P2 ) ≤ max(degX (P1 ), degX (P2 )), pour tous polynômes P1 , P2 ∈ K[X] ce qui
implique alors l’unicité du quotient et du reste de la division euclidienne dans K[X].

Théorème 3.1.3. Soient A un anneau et A, B ∈ A[X] avec B 6= 0 et de coefficient de tête

inversible dans A. Alors il existe un unique couple (Q, R) ∈ A[X]2 tel que A = Q B + R et
deg(R) < deg(B).

Si A est un corps K, alors nous avons le même résultat sans l’hypothèse sur le coefficient
de tête de B.

Théorème 3.1.4. Tout anneau euclidien est principal, i.e., tout idéal d’un anneau euclidien
est engendré par un seul élément.

Démonstration. Soit A un anneau euclidien muni du stathme φ et soit I un idéal non-nul

de A. Soit a ∈ I tel que φ(a) soit minimal, i.e., pour tout x ∈ I \ {0}, φ(x) ≥ φ(a). On
va montrer que I = (a). Soit b ∈ I \ {0}, alors il existe (q, r) ∈ A2 tel que b = q a + r et
r = 0 ou φ(r) < φ(a). Supposons r 6= 0. On a alors r = b − q a ∈ I et φ(r) < φ(a) ce qui
contredit le fait que φ(a) soit minimal. Par conséquent r = 0 et b = q a ∈ (a) ce qui termine
la preuve.
On peut aussi montrer que tout anneau principal est factoriel (admis pour ce cours) de
sorte que l’on a les implications suivantes :

A euclidien ⇒ A principal ⇒ A factoriel.

Certaines propriétés des anneaux factoriels restent vraies dans un anneau principal.

Proposition 3.1.1. Soit A un anneau principal. Alors toute famille d’éléments de A admet
un pgcd et un ppcm.

Démonstration. Soient a, b ∈ A et considérons l’idéal I = (a) + (b) = {x a + y b, x, y ∈ A}.

Comme A est principal, il existe d ∈ A tel que I = (d). On vérifie alors que d est un pgcd
de a et b. De la même façon, on considère l’idéal J = (a) ∩ (b) et on montre que si J = (m),
alors m est un ppcm de a et b.
Une conséquence directe de cette proposition est l’existence de l’égalité de Bézout dans
un anneau principal.

Théorème 3.1.5. Soient A un anneau principal, a, b ∈ A et d = pgcd(a, b) un pgcd de a

et b. Alors il existe (u, v) ∈ A2 tel que u a + v b = d.

En particulier, si a et b sont premiers entres eux, alors il existe (u, v) ∈ A2 tel que
u a + v b = 1.

28
Théorème 3.1.6 (Gauss). Soient A un anneau principal et a, b, c ∈ A. Si c divise le produit
a b et si c est premier avec b, alors c divise a.
Démonstration. Comme c et b sont premiers entres eux, alors il existe (u, v) ∈ A2 tel que
u c + v b = 1. En multipliant cette égalité par a, il vient u a c + v a b = a. Maintenant comme
c divise a b, il existe q ∈ A tel que a b = q c de sorte que l’on obtient (u a + v q) c = a ce qui
montre que c divise a.
Corollaire 3.1.1. Soient A un anneau principal et a, b, p ∈ A avec p irréductible. Si p divise
a b, alors p divise a ou p divise b.
Finalement, on remarque que l’équivalence (a irréductible ⇔ (a) premier) reste vérifiée
dans un anneau principal.

3.2 Calcul du pgcd et des coefficients de Bézout

3.2.1 Dans un anneau euclidien : le cas des entiers
Soit A un anneau effectif euclidien muni d’un stathme φ dans lequel on suppose que l’on sait
effectuer la division euclidienne c’est-à-dire qu’étants donnés a et b dans A, on sait détermi-
ner (q, r) ∈ A2 tel que a = q b + r avec r = 0 ou φ(r) < φ(b). On note alors q = quo(a, b) et
r = a mod b ou r = rem(a, b) le quotient et le reste de la division euclidienne de a par b.

Le problème du calcul de pgcd dans A est alors le suivant : étant donnés a , b ∈ A,

calculer pgcd(a, b), i.e., un pgcd de a et b.

Pour ce faire, on commence par effectuer la division euclidienne de a par b pour obtenir q
et r tels que a = q b+r avec r = 0 ou φ(r) < φ(b). Si r = 0, alors b divise a et pgcd(a, b) = b.
Sinon on a pgcd(a, b) ∼ pgcd(b, r). On effectue alors la division euclidienne de b par r et
ainsi de suite. En notant r0 = a, r1 = b, on définit ainsi la suite des restes successifs :
ri−2 = qi ri−1 + ri , i ≥ 2, ri = 0 ou φ(ri ) < φ(ri−1 ). (3.2)
Comme le stathme φ est à valeur dans N, il existe k ∈ N tel que rk+1 = 0 et rk 6= 0.
On a alors rk = pgcd(a, b). En effet rk+1 = 0 implique rk divise rk−1 et par conséquent
rk = pgcd(rk−1 , rk ) ∼ pgcd(rk−2 , rk−1 ) ∼ · · · ∼ pgcd(r1 , r2 ) ∼ pgcd(r0 , r1 ) = pgcd(a, b).
On obtient alors l’algorithme d’Euclide 7 pour le calcul du pgcd dans un anneau euclidien :
Cet algorithme peut aussi s’écrire à l’aide d’appels récursifs (Algorithme 8).
Comme nous l’avons vu précédemment, en général, il n’y a pas unicité du couple (q, r) dans
la division euclidienne donc ces algorithmes peuvent produire plusieurs suites différentes de
restes successifs. De plus le nombre d’étapes avant de tomber sur un reste nul peut varier :
Exemple 3.2.1. Soit A = Z, a = 30 et b = 18. La suite des restes peut alors s’écrire

 30 = 1 × 18 + 12,
18 = 1 × 12 + 6,
12 = 2 × 6 + 0,


29
Algorithme 7 Algorithme d’Euclide
Entrées : A anneau euclidien, a, b ∈ A.
Sorties : pgcd(a, b) (un pgcd de a et b).
r0 = a, r1 = b, i = 1 ;
Tant que ri 6= 0 Faire
ri+1 = ri−1 mod ri ;
i = i + 1;
Fin tant que
Retourner ri−1 .

Algorithme 8 Algorithme d’Euclide récursif

Entrées : A anneau euclidien, a, b ∈ A.
Sorties : pgcd(a, b).
r = a mod b ;
Si r = 0 Alors
Retourner b ;
Sinon
Calculer pgcd(b, r) à l’aide d’un appel récursif.
Finsi

d’où pgcd(30, 18) = 6. Mais on a aussi


30 = 2 × 18 + (−6),
18 = (−3) × (−6) + 0,
qui nous donne pgcd(30, 18) = −6.
Nous nous intéressons maintenant au calcul des coefficients de Bézout : étant donnés
a , b ∈ A, calculer (u, v) ∈ A2 tel que u a + v b = pgcd(a, b). Pour ceci on va devoir garder
trace des quotients qi calculés à chaque étape de l’algorithme d’Euclide. Définissons la suite
des triplets successifs Ti = (ri , ui , vi ) ∈ A3 définis par T0 = (a, 1, 0), T1 = (b, 0, 1) et
Ti = Ti−2 − qi Ti−1 , i ≥ 2, qi = quo(ri−2 , ri−1 ).
Proposition 3.2.1. Avec les notations ci-dessus, pour tout i, on a ri = ui a + vi b.
Démonstration. Récurrence sur i.
En notant toujours k l’entier tel que rk+1 = 0 et rk 6= 0, on a alors rk = pgcd(a, b) et
uk a + vk b = rk . On obtient alors l’algorithme d’Euclide étendu pour le calcul du pgcd et
des coefficients de Bézout.
Exemple 3.2.2. Reprenons l’exemple (3.2.1) où a = 30 et b = 18. En utilisant la suite de
divisions euclidiennes 
 30 = 1 × 18 + 12,
18 = 1 × 12 + 6,
12 = 2 × 6 + 0,


30
Algorithme 9 Algorithme d’Euclide étendu
Entrées : A anneau euclidien, a, b ∈ A.
Sorties : pgcd(a, b) (un pgcd de a et b) et (u, v) ∈ A2 tel que u a + v b = pgcd(a, b).
r0 = a, u0 = 1, v0 = 0, r1 = b, u1 = 0, v1 = 1, i = 1 ;
Tant que ri 6= 0 Faire
qi+1 = quo(ri−1 , ri ) ;
ri+1 = ri−1 mod ri ;
ui+1 = ui−1 − qi+1 ui ;
vi+1 = vi−1 − qi+1 vi ;
i = i + 1;
Fin tant que
Retourner ri−1 et (ui−1 , vi−1 ).

on obtient le tableau de valeur suivant :

Itération i ri ui vi ri+1 ui+1 vi+1 qi+1
0 30 1 0 18 0 1
1 18 0 1 12 1 −1 1
2 12 1 −1 6 −1 2 1
3 6 −1 2 0 3 −5 2

On a donc (−1) × 30 + 2 × 18 = 6.

3.2.2 Le théorème des restes chinois

Le théorème des restes chinois permet de reconstruire un nombre entier à partir de ses valeurs
modulo certains autres entiers. Il est très utile pour développer des méthodes modulaires en
calcul formel. Le problème s’énonce de la manière suivante : étant donnés quatre nombres
entiers a, b, m et n avec m n ≥ 1, déterminer tous les entiers x ∈ Z tels que :

x ≡ a mod m,
(S) :
x ≡ b mod n.

Théorème 3.2.1. Avec les notations précédentes, si pgcd(m, n) = 1, alors (S) admet une
et une seule solution modulo le produit m n.

Démonstration. On considère l’application ϕ définie par :

ϕ: Z → (Z/m Z) × (Z/n Z),

x 7→ (x mod m, x mod n).

On vérifie alors que ϕ est un morphisme d’anneau. Son noyau ker(ϕ) est donné par

ker(ϕ) = {x ∈ Z | x ≡ 0 mod m et x ≡ 0 mod n} = m n Z,

31
puisque pgcd(m, n) = 1. ϕ induit donc une injection ϕ̃ : Z/(m n Z) → (Z/m Z) × (Z/n Z)
qui est une bijection car Card(Z/(m n Z)) = Card((Z/m Z) × (Z/n Z)) = m n ce qui prouve
le résultat.
D’un point de vue effectif, la solution de (S) peut être calculée à l’aide de l’algorithme
d’Euclide étendu. En effet, en appliquant l’algorithme d’Euclide étendu à m et n premiers
entre eux, on trouve deux entiers u et v tels que u m + v n = 1. On a alors
 
u m ≡ 0 mod m, v n ≡ 1 mod m,
u m ≡ 1 mod n, v n ≡ 0 mod n,
de sorte que x0 = a v n + b u m est une solution particulière de (S). La solution générale de
(S) s’écrit x0 + k m n, k ∈ Z.
Exemple 3.2.3. Prenons a = 2, b = 1, m = 3 et n = 5. En appliquant l’algorithme d’Euclide
étendu à 3 et 5 premiers entre eux, on trouve 2 × 3 + (−1) × 5 = 1. Une solution de (S) est
alors donné par x0 = 2 × (−1) × 5 + 1 × 2 × 3 = −4 qui vérifie bien

x0 ≡ 2 mod 3,
x0 ≡ 1 mod 5.
Le théorème 3.2.1 se généralise au cas où le système (S) possède plus de deux congruences.
Corollaire 3.2.1. Soient m1 , . . . , ms des entiers positifs deux à deux premiers entre Qs eux,
i.e., pgcd(mi , mj ) = 1, pour tout i, j = 1, . . . , s, i 6= j. Soit M = m1 × · · · × ms = i=1 mi
le produit des mi . Alors pour tout s-uplet d’entiers a1 , . . . , as le système de congruences

 x ≡ a1 mod m1 ,

(S) : ..
 .
 x ≡ a mod m ,
s s

admet une et une seule solution modulo M .

Comme précédemment, le système de congruences (S) peut se résoudre en utilisant l’algo-
rithme d’Euclide étendu. Pour i = 1, . . . , s, posons Mi = M/mi qui vérifie pgcd(mi , Mi ) = 1.
En appliquant l’algorithme d’Euclide étendu, on peut Ps calculer un entier wi tel que wi Mi ≡ 1
mod mi . Alors x0 = a1 w1 M1 + · · · + as ws Ms = i=1 ai wi Mi est une solution particulière
de (S). La solution générale de (S) s’écrit x0 + k M, k ∈ Z. Cette méthode est appelée
méthode de Lagrange pour la résolution du système de congruences (S).
Exemple 3.2.4. Considérons le système de congruences

 x ≡ 3 mod 11,
(S) : x ≡ 2 mod 4,
x ≡ 7 mod 15.


Les entiers 11, 4 et 15 sont bien premiers entre eux de sorte que la solution est unique modulo
M = 660. On a M1 = M/11 = 60, M2 = M/4 = 165 et M3 = M/15 = 44.

32
 1. On cherche w1 tel que 60 w1 ≡ 1 mod 11, i.e., 5 w1 ≡ 1 mod 11. L’algorithme d’Eu-
clide étendu donne (−2) × 5 + 1 × 11 = 1 de sorte que w1 = −2 convient ;

2. On cherche w2 tel que 165 w2 ≡ 1 mod 4, i.e., w2 ≡ 1 mod 4. L’algorithme d’Euclide

étendu donne 1 × 1 + 0 × 4 = 1 de sorte que w2 = 1 convient ;

3. On cherche w3 tel que 44 w3 ≡ 1 mod 15, i.e., (−1) w3 ≡ 1 mod 15. L’algorithme
d’Euclide étendu donne (−1) × (−1) + 0 × 15 = 1 de sorte que w3 = −1 convient.

On a donc trouvé la solution particulière x0 = 3×(−2)×60+2×1×165+7×(−1)×44 = −338.

La solution générale est x = −338 + k × 660, k ∈ Z. La solution positive modulo M = 660
est alors −338 + 660 = 322.

Il existe un moyen d’obtenir la solution positive la plus petite en utilisant la méthode de

Garner décrite ci-dessous et qui est basée sur le théorème suivant :

Théorème 3.2.2. Soient m1 , . . . , ms des entiers positifs (pas nécessairement distincts).

Alors l’application

ψ : J0, m1 J× · · · × J0, ms J → J0, m1 · · · ms J,

(r1 , . . . , rs ) 7→ r1 + r2 m1 + r3 m1 m2 + · · · + rs m1 · · · ms−1 ,

est bijective.

Démonstration. Exercice.
Considérons le système de congruences suivant :

 x ≡ a1 mod
 m1 ,
(S) : ..
 .
 x ≡ a mod ms ,
s

où les mi sont deux à deux premiers entre eux et cherchons la solution de (S) appartenant
à J0, m1 · · · ms J. Cette solution s’écrit

x = r1 + r2 m1 + r3 m1 m2 + · · · + rs m1 · · · ms−1 , ri ∈ J0, mi J, i = 1, . . . , s.

Les ri se calculent alors par récurrence comme suit :

• r1 ≡ a1 mod m1 est le reste de la division euclidienne de a1 par m1 ;

• Supposons que l’on connaisse r1 , . . . , rk−1 . On cherche rk tel que x ≡ ak mod mk

donc r1 + r2 m1 + · · · + rk−1 m1 · · · mk−2 + rk m1 · · · mk−1 ≡ ak mod mk . Puisque
pgcd(mk , m1 · · · mk−1 ) = 1, l’algorithme d’Euclide étendu nous donne (uk , vk ) tel que
uk (m1 · · · mk−1 ) + vk mk = 1. On pose alors

rk = uk (ak − (r1 + r2 m1 + · · · + rk−1 m1 · · · mk−2 )) mod mk .

33
Exemple 3.2.5. On considère le système de congruence

x ≡ 1000 mod 1891,
(S) :
x ≡ 501 mod 2499.

Cherchons x sous la forme x = r1 + r2 × 1891 avec r1 ∈ J0, 1891J et r2 ∈ J0, 2499J.

• L’entier positif r1 s’obtient comme reste de la division euclidienne de 1000 par 1891 de
sorte que r1 = 1000 ;

• Calculons r2 . L’algorithme d’Euclide étendu appliqué à 1891 et 2499 nous fournit u

et v tels que u × 1891 + v × 2499 = 1. On trouve u = −1007 et v = 762. D’où
r2 = −1007 (501 − 1000) mod 2499 = 194.

Finalement, on a x = 1000 + 194 × 1891 = 367854.

Remarquons que sur cet exemple, la solution particulière fournit par la méthode de Lagrange
est x0 = 1000 × 762 × 2499 + 501 × (−1007) × 1891 = 950215263.

Le théorème des restes chinois peut se généraliser à tout anneau commutatif unitaire.

Théorème 3.2.3. Soient A un anneau commutatif unitaire, I1 , . . . , Is des idéaux de A tels

que Ij + Ik = A pour tout j, k = 1, . . . , s, j 6= k. Alors on a l’isomorphisme

A/ ∩sj=1 Ij ∼


= (A/I1 ) × · · · × (A/Is ) .

= 1, . . . , s, Ij = nj A, nj ∈ A
Dans le cas particulier où A est principal, alors pour tout j Q
avec les nj premiers entre eux deux à deux et si on note n = si=1 ni , on a :

A/(n A) ∼
= (A/n1 A) × · · · × (A/ns A) .

En particulier, lorsque A = K[X] est l’anneau des polynômes univariés à coefficients

dans un corps K. Si P1 , . . . , Ps sont des polynômes de K[X] premiers entre eux deux à deux,
alors pour tout s-uplet de polynômes (Q1 , . . . , Qs ) ∈ K[X]s , il existe un polynôme P , unique
modulo le produit P1 · · · Ps , tel que

 P ≡ Q1
 mod P1 ,
..
 .
 P ≡ Q mod Ps .
s

Le problème d’interpolation de Lagrange peut alors se voir comme un cas particulier où les
polynômes Pi sont de la forme X − xi avec les xi ∈ K deux à deux distincts et les Qi sont
des constantes yi ∈ K.

34
3.2.3 Dans un anneau de polynômes
Nous nous plaçons maintenant dans le cas d’un anneau de polynômes A[X] à coefficients
dans un anneau effectif A. Étant donnés A, B ∈ A[X], on cherche à déterminer pgcd(A, B).
Lorsque A est un corps K, l’anneau K[X] étant euclidien pour le stathme deg, nous pouvons
appliquer l’algorithme d’Euclide 7 pour calculer pgcd(A, B).

Proposition 3.2.2. L’algorithme d’Euclide 7 calcule un pgcd de deux polynômes A et B

à coefficients dans un corps K en O(deg(A) deg(B)) opérations dans K.

Démonstration. Supposons deg(A) ≥ deg(B). D’après la sous-section 2.2.3, le calcul de P

mod Q peut se faire par la méthode naïve en 2 deg(Q) (deg(P )−deg(Q)+1) opérations dans
P En notant Ri les restes successifs, le coût de l’algorithme d’Euclide est donc majoré par
K.
Pi−1 ) − deg(Ri ) + 1). En majorant tous les deg(Ri ), i ≥ 1 par deg(B)
i≥1 2 deg(Ri ) (deg(R
on obtient 2 deg(B) i≥1 (deg(Ri−1 ) − deg(Ri ) + 1) = 2 deg(B) (deg(R0 ) + deg(B)) puisque
dans le pire cas le nombre d’étapes est égal à deg(B). D’où le résultat puisque R0 = A.

Exemple 3.2.6. Considérons les polynômes A = X 8 + X 6 − 3 X 4 − 3 X 3 + 8 X 2 + 2 X − 5

et B = 3 X 6 + 5 X 4 − 4 X 2 − 9 X + 21 de Z[X] ⊂ Q[X]. La suite des restes obtenue en
appliquant l’algorithme d’Euclide est alors donnée par :

R2 = − 59 X 4 + 19 X 2 − 13 ,




R = − 117 X 2 − 9 X + 441

 ,
 3 25 25


R4 = 233150
19773
X − 102500
6591
,

R5 = − 1288744821




 543589225
,

R6 = 0,


de sorte qu’un pgcd de A et B est donné par R5 = − 1288744821

543589225
. Notons que R5 ∼ 1 donc A
et B sont premiers entre eux.

L’exemple précédent illustre deux problèmes apparaissant lorsque l’on applique l’algo-
rithme d’Euclide à deux polynômes :

1. même si les polynômes de départ sont à coefficients entiers, des dénominateurs appa-
raissent dans les calculs ;

2. la taille des entiers apparaissant dans les restes successifs augmente significativement.

L’apparition des dénominateurs provient du fait que pour chaque division euclidienne nous
devons inverser le coefficient de tête du polynôme par lequel on divise (voir Théorème 3.1.3).
Dans notre exemple, dès la première étape le coefficient de tête de B est 3 donc nous devons
l’inverser pour trouver R2 qui aura donc des dénominateurs dans ses coefficients. Notons
que même si nous partons de deux polynômes unitaires, des dénominateurs peuvent quand
même apparaitre au cours de l’algorithme comme l’illustre l’exemple suivant :

35
Exemple 3.2.7. Considérons les polynômes A = X 4 −13 X 3 +2 X 2 −X−1 et B = X 2 −X−1.
La suite des restes est alors donnée par :

R = −22 X − 10,

 2

41
R3 = − 121


R4 = 0.

Afin d’éviter l’apparition des dénominateurs1 , une idée consiste à utiliser la pseudo-
division euclidienne.

Définition 3.2.1. Soient A et B deux polynômes à coefficients dans un anneau A et b le

coefficient de tête de B. Le pseudo-reste R̃ de A par B noté prem(A, B) est défini par :

bdeg(A)−deg(B)+1 A = Q̃ B + R̃, Q̃, R̃ ∈ A[X], deg(R̃) < deg(B).

On peut donc obtenir un algorithme d’Euclide « sans dénominateurs » en remplaçant les

restes par des pseudo-restes dans l’algorithme d’Euclide.

Algorithme 10 Algorithme d’Euclide par pseudo-restes pour le calcul du pgcd

Entrées : A, B ∈ A.
Sorties : pgcd(A, B) (un pgcd de A et B).
R0 = A, R1 = B, i = 1 ;
Tant que Ri 6= 0 Faire
Ri+1 = prem(Ri−1 , Ri ) ;
i = i + 1;
Fin tant que
Retourner Ri−1 .

Exemple 3.2.8. Reprenons l’exemple 3.2.6 avec la pseudo-division euclidienne. La suite des
pseudo-restes est alors :

R2 = −15 X 4 + 3 X 2 − 9,




R = 15795 X 2 + 30375 X − 59535,


 3


R4 = 1254542875143750 X − 1654608338437500,





 R5 = 12593338795500743100931141992187500.

R6 = 0.


On remarque qu’en appliquant la pseudo-division euclidienne, on fait certes disparaitre les

dénominateurs mais la taille des entiers apparaissant dans les coefficients des pseudo-restes
1
Notons qu’en calcul formel, on évite autant que possible les calculs sur les nombres rationnels pour
lesquels l’addition elle-même requiert des calculs de multiplications et éventuellement des calculs de pgcd.

36
croît très rapidement donc que cela n’est pas très satisfaisant en pratique. Pour diminuer ce
phénomène de croissance des coefficients, une possibilité consiste à diviser à chaque étape le
pseudo-reste par le pgcd de ses coefficients (aussi appelé contenu) pour obtenir une suite de
polynômes primitifs, i.e., dont le contenu vaut 1. On parle alors de suite des pseudo-restes
primitifs.

Exemple 3.2.9. Dans l’exemple 3.2.6, la suite des pseudo-restes primitifs est :
4 2

 R2 = −5 X + X − 3,


R = 13 X 2 + 25 X − 49,


 3


R4 = 4663 X − 6150,





 R5 = 1.

R6 = 0.


Cependant en pratique, cela rajoute un calcul de pgcd des coefficients à chaque étape ce
qui est trop coûteux. Notons finalement qu’il existe un moyen de modifier légèrement cette
approche en prévoyant les facteurs communs qui apparaissent dans les coefficients de sorte à
limiter leur croissance (sans calculer de pgcd dans A) : c’est l’algorithme des sous-résultants
que nous ne détaillerons pas dans ce cours.

Une autre façon d’éviter la croissance des coefficients est d’utiliser une méthode modulaire
c’est-à-dire de réduire modulo un ou plusieurs nombres premiers p et d’appliquer l’algorithme
d’Euclide dans Fp = Z/(p Z). Pour construire un algorithme modulaire pour le calcul du
pgcd, plusieurs choses doivent être prises en considération. Premièrement, si on note D =
pgcd(P, Q), P , Q et D les réductions de P , Q et D modulo un nombre premier p, alors nous
n’avons pas toujours D = pgcd(P , Q). En effet si on note P1 et Q1 les cofacteurs de P et
Q, i.e., P = P1 D et Q = Q1 D, alors on a nécessairement P = P1 D et Q = Q1 D de sorte
que D | pgcd(P , Q) mais pour avoir l’égalité il faut (et il suffit) que P1 et Q1 soient premiers
entre eux dans Z/(p Z) ce qui n’est pas toujours vrai comme l’illustre l’exemple suivant.

Exemple 3.2.10. Soient P = (X 3 + X + 3) (X + 1) et Q = X (X + 1) de sorte que

D = pgcd(P, Q) = X + 1. Pour p = 3, on a P = (X 3 + X) (X + 1) et Q = X (X + 1) de
sorte que pgcd(P , Q) = X (X + 1) 6= D = X + 1. Cela provient du fait que les cofacteurs
X 3 + X + 3 et X ne sont pas premiers entre eux dans Z/(3 Z).

Ceci implique que tout nombre premier ne sera pas « bon » mais qu’il n’y aura qu’un
nombre fini de « mauvais » nombres premiers (avec les notations précédentes, ceux qui
divisent le résultant de P1 et Q1 - cf Sous-section 3.3 suivante). De plus, pour choisir judi-
cieusement ces nombres premiers de sorte à pouvoir reconstruire D à partir Pde D, il nous faut
d
avoir une idée de la taille des coefficients de D que l’on cherche. Pour P = i=1 pi X i ∈ Z[X],
on note ||P ||∞ = max0≤i≤d |pi | la norme infinie de P (aussi appelée hauteur de P ) et
||P ||1 = di=0 |pi | la norme 1 de P . Nous avons alors la borne suivante :
P

37
Lemme 3.2.1 (Borne de Mignotte). Soit P = di=0 pi X i ∈ Z[X] un polynôme primitif,
P
i.e., le pgcd de ses coefficients pi (aussi appelé contenu de P ) vaut 1. Alors si un polynôme
Q ∈ Z[X] de degré n divise P , on a :
√
||Q||∞ ≤ ( d + 1) 2n ||P ||∞ .

Démonstration. Admis pour ce cours.

À partir de ce résultat, deux stratégies s’offrent à nous : soit on choisit un seul « grand »
nombre premier pour calculer directement le pgcd dans Z[X] soit on choisit plusieurs « pe-
tits » nombres premiers, on déroule l’algorithme d’Euclide modulo chacun de ces nombres
premiers p et on reconstruit le pgcd dans Z[X] à l’aide du théorème des restes chinois.
Cette deuxième stratégie nécessite d’appliquer l’algorithme d’Euclide plusieurs fois mais les
coefficients qui interviennent sont plus petits ce qui est un avantage non négligeable pour des
polynômes de départ de degré et de norme infinie élevés. Les algorithmes 11 et 12 montrent
comment chacune de ces stratégies peut être mise en oeuvre pour obtenir un algorithme
modulaire calculant le pgcd de deux polynômes.
Notons que dans  les deux algorithmes la condition utilisée
 pour vérifier que l’on a bien ob-
tenu le pgcd, i.e., ||P̃ ||1 ||S||1 ≤ B et ||Q̃||1 ||S||1 ≤ B est en fait une condition nécessaire
et suffisante pour que pp(S) soit bien le pgcd de P et Q.

Algorithme 11 Algorithme modulaire pour le pgcd utilisant un grand premier

Entrées : P, Q ∈ Z[X] primitifs tels que deg(P ) = n ≥ deg(Q).
Sorties : pgcd(P, Q) (un pgcd de P et Q).
- Soit A un majorant de ||P ||∞ et ||Q||∞ ;
- Calculer le√pgcd b des coefficients dominants de P et Q ;
- Soit B = ( n + 1) 2n A b ;
- Choisir un nombre premier p entre 2 B et 4 B ;
- Appliquer l’algorithme d’Euclide pour calculer le polynôme unitaire R ∈ Z[X] tel que
||R||∞ ≤ p/2 et tel que la réduction modulo p de R est un pgcd des réductions P et Q
de P et Q modulo p ;
- Calculer
 S, P̃ et Q̃ tels que S ≡ b R mod
 p, P̃ S ≡ b P mod p, Q̃ S ≡ b Q mod p ;
Si ||P̃ ||1 ||S||1 ≤ B et ||Q̃||1 ||S||1 ≤ B Alors
- Retourner pp(S) la partie primitive de S, i.e., S divisé par son contenu ;
Sinon
- Choisir un autre nombre premier à l’étape 4 et recommencer.
Finsi

Étant donnés deux polynômes A, B ∈ K[X], l’algorithme d’Euclide étendu 9 peut être
appliqué pour calculer pgcd(A, B) ainsi que deux polynômes U, V ∈ K[X] tels que

U A + V B = pgcd(A, B), deg(U pgcd(A, B)) < deg(B), deg(V pgcd(A, B)) < deg(A).

38
Algorithme 12 Algorithme modulaire pour le pgcd utilisant plusieurs petits premiers
Entrées : P, Q ∈ Z[X] primitif tels que deg(P ) = n ≥ deg(Q).
Sorties : pgcd(P, Q) (un pgcd de P et Q).
- Soit A un majorant de ||P ||∞ et ||Q||∞ ;
- Calculer le pgcd b des coefficients dominants √ de P et Q ;
- Soit k = d2 log2 ((n + 1) b A )e, B = ( n + 1) 2n A b et l = dlog2 (2 B + 1)e
n 2n

- Choisir un ensemble S de 2 l nombres premiers inférieurs à 2 k log(k) ;

- Retirer de S les nombres premiers qui divisent b ;
Pour tout p ∈ S Faire
- Appliquer l’algorithme d’Euclide pour calculer le polynôme unitaire Rp ∈ Z[X] à
coefficients dans {0, . . . , p − 1} et tel que la réduction modulo p de Rp est un pgcd des
réductions P et Q de P et Q modulo p ;
Fin pour
- Soit e = minp∈S (deg(Rp )) et R = {p ∈ S | deg(Rp ) = e} ;
Si Card(R) ≥ l Alors
- Enlever Card(R) − l éléments de R de sorte à n’en garder que l ;
Q des restes chinois pour calculer S, P̃ et Q̃ dans Z[X] de norme
- Appliquer le théorème
infinie inférieure à ( p∈R p)/2 tels que pour tout p ∈ R, S ≡ b Rp mod p, P̃ S ≡ b P
mod  p, Q̃ S ≡ b Q mod p ; 
Si ||P̃ ||1 ||S||1 ≤ B et ||Q̃||1 ||S||1 ≤ B Alors
- Retourner pp(S) la partie primitive de S, i.e., S divisé par son contenu ;
Sinon
- Choisir un autre ensemble S à l’étape 4 et recommencer ;
Finsi
Sinon
- Choisir un autre ensemble S à l’étape 4 et recommencer.
Finsi

39
La condition sur les degrés se montre facilement par récurrence à partir de la proposition 3.2.1
et assure l’unicité des cofacteurs U et V une fois pgcd(A, B) fixé. La complexité de l’algo-
rithme d’Euclide étendu appliqué à deux polynômes A et B avec deg(A) = n ≥ deg(B) reste
quadratique en n. Cet algorithme est très important en calcul formel. Il permet entre autres
d’effectuer des inversions modulaires (et donc des divisions modulaires). En effet, si A et B
sont deux polynômes premiers entre eux et U et V les cofacteurs tels que

U A + V B = 1, deg(U ) < deg(B), deg(V ) < deg(A),

alors V est un inverse de B modulo A : V = 1/B mod A.

Exemple 3.2.11. Soient A = X 2 − X − 1 et B = X 3 − 1. Les polynômes A et B sont

premiers entre eux et on a :

X2 − 1 X −1 X −1 1
− A+ B = 1 ⇐⇒ = 3 mod (X 2 − X − 1).
2 2 2 X −1

Si on note φ le « nombre d’or » défini par A(φ) = 0, on peut par exemple en déduire l’égalité

φ−1 φ−1 (φ − 1)2

= (φ − 1) = .
φ3 − 1 2 2

On peut ainsi utiliser l’algorithme d’Euclide étendu pour « simplifier » des expressions sym-
boliques ce qui est un enjeu crucial en calcul formel.

Remarque 3.2.1. Il existe un algorithme rapide pour le calcul du pgcd (d’entiers ou de

polynômes). Cet algorithme se base sur des multiplications et des divisions euclidiennes ra-
pides et utilise la notion de demi-pgcd. Cet algorithme ne sera pas présenté dans ce cours
(voir TD) mais permet de calculer le pgcd de deux polynômes de degré majoré par n en
O(M(n) log(n)) opérations arithmétiques (où M une fonction telle que le produit de deux
polynômes de degré au plus n puisse se calculer en O(M(n)) avec M(n) ≥ 2 M(n/2)) ce qui
est le meilleur résultat connu à ce jour.

3.3 Résultant de deux polynômes

Dans cette section A désigne un anneau commutatif, unitaire et intègre et K un corps.

3.3.1 Matrice de Sylvester et forme échelonnée

Définition 3.3.1. Soient A = ni=0 ai X i , B = m i
P P
i=0 bi X ∈ A[X] de degrés respectifs n et
m. La matrice de Sylvester de A et B notée Syl(A, B) est la matrice carrée de taille n + m

40
définie par

an an−1 . . . ... a0

 an an−1 ... ... a0 
.. .. ..
 

 . . . 

an an−1 . . . . . . a0
 
Syl(A, B) =  ,
 
bm bm−1 . . . . . . b0 
bm bm−1 ... ... b0
 
 
 .. .. .. 
 . . . 
bm bm−1 . . . . . . b0

où les m premières lignes contiennent les coefficients de A et les n suivantes contiennent les
coefficients de B.

Cette matrice est étroitement liée à l’algorithme d’Euclide étendu car sa transposée
Syl(A, B)T représente la matrice de l’application linéaire

Am−1 [X] × An−1 [X] → An+m−1 [X],

( U , V ) 7→ U A + V B,

où pour k ∈ N, Ak [X] désigne l’ensemble des polynômes de degré inférieur ou égal à k. Les
combinaisons linéaires des lignes de Syl(A, B) donnent alors les coefficients des polynômes
qui peuvent s’écrire sous la forme U A+V B. Si A est un corps K, le pgcd de A et B est alors
caractérisé comme l’élément de plus petit degré qui peut être obtenu ainsi et peut être calculé
à partir d’une forme échelonnée en ligne de Syl(A, B). Cette dernière forme échelonnée en
ligne est donnée par une matrice P Syl(A, B) avec P inversible telle que :

1. Toutes les lignes nulles de P Syl(A, B) sont en-dessous de ses lignes non-nulles ;

2. Le premier coefficient non nul de chaque ligne non nulle de P Syl(A, B) est strictement
à droite du premier coefficient non nul de sa ligne précédente.

Lemme 3.3.1. Soient A et B deux polynômes à coefficients dans un corps K. La dernière

ligne non nulle d’une forme échelonnée en ligne de Syl(A, B) contient les coefficients d’un
pgcd de A et B.

Notons qu’une forme échelonnée en ligne d’une matrice à coefficients dans un corps K peut
s’obtenir en appliquant la réduction de Gauss.

Exemple 3.3.1. Considérons les polynômes A = X 4 − X 3 − 7 X 2 + 2 X + 3 et B = X 3 −

41
4 X 2 + 2 X + 3. On a
 
1 −1 −7 2 3 0 0
 

 0 1 −1 −7
3 0 
 2
0 0 1 −1 −7 2 3 
 

 
Syl(A, B) =  1 −4 2 3 0 0 0 ,
 
 

 0 1 −4 2 3 0 0 
 

 0 0 1 −4 2 3 0 
0 0 0 1 −4 2 3
qui admet comme forme échelonnée en ligne
 
1 −1 −7 2 3 0 0
 
 0 1
 −1 −7 2 
 3 0
 0 0 1 −1 −7 2 3 
 
 
 0 0 0 −14 45 −3 −18  .
 
 
5 12 9
 0 0
 0 0 −7 7 7


 
1 3 
 0 0
 0 0 0 10
− 10 

0 0 0 0 0 0 0
On obtient donc pgcd(A, B) = X − 3. Remarquons aussi que la dimension du noyau de la
matrice donne le degré du pgcd.
Définition 3.3.2. Soient A, B ∈ A[X]. Le résultant de A et B noté Res(A, B) ou encore
ResX (A, B) est le déterminant de la matrice de Sylvester Syl(A, B) de A et B.
Le résultant de A et B appartenant à A[X] est donc un élément de l’anneau A.
Proposition 3.3.1. Soient A et B deux polynômes à coefficients dans un corps K. Alors A
et B sont premiers entre eux si et seulement si Res(A, B) 6= 0.
Démonstration. Le déterminant de la matrice de Sylvester s’obtient comme produit des
éléments diagonaux de sa forme échelonnée en ligne. Par conséquent, Res(A, B) 6= 0 implique
qu’aucun élément diagonal n’est nul donc le pgcd obtenu est sur la dernière ligne donc de
degré 0, i.e., A et B sont premiers entre eux. Inversement si le pgcd de A et B est une
constante non nulle d ∈ K alors le dernier élément diagonal de la forme échelonnée en ligne
de Syl(A, B) est non nul. De plus en multipliant par X k , k = 1, . . . , n + m − 1 l’égalité
U A + V B = d on voit que tout élément diagonal de la forme échelonnée en ligne est non
nul et donc Res(A, B) 6= 0.
La proposition 3.3.1 nous fournit une méthode pour tester si deux polynômes sont pre-
miers entre eux : il suffit de calculer Res(A, B) et de regarder s’il est nul ou non. Notons que
l’on peut montrer que la proposition 3.3.1 se généralise au cas de polynômes à coefficients
dans un anneau factoriel.

42
3.3.2 Propriétés du résultant
Le premier résultat que nous donnons est une formule exprimant le résultant de deux poly-
nômes A et B en fonction des racines de A et B.
Théorème 3.3.1 (Formule de Poisson). Soient A, B ∈ A[X]. Si A et B s’écrivent
A = a (X − α1 ) · · · (X − αn ), B = b (X − β1 ) · · · (X − βm ),
alors on a Y
Res(A, B) = am bn (αi − βj ), (3.3)
i,j
ou encore
Y Y
Res(A, B) = (−1)m n bn A(βj ) = am B(αj ) = (−1)m n Res(B, A).
1≤j≤m 1≤j≤n

Démonstration. Le facteur am bn provient de la multilinéarité du déterminant donc on peut

supposer que les polynômes sont unitaires. Pour montrer l’égalité (3.3) nous allons supposer
que les αi et les βi sont des indéterminées. Q La proposition 3.3.1 appliquée sur le corps K =
Q(α1 , . . . , αn , β1 , . . . , βn ) montre que i,j (αi − βj ) divise Res(A, B). De plus le degré en αi
des m premières lignes de Syl(A, B) est 1 et le degré en αi des n lignes suivantes est 0 de
sorte que le degré du résultant en αi est majoré par m. De la Q même façon le degré en βj de
Res(A, B) est borné par n. Par conséquent Res(A, B) = c i,j (αi − βj ) pour une certaine
constante c. Pour déterminer c, il suffit de considérer le cas où αi = 0 pour tout i = 1, . . . , n.
Dans ce cas la matrice de Sylvester est triangulaire et son déterminant vaut B(0)n ce qui
implique c = 1. Nous avons donc prouvé (3.3). Les égalités suivantes sont des conséquences
directes de (3.3).
Corollaire 3.3.1. Soient A, B, C ∈ A[X]. Alors on a Res(A B, C) = Res(A, C) Res(B, C).
Démonstration. Le résultat est une conséquence de la formule de Poisson que l’on peut
appliquer dans une clôture algébrique du corps des fractions de l’anneau intègre A et du fait
que le résultant est un élément de A.
Proposition 3.3.2. Soient A, B ∈ A[X]. Il existe (U, V ) ∈ A[X]2 tel que
U A + V B = Res(A, B), deg(U ) < deg(B), deg(V ) < deg(A).
Démonstration. Notons Cj , j = 1, . . . , n+m, les colonnes de la matrice de Sylvester Syl(A, B),
n = deg(A)
P et m =−jdeg(B). En ajoutant à la dernière colonne Cn+m de Syl(A, B), l’élément
X n+m n+m
i=1 Cj X , la dernière colonne de la matrice devient

T
X m−1 A(X) X m−2 A(X) . . . A(X) X n−1 B(X) X n−2 B(X) . . . B(X) .
En développant le déterminant par rapport à cette colonne, on obtient donc
Res(A, B) = A(X) X m−1 M1,n+m + X m−2 M2,n+m + · · · + Mm,n+m

+ B(X) X n−1 Mm+1,n+m + X m−2 Mm+2,n+m + · · · + Mm+n,n+m ,

43
où Mi,n+m ∈ A désigne le déterminant de la co-matrice de Syl(A, B) obtenue en supprimant
la ième ligne et la dernière colonne. D’où le résultat.
Corollaire 3.3.2. Soient A, B ∈ A[X]. Alors Res(A, B) = 0 si et seulement si il existe
(U, V ) 6= (0, 0) ∈ A[X]2 tel que

U A + V B = 0, deg(U ) < deg(B), deg(V ) < deg(A).

Démonstration. Supposons Res(A, B) = 0, alors d’après la proposition 3.3.2, il existe un

couple (U, V ) ∈ A[X]2 tel que U A + V B = 0 avec deg(UP ) < deg(B) = m etPdeg(V ) <
deg(A) = n. Vérifions que U 6= 0 ou V 6= 0. En notant U = m−1 i
i=0 ui X et V =
n−1 i
i=0 vi X ,
l’égalité U A + V B = 0 s’écrit Syl(A, B)T w = 0 où w = (um−1 . . . u0 vn−1 . . . v0 )T ∈ An+m .
Or det(Syl(A, B)T ) = 0 donc le système linéaire homogène sSyl(A, B)T w = 0 admet une
solution non nulle. Inversement s’il existe (U, V ) 6= (0, 0) ∈ A[X]2 tel que U A + V B = 0
avec deg(U ) < deg(B) = m et deg(V ) < deg(A) = n, alors le système linéaire homogène de
matrice Syl(A, B)T admet une solution non nulle ce qui implique que son déterminant est
nul.

3.3.3 Calcul du résultant

Étant donnés deux polynômes A, B ∈ K[X] à coefficients dans un corps effectif K, posons
nous le problème du calcul du résultant de A et B. En utilisant la définition, calculer le
résultant revient à calculer le déterminant de la matrice de Sylvester. Le calcul du détermi-
nant d’une matrice pouvant se faire en la même complexité que celui du produit matriciel
(Théorème 2.3.1), les résultats de la sous-section 2.3 montrent que le résultant de A et B
peut se calculer en Õ(n2.81 ) opérations dans K en notant n = max(deg(A), deg(B)). Notons
qu’il existe un algorithme plus efficace (de complexité quadratique) qui calcule le résultant
en O(deg(A) deg(B)) opérations dans K. Nous ne détaillerons pas cet algorithme dans ce
cours mais il est basé sur la division euclidienne et le résultat suivant qui découle de la for-
mule de Poisson (Théorème 3.3.1) : soit A = Q B + R la division euclidienne de A de degré
n par B de degré m. Alors en notant bm le coefficient dominant de B et r = deg(R), on a
Res(A, B) = (−1)n m bn−r
m Res(B, R).

3.3.4 Résultant en plusieurs variables et élimination

L’une des principales applications du résultant en calcul formel est le fait que cela peut
permettre d’éliminer des variables dans l’étude d’un système polynomial. En effet, par défi-
nition le résultant ResX (A, B) de deux polynômes de A[X] est un élément de A. Donc si on
prend le cas particulier où l’anneau A lui même est un anneau de polynômes, par exemple
A = K[Y ] pour un certain corps K, alors le résultant ResX (A, B) de deux polynômes de
K[X, Y ] = K[Y ][X] = A[X] est un élément de A = K[Y ] et nous avons donc éliminé la
variable X. Plus généralement, étant donné un système polynomial défini par des polynômes
de plusieurs variables, on peut alors calculer des résultants pour éliminer successivement des
variables et ramener l’étude de ce système à l’étude de polynômes à une seule variable sur

44
lesquels on sait faire beaucoup plus de choses. Ce point ne sera pas détaillé ici mais nous
illustrerons simplement cet aspect avec le résultat suivant dans le cas de deux polynômes de
deux variables.

Théorème 3.3.2. Soient A, B ∈ K[X, Y ] = K[X][Y ] deux P polynômes bivariés Pà coefficients

dans un corps K algébriquement clos. Si on note A = ni=0 ai (X) Y i et B = m i
i=0 bi (X) Y ,
alors les racines du résultant ResY (A, B) ∈ K[X] sont les abscisses des solutions du système
polynomial {A(X, Y ) = 0, B(X, Y ) = 0} et les racines des termes de têtes an (X) et bm (X).

Démonstration. Exercice à partir du théorème 3.3.1 et de la proposition 3.3.2.

Terminons ce chapitre en illustrant ce résultat sur un exemple très simple.

Exemple 3.3.2. Soient A = X 2 Y + X + 1 et B = X Y − 1. On trouve alors le résultant

ResY (A, B) = −X (2 X + 1). On retrouve donc la racine X = 0 des termes de tête et la
racine X = −1/2 qui correspond à l’abscisse de la solution (−1/2, −2) du système polynomial
{X 2 Y + X + 1 = 0, X Y − 1 = 0}.

45
46
Chapitre 4

Factorisation de polynômes à une

variable

Dans ce chapitre, on note K un corps effectif. On considère un polynôme P ∈ K[X] unitaire

d’une seule variable X. On sait que P s’écrit de manière unique (à permutation près) sous
la forme
P = F1n1 × · · · × Fknk , (4.1)
avec pour i = 1, . . . , k, Fi irréductible, ni ∈ N∗ et les Fi deux à deux distincts. Le problème
de la factorisation s’énonce alors de la manière suivante : étant donné P , calculer les Fi et
les ni pour i = 1, . . . , k.

Deux cas particuliers importants sont celui où le corps K est un corps fini de caractéris-
tique p > 0 et celui où K = Q est le corps des nombres rationnels de caractéristique zéro
(ou, de façon analogue, le cas d’un polynôme à coefficients dans Z). Le cas de la factorisation
sur un corps fini a déjà été traité dans le cours Corps finis donc nous nous concentrerons ici
sur le cas de la factorisation sur un corps de caractéristique nulle. Les algorithmes procèdent
en général en plusieurs étapes la première étape consistant à calculer une factorisation dite
sans carré du polynôme.

4.1 Factorisation sans carré

4.1.1 Définition et existence
Définition 4.1.1. Un polynôme P ∈ K[X] unitaire est dit sans carré (ou sans facteur
multiple) s’il n’existe pas de polynôme F ∈ K[X] \ K tel que F 2 divise P .
Définition 4.1.2. Soit P ∈ K[X] unitaire. La factorisation sans carré de P est donnée par
P = P1 P22 · · · Pss , (4.2)
où pour tout i = 1, . . . , s, Pi est sans carré et les Pi sont premiers entre eux deux à deux. La
partie sans carré de P est alors définie par le produit P1 · · · Ps .

47
 Notons que dans la factorisation sans carré (4.2), les Pi ne sont pas nécessairement irré-
ductibles et certains Pi peuvent être égaux à 1.

Exemple 4.1.1. Soit P = X 8 − 11 X 6 + 36 X 4 − 16 X 2 − 64 ∈ Q[X]. La factorisation sans

carré de P s’écrit alors

P = P1 P22 P33 , P1 = X 2 + 1, P2 = 1, P3 = (X + 2) (X − 2).

Proposition 4.1.1. Tout polynôme P ∈ K[X] unitaire admet une et une seule factorisation
sans carré.

Démonstration. On sait que P s’écrit de manière unique sous la forme P = F1n1 · · · Fknk , avec
pour i = 1, . . . , k, Fi irréductible, ni ∈ N∗ et les Fi deux à deux distincts. On peut alors
écrire P = P1 P22 · · · Pss avec s = max1≤i≤k ni et pour j = 1, . . . , s,

1 si j ∈
/ {n1 , . . . , nk },
Pj = Q
F
ni =j i si j ∈ {n1 , . . . , nk }.

L’unicité découle de l’unicité des Fi (à permutation près).

4.1.2 Calcul de la factorisation sans carré

Intéressons nous maintenant au calcul de factorisation sans carré dans le cas où le corps K
est de caractéristique zéro, e.g., K = Q. Un avantage de la caractéristique nulle est que pour
P ∈ K[X], on a équivalence entre le fait que le polynôme dérivé1 P 0 soit égal à zéro est le
fait que le polynôme P soit constant, i.e., P ∈ K : P 0 = 0 ⇔ P ∈ K.

Proposition 4.1.2. Soient K un corps de caractéristique zéro et P ∈ K[X] un polynôme

unitaire. Alors P est sans carré si et seulement si pgcd(P, P 0 ) = 1, i.e., le polynôme P et
son polynôme dérivée P 0 sont premiers entre eux.

Démonstration. Supposons tout d’abord que P est sans carré et que G = pgcd(P, P 0 ) 6= 1.
Comme P est supposé sans carré on a P = P1 · · · Ps avec les Pi irréductibles, deux à deux
distincts et Pi ∈
/ K car P unitaire. En dérivant, on obtient alors

P 0 = P10 P2 · · · Ps + · · · + P1 · · · Ps−1 Ps0 .

Comme G 6= 1, il existe j ∈ {1, . . . , s} tel que Pj divise P 0 et donc Pj divise ( i6=j Pi ) Pj0 . Or
Q
pgcd(Pj , Pj0 ) = 1 et pgcd(Pj , Pi ) = 1 pour i 6= j car les Pj sont irréductibles et deux à deux
distincts d’où la contradiction. Inversement supposons qu’il existe F ∈ K[X] \ K tel que F 2
divise P . On alors P = F 2 H pour un certain polynôme H ∈ K[X]. En dérivant on obtient
alors P 0 = 2 F 0 F H + F 2 H 0 = F (2 F 0 H + F H 0 ). On a donc F divise P et F divise P 0 d’où
F divise pgcd(P, P 0 ) et donc pgcd(P, P 0 ) 6= 1 puisque F ∈ / K ce qui conclut la preuve.
Pn Pn
1
Si P = i=0 pi X i , avec pi ∈ K, alors P 0 = i=1 i pi X i−1 .

48
Proposition 4.1.3. Soit K un corps
Qs de caractéristique nulle. Soient P1 , . . . , Ps ∈ K[X] \ K
s
des polynômes sans carré, P = i=1 Pi et (c1 , . . . , cs ) ∈ K . Alors si les Pi sont premiers
entre eux deux à deux, on a
s
!
X P 0 Y
pgcd P, ci Pi = Pi .
i=1
P i i t.q. c =0 i

Démonstration. On a
s
! s s
! s  
X P Y X P Y P 0
pgcd P, ci Pi0 = pgcd Pj , ci Pi0 = pgcd Pj , cj P .
i=1
Pi j=1 i=1
Pi j=1
Pj j

D’où le résultat puisque le dernier pgcd vaut 1 si cj 6= 0 et Pj sinon.

À partir de ces résultats, nous pouvons procéder comme suit pour calculer la factorisation
sans carré (4.2) d’un polynôme P ∈ K[X] unitaire. En dérivant la relation (4.2), on obtient

P 0 = P10 P22 · · · Pss + P1 (2 P2 P20 ) P33 · · · Pss + · · · + P1 P22 · · · Ps−1

s−1
(s Ps0 Pss−1 ).

On a alors G = pgcd(P, P 0 ) = P2 P32 · · · Pss−1 et le polynôme H = P/G = P1 · · · Ps est un

polynôme sans carré. On en déduit que le polynôme P1 cherché est P1 = H/pgcd(H, G)
puisque pgcd(H, G) = P2 · · · Ps . On applique ensuite la même méthode à G pour déterminer
P2 et ainsi de suite. On obtient alors l’algorithme 13 parfois appelé algorithme de Musser
pour le calcul de la factorisation sans carré d’un polynôme sur un corps de caractéristique
nulle.

Algorithme 13 Factorisation sans carré d’un polynôme sur un corps de caractéristique 0

Entrées : P ∈ K[X] unitaire où K est un corps de caractéristique zéro.
Sorties : la liste (P1 , . . . , Ps ) telle que (4.2) soit la factorisation sans carré de P .
- Poser i = 1 ;
- Calculer G1 = pgcd(P, P 0 ) et H1 = GP1 ;
Tant que deg(Gi ) ≥ 1 Faire
- Calculer Gi+1 = pgcd(Gi , G0i ) et Hi+1 = GGi+1 i
;
Hi
- Calculer Pi = Hi+1 ;
- i = i + 1;
Fin tant que
- Poser Pi = Gi−1 ;
- Retourner (P1 , . . . , Pi ).

Proposition 4.1.4. Soit M une fonction telle que le produit de deux polynômes de degré
au plus n puisse se calculer en O(M(n)) opérations arithmétiques. Alors l’algorithme 13
est correct et permet de calculer la factorisation sans carré d’un polynôme de degré n à
coefficients dans un corps de caractéristique zéro en O(n M(n) log(n)) opérations dans K.

49
Démonstration. La correction de l’algorithme est laissée à titre d’exercice. Elle se déduit
des explications données ci-dessus. Le nombre d’étape de l’algorithme est majoré par n. Sa
complexité est donc majorée par n calculs de pgcd de polynômes de degré majoré par n
d’où le résultat d’après la remarque 3.2.1.

Illustrons maintenant le déroulement de l’algorithme 13 sur un exemple.

Exemple 4.1.2. Reprenons le polynôme P considéré à l’exemple 4.1.1. L’algorithme 13 se

déroule alors de la manière suivante :

1. i = 1, G1 = pgcd(P, P 0 ) = X 4 − 8 X 2 + 16, H1 = P/G1 = X 4 − 3 X 2 − 4 ;

2. G2 = pgcd(G1 , G01 ) = X 2 − 4, H2 = G1 /G2 = X 2 − 4, P1 = H1 /H2 = X 2 + 1, i = 2 ;

3. G3 = pgcd(G2 , G02 ) = 1, H3 = G2 /G3 = X 2 − 4, P2 = H2 /H3 = 1, i = 3 ;

4. P3 = G2 = X 2 − 4.

On retrouve bien la factorisation sans carré

P = P1 P22 P33 = (X 2 + 1) (X 2 − 4)3 = (X 2 + 1) ((X − 2) (X + 2))3 .

Nous donnons maintenant un autre algorithme 14 appelé algorithme de Yun qui permet de
supprimer le facteur n dans la complexité, i.e., d’obtenir une complexité en O(M(n) log(n))
opérations arithmétiques dans K.

Algorithme 14 Algorithme de Yun pour la factorisation sans carré en caractéristique 0

Entrées : P ∈ K[X] unitaire où K est un corps de caractéristique zéro.
Sorties : la liste (P1 , . . . , Ps ) telle que (4.2) est la factorisation sans carré de P .
- Poser i = 1 ;
0
- Calculer G = pgcd(P, P 0 ), H = G P
et K = PG ;
Tant que deg(H) ≥ 1 Faire
0)
- Calculer R = pgcd(H, K − H 0 ), K = (K−H R
et H = H
R
;
- Poser Pi = R ;
- i = i + 1;
Fin tant que
- Retourner (P1 , . . . , Pi ).

Proposition 4.1.5. Soit M une fonction telle que le produit de deux polynômes de degré au
plus n puisse se calculer en O(M(n)) opérations arithmétiques avec M(n + n0 ) ≥ M(n) +
M(n0 ). Alors l’algorithme 14 est correct et permet de calculer la factorisation sans carré d’un
polynôme de degré n à coefficients dans un corps de caractéristique zéro en O(M(n) log(n))
opérations dans K.

50
Démonstration. Montrons tout d’abord que l’algorithme est correct. Notons Hi et Ki les
polynômes H et K au départ du ième passage dans la boucle Tant que de l’algorithme et
Ri le polynôme R calculé à l’étape i. Soit m la valeur finale de i lorsqu’on sort de la boucle
Tant que et soient Hm et Km les polynômes H et K obtenus à la sortie de la boucle. On
montre par récurrence que pour tout i = 1, . . . , m, on a alors :
s s
Y X Hi 0
Hi = Pj , Ki = (j − i + 1) P,
j=i j=i
Pj j

où les Pi sont les polynômes cherchés vérifiantQs (4.2). Pour i = 1, nous avons H1 = P/G et
0 0 i−1
K1 = P /G. Comme G = pgcd(P, P ) = i=1 Pi , la propriété est vérifiée pour i = 1.
Supposons la propriété vraie pour i ∈ {1, . . . , m − 1}. On a alors, d’après la proposition 4.1.3
s s
!
Y X Hi
Y
Ri = pgcd(Hi , Ki − Hi0 ) = pgcd Pj , (j − i) Pj0 = P j = Pi ,
j=i j=i
Pj j=i

d’où les formules pour Ki+1 et Hi+1 . L’entier m est le premier entier pour lequel deg(Hm ) = 0
donc l’entier s de la factorisation sans carré (4.2) vaut m−1 et on retourne donc bien tous les
facteurs sans carrés P1 , . . . , Ps . Estimons maintenant la complexité de l’algorithme. L’étape
2 nécessite O(M(n) log(n)) opérations dans K. Ensuite, l’étape i de la boucle Tant que
nécessite O(M(deg(Hi )) log(deg(Hi ))) opérations dans K car deg(Ki ) ≤ deg(Hi ) − 1. Le
Pm−1
P de la boucle Tant que est alors majoré par O( i=1 M(deg(H
coût total i )) log(deg(Hi ))) ≤
O(M( m−1 i=1 deg(H ))
Qim−1 log(n)) en utilisant l’inégalité M(n) + M(n 0
) ≤ M(n + n0 ) d’où le
résultat puisque i=1 Hi = P .

Exemple 4.1.3. Reprenons le polynôme P considéré à l’exemple 4.1.1. L’algorithme 14 se

déroule alors de la manière suivante :

1. i = 1, G = pgcd(P, P 0 ) = X 4 − 8 X 2 + 16, H = P/G = X 4 − 3 X 2 − 4, K = P 0 /G =

2 X (4 X 2 − 1) ;

2. R = pgcd(H, K − H 0 ) = X 2 + 1, K = (K − H 0 )/R = 4 X, H = H/R = X 2 − 4,

P1 = R = X 2 + 1, i = 2 ;

3. R = pgcd(H, K −H 0 ) = 1, K = (K −H 0 )/R = 2 X, H = H/R = X 2 −4, P2 = R = 1,

i = 3;

4. R = pgcd(H, K − H 0 ) = X 2 − 4, K = (K − H 0 )/R = 0, H = H/R = 1, P3 = R =

X 2 − 4, i = 4.

On retrouve bien la factorisation sans carré

P = P1 P22 P33 = (X 2 + 1) (X 2 − 4)3 = (X 2 + 1) ((X − 2) (X + 2))3 .

51
 Lorsque le corps K est un corps fini de caractéristique p > 0, alors une difficulté supplé-
mentaire est que si la caractéristique p est plus petite que le degré d’un polynôme P ∈ K[X],
alors P 0 peut être nul sans que P soit nécessairement constant. Dans ce cas, la bonne notion
n’est alors pas la factorisation sans carré mais la factorisation dite séparable que nous n’abor-
derons pas dans ce cours. Notons que la factorisation séparable d’un polynôme de degré n
peut aussi se calculer en O(M(n) log(n)) opérations arithmétiques.

4.2 Factorisation sur Z

4.2.1 Rappels : factorisation sur un corps fini
La factorisation d’un polynôme univarié à coefficients dans un corps fini Fq a été étudiée
dans le cours Corps finis. Les deux principaux algorithmes de factorisation dans ce cas sont
celui de Berlekamp et celui de Cantor et Zassenhaus. Nous résumons les résultats obtenus
dans le théorème suivant :

Théorème 4.2.1. Soit P ∈ Fq [X] unitaire à coefficients dans le corps fini à q éléments Fq .
Alors, on a les résultats suivants :

1. L’algorithme de Berlekamp déterministe (resp. probabiliste) calcule la factorisation

(4.1) de P en O(n3 + n2 q) (resp. O(n3 + n2 log(n) log(q))) opérations dans Fq ;

2. L’algorithme (probabiliste) de Cantor et Zassenhaus calcule la factorisation (4.1) de

P en O(n3 log(q)) opérations dans Fq en utilisant de l’arithmétique classique et en
O(n2 log2 (n) log(log(n)) log(q)) opérations dans Fq en utilisant de l’arithmétique rapide.

On supposera donc dans la section suivante que l’on sait factoriser un polynôme à coef-
ficient dans un corps fini.

4.2.2 Algorithme utilisant un grand nombre premier

Soit P ∈ Z[X] un polynôme primitif. On suppose qu’une factorisation sans carré de P a déjà
été calculée (voir la section précédente) de sorte que le problème de factoriser P est ramené
au problème de factoriser un polynôme sans carré. On supposera donc ici que le polynôme à
factoriser, que l’on note encore P , est primitif et sans carré.

L’existence de la borne de Mignotte (Lemme 3.2.1) sur la taille des coefficients des fac-
teurs d’un polynôme donné nous permet d’utiliser une stratégie modulaire pour calculer la
factorisation (4.1) avec ni = 1, i = 1, . . . , k, d’un polynôme primitif et sans carré à coeffi-
cients dans Z. En notant B la borne de Mignotte, on peut alors choisir un nombre premier
p tel que p > 2 B factoriser le polynôme P dans Z/(p Z) et essayer, par une recherche ex-
haustive, de retrouver un facteur de P sur Z. Cette stratégie mène à l’algorithme 15.

52
Algorithme 15 Factorisation d’un polynôme sur Z en utilisant un grand premier
Entrées : P ∈ Z[X] un polynôme sans-carré et primitif de degré n.
Sorties : P√(si P est irréductible sur Z) ou un facteur non-trivial de P sur Z.
- Soit B = n + 1 2n ||P ||∞ la borne de Mignotte ;
- Choisir un premier p tel que p > 2 B ;
- Factoriser P dans Z/(p Z)[X] : P = F1 · · · Fs avec les Fi irréductibles dans Z/(p Z)[X]
et à coefficients dans {−b p2 c, . . . , b p2 c} ;
Si s = 1 Alors
- P est irréductible dans Z/(p Z)[X] et donc Retourner P qui est irréductible sur Z ;
Sinon
PourQtout I ⊂ {1, . . . , s} tel que I = 6 ∅ Faire
Si i∈I Fi divise P dans Z[X] Alors
Q
- Retourner i∈I Fi qui est un facteur non trivial de P sur Z ;
Finsi
Fin pour
- Retourner P qui est irréductible sur Z.
Finsi

Lorsque l’algorithme 15 renvoie un facteur F de P alors nous pouvons appliquer à nouveau

l’algorithme à F et P/F et ainsi de suite jusqu’à ce que tous les facteurs obtenus soient
irréductibles auquel cas nous avons calculé la factorisation cherchée de P . Notons que dans
le pire cas, nous aurons 2n ensembles I à tester ce qui mène à une complexité exponentielle.
Cependant, si le nombre s de facteurs dans Z/(p Z)[X] est petit, alors cet algorithme a un coût
raisonnable et est utilisable en pratique. Finalement, on remarque qu’un autre inconvénient
de cette méthode est que la borne de Mignotte B peut être grande (si n et ||P ||∞ sont grands)
et dans ce cas, le nombre premier p sera lui aussi grand ce qui impliquera que l’étape de
factorisation dans Z/(p Z)[X] sera très couteuse. L’algorithme proposé dans la sous-section
suivante évite cet inconvénient en utilisant seulement un petit nombre premier p.

4.2.3 Algorithme utilisant la remontée de Hensel

L’algorithme que nous allons présenter maintenant consiste à choisir un petit nombre premier
p convenable, à factoriser P modulo p puis à utiliser la remontée de Hensel pour en déduire
des factorisations de P modulo pk . En choisissant alors k de sorte que pk > 2 B nous pourrons
ainsi calculer des facteurs de P sur Z.

Lemme 4.2.1 (Lemme de Hensel). Soient P ∈ Z[X] et p un nombre premier. On suppose

que P ≡ F1 G1 mod p avec F1 , G1 ∈ Z/(p Z)[X] et F1 et G1 premiers entre eux dans
Z/(p Z)[X]. Alors, pour tout k ∈ N∗ , il existe (Fk , Gk ) ∈ (Z/(pk Z)[X])2 tel que
(
P ≡ Fk Gk mod pk ,
Fk ≡ F1 mod p, Gk ≡ G1 mod p.

53
Démonstration. Raisonnons par récurrence sur k. Par hypothèse, l’énoncé est vrai pour k =
1. Supposons l’énoncé vrai pour k et cherchons (Fk+1 , Gk+1 ) ∈ (Z/(pk+1 Z)[X])2 tel que

P ≡ Fk+1 Gk+1 mod pk+1 , Fk+1 ≡ F1 mod p, Gk+1 ≡ G1 mod p.

On pose Fk+1 = Fk + pk F̂k et Gk+1 = Gk + pk Ĝk , où F̂k et Ĝk sont deux polynômes à
déterminer. On obtient alors

Fk+1 Gk+1 = Fk Gk + pk (F̂k Gk + Fk Ĝk ) + p2 k F̂k Ĝk .

Donc P ≡ Fk+1 Gk+1 mod pk+1 équivaut à P ≡ Fk Gk + pk (F̂k Gk + Fk Ĝk ) mod pk+1 . Or,
par hypothèse de récurrence, on a P = Fk Gk + pk P̂k pour un certain polynôme P̂k , d’où

pk P̂k ≡ pk (F̂k Gk + Fk Ĝk ) mod pk+1 ,

et donc
P̂k ≡ F̂k Gk + Fk Ĝk mod p.
Par hypothèse de récurrence, nous avons aussi Fk ≡ F1 mod p et Gk ≡ G1 mod p ce qui
implique alors
P̂k ≡ F̂k G1 + F1 Ĝk mod p. (4.3)
Comme F1 et G1 sont supposés premiers entre eux dans Z/(p Z)[X], on peut toujours calculer
(F̂k , Ĝk ) vérifiant (4.3) (Cf Chapitre 3) ce qui conclut la preuve.
Notons que la preuve ci-dessus donne une méthode effective pour effectuer la remontée
de Hensel en utilisant l’algorithme d’Euclide étendu (Cf Chapitre 3). Illustrons ceci sur un
exemple :
Exemple 4.2.1. Considérons le polynôme P = X 3 + 8 X − 7, le nombre premier p = 7
et utilisons les notations de la preuve précédente. On a P ≡ F1 G1 mod 7 avec F1 = X
et G1 = X 2 + 1 premiers entre eux dans Z/(7 Z) et P = F1 G1 + 7 (X − 1) de sorte que
P̂1 = X −1. Calculons alors F2 et G2 tels que P ≡ F2 G2 mod 72 , F2 ≡ F1 mod 7, G2 ≡ G1
mod 7. D’après la preuve précédente, on peut prendre F2 = F1 +7 F̂1 et G2 = G1 +7 Ĝ1 où F̂1
et Ĝ1 sont tels que Pˆ1 ≡ F̂1 G1 + F1 Ĝ1 mod 7. En utilisant l’algorithme d’Euclide étendu,
on calcule U et V dans Z/(7 Z) tels que 1 ≡ U F1 + V G1 mod 7. On trouve ici U = 6 X
et V = 1. On pose ensuite F̂1 = Pˆ1 V = X − 1 et F̂1 = Pˆ1 U = 6 X (X − 1) et sorte que
F2 = F1 + 7 F̂1 = 8 X − 7 et G2 = G1 + 7 Ĝ1 = 43 X 2 − 42 X + 1. On peut alors vérifier que
l’on a P ≡ F2 G2 mod 72 , F2 ≡ F1 mod 7, G2 ≡ G1 mod 7.
Dans cette stratégie nous allons devoir choisir un petit nombre premier p. Pour s’assurer
que notre polynôme sans carré à coefficients dans Z reste sans carré dans Z/(p Z), il nous
suffit de supposer que p ne divise pas le discrimant de P , i.e., que P ne divise pas le résultant
Res(P, P 0 ). On obtient l’algorithme 16.

Lorsque l’algorithme 16 renvoie un facteur Gk de P alors nous pouvons appliquer à

nouveau l’algorithme à Gk et P/Gk et ainsi de suite jusqu’à ce que tous les facteurs obtenus

54
Algorithme 16 Factorisation d’un polynôme sur Z en utilisant la remontée de Hensel
Entrées : P ∈ Z[X] un polynôme sans-carré et primitif de degré n.
Sorties : P (si P est irréductible sur Z) ou un facteur non-trivial de P sur Z.
- Soit p un « petit » nombre premier ne divisant pas Res(P, P 0 ) ;
- Factoriser P dans Z/(p Z)[X] : P ≡ F1 · · · Fs mod p ;
Si s = 1 Alors
- P est irréductible dans Z/(p Z)[X] et donc Retourner P qui est irréductible sur Z ;
Sinon
Pour tout I ⊂Q {1, . . . , s} tel queQI =
6 ∅ Faire
- Poser G1 = i∈I Fi et H1 = i∈I / Fi ;
- Utiliser le lemme de Hensel pour calculer Gk et Hk à coefficients dans Z tels que
P ≡ Gk Hk mod pk avec pk > 2 B ;
Si Gk divise P dans Z[X] Alors
- Retourner Gk qui est un facteur non trivial de P sur Z ;
Finsi
Fin pour
- Retourner P qui est irréductible sur Z.
Finsi

soient irréductibles auquel cas nous avons calculé la factorisation cherchée de P . Comme
celle de l’algorithme précédent, la complexité de cet algorithme est exponentielle à cause
de la combinatoire, i.e., du nombre de sous-ensembles de {1, . . . , s} et du fait que dans le
pire des cas on doit tester de manière exhaustive tous les sous-ensembles. Notons qu’il existe
des algorithmes de complexité polynomiale pour factoriser des polynômes sur Z où cette
recherche exhaustive est remplacée par la résolution d’un problème de type sac à dos en
utilisant l’algorithme LLL pour la réduction des réseaux.

55
56
Chapitre 5

Intégration symbolique

Dans ce chapitre nous allons aborder le problème de l’intégration symbolique « sous forme
close », i.e., étant donnée une expression rationnelle formée de fonctions obtenues par com-
position d’exponentielles et de logarithmes, décider si elle admet une primitive de la même
nature et si oui la calculer. Par exemple, les algorithmes que nous allons voir permettent de
montrer que la fonction x 7→ exp(x2 ) n’admet pas de primitive s’écrivant comme expression
rationnelle formée de fonctions obtenues par composition d’exponentielles et de logarithmes.
Par contre la fonction
x (x2 exp(2 x2 ) − ln2 (x + 1))2 + 2 x exp(3 x2 ) (x − (2 x2 + 1) (x + 1) ln(x + 1))

(x + 1) (ln2 (x + 1) − x2 exp(2 x2 ))2

admet pour primitive la fonction

x exp(x2 ) ln(x + 1) 1
1
+ ln ln(x + 1) + x exp(x2 ) − ln ln(x + 1) − x exp(x2 ) .


x−ln(x+1)+ 2 2 2
x exp(2 x ) − ln (x + 1) 2 2

Le chapitre est organisé comme suit. Dans un premier temps nous poserons brièvement
un cadre algébrique pour l’étude de ce problème. Ensuite nous aborderons en détails le cas
particulier où la fonction dont on cherche une primitive est une fraction rationnelle. Enfin,
nous donnerons l’idée de l’algorithme traitant le cas général en se ramenant au cas des
fractions rationnelles.

5.1 Algèbre différentielle et fonctions élémentaires

Définition 5.1.1. Soit A un anneau. Une dérivation de A est une application δ : A → A
vérifiant :

∀f, g ∈ A, δ(f + g) = δ(f ) + δ(g), δ(f g) = δ(f ) g + f δ(g).

Un anneau A muni d’une dérivation δ est appelé anneau différentiel et noté (A, δ). Si A est
un corps, alors on parle de corps différentiel.

57
 L’exemple classique est le corps C(x) des fractions rationnelles à coefficients dans C
muni de la dérivation usuelle δ = d/dx. La fonction exponentielle x 7→ exp(x) est souvent
définie comme une fonction f vérifiant f 0 = f . Cette dernière relation permet alors d’obtenir
récursivement toutes les dérivées de f . Pour modéliser et généraliser cette situation, nous
introduisons la notion d’extension différentielle :

Définition 5.1.2. Soit (K, δ) un corps différentiel. On dit que (K1 , δ1 ) est une extension
différentielle de (K, δ) si K ⊂ K1 et si la restriction de δ1 à K coïncide avec δ.

Par exemple si on prend (K, δ) = (C(x), d/dx), alors l’extension différentielle (K1 , δ1 )
avec K1 = K(X) et δ1 définie par δ1 (X) = X modélise le corps différentiel C(x, exp(x)) muni
de la dérivation usuelle. Dans la suite, par abus de notation, nous noterons souvent de la
même façon, e.g., nous utiliserons le symbole 0 , la dérivation δ sur K et δ1 sur K1 .
Pour modéliser une expression rationnelle formée de fonctions obtenues par composition
d’exponentielles et de logarithmes, nous utilisons le formalisme suivant :

Définition 5.1.3. Soit (K,0 ) un corps différentiel. On dit que (K1 ,0 ) est une extension
élémentaire de (K,0 ) si K1 = K(t1 , . . . , tn ) où, pour i = 1, . . . , n, ti vérifie l’une des trois
conditions suivantes :

1. ti algébrique sur K(t1 , . . . , ti−1 ), i.e., il existe P ∈ K(t1 , . . . , ti−1 )[X] tel que P (ti ) = 0 ;

2. il existe ai ∈ K(t1 , . . . , ti−1 ) tel que t0i = a0i /ai . Le symbole ti modélise alors ln(ai ), le
logarithme de ai ;

3. il existe ai ∈ K(t1 , . . . , ti−1 ) tel que t0i = a0i ti . Le symbole ti modélise alors exp(ai ),
l’exponentielle de ai .

Dans la suite nous considérerons uniquement le cas où (K,0 ) = (C(x), d/dx) et nous
appellerons fonction élémentaire toute fonction appartenant à une extension élémentaire de
(C(x), d/dx).

Exemple 5.1.1. Avec le formalisme de la définition 5.1.3, pour construire la fonction élé-
mentaire exp(x2 ) nous introduisons t1 tel que t01 = 2 x t1 et exp(x2 ) est alors représenté
par t1 . Pour construire ln(ln(x)), nous introduisons t1 tel que t01 = 1/x puis t2 tel que
t02 = t01 /t1 = 1/(x t1 ) et le représentant de ln(ln(x)) est alors t2 . De la même manière, pour
introduire la fonction élémentaire ln(1 + exp(cos(x))), nous procédons comme suit. On intro-
duit tout d’abord t1 tel que t01 = i t1 (t1 modélise alors exp(i x)), puis on introduit t2 tel que
0
t02 = ((t21 + 1)/(2 t1 )) t2 (t2 modélise alors exp(cos(x))) et finalement t3 tel que t03 = t02 /(1+t2 )
qui représente ln(1 + exp(cos(x))).

Calculer la dérivée d’une fonction élémentaire peut se faire simplement en utilisant la

définition successives des éléments ti dans la définition 5.1.3. En outre, la dérivée d’une
fonction élémentaire est une fonction élémentaire. En ce qui concerne l’opération inverse
consistant à calculer une primitive d’une fonction élémentaire, le problème s’avère plus com-
pliqué d’autant plus qu’une fonction élémentaire n’admet pas nécessairement de primitive

58
élémentaire. Le problème de l’intégration symbolique en calcul formel se pose alors de la ma-
nière suivante : étant donnée une fonction élémentaire f , décider si elle admet une primitive
F qui est elle-même une fonction élémentaire et si oui déterminer F . Avant de considérer
le cas général, nous allons tout d’abord nous concentrer sur le cas particulier où la fonction
élémentaire f dont on cherche à calculer une primitive est une fraction rationnelle de C(x).

5.2 Le cas particulier des fractions rationnelles

Étant donnée une fraction rationnelle f /g ∈ C(x), l’objectif de cette section est de calculer
F telle que RF 0 = f /g. Une telle primitive
R
une primitive F = f /g de f /g, i.e., calculer
R
se décompose alors sous la forme F = f /g = qR+ c/d + a/b où q ∈ C[x] est la Rpartie
polynomiale, c/d ∈ C(x) est la partie rationnelle et a/b est la partie logarithmique de f /g.
Le calcul de la primitive se décompose donc en trois étapes :
1. Calcul de la partie polynomiale ;
2. Calcul de la partie rationnelle ;
3. Calcul de la partie logarithmique.
Les trois sous-sections suivantes donnent des algorithmes pour effectuer chacune de ces
étapes.

5.2.1 Calcul de la partie polynomiale

Le calcul de la partie polynomiale ne pose pas de problème particulier et peut être effectué
simplement par l’algorithme 17 :
Algorithme 17 Partie polynomiale de la primitive d’une fraction rationnelle
Entrées : f, g ∈ C[x] deux polynômes de degré au plus n. R R
Sorties : q, h ∈ C[x] tels que deg(q) ≤ n, deg(h) < deg(g) et f /g = q + h/g.
- Appliquer l’algorithme d’Euclide pour calculer le quotient q̃ = di=0 qi xi et le reste r de
P
la division euclidienne f par g ;
Pd qi de i+1
- Calculer q = i=0 i+1 x ;
- Retourner q et h = r.

Proposition 5.2.1. Soit M une fonction telle que le produit de deux polynômes de degré au
plus n puisse se calculer
R en O(M(n)) opérations arithmétiques. L’algorithme 17 calcule la
partie polynomiale de f /g où f et g sont deux polynômes de degré au plus n en O(M(n))
opérations arithmétiques.
Démonstration. Le coût est dominé par celui de la division euclidienne d’où le résultat d’après
la proposition 2.2.6.
R
Nous sommes donc ramené à calculer h/g où h et g sont deux polynômes tels que
deg(h) < deg(g).

59
5.2.2 Calcul de la partie rationnelle
R
Pour calculer la partie rationnelle de h/g nous allons utiliser la décomposition en fractions
partielles qui est une alternative à la décomposition en éléments simples utilisant une fac-
torisation sans carré du dénominateur g au lieu de sa factorisation complète en produits
de facteurs irréductibles. On rappelle que d’après les résultats de la sous-section 4.1.2, la
factorisation sans carré d’un polynôme g de degré n s’écrit g = g1 g22 · · · gss où chaque gi est
un polynôme sans carré (pas nécessairement irréductible) et les gi sont premiers entre eux
deux à deux et peut se calculer en O(M(n) log(n)) opérations arithmétiques. Étant donnée
la factorisation sans carré g = g1 g22 · · · gss , le calcul de la décomposition en fractions partielles
de h/g se déroule en deux étapes. On commence par écrire h/g sous la forme

h p1 p2 ps
= + 2 + ··· + s, ∀i ∈ {1, . . . , s}, deg(pi ) < deg(gii ), (5.1)
g g1 g2 gs

puis on décompose à nouveau chaque pi /gii sous la forme

pi pi,1 pi,2 pi,i
i
= + 2 + ··· + i , ∀i ∈ {1, . . . , s}, ∀j ∈ {1, . . . , i}, deg(pi,j ) < deg(gi ). (5.2)
gi gi gi gi

En multipliant l’équation (5.1) par g = g1 g22 · · · gss , on obtient

h = p1 g22 · · · gss + p2 g1 g33 · · · gss + · · · + ps g1 · · · gs−1

s−1
,

de sorte que
i−1 i+1
∀i ∈ {1, . . . , s}, h ≡ pi g1 · · · gi−1 gi+1 · · · gss mod gii .

Les polynômes gjj et gii étant premiers entre eux pour i 6= j, gjj est inversible modulo gii et son
inverse modulaire gj−j modulo gii peut être calculé à l’aide de l’algorithme d’Euclide étendu
(voir la sous-section 3.2.3). On obtient donc les congruences
−(i−1) −(i+1)
∀i ∈ {1, . . . , s}, pi ≡ h g1−1 · · · gi−1 gi+1 · · · gs−s mod gii .

L’hypothèse deg(pi ) < deg(gii ) nous assure alors que les polynômes pi sont déterminés de
manière unique. Une fois les pi déterminés, pour obtenir l’écriture (5.2), on pose ri = pi
puis on effectue successivement les divisions euclidiennes de rj par gi pour j = i, . . . , 2. On
écrit alors rj = rj−1 gi + pi,j , on définit pi,1 = r1 et on a pi = pi,1 gii−1 + pi,2 gii−2 · · · + pi,i qui
implique (5.2). Finalement on obtient l’algorithme 17 :

Proposition 5.2.2. Soit M une fonction telle que le produit de deux polynômes de degré
au plus n puisse se calculer en O(M(n)) opérations arithmétiques. Soient h et g deux poly-
nômes tels que deg(h) < deg(g) = n et soit g = g1 g22 · · · gss la factorisation sans carré de g.
Alors l’algorithme 18 est correct et calcule la décomposition en fractions partielles de h/g en
O(s M(n) log(n)) opérations arithmétiques.

60
Algorithme 18 Décomposition en fractions partielles
Entrées : h, g ∈ C[x] deux polynômes tels que deg(h) < deg(g) = n et la factorisation
sans carré g = g1 g22 · · · gss de g.
Sorties : des polynômes pi,j ∈ C[x] tels que l’on a (5.1) et (5.2).
Pour tout i = 1, . . . , s Faire
  l’algorithme d’Euclide étendu pour calculer deux polynômes ui et vi tels
- Appliquer
que ui ggi + vi gii = 1 ;
i
- Calculer le reste pi dans la division euclidienne de h ui par gii ;
- Poser ri = pi , puis effectuer successivement les divisions euclidiennes rj = rj−1 gi + pi,j
de rj par gi pour j = i, . . . , 2 et poser pi,1 = r1 ;
Fin pour
- Retourner les pi,j .

Démonstration. La correction de l’algorithme 18 découle directement des explications don-

nées ci-dessus. Étudions sa complexité arithmétique. Le calcul de chaque pi requiert une
division exacte ggi de deux polynômes de degré au plus n, une inversion modulaire en degré
i
au plus n (qui peut se faire à l’aide de l’algorithme d’Euclide étendu - voir la sous-section 3.2.3
- en complexité O(M(n) log(n)) d’après la remarque 3.2.1), une multiplication h ui de poly-
nômes de degré au plus n et une division euclidienne d’un polynôme de degré au plus 2 n par
un polynôme de degré n. La complexité du calcul de tous les pi est donc en O(s M(n) log(n))
opérations arithmétiques. Pour un i fixé, les polynômes rj sont de degrés au plus j deg(gi )
donc le coût du calcul de pi,i , . . . , pi,1 est majoré par O(i2 M(deg(gi ))) opérations arithmé-
tiques. En faisant la somme pour i allant de 1 à s et en utilisant M(n+n0 ) ≥ M(n)+M(n0 ) on
obtient que l’extraction de tous les pi,j peut se faire en O(s M(n)) opérations arithmétiques
d’où le résultat.
Si g = g1 g22 · · · gss désigne la factorisation sans carré de g, nous avons donc réduit le
problème du calcul de primitive de h/g avec deg(h) < deg(g) à celui du calcul de primitive
de i=1 ij=1 pi,j /gij avec pour tout i, j, deg(pi,j ) < deg(gi ). Les polynômes gi sont sans
Ps P
carrés de sorte que pgcd(gi , gi0 ) = 1 (voir Proposition 4.1.2) ce qui entraine qu’ils existent
deux polynômes ui et vi tels que ui gi + vi gi0 = 1. On obtient donc
pi,j ui gi + pi,j vi gi0 −(j − 1) gi0
Z   
−pi,j vi
Z Z Z
pi,j pi,j ui
= = + ,
gij gij gij−1 j−1 gij
puis, en intégrant par parties
−pi,j vi pi,j ui + (pi,j vi )0 /(j − 1)
Z Z
pi,j
= + .
gij (j − 1) gij−1 gij−1
On voit donc que l’on peut ainsi ramener l’intégration d’un pôle d’ordre j à celle d’un
pôle d’ordre j − 1. En appliquant ce procédé de manière itérative, on obtient le procédé de
réduction de Hermite qui permet de calculer la partie rationnelle de la primitive de h/g en
suivant l’algorithme 19.

61
Algorithme 19 Partie rationnelle de la primitive d’une fraction rationnelle : réduction de
Hermite
Entrées : h, g ∈ C[x] deux polynômes Ps Pi tels que deg(h) < deg(g) = n et la décomposition
j
en fractions partielles h/g = i=1 j=1 pi,j /gi où g = g1 g22 · · · gss est la factorisation sans
carré de g. R P RP
Sorties : des polynômes ai , bi , ci et di tels que h/g = i ci /di + i ai /bi et où les bi
2 s−1
divisent g1 g2 · · · gs et les di divisent g2 g3 · · · gs .
- Poser l = 0 ;
Pour tout i = 1, . . . , s Faire
- Appliquer l’algorithme d’Euclide étendu pour calculer ui et vi de degré < deg(gi ) tels
que ui gi + vi gi0 = 1 ;
Pour tout j = i, . . . , 2 Faire
- Effectuer la division euclidienne de vi pi,j par gi : on écrit vi pi,j = q gi + r avec
deg(r) < deg(gi ) ;
- Poser t = ui pi,j + (q gi0 )/(j − 1) ;
- Poser pi,j−1 = pi,j−1 + t + r0 /(j − 1) (normaliser l’expression obtenue pour pi,j−1 ) ;
- Poser l = l + 1 ;
- Poser cl = −r/(j − 1) et dl = gij−1 ;
Fin pour
Fin pour
- Pour i = 1, . . . , s, poser ai = pi,1 et bi = gi ;
- Retourner les polynômes ai , bi , ci et di .

62
Proposition 5.2.3. Soit M une fonction telle que le produit de deux polynômes de degré
au plus n puisse se calculer en O(M(n)) opérations P arithmétiques. Soient h et g deux po-
s Pi
lynômes tels que deg(h) < deg(g) = n et soit h/g = i=1 j=1 pi,j /gij la décomposition en
fraction partielles de h/g où g = g1 g22 · · · gss est la factorisation
R sans carré de g. Alors l’algo-
rithme 19 est correct et calcule la partie rationnelle de h/g en O(M(n) log(n)) opérations
arithmétiques.
Démonstration. La correction de l’algorithme 19 découle directement des explications don-
nées ci-dessus. L’étude de la complexité de l’algorithme est laissée en exercice. Elle utilise
les mêmes arguments que ceux de la preuve de la proposition 5.2.2.
Exemple 5.2.1. Appliquons la réduction de Hermite pour le calcul de la primitive de la
faction rationnelle
1
.
(2 x + 1)2
2

Posons h = 1/4 et g = (x2 + 1/2)2 = g22 avec g2 = x2 + 1/2 sans carré et appliquons
la réduction de Hermite pour calculer la partie rationnelle de h/g. Avec les notations de
l’algorithme 19, nous avons donc h/g = p2,2 /g22 avec p2,2 = h = 1/4. L’algorithme d’Euclide
étendu nous fournit donc u2 = 2 et v2 = −x tels que u2 g2 + v2 g20 = 1. On a alors v2 p2,2 =
q g2 + r avec q = 0 et r = v2 p2,2 = −x/4. On obtient alors p2,1 = u2 p2,2 − 1/4 = 1/4 puis
c1 = x/4 et d1 = g2 . Finalement on pose a2 = p2,1 = 1/4 et b2 = g2 de sorte que
Z Z Z
1 x/4 1/4 x 1
= + = + .
(2 x2 + 1)2 x2 + 1/2 x2 + 1/2 2 (2x2 + 1) 2 (2 x2 + 1)
La forme des polynômes di et bi retournés par l’algorithme 19 montre R que la sortie de
2 s−1
l’algorithme peut s’écrire sous la forme normalisée c/(g2 g3 · · · gs ) + a/(g1 g2 · · · gs ). De
plus, en notant δi , i = 1, . . . , s, le degré de gi et en étudiant le degré des polynômes ci et ai
dans l’algorithme, on voit que le degré de c sera strictement inférieur à δ2 +2 δ3 +· · ·+(s−1) δs
et celui de a sera strictement inférieur à δ1 + δ2 + · · · + δs . On peut alors utiliser une méthode
par coefficients indéterminés pour calculer directement a et b en résolvant un système linéaire.
Cette méthode est dû à Horowitz.
Proposition 5.2.4. Soit MM une fonction telle que le produit de deux matrices carrées
de taille n puisse se calculer en O(MM(n)) opérations arithmétiques. Soient h et g deux
polynômes tels que deg(h) < deg(g) = n et g = g1 g22 · · · gss la factorisation
R sans carré de
g. Alors l’algorithme 20 est correct et calcule la partie rationnelle de h/g en O(MM(n))
opérations arithmétiques.
Démonstration. La correction de l’algorithme découle des explications ci-dessus. Le coût est
clairement dominé par la résolution du système linéaire qui a n équations et n inconnues
d’où le résultat d’après le théorème 2.3.1.
Exemple 5.2.2. Reprenons l’exemple 5.2.1 et appliquons l’algorithme 20 de Horowitz. Nous
avons toujours h = 1/4 et g2 = x2 + 1/2. Avec les notations de l’algorithme 20, on obtient

63
Algorithme 20 Partie rationnelle de la primitive d’une fraction rationnelle : algorithme de
Horowitz
Entrées : h, g ∈ C[x] deux polynômes tels que deg(h) < deg(g) = n et g = g1 g22 · · · gss la
factorisation sans carré de g. R R
Sorties : des polynômes a et c tels que h/g = c/(g2 g32 · · · gss−1 ) + a/(g1 g2 · · · gs ).
- Poser δ = δ1 + δ2 + · · · + δs où pour i = 1, . . . , s, δi = deg(gi ) ;
- Écrire a = a0 + a1 x + · · · aδ−1 xδ−1 et c = c0 + c1 x + · · · + cn−δ−1 xn−δ−1 , où les ai et les
ci sont des coefficients à déterminer ;
- Écrire le système P linéaire obtenu en égalant à zéro tous les coefficients du polynôme
h − c (g1 · · · gs ) + c si=2 (i − 1) g1 g2 · · · gi−1 gi0 gi+1 · · · gs − a g2 g32 · · · gss−1 ;
0

- Résoudre le système linéaire obtenu ;

- Substituer la solution du système linéaire dans a et c et retourner a et c.

δ = 2 de sorte que l’on cherche a et c sous la forme a = a0 + a1 x et c = c0 + c1 x où

a0 , a1 , c0 , c1 sont des constantes à déterminer. On a alors
 
0 0 3 2
 a1  1 c 1 a0
h − c g2 + c g2 − a g2 = −a1 x + (c1 − a0 ) x + 2 c0 − x+ − − ,
2 4 2 2

d’où le système linéaire {−a1 = 0, c1 − a0 = 0, 2 c0 − a1 /2 = 0, 1/4 − c1 /2 − a0 /2 = 0} que

l’on résout pour trouver a0 = 1/4, a1 = 0, c0 = 0, c1 = 1/4. On a donc trouvé a = 1/4 et
c = (1/4) x et on retrouve donc
Z Z Z
1 c a x 1
2 2
= + = 2
+ .
(2 x + 1) g2 g2 2 (2x + 1) 2 (2 x2 + 1)

Notons que l’algorithme 20 de Horowitz est significativement moins bon du point de vue
de la complexité que l’algorithme 19 utilisant la réduction de Hermite (voir le chapitre 2).
Son avantage est qu’il ne nécessite pas le calcul préalable de la décomposition en fractions
partielles et qu’il est très simple à programmer.

5.2.3 Calcul de la partie logarithmique

Nous sommes maintenant ramené au problème du calcul de la partie logarithmique
R de la pri-
mitive d’une fraction rationnelle. Autrement dit, nous voulons calculer a/b où a/b ∈ K(x),
deg(a) < deg(b) et b est un polynôme sans Rcarré. Dans ce cas là, nous allons naturellement
voir apparaître des logarithmes (penser à 1/x) mais nous allons devoir faire face à une
autre difficulté à savoir le fait que les coefficients de la primitive vivent en général dans une
extension algébrique de K. Nous allons utiliser le théorème suivant qui sera admis pour ce
cours.

Théorème 5.2.1. Soit a/b ∈ K(x) une fraction rationnelle telle que a et b sont premiers
entre eux, deg(a) < deg(b) = n et b est un polynôme unitaire et sans carré. Soit L une

64
extension algébrique de K. Alors on a
Z m
a X
= ci ln(vi ),
b i=1

pour des constantes ci ∈ L et des polynômes vi ∈ L[x] premiers entre eux deux à deux si
et seulement si le résultant R(y) = Resx (b, a − b0 y) (parfois appelé résultant de Rothstein-
Trager) se factorise en produit de facteurs linéaires en y dont les ci sont les zéros et les vi
sont les pgcd(b, a − b0 ci ).

Ce théorème mène à l’algorithme 21 de Rothstein-Trager pour le calcul de la partie

logarithmique de la primitive d’une fraction rationnelle.

Algorithme 21 Partie logarithmique de la primitive d’une fraction rationnelle : algorithme

de Rothstein-Trager
Entrées : a, b ∈ K[x] deux polynômes tels que deg(a) < deg(b) = n avec b unitaire et
sans carré.
Sorties : une extension algébrique L de K, des constantes c1 , . . . , cm ∈ L et des polynômes
v1 , . . . , vm ∈ L[x] tels que a/b = m
R P
i=1 ci ln(vi ).
- Calculer le résultant R(y) = Resx (b, a − b0 y) ;
- Factoriser R(y) en produit de facteurs linéaires en introduisant l’extension algébrique
nécessaire L de K ;
- Soient c1 , . . . , cm les zéros deux à deux distincts de R(y) dans L ;
- Pour i = 1, . . . , m, calculer vi = pgcd(b, a − b0 ci ) ;
- Retourner L, les ci et les vi .

Proposition 5.2.5. Soit M une fonction telle que le produit de deux polynômes de degré
au plus n puisse se calculer en O(M(n)) opérations arithmétiques. Soient a, b ∈ K[x] deux
polynômes tels que deg(a) < deg(b) = n avec b unitaire
R et sans carré. Alors l’algorithme
21 est correct et calcule la partie logarithmique de a/b en O(n M(n)2 log(n)) opérations
arithmétiques dans K plus le coût de la factorisation du résultant (qui est de degré au plus
n) à l’étape 2.

Démonstration. La correction de l’algorithme est une conséquence directe du théorème 5.2.1.

Le calcul du résultant de deux polynômes de degré au plus n peut se calculer en O(n2 ) (voir
la sous-section 3.3.3). Le calcul de R(y) nécessite donc O(n2 ) opérations dans K[y] et comme
le degré en y des polynômes qui apparaissent ne dépasse pas n, le résultant R(y) peut se
calculer en O(n2 M(n)) opérations dans K. De même le calcul des pgcd à l’étape 4 coûte
O(M(n) log(n)) opérations dans L soit au plus O(M(n)2 log(n)) opérations dans K puisque
le degré de L sur K ne dépasse pas le degré de R(y) qui est majoré par n. D’où l’estimation
annoncée.

65
Exemple 5.2.3. Reprenons l’exemple 5.2.1 et appliquons l’algorithme de Rothstein-Trager
pour calculer Z
1
.
2 (2 x2 + 1)
On pose donc a = 1/4 et b = x2 + 1/2. On a
1
1 0 2
1
R(y) = Resx (b, a − b0 y) = −2 y 1
4
0 = 2 y2 + .
16
1
0 −2 y 4
√
On introduit maintenant l’extension algébrique Q(i 2) de Q pour factoriser R(y) en facteurs
linéaires sous la forme
√ ! √ !
i 2 i 2
R(y) = 2 y + y− .
8 8
√ √ √
On pose alors c1 = i 8 2 , c2 =√ − i 8 2 et on calcule v1 = pgcd(b, a − b0 c1 ) = x + i
2
2
et
v2 = pgcd(b, a − b0 c2 ) = x − i 2 2 . On obtient alors
Z √ √ ! √ √ !
1 i 2 i 2 i 2 i 2
= c1 ln(v1 ) + c2 ln(v2 ) = ln x + − ln x −
2 (2 x2 + 1) 8 2 8 2

d’où √ √ ! √ !!
Z
1 x i 2 i
2 i 2
= + ln x + − ln x − .
(2 x2 + 1)2 2 (2x2 + 1) 8 2 2

Notons finalement qu’en utilisant un algorithme « rapide » pour calculer le résultant

R(y), nous pouvons gagner un facteur n et obtenir une complexité en O(M(n)2 log(n)) pour
le calcul de la partie logarithmique par l’algorithme de Rothstein-Trager.

5.3 Le cas général des fonctions élémentaires

Étant donnée une fonction élémentaire, nous allons maintenant donner l’idée d’un algorithme
qui décide si elle admet une primitive qui est elle même une fonction élémentaire et si oui la
calcule. L’idée de cet algorithme est basée sur une idée attribuée à Liouville qui consiste à
généraliser le Théorème 5.2.1 aux fonctions élémentaires qui ne sont pas nécessairement des
fractions rationnelles. On obtient alors le résultat suivant qui est admis dans le cadre de ce
cours.

Théorème 5.3.1 (Principe de Risch-Liouville). Soient (K,0 ) une extension élémentaire de

(C(x),0 ), C son corps des constantes, i.e., C = {f ∈ K | f 0 = 0}, et f ∈ K une fonction
élémentaire. S’il existe une extension élémentaire (L,0 ) de (K,0 ) et g ∈ L telle que g 0 = f ,

66
alors il existe v ∈ K, des constantes ci , i = 1, . . . , n, algébriques sur C et des éléments
vi ∈ K(c1 , . . . , cn ), i = 1, . . . , n, tels que
n
0
X vi0
f =v + ci .
i=1
vi

Ce théorème de structure montre que s’il existe une primitive, alors elle sera d’une forme
particulière et on sait à l’avance le type de fonctions qu’il faudra introduire pour l’exprimer.
Ce théorème mène à un algorithme appelé algorithme de Risch qui consiste à reconnaitre
cette forme particulière. Soit f ∈ K = C(t1 , . . . , tn ). Le principe général de l’algorithme
de Risch consiste à procéder de manière récursive pour « éliminer » les ti et se ramèner à
l’intégration d’une fraction rationnelle. Cet algorithme est assez technique et ne sera pas
détaillé dans ce cours.

67
68
Chapitre 6

Équations différentielles linéaires et

séries formelles

Dans ce chapitre, l’objet de notre étude est une équation différentielle linéaire à coefficients
polynomiaux de la forme

an (x) y (n) (x) + an−1 (x) y (n−1) (x) + · · · + a1 (x) y 0 (x) + a0 (x) y(x) = 0, (6.1)

où, pour i = 0, . . . , n, les ai (x) ∈ Q[x] (ou C[x]) sont des polynômes donnés, y est une
fonction inconnue de la variable x et y (i) désigne la ième dérivée de la fonction y par rapport
à x. En notant
 n  n−1  
d d d
L = an (x) + an−1 (x) + · · · + a1 (x) + a0 (x)
dx dx dx

l’opérateur différentiel linéaire associé, l’équation (6.1) s’écrit alors :

L(y) = 0.

L’objectif consiste à donner des méthodes de calcul de certains types de séries solutions
de L(y) = 0 au voisinage d’un point x0 ∈ C.

6.1 Résolution de L(y) = 0 par des séries entières

6.1.1 Structure des solutions
Les solutions de L(y) = 0 forment un espace vectoriel sur Q (ou R ou C). En effet, si c1 , c2
sont des constantes et y1 , y2 des solutions de L(y) = 0, alors on vérifie que L(c1 y1 +c2 y2 ) = 0.
Pour étudier la dimension de cet espace vectoriel, on introduit la notion de wronskien.

69
Définition 6.1.1. Soient y1 , . . . , yr des fonctions d’une variable x. On appelle wronskien de
y1 , . . . , yr et on note Wr(y1 , . . . , yr ) le déterminant

y1 y2 ... yr
y10 y20 ... yr0
Wr(y1 , . . . , yr ) = .. .. .. .
. . .
(r−1) (r−1) (r−1)
y1 y2 . . . yr
Lemme 6.1.1. Des fonctions y1 , . . . , yr d’une variable x sont linéairement indépendantes
sur Q (ou R ou C) si et seulement si Wr(y1 , . . . , yr ) 6= 0.
Démonstration. Des fonctions y1 , . . . , yr sont linéairement dépendantes sur Q si et seulement
si il existe des constantes c1 , . . . , cr ∈ Q non toutes nulles telles que c1 y1 + · · · + cr yr = 0. En
(i) (i)
dérivant cette relation, on obtient c1 y1 +· · ·+cr yr = 0 pour tout i ≥ 1. Le système linéaire
(i) (i)
de r équations {c1 y1 + · · · + cr yr = 0, i = 0, . . . , r − 1} à r inconnues c1 , . . . , cr admet une
solution non nulle si et seulement si son déterminant est égal à zéro d’où le résultat.
Corollaire 6.1.1. Soit L(y) = y (n) + n−1 (i)
P
i=0 ai (x) y (x). Les solutions de L(y) = 0 forment
un espace vectoriel de dimension au plus n.
Démonstration. Supposons que y1 , . . . , yn+1 sont des solutions de L(y) = 0. On a

y1 y2 . . . yn+1
y10 y20 0
. . . yn+1
Wr(y1 , . . . , yn+1 ) = .. .. .. ,
. . .
(n) (n) (n)
y1 y2 . . . yn+1
(n) Pn−1 (i)
et donc, en utilisant les relations yj =− i=0 ai yj , j = 1, . . . , n + 1, on obtient

y1 y2 ... yn+1
y10 y20 ... 0
yn+1
Wr(y1 , . . . , yn+1 ) = .. .. .. .
. . .
Pn−1 (i) Pn−1 (i) Pn−1 (i)
− i=0 ai y1 − i=0 ai y 2 ... − i=0 ai yn+1 (x)

On en déduit donc Wr(y1 , . . . , yn+1 ) = 0 (la dernière ligne du déterminant ci-dessus est une
combinaison linéaire des lignes précédentes) et le lemme 6.1.1 implique alors que y1 , . . . , yn+1
sont linéairement dépendantes d’où le résultat.
Définition 6.1.2. Soient y1 , . . . , yn des solutions de L(y) = 0. Si elles sont linéairement
indépendantes sur C, alors on dit qu’elles forment un système fondamental de solutions de
L(y) = 0.
Dans la suite « résoudre L(y) = 0 » signifiera trouver un système fondamental de so-
lutions. Le terme « solution » peut avoir plusieurs interprétations. Prenons l’exemple de
l’équation L(y) = y 00 + y = 0. On distingue alors :

70
 • les solutions sous forme close (explicite finie) c’est-à-dire dans une classe de fonctions
donnée (sin, cos, exp, . . . ). Ici y(x) = c1 cos(x)+c2 sin(x) est solution pour tout couple
de constantes (c1 , c2 ) ∈ C2 ;

• les solutions sous forme de séries P formelles (on ne se pré-occupe pas de la convergence
+∞
des
P+∞ séries). Ici, si on note y(x) = i
i=0 ci x une telle série formelle, alors on a y 00 (x) =
i 00
P +∞ i
i=0 (i + 1) (i + 2) ci+2 x d’où y(x) + y (x) = i=0 (ci + (i + 1) (i + 2) ci+2 ) x donc les
coefficients ci d’une telle série solution doivent vérifier ci + (i + 1) (i + 2) ci+2 = 0.

• les solutions sous forme de séries analytiques (en 0), i.e., les séries y(x) = ∞ i
P
i=0 yi x
solutions qui convergent dans un disque D(0; R) centré en 0 et de rayon R.

6.1.2 Existence de séries solutions au voisinage d’un point ordinaire

Définition 6.1.3. Soit L(y) = an (x) y (n) (x) + · · · + a1 (x) y 0 (x) + a0 (x) y(x) = 0 où pour
tout i = 0, . . . , n, ai (x) ∈ C[x]. On dit que x0 ∈ C est un point ordinaire de L(y) = 0 si
an (x0 ) 6= 0. Sinon on dit que x0 est un point singulier (ou une singularité) de L(y) = 0.
On dira que l’infini est un point ordinaire (resp. singulier) de L(y) = 0 si 0 est un point
ordinaire (resp. singulier) de l’équation obtenue après le changement de variable X = 1/x.
Théorème 6.1.1 (Théorème de Cauchy). Au voisinage d’un point ordinaire, une équation
différentielle linéaire L(y) = 0 admet un système fondamental de solutions analytiques.
Démonstration. Quitte à faire un changement de variable, on suppose que x0 = 0 est un
point ordinaire de L(y) = 0. Comme an (0) 6= 0, on sait que si y est solution de L(y) = 0,
alors
n−1
(n)
X ai (0) (i)
y (0) = − y (0).
i=0
a n (0)
De plus, on a
P 0
n−1 ai (x)
y (n+1) (x) = − i=0 an (x) y (i) (x)
P 
n−1 ai (x)0 an (x)−ai (x) an (x)0 (i)
Pn−2 ai (x) (i+1) an−1 (x) (n)
= − i=0 an (x)2
y (x) + i=0 an (x) y (x) + an (x)
y (x) ,
Pn−1 ai (x)
puis en utilisant à nouveau la relation y (n) (x) = − i=0 an (x) y (i) (x), on obtient

n−1
(n+1)
X ai,1 (x)
y (x) = y (i) (x), ai,1 (x) ∈ C[x].
i=0
an (x)2

De la même façon, pour tout entier k ∈ N, on obtient des relations de la forme

n−1
(n+k)
X ai,k (x) (i)
y (x) = y (x), ai,k (x) ∈ C[x].
a (x)2k
i=0 n

71
Par conséquent les y (n+k) (0) dépendent tous linéairement de y(0), y 0 (0), . . . , y (n−1) (0). On
construit alors un système fondamental y1 , . . . , yn de solutions sous forme de séries de la
(i)
yj (0) i
manière suivante. Pour j = 1, . . . , n, yj (x) = +∞
P
i=0 i!
x est la série définie par :
(k)
• Pour k = 0, . . . , n − 1, yj (0) = 1 si k = j − 1 et 0 sinon.
(n+k) Pn−1 ai,k (0) (i) aj−1,k (0)
• Pour k ∈ N, yj (0) = i=0 an (0)2k yj (0) = a2nk (0)
.

Les séries y1 , . . . , yn ainsi construites sont solutions

Pnde L(y) = 0 (par construction) et elles
sont linéairement indépendantes. En effet, si y = i=1 ci yi = 0, alors, pour k = 0, . . . , n − 1,
la relation y (k) (0) = 0 implique ck+1 = 0. Elles forment donc un système fondamental
de solutions de L(y) = 0. On admettra que ces séries sont analytiques ce qui termine la
preuve.
Remarque 6.1.1. Les séries y1 , . . . , yn construites dans la preuve précédente sont conver-
gentes (admis) sur un disque D(0, R), où R est le module de la singularité la plus proche.
Le procédé de construction d’un système fondamental de séries solutions d’une équation
différentielle linéaire au voisinage d’un point ordinaire décrit dans la preuve précédente est
souvent appelé méthode de Cauchy. Cette méthode a l’avantage de se généraliser à toutes les
équations différentielles (non nécessairement linéaires) de la forme

(n) P (x, y(x), y 0 (x), . . . , y (n−1) (x))

y (x) = ,
Q(x, y(x), y 0 (x), . . . , y (n−1) (x))


où P et Q sont des polynômes (multivariés). À toute condition initiale x0 , y(x0 ), . . . , y (n−1) (x0 )
telle que Q(x0 , y(x0 ), . . . , y (n−1) (x0 )) 6= 0, on peut associer une unique série formelle solution
de l’équation.

6.1.3 Calcul des séries solutions au voisinage d’un point ordinaire

Il existe des méthodes autres que celle de Cauchy abordée précédemment pour construire
un système fondamental de séries formelles solutions d’une équation différentielle linéaire
L(y) = 0 au voisinage d’un point ordinaire.

Méthode d’Euler-Newton : calcul par approximations successives. L’idée consiste

à chercher chaque solution yi (pour i = 1, . . . , n) sous la forme

xi−1
yi (x) = + xn Σn (x), Σn ∈ C[[x]],
(i − 1)!
où Σn est une série formelle à déterminer. L’équation L(yi ) = 0 donne alors une nouvelle
équation  i−1 
x
L̃(Σn (x)) = −L ,
(i − 1)!

72
pour la série formelle Σn . On pose alors Σn (x) = cn + x Σn+1 (x) et le terme de plus bas degré
de la relation  i−1 
x
L̃(cn + x Σn+1 (x)) = −L
(i − 1)!
donne une équation pour cn que l’on résout. On réitère ensuite ce procédé pour trouver les
coefficients suivants de Σn .

Méthode par relation de récurrence. En remarquant que

(xk )(i) = k (k − 1) · · · (k − i + 1) xk−i ,
Pmi
et en notant ai (x) = j=0 ai,j xj avec ai,j ∈ C les coefficients (polynomiaux) de L, on obtient
n
X mi
n X
X
k k (i)
L(x ) = ai (x) (x ) = k (k − 1) · · · (k − i + 1) ai,j xk−i+j ,
i=0 i=0 j=0

de sorte que
+∞
! +∞ X mi
n X
X X
k
L ck x = k (k − 1) · · · (k − i + 1) ck ai,j xk−i+j .
k=0 k=0 i=0 j=0
P+∞
Il en découle que les coefficients ck d’une série formelle y(x) = k=0 ck xk solution de L(y) = 0
satisfont une relation de récurrence de la forme
M
X
α` (k) ck−` = 0.
`=0

Pour calculer explicitement cette relation de récurrence, on peut procéder comme suit :
1. On décompose l’opérateur différentiel linéaire L sous la forme L = L−A + · · · + LB où
les opérateurs différentiels linéaires Li (i = −A, . . . , B) satisfont
∀k ∈ N, Li (xk ) = αi (k) xk+i ,
et les αi (k) sont des constantes (i.e., ne dépendent pas de la variable x).
2. Si y(x) = +∞ k
P
k=0 ck x est une série formelle solution de L(y) = 0, alors on a
P+∞ k


0 = L(y) = L k=0 ck x
P+∞ k

P+∞ k


= L−A k=0 c k x + · · · + L B k=0 c k x
P+∞ k−A +∞
+ · · · + k=0 αB (k) ck xk+B
P
= k=0 α−A (k) ck x
P+∞
= k=0 (α (k + A) ck+A + · · · + αB (k − B) ck−B ) xk
| −A {z }
(R)

et donc les coefficients ck satisfont la relation de récurrence

(R) : α−A (k + A) ck+A + · · · + αB (k − B) ck−B = 0.

73
 Nous allons maintenant illustrer les trois différentes méthodes abordées ci-dessus sur un
exemple.

Exemple 6.1.1. On considère l’équation différentielle linéaire

L(y) = (1 − x2 ) y 00 (x) − x y(x) = 0,

pour laquelle x0 = 0 est un point ordinaire et on cherche à calculer un système fondamental

de séries formelles solutions au voisinage de 0.

1. Appliquons tout d’abord la méthode de Cauchy. On commence par chercher une pre-
P+∞ y1(i) (0) i
mière solution y1 (x) = i=0 i!
x , avec y1 (0) = 1 et y10 (0) = 0. On a alors les
relations
x 1 + x2 x
y100 (x) = y1 (x), y1000 (x) = y1 (x) + y 0 (x), ...
1 − x2 1−x 2 1 − x2 1

qui impliquent y100 (0) = 0, y1000 (0) = 1, . . . de sorte que y1 (x) = 1 + 61 x3 + O(x4 ). On
peut ainsi calculer autant de termes que l’on veut de la série y1 (x). De même, on
(i)
y2 (0) i
cherche ensuite une seconde série formelle y2 (x) = +∞
P
i=0 i!
x , solution avec cette
fois y2 (0) = 0 et y20 (0) = 1. En utilisant les mêmes relations que précédemment on
obtient alors y200 (0) = 0, y2000 (0) = 0, . . . de sorte que y2 (x) = x + O(x4 ).

2. Appliquons maintenant la méthode d’Euler-Newton. On commence par chercher une

première série y1 (x) solution sous la forme y1 (x) = 1 + x2 Σ2 (x) avec Σ2 ∈ C[[x]] à
déterminer. En développant l’égalité L(1 + x2 Σ2 (x)) = 0, on obtient alors la nouvelle
équation différentielle linéaire

x2 (1 − x2 ) Σ002 (x) + 4 x (1 − x2 ) Σ02 (x) + 2 (1 − x2 ) − x3 Σ2 (x) = x,

pour Σ2 ∈ C[[x]]. En injectant Σ2 (x) = c2 + x Σ3 (x) dans cette dernière équation et en

regardant le terme de plus bas degré on trouve 2 c2 = 0 de sorte que c2 = 0. Puis, en
injectant Σ3 (x) = c3 + x Σ4 (x) (i.e., Σ2 (x) = c3 x + x2 Σ4 (x)) et en regardant le terme
de plus bas degré on trouve 6 c3 = 1 de sorte que c3 = 1/6. On retrouve finalement la
série formelle solution y1 (x) = 1 + 61 x3 + O(x4 ). De la même façon, on cherche ensuite
une seconde solution sous la forme y2 (x) = x + x2 Σ2 (x) et on retrouve la solution
y2 (x) = x + O(x4 ).

3. Appliquons finalement la méthode basée sur laPrelation de récurrence satisfaite par

les coefficients ck d’une série formelle y(x) = +∞ k
k=0 ck x solution de L(y) = 0. On
décompose tout d’abord L sous la forme L = L−2 + L0 + L1 où

L−2 (y) = y 00 , L−2 (xk ) = k (k − 1) xk−2 ,

L0 (y) = −x2 y 00 , L0 (xk ) = −k (k − 1) xk ,

74
 L1 (y) = −x y, L1 (xk ) = −xk+1 .
La relation de récurrence satisfaite par les coefficients ck est alors

(R) : (k + 2) (k + 1) ck+2 − k (k − 1)ck − ck−1 = 0.

Pour k = 0, on obtient c2 = 0, pour k = 1, on obtient 6 c3 − c0 = 0, pour k = 2, on

obtient 12 c4 − 2 c2 − c1 = 0, . . . . En prenant comme conditions initiales c0 = y1 (0) = 1
et c1 = y10 (0) = 0, on obtient c3 = 1/6, c4 = 0 et on retrouve la première solution
y1 (x) = 1 + 61 x3 + O(x4 ). Alternativement, en prenant comme conditions initiales
c0 = y2 (0) = 0 et c1 = y20 (0) = 1, on obtient c3 = 0, c4 = 1/12 et on retrouve la seconde
1
solution y2 (x) = x + 12 x4 + O(x5 ) = x + O(x4 ).

6.2 Résolution de L(y) = 0 par des séries quasi-entières

6.2.1 Séries quasi-entières
Soit L(y) = an (x) y (n) (x) + · · · + a1 (x) y 0 (x) + a0 (x) y(x) = 0 où pour tout i = 0, . . . , n,
ai (x) ∈ C[x]. Nous avons vu dans la section précédente que si 0 est un point ordinaire de
L(y) = 0 (i.e., an (0) 6= 0), alors nous pouvons construire un système fondamental de sé-
ries (formelles) entières solutions dans C[[x]]. En revanche, lorsque 0 est un point singulier
de L(y) = 0 (i.e., an (0) = 0), les séries entières forment une classe trop petite c’est-à-dire
qu’au voisinage d’un point singulier, il n’est en général pas possible de trouver un système
fondamental de solutions sous forme de séries entières. Par exemple, √ l’équation différentielle
0
linéaire 2 x y (x) − y(x) = 0 qui admet pour solution y(x) = c x, où c est une constante
arbitraire n’admet pas de série entière solution au voisinage de zéro. En effet, au voisinage
de zéro nous avons y(x) = x1/2 Σ(x) avec Σ ∈ C[[x]].

Dans la suite, pour λ ∈ C nous noterons xλ la solution de y 0 (x) = λx y(x) telle que y(1) = 1
et log(x) la solution de y 0 (x) = x1 telle que y(1) = 0. Nous avons alors xλ = exp(λ log(x)).
Si k ∈ N, alors nous identifierons xk xλ = xλ+k . Enfin, nous avons les dérivées par rapport à
x (resp. λ) suivantes :
dxλ dxλ
= λ xλ−1 , = log(x) xλ .
dx dλ
Remarque 6.2.1. On peut vérifier que si λ 6∈ N, alors y 0 (x) = λ
x
y(x) n’admet pas de
solution sous forme de série entière (au voisinage de 0).
Définition 6.2.1. On appelle série quasi-entière un objet de la forme xλ Σ(x) où λ ∈ C et
Σ ∈ C[[x]] est une série entière. Si le terme constant de Σ est non nul (i.e., Σ(0) 6= 0), alors
λ est appelé l’ordre de la série quasi-entière xλ Σ(x).
Par exemple les séries de Laurent +∞
P k −k0
P+∞ k
k=−k0 ck x = x k=0 c−k0 +k x Psont des séries
quasi-entières. On appelle série de Puiseux une série de Laurent en x , i.e., +∞
1/r
k=−k0 ck x
k/r
.
λ
Une série quasi-entière x Σ(x) avec λ ∈ Q est donc une série de Puiseux.

75
6.2.2 Équation indicielle et calcul des séries quasi-entières solutions
Soit L = L−A + · · · + LB la décomposition de L telle que les opérateurs différentiels linéaires
Li (i = −A, . . . , B) satisfont

∀k ∈ N, Li (xk ) = αi (k) xk+i .

Définition 6.2.2. Avec les notations précédentes, le polynôme E0 (X) = α−A (X) est appelé
polynôme indiciel de L en 0 et l’équation E0 (X) = 0 est appelée équation indicielle.

En injectant une série quasi-entière y(x) = xλ +∞ k

P
k=0 ck x dans l’équation L(y) = 0, on
montre que :
+∞
! (
X E0 (λ) = 0,
L xλ ck x k = 0 ⇔
k=0
∀k ∈ N, E0 (λ + k + A) ck+A + · · · + αB (λ + k − B) ck−B = 0.
(6.2)

Théorème 6.2.1. Soit L(y) = an (x) y (n) (x) + · · · + a1 (x) y 0 (x) + a0 (x) y(x) = 0 où pour tout
i = 0, . . . , n, ai (x) ∈ C[x].

1. L’équation L(y) = 0 admet une série quasi-entière solution d’ordre λ seulement si

E0 (λ) = 0.

2. Si λ est une racine de E0 (X) = 0 et si, pour tout k ∈ N∗ , E0 (λ + k) 6= 0, alors

l’équation L(y) = 0 admet une série quasi-entière solution d’ordre λ.

Démonstration. Ceci découle de l’équivalence (6.2). Pour 2, si E0 (λ + k) =6 0 pour tout

∗
k ∈ N , alors on peut utiliser la dernière égalité de (6.2), pour calculer successivement
chaque coefficient ck+A .

Les racines de l’équation indicielle E0 (X) = 0 sont appelées exposants de L en 0. Lorsque

0 est un point ordinaire de L(y) = 0, alors avec les notations précédentes, on a A = n et
L−A (xk ) = a k (k − 1) · · · (k − n + 1) xk−n pour une certaine constante non nulle a ∈ C.
L’équation indicielle est alors donnée par E0 (X) = a X (X − 1) · · · (X − n + 1) et nous avons
les n exposants entiers 0, 1, . . . , n − 1 de sorte que l’on retrouve un système fondamental de
séries entières solutions (voir section précédente). Lorsque 0 est un point singulier, alors le
degré du polynôme indiciel E0 (X) est au plus égal à n. Dans le cas où deg(E0 (X)) = n, on
dit que 0 est un point singulier régulier de L(y) = 0.

Nous terminons ce chapitre avec deux exemples qui illustrent le calcul d’un système
fondamental de séries solutions au voisinage d’un point singulier (ici au voisinage de 0).

Exemple 6.2.1. Considérons l’équation d’Euler d’ordre 2 donnée par

L(y) = x2 y 00 (x) + a x y 0 (x) + b y(x) = 0, a, b ∈ C.

76
La décomposition de L s’écrit alors simplement L = L0 car nous avons
L(xk ) = (k (k − 1) + a k + b) xk .
L’équation indicielle est donc E0 (X) = X (X − 1) + a X + b = 0 et nous devons distinguer
deux cas :
1. Cas 1 : E0 (X) a deux racines distinctes λ1 et λ2 . Dans ce cas nous avons deux séries
quasi-entières solutions xλ1 et xλ2 (qui sont linéairement indépendantes).
2. Cas 2 : E0 (X) a une racine double λ. Dans ce cas nous avons une série quasi-entière
solution xλ . Pour construire la deuxième solution nous considérons l’équation L(xk ) =
E0 (k) xk . Comme les dérivations d/dx et d/dk commutent, nous avons
 k
dx d
E0 (k) xk = E0 (k)0 xk + E0 (k) log(x) xk ,


L =
dk dk
ou encore
L log(x) xk = (E0 (k)0 + E0 (k) log(x)) xk .

Comme λ est une racine double de E0 (X), nous avons E0 (λ) = E00 (λ) = 0 ce qui prouve
que log(x) xλ est solution de L(y) = 0. Clairement xλ et log(x) xλ sont linéairement
indépendantes donc nous avons trouvé un système fondamental de solutions.
Exemple 6.2.2. Considérons l’équation
3 0
L(y) = x y 00 (x) + y (x) + y(x) = 0.
2
Nous avons L = L−1 + L0 avec
 
3 3
L−1 (y) = x y + y 0 ,
00
L−1 (x ) = k (k − 1) + k xk−1 ,
k
2 2
L0 (y) = y, L0 (xk ) = xk .
L’équation indicielle est alors
 
3 1
E0 (X) = X (X − 1) + X = X X+ .
2 2
Cette équation indicielle admet donc deux racines distinctes λ1 = 0 et λ2 = −1/2 qui ne
diffèrent pas d’un entier. Le théorème 6.2.1 implique donc que nous avons une série entière
Σ1 ∈ C[[x]] solution et l’autre solution est une série quasi-entière x−1/2 Σ2 (x) avec Σ2 ∈
C[[x]]. La relation de récurrence pour les coefficients ck des séries Σ1 et Σ2 est donnée par
E0 (λ + k + 1) ck+1 + ck = 0.
Pour chacune des séries Σ1 et Σ2 , en se donnant la condition initiale c0 , on peut alors
calculer c1 , c2 , . . . par la formule
−ck
ck+1 = ,
E0 (λ + k + 1)
puisqu’à la fois pour λ = λ1 = 0 et λ = λ2 = −1/2, nous avons E0 (λ + k + 1) 6= 0, pour tout
k ∈ N.

77
Remerciements
Ce cours a été essentiellement préparé à partir :

1. des notes de cours de calcul formel dispensés les années précédentes à l’Université de
Limoges par M. A. Barkatou, J.-A. Weil et d’autres enseignants-chercheurs ;

2. des références [1, 3].

78
Bibliographie

[1] A. Bostan, F. Chyzak, M. Giusti, R. Lebreton, G. Lecerf, B. Salvy et É. Schost Algo-

rithmes Efficaces en Calcul Formel. Notes du cours 2-22 du Master Parisien de Recherche
en Informatique (MPRI), Version du 26/01/2016.

[2] J. von zur Gathen et J. Gerhard. Modern computer algebra. Cambridge University
Press, New York, 1999.

[3] J.-P. Marco, P. Thieullen et J.-A. Weil. Mathématiques L2 (Chapitre 5) - Cours complet
avec 700 tests et exercices corrigés, 838 pages. Pearson Education France, 2007.

[4] M. Mignotte. Mathématiques pour le calcul formel. Presses Universitaires de France,

1989.

79