Accuracy of long double sqrt ()

Question

I noticed a problem accurate to the double double version of sqrt (). The following code demonstrates the problem.

#include <iostream> #include <cmath> #include <cfloat> int main(int argc, char ** argv) { int count=0; long double s=sqrt(3L); std::cout.precision(21); std::cout << "s=" << s << ", s^2=" << s*s << std::endl; while( s*s<3L+LDBL_EPSILON ) { s+=LDBL_EPSILON; std::cout << s << ' ' << s*s << std::endl; ++count; } std::cout << "sqrt(3L) is approximately " << count << " multiples of LDBL_EPSILON away from the correct value." << std::endl; return 0; } 

Compiling and running this with

 >g++ -o sqrt sqrt.cpp && ./sqrt 

gives

 s=1.73205080756887719318, s^2=2.9999999999999996524 1.73205080756887719329 2.99999999999999965284 1.73205080756887719339 2.99999999999999965306 ... (922 lines omitted) 1.73205080756887729347 2.99999999999999999978 1.73205080756887729357 3.00000000000000000022 sqrt(3L) is approximately 926 multiples of LDBL_EPSILON away from the correct value. 

The regular double version of sqrt () gives the double closest to the actual value.

The g ++ version I'm using is

 >g++ -v Using built-in specs. Target: x86_64-linux-gnu Configured with: ../src/configure -v --with-pkgversion='Debian 4.4.5-8' --with-bugurl=file:///usr/share/doc/gcc-4.4/README.Bugs --enable-languages=c,c++,fortran,objc,obj-c++ --prefix=/usr --program-suffix=-4.4 --enable-shared --enable-multiarch --enable-linker-build-id --with-system-zlib --libexecdir=/usr/lib --without-included-gettext --enable-threads=posix --with-gxx-include-dir=/usr/include/c++/4.4 --libdir=/usr/lib --enable-nls --enable-clocale=gnu --enable-libstdcxx-debug --enable-objc-gc --with-arch-32=i586 --with-tune=generic --enable-checking=release --build=x86_64-linux-gnu --host=x86_64-linux-gnu --target=x86_64-linux-gnu Thread model: posix gcc version 4.4.5 (Debian 4.4.5-8) 

Is this a known bug? Should I report this somewhere?

Answer 1

Here you have two problems: firstly, 3L implicitly advances to double not long double , so although you assign the return value to long double , it still uses the low-precision version of sqrt . As an argument, you need static_cast 3 - long double . Secondly, only the double sqrt version is imported into the global namespace, because function overloading is not supported in C, you should use std::sqrt instead.

In this way:

long double s=std::sqrt(static_cast<long double>(3));

Answer 2

Does the value from the "regular" version of double sqrt () experience a greater degree of rounding than a long double? This is what we expect. Perhaps this "grainy" rounding is close to the correct value - closer than a long double sqrt.

To check this, you need to try several values ​​and compare.