Pourquoi une conversion aller-retour via une chaîne n'est-elle pas sûre pour un double?

185

Récemment, j'ai dû sérialiser un double en texte, puis le récupérer. La valeur ne semble pas être équivalente:

double d1 = 0.84551240822557006;
string s = d1.ToString("R");
double d2 = double.Parse(s);
bool s1 = d1 == d2;
// -> s1 is False

Mais selon MSDN: Standard Numeric Format Strings , l'option "R" est censée garantir la sécurité aller-retour.

Le spécificateur de format aller-retour ("R") est utilisé pour garantir qu'une valeur numérique convertie en chaîne sera analysée à nouveau dans la même valeur numérique.

Pourquoi est-ce arrivé?

Philippe Ding
la source
6
J'ai débogué dans mon VS et son retour est vrai ici
Neel
19
Je l'ai reproduit en retournant faux. Question très intéressante.
Jon Skeet
40
.net 4.0 x86 - true, .net 4.0 x64 - false
Ulugbek Umirov
25
Félicitations pour avoir trouvé un bogue aussi impressionnant dans .net.
Aron
14
@Casperah Round trip est spécifiquement destiné à éviter les incohérences en virgule flottante
Gusdor

Réponses:

178

J'ai trouvé le bogue.

.NET effectue les opérations suivantes dans clr\src\vm\comnumber.cpp:

DoubleToNumber(value, DOUBLE_PRECISION, &number);

if (number.scale == (int) SCALE_NAN) {
    gc.refRetVal = gc.numfmt->sNaN;
    goto lExit;
}

if (number.scale == SCALE_INF) {
    gc.refRetVal = (number.sign? gc.numfmt->sNegativeInfinity: gc.numfmt->sPositiveInfinity);
    goto lExit;
}

NumberToDouble(&number, &dTest);

if (dTest == value) {
    gc.refRetVal = NumberToString(&number, 'G', DOUBLE_PRECISION, gc.numfmt);
    goto lExit;
}

DoubleToNumber(value, 17, &number);

DoubleToNumberest assez simple - il appelle simplement _ecvt, qui est dans le runtime C:

void DoubleToNumber(double value, int precision, NUMBER* number)
{
    WRAPPER_CONTRACT
    _ASSERTE(number != NULL);

    number->precision = precision;
    if (((FPDOUBLE*)&value)->exp == 0x7FF) {
        number->scale = (((FPDOUBLE*)&value)->mantLo || ((FPDOUBLE*)&value)->mantHi) ? SCALE_NAN: SCALE_INF;
        number->sign = ((FPDOUBLE*)&value)->sign;
        number->digits[0] = 0;
    }
    else {
        char* src = _ecvt(value, precision, &number->scale, &number->sign);
        wchar* dst = number->digits;
        if (*src != '0') {
            while (*src) *dst++ = *src++;
        }
        *dst = 0;
    }
}

Il s'avère que _ecvtretourne la chaîne 845512408225570.

Remarquez le zéro final?Cela fait toute la différence!
Lorsque le zéro est présent, le résultat est analysé en retour 0.84551240822557006, qui est votre numéro d' origine - il se compare donc égal, et par conséquent, seuls 15 chiffres sont renvoyés.

Cependant, si je tronque la chaîne à ce zéro à 84551240822557 , alors je reviens 0.84551240822556994, ce qui n'est pas votre numéro d'origine, et par conséquent, il renverrait 17 chiffres.

Preuve: exécutez le code 64 bits suivant (dont la plupart j'ai extrait de Microsoft Shared Source CLI 2.0) dans votre débogueur et examinez và la fin de main:

#include <stdlib.h>
#include <string.h>
#include <math.h>

#define min(a, b) (((a) < (b)) ? (a) : (b))

struct NUMBER {
    int precision;
    int scale;
    int sign;
    wchar_t digits[20 + 1];
    NUMBER() : precision(0), scale(0), sign(0) {}
};


#define I64(x) x##LL
static const unsigned long long rgval64Power10[] = {
    // powers of 10
    /*1*/ I64(0xa000000000000000),
    /*2*/ I64(0xc800000000000000),
    /*3*/ I64(0xfa00000000000000),
    /*4*/ I64(0x9c40000000000000),
    /*5*/ I64(0xc350000000000000),
    /*6*/ I64(0xf424000000000000),
    /*7*/ I64(0x9896800000000000),
    /*8*/ I64(0xbebc200000000000),
    /*9*/ I64(0xee6b280000000000),
    /*10*/ I64(0x9502f90000000000),
    /*11*/ I64(0xba43b74000000000),
    /*12*/ I64(0xe8d4a51000000000),
    /*13*/ I64(0x9184e72a00000000),
    /*14*/ I64(0xb5e620f480000000),
    /*15*/ I64(0xe35fa931a0000000),

    // powers of 0.1
    /*1*/ I64(0xcccccccccccccccd),
    /*2*/ I64(0xa3d70a3d70a3d70b),
    /*3*/ I64(0x83126e978d4fdf3c),
    /*4*/ I64(0xd1b71758e219652e),
    /*5*/ I64(0xa7c5ac471b478425),
    /*6*/ I64(0x8637bd05af6c69b7),
    /*7*/ I64(0xd6bf94d5e57a42be),
    /*8*/ I64(0xabcc77118461ceff),
    /*9*/ I64(0x89705f4136b4a599),
    /*10*/ I64(0xdbe6fecebdedd5c2),
    /*11*/ I64(0xafebff0bcb24ab02),
    /*12*/ I64(0x8cbccc096f5088cf),
    /*13*/ I64(0xe12e13424bb40e18),
    /*14*/ I64(0xb424dc35095cd813),
    /*15*/ I64(0x901d7cf73ab0acdc),
};

static const signed char rgexp64Power10[] = {
    // exponents for both powers of 10 and 0.1
    /*1*/ 4,
    /*2*/ 7,
    /*3*/ 10,
    /*4*/ 14,
    /*5*/ 17,
    /*6*/ 20,
    /*7*/ 24,
    /*8*/ 27,
    /*9*/ 30,
    /*10*/ 34,
    /*11*/ 37,
    /*12*/ 40,
    /*13*/ 44,
    /*14*/ 47,
    /*15*/ 50,
};

static const unsigned long long rgval64Power10By16[] = {
    // powers of 10^16
    /*1*/ I64(0x8e1bc9bf04000000),
    /*2*/ I64(0x9dc5ada82b70b59e),
    /*3*/ I64(0xaf298d050e4395d6),
    /*4*/ I64(0xc2781f49ffcfa6d4),
    /*5*/ I64(0xd7e77a8f87daf7fa),
    /*6*/ I64(0xefb3ab16c59b14a0),
    /*7*/ I64(0x850fadc09923329c),
    /*8*/ I64(0x93ba47c980e98cde),
    /*9*/ I64(0xa402b9c5a8d3a6e6),
    /*10*/ I64(0xb616a12b7fe617a8),
    /*11*/ I64(0xca28a291859bbf90),
    /*12*/ I64(0xe070f78d39275566),
    /*13*/ I64(0xf92e0c3537826140),
    /*14*/ I64(0x8a5296ffe33cc92c),
    /*15*/ I64(0x9991a6f3d6bf1762),
    /*16*/ I64(0xaa7eebfb9df9de8a),
    /*17*/ I64(0xbd49d14aa79dbc7e),
    /*18*/ I64(0xd226fc195c6a2f88),
    /*19*/ I64(0xe950df20247c83f8),
    /*20*/ I64(0x81842f29f2cce373),
    /*21*/ I64(0x8fcac257558ee4e2),

    // powers of 0.1^16
    /*1*/ I64(0xe69594bec44de160),
    /*2*/ I64(0xcfb11ead453994c3),
    /*3*/ I64(0xbb127c53b17ec165),
    /*4*/ I64(0xa87fea27a539e9b3),
    /*5*/ I64(0x97c560ba6b0919b5),
    /*6*/ I64(0x88b402f7fd7553ab),
    /*7*/ I64(0xf64335bcf065d3a0),
    /*8*/ I64(0xddd0467c64bce4c4),
    /*9*/ I64(0xc7caba6e7c5382ed),
    /*10*/ I64(0xb3f4e093db73a0b7),
    /*11*/ I64(0xa21727db38cb0053),
    /*12*/ I64(0x91ff83775423cc29),
    /*13*/ I64(0x8380dea93da4bc82),
    /*14*/ I64(0xece53cec4a314f00),
    /*15*/ I64(0xd5605fcdcf32e217),
    /*16*/ I64(0xc0314325637a1978),
    /*17*/ I64(0xad1c8eab5ee43ba2),
    /*18*/ I64(0x9becce62836ac5b0),
    /*19*/ I64(0x8c71dcd9ba0b495c),
    /*20*/ I64(0xfd00b89747823938),
    /*21*/ I64(0xe3e27a444d8d991a),
};

static const signed short rgexp64Power10By16[] = {
    // exponents for both powers of 10^16 and 0.1^16
    /*1*/ 54,
    /*2*/ 107,
    /*3*/ 160,
    /*4*/ 213,
    /*5*/ 266,
    /*6*/ 319,
    /*7*/ 373,
    /*8*/ 426,
    /*9*/ 479,
    /*10*/ 532,
    /*11*/ 585,
    /*12*/ 638,
    /*13*/ 691,
    /*14*/ 745,
    /*15*/ 798,
    /*16*/ 851,
    /*17*/ 904,
    /*18*/ 957,
    /*19*/ 1010,
    /*20*/ 1064,
    /*21*/ 1117,
};

static unsigned DigitsToInt(wchar_t* p, int count)
{
    wchar_t* end = p + count;
    unsigned res = *p - '0';
    for ( p = p + 1; p < end; p++) {
        res = 10 * res + *p - '0';
    }
    return res;
}
#define Mul32x32To64(a, b) ((unsigned long long)((unsigned long)(a)) * (unsigned long long)((unsigned long)(b)))

static unsigned long long Mul64Lossy(unsigned long long a, unsigned long long b, int* pexp)
{
    // it's ok to losse some precision here - Mul64 will be called
    // at most twice during the conversion, so the error won't propagate
    // to any of the 53 significant bits of the result
    unsigned long long val = Mul32x32To64(a >> 32, b >> 32) +
        (Mul32x32To64(a >> 32, b) >> 32) +
        (Mul32x32To64(a, b >> 32) >> 32);

    // normalize
    if ((val & I64(0x8000000000000000)) == 0) { val <<= 1; *pexp -= 1; }

    return val;
}

void NumberToDouble(NUMBER* number, double* value)
{
    unsigned long long val;
    int exp;
    wchar_t* src = number->digits;
    int remaining;
    int total;
    int count;
    int scale;
    int absscale;
    int index;

    total = (int)wcslen(src);
    remaining = total;

    // skip the leading zeros
    while (*src == '0') {
        remaining--;
        src++;
    }

    if (remaining == 0) {
        *value = 0;
        goto done;
    }

    count = min(remaining, 9);
    remaining -= count;
    val = DigitsToInt(src, count);

    if (remaining > 0) {
        count = min(remaining, 9);
        remaining -= count;

        // get the denormalized power of 10
        unsigned long mult = (unsigned long)(rgval64Power10[count-1] >> (64 - rgexp64Power10[count-1]));
        val = Mul32x32To64(val, mult) + DigitsToInt(src+9, count);
    }

    scale = number->scale - (total - remaining);
    absscale = abs(scale);
    if (absscale >= 22 * 16) {
        // overflow / underflow
        *(unsigned long long*)value = (scale > 0) ? I64(0x7FF0000000000000) : 0;
        goto done;
    }

    exp = 64;

    // normalize the mantisa
    if ((val & I64(0xFFFFFFFF00000000)) == 0) { val <<= 32; exp -= 32; }
    if ((val & I64(0xFFFF000000000000)) == 0) { val <<= 16; exp -= 16; }
    if ((val & I64(0xFF00000000000000)) == 0) { val <<= 8; exp -= 8; }
    if ((val & I64(0xF000000000000000)) == 0) { val <<= 4; exp -= 4; }
    if ((val & I64(0xC000000000000000)) == 0) { val <<= 2; exp -= 2; }
    if ((val & I64(0x8000000000000000)) == 0) { val <<= 1; exp -= 1; }

    index = absscale & 15;
    if (index) {
        int multexp = rgexp64Power10[index-1];
        // the exponents are shared between the inverted and regular table
        exp += (scale < 0) ? (-multexp + 1) : multexp;

        unsigned long long multval = rgval64Power10[index + ((scale < 0) ? 15 : 0) - 1];
        val = Mul64Lossy(val, multval, &exp);
    }

    index = absscale >> 4;
    if (index) {
        int multexp = rgexp64Power10By16[index-1];
        // the exponents are shared between the inverted and regular table
        exp += (scale < 0) ? (-multexp + 1) : multexp;

        unsigned long long multval = rgval64Power10By16[index + ((scale < 0) ? 21 : 0) - 1];
        val = Mul64Lossy(val, multval, &exp);
    }

    // round & scale down
    if ((unsigned long)val & (1 << 10))
    {
        // IEEE round to even
        unsigned long long tmp = val + ((1 << 10) - 1) + (((unsigned long)val >> 11) & 1);
        if (tmp < val) {
            // overflow
            tmp = (tmp >> 1) | I64(0x8000000000000000);
            exp += 1;
        }
        val = tmp;
    }
    val >>= 11;

    exp += 0x3FE;

    if (exp <= 0) {
        if (exp <= -52) {
            // underflow
            val = 0;
        }
        else {
            // denormalized
            val >>= (-exp+1);
        }
    }
    else
        if (exp >= 0x7FF) {
            // overflow
            val = I64(0x7FF0000000000000);
        }
        else {
            val = ((unsigned long long)exp << 52) + (val & I64(0x000FFFFFFFFFFFFF));
        }

        *(unsigned long long*)value = val;

done:
        if (number->sign) *(unsigned long long*)value |= I64(0x8000000000000000);
}

int main()
{
    NUMBER number;
    number.precision = 15;
    double v = 0.84551240822557006;
    char *src = _ecvt(v, number.precision, &number.scale, &number.sign);
    int truncate = 0;  // change to 1 if you want to truncate
    if (truncate)
    {
        while (*src && src[strlen(src) - 1] == '0')
        {
            src[strlen(src) - 1] = 0;
        }
    }
    wchar_t* dst = number.digits;
    if (*src != '0') {
        while (*src) *dst++ = *src++;
    }
    *dst++ = 0;
    NumberToDouble(&number, &v);
    return 0;
}
user541686
la source
4
Bonne explication +1. Ce code provient de shared-source-cli-2.0, non? C'est le seul que j'ai trouvé.
Soner Gönül
10
Je dois dire que c'est plutôt pathétique. Les chaînes qui sont mathématiquement égales (comme une avec un zéro à la fin, ou disons 2.1e-1 contre 0,21) devraient toujours donner des résultats identiques, et les chaînes qui sont mathématiquement ordonnées devraient donner des résultats cohérents avec l'ordre.
gnasher729
4
@MrLister: Pourquoi "2.1E-1 ne devrait-il pas être identique à 0.21 juste comme ça"?
user541686
9
@ gnasher729: Je suis plutôt d'accord sur "2.1e-1" et "0.21" ... mais une chaîne avec un zéro à la fin n'est pas exactement égale à un sans - dans le premier, le zéro est un chiffre significatif et ajoute précision.
cHao
4
@cHao: Euh ... cela ajoute de la précision, mais cela n'affecte que la façon dont vous décidez d'arrondir la réponse finale si les sigfigs comptent pour vous, pas la façon dont l'ordinateur doit calculer la réponse finale en premier lieu. Le travail de l'ordinateur est de tout calculer avec la plus grande précision, quelles que soient les précisions de mesure réelles des nombres; c'est le problème du programmeur s'il veut arrondir le résultat final.
user541686
107

Il me semble que ce n'est qu'un bug. Vos attentes sont tout à fait raisonnables. Je l'ai reproduit en utilisant .NET 4.5.1 (x64), en exécutant l'application console suivante qui utilise ma DoubleConverterclasse. DoubleConverter.ToExactStringmontre la valeur exacte représentée par un double:

using System;

class Test
{
    static void Main()
    {
        double d1 = 0.84551240822557006;
        string s = d1.ToString("r");
        double d2 = double.Parse(s);
        Console.WriteLine(s);
        Console.WriteLine(DoubleConverter.ToExactString(d1));
        Console.WriteLine(DoubleConverter.ToExactString(d2));
        Console.WriteLine(d1 == d2);
    }
}

Résultats dans .NET:

0.84551240822557
0.845512408225570055719799711368978023529052734375
0.84551240822556994469749724885332398116588592529296875
False

Résultats en Mono 3.3.0:

0.84551240822557006
0.845512408225570055719799711368978023529052734375
0.845512408225570055719799711368978023529052734375
True

Si vous spécifiez manuellement la chaîne de Mono (qui contient le «006» à la fin), .NET analysera la valeur d'origine. Il semble que le problème réside dans leToString("R") gestion plutôt que dans l'analyse.

Comme indiqué dans d'autres commentaires, il semble que cela soit spécifique à l'exécution sous le CLR x64. Si vous compilez et exécutez le code ci-dessus ciblant x86, tout va bien:

csc /platform:x86 Test.cs DoubleConverter.cs

... vous obtenez les mêmes résultats qu'avec Mono. Il serait intéressant de savoir si le bogue apparaît sous RyuJIT - je ne l'ai pas installé moi-même pour le moment. En particulier, je peux imaginer que cela pourrait être un bogue JIT, ou il est tout à fait possible qu'il existe des implémentations entièrement différentes des composants internes de double.ToStringbasés sur l'architecture.

Je vous suggère de déposer un bug sur http://connect.microsoft.com

Jon Skeet
la source
1
Alors Jon? Pour confirmer, est-ce un bogue dans le JITer, insérant le ToString()? Comme j'ai essayé de remplacer la valeur codée en dur par rand.NextDouble()et il n'y avait aucun problème.
Aron
1
Ouais, c'est définitivement dans la ToString("R")conversion. Essayez de ToString("G32")remarquer qu'il imprime la valeur correcte.
user541686
1
@Aron: Je ne peux pas dire s'il s'agit d'un bogue dans le JITter ou dans une implémentation spécifique à x64 du BCL. Je doute fort que ce soit aussi simple que l'inlining. Tester avec des valeurs aléatoires n'aide pas vraiment beaucoup, IMO ... Je ne suis pas sûr de ce que vous attendez de cela.
Jon Skeet
2
Ce qui se passe, je pense, c'est que le format "aller-retour" produit une valeur qui est 0,498ulp plus grande qu'elle ne devrait l'être, et l'analyse de la logique l'arrondit parfois de manière erronée à cette dernière petite fraction d'un ulp. Je ne suis pas sûr du code que je blâme le plus, car je pense qu'un format «aller-retour» devrait produire une valeur numérique qui se situe dans un quart ULP d'être numériquement correcte; la logique d'analyse qui donne une valeur à moins de 0,75 ulp de ce qui est spécifié est beaucoup plus facile que la logique qui doit produire un résultat à moins de 0,502 ulp de ce qui est spécifié.
supercat
1
Le site Web de Jon Skeet est en panne? Je trouve cela tellement improbable que je ... perds toute foi ici.
Patrick M
2

Récemment, j'essaye de résoudre ce problème . Comme indiqué dans le code , le double.ToString ("R") a la logique suivante:

  1. Essayez de convertir le double en chaîne avec une précision de 15.
  2. Reconvertissez la chaîne en double et comparez-la au double d'origine. S'ils sont identiques, nous renvoyons la chaîne convertie dont la précision est de 15.
  3. Sinon, convertissez le double en chaîne avec une précision de 17.

Dans ce cas, double.ToString ("R") a mal choisi le résultat avec une précision de 15 donc le bogue se produit. Il existe une solution de contournement officielle dans le document MSDN:

Dans certains cas, les valeurs Double formatées avec la chaîne de format numérique standard «R» ne parviennent pas à aller-retour si elles sont compilées à l'aide des commutateurs / platform: x64 ou / platform: anycpu et s'exécutent sur des systèmes 64 bits. Pour contourner ce problème, vous pouvez mettre en forme les valeurs Double à l'aide de la chaîne de format numérique standard «G17». L'exemple suivant utilise la chaîne de format «R» avec une valeur Double qui n'effectue pas un aller-retour et utilise également la chaîne de format «G17» pour effectuer un aller-retour réussi de la valeur d'origine.

Donc, à moins que ce problème ne soit résolu, vous devez utiliser double.ToString ("G17") pour aller-retour.

Mise à jour : il y a maintenant un problème spécifique pour suivre ce bogue.

Jim Ma
la source