Elements 6.0.1
A C++ base framework for the Euclid Software.
Loading...
Searching...
No Matches
Real.h
Go to the documentation of this file.
1
40// Copyright 2005, Google Inc.
41// All rights reserved.
42//
43// Redistribution and use in source and binary forms, with or without
44// modification, are permitted provided that the following conditions are
45// met:
46//
47// * Redistributions of source code must retain the above copyright
48// notice, this list of conditions and the following disclaimer.
49// * Redistributions in binary form must reproduce the above
50// copyright notice, this list of conditions and the following disclaimer
51// in the documentation and/or other materials provided with the
52// distribution.
53// * Neither the name of Google Inc. nor the names of its
54// contributors may be used to endorse or promote products derived from
55// this software without specific prior written permission.
56//
57// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
58// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
59// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
60// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
61// OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
62// SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
63// LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
64// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
65// THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
66// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
67// OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
68//
69// Authors: wan@google.com (Zhanyong Wan), eefacm@gmail.com (Sean Mcafee)
70//
71// The Google C++ Testing Framework (Google Test)
72#ifndef ELEMENTSKERNEL_ELEMENTSKERNEL_REAL_H_
73#define ELEMENTSKERNEL_ELEMENTSKERNEL_REAL_H_
74
75#include <cmath> // for round
76#include <limits> // for numeric_limits
77#include <type_traits> // for is_floating_point
78
79#include "ElementsKernel/Export.h" // ELEMENTS_API
80#include "ElementsKernel/Unused.h" // ELEMENTS_UNUSED
81
83
84namespace Elements {
85
90
91// For testing purposes only. Rather use the isEqual functions for real
92// life comparison
97
98template <std::size_t size>
100public:
101 // This prevents the user from using TypeWithSize<N> with incorrect
102 // values of N.
103 using UInt = void;
104};
105
106// The specialisation for size 4.
107template <>
109public:
110 // unsigned int has size 4 in both gcc and MSVC.
111 //
112 // As base/basictypes.h doesn't compile on Windows, we cannot use
113 // uint32, uint64, and etc here.
114 using Int = int;
115 using UInt = unsigned int;
116};
117
118// The specialisation for size 8.
119template <>
121public:
122 using Int = long long; // NOLINT
123 using UInt = unsigned long long; // NOLINT
124};
125
126template <typename RawType>
129}
130
131template <>
134}
135
136template <>
139}
140
141// This template class represents an IEEE floating-point number
142// (either single-precision or double-precision, depending on the
143// template parameters).
144//
145// The purpose of this class is to do more sophisticated number
146// comparison. (Due to round-off error, etc, it's very unlikely that
147// two floating-points will be equal exactly. Hence a naive
148// comparison by the == operation often doesn't work.)
149//
150// Format of IEEE floating-point:
151//
152// The most-significant bit being the leftmost, an IEEE
153// floating-point looks like
154//
155// sign_bit exponent_bits fraction_bits
156//
157// Here, sign_bit is a single bit that designates the sign of the
158// number.
159//
160// For float, there are 8 exponent bits and 23 fraction bits.
161//
162// For double, there are 11 exponent bits and 52 fraction bits.
163//
164// More details can be found at
165// http://en.wikipedia.org/wiki/IEEE_floating-point_standard.
166//
167// Template parameter:
168//
169// RawType: the raw floating-point type (either float or double)
170template <typename RawType>
172public:
173 // Defines the unsigned integer type that has the same size as the
174 // floating point number.
175 using Bits = typename TypeWithSize<sizeof(RawType)>::UInt;
176
177 // Constants.
178
179 // # of bits in a number.
180 static const std::size_t s_bitcount = 8 * sizeof(RawType);
181
182 // # of fraction bits in a number.
183 static const std::size_t s_fraction_bitcount = std::numeric_limits<RawType>::digits - 1;
184
185 // # of exponent bits in a number.
186 static const std::size_t s_exponent_bitcount = s_bitcount - 1 - s_fraction_bitcount;
187
188 // The mask for the sign bit.
189 static const Bits s_sign_bitmask = static_cast<Bits>(1) << (s_bitcount - 1);
190
191 // The mask for the fraction bits.
192 static const Bits s_fraction_bitmask = ~static_cast<Bits>(0) >> (s_exponent_bitcount + 1);
193
194 // The mask for the exponent bits.
195 static const Bits s_exponent_bitmask = ~(s_sign_bitmask | s_fraction_bitmask);
196
197 // How many ULP's (Units in the Last Place) we want to tolerate when
198 // comparing two numbers. The larger the value, the more error we
199 // allow. A 0 value means that two numbers must be exactly the same
200 // to be considered equal.
201 //
202 // The maximum error of a single floating-point operation is 0.5
203 // units in the last place. On Intel CPU's, all floating-point
204 // calculations are done with 80-bit precision, while double has 64
205 // bits. Therefore, 4 should be enough for ordinary use.
206 //
207 // See the following article for more details on ULP:
208 // http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm.
209 static const std::size_t m_max_ulps = defaultMaxUlps<RawType>();
210
211 // Constructs a FloatingPoint from a raw floating-point number.
212 //
213 // On an Intel CPU, passing a non-normalised NAN (Not a Number)
214 // around may change its bits, although the new value is guaranteed
215 // to be also a NAN. Therefore, don't expect this constructor to
216 // preserve the bits in x when x is a NAN.
217 explicit FloatingPoint(const RawType& x) {
218 m_u.m_value = x;
219 }
220
221 // Static methods
222
223 // Reinterprets a bit pattern as a floating-point number.
224 //
225 // This function is needed to test the AlmostEquals() method.
226 static RawType ReinterpretBits(const Bits& bits) {
227 FloatingPoint fp(0);
228 fp.m_u.m_bits = bits;
229 return fp.m_u.m_value;
230 }
231
232 // Returns the floating-point number that represent positive infinity.
233 static RawType Infinity() {
234 return ReinterpretBits(s_exponent_bitmask);
235 }
236
237 // Non-static methods
238
239 // Returns the bits that represents this number.
240 const Bits& bits() const {
241 return m_u.m_bits;
242 }
243
244 // Returns the exponent bits of this number.
246 return s_exponent_bitmask & m_u.m_bits;
247 }
248
249 // Returns the fraction bits of this number.
251 return s_fraction_bitmask & m_u.m_bits;
252 }
253
254 // Returns the sign bit of this number.
255 Bits signBit() const {
256 return s_sign_bitmask & m_u.m_bits;
257 }
258
259 // Returns true iff this is NAN (not a number).
260 bool isNan() const {
261 // It's a NAN if the exponent bits are all ones and the fraction
262 // bits are not entirely zeros.
263 return (exponentBits() == s_exponent_bitmask) && (fractionBits() != 0);
264 }
265
266 // Returns true iff this number is at most kMaxUlps ULP's away from
267 // rhs. In particular, this function:
268 //
269 // - returns false if either number is (or both are) NAN.
270 // - treats really large numbers as almost equal to infinity.
271 // - thinks +0.0 and -0.0 are 0 DLP's apart.
272 bool AlmostEquals(const FloatingPoint& rhs) const {
273 // The IEEE standard says that any comparison operation involving
274 // a NAN must return false.
275 if (isNan() || rhs.isNan()) {
276 return false;
277 }
278 return distanceBetweenSignAndMagnitudeNumbers(m_u.m_bits, rhs.m_u.m_bits) <= m_max_ulps;
279 }
280
281 // Converts an integer from the sign-and-magnitude representation to
282 // the biased representation. More precisely, let N be 2 to the
283 // power of (kBitCount - 1), an integer x is represented by the
284 // unsigned number x + N.
285 //
286 // For instance,
287 //
288 // -N + 1 (the most negative number representable using
289 // sign-and-magnitude) is represented by 1;
290 // 0 is represented by N; and
291 // N - 1 (the biggest number representable using
292 // sign-and-magnitude) is represented by 2N - 1.
293 //
294 // Read http://en.wikipedia.org/wiki/Signed_number_representations
295 // for more details on signed number representations.
296 static Bits signAndMagnitudeToBiased(const Bits& sam) {
297 if (s_sign_bitmask & sam) {
298 // sam represents a negative number.
299 return ~sam + 1;
300 } else {
301 // sam represents a positive number.
302 return s_sign_bitmask | sam;
303 }
304 }
305
306 // Given two numbers in the sign-and-magnitude representation,
307 // returns the distance between them as an unsigned number.
308 static Bits distanceBetweenSignAndMagnitudeNumbers(const Bits& sam1, const Bits& sam2) {
309 const Bits biased1 = signAndMagnitudeToBiased(sam1);
310 const Bits biased2 = signAndMagnitudeToBiased(sam2);
311 return (biased1 >= biased2) ? (biased1 - biased2) : (biased2 - biased1);
312 }
313
314private:
315 // The data type used to store the actual floating-point number.
317 RawType m_value; // The raw floating-point number.
318 Bits m_bits; // The bits that represent the number.
319 };
320
322};
323
324// Usable AlmostEqual function
325
326template <typename FloatType>
327bool almostEqual2sComplement(ELEMENTS_UNUSED const FloatType& a, ELEMENTS_UNUSED const FloatType& b,
328 ELEMENTS_UNUSED const std::size_t& max_ulps = 0) {
329 return false;
330}
331
332template <typename RawType>
333bool isNan(const RawType& x) {
334
335 using Bits = typename TypeWithSize<sizeof(RawType)>::UInt;
336 Bits x_bits = *reinterpret_cast<const Bits*>(&x);
337
338 Bits x_exp_bits = FloatingPoint<RawType>::s_exponent_bitmask & x_bits;
339 Bits x_frac_bits = FloatingPoint<RawType>::s_fraction_bitmask & x_bits;
340
341 return (x_exp_bits == FloatingPoint<RawType>::s_exponent_bitmask) && (x_frac_bits != 0);
342}
343
344template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
345bool isEqual(const RawType& left, const RawType& right) {
346
347 bool is_equal{false};
348
349 if (not(isNan<RawType>(left) or isNan<RawType>(right))) {
350 using Bits = typename TypeWithSize<sizeof(RawType)>::UInt;
351 Bits l_bits = *reinterpret_cast<const Bits*>(&left);
352 Bits r_bits = *reinterpret_cast<const Bits*>(&right);
353 is_equal = (FloatingPoint<RawType>::distanceBetweenSignAndMagnitudeNumbers(l_bits, r_bits) <= max_ulps);
354 }
355
356 return is_equal;
357}
358
359template <std::size_t max_ulps>
360inline bool isEqual(const float& left, const float& right) {
361 return (isEqual<float, max_ulps>(left, right));
362}
363
364template <std::size_t max_ulps>
365inline bool isEqual(const double& left, const double& right) {
366 return (isEqual<double, max_ulps>(left, right));
367}
368
369template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
370inline bool isNotEqual(const RawType& left, const RawType& right) {
371 return (not isEqual<RawType, max_ulps>(left, right));
372}
373
374template <std::size_t max_ulps>
375inline bool isNotEqual(const float& left, const float& right) {
376 return (isNotEqual<float, max_ulps>(left, right));
377}
378
379template <std::size_t max_ulps>
380inline bool isNotEqual(const double& left, const double& right) {
381 return (isNotEqual<double, max_ulps>(left, right));
382}
383
384template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
385bool isLess(const RawType& left, const RawType& right) {
386 bool is_less{false};
387
388 if (left < right && (not isEqual<RawType, max_ulps>(left, right))) {
389 is_less = true;
390 }
391
392 return is_less;
393}
394
395template <std::size_t max_ulps>
396inline bool isLess(const float& left, const float& right) {
397 return (isLess<float, max_ulps>(left, right));
398}
399
400template <std::size_t max_ulps>
401inline bool isLess(const double& left, const double& right) {
402 return (isLess<double, max_ulps>(left, right));
403}
404
405template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
406bool isGreater(const RawType& left, const RawType& right) {
407 bool is_greater{false};
408
409 if (left > right && (not isEqual<RawType, max_ulps>(left, right))) {
410 is_greater = true;
411 }
412
413 return is_greater;
414}
415
416template <std::size_t max_ulps>
417inline bool isGreater(const float& left, const float& right) {
418 return (isGreater<float, max_ulps>(left, right));
419}
420
421template <std::size_t max_ulps>
422inline bool isGreater(const double& left, const double& right) {
423 return (isGreater<double, max_ulps>(left, right));
424}
425
426template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
427bool isLessOrEqual(const RawType& left, const RawType& right) {
428 bool is_loe{false};
429
430 if (not isGreater<RawType, max_ulps>(left, right)) {
431 is_loe = true;
432 }
433
434 return is_loe;
435}
436
437template <std::size_t max_ulps>
438inline bool isLessOrEqual(const float& left, const float& right) {
439 return (isLessOrEqual<float, max_ulps>(left, right));
440}
441
442template <std::size_t max_ulps>
443inline bool isLessOrEqual(const double& left, const double& right) {
444 return (isLessOrEqual<double, max_ulps>(left, right));
445}
446
447template <typename RawType, std::size_t max_ulps = defaultMaxUlps<RawType>()>
448bool isGreaterOrEqual(const RawType& left, const RawType& right) {
449 bool is_goe{false};
450
451 if (not isLess<RawType, max_ulps>(left, right)) {
452 is_goe = true;
453 }
454
455 return is_goe;
456}
457
458template <std::size_t max_ulps>
459inline bool isGreaterOrEqual(const float& left, const float& right) {
460 return (isGreaterOrEqual<float, max_ulps>(left, right));
461}
462
463template <std::size_t max_ulps>
464inline bool isGreaterOrEqual(const double& left, const double& right) {
465 return (isGreaterOrEqual<double, max_ulps>(left, right));
466}
467
484ELEMENTS_API bool almostEqual2sComplement(const float& left, const float& right,
485 const int& max_ulps = FLT_DEFAULT_MAX_ULPS);
502ELEMENTS_API bool almostEqual2sComplement(const double& left, const double& right,
503 const int& max_ulps = DBL_DEFAULT_MAX_ULPS);
504
518template <typename RawType>
519ELEMENTS_API bool realBitWiseEqual(const RawType& left, const RawType& right) {
520#pragma GCC diagnostic push
521#pragma GCC diagnostic ignored "-Wfloat-equal"
522 return (left == right);
523#pragma GCC diagnostic pop
524}
525
526} // namespace Elements
527
528#endif // ELEMENTSKERNEL_ELEMENTSKERNEL_REAL_H_
529
defines the macros to be used for explicit export of the symbols
Macro to silence unused variables warnings from the compiler.
bool isNan() const
Definition: Real.h:260
FloatingPointUnion m_u
Definition: Real.h:321
Bits fractionBits() const
Definition: Real.h:250
typename TypeWithSize< sizeof(RawType)>::UInt Bits
Definition: Real.h:175
Bits signBit() const
Definition: Real.h:255
static RawType ReinterpretBits(const Bits &bits)
Definition: Real.h:226
FloatingPoint(const RawType &x)
Definition: Real.h:217
static Bits distanceBetweenSignAndMagnitudeNumbers(const Bits &sam1, const Bits &sam2)
Definition: Real.h:308
static Bits signAndMagnitudeToBiased(const Bits &sam)
Definition: Real.h:296
static RawType Infinity()
Definition: Real.h:233
bool AlmostEquals(const FloatingPoint &rhs) const
Definition: Real.h:272
Bits exponentBits() const
Definition: Real.h:245
const Bits & bits() const
Definition: Real.h:240
unsigned int UInt
Definition: Real.h:115
unsigned long long UInt
Definition: Real.h:123
#define ELEMENTS_API
Dummy definitions for the backward compatibility mode.
Definition: Export.h:74
#define ELEMENTS_UNUSED
Definition: Unused.h:39
constexpr std::size_t defaultMaxUlps< float >()
Definition: Real.h:132
bool almostEqual2sComplement(ELEMENTS_UNUSED const FloatType &a, ELEMENTS_UNUSED const FloatType &b, ELEMENTS_UNUSED const std::size_t &max_ulps=0)
Definition: Real.h:327
ELEMENTS_API const double FLT_DEFAULT_TEST_TOLERANCE
Single precision float default test tolerance.
Definition: Real.cpp:35
bool isLess(const RawType &left, const RawType &right)
Definition: Real.h:385
constexpr std::size_t defaultMaxUlps()
Definition: Real.h:127
constexpr std::size_t FLT_DEFAULT_MAX_ULPS
Single precision float default maximum unit in the last place.
Definition: Real.h:87
ELEMENTS_API const double DBL_DEFAULT_TEST_TOLERANCE
Double precision float default test tolerance.
Definition: Real.cpp:36
bool isNan(const RawType &x)
Definition: Real.h:333
constexpr std::size_t DBL_DEFAULT_MAX_ULPS
Double precision float default maximum unit in the last place.
Definition: Real.h:89
ELEMENTS_API bool realBitWiseEqual(const RawType &left, const RawType &right)
This function compares 2 floating point numbers bitwise. These are the strict equivalent of the "=="....
Definition: Real.h:519
bool isGreaterOrEqual(const RawType &left, const RawType &right)
Definition: Real.h:448
bool isLessOrEqual(const RawType &left, const RawType &right)
Definition: Real.h:427
bool isNotEqual(const RawType &left, const RawType &right)
Definition: Real.h:370
bool isGreater(const RawType &left, const RawType &right)
Definition: Real.h:406
constexpr std::size_t defaultMaxUlps< double >()
Definition: Real.h:137
bool isEqual(const RawType &left, const RawType &right)
Definition: Real.h:345