51nod1238 最小公倍数之和 V3 莫比乌斯函数 杜教筛

题意：求\(\sum_{i = 1}^{n}\sum_{j = 1}^{n}lcm(i, j)\).
　　题解：因为是用的莫比乌斯函数求的，所以推导比大部分题解多。。。而且我写式子一般都比较详细，所以可能看上去很多式子，实际上是因为每一步都写了，几乎没有跳过的。所以应该都可以看懂的。
　　末尾的\(e\)函数是指的\(e[1] = 1\),\(e[x] = 0(x != 1)\)这样一个函数
　　\[\sum_{i = 1}^{n}\sum_{j = 1}^{n}lcm(i, j)\]
　　\[\sum_{i = 1}^{n} \sum_{i = 1}^{n} \frac{ij}{gcd(i, j)}\]
　　枚举\(gcd\)
　　\[\sum_{d = 1}^{n} \sum_{i = 1}^{\lfloor {\frac{n}{d}} \rfloor} \sum_{j = 1}^{\lfloor {\frac{n}{d}} \rfloor}[gcd(i, j) == d] \frac{ij}{d}\]
　　因为\((\frac{ijd^2}{d} = ijd)\)，所以：
　　\[\sum_{d = 1}^{n} \sum_{i = 1}^{\lfloor {\frac{n}{d}} \rfloor} \sum_{j = 1}^{\lfloor {\frac{n}{d}} \rfloor}[gcd(i, j) == d] ijd\]
　　\[\sum_{d = 1}^{n}d\sum_{i = 1}^{\lfloor {\frac{n}{d}} \rfloor} \sum_{j = 1}^{\lfloor {\frac{n}{d}} \rfloor}ij[gcd(i, j) == 1]\]
　　\[\sum_{d = 1}^{n} d \sum_{i = 1}^{\lfloor {\frac{n}{d}} \rfloor} \sum_{j = 1}^{\lfloor {\frac{n}{d}} \rfloor} ij \sum_{k | gcd(i, j)} \mu(k)\]
　　枚举k,再枚举k的倍数。
　　\[\sum_{d = 1}^{n}d\sum_{k = 1}^{\lfloor {\frac{n}{d}} \rfloor}\mu(k) \sum_{i = 1}^{\lfloor {\frac{n}{dk}} \rfloor}ik \sum_{j = 1}^{\lfloor {\frac{n}{dk}} \rfloor}jk\]
　　设\(S(n) = \sum_{i = 1}^{n}i\)
　　\[\sum_{d = 1}^{n}d \sum_{k = 1}^{\lfloor {\frac{n}{d}} \rfloor} \mu(k) k ^ 2 S(\frac{n}{dk})\]
　　枚举\(T = dk\)
　　\[\sum_{T = 1}^{n} S(\frac{n}{T})^ 2 \sum_{k | T} \mu(k) k ^ 2 \frac{T}{k}\]
　　\[\sum_{T = 1}^{n} S(\frac{n}{T})^ 2 \sum_{k | T} \mu(k) kT\]
　　\[\sum_{T = 1}^{n} S(\frac{n}{k})^ 2 \cdot T \sum_{k | T} \mu(k)k\]
　　设\(f(T) = T\sum_{k | T} \mu(k) k\)，卷上\(id^2\)，因为\(S(\frac{n}{k})\)可以数论分块，所以我们只需要快速求出区间\([l, r]\)内的\(f\)之和即可,显然求出\(f\)的前缀和即可解决问题
　　\[(f * id^2)(n) = \sum_{i |n}f(i) \frac{n^2}{i^2}=\sum_{i | n}i \sum_{k | i} \mu(k) k \frac{n ^ 2}{i ^ 2}\]
　　\[\sum_{i | n}\sum_{k | i} \mu(k)k\frac{n ^ 2}{i} = n \sum_{i | n}\sum_{k | i} \mu(k) k \frac{n}{i}\]
　　设\[h(i) = \sum_{k | i} \mu(k)k\]，则原式：
　　\[n \sum_{i | n} h(i) \frac{n}{i} = n (h * id)(n)\]
　　\[(f * id ^ 2)(n) = n (h * id)(n)\]
　　\[h(n) = \sum_{k | n}\mu(k)k = (\mu \cdot id) * 1\]
　　\[f * id ^ 2 = n [(\mu \cdot id) * 1 * id] = n[(\mu \cdot id) * id * 1]\]
　　其中\[(\mu \cdot id) * id = \sum_{i | n} \mu(i) i \frac{n}{i} = n \sum_{i | n}\mu(i) = e\]
　　所以
　　\[n[(\mu \cdot id) * id * 1] = n[e * 1] = n\]
　　带入杜教筛的式子：
　　\[g(1)S(n) = \sum_{i = 1}^{n} (f * g)(i) - \sum_{i = 2}^{n}g(i)S(\frac{n}{i})\]
　　\[= \sum_{i = 1}^{n}i - \sum_{i = 2}^{n}i ^ 2 S(\frac{n}{i})\]
　　然后直接上杜教就可以了.
　　其实还有一个问题。。。一开始预处理的前缀和怎么求?
　　要知道前缀和，首先要求出\(f\).
　　因为\(f(T) = T\sum_{k | T} \mu(k) k\)，所以如果我们可以快速求出\(\sum_{k | T}\mu(k)k\),然后就只需要再\(O(n)\)的乘上\(T\)就可以了.
　　我们先预处理出\(\mu(k)\)，然后对于每一个\(k\),枚举它的倍数，统计贡献。那么复杂度为 \(\frac{n}{1} + \frac{n}{2} + ... + \frac{n}{n} = nlogn\)(此处的\(n\)为原题面的\(\frac{2}{3}\)次方，即要预处理的\(f\)个数)

#include<bits/stdc++.h>
using namespace std;
#define R register int
#define LL long long
#define RL register LL
#define AC 3000
#define ac 5000000
#define p 1000000007LL
//#define h(x) ((x <= block) ? sum[x] : S[n / x])

LL n, ans, block;
LL mu[ac], S[AC], sum[ac], inv[AC];
int pri[ac], tot;
bool z[ac], vis[AC];

inline LL h(LL x)
{
    return ((x <= block) ? sum[x] : S[n / x]);
}

inline void up(LL & a, LL b)
{
    a += b;
    if(a >= p) a -= p;
    if(a <= -p) a += p;
}

LL count(LL l, LL r){
    return (r - l + 1) % p * ((r + l) % p) % p * inv[2] % p;
}

void pre()
{
    scanf("%lld", &n), block = pow(n, 0.66666);
    mu[1] = 1;
    for(R i = 2; i <= block; i ++)
    {
        if(!z[i]) pri[++ tot] = i, mu[i] = -1;
        for(R j = 1; j <= tot; j ++)
        {
            int now = pri[j];
            if(i * now > block) break;
            z[i * now] = true;
            if(!(i % now)) break;
            mu[i * now] = - mu[i];
        }
    }
    inv[1] = 1;
    for(R i = 2; i <= 10; i ++) inv[i] = (p - p / i) * inv[p % i] % p;
    for(R i = 1; i <= block; i ++)//枚举mu(i)
        for(R j = 1; j; j ++)//枚举i的倍数
        {
            if(j * i > block) break;
            up(sum[i * j], mu[i] * i % p);
        }
    for(R i = 1; i <= block; i ++) sum[i] = sum[i] * i % p;
    for(R i = 1; i <= block; i ++) up(sum[i], sum[i - 1]);//算出f数组后还要统计前缀和
}

LL get(LL x)
{
    x %= p;
    return x * (x + 1) % p * (2 * x + 1) % p * inv[6] % p;
}

void cal(LL x)
{
    if(x <= block || vis[n / x]) return ;
    LL rnt = count(1, x);
    for(RL i = 2, lim, now; i <= x; i = lim + 1)
    {
        lim = x / (x / i), now = x / i, cal(now);
        up(rnt, - ((get(lim) - get(i - 1)) % p * h(now) % p));
    }
    S[n / x] = rnt, vis[n / x] = true;
}

void work()
{
    for(RL i = 1, lim, now, x; i <= n; i = lim + 1)
    {
        lim = n / (n / i), now = n / i, x = count(1, now);
        up(ans, (x % p * x % p * ((h(lim) - h(i - 1)) % p) % p));
    }
    printf("%lld\n", (ans + p) % p);
}

int main()
{
    //freopen("in.in", "r", stdin);
    pre();
    cal(n);
    work();
//  fclose(stdin);
    return 0;
}

原文地址：https://www.cnblogs.com/ww3113306/p/10166573.html