洛谷-3911 最小公倍数之和

题目描述
对于$A_1,A_2,\cdots ,A_N$
${\sum }_{i=1}^N{\sum }_{j=1}^{N}lcm(A_i,A_j)$
的值。
lcm(a,b) 表示a 和b 的最小公倍数
输入格式
第1 行，1 个整数N。
第2 行，N 个整数$A_1,A_2,\cdots ,A_N$
输出格式
1 个整数，表示所求的值。

输入输出样例
输入 #1
2
2 3

输出 #1 复制
17

说明/提示
• 对于100% 的数据，$1\le N\le 50000;1\le A_i\le 50000$

解释：$F_x={\sum }_{i=1}^n[A_i=x]$
即原始化为
$I={\sum }_{x=1}^p{\sum }_{y=1}^p[x,y]F_x{F}_y$
$={\sum }_{x=1}^p{\sum }_{y=1}^p\frac{xy{F}_x{F}_y}{(x,y)}$
根据经验，考虑枚举 $g=(x,y)$。显然 g 的取值范围为 $[1,p]\cap \mathbb{N}[1,p]\cap N$。
$I={\sum }_{g=1}^p{\sum }_{x=1}^p{\sum }_{y=1}^p\frac{xy{F}_x{F}_y}{g}[(x,y)=g]$
转而枚举 x,y 的倍数。
$I={\sum }_{g=1}^p{\sum }_{x=1}^{\lfloor \frac{p}{g}\rfloor }{\sum }_{y=1}^{\lfloor \frac{p}{g}\rfloor }\frac{gxgy{F}_{gx}{F}_{gy}}{g}[(gx,gy)=g]$
$I={\sum }_{g=1}^p{\sum }_{x=1}^{\lfloor \frac{p}{g}\rfloor }{\sum }_{y=1}^{\lfloor \frac{p}{g}\rfloor }gxy{F}_{gx}{F}_{gy}{\sum }_{d\mid (x,y)}\mu (d)$
$I={\sum }_{g=1}^p{\sum }_{d=1}^{\lfloor \frac{p}{g}\rfloor }{\sum }_{x=1}^{\lfloor \frac{p}{dg}\rfloor }{\sum }_{y=1}^{\lfloor \frac{p}{dg}\rfloor }gdxdy{F}_{dgx}{F}_{dgy}\mu (d)$
$={\sum }_{g=1}^{p}g{\sum }_{d=1}^{\lfloor \frac{p}{g}\rfloor }{d}^2\mu (d){\sum }_{x=1}^{\lfloor \frac{p}{dg}\rfloor }{\sum }_{y=1}^{\lfloor \frac{p}{dg}\rfloor }xy{F}_{dgx}{F}_{dgy}$
$={\sum }_{g=1}^{p}g{\sum }_{d=1}^{\lfloor \frac{p}{g}\rfloor }{d}^2\mu (d)\left({\sum }_{x=1}^{\lfloor \frac{p}{dg}\rfloor }x{F}_{dgx}\right)\left({\sum }_{y=1}^{\lfloor \frac{p}{dg}\rfloor }y{F}_{dgy}\right)$
记 $G_t={\sum }_{x=1}^{\lfloor \frac{p}{t}\rfloor }x{F}_{tx}$
​则原式化为
$I={\sum }_{g=1}^p{\sum }_{d=1}^{\lfloor \frac{p}{g}\rfloor }g{d}^2\mu (d)G_{dg}^2$
$I={\sum }_{T=1}^p{\sum }_{d\mid T}\frac{T}{d}{d}^2\mu (d)G_T^2$
$={\sum }_{T=1}^{p}T{G}_T^2\left({\sum }_{d\mid T}d\mu (d)\right)$
不妨设 $H(x)={\sum }_{d\mid x}d\mu (d)$
${\sum }_{T=1}^{p}T{G}_T^{2}H(T)$
最后对H筛一下就好了

#include<iostream>
#include<cstdio>
typedef long long ll;
const int N = 1e5 + 50;
long long pri[N],mu[N],phi[N],A[N],H[N],T[N],P[N],cnt;
bool isnt_p[N];
long long n;
inline void get_p(){
    isnt_p[1] = true; phi[1] = mu[1] = P[1] = 1;
    bool flag = true;
    for (int i = 1; i < N; ++i){
        if (isnt_p[i] == false){
            pri[++cnt] = i;
            mu[i] = -1;
            phi[i] = i - 1;
            P[i] = 1 - i;
        }
        for (int j = 1; j <= cnt && i * pri[j] < N; ++j){
            isnt_p[i * pri[j]] = true;
            if (i % pri[j]){
                mu[i * pri[j]] = -mu[i];
                P[i * pri[j]] = P[i] * P[pri[j]];
                phi[i * pri[j]] = phi[i] * phi[pri[j]];
            }
            else{
                mu[i * pri[j]] = 0;
                phi[i * pri[j]] = phi[i] * pri[j];
                P[i * pri[j]] = P[i]; 
                break;
            }
        }
    }
}
inline ll read() {ll x; scanf("%lld", &x); return x;}
int main(){
    get_p();
    n = read();
    const int N = 1e5;
    for (int i = 1; i <= n; ++i) A[i] = read(), ++T[A[i]];
    for (int x = 1; x <= N; ++x)
        for (int i = 1; i <= (N / x); ++i)
            H[x] += i * T[x * i];
    ll ans = 0;
    for (int i = 1; i <= N; ++i)
        ans += 1ll * i * H[i] * H[i] * P[i];
    printf("%lld", ans);
}