P1829 [国家集训队]Crash的数字表格 / JZPTAB

DescriptionDescriptionDescription


i=1nj=1mlcm(i,j)\huge \sum_{i=1}^n\sum_{j=1}^m lcm(i,j)i=1∑n​j=1∑m​lcm(i,j)

数据范围:n,m107n,m\leq 10^7n,m≤107

有模数


SolutionSolutionSolution

其实可以用这道题的优化版本预处理f(x)=j=1xlcm(j,x)f(x)=\sum_{j=1}^x lcm(j,x)f(x)=∑j=1x​lcm(j,x)做了

但其实有更好打的方法。。。
莫比乌斯反演。。。

整除分块套一起即可,时间复杂度:O(n+m)O(n+m)O(n+m)


CodeCodeCode

#include<cstdio>
#include<cctype>
#include<algorithm>
#define N 10000000
#define LL long long
#define mod 20101009
using namespace std;
LL mu[N+10],f[N+10],sum[N+10],n,m,k;int prime[N>>3],num,t,x,y;
bool vis[N+10];
inline LL read()
{
	LL f=0,d=1;char c; 
	while(c=getchar(),!isdigit(c)) if(c=='-') d=-1;f=(f<<3)+(f<<1)+c-48;
	while(c=getchar(),isdigit(c)) f=(f<<3)+(f<<1)+c-48;
	return d*f;
}
inline LL solve(int a,int b)
{
	LL res=0;
	for(register int l=1,r=0;l<=a;l=r+1)
	{
		r=min(a/(a/l),b/(b/l));
		res+=(LL)(mu[r]-mu[l-1])*(LL)f[a/l]*f[b/l];
	}
	return res;
}
inline LL Sum(int x,int y){return (1ll*x*(x+1)/2%mod)*(1ll*y*(y+1)/2%mod)%mod;}
inline LL fk1(int x,int y)
{
	int res=0;
	for(register int l=1,r=0;l<=min(x,y);l=r+1) 
	{
		r=min(x/(x/l),y/(y/l));
		res=(res+1ll*(sum[r]-sum[l-1]+mod)*Sum(x/l,y/l)%mod)%mod;
	}
	return res;
}
inline LL fk2(int x,int y)
{
	int res=0;
	for(register int l=1,r=0;l<=min(x,y);l=r+1) 
	{
		r=min(x/(x/l),y/(y/l));
		res=(res+1ll*(r-l+1)*(l+r)/2%mod*fk1(x/l,y/l)%mod)%mod;
	}
	return res;
}
signed main()
{
	n=read();m=read();
	mu[1]=1;k=min(n,m);
	for(register int i=2;i<=k;i++)
	{
		if(vis[i]==0) {prime[++num]=i;mu[i]=-1;}
		for(register int j=1;j<=num&&prime[j]*i<=k;j++)
		{
			vis[prime[j]*i]=true;
			if(i%prime[j]==0) break;
			mu[prime[j]*i]=-mu[i];
		}
	}
	for(register int i=1;i<=k;i++) sum[i]=(sum[i-1]+1ll*i*i%mod*(mu[i]+mod))%mod;
	printf("%lld",fk2(n,m));
}
上一篇:UVa11889 - Benefit


下一篇:CodeForces 55D