diff options
Diffstat (limited to 'src/boost/libs/math/example/float_comparison_example.cpp')
| -rw-r--r-- | src/boost/libs/math/example/float_comparison_example.cpp | 444 |
1 files changed, 444 insertions, 0 deletions
diff --git a/src/boost/libs/math/example/float_comparison_example.cpp b/src/boost/libs/math/example/float_comparison_example.cpp new file mode 100644 index 00000000..6c892fa2 --- /dev/null +++ b/src/boost/libs/math/example/float_comparison_example.cpp @@ -0,0 +1,444 @@ +//!file +//! \brief floating-point comparison from Boost.Test +// Copyright Paul A. Bristow 2015. +// Copyright John Maddock 2015. + +// Use, modification and distribution are subject to the +// Boost Software License, Version 1.0. +// (See accompanying file LICENSE_1_0.txt +// or copy at http://www.boost.org/LICENSE_1_0.txt) + +// Note that this file contains Quickbook mark-up as well as code +// and comments, don't change any of the special comment mark-ups! + +#include <boost/math/special_functions/relative_difference.hpp> +#include <boost/math/special_functions/next.hpp> + +#include <iostream> +#include <limits> // for std::numeric_limits<T>::epsilon(). + +int main() +{ + std::cout << "Compare floats using Boost.Math functions/classes" << std::endl; + + +//[compare_floats_using +/*`Some using statements will ensure that the functions we need are accessible. +*/ + + using namespace boost::math; + +//`or + + using boost::math::relative_difference; + using boost::math::epsilon_difference; + using boost::math::float_next; + using boost::math::float_prior; + +//] [/compare_floats_using] + + +//[compare_floats_example_1 +/*`The following examples display values with all possibly significant digits. +Newer compilers should provide `std::numeric_limits<FPT>::max_digits10` +for this purpose, and here we use `float` precision where `max_digits10` = 9 +to avoid displaying a distracting number of decimal digits. + +[note Older compilers can use this formula to calculate `max_digits10` +from `std::numeric_limits<FPT>::digits10`: +__spaces `int max_digits10 = 2 + std::numeric_limits<FPT>::digits10 * 3010/10000;` +] [/note] + +One can set the display including all trailing zeros +(helpful for this example to show all potentially significant digits), +and also to display `bool` values as words rather than integers: +*/ + std::cout.precision(std::numeric_limits<float>::max_digits10); + std::cout << std::boolalpha << std::showpoint << std::endl; + +//] [/compare_floats_example_1] + +//[compare_floats_example_2] +/*` +When comparing values that are ['quite close] or ['approximately equal], +we could use either `float_distance` or `relative_difference`/`epsilon_difference`, for example +with type `float`, these two values are adjacent to each other: +*/ + + float a = 1; + float b = 1 + std::numeric_limits<float>::epsilon(); + std::cout << "a = " << a << std::endl; + std::cout << "b = " << b << std::endl; + std::cout << "float_distance = " << float_distance(a, b) << std::endl; + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl; + +/*` +Which produces the output: + +[pre +a = 1.00000000 +b = 1.00000012 +float_distance = 1.00000000 +relative_difference = 1.19209290e-007 +epsilon_difference = 1.00000000 +] +*/ + //] [/compare_floats_example_2] + +//[compare_floats_example_3] +/*` +In the example above, it just so happens that the edit distance as measured by `float_distance`, and the +difference measured in units of epsilon were equal. However, due to the way floating point +values are represented, that is not always the case:*/ + + a = 2.0f / 3.0f; // 2/3 inexactly represented as a float + b = float_next(float_next(float_next(a))); // 3 floating point values above a + std::cout << "a = " << a << std::endl; + std::cout << "b = " << b << std::endl; + std::cout << "float_distance = " << float_distance(a, b) << std::endl; + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl; + +/*` +Which produces the output: + +[pre +a = 0.666666687 +b = 0.666666865 +float_distance = 3.00000000 +relative_difference = 2.68220901e-007 +epsilon_difference = 2.25000000 +] + +There is another important difference between `float_distance` and the +`relative_difference/epsilon_difference` functions in that `float_distance` +returns a signed result that reflects which argument is larger in magnitude, +where as `relative_difference/epsilon_difference` simply return an unsigned +value that represents how far apart the values are. For example if we swap +the order of the arguments: +*/ + + std::cout << "float_distance = " << float_distance(b, a) << std::endl; + std::cout << "relative_difference = " << relative_difference(b, a) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(b, a) << std::endl; + + /*` + The output is now: + + [pre + float_distance = -3.00000000 + relative_difference = 2.68220901e-007 + epsilon_difference = 2.25000000 + ] +*/ + //] [/compare_floats_example_3] + +//[compare_floats_example_4] +/*` +Zeros are always treated as equal, as are infinities as long as they have the same sign:*/ + + a = 0; + b = -0; // signed zero + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + a = b = std::numeric_limits<float>::infinity(); + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "relative_difference = " << relative_difference(a, -b) << std::endl; + +/*` +Which produces the output: + +[pre +relative_difference = 0.000000000 +relative_difference = 0.000000000 +relative_difference = 3.40282347e+038 +] +*/ +//] [/compare_floats_example_4] + +//[compare_floats_example_5] +/*` +Note that finite values are always infinitely far away from infinities even if those finite values are very large:*/ + + a = (std::numeric_limits<float>::max)(); + b = std::numeric_limits<float>::infinity(); + std::cout << "a = " << a << std::endl; + std::cout << "b = " << b << std::endl; + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl; + +/*` +Which produces the output: + +[pre +a = 3.40282347e+038 +b = 1.#INF0000 +relative_difference = 3.40282347e+038 +epsilon_difference = 3.40282347e+038 +] +*/ +//] [/compare_floats_example_5] + +//[compare_floats_example_6] +/*` +Finally, all denormalized values and zeros are treated as being effectively equal:*/ + + a = std::numeric_limits<float>::denorm_min(); + b = a * 2; + std::cout << "a = " << a << std::endl; + std::cout << "b = " << b << std::endl; + std::cout << "float_distance = " << float_distance(a, b) << std::endl; + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl; + a = 0; + std::cout << "a = " << a << std::endl; + std::cout << "b = " << b << std::endl; + std::cout << "float_distance = " << float_distance(a, b) << std::endl; + std::cout << "relative_difference = " << relative_difference(a, b) << std::endl; + std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl; + +/*` +Which produces the output: + +[pre +a = 1.40129846e-045 +b = 2.80259693e-045 +float_distance = 1.00000000 +relative_difference = 0.000000000 +epsilon_difference = 0.000000000 +a = 0.000000000 +b = 2.80259693e-045 +float_distance = 2.00000000 +relative_difference = 0.000000000 +epsilon_difference = 0.000000000] + +Notice how, in the above example, two denormalized values that are a factor of 2 apart are +none the less only one representation apart! + +*/ +//] [/compare_floats_example_6] + + +#if 0 +//[old_compare_floats_example_3 +//`The simplest use is to compare two values with a tolerance thus: + + bool is_close = is_close_to(1.F, 1.F + epsilon, epsilon); // One epsilon apart is close enough. + std::cout << "is_close_to(1.F, 1.F + epsilon, epsilon); is " << is_close << std::endl; // true + + is_close = is_close_to(1.F, 1.F + 2 * epsilon, epsilon); // Two epsilon apart isn't close enough. + std::cout << "is_close_to(1.F, 1.F + epsilon, epsilon); is " << is_close << std::endl; // false + +/*` +[note The type FPT of the tolerance and the type of the values [*must match]. + +So `is_close(0.1F, 1., 1.)` will fail to compile because "template parameter 'FPT' is ambiguous". +Always provide the same type, using `static_cast<FPT>` if necessary.] +*/ + + +/*`An instance of class `close_at_tolerance` is more convenient +when multiple tests with the same conditions are planned. +A class that stores a tolerance of three epsilon (and the default ['strong] test) is: +*/ + + close_at_tolerance<float> three_rounds(3 * epsilon); // 'strong' by default. + +//`and we can confirm these settings: + + std::cout << "fraction_tolerance = " + << three_rounds.fraction_tolerance() + << std::endl; // +3.57627869e-007 + std::cout << "strength = " + << (three_rounds.strength() == FPC_STRONG ? "strong" : "weak") + << std::endl; // strong + +//`To start, let us use two values that are truly equal (having identical bit patterns) + + float a = 1.23456789F; + float b = 1.23456789F; + +//`and make a comparison using our 3*epsilon `three_rounds` functor: + + bool close = three_rounds(a, b); + std::cout << "three_rounds(a, b) = " << close << std::endl; // true + +//`Unsurprisingly, the result is true, and the failed fraction is zero. + + std::cout << "failed_fraction = " << three_rounds.failed_fraction() << std::endl; + +/*`To get some nearby values, it is convenient to use the Boost.Math __next_float functions, +for which we need an include + + #include <boost/math/special_functions/next.hpp> + +and some using declarations: +*/ + + using boost::math::float_next; + using boost::math::float_prior; + using boost::math::nextafter; + using boost::math::float_distance; + +//`To add a few __ulp to one value: + b = float_next(a); // Add just one ULP to a. + b = float_next(b); // Add another one ULP. + b = float_next(b); // Add another one ULP. + // 3 epsilon would pass. + b = float_next(b); // Add another one ULP. + +//`and repeat our comparison: + + close = three_rounds(a, b); + std::cout << "three_rounds(a, b) = " << close << std::endl; // false + std::cout << "failed_fraction = " << three_rounds.failed_fraction() + << std::endl; // abs(u-v) / abs(v) = 3.86237957e-007 + +//`We can also 'measure' the number of bits different using the `float_distance` function: + + std::cout << "float_distance = " << float_distance(a, b) << std::endl; // 4 + +/*`Now consider two values that are much further apart +than one might expect from ['computational noise], +perhaps the result of two measurements of some physical property like length +where an uncertainty of a percent or so might be expected. +*/ + float fp1 = 0.01000F; + float fp2 = 0.01001F; // Slightly different. + + float tolerance = 0.0001F; + + close_at_tolerance<float> strong(epsilon); // Default is strong. + bool rs = strong(fp1, fp2); + std::cout << "strong(fp1, fp2) is " << rs << std::endl; + +//`Or we could contrast using the ['weak] criterion: + close_at_tolerance<float> weak(epsilon, FPC_WEAK); // Explicitly weak. + bool rw = weak(fp1, fp2); // + std::cout << "weak(fp1, fp2) is " << rw << std::endl; + +//`We can also construct, setting tolerance and strength, and compare in one statement: + + std::cout << a << " #= " << b << " is " + << close_at_tolerance<float>(epsilon, FPC_STRONG)(a, b) << std::endl; + std::cout << a << " ~= " << b << " is " + << close_at_tolerance<float>(epsilon, FPC_WEAK)(a, b) << std::endl; + +//`but this has little advantage over using function `is_close_to` directly. + +//] [/old_compare_floats_example_3] + + +/*When the floating-point values become very small and near zero, using +//a relative test becomes unhelpful because one is dividing by zero or tiny, + +//Instead, an absolute test is needed, comparing one (or usually both) values with zero, +//using a tolerance. +//This is provided by the `small_with_tolerance` class and `is_small` function. + + namespace boost { + namespace math { + namespace fpc { + + + template<typename FPT> + class small_with_tolerance + { + public: + // Public typedefs. + typedef bool result_type; + + // Constructor. + explicit small_with_tolerance(FPT tolerance); // tolerance >= 0 + + // Functor + bool operator()(FPT value) const; // return true if <= absolute tolerance (near zero). + }; + + template<typename FPT> + bool + is_small(FPT value, FPT tolerance); // return true if value <= absolute tolerance (near zero). + + }}} // namespaces. + +/*` +[note The type FPT of the tolerance and the type of the value [*must match]. + +So `is_small(0.1F, 0.000001)` will fail to compile because "template parameter 'FPT' is ambiguous". +Always provide the same type, using `static_cast<FPT>` if necessary.] + +A few values near zero are tested with varying tolerance below. +*/ +//[compare_floats_small_1 + + float c = 0; + std::cout << "0 is_small " << is_small(c, epsilon) << std::endl; // true + + c = std::numeric_limits<float>::denorm_min(); // 1.40129846e-045 + std::cout << "denorm_ min =" << c << ", is_small is " << is_small(c, epsilon) << std::endl; // true + + c = (std::numeric_limits<float>::min)(); // 1.17549435e-038 + std::cout << "min = " << c << ", is_small is " << is_small(c, epsilon) << std::endl; // true + + c = 1 * epsilon; // 1.19209290e-007 + std::cout << "epsilon = " << c << ", is_small is " << is_small(c, epsilon) << std::endl; // false + + c = 1 * epsilon; // 1.19209290e-007 + std::cout << "2 epsilon = " << c << ", is_small is " << is_small(c, 2 * epsilon) << std::endl; // true + + c = 2 * epsilon; //2.38418579e-007 + std::cout << "4 epsilon = " << c << ", is_small is " << is_small(c, 2 * epsilon) << std::endl; // false + + c = 0.00001F; + std::cout << "0.00001 = " << c << ", is_small is " << is_small(c, 0.0001F) << std::endl; // true + + c = -0.00001F; + std::cout << "0.00001 = " << c << ", is_small is " << is_small(c, 0.0001F) << std::endl; // true + +/*`Using the class `small_with_tolerance` allows storage of the tolerance, +convenient if you make repeated tests with the same tolerance. +*/ + + small_with_tolerance<float>my_test(0.01F); + + std::cout << "my_test(0.001F) is " << my_test(0.001F) << std::endl; // true + std::cout << "my_test(0.001F) is " << my_test(0.01F) << std::endl; // false + + //] [/compare_floats_small_1] +#endif + return 0; +} // int main() + +/* + +Example output is: + +//[compare_floats_output +Compare floats using Boost.Test functions/classes + +float epsilon = 1.19209290e-007 +is_close_to(1.F, 1.F + epsilon, epsilon); is true +is_close_to(1.F, 1.F + epsilon, epsilon); is false +fraction_tolerance = 3.57627869e-007 +strength = strong +three_rounds(a, b) = true +failed_fraction = 0.000000000 +three_rounds(a, b) = false +failed_fraction = 3.86237957e-007 +float_distance = 4.00000000 +strong(fp1, fp2) is false +weak(fp1, fp2) is false +1.23456788 #= 1.23456836 is false +1.23456788 ~= 1.23456836 is false +0 is_small true +denorm_ min =1.40129846e-045, is_small is true +min = 1.17549435e-038, is_small is true +epsilon = 1.19209290e-007, is_small is false +2 epsilon = 1.19209290e-007, is_small is true +4 epsilon = 2.38418579e-007, is_small is false +0.00001 = 9.99999975e-006, is_small is true +0.00001 = -9.99999975e-006, is_small is true +my_test(0.001F) is true + +my_test(0.001F) is false//] [/compare_floats_output] +*/ |