快速幂

一、快速幂

算出数 x 的 k 次方的值。

普通求n次方的算法为o（n）效率较低。

考虑加速幂运算的方法，如果 n = 2^k，可以将其表示为：

x^n = ((x^2)^2)...

只需要k次平方运算即可轻松求得。由此我们联想到，先将n表示为2的幂的和。

n = 2^k1 + 2^k2 + 2^k3 + ...

就有：

x^k = x^(2^k1) + x^(2^k2) + x^(2^k3) + ...

只要在依次求x^2的同时进行计算就可以了，最终得到了o(logn)的算法。

举个例子 求5^22

22 = b(10110)

因此 5^22 = 5^16 + 5^4 + 5^2;

下面给出了快速幂求模的代码。非取模运算摘开mod即可。

typedef long long ll;

ll mod_pow(ll x, ll n, ll mod) {
	ll ans = 1;
	while ( n > 0 ) {
		if ( n&1 ) ans = ans*x%mod;		// 如果二进制最低位为1 则乘上 x^(2^i) 
		x = x*x%mod;					// 将 x 平方 
		n >>= 1;
	}
	return ans; 
}

使用快速幂的好处就是将o（n）优化到了o（logn）

二、矩阵快速幂

有一个很经典的数列叫斐波那契数列，f[n] = f[n-1] + f[n-2];

求第i项的复杂度是 o（n）；不过这个效率太低了。对于n的规模较大的题目无法在规定时间内求出。

斐波那契数列有个通项公式：

试中包含无理数，在计算过程中会有精度问题，尤其是在需要取模的计算中无法简单求出。

而且在其它问题中有一些很难直接求得通项公式。不过这些情况都可以不求出通项，而用矩阵高效的求出第n项的值。

首先，把斐波那契数列的递推式表示成矩阵就是：

记这个01矩阵为A， 则有

因此只要求出A^n即可求出f[n]了，快速求A^n的计算就可以参考前一节的快速幂运算。

在o（logn）时间里求出第n项的值。

typedef long long ll;

// 用二维vector来表示矩阵 
typedef vector<ll> vec;
typedef vector<vec> mat;

// 模 
const int M = 1000000007;

// 计算 A*B
mat mul(mat& A, mat& B) {
	mat C(A.size(), vec(B[0].size()));
	for ( int i = 0 ; i < A.size() ; ++ i ) {
		for ( int k = 0 ; k < B.size() ; ++ k ) {
			for ( int j = 0 ; j < B[0].size() ; ++ j ) {
				C[i][j] = (C[i][j]+A[i][k]*B[k][j]%M)%M;
			}
		}
	}
	return C;
}

// 计算 A^B
mat pow(mat A, ll n) {
	mat B(A.size(), vec(A.size()));
	for ( int i = 0 ; i < A.size() ; ++ i ) {
		B[i][i] = 1;
	}
	while ( n > 0 ) {
		if ( n & 1 ) B = mul(B, A);
		A = mul(A, A);
		n >>= 1;
	}
	return B;
} 

对于本题斐波那契数列的求解代码如下：

// 输入 
ll n;

void solve() {
	mat A(2, vec(2));
	A[0][0] = 1; A[0][1] = 1;
	A[1][0] = 1; A[1][1] = 0;
	A = pow(A, n);
	printf( "%lld\n", A[1][0] );
}

更一般的，对于m项递推式，如果记递推式如下：

则可以把递推式写成如下矩阵形式

通过计算这个矩阵的n次幂，就可以在o（m^3 logn）的时间内计算出第n项的值。

不过，如果递推式中含有常数则稍微复杂一些，需变形如下：