矩阵乘法求线性递推式的第n项

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <cstring>
#include <string.h>
#include <vector>
#include <queue>      //poj2118  矩阵乘法求线性递推式的第n项
#define mod 10000
using namespace std;
struct matrix
{
int c[105][105];
void init()
{
for(int j=0;  j<105; ++j)
{
for(int t=0; t<105; ++t)
c[j][t]=0;
}
}
}p;
int k, a[105], b[105];
matrix mul(matrix a1, matrix b1)
{
matrix q;
q.init();
int j, t, g;
for(j=0; j<k; ++j)
{
for(t=0; t<k; ++t)
{
for(g=0; g<k; ++g)
{
q.c[j][t]+=a1.c[j][g]*b1.c[g][t];
q.c[j][t]%=mod;
}
}
}
return q;
}
/*
matrix f(int x)  //二分法求矩阵的多次幂
{
if(x==1)return p;
matrix s=f(x/2);
s=mul(s, s);
if(x&1)
s=mul(s, p);
return s;
}
*/
matrix f(int x)   //矩阵多次幂可以防止爆栈
{
matrix q, s=p;
int t, j;
for(j=0; j<k; ++j)    //构建单位矩阵
{
for(t=0; t<k; ++t)
{
if(j==t)
q.c[j][t]=1;
else q.c[j][t]=0;
}
}
while(x)
{
if(x&1)
q=mul(s, q);
x=x>>1;
s=mul(s, s);
}
return q;
}

int main()
{
int t, j, i, ans;
while(scanf("%d", &k)&&k)
{
p.init();
for(t=0; t<k; ++t)
scanf("%d", a+t);
for(t=1; t<=k; ++t)
scanf("%d", b+t);
scanf("%d", &i);
if(i<k)
{
printf("%d\n", a[i]);
continue;
}
for(t=k, j=0; t>0; --t, ++j)   //构建初始矩阵最后一行的系数,注意这里系数的顺序
p.c[k-1][j]=b[t];
for(j=0; j<k-1; ++j)    //右上角的(k-1)*(k-1)的小矩阵中的主对角线上填1，矩阵第n行填对应的系数，其它地方都填0
p.c[j][j+1]=1;
i=i-k+1;
matrix result=f(i);
for(t=ans=0; t<k; ++t)  //数组a相当于另一个k行一列矩阵
{
ans+=result.c[k-1][t]*a[t];
ans%=mod;
}
printf("%d\n", ans);
}
return 0;
}