費式數列O(logN)

Created at 2019-02-23 Updated at 2019-02-28 Category 學習紀錄 Tag 學習

費式數列$O(logN)$解法(使用python)

前言

一般面對費式數列的題目時,最直觀的想法就是使用遞迴,複雜度$O(2^n)$

1
2
3
4
5
def fibo_recursive(num):
if num < 2:
return num
else:
return fibo_recursive(num-1) + fibo_recursive(num-2)

這東西光num=40就跑得有點久了,於是有一種$O(N)$的解法

1
2
3
4
5
6
7
8
9
10
11
def fibo_for_loop(num):
a = 0
b = 1
if num < 2:
return num
else:
for i in range (1, num):
temp = b
b = a + b
a = temp
return b

這方法已經滿快了但是…

以下就來推導$O(logN)$

推導

我們都知道費式數列的規則:

其中

我們現在有個問題是,我們能否推出函數$F(N) = ?$

我們來一步一步地慢慢推。

首先我們把

寫成

為甚麼要多寫個呢?
因為這樣我們就能以矩陣來表示這個式子

為了接下來方便運算我們加上

當中的

是為了使等號成立。
來到這裡我們試著帶入數字看看

這裡運算可以知道等試是成立的

再繼續帶

我們可以發現我們上面已經有

所以我們把這條式子帶入將會變成

再帶下個則會變成

以此類推我們能得到

有沒有感受到好像發現了甚麼東西?我們成功的推導到這裡,也就是說只要把$N$帶入就能算出答案了。但問題又來了,如果$N=50$,矩陣50次方怎麼算?接下來就要搬出線性代數所教的”對角化”,至於對角化詳細的內容就請自行去翻課本或查資料,這裡就不多講了。簡單來說就是一個矩陣$A$能夠被寫成

其中

我們利用對角化求出我們的$P$與$D$分別是

接著將其帶入

由於我們只需要找$F(N)$所以只需要算

經過整理會變成

最後將$\lambda_1\lambda_2$帶入就會得到

這樣就成功推導出來了。
至於為甚麼會說這個方法是$O(logN)$呢?
因為以上的式子所要計算的是次方那個部分,而次方的問題可以用快速冪或矩陣快速冪去處理,因此時間複雜度是$O(logN)$

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
def fast_pow(a, exp):

base = a
if exp == 0:
return 0
else:
result = 1
while exp:
if exp & 1:
result *= base
base *= base
exp >>= 1
return result
def fibo_formula(num):
return (math.sqrt(5)/5)*(fast_pow((1+math.sqrt(5))/2, num) - fast_pow((1-math.sqrt(5))/2, num))

然而如果只接用推導出來的公式用快速冪去算,因為有$\sqrt5$的關係,$F(62)$開始會有誤差,因此必須使用矩陣快速冪的方法。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
def matrixMultiply(matrix, base):
row = len(matrix)
col = len(base[0])
matrix_ans = [[0 for i in range (col)] for i in range (row)]

for i in range (0, row):
for j in range (0, col):
for k in range (0, row):
matrix_ans[i][j] = matrix_ans[i][j] + matrix[i][k] * base[k][j]

return matrix_ans

def matrixPow(A, exp):
res_matrix = [[1, 0], [0, 1]]
while exp:
if exp & 1:
res_matrix = matrixMultiply(res_matrix, A)
A = matrixMultiply(A, A)
exp >>= 1
return res_matrix

def get_fibo_number(num):
if num < 2:
return num
else:
A = [[1, 1], [1, 0]]
base = [[1], [0]]
ans = matrixMultiply(matrixPow(A, num-1), base)
return ans[0][0]

Site by ICEJJ using Hexo & Random

Hide