[NOI2015] 寿司晚宴

题意：

有两人从$[2,n]$的正整数里分别选出两个集合，求这两个集合互质的方案数，对p取模。

定义集合$A,B$互质为$forall ain A,bin B,(a,b)=1$。

$nleq 500,pleq 10^9$。

题解：

显然得状压dp，但是500以内有100个质数，不能直接做。

不过数论题一般都可以对质数根号分治一下，把$>sqrt{n}$的大质数拿出来单独做。

注意到如果按大质数排序，那么对于一种合法方案，每个大质数相同的连通块内都只能有一个人选数。

于是我们就可以考虑把每个数拆成两部分：所有$leq sqrt{n}$的质数（一共8个）状态的压缩，和一个大质数（没有就是1）。

考虑dp，由于没法按正常dp那样算完两个状态再回头加，于是令$dp(i,j,k)$表示考虑前i个数，第1/2个人的小质数状态是j/k的方案数。

如果大质数是1，那么直接转移，$dp(i,j|s_i ,k)+=dp(i-1,j,k),dp(i,j,k|s_i )+=dp(i-1,j,k)$。

如果大质数不是1，那么需要将两个人分开转移。

我们设辅助数组$g(j,k),h(j,k)$表示在当前块内只转移第1/2个人，且两人状态是j/k的方案数。

那么直接$g(j|s_i ,k)+=g(j,k),h(j,k|s_i )+=h(j,k)$，最后令$dp(i,j,k)=dp(las_i ,j,k)+g(j,k)+h(j,k)$即可。

这个dp本质上是个背包，所以如果压了一位状态则需要倒着枚举背包大小。

复杂度$O(2^{16}n)$。

套路：

• 数论题$ ightarrow$把$>sqrt{n}$的大质数拿出来单独做。

代码：

#include<bits/stdc++.h>
#define maxn 505
#define maxm 500005
#define inf 0x7fffffff
#define ll long long
#define rint register int
#define debug(x) cerr<<#x<<": "<<x<<endl
#define fgx cerr<<"--------------"<<endl
#define dgx cerr<<"=============="<<endl

using namespace std;
int mod,p[maxn],ish[maxn],las[maxn],f[maxn][256][256];
int g[256][256],h[256][256];
struct element{
    int sta,pri;
    bool operator<(const element b)const{return pri<b.pri;}
}A[maxn];

inline int read(){
    int x=0,f=1; char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c=='-') f=-1;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    return x*f;
}

inline int mo(int x){return x>=mod?x-mod:x;}
inline void init(int n){
    for(int i=2;i<=n;i++){
        if(!ish[i]) p[++p[0]]=i,las[i]=i;
        for(int j=1;j<=p[0];j++){
            if(p[j]*i>n) break;
            ish[p[j]*i]=1,las[p[j]*i]=p[j];
            if(i%p[j]==0) break;
        }
    }
}

int main(){
    int n=read(); mod=read(),init(n);
    for(int i=2;i<=n;i++){
        int x=i;
        for(int j=1;j<=min(p[0],8);j++)
            if(i%p[j]==0){
                A[i].sta|=(1<<(j-1));
                while(x%p[j]==0) x/=p[j];
            }
        A[i].pri=x;
    }
    sort(A+2,A+1+n),f[1][0][0]=1;
    int ans=0;
    for(int i=2,j;i<=n;i=j+1){
        j=i; while(A[j].pri==A[j+1].pri) j++;
        int x=A[i].pri; 
        if(x==1){
            for(int k=i;k<=j;k++){
                int s=A[k].sta;
                for(int s1=255;s1>=0;s1--)
                    for(int s2=255;s2>=0;s2--){
                        f[k][s1][s2]=mo(f[k][s1][s2]+f[k-1][s1][s2]);
                        f[k][s1|s][s2]=mo(f[k][s1|s][s2]+f[k-1][s1][s2]); 
                        f[k][s1][s2|s]=mo(f[k][s1][s2|s]+f[k-1][s1][s2]);
                    }
            }
        }
        else{
            memset(g,0,sizeof(g));
            memset(h,0,sizeof(h));
            for(int s1=255;s1>=0;s1--)
                for(int s2=255;s2>=0;s2--)
                    g[s1][s2]=h[s1][s2]=f[i-1][s1][s2];
            for(int k=i;k<=j;k++){
                int s=A[k].sta;
                for(int s1=255;s1>=0;s1--)
                    for(int s2=255;s2>=0;s2--){
                        g[s1|s][s2]=mo(g[s1|s][s2]+g[s1][s2]);
                        h[s1][s2|s]=mo(h[s1][s2|s]+h[s1][s2]);
                    }
            }
            for(int s1=255;s1>=0;s1--)
                for(int s2=255;s2>=0;s2--){
                    f[j][s1][s2]=mo(mo(g[s1][s2]+h[s1][s2])+mod-f[i-1][s1][s2]);
                }
        }
    }
    for(int s1=0;s1<256;s1++)
        for(int s2=0;s2<256;s2++)
            if((s1&s2)==0)
                ans+=f[n][s1][s2],ans%=mod;
    printf("%d
",ans);
    return 0;
}