「POJ 3070」Fibonacci - 矩阵乘法

发表于 更新于 分类于 ， 阅读次数： Disqus：

题意描述

POJ 链接

给定一个数 \(n\)，求斐波那契数列第 \(n\) 项的后四位。

数据范围：\(1 \leq n \leq 10^9\)

解题思路

我们可以使用矩阵来加速地推。

我们知道斐波那契的通项公式为：

\[ a_i = a_{i - 1} + a_{i - 2} \]

故我们可以设矩阵 \(F(n)\)：

\[ F(n) = \begin{bmatrix} a_n & a_{n + 1} \end{bmatrix} \]

同时我们可以设：

\[ F(n) = F(n - 1) \times X \]

即为：

\[ \begin{bmatrix} a_n & a_{n + 1} \end{bmatrix} = \begin{bmatrix} a_{n - 1} & a_n \end{bmatrix} \times X \]

我们可以推出：

\[ \begin{bmatrix} a_n & a_{n + 1} \end{bmatrix} = \begin{bmatrix} a_{n - 1} & a_n \end{bmatrix} \times \begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix} \]

所以我们最终可以推出下列式子：

\[ F(n) = F(1) \times \begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix} ^ {n - 1} = \begin{bmatrix} 1 & 1 \end{bmatrix} \times \begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix} ^ {n - 1} \]

最终利用快速幂就可以在 \(O(\log n)\) 的时间复杂度中求出答案。

代码演示

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
#include <cstdio>
#include <cstring>
#include <algorithm>

int fib[1 + 1][2 + 1], x[2 + 1][2 + 1];

void init() {
    fib[1][1] = fib[1][2] = 1;
    x[1][1] = 0;
    x[1][2] = x[2][1] = x[2][2] = 1;
}

void mul(int a[][2 + 1], int ax, int ay, int b[][2 + 1], int bx, int by, int ans[][2 + 1], int mod) {
    int s1[2 + 1][2 + 1], s2[2 + 1][2 + 1];
    memcpy(s1, a, sizeof(int) * (ax + 1) * (ay + 1));
    memcpy(s2, b, sizeof(int) * (bx + 1) * (by + 1));

    for (int i = 1; i <= ax; i++) {
        for (int j = 1; j <= by; j++) {
            int sum = 0;
            for (int k = 1; k <= ay; k++) sum = (sum + s1[i][k] * s2[k][j]) % mod;
            ans[i][j] = sum;
        }
    }
}

void pow(int a[][2 + 1], int ax, int ay, int p, int ans[][2 + 1], int mod) {
    int s1[2 + 1][2 + 1];
    memcpy(s1, a, sizeof(int) * (2 + 1) * (2 + 1));
    memset(ans, 0, sizeof(int) * (2 + 1) * (2 + 1));
    for (int i = 1; i <= std::min(ax, ay); i++) ans[i][i] = 1;


    for (; p; p >>= 1) {
        if (p & 1) mul(ans, 2, 2, s1, 2, 2, ans, mod);
        mul(s1, 2, 2, s1, 2, 2, s1, mod);
    }
}

int main() {
    int n;

    while (scanf("%d", &n) != EOF) {
        if (n == -1) break;
        if (n == 0) {
            printf("0\n");
            continue;
        }
        if (n == 1 || n == 2) {
            printf("1\n");
            continue;
        }
        if (n == 3) {
            printf("2\n");
            continue;
        }

        init();
        pow(x, 2, 2, n - 1, x, 10000);
        mul(fib, 1, 2, x, 2, 2, fib, 10000);

        printf("%d\n", fib[1][1]);
    }

    return 0;
}