poj3904&&SP4191 Sky Code

guapi2333

于 2021-06-12 17:19:22 发布

阅读量75

点赞数

分类专栏：数论

本文链接：https://blog.csdn.net/guapi2333/article/details/117851602

版权

数论专栏收录该内容

25 篇文章 3 订阅

订阅专栏

这篇博客介绍了如何结合容斥原理和Möbius函数来解决一类数论问题，具体涉及计算给定数组中满足最大公约数为1的四元组数量。首先通过埃氏筛法计算Möbius函数，然后通过组合数学公式得出答案。代码实现中展示了C++的解决方案，对于每个数i，计算其在数组中出现次数scnt[i]，并应用容斥原理求解最终答案。

摘要由CSDN通过智能技术生成

Label

容斥原理+Möbius函数解决数论相关问题

Description

多组询问：给定 $n$ 个数 $a_i(n,a_i\leq 10^4)$ ，试求 $\sum_{a,b,c,d\in a}[gcd(a,b,c,d)==1]$

Solution

有了解决p3455的基础，我们不难推导此题的解法。

设 $d_i$ 为 $a$ 中被 $d$ 整除的数的个数，我们利用类似埃氏筛的方法求出 $d_1\sim d_n$ ，接下来直接套用p3455的容斥思路与公式可得：
$ans=\sum_{i=1}^{max\{a_i\}}\mu(i)C_{d_i}^{4}$

Code

#include<cstdio>
#include<cstring>
#include<iostream>
#define ri register int
#define ll long long
using namespace std;

const int MAXN=1e4;
int n,miu[MAXN+20],prime[MAXN],cnt,a[MAXN+20],maxn,scnt[MAXN+20];
ll ans;
bool notprime[MAXN+20],vis[MAXN+20];

void Mobius()
{
	miu[1]=1,notprime[1]=true;
	for(ri i=2;i<=MAXN;++i)
	{
		if(!notprime[i]) prime[++cnt]=i,miu[i]=-1;
		for(ri j=1;j<=cnt&&i*prime[j]<=MAXN;++j)
		{
			notprime[i*prime[j]]=true;
			if(i%prime[j]==0) break;
			else miu[i*prime[j]]=-miu[i];
		}
	}	
}

int main()
{
	Mobius();	
	while(~scanf("%d",&n))
	{
		memset(scnt,0,sizeof(scnt)); 
		memset(vis,0,sizeof(vis));
		maxn=0,ans=0LL;
		for(ri i=1;i<=n;++i)
		{
			scanf("%d",&a[i]);
			maxn=max(maxn,a[i]); vis[a[i]]=1;
		}
		scnt[1]=n;
		for(ri i=2;i<=maxn;++i)
			for(ri j=1;i*j<=maxn;++j)
				if(vis[i*j]) ++scnt[i];
		for(ri i=1;i<=maxn;++i)
			if(scnt[i]>=4) ans+=(ll)miu[i]*(ll)(scnt[i]-3)*(ll)(scnt[i]-2)*(ll)(scnt[i]-1)*(ll)scnt[i]/24LL;
		cout<<ans<<'\n';
	}
	return 0;
}