[NOI2015]寿司晚宴
这是一篇需要龟速乘的思博题解
我们考虑一下\(n\leq 30\)的睿智暴力，显然质因数个数少得一批，互质的条件又等价于没有公共的质因子，所以我们直接状压质因子，\(dp[i][S_1][S_2]\)表示前\(i\)个数考虑完了，第一个人选择的质因子状态为\(S_1\)，第二个人选择的质因子状态为\(S_2\)，转移的话就考虑当前这个数是给第一个人，第二个人还是两个人都不要就好了，显然可以滚动数组优化，又因为\(S_1\)和\(S_2\)不能有交，于是可以直接枚举补集
复杂度大概是\(O(n3^{\pi(n)})\)，\(\pi(n)\)表示不大于\(n\)的质数个数
我们发现\(n\leq 500\)时这个暴力不适用的原因是质因子个数太多了，但是考虑到一个非常小学生的结论，一个数\(n\)最多只会有一个超过\(\sqrt{n}\)的质因子，于是我们可以考虑单独讨论这个超过\(\sqrt{n}\)的质因子
状态还是暴力\(dp\)的状态，但是我们只压不超过\(\sqrt{n}\)的质因子，我们把所有数分解之后按照剩下的那个质因子排序，对于相同的质因子我们一起考虑，对于一个具有相同超过根号质因子的区间，我们发现这里面的数不能同时分给两个人，这样这个两个人就会具有一个相同的质因子了
对于一段最大值因子相同的区间，我们再做一个\(dp\)，设\(f[S]\)表示在这个区间里选择的数小于根号的质因子状态为\(f\)的方案数，做完这个\(dp\)后我们就考虑把这些\(f\)分给第一个人第二个人还是都不选，这里的转移还是一样，也是需要保证小于\(\sqrt{n}\)的质因子没有交，这里显然需要一个乘法，模数又是\(10^{10}\)级别的，于是需要龟速乘一下
复杂度大概是\(O(n2^{\pi(\sqrt{n})}+(\pi(n)-\pi(\sqrt{n}))5^{\pi(\sqrt{n})})\)
有点长的代码
#include <bits/stdc++.h>
#define re register
#define LL long long
int n, len, o, l, is[505], p[505];
LL dp[2][257][257], f[505][257], mod;
struct Number {int res, s;}a[505];
inline int cmp(const Number &A, const Number &B) { return A.res < B.res; }
inline void solve(int L, int R) {
memset(f[L-1], 0, sizeof(f[L-1]));
f[L - 1][0] = 1;
for (re int i = L; i <= R; i++)
for (re int j = 0; j <= len; j++) {
f[i][j | a[i].s] = (f[i][j | a[i].s] + f[i - 1][j]) % mod;
f[i][j] = (f[i][j] + f[i - 1][j]) % mod;
}
f[R][0] = (f[R][0] - 1 + mod) % mod;
}
inline LL mul(LL a, LL b) {
LL S = 0;
for (; b; b >>= 1ll, a = (a + a) % mod)
if (b & 1ll) S = (S + a) % mod;
return S;
}
int main() {
scanf("%d%lld", &n, &mod);
is[1] = 1;
for (re int i = 2; i <= n; i++) {
if (!is[i]) p[++p[0]] = i;
for (re int j = 1; j <= p[0] && p[j] * i <= n; j++) {
is[p[j] * i] = 1;
if (i % p[j] == 0) break;
}
}
for (re int i = 2; i <= n; i++) {
int x = i, s = 0;
for (re int j = 1; j <= p[0] && p[j] * p[j] <= n; j++) {
if (x % p[j]) continue;
s |= 1 << (j - 1);
while (x % p[j] == 0) x /= p[j];
}
a[i - 1].res = x, a[i - 1].s = s;
}
for (re int i = 1; i <= p[0]; i++)
if (p[i] * p[i] <= n) ++len;
len = (1 << len) - 1;
std::sort(a + 1, a + n, cmp);
dp[0][0][0] = 1;
for (re int i = 1; i < n; i++) {
if (a[i].res > 1) {l = i;break;}
memset(dp[o ^ 1], 0, sizeof(dp[o ^ 1]));
for (re int j = 0; j <= len; j++) {
int S = len ^ j;
for (re int t = S; t; t = (t - 1) & S) {
if (!dp[o][j][t]) continue;
dp[o ^ 1][j][t] = (dp[o ^ 1][j][t] + dp[o][j][t]) % mod;
if (((j | a[i].s) & t) == 0)
dp[o ^ 1][j | a[i].s][t] = (dp[o ^ 1][j | a[i].s][t] + dp[o][j][t]) % mod;
if ((j & (t | a[i].s)) == 0)
dp[o ^ 1][j][t | a[i].s] = (dp[o ^ 1][j][t | a[i].s] + dp[o][j][t]) % mod;
}
dp[o ^ 1][j][0] = (dp[o ^ 1][j][0] + dp[o][j][0]) % mod;
dp[o ^ 1][j | a[i].s][0] = (dp[o ^ 1][j | a[i].s][0] + dp[o][j][0]) % mod;
if ((a[i].s & j) == 0) dp[o ^ 1][j][a[i].s] = (dp[o ^ 1][j][a[i].s] + dp[o][j][0]) % mod;
}
o ^= 1;
}
for (re int i = l + 1; i <= n; i++) {
if (a[i].res == a[i - 1].res) continue;
solve(l, i - 1);
memset(dp[o ^ 1], 0, sizeof(dp[o ^ 1]));
for(re int k = 0; k <= len; k++) {
if(!f[i - 1][k]) continue;
for (re int j = 0; j <= len; j++) {
int S = len ^ j;
for (re int t = S; t; t = (t - 1) & S) {
if (!dp[o][j][t]) continue;
if (((j | k) & t) == 0)
dp[o ^ 1][j | k][t] = (dp[o ^ 1][j | k][t] + mul(dp[o][j][t], f[i - 1][k])) % mod;
if ((j & (t | k)) == 0)
dp[o ^ 1][j][t | k] = (dp[o ^ 1][j][t | k] + mul(dp[o][j][t], f[i - 1][k])) % mod;
}
dp[o ^ 1][j | k][0] = (dp[o ^ 1][j | k][0] + mul(dp[o][j][0], f[i - 1][k])) % mod;
if ((j & k) == 0) dp[o ^ 1][j][k] = (dp[o ^ 1][j][k] + mul(dp[o][j][0], f[i - 1][k])) % mod;
}
}
for (re int j = 0; j <= len; j++) {
int S = len ^ j;
for (re int t = S; t; t = (t - 1) & S) {
if (!dp[o][j][t]) continue;
dp[o ^ 1][j][t] = (dp[o ^ 1][j][t] + dp[o][j][t]) % mod;
}
dp[o ^ 1][j][0] = (dp[o ^ 1][j][0] + dp[o][j][0]) % mod;
}
l = i;o ^= 1;
}
LL ans = 0;
for (re int i = 0; i <= len; i++)
for (re int j = 0; j <= len; j++)
if ((i & j) == 0) ans = (ans + dp[o][i][j]) % mod;
printf("%lld\n", ans);
return 0;
}
