(本题解选自洛谷并增加了一些内容)

暴力 $dp$ 好有道理啊。

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

显然是

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

贝祖定理：设 $a,b$ 是不全为 $0$ 的整数，则对任意整数 $x,y$，都有 $\gcd (a,b)∣ax+by$，且存在整数 $x,y$，使得 $ax+by= \gcd (a,b)$。扩展成多变量的形式仍然成立。（详情可见这里）

从这个定理出发，我们就可以得到下述定理：

对于关于 $x_1,x_2,\dots,x_n$ 的不定方程 $v_1x_1+v_2x_2+\dots+v_nx_n=w$ 有解当且仅当 $\gcd(v_1,v_2,v_3,\dots,v_n)|w$

在此处，我们先将方程转化成 $v_1x_1+v_2x_2+...+v_nx_n-w=k\times P$，移项得 $v_1x_1+v_2x_2+...+v_nx_n+(-k)P=w$，然后把 $P$ 也当成未知数就行了。

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

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

实际上，你要真暴力转移照样超时，所以 dp 还是参考另一篇的代码吧。

但是我们来反演一波吧。

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

非常显然存在：

反演可得：

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

该性质也是 dp 优化的重要性质。

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

最后我们的答案是。

我们累加一遍就好了。

复杂度很玄学，毕竟有一个 $\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;
}

题解

我们首先考虑 $n=1$ 的情况，通过打表可以发现，对于体积为 $v$ 的物品，在可以取无限次的条件下，能够取到的体积对 $P$ 取膜的值可以为

$gcd(v,P),2*gcd(v,P),3*gcd(v,P)...$

所以我们可以把这个物品的体积看做 $gcd(v,P)$，这对答案不会有任何影响。

由此我们得到启发，$O(\sqrt{P})$ 预处理出 $P$ 的每个约数，存储每个约数出现的次数。由于同一个约数多次选取对答案没有任何影响，容易得出，如果 $P$ 的一个约数出现次数为 $x$，则有 $2^x-1$ 种方案选取这个约数。

继续打表我们又可以发现，如果选取了 $P$ 的两个约数 $a_{1},a_{2}$，能够取到的值膜 $P$ 的值可以为

$gcd(a_{1},a_{2}),2*gcd(a_{1},a_{2}),3*gcd(a_{1},a_{2})...$

我们可以口胡证明，当选取约数个数大于 2 时，上式依然成立。

由此我们可以写出DP方程，设 $F[i][j]$ 表示选到 $P$ 的第 $i$ 个约数，选出的数的 gcd 值为 $P$ 的第 $j$ 个约数的方案数，则转移方程为：

$F[i][j]=F[i-1][j]+(1+\sum{_{gcd(a[k],a[i])==a[j]}F[i-1][k]})*(2^{s[i]}-1)$

其中 $a[i]$ 表示 $P$ 的第 $i$ 个约数，$s[i]$ 表示第 $i$ 个约数的出现次数。

该部分时间复杂度为 $O(M^2logM)$，空间复杂度为 $O(M^2)$ 或 $O(2M)$ （滚动数组），其中 $M$ 为 $P$ 的约数个数。

这个时候，对于每一个询问 $w_{i}$，我们容易得出答案就是

$\sum{_{a[i]|w_{i}}F[n][i]}$

但是询问数量可以达到 $10^6$ ，$M$ 最大可以达到 $10^3$ 以上，$O(qM)$ 的暴力统计依然会TLE。

事实上我们可以发现，如果一个数既是 $P$ 的约数，又是 $w_{i}$ 的约数，那么它一定是 $gcd(P,w_{i})$ 的约数，因此我们只需要统计

$\sum{_{a[i]|gcd(P,w_{i})}F[n][i]}$

而这可以在DP之后就预处理出来：

$G[i]=\sum{_{a[j]|a[i]}F[n][j]}$

因此我们可以 $O(1)$ 回答询问，这样总复杂度就是 $O(\sqrt{P}+M^2logM+q)$，这道题就可以AC啦

代码如下：

#include<bits/stdc++.h>
using namespace std;

const int N=1000010,M=10050,mod=1000000007;
int n,q,P,v[N],num[M],tot[M],cnt=0,f[2][M],g[M],now=0,sum[N];

inline void init(){
    sum[1]=2;
    for(register int i=2;i<=n;i++)sum[i]=sum[i-1]*2%mod;
    for(register int i=1;i<=n;i++)sum[i]=(sum[i]+mod-1)%mod;
    
    for(register int i=1;i<=sqrt(P);i++)if(P%i==0)num[++cnt]=i;
    for(register int i=cnt;i>1;i--)if(P/num[i]!=num[i])num[++cnt]=P/num[i];
    for(register int i=1;i<=n;i++){
        int pos=lower_bound(num+1,num+cnt+1,v[i])-num;
        tot[pos]++;
    }
    
    for(register int i=1;i<=cnt;i++)if(tot[i]){
        now^=1;
        for(register int j=1;j<=cnt;j++)f[now][j]=f[now^1][j];
        for(register int j=1;j<=cnt;j++)if(f[now^1][j]){
            int nxt=__gcd(num[j],num[i]);
            int pos=lower_bound(num+1,num+cnt+1,nxt)-num;
            (f[now][pos]+=1LL*f[now^1][j]*sum[tot[i]]%mod)%=mod;
        }
        (f[now][i]+=sum[tot[i]])%=mod;
    }
    for(register int i=1;i<=cnt;i++){
    	for(register int j=1;j<=i;j++)if(num[i]%num[j]==0){
            (g[i]+=f[now][j])%=mod;
        }
    }
}

int main(){
    scanf("%d%d%d",&n,&q,&P);
    for(register int i=1;i<=n;i++)scanf("%d",&v[i]),v[i]=__gcd(v[i],P);
    init();
    for(register int i=1;i<=q;i++){
        int x,ans=0;scanf("%d",&x);x=__gcd(x,P);
        int pos=lower_bound(num+1,num+cnt+1,x)-num;
        printf("%d\n",g[pos]);
    }
    return 0;
}