Skip to content

Commit c52d324

Browse files
koba-e964koba-e964
committed
Add yukicoder/3572-ex.rs yukicoder/3572.rs yukicoder/3573.rs
1 parent 3bc869b commit c52d324

3 files changed

Lines changed: 494 additions & 0 deletions

File tree

‎yukicoder/3572-ex.rs‎

Lines changed: 153 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,153 @@
1+
/// Verified by https://atcoder.jp/contests/abc198/submissions/21774342
2+
mod mod_int {
3+
use std::ops::*;
4+
pub trait Mod: Copy { fn m() -> i64; }
5+
#[derive(Copy, Clone, Hash, PartialEq, Eq, PartialOrd, Ord)]
6+
pub struct ModInt<M> { pub x: i64, phantom: ::std::marker::PhantomData<M> }
7+
impl<M: Mod> ModInt<M> {
8+
// x >= 0
9+
pub fn new(x: i64) -> Self { ModInt::new_internal(x % M::m()) }
10+
fn new_internal(x: i64) -> Self {
11+
ModInt { x: x, phantom: ::std::marker::PhantomData }
12+
}
13+
pub fn pow(self, mut e: i64) -> Self {
14+
debug_assert!(e >= 0);
15+
let mut sum = ModInt::new_internal(1);
16+
let mut cur = self;
17+
while e > 0 {
18+
if e % 2 != 0 { sum *= cur; }
19+
cur *= cur;
20+
e /= 2;
21+
}
22+
sum
23+
}
24+
#[allow(dead_code)]
25+
pub fn inv(self) -> Self { self.pow(M::m() - 2) }
26+
}
27+
impl<M: Mod> Default for ModInt<M> {
28+
fn default() -> Self { Self::new_internal(0) }
29+
}
30+
impl<M: Mod, T: Into<ModInt<M>>> Add<T> for ModInt<M> {
31+
type Output = Self;
32+
fn add(self, other: T) -> Self {
33+
let other = other.into();
34+
let mut sum = self.x + other.x;
35+
if sum >= M::m() { sum -= M::m(); }
36+
ModInt::new_internal(sum)
37+
}
38+
}
39+
impl<M: Mod, T: Into<ModInt<M>>> Sub<T> for ModInt<M> {
40+
type Output = Self;
41+
fn sub(self, other: T) -> Self {
42+
let other = other.into();
43+
let mut sum = self.x - other.x;
44+
if sum < 0 { sum += M::m(); }
45+
ModInt::new_internal(sum)
46+
}
47+
}
48+
impl<M: Mod, T: Into<ModInt<M>>> Mul<T> for ModInt<M> {
49+
type Output = Self;
50+
fn mul(self, other: T) -> Self { ModInt::new(self.x * other.into().x % M::m()) }
51+
}
52+
impl<M: Mod, T: Into<ModInt<M>>> AddAssign<T> for ModInt<M> {
53+
fn add_assign(&mut self, other: T) { *self = *self + other; }
54+
}
55+
impl<M: Mod, T: Into<ModInt<M>>> SubAssign<T> for ModInt<M> {
56+
fn sub_assign(&mut self, other: T) { *self = *self - other; }
57+
}
58+
impl<M: Mod, T: Into<ModInt<M>>> MulAssign<T> for ModInt<M> {
59+
fn mul_assign(&mut self, other: T) { *self = *self * other; }
60+
}
61+
impl<M: Mod> Neg for ModInt<M> {
62+
type Output = Self;
63+
fn neg(self) -> Self { ModInt::new(0) - self }
64+
}
65+
impl<M> ::std::fmt::Display for ModInt<M> {
66+
fn fmt(&self, f: &mut ::std::fmt::Formatter) -> ::std::fmt::Result {
67+
self.x.fmt(f)
68+
}
69+
}
70+
impl<M: Mod> ::std::fmt::Debug for ModInt<M> {
71+
fn fmt(&self, f: &mut ::std::fmt::Formatter) -> ::std::fmt::Result {
72+
let (mut a, mut b, _) = red(self.x, M::m());
73+
if b < 0 {
74+
a = -a;
75+
b = -b;
76+
}
77+
write!(f, "{}/{}", a, b)
78+
}
79+
}
80+
impl<M: Mod> From<i64> for ModInt<M> {
81+
fn from(x: i64) -> Self { Self::new(x) }
82+
}
83+
// Finds the simplest fraction x/y congruent to r mod p.
84+
// The return value (x, y, z) satisfies x = y * r + z * p.
85+
fn red(r: i64, p: i64) -> (i64, i64, i64) {
86+
if r.abs() <= 10000 {
87+
return (r, 1, 0);
88+
}
89+
let mut nxt_r = p % r;
90+
let mut q = p / r;
91+
if 2 * nxt_r >= r {
92+
nxt_r -= r;
93+
q += 1;
94+
}
95+
if 2 * nxt_r <= -r {
96+
nxt_r += r;
97+
q -= 1;
98+
}
99+
let (x, z, y) = red(nxt_r, r);
100+
(x, y - q * z, z)
101+
}
102+
} // mod mod_int
103+
104+
macro_rules! define_mod {
105+
($struct_name: ident, $modulo: expr) => {
106+
#[derive(Copy, Clone, PartialEq, Eq, PartialOrd, Ord, Hash)]
107+
pub struct $struct_name {}
108+
impl mod_int::Mod for $struct_name { fn m() -> i64 { $modulo } }
109+
}
110+
}
111+
const MOD: i64 = 998_244_353;
112+
define_mod!(P, MOD);
113+
type MInt = mod_int::ModInt<P>;
114+
115+
fn c_part(d: usize, n: usize) -> Vec<MInt> {
116+
let mut dp = vec![MInt::new(0); n + 1];
117+
let mut parts = vec![1];
118+
if d % 2 == 0 {
119+
parts.push(1);
120+
for _ in 0..d / 2 - 1 {
121+
parts.push(2);
122+
}
123+
} else {
124+
for _ in 0..d / 2 {
125+
parts.push(2);
126+
}
127+
}
128+
dp[d] += 1;
129+
for &part in &parts {
130+
let mut ep = dp.clone();
131+
for i in d + 1..n + 1 {
132+
if i >= part {
133+
let val = ep[i - part];
134+
ep[i] += val;
135+
}
136+
}
137+
dp = ep;
138+
}
139+
dp
140+
}
141+
142+
fn main() {
143+
let n = 11;
144+
let mut dp = vec![MInt::new(0); n + 1];
145+
for i in 1..n + 1 {
146+
let sub = c_part(i, n);
147+
eprintln!("{i} => {sub:?}");
148+
for j in 0..n + 1 {
149+
dp[j] += sub[j];
150+
}
151+
}
152+
eprintln!("{dp:?}");
153+
}

‎yukicoder/3572.rs‎

Lines changed: 191 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,191 @@
1+
fn getline() -> String {
2+
let mut ret = String::new();
3+
std::io::stdin().read_line(&mut ret).unwrap();
4+
ret
5+
}
6+
7+
/// Verified by https://atcoder.jp/contests/abc198/submissions/21774342
8+
mod mod_int {
9+
use std::ops::*;
10+
pub trait Mod: Copy { fn m() -> i64; }
11+
#[derive(Copy, Clone, Hash, PartialEq, Eq, PartialOrd, Ord)]
12+
pub struct ModInt<M> { pub x: i64, phantom: ::std::marker::PhantomData<M> }
13+
impl<M: Mod> ModInt<M> {
14+
// x >= 0
15+
pub fn new(x: i64) -> Self { ModInt::new_internal(x % M::m()) }
16+
fn new_internal(x: i64) -> Self {
17+
ModInt { x: x, phantom: ::std::marker::PhantomData }
18+
}
19+
pub fn pow(self, mut e: i64) -> Self {
20+
debug_assert!(e >= 0);
21+
let mut sum = ModInt::new_internal(1);
22+
let mut cur = self;
23+
while e > 0 {
24+
if e % 2 != 0 { sum *= cur; }
25+
cur *= cur;
26+
e /= 2;
27+
}
28+
sum
29+
}
30+
#[allow(dead_code)]
31+
pub fn inv(self) -> Self { self.pow(M::m() - 2) }
32+
}
33+
impl<M: Mod> Default for ModInt<M> {
34+
fn default() -> Self { Self::new_internal(0) }
35+
}
36+
impl<M: Mod, T: Into<ModInt<M>>> Add<T> for ModInt<M> {
37+
type Output = Self;
38+
fn add(self, other: T) -> Self {
39+
let other = other.into();
40+
let mut sum = self.x + other.x;
41+
if sum >= M::m() { sum -= M::m(); }
42+
ModInt::new_internal(sum)
43+
}
44+
}
45+
impl<M: Mod, T: Into<ModInt<M>>> Sub<T> for ModInt<M> {
46+
type Output = Self;
47+
fn sub(self, other: T) -> Self {
48+
let other = other.into();
49+
let mut sum = self.x - other.x;
50+
if sum < 0 { sum += M::m(); }
51+
ModInt::new_internal(sum)
52+
}
53+
}
54+
impl<M: Mod, T: Into<ModInt<M>>> Mul<T> for ModInt<M> {
55+
type Output = Self;
56+
fn mul(self, other: T) -> Self { ModInt::new(self.x * other.into().x % M::m()) }
57+
}
58+
impl<M: Mod, T: Into<ModInt<M>>> AddAssign<T> for ModInt<M> {
59+
fn add_assign(&mut self, other: T) { *self = *self + other; }
60+
}
61+
impl<M: Mod, T: Into<ModInt<M>>> SubAssign<T> for ModInt<M> {
62+
fn sub_assign(&mut self, other: T) { *self = *self - other; }
63+
}
64+
impl<M: Mod, T: Into<ModInt<M>>> MulAssign<T> for ModInt<M> {
65+
fn mul_assign(&mut self, other: T) { *self = *self * other; }
66+
}
67+
impl<M: Mod> Neg for ModInt<M> {
68+
type Output = Self;
69+
fn neg(self) -> Self { ModInt::new(0) - self }
70+
}
71+
impl<M> ::std::fmt::Display for ModInt<M> {
72+
fn fmt(&self, f: &mut ::std::fmt::Formatter) -> ::std::fmt::Result {
73+
self.x.fmt(f)
74+
}
75+
}
76+
impl<M: Mod> ::std::fmt::Debug for ModInt<M> {
77+
fn fmt(&self, f: &mut ::std::fmt::Formatter) -> ::std::fmt::Result {
78+
let (mut a, mut b, _) = red(self.x, M::m());
79+
if b < 0 {
80+
a = -a;
81+
b = -b;
82+
}
83+
write!(f, "{}/{}", a, b)
84+
}
85+
}
86+
impl<M: Mod> From<i64> for ModInt<M> {
87+
fn from(x: i64) -> Self { Self::new(x) }
88+
}
89+
// Finds the simplest fraction x/y congruent to r mod p.
90+
// The return value (x, y, z) satisfies x = y * r + z * p.
91+
fn red(r: i64, p: i64) -> (i64, i64, i64) {
92+
if r.abs() <= 10000 {
93+
return (r, 1, 0);
94+
}
95+
let mut nxt_r = p % r;
96+
let mut q = p / r;
97+
if 2 * nxt_r >= r {
98+
nxt_r -= r;
99+
q += 1;
100+
}
101+
if 2 * nxt_r <= -r {
102+
nxt_r += r;
103+
q -= 1;
104+
}
105+
let (x, z, y) = red(nxt_r, r);
106+
(x, y - q * z, z)
107+
}
108+
} // mod mod_int
109+
110+
macro_rules! define_mod {
111+
($struct_name: ident, $modulo: expr) => {
112+
#[derive(Copy, Clone, PartialEq, Eq, PartialOrd, Ord, Hash)]
113+
pub struct $struct_name {}
114+
impl mod_int::Mod for $struct_name { fn m() -> i64 { $modulo } }
115+
}
116+
}
117+
const MOD: i64 = 998_244_353;
118+
define_mod!(P, MOD);
119+
type MInt = mod_int::ModInt<P>;
120+
121+
fn convolution(a: &[MInt], b: &[MInt]) -> Vec<MInt> {
122+
if a.is_empty() || b.is_empty() {
123+
return vec![];
124+
}
125+
let n = a.len() - 1;
126+
let m = b.len() - 1;
127+
let mut ans = vec![MInt::new(0); n + m + 1];
128+
for i in 0..n + 1 {
129+
for j in 0..m + 1 {
130+
ans[i + j] += a[i] * b[j];
131+
}
132+
}
133+
ans
134+
}
135+
136+
// Finds [x^n] p(x)/q(x)
137+
// Ref: https://qiita.com/ryuhe1/items/da5acbcce4ac1911f47a
138+
// Verified by: https://atcoder.jp/contests/tdpc/submissions/24583334
139+
// Depends on: MInt.rs
140+
fn bostan_mori(p: &[MInt], q: &[MInt], mut n: i64) -> MInt {
141+
if p.is_empty() {
142+
return 0.into();
143+
}
144+
assert!(p.len() < q.len());
145+
let mut p = p.to_vec();
146+
let mut q = q.to_vec();
147+
while n > 0 {
148+
let mut qn = q.clone();
149+
for i in 0..qn.len() {
150+
if i % 2 == 1 {
151+
qn[i] = -qn[i];
152+
}
153+
}
154+
let num = convolution(&p, &qn);
155+
let den = convolution(&q, &qn);
156+
let mut nxt_p = vec![MInt::new(0); q.len() - 1];
157+
let mut nxt_q = vec![MInt::new(0); q.len()];
158+
for i in 0..q.len() - 1 {
159+
let to = 2 * i + (n % 2) as usize;
160+
if to < num.len() {
161+
nxt_p[i] = num[to];
162+
}
163+
}
164+
for i in 0..q.len() {
165+
nxt_q[i] = den[2 * i];
166+
}
167+
p = nxt_p;
168+
q = nxt_q;
169+
n /= 2;
170+
}
171+
p[0] * q[0].inv()
172+
}
173+
174+
// https://yukicoder.me/problems/no/3572 (3)
175+
// Solved with hints
176+
// 根の集合は、あるdに対して1のd乗根全体、あるいは0、あるいはそれらのunionでなければならない。
177+
// 元々の問題の答えを a(n) とし、0を考えないときの答えを b(n) とすると、 a(n) = \sum_{1<=k<=n}b(k) + 1 が成立する。
178+
// b(n) は d ごとに分けると理解しやすい。b(n) のうち、根の集合が1のd乗根全体であるものを c(d,n) と呼ぶ。
179+
// x^d-1 のQ上の因数分解の次数の多重集合を f(d) とすると、 c(d,n) は c(d,d) = 1 から開始して、 f(d) でナップサック数え上げをして得られる数列である。
180+
// -> ChatGPT に聞いたら問題文を誤読していることがわかった。係数は実数なので、f(d) としてもR上の因数分解の次数を見るべきである。
181+
// こうなると f(d) の要素は 1,2 のみであり、 b(n) の母関数は簡単に計算でき (x+x^2)/(1-x-2x^2+2x^3) である。
182+
// これにより a(n) の母関数も計算でき、 (1-x^2+2x^3)/(1-2x-x^2+4x^3-2x^4) である。
183+
fn main() {
184+
let n = getline().trim().parse::<i64>().unwrap();
185+
let p = vec![MInt::new(1), MInt::new(0), -MInt::new(1), MInt::new(2)];
186+
let q = vec![MInt::new(1), -MInt::new(2), -MInt::new(1), MInt::new(4), -MInt::new(2)];
187+
for k in 1..10 {
188+
eprintln!("{k} => {}", bostan_mori(&p, &q, k));
189+
}
190+
println!("{}", bostan_mori(&p, &q, n));
191+
}

0 commit comments

Comments
 (0)