Incorrect sqrt result (is there a maximal precision?)

[jakobkogler]

The following code will print about 21'000 digits (which is what you expect given precision 70'000 and given that one decimal digits is about 3.3 binary digits).

from bigfloat import *
print(sqrt(2, precision(70000)))

However, from those digits only about the first 5'000 digits are correct. And the last 5'000 digits are all zeros.
Increasing the precision will not change the result, it will just add more zeros at the end.

This clearly looks like a bug to me. I wrote the same thing in C++ (using MPFR C++), and there I get the correct digits. The first 20'000 of them are correct (comparing to https://oeis.org/A002193/b002193.txt).

#include <iostream>
#include "mpreal.h"

int main() {
    using mpfr::mpreal;    
    mpreal::set_default_prec(70'000);
    std::cout.precision(20'000);
    std::cout << mpfr::sqrt(2) << std::endl;
}

[mdickinson]

Thanks for the report. I confirm that I see the same results as you.

I think what you're seeing is due to the default context having an emin of -16493 (to match IEEE 754 quadruple precision). This means that because the default context also has subnormalize=True (again to match IEEE 754), the result ends up being quantized to be an integer multiple of some power of two 2**-16493. Can you try decreasing the emin to something like -999999999 and see what happens?

EDIT: subnormalize=True isn't relevant here; edited the comment to fix.

[mdickinson]

It's possible that we want to change the way that the context treats emin here - that it should be automatically adjusted when the precision is increased. The emin exposed in the context matches MPFR's concept of emin, but not the traditional IEEE 754 concept of the minimum exponent - it may be that we want to adjust the way the context interprets emin to match IEEE 754 instead.

[mdickinson]

Marking as a bug: the code is technically behaving as intended, but that behaviour is surprising; I think this is a design bug that we should try and fix.

[jakobkogler]

I can confirm that sqrt(2, Context(precision=70000, emin=-9999999)), does indeed compute the correct result.

And yes, that result is surprising.
Especially as a human you don't usually look at last digits of the digits expansion. You just assume that it is correct, if you use a library that claims that it can do high precision arithmetic. So forgetting about setting the emin correctly it's a very simple error to make.