「NOI 2015」寿司晚宴

题目链接：https://www.lydsy.com/JudgeOnline/problem.php?id=4197

---

Description

有 $n-1$ 个寿司，第 $i$ 个寿司的美味度为 $i+1$。小 $\text{G}$ 和小 $\text{W}$ 每人选择一些寿司来品尝。规定一种方案为不和谐的当且仅当：小 $\text{G}$ 和小 $\text{W}$ 品尝的寿司中分别存在美味度为 $x$ 和 $y$ 的寿司，且 $x$ 和 $y$ 不互质。求一共有多少种方案是和谐的。

$2\leqslant n\leqslant 500$

---

Solution

问题的本质为：两个人在 $2\sim n$ 中选取数字，第一个人和第二个人选择的数字的质因子集合分别为 $S_1$ 和 $S_2$，满足 $S_1\cap S_2=\varnothing$ 的选择方案数。

$30\ \text{pts}\quad n\leqslant 30$

  注意到 $30$ 以内的质数只有 $10$ 个，考虑使用状压 $\text{DP}$ 计算答案。可以将每个数分解质因数并压成二进制数 $S$（其中 $S$ 的第 $i$ 位表示是否有第 $i$ 个质因子）。
  设 $f[i][S_1][S_2]$ 表示考虑到第 $i$ 个数，第一个人选择的数的质因子集合为 $S_1$，第二个人选择的数的质因子为 $S_2$，有如下状态转移方程：
  $f[i][S_1|k][S_2]+=f[i-1][S_1][S_2]$（$k\cap S_2=\varnothing$）
  $f[i][S_1][S_2|k]+=f[i-1][S_1][S_2]$（$k\cap S_1=\varnothing$）
  其中第一维可以压缩，答案为 $\displaystyle\sum_{S_1\cap S_2=\varnothing}f[n][S_1][S_2]$

$100\ \text{pts}\quad n\leqslant 500$

  按照之前的状态压缩显然不可行，考虑如何将状态继续压缩。
  注意每个数字 $x$ 至多只有一个大于 $\sqrt{x}$ 的质因子，因此我们可以将不大于 $\sqrt{x}$ 的质因子进行 $\text{DP}$，单独考虑大于 $\sqrt{x}$ 的质因子。
  小于 $\sqrt{x}_{\max}=\sqrt{500}\approx 22.4$ 的质因子只有 $8$ 个，于是我们对这些质因子按照一般的状压 $\text{DP}$ 进行计算。
  考虑大于 $\sqrt{x}$ 的质因子 $p$，如果其中一个人选择了含有质因子 $p$ 的某个数，那么另一个人就不能选择任何含有质因子 $p$ 的数字。设 $g[k][p][S_1][S_2]$ 表示第 $k$ 个人选择了质因子 $p$，第一个人和第二个人选择的小于 $\sqrt{n}$ 的质因子集合分别为 $S_1$ 和 $S_2$，转移显然。
  在计算完 $g$ 数组后，再更新 $f$ 的答案。考虑如何将两个数组合并：（稍有常识的人就能发现）
  $f[S_1][S_2]=g[0][S_1][S_2]+g[1][S_1][S_2]-f[S_1][S_2]-f[S_1][S_2]$
最后减去 $f[S_1][S_2]$ 是因为两个人都不选这个集合的数的情况被统计了 $2$ 次。

  时间复杂度：$\Theta(n \cdot 2^{16})$

---

Code

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

const int N=505;
const int P=8;
const int prime[P]={2,3,5,7,11,13,17,19};
int n,p,f[1<<P][1<<P],g[2][1<<P][1<<P];
pair<int,int> s[N];

int main() {
    scanf("%d%d",&n,&p);
    for(int i=2;i<=n;++i) {
        int tmp=i;
        for(int j=0;j<8;++j) while(tmp%prime[j]==0) s[i].second|=(1<<j),tmp/=prime[j];
        s[i].first=tmp;
    }
    std::sort(s+1,s+n+1);
    f[0][0]=1;
    for(int i=2;i<=n;++i) {
        if(i==2||s[i].first==1||s[i].first!=s[i-1].first) {
            memcpy(g[0],f,sizeof(f));
            memcpy(g[1],f,sizeof(f));
        }
        for(int s1=(1<<P)-1;s1>=0;--s1)
            for(int s2=(1<<P)-1;s2>=0;--s2) {
                if(s1&s2) continue;
                if((s2&s[i].second)==0) (g[0][s1|s[i].second][s2]+=g[0][s1][s2])%=p;
                if((s1&s[i].second)==0) (g[1][s1][s2|s[i].second]+=g[1][s1][s2])%=p;
            }
        if(i==n||s[i].first==1||s[i].first!=s[i+1].first) {
            for(int s1=(1<<P)-1;s1>=0;--s1)
                for(int s2=(1<<P)-1;s2>=0;--s2)
                    if((s1&s2)==0) f[s1][s2]=((g[0][s1][s2]+g[1][s1][s2]-f[s1][s2])%p+p)%p;
        }
    }
    int ans=0;
    for(int i=(1<<P)-1;i>=0;--i)
        for(int j=(1<<P)-1;j>=0;--j)
            if((i&j)==0) (ans+=f[i][j])%=p;
    printf("%d\n",ans);
    return 0;
}