SDOI2015 约数个数和

SDOI2015 约数个数和

标签: 数学 莫比乌斯反演


题目描述

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

输入输出格式

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

输出格式:
T行,每行一个整数,表示你所求的答案。

输入输出样例

输入样例#1:
2
7 4
5 6

输出样例#1:
110
121

数据范围

1<=N, M<=50000
1<=T<=50000

题解

想要做这道题,我们必须对这个约数个数函数有一个了解,存在这么一个公式:
$$d(ij)=\sum_{x|i}\sum_{y|j}[gcd(x, y) == 1]$$
这里有一个比较好的证明
然后原式为:
$$\sum_{i=1}^{N}\sum_{j=1}^{M}\sum_{x|i}\sum_{y|j}[gcd(x,y)==1]$$
考虑更改枚举项,我们直接枚举x, y的倍数:
$$ \sum_{x=1}^{N}\sum_{y=1}^{M}\sum_{i=1}^{\frac{N}{x}}\sum_{j=1}^{\frac{M}{y}}[gcd(x,y)==1]$$
化简
$$ \sum_{x=1}^{N}\sum_{y=1}^{M}\left \lfloor \frac{N}{x} \right \rfloor\left \lfloor \frac{M}{y} \right \rfloor[gcd(x,y)==1]$$
看起来似乎跟莫比乌斯有点关系了
于是我们不妨直接设
$$ f(x)=\sum_{i=1}^{N}\sum_{j=1}^{M}\left \lfloor \frac{N}{i} \right \rfloor\left \lfloor \frac{M}{j} \right \rfloor[gcd(i,j)==x]$$
$$F(x)=\sum_{x|d}f(d)$$
一般都是从$F(x)$入手
$$ F(x)=\sum_{i=1}^{M}\sum_{j=1}^{N}\left \lfloor \frac{N}{i} \right \rfloor\left \lfloor \frac{M}{j} \right \rfloor[x|gcd(i,j)]$$
化简到这里,就很容易联想到更改枚举项,我们枚举x的倍数:
$$ F(x)=\sum_{i=1}^{\frac{N}{x}}\sum_{j=1}^{\frac{M}{x}}\left \lfloor \frac{N}{i} \right \rfloor\left \lfloor \frac{M}{j} \right \rfloor\frac{1}{x^{2}}$$
那么我们预处理一个$g(x)=\sum_{i=1}^{n}\left \lfloor \frac{n}{i} \right \rfloor$
那么$F(x)$就可以在O(1)算出来了
那么关于$g(x)$如何预处理,它其实就是前n个数的约数个数的前缀和,可以线性筛(具体百度)
做了这么多,回到正题:
$$Ans=f(1)$$
$$Ans =\sum_{d=1}^{min(n,m)}u(d)F(d)$$(莫比乌斯反演)
$$Ans =\sum_{x=1}^{min(n,m)}u(x)\sum_{i=1}^{\frac{N}{x}}\sum_{j=1}^{\frac{M}{x}}\left \lfloor \frac{N}{ix} \right \rfloor\left \lfloor \frac{M}{jx} \right \rfloor$$
由于前面处理了一个$g(x)$所以:
$$Ans=\sum_{i=1}^{min(n,m)}u(x)g(\frac{n}{x})g(\frac{m}{x})$$
化简到这里于是就可以完结撒花啦✿✿ヽ(°▽°)ノ✿
最后记着用整除分块优化~~

代码

#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<algorithm>
#define N 60010
#define ll long long
using namespace std;

int t, n, m, cnt = 0;
int u[N], pri[N], fac[N], d[N], sum[N];
ll g[N], ans;
bool vis[N];

void get()
{
    int t = N - 10;
    u[1] = fac[1] = 1;
    for(int i = 2; i <= t; i ++)
    {
        if(!vis[i]) 
        {
            pri[++ cnt] = i; u[i] = -1;
            fac[i] = 2; d[i] = 1;
        }
        for(int j = 1; j <= cnt && i * pri[j] <= t; j ++)
        {
            vis[i * pri[j]] = true;
            if(i % pri[j]) 
            {
                u[i * pri[j]] = -u[i];
                fac[i * pri[j]] = fac[i] * 2;
                d[i * pri[j]] = 1;
            }
            else 
            {
                fac[i * pri[j]] = fac[i] / (d[i] + 1) * (d[i] + 2);
                d[i * pri[j]] = d[i] + 1;
                break;
            }
        }
    } 
    for(int i = 1; i <= t; i ++) 
    {
        g[i] = g[i - 1] + fac[i];
        sum[i] = sum[i - 1] + u[i];
    }
}

int main()
{
    get();
    scanf("%d", &t);
    while(t --)
    {
        scanf("%d%d", &n, &m);
        ans = 0;
        int d = min(n, m);
        for(int l = 1, r; l <= d; l = r + 1)
        {
            r = min(n / (n / l), m / (m / l));
            ans += 1ll * (sum[r] - sum[l - 1]) * g[n / l] * g[m / l];
        }
        printf("%lld\n", ans);
    }
    return 0;
} 

心得

感觉逐渐摸索出了一些更改枚举项的技巧,不过设$f(x)$函数思维还是比较僵化,没有想到直接把那么一大坨式子直接设为$f(x)$。看来以后化简不下去了就直接设函数也是一种方法(然后发现是不会化简…)


 上一篇
bzoj2154 Crash的数字表格 bzoj2154 Crash的数字表格
bzoj2154 Crash的数字表格标签: 数学 莫比乌斯反演 题目描述今天的数学课上,Crash小朋友学习了最小公倍数(Least Common Multiple)。对于两个正整数a和b,LCM(a, b)表示能同时被a和b整除的最小
2019-02-13
下一篇 
POI200 7ZAP-Queries POI200 7ZAP-Queries
[POI2007]ZAP-Queries标签: 数学 莫比乌斯反演 题目描述FGD正在破解一段密码,他需要回答很多类似的问题:对于给定的整数a,b和d,有多少正整数对x,y,满足1<=x<=a,1<=y<=b,并且
2019-02-10