增长的阶
练习1.15 b) 在求值( sine a )时,由过程$sine$所产生的计算过程使用的空间和增长的阶是什么?
(define(cube x) ( * x x x ))
(define( p x ) ( - ( * 3 x ) ( * 4 ( cube x ) ) ))
(
define( sine _angle )
(
if ( not ( > ( abs _angle ) 0.1 ) )
_angle
( p ( sine ( / _angle 3.0 ) ) )
)
)
根据主定理分析, $T(n) = aT( \frac{b}{n} )+ f(n) $ 其中$a = 1, b = 3, f(n) = \Theta ( 1 )$
符合主定理的第二条,若$ f(n) = \Theta( n ^ \log_b a ) $, 那么$T( n ) = \Theta( {n ^ \log_b a} log n )$
$ \log_b a = \log_3 1 = 0$ $\therefore $ $T( n ) = \Theta( {n ^ \log_b a} log n ) = \Theta( log n ) $
因此计算sine增长的阶是$\Theta(log n)$
练习1.26 将原( expmod )的算法更改了一下,时间复杂度发生了改变,由O(log n)变成了O( n ),试分析两个算法的时间复杂度.原算法在下文的快速幂取模.以下是更改后的代码.
(define (expmod base exp m)
(
cond ( ( = exp 0 ) 1 )
( (even? exp)
( remainder ( * ( expmod base ( / exp 2 ) m )
( expmod base ( / exp 2 ) m )
) m ))
(
else ( remainder ( * base ( expmod base ( - exp 1 ) m ) ) m )
)
)
)
[latexpage]
原算法的递推关系式是$T( n ) = T( \frac{n}{2} ) + \Theta( 1 )$ , 根据主定理计算出算法的时间复杂度为$ \Theta( log n ) $
更改后的递推关系式是$T( n ) = 2T( \frac{n}{2} ) + \Theta( 1 )$, 根据主定理计算出算法的时间复杂度为$ \Theta( n ) $
看上去是相同的计算实际上square(x)将计算减少了很多
幂运算
朴素幂运算
$ b^n = b * b ^ {n - 1} $
$ b ^ 0 = 1 $
可以直接翻译为以下过程
(
define ( expt b n )
( if ( = n 0 )
1
( * b ( expt b ( - n 1 ) ) )
)
)
这是一个线性的递归计算过程,需要$\Theta (n)$的空间( 因为递归运算的堆栈需要空间 )以及
$ \Theta (n) $的时间.
一个将空间优化为$ \Theta(1) $ 的迭代版本
(
define ( expt2 b n )
( expt_iter b n 1 )
)
(
define ( expt_iter b counter product )
(
if ( = counter 0 )
product
(
expt_iter b ( - counter 1 ) ( * b product )
)
)
)
快速幂运算
在计算$b ^ 8$时可以通过连续求平方,以更少的步骤完成幂运算.例如,在朴素幂运算下$ b ^ 8 $:
$ b ^ 8 = b * ( b * ( b * ( b * ( b * ( b * ( b * b ) ) ) ) ) ) $
用连续求平方:
$ b ^ 2 = b * b $
$ b ^ 4 = b ^ 2 * b ^ 2 $
$ b ^ 8 = b ^ 4 * b ^ 4 $
连续求平方发的一般乘幂计算:
$ b ^ n = ( b ^ {n / 2} ) $, 若n是偶数
$ b ^ n = b * b ^ { n - 1 } $, 若n是奇数
将这个方法定义成以下过程:
(
define ( fast_expt b n )
(
cond ( ( = n 0 ) 1 )
( ( is_even n ) ( square ( fast_expt b ( / n 2 ) )) )
( else ( * b ( fast_expt b ( - n 1 ) ) ) )
)
)
(define ( square n )( * n n ))
(define ( is_even n )( = ( remainder n 2 ) 0))
时间复杂度为$\Theta ( log n )$
练习1.16 利用关系$ ( b ^ { n / 2 } ) ^ 2 = ( b ^ 2 ) ^ { n / 2 } $将上面的快速幂算法改成按照迭代方式的求幂计算.
(
define ( fast_expt b n )
( fast_expt_iter b n 1 )
)
(
define ( fast_expt_iter b counter product )
(
cond ( ( = counter 0 ) product )
( ( is_even counter ) ( fast_expt_iter ( square b ) ( / counter 2 ) product ))
( else ( fast_expt_iter b ( - counter 1 ) ( * product b ) ) )
)
)
除了基数b和指数n之外,还有一个附加的状态变量a(上面的代码是product),并定义好状态变换,使得从一个状态到另一个状态时乘积$ a b ^ n $不变.在计算过程开始时a取1,并用计算过程结束时a的值作为回答.一般说,定义一个不变量,要求它在状态之间保持不变,这一技术是思考迭代算法设计问题时的一种强有力的方法.
俄罗斯农夫算法
练习1.18 利用反复加法的方式求出乘积.
(
define ( * a b )
(
if ( = b 0 )
0
( + a ( * a ( - b 1 ) ) )
)
)
这是一个线性算法,在练习1.17中要求将其优化为对数复杂度的算法,并且使用double和halve两个方法;在练习1.18中要求再改为迭代计算过程.
练习1.17
( define (double x)( * x 2 ) )
( define ( halve y ) ( / y 2 ) )
( define (is_even n) ( = ( remainder n 2 ) 0 ))
(
define ( fast_mul1 a b )
(
cond ( ( = b 0 ) 0)
( ( is_even b ) ( double ( fast_mul1 a ( halve b ) ) ) )
( else ( + a (fast_mul1 a ( - b 1 ) ) ) )
)
)
练习1.18
(
define ( fast_mul a b )
( fast_mul_iter a b 0 )
)
;Russian farmer algorithm
(
define ( fast_mul_iter a b product )
(
cond ( ( = b 0 ) product)
( ( is_even b ) ( fast_mul_iter ( double a ) ( halve b ) product ) )
( else ( fast_mul_iter a ( - b 1 ) ( + product a ) ) )
)
)
这里的$product$就是附加的状态变量,也是最后要输出的答案.
有一个乘法表达式:$ans = a * b$
当b为奇数时 $ ans = a * b_1 + c_1, b_1 = ( b - 1 ), c_1 = a$并将$c_1$累加到$product$上
当b为偶数时$ ans = (a * 2) * ( b / 2 ) $, 这时b的状态就有可能会改变
当b为0时,就可以输出$product$了
这一算法又被称为”俄罗斯农夫算法”
斐波那契数列
朴素斐波那契函数
根据斐波那契数列的定义:
[latexpage]
$Fib( x )= \begin{cases} 0, &x = 0 \cr 1, & x =1 \cr Fib( x - 1 ) + Fib( x - 2 ) , & x > 1 \end{cases}$
将这个过程翻译为一个计算斐波那契数的递归过程
(
define ( fib n )
(
cond (( = n 0 ) 0 )
( ( = n 1 ) 1 )
( else ( + ( fib ( - n 1 ) )
( fib ( - n 2 ) ))
)
)
)
但是这个过程有太多的冗余计算,时间复杂度可以说是一个指数级别的.
$ Fib(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5} $
迭代的斐波那契函数
用一对整数a和b将它们分别初始化为$ Fib(1) = 1 $ 和$ Fib(0) = 0 $,然后反复同时使用下面的变换规则:
$ a \leftarrow a + b $
$ b \leftarrow a $
不难证明在$n$次这些变换后,a和b将分别等于Fib( n + 1 ) 和 Fib( n ), 将其翻译为下面过程:
( define ( fib2 n )
( fib_iter 1 0 n)
)
( define ( fib_iter a b count )
(if ( = count 0 )
b
( fib_iter ( + a b ) a ( - count 1 ) )
)
)
另外两种优化思路是(不再贴出代码):
1.使用记忆化搜索
2.递推方法,用数组记录fib[0] = 0, fib[1] = 1, 并递推fib[i] = fib[i - 1] + fib[i - 2].
[latexpage]
练习1.13 证明$Fib( n ) $是最接$\frac{\phi ^ n }{\sqrt 5}$的整数，其中$\phi = \frac{1 + \sqrt 5 }{2}$，证明$ Fib(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5} $
先证明$ Fib(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5} $
令$A(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5}$
$Fib(1) = 1$ $A(1) = 1$ $ Fib(1) = A(1) $
$Fib(2) = 1$ $A(1) = 1$ $ Fib(2) = A(2) $
设当$k = n $时有$ Fib(n) = A(n) $，则当$k = n + 1$时有
$Fib(n + 1) = Fib( n ) + Fib(n - 1)$
$Fib( n ) = A(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5}$
$Fib( n - 1 ) = A( n - 1 ) = \frac{ \phi ^ {n - 1} - \gamma ^ {n - 1}}{\sqrt 5}$
则有 $Fib( n ) + Fib(n - 1) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5} + \frac{ \phi ^ {n - 1} - \gamma ^ {n - 1}}{\sqrt 5} $ $= \frac{\phi ^ n( 1 + \frac{1}{\phi} ) - \gamma ^ n( 1 + \frac{1}{\gamma } )}{\sqrt 5}$
其中 $ 1 + \frac{1}{\phi ^ n} = \frac{3 + \sqrt 5}{1 + \sqrt 5} = \frac{(3 + \sqrt 5)(1 - \sqrt5)}{(1 + \sqrt5)( 1 - \sqrt 5 )} = \frac{1 + \sqrt 5}{2} = \phi$
$ 1 + \frac{1}{\gamma ^ n} = \frac{3 - \sqrt 5}{1 - \sqrt 5} = \frac{(3 - \sqrt 5)(1 + \sqrt5)}{(1 - \sqrt5)( 1 + \sqrt 5 )} = \frac{1 - \sqrt 5}{2} = \gamma$
则 $Fib( n ) + Fib(n - 1) = \frac{\phi ^ {n + 1} - \gamma ^ { n + 1 }}{\sqrt 5} = A(n+1)$
$ \therefore Fib(n + 1) = A(n+1)$
$ \therefore Fib(n) = \frac{ \phi ^ n - \gamma ^ n}{\sqrt 5} $
这个证明也可以用特征根方程来证明,这里不再赘述.
$| {Fib(n) - \frac{\phi}{\sqrt 5} | = | \frac{\gamma ^ n}{\sqrt 5} |$
其中 $| \gamma | = | \frac{1 - \sqrt 5}{2} | < 1$
$ | \gamma ^ 2 | < 1 $
$ | \frac{\gamma ^ 2}{\sqrt 5} | <, 0.5$
$ \therefore $ $ Fib( n )$ 是最接 $\frac{\phi ^ n }{\sqrt 5}$ 的整数
练习1.19
[latexpage]
主要是求出$p_1$和$q_1$ ,题目说的操作两次$T_{pq}$相当于操作一次$T_{p’q’}$
其中$T_{pq}$是对于对偶$(a,b)$按照如下的规则变化
$a \leftarrow bq + aq + ap$
$ b \leftarrow bp + aq $
对$T_{pq}$操作两次有
$a_1 = bq + a( p + q )$
$ b_1 = bp + aq $
$ a_2 = b_1 q + a_1 ( p + q ) = ( b1 + aq ) q + ( bq + a ( p + q )( p + q ) )$
$= b( 2pq + q^2 ) + q( p ^2 + 2 p ^2 + 2pq )$
其中$ q_1 = 2pq + q^2 $ $ p_1 = p ^ 2 + q ^2 $
$ b_2 = b_1 q + a_1 q = ( bp + aq ) p + ( bq + a( p + q ) ) q = b( p^2 + q^2 ) + a( 2pq + p^2 )$
这里$b_2$求出的$q_1$和$p_1$可以作为一个验证.
$\therefore $ $ q_1 = 2pq + q^2 $ $ p_1 = p ^ 2 + q ^2 $
然后贴出代码
(define (fib n )
( fib_iter 1 0 0 1 n )
)
(define (even? n ) ( = ( remainder n 2 ) 0 ) )
(define (square x ) ( * x x ) )
(define (fib_iter a b p q count )
(
cond( ( = count 0 ) b )
( ( even? count )
( fib_iter a
b
( + ( square p) (square q) )
( + ( * 2 p q ) ( square q ) )
( / count 2 )
)
)
(
else ( fib_iter ( + ( * b q ) ( * a q ) ( * a p ) )
( + ( * b p ) ( * a q ) )
p
q
( - count 1 )
)
)
)
)
最大公约数的Lame定理
如果欧几里得算法需要用k步计算出一对整数的GCD,那么这对数中较小的那个数必然大于或者等于第k个斐波那契数.
[latexpage]考虑数对序列$( a_k, b_k )$,其中$a_k \geq b_k$, 假设欧几里得算法在第k步结束.
如果$ ( a_{k+1}, b_{k+1} ) \rightarrow ( a_{k}, b_{k} ) \rightarrow ( a_{k-1}, b_{k-1} )$是归约序列中连续的三个数对,必然有$ b_{k + 1} \geq b_k + b_{k - 1} $.这里的每个归约步骤都是通过应用变换$ a_{k-1} = b_k $, $ b_{ k - 1 } = a_{k}$除以$ b_k $的余数.
$ a_k = q b_k + b{ k - 1 } $ 其中$ q \geq 1 $所以有 $ a_k = q b_k + b_{k - 1} \geq b_k + b_{ k - 1 }$,在前一个归约步骤中有$ b_{ k + 1 } = a_k $, 因此 $b_{k + 1} = a_k = q b_k + b_{ k - 1 } \geq b_k + b_{ k - 1 }$
要证明$ b_k \geq Fib(k) $, 一定有$ b_1 \geq Fib(1) = 1 $
设$ b_{ k } \geq Fib(k) $
$ b_{ k - 1 } \geq Fib(k - 1) $.
则 $b_{ k } +b_{ k- 1} \geq Fib(k) + Fib(k - 1)$
$ \because b_{k + 1} \geq b_k + b_{ k - 1 } \geq Fib(k) + Fib(k - 1) $
$ \therefore b_{k + 1} \geq Fib( k + 1 )$
$ \therefore b_k \geq Fib(k)$
素数检测
[latexpage]提供两种素数检测方法,第一种具有$\Theta ( \sqrt n )$的增长阶, 另一个是运用费马小定理的概率算法,具有$\Theta ( log n )$的增长阶.
传统素数检测方法
[latexpage]通过枚举2到$ \sqrt(n) $的所有整数来验证n是否为素数.时间复杂度为$\Theta ( \sqrt n )$
(define ( smallest_divisor n ) ( find_divisor n 2 ))
(define ( square x ) ( * x x ) )
( define ( find_divisor n test_divisor )
(cond ( ( > ( square test_divisor ) n ) n )
( ( divides? test_divisor n ) test_divisor )
( else ( find_divisor n ( + test_divisor 1 ) ) )
)
)
(define ( divides? a b ) ( = ( remainder b a ) 0 ))
(define ( prime? n )
( = n ( smallest_divisor n ) )
)
费马检查
[latexpage]费马小定理:如果n是一个素数,a是小于n的任意正整数,那么a的n次方与a模n同余.(更多详细内容可见费马小定理)对于给定的整数n, 随机任取一个$a < n$ 计算出$ a^n $取模n的余数.如果结果不等于a,那么n就肯定不是素数;如果结果等于a,那么n就有很大可能是素数.检查的a越多,那么n是素数的可能性也就越大,也存在一些数据能够欺骗费马检查,所以综上这是一个概率算法.
要完成费马检查需要先写一个快速幂取模的过程.
快速幂取模
快速幂取模基于以下定理
[latexpage]$(a * b) \% c= ( a \% c ) * ( b \% c ) \% c $
[latexpage]可以将一个大数分成几部分来求模以达到快速运算的目的,这个过程类似之前的fast-exp过程.其增长的阶是$ \Theta( log n ) $.
费马检查的代码
(define (square x ) ( * x x ) )
(define (expmod base exp m)
(
cond ( ( = exp 0 ) 1 )
( (even? exp)
( remainder ( square ( expmod base ( / exp 2 ) m ) ) m ))
(
else ( remainder ( * base ( expmod base ( - exp 1 ) m ) ) m )
)
)
)
(define (fermat_test n )
(define (try_it a )
( = ( expmod a n n ) a )
)
( try_it ( + 1 ( random( - n 1 ) ) ) )
)
(define (fast_prime? n times)
(
cond( ( = times 0 ) true )
( ( fermat_test n ) ( fast_prime? n ( - times 1 ) ) )
( else false )
)
)
能够欺骗费马检查的数成为_Carmichael_数,其中最小的只有561,在练习1.28中提供了一种优化过的不会被欺骗的费马检查.以下是对561进行费马检查的中间生成数据,看看它是如何欺骗费马检查的.
exp = 1 * = 3 mod = 3
exp = 2 squared = 9 mod = 9
exp = 4 squared = 81 mod = 81
exp = 8 squared = 6561 mod = 390
exp = 16 squared = 152100 mod = 69
exp = 17 * = 207 mod = 207
exp = 34 squared = 42849 mod = 213
exp = 35 * = 639 mod = 78
exp = 70 squared = 6084 mod = 474
exp = 140 squared = 224676 mod = 276
exp = 280 squared = 76176 mod = 441
exp = 560 squared = 194481 mod = 375
exp = 561 * = 1125 mod = 3
Miller-Rabin检查
[latexpage]费马检查的一种不会被欺骗的变形成为Miller-Rabin检查,来源与费马小定理的一个变形. 如果n是素数,a是小于n的整数,则a的( n - 1 )次幂与1模n同余.这个算法需要查看是否遇到了”1取模的非平凡平方根”,也就是说,是不是存在不等于1或者n-1的数使平方取模等于1.如果1取模的非平凡平方根存在,那么n就不是素数.并且,如果n是非素数的奇数,至少有一半的数a < n 按照这种方式计算$a^{ n - 1 }$,将会遇到1取模n的非平凡平方根.
[latexpage] Miller-Rabin检查基于的定理: 如果p是素数，x是小于p的正整数，且$x^2$ mod p = 1，那么要么x=1，要么x=p-1。这是显然的，因为$x^2$ mod p = 1相当于p能整除$x^2$-1，也即p能整除(x+1)(x-1)。由于p是素数，那么只可能是x-1能被p整除(此时x=1)或x+1能被p整除(此时x=p-1)。
练习1.28 编写Miller-Rabin检查
(define (square x ) ( * x x ) )
(define (Miller_Rabin base exp m )
(
cond( ( = exp 0 ) 1)
( ( square_check base m ) 0)
( ( even? exp ) ( remainder ( square ( Miller_Rabin base ( / exp 2 ) m ) ) m ) )
(else
(remainder ( Miller_Rabin base ( - exp 1 ) m ) m)
)
)
)
(define (square_check a n )
(
and ( not ( = a 1 ) )
( not ( = a ( - n 1 ) ) )
( = 1 ( remainder ( square a ) n ) )
)
)
(define (Miller_Rabin_test n )
(define (try_it a )
( = ( Miller_Rabin a n n ) a )
)
( try_it ( + 1 ( random( - n 1 ) ) ) )
)
(define (fast_prime? n times)
(
cond( ( = times 0 ) true )
( ( Miller_Rabin_test n ) ( fast_prime? n ( - times 1 ) ) )
( else false )
)
)
Miller_Rabin方法参考的文章:http://www.matrix67.com/blog/archives/234