https://www.luogu.com.cn/problem/P3702

稍微有点套路......

先不管那个质数的条件

考虑\(DP\)，令\(f_{i,j}\)表示\(i\)个数和\(\%p\)为\(j\)的方案数。

有转移

\[f_{i,c}=\sum\limits_{(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}=\sum\limits_{(a+b) \% p=c} f_{i,a}f_{i,b} \]

进一步的如果我们知道\(f_i,f_j\)，可以推到

\[f_{i+j,c}=\sum\limits_{(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\n", (ans1[0] - ans2[0] + mod) % mod);
  return 0;
}