Skip to content

数论分块 学习笔记

数论分块概述

数论分块可以快速计算一些含有除法向下取整的和式(即形如 i=1nf(i)g(ni) 的和式当可以在 O(1) 内计算 f(r)f(l) 或已经预处理出 f 的前缀和时,数论分块就可以在 O(n) 的时间内计算上述和式的值。

它主要利用了富比尼定理(Fubini's theoremni 相同的数打包同时计算。

NOTE

富比尼定理

又称「算两次以意大利数学家圭多・富比尼(Guido Fubini)命名。

富比尼定理的积分形式:只要二重积分 |f(x,y)|dxdy 有界,则可以逐次计算二重积分,并且可以交换逐次积分的顺序。

积分号也是特殊的求和号,因此在一般求和中,富比尼定理往往呈现为更换计数顺序,即交换两个求和号。

组合数学中的富比尼定理表现为,用两种不同的方法计算同一个量,从而建立相等关系。

引理 1

a,b,cZ,abc=abc

证:

ab=ab+r(0r<1)abc=ab1c=1c(ab+r)=abc+rc=abc

证毕

引理 2

nN+,|{nddN+,dn}|2n

|V| 表示集合 V 的元素个数

证:

对于 dnndn 种取值

对于 d>n,有 ndn,也只有 n 种取值

综上,得证

数论分块结论

对于常数 n,使得式子

ni=nj

成立且满足 ijnj 值最大为 nni,即值 ni 所在块的右端点为 nni

证明

​ 令 k=ni,可以知道 kni

nknni=i=ij=max满足条件的所有 i=imax=nk=nni

实现过程

数论分块的过程大概如下:考虑和式

i=1nf(i)ni

那么由于我们可以知道 ni 的值成一个块状分布(就是同样的值都聚集在连续的块中那么就可以用数论分块加速计算,降低时间复杂度。

利用上述结论,我们先求出 f(i)前缀和(记作 s(i)=j=1if(j)然后每次以 [l,r]=[l,nni] 为一块,分块求出贡献累加到结果中即可。

伪代码如下:

1Calculate s(i), the prefix sum of f(i).2l13r04result05while ln do:6rnn/l7resultresult+[s(r)s(l1)]×nl8lr+19end while 

最终得到的 result 即为所求的和式。

向上取整的数论分块

向上取整与前文所述的向下取整十分类似,但略有区别:

对于常数 n,使得式子

ni=nj

成立且满足 ijnj 值最大为 n1n1i,即值 ni 所在块的右端点为 n1n1i

证明

ni=n1i+1,可以发现 n 的上取整分块与 n1 的下取整分块是一样的。

WARNING

i=n 时,上式会出现分母为 0 的错误,需要特殊处理。


N 维数论分块

求含有 a1ia2iani 的和式时,数论分块右端点的表达式从一维的 ni 变为 minj=1n{aji},即对于每一个块的右端点取最小(最接近左端点)的那个作为整体的右端点。可以借助下图理解:

一般我们用的较多的是二维形式,此时可将代码中 r = n / (n / i) 替换成 r = min(n / (n / i), m / (m / i))


应用例 1:Chain Reaction

题意

有一排 n 个怪兽,每个怪兽初始血量为 ai,一次攻击会使一段连续存活的怪兽血量减 k,血量不大于 0 视作死亡,对于所有 k 求出击杀所有怪兽所需攻击次数,n,ai105

思路

太、太、太、太聪明了

对于某个 k,使用差分数组的技巧,令 a0=0,答案就是 i=1nmax(0,aikai1k)

而如果有 aiai1 ,则有 aikai1k0 。所以可以分别计算 aikai1k 在一段 k 上的贡献,然后差分累加出 Ansi 即可。

复杂度 O(ai)。也存在其他解法。

AC 代码

c++
#include <algorithm>
#include <cmath>
#include <ctime>
#include <iostream>
#include <map>
#include <queue>
#include <random>
#include <set>
#include <string>
#include <vector>

using namespace std;

mt19937_64 rnd(time(0));

typedef long long ll;
typedef uint64_t ull;

const ll INF = 1e18, Mod = 1e9 + 7;

ll n, m;
ll h[100010];
ll diff[100010];

int main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    cout.tie(0);

    cin >> n;

    for (ll i = 1; i <= n; i++)
        cin >> h[i], m = max(m, h[i]);

    // cout<<m<<endl;

    for (ll i = 1; i <= n; i++) {
        if (h[i] > h[i - 1]) {
            ll l = 1, r = 0;
            while (l <= h[i]) {
                if (l == h[i]) {
                    diff[l] += 1;
                    break;
                }
                r = min(m, (h[i] - 1) / ((h[i] - 1) / l));

                // cout<<l<<' '<<r<<endl;
                diff[l] += (h[i] + l - 1) / l;
                diff[r + 1] -= (h[i] + l - 1) / l;

                l = r + 1;
            }
            l = 1, r = 0;
            // cout<<l<<' '<<r<<endl;
            while (l <= h[i - 1]) {
                if (l == h[i - 1]) {
                    diff[l] -= 1;
                    break;
                }
                r = min(m, (h[i - 1] - 1) / ((h[i - 1] - 1) / l));

                // cout<<l<<' '<<r<<endl;
                diff[l] -= (h[i - 1] + l - 1) / l;
                diff[r + 1] += (h[i - 1] + l - 1) / l;

                l = r + 1;
            }
        }
    }

    ll ans = 0;
    for (ll i = 1; i <= m; i++) {
        ans += diff[i];
        cout << ans << ' ';
    }
    cout << endl;
    return 0;
}

应用例 2:Monster

思路

枚举 加 k 攻击力 的模块的个数 然后用数论分块来在区间内跳跃遍历答案。

AC 代码

c++
#include <algorithm>
#include <cmath>
#include <ctime>
#include <iostream>
#include <map>
#include <queue>
#include <random>
#include <set>
#include <string>
#include <vector>
#define OOO cout << ">>>>"; // 调试标识

using namespace std;

mt19937_64 rnd(time(0));

typedef long long ll;
typedef uint64_t ull;

const ll INF = 1e18, Mod = 1e9 + 7;

ll tt;

int main() {
    ios::sync_with_stdio(0);
    cin.tie(0);
    cout.tie(0);

    cin >> tt;

    while (tt--) {
        ll x, y, z, k, ans = INF;
        cin >> x >> y >> z >> k;
        for (ll r = 0;; r++) {
            ll h = z - r * (r + 1) * k / 2, dmg = r * k,
               cost = r * (k * x + y);
            if (h <= 0) {
                ans = min(ans, cost);
                break;
            } else
                for (ll l = max(dmg, 1ll); l <= min(dmg + k, h);) {
                    ans = min(ans, cost + (l - dmg) * x +
                                       (ll)ceil(1.0 * h / l) * y);
                    l = l ^ h ? ((h - 1) / ((h - 1) / l) + 1) : (h + 1);
                }
        }
        cout << ans << endl;
    }
    cout.flush();
    return 0;
}

其中保留了部分调试输出,可以使用它们来查看具体的处理过程,以深化对算法的理解。

数论分块扩展

以计算含有 nd 的和式为例。考虑对于一个正整数 n,如何求出集合

S={nddN+,dn}

的所有值,以及对每一种值求出哪些 d 会使其取到这个值。可以发现:

  1. 因为 nd 是单调不增的,所以对于所有 vS,使得 nd=vd 必然是一段区间。
  2. 对于任意正整数 tn,我们对 t>tvS 分别分析,可以发现 t+n/t2|S|,取 t=n3 得到 |S| 的一个上界为 O(n3)

这些结论与数论分块所需的引理相似,因此猜测可以写为数论分块形式。

结论是:使得式子

np=nq

成立的最大的 q 满足 pqn

nn/p2

证明

​ 令 v=np=nq,那么

vnqv2nqqnv2qnv2

同理 pn/v2。同时

nn/v2nn/v2=v=v

又由 pn/v2 以及单调性可推出

v=npnn/v2

所以

nn/v2=v

进而 q=n/v2 是最大的使得 n/p=n/q 成立的 q


故原问题可以写为数论分块形式,代码与数论分块形式并无二异。

两个更加通用的结论

  1. 对于某个不超过 nα/β 的正整数 i,使式子 nαiβ=nαjβ 成立的最大的 jnα/βnα/iβ1/β
  2. 集合 {nαdβ:d=1,2,,n} 的大小不超过 min{n,2nα/(1+β)}

习题

  1. CQOI2007 余数求和(需要一点转化和特判)

  2. UVa11526 H(n)(几乎可以当做模板题)

  3. POI2007 ZAP-Queries(数论分块一般配合 莫比乌斯反演 用以进一步降低复杂度;本题需要用到 [n=1]=d|nμ(d) 这一条莫反结论)