Ограниченная точность с плавающей запятой и проблема генерирования бесконечно гармонического сигнала - PullRequest
0 голосов
/ 22 декабря 2018

Предположим, нам нужно создать очень длинный гармонический сигнал, в идеале бесконечно длинный.На первый взгляд решение кажется тривиальным:

Sample1:

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

К сожалению, это нерабочее решение.Для t >> dt из-за ограниченной точности с плавающей точкой будут получены неправильные значения.К счастью, можно вспомнить, что sin(2*PI* n + x) = sin(x), где n - произвольное целочисленное значение, поэтому при модификации примера несложно получить «бесконечный» аналог

Sample2:

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

Для одного физического моделирования мне нужно было получить бесконечный сигнал, который представляет собой сумму гармонических сигналов , например:

Sample3:

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

В этой форме код некорректно работает для больших t, по аналогичным причинам для Sample1.Вопрос - как получить «бесконечную» реализацию алгоритма Sample3?Я предполагаю, что решение должно выглядеть как Sample2.Очень важное примечание - вообще говоря, w [i] является произвольным, а не гармониками , то есть все частоты не кратны некоторой базовой частоте, поэтому я не могу найти общую tau.Использование типов с большей точностью (double, long double) не допускается.

Спасибо за ваш совет!

1 Ответ

0 голосов
/ 22 декабря 2018

Вы можете выбрать произвольный tau и сохранить напоминания фазы для каждого мода при вычитании его из 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);
       }
   }
}
...