synaps/arithm/epsilon.h

/********************************************************************
 *   This file is part of the source code of the SYNAPS library.
 *   Author(s): J.P. Pavone, GALAAD, INRIA
 *   $Id: epsilon.H,v 1.1 2005/07/11 09:15:23 mourrain Exp $
 ********************************************************************/
#ifndef SYNAPS_NUMERICS_EPSILON_HPP
#define SYNAPS_NUMERICS_EPSILON_HPP
#include <synaps/init.h>
#include <synaps/arithm/rounding.h>
//--------------------------------------------------------------------
__BEGIN_NAMESPACE_SYNAPS
//--------------------------------------------------------------------
namespace numerics
{
  /* epsilon<T>::result is a numerical zero */
  struct false_t{};
  struct true_t {};
  template < typename T > struct inexact { typedef false_t result_t; };
  template <> struct inexact<long double>{ typedef true_t  result_t; };
  template < typename T > struct epsilon { 
    static T result;  
    static T _init_(const true_t&) {
      volatile T  half(0.5), curr(1.0), one(1.0), prev, eps(1.0);
      do 
        { 
          prev = curr;
          eps *= half;
          curr = one+eps;
        } 
      while ((curr!=1.0)&&(curr!=prev));
      return eps;
    };
    static T _init_(const false_t&) { return T(0); };
    static T _init_() { return _init_(typename inexact< typename LongVersion<T>::result_t >::result_t()); };
  };
  template < typename T > T epsilon<T>::result = epsilon<T>::_init_();
};
//--------------------------------------------------------------------
__END_NAMESPACE_SYNAPS
/********************************************************************/
#endif //