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的方案数，用所有的减去没有素数的，写出方程发现可以矩阵转移，然后乱敲

---

1 #include
       
         2 using namespace std; 3 template 
        
          inline void read(_T &_x) { 4     int _t; bool flag = false; 5     while ((_t = getchar()) != '-' && (_t < '0' || _t > '9')) ; 6     if (_t == '-') _t = getchar(), flag = true; _x = _t - '0'; 7     while ((_t = getchar()) >= '0' && _t <= '9') _x = _x * 10 + _t - '0'; 8     if (flag) _x = -_x; 9 }10 using namespace std;11 typedef long long LL;12 const int mod = 20170408;13 const int maxv = 200;14 struct Mat {15     int n, m, a[maxv][maxv];16     Mat() {}17     Mat(int x, int y):n(x), m(y) {18         for (register int i = 0, j; i < n; ++i)19             for (j = 0; j < m; ++j)20                 a[i][j] = 0;21     }22     inline void init() {
     for (int i = 0; i < n; ++i) a[i][i] = 1; }23     inline Mat operator * (Mat B) {24         Mat C(n, B.m);25         for (register int i = 0, j, k; i < n; ++i)26             for (j = 0; j < B.m; ++j)27                 for (k = 0; k < m; ++k) {28                     C.a[i][j] += (int)((LL)a[i][k] * B.a[k][j] % mod);29                     if (C.a[i][j] >= mod) C.a[i][j] -= mod;30                 }31         return C;32     }33     inline Mat operator ^ (int t) {34         Mat res(n, m), tmp = *this; res.init();35         while (t) {36             if (t & 1) res = res * tmp;37             tmp = tmp * tmp, t >>= 1;38         }39         return res;40     }41 };42 const int maxp = 210;43 const int maxm = 20000010;44 int n, m, p;45 int cnt_p[maxp], cnt_n[maxp];46 bool vis[maxm];47 int prime[maxm / 10], pcnt;48 inline void Init() {49     cnt_p[1 % p] = 1;50     for (register int i = 2, j; i <= m; ++i) {51         if (!vis[i]) {52             prime[++pcnt] = i;53         } else {54             ++cnt_p[i % p];55         }56         for (j = 1; j <= pcnt && i * prime[j] <= m; ++j) {57             vis[i * prime[j]] = true;58             if (i % prime[j] == 0) break;59         }60     }61     int tmpa = m / p, tmpb = m % p;62     for (register int i = 0; i < p; ++i) {63         cnt_n[i] = tmpa;64         cnt_n[i] += (i && i <= tmpb);65     }66 }67 inline int getres(Mat &a) {68     Mat x(p, 1);69     x.a[0][0] = 1;70     return ((a ^ n) * x).a[0][0];71 }72 int main() {73     //freopen();74     //freopen();75     read(n), read(m), read(p);76     Init();77     Mat a(p, p), b(p, p);78     for (int i = 0, j, to; i < p; ++i) {79         for (j = 0; j < p; ++j) {80             to = i + j;81             if (to >= p) to -= p;82             a.a[i][to] = cnt_n[j];83             b.a[i][to] = cnt_p[j];84         }85     }86     int ans = getres(a) - getres(b);87     if (ans < 0) ans += mod;88     cout << ans << endl;89     return 0;90 }