All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
Classes | Functions
Gaudi::Math Namespace Reference

Classes

class  Lomont
 The equality comparison of double numbers using as the metric the maximal number of Units in the Last Place (ULP). More...
 
class  Lomont< double >
 the specialization for double numbers More...
 
class  Lomont< float >
 the specialization for float numbers More...
 

Functions

GAUDI_API bool lomont_compare_float (const float af, const float bf, const unsigned short maxULPs)
 equality comparison of float numbers using as the metric the maximal number of Units in the Last Place (ULP). More...
 
GAUDI_API bool lomont_compare_double (const double af, const double bf, const unsigned int maxULPs)
 equality comparison of double numbers using as the metric the maximal number of Units in the Last Place (ULP). More...
 
GAUDI_API float next_float (const float af, const short ulps)
 Get the floating number that representation is different with respect to the argument for the certain number of "Units in the Last Position". More...
 
GAUDI_API double next_double (const double af, const short ulps)
 Get the floating number that representation is different with respect to the argument for the certain number of "Units in the Last Position". More...
 

Function Documentation

bool Gaudi::Math::lomont_compare_double ( const double  af,
const double  bf,
const unsigned int  maxULPs 
)

equality comparison of double numbers using as the metric the maximal number of Units in the Last Place (ULP).

It is a slightly modified version of very efficient implementation of the initial Bruce Dawson's algorithm by Chris Lomont.

See also
www.lomont.org
http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm

C.Lomont claims the algorithm is factor 2-10 more efficient with respect to classical Knuth's algorithm for comparison of the floating number using the relative precision.

The effective relative difference depends on the choice of maxULPS:

  • For the case of maxULPs=1, (of course it is totally unphysical case!!!) the effective relative precision r = |a-b|/(|a|+|b|)is ~6e-16 for |a|,|b|>1.e-304, and then it quickly goes to ~1
const double a = ... ;
const double b = ... ;
const bool equal = Gaudi::Math::lomont_compare_double ( a , b ) ;
Parameters
afthe first number
bfthe second number
maxULPsthe maximal metric deviation in the terms of maximal number of units in the last place
Author
Vanya BELYAEV Ivan..nosp@m.Bely.nosp@m.aev@n.nosp@m.ikhe.nosp@m.f.nl
Date
2008-11-08

the final check

Definition at line 291 of file Lomont.cpp.

294 {
295  // ==========================================================================
297  BOOST_STATIC_ASSERT( std::numeric_limits<double> ::is_specialized &&
298  std::numeric_limits<Long> ::is_specialized &&
299  std::numeric_limits<ULong> ::is_specialized &&
300  boost::integer_traits<ULong> ::is_specialized &&
301  boost::integer_traits<Long> ::is_specialized &&
302  sizeof(double)==sizeof(Long) &&
303  sizeof(double)==sizeof(ULong) &&
304  64 == std::numeric_limits<ULong>::digits ) ;
305  // ==========================================================================
306  Cast_D caster ;
307 
308  //Long ai = *reinterpret_cast<const Long*>( &af ) ;
309  //Long bi = *reinterpret_cast<const Long*>( &bf ) ;
310  Long ai = caster.d2l ( af ) ;
311  Long bi = caster.d2l ( bf ) ;
312 
313  Long test = (((ULong)(ai^bi))>>63)-1;
314 
315  // assert ( (0==test) || ( boost::integer_traits<ULong>::const_max == test ) ) ;
316 
317  Long diff = ((( boost::integer_traits<Long>::const_min - ai ) & (~test)) | ( ai& test )) - bi ;
318 
319  Long maxDiff_ = maxULPs ;
320 
321  Long v1 = maxDiff_ + diff ;
322  Long v2 = maxDiff_ - diff ;
323 
324  return 0<=(v1|v2) ;
325 }
validate_patch.test
def test
Definition: validate_patch.py:166
bool Gaudi::Math::lomont_compare_float ( const float  af,
const float  bf,
const unsigned short  maxULPs 
)

equality comparison of float numbers using as the metric the maximal number of Units in the Last Place (ULP).

It is a slightly modified version of very efficient implementation of the initial Bruce Dawson's algorithm by Chris Lomont.

See also
www.lomont.org
http://www.cygnus-software.com/papers/comparingfloats/comparingfloats.htm

C.Lomont claims the algorithm is factor 2-10 more efficient with respect to classical Knuth's algorithm for comparison of the floating number using the relative precision.

The effective relative difference depends on the choice of maxULPS:

  • For the case of maxULPs=1, (of course it is totally unphysical case!!!) the effective relative precision r = |a-b|/(|a|+|b|)is between 3.5e-8 and 5.5e-8 for |a|,|b|>1.e-37, and then it quickly goes to ~1
  • For the case of maxULPS=10 the effective relative precision is between 3e-8 and 6e-7 for |a|,|b|>1.e-37, and then it quickly goes to ~1
  • For the case of maxULPS=100 the effective relative precision is around ~6e-6 for |a|,|b|>1.e-37, and then it quickly goes to ~1
  • For the case of maxULPS=1000 the effective relative precision is around ~6e-5 for |a|,|b|>1.e-37, and then it quickly goes to ~1
const float a = ... ;
const float b = ... ;
const bool equal = Gaudi::Math::lomont_compare_float ( a , b ) ;
Parameters
afthe first number
bfthe second number
maxULPsthe maximal metric deviation in the terms of maximal number of units in the last place
Author
Vanya BELYAEV Ivan..nosp@m.Bely.nosp@m.aev@n.nosp@m.ikhe.nosp@m.f.nl
Date
2008-11-08

Definition at line 195 of file Lomont.cpp.

198 {
199  // ==========================================================================
200  // prerequisites:
201  BOOST_STATIC_ASSERT( std::numeric_limits<float> ::is_specialized &&
202  boost::integer_traits<int> ::is_specialized &&
203  boost::integer_traits<unsigned int> ::is_specialized &&
204  sizeof(float)==sizeof(int) &&
205  sizeof(float)==sizeof(unsigned int) &&
206  32 == boost::integer_traits<unsigned int>::digits ) ;
207  // ==========================================================================
208 
209  Cast_F caster ;
210 
211  //int ai = *reinterpret_cast<const int*>( &af ) ;
212  //int bi = *reinterpret_cast<const int*>( &bf ) ;
213  int ai = caster.f2i ( af ) ;
214  int bi = caster.f2i ( bf ) ;
215 
216  int test = (((unsigned int)(ai^bi))>>31)-1;
217 
218  // assert ( (0==test) || ( boost::integer_traits<unsigned int>::const_max == test ) ) ;
219 
220  int diff = ((( boost::integer_traits<int>::const_min - ai ) & (~test)) | ( ai& test )) - bi ;
221 
222  int maxDiff_ = maxULPs ;
223 
224  int v1 = maxDiff_ + diff ;
225  int v2 = maxDiff_ - diff ;
226 
227  return 0<=(v1|v2) ;
228 }
validate_patch.test
def test
Definition: validate_patch.py:166
double Gaudi::Math::next_double ( const double  af,
const short  ulps 
)

Get the floating number that representation is different with respect to the argument for the certain number of "Units in the Last Position".

For ulps=1, it is just next float number, for ulps=-1 is is the previous one.

This routine is very convenient to test the parameter maxULPS for the routine Gaudi::Math::lomont_compare_double

See also
Gaudi::Math::lomont_compare_double
Parameters
afthe reference number
ulpsthe bias
Returns
the biased float number (on distance "ulps")
Author
Vanya BELYAEV Ivan..nosp@m.Bely.nosp@m.aev@n.nosp@m.ikhe.nosp@m.f.nl
Date
2008-11-08

the final check

Definition at line 344 of file Lomont.cpp.

345 {
346  // ==========================================================================
348  BOOST_STATIC_ASSERT( std::numeric_limits<double> ::is_specialized &&
349  std::numeric_limits<Long> ::is_specialized &&
350  sizeof(double)==sizeof(Long) ) ;
351  // ==========================================================================
352  Cast_D caster ;
353 
354  //Long al = *reinterpret_cast<const Long*>( &ad ) ;
355  Long al = caster.d2l ( ad ) ;
356  al += ulps ;
357  //return *reinterpret_cast<double*>(&al) ;
358  return caster.l2d ( al ) ;
359 }
float Gaudi::Math::next_float ( const float  af,
const short  ulps 
)

Get the floating number that representation is different with respect to the argument for the certain number of "Units in the Last Position".

For ulps=1, it is just next float number, for ulps=-1 is is the previous one.

This routine is very convenient to test the parameter maxULPS for the routine Gaudi::Math::lomont_compare_float

See also
Gaudi:Math::lomont_compare_float
Parameters
afthe reference number
ulpsthe bias
Returns
the biased float number (on distance "ulps")
Author
Vanya BELYAEV Ivan..nosp@m.Bely.nosp@m.aev@n.nosp@m.ikhe.nosp@m.f.nl
Date
2008-11-08

the final check

Definition at line 247 of file Lomont.cpp.

248 {
250  BOOST_STATIC_ASSERT( std::numeric_limits<float> ::is_specialized &&
251  std::numeric_limits<int> ::is_specialized &&
252  sizeof(float)==sizeof(int) ) ;
253  // ==========================================================================
254  Cast_F caster ;
255 
256  // int ai = *reinterpret_cast<const int*>( &af ) ;
257  int ai = caster.f2i ( af ) ;
258  ai += ulps ;
259  // return *reinterpret_cast<float*>(&ai) ;
260  return caster.i2f ( ai ) ;
261 }