九度OJ题目1081：递推数列-快速幂

矩阵快速幂模版。

题目描述：

给定a0,a1,以及an=p*a(n-1) + q*a(n-2)中的p,q。这里n >= 2。 求第k个数对10000的模。

输入：

输入包括5个整数：a0、a1、p、q、k。

输出：

第k个数a(k)对10000的模。

样例输入：

20 1 1 14 5

样例输出：

8359

思路：

求斐波那契数列的时候，递推公式为

f(n) = f(n-1)+f(n-2)

如果n很大的话，需要用到矩阵。

[f(2)f(1)]=[1110]∗[f(1)f(0)]

递推公式得出：

[f(n)f(n−1)]=[1110]n−1[f(1)f(0)]

本题目和斐波那契数列类似的：

[f(n)f(n−1)]=[p1q0]n−1[f(1)f(0)]

注意特判k为1和0时候的值

#include<stdio.h>
#include<math.h>
#include<bitset>
#include<iostream>
using namespace std;
const int MOD =10000;
typedef int Matrix[3][3];
Matrix mp[32];
//返回MOD结果
void mod(int &n)
{
n %= MOD;
}
void mat_mul(const Matrix &A,const Matrix &B,Matrix &R,const int n)
{
for(int i = 0 ; i < n;i++)
{
for(int j=0;j<n;j++)
{
int sum = 0;
for(int k = 0 ; k < n;k++)
{
sum += A[i][k]*B[k][j] ;
mod(sum);
}
R[i][j] = sum;
}
}
}
void mat_assign(Matrix &A, const Matrix &B, const int n)
{
for(int i = 0 ; i <n;i++)
for(int j = 0;j<n;j++)
A[i][j] = B[i][j];
}
void mat_assign(Matrix &A,int v, const int n)
{
for(int i =0;i<n;i++)
for(int j = 0 ; j <n;j++)
if(i==j) A[i][j] =v;
else A[i][j] = 0;
}
void mat_pow_2(const Matrix &A , const int p, const int n)
{
mat_assign(mp[0],A,n);
for(int i =1 ; i <= p;i++)
mat_mul(mp[i-1],mp[i-1],mp[i],n);
}
//计算A的p次
void mat_pow(const Matrix &A,const int p, Matrix &R,const int n)
{
int max_p = (int)log2(p);
mat_pow_2(A,max_p,n);
bitset<32>bin(p);
mat_assign(R,1,n);
for(int i = 0 ; i <= max_p;i++)
{
if(bin[i]==1)
{
Matrix TT;
mat_mul(R,mp[i],TT,n);
mat_assign(R,TT,n);
}
}
}
int main()
{
int a0,a1,p,q,k;
while(~scanf("%d%d%d%d%d",&a0,&a1,&p,&q,&k))
{
if(k==1)
{
cout<<a1<<endl;
}
else if(k==0)
{
cout<<a0<<endl;
}
else
{
Matrix TT;
Matrix T;
T[0][0] = p;
T[0][1] = q;
T[1][0] = 1;
T[1][1] = 0;
mat_pow(T,k-1,TT,2);
//          cout<<TT[0][0]<<" "<<TT[0][1]<<endl;
//          cout<<TT[1][0]<<" "<<TT[1][1]<<endl;
int ans = (TT[0][0]*a1)%MOD+(TT[0][1]*a0)%MOD;
cout<<ans%MOD<<endl;
}

}

}

但是这是一道清华机试题，机试题当然是没有模板的，因此考虑些一个更简单的代码，直接就是二分算矩阵的幂

#include<stdio.h>
#include<math.h>
#include<bitset>
#include<iostream>
using namespace std;
const int MOD = 10000;
struct Matrix
{
int a[5][5];
};
Matrix mul(const Matrix &A, Matrix &B)
{
Matrix tmp;
tmp.a[0][0] = (A.a[0][0]*B.a[0][0]+A.a[0][1]*B.a[1][0])%MOD;
tmp.a[0][1] = (A.a[0][0]*B.a[0][1]+A.a[0][1]*B.a[1][1])%MOD;
tmp.a[1][0] = (A.a[1][0]*B.a[0][0]+A.a[1][1]*B.a[1][0])%MOD;
tmp.a[1][1] = (A.a[1][0]*B.a[0][1]+A.a[1][1]*B.a[1][1])%MOD;
return tmp;
}
int main()
{
int a0,a1,p,q,k;
while(~scanf("%d%d%d%d%d",&a0,&a1,&p,&q,&k))
{
if(k==1) cout<<a1<<endl;
else if(k==0) cout<<a0<<endl;
else
{
Matrix TT;
Matrix T;
T.a[0][0] = p;
T.a[0][1] = q;
T.a[1][0] = 1;
T.a[1][1] = 0;
TT.a[0][0]= 1;
TT.a[0][1]= 0;
TT.a[1][0]= 0;
TT.a[1][1]= 1;

k--;
while(k>=1)
{
if((k&1)==1)
TT = mul(TT,T);
T = mul(T,T);
k = k>>1;
}
//          cout<<TT[0][0]<<" "<<TT[0][1]<<endl;
//          cout<<TT[1][0]<<" "<<TT[1][1]<<endl;
int ans = (TT.a[0][0]*a1)%MOD+(TT.a[0][1]*a0)%MOD;
cout<<ans%MOD<<endl;
}

}
}