realroot_doc 0.1.1
/Users/mourrain/Devel/mmx/realroot/include/realroot/rounding_mode.hpp
Go to the documentation of this file.
00001 /********************************************************************
00002  *   This file is part of the source code of the realroot library.
00003  *   Author(s): J.P. Pavone, GALAAD, INRIA
00004  *   $Id: fpu.h,v 1.1 2005/07/11 09:15:23 mourrain Exp $
00005  ********************************************************************/
00006 #ifndef mmx_numerics_rounding_mode_hpp
00007 #define mmx_numerics_rounding_mode_hpp
00008 //--------------------------------------------------------------------
00009 #include <fenv.h>
00010 //--------------------------------------------------------------------
00011 namespace mmx {
00012 //--------------------------------------------------------------------
00013 namespace numerics
00014 {
00015   inline int get_cw()  { return fegetround(); };
00016   inline int get_rnd() { return fegetround(); };
00017   inline void set_cw( int cw ) { fesetround(cw);};
00018   inline static int rnd_up() { return FE_UPWARD;};
00019   inline static int rnd_dw() { return FE_DOWNWARD;};
00020   inline static int rnd_nr() { return FE_TONEAREST;};
00021   inline static int rnd_z()  { return FE_TOWARDZERO;};
00022   template< typename T > struct LongVersion { typedef T result_t; };
00023   template<> struct LongVersion<float>  { typedef long double result_t; };
00024   template<> struct LongVersion<double> { typedef long double result_t; };
00025   template< class T > struct fpu_rounding { 
00026     static const int id = 0; 
00027     typedef int rnd_t;
00028     inline static rnd_t  getrnd() { return 0; };
00029     inline static void   setrnd( rnd_t r ) { };
00030     inline static int rnd_up() { return 0; };
00031     inline static int rnd_dw() { return 0; };
00032     inline static int rnd_nr() { return 0; };
00033     inline static int rnd_z()  { return 0; };
00034     fpu_rounding( int ){}; 
00035   };
00036   template<> struct fpu_rounding<long double> 
00037   { 
00038     typedef int rnd_t;
00039     inline static rnd_t  getrnd() { return fegetround(); };
00040     inline static void   setrnd( rnd_t r ) { fesetround(r); };
00041     inline static int rnd_up() { return FE_UPWARD;};
00042     inline static int rnd_dw() { return FE_DOWNWARD;};
00043     inline static int rnd_nr() { return FE_TONEAREST;};
00044     inline static int rnd_z()  { return FE_TOWARDZERO;};
00045     static const int id = 1;
00046     rnd_t m_prev;
00047     bool  m_chg;
00048     fpu_rounding( rnd_t rnd ) 
00049     {
00050       m_chg = false;
00051       m_prev = getrnd();
00052       if ( m_prev != rnd ) 
00053         { 
00054           m_chg = true;
00055           setrnd(rnd); 
00056         };
00057     };
00058     ~fpu_rounding() { 
00059       if ( m_chg ) setrnd(m_prev);
00060     };
00061   };
00062  
00063   
00064   template< typename T >
00065   struct rounding : fpu_rounding< typename LongVersion< T >::result_t >  
00066   { rounding( typename fpu_rounding<typename LongVersion<T>::result_t >::rnd_t rnd ) : fpu_rounding< typename LongVersion< T >::result_t >( rnd ) {}; };
00067   
00068   template< typename T > inline T div_dw(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_dw()); return a/b; };
00069   template< typename T > inline T div_up(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_up()); return a/b; };
00070   template< typename T > inline T add_dw(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_dw()); return a+b; };
00071   template< typename T > inline T add_up(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_up()); return a+b; };
00072   template< typename T > inline T mul_dw(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_dw()); return a*b; };
00073   template< typename T > inline T mul_up(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_up()); return a*b; };
00074   template< typename T > inline T sub_dw(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_dw()); return a-b; };
00075   template< typename T > inline T sub_up(const T& a, const T& b) { rounding<T> r(rounding<T>::rnd_up()); return a-b; };
00076 
00077   namespace rdw
00078   {
00079     template< typename T > inline  
00080     T upadd( const T& a, const T& b ) { T tmp = -a; tmp -= b; return -tmp;  }; 
00081     template< typename T > inline  
00082     T dwadd( const T& a, const T& b ) { return    a+b;    };
00083     template< typename T > inline  
00084     T upsub( const T& a, const T& b ) { T tmp =  b; tmp -= a; return -tmp;  };
00085     template< typename T > inline  
00086     T dwsub( const T& a, const T& b ) { return    a-b;    };
00087     template< typename T > inline  
00088     T upmul( const T& a, const T& b ) { T tmp = -a; tmp *= b; return -tmp;  };
00089     template< typename T > inline  
00090     T dwmul( const T& a, const T& b ) { return    a*b;    };
00091     template< typename T > inline  
00092     T updiv( const T& a, const T& b ) { T tmp = -a; tmp /= b; return -tmp;  };
00093     template< typename T > inline  
00094     T dwdiv( const T& a, const T& b ) { return    a/b;    };
00095     template< typename T > inline  
00096     T updet( const T& a, const T& b, const T& c, const T& d ) 
00097     { return upsub(upmul(a,d),dwmul(b,c)); };
00098     template< typename T > inline
00099     T dwdet( const T& a, const T& b, const T& c, const T& d ) 
00100     { return dwsub(dwmul(a,d),upmul(b,c)); };  
00101   };
00102  
00103   namespace rup
00104   {
00105     template< typename T > inline  
00106     T upadd( const T& a, const T& b ) { return    a+b;    };
00107     template< typename T > inline  
00108     T dwadd( const T& a, const T& b ) { T tmp = -a; tmp -= b; return -tmp;  }; 
00109     template< typename T > inline  
00110     T upsub( const T& a, const T& b ) { return    a-b;    };
00111     template< typename T > inline  
00112     T dwsub( const T& a, const T& b ) { T tmp =  b; tmp -= a; return -tmp;  };
00113     template< typename T > inline  
00114     T upmul( const T& a, const T& b ) { return    a*b;    };
00115     template< typename T > inline  
00116     T dwmul( const T& a, const T& b ) { T tmp = -a; tmp *= b; return -tmp;  };
00117     template< typename T > inline  
00118     T updiv( const T& a, const T& b ) { return    a/b;    };
00119     template< typename T > inline  
00120     T dwdiv( const T& a, const T& b ) { T tmp = -a; tmp /= b; return -tmp;  };
00121     template< typename T > inline  
00122     T updet( const T& a, const T& b, const T& c, const T& d ) 
00123     { return upsub(upmul(a,d),dwmul(b,c)); };
00124     template< typename T > inline
00125     T dwdet( const T& a, const T& b, const T& c, const T& d ) 
00126     { return dwsub(dwmul(a,d),upmul(b,c)); };
00127   }
00128 }; //namespace numerics
00129 //======================================================================
00130 } //namespace mmx
00131 //======================================================================
00132 #endif //