Résultat à virgule flottante différent avec l'optimisation activée - bug du compilateur?

Question

Le code ci-dessous fonctionne sur Visual Studio 2008 avec et sans optimisation. Mais il ne fonctionne que sur g++ sans optimisation (O0).

#include <cstdlib>
#include <iostream>
#include <cmath>

double round(double v, double digit)
{
    double pow = std::pow(10.0, digit);
    double t = v * pow;
    //std::cout << "t:" << t << std::endl;
    double r = std::floor(t + 0.5);
    //std::cout << "r:" << r << std::endl;
    return r / pow;
}

int main(int argc, char *argv[])
{
    std::cout << round(4.45, 1) << std::endl;
    std::cout << round(4.55, 1) << std::endl;
}

Le résultat devrait être:

4.5
4.6

Mais g++ avec optimisation (O1 - O3) donnera comme résultat:

4.5
4.5

Si j'ajoute le mot clé volatile avant t, cela fonctionne, y aurait-il un bug lié à l'optimisation?

Testé sur g++ 4.1.2 et 4.4.4.

Voici le résultat sur ideone: http://ideone.com/Rz937

Et l'option que j'ai testée sur g++ est simple:

g++ -O2 round.cpp

Le résultat le plus intéressant, même si j'active l'option /fp:fast sur Visual Studio 2008, le résultat reste correct.

Question supplémentaire:

Je me demandais, devrais-je toujours activer l'option -ffloat-store?

Parce que la version de g++ que j'ai testée est fournie avec CentOS/Red Hat Linux 5 et CentOS/Redhat 6.

J'ai compilé beaucoup de mes programmes sur ces plateformes, et je crains que cela ne provoque des bugs inattendus dans mes programmes. Il semble un peu difficile d'enquêter sur tout mon code C++ et les bibliothèques utilisées pour savoir s'ils ont de tels problèmes. Des suggestions?

Est-ce que quelqu'un s'intéresse au fait que même avec /fp:fast activé, Visual Studio 2008 fonctionne toujours? Il semble que Visual Studio 2008 soit plus fiable sur ce problème que g++?

Answer 1

J'ai creusé davantage dans ce problème et je peux apporter plus de précisions. Tout d'abord, les représentations exactes de 4,45 et 4,55 selon gcc sur x84_64 sont les suivantes (avec libquadmath pour afficher la dernière précision):

float 32:   4,44999980926513671875
double 64:  4,45000000000000017763568394002504646778106689453125
doublex 80: 4,449999999999999999826527652402319290558807551860809326171875
quad 128:   4,45000000000000000000000000000000015407439555097886824447823540679418548304813185723105561919510364532470703125

float 32:   4,55000019073486328125
double 64:  4,54999999999999982236431605997495353221893310546875
doublex 80: 4,550000000000000000173472347597680709441192448139190673828125
quad 128:   4,54999999999999999999999999999999984592560444902113175552176459320581451695186814276894438080489635467529296875

Comme Maxim l'a dit ci-dessus, le problème est dû à la taille de 80 bits des registres FPU.

Mais pourquoi le problème ne se produit-il jamais sur Windows? Sur IA-32, le FPU x87 était configuré pour utiliser une précision interne pour la mantisse de 53 bits (équivalente à une taille totale de 64 bits: double). Pour Linux et Mac OS, la précision par défaut de 64 bits était utilisée (équivalente à une taille totale de 80 bits: long double). Donc le problème pourrait être possible, ou pas, sur ces différentes plateformes en changeant le mot de contrôle du FPU (en supposant que la séquence d'instructions déclencherait le bug). Le problème a été signalé à gcc comme bug 323 (lire au moins le commentaire 92!).

Pour afficher la précision de la mantisse sur Windows, vous pouvez compiler ceci en 32 bits avec VC++:

#include "stdafx.h"
#include   
#include   

int main(void)
{
    char t[] = { 64, 53, 24, -1 };
    unsigned int cw = _control87(0, 0);
    printf("mantisse est de %d bits\n", t[(cw >> 16) & 3]);
}

et sur Linux/Cygwin:

#include 

int main(int argc, char **argv)
{
    char t[] = { 24, -1, 53, 64 };
    unsigned int cw = 0;
    __asm__ __volatile__ ("fnstcw %0" : "=m" (*&cw));
    printf("mantisse est de %d bits\n", t[(cw >> 8) & 3]);
}

Notez qu'avec gcc vous pouvez définir la précision du FPU avec -mpc32/64/80, bien que cela soit ignoré sous Cygwin. Mais gardez à l'esprit que cela modifiera la taille de la mantisse, mais pas celle de l'exposant, laissant la porte ouverte à d'autres types de comportements différents.

Sur l'architecture x86_64, SSE est utilisé comme l'a dit tmandry, donc le problème ne se produira pas à moins que vous ne forciez l'ancien FPU x87 pour les calculs en FP avec -mfpmath=387, ou à moins que vous ne compiliez en mode 32 bits avec -m32 (vous aurez besoin du package multilib). J'ai réussi à reproduire le problème sur Linux avec différentes combinaisons de drapeaux et de versions de gcc:

g++-5 -m32 floating.cpp -O1
g++-8 -mfpmath=387 floating.cpp -O1

J'ai essayé quelques combinaisons sur Windows ou Cygwin avec VC++/gcc/tcc mais le bug n'est jamais apparu. Je suppose que la séquence d'instructions générée n'est pas la même.

Enfin, notez qu'un moyen exotique d'éviter ce problème avec 4,45 ou 4,55 serait d'utiliser _Decimal32/64/128, mais le support est vraiment rare... J'ai passé beaucoup de temps juste pour pouvoir faire un printf avec libdfp !

Answer 2

Personnellement, j'ai rencontré le même problème en passant de gcc à VS. Dans la plupart des cas, je pense qu'il est préférable d'éviter l'optimisation. La seule fois où cela vaut la peine, c'est lorsque vous travaillez avec des méthodes numériques impliquant de grandes quantités de données en virgule flottante. Même après avoir désassemblé, je suis souvent déçu par les choix des compilateurs. Très souvent, il est juste plus facile d'utiliser des intrinsèques du compilateur ou d'écrire directement l'assemblage.