一道水题---lucas的数论

题目很简洁:

求

$$\sum_{i=1}^n \sum_{j=1}^nf(i\times j)$$
其中 $f(x)$表示 $x$ 的约数个数.答案很大,mod 1000000007后输出.
数据范围：$1\leq n\leq 10^{9}$

解法：

$d|(i\times j)$ 等价于$\frac {d}{gcd(i,d)}|j$
设$e=gcd(i,d),\ i=k_1\cdot e,\ d=k_2\cdot e$
则只需且必须：$k_1|i且k_2|j，gcd(k_1,k_2)=1$
所以

$$ans=\sum_{i=1}^n \sum_{j=1}^nf(i\times j)$$
$$=\sum_{i=1}^n \sum_{j=1}^n\sum_{k_1|i}^n\sum_{k_2|j}^n[gcd(k_1,k_2)=1]$$
$$=\sum_{k_1=1}^n \sum_{k_2=1}^n[gcd(k_1,k_2)=1]\sum_{k_1|i}^n\sum_{k_2|j}^n 1$$
$$=\sum_{k_1=1}^n \sum_{k_2=1}^n[gcd(k_1,k_2)=1]\cdot[\frac{n}{k_1}]\cdot[\frac{n}{k_2}]$$
到了这一步，就用本人比较喜欢的方式进行莫比乌斯反演.
令 $$f(x)=\sum_{k_1=1}^n \sum_{k_2=1}^n[gcd(k_1,k_2)=x]\cdot[\frac{n}{k_1}]\cdot[\frac{n}{k_2}]$$
$$F(x)=\sum_{k_1=1}^n \sum_{k_2=1}^n[x|gcd(k_1,k_2)]\cdot[\frac{n}{k_1}]\cdot[\frac{n}{k_2}]$$
则由f,F的含义可得 $$F[i]=\sum_{i|d}^nf(d)$$
反演得
$$f(i)=\sum_{i|d}^n\mu(\frac{d}{i})\cdot F(d)$$ 特别地 $$f(1)=\sum_{d=1}^n\mu(d)\cdot F(d)$$
又 $$F(d)=\sum_{d|k_1}^n \sum_{d|k_2}^n\cdot[\frac{n}{k_1}]\cdot[\frac{n}{k_2}]$$ $$=\sum_{k_1=1}^{\frac{n}{d}}\sum_{k_2=1}^{\frac{n}{d}}[\frac{n}{k_1d}]\cdot[\frac{n}{k_2d}]=(\sum_{k=1}^{\frac{n}{d}}[\frac{n}{kd}])^{2}$$
所以 $$ans=f(1)=\sum_{d=1}^n\mu(d)\cdot(\sum_{k=1}^{\frac{n}{d}}[\frac{n}{kd}])^{2}$$
考虑到 $\frac{n}{d}$的不同取值只有不到 $2\sqrt{n}$种，故外层可以在O( $\sqrt{n}$）内枚举,内层同样的道理,也可以快速处理,复杂度不好分析,大约O( $n^{\frac{3}{4}}$).
至于 $\mu(d)$怎么办,用杜教筛求快速前缀和即可,杜教筛复杂度O( $n^{\frac{2}{3}}$).

参考代码

#include<cstdio>
#include<map>
using namespace std; 
#define MAXN 3000010
#define mod 1000000007
#define LL long long
int lim=1000000;
int n,pri[MAXN],tot,mu[MAXN*10];
bool vis[MAXN*10];
void init()
{
    mu[1]=1;
    for(int i=2;i<=lim;i++)
    {
        if(!vis[i])
        {
            pri[++tot]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=tot&&i*pri[j]<=lim;j++)
        {
            vis[i*pri[j]]=1;
            if(i%pri[j]==0)
            {
                mu[i*pri[j]]=0;
                break;
            }
            else
            {
                mu[i*pri[j]]=-mu[i];
            }
        }
    }
    for(int i=1;i<=lim;i++)
    {
        mu[i]+=mu[i-1];
    }
    return ;
}
map<LL,LL>MU;
LL f(LL x)
{
    if(x<=lim)
    {
        return mu[x];
    }
    if(MU.count(x))
    {
        return MU[x];
    }
    LL ans=0,e,ee;
    for(e=2;e<=x;e=ee+1)
    {
        ee=(x/(x/e));
        ans+=f(x/e)*(ee-e+1);
    }
    return MU[x]=1-ans;
}
int ans=0;
inline int cal(int x)
{
    int ff,fff,an=0;
    for(ff=1;ff<=x;ff=fff+1)
    {
        fff=x/(x/ff);
        an+=1ll*(fff-ff+1)*(x/ff)%mod;
        an-=an>=mod?mod:0;
    }
    return an;
}
int main()
{
    scanf("%d",&n);
    init();
    int e,ee,q;
    for(e=1;e<=n;e=ee+1)
    {
        ee=n/(n/e);
        q=cal(n/e);
        ans+=1ll*q*q%mod*(f(ee)-f(e-1))%mod;
        ans+=ans<0?mod:0;
        ans-=ans>=mod?mod:0;
    }
    printf("%d\n",ans);  
    return 0;
}