diff --git a/code/mathematics/divisors.cpp b/code/mathematics/divisors.cpp new file mode 100644 index 0000000..33d5a0a --- /dev/null +++ b/code/mathematics/divisors.cpp @@ -0,0 +1,15 @@ +#include "factor.cpp" +template +vector divisors(T n) { + vector res{{ 1 }}; + for (auto p : factor(n)) { + auto offset = res.size(); + for (size_t i = 0; i < offset; i++) { + for (T j = 1, m = p.fs; j <= p.sc; j++, m *= p.fs) { + res.push_back(m * res[i]); + } + } + } + return res; +} +// vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/divisors.test.cpp b/code/mathematics/divisors.test.cpp new file mode 100644 index 0000000..87869e4 --- /dev/null +++ b/code/mathematics/divisors.test.cpp @@ -0,0 +1,4 @@ +// Field testing: Kattis positivedivisors +void test() { +} +// vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/egcd.cpp b/code/mathematics/egcd.cpp index 85a3a45..08983b2 100644 --- a/code/mathematics/egcd.cpp +++ b/code/mathematics/egcd.cpp @@ -1,5 +1,6 @@ -ll egcd(ll a, ll b, ll& x, ll& y) { +template +T egcd(T a, T b, T& x, T& y) { if (b == 0) { x = 1; y = 0; return a; } - ll d = egcd(b, a % b, x, y); + T d = egcd(b, a % b, x, y); x -= a / b * y; swap(x, y); return d; } // vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/factor.cpp b/code/mathematics/factor.cpp new file mode 100644 index 0000000..3cc8d0c --- /dev/null +++ b/code/mathematics/factor.cpp @@ -0,0 +1,27 @@ +#include "miller_rabin.cpp" +#include "pollard_rho.cpp" +template +map factor(T n) { + if (n == 1) return {}; + if (is_probable_prime(n, 20)) return { { n, 1 } }; + map res; + while (n % 2 == 0) ++res[2], n /= 2; + for (T i = 3; i*i*i <= n; i+=2) { + while (n % i == 0) ++res[i], n /= i; + } + if (is_probable_prime(n, 20)) { + res[n] += 1; + return res; + } + for (T seed_value : {2,3,5,7,11,13,1031}) { + T a_factor = rho(n, seed_value); + if (a_factor != 1) { + ++res[a_factor]; + ++res[n / a_factor]; + return res; + } + } + if (n != 1) ++res[n]; + return res; +} +// vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/factor.test.cpp b/code/mathematics/factor.test.cpp new file mode 100644 index 0000000..87869e4 --- /dev/null +++ b/code/mathematics/factor.test.cpp @@ -0,0 +1,4 @@ +// Field testing: Kattis positivedivisors +void test() { +} +// vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/gcd.cpp b/code/mathematics/gcd.cpp index 9436ef9..0745956 100644 --- a/code/mathematics/gcd.cpp +++ b/code/mathematics/gcd.cpp @@ -1,2 +1,3 @@ -ll gcd(ll a, ll b) { return b == 0 ? a : gcd(b, a % b); } +template +T gcd(T a, T b) { return b == 0 ? a : gcd(b, a % b); } // vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/miller_rabin.cpp b/code/mathematics/miller_rabin.cpp index f0959f4..444de61 100644 --- a/code/mathematics/miller_rabin.cpp +++ b/code/mathematics/miller_rabin.cpp @@ -1,12 +1,13 @@ #include "mod_pow.cpp" -bool is_probable_prime(ll n, int k) { +template +bool is_probable_prime(T n, int k) { if (~n & 1) return n == 2; if (n <= 3) return n == 3; - int s = 0; ll d = n - 1; + int s = 0; T d = n - 1; while (~d & 1) d >>= 1, s++; while (k--) { - ll a = (n - 3) * rand() / RAND_MAX + 2; - ll x = mod_pow(a, d, n); + T a = (n - 3) * rand() / RAND_MAX + 2; + T x = mod_pow(a, d, n); if (x == 1 || x == n - 1) continue; bool ok = false; rep(i,0,s-1) { diff --git a/code/mathematics/mod_inv.cpp b/code/mathematics/mod_inv.cpp index beac770..df2222e 100644 --- a/code/mathematics/mod_inv.cpp +++ b/code/mathematics/mod_inv.cpp @@ -1,5 +1,6 @@ #include "egcd.cpp" -ll mod_inv(ll a, ll m) { - ll x, y, d = egcd(a, m, x, y); +template +T mod_inv(T a, T m) { + T x, y, d = egcd(a, m, x, y); return d == 1 ? smod(x,m) : -1; } // vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/pollard_rho.cpp b/code/mathematics/pollard_rho.cpp index cf5d248..5823df6 100644 --- a/code/mathematics/pollard_rho.cpp +++ b/code/mathematics/pollard_rho.cpp @@ -1,19 +1,18 @@ -// public static int[] seeds = new int[] {2,3,5,7,11,13,1031}; -// public static BigInteger rho(BigInteger n, -// BigInteger seed) { -// int i = 0, -// k = 2; -// BigInteger x = seed, -// y = seed; -// while (i < 1000000) { -// i++; -// x = (x.multiply(x).add(n) -// .subtract(BigInteger.ONE)).mod(n); -// BigInteger d = y.subtract(x).abs().gcd(n); -// if (!d.equals(BigInteger.ONE) && !d.equals(n)) { -// return d; } -// if (i == k) { -// y = x; -// k = k*2; } } -// return BigInteger.ONE; } +#include "gcd.cpp" +template +T rho(T n, T seed) { + T i = 0, k = 2; + T x = seed, y = seed; + while (i < 1000000) { + i++; + x = (x*x + n - 1) % n; + T d = gcd(abs(y - x), n); + if (d != 1 && d != n) return d; + if (i == k) { + y = x; + k = k*2; + } + } + return 1; +} // vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/pollard_rho.java b/code/mathematics/pollard_rho.java new file mode 100644 index 0000000..65ee9a4 --- /dev/null +++ b/code/mathematics/pollard_rho.java @@ -0,0 +1,18 @@ +// public static int[] seeds = new int[] {2,3,5,7,11,13,1031}; +// public static BigInteger rho(BigInteger n, +// BigInteger seed) { +// int i = 0, +// k = 2; +// BigInteger x = seed, +// y = seed; +// while (i < 1000000) { +// i++; +// x = (x.multiply(x).add(n) +// .subtract(BigInteger.ONE)).mod(n); +// BigInteger d = y.subtract(x).abs().gcd(n); +// if (!d.equals(BigInteger.ONE) && !d.equals(n)) { +// return d; } +// if (i == k) { +// y = x; +// k = k*2; } } +// vim: cc=60 ts=2 sts=2 sw=2: diff --git a/code/mathematics/pollard_rho.test.cpp b/code/mathematics/pollard_rho.test.cpp index db33444..87869e4 100644 --- a/code/mathematics/pollard_rho.test.cpp +++ b/code/mathematics/pollard_rho.test.cpp @@ -1,4 +1,4 @@ +// Field testing: Kattis positivedivisors void test() { - } // vim: cc=60 ts=2 sts=2 sw=2: