好题呀。

首先两者所拿的集合必然是不相交的。首先我们考虑30分的暴力做法，30以内质数很少，我们可以进行一波状压，用\(dp[i][j]\)表示第一个人取集合i，第二个人取集合j的方案数。我们记录当前数字的质因子集合为k，之后就可以进行转移了。

但是数据范围到了500，就会导致无法状压那么多的质数，我们考虑500以内的数，最多只可能有1个大于22的质因子。

这就意味着我们完全可以把它单独提出来处理，之后剩下的8个质数我们继续进行状压。

具体的操作方法就是，我们首先把所有的数按照最大质因子进行排序，相同的排在一起，因为相同的两者是不能同时选的。对于每个最大质因子，我们用\(f1[][]\)表示这个最大质因子让第一个人选取，\(f2[][]\)表示这个最大质因子让第二个人选取。

之后在小范围内进行状压DP即可。

记录当前元素的质因子集合为k，转移方程是:

\(f1[S_1 | k][S_2] += f1[S_1][S_2]\)\((k \& S_2 == 0)\)

\(f2[S_1][S_2 | k] += f2[S_1][S_2]\)\((k \& S_1 == 0)\)

之后合并答案的时候，\(dp[i][j] = f1[i][j] + f2[i][j] - dp[i][j]\)，其中减去一次是考虑到两个人都没有选择这个最大质因子，这种情况被计算了两次。

时间复杂度\(O(n2^{16})\)

#include
    
     #define rep(i,a,n) for(int i = a;i <= n;i++)#define per(i,n,a) for(int i = n;i >= a;i--)#define enter putchar('\n')#define fr friend inline#define y1 pojusing namespace std;typedef long long ll;const int M = 1005;const int INF = 1000000009;const double eps = 1e-7;int read(){   int ans = 0,op = 1;char ch = getchar();   while(ch < '0' || ch > '9') {if(ch == '-') op = -1;ch = getchar();}   while(ch >= '0' && ch <= '9') ans = ans * 10 + ch - '0',ch = getchar();   return ans * op;}int n,P,p[9] = {0,2,3,5,7,11,13,17,19};int tot,dp[300][300],f1[300][300],f2[300][300],ans;bool np[M];struct node{   int val,mprime,S;   bool operator < (const node &g){return mprime < g.mprime;}}a[M];int inc(int a,int b){return (a+b) % P;}void get(node &x){   int cur = 0,t = x.val;   rep(i,1,8)   {      if(t % p[i]) continue;      cur |= (1 << (i-1));      while(!(t % p[i])) t /= p[i];   }   x.S = cur,x.mprime = (t > 1) ? t : -1;}int main(){   n = read(),P = read();   rep(i,2,n) a[i-1].val = i,get(a[i-1]);   sort(a+1,a+n),dp[0][0] = 1;   rep(i,1,n-1)   {      if(i == 1 || a[i].mprime != a[i-1].mprime || a[i].mprime == -1) memcpy(f1,dp,sizeof(f1)),memcpy(f2,dp,sizeof(f2));      per(j,255,0)      per(k,255,0)      {     if(j & k) continue;     if(!(a[i].S & j)) f2[j][k | a[i].S] = inc(f2[j][k | a[i].S],f2[j][k]);     if(!(a[i].S & k)) f1[j | a[i].S][k] = inc(f1[j | a[i].S][k],f1[j][k]);      }      if(i == n-1 || a[i].mprime != a[i+1].mprime || a[i].mprime == -1)      {     rep(j,0,255)        rep(k,0,255) if(!(j & k)) dp[j][k] = inc(f1[j][k],inc(f2[j][k],P - dp[j][k]));      }   }   rep(i,0,255)      rep(j,0,255) if(!(i&j) && dp[i][j]) ans = inc(ans,dp[i][j]);   printf("%d\n",ans);   return 0;}