All entries for Wednesday 04 April 2018
April 04, 2018
Valgrind and Extended Precision
Valgrind is an invaluable tool for checking for memory errors in compiled code. Briefly, it runs your code inside its own runtime which includes checked versions of system functions like malloc, so can detect undefined variables, memory leaks and much more. Using this on one of my C++ codes, I found something rather odd. When I ran normally, I had no problems, but inside valgrind I got unexpected NaN results. The actual code, in Minimum Working Example form was
#include <boost/math/special_functions.hpp> #include <iostream> int main(){ double arg = -5.58377179584844451; std::cout<<"Bessel call "<<boost::math::cyl_bessel_j(5+1, (float) arg)<<'\n'; std::cout<<"Bessel call with double "<<boost::math::cyl_bessel_j(5+1, arg)<<'\n'; }
Running normally, I got the expected answer of 0.196642 for both calls. Run inside valgrind, the float version was unchanged, but the double version became NaN! This was very strange. Until I constructed the example above, I assumed the bug was elsewhere within my code, and not with the Bessel functions themselves.
After a lot of head-scratching and searching, I found a link to the valgrind manualand the answer:
Precision: There is no support for 80 bit arithmetic. Internally, Valgrind represents all such "long double" numbers in 64 bits, and so there may be some differences in results. Whether or not this is critical remains to be seen.
After a bit more searching and then checking of DBL_MIN and LDBL_MIN from limits.h I had an answer:
#include <limits> #include <iostream> int main(){ std::cout<<__LDBL_MIN__<<'\n'; std::cout<<__DBL_MIN__<<'\n'; std::cout<<"Minimum value is "<<std::numeric_limits<long double>::min()<<'\n'; }
Normally one gets a result like:
3.3621e-4932 2.22507e-308 Minimum value is 3.3621e-4932
In valgrind once again:
0 2.22507e-308 Minimum value is 0
So the limited 80-bit support referred to was, on my platform at least, breaking __LDBL_MIN__ which should be the smallest representable non zero value. Since the real problem is happening off inside boost::bessel somewhere, there's no simple fix. Most likely, the problem is that boost uses __LDBL_MIN__ to soften a division, preventing a divide by zero which later results in a NaN. The easier workaround requires that you detect running inside valgrind and adjust calculation. This is easy using the supplied macro:
#include "valgrind.h" #ifdef RUNNING_ON_VALGRIND if(RUNNING_ON_VALGRIND){ //Alternative path here } #endif
Note valgrind.h is specially licensed, differently to the rest of Valgrind, so that including it in your code does not require you open-source your code. Inside the special path, I used a float version of the bessel call, since inside valgrind I am interested in general output, but don't need the full precision of a double.
Another alternative overrides __LDBL_MIN__ before including the boost header. I chose to set it equal to __DBL_MIN__ Note this only works with header-only libraries and is a bit risky. The code becomes something like:
#undef __LDBL_MIN__ #define __LDBL_MIN__ __DBL_MIN__ #include <boost/math/special_functions.hpp> #include <limits> #include <iostream> int main(){ std::cout<<__LDBL_MIN__<<'\n'; std::cout<<__DBL_MIN__<<'\n'; std::cout<<"Minimum value is "<<std::numeric_limits<long double>::min()<<'\n'; double arg = -5.58377179584844451; std::cout<<"Bessel call with double "<<boost::math::cyl_bessel_j(5+1, arg)<<'\n'; }
Sometimes, the bug really is in the compiler or tool!