floating-point - 反复规范化 IEEE 浮点向量会使其发生变异吗？

Question

如果我采用(非零)浮点向量((x, y, z) 向量)，并将其归一化为单位长度，第二次归一化保证返回相同结果？

Answer 1

我不知道文献中有相关结果。快速测试表明，对一个 3D 向量进行两次归一化经常会导致归一化向量与重新归一化向量之间的微小差异，即使在小心地准确执行归一化时，例如，通过更高精度的算术执行也是如此。我使用下面的 ISO-C99 程序进行此快速测试，使用我的编译器针对 x64 平台的“严格”浮点设置 (icl/fp:strict) 对其进行编译。

#include <stdio.h>
#include <stdlib.h>
#include <stdint.h>
#include <string.h>
#include <math.h>

#define USE_HYPOT  (1)

// Fixes via: Greg Rose, KISS: A Bit Too Simple. http://eprint.iacr.org/2011/007
static uint32_t z=362436069,w=521288629,jsr=362436069,jcong=123456789;
#define znew (z=36969*(z&0xffff)+(z>>16))
#define wnew (w=18000*(w&0xffff)+(w>>16))
#define MWC  ((znew<<16)+wnew)
#define SHR3 (jsr^=(jsr<<13),jsr^=(jsr>>17),jsr^=(jsr<<5)) /* 2^32-1 */
#define CONG (jcong=69069*jcong+13579)                     /* 2^32 */
#define KISS ((MWC^CONG)+SHR3)

float uint32_as_float (uint32_t a)
{
    float r;
    memcpy (&r, &a, sizeof(r));
    return r;
}

void normalize_3d (float *a, float *b, float *c)
{
#if USE_HYPOT
    double l = hypot (hypot ((double)a[0], (double)b[0]), (double)c[0]);
#else // USE_HYPOT
    double l = sqrt ((double)a[0]*a[0] + (double)b[0]*b[0] + (double)c[0]*c[0]);
#endif // USE_HYPOT
    *a = (float)((double)a[0] / l);
    *b = (float)((double)b[0] / l);
    *c = (float)((double)c[0] / l);
}

int main (void)
{
    float a, aa, aaa, b, bb, bbb, c, cc, ccc;
    do {
        /* generate random vector */
        do {
            a = uint32_as_float (KISS & ~0x80000000u);
        } while (isnanf (a) || (a > 0x1.0p126f) || (a < 0x1.0p-126f));
        do {
            b = uint32_as_float (KISS & ~0x80000000u);
        } while (isnanf (b) || (b > 0x1.0p126f) || (b < 0x1.0p-126f));
        do {
            c = uint32_as_float (KISS & ~0x80000000u);
        } while (isnanf (c) || (c > 0x1.0p126f) || (c < 0x1.0p-126f));

        /* normalize vector once */
        aa = a; bb = b; cc = c;
        normalize_3d (&aa, &bb, &cc);

        /* re-normalize normalized vector */
        aaa = aa; bbb = bb; ccc = cc;
        normalize_3d (&aaa, &bbb, &ccc);

        /* check whether normalized vector is equal to re-normalized one */
        if ((aa != aaa) || (bb != bbb) || (cc != ccc)) {
            printf ("norm = (%15.6a, %15.6a, %15.6a)  re-norm = (%15.6a, %15.6a, %15.6a)\n", aa, bb, cc, aaa, bbb, ccc);
        }
    } while (1);
    return EXIT_SUCCESS;
}