M斐波那契数列

发表于 | 更新于 | 分类于 | 评论数： | 阅读次数：

题目传送门

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

题目解决

题目分析:

  1. 注意到n的范围最大到10^9,所以简单的进行递推时间和空间上都无法处理到限制范围内.
  2. 这道题涉及到取余运算,所以我们应该知道取余运算的其中一个性质:$$(ab) \mod p = ((a \mod p) (b \mod p)) \mod p$$
  3. 注意到该递推式和斐波那契数列有相似之处.

知识点:

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)的算法:

1
2
3
4
ans = 1;    
for(i = 1; i <= n; i++){
    ans *= x;
}

A:   O(logN)的算法:

  1. 将N转换成二进制表示形式 $$ N = 2 ^ {k_{1}} + 2 ^ {k_{2}} + 2 ^ {k_{3}} + …. + 2 ^ {k_{n}}$$
  2. 要求解$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}}}$$
  3. 所以只要在依次求$x ^ {2^i} $ 的同时进行计算就好了,最终得到了O(logN)的计算幂运算的算法

举个例子:计算 $ x ^ {22} $

  1. 十进制下的22的二进制表示 : $ 22_{10} = 10110_{2}$
  2. $ N = 22 = 2^{4} + 2^{2} + 2^{1}$
  3. $ x^{22} = x^{2^{4} + 2^{2} + 2^{1}} = x^{16}x^{4}x^{2}$
  4. 可以看到如果是O(N)的算法,当前例子需要计算22次,而O(logN)的算法只需要计算5次

\begin{array}{|rrrrrrrr|}
\hline
& \color{red}{CODE:}&
\hline
\end{array}

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
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}

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
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
#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("%lld\n", powN(a, An_1.arr[0][0]) * powN(b, An_1.arr[1][0]) % MOD1);
    }
    return 0;
}

思考:

Q: 为什么要使用费马小定理?
A: 因为 $ x \mod p-1 <= x $ 所以可以减少计算整数幂时候的计算次数,不使用费马小定理会超时.

作者 @何军
2018 年 03月 24日

emm.坚持原创技术分享，您的支持将鼓励我这个小菜鸡的创作！