float_next.qbk 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292
  1. [section:next_float Floating-Point Representation Distance (ULP),
  2. and Finding Adjacent Floating-Point Values]
  3. [@http://en.wikipedia.org/wiki/Unit_in_the_last_place Unit of Least Precision or Unit in the Last Place]
  4. is the gap between two different, but as close as possible, floating-point numbers.
  5. Most decimal values, for example 0.1, cannot be exactly represented as floating-point values,
  6. but will be stored as the
  7. [@http://en.wikipedia.org/wiki/Floating_point#Representable_numbers.2C_conversion_and_rounding
  8. closest representable floating-point].
  9. Functions are provided for finding adjacent greater and lesser floating-point values,
  10. and estimating the number of gaps between any two floating-point values.
  11. The floating-point type (FPT) must have has a fixed number of bits in the representation.
  12. The number of bits may set at runtime, but must be the same for all numbers.
  13. For example, __NTL_quad_float type (fixed 128-bit representation),
  14. __NTL_RR type (arbitrary but fixed decimal digits, default 150) or
  15. __multiprecision __cpp_dec_float and__cpp_bin_float are fixed at runtime,
  16. but [*not] a type that extends the representation to provide an exact representation
  17. for any number, for example [@http://keithbriggs.info/xrc.html XRC eXact Real in C].
  18. The accuracy of mathematical functions can be assessed and displayed in terms of __ULP,
  19. often as a ulps plot or by binning the differences as a histogram.
  20. Samples are evaluated using the implementation under test and compared with 'known good'
  21. representation obtained using a more accurate method. Other implementations, often using
  22. arbitrary precision arithmetic, for example __WolframAlpha are one source of references
  23. values. The other method, used widely in Boost.Math special functions, it to carry out
  24. the same algorithm, but using a higher precision type, typically using Boost.Multiprecision
  25. types like `cpp_bin_float_quad` for 128-bit (about 35 decimal digit precision), or
  26. `cpp_bin_float_50` (for 50 decimal digit precision).
  27. When converted to a particular machine representation, say `double`, say using a `static_cast`,
  28. the value is the nearest representation possible for the `double` type. This value
  29. cannot be 'wrong' by more than half a __ulp, and can be obtained using the Boost.Math function `ulp`.
  30. (Unless the algorithm is fundamentally flawed, something that should be revealed by 'sanity'
  31. checks using some independent sources).
  32. See some discussion and example plots by Cleve Moler of Mathworks
  33. [@https://blogs.mathworks.com/cleve/2017/01/23/ulps-plots-reveal-math-function-accurary/
  34. ulps plots reveal math-function accuracy].
  35. [section:nextafter Finding the Next Representable Value in a Specific Direction (nextafter)]
  36. [h4 Synopsis]
  37. ``
  38. #include <boost/math/special_functions/next.hpp>
  39. ``
  40. namespace boost{ namespace math{
  41. template <class FPT>
  42. FPT nextafter(FPT val, FPT direction);
  43. }} // namespaces
  44. [h4 Description - nextafter]
  45. This is an implementation of the `nextafter` function included in the C99 standard.
  46. (It is also effectively an implementation of the C99 `nexttoward` legacy function
  47. which differs only having a `long double` direction,
  48. and can generally serve in its place if required).
  49. [note The C99 functions must use suffixes f and l to distinguish `float` and `long double` versions.
  50. C++ uses the template mechanism instead.]
  51. Returns the next representable value after /x/ in the direction of /y/. If
  52. `x == y` then returns /x/. If /x/ is non-finite then returns the result of
  53. a __domain_error. If there is no such value in the direction of /y/ then
  54. returns an __overflow_error.
  55. [warning The template parameter FTP must be a floating-point type.
  56. An integer type, for example, will produce an unhelpful error message.]
  57. [tip Nearly always, you just want the next or prior representable value,
  58. so instead use `float_next` or `float_prior` below.]
  59. [h4 Examples - nextafter]
  60. The two representations using a 32-bit float either side of unity are:
  61. ``
  62. The nearest (exact) representation of 1.F is 1.00000000
  63. nextafter(1.F, 999) is 1.00000012
  64. nextafter(1/f, -999) is 0.99999994
  65. The nearest (not exact) representation of 0.1F is 0.100000001
  66. nextafter(0.1F, 10) is 0.100000009
  67. nextafter(0.1F, 10) is 0.099999994
  68. ``
  69. [endsect] [/section:nextafter Finding the Next Representable Value in a Specific Direction (nextafter)]
  70. [section:float_next Finding the Next Greater Representable Value (float_next)]
  71. [h4 Synopsis]
  72. ``
  73. #include <boost/math/special_functions/next.hpp>
  74. ``
  75. namespace boost{ namespace math{
  76. template <class FPT>
  77. FPT float_next(FPT val);
  78. }} // namespaces
  79. [h4 Description - float_next]
  80. Returns the next representable value which is greater than /x/.
  81. If /x/ is non-finite then returns the result of
  82. a __domain_error. If there is no such value greater than /x/ then
  83. returns an __overflow_error.
  84. Has the same effect as
  85. nextafter(val, (std::numeric_limits<FPT>::max)());
  86. [endsect] [/section:float_next Finding the Next Greater Representable Value (float_prior)]
  87. [section:float_prior Finding the Next Smaller Representable Value (float_prior)]
  88. [h4 Synopsis]
  89. ``
  90. #include <boost/math/special_functions/next.hpp>
  91. ``
  92. namespace boost{ namespace math{
  93. template <class FPT>
  94. FPT float_prior(FPT val);
  95. }} // namespaces
  96. [h4 Description - float_prior]
  97. Returns the next representable value which is less than /x/.
  98. If /x/ is non-finite then returns the result of
  99. a __domain_error. If there is no such value less than /x/ then
  100. returns an __overflow_error.
  101. Has the same effect as
  102. nextafter(val, -(std::numeric_limits<FPT>::max)()); // Note most negative value -max.
  103. [endsect] [/section:float_prior Finding the Next Smaller Representable Value (float_prior)]
  104. [section:float_distance Calculating the Representation Distance
  105. Between Two floating-point Values (ULP) float_distance]
  106. Function float_distance finds the number of gaps/bits/ULP between any two floating-point values.
  107. If the significands of floating-point numbers are viewed as integers,
  108. then their difference is the number of ULP/gaps/bits different.
  109. [h4 Synopsis]
  110. ``
  111. #include <boost/math/special_functions/next.hpp>
  112. ``
  113. namespace boost{ namespace math{
  114. template <class FPT>
  115. FPT float_distance(FPT a, FPT b);
  116. }} // namespaces
  117. [h4 Description - float_distance]
  118. Returns the distance between /a/ and /b/: the result is always
  119. a signed integer value (stored in floating-point type FPT)
  120. representing the number of distinct representations between /a/ and /b/.
  121. Note that
  122. * `float_distance(a, a)` always returns 0.
  123. * `float_distance(float_next(a), a)` always returns -1.
  124. * `float_distance(float_prior(a), a)` always returns 1.
  125. The function `float_distance` is equivalent to calculating the number
  126. of ULP (Units in the Last Place) between /a/ and /b/ except that it
  127. returns a signed value indicating whether `a > b` or not.
  128. If the distance is too great then it may not be able
  129. to be represented as an exact integer by type FPT,
  130. but in practice this is unlikely to be a issue.
  131. [endsect] [/section:float_distance Calculating the Representation Distance
  132. Between Two floating-point Values (ULP) float_distance]
  133. [section:float_advance Advancing a floating-point Value by a Specific
  134. Representation Distance (ULP) float_advance]
  135. Function `float_advance` advances a floating-point number by a specified number
  136. of ULP.
  137. [h4 Synopsis]
  138. ``
  139. #include <boost/math/special_functions/next.hpp>
  140. ``
  141. namespace boost{ namespace math{
  142. template <class FPT>
  143. FPT float_advance(FPT val, int distance);
  144. }} // namespaces
  145. [h4 Description - float_advance]
  146. Returns a floating-point number /r/ such that `float_distance(val, r) == distance`.
  147. [endsect] [/section:float_advance]
  148. [section:ulp Obtaining the Size of a Unit In the Last Place - ULP]
  149. Function `ulp` gives the size of a unit-in-the-last-place for a specified floating-point value.
  150. [h4 Synopsis]
  151. ``
  152. #include <boost/math/special_functions/ulp.hpp>
  153. ``
  154. namespace boost{ namespace math{
  155. template <class FPT>
  156. FPT ulp(const FPT& x);
  157. template <class FPT, class Policy>
  158. FPT ulp(const FPT& x, const Policy&);
  159. }} // namespaces
  160. [h4 Description - ulp]
  161. Returns one [@http://en.wikipedia.org/wiki/Unit_in_the_last_place unit in the last place] of ['x].
  162. Corner cases are handled as follows:
  163. * If the argument is a NaN, then raises a __domain_error.
  164. * If the argument is an infinity, then raises an __overflow_error.
  165. * If the argument is zero then returns the smallest representable value: for example for type
  166. `double` this would be either `std::numeric_limits<double>::min()` or `std::numeric_limits<double>::denorm_min()`
  167. depending whether denormals are supported (which have the values `2.2250738585072014e-308` and `4.9406564584124654e-324` respectively).
  168. * If the result is too small to represent, then returns the smallest representable value.
  169. * Always returns a positive value such that `ulp(x) == ulp(-x)`.
  170. [*Important:] The behavior of this function is aligned to that of [@http://docs.oracle.com/javase/7/docs/api/java/lang/Math.html#ulp%28double%29
  171. Java's ulp function], please note
  172. however that this function should only ever be used for rough and ready calculations as there are enough
  173. corner cases to trap even careful programmers. In particular:
  174. * The function is asymetrical, which is to say, given `u = ulp(x)` if `x > 0` then `x + u` is the
  175. next floating-point value, but `x - u` is not necessarily the previous value. Similarly, if
  176. `x < 0` then `x - u` is the previous floating-point value, but `x + u` is not necessarily the next
  177. value. The corner cases occur at power of 2 boundaries.
  178. * When the argument becomes very small, it may be that there is no floating-point value that
  179. represents one ULP. Whether this is the case or not depends not only on whether the hardware
  180. may ['sometimes] support denormals (as signalled by `std::numeric_limits<FPT>::has_denorm`), but also whether these are
  181. currently enabled at runtime (for example on SSE hardware, the DAZ or FTZ flags will disable denormal support).
  182. In this situation, the `ulp` function may return a value that is many orders of magnitude too large.
  183. In light of the issues above, we recomend that:
  184. * To move between adjacent floating-point values always use __float_next, __float_prior or __nextafter (`std::nextafter`
  185. is another candidate, but our experience is that this also often breaks depending which optimizations and
  186. hardware flags are in effect).
  187. * To move several floating-point values away use __float_advance.
  188. * To calculate the edit distance between two floats use __float_distance.
  189. There is none the less, one important use case for this function:
  190. If it is known that the true result of some function is x[sub t] and the calculated result
  191. is x[sub c], then the error measured in ulp is simply [^fabs(x[sub t] - x[sub c]) / ulp(x[sub t])].
  192. [endsect] [/section ulp]
  193. [endsect] [/ section:next_float Floating-Point Representation Distance (ULP),
  194. and Finding Adjacent Floating-Point Values]
  195. [/
  196. Copyright 2008 John Maddock and Paul A. Bristow.
  197. Distributed under the Boost Software License, Version 1.0.
  198. (See accompanying file LICENSE_1_0.txt or copy at
  199. http://www.boost.org/LICENSE_1_0.txt).
  200. ]