# POJ 3233 二分求等比数列 矩阵快速幂

xiaoxiao2022-05-22  9

（3）当时，那么有

#include <cstdio> #include <iostream> #include <cstring> #define sf scanf #define pf printf using namespace std; const int maxn = 35; int n; __int64 mod; struct Mat{ __int64 C[maxn][maxn]; Mat(){memset(C,0,sizeof C);for(int i = 0;i < n;++i) C[i][i] = 1;} void out(){ for(int i = 0;i < n;++i){ pf("%I64d",C[i][0]); for(int j = 1;j < n;++j) pf(" %I64d",C[i][j]); pf("\n"); } } }; Mat operator*(const Mat& A,const Mat& B){ Mat ret; for(int i = 0;i < n;++i){ for(int j = 0;j < n;++j){ ret.C[i][j] = 0; for(int k = 0;k < n;++k){ ret.C[i][j] += A.C[i][k]*B.C[k][j]; ret.C[i][j] %= mod; } } } return ret; } Mat operator+(const Mat& A,const Mat& B){ Mat ret; for(int i = 0;i < n;++i) for(int j = 0;j < n;++j) {ret.C[i][j] = A.C[i][j] + B.C[i][j];ret.C[i][j] %= mod;} return ret; } Mat operator^(Mat A,int k){ Mat ret; while(k){ if(k & 1) ret = ret * A; A = A * A; k = k >> 1; } return ret; } Mat getSum(const Mat A,int k){ if(k == 1){ return A; } else{ Mat ret = getSum(A,k / 2); if(k & 1){ Mat tmp = A ^ (k / 2 + 1); return ret + (ret * tmp) + tmp; } else{ return ret + ret * (A ^ (k / 2)); } } } int main(){ int k; while( ~sf("%d%d%I64d",&n,&k,&mod)){ Mat st; for(int i = 0;i < n;++i) for(int j = 0;j < n;++j) {sf("%I64d",&st.C[i][j]);st.C[i][j] %= mod;} getSum(st,k).out(); } return 0; }