Description
已知 f(x)=∑d|xμ(d)∗d
现在请求出下面式子的值
∑i=1n∑j=1nf(gcd(i,j))∗f(lcm(i,j))
由于值可能过大所以请对 109+7 取模
Solution
首先, f(x) 是一个积性函数,这个由莫比乌斯反演的推理可知,
那么,把原来的式子换一个方式写:
设
i=∏Paii,j=∏Pbii
∑i=1n∑j=1n(∏f(Pmin(ai,bi)i))∗(∏f(Pmax(ai,bi)i))
即:
∑i=1n∑j=1nf(i)∗f(j)
即:
(∑i=1nf(i))2
设: ans=∑ni=1f(i)
继续:
ans=∑i=1n∑d|iμ(d)∗d
ans=∑i=1n⌊ni⌋μ(i)∗i
好了,剩下的就是杜教筛版子了
Code
#include <cstdio>
#include <map>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define Gsum(l,r) ((((LL)(l+r)*(r-l+1))>>1)%mo)
using namespace std;
typedef long long LL;
const int N=200500,mo=1e9+7;
int read(int &n)
{
char ch=' ';int q=0,w=1;
for(;(ch!='-')&&((ch<'0')||(ch>'9'));ch=getchar());
if(ch=='-')w=-1,ch=getchar();
for(;ch>='0' && ch<='9';ch=getchar())q=q*10+ch-48;n=q*w;return n;
}
int n;
LL ans;
int pr[N/2],mu[N+1];
bool prz[N+1];
map<int,int>Jmu;
LL Gmu(int n)
{
if(n<=N)return mu[n];
if(Jmu[n]!=2333&&Jmu[n])return Jmu[n];
LL ans=0;
for(LL nx,i=2;i<=n;i=nx+1)
{
nx=n/(n/i);
ans=(ans+Gmu(n/i)*Gsum(i,nx))%mo;
}
ans=(1-ans+mo);
if(ans>=mo)ans-=mo;
return Jmu[n]=ans;
}
int main()
{
int q,w;
read(n);
mu[1]=1;
Jmu[0]=Jmu[1]=2333;
fo(i,2,N)
{
Jmu[i]=2333;
if(!prz[i])pr[++pr[0]]=i,mu[i]=-1;
fo(j,1,pr[0])
{
int t=i*pr[j];
if(t>N)break;
prz[t]=1;
if(i%pr[j]==0)break;
mu[t]=-mu[i];
}
}
fo(i,2,N)mu[i]=((LL)mu[i]*i+mu[i-1])%mo;
ans=0;
for(int i=1,nx;i<=n;i=nx+1)
{
nx=n/(n/i);
ans=(ans+(LL)(n/i)*(Gmu(nx)-Gmu(i-1)))%mo;
}
ans=ans*ans%mo;
ans=(ans+mo)%mo;
printf("%lld\n",ans);
return 0;
}