[bzoj2154]Crash的数字表格

[bzoj2154]Crash的数字表格

标签： 莫比乌斯反演 积性函数 线性筛

---

题目链接

题解

我的貌似和网上大多数的题解都不一样......

首先设(f(d)为sumsum ij(gcd(i,j)==d))
那么(F(d)为sumsum ij(d|gcd(i,j)))
由莫比乌斯反演定理可知，(f(i)=sum mu(d/i)F(d))。

[ egin{align} sum_{i=1}^n frac{f(i)}i & =sum_{i=1}^n frac{sum_{i|d}mu(d/i)F(d)}{i} \ & = sum_{i=1}^n F(i)sum_{d|i}frac{mu(i/d)}{d} \ & = sum_{i=1}^n F(i)frac{sum_{d|i} mu(d)d}{i} end{align} ]

这就非常舒服了。
现在问题就是如何快速求出$ sum_{d|i} mu(d)d ( 这是一个积性函数，并且有很好的性质： 对于质数p， ) p|n (时，) f(pn)=f(n) ( 否则，) f(pn)=f(p)f(n) $

Code

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<cmath>
#include<queue>
#include<stack>
#include<set>
#include<map>
using namespace std;
#define ll long long
#define REP(i,a,b) for(register int i=(a),_end_=(b);i<=_end_;i++)
#define DREP(i,a,b) for(int i=(a),_end_=(b);i>=_end_;i--)
#define EREP(i,a) for(int i=start[(a)];i;i=e[i].next)
inline int read()
{
	int sum=0,p=1;char ch=getchar();
	while(!(('0'<=ch && ch<='9') || ch=='-'))ch=getchar();
	if(ch=='-')p=-1,ch=getchar();
	while('0'<=ch && ch<='9')sum=sum*10+ch-48,ch=getchar();
	return sum*p;
}

const int maxn=1e7+20;
const int mod=20101009;
const int inv4=15075757;

int prime[maxn],tot,g[maxn],mark[maxn],mu[maxn];

int n,m;
void prepare()
{
	g[1]=mu[1]=1;
	REP(i,2,n)
	{
		if(!mark[i])
		{
			prime[++tot]=i;
			mu[i]=-1;g[i]=-i+1;
		}
		for(register int j=1;j<=tot && i*prime[j]<=n;j++)
		{
			mark[i*prime[j]]=1;
			if(i%prime[j])mu[i]=-1,g[i*prime[j]]=(ll)g[i]*g[prime[j]]%mod;
			else mu[i]=0,g[i*prime[j]]=g[i];
		}
	}
}


void init()
{
	n=read();m=read();
	if(n>m)swap(n,m);
}

int ans;

void doing()
{
	REP(i,1,n)
	{
		register int x=n/i,y=m/i;
		ans=(ans+((ll)(x)*(x+1)%mod)*((ll)(y)*(y+1)%mod)%mod*((ll)i*inv4%mod*g[i]%mod)%mod)%mod;
	}
	printf("%d
",(ans+mod)%mod);
}

int main()
{
	freopen("number.in","r",stdin);
	freopen("number.out","w",stdout);
	init();
	prepare();
	doing();
	return 0;
}