算法導論(第三版)參考答案:練習4.2-1,練習4.2-2,練習4.2-3,練習4.2-4,練習4.2-5,練習4.2-6,練習4.2-7
Exercise 4.2-1
Use Strassen’s algorithm to compute the matrix product
>(>1>735>)>(>6>482>)>
Show your work.
步驟1:
A11=(1),A12=(2),B21=(7),B22=(5)B11=(6),B12=(8),B21=(4),B22=(2)
步驟2;
S1=B12−B22=6S2=A11+A12=4S3=A21+A22=12S4=B21−B11=−2S5=A11+A22=6S6=B11+B22=8S7=A12−A22=−2S8=B21+B22=6S9=A11+A21=−6S10=B11+B12=14
步驟3:
P1P2P3P4P5P6P7=A11⋅S1=1⋅6=6=S2⋅B22=4⋅2=8=S3⋅B11=6⋅12=72=A22⋅S4=6⋅12=72=S5⋅S6=6⋅8=48=S7⋅S8=(−2)⋅6=−12=S9⋅S10=(−6)⋅14=−84
步驟4:
C11C12C21C22=P5+P4−P2+P6=48+(−10)−8+(−12)=18=P1+P2=6+8=14=P3+P4=72+(−10)=62=P5+P1−P3−P7=48+6−72−(−84)=66
結果為:
C=(C11C21C12C22)=(18621466)
Exercise 4.2-2
Write pseudocode for Strassen’s algorithm
Pseudocode
STRASSEN-MATRIX-MULTIPLY-RECURSIVE(A, B)
n = A.rows
let C be a new n×n matrix
if n ==
c11 = a11 * b11
else partition A, B, and C as in equations ()
calculate addition: S1-S10
P1 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(A11, S1)
p2 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(S2, B22)
P3 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(S3, B11)
P4 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(A22, S4)
P5 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(S5, S6)
P6 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(S7, S8)
P7 = STRASSEN-MATRIX-MULTIPLY-RECURSIVE(S9, S10)
C11 = P5 + P4 - P2 + P6
C12 = P1 + P2
C21 = P3 + P4
C22 = P5 + P1 - P3 - P7
return C
C code
#include <alloca.h>
#include <stdio.h>
// The matrix representation. One structure to fit both a matrix and a
// submatrix.
typedef struct {
int x;
int y;
int size;
int original_size;
int *data;
} matrix;
// Functions to index matrices
int get(matrix m, int x, int y) {
return m.data[m.original_size * (m.x + x) + m.y + y];
};
void put(matrix m, int x, int y, int value) {
m.data[m.original_size * (m.x + x) + m.y + y] = value;
};
// Matrix building
matrix create_matrix(int size, int *data) {
matrix result;
result.x = ;
result.y = ;
result.size = size;
result.original_size = size;
result.data = data;
return result;
}
matrix submatrix(matrix A, int x, int y, int size) {
matrix result;
result.x = A.x + x;
result.y = A.y + y;
result.size = size;
result.original_size = A.original_size;
result.data = A.data;
return result;
}
#define INIT_ON_STACK(m_, size_) \
m_.x = ; \
m_.y = ; \
m_.size = size_; \
m_.original_size = size_; \
m_.data = alloca(size_ * size_ * sizeof(int));
// Adding and subtracting matrices
void plus(matrix C, matrix A, matrix B) {
for (int i = ; i < C.size; i++) {
for (int j = ; j < C.size; j++) {
put(C, i, j, get(A, i, j) + get(B, i, j));
}
}
}
void minus(matrix C, matrix A, matrix B) {
for (int i = ; i < C.size; i++) {
for (int j = ; j < C.size; j++) {
put(C, i, j, get(A, i, j) - get(B, i, j));
}
}
}
void add(matrix T, matrix S) {
for (int i = ; i < T.size; i++) {
for (int j = ; j < T.size; j++) {
put(T, i, j, get(T, i, j) + get(S, i, j));
}
}
}
void sub(matrix T, matrix S) {
for (int i = ; i < T.size; i++) {
for (int j = ; j < T.size; j++) {
put(T, i, j, get(T, i, j) - get(S, i, j));
}
}
}
void zero(matrix m) {
for (int i = ; i < m.size; i++) {
for (int j = ; j < m.size; j++) {
put(m, i, j, );
}
}
}
// A function to print matrices
void print_matrix(matrix m) {
printf("%dx%d (+%d+%d) (%d)\n", m.size, m.size, m.x, m.y, m.original_size);
printf("==============\n");
for (int i = ; i < m.size; i++) {
for (int j = ; j < m.size; j++) {
printf("%4d", get(m, i, j));
}
printf("\n");
}
printf("\n");
}
// Strassen's algorithm
void strassen(matrix C, matrix A, matrix B) {
int size = A.size,
half = size / ;
if (A.size == ) {
put(C, , , get(A, , ) * get(B, , ));
} else {
matrix s1, s2, s3, s4, s5, s6, s7, s8, s9, s10;
matrix p1, p2, p3, p4, p5, p6, p7;
INIT_ON_STACK(s1, half);
INIT_ON_STACK(s2, half);
INIT_ON_STACK(s3, half);
INIT_ON_STACK(s4, half);
INIT_ON_STACK(s5, half);
INIT_ON_STACK(s6, half);
INIT_ON_STACK(s7, half);
INIT_ON_STACK(s8, half);
INIT_ON_STACK(s9, half);
INIT_ON_STACK(s10, half);
INIT_ON_STACK(p1, half);
INIT_ON_STACK(p2, half);
INIT_ON_STACK(p3, half);
INIT_ON_STACK(p4, half);
INIT_ON_STACK(p5, half);
INIT_ON_STACK(p6, half);
INIT_ON_STACK(p7, half);
matrix a11 = submatrix(A, , , half);
matrix a12 = submatrix(A, , half, half);
matrix a21 = submatrix(A, half, , half);
matrix a22 = submatrix(A, half, half, half);
matrix b11 = submatrix(B, , , half);
matrix b12 = submatrix(B, , half, half);
matrix b21 = submatrix(B, half, , half);
matrix b22 = submatrix(B, half, half, half);
matrix c11 = submatrix(C, , , half);
matrix c12 = submatrix(C, , half, half);
matrix c21 = submatrix(C, half, , half);
matrix c22 = submatrix(C, half, half, half);
minus(s1, b12, b22);
plus(s2, a11, a12);
plus(s3, a21, a22);
minus(s4, b21, b11);
plus(s5, a11, a22);
plus(s6, b11, b22);
minus(s7, a12, a22);
plus(s8, b21, b22);
minus(s9, a11, a21);
plus(s10, b11, b12);
strassen(p1, a11, s1);
strassen(p2, s2, b22);
strassen(p3, s3, b11);
strassen(p4, a22, s4);
strassen(p5, s5, s6);
strassen(p6, s7, s8);
strassen(p7, s9, s10);
zero(c11);
zero(c12);
zero(c21);
zero(c22);
add(c11, p5);
add(c11, p4);
sub(c11, p2);
add(c11, p6);
add(c12, p1);
add(c12, p2);
add(c21, p3);
add(c21, p4);
add(c22, p5);
add(c22, p1);
sub(c22, p3);
sub(c22, p7);
}
}
Exercise 4.2-3
How would you modify Strassen’s algorithm to multiply n×n matrices in which n is not an exact power of 2? Show that the resulting algorithm runs in time Θ(nlg7) .
可以将矩陣填充為 n×n 的矩陣,在 Θ(nlg7) 内求解完成後,再将填充元素剝離掉。
Exercise 4.2-4
What is the largest k such that if you can multiply 3×3 matrices using k multiplications (not assuming commutativity of multiplication), then you can multiply n×n matrices is time o(nlg7) ? What would the running time of this algorithm be?
根據題意,可得如下遞歸式:
T(n)=kT(n/3)+Θ(n)
根據主方法可知運作時間: Θ(nlog3k) 。比較 Θ(nlog3k) 和 Θ(nlg7)
nlog3k<nlg7⇓log3k<lg7⇓k<3lg7≈21.84986
Exercise 4.2-5
V. Pan has discovered a way of multiplying 68×68 matrices using 132,464 multiplications, a way of multiplying 70×70 matrices using 143,640 multiplications, and a way of multiplying 72×72 matrices using 155,424 multiplications. Which method yields the best asymptotic running time when used in a divide-and-conquer matrix-multiplication algorithm? How does it compare to Strassen’s algorithm?
log68132464≈2.795128log70143640≈2.795122log72155424≈2.795147
70×70 的矩陣乘法最快, lg7≈2.81 ,比Strassen算法好。
Exercise 4.2-6
How quickly can you multiply a kn×n matrix by an n×kn matrix, using Strassen’s algorithm as a subroutine? Answer the same question with the order of the input matrices reversed.
(kn×n)⋅(n×kn) 需要 k2 個 (n×n) 矩陣相乘,即 T=k2⋅Θ(nlg7) 。
(n×kn)⋅(kn×n) 隻需要 k 個 (n×n) 矩陣相乘,加上 k−1 個 (n×n) 矩陣相加,即 T=k⋅Θ(nlg7)+(k−1)Θ(n2)=k⋅Θ(nlg7) 。
Exercise 4.2-7
Show how to multiply the complex numbers a+bi and c+di using only three multiplications of real numbers. The algorithm should take a , b, c and d as input and produce the real component ac−bd and the imaginary component ad+bc separately.
(a+bi)⋅(c+di)=(ac−bd)+(ad+bc)i=ac−bd+[(a+b)⋅(c+d)−ac−bd]i