float_comparison_example.cpp 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444
  1. //!file
  2. //! \brief floating-point comparison from Boost.Test
  3. // Copyright Paul A. Bristow 2015.
  4. // Copyright John Maddock 2015.
  5. // Use, modification and distribution are subject to the
  6. // Boost Software License, Version 1.0.
  7. // (See accompanying file LICENSE_1_0.txt
  8. // or copy at http://www.boost.org/LICENSE_1_0.txt)
  9. // Note that this file contains Quickbook mark-up as well as code
  10. // and comments, don't change any of the special comment mark-ups!
  11. #include <boost/math/special_functions/relative_difference.hpp>
  12. #include <boost/math/special_functions/next.hpp>
  13. #include <iostream>
  14. #include <limits> // for std::numeric_limits<T>::epsilon().
  15. int main()
  16. {
  17. std::cout << "Compare floats using Boost.Math functions/classes" << std::endl;
  18. //[compare_floats_using
  19. /*`Some using statements will ensure that the functions we need are accessible.
  20. */
  21. using namespace boost::math;
  22. //`or
  23. using boost::math::relative_difference;
  24. using boost::math::epsilon_difference;
  25. using boost::math::float_next;
  26. using boost::math::float_prior;
  27. //] [/compare_floats_using]
  28. //[compare_floats_example_1
  29. /*`The following examples display values with all possibly significant digits.
  30. Newer compilers should provide `std::numeric_limits<FPT>::max_digits10`
  31. for this purpose, and here we use `float` precision where `max_digits10` = 9
  32. to avoid displaying a distracting number of decimal digits.
  33. [note Older compilers can use this formula to calculate `max_digits10`
  34. from `std::numeric_limits<FPT>::digits10`:
  35. __spaces `int max_digits10 = 2 + std::numeric_limits<FPT>::digits10 * 3010/10000;`
  36. ] [/note]
  37. One can set the display including all trailing zeros
  38. (helpful for this example to show all potentially significant digits),
  39. and also to display `bool` values as words rather than integers:
  40. */
  41. std::cout.precision(std::numeric_limits<float>::max_digits10);
  42. std::cout << std::boolalpha << std::showpoint << std::endl;
  43. //] [/compare_floats_example_1]
  44. //[compare_floats_example_2]
  45. /*`
  46. When comparing values that are ['quite close] or ['approximately equal],
  47. we could use either `float_distance` or `relative_difference`/`epsilon_difference`, for example
  48. with type `float`, these two values are adjacent to each other:
  49. */
  50. float a = 1;
  51. float b = 1 + std::numeric_limits<float>::epsilon();
  52. std::cout << "a = " << a << std::endl;
  53. std::cout << "b = " << b << std::endl;
  54. std::cout << "float_distance = " << float_distance(a, b) << std::endl;
  55. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  56. std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl;
  57. /*`
  58. Which produces the output:
  59. [pre
  60. a = 1.00000000
  61. b = 1.00000012
  62. float_distance = 1.00000000
  63. relative_difference = 1.19209290e-007
  64. epsilon_difference = 1.00000000
  65. ]
  66. */
  67. //] [/compare_floats_example_2]
  68. //[compare_floats_example_3]
  69. /*`
  70. In the example above, it just so happens that the edit distance as measured by `float_distance`, and the
  71. difference measured in units of epsilon were equal. However, due to the way floating point
  72. values are represented, that is not always the case:*/
  73. a = 2.0f / 3.0f; // 2/3 inexactly represented as a float
  74. b = float_next(float_next(float_next(a))); // 3 floating point values above a
  75. std::cout << "a = " << a << std::endl;
  76. std::cout << "b = " << b << std::endl;
  77. std::cout << "float_distance = " << float_distance(a, b) << std::endl;
  78. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  79. std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl;
  80. /*`
  81. Which produces the output:
  82. [pre
  83. a = 0.666666687
  84. b = 0.666666865
  85. float_distance = 3.00000000
  86. relative_difference = 2.68220901e-007
  87. epsilon_difference = 2.25000000
  88. ]
  89. There is another important difference between `float_distance` and the
  90. `relative_difference/epsilon_difference` functions in that `float_distance`
  91. returns a signed result that reflects which argument is larger in magnitude,
  92. where as `relative_difference/epsilon_difference` simply return an unsigned
  93. value that represents how far apart the values are. For example if we swap
  94. the order of the arguments:
  95. */
  96. std::cout << "float_distance = " << float_distance(b, a) << std::endl;
  97. std::cout << "relative_difference = " << relative_difference(b, a) << std::endl;
  98. std::cout << "epsilon_difference = " << epsilon_difference(b, a) << std::endl;
  99. /*`
  100. The output is now:
  101. [pre
  102. float_distance = -3.00000000
  103. relative_difference = 2.68220901e-007
  104. epsilon_difference = 2.25000000
  105. ]
  106. */
  107. //] [/compare_floats_example_3]
  108. //[compare_floats_example_4]
  109. /*`
  110. Zeros are always treated as equal, as are infinities as long as they have the same sign:*/
  111. a = 0;
  112. b = -0; // signed zero
  113. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  114. a = b = std::numeric_limits<float>::infinity();
  115. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  116. std::cout << "relative_difference = " << relative_difference(a, -b) << std::endl;
  117. /*`
  118. Which produces the output:
  119. [pre
  120. relative_difference = 0.000000000
  121. relative_difference = 0.000000000
  122. relative_difference = 3.40282347e+038
  123. ]
  124. */
  125. //] [/compare_floats_example_4]
  126. //[compare_floats_example_5]
  127. /*`
  128. Note that finite values are always infinitely far away from infinities even if those finite values are very large:*/
  129. a = (std::numeric_limits<float>::max)();
  130. b = std::numeric_limits<float>::infinity();
  131. std::cout << "a = " << a << std::endl;
  132. std::cout << "b = " << b << std::endl;
  133. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  134. std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl;
  135. /*`
  136. Which produces the output:
  137. [pre
  138. a = 3.40282347e+038
  139. b = 1.#INF0000
  140. relative_difference = 3.40282347e+038
  141. epsilon_difference = 3.40282347e+038
  142. ]
  143. */
  144. //] [/compare_floats_example_5]
  145. //[compare_floats_example_6]
  146. /*`
  147. Finally, all denormalized values and zeros are treated as being effectively equal:*/
  148. a = std::numeric_limits<float>::denorm_min();
  149. b = a * 2;
  150. std::cout << "a = " << a << std::endl;
  151. std::cout << "b = " << b << std::endl;
  152. std::cout << "float_distance = " << float_distance(a, b) << std::endl;
  153. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  154. std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl;
  155. a = 0;
  156. std::cout << "a = " << a << std::endl;
  157. std::cout << "b = " << b << std::endl;
  158. std::cout << "float_distance = " << float_distance(a, b) << std::endl;
  159. std::cout << "relative_difference = " << relative_difference(a, b) << std::endl;
  160. std::cout << "epsilon_difference = " << epsilon_difference(a, b) << std::endl;
  161. /*`
  162. Which produces the output:
  163. [pre
  164. a = 1.40129846e-045
  165. b = 2.80259693e-045
  166. float_distance = 1.00000000
  167. relative_difference = 0.000000000
  168. epsilon_difference = 0.000000000
  169. a = 0.000000000
  170. b = 2.80259693e-045
  171. float_distance = 2.00000000
  172. relative_difference = 0.000000000
  173. epsilon_difference = 0.000000000]
  174. Notice how, in the above example, two denormalized values that are a factor of 2 apart are
  175. none the less only one representation apart!
  176. */
  177. //] [/compare_floats_example_6]
  178. #if 0
  179. //[old_compare_floats_example_3
  180. //`The simplest use is to compare two values with a tolerance thus:
  181. bool is_close = is_close_to(1.F, 1.F + epsilon, epsilon); // One epsilon apart is close enough.
  182. std::cout << "is_close_to(1.F, 1.F + epsilon, epsilon); is " << is_close << std::endl; // true
  183. is_close = is_close_to(1.F, 1.F + 2 * epsilon, epsilon); // Two epsilon apart isn't close enough.
  184. std::cout << "is_close_to(1.F, 1.F + epsilon, epsilon); is " << is_close << std::endl; // false
  185. /*`
  186. [note The type FPT of the tolerance and the type of the values [*must match].
  187. So `is_close(0.1F, 1., 1.)` will fail to compile because "template parameter 'FPT' is ambiguous".
  188. Always provide the same type, using `static_cast<FPT>` if necessary.]
  189. */
  190. /*`An instance of class `close_at_tolerance` is more convenient
  191. when multiple tests with the same conditions are planned.
  192. A class that stores a tolerance of three epsilon (and the default ['strong] test) is:
  193. */
  194. close_at_tolerance<float> three_rounds(3 * epsilon); // 'strong' by default.
  195. //`and we can confirm these settings:
  196. std::cout << "fraction_tolerance = "
  197. << three_rounds.fraction_tolerance()
  198. << std::endl; // +3.57627869e-007
  199. std::cout << "strength = "
  200. << (three_rounds.strength() == FPC_STRONG ? "strong" : "weak")
  201. << std::endl; // strong
  202. //`To start, let us use two values that are truly equal (having identical bit patterns)
  203. float a = 1.23456789F;
  204. float b = 1.23456789F;
  205. //`and make a comparison using our 3*epsilon `three_rounds` functor:
  206. bool close = three_rounds(a, b);
  207. std::cout << "three_rounds(a, b) = " << close << std::endl; // true
  208. //`Unsurprisingly, the result is true, and the failed fraction is zero.
  209. std::cout << "failed_fraction = " << three_rounds.failed_fraction() << std::endl;
  210. /*`To get some nearby values, it is convenient to use the Boost.Math __next_float functions,
  211. for which we need an include
  212. #include <boost/math/special_functions/next.hpp>
  213. and some using declarations:
  214. */
  215. using boost::math::float_next;
  216. using boost::math::float_prior;
  217. using boost::math::nextafter;
  218. using boost::math::float_distance;
  219. //`To add a few __ulp to one value:
  220. b = float_next(a); // Add just one ULP to a.
  221. b = float_next(b); // Add another one ULP.
  222. b = float_next(b); // Add another one ULP.
  223. // 3 epsilon would pass.
  224. b = float_next(b); // Add another one ULP.
  225. //`and repeat our comparison:
  226. close = three_rounds(a, b);
  227. std::cout << "three_rounds(a, b) = " << close << std::endl; // false
  228. std::cout << "failed_fraction = " << three_rounds.failed_fraction()
  229. << std::endl; // abs(u-v) / abs(v) = 3.86237957e-007
  230. //`We can also 'measure' the number of bits different using the `float_distance` function:
  231. std::cout << "float_distance = " << float_distance(a, b) << std::endl; // 4
  232. /*`Now consider two values that are much further apart
  233. than one might expect from ['computational noise],
  234. perhaps the result of two measurements of some physical property like length
  235. where an uncertainty of a percent or so might be expected.
  236. */
  237. float fp1 = 0.01000F;
  238. float fp2 = 0.01001F; // Slightly different.
  239. float tolerance = 0.0001F;
  240. close_at_tolerance<float> strong(epsilon); // Default is strong.
  241. bool rs = strong(fp1, fp2);
  242. std::cout << "strong(fp1, fp2) is " << rs << std::endl;
  243. //`Or we could contrast using the ['weak] criterion:
  244. close_at_tolerance<float> weak(epsilon, FPC_WEAK); // Explicitly weak.
  245. bool rw = weak(fp1, fp2); //
  246. std::cout << "weak(fp1, fp2) is " << rw << std::endl;
  247. //`We can also construct, setting tolerance and strength, and compare in one statement:
  248. std::cout << a << " #= " << b << " is "
  249. << close_at_tolerance<float>(epsilon, FPC_STRONG)(a, b) << std::endl;
  250. std::cout << a << " ~= " << b << " is "
  251. << close_at_tolerance<float>(epsilon, FPC_WEAK)(a, b) << std::endl;
  252. //`but this has little advantage over using function `is_close_to` directly.
  253. //] [/old_compare_floats_example_3]
  254. /*When the floating-point values become very small and near zero, using
  255. //a relative test becomes unhelpful because one is dividing by zero or tiny,
  256. //Instead, an absolute test is needed, comparing one (or usually both) values with zero,
  257. //using a tolerance.
  258. //This is provided by the `small_with_tolerance` class and `is_small` function.
  259. namespace boost {
  260. namespace math {
  261. namespace fpc {
  262. template<typename FPT>
  263. class small_with_tolerance
  264. {
  265. public:
  266. // Public typedefs.
  267. typedef bool result_type;
  268. // Constructor.
  269. explicit small_with_tolerance(FPT tolerance); // tolerance >= 0
  270. // Functor
  271. bool operator()(FPT value) const; // return true if <= absolute tolerance (near zero).
  272. };
  273. template<typename FPT>
  274. bool
  275. is_small(FPT value, FPT tolerance); // return true if value <= absolute tolerance (near zero).
  276. }}} // namespaces.
  277. /*`
  278. [note The type FPT of the tolerance and the type of the value [*must match].
  279. So `is_small(0.1F, 0.000001)` will fail to compile because "template parameter 'FPT' is ambiguous".
  280. Always provide the same type, using `static_cast<FPT>` if necessary.]
  281. A few values near zero are tested with varying tolerance below.
  282. */
  283. //[compare_floats_small_1
  284. float c = 0;
  285. std::cout << "0 is_small " << is_small(c, epsilon) << std::endl; // true
  286. c = std::numeric_limits<float>::denorm_min(); // 1.40129846e-045
  287. std::cout << "denorm_ min =" << c << ", is_small is " << is_small(c, epsilon) << std::endl; // true
  288. c = (std::numeric_limits<float>::min)(); // 1.17549435e-038
  289. std::cout << "min = " << c << ", is_small is " << is_small(c, epsilon) << std::endl; // true
  290. c = 1 * epsilon; // 1.19209290e-007
  291. std::cout << "epsilon = " << c << ", is_small is " << is_small(c, epsilon) << std::endl; // false
  292. c = 1 * epsilon; // 1.19209290e-007
  293. std::cout << "2 epsilon = " << c << ", is_small is " << is_small(c, 2 * epsilon) << std::endl; // true
  294. c = 2 * epsilon; //2.38418579e-007
  295. std::cout << "4 epsilon = " << c << ", is_small is " << is_small(c, 2 * epsilon) << std::endl; // false
  296. c = 0.00001F;
  297. std::cout << "0.00001 = " << c << ", is_small is " << is_small(c, 0.0001F) << std::endl; // true
  298. c = -0.00001F;
  299. std::cout << "0.00001 = " << c << ", is_small is " << is_small(c, 0.0001F) << std::endl; // true
  300. /*`Using the class `small_with_tolerance` allows storage of the tolerance,
  301. convenient if you make repeated tests with the same tolerance.
  302. */
  303. small_with_tolerance<float>my_test(0.01F);
  304. std::cout << "my_test(0.001F) is " << my_test(0.001F) << std::endl; // true
  305. std::cout << "my_test(0.001F) is " << my_test(0.01F) << std::endl; // false
  306. //] [/compare_floats_small_1]
  307. #endif
  308. return 0;
  309. } // int main()
  310. /*
  311. Example output is:
  312. //[compare_floats_output
  313. Compare floats using Boost.Test functions/classes
  314. float epsilon = 1.19209290e-007
  315. is_close_to(1.F, 1.F + epsilon, epsilon); is true
  316. is_close_to(1.F, 1.F + epsilon, epsilon); is false
  317. fraction_tolerance = 3.57627869e-007
  318. strength = strong
  319. three_rounds(a, b) = true
  320. failed_fraction = 0.000000000
  321. three_rounds(a, b) = false
  322. failed_fraction = 3.86237957e-007
  323. float_distance = 4.00000000
  324. strong(fp1, fp2) is false
  325. weak(fp1, fp2) is false
  326. 1.23456788 #= 1.23456836 is false
  327. 1.23456788 ~= 1.23456836 is false
  328. 0 is_small true
  329. denorm_ min =1.40129846e-045, is_small is true
  330. min = 1.17549435e-038, is_small is true
  331. epsilon = 1.19209290e-007, is_small is false
  332. 2 epsilon = 1.19209290e-007, is_small is true
  333. 4 epsilon = 2.38418579e-007, is_small is false
  334. 0.00001 = 9.99999975e-006, is_small is true
  335. 0.00001 = -9.99999975e-006, is_small is true
  336. my_test(0.001F) is true
  337. my_test(0.001F) is false//] [/compare_floats_output]
  338. */