欧拉函数的应用：快速求解 1~n 中两两数字的最小公倍数

今天遇到一个十分 Dark 的题目，让你求：

$$\sum_{i=1}^{n-1} \sum_{j=i+1}^{n} lcm(i,j)$$

一共 $T$ 组数据，数据范围是：$T \leqslant 2 \ast 10^5, n \leqslant 3\ast 10^6$……

题目链接：LightOJ 1375 LCM Extreme

这题据说可以用莫比乌斯反演做，但是看了一大堆似乎还是没怎么懂……但是翻大佬的博客找到一个精妙的做法～

首先把上面那个公式分开，先求：

$$\sum_{i=1}^{n} lcm(n,i)$$

显然有：

$$lcm(n,i)=\frac {ni} {gcd(n,i)}$$

设 $g=gcd(n,i)$，则 $\displaystyle \frac n g$ 与 $\displaystyle \frac i g$ 互质。
原式 = $\displaystyle g\ast \frac {n} {g} \ast \frac {i} {g}$ ，并且显然 $\displaystyle \frac i g \leqslant \frac n g$
接下来要用到一个公式：小于 $x$ 的与 $x$ 互质的数之和为 $\displaystyle \frac {x \ast \varphi(x)} {2}$ ；
所以小于 $\displaystyle \frac n g$ 的且与其互质之数之和：

$$\frac {(\frac n g)\ast \varphi(\frac n g)} {2}$$

令 $\displaystyle d=\frac n g$，则 $d$ 是 $n$ 的因子。方便后面枚举。得到：

$$\sum_{i=1}^{n-1} \frac i g=\sum_{d|n} \frac {d \ast \varphi (d)} {2}$$

因为我们要让左边得到 $lcm(n,i)$，根据 $\displaystyle lcm(n,i)= \displaystyle n\ast \frac {n} {g} \ast \frac {i} {g}$，我们可以将左边的 $\displaystyle \frac n g \ast g$ 提取出来得到：

$$\sum_{i=1}^{n-1} lcm(n,i)=\sum_{d|n} n \ast \frac {d \ast \varphi (d)} {2}$$

哇！！！就要成功了！接下来考虑 $i=n$ 的情况，因为 $i=n$ 的时候 $\displaystyle \frac n g=1$，需要特判：

$$\sum_{i=1}^{n} lcm(n,i)=n+ \sum_{d|n} n \ast \frac {d \ast \varphi(d)} 2$$

最后我们只需要用欧拉筛进行一次 O(N) 求欧拉函数就可以方便地求出 $\displaystyle \sum_{i=1}^{n} lcm(n,i)$ 了～枚举因子的复杂度是调和级数，那么这个程序复杂度由 $\Theta (n)$ 降低到 $\log(n)$。

现在回到原问题，求出 $\displaystyle \sum_{i=1}^{n-1} \sum_{j=i+1}^{n} lcm(i,j)$，只要枚举 N 求解，前缀和累加一下就好了！！

还有说下这题题目让我们对 $2^{32}$ 取模，如果直接开 long long 再模可能会爆内存。处理技巧是直接开 unsigned int，然后自然溢出（也就是不用模的）。最后输出不能用 %d，要用 %u。

复杂度 $\Theta(n \log n)$。

贴上代码（要注意一点细节，不然还是可能爆内存……）：

#include<cstdio>
#include<cstring>
#include<iostream>
using namespace std;
const int maxn=3000005,N=3000000;
int T,n,prime[maxn];
bool vis[maxn];
unsigned long long f[maxn];
unsigned int phi[maxn];
inline int read(){
    int ret=0,f=1;char ch=getchar();
    while (ch<'0'||ch>'9') {if (ch=='-') f=-1;ch=getchar();}
    while (ch>='0'&&ch<='9') ret=ret*10+ch-'0',ch=getchar();
    return  ret*f;
}
inline void BuildPhi(){ // 线性筛法求欧拉函数，黑科技！
    phi[1]=1;
    memset(vis,1,sizeof(vis));
    vis[1]=false;
    for (int i=2;i<=N;i++){
        if (vis[i]){
            phi[i]=i-1;
            prime[++prime[0]]=i;
        }
        for (int j=1;j<=prime[0];j++){
            if (i*prime[j]>N) break;
            vis[i*prime[j]]=false;
            if (i%prime[j]) phi[i*prime[j]]=(prime[j]-1)*phi[i];
            else {phi[i*prime[j]]=prime[j]*phi[i];break;}
        }
    }
}
inline void BuildSum(){
    for (int i=1;i<=N;i++){
        for (int j=i;j<=N;j+=i) f[j]+=(unsigned long long)phi[i]*i/2*j; // i 作为因子
        f[i]=f[i-1]+f[i]; // i 之前的一定都求好了，干脆累加起来
    }
}
int main(){
    T=read();
    BuildPhi();
    BuildSum();
    for (int t=1;t<=T;t++){
        n=read();
        printf("Case %d: %llu\n",t,f[n]);
    }
    return 0;
}