HAOI2011 - Problem b

【LuoguP2522】【BZOJ2301】【HAOI2011】Problem b
题目链接1: https://www.luogu.org/problem/P2522
题目链接2: https://www.lydsy.com/JudgeOnline/problem.php?id=2301

前言

机房很久前就人均会反演，今天我终于看懂是啥了。

然后把整除分块理解透彻！

先做道最简单的题练练手，感谢莫比乌斯反演 - OI WIKI和大佬CTZ的帮助。

分析

本题题目翻译成公式：

$$ \sum_{x=a}^{b}\sum_{y=c}^{d} [gcd(x,y)=k] $$

可以按照二维前缀和的思想分治一下，先设

$$ F(n,m) = \sum_{i=1}^{n}\sum_{j=1}^{m} [gcd(i,j)=k] $$

那么

$$ ans = F(b,d) - F(b,c-1) - F(a-1,d) + F(a-1,c-1) $$

然后对F(n,m)进行化简，首先两边同除以K，至于可行性，记录一下大佬的亲自讲解

如果原式中$gcd(i,j)=k$，那么$i,j$一定是k的倍数，现在设

$$ i = i’ * k , j = j’ * k $$

那么一定$gcd(i’,j’)=1$，所以得出

$$ F(n,m) = \sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{k} \rfloor} [gcd(i,j)=1] $$

这就可以莫比乌斯反演了，首先记住一个通用公式

$$ [gcd(i,j)=1] = \sum_{d \mid gcd(i,j)}\mu(d) $$

然后代入F(n,m)中

$$ F(n,m) = \sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{k} \rfloor} \sum_{d \mid gcd(i,j)}\mu(d) $$

先默认$n < m$，那么可知$d \in [1,\lfloor \frac{n}{k} \rfloor]$

然后改变一下位置

$$ F(n,m) = \sum_{d=1}^{\lfloor \frac{n}{k} \rfloor}\mu(d)\sum_{i=1}^{\lfloor \frac{n}{k} \rfloor}[i \mid p]\sum_{j=1}^{\lfloor \frac{m}{k} \rfloor}[j \mid p] $$

注意，右边是一个整体，后面的$[i\mid d],[j\mid d]$取1或0，也就是表示能否整除。

这样还能进一步化简，根据常识易得

$$ \sum_{i=1}^{n} i \mid p = \lfloor \frac{n}{p} \rfloor$$

从而推出

$$ F(n,m) = \sum_{d=1}^{\lfloor \frac{n}{k} \rfloor}\mu(d)\lfloor \frac{n}{kd} \rfloor{\lfloor \frac{m}{kd} \rfloor} $$

式子终于推完了，然后就是整除分块。

其实我一开始没看懂整除分块是因为不知道有啥用，其实是用来$O(\sqrt{n})$求这个式子

$$ \sum_{i=1}^{n} \lfloor \frac{i}{p} \rfloor $$

这样就好说了，我们知道某些$\lfloor \frac{i}{p} \rfloor$的值是相等的，并且这些项在一个连续区间$[l,r]$内。

所以枚举时找一个右端点$r$，直接计算$[i,r]$这个区间的结果，这个$r= \lfloor \frac{n}{\lfloor \frac{n}{i} \rfloor} \rfloor$。

具体证明不难，数论分块 - OI WIKI有略证。

然后就写代码，还需要写莫比乌斯函数筛、前缀和。

这个莫比乌斯函数筛我有待研究。

代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43

#include <iostream>
#include <stdio.h>
#define MAXA 50001
using namespace std;
int miu[50005], notzs[50005], zs[50005], sum[50005];
void shai() {
	notzs[1]= 1, miu[1]= 1;
	for(int i= 2; i <= MAXA; i++) {
		if(!notzs[i]) miu[i]= -1, zs[++zs[0]]= i;
		for(int j= 1; j <= zs[0] && i * zs[j] <= MAXA; j++) {
			notzs[i * zs[j]]= 1;
			if(i % zs[j] == 0) {
				miu[i * zs[j]]= 0;
				break;
			}
			miu[i * zs[j]]= -miu[i];
		}
	}
	for(int i= 1; i <= MAXA; i++) sum[i]= sum[i - 1] + miu[i];
	return;
}
int query(int n, int m, int k) {
	int res= 0;
	n/= k, m/= k;
	if(n > m) swap(n, m);
	for(int i= 1, j; i <= n; i= j + 1) {
		j= min(n / (n / i), m / (m / i));
		res+= (sum[j] - sum[i - 1]) * (n / i) * (m / i);
	}
	return res;
}
#define query(x, y) query(x, y, k)
int t;
int main() {
	shai();
	cin >> t;
	int a, b, c, d, k;
	while(t--) {
		cin >> a >> b >> c >> d >> k;
		cout << query(b, d) - query(b, c - 1) - query(a - 1, d) + query(a - 1, c - 1) << endl;
	}
	return 0;
}

后记

然后你还可以把【LuoguP3455】【POI2007】ZAP-Queries 一块通过了，只需要改一下输入输出。