Алгоритм умножения с фиксированной запятой - PullRequest
Новая
       19

Алгоритм умножения с фиксированной запятой

8 голосов
/ 04 марта 2011

Я пытаюсь изменить масштаб временной метки (только дробная часть секунд) с наносекунд (единицы 10 ^ -9 секунд) до нижней половины метки времени NTP (единицы 2 ^ -32 секунд). Фактически это означает умножение на 4,2949673. Но мне нужно делать это без математики с плавающей точкой и без использования целых чисел больше 32 бит (на самом деле, я пишу это для 8-битного микроконтроллера, поэтому даже 32-битная математика стоит дорого, особенно деления).

Я придумал несколько алгоритмов, которые работают достаточно хорошо, но у меня нет реальных основ численных методов, поэтому я был бы признателен за любые предложения о том, как их улучшить, или любые другие алгоритмы, которые будут точнее и / или быстрее.

Алгоритм 1 

uint32_t intts = (ns >> 16) * 281474 + (ns << 16) / 15259 + ns / 67078;

Первые две константы были выбраны для того, чтобы немного превышать, а не превышать, правильную цифру, и конечный коэффициент 67078 был определен эмпирически, чтобы исправить это. Производит результаты в пределах +/- 4 единиц NTP от правильного значения, которое +/- 1 нс - приемлемо, но остаточное значение изменяется с ns. Я думаю, я мог бы добавить еще один термин.

Алгоритм 2 

uint32_t ns2 = (2 * ns) + 1;
uint32_t intts = (ns2 << 1)
  + (ns2 >> 3) + (ns2 >> 6) + (ns2 >> 8) + (ns2 >> 9) + (ns2 >> 10)
  + (ns2 >> 16) + (ns2 >> 18) + (ns2 >> 19) + (ns2 >> 20) + (ns2 >> 21)
  + (ns2 >> 22) + (ns2 >> 24) + (ns2 >> 30) + 3;

Основано на двоичном расширении 4.2949673 (фактически на двоичном расширении 2.14748365, поскольку я начинаю с удвоения и добавления одного для выполнения округления). Возможно, быстрее, чем алгоритм 1 (я еще не получил тесты). +3 было определено эмпирически, чтобы исключить недолёт от усечения всех этих младших битов, но он не выполняет наилучшую возможную работу.

Ответы [ 2 ]

8 голосов
/ 04 марта 2011 
uint32_t convert(uint32_t x) {
    const uint32_t chi = 0x4b82;
    const uint32_t clo = 0xfa09;
    const uint32_t round = 0x9525;
    const uint32_t xhi = x >> 16;
    const uint32_t xlo = x & 0xffff;
    uint32_t lowTerm = xlo*clo;
    uint32_t crossTerms = xhi*clo + xlo*chi;
    uint32_t rounded = crossTerms + (lowTerm >> 16) + round >> 16;
    uint32_t highTerm = xhi*chi;
    return (x << 2) + highTerm + rounded;
}

Базовое умножение с фиксированной запятой, моделирование продукта 32x32 -> 64 с использованием четырех продуктов 16x16 -> 32. Константа round была выбрана для минимизации ошибки с помощью простого двоичного поиска. Это выражение хорошо для +/- 0,6 NTP во всем допустимом диапазоне.

Ведущий 4 в масштабном коэффициенте обрабатывается в смену. Компиляторы обычно могут генерировать довольно приличный код для такого типа вещей, но его часто можно упростить с помощью рукописной сборки, если вам нужно.

Если вам не нужна такая большая точность, вы можете избавиться от lowTerm и round и получить ответ, который хорош для +/- 1,15 NTP: 

uint32_t convert(uint32_t x) {
    const uint32_t chi = 0x4b82;
    const uint32_t clo = 0xfa09;
    const uint32_t xhi = x >> 16;
    const uint32_t xlo = x & 0xffff;
    uint32_t crossTerms = xhi*clo + xlo*chi;
    uint32_t highTerm = xhi*chi;
    return (x << 2) + highTerm + (crossTerms >> 16) + 1;
}
1 голос
/ 04 марта 2011

Я мог бы заявить об очевидном ... но вы гуглили интернет-библиотеки для математических библиотек с фиксированной точкой?Их много.Вот хороший пример с реализациями C ++ и x86 в архивах Flipcode:

http://www.flipcode.com/archives/Fixed_Point_Routines.shtml

...