给出M个不同的数，每个数有无穷多个。从中选N个组成环，问最大公约数为1的方案有多少种。
上一题是这题的基础，不同的是每个数可以选多个，tot[i]表示m个数中i的倍数有多少个，那么用构成长度为l最大公约数为i的倍数的方案有tot[i]^l种，之后再用容斥求一下最大公约数为1的方案即可。
另外要考虑循环同构，套一个Burnside就行了，坑的就是N是10007倍数的情况，WA了N次。。要将MOD*n。。
1 #include <stdio.h> 2 #include <string.h> 3 #include <algorithm> 4 #define MAXN 20010 5 int MOD; 6 using namespace std; 7 typedef long long LL; 8 int cas, n, m, tx; 9 int pri[MAXN], notp[MAXN], flag[MAXN], prs; 10 int tot[MAXN], maxm; 11 int fac[MAXN][305], facs[MAXN], phi[MAXN * 10]; 12 void init(){ 13 for (int i = 2; i < MAXN; i++) if (!notp[i]) { 14 flag[i] = 1; 15 for (int j = i*2; j < MAXN; j += i) { 16 if (!notp[j]) notp[j] = flag[j] = 1; 17 else if (flag[j]) flag[j]++; 18 if (j%(i*i) == 0) flag[j] = 0; 19 } 20 } 21 for (int i = 2; i < MAXN; i++) 22 if (!notp[i]) pri[prs++] = i; 23 for (int i = 1; i < MAXN; i++) 24 for (int j = 1; j * j <= i; j++) if (i % j == 0) { 25 fac[i][facs[i]++] = j; 26 if (j * j != i)fac[i][facs[i]++] = i/j; 27 } 28 const int maxn = 100000; 29 for (int i = 1; i <= maxn; i ++) phi[i] = i; 30 for (int i = 2; i <= maxn; i ++) if(phi[i] == i) { 31 for(int j = i; j <= maxn; j += i) 32 phi[j] = phi[j] / i * (i - 1); 33 } 34 } 35 int powmod(int x, int d){ 36 int ans = 1, tmp = x; 37 for (; d; d >>= 1, tmp = (LL)tmp * tmp % MOD) 38 if (d&1) ans = (LL)ans * tmp % MOD; 39 return ans; 40 } 41 int cal(int l){ 42 int ans = powmod(tot[1], l); 43 for (int i = 2; i <= maxm; i++) { 44 if (flag[i] == 0) continue; 45 if (flag[i] & 1) ans = (ans - powmod(tot[i], l)) % MOD; 46 else ans = (ans + powmod(tot[i], l)) % MOD; 47 } 48 return (ans % MOD + MOD) % MOD; 49 } 50 void solve(){ 51 int ans = 0; 52 for(int i = 1; i * i <= n; i ++) if(n % i == 0) { 53 ans = (ans + (LL)cal(i) * phi[n / i]) % MOD; 54 if(i * i != n) 55 ans = (ans + (LL)cal(n / i) * phi[i]) % MOD; 56 } 57 printf("%d\n", (ans + MOD) % MOD / n); 58 } 59 int main(){ 60 //freopen("test.in","r",stdin); 61 init(); 62 scanf("%d", &cas); 63 while (cas--) { 64 scanf("%d%d", &m, &n); 65 MOD = 10007 * n; 66 memset(tot, 0, sizeof tot); 67 maxm = 0; 68 for (int i = 1; i <= m; i++) { 69 scanf("%d", &tx); 70 maxm = max(tx, maxm); 71 for (int j = 0; j < facs[tx]; j++) 72 tot[fac[tx][j]]++; 73 } 74 solve(); 75 } 76 return 0; 77 }