sin(pi)没有返回精确的零? [英] sin(pi) not returning exact zero?

问题描述

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

int main(int argc, char* argv[])
{
    double pi = acos(-1);

    for(int i = 0; i<65;i++)
    {
        double resultD = sin( pi * (2.0 *i / 128.0 ));
        float resultF  = sin( pi * (2.0f*i / 128.0f));

        printf("sine[%d]\n %.25f\n %.25f \n", i,resultD,resultF);
    }
}

在上面的代码中，我希望最后一个值干净为零.但是，它类似于0.0000000000000001224606354. 另一方面，在i=0处的值是干净的0.000000000000000，在i=32处的值是干净的1.000000000000000.

In the above code, I would expect to get last value to be clean zero. However, it is something like 0.0000000000000001224606354. On the other hand, the values at i=0 is clean 0.000000000000000, and the value at i=32 is clean 1.000000000000000.

我正在将浮点数移植到定点代码，这种微小的差异会逐步升级.
有没有办法强制i=64处的值精确为零?

I am porting floating-point to fixed-point code and this small difference escalates.
Is there a way to force the value at i=64 to be exact zero?

推荐答案

sin(pi)不会返回0，因为pi不是π.

sin(pi) does not return 0 because pi is not π.

π是不合理的数字.即使具有更高的精度和最佳的acos()/sin()例程，也只能近似π的值.从acos(-1)返回的值是 machine-pi ，它是最接近π的double.

π is an irrational number. Even with extended precision and the best acos()/sin() routines, the value of π can only be approximated. The value returned from acos(-1) is machine-pi, the closest double to π.

printf("machine pi %.20f...\n", pi);
// machine pi 3.14159265358979311600...
// math pi    3.14159265358979323846...
//                             ^^^^^...

当代码计算pi * (2.0*64/128.0)时，可以预期乘积是精确的，没有舍入误差，但是pi的值此处是 machine-pi 而不是π. machine-pi 的正弦确实约为1.22 ... e-16

When code calculates pi * (2.0*64/128.0), the product can be expected to be exact with no rounding error, yet the value of pi here is machine-pi and not π. The sine of machine-pi is indeed about 1.22...e-16

我正在将浮点代码移植到定点代码，这种微小的差异会逐步升级.

I am porting floating-point to fixed-point code and this small difference escalates.

对于三角函数，如果原始参数是度(或在OP的情况下为圆的1/128)或某个精确周期，则与无理2 *π弧度不同:请首先执行范围缩小在原始单位中，则可能随后使用trig身份. 最后转换为弧度并应用trig函数.

For trigonometric functions, if the original argument is degrees (or in OP's case 1/128th of a circle) or some exact period, unlike the irrational 2*π radians: perform the range reductions first in original units, then perhaps then use trig identities. Lastly convert to radians and apply the trig function.

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

#ifndef M_PI
#define M_PI 3.1415926535897932384626433832795
#endif

// Binary angle measurement
#define BAM 128

static double bam128tor(int bam) {
  return M_PI * 2 / BAM * bam;
}

double sin_BAM128(int bam) {
  // mod by the integral number of units in a circle
  bam %= BAM;
  if (bam < 0) bam += BAM;

  unsigned quadrant = 0;
  if (bam >= BAM / 2) {
    bam -= BAM / 2;
    quadrant += 2;
  }
  if (bam >= BAM / 4) {
    bam -= BAM / 4;
    quadrant += 1;
  }
  if (bam >= BAM / 8) {
    bam = -(BAM / 4 - bam);
    quadrant = (quadrant + 1) % 4;
  }
  printf("%2u %4d ", quadrant, bam);
  switch (quadrant % 4) {
    case 0:
      return sin(bam128tor(bam));
    case 1:
      return cos(bam128tor(bam));
    case 2:
      return sin(bam128tor(-bam)) * 1.0;
    default:
      return -cos(bam128tor(bam));
  }
}

int main(void) {
  double pi = acos(-1);
  for (int i = 0; i <= BAM; i += 4) {
    double resultD = sin(pi * (2.0 * i / 128.0));
    //float resultF = sin(pi * (2.0f * i / 128.0f));
    printf("sine[%3d]\t% .20f\t% .20f\n", i, resultD, sin_BAM128(i));
  }
}

输出

 0    0 sine[  0]    0.00000000000000000000  0.00000000000000000000
 0    4 sine[  4]    0.19509032201612824808  0.19509032201612824808
 0    8 sine[  8]    0.38268343236508978178  0.38268343236508978178
 0   12 sine[ 12]    0.55557023301960217765  0.55557023301960217765
 1  -16 sine[ 16]    0.70710678118654746172  0.70710678118654757274
 1  -12 sine[ 20]    0.83146961230254523567  0.83146961230254523567
 1   -8 sine[ 24]    0.92387953251128673848  0.92387953251128673848
 1   -4 sine[ 28]    0.98078528040323043058  0.98078528040323043058
 1    0 sine[ 32]    1.00000000000000000000  1.00000000000000000000
 1    4 sine[ 36]    0.98078528040323043058  0.98078528040323043058
 1    8 sine[ 40]    0.92387953251128673848  0.92387953251128673848
 1   12 sine[ 44]    0.83146961230254545772  0.83146961230254523567
 2  -16 sine[ 48]    0.70710678118654757274  0.70710678118654746172
 2  -12 sine[ 52]    0.55557023301960206663  0.55557023301960217765
 2   -8 sine[ 56]    0.38268343236508989280  0.38268343236508978178
 2   -4 sine[ 60]    0.19509032201612860891  0.19509032201612824808
 2    0 sine[ 64]    0.00000000000000012246  0.00000000000000000000 exact
 2    4 sine[ 68]   -0.19509032201612835911 -0.19509032201612824808
 2    8 sine[ 72]   -0.38268343236508967076 -0.38268343236508978178
 2   12 sine[ 76]   -0.55557023301960195560 -0.55557023301960217765
 3  -16 sine[ 80]   -0.70710678118654746172 -0.70710678118654757274
 3  -12 sine[ 84]   -0.83146961230254523567 -0.83146961230254523567
 3   -8 sine[ 88]   -0.92387953251128651644 -0.92387953251128673848
 3   -4 sine[ 92]   -0.98078528040323031956 -0.98078528040323043058
 3    0 sine[ 96]   -1.00000000000000000000 -1.00000000000000000000
 3    4 sine[100]   -0.98078528040323043058 -0.98078528040323043058
 3    8 sine[104]   -0.92387953251128662746 -0.92387953251128673848
 3   12 sine[108]   -0.83146961230254545772 -0.83146961230254523567
 0  -16 sine[112]   -0.70710678118654768376 -0.70710678118654746172
 0  -12 sine[116]   -0.55557023301960217765 -0.55557023301960217765
 0   -8 sine[120]   -0.38268343236509039240 -0.38268343236508978178
 0   -4 sine[124]   -0.19509032201612871993 -0.19509032201612824808
 0    0 sine[128]   -0.00000000000000024493  0.00000000000000000000 exact

这篇关于sin(pi)没有返回精确的零?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！