标签：约数 lfloor GMOJ4496 frac gcd prm sum rfloor 互补

题目

题目链接：https://gmoj.net/senior/#main/show/4496
求

\[\sum^{n}_{i=1}\sum_{j|i}\gcd(j,\frac{i}{j}) \]

\(n\leq 10^{11}\)。

思路

\[\sum^{n}_{i=1}\sum_{j|i}\gcd(j,\frac{i}{j}) \]

\[=\sum^{n}_{i=1}\sum^{\lfloor\frac{n}{i}\rfloor}_{j=1}\gcd(i,j) \]

枚举 \(\gcd=d\)

\[\sum^{n}_{d=1}d\sum^{\lfloor\frac{n}{d^2}\rfloor}_{i=1}\sum^{\lfloor\frac{n}{d^2i}\rfloor}_{j=1}[\gcd(i,j)=1] \]

令 \(g(n)=\sum^{n}_{i=1}\sum^{\lfloor\frac{n}{i}\rfloor}_{j=1}[\gcd(i,j)=1]\)，则原式 \(=\sum^{n}_{d=1}g(\lfloor\frac{n}{d^2}\rfloor)\)。

\[g(n)=\sum^{\sqrt{n}}_{d=1}\mu(d)\sum^{\lfloor\frac{n}{d^2}\rfloor}_{i=1}\lfloor\frac{n}{d^2i}\rfloor \]

整除分块即可。记忆化一下可以卡进 \(1s\)。

代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

const int N=350010,M=10000000;
int m,prm[N],mu[N];
bool v[N];
ll n,ans,f[M];

void findprm(int n)
{
	mu[1]=1;
	for (int i=2;i<=n;i++)
	{
		if (!v[i]) prm[++m]=i,mu[i]=-1;
		for (int j=1;j<=m;j++)
		{
			if (i>n/prm[j]) break;
			v[i*prm[j]]=1; mu[i*prm[j]]=-mu[i];
			if (!(i%prm[j])) { mu[i*prm[j]]=0; break; }
		}
	}
}

ll calc2(ll n)
{
	if (n<M && f[n]) return f[n];
	ll res=0;
	for (ll l=1,r;l<=n;l=r+1)
	{
		r=n/(n/l);
		res+=(r-l+1)*(n/l);
	}
	if (n<M) f[n]=res;
	return res;
}

ll calc1(ll n)
{
	ll res=0;
	for (ll i=1;i*i<=n;i++)
		res+=mu[i]*calc2(n/i/i);
	return res;
}

int main()
{
	freopen("gcd.in","r",stdin);
	freopen("gcd.out","w",stdout);
	findprm(N-1);
	scanf("%lld",&n);
	for (ll i=1;i*i<=n;i++)
		ans+=i*calc1(n/i/i);
	cout<<ans;
	return 0;
}

标签：约数,lfloor,GMOJ4496,frac,gcd,prm,sum,rfloor,互补
来源： https://www.cnblogs.com/stoorz/p/15006209.html

本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享；
2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关；
3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关；
4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除；
5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。