Hellc

[BZOJ 2226] LCMSum - 数论函数

【题目描述】

给定 $n$，求： $$ \sum_{i = 1}^n \text{lcm}(i, n) $$ $1 \leq n \leq 10^6$

$T$ 组询问，$1 \leq T \leq 3 \times 10^5$

【题目链接】

【解题思路】

$$ \begin{align*} \text{ans}(n) & =\sum_{i = 1}^n \text{lcm}(i, n) \\ &= \sum_{i = 1}^n\frac {i \cdot n} {\gcd(i, n)} \\ &= n\sum_{i = 1}^n\frac {i} {\gcd(i, n)} \\ &= n\sum_{d | n}\sum_{i = 1}^n \frac {i} {d} [\gcd(i, n) = d] \\ &= n\sum_{d | n}\sum_{i = 1}^{\frac n d} i [\gcd(i, \frac n d) = 1] \\ \end{align*} $$

设 $$ f(x) = \sum_{i = 1}^x i [\gcd(i, x) = 1] $$ 其含义是，$1 \sim x$ 中与 $x$ 互质的数的和，则 $$ \begin{align*} \text{ans}(n) &= n\sum_{d | n} f(\frac n d) \\ &= n\sum_{d | n} f(d) \end{align*} $$ 如何快速计算 $f(x)$ 呢？有结论： $$ f(x) = \begin{cases} 1 & x = 1 \\ \frac {x \cdot \varphi(x)} {2} & x \geq 2, \end{cases} $$ 证明这个结论比较简单。

$x = 1, 2$ 时易得其成立。

对于 $x > 2$ ，与 $x$ 互质的数总是成对出现，就是说如果 $i$ 与 $x$ 互质，那么 $x - i$ 也与 $x$ 互质（因为 $\gcd(i, x) = \gcd(x - i, x)$ ），而它们的和是 $x$ ，这样的数对有 $\frac {\varphi(x)} 2$ 个，总和也就计算出来了。

式子推到这儿，我们已经可以暴力了，线性筛预处理出 $\varphi$ ，然后枚举约数就可以 $O(\sqrt n)$ 回答每次询问了。

总复杂度 $O(T\sqrt n)$，总计算量大概在 $3 \times 10^8$左右，在 SPOJ上虚极了呢，然而在BZOJ这样是可以过的。

我们还是说说复杂度更正确的做法吧， $\sum_{d \mid n}f(d)$ 是积性的，所以我们可以先用线性筛 $O(n)$ 将其预处理出来，然后就可以$O(1)$ 查询了。

为方便，先不管 $\frac 1 2$，以及$x = 1$时的特殊情况，设 $H(x) = x \varphi(x)$，我们要筛的是 $F(n) = \sum_{d | n} H(d)$，即 $F = H \times u$。

我们考虑 $F(p^k)$： $$ \begin{align*} F(p^k) &= \sum_{i = 0}^k p^i \varphi(p^i) \\ &= 1 + \sum_{i = 1}^k p^i \varphi(p^i) \\ &= 1 + \sum_{i = 1}^k p^i (p - 1)p^{i - 1} \\ &= 1 + (p - 1)\sum_{i = 1}^k p^{2i - 1} \\ &= 1 + (p - 1)p\frac {1 - p^{2k}} {1 - p^2} \\ &= \frac {p^{2k + 1} + 1} {p + 1} \\ \end{align*} $$ 其中用到了一次「等比数列求和公式」。

这样我们就可以线性筛出 $F$ 了，然后考虑如何计算答案。 $$ \begin{align*} \text{ans}(n) & = n\sum_{d \mid n} f(d) \\ & = \frac 1 2 n \sum_{d \mid n} 2f(d) \\ & = \frac 1 2 n (F(n) + 1) \\ \end{align*} $$ 好啦～

总结：

在推式子时要善于观察式子的意义，以便利用结论，不要像我一开始那样虽然知道这个结论，但没有观察到，还强行把最后的 $\gcd = 1$ 用 $\mu$ 展开了...

线性筛预可以快速预处理积性函数，关键在于求 $f(p^k)$ 。

【AC代码】

#include <cstdio>

typedef long long int64;

const int64 MAX_N = 1000000;

inline int64 calc(int64 p, int64 pow){
    return (pow * pow * p + 1) / (p + 1);
}

bool isNotPrime[MAX_N + 1];
int primes[MAX_N + 1], m;
int min[MAX_N + 1];
int64 f[MAX_N + 1];

inline void seive(int64 n = MAX_N){
    isNotPrime[1] = true, min[1] = 0, f[1] = 1;
    for(int i = 2; i <= n; i++){
        if(!isNotPrime[i]){
            primes[m++] = i;
            min[i] = i;
            f[i] = (int64)i * i - i + 1;
        }

        for(int j = 0, p, x; j < m; j++){
            if((x = i * (p = primes[j])) > n) break;

            isNotPrime[x] = true;
            if(i % p == 0) min[x] = min[i] * p; else min[x] = p;

            if(x == min[x]){
                f[x] = calc(p, x);
            } else{
                f[x] = f[ min[x] ] * f[ x / min[x] ];
            }

            if(i % p == 0) break;
        }
    }
}

inline int64 solve(int n){
    return n * (f[n] + 1) >> 1;
}

int main(){
    seive();

    int T_T;

    scanf("%d", &T_T);
    while(T_T--){
        int n;
        
        scanf("%d", &n);

        printf("%lld\n", solve(n));
    }

    return 0;
}