洛谷P3327 [SDOI2015]约数个数和【莫比乌斯反演】

最新推荐文章于 2022-05-12 09:35:21 发布

weixin_30527423

最新推荐文章于 2022-05-12 09:35:21 发布

阅读量77

点赞数

CC 4.0 BY-SA版权

原文链接：http://www.cnblogs.com/Mychael/p/8361823.html

本文介绍了一种解决特定数学问题的方法，该问题是求给定范围内所有数的约数个数之和，并通过莫比乌斯反演实现了高效算法。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题目

设d(x)为x的约数个数，给定N、M，求\(\sum_{i = 1}^{N} \sum_{j = 1}^{M} d(ij)\)

输入格式

输入文件包含多组测试数据。第一行，一个整数T，表示测试数据的组数。接下来的T行，每行两个整数N、M。

输出格式

T行，每行一个整数，表示你所求的答案。

输入样例

2
7 4
5 6

输出样例

110
121

提示

1<=N, M<=50000

1<=T<=50000

题解

好神的题【~~是我太弱吧~~】
首先上来就伤结论。。

题目所求

\[ans = \sum_{i = 1}^{N} \sum_{j = 1}^{M} d(ij)\]

有一个这样的结论：

\[d(ij) = \sum_{x|i}\sum_{y|j} [gcd(x,y) == 1]\]

那么就转化为了：

\[ans =\sum_{i = 1}^{N} \sum_{j = 1}^{M} \sum_{x|i}\sum_{y|j} [gcd(x,y) == 1]\]

我们考虑对于每一对互质的x、y，x会被枚举\(\lfloor \frac{N}{x} \rfloor\)次，y会被枚举\(\lfloor \frac{M}{y} \rfloor\)次
所以有

\[ans =\sum_{i = 1}^{N} \sum_{j = 1}^{M} \lfloor \frac{N}{i} \rfloor \lfloor \frac{M}{j} \rfloor [gcd(i,j) == 1]\]

那么可以进行莫比乌斯反演了
令
\[f(n) = \sum_{i = 1}^{N} \sum_{j = 1}^{M} \lfloor \frac{N}{i} \rfloor \lfloor \frac{M}{j} \rfloor [gcd(i,j) == n]\]

\[F(n) = \sum_{i = 1}^{N} \sum_{j = 1}^{M} \lfloor \frac{N}{i} \rfloor \lfloor \frac{M}{j} \rfloor [n | gcd(i,j)]\]

那么有

\[ \begin{aligned} F(d) &= (\sum_{i = 1}^{N}\lfloor \frac{N}{i} \rfloor) * (\sum_{j = 1}^{M} \lfloor \frac{M}{j} \rfloor) [d | gcd(i,j)] \\ &= (\sum_{i = 1}^{\lfloor \frac{N}{d} \rfloor}\lfloor \frac{N}{id} \rfloor) * (\sum_{j = 1}^{\lfloor \frac{M}{d} \rfloor} \lfloor \frac{M}{jd} \rfloor) \\ &= (\sum_{i = 1}^{\lfloor \frac{N}{d} \rfloor}\lfloor \frac{\lfloor \frac{N}{d} \rfloor}{i} \rfloor) * (\sum_{j = 1}^{\lfloor \frac{M}{d} \rfloor} \lfloor \frac{\lfloor \frac{M}{d} \rfloor}{j} \rfloor) \end{aligned} \]

其中\(\sum_{i = 1}^{N}\lfloor \frac{N}{i} \rfloor\)可以\(O(n\sqrt{n})\)预处理出，我们记为\(sum(n)\)
那么
\[F(n) = sum(\lfloor \frac{N}{n} \rfloor) * sum(\lfloor \frac{M}{n} \rfloor)\]

\[ans = f(1) = \sum_{d = 1}^{N} \mu(d) * F(d) = \sum_{d} \mu(d) sum(\lfloor \frac{N}{d} \rfloor) * sum(\lfloor \frac{M}{d} \rfloor)\]

分块计算
复杂度\(O(T\sqrt{N} + N\sqrt{N})\)

#include<iostream>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define LL long long int
#define REP(i,n) for (int i = 1; i <= (n); i++)
#define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
#define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<' '; puts("");
using namespace std;
const int maxn = 50005,maxm = 100005,INF = 1000000000;
inline int read(){
    int out = 0,flag = 1; char c = getchar();
    while (c < 48 || c > 57) {if (c == '-') flag = -1; c = getchar();}
    while (c >= 48 && c <= 57) {out = (out << 3) + (out << 1) + c - '0'; c = getchar();}
    return out * flag;
}
int prime[maxn],primei,mu[maxn],isn[maxn];
LL sum[maxn];
void init(){
    mu[1] = 1;
    for (int i = 2; i < maxn; i++){
        if (!isn[i]) prime[++primei] = i,mu[i] = -1;
        for (int j = 1; j <= primei && i * prime[j] < maxn; j++){
            isn[i * prime[j]] = true;
            if (i % prime[j] == 0){mu[i * prime[j]] = 0; break;}
            mu[i * prime[j]] = -mu[i];
        }
    }
    for (int i = 1; i < maxn; i++) mu[i] += mu[i - 1];
    for (int n = 1,nxt; n <= 50000; n++){
        for (int i = 1; i <= n; i = nxt + 1){
            nxt = n / (n / i);
            sum[n] += (LL)(nxt - i + 1) * (n / i);
        }
    }
}
int main(){
    init();
    int T = read(),n,m;
    while (T--){
        n = read(); m = read();
        if (n > m) swap(n,m);
        LL ans = 0; int nxt;
        for (int i = 1; i <= n; i = nxt + 1){
            nxt = min(n / (n / i),m / (m / i));
            ans += sum[n / i] * sum[m / i] * (mu[nxt] - mu[i - 1]);
        }
        printf("%lld\n",ans);
    }
    return 0;
}

转载于:https://www.cnblogs.com/Mychael/p/8361823.html