This blog post is only available in Simplified Chinse.
斐波那契（Fibonacci）数列的递推式是： 。根据这个递推式，我们可以在 复杂度内求出第 n 项，但是当 n 很大时，这种方法就显得很慢。其实利用矩阵快速幂，我们可以在 内求出第 n 项。
矩阵
定义
一个 m×n 的矩阵（matrix）是一个由 m 行（row）n 列（column）元素排列成的矩形阵列。矩阵里的元素可以是数字、符号或数学式。例如，以下就是一个 2×3 的矩阵：
1 & 2 & 3 \\\\ 4 & 6 &5 \end{bmatrix}$$ 矩阵是线性代数的知识，更多这里就不介绍了。 ### 矩阵的加减法 矩阵的加法减法很简单： $$\begin{bmatrix} 1 & 2 & 3 \\\\ 4 & 6 & 5 \end{bmatrix} + \begin{bmatrix} 2 & 1 & 0 \\\\ 2 & 0 & 1 \end{bmatrix} = \begin{bmatrix} 3 & 3 & 3 \\\\ 6 & 6 & 6 \end{bmatrix}$$ 减法也同理。 ### 矩阵乘以矩阵 矩阵之间的相乘就比较复 (bian) 杂 (tai) 了，根据 [Wikipedia](https://zh.wikipedia.org/wiki/%E7%9F%A9%E9%99%A3%E4%B9%98%E6%B3%95#%E7%94%B1%E5%AE%9A%E7%BE%A9%E7%9B%B4%E6%8E%A5%E8%A8%88%E7%AE%97 "Wikipedia")： 矩阵相乘最重要的方法是一般矩阵乘积。它只有在第一个矩阵的列数（column）和第二个矩阵的行数（row）相同时才有定义。一般单指矩阵乘积时，指的便是一般矩阵乘积。若 $A$ 为 $ m\times n $ 矩阵，$B$ 为 $n\times p$ 矩阵，则他们的乘积 $AB$ (有时记做 $A \cdot B$）会是一个 $m\times p$ 矩阵。其乘积矩阵的元素如下面式子得出： $$AB_{ij} = \sum_{r=1}^{n}a_{ir}b_{rj}= a_{i1}b_{1j}+a_{i2}b_{2j}+\dots+a_{in}b_{nj}$$ 看个例子就懂了： $$\begin{bmatrix} 1 & 2 \\\\ 4 & 3 \end{bmatrix} + \begin{bmatrix} 2 & 1 \\\\ 2 & 0 \end{bmatrix} = \begin{bmatrix} 3 & 3 \\\\ 6 & 6 \end{bmatrix}$$ 第一个矩阵第一行的每个数字（1 和 2）各自乘以第二个矩阵第一列对应位置的数字（2 和 2），然后将乘积相加（ 1 x 2 + 2 x 2），得到结果矩阵左上角的值 6。同理，其他的几个值也是这么算的。维基百科上有一张形象的图：  ## 矩阵快速幂求斐波那契数列 其实 [POJ上这道模板题](http://poj.org/problem?id=3070 "POJ上这道模板题") 已经告诉你怎么做了： $$\begin{bmatrix} F_{n+1} & F_n \\\\ F_n & F_{n-1} \end{bmatrix} = \begin{bmatrix} 1 & 1 \\\\ 1 & 0 \end{bmatrix}^n = \underbrace{ \begin{bmatrix} 1 & 1 \\\\ 1 & 0 \end{bmatrix} \begin{bmatrix} 1 & 1 \\\\ 1 & 0 \end{bmatrix}\cdots \begin{bmatrix} 1 & 1 \\\\ 1 & 0 \end{bmatrix} }_ {\text{n times}}$$ 也就是说矩阵 $$\begin{bmatrix} 1 & 1 \\\\ 1 & 0 \end{bmatrix} $$ 的 n 次幂里面四个数就分别是 $F_{n+1}$、$F_n$、$F_n$、$F_{n-1}$ 。 ## 例题 ### 模板题 [洛谷 3390 【模板】矩阵快速幂](https://www.luogu.org/problemnew/show/P3390 "洛谷 3390 【模板】矩阵快速幂") [POJ3070 - Fibonacci](http://poj.org/problem?id=3070 "POJ3070 - Fibonacci") ### 综合 [CodeForces 551D. GukiZ and Binary Operations](https://codeforces.com/problemset/problem/551/D "CodeForces 551D. GukiZ and Binary Operations") ## 代码 以下是洛谷上的模板题代码： ```cpp #include<cstdio> #include<cstring> #include<iostream> using namespace std; const int maxn=105,tt=1000000007;; int n; long long m; inline int read(){ int ret=0,f=1;char ch=getchar(); while (ch<'0'||ch>'9') {if (ch=='-') f=-1;ch=getchar();} while (ch>='0'&&ch<='9') ret=ret*10+ch-'0',ch=getchar(); return ret*f; } struct Matrix{ int a[maxn][maxn]; void ReadMatrix(){ memset(a,0,sizeof(a)); for (int i=0;i<n;i++) for (int j=0;j<n;j++) a[i][j]=read(); } void Init(){ memset(a,0,sizeof(a)); for (int i=0;i<n;i++) a[i][i]=1; } void Clear(){ memset(a,0,sizeof(a)); } Matrix operator *(Matrix b){ Matrix c; c.Clear(); for (int i=0;i<n;i++) for (int j=0;j<n;j++) for (int k=0;k<n;k++) c.a[i][j]=((long long)c.a[i][j]+(long long)a[i][k]*b.a[k][j])%tt; return c; } Matrix operator ^(long long b){ Matrix ret; ret.Init(); Matrix w; w.Clear(); for (int i=0;i<n;i++) for (int j=0;j<n;j++) w.a[i][j]=a[i][j]; while (b){ if (b%2) ret=ret*w; b=b/2;w=w*w; } return ret; } void Write(){ for (int i=0;i<n;i++){ for (int j=0;j<n;j++) printf("%d ",a[i][j]); printf("\n"); } } }; int main(){ n=read();scanf("%lld",&m); Matrix a; a.ReadMatrix(); // printf("Read part finished.\n"); a=a^m; a.Write(); return 0; } ``` 下面是 POJ 上模板题代码： ```cpp #include<cstdio> #include<cstring> #include<iostream> using namespace std; const int tt=10000; int n; struct Matrix{ int a[3][3]; void Init(){ memset(a,0,sizeof(a)); a[0][0]=a[1][0]=a[0][1]=1; a[1][1]=0; } void Clear(){ memset(a,0,sizeof(a)); } Matrix operator *(Matrix b){ Matrix c; c.Clear(); for (int i=0;i<2;i++) for (int j=0;j<2;j++) for (int k=0;k<2;k++) c.a[i][j]=((long long)c.a[i][j]+(long long)a[i][k]*b.a[k][j])%tt; return c; } Matrix operator ^(long long b){ Matrix ret; ret.Init(); Matrix w; w.Clear(); for (int i=0;i<2;i++) for (int j=0;j<2;j++) w.a[i][j]=a[i][j]; while (b){ if (b%2) ret=ret*w; b=b/2;w=w*w; } return ret; } void Write(){ for (int i=0;i<2;i++){ for (int j=0;j<2;j++) printf("%d ",a[i][j]); printf("\n"); } } }; int main(){ scanf("%d",&n);n--; while (n!=-2){ if (n==-1) printf("0\n"); else{ Matrix a; a.Init(); a=a^n; printf("%d\n",a.a[1][0]); } scanf("%d",&n);n--; } return 0; } ``` ## 参考 [矩阵 - 维基百科，自由的百科全书](https://zh.wikipedia.org/wiki/%E7%9F%A9%E9%98%B5#%E7%9F%A9%E9%98%B5%E7%9A%84%E8%AE%A1%E7%AE%97 "矩阵 - 维基百科，自由的百科全书") [矩阵乘法 - 维基百科，自由的百科全书](https://zh.wikipedia.org/wiki/File:Matrix_multiplication_diagram.PNG "矩阵乘法 - 维基百科，自由的百科全书") [理解矩阵乘法 - 阮一峰的网络日志](http://www.ruanyifeng.com/blog/2015/09/matrix-multiplication.html "理解矩阵乘法 - 阮一峰的网络日志")