Каков первый дубль, который отклоняется от соответствующего значения long на дельту?

Question

Я хочу знать первое двойное значение от 0d вверх, которое отклоняется на длину "того же значения" на некоторую дельту, скажем, 1e-8.Но здесь у меня ничего не получается.Я пытаюсь сделать это на C, хотя обычно использую управляемые языки, на всякий случай.Пожалуйста, помогите.

#include <stdio.h>
#include <limits.h>
#define DELTA 1e-8

int main() {
    double d = 0; // checked, the literal is fine
    long i;
    for (i = 0L; i < LONG_MAX; i++) {
         d=i; // gcc does the cast right, i checked
         if (d-i > DELTA || d-i < -DELTA) {
              printf("%f", d);
              break;
         }
    }
}

Я предполагаю, что проблема в том, что d-i преобразует i в double и, следовательно, d == i, и тогда разница всегда равна 0.Как еще я могу определить это правильно - я бы предпочел забавное приведение C вместо сравнения строк, что заняло бы целую вечность.

ОТВЕТ:все именно так, как мы ожидали.2 ^ 53 + 1 = 9007199254740993 - это первая точка различия в соответствии со стандартными инструментами C / UNIX / POSIX.Большое спасибо Паксу за его программу.И я думаю, математика снова побеждает.

Answer 1

Удвоения в IEE754 имеют точность 52 бита, что означает, что они могут хранить числа с точностью до (как минимум) 2⁵¹.

Если ваши длинные значения 32-разрядные, они будут иметь только (положительный) диапазон от 0 до 2³¹ таким образом, не существует 32-битной длины, которая не могла бы быть представлена точно как double .Для 64-разрядной версии это будет (примерно) 2⁵² так что я бы начал примерно с этого, а не с нуля.

Вы можете использовать следующую программу, чтобы определить, где начинают возникать сбои.В более ранней версии я полагался на тот факт, что последняя цифра в числе, которое непрерывно удваивается, следует за последовательностью {2,4,8,6}.Однако в конце концов я решил использовать известный надежный инструмент (bc) для проверки всего числа, а не только последней цифры.

Имейте в виду , что это мочь быть затронутым действиями sprintf() вместо реальной точности удвоений (лично я так не думаю, поскольку у него не было проблем с определенными числами до 2¹⁴³).

Это и есть программа:

#include <stdio.h>
#include <string.h>

int main() {
    FILE *fin;
    double d = 1.0; // 2^n-1 to avoid exact powers of 2.
    int i = 1;
    char ds[1000];
    char tst[1000];

    // Loop forever, rely on break to finish.
    while (1) {
        // Get C version of the double.
        sprintf (ds, "%.0f", d);

        // Get bc version of the double.
        sprintf (tst, "echo '2^%d - 1' | bc >tmpfile", i);
        system(tst);
        fin = fopen ("tmpfile", "r");
        fgets (tst, sizeof (tst), fin);
        fclose (fin);
        tst[strlen (tst) - 1] = '\0';

        // Check them.
        if (strcmp (ds, tst) != 0) {
            printf( "2^%d - 1 <-- bc failure\n", i);
            printf( "   got       [%s]\n", ds);
            printf( "   expected  [%s]\n", tst);
            break;
        }

        // Output for status then move to next.
        printf( "2^%d - 1 = %s\n", i, ds);
        d = (d + 1) * 2 - 1;  // Again, 2^n - 1.
        i++;
    }
}

Это продолжается до тех пор, пока:

2^51 - 1 = 2251799813685247
2^52 - 1 = 4503599627370495
2^53 - 1 = 9007199254740991
2^54 - 1 <-- bc failure
   got       [18014398509481984]
   expected  [18014398509481983]

примерно там, где я ожидал, что это потерпит неудачу.

В качестве отступления, я изначально использовал числа вида 2ⁿ но это подтолкнуло меня к:

2^136 = 87112285931760246646623899502532662132736
2^137 = 174224571863520493293247799005065324265472
2^138 = 348449143727040986586495598010130648530944
2^139 = 696898287454081973172991196020261297061888
2^140 = 1393796574908163946345982392040522594123776
2^141 = 2787593149816327892691964784081045188247552
2^142 = 5575186299632655785383929568162090376495104
2^143 <-- bc failure
   got       [11150372599265311570767859136324180752990210]
   expected  [11150372599265311570767859136324180752990208]

с размером двойного файла, равным 8 байтам (проверено с помощью sizeof).Оказалось , что эти числа были двоичной формы "1000..." который гораздо дольше может быть представлен двойниками.Вот тогда-то я и перешел на использование 2ⁿ-1, чтобы получить лучший битовый шаблон:все до единого биты.

Answer 2

Первый лонг, который будет "неправильным" при приведении к double, не будет отключен на 1e-8, он будет отключен на 1.До тех пор, пока double может соответствовать long по своему значению, он будет представлять его точно.

Я забыл точно, сколько битов имеет double для сравнения точности со смещением, но это подскажет вам максимальный размер, который он может представлять.Первое значение long, которое будет неправильным, должно иметь двоичную форму 10000 ..., так что вы можете найти его намного быстрее, начав с 1 и сдвинув влево.

В Википедии указано 52 бита в значимом значении, не считая неявного начального 1.Это должно означать, что первое значение long, которое будет приведено к другому значению, равно 2^ 53.

Answer 3

Хотя я не решаюсь упоминать Fortran 95 и его преемников в этом обсуждении, я упомяну, что Fortran начиная со стандарта 1990 года предлагает встроенную функцию интервала, которая сообщает вам, в чем разница между представимыми реалами относительно данного РЕАЛЬНОГО.Вы могли бы выполнить бинарный поиск по этому вопросу, останавливаясь, когда ИНТЕРВАЛ (X) > ДЕЛЬТА.Для компиляторов, которые используют ту же модель с плавающей запятой, что и та, которая вас интересует (вероятно, это стандарт IEEE754), вы должны получить те же результаты.

Answer 4

С другой стороны, я думал, что удвоения могут точно представлять все целые числа (в пределах их границ).

Если это не так, то вам захочется привести i и d к чему-то с БОЛЬШЕЙ точностью, чем любой из них.Возможно, длинный дубль сработает.