Luogu

BZOJ

分析

本思路来源于这篇题解


$\mathrm{DP}$ + 矩阵快速幂。

显然答案为和为 $p$ 的倍数的方案数 - 和为 $p$ 的倍数且不含质数的方案数

先算前面那一部分。

设 $f[i][j]$ 表示选 $i$ 个数,和 $\%p$ 为 $j$ 的方案数。显然有 $f[i][j]=\sum\limits_{k=0}^{p-1}f[i-1][k]\times cnt[(j-k+p)\%p]$ 。这里的 $cnt[i]$ 表示 $1\sim m$ 中 $\%p$ 为 $i$ 的个数。

再算后面那一部分。

设 $g[i][j]$ 表示选 $i$ 个合数,和 $\%p$ 为 $j$ 的方案数。显然有 $g[i][j]=\sum\limits_{k=0}^{p-1}g[i-1][k]\times compo[(j-k+p)\%p]$ 。这里的 $compo[i]$ 表示 $1\sim m$ 中的合数 $\%p$ 为 $i$ 的个数。

答案为 $f[n][0]-g[n][0]$ 。这样就可以获得 $20$ 分的好成绩。


可以发现, $f$ 和 $g$ 可以用矩阵快速幂来转移。

于是就获得了 $100$ 分。

另外,这个转移矩阵是循环矩阵,所以其实可以做到 $O(n^2)$ 的矩阵乘法。 然而没有必要


具体实现及细节见代码。

另外,不知道为什么在 $\mathrm{BZOJ}$ 上会 $\mathrm{TLE}$ 。

代码

// =================================
//   author: M_sea
//   website: http://m-sea-blog.com/
// =================================
#include <algorithm>
#include <iostream>
#include <cstdlib>
#include <cstring>
#include <cstdio>
#include <cmath>
#define re register
typedef long long ll;
using namespace std;

inline int read() {
    int X=0,w=1; char c=getchar();
    while (c<'0'||c>'9') { if (c=='-') w=-1; c=getchar(); }
    while (c>='0'&&c<='9') X=X*10+c-'0',c=getchar();
    return X*w;
}

const int N=100+10;
const int M=20000000+10;
const int mod=20170408;

int cnt[M],compo[M];
int pr[M],v[M],tot=0;

inline void sieve(int n) {
    v[1]=1;
    for (re int i=2;i<=n;++i) {
        if (!v[i]) pr[++tot]=i;
        for (re int j=1;j<=tot&&i*pr[j]<=n;++j) {
            v[i*pr[j]]=1;
            if (i%pr[j]==0) break;
        }
    }
}

struct Matrix {
    int n,m;
    int s[N][N];

    Matrix(int n_=0,int m_=0):m(m_),n(n_) { }
    int* operator [](int i) { return s[i]; }
} P,Q,V,W;

Matrix operator *(Matrix a,Matrix b) {
    Matrix c(a.n,b.m);
    for (re int i=1;i<=c.n;++i)
        for (re int j=1;j<=c.m;++j) {
            ll s=0;
            for (re int k=1;k<=a.m;++k) s+=1ll*a[i][k]*b[k][j];
            c[i][j]=s%mod;
        }
    return c;
}

inline Matrix qpow(Matrix a,int b) {
    Matrix c(a.n,a.m);
    for (re int i=1;i<=c.n;++i) c[i][i]=1;
    for (;b;b>>=1,a=a*a)
        if (b&1) c=c*a;
    return c;
}

int main() {
    int n=read(),m=read(),p=read();

    for (re int i=0;i<p;++i) cnt[i]=m/p;
    for (re int i=1;i<=m%p;++i) ++cnt[i];

    sieve(m);
    for (re int i=1;i<=m;++i)
        if (v[i]) ++compo[i%p];

    P.n=P.m=p;
    P[1][1]=cnt[0];
    for (re int i=2;i<=p;++i) P[1][i]=cnt[p-i+1];
    for (re int i=2;i<=p;++i) {
        for (re int j=2;j<=p;++j) P[i][j]=P[i-1][j-1];
        P[i][1]=P[i-1][p];
    }

    Q.n=Q.m=p;
    Q[1][1]=compo[0];
    for (re int i=2;i<=p;++i) Q[1][i]=compo[p-i+1];
    for (re int i=2;i<=p;++i) {
        for (re int j=2;j<=p;++j) Q[i][j]=Q[i-1][j-1];
        Q[i][1]=Q[i-1][p];
    }

    V.n=p,V.m=1;
    for (re int i=1;i<=p;++i) V[i][1]=cnt[i-1];

    W.n=p,W.m=1;
    for (re int i=1;i<=p;++i) W[i][1]=compo[i-1]%mod;

    P=qpow(P,n-1)*V,Q=qpow(Q,n-1)*W;
    printf("%d\n",(P[1][1]-Q[1][1]+mod)%mod);
    return 0;
}
最后修改:2019 年 09 月 26 日 01 : 17 PM