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; }