L’écriture d’un entier A en base N ressemble à celle d’un polynôme :

A = a + bN + cN ² +… où a, b, c… sont des nombres entiers compris entre 0 et N – 1. Pour additionner deux nombres à n chiffres, en utilisant l’algorithme classique, enseigné à l’école élémentaire, on effectue n additions de nombres à un chiffre, avec, éventuellement, n additions correspondant aux retenues. On vérifie effectivement qu’en utilisant un ordinateur le temps de calcul pour de très grands nombres est proportionnel à n. Plus précisément, il est majoré par une constante multipliée par n, ce que l’on note O(n) (« grand o de n »). On dit que sa complexité – ce qui est un modèle mathématique du temps de calcul – est linéaire.

En utilisant l’algorithme classique, enseigné dans les écoles, la multiplication n’est pas linéaire mais quadratique, c’est-à-dire que sa complexité est un O(n²) car on effectue alors n ² multiplications de nombres à un chiffre, puis environ n additions. On vérifie effectivement qu’en utilisant un ordinateur le temps de calcul pour des grands nombres est proportionnel à n ².

Les ordinateurs sont devenus si puissants que ces mesures ne deviennent visibles que pour des nombres de plusieurs milliers de chiffres. Mais si vous doublez le nombre de chiffres, vous multipliez le temps de calcul par 4. Dans le cas de l’addition, vous ne faites que doubler le temps de calcul ; là est toute la différence entre complexité linéaire et complexité quadratique.

L’algorithme de Strassen

En 1971, le mathématicien allemand Volker Strassen (né en 1936) a trouvé un algorithme dont la complexité est en O(n (ln n) (ln ln n)), où ln désigne le logarithme népérien (voir Tangente 177, 2017). Strassen a en outre conjecturé qu’il existait un algorithme quasi linéaire, c’est-à-dire de complexité O(n ln n).

L’algorithme de Strassen est fondé sur la notion de transformée de Fourier discrète. À un nombre A de n chiffres, la transformée de Fourier discrète T associe le n-uplet de nombres complexes T(A) suivants :

(A(1), A(ω), A(ω ²)… A(ω ^n–1)), où A(X) = a + bX + cX ² + … 

et  \( \omega = \cos \dfrac{2\pi}{n} + i \sin \dfrac{2\pi}{n} \) (voir FOCUS sur les racines de l’unité).

Ce calcul requiert a priori n² multiplications, son temps de calcul est donc quadratique.

Une méthode plus subtile permet d’effectuer ce calcul dans un temps quasi linéaire, ce qui est beaucoup plus rapide. Par exemple, en termes d’ordres de grandeur, cela fait passer d’un million d’opérations à sept mille opérations, ou encore d’un milliard d’opérations à trois cent mille opérations. Avant d’examiner cette méthode, voyons l’intérêt de la transformation de Fourier discrète T.

Le principal avantage de T est d’être inversible, c’est-à-dire que la connaissance du n-uplet T(A) permet de reconstituer l’entier A. Le calcul consiste simplement à résoudre le système d’équations suivant :
a + b + c + … = A(1),

a + b ω + c ω ² + … = A(ω),

…

a + b ω ^n–1 + c ω ^2 (n–1) + … = A(ω ^n–1),

où les inconnues sont a, b, c… Comme la somme des racines n^ièmes de l’unité est nulle (voir FOCUS), il suffit d’ajouter toutes les équations pour obtenir a :

\( a = \dfrac{\text{A} (1) + \text{A}(\omega) + \text{A}(\omega^2) + … + \text{A}(\omega^{n-1}) } {n}.\)  

De même, on obtient b en multipliant la dernière équation par ω, l’avant-dernière par ω ², et ainsi de suite. Par exemple :

\( b = \dfrac{\text{A} (1) + \omega^{n-1} \text{A}(\omega) + \omega^{n-2} \text{A} (\omega^2) + … + \omega \text{A} (\omega ^{n-1}) } {n}.\)   

 Le passage de T(A) à A est de même nature que celui de A à T(A), la seule différence est la division finale par n. La méthode rapide permet de l’effectuer dans un temps quasi linéaire.

Gagner du temps dans les calculs

Étant donnés deux entiers A et B, on peut calculer leurs transformées de Fourier discrètes T(A) et T(B) dans un temps quasi linéaire. On les multiplie alors pour obtenir C(1) = A(1) B(1), C(ω) = A(ω) B(ω), etc., ce qui demande n multiplications. Il est alors facile de reconstituer C. Le calcul montre que l’on manque la quasi-linéarité de peu, puisque l’on obtient une complexité en O(n ( ln n) ( ln ln n)). On obtient donc le produit C = A × B en un temps proportionnel à n × ( ln n) × (ln ( ln n)) au lieu de n ² avec l’algorithme traditionnel appris à l’école. L’intérêt principal de la transformée de Fourier discrète se situe là : gagner du temps dans les calculs arithmétiques, qui sont à la base de tous les autres.

Soit 2m le nombre de chiffres de A (au besoin, on ajoute un zéro en tête), et séparons les termes d’ordre pair de ceux d’ordre impair. En considérant les polynômes associés, on obtient A(X) = A₀( X ² ) + X A ₁( X ²) où A ₀ et A ₁ sont de degrés m. 

Par exemple, le polynôme 7X⁵ – X ⁴ + 2X³ + 11X² – X + 4 peut s’écrire
(–Y² + 11Y + 4) + X( 7Y² + 2Y – 1) avec Y = X ².

 La formule  

\( \omega = \cos \dfrac{2\pi}{n} + i \sin \dfrac{2\pi}{n}, \)

implique  

\( \omega ^2 = \cos \dfrac{2 \pi}{m} + i \sin \dfrac{2\pi}{m}, \)

donc : A(1) = A₀(1) + A₁(1),

           A(ω) = A₀( ω ²) + ω A₁( ω ²),

           …

A(ω ^m–1) = A₀(ω ^2 (m–1)) + ω ^m–1 A₁(ω ^{2 ( m–1)}).

Ces formules montrent que, pour calculer une transformée de Fourier discrète d’un nombre à 2m chiffres, il suffit d’en calculer deux de nombres à m chiffres puis d’effectuer m multiplications et m additions. Si l’on note S(m) le nombre de multiplications nécessaires pour calculer la transformée de Fourier discrète d’un nombre de m chiffres, la fonction S vérifie S(2m) = 2 S(m) + m. En résolvant cette relation de récurrence, on s’assure que S(m) est de l’ordre de grandeur de m ln m (voir FOCUS). Le calcul se fait donc de façon récursive : celui de T(A) appelle ceux de T(A₀) et T(A₁), et ainsi de suite, jusqu’au cas évident où les nombres n’ont qu’un seul chiffre. On obtient ainsi T(A) et on en déduit A. Pour obtenir un algorithme de multiplication quasi linéaire, David Harvey (de l’université de Nouvelle-Galles du Sud, Australie) et Joris van der Hoeven (de l’École polytechnique et du CNRS) utilisent des anneaux de polynômes à plusieurs variables (au lieu de la seule variable X pour nos polynômes A, B, C, A₀ et A₁) dans lesquels la multiplication est particulièrement efficace. L’utilité du nouvel algorithme est d’accélérer tous les gros calculs, ceux des décimales de π comme ceux des prévisions météorologiques.