is this a bug in sqrt function

Question

I create an application to compute the primes numbers in 64-bit range so when i tried to compute the square root of 64-bit number using sqrt function from math.h i found the answer is not accurate for example when the input is ~0ull the answer should be ~0u but the one I get is 0x100000000 which is not right, so i decided to create my own version using assembly x86 language to see if this is a bug, here is my function:

inline unsigned prime_isqrt(unsigned long long value)
{
    const unsigned one = 1;
    const unsigned two = 2;

    __asm
    {
        test dword ptr [value+4], 0x80000000
        jz ZERO
        mov eax, dword ptr [value]
        mov ecx, dword ptr [value + 4]

        shrd eax, ecx, 1
        shr  ecx, 1
        mov  dword ptr [value],eax 
        mov  dword ptr [value+4],ecx

        fild value
        fimul two
        fiadd one
        jmp REST
ZERO: 
        fild value
REST: 
        fsqrt
        fisttp value
        mov eax, dword ptr [value]
    }
}

the input is an odd number to get its square root. When i test my function with same input the result was the same.

What i don't get is why those functions round the result or to be specific why sqrt instruction round the result?

Answer 1

sqrt doesn't round anything - you do when you convert your integer into a double. A double can't represent all the numbers that a 64-bit integer can without loss of precision. Specifically starting at 2⁵³, there are multiple integers that will be represented as the same double value.

So if you convert an integer above 2⁵³ to double, you lose some of the least significant bits, which is why (double)(~0ull) is 18446744073709552000.0, not 18446744073709551615.0 (or to be more precise the latter is actually equal to the former because they represent the same double number).

Answer 2

You're not very clear about the C++ function that you're calling. sqrt is an overloaded name. You probably wanted sqrt(double(~0ull)). There's no sqrt overload which takes an unsigned long long.