Description (générale)

Un ordinateur n'a qu'une mémoire finie, et ne peux représenter une infinité de nombres. Si vous ajoutez 0.1 + 0.1 en C++, vous n'obtiendrez probablement pas le résulat espéré 0.2. Le problème est que ce que vous obtiendrez dépend de nombreux paramètres. Avec une configuration, vous pouvez obtenir 0.1995, et avec une autre, 0.1999. Pour la grande majorité des programmes ces petites différences n'ont pas d'importance. Mais si vous désirez reproduire les mêmes résultats 2 fois de suite, pour une expérience scientifique par exemple, et ceci sur différentes machines, voire même sur une unique machine mais avec différentes options, alors vous devez être plus prudents. Beaucoup plus prudents, car ces petites différences ont parfois tendance à s'accumuler assez vite.

La bibliothèque de calcul flottants reproductible qui est le sujet de cette page vous permet de contrôler comment sont effectués les calculs en C++. Le but est de rendre vos programmes fiables et reproductibles.

Description (détaillée)

Les résultats de calculs en virgule flottante dépendent beaucoup de l'implémentation matérielle du processeur (FPU), du compilateur et de ses optimisations, et de la bibliothèque mathématique du système (libm). Une expérience est souvent reproductible seulement sur la même machine avec la même bibliothèque mathématique, le même compilateur utilisant les mêmes options. Et même là encore, le standard C++ ne grarantit PAS des résultats reproductibles. Exemple:

    double x = 1.0;  x /= 10.0;
    double y = x;
    ...
    if (y == x) {
        // CECI N'EST PAS TOUJOURS VRAI !
    }

Un problème similaire et plus général est la génération de nombres aléatoires, souvents nécessaires pour rendre une expérience reproductible:

    set_random_seed(42);
    double x = random_number();

Ce code peut ne pas retourner le même x en fonction du processeur (FPU), du système, du compilateur, etc.

Ces problèmes sont liés à:

• Certains processeurs (comme x87 sous Linux) gardent par défaut une précision interne (80 bits) plus grande que la taille mémoire du type C++. Dans le premier exemple, si x est sur la pile mais y dans un registre, la comparaison va échouer. Pire, le fait que x ou y soit ou non dans un registre dépend de ce qu'il y a dans la section "...". Ce problème peut être résolu en restreignant la précision interne du processeur à la taille mémoire du type C++. Malheureusement, sur x87, certaines opérations (commes les fonctions transcendentales) sont toujours calculées à la précision maximum.
• La façon dont votre processeur implémente la norme IEEE754, et en particulier, les opérations sur les nombres dénormalisés. L'exemple arithmetic_test fournit avec le code source produit systématiquement des résultats différents entre les processeurs SSE et x87. Ce problème n'est PAS résolu de manière générale en restreignant la précision interne du processeur, bien qu'il y ait des relations entre ces deux problèmes (en particulier, un flottant dénormalisé simple précision peut parfois être représenté par un flottant normalisé en double précision).
• La façon dont votre compilateur interagit avec la bibliothèque mathématique du système (libm), et surtout en fonction des paramètres d'optimizations. En particulier, g++-3.4 et g++-4.0 peuvent produire des résultats différents sur le même système et avec la même bibliothèque. Ce problème est partiellement résolu en sélectionnant des options adéquates du compilateur.
• Finalement, le standard IEEE754 laisse aussi une certaine marge pour l'implémentation des types NaN. En particulier, la sérialization binaire de résultats numérique NaN peut différer.Ce problème est résolu en comparant les significations logiques des types NaN et non leur représentation binaire.

D'autres points à considérer:

• Le processeur SSE a une option pour accélérer les calculs en approximant les nombres dénormaux par 0. Pour certaines application ceci est carrément faux, mais pour d'autres applications comme le traitement du signal (DSP) en temps réel, c'est actactement ce qu'il faut (certains rapports font état d'une dégradation par 30 des performances lors de multiplications de nombres dénormalisés par rapport à des nombre normalisés). En raison du mode par défaut "arrondi à plus proche", les nombres dénormalisés tendent à persister plutôt qu'à s'annuler vers 0. La désactivation des nombres dénormalisés est prévue dans le processeur SSE, mais malheureusement ceci n'est pas reproductible sur 87 en natif. Dans ce cas, il est possible de verifier l'état d'un résultat après chaque opération et affectation, et d'émuler l'approximation par 0 à l'aide d'une classe C++ encapsulant le type de base. Note: ceci résoud aussi le problème ci-dessus de non reproductibilité des nombres dénormaux entre SSE et x87.
• Aucune dépendence externe. Plus il y a de dépendences, et plus il y a de risques d'une mauvaise version ou configuration. Cette bibliothèque est indépendente ; elle fournit l'intégralité de la libm sans reposer sur des liens externes. De cette façon, elle peut être incluse telle quelle dans un projet, avec configuration et spécialisation minimale.

La bibliothèque

Téléchargez la dernière version 0.3 sous la forme d'une archive tar.bz2 ou bien comme tar.gz.

Fonctionalités

• Support pour x87 (vieux et nouveaux PCs) et SSE (nouveaux PCs seulement). Support d'une implémentation IEEE754 logicielle.
• Support pour les types Simple, Double and Extended correspondant aux types C++ natif (float, double, long double), mais qui en plus gèrent la précision interne du processeur et les nombres dénormaux. Ces types sont selon la configuration ou bien des alias vers les types natifs, ou bien des classes C++ les encapsulant.
• Jeu complet de fonctions de la libm, et non seulement une sous-partie, en tout cas pour les precisions Simple et Double.
• Génération de nombres aléatoires, grâce à l'excellent Mersenne Twister créé par Takuji Nishimura et Makoto Matsumoto.
• Pas de dépéndences externes, ce qui rend la bibliothèque facile à intégrer dans un projet.
• Support pour supprimer les nombres dénormaux, y compris avec x87. Ceci augmente les performances du processeur SSE, mais ralentit x87 (sauf si vous utilisez vraiment beaucoup de dénormaux, auxquel cas le coût de l'encapsulation est compensé par le gain sur la vitesse des opérations).
• Bonus: Gestion des modes d'arrondi de des NaN avec des fonctions C99, et support de la libm y compris avec l'implémentation logicielle (softfloat) des nombres en virgule flottante (vive le C++!).

Grille de configuration

Une case de cette grille doit être sélectionnée lors de la configuration. Tous les types présents dans cette case sont alors disponible dans le programme.
Les types marqués par * sont en fait des alias vers les types natifs float / double / long double, et le bon fonctionnement est assuré par support matériel (FPU).
Les autres types sont des classes d'encapsulation qui se comportent comme des types de base.

Excepté pour la représentation binaire des valeurs NaN:

• Des résultats identiques sont attendus pour la même précision des cases "SSE avec dénormaux" et "Logiciel avec dénormaux".
• Des résultats identiques sont attendus pour la même précision des cases "SSE sans dénormaux" et "x87 sans dénormaux".
• Des résultats identiques sont attendus pour "Extended, x87 avec dénormaux" et "Extended, Logiciel avec dénormaux".
• Des résultats presque identiques sont attendus pour la même précision des cases "SSE avec dénormaux" et "x87 avec dénormaux". Les différences interviennent de façon peu fréquente, et seulement pour les calculs faisant intervenir des nombres dénormaux.
• Toutes les autres combinaisons donnent des résultats différents.

Critères de comparaison

• La meilleure performance est obtenue par la combinaison "SSE sans dénormaux / Simple". Ce qui importe le plus pour la performance est d'abord l'utilisation de types natifs, puis la taille mémoire, et enfin la présence ou non de dénormaux (sauf si vous utilisez beaucoup de dénormaux, auxquel cas leur désactivation devient le critère le plus important).
• La meilleure précision est obtenue par les types Extended en modes avec dénormaux. Ce qui importe le plus pour la précision est d'abord la taille mémoire, puis l'utilisation de dénormaux.
• La meilleure conformance au standard IEEE754 est obtenue par l'implémentation Logiciel (et les résultats équivalents avec SSE/x87). La présence de nombre dénormaux est requise par le standard.

Documentation

Remerciements

Ce code repose en grande partie sur la bibliothèque mathématique GNU Libm, qui utilise elle-même des portions de la netlib fdlibm de Sun, ainsi que de GNU MP, et de la "IBM's accurate portable mathematical library".

Ce code utilise également SoftFloat, la bibliothèque de référence pour l'implémentation logicielle des calculs en virgule flottante.

Le générateur de nombres aléatoires est le Mersene Twister, créé par Takuji Nishimura et Makoto Matsumoto, et adapté en C++ pour ce project. Merci de lire la license (style BSD) dans le fichier Random.cpp si vous avez l'intention de créer des paquets binaires de la présente bibliothèque Streflop (sans oublier d'être conforme à la LGPL).

Merci de lire l'historique et les informations pour le "copyright" dans les fichiers README.txt des sous-répertoires libm et softfloat, ainsi que le fichier LGPL.txt dans le répertoire principal. Faites ce que vous voulez avec cette bibliothèque, mais à vos propres risques, et conformément à la LGPL (voire le fichier correspondant).

Merci à Tobi Vollebregt pour le retour, le portage vers Win32, les tests, et des patches, particulièrement à propos des modes d'arrondi.

Usage

Inclure "streflop.h" à la place de <math.h>, et lier avec streflop.a au lieu de libm. Vous pouvez sans risque utiliser libm dans une autre partie du même programme, la confusion n'est pas possible lors de l'édition des liens. Cependant il est important d'inclure le bon fichier d'en-tête.

Utilisez l'espace de nommage streflop, et les types Simple, Double et Extended selon les besoins, au lieu de float, double, et long double. Les types streflop peuvent selon les cas n'être de toutes façon que de simples alias vers les types de base (voire la grille de configuration).

Vous devez également appeler streflop_init<FloatType> où FloatType = Simple, Double, ou Extended, avant d'utiliser ce type. Vous ne devez alors utiliser que ce type (ex: Simple) jusqu'au prochain appel de streflop_init. En d'autres termes, séparez votre code en blocs logiques n'utilisant qu'un seul type à la fois, et passez d'un bloc à l'autre grâce à streflop_init. Dans le cas le plus simple, choisissez un type unique et appelez streflop_init en début de programme. Le but de ces fonctions d'initialisation est de configurer correctement le processeur (FPU), ce qui est particulièrement important pour les alias vers les types natifs.

Editez Makefile.common pour choisir une case dand la grille de configuration. Compilez, et liez avec streflop.a au lieu de libm.

Merci de lire le fichier README.txt et l'exemple arithmeticTest.cpp qui sont fournis avec le code source.

Notes

Prenez garde aux options d'optimisation trop aggressives! En particulier, ce code repose sur reinterpret_cast et des unions, et le compilateur ne doit pas présupposer que la mémoire est strictement réservée pour une unique variable (strict aliasing). Ceci est malheureusement le cas pour g++ dans les niveaux d'optimisation 2 et 3. De même, le compilateur ne doit pas présupposer que les NaN peuvent être ignorés, ou que le processeur est dans mode d'arrondi constant. Ex: -O3 -fno-strict-aliasing -frounding-math -fsignaling-nans.

Vous devez si nécessaire spécifier les options correctes de sélection du processeur (FPU), comme par exemple -mfpmath=sse -msse -msse2. Le -msse2 est important, il y a des cas où g++ refuse d'utiliser SSE (and retombe silencieusement sur x87) lorsque using -msse est spécifié sans -msse2. Ceci signifie aussi que vous ne pouvez pas utiliser le mode SSE de cette bibliothèque sur des processeurs où seul SSE1 est présent (sans SSE2, comme certains athlon-xp par exemple).

La libm du système produira quasiment certainement des résultats différents selon le FPU, le compilateur, les options, etc. L'idée de cette bibliothèque est d'augmenter la reproductibilité de vos expériences comparé à la libm du système, en tenant compte des limitations dues aux dénormaux, etc. Si vous voulez des résultats guarantis (mais bien plus lents) indépendamment du type de machine, du FPU, sans se soucier des dénormaux ou d'autres choses encore, alors veuillez utiliser une implémentation purement logicielle comme GNU MP. Si vous désirez utilisez votre matériel tout en gardant une certaine reproductibilité, alors cette bibliothèque est pour vous. Si elle ne remplit pas complètement vos besoins, alors améliorez-là: après tout, ceci est un logiciel libre :)

Les fonctionnalités C99 suivantes sont implémentées en ce qui concerne les modes d'arrondis et la gestion des NaN (y compris pour le mode Logiciel): fe(get|set)round, fe(get|set)env, et feholdexcept. Vous pouvez les utilisez pour changer les mode d'arrondi ou générer des exceptions sur les NaNs. Ces fonctions devraient fonctionner normalement, tout du moins pour autant que le FPU fonctionne normalement*, mais elles n'ont pas été testées extensivement. * en particulier certains rapports font état que x87 échoue parfois à détecter certaines conditions faisant intervenir des nombres dénormaux.

Il y a bien sûr la possibilité que des bugs inconnus soient présents. Comme cette bibliothèque est basée sur GNU libm 2.4, tout bug trouvé dans cette version a des chance de se retrouver aussi dans streflop.

Le support pour le mode Extended est incomplet. Il manque des fonctions, en particulier concenant la trigonométrie. L'implémentation ldbl-96 de la libm ne contient pas ces fonctions. Comme streflop assure la séparation stricte des précisions Extended et Double, ces fonctions manquantes ont été ici implémentées en passant temporairement en mode Double grâce à streflop_init<Double>, puis en appelant la fonction Double correspondante pour obtenir un résultat temporaire sur la pile, puis en revenant en mode Extended par streflop_init<Extended>, pour enfin convertir le résultat en nombre Extended.