pub struct MontgomeryU64 { // r = 2 ** 64 modulo: u64, // odd and (m < 2 ** 63) modulo_prime: u64, // modulo * modulo_prime = -1 mod r r: u64, // r mod m r2: u64, // r^{2} mod m r3: u64, // r^{3} mod m } impl MontgomeryU64 { /// Returns /// (x, y, g) /// s.t. a * x + b * y = g /// g = gcd(a, b) fn extended_gcd(a: i128, b: i128) -> (i128, i128, u128) { let (mut xs, mut ys, mut s) = (1, 0, a); let (mut xt, mut yt, mut t) = (0, 1, b); while s % t != 0 { let q = s / t; let (xu, yu, u) = (xs - q * xt, ys - q * yt, s - q * t); (xs, ys, xt, yt) = (xt, yt, xu, yu); (s, t) = (t, u); } if t < 0 { (-xt, -yt, (-t) as u128) } else { (xt, yt, t as u128) } } pub fn new(modulo: u64) -> Self { assert!(modulo % 2 == 1); assert!(modulo < 1 << 63); let modulo_prime = { let (_, b, _) = Self::extended_gcd(1_i128 << 64, modulo as i128); if b <= 0 { (-b) as u64 } else { (-b + (1_i128 << 64)) as u64 } }; let r = ((1_u128 << 64) % modulo as u128) as u64; let r2 = ((r as u128 * r as u128) % modulo as u128) as u64; let r3 = ((r2 as u128 * r as u128) % modulo as u128) as u64; Self { modulo, modulo_prime, r, r2, r3, } } /// Returns: /// a * r mod N pub fn to_montgomery_form(&self, a: u64) -> u64 { self.mul(a, self.r2) } /// Montgomery reduction /// /// Returns: /// t * r^{-1} mod N pub fn reduction(&self, t: u128) -> u64 { let t = { ((t + (((t as u64).wrapping_mul(self.modulo_prime)) as u128 * self.modulo as u128)) >> 64) as u64 }; if t < self.modulo { t } else { t - self.modulo } } pub fn modulo(&self) -> u64 { self.modulo } pub fn add(&self, ar: u64, br: u64) -> u64 { let t = ar + br; if t < self.modulo { t } else { t - self.modulo } } pub fn sub(&self, ar: u64, br: u64) -> u64 { let (t, f) = ar.overflowing_sub(br); if !f { t } else { t.wrapping_add(self.modulo) } } pub fn mul(&self, ar: u64, br: u64) -> u64 { self.reduction(ar as u128 * br as u128) } pub fn inv(&self, ar: u64) -> Option { let (x, _, g) = Self::extended_gcd(ar as i128, self.modulo as i128); if g != 1 { None } else { Some(if x < 0 { self.mul((x as i64 + self.modulo as i64) as u64, self.r3) } else { self.mul(x as u64, self.r3) }) } } pub fn pow(&self, ar: u64, mut n: usize) -> u64 { let mut res = self.r; let mut x = ar; while n > 0 { if n % 2 == 1 { res = self.mul(res, x); } x = self.mul(x, x); n /= 2; } res } } /// Returns: /// if n is prime number: /// true /// else: /// false /// /// Algorithm: /// Miller-Rabin /// /// References: /// - [Deterministic variants of the Miller-Rabin primality test. Miller-Rabin SPRP bases records](https://miller-rabin.appspot.com/) /// - [64bit数の素数判定](https://zenn.dev/mizar/articles/791698ea860581) pub fn is_prime(n: u64) -> bool { if n == 0 || n == 1 { return false; } if n == 2 { return true; } if n % 2 == 0 { return false; } let s = (n - 1).trailing_zeros(); let d = (n - 1) >> s; let maybe_prime = |a| { let a = a % n; if a == 0 { return true; } let mont = MontgomeryU64::new(n); let mut ad = mont.pow(mont.to_montgomery_form(a), d as usize); if mont.reduction(ad as u128) == 1 || mont.reduction(ad as u128) == n - 1 { return true; } for _ in 1..s { ad = mont.pow(ad, 2); if mont.reduction(ad as u128) == n - 1 { return true; } } false }; [2, 325, 9375, 28178, 450775, 9780504, 1795265022] .into_iter() .all(maybe_prime) } use proconio::input; fn main() { input! { q: u64, } for _ in 0..q { input! { n: u64, } println!("{} {}", n, if is_prime(n) { 1 } else { 0 }); } }