https://www.luogu.com.cn/problem/P3702
稍微有点套路......
先不管那个质数的条件
考虑(DP),令(f_{i,j})表示(i)个数和(\%p)为(j)的方案数。
有转移
[f_{i,c}=sumlimits_{(a+b) \% p = c} f_{i-1,a}f_{i-1,b}
]
但这样转移(i)就有(10^9)个...
注意到如果我们没必要从(i-1)转移到(i)
如果知道了(f_i),那么我们可以推到(f_{2i}),即
[f_{2i,c}=sumlimits_{(a+b) \% p=c} f_{i,a}f_{i,b}
]
进一步的如果我们知道(f_i,f_j),可以推到
[f_{i+j,c}=sumlimits_{(a+b) \% p = c}f_{i,a}f_{j,b}
]
可以像快速幂那样倍增转移,边界就是(f_{1,i}=1...m)中(\% p = i)的数的个数。
然后发现这样选出来的不一定合法,减掉没有素数构成的序列个数就好了(这个可以和上面类似的(DP))。
所以这个(p)可以搞到(10^5)的样子然后FFT转移的吧......
但这里暴力卷积就好了....
复杂度(O(p^2 log n))
#include <bits/stdc++.h>
using namespace std;
const int N = 210, mod = 20170408;
int p;
struct poly {
int a[N];
poly() { memset(a, 0, sizeof(a)); }
inline int &operator [] (int i) { return a[i]; }
inline int operator [] (int i) const { return a[i]; }
poly operator * (const poly &o) const {
static int A[N];
for (int i = 0; i <= p * 2; i++) A[i] = 0;
for (int i = 0; i < p; i++)
for (int j = 0; j < p; j++) {
int k = (i+j) % p; A[k] += 1ll * a[i] * o[j] % mod;
if (A[k] >= mod) A[k] -= mod;
}
poly ans;
for (int i = 0; i < p; i++) ans[i] = A[i];
return ans;
}
};
poly qpow(poly A, int n) {
poly ans; ans[0] = 1;
for (; n; n >>= 1, A = A * A)
if (n & 1) ans = ans * A;
return ans;
}
vector<int> ps;
bool vis[20000010];
void sieve(int n) {
for (int i = 2; i <= n; i++) {
if (!vis[i]) ps.push_back(i);
for (auto j : ps) {
if (i * j > n) break;
vis[i*j] = 1; if (i % j == 0) break;
}
}
}
int n, m;
poly A, B, ans1, ans2;
int main() {
scanf("%d%d%d", &n, &m, &p);
sieve(m);
for (int i = 1; i <= m; i++) ++A[i % p], ++B[i % p];
for (auto pr : ps) --B[pr % p];
ans1 = qpow(A, n), ans2 = qpow(B, n);
printf("%d
", (ans1[0] - ans2[0] + mod) % mod);
return 0;
}