[POJ]3233. Matrix Power Series

Analysis

    找了几道矩阵乘法的题目练一练，这道是最水的。大意是给定N*N矩阵A，求A+A^1+A^2+....+A^k。那么使用二分的方法，将要求的矩阵变形为

A+A^1+A^2+....+A^(k div 2) + A^(k div 2)*(A+A^1+A^2+....+A^(k div 2))

就可以在log的时间内出解了。

Accepted Code

type    matrix=array[1..30,1..30] of longint;var    a,ans,tmp:matrix;    n,k,m,i,j:longint;function matrixpower(a:matrix;var b:matrix;p:longint):matrix;var    tmp,tmp2,b2,b3:matrix;    i,j,k:longint;begin    if p=1 then    begin        b:=a;        matrixpower:=a;    end    else    begin        tmp:=matrixpower(a,b,p shr 1);        if p and 1=0 then        begin            matrixpower:=tmp;            for i:=1 to n do                for j:=1 to n do                    for k:=1 to n do                        matrixpower[i,j]:=(matrixpower[i,j]+tmp[i,k]*b[k,j]) mod m;            for i:=1 to n do                for j:=1 to n do                    b2[i,j]:=0;            for i:=1 to n do                for j:=1 to n do                    for k:=1 to n do                        b2[i,j]:=(b2[i,j]+b[i,k]*b[k,j]) mod m;            b:=b2;        end        else        begin            tmp2:=tmp;            for i:=1 to n do                for j:=1 to n do                    for k:=1 to n do                        tmp2[i,j]:=(tmp2[i,j]+tmp[i,k]*b[k,j]) mod m;            for i:=1 to n do                for j:=1 to n do                begin                    b2[i,j]:=0;                    b3[i,j]:=0;                    matrixpower[i,j]:=0;                end;            for i:=1 to n do                for j:=1 to n do                    for k:=1 to n do                        b2[i,j]:=(b2[i,j]+b[i,k]*b[k,j]) mod m;            for i:=1 to n do                for j:=1 to n do                    for k:=1 to n do                        matrixpower[i,j]:=(matrixpower[i,j]+tmp2[i,k]*a[k,j]) mod m;            for i:=1 to n do                for j:=1 to n do                    matrixpower[i,j]:=(matrixpower[i,j]+a[i,j]) mod m;            for i:=1 to n do                for j:=1 to n do                begin                    for k:=1 to n do                        b3[i,j]:=(b3[i,j]+b2[i,k]*a[k,j]) mod m;                    b3[i,j]:=b3[i,j] mod m;                end;            b:=b3;        end;    end;end;begin    readln(n,k,m);    for i:=1 to n do        for j:=1 to n do            read(a[i,j]);    ans:=matrixpower(a,tmp,k);    for i:=1 to n do    begin        for j:=1 to n-1 do            write(ans[i,j],' ');        writeln(ans[i,n]);    end;end.