これがやりたかっただけだろシリーズ。素因数分解のためにポラードのロー法を貼った。
https://atcoder.jp/contests/caddi2018/tasks/caddi2018_a
問題
N個の正の整数a_1, ..., a_Nについて、それらの積はPである。
このとき、gcd(a_1,...,a_N)としてあり得る値の最大値を求めよ。
1 <= N, P <= 10^12
解法
Pを素因数分解してp_i^e_iの積になったとする。p_i^floor(e_i/N)の積が答えである。
素因数分解は普通にやってもO(sqrt(P))くらいでできるが、ポラードのロー法などの速いアルゴリズムならO(P^{1/4})くらいでできるはず。
実装上の注意点
- e_iはそこまで大きくならないので、p_i^floor(e_i/n)の積は単純に掛けていくだけでよい。
提出: Submission #3837935 - CADDi 2018 (Rust)
mod pollard_rho { use std::collections::HashMap; /// binary gcd pub fn gcd(mut x: i64, mut y: i64) -> i64 { if y == 0 { return x; } if x == 0 { return y; } let mut sh = 0; while ((x | y) & 1) == 0 { x >>= 1; y >>= 1; sh += 1; } while (x & 1) == 0 { x >>= 1; } while y != 0 { while (y & 1) == 0 { y >>= 1; } if x > y { let t = x; x = y; y = t; } y -= x; } x << sh } fn add_mod(x: i64, y: i64, n: i64) -> i64 { let z = x + y; if z >= n { z - n } else { z } } fn mul_mod(x: i64, mut y: i64, n: i64) -> i64 { assert!(x >= 0); assert!(x < n); let mut sum = 0; let mut cur = x; while y > 0 { if (y & 1) == 1 { sum = add_mod(sum, cur, n); } cur = add_mod(cur, cur, n); y >>= 1; } sum } fn mod_pow(x: i64, mut e: i64, n: i64) -> i64 { let mut prod = if n == 1 { 0 } else { 1 }; let mut cur = x % n; while e > 0 { if (e & 1) == 1 { prod = mul_mod(prod, cur, n); } cur = mul_mod(cur, cur, n); e >>= 1; } prod } pub fn is_prime(n: i64) -> bool { if n <= 1 { return false; } let small = [2, 3, 5, 7, 11, 13]; if small.iter().any(|&u| u == n) { return true; } if small.iter().any(|&u| n % u == 0) { return false; } let mut d = n - 1; let mut e = 0; while (d & 1) == 0 { d >>= 1; e += 1; } let a = [2, 325, 9375, 28178, 450775, 9780504, 1795265022]; a.iter().all(|&a| { if a >= n { return true; } let mut x = mod_pow(a, d, n); if x == 1 { return true; } for _ in 0 .. e { if x == n - 1 { return true; } x = mul_mod(x, x, n); if x == 1 { return false; } } x == 1 }) } fn pollard_rho(n: i64, c: &mut i64) -> i64 { if n % 2 == 0 { return 2; } loop { let mut x: i64 = 2; let mut y = 2; let mut d = 1; let cc = *c; let f = |i| add_mod(mul_mod(i, i, n), cc, n); while d == 1 { x = f(x); y = f(f(y)); d = gcd((x - y).abs(), n); } if d == n { *c += 1; continue; } return d; } } /// Outputs (p, e) in p's ascending order. pub fn factorize(x: i64) -> Vec<(i64, usize)> { if x <= 1 { return Vec::new(); } let mut hm = HashMap::new(); let mut pool = vec![x]; let mut c = 1; while let Some(u) = pool.pop() { if is_prime(u) { *hm.entry(u).or_insert(0) += 1; continue; } let p = pollard_rho(u, &mut c); pool.push(p); pool.push(u / p); } let mut v: Vec<_> = hm.into_iter().collect(); v.sort(); v } } // mod pollard_rho use pollard_rho::*; fn solve() { let out = std::io::stdout(); let mut out = BufWriter::new(out.lock()); macro_rules! puts { ($format:expr) => (write!(out,$format).unwrap()); ($format:expr, $($args:expr),+) => (write!(out,$format,$($args),*).unwrap()) } input! { n: i64, p: i64, } let fac = factorize(p); let mut ans = 1; for (a, e) in fac { let e = e as i64 / n; for _ in 0 .. e { ans *= a; } } puts!("{}\n", ans); }
まとめ
ポラードのロー法で殴ったので記念に書いただけです。ライブラリとしてでも使ってください。