并不对劲的bzoj4804:欧拉心算

题目大意

\(t\)(\(t\leq5000\))组询问，每次询问给出\(n\)(\(n\leq10^7\))，求：

\[\sum_{i=1}^{n}\sum_{j=1}^{n}\phi(gcd(i,j)) \]

题解

枚举gcd，原式变为：

\[\sum_{k=1}^{n}\phi(k)\sum_{i=1}^{n}\sum_{j=1}^{n}[gcd(i,j)=k] \]

\[\sum_{k=1}^{n}\phi(k)\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}[gcd(i,j)=1] \]

发现\(\sum_{j=1}^{i}[gcd(i,j)=1] = \phi(i)\)(1)
那么将\(\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}[gcd(i,j)=1]\)中\(i>j\)和\(i<j\)分开考虑，相当于是把(1)式算了两遍
但是\(i=j=1\)算重（chong二声）了，所以是两个(1)式-1
即\(\sum_{i=1}^{\lfloor\frac{n}{k}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{k}\rfloor}[gcd(i,j)=1] = (\sum_{i=1}{\lfloor\frac{n}{k}\rfloor}2*\phi(i))-1\)
那么原式=\(\sum_{k=1}^{n}\phi(k)( (\sum_{i=1}{\lfloor\frac{n}{k}\rfloor}2*\phi(i))-1)\)
直接整除分块就行了

代码

#include<algorithm>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<ctime>
#include<iomanip>
#include<iostream>
#include<map>
#include<queue>
#include<set>
#include<stack>
#include<vector>
#define rep(i,x,y) for(register int i=(x);i<=(y);++i)
#define dwn(i,x,y) for(register int i=(x);i>=(y);--i)
#define maxn 10000001
#define LL long long
#define lim (maxn-1)
using namespace std;
int read()
{
	int x=0,f=1;char ch=getchar();
	while(!isdigit(ch)&&ch!='-')ch=getchar();
	if(ch=='-')f=-1,ch=getchar();
	while(isdigit(ch))x=(x<<1)+(x<<3)+ch-'0',ch=getchar();
	return x*f;
}
void write(LL x)
{
	if(x==0){putchar('0'),putchar('\n');return;}
	int f=0;char ch[20];
	if(x<0)putchar('-'),x=-x;
	while(x)ch[++f]=x%10+'0',x/=10;
	while(f)putchar(ch[f--]);
	putchar('\n');
	return;
}
int no[maxn],p[maxn],cnt,t,n;
LL phi[maxn],f[maxn];
int main()
{
	no[1]=phi[1]=1;
	rep(i,2,lim)
	{
		if(!no[i])phi[i]=i-1,p[++cnt]=i;
		for(int j=1;j<=cnt&&i*p[j]<=lim;j++)
		{
			no[i*p[j]]=1;
			if(i%p[j]==0){phi[i*p[j]]=phi[i]*p[j];break;}
			phi[i*p[j]]=phi[i]*phi[p[j]];
		}
	}
	rep(i,1,lim)phi[i]+=phi[i-1];
	rep(i,1,lim)f[i]=phi[i]*2ll-1ll;
	t=read();
	while(t--)
	{
		n=read();LL ans=0;
		for(int l=1,r=0;l<=n;l=r+1)
		{
			r=n/(n/l);
			ans+=(phi[r]-phi[l-1])*f[n/l];
		}
		write(ans);
	}
	return 0;
}