C 乘法或加法浮点结果 NaN

Question

我在程序中使用 libresample。一段时间后(大约 50 分钟)，它在一个工作站的 lib 函数 lrsFilterUD() 中崩溃。

float lrsFilterUD(float Imp[],  /* impulse response */
              float ImpD[], /* impulse response deltas */
              UWORD Nwing,  /* len of one wing of filter */
              BOOL Interp,  /* Interpolate coefs using deltas? */
              float *Xp,    /* Current sample */
              double Ph,    /* Phase */
              int Inc,    /* increment (1 for right wing or -1 for left) */
              double dhb)
{
   float a;
   float *Hp, *Hdp, *End;
   float v, t;
   double Ho;

   v = 0.0; /* The output value */
   Ho = Ph*dhb;
   End = &Imp[Nwing];
   if (Inc == 1)        /* If doing right wing...              */
   {                      /* ...drop extra coeff, so when Ph is  */
      End--;            /*    0.5, we don't do too many mult's */
      if (Ph == 0)      /* If the phase is zero...           */
         Ho += dhb;     /* ...then we've already skipped the */
   }                         /*    first sample, so we must also  */
                        /*    skip ahead in Imp[] and ImpD[] */

   if (Interp)
      while ((Hp = &Imp[(int)Ho]) < End) {
         t = *Hp;       /* Get IR sample */
         Hdp = &ImpD[(int)Ho];  /* get interp bits from diff table*/
         a = Ho - floor(Ho);      /* a is logically between 0 and 1 */
         t += (*Hdp)*a; /* t is now interp'd filter coeff */
         t *= *Xp;      /* Mult coeff by input sample */
         v += t;            /* The filter output */
         Ho += dhb;     /* IR step */
         Xp += Inc;     /* Input signal step. NO CHECK ON BOUNDS */
      }
   else 
      while ((Hp = &Imp[(int)Ho]) < End) {
         dprintf("while begin: Hp = %p, *Hp = %a, (int)Ho = %d, Imp[(int)Ho] = %a, &Imp[(int)Ho] = %p", Hp, *Hp, (int)Ho, Imp[(int)Ho], &Imp[(int)Ho]);
         t = *Hp;       /* Get IR sample */
         dprintf("before t = %a, *Xp = %a, Xp = %p", t, *Xp, Xp);
         t *= *Xp;      /* Mult coeff by input sample */
         dprintf("after2 t = %a, v = %a", t, v);
         v += t;            /* The filter output */
         dprintf("v = %a", v);
         Ho += dhb;     /* IR step */
         Xp += Inc;     /* Input signal step. NO CHECK ON BOUNDS */
      }

   return v;
}

我记录了 t、*Xp、Xp 乘法前后的值:

while begin: Hp = 0xaf5daa8, *Hp = -0.009034, (int)Ho = 16384, Imp[(int)Ho] = -0.009034, &Imp[(int)Ho] = 0xaf5daa8
before multiplication t = -0.009034, *Xp = 0.000000, Xp = 0xaebe9b8
after multiplication t = nan

这段代码运行了很多次，崩溃前有相同的t和Xp值:

before multiplication t = -0.009034, *Xp = 0.000000, Xp = 0xaebe9c8
after multiplication t = -0.000000, v = 282.423676

或者另一种情况:

before addition t = -460.799988, v = 0.000000
after addition v = nan

什么可能导致 nan？这是在 Linux 上使用 gcc 4.1.2 编译的。

更新:将变量打印为 %a。结果:

//t = 0x1.2806bap+2
//Hp = 0xb3bb870
t = *Hp;
//t = nan

更新二:如果用icpc编译代码就没有这个问题。那么是否存在特定于编译器的问题？

Answer 1

显然，-0.009034•0.000000 不应产生 NaN。因此，问题中呈现的代码和数据不是实际计算的准确表示，或者计算实现有缺陷。

如果我们假设硬件和基本计算实现没有缺陷，那么调查的一些可能性包括:

• t 和 *Xp 的记录未能立即记录 t 和 *Xp 的正确值在乘法之前或乘法之后的 t 的正确值。
• t 或*Xp 的值显示不正确。例如，用于显示 *Xp 的格式显示“0.000000”，即使 *Xp 具有其他值，例如 NaN。
• Xp 指向不合适的地方，导致 *Xp 不可靠(例如，被外部操作更改)。
• 显示的代码不准确。例如，它来自已更改的旧源代码，或者它是新源代码但正在执行以前编译的代码。

注意:当使用浮点对象进行调试时，您应该不使用“%f”等格式打印，尤其是不要使用默认的位数。您应该使用“%a”打印，它使用十六进制表示打印浮点值的精确值。在许多情况下，您也可以使用“%.99g”，前提是您的 C 实现能够很好​​地将浮点值转换为十进制。