[HAOI2018]奇怪的背包

时间：2023-03-10 06:29:19
[HAOI2018]奇怪的背包

题目

暴力\(dp\)好有道理啊

于是我们来个反演吧

考虑一个体积序列\(\{v_1,v_2,...v_n\}\)能凑成\(w\)的条件

显然是

\[v_1x_1+v_2x_2+...+v_nx_n\equiv w(mod\ P)
\]

根据贝祖定理，我们知道上面的同余方程有解的条件是

\[gcd(v_1,v_2...v_n,P)|w
\]

现在题目转化成了求有多少个子集满足\(gcd(v_1,v_2..v_n,P)|w\)了

显然一个暴力\(dp\)记录一下当前\(gcd\)转移就好了，由于\(gcd\)必然是\(P\)的约数，\(P\)的约数是\(\sqrt{P}\)级别，我们我们只存这些有用的状态，之后每次暴力转移就好了，复杂度\(O(n\sqrt{P})\)就没了

但是我们来反演一波吧

设\(F(n)\)表示有多少个子集的\(gcd\)是\(n\)的倍数，\(f(n)\)表示有多少个子集的\(gcd\)是\(n\)

非常显然存在

\[F(n)=\sum_{n|d}f(d)
\]

反演可得

\[f(n)=\sum_{n|d}F(d)\mu(\frac{d}{n})
\]

\(F(n)\)我们随便搞一搞就好了，若有\(s\)个数都有\(n\)这个约数，那么\(F(n)=2^s-1\)

之后反演得到\(f\)，这里直接枚举\(P\)的约数，复杂度是\(O(\sigma^2(P))\)，但是我们还需要求一个\(\mu\)，配合线筛预处理还是挺快的

最后我们的答案是

\[Ans(n)=\sum_{d|n}f(d)
\]

我们累加一遍就好了

复杂度很玄学，毕竟有一个\(\mu\)需要求，可能是\(O(\sigma^2(P)\sqrt{P})\)

代码

#include<algorithm>

#include<iostream>

#include<cstring>

#include<cstdio>

#include<cmath>

#define re register

#define LL long long

#define max(a,b) ((a)>(b)?(a):(b))

#define min(a,b) ((a)<(b)?(a):(b))

const int maxn=1e6+5;

const int mod=1e9+7;

const int M=2e5+5;

inline int read() {

    char c=getchar();int x=0;while(c<'0'||c>'9') c=getchar();

    while(c>='0'&&c<='9') x=(x<<3)+(x<<1)+c-48,c=getchar();return x;

}

int d[M],F[M],ans[M],f[M],g[M];

int pw[maxn],a[maxn];

int n,m,P,Sqr,tot;

int is[M],p[M>>1],mu[M],id2[M],id1[M];

inline int find(int x) {return (x<=Sqr)?id1[x]:id2[P/x];}

inline int gcd(int a,int b) {return !b?a:gcd(b,a%b);}

inline void Init() {

    is[1]=1;mu[1]=1;

    for(re int i=2;i<=Sqr;i++) {

        if(!is[i]) p[++p[0]]=i,mu[i]=-1;

        for(re int j=1;j<=p[0]&&p[j]*i<=Sqr;j++) {

            is[p[j]*i]=1;

            if(i%p[j]==0) break;

            mu[p[j]*i]=-1*mu[i];

        }

    }

    for(re int i=1;i<=tot;i++)

    if(d[i]<=Sqr) id1[d[i]]=i;

        else id2[P/d[i]]=i;

}

inline int getMu(int x) {

    if(x<=Sqr) return mu[x];

    int now=1,t=x;

    for(re int i=1;p[i]*p[i]<=t&&i<=p[0];i++)

    if(x%p[i]==0) {

        x/=p[i];

        if(x%p[i]==0) return 0;

        now=-1*now;

        if(x==1) break;

    }

    if(x!=1) now=-1*now;

    return now;

}

int main() {

    n=read(),m=read(),P=read();

    for(re int i=1;i*i<=P;i++)

    if(P%i==0) {

        d[++tot]=i;

        if(i!=P/i) d[++tot]=P/i;

    }

    pw[0]=1;

    for(re int i=1;i<=n;i++) pw[i]=(pw[i-1]+pw[i-1])%mod;

    for(re int i=1;i<=n;i++) a[i]=read();

    std::sort(d+1,d+tot+1);

    Sqr=std::sqrt(P);Init();

    for(re int i=1;i<=n;i++) {

        int k=gcd(a[i],P);

        g[find(k)]++;

    }

    for(re int i=1;i<=tot;i++)

        for(re int j=i;j<=tot;j++) {

            if(d[j]%d[i]) continue;

            F[i]+=g[j];

        }

    for(re int i=1;i<=tot;i++) F[i]=pw[F[i]]-1;

    for(re int i=1;i<=tot;i++)

        for(re int j=i;j<=tot;j++) {

            if(d[j]%d[i]) continue;

            f[i]+=F[j]*getMu(d[j]/d[i]);

            f[i]%=mod;

            f[i]=(f[i]+mod)%mod;

        }

    for(re int i=1;i<=tot;i++)

        for(re int j=i;j<=tot;j++) {

            if(d[j]%d[i]) continue;

            ans[j]=(ans[j]+f[i])%mod;

        }

    for(re int i=1;i<=m;i++) {

        int x=read();

        int k=gcd(x,P);

        printf("%d\n",ans[find(k)]);

    }

    return 0;

}

相关文章