BZOJ 3994 约数个数和(莫比乌斯反演)

题目链接：
BZOJ 3994 约数个数和
题意：
$求\sum_{i=1}^{n}\sum_{j=1}^{m}d(i*j),定义d(i)为i的约数个数.n,m\in [1, 50000]$
分析：

$$ans =\sum_{gcd(i,j)=1}\lfloor{\frac{n}{i}}\rfloor \lfloor{\frac{m}{j}}\rfloor = \sum_{i=1}^{n}\lfloor{\frac{n}{i}}\rfloor\sum_{j=1}^{m} \lfloor{\frac{m}{j}}\rfloor$$
$定义g(n)=\sum_{i=1}^{n}\lfloor{\frac{n}{i}}\rfloor=\sum_{i= 1}^{n}d(i)$,只需要预处理出 $g(n)$就可以在 $O(\sqrt n)$时间范围内解决问题。如果选择分步加速的话，预处理的复杂度是 $O(n\sqrt n)$,但是其实我们考虑每个数 $i$的约数个数，然后$g(n)$就是前缀和了。在线性筛时每个合数是被最小质因子筛掉的,我们只需要记录这个每个数 $m$最小质因子的幂次$num[m],d[m]$记录 $m$的约数个数,显然$d(m)$是积性函数.对于 $m=i*prime[j]$，如果 $i\ \% \ prime[j]！=0$,根据积性函数性质 $num[m]=1,d[m]=d[i]*d[prime[j]],$否则 $num[m]=num[i]+1,$因为 $m$的约数个数是:$(e_1+1)*(e_2+1)*...*(e_k+1),e_i$是质因子分解后各质因子的幂次，我们考虑 $m$从$i$的转移过程， $m$只比$i$在 $prime[j]$的幂次上多 $1$,所以$d[m]=\frac{d[i]}{num[i]+1]}*(num[i]+2)$.

亲测$g$如果分步加速预处理$6000+MS,$用线性筛$2000+MS$

#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <bitset>
#include <string>
using namespace std;
typedef long long ll;
const int MAX_N = 50010;

int prime_cnt, prime[MAX_N], mu[MAX_N], num[MAX_N];
ll g[MAX_N], d[MAX_N], sum[MAX_N];
bitset<MAX_N> bs;

void GetMu()
{
    bs.set();
    prime_cnt = 0;
    mu[1] = d[1] = 1;
    for(int i = 2; i < MAX_N; ++i) {
        if(bs[i]) {
            prime[prime_cnt++] = i;
            mu[i] = -1;
            num[i] = 1;
            d[i] = 2;
        }
        for(int j = 0; j < prime_cnt && i * prime[j]  < MAX_N; ++j) {
            bs[i * prime[j]] = 0;
            if(i % prime[j]) {
                mu[i * prime[j]] = -mu[i];
                num[i * prime[j]] = 1;
                d[i * prime[j]] = d[i] * d[prime[j]];
            }else {
                mu[i * prime[j]] = 0;
                num[i * prime[j]] = num[i] + 1;
                d[i * prime[j]] = d[i] / (num[i] + 1) * (num[i] + 2);
                break;
            }
        }
    }
    for(int i = 1; i < MAX_N; ++i) {
        sum[i] = sum[i - 1] + mu[i];
        g[i] = g[i - 1] + d[i];
    }
}

inline ll solve(int n, int m)
{
    ll res = 0;
    int top = min(n, m), last;
    for(int i = 1; i <= top; i = last + 1) {
        last = min(n / (n / i), m / (m / i));
        res += (sum[last] - sum[i - 1]) * g[n / i] * g[m / i];
    }
    return res;
}

int main()
{
    GetMu();
    int T, n, m;
    scanf("%d", &T);
    while(T--) {
        scanf("%d%d", &n, &m);
        printf("%lld\n", solve(n, m));
    }
    return 0;
}