斐波那契数列的矩阵快速幂算法

我们在学习C++算法的时候会遇到一个十分基本的问题：

我们都知道斐波那契数列：1，1，2，3，5，8，13，⋯，现在已知第一个1记为第一项，第二个1记为第二项，2记为第三项，以此类推，现在我向里面输入一个n，请你告诉我数列中第n项的数值是多少？

这个题很多人第一眼会感觉十分简单，确实，它确实不难，我们都十分熟悉斐波那契数列的递推公式：a_n = a_n − 1 + a_n − 2(n ≥ 2)

这天然就给了我们递归的土壤：

#include<bits/stdc++.h>
using namespace std;
int Fibonacci(int n){
    if(n == 1 || n == 2)
        return 1;
    else
    	return Fibonacci (n - 1) + Fibonacci (n - 2);
}
int main()
{
    int n;
    cin >> n;
    cout << Fibonacci(n) <<endl;
    return 0;
}

基本上系统学习过C++的同学都可以写出这段代码，但是这里存在一个非常明显的不足：时间复杂度太高，粗略估计：O(2ⁿ)，在这个条件下，基本上n ≥ 50的时候计算机就要算很久才能出结果的（我之前记得计算机每秒的运算次数大概在10⁹次左右，所以原本以为在40的时候差不多就要卡住了，但是实际操作下来的话50左右才会被卡住，只能说是计算能力又进步了，我out了​​）

当然，进一步想，我们的斐波那契数列计算太过于重复了，每一次递归，都要再计算一遍很多之前已经计算出结果的函数，那如果我们用一个数组，将我们计算过的数和相数都放在里面，不就大大简便了呢？所以这就是动态规划：

#include<bits/stdc++.h>
using namespace std;
int main()
{
	int n;
	cin >> n;
	long long *dp = new long long [n + 1];
	dp[1] = 1;
	dp[0] = 0;
	for(int i = 2; i <= n; i++)
		dp[i] = dp[i - 1] + dp[i - 2];
	cout << dp[n] <<endl;
	return 0;
}

（我的码风不是特别规范，大家见谅，这里面的算是最基础的动态规划，比较简单，我就不详细说明了）

通过这个算法，我们已经将时间复杂度降到了O(n)的级别，基本上在long long数值范围内的斐波那契数列值我们都可以计算出来了（超出long long数值范围的部分属于数字处理了，在这篇博客里面暂时不讨论，当然不是笔者不会的原因）

不过正经讲，在笔者做过的信息学竞赛题目中，这类结果一般是让你取模的，所以说更高效的算法也是有必要的。

在算法结构中，有一种十分巧妙的方法叫做快速幂。如果我们要计算xⁿ，这里面最为直接的方法就是将n个x相乘，时间复杂度为O(n)，但这个方法可以进一步简化，将时间复杂度降到O(log₂n)，代码如下：

#include<bits/stdc++.h>
using namespace std;
int main()
{
	int x, n;
	cin >> x >> n;
	int res = 1;
	while(n > 0){
		if(n & 1)res *= x;
		x = x * x;
		n >> = 1;
	}
	cout << res << endl;
	return 0;
}

这种方法其实非常好理解，就是将n转换为二进制，每一次只取n的末位来看，如果是1，那么就将结果乘上x，如果是0，那就不用管，同时，由于数位不断地左移（注意，n >  >  = 1是将n右移，但是从位数上来说，就是左移），x需要自增，也就是变为x²。

这个算法的时间复杂度，就是只有O(log₂n)。

我们知道了这个方法，对于我们解决这个问题有什么帮助呢？

对于斐波那契数列，特殊矩阵的幂可以计算出斐波那契额数列的结果：

$$ \begin{align*} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} &= \begin{bmatrix} 2 & 1 \\ 1 & 0 \end{bmatrix} \\ \begin{bmatrix} 2 & 1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} &= \begin{bmatrix} 3 & 2 \\ 2 & 1 \end{bmatrix} \\ \begin{bmatrix} 3 & 2 \\ 2 & 1 \end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} &= \begin{bmatrix} 5 & 3 \\ 3 & 2 \end{bmatrix}\end{align*} $$

$$ \begin{matrix} \vdots \\ \end{matrix} $$

$$ \begin{align*}\begin{bmatrix} a_n & a_{n-1} \\ a_{n-1} & a_{n-2} \end{bmatrix} \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} &= \begin{bmatrix} a_n + a_{n-1} & a_n \\ a_{n-1} + a_{n-2} & a_{n-1} \end{bmatrix} = \begin{bmatrix} a_{n+1} & a_n \\ a_n & a_{n-1} \end{bmatrix} \end{align*} $$

从中我们可以很轻易地得到： $$ \left[ \begin{matrix} 1&1\\ 1&0 \end{matrix} \right]^n= \left[ \begin{matrix} a_{n+1}&a_n\\ a_n&a_{n-1} \end{matrix}\\ \right] $$

其中a_n是斐波那契数列中的第n项

其实通过这些铺垫，我们通过矩阵快速幂来解决斐波那契问题的答案已经呼之欲出了，最后的一道拦路虎就是如何运用C++来计算矩阵的乘法，这个问题我觉得没有单独拎出来讨论的必要（绝对不是因为笔者懒了不想写了,现在已经凌晨1：30了，好困）

接下来是完整的代码（幸好之前已经写过了，可以早点睡觉）

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef vector<ll> vec;
typedef vector<vec> mac;//mac就是用来表示矩阵的
//计算矩阵的相乘
mac calmartix(mac a,mac b){
	mac ans(2, vec(2) );//只需要一个2*2的矩阵就可以了
	int cnt;
	for(int i = 0;i < 2; i++)
		for(int j = 0;j < 2; j++){
			cnt = 0;
			for(int k = 0;k < 2;k++)
				cnt += a[i][k] * b[k][j];//矩阵的计算
			ans[i][j] = cnt;
		}
	return ans;
}
//快速幂
mac quickcal(mac a,ll n){
	mac res(2, vec(2));
	res[1][1] = 1;
	res[0][0] = 1;
	res[1][0] = 0;
	res[0][1] = 0;//这时候的res是一个单位矩阵
	while(n > 0){
		if(n & 1) res = calmartix(res , a);
		a = calmartix(a , a);
		n >>= 1;
	}//矩阵的快速幂
	return res;
} 
int main()
{
	int n;
	cin >> n;
	mac ans(2,vec(2));
	mac a(2,vec(2));
	a[0][0] = 1;
	a[0][1] = 1;
	a[1][0] = 1;
	a[1][1] = 0;//这里的a就是我们之前讨论的特殊矩阵
	ans = quickcal(a ,n);
	cout << ans[0][1] << endl; //这里你也可以cout << ans[1][0] << endl;都一样看你心情
	return 0;
}

矩阵的计算次数是一个常数，我们在考虑时间复杂度的时候可以忽略，这个算法的时间复杂度就是O(log₂n)

有兴趣的读者可以自己动手亲自写一写~