POJ 3070 Fibonacci数列 矩阵乘法及乘幂求法

可以先推算一下，会发现m0,m1,m2，m3等为方阵的乘幂，其第1行第2列就是Fibonacci数列的项，关键是如何求矩阵的乘幂呢？

分奇数幂和偶数幂两种情况

A(2)^n因为方阵的乘法有结合律，所以A(2)^n=A(2)^(n/2)*A(2)^(n/2)，不妨设n是偶数
所以求A(n)就可以化成求A(n/2)并作一次乘法，所以递归方程是：A(n)=A(n/2)^2

但是如果n为奇数，可以A(n)=A((n-1)/2) * A((n+1)/2或者A(n)=A((n-1)/2)^2 * A(1)

#include <iostream>
using namespace std;
struct matrix{
	int a[2][2];
};
matrix m0,m1,m2;
void initial(){
	m0.a[0][0]=1;m0.a[0][1]=0;m0.a[1][0]=0;m0.a[1][1]=1;
	m1.a[0][0]=1;m1.a[0][1]=1;m1.a[1][0]=1;m1.a[1][1]=0;
	m2.a[0][0]=2;m2.a[0][1]=1;m2.a[1][0]=1;m2.a[1][1]=1;
}
matrix mul(matrix x,matrix y){
	matrix res;
	res.a[0][0] = (x.a[0][0]*y.a[0][0] + x.a[0][1]*y.a[1][0]) % 10000;
	res.a[0][1] = (x.a[0][0]*y.a[0][1] + x.a[0][1]*y.a[1][1]) % 10000;
	res.a[1][0] = (x.a[1][0]*y.a[0][0] + x.a[1][1]*y.a[1][0]) % 10000;
	res.a[1][1] = (x.a[1][0]*y.a[0][1] + x.a[1][1]*y.a[1][1]) % 10000;
	return res;
}
matrix powm(int n){
	if (n == 0)
		return m0;
	if (n == 1)
		return m1;
	if (n == 2)
		return m2;
	int m = n/2;//注意若n为奇数的情况
	matrix temp = powm(m);
	temp = mul(temp,temp);
	if (n%2 == 1)
	{
		temp = mul(temp,m1);
	}
	return temp;
}
int main(){
	int n;
	matrix res;
	initial();
	while (cin>>n,n != -1)
	{
		res =  powm(n);
		cout<<res.a[0][1]<<endl;
	}
	return 0;	
}

关于矩阵的乘幂运算也可以如下运算，用于POJ3744中：

#include <iostream>  
#include <cmath>  
#include <algorithm>  
using namespace std;  
const double esp = 1e-15;  
struct Matrix{  
    double matrix[2][2];  
    Matrix operator*(const Matrix & m){//重载运算符*，矩阵乘法  
        Matrix ans = {0,0,0,0};  
        for (int i = 0;i < 2;i++)  
            for (int k =0 ;k < 2 ;k ++)  
                if (matrix[i][k]>esp)  
                    for (int j = 0;j < 2;j++)  
                    {  
                        ans.matrix[i][j]+=matrix[i][k]*m.matrix[k][j];  
                    }  
                return ans;  
    }  
}unit = {1,0,0,1};  
Matrix pow(Matrix a,int n){//矩阵求幂  
    Matrix p = a,ans = unit;  
    while(n){  
        if (n&1)  
        {  
            ans = ans*p;  
        }  
        p = p*p;  
        n>>=1;  
    }  
    return ans;  
}