首页 > 解决方案 > 有限浮点精度和无限谐波信号生成问题

问题描述

假设我们需要生成一个很长的谐波信号,理想情况下是无限长的。乍一看,解决方案似乎微不足道:

样品1:

   float t = 0;
   while (runned)
   {
      float v = sinf(w * t);
      t += dt;
   }

不幸的是,这是一个无效的解决方案。因为t >> dt由于浮点精度有限,将获得不正确的值。幸运的是,我们可以记住sin(2*PI* n + x) = sin(x)其中 n - 任意整数值,因此修改示例并不难获得“无限”模拟

样品2:

   float t = 0;
   float tau = 2 * M_PI / w;
   while (runned)
   {
      float v = sinf(w * t);
      t += dt;
      if (t > tau) t -= tau;
   }

对于一个物理模拟,我需要得到一个无限信号,它是谐波信号的总和,如下所示:

样品3:

   float getSignal(float x)
   {
      float ret = 0;
      for (int i = 0; i < modNum; i++)
         ret += sin(w[i] * x);
      return ret;
   }

   float t = 0;
   while (runned)
   {
      float v = getSignal(t);
      t += dt;
   }

在这种形式中,代码对于 large 无法正常工作t,原因与 Sample1 类似。问题是 - 如何获得 Sample3 算法的“无限”实现?我假设解决方案应该看起来像 Sample2。一个非常重要的注意事项 - 一般来说,w[i] 是任意的而不是谐波,也就是说,所有频率都不是某些基频的倍数,所以我找不到 common tau。不允许使用精度更高的类型(双精度、长双精度)。

谢谢你的建议!

标签: c++floating-pointtrigonometry

解决方案


您可以选择任意一个tau并在从中减去每个 mod 时存储阶段提醒t(正如@Damien 在评论中建议的那样)。

此外,将时间表示为t = dt * it整数it可以提高数值稳定性(我认为)。

也许是这样的:

int ndt = 1000;       // accumulate phase every 1000 steps for example
float tau = dt * ndt;

std::vector<float> phases(modNum, 0.0f);

int it = 0;
float t = 0.0f;
while (runned)
{
   t = dt * it;

   float v = 0.0f;
   for (int i = 0; i < modNum; i++)
   {
       v += sinf(w[i] * t + phases[i]);
   }

   if (++it >= ndt)
   {
       it = 0;
       for (int i = 0; i < modNum; ++i)
       {
           phases[i] = fmod(w[i] * tau + phases[i], 2 * M_PI);
       }
   }
}

推荐阅读