P5366 [SNOI2017]遗失的答案 题解

题目描述

给定 \(n,g,l\) ， \([1,n]\) 中的每个数可选可不选，但是要求所有选择的数的 \(\gcd\) 为 \(g\) ， \(\text{lcm}\) 为 \(l\) 。

\(q\) 次询问，每次给定 \(x\) ，要求 \(x\) 必选，求方案数，对 \(10^9+7\) 取模。

数据范围

• \(1\le q\le 10^5\) 。
• \(1\le n,g,l,x\le10^8\) 。

时间限制 \(\texttt{1s}\) ，空间限制 \(\texttt{250MB}\) 。

分析

特判掉 \(g\) 不整除 \(l\) 的情况，否则可以把 \(g\) 视为 \(1\) 。

注意到 \(\max\limits_{1\le n\le 10^8}\omega(n)=8,\max\limits_{1\le n\le 10^8}d(n)=768\) ，因此 \(q\) 很大没有意义，只有不超过 \(768\) 个数可能被选。

先翻译条件， \(\forall p\mid l\) ，存在 \(a,b\) 被选满足 \(v_p(a)=0,v_p(b)=v_p(l)\) 。

记 \(m=\omega(l)\le 8\) ，将 \(l\) 的每个因子 \(x\) 看成一个长为 \(2m\) 的二进制数，记录 \([v_p(x)=0]\) 和 \([v_p(x)=v_p(l)]\) 的信息。

至此，问题转化为，有 \(d(l)\) 个 \(2m\) 位二进制数 \(a_1,\cdots,a_{d(l)}\) ，对每个数，求有多少种方案钦定它被选，并且选择的数按位或等于 \(2^{2m}-1\) 。

接下来有两条路可走。

FWT

如果没有钦定必选的限制，显然答案为 \([z^{2^{2m}-1}]\prod(1+z^{a_i})\) ，这里乘法定义为或卷积。

对于钦定 \(x\) 必选的限制，维护 \(1+z^{a_i}\) 的前缀后缀积，再单独乘上 \(z^{a_i}\) 即可。

时间复杂度 \(\mathcal O(d(l)\cdot m\cdot 2^{2m})\) 。

#include<bits/stdc++.h>
#define poly array<int,1<<16>
using namespace std;
const int mod=1e9+7;
int g,k,l,m,n,q,x;
int p[10],v[10],id[1<<16],res[1<<16];
poly a[800],pre[800],suf[800];
void get_prime(int n)
{
    for(int i=2;i*i<=n;i++)
    {
        if(n%i) continue;
        while(n%i==0) n/=i,v[m]++;
        p[m++]=i;
    }
    if(n!=1) v[m]=1,p[m++]=n;
}
int get_num(int x)
{
    if(x>n) return -1;
    int s=0;
    for(int i=0,cnt=0;i<m;i++)
    {
        for(cnt=0;x%p[i]==0;cnt++,x/=p[i]) ;
        if(!cnt) s|=1<<i;
        if(cnt==v[i]) s|=1<<(m+i); 
    }
    return s;
}
int add(int x,int y)
{
    if((x+=y)>=mod) x-=mod;
    return x;
}
int dec(int x,int y)
{
    if((x-=y)<0) x+=mod;
    return x;
}
void fwt(poly &f,int n,int op)
{
    for(int k=2,m=1;k<=n;k<<=1,m<<=1)
        for(int i=0;i<n;i+=k)
            for(int j=i;j<i+m;j++)
                f[j+m]=op==1?add(f[j+m],f[j]):dec(f[j+m],f[j]);
}
poly operator*(poly f,poly g)
{
    for(int i=0;i<1<<(2*m);i++) f[i]=1ll*f[i]*g[i]%mod;
    return f;
}
void add(int x)
{
    if(!~(x=get_num(x))) return ;
    id[x]=++k,a[k][0]++,a[k][x]++;
    fwt(a[k],1<<(2*m),1),pre[k]=pre[k-1]*a[k];
}
int main()
{
    scanf("%d%d%d%d",&n,&g,&l,&q),n/=g;
    if(l%g)
    {
        while(q--) printf("0\n");
        return 0;
    }
    get_prime(l/=g),pre[0][0]=1,fwt(pre[0],1<<(2*m),1);
    for(int i=1;i*i<=l;i++)
    {
        if(l%i) continue;
        add(i);
        if(i*i!=l) add(l/i);
    }
    suf[k+1]=pre[0];
    for(int i=k;i>=1;i--) suf[i]=suf[i+1]*a[i];
    memset(res,-1,sizeof(res));
    while(q--)
    {
        scanf("%d",&x);
        if(x%g||l%(x/=g)||!~(x=get_num(x))) printf("0\n");
        else if(res[x]!=-1) printf("%d\n",res[x]);
        else
        {
            poly tmp={0};tmp[x]=1,fwt(tmp,1<<(2*m),1);
            tmp=pre[id[x]-1]*tmp*suf[id[x]+1],fwt(tmp,1<<(2*m),-1);
            printf("%d\n",res[x]=tmp[(1<<(2*m))-1]);
        }
    }
    return 0;
}

子集反演

记 \(cnt_S\) 为数组 \(a\) 中 \(S\) 的出现次数， \(f_S\) 为选若干个数按位与为 \(S\) 的方案数， \(g_S\) 表示选若干个数按位与为 \(S\) 的子集的方案数。

子集反演：

\[g_S=\sum_{T\subseteq S}f_T\Rightarrow f_S=\sum_{T\subseteq S}(-1)^{|S|-|T|}g_T\\ \]

如果不考虑 \(x\) 的限制，显然 \(g_S=2^{\sum\limits_{T\subseteq S}cnt_T}\) ，并且 \(f_{2^{2m}-1}\) 即为所求。

现在强制 \(x\) 被选，首先不是 \(x\) 的超集的集合不再产生贡献，其次对于所有 \(x\) 的超集 \(S\) ， \(g_S\) 要除以 \(2\) 。

枚举 \(x\) 的超集暴力统计答案即可。

时间复杂度 \(\mathcal O(3^m)\) 。

#include<bits/stdc++.h>
using namespace std;
const int maxn=1<<16,mod=1e9+7;
int g,l,m,n,q,x,all;
int p[maxn],v[maxn];
int pw[maxn],sz[maxn],ans[maxn];
void get_prime(int n)
{
    for(int i=2;i*i<=n;i++)
    {
        if(n%i) continue;
        while(n%i==0) v[m]++,n/=i;
        p[m++]=i;
    }
    if(n!=1) v[m]=1,p[m++]=n;
}
int get_num(int x)
{
    int res=0;
    for(int i=0;i<m;i++)
    {
        int cnt=0;
        while(x%p[i]==0) cnt++,x/=p[i];
        if(!cnt) res|=1<<i;
        if(cnt==v[i]) res|=1<<(m+i); 
    }
    return res;
}
int solve(int x)
{
    if(x>n||x%g) return 0;///注意不要漏判x>n的情况
    x/=g;
    if(l%x) return 0;
    return ans[get_num(x)];
}
int main()
{
    scanf("%d%d%d%d",&n,&g,&l,&q);
    if(l%g)
    {
        while(q--) printf("0\n");
        return 0;
    }
    l/=g,get_prime(l);
    for(int i=1;i<=n/g&&i*i<=l;i++)
    {
        if(l%i) continue;
        sz[get_num(i)]++;
        if(i*i!=l&&l/i<=n/g) sz[get_num(l/i)]++;
    }
    pw[0]=1,all=1<<(2*m);
    for(int i=1;i<=768;i++) pw[i]=2*pw[i-1]%mod;
    for(int i=0;i<2*m;i++)
        for(int s=0;s<all;s++)
            if(s>>i&1) sz[s]+=sz[s^(1<<i)];
    for(int s=0;s<all;s++)
        for(int i=s;;i=(i+1)|s)
        {
            ans[s]=(ans[s]+(__builtin_parity(i)?mod-1ll:1)*pw[sz[i]-1])%mod;
            if(i==all-1) break;
        }
    while(q--) scanf("%d",&x),printf("%d\n",solve(x));
    return 0;
}