試圖寫一個代碼尋找機器epsilon

Question

我想找出C中的各種浮點格式（即float，double和long double）的精度級別。以下是我目前使用的代碼：

#include <stdio.h> 
#define N 100000 

int main(void) 
{ 
    float max = 1.0, min = 0.0, test; 
    int i;        /* Counter for the conditional loop */ 

    for (i = 0; i < N; i++) { 
     test = (max + min)/2.0; 
     if((1.0 + test) != 1.0)   /* If too high, set max to test and try again */ 
    max = test; 
    if((1.0 + test) == 1.0)  /* If too low, set min to test and try again */ 
     min = test; 
    } 
    printf("The epsilon machine is %.50lf\n", max); 
    return 0; 
} 

這給出了大約2^-64的預期值。然而，當我將減速度變成雙打或「長雙打」時，我得到的答案相同，我應該得到一個較小的值，但我不知道。任何人有任何想法？

Answer 1

猜測，爲什麼你得到了相同的答案：

if((1.0 + test) != 1.0) 

這裏1.0是雙恆定的，所以它的推廣您的浮到一張雙人牀和執行另外一個雙。你可能想要在這裏聲明一個臨時的浮點數來執行加法，或者使這些浮點數字常量（1.0f IIRC）。

您可能也會陷入臨時浮點超精度問題中，可能需要強制它將中間體存儲在內存中以降低到正確的精度。

---

下面是重做範圍搜索方法的快速步驟，但計算正確類型的測試。不過，我得到的答案稍大。

#include <stdio.h> 
#define N 100000 
#define TYPE float 

int main(void) 
{ 
    TYPE max = 1.0, min = 0.0, test; 
    int i; 

    for (i = 0; i < N; i++) 
    { 
     TYPE one_plus_test; 

     test = (max + min)/((TYPE)2.0); 
     one_plus_test = ((TYPE)1.0) + test; 
     if (one_plus_test == ((TYPE)1.0)) 
     { 
     min = test; 
     } 
     else 
     { 
     max = test; 
     } 
    } 
    printf("The epsilon machine is %.50lf\n", max); 
    return 0; 
} 

Answer 2

我不知道你的算法應該如何工作。這一個（C++）提供了正確的答案：

#include <iostream> 

template<typename T> 
int epsilon() { 
    int pow = 0; 
    T eps = 1; 
    while (eps + 1 != 1) { 
     eps /= 2; 
     --pow; 
    } 
    return pow + 1; 
} 

int main() { 
    std::cout << "Epsilon for float: 2^" << epsilon<float>() << '\n'; 
    std::cout << "Epsilon for double: 2^" << epsilon<double>() << '\n'; 
} 

此計算的最小值，使得當添加到1，仍然是可區分從1

輸出：

Epsilon for float: 2^-23 
Epsilon for double: 2^-52 

Answer 3

甲這種代碼的問題是編譯器會將浮點變量加載到微處理器的浮點寄存器中。如果您的微處理器僅具有雙精度浮點寄存器，則精度將與float和double相同。

您需要找到一種方法來強制編譯器將浮點值存儲回每次兩次計算之間的內存（到正確類型的變量中）。這樣它不得不放棄寄存器的附加精度。但是今天的編譯器在優化你的代碼方面很聰明。所以這可能很難實現。

Answer 4

這取決於你的意思是「精確度」。

浮點數有「常規」（正常）值，但也有特殊的子正常數。如果要找出不同的限制，C標準具有預定義的常量：

#include <math.h> 
#include <stdio.h> 
#include <float.h> 

int main(void) 
{ 
    printf("%30s: %g\n", "FLT_EPSILON", FLT_EPSILON); 
    printf("%30s: %g\n", "FLT_MIN", FLT_MIN); 
    printf("%30s: %g\n", "nextafterf(0.0, 1.0)", nextafterf(0.0, 1.0)); 
    printf("%30s: %g\n", "nextafterf(1.0, 2.0)-1", (nextafterf(1.0, 2.0) - 1.0f)); 
    puts(""); 
    printf("%30s: %g\n", "DBL_EPSILON", DBL_EPSILON); 
    printf("%30s: %g\n", "DBL_MIN", DBL_MIN); 
    printf("%30s: %g\n", "nextafter(0.0, 1.0)", nextafter(0.0, 1.0)); 
    printf("%30s: %g\n", "nextafter(1.0, 2.0)-1", (nextafter(1.0, 2.0) - 1.0)); 
    puts(""); 
    printf("%30s: %Lg\n", "LDBL_EPSILON", LDBL_EPSILON); 
    printf("%30s: %Lg\n", "LDBL_MIN", LDBL_MIN); 
    printf("%30s: %Lg\n", "nextafterl(0.0, 1.0)", nextafterl(0.0, 1.0)); 
    printf("%30s: %Lg\n", "nextafterl(1.0, 2.0)-1", (nextafterl(1.0, 2.0) - 1.0)); 
    return 0; 
} 

上述程序打印紙4點的值對於每種類型：

• 和至少值1之間的差大於1在該類型（TYPE _EPSILON），在給定類型（TYPE _MIN）
• 最小正歸一化值。這不包括subnormal numbers,
• 給定類型的最小正值（nextafter * (0 ... )）。這包括次正態數字，
• 大於1的最小數字。這與類型相同_EPSILON，但以不同方式計算。

根據「精度」的含義，以上任何一項或任何一項都不會對您有用。

下面是上述程序的我的電腦上輸出：

   FLT_EPSILON: 1.19209e-07 
        FLT_MIN: 1.17549e-38 
     nextafterf(0.0, 1.0): 1.4013e-45 
    nextafterf(1.0, 2.0)-1: 1.19209e-07 

       DBL_EPSILON: 2.22045e-16 
        DBL_MIN: 2.22507e-308 
     nextafter(0.0, 1.0): 4.94066e-324 
    nextafter(1.0, 2.0)-1: 2.22045e-16 

       LDBL_EPSILON: 1.0842e-19 
        LDBL_MIN: 3.3621e-4932 
     nextafterl(0.0, 1.0): 3.6452e-4951 
    nextafterl(1.0, 2.0)-1: 1.0842e-19 

Answer 5

我想補充一點，你可以通過使用long double得到一個浮點運算精度最高。

要應用到@ RUP的解決方案，只是改變TYPE到long double和printf聲明：

printf("The epsilon machine is %.50Lf\n", max); 

這是小量使用float我的機器上：

0.00000005960465188081798260100185871124267578125000 

而且使用long double：

0.00000000000000000005421010862427522170625011179761 

區別非常顯着。

Answer 6

IEEE 754浮點格式具有當重新解釋爲具有相同寬度的二進制補碼整數時，它們在正值上單調遞增且在負值上單調遞減（參見32位浮點數的二進制表示）的屬性。他們也有財產0 < | f（x）| <∞，| f（x + 1） - f（x）| ≥| f（x）-f（x-1）| （其中f（x）是x的上述整數重新解釋）。在允許類型雙擊並且總是使用IEEE 754-1985的語言中，我們可以利用這個來計算一個固定時間內的機器epsilon。例如，在C：

typedef union { 
    long long i64; 
    double d64; 
} dbl_64; 

double machine_eps (double value) 
{ 
    dbl_64 s; 
    s.d64 = value; 
    s.i64++; 
    return s.d64 - value; 
} 

從https://en.wikipedia.org/wiki/Machine_epsilon