Somme pythagoricienne

En mathématiques la somme pythagoricienne de deux nombres a et b est le nombre √a²+b², que l'on peut voir comme la longueur de l'hypoténuse d'un triangle rectangle de côtés de longueurs a et b par le théorème de Pythagore. L'opération associée, l'addition pythagoricienne, est commutative et associative. La somme pythagoricienne peut être calculée efficacement sur machine par un algorithme donné par Moler et Morisson, qui s'avère être une application de la méthode de Halley.

Définition et propriétés

Si on note a ⊕ b = √a²+b², l'opération ainsi définie est commutative et associative, et :

a₁ ⊕ a₂ ⊕ … ⊕ a_n = √a₁² + a₂² + … + a_n².

Calcul numérique

La somme pythagoricienne des longueurs a et b des deux côtés de l'angle droit d'un triangle rectangle est la longueur de l'hypoténuse d'après le théorème de Pythagore. Le calcul de celle-ci donne la norme euclidienne en dimension 2, et plus généralement dans un espace de dimension quelconque. Elle permet ainsi de calculer le module d'un nombre complexe, ou de passer d'une représentation cartésienne à une représentation polaire.

Le calcul peut procèder par élévation au carré de a et b puis racine carrée de la somme, mais si a ou b prennent des valeurs très élevées, le calcul intermédiaire de a² ou de b² sur ordinateur peut créer une erreur de dépassement ou de soupassement.

Typiquement, le nombre le plus grand représentable en simple précision est de l'ordre de 3 × 10³⁸ ; ainsi, si a ou b est supérieur à 2 × 10¹⁹, le logiciel déclarera le résultat « infini » (erreur de dépassement), même si le résultat final est représentable (par exemple si l'autre longueur est petite en valeur absolue) ; par exemple, ${\sqrt {(2\times 10^{19})^{2}+1^{2}}}\simeq 2\times 10^{19}$. Avec une représentation en double précision, le nombre le plus grand représentable est environ 1,8 × 10³⁰⁸, donc le résultat sera déclaré « infini » si a ou b est supérieur à 1,3 × 10¹⁵⁴.

À l'inverse, la valeur absolue la plus petite représentable est environ 1,1 × 10⁻³⁸ en simple précision et 2,2 × 10⁻³⁰⁸ en double précision. Ainsi, en simple précision, le calcul pour a = b = 2 × 10⁻³⁸ donnera 0 (erreur de soupassement), alors que la valeur attendue est représentable et vaut environ 1,4 × 10⁻³⁸.

La méthode de Moler et Morrison, une méthode itérative cubique (on gagne à peu près 3 chiffres significatifs à chaque itération), dérivée de la méthode de Halley, évite ce problème.

Cette méthode consiste en se ramenant d'abord à a ≥ b, puis à remplacer a et b, par deux nombres p et q positifs avec p > a ≥ b > q ayant même somme pythagoricienne, c'est-à-dire que √a²+b² = √p²+q². À chaque étape, l'algorithme diminue la valeur de q et augmente la valeur de p tout en conservant la somme pythagoricienne. Lorsque la valeur de q devient très petite, on a donc √a²+b² ≈ p. Concrètement, l'algorithme s'écrit en pseudo-code :

fonction [résultat] = somme_pythagoricienne(a, b)

// initialisation
p ← max(|a|, |b|)
q ← min(|a|, |b|)

// Calcul
tant que (q a une valeur significative) faire
  r ← (q/p)^2
  s ← r/(4 + r)
  p ← p + 2*s*p
  q ← s*q
fin de boucle

résultat ← p

fin de fonction

Pour un calcul de norme euclidienne dans un espace de dimension supérieure à 2, il suffit de calculer la somme pythagoricienne de manière successive. Pour un vecteur X de dimension n :

fonction [résultat] = norme2(X)

// initialisation
n ← dimension(X)
s ← 0

// Calcul
pour i allant de 1 à n
  s ← somme_pythagoricienne(s, X(i))
fin de boucle

résultat ← s

fin de fonction

Bibliographie

• (en) Cleve Moler et Donald Morrison, « Replacing Square Roots by Pythagorean Sums », IBM Journal of Research and Development, vol. 27, n^o 6,‎ 1983, p. 577–581 (DOI 10.1147/rd.276.0577, CiteSeer^x 10.1.1.90.5651, lire en ligne).
• (en) Augustin A. Dubrulle, « A Class of Numerical Methods for the Computation of Pythagorean Sums », IBM Journal of Research and Development, vol. 27, n^o 6,‎ 1983, p. 582–589 (DOI 10.1147/rd.276.0582, CiteSeer^x 10.1.1.94.3443, lire en ligne).