分享

黄金数列的急速运算

 imelee 2014-01-18

算法的时间复杂度

Edit:
用了代码高亮脚本以后超过字数上限
高亮后的代码会放在回复里……

========

首先谢谢小年给开了日志~~


然后抄送某个多线程反而比单线程慢的小破孩。。。


看了下发现代码缩进神马的都不见了。。。。纳斯达克老师新的脚本里貌似没有更新代码支持部分呢……

最后,第一次写科普文,求各种批评指点,尤其是在非技术的方面。。

谢谢啦~~
=================


在信息化的时代里每个人都应该懂一点编程。重复性的机械劳动,都应该交给机器去做,从而把人的精力放在更有挑战性的任务上。写代码第一个要避免的就是syntex error和类似这里(http://www./post/65643/)所写非所想的错误;第二个就是要注意算法的复杂度。一个程序没有bug,但是因为算法不够优化,最后要很长时间才能得出结果,发生这种事呢,大家都不想的。这里就从斐波那契数列入手,来讲讲算法的时间复杂度。


0. 准备工作: Landau符号

Landau符号其实是由德国数论学家保罗·巴赫曼(Paul Bachmann)在其1892年的著作《解析数论》首先引入,由另一位德国数论学家艾德蒙·朗道(Edmund Landau)推广。Landau符号的作用在于用简单的函数来描述复杂函数行为,给出一个上或下(确)界。在计算算法复杂度时一般只用到大O符号,Landau符号体系中的小o符号、Θ符号等等比较不常用。这里的O,最初是用大写希腊字母,但现在都用大写英语字母O;小o符号也是用小写英语字母o,Θ符号则维持大写希腊字母Θ。

f (n) = Ο(g (n)) 表示存在一个常数C,使得在当n趋于正无穷时总有 f (n) ≤ C * g(n)。简单来说,就是f(n)在n趋于正无穷时最大也就跟g(n)差不多大。虽然对g(n)没有规定,但是一般都是取尽可能简单的函数。例如,O(2n^2+n +1) = O (3n^2+n+3) = O (7n^2 + n) = O ( n^2 ) ,一般都只用O(n^2)表示就可以了。注意到大O符号里隐藏着一个常数C,所以g(n)里一般不加系数。如果把f(n)当做一棵树,那么O(g(n))所表达的就是树干,只关心其中的主干,其他的细枝末节全都抛弃不管。

在时间复杂度计算中常见的级别有O(1), O(log n) , O(n), O(n * log n) , O(n^k), O(a^n), O(n!),其大小逐级上升:



注意到所有的对数只不过相差一个常数,所以这里都用了常用对数。另外一般程序只处理32位的数据,因此最大整数是2^32-1,大约等于10^9。因此log n可以认为是近似等于10的常数,也就是说O(log n)的近似等于O(1),O(n * log n)近似等于O(n),这点也在上表中有所反应。

全国有10^9个人,一个人大约有10^14个细胞,一摩尔基本粒子中有10^24个(大致相当于全国人民身上的细胞总数),宇宙中大概有10^80个基本粒子,所以上表中的某些数字,你懂的……

1 算法的时间复杂度

好了,了解了Landou或者叫大O符号,我们就可以来讨论算法的时间复杂度了。一个算法,对于一个规模是n的输入,需要T(n)的时间进行运算。T(n)的具体值,一方面取决于算法,另一方面和计算用的工具也有关。但我们关心的是,当n扩大100倍以后,T(n)是几乎不变,还是扩大了100倍、10000倍或者怎样。也就是说,我们关心的是规模扩大后,运算时间增长得有多快,是一个相对的结果,而不是绝对结果。

建立一个空列表,需要的时间总是固定的,就是O(1);如果有一句for或者while的语句,把同一行代码运算了n次,那么就是O(n)级的复杂度;如果for下面还有for,各自要算n次,那么就是n^2的复杂度,以此类推。

复杂度举例:

* O(1) 常数级复杂度,也就是说程序运行的时间与需要处理的数据大小无关。通常把比较大小、加减乘除等简单的运算都当做常数级复杂度。 值得注意的是,在处理大数(二进制下数据长度超过32位或者十进制下超过8位)时,将加减乘除等运算当做常数复杂度不再适用。

* O(log n) 将一个10进制整数转化为2进制整数

* O(n):判断一个元素是否属于一个大小为n的集合/列表;找出n个数中的最大值;

* O(n * log n) 快速排序法

* O(n^2) 最直白的两两比较然后排序方法,需要n*(n-1)/2次比较,所以是O(n^2)级。

* O(2^n) 列出所有长度为n的0,1字符串之类的穷举计算

* O(n!) 列出n个元素的全排列之类的穷举计算

一般来说多项式级的复杂度是可以接受的,很多问题都有多项式级的解——也就是说,这样的问题,对于一个规模是n的输入,在n^k的时间内得到结果,称为P问题。有些问题要复杂些,没有多项式时间的解,但是可以在多项式时间里验证某个猜测是不是正确。比如问4294967297是不是质数?如果要直接入手的话,那么要把小于4294967297的平方根的所有素数都拿出来,看看能不能整除。还好欧拉告诉我们,这个数等于641和6700417的乘积,不是素数,很好验证的,顺便麻烦转告费马他的猜想不成立。大数分解、Hamilton回路之类的问题,都是可以多项式时间内验证一个“解”是否正确,这类问题叫做NP问题。(据信程序猿找妹子也是一个NP问题。。)

P问题显然都是NP问题——既然你能在多项式时间里得到答案,那么只要比较一下答案和猜测是不是一样就可以验证一个解。反过来,多数人相信,NP问题不都是P问题。几个月前看到一篇论文说NP不等于P,但没有看到后续报道,也读不懂原文,摊手……如果NP问题都是P问题,那么大数分解不再是难题,从而基于大数分解的RSA加密系统顿时崩溃;同样md5算法也不再有效,甚至现在所有的多项式复杂度的加密算法都没用了——不过这样一来程序猿也可以方便地找到妹子了,说起来并不是一无是处嘛……

2 斐波那契数列

早在1150年印度科学家Gopala就研究过斐波那契数列,但正如Landau抢走了Landau符号的命名一样,斐波那契也获得了不属于他的数列的冠名权。斐波那契数列是和兔子紧紧联系在一起的。 斐波那契的兔子是这样的一个模型:一开始有一对小兔子。小兔子需要一个月以后才能长成大兔子从而拥有生育能力;所有拥有生育能力的大兔子每个月都会生一对小兔子。所有的兔子都长生不老,问第n个月一共有多少兔子呢?

维基百科词条中有这样一张直观的图:



让我们列一个表来看看每个月大兔子小兔子各有多少只:



从表中我们可以找到这样的规律:

F(n)=F(n-1)+F(n-2)

这就是斐波那契数列的递推公式,可以用数学归纳法证明。

利用一些数学工具我们可以得到斐波那契数列的通项公式:



其中小括号里的两个数是方程 x^2-x-1=0的解,跟著名的黄金分割比有密切的关系。

用上面的Landau符号,有F(n) = O(^n),其中是黄金分割比:

 

所以说,兔子的增长速度还是非常快的……




斐波那契数列还可以通过下面的矩阵乘方来实现,这点会在后面的算法中用到。对矩阵不熟悉的读者,暂时可以认为这里2阶矩阵的乘法是两组4元数组之间的运算,结果仍然是一组4元数组;对于相同的矩阵(4元数组),这个运算满足交换律和结合律。

 



3 三种程序计算斐波那契数列第n项

第一种算法,too simple, sometimes naif,直接从递推公式下手:

def fibo_1(n):
'''
简单递归方法求斐波那契数列第n项。
输入变量n应当为正整数。(对非法输入没有警告)
'''
if n < 3:
res=1 #设定初始状态:F(1)=F(2)=1
else:
res=fibo_1(n-1)+fibo_1(n-2) #递归公式
return res

这样计算的复杂度是怎样能?让我们看看递归次数:



可以看到计算F(n),需要的递归次数是F(n-2)。

指数级的复杂度啊……

不行不行,我要提高自身修养,看看第二个算法

def fibo_2(n):
'''
通过递归方法求斐波那契数列第n项
变量n应当为正整数。(对非法输入没有警告)
记录最后两项计算结果
'''
a=1
b=1 #设定初始状态,a=F(1),b=F(2)
for i in range(2,n):
a,b=a+b,a # 由数学归纳法易得第i步中a=F(i),b=F(i-1)
return a #循环结束时i=n,a=F(n)

分析:这次好多了,for语句,循环n-2次,O(n)级的复杂度。 

最后来个利用矩阵乘方的算法,只需要O(log n)的复杂度,可以跟美国的那个谁谈笑风生了:

def fibo_3(n):
'''
将求斐波那契数列转化为求2阶矩阵n次幂问题。
'''
def prod(m1,m2):
'''
两个矩阵的乘积
'''
a1,b1,c1,d1=m1
a2,b2,c2,d2=m2
return (a1*a2+b1*c2,a1*b2+c1*d2,a2*c1+c2*d1,b2*c1+d1*d2)
def pui(m,n):
'''
求一个矩阵m的n次幂。
如n=2k+1,则m^n=(m^k)^2*m;
如n=2k,则m^n=(m^k)^2。
'''
if n == 1:
return m
else:
res=pui(m,n//2) #上述的k=n//2,//表示带余除法的商,例如 5//2 = 2
res=prod(res,res) #得到上述的(m^k)^2
if n%2==1: # n%2表示n被2除的余数
res=prod(res,m) #如果n是奇数,那么结果还要再乘以m
return res
m=(1,1,1,0) #初始化矩阵
(a,b,c,d)=pui(m,n-1) # n次幂的第一项是F(n+1)
return a


这里用到的乘方算法,对于奇数2k+1,只要先算k次幂,然后自乘,再乘以最初的矩阵;对于偶数2k,只需先算k次幂,然后自乘即可。而计算k次幂时仍然可以采用上面的算法。如果n在2进制表达中一共有m位,并且其中有l个1,那么一共需要进行m+l次乘法,小于2*m=O(log n)。

多说无益,让我们来看看三种算法处理相同规模的输入各自要多少时间:

def timing(f,n):
'''计算f(n)需要多少时间'''
import time #载入时间模块
t1=time.time() #读取当前时刻(计算开始时刻)
res=f(n) #计算f(n)
t2=time.time() #再次读取当前时刻(计算结束时刻)
return t2-t1 #两次时刻差值既是计算f(n)所需要的时间,单位为秒



第一种天然呆的算法果然弱爆了,算F(30)就要0.6秒;第二种算法算十万要3秒多;第三种算法算一百万也不过3秒不到。上面的算法有一个问题,就是明明第二种算法是线性的,第三种算法是对数级的,但是当我尝试多加一个0时,所需要的时间完全不能忍受,不得不强制停止,大大超过按照前几次测试推算的结果。为什么呢?答案就在于,斐波那契的兔子繁殖特别快,中间过程的数值很快超过了2^32,python不得不启用了大数运算。而数字越大,进行乘法、加法运算所需要的时间也越长,不能再当做常数时间。那么让我们来看看下面这个简化,只求计算结果的末7位,也就是说利用求余数的方法把所有的数字都限制在10^8以内:

def fibo_1_bis(n):
'''
简单递归方法求斐波那契数列第n项。
输入变量n应当为正整数。(对非法输入没有警告)
'''
if n < 3:
res=1 #设定初始状态:F(1)=F(2)=1
else:
res=(fibo_1(n-1)+fibo_1(n-2))%100000000 #递归公式
return res


def fibo_2_bis(n):
'''
通过递归方法求斐波那契数列第n项
变量n应当为正整数。(对非法输入没有警告)
记录最后两项计算结果
'''
a=1
b=1 #设定初始状态,a=F(1),b=F(2)
for i in range(2,n):
a,b=(a+b)%100000000,a # 由数学归纳法易得第i步中a=F(i),b=F(i-1)
return a #循环结束时i=n,a=F(n)


def fibo_3_bis(n):
'''
将求斐波那契数列转化为求2阶矩阵n次幂问题。
'''
def prod(m1,m2):
'''
两个矩阵的乘积
'''
a1,b1,c1,d1=m1
a2,b2,c2,d2=m2
return ((a1*a2+b1*c2)%100000000,(a1*b2+c1*d2)%100000000,(a2*c1+c2*d1)%100000000,(b2*c1+d1*d2)%100000000)
def pui(m,n):
'''
求一个矩阵m的n次幂。
如n=2k+1,则m^n=(m^k)^2*m;
如n=2k,则m^n=(m^k)^2。
'''
if n == 1:
return m
else:
res=pui(m,n//2) #上述的k=n//2,//表示带余除法的商,例如 5//2 = 2
res=prod(res,res) #得到上述的(m^k)^2
if n%2==1: # n%2表示n被2除的余数
res=prod(res,m) #如果n是奇数,那么结果还要再乘以m
return res
m=(1,1,1,0) #初始化矩阵
(a,b,c,d)=pui(m,n-1) # n次幂的第一项是F(n+1)
return a



再来看看需要多少时间:





可以看到第一种算法需要的时间基本满足斐波那契的递推公式,和前面的预测相符;第二种算法需要的时间也和预测的线性增长相符;第三种算法,不管我输入多少个0,都是瞬间给出答案,O(log n)级的算法真心屌爆了有没有!!!


附python文件fibo.py下载地址:

http:///file/aqyjof34

参考了网上的思路,写了个Java版的: 

Java代码  收藏代码
  1. public class Fibonacci {  
  2.   
  3.     final  static int[] A={1,1,1,0};  
  4.       
  5.       
  6.     public static void main(String[] args) {  
  7.         int n=7;  
  8.         for(int i=0;i<=n;i++){  
  9.             int f=fibonacci(i);  
  10.             System.out.print(f+" ");  
  11.         }  
  12.     }  
  13.   
  14.     static int fibonacci(int n){  
  15.         if(n==0)return 0;  
  16.         if(n==1)return 1;  
  17.         if(n>1){  
  18.             int[] re=power(n-1);  
  19.             return re[0];//矩阵乘积的第00项为所求  
  20.         }  
  21.         return -1;  
  22.     }  
  23.       
  24.     //a^n=a^(n/2)*a^(n/2)   or   a^n=a^(n/2)*a^(n/2)*a  
  25.     static int[] power(int n){  
  26.         int[] a=new int[4];  
  27.         if(n==1){  
  28.             a=A;  
  29.         }  
  30.         else if(n%2==0){  
  31.             a=matixMultiply(power(n/2),power(n/2));  
  32.         }else if(n%2==1){  
  33.             int[] temp=matixMultiply(power(n/2),power(n/2));  
  34.             a=matixMultiply(A,temp);  
  35.         }  
  36.         return a;  
  37.     }  
  38.       
  39.     //矩阵乘法  
  40.     // return A*B  
  41.     static int[] matixMultiply(int[] a,int [] b){  
  42.         int[] re=new int[4];  
  43.         re[0]=a[0]*b[0]+a[1]*b[2];  
  44.         re[1]=a[0]*b[1]+a[1]*b[3];  
  45.         re[2]=a[2]*b[0]+a[3]*b[2];  
  46.         re[3]=a[2]*b[1]+a[3]*b[3];  
  47.         return re;  
  48.     }  
  49. }  

    本站是提供个人知识管理的网络存储空间,所有内容均由用户发布,不代表本站观点。请注意甄别内容中的联系方式、诱导购买等信息,谨防诈骗。如发现有害或侵权内容,请点击一键举报。
    转藏 分享 献花(0

    0条评论

    发表

    请遵守用户 评论公约

    类似文章 更多