Julia: Loss of accuracy of BigFloat on Win10-64bit plattform

[egarten]

I found an issue regarding BigFloats using julia0.6.0 on a Win10 system. This error does not appear on linux or macOS systems:

# Set precision to 7 million bits, so I expect about 2 million valid digits
setprecision(7000000)
# Get value of pi with high precision
bigpi=BigFloat(pi)
# Calculate difference of sin(pi/2)-1. Should be zero in exact arithmetic
sin(bigpi/BigFloat(2))-1
# Calualate difference of sin(pi/6)-1/2. Should be zero in exact arithmetic
sin(bigpi/BigFloat(6))-1/BigFloat(2)

For the first difference I get the value of zero. That means the result of sin(pi/2) is valid up to 2 million digits.

For the seconds difference I get an value of about 4.76993...e-306923. That means sin(pi/6) is only calculated up to about 300,000 digits and not 2 million digits.

Some more details have already been discussed:
https://discourse.julialang.org/t/accuracy-of-bigfloat-calculating-sin-with-high-accuracy/4604

[simonbyrne]

Anyone want to try it on 32-bit linux? Would be good to know if it's a Clong issue.

FWIW, there are some accuracy bounds of the implementation in http://www.mpfr.org/algorithms.pdf

[stevengj]

Note that the problem seems to be specifically in the sin function. If I do big(pi)/6 in this precision on MacOS (where sin is accurate), the last few digits are 399664339975e-01. On Win64, I get exactly the same answer, so it is computing π/6 correctly.

[vtjnash]

For any library which uses libm where we care about numeric accuracy, we should be statically linking a copy of libopenlibm on platforms where the system libm is not accurate. We usually do this on Windows (since the old win32 crt math library was not very accurate) for many libraries. Do we know if MPFR is using libm functions though?

[simonbyrne]

I wouldn't think so.

[simonbyrne]

I see the same issue on a 32-bit linux VM.

[simonbyrne]

Okay, this is an issue with MPFR. The code here on 32-bit linux gives:

Diff is 4.7699306332221186288894962846258624262005363580770e-306923

[stevengj]

@simonbyrne, thanks for looking into this. Can you send a bug report to [email protected] as described here?

[simonbyrne]

Will do.

[simonbyrne]

See https://sympa.inria.fr/sympa/arc/mpfr/2017-07/msg00002.html

[tkelman]

looks like there's a patch already that we could try applying?

[stevengj]

Great, looks like they have a fix for MPFR 3.1. We could apply a patch before then if need be.

[stevengj]

MPFR 3.1.5 includes the fix, so once we upgrade then we can delete our patch.

[simonbyrne]

I just tried this on 0.6.1, which is supposed to have the fix in, but still see the issue.

[simonbyrne]

Ah, it was fixed in the 3.1.6 release, not 3.1.5 (which was released over a year ago). The link @stevengj pointed to are the bugs _since_ the 3.1.5 release.

[StefanKarpinski]

Should we have a Julia-side regression test for this? Defense in depth and all that.

[KristofferC]

Don't we have that? https://github.com/JuliaLang/julia/commit/35011da9a9d2ebcbcfbc07d398615e2fb4260ed1#diff-26b791854ff5dbe1c3ed1392119bcf42R892