蓝书（算法竞赛进阶指南）刷题记录——POJ3904 Sky Code（容斥原理）

最新推荐文章于 2024-06-13 20:53:49 发布

原创最新推荐文章于 2024-06-13 20:53:49 发布 · 227 阅读

1 ·

CC 4.0 BY-SA版权

本文深入解析了POJ3904数论题，通过转化问题求解四元组数量，利用容斥原理和Mobius函数简化计算，提供了完整的代码实现。

题目：POJ3904.
题目大意：给定 $n$ 个数 $a_i$ ，求有多少个四元组 $(w, x, y, z)$ 满足 $gcd(a_w,a_x,a_y,a_z)=1$ .
$1\leq n,a_i\leq 10^4$ .

~~我容斥果然还是太菜了，这种题都要看题解…~~

首先 $g c d (a, b, c, d) = 1$ 不代表两两互质.

考虑把问题转化为用所有四元组的数量 $C_{n}^{4}$ 减掉 $g c d$ 不为 $1$ 的四元组数量.

考虑如何求出后者，我们统计 $c n t [i]$ 表示有多少个数是 $i$ 的倍数，那么有 $C_{cnt[i]}^{4}$ 个四元组的 $g c d$ 是 $i$ 的倍数.

想到这一步就开始考虑容斥，很容易发现 $i$ 的容斥系数即为Mobius函数 $\mu(i)$ ，线性筛预处理即可，时间复杂度 $O(n\log n)$ .

代码如下：

#include<iostream>
#include<cstdio>
#include<algorithm>
  using namespace std;
  
#define Abigail inline void
typedef long long LL;

const int N=10000;

LL C(LL n){return n<4?0:n*(n-1)*(n-2)*(n-3)/24;}

int pr[N+9],tp,b[N+9];
LL mu[N+9];

void sieve(int n){
  for (int i=2;i<=n;++i) b[i]=1;
  mu[1]=1;
  for (int i=2;i<=n;++i){
  	if (b[i]) pr[++tp]=i,mu[i]=-1;
  	for (int j=1;j<=tp&&i*pr[j]<=n;++j){      //把n写成tp调了好久... 
	  b[i*pr[j]]=0;
	  if (i%pr[j]==0){mu[i*pr[j]]=0;break;}
	  mu[i*pr[j]]=-mu[i];
	}
  }
}

int n,mx;
LL c[N+9],cnt[N+9],ans;

Abigail start(){
  sieve(N);
}

Abigail into(){
  int a;
  for (int i=1;i<=n;++i){
    scanf("%d",&a);
    mx=max(mx,a);
    ++cnt[a];
  }
}

Abigail work(){
  for (int i=1;i<=mx;++i)
    for (int j=i;j<=mx;j+=i)
      c[i]+=cnt[j];
  for (int i=1;i<=mx;++i)
    ans+=mu[i]*C(c[i]);
}

Abigail outo(){
  printf("%lld\n",ans);
  for (int i=1;i<=mx;++i)
    cnt[i]=c[i]=0;
  mx=ans=0;
}

int main(){
  start();
  while (~scanf("%d",&n)){
    into();
    work();
    outo();
  }
  return 0;
}