Alice想要得到一个长度为n的序列,序列中的数都是不超过m的正整数,而且这n个数的和是p的倍数。
Alice还希望,这n个数中,至少有一个数是质数。
Alice想知道,有多少个序列满足她的要求。
对100%的数据,1≤n≤10^9,1≤m≤2×10^7,1≤p≤100
一个比较显然的想法就是总方案减去一个素数都没有的方案
考虑怎么计算总方案:
设 $f_{i,j}$ 表示已经填了 $i$ 个数,当前和模 $p$ 的值为 $j$ 的方案数,显然有转移:
$$
f_{i,j}=\sum_{x=0}^{p-1}f_{i-1,j-x}g_x
$$
其中 $g_x$ 表示 $1 \sim m$ 中有多少个数字模 $p$ 等于 $x$,这个可以直接 $O(m)$ 的扫一遍后的得出
然后就只需要做一个素数都没有的方案了
发现这个只需要修改一下 $g_x$ 就行了,设 $g'_x$ 表示 $1 \sim m$ 中有多少个数字模 $p$ 等于 $x$ 且不是素数
这个直接线性筛素数后修改 $g$ 数组就行了
然后套个矩阵快速幂优化转移就行了


1 // luogu-judger-enable-o2 2 #include <bits/stdc++.h> 3 using namespace std; 4 typedef long long ll; 5 6 const int mod = 20170408, N = 100; 7 8 ll n, m, p, ans; 9 10 struct Mat { 11 int a[N][N]; 12 Mat() { memset(a, 0, sizeof a); } 13 int *operator [] (int x) { return a[x]; } 14 Mat operator * (Mat b) { 15 Mat c; 16 for(int i = 0 ; i < p ; ++ i) 17 for(int j = 0 ; j < p ; ++ j) 18 for(int k = 0 ; k < p ; ++ k) 19 c[i][j] = (c[i][j] + 1ll * a[i][k] * b[k][j] % mod) % mod; 20 return c; 21 } 22 }; 23 24 Mat pw(Mat a, ll b) { 25 Mat r; 26 27 for(int i = 0 ; i < p ; ++ i) r[i][i] = 1; 28 29 for( ; b ; b >>= 1, a = a * a) 30 if(b & 1) 31 r = r * a; 32 return r; 33 } 34 35 const int M = 2e7 + 10; 36 ll pri[1270607 + 10], tot, cnt[N]; 37 bool vis[M]; 38 39 int main() { 40 cin >> n >> m >> p; 41 42 for(int i = 2 ; i <= m ; ++ i) { 43 if(!vis[i]) { 44 pri[++ tot] = i; 45 } 46 for(int j = 1 ; j <= tot && i * pri[j] <= m ; ++ j) { 47 vis[i * pri[j]] = 1; 48 if(i % pri[j] == 0) break; 49 } 50 } 51 52 for(int i = 1 ; i <= m ; ++ i) 53 ++ cnt[i % p]; 54 55 Mat base, coe; 56 base[0][0] = 1; 57 58 for(int i = 0 ; i < p ; ++ i) { 59 for(int j = 0 ; j < p ; ++ j) { 60 coe[i][((i - j) % p + p) % p] = cnt[j]; 61 } 62 } 63 ans = (base * pw(coe, n))[0][0]; 64 65 for(int i = 1 ; i <= tot ; ++ i) cnt[pri[i] % p] --; 66 67 for(int i = 0 ; i < p ; ++ i) { 68 for(int j = 0 ; j < p ; ++ j) { 69 coe[i][((i - j) % p + p) % p] = cnt[j]; 70 } 71 } 72 ans -= (base * pw(coe, n))[0][0]; 73 74 // cout << "tot = " << tot << endl; 75 76 cout << (ans % mod + mod) % mod << endl; 77 }