Effective implementation of algorithms (Master Thesis)
Effective and error-free implementation of algorithms
|
00001 #ifndef H_MATH_FACTORIZE_ORACLE_BRENT 00002 #define H_MATH_FACTORIZE_ORACLE_BRENT 00003 00004 #include "utils/preconditions/preconditions.h" 00005 #include <vector> 00006 #include <utility> 00007 #include <limits> 00008 #include <stdexcept> 00009 #include "math/powermod/powermod.h" 00010 #include "math/gcd/gcd.h" 00011 #include "utils/rand/rand.h" 00012 #include "utils/branch_predict/branch_predict.h" 00013 00014 namespace math { 00015 namespace factorize { 00016 00017 template <class Powermod> 00018 class OracleBrent_ { 00019 public: 00020 template <typename T> 00021 static T findFactor(T number) 00022 { 00023 Preconditions::check(number > 1, "Too small or negative to factorize"); 00024 Preconditions::check(!math::primes::PrimesFast::isPrime(number), 00025 "nothing to factorize!"); 00026 if (number == 4) { 00027 // This is special case, as pollard/brent will fail for any combination of 00028 // parameters 00029 return 2; 00030 } 00031 00032 const int MAX_ITERATIONS = 1000; 00033 Rand r(47); 00034 for (int i = 0; i < MAX_ITERATIONS; i++) { 00035 T start = r.rand(0, std::min(number, (T) 10000)); 00036 T a = r.rand(0, std::min(number - 1, (T) 10000)); 00037 T res = brent(number, start, a); 00038 if (res != 0) { 00039 assert(res > 1); 00040 assert(number % res == 0); 00041 return res; 00042 } 00043 } 00044 00045 throw std::runtime_error("Can't factorize number"); 00046 } 00047 00048 private: 00049 template <typename T> 00050 static T advance(T x, T n, T a) { 00051 Preconditions::checkRange(x, n); 00052 x = (Powermod::multmod(x, x, n) + a) % n; 00053 if (x < 0) x+= n; // fix on negative value overflow 00054 return x; 00055 } 00056 00057 template <typename T> 00058 static T brent(T number, T start, T a) { 00059 typedef int IterationType; 00060 T tortoise = start; 00061 T hare = start; 00062 T d = 1; 00063 IterationType power_of_two = 1; 00064 IterationType cycle = 0; 00065 00066 while (d == 1) { 00067 if (UNLIKELY(cycle == power_of_two)) { 00068 tortoise = hare; 00069 cycle = 0; 00070 if (power_of_two < std::numeric_limits<IterationType>::max() / 2) { 00071 power_of_two *= 2; 00072 } 00073 } 00074 cycle++; 00075 hare = advance(hare, number, a); 00076 T diff = (tortoise < hare) ? hare - tortoise : tortoise - hare; 00077 d = math::gcd::gcd(diff, number); 00078 } 00079 if (d == number) { 00080 return 0; // can't factorize 00081 } 00082 return d; 00083 } 00084 00085 }; 00086 00087 typedef OracleBrent_<math::powermod::PowermodExtended> OracleBrent; 00088 00089 } // namespace factorize 00090 } // namespace math 00091 #endif