试图编写代码来查找机器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

这取决于你所说的"精确等级".

浮点数具有"常规"(正常)值,但也有特殊的次正常数字.如果要查找不同的限制,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个值:

• 该类型(TYPE_EPSILON)中1和最小值之间的差值大于1 ,
• 给定类型(TYPE_MIN)中的最小正标准化值.这不包括次正规数,
• 给定类型中的最小正值(nextafter*(0... )).这包括次正规数,
• 最小数量大于1.这与TYPE相同_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 2

猜猜为什么你会得到相同的答案：

if( (1.0 + test) != 1.0)

这里 1.0 是一个双精度常量，因此它将您的浮点数提升为双精度，并以双精度执行加法。您可能想在此处声明一个临时浮点数来执行加法，或者使这些浮点数常量 ( 1.0fIIRC)。

您还可能陷入临时浮点精度过高的问题，并且可能需要强制它将中间值存储在内存中以降低到正确的精度。

这是重做范围搜索方法但以正确的类型计算测试的快速方法。不过，我得到的答案有点太大了。

#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;
}