[bzoj4176]Lucas的数论

Orz w_yqts
Orz xudyh
Orz popoqqq
求∑nx=1∑ny=1d(ij)
n<=109

首先证明d(nm)=∑i|n∑j|m[gcd(i,j)=1]
若有n=n′∗pk1,m=m′∗pk2
则素数p对左式d(nm)的贡献为k1+k2+1
对于右式则有(pk1,1),(pk1−1,1),...,(p,1),(1,1),(1,p),...,(1,pk2−1),(1,pk2)
贡献也是k1+k2+1
所以等号成立
直接把这个式子代进去反演一下就好了……?

μ的部分用杜教筛做…

#include <bits/stdc++.h>using namespace std;#define p 1000000007#define ll long long#define N 2500000map <int,int> f;int pn,pr[N],flag[N],miu[N];int n;inline int sqr(int x){    return (ll)(x)*x%p;}inline int calc(int n){    int res=0;    for (int i=1,pos;i<=n;i=pos+1)    {        pos=n/(n/i);        res+=(ll)(pos-i+1)*(n/i)%p;        res%=p;    }    return res;}inline int calcmiu(int n){    if (n<N) return miu[n];    if (f[n]) return f[n];    int res=1;    for (int i=2,pos;i<=n;i=pos+1)    {        pos=n/(n/i);        res-=calcmiu(n/i)*(pos-i+1);    }    return f[n]=res;}void init(){    miu[1]=flag[1]=1;    for (int i=2;i<N;++i)    {        if (!flag[i]) pr[++pn]=i,miu[i]=-1;        for (int j=1;j<=pn && i*pr[j]<N;++j)        {            flag[i*pr[j]]=1;            if (i%pr[j]==0) {miu[i*pr[j]]=0;break;}            miu[i*pr[j]]=-miu[i];        }    }    for (int i=1;i<N;++i) miu[i]+=miu[i-1];}inline int mobius(int n){    int res=0;    for (int i=1,pos;i<=n;i=pos+1)    {        pos=n/(n/i);        res+=(ll)(calcmiu(pos)-calcmiu(i-1))*sqr(calc(n/i))%p;        res%=p;    }    return (res+p)%p;}main(){    init();    cin>>n;    //for (int i=1;i<=n;++i) cout<<miu[i]<<' ';cout<<endl;    cout<<mobius(n)<<endl;}