POJ 3070 Fibonacci （矩阵快速幂）

本文链接：https://blog.youkuaiyun.com/shao824714565/article/details/78034845

本文详细解析了如何运用矩阵快速幂的方法高效求解POJ 3070题目的Fibonacci数列问题，通过矩阵运算大大减少了时间复杂度，提高了算法效率。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

In the Fibonacci integer sequence, F₀ = 0, F₁ = 1, and F_n = F_n_{− 1} + F_n_{− 2} for n ≥ 2. For example, the first ten terms of the Fibonacci sequence are:

0, 1, 1, 2, 3, 5, 8, 13, 21, 34, …

An alternative formula for the Fibonacci sequence is

Given an integer n, your goal is to compute the last 4 digits of F_n.

Input

The input test file will contain multiple test cases. Each test case consists of a single line containing n (where 0 ≤ n ≤ 1,000,000,000). The end-of-file is denoted by a single line containing the number −1.

Output

For each test case, print the last four digits of F_n. If the last four digits of F_nare all zeros, print ‘0’; otherwise, omit any leading zeros (i.e., print F_n mod 10000).

Sample Input

Sample Output

如果不是提前知道斐波那契数列可以用矩阵求解，真的不好把他们两个扯上关系。

第一部分：矩阵的基础知识

1.结合性（AB)C=A(BC）．

2.对加法的分配性（A+B)C=AC+BC，C(A+B）=CA+CB ．

3.对数乘的结合性 k(AB）=（kA)B =A(kB）．

4.关于转置 (AB)'=B'A'．

一个矩阵就是一个二维数组，为了方便声明多个矩阵，我们一般会将矩阵封装一个类或定义一个矩阵的结构体，我采用的是后者。（弱鸡的我也直只会用结构体实现）

第二部分：矩阵相乘

若A为n×k矩阵，B为k×m矩阵，则它们的乘积AB(有时记做A·B)将是一个n×m矩阵。前一个矩阵的列数应该等于后一个矩阵的行数，得出的矩阵行数等于前一个矩阵的行数，列数等于后一个矩阵的行数。

其乘积矩阵AB的第i行第j列的元素为：

举例：A、B均为3*3的矩阵：C=A*B，下面的代码会涉及到两种运算顺序，第一种就是直接一步到位求，第二种就是每次求一列，比如第一次，C00+=a00*b00,C01+=a00*b01……第二次C00+=a00*b10，C01+=a01*b11……以此类推。。。

C00 = a00*b00 + a01*b10 + a02*b20

C01 = a00*b01 + a01*b11 + a02*b21

C02 = a00*b02 + a01*b12 + a02*b22

C10 = a10*b00 + a11*b10 + a12*b20

C11 = a10*b00 + a11*b11 + a12*b21

C12 = a10*b02 + a11*b12 + a12*b22

C20 = a20*b00 + a21*b10 + a22*b20

C21 = a20*b01 + a21*b11 + a22*b21

C22 = a20*b02 + a21*b12 + a22*b22

先来一个矩阵相乘的代码：

struct matrix
{
    int a[2][2];
};
matrix Multiply_Matrix(matrix x,matrix y)//两个矩阵相乘
{
    matrix z;//存储新的矩阵
    memset(z.a,0,sizeof(z.a));
    for(int i=0; i<2; i++)
        for(int j=0; j<2; j++)
            for(int k=0; k<2; k++)
            {
                z.a[i][j]+=x.a[i][k]*y.a[k][j];
                z.a[i][j]%=MOD;
            }
    return z;
}

第三部分：矩阵快速幂 //其实和普通快速幂类似，只不过这里需要得到的是一个矩阵

神马是幂？【很多时候会被高大上的名字吓到。。。导致学习效率降低。。。其实没辣么可怕，很简单！！！】

幂又称乘方。表示一个数字乘若干次的形式,如n个a相乘的幂为a^n ，或称a^n为a的n次幂。a称为幂的底数，n称为幂的指数。——引自.度娘百科

这类题，指数都是很大很大很大很大很大很大很大的。。。霸王硬上弓的话，很容易超时的 T_T 。。。所以得快速幂→_→

学过之后发现，其实矩阵快速幂的核心思想跟以前学过的快速幂取模非常非常相似，只是矩阵乘法需要另外写个函数，就是上面那个代码。。。

快速幂的思路就是：

设A为矩阵，求A的N次方，N很大，1000000左右吧。。。

先看小一点的，A的9次方

A^9

= A*A*A*A*A*A*A*A*A 【一个一个乘，要乘9次】

= A*(A*A)*(A*A)*(A*A)*(A*A)【保持格式的上下统一，所以加上这句】

= A*(A^2)^4 【A平方后，再四次方，还要乘上剩下的一个A，要乘6次】

= A*((A^2)^2)^2【A平方后，再平方，再平方，还要乘上剩下的一个A，要乘4次】

也算是一种二分思想的应用吧，1000000次幂，暴力要乘1000000次，快速幂就只要（log2底1000000的对数）次，大约20次。。。这。。。我没错吧。。。

单位矩阵： n*n的矩阵 mat ( i , i )=1; 任何一个矩阵乘以单位矩阵就是它本身 n*单位矩阵=n，可以把单位矩阵等价为整数1。（单位矩阵用在矩阵快速幂中）

例如下图就是一个7*7的单位矩阵：

下面来实现一个矩阵快速幂：

通过位运算，减少计算次数

matrix Quick_Matrix(matrix B,int n)//矩阵快速幂
{
    int i,j;
    matrix ans;
    for(i=0; i<2; i++)//初始化ans
        for(j=0; j<2; j++)
            ans.a[i][j]=(i==j);
    while(n)
    {
        if(n&1)
            ans=Multiply_Matrix(ans,B);
        n>>=1;
        B=Multiply_Matrix(B,B);
    }
    return ans;
}

但是如何把斐波那契数列和矩阵联系起来呢？

用递归很巧妙的就联系起来了

附上代码：

#include<stdio.h>
#include<string.h>
#include<string>
#include<algorithm>
using namespace std;
#define MOD 10000
struct matrix
{
    int a[2][2];
};
matrix Multiply_Matrix(matrix x,matrix y)//两个矩阵相乘
{
    matrix z;//存储新的矩阵
    memset(z.a,0,sizeof(z.a));
    for(int i=0; i<2; i++)
        for(int j=0; j<2; j++)
            for(int k=0; k<2; k++)
            {
                z.a[i][j]+=x.a[i][k]*y.a[k][j];
                z.a[i][j]%=MOD;
            }
    return z;
}
matrix Quick_Matrix(matrix B,int n)//矩阵快速幂
{
    int i,j;
    matrix ans;
    for(i=0; i<2; i++)//初始化ans
        for(j=0; j<2; j++)
            ans.a[i][j]=(i==j);
    while(n)
    {
        if(n&1)
            ans=Multiply_Matrix(ans,B);
        n>>=1;
        B=Multiply_Matrix(B,B);
    }
    return ans;
}
int main()
{
    int i,j,n;
    matrix A,C;
    while(~scanf("%d",&n)&&n!=-1)
    {
        memset(C.a,0,sizeof(C.a));
        A.a[0][0]=A.a[0][1]=A.a[1][0]=1;
        A.a[1][1]=0;
        if(n==0)
        {
            printf("0\n");
            continue;
        }
        C=Quick_Matrix(A,n);
        printf("%d\n",C.a[0][1]);
    }
    return 0;
}