参考:http://blog.csdn.net/wzf_2000/article/details/54630931
有这样一个显然的结论:当( |mu(n)|==1 )时,( phi(nk)=phi(k)sum_{d|gcd(n,k)}phi(frac{n}{d}) )然后看n的范围比较友好就先不去管它,先看后面的:
[if |mu(i)|==1
]
[sum_{k=1}^{i}sum_{d|i,d|k}phi(frac{n}{d})phi(k)
]
[=sum_{d|i}phi(frac{n}{d})sum_{k=1}^{left lfloor frac{m}{d}
ight
floor}phi(dk)
]
发现这形成了一个子问题的形式,于是可以用杜教筛。
对于其他的部分,k是i的因数中最大的且( |mu(k)|==1 )的数:
[if |mu(i)|==0
]
[sum_{j=1}^{m}phi(ij)=frac{isum_{j=1}^{m}phi(kj)}{k}
]
时间复杂度不会算
#include<iostream>
#include<cstdio>
#include<map>
#include<cstring>
using namespace std;
const int N=1000005,m=1000000,mod=1e9+7;
int phi[N],mi[N],q[N],tot,n,k,s[N],ans[N];
bool v[N];
map<long long,int>mp;
int S(int n,int l)
{
if(l<=1)
return phi[n*l];
if(n==1)
{
if(l<=m)
return s[l];
if(ans[k/l]!=-1)
return ans[k/l];
long long re=(long long)l*(l+1)/2%mod;
for(int i=2,la;i<=l;i=la+1)
{
la=l/(l/i);
if(l/i<=m)
re=(re-(long long)s[l/i]*(la-i+1)%mod)%mod;
else
re=(re-(long long)S(1,l/i)*(la-i+1)%mod)%mod;
}
return ans[k/l]=(re%mod+mod)%mod;
}
if(mp[(long long)n*mod+l])
return mp[(long long)n*mod+l];
long long re=0ll;
for(int i=1;i*i<=n;i++)
if(n%i==0)
{
re=(re+(long long)phi[n/i]*S(i,l/i)%mod)%mod;
if(i*i!=n)
re=re+(long long)phi[i]*S(n/i,l/(n/i))%mod;
}
return mp[(long long)n*mod+l]=(re%mod+mod)%mod;
}
int main()
{
memset(ans,-1,sizeof(ans));
mi[1]=phi[1]=1;
for(int i=2;i<=m;i++)
{
if(!v[i])
{
q[++tot]=i;
phi[i]=i-1;
mi[i]=i;
}
for(int j=1;j<=tot&&i*q[j]<=m;j++)
{
int k=i*q[j];
v[k]=1;
if(i%q[j]==0)
{
phi[k]=phi[i]*q[j];
mi[k]=mi[i];
break;
}
phi[k]=phi[i]*(q[j]-1);
mi[k]=mi[i]*q[j];
}
}
for(int i=1;i<=m;i++)
s[i]=(s[i-1]+phi[i])%mod;
scanf("%lld%lld",&n,&k);
if(n>k)
swap(n,k);
long long ans=0ll;
for(int i=1;i<=n;i++)
ans=(ans+((long long)i/mi[i]*S(mi[i],k)%mod))%mod;
printf("%lld
",(ans%mod+mod)%mod);
return 0;
}