不多说贴代码
牛顿迭代
#include<stdio.h>
#include<string.h>
int n,m,k;
struct matrix{
int a[33][33];
}a,b,c;
//s13=s6+a7+a7*s6;
//s12=s6+a6*s6;
//A+A^2+A^3+A^4+A^5+A^6=A+A^2+A^3 +A^3*(A+A^2+A^3)
//A+A^2+A^3+A^4= A+A^2 + A^2*(A+A^2) //偶
//A+A^2+A^3= A+ A^2 +A^2*(A) //积
//A+A^2 = A+A*A
//A
matrix cheng(matrix x,matrix y){ //x*y;
memset(c.a,0,sizeof(c));
int i,k,j;
for(i=0;i<n;i++)
for(j=0;j<n;j++)
for(k=0;k<n;k++){
c.a[i][j]+=(x.a[i][k]*y.a[k][j]);
if(c.a[i][j]>=m) c.a[i][j]%=m;
}
return c;
}
matrix jia(matrix x,matrix y){ //x+y;
memset(c.a,0,sizeof(c));
int i,j;
for(i=0;i<n;i++)
for(j=0;j<n;j++){
c.a[i][j]=(x.a[i][j] + y.a[i][j]);
if(c.a[i][j]>=m) c.a[i][j]%=m;
}
return c;
}
//s13=s6+a7+a7*s6;
//s12=s6+a6*s6;
//A+A^2+A^3+A^4+A^5+A^6+A^7=s3+a4+a4*s3
//A+A^2+A^3+A^4+A^5+A^6=A+A^2 +A^3 +A^3*(A+A^2+A^3)
//A+A^2+A^3+A^4+A^5=s2+a3+a3*s2//奇
//A+A^2+A^3+A^4= A+A^2 + A^2*(A+A^2) //偶
//A+A^2+A^3= A+ A^2 +A^2*(A) //积
//A+A^2 = A+A*A
//A
matrix ans(matrix e, int x){
matrix aa,bb,cc;
int i,j;
if(x==1){
return a;
}
if(x&1){//奇数
aa=ans(e,x>>1);
cc=cheng(a,e);
bb=jia(aa,jia(cc,cheng(cc,aa)));
a=cheng(cc,a);
}
else{//偶数
aa=ans(e,x>>1);
bb=jia(aa,cheng(a,aa));
a=cheng(a,a);
}
/*
printf("\n\nk=%d\n",x);
for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",bb.a[i][j]," \n"[j==n-1]);
*/
return bb;
}
int main(){
int i,j;
while(~scanf("%d%d%d",&n,&k,&m)){
for(i=0;i<n;i++)
for(j=0;j<n;j++){
scanf("%d",&(a.a[i][j]));
if(a.a[i][j]>=m) a.a[i][j]%=m;
}
b=ans(a,k);
for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",b.a[i][j]," \n"[j==n-1]);
}
return 0;
}
矩阵套矩阵
#include<stdio.h>
#include<string.h>
int n,m,k;
struct matrix{
int a[33][33];
}a,b,c;
struct node{
matrix a[2][2];
};
/*
s[k]=A+A^2+...+A^k
矩阵:
|s[0] A|*|I O|^k=|s[k] A^(k+1)| (I为单位矩阵,O为0矩阵,A为输入的矩阵)
|O O| |I A| |O O |
*/
matrix cheng(matrix x, matrix y ){
matrix c;
memset(c.a,0,sizeof(c.a));
int i,j,k;
for(i=0;i<n;i++)
for(j=0;j<n;j++)
for(k=0;k<n;k++){
c.a[i][j]+=(x.a[i][k] * y.a[k][j]);
if(c.a[i][j]>=m) c.a[i][j]%=m;
}
return c;
}
matrix jia(matrix x,matrix y){
matrix c;
memset(c.a,0,sizeof(c.a));
int i,j;
for(i=0;i<n;i++)
for(j=0;j<n;j++){
c.a[i][j]=(x.a[i][j] + y.a[i][j]);
if(c.a[i][j]>=m) c.a[i][j]%=m;
}
return c;
}
node cheng_matrix(node x,node y){
node c;
int i,j,k;
for(i=0;i<2;i++)
for(j=0;j<2;j++)
memset(c.a[i][j].a,0,sizeof(c.a[i][j].a));
for(i=0;i<2;i++)
for(j=0;j<2;j++)
for(k=0;k<2;k++){
c.a[i][j]=jia( c.a[i][j] , cheng(x.a[i][k], y.a[k][j] ));
}
return c;
}
int main(){
int i,j,k1,k2;
node x,y;
while(~scanf("%d%d%d",&n,&k,&m)){
memset(a.a,0,sizeof(a.a));
for(i=0;i<n;i++){
b.a[i][i]=1;
for(j=0;j<n;j++){
scanf("%d",&(a.a[i][j]));
if(a.a[i][j]>=m) a.a[i][j]%=m;
}
}
/*
s[k]=A+A^2+...+A^k
矩阵:
|s[0] A|*|I O|^k=|s[k] A^(k+1)| (I为单位矩阵,O为0矩阵,A为输入的矩阵)
|O O| |I A| |O O |
*/
x.a[0][0]=c;
x.a[1][0]=c;
x.a[1][1]=c;
y.a[0][1]=c;
x.a[0][1]=a;
y.a[1][1]=a;
y.a[0][0]=b;
y.a[1][0]=b;
/* for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",x.a[0][0].a[i][j]," \n"[j==n-1]);
printf("\n\n\n\n\n\n\n");
*/
while(k){
if(k&1) x=cheng_matrix(x,y);
y=cheng_matrix(y,y);
k>>=1;
/*
for(k1=0;k1<2;k1++)
for(k2=0;k2<2;k2++){
printf("X%d%d:\n",k1,k2);
for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",x.a[k1][k2].a[i][j]," \n"[j==n-1]);
}
for(k1=0;k1<2;k1++)
for(k2=0;k2<2;k2++){
printf("y%d%d:\n",k1,k2);
for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",y.a[k1][k2].a[i][j]," \n"[j==n-1]);
*/
}
for(i=0;i<n;i++)
for(j=0;j<n;j++)
printf("%d%c",x.a[0][0].a[i][j]," \n"[j==n-1]);
}
return 0;
}