算法导论 第二部分:排序和顺序统计量

堆排序

堆

(二叉)堆是一个数组,它可以被看成一个近似的完全二叉树。

1 2 3 4 5 6 7 8 9 10
16 14 10 8 7 9 3 2 4 1
graph TD
    1((16)) --> 2((14))
    1 --> 3((10))
    2 --> 4((8))
    2 --> 5((7))
    3 --> 6((9))
    3 --> 7((3))
    4 --> 8((2))
    4 --> 9((4))
    5 --> 10((1))

对于一个数组 , 是数组的元素个数, 中存放的是堆的有效元素,这里 ,我们很容易得到对于结点 的父结点、左孩子和右孩子的下标。
注:与书上不同,我们使用从 开始的堆(0-based)

1
2
3
4
5
6
7
8
9
10
11
12
13
inline static int PARENT(int i) {
    return (i - 1) / 2;
}
inline static int LEFT(int i) {
    return 2 * i + 1;
}
inline static int RIGHT(int i) {
    return 2 * i + 2;
}
/**
i=2时,父结点应该是0
i=0时,左结点应该是1,右结点应该是2
**/

二叉堆可以分为两种:最大堆和最小堆。在这两种堆中结点都要满足堆的性质。
在最大堆中,最大堆的性质是指除了根结点以外的所有结点 都要满足

与此相反,最小堆性质是指除了根结点以外的所有结点 都要满足

如果把堆看成一棵树,那么我们定义一个堆中结点的高度是该结点到叶结点最长简单路径上边的数目。堆的高度就是根结点的高度。

维护堆的性质

根据最大堆性质,我们需要确保当前结点 的值是 、 和 的最大者。同时,对于发生了交换的情况,我们还需要递归检查交换的子结点下是否还满足堆的性质。
我们给出代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
template<typename T> void MAX_HEAPIFY(vector <T>&a, int i, int heap_size) {
    int l = LEFT(i);
    int r = RIGHT(i);
    int max = i;
    if (l < heap_size && a[l] > a[i]) {
        max = l;
    }
    if (r < heap_size && a[r] > a[max]) {
        max = r;
    }
    if (max != i) {
        swap(a[max], a[i]);
        MAX_HEAPIFY(a, max, heap_size);
    }
}

有递归式

根据主定理,有

即对于高度为 的结点来说,它的时间复杂度是

建堆

只需要从叶结点的父结点开始从下往上维护最大堆即可

1
2
3
4
5
6
7
8
template<typename T> void BUILD_MAX_HEAP(vector <T>&a) {
    for (int i = a.size() / 2; i >= 0; i--) {
        MAX_HEAPIFY(a, i, a.size());
    }
}
/** 来自ChatGPT:
0-based 最后一个非叶子是 `heap_size/2 - 1`,你用了 `heap_size/2`,会多 heapify 一次(不致命但没必要)
**/

堆排序算法

首先把数组建成一个最大堆,由于最大的元素始终放在 ,因此只需从 开始往前遍历交换元素,并逐渐缩小堆的大小就可以实现排序。

1
2
3
4
5
6
7
8
9
template<typename T> void HEAP_SORT(vector <T>&a) {
    BUILD_MAX_HEAP(a);
    int heap_size = a.size();
    for (int i = heap_size - 1; i > 0; i--) {
        swap(a[i], a[0]);
        heap_size -= 1;
        MAX_HEAPIFY(a, 0, heap_size);
    }
}

优先队列

优先队列(priority queue) 是一种用来维护由一组元素构成的集合 的数据结构,其中每一个元素都会有一个相关的值,被称为 关键字(key)。一个 最大优先队列 支持以下操作:

  • 把元素 插入到 中
  • 返回 中具有最大关键字的元素
  • 去掉并返回 中具有最大关键字的元素
  • 将元素 的关键字增加到 ,这里 不小于 的原关键字值
    最大优先队列的用途有很多,比如调度器要选择优先级最大的任务来执行。
    相应地,最小优先队列支持的操作包括 、、和
    。最小优先队列可以被用于基于事件驱动的模拟器。队列中保存要模拟的事件,每个事件都有一个发生时间作为关键字,事件必须按照发生的时间顺序进行模拟。
    显然优先队列可以用堆来实现,实现时,我们需要在堆中的每个元素中存储对应对象的句柄(handle),通常是数组的下标,但具体依赖于应用场景。
    在 时间内实现 操作
1
2
3
template<typename T> T HEAP_MAXMIUM (vector <T>&a) {
    return a[0];
}

在 时间内实现 操作

1
2
3
4
5
6
7
8
9
10
template<typename T> T HEAP_EXTRACT_MAX(vector <T>&a) {
    int heap_size = a.size();
    if (heap_size < 1) { throw runtime_error("heap underflow"); }
    T max = a[0];
    a[0] = a[heap_size - 1];
    heap_size -= 1;
    a.pop_back();
    MAX_HEAPIFY(a, 0, heap_size);
    return max;
}

在 时间内实现 操作。在更改元素的关键字后,可能会破坏最大堆性质,所以我们像插入排序那样,将 与 层层比较,找到新的合适的位置。

1
2
3
4
5
6
7
8
template<typename T> void HEAP_INCREASE_KEY(vector <T>&a, int i, T key) {
    if (key < a[i]) { throw runtime_error("New key is smaller than the old key"); }
    a[i] = key;
    while (i > 0 && a[i] > a[PARENT(i)]) {
        swap(a[i], a[PARENT(i)]);
        i = PARENT(i);
    }
}

在 时间内实现 操作

1
2
3
4
5
6
7
8
9
template<typename T> void MAX_HEAP_INSERT(vector <T>&a, T key, int& heap_size) {
    heap_size += 1;
    a.push_back(key);
    int i = heap_size - 1;
    while (i > 0 && a[i] > a[PARENT(i)]) {
        swap(a[i], a[PARENT(i)]);
        i = PARENT(i);
    }
}

总之,在一个包含 个元素的堆中,所有优先队列的操作都可以在 的时间内完成。

快速排序

快速排序的描述

快速排序也是一个分治算法,相较于归并排序,快速排序拥有更小的常数。
分解: 将数组分为 和 两部分,使得前面的数组中的每一个元素都小于等于 ,而 也小于等于 中的每一个元素。计算下标 是分解的一部分。
解决: 递归调用快速排序,对子数组 和 进行排序。
合并: 因为子数组都是原址排序的,不需要合并操作。
下面是实现:

1
2
3
4
5
6
7
template <typename T> void QUICKSORT(vector <T>&a, int p, int r) {
    if (p < r) {
        int q = PARTITION(a, p, r);
        QUICKSORT(a, p, q-1);
        QUICKSORT(a, q+1, r);
    }
}

数组划分

1
2
3
4
5
6
7
8
9
10
11
12
template <typename T> int PARTITION(vector <T>&a, int p, int r) {
    T x = a[r];
    int i = p - 1;
    for (int j = p; j <= r-1; j++) {
        if (a[j] <= x) {
            i += 1;
            swap(a[i], a[j]);
        }
    }
    swap(a[i+1], a[r]);
    return i + 1;
}

划分函数 先是选择了一个 主元(pivot) ,并围绕它来划分数组 。 的时候,将 递增 并将 和 交换,保证 中的元素总是不大于 (相当于把新的元素并入了),而且 的元素总是大于 。 将从 到 遍历,最后得到的 和 数组就是分解所要求的。最后交换 与 使得主元 位于中间。

p … i … j … r
x(主元)
无限制

作为循环不变量,循环每一轮都满足对于任意数组下标 ,有:

  1. 若 ,则
  2. 若 ,则
  3. 若 ,则

我们拿循环不变式来证明一下:
初始化: 在循环第一次开始前, , 之间没有元素,明显满足循环不变量。
保持: 考虑两种情况:当 时,我们只对 ,之后条件2仍然成立,其他条件不变。当 时,将 并交换了 和 ,并且 ,很明显三个条件仍然满足。
终止: 终止时,,并将 交换,并返回主元的下标。
实际上, 严格小于 的每一个元素。

快速排序的性能

最坏情况划分

当 将数组分为了 和 个元素时,是快速排序的最坏情况。假设每次划分都出现这种划分,划分的时间复杂度为 ,递归式可以表示为 ,解得

最好情况划分

假设划分的比例为 ,则递归式为

写出递归树可以看出,递归在 处停止,每层代价为 ,因此总代价为 。
只要划分是常数比例的,算法运行时间总是

对于平均情况的直观观察

当出现一次差的划分 和一次好的划分 时,,差的划分的代价可以吸收到好的划分的代价,仍然是 。可以做出猜测平均情况也是

快速排序的随机化版本

随机选取一个元素与交换

1
2
3
4
5
6
7
template<typename T> int RAND_PARTITION(vector<T>&a, int p, int r) {
    std::random_device rd;
    std::mt19937 gen(rd());
    std::uniform_int_distribution<> distrib(p, r);
    swap(a[distrib(gen)], a[r]);
    return PARTITION(a, p, r);
}

同时修改原算法

1
2
3
4
5
6
7
template <typename T> void RAND_QUICKSORT(vector <T>&a, int p, int r) {
    if (p < r) {
        int q = RAND_PARTITION(a, p, r);
        RAND_QUICKSORT(a, p, q-1);
        RAND_QUICKSORT(a, q+1, r);
    }
}

快速排序分析

最坏情况分析

设最坏的情况为,有递归式

子问题规模加和为 ,则 是从 变化到 ,不妨猜测 成立,则

式子 在端点取到最大值,且这个式子二阶导数是正的,表达式的上界为 代入 中,我们有

故
还可以证明,所以运行时间是

期望运行时间

快速排序运行期间,最多会调用 次 函数,每次运行该函数就需要执行若干次for循环,该函数运行时间取决于判断 a[j] <= x 被执行的次数。那么很明显有下面的引理。
引理7.1: 快速排序在包含 个元素数组上运行时,设在 函数循环中判断被执行的次数为 ,则整个算法运行时间为
我们令 中第 小的元素为 ,并记 为包含 和他们之间的元素组成的集合。
我们使用指示器随机变量

注意到对于每个 ,算法至多比较一次。总比较次数为

两边取期望

对于概率那一项,我们先来看什么时候不会发生比较。我们假设元素的值都是互异的,则当 的主元被选中后, 和 以后不会再发生比较。当主元 被选为主元,则它将与除了它自己的其他所有元素比较,对于 也是。因此,当且仅当 中将被选中主元的元素是 或 时 和 会发生比较。
让我们来计算概率,因为主元是随机选择的,在 中没有被选中主元时,每个元素都等可能且独立地被选中,概率为 ,于是

代回,有

令 ,

因此在输入元素互异时,快速排序的平均运行时间为

线性时间排序

排序算法的下界

上面介绍的算法都是比较排序,即在排序的最终结果中,各元素的次序依赖于它们之间的比较。在这种算法中,我们只能使用元素之间的相互比较来获得次序信息。我们假设这里的比较操作全是 且输入元素全是互异的。

决策树模型

比较排序可以抽象为一颗决策树,其中每个内部结点以 标记,其中 , 是输入元素个数。每个结点表示比较 与 。以 为例,左子树表示确定了 之后的比较,右子树表示确定了 之后的比较,算法沿着决策树运行,每到达一个叶结点,就代表确定了一个顺序 。对于正确的比较排序算法, 个元素的 种排列都是可达的。

最坏情况的下界

比较排序算法中最坏情况比较次数就等于其决策树的高度。当决策树中每种排列都以可达叶节点的形式出现,该决策树高度的下界就是算法运行时间的下界。
定理8.1: 在最坏的情况下,任何比较排序算法都需要做 次比较。
证明: 考虑一棵高度为 有 个可达叶结点的决策树,它对应一个由 个元素所做的比较排序,总共有 种排序,因此

两边取对数,(公式3.19)

因此,堆排序和归并排序都是渐进最优的比较排序算法。

计数排序

注意这种排序算法只适用于排序 到 的整数数组。要排序负数,还需要加上一个偏移。
给出代码(注意结果数组应能够放下排序后的数组)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
template<typename T> void COUNTING_SORT(const vector<T>&a, vector<T>&res) {
    static_assert(std::is_integral_v<T>, "Must be int... to sort");
    if (a.empty()) { throw std::runtime_error("Cannot sort an empty array"); }
    res.resize(a.size());
    T max = *std::max_element(a.begin(),a.end());
    const size_t k = max >= 0 ? max : throw std::runtime_error("cannot sort max < 0") ;
    vector<size_t> Temp(k+1,0);

    for (T elem : a) {
        elem >= 0 ? Temp[elem] = Temp[elem] + 1 : throw std::runtime_error("cannot sort elem < 0");
    }
    for (size_t i = 1; i <= k; i++) {
        Temp[i] = Temp[i-1] + Temp[i];
    }
    for (size_t i = a.size(); i-- > 0;) {
        res[Temp[a[i]] - 1] = a[i];
        Temp[a[i]] = Temp[a[i]] - 1;
    }
}

我们首先确定输入数据的最大值,遍历输入数组并以元素的值为下标,统计每个元素的个数。然后再遍历计数数组,计算出不大于 的数字的个数,然后以计数数组的值为下标,向结果数组写回。与书上不同的是,我们这里使用 来把下标转化为 的数组。
计数排序总运行时间是 ,其中 ,当 时一般采用这种排序,运行时间是 。与其它排序不同,计数排序是稳定的。

基数排序

按照最低有效位挨个比较每位的大小,用用时稳定的计数排序比较每位大小。

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
template<typename T> void COUNT_SORT_EXP(std::vector<T>&a, const size_t& exp, std::vector<T>&out) {
    const int base = 10;
    out.resize(a.size());
    std::vector<size_t>cnt(base,0);
    for (size_t item : a) {
        int digit = (item / exp) % base;
        cnt[digit]++;
    }
    for (size_t i = 1; i < base; i++) {
        cnt[i] = cnt[i-1] + cnt[i];
    }
    for (size_t i = a.size(); i-- > 0;) {
        size_t x = a[i];
        int digit = (x / exp) % base;
        size_t pos = cnt[digit] - 1;
        cnt[digit] -= 1;
        out[pos] = x;
    }
}

template<typename T> void RADIX_SORT(std::vector<T>&a, const size_t& d) {
    std::vector<T>temp(a.size(),0);
    size_t exp = 1;
    for (int i = 1; i <= d; i++) {
        COUNT_SORT_EXP(a, exp, temp);
        a.swap(temp);
        exp = exp * 10;
    }
}

对于 个 位数,每一位数有 个取值,使用的稳定的排序算法耗时为 ,则总耗时为 。
引理8.4: 给定 个 位数和任何正整数 ,如果稳定的排序算法耗时为 ,那么总耗时为 。
就是将 位数分为 个数,其中每个不超过 。

桶排序

桶排序仅限于对输入数据在 范围内。桶排序将 放在 的列表中,并用插入排序对每个 进行排序,最后重新拼成排好序的数组。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
template<typename T> void BUCKET_SORT(std::vector<T>&a) {
    static_assert(std::is_floating_point_v<T>);
    size_t n = a.size();
    std::vector<std::vector<T>>b(n, std::vector<T>{});
    for (size_t i = 0; i < n; i++) {
        __int128_t index = n*a[i] >= n ? n-1 : n*a[i];
        b[index].push_back(a[i]);
    }
    for (size_t i = 0; i < n; i++) {
        InsertionSort(b[i]);
    }
    a.clear();
    for (size_t i = 0; i < n; i++) {
        for (size_t j = 0; j < b[i].size(); j++) {
            a.push_back(b[i][j]);
        }
    }
}

桶排序的时间代价为

先算平均时间,对两边取期望,有

断言

我们定义对于所有的 ,有

因此

展开,有

其中

于是

当 时, 是独立的,于是

代回公式,得

所以桶排序的期望运行时间为

中位数和顺序统计量

在由 个元素组成的集合中,第 个顺序统计量(order statistic) 是该集合中第 小的元素。当 为奇数的时候,中位数在 处,当 是偶数时,中位数在 和 处。不考虑 的奇偶性时,中位数是 (下中位数)和 (上中位数)。我们讨论的中位数全是下中位数。
选择问题:

  • 输入: 一个包含 个互异的数的集合 和一个整数 。
  • 输出: 元素 ,且 中恰好有 个其他元素小于它。

最大值和最小值

最优的算法要进行 次比较找到最小值

1
2
3
4
5
6
7
8
9
template<typename T> T MINIMUM(const std::vector<T>&a) {
    T min = a[0];
    for (size_t i = 2; i < a.size(); i++) {
        if (min > a[i]) {
            min = a[i];
        }
    }
    return min;
}

为了同时找到最大值,我们只需要比较 次即可,方法是每次从数组读入两个,相互比较再和已知的最大值和最小值比较。

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
template<typename T> std::vector<T> MAX_AND_MINIMUM(const std::vector<T>&a) {
    T min = a[0];
    T max = a[1];
    for (size_t i = 1; i < a.size(); i += 2) {
        if (a[i] > a[i-1]) {
            if (a[i] > max) {
                max = a[i];
            }
            if (a[i-1] < min) {
                min = a[i-1];
            }
        }
        else {
            if (a[i] < min) {
                min = a[i];
            }
            if (a[i-1] > max) {
                max = a[i-1];
            }
        }
    }
    std::vector<T>res;
    res.push_back(min);
    res.push_back(max);
    return res;
}

期望为线性时间的选择算法

我们借用一下快速排序里的 函数,还记得快速排序的划分函数吗,我们随机选择一个主元,然后将数组划分为小于等于主元和大于主元的两个子数组。

1
2
3
4
5
6
7
8
9
10
11
12
template<typename T> T RAND_SELECT(std::vector<T>&a, int p, int r, int i) {
    if (p == r) return a[p];
    size_t q = RAND_PARTITION(a, p, r);
    size_t k = q - p + 1;
    if (i == k) {
        return a[q];
    }
    else if (i < k) {
        return RAND_SELECT(a, p, q-1, i);
    }
    else return RAND_SELECT (a, q+1, r, i-k);
}

这个算法的期望运行时间是 。
因为是随机选择主元,所以对于每一个 ,子数组 有 个元素的概率是 。我们定义指示器随机变量为

因为元素互异,所以

因为调用这个函数时,并不知道主元 是不是正确答案,算法要在 或 上递归,取决于第 小的元素相对与 的位置。为了得到上界,我们总是假设该元素落在 的右边。当 时,我们处理的两个子数组的规模为 和 ,有递归式

两边取期望

下面来考虑 ,我们有

若 是偶数,那么从 到 的每一项都会在总和中出现两次,若 是奇数,那么除了 只出现一次外,其他的项也会出现两次,有

我们将使用替代法来求得 。设满足这个递归式初始条件的某个常数 ,有 。假设对小于某个常数的 ,有 。还要选出最后一项 。这个归纳假设可以得到

为了完成证明,我们还需要证明对于足够大的 ,这个式子最多是 ,即 ,化简得 。只要我们的常数 满足 ,即 ,就能两边同除 ,得

因此若假设对所有 ,都有 ,那么就有 。

最坏情况为线性时间的选择算法

描述:

  1. 将输入的数组分为 组,每组个元素,最后一组有 个元素
  2. 然后,用插入排序寻找它们的中位数
  3. 再递归调用 找出这 个中位数的中位数 ,这里我们寻找的是下中位数。
  4. 用快速排序的 函数的修改版,以 为主元,划分输入数组,设划分后有 个数不大于 ,有 个数在高区。
  5. 如果 ,那么 就是结果,如果 ,那么在高区递归查找第 小的元素,如果 ,那么在低区递归查找第 小的元素。
    给出代码:
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
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
template<typename T>void insertion_sort(std::vector<T>&arr, size_t l, size_t r) {
    size_t len = r - l + 1;
    if (len <= 1) return;
    for (size_t j = l + 1; j < r + 1; j++) {
        T key = arr[j];
        size_t i = j - 1;
        while(i >= l && arr[i] > key ) {
            arr[i+1] = arr[i]; i--;
        }
        arr[i+1] = key;
    }
}
// 三路划分:< pivot, == pivot, > pivot
// 返回等于 pivot 的区间 [eqL, eqR](闭区间)
template <typename T> std::pair<size_t, size_t> partition3(std::vector<T>& a, size_t l, size_t r, const T& pivot) {
    size_t lt = l;      // a[l..lt-1] < pivot
    size_t i  = l;      // current
    size_t gt = r;      // a[gt+1..r] > pivot

    while (i <= gt) {
        if (a[i] < pivot) {
            std::swap(a[lt], a[i]);
            lt++;
            i++;
        } else if (a[i] > pivot) {
            std::swap(a[i], a[gt]);
            gt--;
        } else {
            ++i;
        }
    }
    // [l, lt-1] < pivot, [lt, gt] == pivot, [gt+1, r] > pivot
    return {lt, gt};
}
// BFPRT: 在 a[l..r] 中找第 k 小
template <typename T> T bfprt_select(std::vector<T>& a, size_t l, size_t r, size_t k) {
    size_t n = r - l + 1;
    // 小规模直接排序返回
    if (n <= 5) {
        insertion_sort(a, l, r);
        return a[l + k];
    }
    // 1) 分组5个:排序每组,把每组中位数搬到前面
    size_t m = 0; // medians count
    for (size_t start = l; start <= r; start += 5) {
        size_t end = std::min(start + 4, r);
        insertion_sort(a, start, end);
        // 每组中位数索引(下中位数)
        size_t len = end - start + 1;
        size_t med_idx = start + (len - 1) / 2;
        // 把中位数放到 a[l + m]
        std::swap(a[l + m], a[med_idx]);
        ++m;
    }
    // 2) 递归选“中位数的中位数”作为 pivot(m 个中位数里找第 m/2 小)
    T pivot = bfprt_select(a, l, l + m - 1, m / 2);
    // 3) 围绕 pivot 三路划分
    auto [eqL, eqR] = partition3(a, l, r, pivot);
    // 4) 决定往哪边递归(只递归一边)
    size_t leftSize = (eqL > l) ? (eqL - l) : 0;
    size_t midSize  = eqR - eqL + 1;
    if (k < leftSize) {
        return bfprt_select(a, l, eqL - 1, k);
    } else if (k < leftSize + midSize) {
        return pivot; // 在 ==pivot 区间内
    } else {
        return bfprt_select(a, eqR + 1, r, k - leftSize - midSize);
    }
}
template <typename T> T SELECT(std::vector<T>& a, size_t k) {
    return bfprt_select(a, 0, a.size() - 1, k - 1); // 把第k小转换成1-based,这里要减1
}

先确定大于主元划分 的元素个数的下界,在第二步找出的中位数中,至少有一半大于等于 ,因此除了 不能被整除的时候的最后一组和包含 的那一组外,至少有一半的组中,有 个元素大于 。不算这两个组,大于 的元素的个数为

同理,至少有 个元素小于 ,因此在第五步, 最多作用 个元素。
很明显步骤 需要 的时间,步骤三需要 的时间,步骤五至多需要 的时间。我们假设 是单调递增的,而且假设对于任何少于 的规模需要 的时间。可以得到下面的递归式:

我们再次使用替换法,要选择某个常数 和所有的 ,有 。
首先对于 ,这个假设显然成立。
而对于 的情况,我们还要选择常数 使得 的上界是 。代入递归式,得

上式最多是 ,下式成立:

当 时,等价于 ,我们假设了 ,所以任何满足 就能够满足 。注意,我们可以任意的选择大于的数字,并没有特殊之处,只需要根据假设选择合适的就可以。
由此可见, 最坏情况下的运行时间是。
本章讲的两个选择算法,也是通过比较元素来获得它们之间的相对次序的,然而我们并没有对它们进行排序,所以它们不受 的制约。