結果

問題 No.2578 Jewelry Store
ユーザー akakimidori
提出日時 2023-12-06 00:20:36
言語 Rust
(1.83.0 + proconio)
結果
AC  
実行時間 571 ms / 3,500 ms
コード長 26,346 bytes
コンパイル時間 17,120 ms
コンパイル使用メモリ 379,948 KB
実行使用メモリ 11,008 KB
最終ジャッジ日時 2024-09-27 00:41:39
合計ジャッジ時間 23,945 ms
ジャッジサーバーID
(参考情報)
judge1 / judge5
このコードへのチャレンジ
(要ログイン)
ファイルパターン 結果
other AC * 54
権限があれば一括ダウンロードができます
コンパイルメッセージ
warning: type alias `Map` is never used
   --> src/main.rs:940:6
    |
940 | type Map<K, V> = BTreeMap<K, V>;
    |      ^^^
    |
    = note: `#[warn(dead_code)]` on by default

warning: type alias `Set` is never used
   --> src/main.rs:941:6
    |
941 | type Set<T> = BTreeSet<T>;
    |      ^^^

warning: type alias `Deque` is never used
   --> src/main.rs:942:6
    |
942 | type Deque<T> = VecDeque<T>;
    |      ^^^^^

warning: variable `D` should have a snake case name
   --> src/main.rs:782:9
    |
782 |     let D = if n % 4 == 1 { 2 * n } else { n };
    |         ^ help: convert the identifier to snake case: `d`
    |
    = note: `#[warn(non_snake_case)]` on by default

warning: variable `sqrt_D` should have a snake case name
   --> src/main.rs:783:10
    |
783 |     let (sqrt_D, mut Q) = D.isqrt_rem();
    |          ^^^^^^ help: convert the identifier to snake case: `sqrt_d`

warning: variable `Q` should have a snake case name
   --> src/main.rs:783:22
    |
783 |     let (sqrt_D, mut Q) = D.isqrt_rem();
    |                      ^ help: convert the identifier to snake case: `q`

warning: variable `Q_hat` should have a snake case name
   --> src/main.rs:784:13
    |
784 |     let mut Q_hat = 1_u64;
    |             ^^^^^ help: convert the identifier to snake case: `q_hat`

warning: variable `P` should have a snake case name
   --> src/main.rs:785:13
    |
785 |     let mut P = sqrt_D;
    |             ^ help: convert the identifier to snake case: `p`

warning: variable `L` should have a snake case name
   --> src/main.rs:786:9
    |
786 |     let L = (2.0 * (2.0 * (D as f64).sqrt()).sqrt()) as u64;
    |         ^ help: convert the identifier to snake case: `l`

warning: variable `B` should have a snake case name
   --> src/main.rs:787:9
    |
787 |     let B = 2 * L;
    |         ^ help: convert the identifier to snake case: `b`

warning: variable `P_prime` should have a snake case name
   --> src/main.rs:798:13
    |
798 |         let P_prime = q * Q -

ソースコード

diff #
プレゼンテーションモードにする

// ---------- begin bitwise transform ----------
pub fn bitwise_transform<T, F>(a: &mut [T], mut f: F)
where
F: FnMut(&mut T, &mut T)
{
let n = a.len().trailing_zeros() as usize;
assert!(a.len() == 1 << n);
for i in 0..n {
for a in a.chunks_exact_mut(2 << i) {
let (x, y) = a.split_at_mut(1 << i);
for (x, y) in x.iter_mut().zip(y) {
f(x, y);
}
}
}
}
// ---------- end bitwise transform ----------
use std::ops::*;
// ---------- begin trait ----------
pub trait Zero: Sized + Add<Self, Output = Self> {
fn zero() -> Self;
fn is_zero(&self) -> bool;
}
pub trait One: Sized + Mul<Self, Output = Self> {
fn one() -> Self;
fn is_one(&self) -> bool;
}
pub trait Ring: Zero + One + Sub<Output = Self> {}
pub trait Field: Ring + Div<Output = Self> {}
// ---------- end trait ----------
// ---------- begin modint ----------
pub const fn pow_mod(mut r: u32, mut n: u32, m: u32) -> u32 {
let mut t = 1;
while n > 0 {
if n & 1 == 1 {
t = (t as u64 * r as u64 % m as u64) as u32;
}
r = (r as u64 * r as u64 % m as u64) as u32;
n >>= 1;
}
t
}
pub const fn primitive_root(p: u32) -> u32 {
let mut m = p - 1;
let mut f = [1; 30];
let mut k = 0;
let mut d = 2;
while d * d <= m {
if m % d == 0 {
f[k] = d;
k += 1;
}
while m % d == 0 {
m /= d;
}
d += 1;
}
if m > 1 {
f[k] = m;
k += 1;
}
let mut g = 1;
while g < p {
let mut ok = true;
let mut i = 0;
while i < k {
ok &= pow_mod(g, (p - 1) / f[i], p) > 1;
i += 1;
}
if ok {
break;
}
g += 1;
}
g
}
pub const fn is_prime(n: u32) -> bool {
if n <= 1 {
return false;
}
let mut d = 2;
while d * d <= n {
if n % d == 0 {
return false;
}
d += 1;
}
true
}
#[derive(Clone, Copy, PartialEq, Eq)]
pub struct ModInt<const M: u32>(u32);
impl<const M: u32> ModInt<{ M }> {
const REM: u32 = {
let mut t = 1u32;
let mut s = !M + 1;
let mut n = !0u32 >> 2;
while n > 0 {
if n & 1 == 1 {
t = t.wrapping_mul(s);
}
s = s.wrapping_mul(s);
n >>= 1;
}
t
};
const INI: u64 = ((1u128 << 64) % M as u128) as u64;
const IS_PRIME: () = assert!(is_prime(M));
const PRIMITIVE_ROOT: u32 = primitive_root(M);
const ORDER: usize = 1 << (M - 1).trailing_zeros();
const fn reduce(x: u64) -> u32 {
let _ = Self::IS_PRIME;
let b = (x as u32 * Self::REM) as u64;
let t = x + b * M as u64;
let mut c = (t >> 32) as u32;
if c >= M {
c -= M;
}
c as u32
}
const fn multiply(a: u32, b: u32) -> u32 {
Self::reduce(a as u64 * b as u64)
}
pub const fn new(v: u32) -> Self {
assert!(v < M);
Self(Self::reduce(v as u64 * Self::INI))
}
pub const fn const_mul(&self, rhs: Self) -> Self {
Self(Self::multiply(self.0, rhs.0))
}
pub const fn pow(&self, mut n: u64) -> Self {
let mut t = Self::new(1);
let mut r = *self;
while n > 0 {
if n & 1 == 1 {
t = t.const_mul(r);
}
r = r.const_mul(r);
n >>= 1;
}
t
}
pub const fn inv(&self) -> Self {
assert!(self.0 != 0);
self.pow(M as u64 - 2)
}
pub const fn get(&self) -> u32 {
Self::reduce(self.0 as u64)
}
pub const fn zero() -> Self {
Self::new(0)
}
pub const fn one() -> Self {
Self::new(1)
}
}
impl<const M: u32> Add for ModInt<{ M }> {
type Output = Self;
fn add(self, rhs: Self) -> Self::Output {
let mut v = self.0 + rhs.0;
if v >= M {
v -= M;
}
Self(v)
}
}
impl<const M: u32> Sub for ModInt<{ M }> {
type Output = Self;
fn sub(self, rhs: Self) -> Self::Output {
let mut v = self.0 - rhs.0;
if self.0 < rhs.0 {
v += M;
}
Self(v)
}
}
impl<const M: u32> Mul for ModInt<{ M }> {
type Output = Self;
fn mul(self, rhs: Self) -> Self::Output {
self.const_mul(rhs)
}
}
impl<const M: u32> Div for ModInt<{ M }> {
type Output = Self;
fn div(self, rhs: Self) -> Self::Output {
self * rhs.inv()
}
}
impl<const M: u32> AddAssign for ModInt<{ M }> {
fn add_assign(&mut self, rhs: Self) {
*self = *self + rhs;
}
}
impl<const M: u32> SubAssign for ModInt<{ M }> {
fn sub_assign(&mut self, rhs: Self) {
*self = *self - rhs;
}
}
impl<const M: u32> MulAssign for ModInt<{ M }> {
fn mul_assign(&mut self, rhs: Self) {
*self = *self * rhs;
}
}
impl<const M: u32> DivAssign for ModInt<{ M }> {
fn div_assign(&mut self, rhs: Self) {
*self = *self / rhs;
}
}
impl<const M: u32> Neg for ModInt<{ M }> {
type Output = Self;
fn neg(self) -> Self::Output {
if self.0 == 0 {
self
} else {
Self(M - self.0)
}
}
}
impl<const M: u32> std::fmt::Display for ModInt<{ M }> {
fn fmt<'a>(&self, f: &mut std::fmt::Formatter<'a>) -> std::fmt::Result {
write!(f, "{}", self.get())
}
}
impl<const M: u32> std::fmt::Debug for ModInt<{ M }> {
fn fmt<'a>(&self, f: &mut std::fmt::Formatter<'a>) -> std::fmt::Result {
write!(f, "{}", self.get())
}
}
impl<const M: u32> std::str::FromStr for ModInt<{ M }> {
type Err = std::num::ParseIntError;
fn from_str(s: &str) -> Result<Self, Self::Err> {
let val = s.parse::<u32>()?;
Ok(ModInt::new(val))
}
}
impl<const M: u32> From<usize> for ModInt<{ M }> {
fn from(val: usize) -> ModInt<{ M }> {
ModInt::new((val % M as usize) as u32)
}
}
// ---------- end modint ----------
// ---------- begin precalc ----------
pub struct Precalc<const MOD: u32> {
fact: Vec<ModInt<MOD>>,
ifact: Vec<ModInt<MOD>>,
inv: Vec<ModInt<MOD>>,
}
impl<const MOD: u32> Precalc<MOD> {
pub fn new(size: usize) -> Self {
let mut fact = vec![ModInt::one(); size + 1];
let mut ifact = vec![ModInt::one(); size + 1];
let mut inv = vec![ModInt::one(); size + 1];
for i in 2..=size {
fact[i] = fact[i - 1] * ModInt::from(i);
}
ifact[size] = fact[size].inv();
for i in (2..=size).rev() {
inv[i] = ifact[i] * fact[i - 1];
ifact[i - 1] = ifact[i] * ModInt::from(i);
}
Self { fact, ifact, inv }
}
pub fn fact(&self, n: usize) -> ModInt<MOD> {
self.fact[n]
}
pub fn ifact(&self, n: usize) -> ModInt<MOD> {
self.ifact[n]
}
pub fn inv(&self, n: usize) -> ModInt<MOD> {
assert!(0 < n);
self.inv[n]
}
pub fn perm(&self, n: usize, k: usize) -> ModInt<MOD> {
if k > n {
return ModInt::zero();
}
self.fact[n] * self.ifact[n - k]
}
pub fn binom(&self, n: usize, k: usize) -> ModInt<MOD> {
if n < k {
return ModInt::zero();
}
self.fact[n] * self.ifact[k] * self.ifact[n - k]
}
}
// ---------- end precalc ----------
impl<const M: u32> Zero for ModInt<{ M }> {
fn zero() -> Self {
Self::zero()
}
fn is_zero(&self) -> bool {
self.0 == 0
}
}
impl<const M: u32> One for ModInt<{ M }> {
fn one() -> Self {
Self::one()
}
fn is_one(&self) -> bool {
self.get() == 1
}
}
impl<const M: u32> Ring for ModInt<{ M }> {}
impl<const M: u32> Field for ModInt<{ M }> {}
// ---------- begin array op ----------
struct NTTPrecalc<const M: u32> {
sum_e: [ModInt<{ M }>; 30],
sum_ie: [ModInt<{ M }>; 30],
}
impl<const M: u32> NTTPrecalc<{ M }> {
const fn new() -> Self {
let cnt2 = (M - 1).trailing_zeros() as usize;
let root = ModInt::new(ModInt::<{ M }>::PRIMITIVE_ROOT);
let zeta = root.pow((M - 1) as u64 >> cnt2);
let mut es = [ModInt::zero(); 30];
let mut ies = [ModInt::zero(); 30];
let mut sum_e = [ModInt::zero(); 30];
let mut sum_ie = [ModInt::zero(); 30];
let mut e = zeta;
let mut ie = e.inv();
let mut i = cnt2;
while i >= 2 {
es[i - 2] = e;
ies[i - 2] = ie;
e = e.const_mul(e);
ie = ie.const_mul(ie);
i -= 1;
}
let mut now = ModInt::one();
let mut inow = ModInt::one();
let mut i = 0;
while i < cnt2 - 1 {
sum_e[i] = es[i].const_mul(now);
sum_ie[i] = ies[i].const_mul(inow);
now = ies[i].const_mul(now);
inow = es[i].const_mul(inow);
i += 1;
}
Self { sum_e, sum_ie }
}
}
struct NTTPrecalcHelper<const MOD: u32>;
impl<const MOD: u32> NTTPrecalcHelper<MOD> {
const A: NTTPrecalc<MOD> = NTTPrecalc::new();
}
pub trait ArrayAdd {
type Item;
fn add(&self, rhs: &[Self::Item]) -> Vec<Self::Item>;
}
impl<T> ArrayAdd for [T]
where
T: Zero + Copy,
{
type Item = T;
fn add(&self, rhs: &[Self::Item]) -> Vec<Self::Item> {
let mut c = vec![T::zero(); self.len().max(rhs.len())];
c[..self.len()].copy_from_slice(self);
c.add_assign(rhs);
c
}
}
pub trait ArrayAddAssign {
type Item;
fn add_assign(&mut self, rhs: &[Self::Item]);
}
impl<T> ArrayAddAssign for [T]
where
T: Add<Output = T> + Copy,
{
type Item = T;
fn add_assign(&mut self, rhs: &[Self::Item]) {
assert!(self.len() >= rhs.len());
self.iter_mut().zip(rhs).for_each(|(x, a)| *x = *x + *a);
}
}
impl<T> ArrayAddAssign for Vec<T>
where
T: Zero + Add<Output = T> + Copy,
{
type Item = T;
fn add_assign(&mut self, rhs: &[Self::Item]) {
if self.len() < rhs.len() {
self.resize(rhs.len(), T::zero());
}
self.as_mut_slice().add_assign(rhs);
}
}
pub trait ArraySub {
type Item;
fn sub(&self, rhs: &[Self::Item]) -> Vec<Self::Item>;
}
impl<T> ArraySub for [T]
where
T: Zero + Sub<Output = T> + Copy,
{
type Item = T;
fn sub(&self, rhs: &[Self::Item]) -> Vec<Self::Item> {
let mut c = vec![T::zero(); self.len().max(rhs.len())];
c[..self.len()].copy_from_slice(self);
c.sub_assign(rhs);
c
}
}
pub trait ArraySubAssign {
type Item;
fn sub_assign(&mut self, rhs: &[Self::Item]);
}
impl<T> ArraySubAssign for [T]
where
T: Sub<Output = T> + Copy,
{
type Item = T;
fn sub_assign(&mut self, rhs: &[Self::Item]) {
assert!(self.len() >= rhs.len());
self.iter_mut().zip(rhs).for_each(|(x, a)| *x = *x - *a);
}
}
impl<T> ArraySubAssign for Vec<T>
where
T: Zero + Sub<Output = T> + Copy,
{
type Item = T;
fn sub_assign(&mut self, rhs: &[Self::Item]) {
if self.len() < rhs.len() {
self.resize(rhs.len(), T::zero());
}
self.as_mut_slice().sub_assign(rhs);
}
}
pub trait ArrayDot {
type Item;
fn dot(&self, rhs: &[Self::Item]) -> Vec<Self::Item>;
}
impl<T> ArrayDot for [T]
where
T: Mul<Output = T> + Copy,
{
type Item = T;
fn dot(&self, rhs: &[Self::Item]) -> Vec<Self::Item> {
assert!(self.len() == rhs.len());
self.iter().zip(rhs).map(|p| *p.0 * *p.1).collect()
}
}
pub trait ArrayDotAssign {
type Item;
fn dot_assign(&mut self, rhs: &[Self::Item]);
}
impl<T> ArrayDotAssign for [T]
where
T: MulAssign + Copy,
{
type Item = T;
fn dot_assign(&mut self, rhs: &[Self::Item]) {
assert!(self.len() == rhs.len());
self.iter_mut().zip(rhs).for_each(|(x, a)| *x *= *a);
}
}
pub trait ArrayMul {
type Item;
fn mul(&self, rhs: &[Self::Item]) -> Vec<Self::Item>;
}
impl<T> ArrayMul for [T]
where
T: Zero + One + Copy,
{
type Item = T;
fn mul(&self, rhs: &[Self::Item]) -> Vec<Self::Item> {
if self.is_empty() || rhs.is_empty() {
return vec![];
}
let mut res = vec![T::zero(); self.len() + rhs.len() - 1];
for (i, a) in self.iter().enumerate() {
for (res, b) in res[i..].iter_mut().zip(rhs.iter()) {
*res = *res + *a * *b;
}
}
res
}
}
// transform len=1NTT
pub trait ArrayConvolution {
type Item;
fn transform(&mut self, len: usize);
fn inverse_transform(&mut self, len: usize);
fn convolution(&self, rhs: &[Self::Item]) -> Vec<Self::Item>;
}
impl<const M: u32> ArrayConvolution for [ModInt<{ M }>] {
type Item = ModInt<{ M }>;
fn transform(&mut self, len: usize) {
let f = self;
let n = f.len();
let k = (n / len).trailing_zeros() as usize;
assert!(len << k == n);
assert!(k <= ModInt::<{ M }>::ORDER);
let pre = &NTTPrecalcHelper::<{ M }>::A;
for ph in 1..=k {
let p = len << (k - ph);
let mut now = ModInt::one();
for (i, f) in f.chunks_exact_mut(2 * p).enumerate() {
let (x, y) = f.split_at_mut(p);
for (x, y) in x.iter_mut().zip(y.iter_mut()) {
let l = *x;
let r = *y * now;
*x = l + r;
*y = l - r;
}
now *= pre.sum_e[(!i).trailing_zeros() as usize];
}
}
}
fn inverse_transform(&mut self, len: usize) {
let f = self;
let n = f.len();
let k = (n / len).trailing_zeros() as usize;
assert!(len << k == n);
assert!(k <= ModInt::<{ M }>::ORDER);
let pre = &NTTPrecalcHelper::<{ M }>::A;
for ph in (1..=k).rev() {
let p = len << (k - ph);
let mut inow = ModInt::one();
for (i, f) in f.chunks_exact_mut(2 * p).enumerate() {
let (x, y) = f.split_at_mut(p);
for (x, y) in x.iter_mut().zip(y.iter_mut()) {
let l = *x;
let r = *y;
*x = l + r;
*y = (l - r) * inow;
}
inow *= pre.sum_ie[(!i).trailing_zeros() as usize];
}
}
let ik = ModInt::new(2).inv().pow(k as u64);
for f in f.iter_mut() {
*f *= ik;
}
}
fn convolution(&self, rhs: &[Self::Item]) -> Vec<Self::Item> {
if self.len().min(rhs.len()) <= 32 {
return self.mul(rhs);
}
const PARAM: usize = 10;
let size = self.len() + rhs.len() - 1;
let mut k = 0;
while (size + (1 << k) - 1) >> k > PARAM {
k += 1;
}
let len = (size + (1 << k) - 1) >> k;
let mut f = vec![ModInt::zero(); len << k];
let mut g = vec![ModInt::zero(); len << k];
f[..self.len()].copy_from_slice(self);
g[..rhs.len()].copy_from_slice(rhs);
f.transform(len);
g.transform(len);
let mut buf = [ModInt::zero(); 2 * PARAM - 1];
let buf = &mut buf[..(2 * len - 1)];
let pre = &NTTPrecalcHelper::<{ M }>::A;
let mut now = ModInt::one();
for (i, (f, g)) in f
.chunks_exact_mut(2 * len)
.zip(g.chunks_exact(2 * len))
.enumerate()
{
let mut r = now;
for (f, g) in f.chunks_exact_mut(len).zip(g.chunks_exact(len)) {
buf.fill(ModInt::zero());
for (i, f) in f.iter().enumerate() {
for (buf, g) in buf[i..].iter_mut().zip(g.iter()) {
*buf = *buf + *f * *g;
}
}
f.copy_from_slice(&buf[..len]);
for (f, buf) in f.iter_mut().zip(buf[len..].iter()) {
*f = *f + r * *buf;
}
r = -r;
}
now *= pre.sum_e[(!i).trailing_zeros() as usize];
}
f.inverse_transform(len);
f.truncate(self.len() + rhs.len() - 1);
f
}
}
// ---------- end array op ----------
pub trait ModularArithmetic: Sized {
fn add_mod(self, other: Self, m: Self) -> Self;
fn sub_mod(self, other: Self, m: Self) -> Self;
fn mul_mod(self, other: Self, m: Self) -> Self;
fn div_mod(self, other: Self, m: Self) -> Option<Self>;
fn inv_mod(self, m: Self) -> Option<Self>;
fn pow_mod(self, exp: Self, m: Self) -> Self;
}
impl ModularArithmetic for u64 {
fn add_mod(self, other: u64, m: u64) -> u64 {
let value = self + other;
if value >= m {
value - m
} else {
value
}
}
fn sub_mod(self, other: u64, m: u64) -> u64 {
if self >= other {
self - other
} else {
self + m - other
}
}
fn mul_mod(self, other: u64, m: u64) -> u64 {
//(self * other) % m
((self as u128 * other as u128) % m as u128) as u64
}
fn inv_mod(self, m: u64) -> Option<u64> {
1.div_mod(self, m)
}
fn div_mod(self, other: u64, m: u64) -> Option<u64> {
let mut x = (m, 0u64);
let mut y = (other, 1u64);
while y.0 != 0 {
x.1 = x.1.wrapping_sub((x.0 / y.0).wrapping_mul(y.1));
x.0 %= y.0;
std::mem::swap(&mut x, &mut y);
}
if self % x.0 != 0 {
None
} else {
Some((self / x.0).mul_mod(if x.1 >= m { x.1.wrapping_add(m) } else { x.1 }, m))
}
}
fn pow_mod(self, mut exp: Self, m: Self) -> Self {
let mut x = 1;
let mut a = self;
while exp > 0 {
if exp % 2 == 1 {
x = x.mul_mod(a, m);
}
a = a.mul_mod(a, m);
exp >>= 1;
}
x
}
}
pub trait MillerRabin {
fn is_prime(self) -> bool;
}
impl MillerRabin for u64 {
fn is_prime(self: u64) -> bool {
if self <= 3 {
return self == 2 || self == 3;
}
if self % 2 == 0 {
return false;
}
let r = (self - 1).trailing_zeros();
let d = (self - 1) >> r;
for &a in &[2, 325, 9375, 28178, 450775, 9780504, 1795265022] {
let a = a % self;
if a == 0 {
continue;
}
let mut x = a.pow_mod(d, self);
if x == 1 || x == self - 1 {
continue;
}
let mut i = r;
while i > 0 {
x = x.mul_mod(x, self);
if x == 1 {
return false;
} else if x == self - 1 {
break;
}
i -= 1;
}
if i == 0 {
return false;
}
}
return true;
}
}
pub trait IntegerSquareRootWithRemainder {
type Output;
/// Find ($s$, $q$) such that $n = s^2 + q$ with $n \in [s^2, (s+1)^2)$.
fn isqrt_rem(self) -> (Self::Output, Self::Output);
}
impl IntegerSquareRootWithRemainder for u64 {
type Output = u64;
fn isqrt_rem(self) -> (Self::Output, Self::Output) {
let mut one = 1 << 62;
while one > self {
one >>= 2;
}
let mut rem = self;
let mut res = 0;
while one > 0 {
if rem >= res + one {
rem -= res + one;
res = res + 2 * one;
}
res >>= 1;
one >>= 2;
}
(res, rem)
}
}
pub fn binary_gcd(mut a: u64, mut b: u64) -> u64 {
if a == 0 || b == 0 {
return a + b;
}
let k = (a | b).trailing_zeros();
a >>= a.trailing_zeros();
b >>= b.trailing_zeros();
while a != b {
if a < b {
std::mem::swap(&mut a, &mut b);
}
a -= b;
a >>= a.trailing_zeros();
}
a << k
}
fn squfof_internal(n: u64) -> Option<u64> {
let (sqrt_n, rem) = n.isqrt_rem();
if rem == 0 {
return Some(sqrt_n);
}
let D = if n % 4 == 1 { 2 * n } else { n };
let (sqrt_D, mut Q) = D.isqrt_rem();
let mut Q_hat = 1_u64;
let mut P = sqrt_D;
let L = (2.0 * (2.0 * (D as f64).sqrt()).sqrt()) as u64;
let B = 2 * L;
let mut queue = VecDeque::new();
// println!("Step 1: {} {} {}", P, Q, L);
let mut i = 1;
let r = loop {
// eprintln!("Step 2 {}: {} {}", i, P, Q);
// 2a
let q = (sqrt_D + P) / Q;
let P_prime = q * Q - P;
// 2b
if Q <= L {
if Q % 2 == 0 {
queue.push_back((Q / 2, P % (Q / 2)));
} else if 2 * Q <= L {
queue.push_back((Q, P % Q));
}
}
// 2c
let t = Q_hat.wrapping_add(q.wrapping_mul(P.wrapping_sub(P_prime)));
Q_hat = Q;
Q = t;
P = P_prime;
// 2d
let (r, rem) = Q.isqrt_rem();
// maybe a typo in the paper (odd / even)
if i % 2 == 1 && rem == 0 {
let mut idx = 0;
while idx < queue.len() {
let (r_, t) = queue[idx];
if r == r_ && P % r == t % r {
break;
}
idx += 1;
}
if idx == queue.len() {
break r;
} else if r > 1 {
queue = queue.split_off(idx + 1);
} else {
return None;
//unreachable!()
}
}
i += 1;
if i >= B {
return None;
}
};
// 3
// eprintln!("Step 3: Q_hat = {}, P = {}, r = {}", Q_hat, P, r);
Q_hat = r;
P = P + r * ((sqrt_D - P) / r);
Q = (D - P * P) / Q_hat;
assert!((D - P * P) % Q_hat == 0);
// 4
loop {
// eprintln!("{} {}", P, Q);
let q = (sqrt_D + P) / Q;
let P_prime = q * Q - P;
if P == P_prime {
break;
}
let t = Q_hat.wrapping_add(q.wrapping_mul(P.wrapping_sub(P_prime)));
Q_hat = Q;
Q = t;
P = P_prime;
}
// step 5
if Q % 2 == 0 {
return Some(Q / 2);
} else {
return Some(Q);
}
}
pub fn squfof(n: u64) -> u64 {
for m in 2.. {
if m > 1 && n % m == 0 {
return m;
}
if let Some(k) = squfof_internal(m * n) {
let g = binary_gcd(k, n);
if g != 1 && g != n {
return g;
}
}
}
unreachable!();
}
pub trait Factorize: MillerRabin {
type Output;
fn factorize(self) -> HashMap<Self::Output, u64>;
}
impl Factorize for u64 {
type Output = u64;
fn factorize(self) -> HashMap<Self::Output, u64> {
let mut result = HashMap::new();
let mut stack = vec![self];
while let Some(n) = stack.pop() {
if n == 1 {
continue;
} else if n.is_prime() {
*result.entry(n).or_insert(0) += 1;
} else {
let x = squfof(n);
stack.push(x);
stack.push(n / x);
}
}
result
}
}
// ---------- begin scannner ----------
#[allow(dead_code)]
mod scanner {
use std::str::FromStr;
pub struct Scanner<'a> {
it: std::str::SplitWhitespace<'a>,
}
impl<'a> Scanner<'a> {
pub fn new(s: &'a String) -> Scanner<'a> {
Scanner {
it: s.split_whitespace(),
}
}
pub fn next<T: FromStr>(&mut self) -> T {
self.it.next().unwrap().parse::<T>().ok().unwrap()
}
pub fn next_bytes(&mut self) -> Vec<u8> {
self.it.next().unwrap().bytes().collect()
}
pub fn next_chars(&mut self) -> Vec<char> {
self.it.next().unwrap().chars().collect()
}
pub fn next_vec<T: FromStr>(&mut self, len: usize) -> Vec<T> {
(0..len).map(|_| self.next()).collect()
}
}
}
// ---------- end scannner ----------
use std::collections::*;
use std::io::Write;
type Map<K, V> = BTreeMap<K, V>;
type Set<T> = BTreeSet<T>;
type Deque<T> = VecDeque<T>;
fn main() {
use std::io::Read;
let mut s = String::new();
std::io::stdin().read_to_string(&mut s).unwrap();
let mut sc = scanner::Scanner::new(&s);
let out = std::io::stdout();
let mut out = std::io::BufWriter::new(out.lock());
run(&mut sc, &mut out);
}
type M = ModInt<998244353>;
fn run<W: Write>(sc: &mut scanner::Scanner, out: &mut std::io::BufWriter<W>) {
let t: u32 = sc.next();
let m: u64 = sc.next();
let p = m.factorize();
let len = p.len();
for _ in 0..t {
let n: usize = sc.next();
let b = sc.next::<M>();
let c = sc.next::<M>();
let d = sc.next::<M>();
let a: Vec<u64> = sc.next_vec(n);
let mut val = b;
let mut dp = vec![M::one(); 1 << len];
for &a in a.iter() {
if m % a == 0 {
let mut a = a;
let mut bit = 0;
for (i, (&p, &k)) in p.iter().enumerate() {
let mut c = 0;
while a % p == 0 {
a /= p;
c += 1;
}
if c == k {
bit |= 1 << i;
}
}
dp[bit] *= M::one() + val;
}
val = val * c + d;
}
bitwise_transform(&mut dp, |a, b| *b *= *a);
let mut ans = M::zero();
for (i, dp) in dp.iter().rev().enumerate() {
if i.count_ones() % 2 == 0 {
ans += *dp - M::one();
} else {
ans -= *dp - M::one();
}
}
writeln!(out, "{}", ans).ok();
}
/*
let mut v = 1u128;
for (i, p) in [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47].iter().enumerate() {
v *= p;
println!("{}: {}", v, i + 1);
}
println!("{}", 10u64.pow(18));
*/
}
הההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההההה
XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
0