文章目录

1 前言
2 矩阵乘法
- 2.1 朴素解法
  - 2.1.1 分解
  - 2.1.2 解决
- 2.2 Strassen算法
  - 2.2.1 分解
  - 2.2.2 解决
3 分治法分析
4 后记

前言

这一章将介绍重要的算法设计策略——分治法（divide and conquer）。本章基于《Introduction to Algorithms》4th 第四章，包含分治法解决问题（以矩阵乘法为例）和分析分治法算法的运行时间——代入法、递归树法以及主方法。

矩阵乘法

这一章举例是以N x N的两个矩阵A, B相乘，得到N x N矩阵C。

朴素解法

在线性代数中我们学到，两个矩阵相乘，C的第i行第j列的数C_{ij =}A的第i行 x B的第j列（逐位相乘后逐位相加）。

抽象成数学语言，设C_ijk = A的第i行第k个数 x B的第j列第k个数，于是有：

我们用伪代码很容易实现：

MATRIX-MULTIPLY(A, B, C, n)
    for i = 1 to n
        for j = 1 to n
            for k = 1 to n
                c_ij = c_ij + a_ik * b_kj

很明显，它有三层循环，时间复杂度为Θ(n³)。

我们可以考虑使用分治法改造它，使用分治法我们要做到：

分解(Divide)：分治法中的分，将原问题分解为若干子问题。
解决(Conquer)：分治法中的治，求出子问题的解。
合并(Combine)：将子问题的解合并为原问题的解。

分解

我们可以把矩阵分为4个n/2 X n/2的子矩阵（n是2的幂），分别为A_11, A_12, A_21, A_22和B_11, B_12, B_21, B_22，同时也会生成结果为4个n/2 X n/2的子矩阵：C_11, C_12, C_21, C_22

此时我们的公式就可以从C = A X B改写为如下形式：

接下来我们要做的就是继续分割每个矩阵到最小的部分。每次都分成4个部分，直到无法再分。

解决

我们已经分割成子问题去解决了，接下来我们只需要去解决子矩阵之间的乘法即可。

MATRIX-MULTIPLY-RECURSIVE(A, B, C, n)
    if n == 1
        c_11 = c_11 + a_11 · b_11
        return
    partition A, B, C into n/2 × n/2 submatrices A_11, A_12, A_21, A_22; B_11, B_12, B_21, B_22; C_11, C_12, C_21, C_22;
    MATRIX-MULTIPLY-RECURSIVE(A_11, B_11, C_11, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_11, B_12, C_12, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_21, B_11, C_21, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_21, B_12, C_22, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_12, B_21, C_11, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_12, B_22, C_12, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_22, B_21, C_21, n/2)
    MATRIX-MULTIPLY-RECURSIVE(A_22, B_22, C_22, n/2)

如之前图，我们执行了8次乘法运算得到了C矩阵，运用了8个递归表达式。设运行时间为T(n)，故8次递归的运行时间为8T(n/2)，再加上其他语句，我们可以得到：T(n) = 8T(n/2) + Θ(1)

其时间复杂度为Θ(n^log8) = Θ(n³)但是这样的时间复杂度依旧很高，我们还可以进一步优化。

Strassen算法

1969年，Volker Strassen发明了算法，将矩阵乘法的时间复杂度降低到了Θ(n^log7)，算法在运用分治法解决问题过程中，减少了递归树的茂密程度。Strassen算法将递归树的 8 个分枝减少为 7 个分枝，即只进行了7次递归即可解决。

从本质上来看，Strassen算法把原本计算时间长的8的乘法，转成了较多加减法的乘法，减少乘法的数量。这个分解思想来自于《On the Complexity of Matrix Multiplication》，这里只是介绍其过程，只需要知道这个算法将矩阵乘法的时间复杂度降低到了Θ(n^log7)。

分解

和之前一样，我们先把矩阵分为4个n/2 X n/2的子矩阵，得到：

接下来，创建10个 $n / 2 \times n / 2$ 的矩阵 $S_{1}, S_{2}, \dots, S_{10}$ 保存上一步的子矩阵的和或差。定义为：

这 10 个矩阵分为两组，每组 5 个，都是 $A$ 的子矩阵或者 $B$ 的子矩阵内部之间相互加减。这个过程涉及矩阵加减法，时间复杂度为Θ(n²)。

解决

接下来我们递归计算7个矩阵乘积 $P_{1}, P_{2}, \dots, P_{7}$

最后可以直接计算出 $C_{11}, C_{12}, C_{21}, C_{22}$

Strassen算法用了 7 次矩阵乘法 18 次矩阵加减法。

无论如何，乘法运算的时间复杂度都高于加法运算的，所以通过增加常数个数的加法运算而减少乘法运算的策略都是非常划算的。

它的伪代码为：

STRASSEN(A, B, n)
    let C be a new n × n matrix
    if n == 1
        c[1, 1] = a[1, 1] · b[1, 1]
        return
    A[1, 1] = A[1 : n / 2][1 : n / 2]
    A[1, 2] = A[1 : n / 2][n / 2 + 1 : n]
    A[2, 1] = A[n / 2 + 1 : n][1 : n / 2]
    A[2, 2] = A[n / 2 + 1 : n][n / 2 + 1 : n]
    B[1, 1] = B[1 : n / 2][1 : n / 2]
    B[1, 2] = B[1 : n / 2][n / 2 + 1 : n]
    B[2, 1] = B[n / 2 + 1 : n][1 : n / 2]
    B[2, 2] = B[n / 2 + 1 : n][n / 2 + 1 : n]
    S[1] = B[1, 2] - B[2, 2]
    S[2] = A[1, 1] + A[1, 2]
    S[3] = A[2, 1] + A[2, 2]
    S[4] = B[2, 1] - B[1, 1]
    S[5] = A[1, 1] + A[2, 2]
    S[6] = B[1, 1] + B[2, 2]
    S[7] = A[1, 2] - A[2, 2]
    S[8] = B[2, 1] + B[2, 2]
    S[9] = A[1, 1] - A[2, 1]
    S[10] = B[1, 1] + B[1, 2]
    P[1] = STRASSEN(A[1, 1], S[1], n / 2)
    P[2] = STRASSEN(S[2], B[2, 2], n / 2)
    P[3] = STRASSEN(S[3], B[1, 1], n / 2)
    P[4] = STRASSEN(A[2, 2], S[4], n / 2)
    P[5] = STRASSEN(S[5], S[6], n / 2)
    P[6] = STRASSEN(S[7], S[8], n / 2)
    P[7] = STRASSEN(S[9], S[10], n / 2)
    C[1 : n / 2][1 : n / 2] = P[5] + P[4] - P[2] + P[6]
    C[1 : n / 2][n / 2 + 1 : n] = P[1] + P[2]
    C[n / 2 + 1 : n][1 : n / 2] = P[3] + P[4]
    C[n / 2 + 1 : n][n / 2 + 1 : n] = P[5] + P[1] - P[3] - P[7]
    return C

分治法分析

求解分治法的运行时间没有这么直观，我们通过递归知道了分治法的时间表达式后，可以通过代入法、递归树法以及主方法来求解分治法的运行时间。

我们在分析中，存在递归式，一个递归式就是一个等式（方程式），它通过通常更小的函数参数值来描述一个函数。即是带有递归时间表达式。

代入法

无渐进符号递归式

代入法(substitution method)分两步：

使用常数符号猜测解的形式。
使用数学归纳法求出解中常数，并证明解是正确的。

一般来说，我们的猜测有：O(nlgn)、O(n^lgn)等等。接下来通过一个例子来解释代入法的过程。

例如eg.1：

我们猜测T(n) = nlgn + n ①，因为2T(n/2)常常是nlogn的复杂度。进行数学归纳法的证明：

当n=1时，T(n) = 1lg1 +1 = 1，符合T(n)特征。
当n>1时，当n=n/2时，T(n/2) = (n/2) lg (n/2) + n/2 ②
把②带入原式T(n) = 2T(n/2) + n中，得T(n) = n lg(n/2)+ 2n。继续化简，得T(n) = n (lgn-lg2)+ 2n => T(n) = n (lgn-1)+ 2n => T(n) = nlgn + n 得到式子①，则证明成功。

此时T(n) = nlgn + n，可以使用渐进符号表示为T(n) = Θ(nlgn).

带有渐进符号的递归式

这种是递归式中不包含渐进符号时我们做的，当然更加常见的是带有渐进符号的递归式，我们一般需要进行下面步骤：

用常数组合代替渐进符号
分别分析上下界，往往使用不同常数组合代替。

例如eg.2：T(n) = 2T(n/2) + Θ(n).

首先我们如果想分析它的上界，即T(n) = 2T(n/2) + O(n)。我们用cn代替O(n)，并且写成不等式的形式：T(n) ≤ 2T(n/2) + cn ，c为正常数。

根据之前求到的，我们可以猜T(n) ≤ dn lgn，d为其他正常数。和之前一样，对它进行数学归纳法证明：

当n=1时，T(n) ≤ dlg1= 0，是个常数，符合T(n)特征。
当n>1时，当n=n/2时，T(n/2) ≤ (dn/2) lg(n/2) ③
=> T(n) ≤ 2T(n/2) + cn
= dn lg(n/2) + cn
= dn (lgn-1)+ cn
= dn lgn - dn + cn
≤ dn lgn （当且仅当 -dn + cn ≤ 0 => c ≤ d时成立，常数无所谓）
此时可以得到T(n) ≤ dn lgn，写作T(n) = O(nlgn)

同理我们还可以证明其下界，用T(n) ≥ 2T(n/2) + cn表示，猜测T(n) ≥ dn lg n，得到T(n) = Ω(nlgn)。

因此，我们证明了T(n) = Θ(nlgn)。

常数递归式减去低阶项

当然，我们还可能遇到特殊的渐进符号，例如Θ(1)。此时为了配平，右边减去低阶项，一阶项的低阶项为常数项。

例如eg.3：T(n) = 2T(n/2) + Θ(1).

我们猜它：T(n) = O(n)，则此时写成T(n) ≤ cn – d ，而不是T(n) ≤ cn。因为在验证中，若不等式写为T(n) ≤ cn，则会出现下面过程：

T(n) ≤ 2(c(n/2)) + Θ(1)
= cn + Θ(1)
≤ cn （X）

此时这个不等式不一定成立，因为Θ(1)不一定<0，无法得到一个常数范围使其成立，证明不成立。我们应该使用T(n) ≤ cn – d，作为配平。

T(n) ≤ cn – d, where d ≥ 0 is a constant:
T(n) ≤ 2(c(n/2) – d) + Θ(1)
= cn – 2d + Θ(1)
≤ cn – d – (d – Θ(1))
≤ cn – d

此时我们只需要d – Θ(1)>0 => d>Θ(1)即可成立。此时有一个范围使其成立，证明成立。

常见错误

我们证明的时候可能会出现一些常见的问题，例如：T(n) = 2T(n/2) + Θ(n)，猜测T(n) = O(n)

错误写法，直接带入渐进符号：T(n) ≤ 2O(n/2) + Θ(n) = 2O(n) + Θ(n) = O(n)
因为此时O(n)所隐藏的系数不断变化，也就是表达式不一定成立。
错误写法，不替换渐进符号：T(n) ≤ 2c(n/2) + Θ(n) = cn + Θ(n) = O(n)
我们无法确定Θ(n)所隐藏的系数一定大于 c，表达式不一定成立。因此我们需要使用常数组合表达渐进符号。

递归树法

在递归树中，每个结点表示一个单一子问题的代价（需要时间），子问题对应某次递归函数的调用。我们将树中每层中的代价求和，得到每层代价，然后将所有层的代价求和，得到所有层次递归调用的总代价。

例如分析：T(n) = 3T(n/4) + Θ(n²). 此时Θ(n²)是个二次匿名函数，我们可以用cn²表示，此时c为常数且c>0。故我们可以分析：T(n) = 3T(n/4) + cn²。因为本质上递归中执行的为若干个cn²（递归的终点是执行最后一次cn²），故我们通过分析cn²得到递归树。

因为我们递归有3个，故从T(n)处可以分成三个分支：

T(n/4)又可以通过T(n) = 3T(n/4) + cn²得到，T(n/4) = 3T(n/16) + cn². 故每个T(n/4)又可以有三个分支：

同理我们可以继续分T(n/16)，直到最后n分到1为止，此时递归树为：

我们发现到了最后的时候，全都是Θ(1)的运行时间。我们分析每一层，第0层运行时间为cn²，第1层加起来为(3/16)cn²，第2层则为(3/16)²cn²。此时它们与层数呈等比数列，第i层运行时间为(3/16)ⁱcn²④。每个子结点规模为父母结点的1/4，最深处我们定义为d，此时最深处子节点规模为1，n/(4)^d = 1 => d = log₄ⁿ。在树中，高度与深度的数值一样，故高度为h = log₄ⁿ⑤。

我们知道了高度，就可以求高度为h_i处的结点的位置，那一层的结点数为3^h_i。例如第0层，结点数为3⁰ = 1，第一层为3¹ = 3，第二层为3² = 9. 在最后一层，结点数为3^log₄n 。

此时根据数学知识得，我们有一个指数公式：a^log_bc = c^log_ba。

所以可以得到3^log₄n = n^{log₄ 3}。最后一层全都是Θ(1)常数，故最后一层运行时间为：Θ(1) * n^{log₄ 3}= Θ(n^{log₄ 3}) ⑥。

我们把所有层的运行时间加起来，除了最后一层外，其他层的和是一个等比数列求和。

此时我们可以求它的上下界。例如求其上界，继续对它化简得：

这里运用了放缩，使用一个无穷范围以消除指数。这样我们就求到了它的上界。

当然有些时候我们可能遇到树叶结点参差不齐的情况，例如：T(n) = T(n/3) + T(2n/3) + Θ(n).

这棵递归树是不平衡的，很明显一个父结点的两个子结点的代价不一样，左子结点的数据规模是父结点的 1/3，右子结点的数据规模是父结点的 2/3，从根结点到不同的叶结点的路径长度也不一样。

同样的，我们改为T(n) = T(n/3) + T(2n/3) + cn进行分析。

我们发现，从根节点开始的路径在2/3的分支上更长，因为我们假设均为1/3分，或均为2/3分，得到节点数为2^log₃n和2^log_3/2n，此时根据数学知识2^log_3/2n更加大，故此时2^log_3/2n的分支更长。

注：(2/3)ⁱn=1 => i = log_3/2 n，完全二叉树的情况下，节点数为2ⁱ

故我们得到这样的递归树：

此时最深的深度为log_3/2n。因为有叶子节点的存在，导致每一层的代价和≤cn（没到叶子节点的每层代价和为cn，到了叶子节点的小于cn）。此时我们计算总的代价和，可以用O(cnlog_3/2n) = O(nlgn)表示（因为假设最多的情况为完全满二叉树，高度为log_3/2 n，实际上达不到，作为一个粗略的上界）。

因为是一个粗略的上界，我们需要使用代入进行验证。即证明：T(n)≤dn lgn，d>0是个常数。

T(n) ≤ T(n/3) + T(2n/3) + Θ(n)

≤ d(n/3) lg(n/3) + d(2n/3) lg(2n/3) +cn

= d(n/3)(lgn-lg3) + d(2n/3)(lgn-lg(3/2)) +cn

= dn lgn - (d(n/3)lg3 + d(2n/3)lg(3/2)) + cn

= dn lgn - dn((2lg3)/3 + (2lg3/2)/3) +cn

≤ dn lgn 当且仅当 c/(lg3-(2/3)) ≤ d

此时得到的一个粗略上界即可，无需对递归树的每层代价进行精确计算。

同理我们还可以对它进行下界的估计，发现T(n)≥dn lgn 在d>0时也成立，故T(n) = Θ(nlogn)。

Master method

主方法主要来自于递归树法的总结，适用于形如T(n) = aT(n/b) + f(n)的递归式。此时a ≥ 1, b > 1，把f(n)称为驱动函数（driving function），它是一个渐进非负函数。

同时，我们的递归式还可以推导为：T(n) = a'T(⌊n/b⌋) + a''T(⌈n/b⌉) + f(n)，a', a'' ≥ 0 and a' + a'' = a.

则分为下面三种情况：

如果存在ε>0，f(n) = O(n^{log_b a-ε})，则T(n) = Θ(n^{log_b a})；
如果存在k≥0，f(n) = Θ(n^{log_b a} lg^k n)，则T(n) = Θ(n^{log_b a} lg^k+1 n)；
如果存在ε>0，f(n) = Ω(n^{log_b a+ε})，且对于某个常数c<1和所有足够大的n有af(n/b) ≤ cf(n)则T(n) = Θ(f(n)).

主方法的含义是：f(n) 与 n^log_ba 两个函数中较大的那个决定递归式的解。其中f(n)通常是使用常数c的组合表示，在递归树法中分析代价所用，而n^log_ba常常指递归树中最后一层节点的代价（Θ(n^log_ba)）。

当f(n) < n^log_ba时，递归式的解T(n) = Θ(n^{log_b a})；这种情况下代价主要来源于叶子节点。
当f(n) ≈ n^log_ba时，递归式的解乘上一个对数因子，T(n) = Θ(n^{log_b a} lg^k+1 n) = Θ(f(n) lg^k+1 n)；这种情况下每一层的代价基本相等为n^{log_b a} lg^k ，而其高度一般为Θ(lgn)。
当f(n) > n^log_ba时，递归式的解T(n) = Θ(f(n))。这种情况下代价主要来源于根节点。

此时的大小比较时多项式的渐进大小比较，而不是常数项与低次项的比较。实际上比的是两个函数渐近界的大小，即哪个函数增长得快，它们至少要相差nε个因子。

下面分情况分析如何使用主方法。

eg.1：T(n) = 2T(n/4) + 1.

分析可知：f(n) = 1， n^log_ba= n^{log4 2} = n^1/2。很明显，f(n) < n^log_ba，可以使用主方法中的情况1，得到T(n) = Θ(n^1/2)

eg.2：T(n) = 2T(n/4) + √nlg²n

分析可知：f(n) = √nlg²9， n^log_ba= n^{log4 2} = n^1/2。此时在增长率上，f(n) = Θ(n^1/2lg²n) = Θ(n^log_balg²n)，它们增长率相当，适用于情况2，得到T(n) = Θ(n^1/2lg³n)。

eg.3：T(n) = 2T(n/4) + n².

分析可知：f(n) = n²， n^log_ba= n^{log4 2} = n^1/2。很明显，f(n) > n^log_ba，且此时满足2f(n/4) = (1/8)n²≤ cf(n)，此时0<c<1。可以使用主方法中的情况3，得到T(n) = Θ(n²)。