基准时间限制:5 秒 空间限制:262144 KB 分值: 640 难度:8级算法题
给出一个数N,输出小于等于N的所有数,两两之间的最大公约数之和。
相当于计算这段程序(程序中的gcd(i,j)表示i与j的最大公约数):
由于结果很大,输出Mod 1000000007的结果。
G=0;
for(i=1;i<=N;i++)
for(j=1;j<=N;j++)
{
G = (G + gcd(i,j)) % 1000000007;
}
Input
输入一个数N。(2 <= N <= 10^10)
Output
输出G Mod 1000000007的结果。
Input示例
100
Output示例
31080
这个其实要比最小公倍数的版本要好做,之前我也推过,很容易我们可以得到:
ans=∑d=1n⌊nd⌋2ϕ(d)
a
n
s
=
∑
d
=
1
n
⌊
n
d
⌋
2
ϕ
(
d
)
现在只是n太大了,那么在这里我们仍然可以分块来处理,只是需要处理欧拉函数的前缀和,那这个用杜教筛很容易就解决了。有:
sum(n)=n2+n2−∑i=2nsum(⌊nk⌋)
s
u
m
(
n
)
=
n
2
+
n
2
−
∑
i
=
2
n
s
u
m
(
⌊
n
k
⌋
)
代码:
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<map>
#include<cmath>
#define maxx 5000000
#define mod 1000000007
#define ll long long
using namespace std;
bool isP[maxx];
int prime[400000];
int cnt;
ll phi[maxx];
ll inv2=500000004;
void init()
{
phi[1]=1;
for(int i=2;i<maxx;i++)
{
if(!isP[i])
{
prime[cnt++]=i;
phi[i]=i-1;
}
for(int j=0;j<cnt&&(ll)i*prime[j]<maxx;j++)
{
isP[i*prime[j]]=true;
if(i%prime[j])
phi[i*prime[j]]=phi[i]*(prime[j]-1);
else
{
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
}
}
for(int i=1;i<maxx;i++)
phi[i]=(phi[i]+phi[i-1])%mod;
}
map<ll,ll>M;
ll work(ll x)
{
if(x<maxx)return phi[x];
if(M[x])return M[x];
ll t=x%mod;
ll ans=t*(t+1)%mod*inv2%mod;
for(ll i=2,last;i<=x;i=last+1)
{
last=x/(x/i);
ans=(ans-(last-i+1)%mod*work(x/i)%mod)%mod;
}
if(ans<0)ans=(ans+mod)%mod;
M[x]=ans;
return ans;
}
int main()
{
init();
ll ans=0;
ll n;
cin>>n;
for(ll i=1,last;i<=n;i=last+1)
{
last=n/(n/i);
ll now=(n/i)%mod;
now=now*now%mod;
ll base=(work(last)-work(i-1)+mod)%mod;
ans=(ans+now*base%mod)%mod;
}
cout<<ans<<endl;
return 0;
}