这是一篇需要龟速乘的思博题解
我们考虑一下(nleq 30)的睿智暴力，显然质因数个数少得一批，互质的条件又等价于没有公共的质因子，所以我们直接状压质因子，(dp[i][S_1][S_2])表示前(i)个数考虑完了，第一个人选择的质因子状态为(S_1)，第二个人选择的质因子状态为(S_2)，转移的话就考虑当前这个数是给第一个人，第二个人还是两个人都不要就好了，显然可以滚动数组优化，又因为(S_1)和(S_2)不能有交，于是可以直接枚举补集
复杂度大概是(O(n3^{pi(n)}))，(pi(n))表示不大于(n)的质数个数
我们发现(nleq 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
", ans);
return 0;
}