c# - 浮点范围缩减

Question

我正在使用 Mono 在 C# 中实现一些 32 位浮点三角函数，希望使用 Mono.Simd。我目前只缺少固体范围减少。我现在有点卡住了，因为显然 Mono's SIMD extensions不包括 float 和整数之间的转换，这意味着我无法使用通常的方法进行舍入/截断。但是，我可以在整数和 float 之间按位转换。

这样的事情能做吗？如果需要，我可以上下缩放域，但理想情况下，范围缩小应该导致域为 [0, 2 pi] 或 [-pi, pi]。我有一种预感，如果域是 2 的幂，就可以用指数做一些 IEEE 魔法，但我真的不确定如何去做。

编辑:好吧，我试过弄乱这个 C 代码，感觉就像我在某事的边缘(它不起作用但小数部分总是正确的，至少在十进制/base10 中......)。核心原则似乎是获取你的域和输入指数之间的指数差异，并用移位的尾数和调整后的指数组成一个新的 float 。但它不适用于负数，我不知道如何处理2 的非幂(或任何小数 - 事实上，除 2 以外的任何东西都不起作用!)。

// here's another more correct attempt:
float fmodulus(float val, int domain)
{
    const int mantissaMask = 0x7FFFFF;
    const int exponentMask = 0x7F800000;

    int ival = *(int*)&val;

    int mantissa = ival & mantissaMask;
    int rawExponent = ival & exponentMask;
    int exponent = (rawExponent >> 23) - (129 - domain);
    // powers over one:
    int p = exponent;

    mantissa <<= p;
    rawExponent = exponent >> p;
    rawExponent += 127;
    rawExponent <<= 23;

    int newVal = rawExponent & exponentMask;
    newVal |= mantissa & mantissaMask;

    float ret = *(float*)&newVal;
    
    return ret;
}

float range_reduce(float value, int range )
{
    const int mantissaMask = 0x7FFFFF;
    const int exponentMask = 0x7F800000;

    int ival = *(int*)&value;
    // grab exponent:
    unsigned exponent = (ival & exponentMask) >> 23;
    // grab mantissa:
    unsigned mantissa = ival & mantissaMask;

    // remove bias, and see how much the exponent is over range/domain
    unsigned char erange = (unsigned char)(exponent - (125 + range));
    // check if sign bit is set - that is, the exponent is under our range
    if (erange & 0x80)
    {
        // don't do anything then.
        erange = 0;
    }

    // shift mantissa (and chop off bits) by the reduced amount
    int inewVal = (mantissa << (erange)) & mantissaMask;
    // add exponent, and subtract the amount we reduced the argument with
    inewVal |= ((exponent - erange) << 23) & exponentMask;

    // reinterpret
    float newValue = *(float*)&inewVal;
    return newValue;
    //return newValue - ((erange) & 0x1 ? 1.0f : 0.0f);
}

int main()
{
    float val = 2.687f;
    int ival = *(int*)&val;
    float correct = fmod(val, 2);
    float own = range_reduce(val, 2);

    getc(stdin);
}

编辑 2:

好吧，我真的想从 IEEE 二进制系统的角度来理解这一点。如果我们这样写模运算:

output = input % 2

[exponent] + [mantissa_bit_n_times_exponent]

3.5     = [2] + [1 + 0.5]                   ->[1] + [0.5]       = 1.5
4.5     = [4] + [0 + 0 + 0.5]               ->[0.5] + [0]       = 0.5
5.5     = [4] + [0 + 1 + 0.5]               ->[1] + [0.5]       = 1.5
2.5     = [2] + [0 + 0.5]                   ->[0.5] + [0]       = 0.5
2.25    = [2] + [0 + 0 + 0.25]              ->[0.25]            = 0.25
2.375   = [2] + [0 + 0 + 0.25 + 0.125]      ->[0.25] + [0.125]  = 0.375
13.5    = [8] + [4 + 0 + 1 + 0.5]           ->[1] + [0.5]       = 1.5
56.5    = [32] + [16 + 8 + 0 + 0 + 0 + 0.5] ->[0.5]             = 0.5

我们可以看到所有情况下的输出都是一个新数，没有原始指数，尾数移动了一个量(基于指数和第一个指数后尾数的第一个非零位- 尾数的位被忽略)到指数中。但我不确定这是否是正确的方法，它在纸面上效果很好。

编辑3:我卡在单声道版本 2.0.50727.1433

Answer 1

您可以将问题简化为取 float mod 1。为简化此操作，您可以使用位运算计算 float 的底数，然后使用浮点减法。以下是这些操作的(不安全)C# 代码:

// domain is assumed to be positive
// returns value in [0,domain)
public float fmodulus(float val, float domain)
{
    if (val < 0)
    {
        float negative = fmodulus(-val, domain);
        if (domain - negative == domain)
            return 0;
        else
            return domain-negative;
    }

    if (val < domain)
        return val; // this avoids losing accuracy

    return fmodOne(val / domain) * domain;
}

// assumes val >= 1, so val is positive and the exponent is at least 0 
unsafe public float fmodOne(float val)
{
    int iVal = *(int*)&val;
    int uncenteredExponent = iVal >> 23;
    int exponent = uncenteredExponent - 127; // 127 corresponds to 2^0 times the mantissa
    if (exponent >= 23) 
        return 0; // not enough precision to distinguish val from an integer

    int unneededBits = 23 - exponent; // between 0 and 23
    int iFloorVal = (iVal >> unneededBits) << unneededBits; // equivalent to using a mask to zero the bottom bits of the mantissa
    float floorVal = *(float*)&iFloorVal; // convert the bit pattern back to a float

    return val-floorVal;
}

例如，fmodulus(100.1f, 1) 是 0.09999847。 100.1f的位模式是

0 10000101 10010000011001100110011

floorVal(100f)的位模式是

0 10000101 10010000000000000000000

浮点减法给出接近 0.1f 的值:

0 01111011 10011001100110000000000

实际上，我很惊讶最后 8 位被清零了。我认为只有 0.1f 的最后 6 位应该被替换为 0。也许可以比依赖浮点减法做得更好。