Effective implementation of algorithms (Master Thesis)
Effective and error-free implementation of algorithms
|
00001 #ifndef H_MATH_POWERMOD_MULTMOD_EXTENDED 00002 #define H_MATH_POWERMOD_MULTMOD_EXTENDED 00003 00007 #include "utils/preconditions/preconditions.h" 00008 #include "utils/static_assert/static_assert.h" 00009 #include <limits.h> 00010 #include "utils/macros/unused.h" 00011 00012 namespace math { 00013 namespace powermod { 00014 00024 template <int shift> 00025 class MultmodExtended { 00026 STATIC_ASSERT(shift > 0, ""); 00027 00028 public: 00029 template <typename T> 00030 static T max_argument(T UNUSED(x)) { 00031 return std::numeric_limits<T>::max() >> (shift); 00032 } 00036 template <typename T> 00037 static T multmod(T a, T b, T modulo) { 00038 STATIC_ASSERT(std::numeric_limits<T>::is_integer, "") 00039 // TODO(check this) 00040 STATIC_ASSERT(shift * 2 < std::numeric_limits<T>::digits, 00041 ""); 00042 const T max = max_argument((T) 0); 00043 Preconditions::checkRange(modulo, (T) 1, max); 00044 Preconditions::checkRange(a, modulo); 00045 Preconditions::checkRange(b, modulo); 00046 const T step = ((T)1) << shift; 00047 const T mask = step - 1; 00048 00049 T result = 0; 00050 while (b != 0) { 00051 result = (result + a * (b & mask)) % modulo; 00052 a = (a << shift) % modulo; 00053 b >>= shift; 00054 } 00055 00056 return result; 00057 } 00058 }; 00059 00067 class MultmodExtendedOpt { 00068 public: 00069 template <typename T> 00070 static T max_argument(T UNUSED(x)) { 00071 return std::numeric_limits<T>::max() >> (1); 00072 } 00073 00074 template<typename T> 00075 static T multmod(T a, T b, T modulo) { 00076 // We will be using shift-by-one strategy 00077 const T max = max_argument((T) 0); 00078 Preconditions::checkRange(modulo, (T) 1, max); 00079 Preconditions::checkRange(a, modulo); 00080 Preconditions::checkRange(b, modulo); 00081 T result = 0; 00082 while (b != 0) { 00083 if ((b & 1) == 1) { 00084 result += a; 00085 } 00086 b >>= 1; 00087 a <<= 1; 00088 // Note that both result and a may be at most 00089 // 2*modulo-1 and we may optimize the slow remainder calculation 00090 if (result > modulo) result -= modulo; 00091 if (a > modulo) a -= modulo; 00092 } 00093 00094 return result; 00095 } 00096 }; 00097 00098 } // namespace powermod 00099 } // namespace math 00100 #endif