Is it possible to write Quake's fast InvSqrt() function in C#?

Question

This is just to satisfy my own curiosity.

Is there an implementation of this:

float InvSqrt (float x)
{
   float xhalf = 0.5f*x;
   int i = *(int*)&x;
   i = 0x5f3759df - (i>>1);
   x = *(float*)&i;
   x = x*(1.5f - xhalf*x*x);
   return x;
}

in C#? If it exists, post the code.

I guess I should have mentioned I was looking for a "safe" implementation... Either way, the BitConverter code solves the problem. The union idea is interesting. I'll test it and post my results.

Edit: As expected, the unsafe method is the quickest, followed by using a union (inside the function), followed by the BitConverter. The functions were executed 10000000 times, and the I used the System.Diagnostics.Stopwatch class for timing. The results of the calculations are show in brackets.

Input: 79.67
BitConverter Method: 00:00:01.2809018 (0.1120187)
Union Method: 00:00:00.6838758 (0.1120187)
Unsafe Method: 00:00:00.3376401 (0.1120187)

For completeness, I tested the built-in Math.Pow method, and the "naive" method (1/Sqrt(x)).

Math.Pow(x, -0.5): 00:00:01.7133228 (0.112034710535584)
1 / Math.Sqrt(x): 00:00:00.3757084 (0.1120347)

The difference between 1 / Math.Sqrt() is so small that I don't think one needs to resort to the Unsafe Fast InvSqrt() method in C# (or any other unsafe method). Unless one really needs to squeeze out that last bit of juice from the CPU... 1/Math.Sqrt() is also much more accurate.

Downcast

For completeness, you should run a test using "1/math.Sqrt()" also. – Pyrrhonism

I did, and another scenario. I'll update the benchmarks as soon as I can verify the results. – Downcast

I think you should have run your tests with a larger set, so that they would take a few seconds to complete. – Vibrate

Answer 1

You should be able to use the StructLayout and FieldOffset attributes to fake a union for plain old data like floats and ints.

[StructLayout(LayoutKind.Explicit, Size=4)]
private struct IntFloat {
    [FieldOffset(0)]
    public float floatValue;

    [FieldOffset(0)]
    public int intValue;

    // redundant assignment to avoid any complaints about uninitialized members
    IntFloat(int x) {
        floatValue = 0;
        intValue = x;
    }

    IntFloat(float x) { 
        intValue = 0;
        floatValue = x;
    }

    public static explicit operator float (IntFloat x) {
        return x.floatValue;
    }

    public static explicit operator int (IntFloat x) { 
        return x.intValue;
    }

    public static explicit operator IntFloat (int i) {
        return new IntFloat(i);
    }
    public static explicit operator IntFloat (float f) { 
        return new IntFloat(f);
    }
}

Then translating InvSqrt is easy.

Demi

Just a word of warning that I just implemented (only for fun) FastInvSqrt using this trick and it was significantly slower than using the more direct implementation with pointers. Of course these days even with pointers FastInvSqrt is slower than 1/Math.Sqrt(x)... – Yuille

Not terribly surprised. Thanks for the feedback! =) – Demi

As a side note: github.com/ekmett/approximate/blob/master/cbits/fast.c provides approximate log, exp, pow, etc. those still at least seem to be a win in terms of clock cycles, but probably not in C#. =) – Demi

Confirmed to be slower, although not very much, nice code, good to know the StructLayout stuff. Also note that repeating the line that comes before the return gives you extra precision (and also takes more time). – Parachronism

Answer 2

Use BitConverter if you want to avoid unsafe code.

float InvSqrt(float x)
{
    float xhalf = 0.5f * x;
    int i = BitConverter.SingleToInt32Bits(x);
    i = 0x5f3759df - (i >> 1);
    x = BitConverter.Int32BitsToSingle(i);
    x = x * (1.5f - xhalf * x * x);
    return x;
}

The code above uses new methods introduced in .NET Core 2.0. For .NET Framework, you have to fall back to the following (which performs allocations):

float InvSqrt(float x)
{
    float xhalf = 0.5f * x;
    int i = BitConverter.ToInt32(BitConverter.GetBytes(x), 0);
    i = 0x5f3759df - (i >> 1);
    x = BitConverter.ToSingle(BitConverter.GetBytes(i), 0);
    x = x * (1.5f - xhalf * x * x);
    return x;
}

Otherwise, the C# code is exactly the same as the C code you gave, except that the method needs to be marked as unsafe:

unsafe float InvSqrt(float x) { ... }

Flambeau

But that does leave open the question: "Is it still fast using BitConverter?" – Pyrrhonism

The BitConverter method does cause two array allocations, which could be a performance issue. (Unfortunately, there's no BitConverter.SingleToInt32Bits method, analogous to BitConverter.DoubleToInt64Bits, which should be both fast and inlineable.) – Flambeau

In .NET Core BitConverter allows reading and writing to Span<T>. No allocations then need to be made. – Steinbach

@Steinbach .NET Core also adds the missing BitConverter.SingleToInt32Bits method I wished for in my eleven-year-old comment. I've updated the answer with that, as it's much simpler than using the Span-related methods. – Flambeau

Answer 3

Definitely possible in unsafe mode. Note that even though in the Quake 3 source code the constant 0x5f3759df was used, numerical research showed that the constant 0x5f375a86 actually yields better results for Newton Approximations.

Headwaiter

Answer 4

I don't see why it wouldn't be possible using the unsafe compiler option.

Necessaries

If somebody programs in unity he has this reason: Some devices and web player don't support unsafe code due to safety reasons. – Randolf