fork download
copy 
  1. #include <bits/stdc++.h>
  2.  
  3. #define rep(i, a, b) for (int i = a; i < (b); ++i)
  4. #define trav(a, x) for (auto &a : x)
  5. #define all(x) begin(x), end(x)
  6. #define sz(x) (int)(x).size()
  7. #define all(x) begin(x), end(x)
  8. typedef long long ll;
  9. using namespace std;
  10. typedef vector<int> vi;
  11.  
  12. typedef unsigned long long ull;
  13.  
  14. struct timeit {
  15.     decltype(chrono::high_resolution_clock::now()) begin;
  16.     const string label;
  17.     timeit(string label = "???") : label(label) { begin = chrono::high_resolution_clock::now(); }
  18.     ~timeit() {
  19.         auto end = chrono::high_resolution_clock::now();
  20.         auto duration = chrono::duration_cast<chrono::milliseconds>(end - begin).count();
  21.         cout << duration << "ms elapsed [" << label << "]" << endl;
  22.     }
  23. };
  24. mt19937_64 rng(5);
  25. uniform_int_distribution<ull> uni(0, 1ull << 62);
  26.  
  27. int gcd_count = 0;
  28. ull gcd(ull u, ull v) {
  29.     gcd_count++;
  30.     return __gcd(u, v);
  31. }
  32. // ull gcd(ull u, ull v) {
  33. //     gcd_count++;
  34. //     int shift;
  35. //     if (u == 0)
  36. //         return v;
  37. //     if (v == 0)
  38. //         return u;
  39. //     shift = __builtin_ctzl(u | v);
  40. //     u >>= __builtin_ctzl(u);
  41. //     do {
  42. //         v >>= __builtin_ctzl(v);
  43. //         if (u > v)
  44. //             swap(u, v);
  45. //         v = v - u;
  46. //     } while (v != 0);
  47. //     return u << shift;
  48. // }
  49.  
  50. typedef long double ld;
  51.  
  52. int mul_cnt = 0;
  53. ull mod_mul(ull a, ull b, ull M) {
  54.     mul_cnt++;
  55.     ll ret = a * b - M * ull(ld(a) * ld(b) / ld(M));
  56.     return ret + M * (ret < 0) - M * (ret >= (ll)M);
  57. }
  58. ull mod_pow(ull b, ull e, ull mod) {
  59.     ull ans = 1;
  60.     for (; e; b = mod_mul(b, b, mod), e /= 2)
  61.         if (e & 1)
  62.             ans = mod_mul(ans, b, mod);
  63.     return ans;
  64. }
  65. ull mod_add(ull x, ull y, ull p) { return (x += y) < p ? x : x - p; }
  66. bool isPrime(ll n) {
  67.     if (n == 1)
  68.         return 0;
  69.     if (n < 2 || n % 6 % 4 != 1)
  70.         return n - 2 < 2;
  71.     ll A[] = {2, 13, 23, 1662803}, s = __builtin_ctzll(n - 1), d = n >> s;
  72.     for (auto a : A) { // ^ count trailing zeroes
  73.         ll p = mod_pow(a, d, n), i = s;
  74.         while (p != 1 && p != n - 1 && a % n && i--)
  75.             p = mod_mul(p, p, n);
  76.         if (p != n - 1 && i != s)
  77.             return 0;
  78.     }
  79.     return 1;
  80. }
  81. namespace dacin {
  82. namespace Factor {
  83. template <typename T, typename SFINAE = void> struct bigger_type {};
  84. template <typename T> using bigger_type_t = ull;
  85. template <typename T>
  86. struct bigger_type<T, typename enable_if<is_integral<T>::value && is_signed<T>::value && sizeof(T) == 4>::type> {
  87.     using type = long long;
  88. };
  89. template <typename T>
  90. struct bigger_type<T, typename enable_if<is_integral<T>::value && !is_signed<T>::value && sizeof(T) == 4>::type> {
  91.     using type = unsigned long long;
  92. };
  93. /*template<typename T> struct bigger_type<T, typename enable_if<is_integral<T>::value && is_signed<T>::value &&
  94. sizeof(T) == 8>::type>{using type = __int128;}; template<typename T> struct bigger_type<T, typename
  95. enable_if<is_integral<T>::value && !is_signed<T>::value && sizeof(T) == 8>::type>{using type = unsigned __int128;};*/
  96.  
  97. template <typename int_t = unsigned long long> struct Mod_Int {
  98.     static inline int_t add(int_t const &a, int_t const &b, int_t const &mod) {
  99.         int_t ret = a + b;
  100.         if (ret >= mod)
  101.             ret -= mod;
  102.         return ret;
  103.     }
  104.     static inline int_t sub(int_t const &a, int_t const &b, int_t const &mod) { return add(a, mod - b); }
  105.     static inline int_t mul(int_t const &a, int_t const &b, int_t const &mod) {
  106.         return a * static_cast<bigger_type_t<int_t>>(b) % mod;
  107.     }
  108.     static inline int_t pow(int_t const &a, int_t const &b, int_t const &mod) {
  109.         int_t ret = 1;
  110.         int_t base = a;
  111.         for (int i = 0; b >> i; ++i) {
  112.             if ((b >> i) & 1)
  113.                 ret = mul(ret, base, mod);
  114.             base = mul(base, base, mod);
  115.         }
  116.         return ret;
  117.     }
  118. };
  119.  
  120. template <typename T> bool is_prime(T x) {
  121.     static_assert(is_integral<T>::value);
  122.     if (x < T(4))
  123.         return x > T(1);
  124.     for (T i = 2; i * i <= x; ++i) {
  125.         if (x % i == 0)
  126.             return false;
  127.     }
  128.     return true;
  129. }
  130.  
  131. template <typename T> bool miller_rabin_single(T const &x, T base) {
  132.     static_assert(is_integral<T>::value);
  133.     if (x < T(4))
  134.         return x > T(1);
  135.     if (x % 2 == 0)
  136.         return false;
  137.     base %= x;
  138.     if (base == 0)
  139.         return true;
  140.  
  141.     T xm1 = x - 1;
  142.     int j = 1;
  143.     T d = xm1 / 2;
  144.     while (d % 2 == 0) { // could use __builtin_ctz
  145.         d /= 2;
  146.         ++j;
  147.     }
  148.     T t = Mod_Int<T>::pow(base, d, x);
  149.     if (t == T(1) || t == T(xm1))
  150.         return true;
  151.     for (int k = 1; k < j; ++k) {
  152.         t = Mod_Int<T>::mul(t, t, x);
  153.         if (t == xm1)
  154.             return true;
  155.         if (t <= 1)
  156.             break;
  157.     }
  158.     return false;
  159. }
  160.  
  161. template <typename T> bool miller_rabin_multi(T const &) { return true; }
  162. template <typename T, typename... S> bool miller_rabin_multi(T const &x, T const &base, S const &... bases) {
  163.     static_assert(is_integral<T>::value);
  164.     if (!miller_rabin_single(x, base))
  165.         return false;
  166.     return miller_rabin_multi(x, bases...);
  167. }
  168.  
  169. template <typename T> bool miller_rabin(T const &x) {
  170.     static_assert(is_integral<T>::value);
  171.     if (x < 316349281)
  172.         return miller_rabin_multi(x, T(11000544), T(31481107));
  173.     if (x < 4759123141ull)
  174.         return miller_rabin_multi(x, T(2), T(7), T(61));
  175.     return miller_rabin_multi(x, T(2), T(325), T(9375), T(28178), T(450775), T(9780504), T(1795265022));
  176. }
  177.  
  178. template <typename T> T isqrt(T const &x) {
  179.     static_assert(is_integral<T>::value);
  180.     assert(x >= T(0));
  181.     T ret = static_cast<T>(sqrtl(x));
  182.     while (ret > 0 && ret * ret > x)
  183.         --ret;
  184.     while (x - ret * ret > 2 * ret)
  185.         ++ret;
  186.     return ret;
  187. }
  188. template <typename T> T icbrt(T const &x) {
  189.     static_assert(is_integral<T>::value);
  190.     assert(x >= T(0));
  191.     T ret = static_cast<T>(cbrt(x));
  192.     while (ret > 0 && ret * ret * ret > x)
  193.         --ret;
  194.     while (x - ret * ret * ret > 3 * ret * (ret + 1))
  195.         ++ret;
  196.     return ret;
  197. }
  198. vector<uint16_t> saved;
  199. // fast prime factorization from
  200. // https://g...content-available-to-author-only...b.com/radii/msieve
  201. uint64_t squfof_iter_better(uint64_t const &x, uint64_t const &k, uint64_t const &it_max, uint32_t cutoff_div) {
  202.     if (__gcd((uint64_t)k, x) != 1)
  203.         return __gcd((uint64_t)k, x);
  204.     // cerr << "try: " << x << " " << k << "\n";
  205.     saved.clear();
  206.     uint64_t scaledn = k * x;
  207.     if (scaledn >> 62)
  208.         return 1;
  209.     uint32_t sqrtn = isqrt(scaledn);
  210.     uint32_t cutoff = isqrt(2 * sqrtn) / cutoff_div;
  211.     uint32_t q0 = 1;
  212.     uint32_t p1 = sqrtn;
  213.     uint32_t q1 = scaledn - p1 * p1;
  214.  
  215.     if (q1 == 0) {
  216.         uint64_t factor = __gcd(x, (uint64_t)p1);
  217.         return factor == x ? 1 : factor;
  218.     }
  219.  
  220.     uint32_t multiplier = 2 * k;
  221.     uint32_t coarse_cutoff = cutoff * multiplier;
  222.     // cerr << "at: " << multiplier << "\n";
  223.  
  224.     uint32_t i, j;
  225.     uint32_t p0 = 0;
  226.     uint32_t sqrtq = 0;
  227.  
  228.     for (i = 0; i < it_max; ++i) {
  229.         uint32_t q, bits, tmp;
  230.  
  231.         tmp = sqrtn + p1 - q1;
  232.         q = 1;
  233.         if (tmp >= q1)
  234.             q += tmp / q1;
  235.  
  236.         p0 = q * q1 - p1;
  237.         q0 = q0 + (p1 - p0) * q;
  238.  
  239.         if (q1 < coarse_cutoff) {
  240.             tmp = q1 / __gcd(q1, multiplier);
  241.  
  242.             if (tmp < cutoff) {
  243.                 saved.push_back((uint16_t)tmp);
  244.             }
  245.         }
  246.  
  247.         bits = 0;
  248.         tmp = q0;
  249.         while (!(tmp & 1)) {
  250.             bits++;
  251.             tmp >>= 1;
  252.         }
  253.         if (!(bits & 1) && ((tmp & 7) == 1)) {
  254.  
  255.             sqrtq = (uint32_t)isqrt(q0);
  256.  
  257.             if (sqrtq * sqrtq == q0) {
  258.                 for (j = 0; j < saved.size(); ++j) {
  259.                     if (saved[j] == sqrtq)
  260.                         break;
  261.                 }
  262.                 if (j == saved.size())
  263.                     break;
  264.                 // else cerr << "skip " << i << "\n";;
  265.             }
  266.         }
  267.         tmp = sqrtn + p0 - q0;
  268.         q = 1;
  269.         if (tmp >= q0)
  270.             q += tmp / q0;
  271.  
  272.         p1 = q * q0 - p0;
  273.         q1 = q1 + (p0 - p1) * q;
  274.  
  275.         if (q0 < coarse_cutoff) {
  276.             tmp = q0 / __gcd(q0, multiplier);
  277.  
  278.             if (tmp < cutoff) {
  279.                 saved.push_back((uint16_t)tmp);
  280.             }
  281.         }
  282.     }
  283.  
  284.     if (sqrtq == 1) {
  285.         return 1;
  286.     }
  287.     if (i == it_max) {
  288.         return 1;
  289.     }
  290.  
  291.     q0 = sqrtq;
  292.     p1 = p0 + sqrtq * ((sqrtn - p0) / sqrtq);
  293.     q1 = (scaledn - (uint64_t)p1 * (uint64_t)p1) / (uint64_t)q0;
  294.  
  295.     for (j = 0; j < it_max; ++j) {
  296.         uint32_t q, tmp;
  297.  
  298.         tmp = sqrtn + p1 - q1;
  299.         q = 1;
  300.         if (tmp >= q1)
  301.             q += tmp / q1;
  302.  
  303.         p0 = q * q1 - p1;
  304.         q0 = q0 + (p1 - p0) * q;
  305.  
  306.         if (p0 == p1) {
  307.             q0 = q1;
  308.             break;
  309.         }
  310.  
  311.         tmp = sqrtn + p0 - q0;
  312.         q = 1;
  313.         if (tmp >= q0)
  314.             q += tmp / q0;
  315.  
  316.         p1 = q * q0 - p0;
  317.         q1 = q1 + (p0 - p1) * q;
  318.  
  319.         if (p0 == p1)
  320.             break;
  321.     }
  322.     if (j == it_max) {
  323.         cerr << "RNG\n";
  324.         return 1;
  325.     } // random fail
  326.  
  327.     uint64_t factor = __gcd((uint64_t)q0, x);
  328.     if (factor == x)
  329.         factor = 1;
  330.     return factor;
  331. }
  332. uint64_t squfof(uint64_t const &x) {
  333.     static array<uint32_t, 16> multipliers{
  334.         1,     3,      5,      7,         11,         3 * 5,      3 * 7,      3 * 11,
  335.         5 * 7, 5 * 11, 7 * 11, 3 * 5 * 7, 3 * 5 * 11, 3 * 7 * 11, 5 * 7 * 11, 3 * 5 * 7 * 11};
  336.  
  337.     uint64_t cbrt_x = icbrt(x);
  338.     if (cbrt_x * cbrt_x * cbrt_x == x)
  339.         return cbrt_x;
  340.  
  341.     uint32_t iter_lim = isqrt(isqrt(x)) + 10;
  342.     // uint32_t iter_lim = 300;
  343.     for (uint32_t iter_fact = 1; iter_fact < 20000; iter_fact *= 4) {
  344.         for (uint32_t const &k : multipliers) {
  345.             if (numeric_limits<uint64_t>::max() / k <= x)
  346.                 continue; // would overflow
  347.             uint32_t const it_max = iter_fact * iter_lim;
  348.             uint64_t factor = squfof_iter_better(x, k, it_max, 1);
  349.             if (factor == 1 || factor == x)
  350.                 continue;
  351.             return factor;
  352.         }
  353.     }
  354.     cerr << "failed to factor: " << x << "\n";
  355.     assert(0);
  356.     return 1;
  357. }
  358.  
  359. template <typename T> vector<T> factorize(T x) {
  360.     static_assert(is_integral<T>::value);
  361.     vector<T> ret;
  362.     const uint32_t trial_limit = 1000;
  363.     auto trial = [&](int i) {
  364.         while (x % i == 0) {
  365.             x /= i;
  366.             ret.push_back(i);
  367.         }
  368.     };
  369.     trial(2);
  370.     trial(3);
  371.     for (uint32_t i = 5, j = 2; i < trial_limit && i * i <= x; i += j, j = 6 - j) {
  372.         trial(i);
  373.     }
  374.     if (x > 1) {
  375.         static stack<T> s;
  376.         s.push(x);
  377.         while (!s.empty()) {
  378.             x = s.top();
  379.             s.pop();
  380.             if (!miller_rabin(x)) {
  381.                 T factor = squfof(x);
  382.                 if (factor == 1 || factor == x) {
  383.                     assert(0);
  384.                     return ret;
  385.                 }
  386.                 s.push(factor);
  387.                 s.push(x / factor);
  388.             } else {
  389.                 ret.push_back(x);
  390.             }
  391.         }
  392.     }
  393.     sort(ret.begin(), ret.end());
  394.     return ret;
  395. }
  396. } // namespace Factor
  397.  
  398. } // namespace dacin
  399.  
  400. const int maxp = 1e6, maxv = 25;
  401. namespace cancer {
  402.  
  403. ll pollard(ll n) {
  404.     static ll seq[maxp];
  405.     while (1) {
  406.         ll x = rand() % n, y = x, c = rand() % n;
  407.         ll *px = seq, *py = seq, tim = 0, prd = 1;
  408.         while (1) {
  409.             *py++ = y = mod_add(mod_mul(y, y, n), c, n);
  410.             *py++ = y = mod_add(mod_mul(y, y, n), c, n);
  411.             if ((x = *px++) == y)
  412.                 break;
  413.             ll tmp = prd;
  414.             prd = mod_mul(prd, abs(y - x), n);
  415.             if (!prd)
  416.                 return gcd(tmp, n);
  417.             if ((++tim) == maxv) {
  418.                 if ((prd = gcd(prd, n)) > 1 && prd < n)
  419.                     return prd;
  420.                 tim = 0;
  421.             }
  422.         }
  423.         if (tim && (prd = gcd(prd, n)) > 1 && prd < n)
  424.             return prd;
  425.     }
  426. }
  427. } // namespace cancer
  428. namespace kactl {
  429. ull pollard(ull n) {
  430.     auto f = [n](ull x) { return (mod_mul(x, x, n) + 1) % n; };
  431.     if (!(n & 1))
  432.         return 2;
  433.     for (ull i = 2;; i++) {
  434.         ull x = i, y = f(x), p;
  435.         while ((p = gcd(n + y - x, n)) == 1)
  436.             x = f(x), y = f(f(y));
  437.         if (p != n)
  438.             return p;
  439.     }
  440. }
  441. } // namespace kactl
  442. namespace brent {
  443. ll pollard(ll n) {
  444.     if (!(n & 1))
  445.         return 2;
  446.     if (n == 25)
  447.         return 5;
  448.     auto f = [n](ull x) { return mod_add(mod_mul(x, x, n), 1, n); };
  449.     for (ll i = 2;; i = rand() % n) {
  450.         ll power = 1, lam = 1;
  451.         ll x = i, y = f(x);
  452.         ll tim = 0, prd = 1;
  453.         while (1) {
  454.             if (power == lam) {
  455.                 x = y;
  456.                 power *= 2;
  457.                 lam = 0;
  458.             }
  459.             y = f(y);
  460.             lam++;
  461.             if (x == y)
  462.                 break;
  463.             prd = mod_mul(prd, n + y - x, n);
  464.             if (!prd) {
  465.                 return gcd(n + y - x, n);
  466.             }
  467.             if ((++tim) == maxv) {
  468.                 if ((prd = gcd(prd, n)) > 1 && prd < n)
  469.                     return prd;
  470.                 tim = 0;
  471.             }
  472.         }
  473.         if ((prd = gcd(prd, n)) > 1 && prd < n)
  474.             return prd;
  475.     }
  476. }
  477. } // namespace brent
  478. namespace china {
  479. ull n, niv_n;
  480. void setn(ull n_) {
  481.     n = n_;
  482.     niv_n = 1;
  483.     ull x = n;
  484.     for (int i = 1; i <= 62; i++) {
  485.         niv_n *= x;
  486.         x *= x;
  487.     }
  488. }
  489.  
  490. ull rho(ull seed, ull n, ull inc) {
  491.     static const int step = 1 << 9;
  492.  
  493.     setn(n);
  494.     auto f = [n](ull x, ull inc) { return mod_mul(x, x, n) + inc; };
  495.  
  496.     auto sub = [&](ull x, ull y) { return x > y ? x - y : y - x; };
  497.  
  498.     ull y = f(seed, inc), a = f(seed, inc);
  499.  
  500.     for (int l = 1;; l <<= 1) {
  501.         ull x = y;
  502.         for (int i = 1; i <= l; i++)
  503.             y = f(y, inc);
  504.         for (int k = 0; k < l; k += step) {
  505.             int d = min(step, l - k);
  506.             ull g = seed, y0 = y;
  507.  
  508.             for (int i = 1; i <= d; i++) {
  509.                 y = f(y, inc);
  510.                 g = mod_mul(g, sub(x, y), n);
  511.             }
  512.  
  513.             g = gcd(g, n);
  514.  
  515.             if (g == 1)
  516.                 continue;
  517.             if (g != n)
  518.                 return g;
  519.  
  520.             ull y = y0;
  521.  
  522.             while (gcd(sub(x, y), n) == 1)
  523.                 y = f(y, inc);
  524.  
  525.             return gcd(sub(x, y), n) % n;
  526.         }
  527.     }
  528.     return 0;
  529. }
  530. int const num_tries = 10;
  531.  
  532. ull pollard(ull x) {
  533.     if (x % 2 == 0)
  534.         return 2;
  535.     if (x % 3 == 0)
  536.         return 3;
  537.  
  538.     for (ull i = 2; i < num_tries; i++) {
  539.         ull ans = rho(rand(), x, i);
  540.         if (ans != 0 and ans != x)
  541.             return ans;
  542.     }
  543.     return 0;
  544. }
  545. } // namespace china
  546.  
  547. namespace pajene {
  548. ull pollard(ull n) {
  549.     auto f = [n](ull x) { return mod_mul(x, x, n) + 1; };
  550.     if (!(n & 1))
  551.         return 2;
  552.     ull x = 0, y = 0, tim = 0, prd = 1, i = 1, tmp;
  553.     for (; prd; y = f(f(y)), x = f(x)) {
  554.         if (x == y || ++tim == 40) {
  555.             tim = 0;
  556.             if ((prd = gcd(prd, n)) > 1)
  557.                 return prd;
  558.             while (x == y)
  559.                 x = ++i, y = f(x);
  560.         }
  561.         tmp = prd, prd = mod_mul(prd, n + y - x, n);
  562.     }
  563.     return gcd(tmp, n);
  564. }
  565. } // namespace pajene
  566.  
  567. namespace newkactl {
  568. int maxv = 40;
  569. ull pollard(ull n) {
  570.     auto f = [n](ull x) { return mod_mul(x, x, n) + 1; };
  571.     if (!(n & 1))
  572.         return 2;
  573.     for (ull i = 2;; i = rand() % n) {
  574.         ull x = i, y = f(x), tim = 0, prd = 1;
  575.         while ((x = f(x)) != (y = f(f(y)))) {
  576.             if (!(prd = mod_mul(prd, n + y - x, n)))
  577.                 return gcd(n + y - x, n);
  578.             if ((++tim) == maxv) {
  579.                 prd = gcd(prd, n), tim = 0;
  580.                 if (prd > 1)
  581.                     return prd;
  582.             }
  583.         }
  584.         if ((prd = gcd(prd, n)) > 1)
  585.             return prd;
  586.     }
  587. }
  588. } // namespace newkactl
  589. int ptot, pr[maxp], d[maxp];
  590. vector<ull> factor(ull n, const string &version) {
  591.     if (n == 1)
  592.         return {};
  593.     if (isPrime(n))
  594.         return {n};
  595.     ull x = 0;
  596.     if (version == "cancer") {
  597.         x = cancer::pollard(n);
  598.     } else if (version == "pajene") {
  599.         x = pajene::pollard(n);
  600.     } else if (version == "newkactl") {
  601.         x = newkactl::pollard(n);
  602.     } else if (version == "brent") {
  603.         x = brent::pollard(n);
  604.     } else if (version == "china") {
  605.         x = china::pollard(n);
  606.     } else {
  607.         x = kactl::pollard(n);
  608.     }
  609.  
  610.     auto l = factor(x, version), r = factor(n / x, version);
  611.     l.insert(l.end(), all(r));
  612.     return l;
  613. }
  614.  
  615. void sieve() {
  616.     for (int p = 2; p < maxp; p++) {
  617.         if (!d[p])
  618.             pr[ptot++] = d[p] = p;
  619.         for (int j = 0, k; (k = p * pr[j]) < maxp; ++j) {
  620.             d[k] = pr[j];
  621.             if (d[p] == pr[j])
  622.                 break;
  623.         }
  624.     }
  625. }
  626.  
  627. const int ITERS = 1e2;
  628. int main() {
  629.     sieve();
  630.     vector<ull> primes;
  631.     const ull MN_SIZE = 1e9 - 1e5;
  632.     const ull MX_SIZE = 1e9;
  633.     for (int i = MN_SIZE; i < MX_SIZE; i++) {
  634.         if (isPrime(i)) {
  635.             primes.push_back(i);
  636.         }
  637.     }
  638.     vector<ull> vals;
  639.     for (int i = 0; i < ITERS; i++)
  640.         vals.push_back(uni(rng));
  641.  
  642.     vector<ull> spvals;
  643.     for (int i = 0; i < ITERS; i++) {
  644.         ull a = primes[uni(rng) % primes.size()];
  645.         ull b = primes[uni(rng) % primes.size()];
  646.         spvals.push_back(a * b);
  647.     }
  648.  
  649.     vector<ull> seq;
  650.     for (int i = 0; i < ITERS; i++) {
  651.         seq.push_back(i + 2);
  652.     }
  653.     vals = spvals;
  654.     // vals = seq;
  655.  
  656.     {
  657.         gcd_count = 0;
  658.         mul_cnt = 0;
  659.         timeit x("pajene rho");
  660.         int ans = 0;
  661.         for (int i = 0; i < ITERS; i++) {
  662.             ull x = vals[i];
  663.             // cout << x << endl;
  664.             auto res = factor(x, "pajene");
  665.             for (auto x : res)
  666.                 ans += x;
  667.         }
  668.         cout << ans << endl;
  669.         cout << gcd_count << ' ' << mul_cnt << endl;
  670.     }
  671.     // {
  672.     //     cout << endl;
  673.     //     gcd_count = 0;
  674.     //     mul_cnt = 0;
  675.     //     timeit x("brent rho");
  676.     //     int ans = 0;
  677.     //     for (int i = 0; i < ITERS; i++) {
  678.     //         ull x = vals[i];
  679.     //         // cout << x << endl;
  680.     //         auto res = factor(x, "brent");
  681.     //         ans += res[0];
  682.     //     }
  683.     //     cout << ans << endl;
  684.     //     cout << gcd_count << ' ' << mul_cnt << endl;
  685.     // }
  686.     // {
  687.     //     cout << endl;
  688.     //     gcd_count = 0;
  689.     //     mul_cnt = 0;
  690.     //     timeit x("new kactl rho");
  691.     //     int ans = 0;
  692.     //     for (int i = 0; i < ITERS; i++) {
  693.     //         ull x = vals[i];
  694.     //         // cout << x << endl;
  695.     //         auto res = factor(x, "newkactl");
  696.     //         for (auto x : res)
  697.     //             ans += x;
  698.     //     }
  699.     //     cout << ans << endl;
  700.     //     cout << gcd_count << ' ' << mul_cnt << endl;
  701.     // }
  702.     // {
  703.     //     gcd_count = 0;
  704.     //     timeit x("cancer rho");
  705.     //     int ans = 0;
  706.     //     for (int i = 0; i < ITERS; i++) {
  707.     //         ull x = vals[i];
  708.     //         auto res = factor(x, "cancer");
  709.     //         for (auto x : res)
  710.     //             ans += x;
  711.     //     }
  712.     //     cout << ans << endl;
  713.     //     cout << gcd_count << endl;
  714.     // }
  715.     // {
  716.     //     gcd_count = 0;
  717.     //     timeit x("new kactl cache rho");
  718.     //     int ans = 0;
  719.     //     for (int i = 0; i < ITERS; i++) {
  720.     //         ull x = vals[i];
  721.     //         auto res = factor(x, "newkactlcache");
  722.     //         for (auto x : res)
  723.     //             ans += x;
  724.     //     }
  725.     //     cout << ans << endl;
  726.     //     cout << gcd_count << endl;
  727.     // }
  728.     // {
  729.     //     gcd_count = 0;
  730.     //     timeit x("china rho");
  731.     //     int ans = 0;
  732.     //     for (int i = 0; i < ITERS; i++) {
  733.     //         ull x = vals[i];
  734.     //         auto res = factor(x, "china");
  735.     //         for (auto x : res)
  736.     //             ans += x;
  737.     //     }
  738.     //     cout << ans << endl;
  739.     //     cout << gcd_count << endl;
  740.     // }
  741.     {
  742.         cout << endl;
  743.         gcd_count = 0;
  744.         mul_cnt = 0;
  745.         timeit x("kactl rho");
  746.         int ans = 0;
  747.         for (int i = 0; i < ITERS; i++) {
  748.             ull x = vals[i];
  749.             auto res = factor(x, "kactl");
  750.             for (auto x : res)
  751.                 ans += x;
  752.         }
  753.         cout << ans << endl;
  754.         cout << gcd_count << ' ' << mul_cnt << endl;
  755.     }
  756.     // {
  757.     //     cout << endl;
  758.     //     timeit x("dacin factorize");
  759.     //     int ans = 0;
  760.     //     for (int i = 0; i < ITERS; i++) {
  761.     //         ull x = vals[i];
  762.     //         auto res = dacin::Factor::factorize(x);
  763.     //         for (auto x : res)
  764.     //             ans += x;
  765.     //     }
  766.     //     cout << ans << endl;
  767.     // }
  768. }
Success #stdin #stdout 1.32s 10920KB
comments (?)
stdin
copy 
23451235245738996677668843235
stdout
copy 
-1873832444
58314 9359255
143ms elapsed [pajene rho]

-1873832444
2326645 7024650
1143ms elapsed [kactl rho]
https://ideone.com/Rsw8wY
language:
C++ (gcc 8.3)
created:
6 days ago
visibility:
public

Share or Embed source code

Discover > Sphere Engine API

The brand new service which powers Ideone!

Discover > IDE Widget

Widget for compiling and running the source code in a web browser!