realroot_doc 0.1.1
|
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 //