Introduction to Algorithms - Third Edition
Part I. Foundations
Chapter 4. Divide-and-Conquer
当子问题足够大以进行递归求解时,称其为递归情况(recursive case)。一旦子问题变得很小,以至不再递归的程度,就说递归“结束(bottom out)”了,已经回到了基本情况(base case)。
递归式
递归式与分治法范式提供了一种自然的方式来表征分治算法的运行时间。
递归式(recurrence)是一组等式或不等式,其描述的函数是用在更小的输入下该函数的值来定义的。
本章介绍三种解递归式的方法,即找出解的渐近
Θ
\Theta
Θ 或
O
O
O 界的方法。
代换法(substitution method):先猜测一个界,再使用数学归纳法证明猜测的正确性。
递归树方法(recursion-tree method):将递归式转化为树,树中的结点代表不同递归层次付出的代价。使用边界求和的技巧来解递归。
主方法(master method):给出递归形式
T
(
n
)
=
a
T
(
n
/
b
)
+
f
(
n
)
T(n) = aT(n/b) + f(n)
T(n)=aT(n/b)+f(n) (4.2)
的界,其中
a
≥
1
a\ge 1
a≥1,
b
>
1
b>1
b>1,
f
(
n
)
f(n)
f(n) 是给定的函数。
有时,递归式是不等式。例, T ( n ) ≤ 2 T ( n ) + Θ ( n ) T(n) \le 2T(n) + \Theta(n) T(n)≤2T(n)+Θ(n),此时,递归式描述了 T ( n ) T(n) T(n) 的上界,其解使用 O O O 形式。相似地, T ( n ) ≥ 2 T ( n ) + Θ ( n ) T(n) \ge 2T(n) + \Theta(n) T(n)≥2T(n)+Θ(n),此递归式给出 T ( n ) T(n) T(n) 的下界,其解使用 Ω \Omega Ω 形式。
递归式中的技术细节
在描述并解递归式时,常忽略上取整、下取整和边界条件。进行分析时先假设没有这些细节,而后再确定它们重要与否。但我们需要知道它们在什么情况下是重要的。经验和一些定理说明:这些细节不会影响许多递归式(表征分治算法)的渐近界(见定理 4.1)。但是,在本章中,我们将介绍其中一些细节,并说明递归解决方法的要点。
4.1 最大子数组问题(The maximum-subarray problem)
问题:允许购买一次股票,然后在以后的某个日期将其出售,目标是使利润最大化。
下表显示了 17 天股票的价格。
Day | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 | 14 | 15 | 16 |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Price | 100 | 113 | 110 | 85 | 105 | 102 | 86 | 63 | 81 | 101 | 94 | 106 | 101 | 79 | 94 | 90 | 97 |
Change | 13 | -3 | -25 | 20 | -3 | -16 | -23 | 18 | 20 | -7 | 12 | -5 | -12 | -5 | -4 | 7 |
图 4.1 表格的最后一行给出了较前一日的价格变化。
暴力解决方案(A brute-force solution)
尝试所有可能的购买和出售日期对,购买日期要早于出售日期。 n n n 天的时间段有 ( n 2 ) \binom{n}{2} (2n) 个这样的日期对。 ( n 2 ) \binom{n}{2} (2n) 是 Θ ( n 2 ) \Theta(n^2) Θ(n2),这种方法将花费 Ω ( n 2 ) \Omega(n^2) Ω(n2) 时间。
转换
将图 4.1 中的表格的最后一行视为数组
A
A
A,找到
A
A
A 的非空连续子数组,使其值之和最大。 我们将此连续子数组称为最大子数组(maximum subarray)。
股票价格的变化是最大子数组的问题。子数组
A
[
8..11
]
A[8..11]
A[8..11] 的和为 43,在数组
A
A
A 的任何连续子数组中具有最大的和。因此,在第 7 日购买股票,在第 11 日售出,可以赚取最大利益为每股 43 美元。
分而治之解决方案
程序 FIND-MAX-CROSSING-SUBARRAY 找出子数组 A [ l o w . . h i g h ] A[low..high] A[low..high] 中跨越中点(称为 m i d mid mid)的最大子数组,运行时间可表示为 Θ ( n ) \Theta(n) Θ(n),其中 n n n 是子数组 A [ l o w . . h i g h ] A[low..high] A[low..high] 的元素个数。
FIND-MAX-CROSSING-SUBARRAY ( A , l o w , m i d , h i g h ) (A,low,mid,high) (A,low,mid,high)
1 l e f t _ s u m = − ∞ left\_sum = -\infty left_sum=−∞
2 s u m = 0 sum=0 sum=0
3 for i = m i d i=mid i=mid downto l o w low low
4 s u m = s u m + A [ i ] sum = sum + A[i] sum=sum+A[i]
5 if s u m > l e f t _ s u m sum > left\_sum sum>left_sum
6 l e f t _ s u m = s u m left\_sum = sum left_sum=sum
7 m a x _ l e f t = i max\_left = i max_left=i
8 r i g h t _ s u m = − ∞ right\_sum = -\infty right_sum=−∞
9 s u m = 0 sum = 0 sum=0
10 for j = m i d + 1 j=mid+1 j=mid+1 to h i g h high high
11 s u m = s u m + A [ j ] sum = sum + A[j] sum=sum+A[j]
12 if s u m > r i g h t _ s u m sum > right\_sum sum>right_sum
13 r i g h t _ s u m = s u m right\_sum = sum right_sum=sum
14 m a x _ r i g h t = j max\_right = j max_right=j
15 return ( m a x _ l e f t , m a x _ r i g h t , l e f t _ s u m + r i g h t _ s u m ) (max\_left, max\_right, left\_sum+right\_sum) (max_left,max_right,left_sum+right_sum)
程序 FIND-MAX-CROSSING-SUBARRAY 通过分而治之策略解决最大子数组问题。
FIND-MAXIMUM-SUBARRAY ( A , l o w , h i g h ) (A, low, high) (A,low,high)
1 if h i g h = = l o w high == low high==low
2 return ( l o w , h i g h , A [ l o w ] ) (low, high, A[low]) (low,high,A[low]) // 基本情况:只包含一个元素
3 else
4 m i d = ⌊ l o w + h i g h ⌋ / 2 mid = \lfloor low+high\rfloor / 2 mid=⌊low+high⌋/2
5 ( l e f t _ l o w , l e f t _ h i g h , l e f t _ s u m ) = (left\_low,left\_high,left\_sum) = (left_low,left_high,left_sum)= FIND-MAXIMUM-SUBARRAY ( A , l o w , m i d ) (A, low, mid) (A,low,mid)
6 ( r i g h t _ l o w , r i g h t _ h i g h , r i g h t _ s u m ) = (right\_low,right\_high,right\_sum) = (right_low,right_high,right_sum)= FIND-MAXIMUM-SUBARRAY ( A , m i d + 1 , h i g h ) (A, mid+1, high) (A,mid+1,high)
7 ( c r o s s _ l o w , c r o s s _ h i g h , c r o s s _ s u m ) = (cross\_low,cross\_high,cross\_sum) = (cross_low,cross_high,cross_sum)= FIND-MAX-CROSSING-SUBARRAY ( A , l o w , m i d , h i g h ) (A,low,mid,high) (A,low,mid,high)
8 if l e f t _ s u m ≥ r i g h t _ s u m left\_sum \ge right\_sum left_sum≥right_sum and l e f t _ s u m ≥ c r o s s _ s u m left\_sum \ge cross\_sum left_sum≥cross_sum
9 return ( l e f t _ l o w , l e f t _ h i g h , l e f t _ s u m ) (left\_low,left\_high,left\_sum) (left_low,left_high,left_sum)
10 elseif r i g h t _ s u m ≥ l e f t _ s u m right\_sum \ge left\_sum right_sum≥left_sum and r i g h t _ s u m ≥ c r o s s _ s u m right\_sum \ge cross\_sum right_sum≥cross_sum
11 return ( r i g h t _ l o w , r i g h t _ h i g h , r i g h t _ s u m ) (right\_low,right\_high,right\_sum) (right_low,right_high,right_sum)
12 else
13 return ( c r o s s _ l o w , c r o s s _ h i g h , c r o s s _ s u m ) (cross\_low,cross\_high,cross\_sum) (cross_low,cross_high,cross_sum)
Python 代码实现如下:
def FindMaxCrossingSubarray(A, low, mid, high):
left_sum = float('-inf')
ssum = 0
for i in range(mid, low-1, -1):
ssum += A[i]
if ssum > left_sum:
left_sum = ssum
max_left = i
right_sum = float('-inf')
ssum = 0
for j in range(mid+1, high+1):
ssum += A[j]
if ssum > right_sum:
right_sum = ssum
max_right = j
return (max_left, max_right, left_sum + right_sum)
def FindMaximumSubarray(A, low, high):
if high == low:
return (low, high, A[low])
else:
mid = (low + high) / 2
(left_low, left_high, left_sum) = FindMaximumSubarray(A, low, mid)
(right_low, right_high, right_sum) = FindMaximumSubarray(A, mid+1, high)
(cross_low, cross_high, cross_sum) = FindMaxCrossingSubarray(A, low, mid, high)
if left_sum >= right_sum and left_sum >= cross_sum:
return (left_low, left_high, left_sum)
elif right_sum >= left_sum and right_sum >= cross_sum:
return (right_low, right_high, right_sum)
else:
return (cross_low, cross_high, cross_sum)
if __name__ == "__main__":
A = [13,-3,-25,20,-3,-16,-23,18,20,-7,12,-5,-12,-5,-4,7]
print FindMaximumSubarray(A, 0, len(A)-1) # 打印 (7, 10, 43)
分析分治算法
在元素个数为
n
n
n 的子数组中,FIND-MAXIMUM-SUBARRAY 的运行时间
T
(
n
)
T(n)
T(n) 可以用递归式表示为:
T
(
n
)
=
{
Θ
(
1
)
if
n
=
1
2
T
(
n
/
2
)
+
Θ
(
n
)
if
n
>
1
T(n) = \begin{cases} \Theta(1) &\text{if } n=1 \\ 2T(n/2) + \Theta(n) &\text{if } n >1 \end{cases}
T(n)={Θ(1)2T(n/2)+Θ(n)if n=1if n>1 (4.7)
这个递归式的解为
T
(
n
)
=
Θ
(
n
lg
n
)
T(n) = \Theta(n\lg n)
T(n)=Θ(nlgn)。
练习
4.1-1 当
A
A
A 中的所有元素都是负数时,FIND-MAXIMUM-SUBARRAY 返回什么?
解: 返回
A
A
A 中最大的元素及其下标。
4.1-2 写出使用暴力方法解决最大子数组问题的伪代码。运行时间为
Θ
(
n
2
)
\Theta(n^2)
Θ(n2)。
解: FIND-MAX-SUBARRAY(
A
A
A)
1
m
a
x
_
s
u
m
=
−
∞
max\_sum = -\infty
max_sum=−∞
2 for
i
=
1
i=1
i=1 to
A
.
l
e
n
g
t
h
A.length
A.length
3
s
u
m
=
0
sum = 0
sum=0
4 for
j
=
i
j=i
j=i to
A
.
l
e
n
g
t
h
A.length
A.length
5
s
u
m
=
s
u
m
+
A
[
j
]
sum = sum + A[j]
sum=sum+A[j]
6 if
s
u
m
>
m
a
x
_
s
u
m
sum > max\_sum
sum>max_sum
7
m
a
x
_
s
u
m
=
s
u
m
max\_sum = sum
max_sum=sum
8
l
o
w
=
i
low = i
low=i
9
h
i
g
h
=
j
high = j
high=j
10 return
(
l
o
w
,
h
i
g
h
,
m
a
x
_
s
u
m
)
(low, high, max\_sum)
(low,high,max_sum)
4.1-3 在计算机上分别使用暴力方法和递归算法实现最大子数组问题。什么问题规模
n
0
n_0
n0 给出递归算法优于暴力算法的交叉点?当问题规模小于
n
0
n_0
n0 时,将递归算法的基本情况改为使用暴力算法。这会改变交叉点吗?
解: 令
Θ
(
n
lg
n
)
≥
Θ
(
n
2
)
\Theta(n\lg n) \ge \Theta(n^2)
Θ(nlgn)≥Θ(n2) 的最小正解为
m
0
m_0
m0,则有
n
0
≥
m
0
n_0 \ge m_0
n0≥m0。
当子问题规模小于
n
0
n_0
n0 时,使用暴力算法使得子问题花费时间减少,递归算法的总运行时间减少。会改变交叉点,交叉点会变小。
4.1-4 假设改变最大子数组的定义,允许结果为空子数组,其中空子数组的值的和为 0。如何更改不允许空子数组的任何算法,以允许结果为空子数组?
解: 方法一:将求出的最大子数组的和与 0 相比较,取较大的值。
方法二:将初始的
m
a
x
_
s
u
m
=
−
∞
max\_sum=-\infty
max_sum=−∞ 改为
m
a
x
_
s
u
m
=
0
max\_sum = 0
max_sum=0,
l
e
f
t
_
s
u
m
left\_sum
left_sum 等作同样处理。将标明子数组范围的索引也初始化为 -1(或其他可行的数字)。
4.1-5 针对最大子数组问题,使用下述思想开发一个非递归的,线性时间的算法。从数组的左端开始,向右移动,记录到目前为止的最大子数组。知道
A
[
1..
j
]
A[1..j]
A[1..j] 的一个最大子数组,扩展一下,找到以索引
j
+
1
j+1
j+1 结尾的子数组的最大子数组,可通过以下观察:
A
[
1..
j
+
1
]
A[1..j+1]
A[1..j+1] 的最大子数组要么是
A
[
1..
j
]
A[1..j]
A[1..j] 的最大子数组,要么是子数组
A
[
i
.
.
j
+
1
]
A[i..j+1]
A[i..j+1],其中
1
≤
i
≤
j
+
1
1 \le i \le j+1
1≤i≤j+1。在知道索引
j
j
j 结尾的最大子数组的基础上,在常数时间内,确定
A
[
i
.
.
j
+
1
]
A[i..j+1]
A[i..j+1] 形式的最大子数组。
解: LINEAR-FIND-MAX-SUBARRAY(
A
A
A)
1
m
a
x
_
s
u
m
=
A
[
1
]
max\_sum = A[1]
max_sum=A[1]
2
s
u
m
=
A
[
1
]
sum = A[1]
sum=A[1]
3
l
o
w
=
1
low = 1
low=1
4
h
i
g
h
=
1
high = 1
high=1
5
s
t
a
r
t
=
1
start = 1
start=1
6 for
i
=
2
i=2
i=2 to
A
.
l
e
n
g
t
h
A.length
A.length
7
s
u
m
=
s
u
m
+
A
[
i
]
sum = sum + A[i]
sum=sum+A[i]
8 if
s
u
m
>
m
a
x
_
s
u
m
sum > max\_sum
sum>max_sum
9
m
a
x
_
s
u
m
=
s
u
m
max\_sum = sum
max_sum=sum
10
l
o
w
=
s
t
a
r
t
low = start
low=start
11
h
i
g
h
=
i
high = i
high=i
12 if
s
u
m
<
0
sum < 0
sum<0
13
s
u
m
=
0
sum = 0
sum=0
14
s
t
a
r
t
=
i
start = i
start=i
15 return
(
l
o
w
,
h
i
g
h
,
m
a
x
_
s
u
m
)
(low, high, max\_sum)
(low,high,max_sum)
4.2 斯特拉森算法(Strassen’s algorithm)求解矩阵相乘
如果
A
=
(
a
i
j
)
A=(a_{ij})
A=(aij) 和
B
=
(
b
i
j
)
B = (b_{ij})
B=(bij) 是
n
×
n
n \times n
n×n 矩阵,那么在乘积
C
=
A
⋅
B
C = A\cdot B
C=A⋅B 中,定义其元素
c
i
j
c_{ij}
cij 为:
C
i
j
=
∑
k
=
1
n
a
i
k
⋅
b
k
j
C_{ij} = \sum_{k=1}^{n}a_{ik}\cdot b_{kj}
Cij=∑k=1naik⋅bkj, (4.8)
其中,
i
,
j
=
1
,
2
,
…
,
n
i,j = 1,2,\dots,n
i,j=1,2,…,n。
假设每个矩阵有一个属性
r
o
w
s
rows
rows,给出矩阵的行数。
SQUARE-MATRIX-MULTIPLY ( A , B ) (A, B) (A,B)
1 n = A . r o w s n = A.rows n=A.rows
2 let C C C be a new n × n n \times n n×n matrix
3 for i = 1 i=1 i=1 to n n n
4 for j = 1 j=1 j=1 to n n n
5 c i j = 0 c_{ij} = 0 cij=0
6 for k = 1 k = 1 k=1 to n n n
7 c i j = c i j + a i k ⋅ b k j c_{ij} = c_{ij} + a_{ik} \cdot b_{kj} cij=cij+aik⋅bkj
8 return C C C
程序 SQUARE-MATRIX-MULTIPLY 的运行时间为
Θ
(
n
3
)
\Theta(n^3)
Θ(n3)。
针对
n
×
n
n \times n
n×n 矩阵相乘的 Strassen 递归算法的运行时间为
Θ
(
n
lg
7
)
\Theta(n^{\lg 7})
Θ(nlg7),所以其运行时间为
o
(
n
2.81
)
o(n^{2.81})
o(n2.81)。
一个简单的分治算法
假设
n
n
n 是 2 的幂次,把
A
A
A,
B
B
B 和
C
C
C 都分成四个
n
/
2
×
n
/
2
n/2 \times n/2
n/2×n/2 矩阵。则等式
C
=
A
⋅
B
C = A \cdot B
C=A⋅B 可以写成:
(
C
11
C
12
C
21
C
22
)
=
(
A
11
A
12
A
21
A
22
)
⋅
(
B
11
B
12
B
21
B
22
)
\begin{pmatrix} C_{11} & C_{12} \\ C_{21} & C_{22} \end{pmatrix} = \begin{pmatrix} A_{11} & A_{12} \\ A_{21} & A_{22} \end{pmatrix} \cdot \begin{pmatrix} B_{11} & B_{12} \\ B_{21} & B_{22} \end{pmatrix}
(C11C21C12C22)=(A11A21A12A22)⋅(B11B21B12B22) (4.10)
等式(4.10)对应四个等式:
C
11
=
A
11
⋅
B
11
+
A
12
⋅
B
21
C_{11} = A_{11}\cdot B_{11} + A_{12} \cdot B_{21}
C11=A11⋅B11+A12⋅B21, (4.11)
C
12
=
A
11
⋅
B
12
+
A
12
⋅
B
22
C_{12} = A_{11}\cdot B_{12} + A_{12} \cdot B_{22}
C12=A11⋅B12+A12⋅B22, (4.12)
C
21
=
A
21
⋅
B
11
+
A
22
⋅
B
21
C_{21} = A_{21}\cdot B_{11} + A_{22} \cdot B_{21}
C21=A21⋅B11+A22⋅B21, (4.13)
C
22
=
A
21
⋅
B
12
+
A
22
⋅
B
22
C_{22} = A_{21}\cdot B_{12} + A_{22} \cdot B_{22}
C22=A21⋅B12+A22⋅B22。 (4.14)
SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A , B ) (A,B) (A,B)
1 n = A . r o w s n = A.rows n=A.rows
2 let C C C be a new n × n n\times n n×n matrix
3 if n = = 1 n==1 n==1
4 c 11 = a 11 ⋅ b 11 c_{11} = a_{11} \cdot b_{11} c11=a11⋅b11
5 else partition A A A, B B B and C C C into n / 2 × n / 2 n/2 \times n/2 n/2×n/2 matrices
6 C 11 = C_{11} = C11= SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 11 , B 11 ) (A_{11},B_{11}) (A11,B11)
+ + + SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 12 , B 21 ) (A_{12},B_{21}) (A12,B21)
7 C 12 = C_{12} = C12= SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 11 , B 12 ) (A_{11},B_{12}) (A11,B12)
+ + + SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 12 , B 22 ) (A_{12},B_{22}) (A12,B22)
8 C 21 = C_{21} = C21= SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 21 , B 11 ) (A_{21},B_{11}) (A21,B11)
+ + + SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 22 , B 21 ) (A_{22},B_{21}) (A22,B21)
9 C 22 = C_{22} = C22= SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 21 , B 12 ) (A_{21},B_{12}) (A21,B12)
+ + + SQUARE-MATRIX-MULTIPLY-RECURSIVE ( A 22 , B 22 ) (A_{22},B_{22}) (A22,B22)
10 return C C C
在第 5 行如何划分矩阵?① 创建 12 个新的
n
/
2
×
n
/
2
n/2 \times n/2
n/2×n/2 矩阵,需要花费
Θ
(
n
2
)
\Theta(n^2)
Θ(n2) 时间来复制元素;② 使用索引计算来指定子矩阵,第 5 行运行时间则只有
Θ
(
1
)
\Theta(1)
Θ(1)。
第 6 - 9 行矩阵相加所用时间为
Θ
(
n
2
)
\Theta(n^2)
Θ(n2)。
程序 SQUARE-MATRIX-MULTIPLY-RECURSIVE 运行时间的递归式表示为:
T
(
n
)
=
{
Θ
(
1
)
if
n
=
1
,
8
T
(
n
/
2
)
+
Θ
(
n
2
)
if
n
>
1.
T(n) = \begin{cases} \Theta(1) &\text{if } n=1, \\ 8T(n/2) + \Theta(n^2) &\text{if } n >1. \end{cases}
T(n)={Θ(1)8T(n/2)+Θ(n2)if n=1,if n>1. (4.17)
递归式的解为
T
(
n
)
=
Θ
(
n
3
)
T(n) = \Theta(n^3)
T(n)=Θ(n3)。
记住:虽然渐近符号包含常数乘法的因子,但递归式,比如
T
(
n
/
2
)
T(n/2)
T(n/2) 不包含。
Strassen 方法
四个步骤:
- 将输入矩阵 A A A, B B B 和输出矩阵 C C C 分解成 n / 2 × n / 2 n/2 \times n/2 n/2×n/2 子矩阵。使用索引计算花费 Θ ( 1 ) \Theta(1) Θ(1) 时间。
- 创建 10 个 n / 2 × n / 2 n/2 \times n/2 n/2×n/2 矩阵 S 1 , S 2 , … , S 10 S_1,S_2,\dots,S_{10} S1,S2,…,S10,它们是步骤 1 中创建的两个矩阵的和或者差。创建 10 个矩阵花费时间 Θ ( n 2 ) \Theta(n^2) Θ(n2)。
- 使用步骤 1 创建的子矩阵和步骤 2 中的 10 个矩阵,递归计算 7 个矩阵积 P 1 , P 2 , … , P 7 P_1,P_2,\dots,P_7 P1,P2,…,P7。矩阵 P i P_i Pi 是 n / 2 × n / 2 n/2 \times n/2 n/2×n/2 的。
- 通过加减矩阵 P i P_i Pi 的各种组合,计算出结果矩阵 C C C 所需的子矩阵 C 11 C_{11} C11, C 12 C_{12} C12, C 21 C_{21} C21, C 22 C_{22} C22。计算 4 个矩阵花费时间 Θ ( n 2 ) \Theta(n^2) Θ(n2)。
Strassen 算法运行时间
T
(
n
)
T(n)
T(n) 的递归式表示为:
T
(
n
)
=
{
Θ
(
1
)
if
n
=
1
,
7
T
(
n
/
2
)
+
Θ
(
n
2
)
if
n
>
1.
T(n) = \begin{cases} \Theta(1) &\text{if } n=1, \\ 7T(n/2) + \Theta(n^2) &\text{if } n >1. \end{cases}
T(n)={Θ(1)7T(n/2)+Θ(n2)if n=1,if n>1. (4.18)
递归式(4.18)的解为
T
(
n
)
=
Θ
(
n
lg
7
)
T(n) = \Theta(n^{\lg7})
T(n)=Θ(nlg7)。
在步骤 2 中,创建了下面 10 个矩阵:
S
1
=
B
12
−
B
22
S_1 = B_{12} - B_{22}
S1=B12−B22,
S
2
=
A
11
+
A
12
S_2 = A_{11} + A_{12}
S2=A11+A12,
S
3
=
A
21
+
A
22
S_3 = A_{21} + A_{22}
S3=A21+A22,
S
4
=
B
21
−
B
11
S_4 = B_{21} - B_{11}
S4=B21−B11,
S
5
=
A
11
+
A
22
S_5 = A_{11} + A_{22}
S5=A11+A22,
S
6
=
B
11
+
B
22
S_6 = B_{11} + B_{22}
S6=B11+B22,
S
7
=
A
12
−
A
22
S_7 = A_{12} - A_{22}
S7=A12−A22,
S
8
=
B
21
+
B
22
S_8 = B_{21} + B_{22}
S8=B21+B22,
S
9
=
A
11
−
A
21
S_9 = A_{11} - A_{21}
S9=A11−A21,
S
10
=
B
11
+
B
12
S_{10} = B_{11} + B_{12}
S10=B11+B12。
在步骤 3 中,
P
1
=
A
11
⋅
S
1
P_1 = A_{11}\cdot S_1
P1=A11⋅S1,
P
2
=
S
2
⋅
B
22
P_2 = S_2 \cdot B_{22}
P2=S2⋅B22,
P
3
=
S
3
⋅
B
11
P_3 = S_3 \cdot B_{11}
P3=S3⋅B11,
P
4
=
A
22
⋅
S
4
P_4 = A_{22} \cdot S_4
P4=A22⋅S4,
P
5
=
S
5
⋅
S
6
P_5 = S_5 \cdot S_6
P5=S5⋅S6,
P
6
=
S
7
⋅
S
8
P_6 = S_7 \cdot S_8
P6=S7⋅S8,
P
7
=
S
9
⋅
S
10
P_7 = S_9 \cdot S_{10}
P7=S9⋅S10。
在步骤 4 中,
C
11
=
P
5
+
P
4
−
P
2
+
P
6
C_{11} = P_5 + P_4 - P_2 + P_6
C11=P5+P4−P2+P6,
C
12
=
P
1
+
P
2
C_{12} = P_{1} + P_{2}
C12=P1+P2,
C
21
=
P
3
+
P
4
C_{21} = P_{3} + P_4
C21=P3+P4,
C
22
=
P
5
+
P
1
−
P
3
−
P
7
C_{22} = P_5 + P_1 - P_3 - P_7
C22=P5+P1−P3−P7。
练习
4.2-1 使用 Strassen 算法计算矩阵积
(
1
3
7
5
)
(
6
8
4
2
)
\begin{pmatrix} 1 & 3 \\ 7 & 5 \end{pmatrix} \begin{pmatrix} 6 & 8 \\ 4 & 2 \end{pmatrix}
(1735)(6482)。
解:
A
11
=
(
1
)
A_{11} = (1)
A11=(1),
A
12
=
(
3
)
A_{12} = (3)
A12=(3),
A
21
=
(
7
)
A_{21} = (7)
A21=(7),
A
22
=
(
5
)
A_{22} = (5)
A22=(5),
B
11
=
(
6
)
B_{11} = (6)
B11=(6),
B
12
=
(
8
)
B_{12} = (8)
B12=(8),
B
21
=
(
4
)
B_{21} = (4)
B21=(4),
B
22
=
(
2
)
B_{22} = (2)
B22=(2)。
S
1
=
(
6
)
S_1 = (6)
S1=(6),
S
2
=
(
4
)
S_2 = (4)
S2=(4),
S
3
=
(
12
)
S_3 = (12)
S3=(12),
S
4
=
(
−
2
)
S_4 = (-2)
S4=(−2),
S
5
=
(
6
)
S_5 = (6)
S5=(6),
S
6
=
(
8
)
S_6 = (8)
S6=(8),
S
7
=
(
−
2
)
S_7 = (-2)
S7=(−2),
S
8
=
(
6
)
S_8 = (6)
S8=(6),
S
9
=
(
−
6
)
S_9 = (-6)
S9=(−6),
S
10
=
(
14
)
S_{10} = (14)
S10=(14)。
P
1
=
(
6
)
P_1 = (6)
P1=(6),
P
2
=
(
8
)
P_2 = (8)
P2=(8),
P
3
=
(
72
)
P_3 = (72)
P3=(72),
P
4
=
(
−
10
)
P_4 = (-10)
P4=(−10),
P
5
=
(
48
)
P_5 =(48)
P5=(48),
P
6
=
(
−
12
)
P_6 = (-12)
P6=(−12),
P
7
=
(
−
84
)
P_7 = (-84)
P7=(−84)。
C
11
=
(
18
)
C_{11} = (18)
C11=(18),
C
12
=
(
14
)
C_{12} = (14)
C12=(14),
C
21
=
(
62
)
C_{21} = (62)
C21=(62),
C
22
=
(
66
)
C_{22} = (66)
C22=(66)。
所以,所求结果为
C
=
(
18
14
62
66
)
C = \begin{pmatrix} 18 & 14 \\ 62 & 66 \end{pmatrix}
C=(18621466)。
4.2-2 写出 Strassen 算法的伪代码。
解: SQUARE-MATRIX-MULTIPLY-STRASSEN
(
A
,
B
)
(A, B)
(A,B)
n
=
A
.
r
o
w
s
n = A.rows
n=A.rows
let
C
C
C be a new
n
×
n
n\times n
n×n matrix
if
n
=
=
1
n==1
n==1
c
11
=
a
11
⋅
b
11
c_{11} = a_{11} \cdot b_{11}
c11=a11⋅b11
else partition
A
A
A,
B
B
B and
C
C
C into
n
/
2
×
n
/
2
n/2 \times n/2
n/2×n/2 matrices
S
1
=
B
12
−
B
22
S_1 = B_{12} - B_{22}
S1=B12−B22
S
2
=
A
11
+
A
12
S_2 = A_{11} + A_{12}
S2=A11+A12
S
3
=
A
21
+
A
22
S_3 = A_{21} + A_{22}
S3=A21+A22
S
4
=
B
21
−
B
11
S_4 = B_{21} - B_{11}
S4=B21−B11
S
5
=
A
11
+
A
22
S_5 = A_{11} + A_{22}
S5=A11+A22
S
6
=
B
11
+
B
22
S_6 = B_{11} + B_{22}
S6=B11+B22
S
7
=
A
12
−
A
22
S_7 = A_{12} - A_{22}
S7=A12−A22
S
8
=
B
21
+
B
22
S_8 = B_{21} + B_{22}
S8=B21+B22
S
9
=
A
11
−
A
21
S_9 = A_{11} - A_{21}
S9=A11−A21
S
10
=
B
11
+
B
12
S_{10} = B_{11} + B_{12}
S10=B11+B12
P
1
=
P_1 =
P1= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
A
11
,
S
1
)
(A_{11}, S_1)
(A11,S1)
P
2
=
P_2 =
P2= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
S
2
,
B
22
)
(S_2, B_{22})
(S2,B22)
P
3
=
P_3 =
P3= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
S
3
,
B
11
)
(S_3, B_{11})
(S3,B11)
P
4
=
P_4 =
P4= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
A
22
,
S
4
(A_{22},S_4
(A22,S4)
P
5
=
P_5 =
P5= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
S
5
,
S
6
)
(S_5,S_6)
(S5,S6)
P
6
=
P_6 =
P6= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
S
7
,
S
8
)
(S_7, S_8)
(S7,S8)
P
7
=
P_7 =
P7= SQUARE-MATRIX-MULTIPLY-STRASSEN
(
S
9
,
S
10
)
(S_9, S_{10})
(S9,S10)
C
11
=
P
5
+
P
4
−
P
2
+
P
6
C_{11} = P_5 + P_4 - P_2 + P_6
C11=P5+P4−P2+P6
C
12
=
P
1
+
P
2
C_{12} = P_{1} + P_{2}
C12=P1+P2
C
21
=
P
3
+
P
4
C_{21} = P_{3} + P_4
C21=P3+P4
C
22
=
P
5
+
P
1
−
P
3
−
P
7
C_{22} = P_5 + P_1 - P_3 - P_7
C22=P5+P1−P3−P7
return
C
C
C
4.2-3 两个
n
×
n
n \times n
n×n 矩阵相乘,其中
n
n
n 不是 2 的幂次,如何修改 Strassen 算法?证明得到的算法的运行时间是
Θ
(
n
lg
7
)
\Theta(n ^{\lg7})
Θ(nlg7)。
解: 将输入矩阵
A
A
A,
B
B
B 和输出矩阵
C
C
C 分解成
⌊
n
/
2
⌋
×
⌊
n
/
2
⌋
\lfloor n/2\rfloor \times \lfloor n/2\rfloor
⌊n/2⌋×⌊n/2⌋ 子矩阵、
⌊
n
/
2
⌋
×
⌈
n
/
2
⌉
\lfloor n/2\rfloor \times \lceil n/2\rceil
⌊n/2⌋×⌈n/2⌉ 子矩阵、
⌈
n
/
2
⌉
×
⌊
n
/
2
⌋
\lceil n/2\rceil \times \lfloor n/2\rfloor
⌈n/2⌉×⌊n/2⌋ 子矩阵和
⌈
n
/
2
⌉
×
⌈
n
/
2
⌉
\lceil n/2\rceil \times \lceil n/2\rceil
⌈n/2⌉×⌈n/2⌉ 子矩阵。
证明:使用主方法证明(见章节 4.5)。
描述 Strassen 算法运行时间的递归式为递归式 (4.18),
T
(
n
)
=
7
T
(
n
/
2
)
+
Θ
(
n
2
)
T(n) = 7T(n/2) + \Theta(n^2)
T(n)=7T(n/2)+Θ(n2)。
有
a
=
7
a = 7
a=7,
b
=
2
b = 2
b=2,
f
(
n
)
=
Θ
(
n
2
)
f(n) = \Theta(n^2)
f(n)=Θ(n2)。
所以
n
log
b
a
=
n
log
2
7
=
n
lg
7
n^{\log_b a} = n^{\log_2 7} = n^{\lg 7}
nlogba=nlog27=nlg7,又
2.80
<
lg
7
<
2.81
2.80 < \lg7 <2.81
2.80<lg7<2.81,
所以,存在
ϵ
=
0.8
\epsilon = 0.8
ϵ=0.8,有
f
(
n
)
=
O
(
n
lg
7
−
ϵ
)
f(n) = O(n^{\lg7 - \epsilon})
f(n)=O(nlg7−ϵ),符合主定理的情况①。
所以,
T
(
n
)
=
Θ
(
n
lg
7
)
T(n) = \Theta(n^{\lg7})
T(n)=Θ(nlg7)。
4.2-4 如果可以使用 k k k 次乘法操作完成 3 × 3 3 \times 3 3×3 矩阵相乘(不考虑乘法的可交换性),那么可以在 o ( n lg 7 ) o(n^{\lg 7}) o(nlg7) 时间内完成 n × n n\times n n×n 矩阵相乘,在这种情况下, k k k 的最大值是多少?这个算法的运行时间是多少?
解: 如果可以使用
k
k
k 次乘法操作完成
3
×
3
3 \times 3
3×3 矩阵相乘,那么可以通过递归地进行
n
/
3
×
n
/
3
n/3 \times n/3
n/3×n/3 矩阵相乘,完成
n
×
n
n\times n
n×n 矩阵相乘,其运行时间是
T
(
n
)
=
k
T
(
n
/
3
)
+
Θ
(
n
2
)
T(n) = kT(n/3) + \Theta(n^2)
T(n)=kT(n/3)+Θ(n2)。
例如,等式
C
=
A
⋅
B
C = A \cdot B
C=A⋅B 可以写成
(
C
11
C
12
C
13
C
21
C
22
C
23
C
31
C
32
C
33
)
=
(
A
11
A
12
A
13
A
21
A
22
A
23
A
31
A
32
A
33
)
⋅
(
B
11
B
12
B
13
B
21
B
22
B
23
B
31
B
32
B
33
)
\begin{pmatrix} C_{11} & C_{12} & C_{13} \\ C_{21} & C_{22} & C_{23} \\ C_{31} & C_{32} & C_{33} \end{pmatrix} = \begin{pmatrix} A_{11} & A_{12} & A_{13} \\ A_{21} & A_{22} & A_{23} \\ A_{31} & A_{32} & A_{33} \end{pmatrix} \cdot \begin{pmatrix} B_{11} & B_{12} & B_{13} \\ B_{21} & B_{22} & B_{23} \\ B_{31} & B_{32} & B_{33} \end{pmatrix}
⎝⎛C11C21C31C12C22C32C13C23C33⎠⎞=⎝⎛A11A21A31A12A22A32A13A23A33⎠⎞⋅⎝⎛B11B21B31B12B22B32B13B23B33⎠⎞
使用主方法解这个递归式,考虑
n
log
3
k
n^{\log_3 k}
nlog3k 和
n
2
n^2
n2 的比值:
如果
log
3
k
=
2
\log_3 k = 2
log3k=2,符合情况②,
T
(
n
)
=
Θ
(
n
2
lg
n
)
T(n) = \Theta(n^2\lg n)
T(n)=Θ(n2lgn)。这种情况下,
k
=
9
k=9
k=9,且
T
(
n
)
=
o
(
n
lg
7
)
T(n) = o(n^{\lg7})
T(n)=o(nlg7)。
如果
log
3
k
<
2
\log_3k < 2
log3k<2,符合情况③,
T
(
n
)
=
Θ
(
n
2
)
T(n) = \Theta(n^2)
T(n)=Θ(n2)。这种情况下,
k
<
9
k < 9
k<9,且
T
(
n
)
=
o
(
n
2
)
T(n) = o(n^2)
T(n)=o(n2)。
如果
log
3
k
>
2
\log_3k > 2
log3k>2,符合情况①,
T
(
n
)
=
Θ
(
n
log
3
k
)
T(n) = \Theta(n^{\log_3 k})
T(n)=Θ(nlog3k)。这种情况下,
k
>
9
k > 9
k>9。当
log
3
k
<
lg
7
\log_3k < \lg7
log3k<lg7 时(即
k
<
3
lg
7
≈
21.85
k < 3^{\lg7} \approx 21.85
k<3lg7≈21.85),
T
(
n
)
=
o
(
n
2
)
T(n) = o(n^2)
T(n)=o(n2)。
k
k
k 的最大值为 21。
因此, k = 21 k=21 k=21,运行时间为 Θ ( n log 3 k ) = Θ ( n log 3 21 ) = O ( n 2.80 ) \Theta(n^{\log_3k}) = \Theta(n^{\log_3 21}) = O(n^{2.80}) Θ(nlog3k)=Θ(nlog321)=O(n2.80)。(因为 log 3 21 ≈ 2.77 \log_3 21 \approx 2.77 log321≈2.77)
4.2-5 V. Pan 发现了一个使用 132,464 次乘法操作完成
68
×
68
68 \times 68
68×68 矩阵相乘的方式,一个使用 143,640 次乘法操作完成
70
×
70
70 \times 70
70×70 矩阵相乘的方式,一个使用 155,424 次乘法操作完成
72
×
72
72 \times 72
72×72 矩阵相乘的方式。当使用在分治的矩阵相乘算法中,哪个方法产生最好的渐近运行时间?与 Strassen 算法相比如何?
解: 这 3 种方法的运行时间的递归式分别为:
T
(
n
)
=
132464
T
(
n
/
68
)
+
Θ
(
n
2
)
T(n) = 132464T(n/68) + \Theta(n^2)
T(n)=132464T(n/68)+Θ(n2);
H
(
n
)
=
143640
S
(
n
/
70
)
+
Θ
(
n
2
)
H(n) = 143640S(n/70) + \Theta(n^2)
H(n)=143640S(n/70)+Θ(n2);
Q
(
n
)
=
155424
Q
(
n
)
+
Θ
(
n
2
)
Q(n) = 155424Q(n) + \Theta(n^2)
Q(n)=155424Q(n)+Θ(n2)。
使用主方法解递归式,得
T
(
n
)
=
Θ
(
n
log
68
132464
)
T(n) = \Theta(n^{\log_{68} 132464})
T(n)=Θ(nlog68132464),(
log
68
132464
≈
2.79512848
\log_{68} 132464 \approx 2.79512848
log68132464≈2.79512848)
H
(
n
)
=
Θ
(
n
log
70
143640
)
H(n) = \Theta(n^{\log_{70} 143640})
H(n)=Θ(nlog70143640),(
log
70
143640
≈
2.79512268
\log_{70} 143640 \approx 2.79512268
log70143640≈2.79512268)
Q
(
n
)
=
Θ
(
n
log
72
155424
)
Q(n) = \Theta(n^{\log_{72} 155424})
Q(n)=Θ(nlog72155424)。(
log
68
132464
≈
2.79514739
\log_{68} 132464 \approx 2.79514739
log68132464≈2.79514739)
所以,第 2 个方法产生最好的渐近运行时间。比 Strassen 算法的渐近运行时间好。
4.2-6 使用 Strassen 算法作为子程序,一个
k
n
×
n
kn \times n
kn×n 矩阵与一个
n
×
k
n
n \times kn
n×kn 矩阵相乘有多快?将输入矩阵交换位置,回答相同的问题。
解: 本题中两个矩阵相乘示例如下,
(
C
11
C
12
…
C
1
k
C
21
C
22
…
C
2
k
…
C
k
1
C
k
2
…
C
k
k
)
=
(
A
11
A
21
…
A
k
1
)
⋅
(
B
11
B
12
…
B
1
k
)
\begin{pmatrix} C_{11} & C_{12} & \dots & C_{1k} \\ C_{21} & C_{22} & \dots & C_{2k} \\ \dots \\C_{k1} & C_{k2} & \dots & C_{kk} \end{pmatrix} = \begin{pmatrix} A_{11} \\ A_{21} \\ \dots \\ A_{k1} \end{pmatrix} \cdot \begin{pmatrix} B_{11} & B_{12} & \dots & B_{1k} \end{pmatrix}
⎝⎜⎜⎛C11C21…Ck1C12C22Ck2………C1kC2kCkk⎠⎟⎟⎞=⎝⎜⎜⎛A11A21…Ak1⎠⎟⎟⎞⋅(B11B12…B1k)
上面的
A
i
j
A_{ij}
Aij、
B
i
j
B_{ij}
Bij 和
C
i
j
C_{ij}
Cij 都是
n
×
n
n \times n
n×n 矩阵。
C
11
=
A
11
⋅
B
11
,
C
12
=
A
11
⋅
B
12
C_{11} = A_{11} \cdot B_{11},C_{12} = A_{11} \cdot B_{12}
C11=A11⋅B11,C12=A11⋅B12,
…
\dots
…,
C
1
k
=
A
11
⋅
B
1
k
C_{1k} = A_{11} \cdot B_{1k}
C1k=A11⋅B1k,
C
21
=
A
21
⋅
B
11
C_{21} = A_{21} \cdot B_{11}
C21=A21⋅B11,
C
22
=
A
21
⋅
B
12
C_{22} = A_{21} \cdot B_{12}
C22=A21⋅B12,
…
\dots
…,
C
2
k
=
A
21
⋅
B
1
k
C_{2k} = A_{21} \cdot B_{1k}
C2k=A21⋅B1k,
……
C
k
1
=
A
k
1
⋅
B
11
C_{k1} = A_{k1} \cdot B_{11}
Ck1=Ak1⋅B11,
C
k
2
=
A
k
1
⋅
B
12
C_{k2} = A_{k1} \cdot B_{12}
Ck2=Ak1⋅B12,
…
\dots
…,
C
k
k
=
A
k
1
⋅
B
1
k
C_{kk} = A_{k1} \cdot B_{1k}
Ckk=Ak1⋅B1k。
因为使用 Strassen 算法作为子程序,所以每次进行
C
i
j
C_{ij}
Cij 的计算的运行时间为
Θ
(
n
lg
7
)
\Theta(n^{\lg7})
Θ(nlg7),一共计算
k
2
k^2
k2 次。
所以本题矩阵相乘的运行时间为
Θ
(
k
2
n
lg
7
)
\Theta(k^2 n^{\lg7})
Θ(k2nlg7)。
将输入矩阵交换位置,运行时间相同,是 Θ ( k 2 n lg 7 ) \Theta(k^2 n^{\lg7}) Θ(k2nlg7)。
4.2-7 设计算法,仅使用 3 次实数乘法,实现复数
a
+
b
i
a+bi
a+bi 和
c
+
d
i
c+di
c+di 的相乘。这个算法接受
a
a
a,
b
b
b,
c
c
c 和
d
d
d 作为输入,分别生成实数部分
a
c
−
b
d
ac-bd
ac−bd 和虚数部分
a
d
+
b
c
ad+bc
ad+bc。
解: MULTIPLY-COMPLEX(
a
,
b
,
c
,
d
a,b,c,d
a,b,c,d)
m
u
l
t
1
=
(
a
+
b
)
⋅
(
c
+
d
)
mult1 = (a+b)\cdot(c+d)
mult1=(a+b)⋅(c+d) // 结果为 ac + ad + bc + bd
m
u
l
t
2
=
a
⋅
c
mult2 = a\cdot c
mult2=a⋅c
m
u
l
t
3
=
b
⋅
d
mult3 = b\cdot d
mult3=b⋅d
a
n
s
1
=
m
u
l
t
2
−
m
u
l
t
3
ans1 = mult2 - mult3
ans1=mult2−mult3
a
n
s
2
=
m
u
l
t
1
−
m
u
l
t
2
−
m
u
l
t
3
ans2 = mult1 - mult2 - mult3
ans2=mult1−mult2−mult3
return
(
a
n
s
1
,
a
n
s
2
)
(ans1,ans2)
(ans1,ans2)
学习笔记目录:【算法导论】目录