4818: [Sdoi2017]序列计数
Time Limit: 30 Sec Memory Limit: 128 MB
Description
Alice想要得到一个长度为n的序列，序列中的数都是不超过m的正整数，而且这n个数的和是p的倍数。Alice还希望
，这n个数中，至少有一个数是质数。Alice想知道，有多少个序列满足她的要求。
Input
一行三个数，n,m,p。
1<=n<=10^9,1<=m<=2×10^7,1<=p<=100
Output
一行一个数，满足Alice的要求的序列数量，答案对20170408取模。
Sample Input
3 5 3
Sample Output
33
月考后日常颓水题，好久没有1A了，可能SDOI比较良心
SBDP，D[i]表示和为i的方案数，用所有的减去没有素数的，写出方程发现可以矩阵转移，然后乱敲
#include<bits/stdc++.h>
using namespace std;
template <class _T> inline void read(_T &_x) {
int _t; bool flag = false;
while ((_t = getchar()) != '-' && (_t < '' || _t > '')) ;
if (_t == '-') _t = getchar(), flag = true; _x = _t - '';
while ((_t = getchar()) >= '' && _t <= '') _x = _x * + _t - '';
if (flag) _x = -_x;
}
using namespace std;
typedef long long LL;
const int mod = ;
const int maxv = ;
struct Mat {
int n, m, a[maxv][maxv];
Mat() {}
Mat(int x, int y):n(x), m(y) {
for (register int i = , j; i < n; ++i)
for (j = ; j < m; ++j)
a[i][j] = ;
}
inline void init() {for (int i = ; i < n; ++i) a[i][i] = ; }
inline Mat operator * (Mat B) {
Mat C(n, B.m);
for (register int i = , j, k; i < n; ++i)
for (j = ; j < B.m; ++j)
for (k = ; k < m; ++k) {
C.a[i][j] += (int)((LL)a[i][k] * B.a[k][j] % mod);
if (C.a[i][j] >= mod) C.a[i][j] -= mod;
}
return C;
}
inline Mat operator ^ (int t) {
Mat res(n, m), tmp = *this; res.init();
while (t) {
if (t & ) res = res * tmp;
tmp = tmp * tmp, t >>= ;
}
return res;
}
};
const int maxp = ;
const int maxm = ;
int n, m, p;
int cnt_p[maxp], cnt_n[maxp];
bool vis[maxm];
int prime[maxm / ], pcnt;
inline void Init() {
cnt_p[ % p] = ;
for (register int i = , j; i <= m; ++i) {
if (!vis[i]) {
prime[++pcnt] = i;
} else {
++cnt_p[i % p];
}
for (j = ; j <= pcnt && i * prime[j] <= m; ++j) {
vis[i * prime[j]] = true;
if (i % prime[j] == ) break;
}
}
int tmpa = m / p, tmpb = m % p;
for (register int i = ; i < p; ++i) {
cnt_n[i] = tmpa;
cnt_n[i] += (i && i <= tmpb);
}
}
inline int getres(Mat &a) {
Mat x(p, );
x.a[][] = ;
return ((a ^ n) * x).a[][];
}
int main() {
//freopen();
//freopen();
read(n), read(m), read(p);
Init();
Mat a(p, p), b(p, p);
for (int i = , j, to; i < p; ++i) {
for (j = ; j < p; ++j) {
to = i + j;
if (to >= p) to -= p;
a.a[i][to] = cnt_n[j];
b.a[i][to] = cnt_p[j];
}
}
int ans = getres(a) - getres(b);
if (ans < ) ans += mod;
cout << ans << endl;
return ;
}