float_comparison.qbk 6.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120
  1. [section:float_comparison Floating-point Comparison]
  2. [import ../../example/float_comparison_example.cpp]
  3. Comparison of floating-point values has always been a source of endless difficulty and confusion.
  4. Unlike integral values that are exact, all floating-point operations
  5. will potentially produce an inexact result that will be rounded to the nearest
  6. available binary representation. Even apparently inocuous operations such as assigning
  7. 0.1 to a double produces an inexact result (as this decimal number has no
  8. exact binary representation).
  9. Floating-point computations also involve rounding so that some 'computational noise' is added,
  10. and hence results are also not exact (although repeatable, at least under identical platforms and compile options).
  11. Sadly, this conflicts with the expectation of most users, as many articles and innumerable cries for help show all too well.
  12. Some background reading is:
  13. * Knuth D.E. The art of computer programming, vol II, section 4.2, especially Floating-Point Comparison 4.2.2, pages 198-220.
  14. * [@http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html David Goldberg, "What Every Computer Scientist Should Know About Floating-Point Arithmetic"]
  15. * [@http://adtmag.com/articles/2000/03/16/comparing-floatshow-to-determine-if-floating-quantities-are-close-enough-once-a-tolerance-has-been-r.aspx
  16. Alberto Squassabia, Comparing floats listing]
  17. * [@https://code.google.com/p/googletest/wiki/AdvancedGuide#Floating-Point_Comparison Google Floating-Point_Comparison guide]
  18. * [@https://www.boost.org/doc/libs/release/libs/test/doc/html/boost_test/testing_tools/extended_comparison/floating_point.html Boost.Test Floating-Point_Comparison]
  19. Boost provides a number of ways to compare floating-point values to see if they are tolerably close enough to each other,
  20. but first we must decide what kind of comparison we require:
  21. * Absolute difference/error: the absolute difference between two values ['a] and ['b] is simply `fabs(a-b)`.
  22. This is the only meaningful comparison to make if we know that the result may have cancellation error (see below).
  23. * The edit distance between the two values: i.e. how many (binary) floating-point values are between two values ['a] and ['b]?
  24. This is provided by the function __float_distance, but is probably only useful when you know that the distance should be very small.
  25. This function is somewhat difficult to compute, and doesn't scale to values that are very far apart. In other words, use with care.
  26. * The relative distance/error between two values. This is quick and easy to compute, and is generally the method of choice when
  27. checking that your results are "tolerably close" to one another. However, it is not as exact as the edit distance when dealing
  28. with small differences, and due to the way floating-point values are encoded can "wobble" by a factor of 2 compared to the "true"
  29. edit distance. This is the method documented below: if `float_distance` is a surgeon's scalpel, then `relative_difference` is more like a Swiss army knife: both have important but different use cases.
  30. [h5:fp_relative Relative Comparison of Floating-point Values]
  31. `#include <boost/math/special_functions/relative_difference.hpp>`
  32. template <class T, class U>
  33. ``__sf_result`` relative_difference(T a, U b);
  34. template <class T, class U>
  35. ``__sf_result`` epsilon_difference(T a, U b);
  36. The function `relative_difference` returns the relative distance/error ['E] between two values as defined by:
  37. [expression E = fabs((a - b) / min(a,b))]
  38. The function `epsilon_difference` is a convenience function that returns `relative_difference(a, b) / eps` where
  39. `eps` is the machine epsilon for the result type.
  40. The following special cases are handled as follows:
  41. * If either of ['a] or ['b] is a NaN, then returns the largest representable value for T: for example for type `double`, this
  42. is `std::numeric_limits<double>::max()` which is the same as `DBL_MAX` or `1.7976931348623157e+308`.
  43. * If ['a] and ['b] differ in sign then returns the largest representable value for T.
  44. * If both ['a] and ['b] are both infinities (of the same sign), then returns zero.
  45. * If just one of ['a] and ['b] is an infinity, then returns the largest representable value for T.
  46. * If both ['a] and ['b] are zero then returns zero.
  47. * If just one of ['a] or ['b] is a zero or a denormalized value, then it is treated as if it were the
  48. smallest (non-denormalized) value representable in T for the purposes of the above calculation.
  49. These rules were primarily designed to assist with our own test suite, they are designed to be robust enough
  50. that the function can in most cases be used blindly, including in cases where the expected result is actually
  51. too small to represent in type T and underflows to zero.
  52. [h5 Examples]
  53. [compare_floats_using]
  54. [compare_floats_example_1]
  55. [compare_floats_example_2]
  56. [compare_floats_example_3]
  57. [compare_floats_example_4]
  58. [compare_floats_example_5]
  59. [compare_floats_example_6]
  60. All the above examples are contained in [@../../example/float_comparison_example.cpp float_comparison_example.cpp].
  61. [h5:small Handling Absolute Errors]
  62. Imagine we're testing the following function:
  63. double myspecial(double x)
  64. {
  65. return sin(x) - sin(4 * x);
  66. }
  67. This function has multiple roots, some of which are quite predicable in that both
  68. `sin(x)` and `sin(4x)` are zero together. Others occur because the values returned
  69. from those two functions precisely cancel out. At such points the relative difference
  70. between the true value of the function and the actual value returned may be ['arbitrarily
  71. large] due to [@http://en.wikipedia.org/wiki/Loss_of_significance cancellation error].
  72. In such a case, testing the function above by requiring that the values returned by
  73. `relative_error` or `epsilon_error` are below some threshold is pointless: the best
  74. we can do is to verify that the ['absolute difference] between the true
  75. and calculated values is below some threshold.
  76. Of course, determining what that threshold should be is often tricky,
  77. but a good starting point would be machine epsilon multiplied by the largest
  78. of the values being summed. In the example above, the largest value returned
  79. by `sin(whatever)` is 1, so simply using machine epsilon as the target for
  80. maximum absolute difference might be a good start (though in practice we may need
  81. a slightly higher value - some trial and error will be necessary).
  82. [endsect] [/section:float_comparison Floating-point comparison]
  83. [/
  84. Copyright 2015 John Maddock and Paul A. Bristow.
  85. Distributed under the Boost Software License, Version 1.0.
  86. (See accompanying file LICENSE_1_0.txt or copy at
  87. http://www.boost.org/LICENSE_1_0.txt).
  88. ]