定点乘法答案

【问题标题】：Fixed point multiplication定点乘法
【发布时间】：2012-12-04 14:39:36
【问题描述】：

我需要根据一个非常数因子将一个值从一个单位转换为另一个单位。输入值范围是0到1073676289，取值范围是0到1155625。转换可以这样描述：

output = input * (range / 1073676289)

我自己最初的定点实现感觉有点笨拙：

// Input values (examples)
unsigned int input = 536838144;  // min 0, max 1073676289
unsigned int range = 1155625;    // min 0, max 1155625

// Conversion
unsigned int tmp = (input >> 16) * ((range) >> 3u);
unsigned int output = (tmp / ((1073676289) >> 16u)) << 3u;

我的代码可以改进得更简单或更准确吗？

【问题讨论】：

您的预期输出范围是多少？（是不是 (ab)b

标签： c fixed-point

【解决方案1】：

这将为您提供没有浮点值的最佳精度，结果将四舍五入到最接近的整数值：

output = (input * (long long) range + 536838144) / 1073676289;

【讨论】：

@banuj: 536838144 是 1073676289 的一半；并导致最终结果四舍五入到最接近的整数而不是截断。例如，1073676288/1073676289 = 0 但(1073676288+536838144)/1073676289 = 1。
@banuj：正如布伦丹所说，它是为了四舍五入。更准确地说，它是地板（分母/2）。整数的另外两个“边缘”示例：(536838144 + 536838144) / 1073676289 = 0 但(536838145 + 536838144) / 1073676289 = 1。

【解决方案2】：

问题是input * range 会溢出一个 32 位整数。通过使用 64 位整数来解决此问题。

uint64_least_t tmp;
tmp = input;
tmp = tmp * range;
tmp = tmp / 1073676289ul;
output = temp;

【讨论】：

【解决方案3】：

快速访问 Google 让我注意到了 http://sourceforge.net/projects/fixedptc/

这是一个头文件中的 c 库，用于管理 32 位或 64 位整数中的定点数学。

用下面的代码做一点实验：

#include <stdio.h>
#include <stdint.h>

#define FIXEDPT_BITS        64

#include "fixedptc.h"

int main(int argc, char ** argv)
{
    unsigned int input = 536838144;  // min 0, max 1073676289
    unsigned int range = 1155625;    // min 0, max 1155625

    // Conversion
    unsigned int tmp = (input >> 16) * ((range) >> 3u);
    unsigned int output = (tmp / ((1073676289) >> 16u)) << 3u;

    double output2 = (double)input * ((double)range / 1073676289.0);

    uint32_t output3 = fixedpt_toint(fixedpt_xmul(fixedpt_fromint(input), fixedpt_xdiv(fixedpt_fromint(range), fixedpt_fromint(1073676289))));

    printf("baseline = %g, better = %d, library = %d\n", output2, output, output3);

    return 0;
}

得到以下结果：

baseline = 577812, better = 577776, library = 577812

显示出比您的代码更好的精度（匹配浮点数）。在引擎盖下它并没有做任何非常复杂的事情（并且在 32 位中根本不起作用）

/* Multiplies two fixedpt numbers, returns the result. */
static inline fixedpt
fixedpt_mul(fixedpt A, fixedpt B)
{
    return (((fixedptd)A * (fixedptd)B) >> FIXEDPT_FBITS);
}


/* Divides two fixedpt numbers, returns the result. */
static inline fixedpt
fixedpt_div(fixedpt A, fixedpt B)
{
    return (((fixedptd)A << FIXEDPT_FBITS) / (fixedptd)B);
}

但它确实表明您可以获得所需的精度。你只需要 64 位就可以做到这一点

【讨论】：

【解决方案4】：

你不会比output = input * (range / 1073676289)更简单

如以下 cmets 所述，如果您被限制为整数运算，那么对于 range < 1073676289：range / 1073676289 == 0，您最好使用：

output = range < 1073676289 ? 0 : input

如果这不是您想要的，而您实际上想要精确，那么

output = (input * range) / 1073676289

将是要走的路。

如果你需要做很多这些，那么我建议你使用double 并让你的编译器矢量化你的操作。精度也可以。

【讨论】：

因为range 总是小于 1073676289 并且您使用的是整数除法； output = input * (range / 1073676289) 与 output = input * 0 相同。
如果我不清楚，我很抱歉，但我不能使用浮点值（也不能使用双精度值）。