【数学】数论分块（整除分块）

Description

数论分块，通常用于快速求解形如 $\sum _{i=1}^{n}f(i)\cdot g\left(\lfloor \frac{n}{i}\rfloor \right)$ 的和式，所以通常被称为 整除分块，当能用 $O(1)$ 计算出 $\sum _{i=l}^{r}f(i)$ 时，数论分块便能用 $O(\sqrt{n})$ 的时间计算出上式的值、数论分块经常搭配 莫比乌斯反演 一起使用。

Conclusion

$\forall n,l\in {\mathbb{N}}^{\ast },l\le n$，使得
$$\lfloor \frac{n}{l}\rfloor =\lfloor \frac{n}{r}\rfloor$$
成立的最大的满足 $l\le r\le n$ 的 $r$ 的值为 $\lfloor \frac{n}{\lfloor \frac{n}{l}\rfloor }\rfloor$。

也就是说，若某个块内所有数的值为 $\lfloor \frac{n}{l}\rfloor =k$，那么这个块的右端点就是 $r=\lfloor \frac{n}{k}\rfloor$。

Proof

对于这个块中的任意一个数 $x$，应当满足 $\lfloor \frac{n}{x}\rfloor =k$，即 $n=xk+r(0\le r<x)$。

当 $n,k$ 已知时，只要确定 $x$，就有一个 $r$ 与之对应。
$$\begin{gathered} n=xk+r \\ n\ge xk \\ x\le \frac{n}{k} \end{gathered}$$
即 $x$ 的最大值为 $\lfloor \frac{n}{k}\rfloor$。

即 $r=\lfloor \frac{n}{\lfloor \frac{n}{l}\rfloor }\rfloor$。

证毕。

Method

一个块内的所有数都相等，所以每块每块地求和即可。

在找到上一个块的右端点后，加一就可以得到下一个块的左端点。

UVA11526 H(n)

模板题，求 $\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor$。

//18 = 9 + 9 = 18.
#include <iostream>
#include <cstdio>
#define Debug(x) cout << #x << "=" << x << endl
#define int long long
using namespace std;

int H(int n)
{
	int res = 0;
	for (int l = 1, r; l <= n; l = r + 1)
	{
		int k = n / l;
		r = n / k;
		res += k * (r - l + 1);
	}
	return res;
}

signed main()
{
	int t;
	scanf("%lld", &t);
	while (t--)
	{
		int n;
		scanf("%lld", &n);
		printf("%lld\n", H(n));
	}
	return 0;
}

Complexity

Lemma

$\forall n,i\in {\mathbb{N}}^{\ast },i\le n$，$\lfloor \frac{n}{i}\rfloor$ 的不同值最多有 $2\lfloor \sqrt{n}\rfloor$ 个，即最多有 $2\lfloor \sqrt{n}\rfloor$ 个块。

Proof

$\forall i\le \lfloor \sqrt{n}\rfloor$，$i$ 只有 $\lfloor \sqrt{n}\rfloor$ 种取值，则 $\lfloor \frac{n}{i}\rfloor$ 只有至多 $\lfloor \sqrt{n}\rfloor$ 种取值；

$\forall i>\lfloor \sqrt{n}\rfloor$，有 $i\ge \lfloor \sqrt{n}\rfloor +1>\sqrt{n}$，$\lfloor \frac{n}{i}\rfloor \le \lfloor \frac{n}{\sqrt{n}}\rfloor =\lfloor \sqrt{n}\rfloor$，也只有至多 $\lfloor \sqrt{n}\rfloor$ 种取值。

所以最多只有 $2\lfloor \sqrt{n}\rfloor$ 种取值。

证毕。

---

由引理可知最多有 $2\lfloor \sqrt{n}\rfloor$ 个块，即 $\mathrm{for}$ 循环最多会执行 $2\lfloor \sqrt{n}\rfloor$ 次，所以时间复杂度为 $O(\sqrt{n})$。

Extension

$N$ 维数论分块。

求形如 $\sum _{i=1}^n\lfloor \frac{a_1}{i}\rfloor \lfloor \frac{a_2}{i}\rfloor \cdots \lfloor \frac{a_m}{i}\rfloor$ 的和式的值。

$r=\underset{i=1}{\overset{m}{\min }}\left\{\lfloor \frac{a_i}{\lfloor \frac{a_i}{l}\rfloor }\rfloor \right\}$ 即可，即对于每一个块的右端点取最小（最接近左端点）的那个作为整体的右端点。

较常用的是 $2$ 维数论分块。求 $\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor$ 时 $r=\min \left(\lfloor \frac{n}{\lfloor \frac{n}{l}\rfloor }\rfloor ,\lfloor \frac{m}{\lfloor \frac{m}{l}\rfloor }\rfloor \right)$，也就是在代码中 r = min(n / (n / l), m / (m / l));。

Problems

A

P2261 [CQOI2007]余数求和

Description

计算 $\sum _{i=1}^{n}k\mathrm{m}\mathrm{o}\mathrm{d}i$。

Solution

$$\begin{array}{rl}\sum _{i=1}^{n}k\mathrm{m}\mathrm{o}\mathrm{d}i& =\sum _{i=1}^{n}k-\lfloor \frac{k}{i}\rfloor i\\ & =nk-\sum _{i=1}^n\lfloor \frac{k}{i}\rfloor i\end{array}$$

观察后面一项，对于左端点为 $l$，右端点为 $r$ 的块，贡献就是

$$\begin{array}{rl}\sum _{i=l}^r\lfloor \frac{k}{i}\rfloor i& =\sum _{i=l}^r\lfloor \frac{k}{l}\rfloor i\\ & =\lfloor \frac{k}{l}\rfloor \sum _{i=1}^{r}i\\ & =\lfloor \frac{k}{l}\rfloor \frac{(l+r)(r-l+1)}{2}\end{array}$$
但是这次 $i$ 要循环到 $n$ 而非 $k$，当 $n>k$ 时 $\lfloor \frac{k}{l}\rfloor$ 有可能为 $0$，这样 $r=\lfloor \frac{k}{\lfloor \frac{k}{l}\rfloor }\rfloor$ 就无意义了。

发现当 $i>k$ 时 $\lfloor \frac{k}{i}\rfloor i=0$，所以 $i$ 循环到 $\min (n,k)$ 即可。

在取右端点时 $r=\min \left(n,\lfloor \frac{k}{\lfloor \frac{k}{l}\rfloor }\rfloor \right)$。

Code

//18 = 9 + 9 = 18.
#include <iostream>
#include <cstdio>
#define Debug(x) cout << #x << "=" << x << endl
typedef long long ll;
using namespace std;

ll block(int n, int k)
{
	ll res = 0;
	for (int l = 1, r; l <= min(n, k); l = r + 1)
	{
		r = min(n, k / (k / l));
		res += (ll)(k / l) * ((ll)(l + r) * (r - l + 1) >> 1);
	}
	return res;
}

int main()
{
	int n, k;
	scanf("%d%d", &n, &k);
	printf("%lld\n", (ll)n * k - block(n, k));
	return 0;
}

B

P2260 [清华集训2012]模积和

Description

求
$$\left[\sum _{i=1}^n\sum _{j=1}^m(n\mathrm{m}\mathrm{o}\mathrm{d}i)(m\mathrm{m}\mathrm{o}\mathrm{d}j),i\ne j\right]\mathrm{m}\mathrm{o}\mathrm{d}19940417$$

Solution

假设 $n\le m$（否则交换）。
$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^m(n\mathrm{m}\mathrm{o}\mathrm{d}i)(m\mathrm{m}\mathrm{o}\mathrm{d}j)-\sum _{i=1}^n(n\mathrm{m}\mathrm{o}\mathrm{d}i)(m\mathrm{m}\mathrm{o}\mathrm{d}i)\\ & =\sum _{i=1}^{n}n\mathrm{m}\mathrm{o}\mathrm{d}i\cdot \sum _{j=1}^{m}m\mathrm{m}\mathrm{o}\mathrm{d}j-\sum _{i=1}^n(n\mathrm{m}\mathrm{o}\mathrm{d}i)(m\mathrm{m}\mathrm{o}\mathrm{d}i)\end{array}$$
前面两项都形如 $\sum _{i=1}^{k}k\mathrm{m}\mathrm{o}\mathrm{d}i$，用上一题的思路求解。
$$\begin{array}{rl}\sum _{i=1}^n(n\mathrm{m}\mathrm{o}\mathrm{d}i)(m\mathrm{m}\mathrm{o}\mathrm{d}i)& =\sum _{i=1}^n(n-\lfloor \frac{n}{i}\rfloor i)(m-\lfloor \frac{m}{i}\rfloor i)\\ & =\sum _{i=1}^n(nm-\lfloor \frac{n}{i}\rfloor i\cdot m-\lfloor \frac{m}{i}\rfloor i\cdot n+\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor i^2)\\ & =n^{2}m-m\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor i-n\sum _{i=1}^n\lfloor \frac{m}{i}\rfloor i+\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor i^2\end{array}$$
第 $2,3$ 项都是模板，第 $4$ 项就是一个扩展版中的 $2$ 维数论分块，对于左端点为 $l$，右端点为 $r$ 的块，贡献就是
$$\begin{array}{rl}\sum _{i=l}^r\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor i^2& =\lfloor \frac{n}{l}\rfloor \lfloor \frac{m}{l}\rfloor \sum _{i=l}^r{i}^2\\ & =\lfloor \frac{n}{l}\rfloor \lfloor \frac{m}{l}\rfloor (\sum _{i=1}^r{i}^2-\sum _{i=1}^{l-1}{i}^2)\end{array}$$
平方和有公式 $\sum _{i=1}^n{i}^2=\frac{n(n+1)(2n+1)}{6}$。

Code

函数 $\mathrm{b}\mathrm{l}\mathrm{o}\mathrm{c}\mathrm{k}1(n,m)$ 求的是 $\sum _{i=1}^n\lfloor \frac{m}{i}\rfloor i$（其中保证 $n\le m$）。

函数 $\mathrm{b}\mathrm{l}\mathrm{o}\mathrm{c}\mathrm{k}2(n,m)$ 求的是 $\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor \lfloor \frac{m}{i}\rfloor i^2$。

//18 = 9 + 9 = 18.
#include <iostream>
#include <cstdio>
#define Debug(x) cout << #x << "=" << x << endl
#define int long long
using namespace std;

const int MOD = 19940417;

int x, y;

void exgcd(int a, int b)
{
	if (!b)
	{
		x = 1, y = 0;
		return;
	}
	exgcd(b, a % b);
	int tmp = x;
	x = y;
	y = tmp - a / b * y;
}

int inv(int a)
{
	exgcd(a, MOD);
	x = (x % MOD + MOD) % MOD;
	return x;
}

int block1(int n, int m)
{
	int res = 0;
	for (int l = 1, r; l <= n; l = r + 1)
	{
		r = min(n, m / (m / l));
		res = (res + (m / l) * ((l + r) * (r - l + 1) / 2 % MOD) % MOD) % MOD;
	}
	return res;
}

int part1(int n, int m)
{
	return ((n * n % MOD - block1(n, n)) * (m * m % MOD - block1(m, m)) % MOD + MOD) % MOD;
}

int sum(int n)
{
	return n * (n + 1) % MOD * (2 * n + 1) % MOD * inv(6) % MOD;
}

int block2(int n, int m)
{
	int res = 0;
	for (int l = 1, r; l <= n; l = r + 1)
	{
		r = min(n / (n / l), m / (m / l));
		res = (res + (n / l) * (m / l) % MOD * (sum(r) - sum(l - 1)) % MOD) % MOD;
	}
	return res;
}

int part2(int n, int m)
{
	int a = n * n % MOD * m % MOD, b = block1(n, n) * m % MOD, c = block1(n, m) * n % MOD, d = block2(n, m);
	return ((a - b - c + d) % MOD + MOD) % MOD;
}

signed main()
{
	int n, m;
	scanf("%lld%lld", &n, &m);
	if (n > m)
	{
		swap(n, m);
	}
	printf("%lld\n", ((part1(n, m) - part2(n, m)) % MOD + MOD) % MOD);
	return 0;
}

Reference

• [1] OI Wiki：数论分块
• [2] *ACoder*：总结与思考：数论分块