题目大意
给出一段序列，求其中最大公约数为1的四元组的个数。
思路
我们要用到反演、正难则反的思想。对于每一个大于1的数字(x)，求出最大公约数为(x)的四元组的个数(g(x))，然后用排列中所有四元组的组合个数减去(sum g(x))即可。
直接求(g(x))没有什么思路，但是求公约数中存在(x)的四元组的个数(f(x))会比较容易。枚举约数中存在x的数列元素的个数(n)，则有
[f(x)=C_n^4
]
那么怎么把(f(x))变为(g(x))呢？这要用到莫比乌斯反演。
莫比乌斯反演
莫比乌斯函数
[mu(x)=
egin{cases}
1 & ext{若$x$=1}\
0 & ext{若对$x$质因数分解得到的每个质数的次数中存在大于1的}\
(-1)^k & ext{$k$为$x$的质因数个数}
end{cases}
]
莫比乌斯反演公式
若
[f(x)=sum_{k=1}^{lfloor frac{N}{x}
floor}g(kx) ag{1}
]
则
[g(x)=sum_{k=1}^{lfloor frac{N}{x}
floor}f(kx)mu(k)
]
我们发现这道题若把序列中的数字最大值作为(N),(f(x),g(x))恰好满足该关系(1)。于是我们跟着公式求(g(x))即可。
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;
#define ll long long
const int MAX_N = 10010, MAX_R = 5, MAX_PRIME_CNT = MAX_N;
ll C[MAX_N][MAX_R];
int Num[MAX_N], Mu[MAX_N];
ll F[MAX_N];
void GetMu(int *mu, int n)
{
static bool NotPrime[MAX_N];
static int prime[MAX_PRIME_CNT];
memset(NotPrime, false, sizeof(NotPrime));
int primeCnt = 0;
mu[1] = 1;
for (int i = 2; i <= n; i++)
{
if (!NotPrime[i])
{
prime[primeCnt++] = i;
mu[i] = -1;
}
for (int j = 0; j < primeCnt; j++)
{
if (i*prime[j] > n)
break;
NotPrime[i*prime[j]] = true;
if (i%prime[j] == 0)
{
mu[i*prime[j]] = 0;
break;
}
else
mu[i*prime[j]] = -mu[i];
}
}
}
void GetC(int r, int n)
{
memset(C, 0, sizeof(C));
for (int i = 1; i <= n; i++)
{
C[i][0] = 1;
for (int j = 1; j <= min(i, r); j++)
{
if (i == j)
C[i][j] = 1;
else
C[i][j] = C[i - 1][j - 1] + C[i - 1][j];
}
}
}
ll Proceed(int maxN, int n)
{
memset(F, 0, sizeof(F));
for (int i = 2; i <= maxN; i++)
{
int cnt = 0;
for (int j = 1; j <= maxN / i; j++)
cnt += Num[i * j];
F[i] = C[cnt][4];
}
ll ans = 0;
for (int i = 2; i <= maxN; i++)
for (int j = 1; j <= maxN / i; j++)
ans += F[i * j] * Mu[j];
return C[n][4] - ans;
}
int main()
{
GetMu(Mu, 10000);
GetC(4, 10000);
int n, maxN = 0, x;
while (~scanf("%d", &n))
{
memset(Num, 0, sizeof(Num));
for (int i = 1; i <= n; i++)
{
scanf("%d", &x);
Num[x]++;
maxN = max(maxN, x);
}
printf("%lld
", Proceed(maxN, n));
}
return 0;
}