M斐波那契數列
- 2020 年 2 月 18 日
- 筆記
Problem Description
M斐波那契數列F[n]是一種整數數列，它的定義如下：
F[0] = a F1 = b F[n] = F[n-1] * F[n-2] ( n > 1 )
現在給出a, b, n，你能求出F[n]的值嗎？
Input
輸入包含多組測試數據； 每組數據佔一行，包含3個整數a, b, n（ 0 <= a, b, n <= 10^9 ）
Output
對每組測試數據請輸出一個整數F[n]，由於F[n]可能很大，你只需輸出F[n]對1000000007取模後的值即可，每組數據輸出一行。
Sample Input
0 1 0 6 10 2
Sample Output
0 60
題目解決
題目分析:
- 注意到n的範圍最大到10^9,所以簡單的進行遞推時間和空間上都無法處理到限制範圍內.
- 這道題涉及到取余運算,所以我們應該知道取余運算的其中一個性質:$$(ab) mod p = ((a mod p) (b mod p)) mod p$$
- 注意到該遞推式和斐波那契數列有相似之處.
知識點:
1.斐波那契數列的矩陣運算
1.斐波那契數列
$fib(0) = 0$ $fib(1) = 1$ $fib(n) = fib(n-1)+fib(n-2)$ | n | 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 | 13 | 14 | … | … | n-1 | n | | —–: | —–: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | | fib(n) | 0 | 1 | 1 | 2 | 3 | 5 | 8 | 13 | 21 | 34 | 55 | 89 | 144 | 233 | 377 | … | … | fib(n-1) | fib(n) |
2.將斐波那契數列寫到矩陣中
1.考慮按照以下表達式構造n個矩陣:
$$ A_{n} = begin{bmatrix} fib(n) & fib(n+1) fib(n+1) & fib(n+2) end{bmatrix} $$
2.得到如下矩陣
$A_{0} = begin{bmatrix}0 & 1 1 & 1 end{bmatrix} $ $A_{1} = begin{bmatrix}1 & 1 1 & 2 end{bmatrix} $ $ A_{2} = begin{bmatrix}1 & 2 2 & 3 end{bmatrix} $ $A_{3} = begin{bmatrix}2 & 3 3 & 5 end{bmatrix} $ $A_{4} = begin{bmatrix}3 & 5 5 & 8 end{bmatrix} $ … … …
$A_{n} = begin{bmatrix}fib(n) & fib(n+1) fib(n+1) & fib(n+2) end{bmatrix} $
3.矩陣乘法運算
$A_{1} = begin{bmatrix}1 & 1 1 & 2 end{bmatrix} = A_{0}A_{0}$ $ A_{2} = begin{bmatrix}1 & 2 2 & 3 end{bmatrix} = A_{0}A_{0}A_{0}$ $A_{3} = begin{bmatrix}2 & 3 3 & 5 end{bmatrix} = A_{0}A_{0}A_{0}A_{0} $
$ A_{4} = begin{bmatrix}3 & 5 5 & 8 end{bmatrix} = A_{0}A_{0}A_{0}A_{0}A_{0}$ $A_{n} = begin{bmatrix}fib(n) & fib(n+1) fib(n+1) & fib(n+2) end{bmatrix} = A_{0} ^ {n+1}$
begin{array}{|rrrrrrrr|} hline color{red}{綜上所述,要想求解f(n),求得A(n)或者A(n-2)即可} & hline end{array}
2.快速冪:
Q: 求解 $ x ^ N $ A: O(N)的演算法:
ans = 1; for(i = 1; i <= n; i++){ ans *= x; }
A: O(logN)的演算法:
- 將N轉換成二進位表示形式 $$ N = 2 ^ {k_{1}} + 2 ^ {k_{2}} + 2 ^ {k_{3}} + …. + 2 ^ {k_{n}}$$
- 要求解$x ^ N $,相當與計算$$ x ^ N = x ^ {2 ^ {k_{1}} + 2 ^ {k_{2}} + 2 ^ {k_{3}} + …. + 2 ^ {k_{n}}} = x^{2 ^ {k_{1}}}x^{2 ^ {k_{2}}}x^{2 ^ {k_{3}}}…x^{2 ^ {k_{n}}}$$
- 所以只要在依次求$x ^ {2^i} $ 的同時進行計算就好了,最終得到了O(logN)的計算冪運算的演算法
舉個例子:計算 $ x ^ {22} $
- 十進位下的22的二進位表示 : $ 22_{10} = 10110_{2}$
- $ N = 22 = 2^{4} + 2^{2} + 2^{1}$
- $ x^{22} = x^{2^{4} + 2^{2} + 2^{1}} = x^{16}x^{4}x^{2}$
- 可以看到如果是O(N)的演算法,當前例子需要計算22次,而O(logN)的演算法只需要計算5次
begin{array}{|rrrrrrrr|} hline & color{red}{CODE:}& hline end{array}
typedef long long LL; // 計算x^n,複雜度O(logN) LL powN(LL x, LL n){ LL r = 1; // 計算結果 while(n){ if(n&1) r *= x; // 如果n的二進位表示i(i >= 0)的位置是1,則乘上x^(2^i) x *= x; // 將x平方 n >>= 1; // 右移一位 } return r; // 返回結果 } // 計算x^n mod p, ( a * b ) % p = ( ( a % p ) * ( b * p ) ) % p 複雜度O(logN) LL mod_powN(LL x, LL n, LL p){ LL r = 1; // 計算結果 while(n){ if(n&1) r = r * x % p; // 如果n的二進位表示i(i >= 0)的位置是1,則乘上x^(2^i) x = x * x % p; // 將x平方 n >>= 1; // 右移一位 } return r; // 返回結果 }
3.費馬小定理:
費馬小定理是數論中的一個定理：假如a是一個整數，p是一個質數，那麼 $a^{p}-a$是p的倍數，可以表示為 $$ a ^ p equiv a (mod p) $$ 如果a不是p的倍數，這個定理也可以寫成 $$ a ^ {p-1} equiv 1 (mod p)$$ 費馬小定理是歐拉定理的一個特殊情況：如果n和a的最大公因數是1，那麼 $$ a^{φ(n)} equiv 1 (mod n) $$ 這裡φ(n)是歐拉函數。歐拉函數的值是所有小於或等於n的正整數中與n互質的數的個數。假如n是一個素數（質數），則φ(n) = n-1，即費馬小定理 注: $equiv$是同餘符號 $ a mod p = b mod p$ 可表示為 $ a equiv b (mod p)$
推導以下表達式: 當p是質數時: $$ a ^ n equiv a ^ {(n mod (p-1))} (mod p) $$ 成立
證: $$ n = k(p-1) + (n mod (p-1)) (k為整數) $$ $$ a ^ n equiv a ^ {k(p-1) + (n mod (p-1))} equiv {a ^ {p-1}} ^ {k} a ^ {(n mod (p-1))} (mod p)$$ 由費馬小定理可知: $$ {a ^ {p-1}} ^ {k} equiv a ^ {p-1} equiv 1(mod p)$$ 得證: $$ a ^ n equiv a ^ {(n mod (p-1))} (mod p) $$
回到題目:
1.分析
我們注意到雖然題目中的M斐波那契數列和斐波那契數列有異曲同工之處,但是一個計算的是乘法,一個是加法,但加法和乘法之間有著密切的聯繫.
2.算一算
先計算前幾項的值: $$ f(0) = a = a ^ 1 b ^ 0 $$ $$ f(1) = b = a ^ 0 b ^ 1 $$ $$ f(2) = a b = a ^ 1 b ^ 1 $$ $$ f(3) = a b b = a ^ 1 b ^ 2 $$ $$ f(4) = a a b b b = a ^ 2 b ^ 3 $$ $$ f(4) = a a a b b b b b= a ^ 3 b ^ 5 $$ $$ ……………………………………………………..$$ $$ f(n) = f(n-1)f(n-2) = a ^ x b ^ y $$ |n| 0 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 11 | 12 | 13 | … | … | n-1 | n | |——–| —–: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | :—-: | |x| 1 | 0 | 1 | 1 | 2 | 3 | 5 | 8 | 13 | 21 | 34 | 55 | 89 | … | … | fib(n-2) | fib(n-1) | |y| 0 | 1 | 1 | 2 | 3 | 5 | 8 | 13 | 21 | 34 | 55 | 89 | 144 | … | … | fib(n-1) | fib(n) |
不難得出: $$ f(0) = a $$ $$ f(1) = b $$ $$ f(n) = f(n-1)f(n-2) = a ^ {fib(n-1)} b ^ {fib(n)} (n >= 2) $$
通過上面費馬小定理證明的結論,當p是質數時: $$ a ^ n equiv a ^ {(n mod (p-1))} (mod p) $$ 所以我們實際上需要計算的就是下面的表達式(p = 1000000007): $$ f(n) mod p= f(n-1)f(n-2) mod p = a ^ {fib(n-1) mod p-1} b ^ {fib(n) mod p-1} mod p(n >= 2) $$
1.計算類似於計算整數的快速冪演算法計算矩陣的冪 $A_{n-1} = A_{0} ^ {n}$ $A_{0} = begin{bmatrix}0 & 1 1 & 1 end{bmatrix}$ $A_{n-1} = begin{bmatrix}fib(n-1) & fib(n) fib(n) & fib(n+1) end{bmatrix} = A_{0} ^ {n}$
計算矩陣的冪時,我們實際上需要得到的結果是: $$ A_{n-1} mod p-1 = begin{bmatrix}fib(n-1) mod p-1 & fib(n) mod p-1 fib(n) mod p-1 & fib(n+1) mod p-1 end{bmatrix} = (A_{0} mod p-1) ^ {n} mod p-1 $$
如果理解了上面的快速冪演算法的話,這實現起來將會非常簡單
2.得到$fib(n-1) mod p-1$ 和 $fib(n) mod p-1$ 之後,設這兩個數為x,y,用整數的快速冪演算法計算: $$a ^ x b^y mod p = (a ^ x mod p b^y mod p) mod p$$
begin{array}{|rrrrrrrr|} hline & color{red}{CODE:}& hline end{array}
#include <bits/stdc++.h> using namespace std; typedef long long LL; const int MOD1 = 1e9+7; // p const int MOD2 = 1e9+6; // p-1 LL a,b,n; typedef struct Matrix{ // 矩陣結構體 LL arr[2][2]; }Matrix; Matrix unit = { // 單位矩陣 1,0, 0,1 }; Matrix A0 = { // A0矩陣 0,1, 1,1 }; // 矩陣a*b Matrix matrixMulti(Matrix a, Matrix b){ int i,j,k; Matrix tmp; for(i = 0; i < 2; i++){ for(j = 0; j < 2; j++){ tmp.arr[i][j] = 0; for(k = 0; k < 2; k++){ tmp.arr[i][j] += (a.arr[i][k]*b.arr[k][j])%MOD2; tmp.arr[i][j] %= MOD2; } } } return tmp; } // 矩陣快速冪,計算的結果對MOD2取余 Matrix matrixPow(Matrix a, LL n){ Matrix r = unit, base = a; while(n){ if(n&1) r = matrixMulti(r,base); base = matrixMulti(base, base); n >>= 1; } return r; } // 整數快速冪,計算的結果對MOD1取余 LL powN(LL a, LL n){ LL r = 1, base = a%MOD1; while(n){ if(n&1) r = r*base%MOD1; base = base*base%MOD1; n >>= 1; } return r; } int main(){ //freopen("in.txt", "r", stdin); Matrix An_1; while(~scanf("%lld%lld%lld", &a, &b, &n)){ // 輸入a,b,n An_1 = matrixPow(A0,n); // 計算An_1 printf("%lldn", powN(a, An_1.arr[0][0]) * powN(b, An_1.arr[1][0]) % MOD1); } return 0; }
思考:
Q: 為什麼要使用費馬小定理? A: 因為 $ x mod p-1 <= x $ 所以可以減少計算整數冪時候的計算次數,不使用費馬小定理會超時.
