算法 : 第4版 / (美 塞奇威客 (Sedgewick, R. , (美 韦恩 (Wayne, K. 著 ; 谢路云译. -- 北京 : 人民邮电出版社, 2012.10 (2021.5重印（以下简称原书或书）中的练习题 1.1.27 的三种解法（C++ 实现），并对包含原书题中的递归方法在内的四种解法的执行时间进行了统计和对比。
◆ 要求
b(n,k,p = 1.0 ( n == 0 && k == 0
b(n,k,p = 0.0 ( n < 0 || k < 0
b(n,k,p = (1.0-p * b(n-1,k,p + p * b(n-1,k-1,p
◆ 解一
依然采用递归的方式，但使用二维数组保存中间结果。
static long double binomial1(int N, int K, long double p // #1
{
long double x;
long double** b = new long double*[N+1]; // #2
long double* data = new long double[(N+1*(K+1];
...
x = binomial1_impl(b, N, K, p; // #3
...
return x;
}
static long double binomial1_impl(long double** b, int N, int K, long double p
{
if (N == 0 && K == 0 return 1.0L;
if (N < 0 || K < 0 return 0.0L;
if (b[N][K] == -1 {
... // #4
b[N][K] = (1.0L-p * binomial1_impl(b, N-1, K, p + p * binomial1_impl(b, N-1, K-1, p;
}
return b[N][K];
}
保持对外的接口不变（#1），创建一个二维数组 b[0..N][0..K] 保存中间计算结果（#2），并将其传给算法实现（#3）。算法虽然还是用递归调用（#4），但由于中间结果保存在全局的二维数组中，不用频繁地压栈和弹栈去获取中间数据。此解法网络上也见于 [github] reneargento/algorithms-sedgewick-wayne 和 [github] aistrate/AlgorithmsSedgewick 。
◆ 解二
以 N = K = 4 为例，
0 1 2 3 4
0 * * * * * <-- * 代表待计算的单元
1 * * * * *
2 * * * * *
3 * * * * *
4 * * * * ? <-- 最终计算结果的单元 (4,4
仔细考察递归关系式的特点，b(-1,*,p = 0.0, b(*,-1,p = 0.0。由
b(0,1,p = (1.0-p * b(-1,1,p + p * b(-1,0,p
= (1.0-p * 0.0 + p * 0.0
= 0.0
b(0,2,p = (1.0-p * b(-1,1,p + p * b(-1,1,p
= (1.0-p * 0.0 + p * 0.0
= 0.0
...
可推论出，二维数组中的第 0 行中的所有单元（不含b[0][0]）均为 0.0；由
b(1,0,p = (1.0-p * b(0,0,p + p * b(0,-1,p
= (1.0-p * 1.0 + p * 0.0
= 1.0-p
b(2,0,p = (1.0-p * b(1,0,p + p * b(1,-1,p
= (1.0-p * (1.0-p + p * 0.0
= (1.0-p^2
...
可推论出，二维数组中的第 0 列的单元为 (1.0-p^y。
第一种情况，N < K：（如 N = 4, K = 6）
0 1 2 3 4 5 6
0 - - G # # # # <-- G 点所在单元为 0.0
1 - - - * * * * <-- '-' 代表不用计算的单元
2 - - - - * * *
3 - - - - - * *
4 - - - - - - ? <-- 最终结果的存储单元
G 点为 b(0,K-N。按照递推关系式容易推导出，'#' 和 '*' 区域均为 0.0，所以最终结果即 0.0。
0 1 2 3 4 5 6
0 G # # # # # # <-- G 点所在单元为 1.0
1 - * * * * * *
2 - - * * * * *
3 - - - * * * *
4 - - - - * * *
5 - - - - - * *
6 - - - - - - ?
G 点为 b(0,0。按照递推关系式容易推导出，数组中 n = k 的单元为 p^n。所以最终结果即 p^N。
0 1 2 3 4
0 # # # # #
1 # # # # #
2 G # # # # <-- G 点所在单元为 (1.0-p^2
3 - * * * *
4 - - * * *
5 - - - * *
6 - - - - ?
G 点为 b(N-K,0。可先计算 '#' 区域中的单元，再计算 '*' 区域中的单元，得出最终结果。处理'#'区域时，为避免大量的数组下标越界判断，可以考虑先计算 0 行和 0 列的所有单元。
static long double binomial2(int N, int K, long double p
{
long double x;
if (N < K { // #1
x = 0.0L;
} else if (N == K { // #2
x = powl(p, N;
} else { // #3
...
b[0][0] = 1.0L;
// process '#' area // #4
// calcuate [1..N-K][0]
for (i = 1; i <= N-K; ++i
b[i][0] = powl(1.0L-p, i;
// calcuate [0][1..K]
for (j = 1; j <= K; ++j
b[0][j] = 0.0L;
// calcuate [1..N-K][1..K]
for (i = 1; i <= N-K; ++i
for (j = 1; j <= K; ++j
b[i][j] = (1.0L-p * b[i-1][j] + p * b[i-1][j-1];
// process '*' area // #5
for (i = N-K+1; i <= N; ++i
for (j = i-(N-K; j <= K; ++j
b[i][j] = (1.0L-p * b[i-1][j] + p * b[i-1][j-1];
x = b[N][K]; // #6
...
}
return x;
}
三条分支（#1、#2、#3）分别对应前述三种情况。在第三种情况下，再先处理 '#' 区域（#4），然后采用递推求值的方式处理 '*' 区域（#5），最后得到结果（#6）。
◆ 解三
0 1 2 3 4
0 1.0
1 1.0-p p
2 (1.0-p^2 2*(1.0-p*p p^2
3 (1.0-p^3 3*[(1.0-p^2]*p 3*(1-p*(p^2 p^3
4 (1.0-p^4 4*[(1.0-p^3]*p 6*[(1.0-p^2]*(p^2 4*(1.0-p*(p^3 p^4
可以发现，从第 0 行到第 N 行的非零单元即“杨辉三角形”，第 n 行中的非零单元之和构成 [(1.0-p + p]^k 的展开式。因此，解二中的第三种情况，可结合利用通项公式 C(N,K*[(1.0-p^(N-K]*(p^K 来解决。
static long double binomial3(int N, int K, long double p
{
long double x;
if ...
} else {
x = combination(N, K * powl(1.0L-p, N-K * powl(p, K;
}
return x;
}
◆ 测试
编译并执行程序，
$ g++ -std=c++11 -o 27.out 27.cpp
$ ./27.out 10 5 0.25
为易于显示两者之间的差异，笔者选择了硬件配置偏低的测试环境。
- 硬件配置：Raspberry Pi 3 Model B
- Quad Core 1.2GHz 64bit
- 1G RAM
- 16G MicroSD
- 100 Base Ethernet
- 软件配置：Raspbian Stretch
- g++ (Raspbian 6.3.0-18+rpi1+deb9u1 6.3.0 20170516
结果如下图所示，
◆ 最后
完整示例代码和测试结果，请参考 [gitee] cnblogs/17328989 。
写作过程中，笔者参考了 [github] reneargento/algorithms-sedgewick-wayne 和 [github] aistrate/AlgorithmsSedgewick 的实现方式。致 reneargento 和 aistrate 两位。