矩阵的快速取幂法——求类Fibonacci数列（3）

AOJ 395 

Advanced Fibonacci

Time Limit: 1000 ms       Memory Limit: 64 MB    Total Submission: 16       Accepted: 5 

Description

如果你感觉到普通Fibonacci数列是简约的美，那么我们来试试通俗的Fibonacci，它是这样定义的

F0=x
F1=y
Fn=p*Fn-1+q*Fn-2，
同样，这次我们要计算这个数列的前10^9次项
Input输入包括多组数据，每组数据包括一行格式为 x,y,p,q,n，其中1 =< x,y,p,q =< 1000,n<=10^9,输入以0 0 0 0 0结束
Output每行输出给定Fibonacci的第N项，由于数据比较庞大，给出结果 mod 10003的值

Sample Input

1 1 1 1 1

1 1 1 1 3
1 2 1 3 2
1 2 1 3 3
1 2 1 3 4
0 0 0 0 0
Sample Output

1[EOL]
2[EOL]
2[EOL]
5[EOL]
11[EOL][EOF]

代码：

#include<iostream>
#include<cstring>
using namespace std;
int sq[35][2][2];
void mlt(int p,int q)
{ 
    int i,j,k,l; 
    memset(sq,0,sizeof(sq));
    sq[0][0][0]=p,sq[0][0][1]=q; 
    sq[0][1][0]=1,sq[0][1][1]=0; 
    for(l=1;l<32;l++) 
        for(i=0;i<2;i++) 
            for(j=0;j<2;j++) 
                for(k=0;k<2;k++) 
                    sq[l][i][j]=(sq[l][i][j]+sq[l-1][i][k]*sq[l-1][k][j])%10003; 
}
int main()
{ 
    int n,ans,i,j,k,l,x,y,p,q; 
    int m[2][2]={0}; 
    while(cin>>x>>y>>p>>q>>n&&!(x==0&&y==0&&p==0&&q==0&&n==0)) 
    { 
        if(n==1) 
        {
            cout<<x<<endl;
            continue;
        } 
        if(n==2) 
        {
            cout<<y<<endl;
            continue;
        } 
        int t[2][2]={{1,0},{0,1}}; 
        n=n-2;l=0; 
        mlt(p,q); 
        while(n) 
        { 
            if(n&1) 
            { 
                memset(m,0,sizeof(m)); 
                for(i=0;i<2;i++) 
                    for(j=0;j<2;j++) 
                        for(k=0;k<2;k++) 
                            m[i][j]=(m[i][j]+sq[l][i][k]*t[k][j])%10003; 
                for(i=0;i<2;i++) 
                    for(j=0;j<2;j++) 
                        t[i][j]=m[i][j]; 
            } 
            l++;
            n=n>>1; 
        } 
        ans=(t[0][0]*y+t[0][1]*x)%10003;
        cout<<ans<<endl; 
    } 
    return 0;
}