计算范数/内积的C算法

html5 • • 问答

我需要检查 R^2 中的一个点是否位于半径 r 相对较大的圆（最多 10^5）中。显然，我通常只会将内积与 r^2 进行比较，但这是在嵌入式环境中，这不适用于足够大的 int32_t 值，因为正交会溢出类型（最多 32 位类型）。

可能的解决方案：

我可以从两个 32 位整数中手动拼凑出一个 64 位产品（可能我最终会这样做）。

我可以将所有内容除以 10（或任何值），然后进行通常的内积比较，但我失去了精度。

我可以尝试检查圆圈中内接的 n 边形内部，但这是很多计算、表格等，而且我仍然没有精确度。

是否有通常用于此类事情的算法？

回答

恐怕计算 64 位结果是最简单的解决方案。检查您的编译器是否可以为此生成有效的内联代码：

int check_distance(int x, int y, int r) {
    return (long long)x * x + (long long)y * y <= (long long)r * r;
}

如果生成的代码看起来很慢，您可以添加一个测试来检查是否需要 64 位操作。假设x，y和r是正数，这里是一个使用无符号算术和精确宽度类型的解决方案<stdint.h>：

int check_distance(uint32_t x, uint32_t y, uint32_t r) {
    if (x <= 46340 && y <= 46340 && r <= 0xffff) {
        /* 32-bit unsigned expression does not overflow */
        return x * x + y * y <= r * r;
    } else {
        return (uint64_t)x * x + (uint64_t)y * y <= (uint64_t)r * r;
    }
}

请注意常量 46340 是floor(sqrt(pow(2, 31)))： 如果x和y都大于此值，x*x + y*y则将超过 2 32。

这是一个具有更快测试的替代方法，但对于稍小的值，它将回退到 64 位操作：

int check_distance(uint32_t x, uint32_t y, uint32_t r) {
    if ((x | y | r) <= 0x7fff) {
        /* 32-bit unsigned expression does not overflow */
        return x * x + y * y <= r * r;
    } else {
        return (uint64_t)x * x + (uint64_t)y * y <= (uint64_t)r * r;
    }
}

那么如果你真的不想使用编译器的 64 位算法，你可以显式地编写计算。考虑 的范围x，y并r指定为<= 100000，将值右移 2 位保持x并y低于 46340：

int check_distance(uint32_t x, uint32_t y, uint32_t r) {
    if (x <= 46340 && y1 <= 46340 && r1 <= 0xffff) {
        /* 32-bit unsigned expression does not overflow */
        return x * x + y * y <= r * r;
    } else {
        /* shift all values right 2 bits to keep them below 46340 */
        uint32_t x0 = x & 3;
        uint32_t y0 = y & 3;
        uint32_t r0 = r & 3;
        uint32_t x1 = x >> 2;
        uint32_t y1 = y >> 2;
        uint32_t r1 = r >> 2;
        uint32_t x2_lo = x0 * (x0 + x1 * 8);
        uint32_t y2_lo = y0 * (y0 + y1 * 8);
        uint32_t d2_lo = x2_lo + y2_lo;
        uint32_t d2_hi = x1 * x1 + y1 * y1 + (d2_lo >> 4);
        uint32_t r2_lo = r0 * (r0 + r1 * 8);
        uint32_t r2_hi = r1 * r1 + (r2_lo >> 4);
        return d2_hi < r2_hi || (d2_hi == r2_hi && (d2_lo & 15) <= (r2_lo & 15));
    }
}

最后，将值移位 5 位允许最多 1000000 的数字：

int check_distance(uint32_t x, uint32_t y, uint32_t r) {
    if (x <= 46340 && y1 <= 46340 && r1 <= 0xffff) {
        /* 32-bit unsigned expression does not overflow */
        return x * x + y * y <= r * r;
    } else {
        /* shift all values right 5 bits to keep them below 46340 */
        uint32_t x0 = x & 31;
        uint32_t y0 = y & 31;
        uint32_t r0 = r & 31;
        uint32_t x1 = x >> 5;
        uint32_t y1 = y >> 5;
        uint32_t r1 = r >> 5;
        uint32_t x2_lo = x0 * (x0 + x1 * 64);
        uint32_t y2_lo = y0 * (y0 + y1 * 64);
        uint32_t d2_lo = x2_lo + y2_lo;
        uint32_t d2_hi = x1 * x1 + y1 * y1 + (d2_lo >> 10);
        uint32_t r2_lo = r0 * (r0 + r1 * 64);
        uint32_t r2_hi = r1 * r1 + (r2_lo >> 10);
        return d2_hi < r2_hi || (d2_hi == r2_hi && (d2_lo & 1023) <= (r2_lo & 1023));
    }
}

以上所有版本都会为指定的范围产生精确的结果。如果您不需要精确的结果，您可以移动这些值以将它们置于适当的范围内并执行 32 位计算：

int check_distance(uint32_t x, uint32_t y, uint32_t r) {
    while (x > 46340 || y > 46340 || r > 0xffff) {
        x >>= 1;
        y >>= 1;
        r >>= 1;
    }
    /* 32-bit unsigned expression no longer overflows */
    return x * x + y * y <= r * r;
}

以上是计算范数/内积的C算法的全部内容。
THE END
分享
二维码
< <上一篇
下一篇>>