poj 3233 Matrix Power Series（矩阵快速幂+二分求和）_matrix power series vjudge-优快云博客

本文介绍了一种利用矩阵快速幂技巧高效求解特定形式矩阵和的方法。通过将原问题转化为递归形式，采用分治策略，可以显著降低计算复杂度。适用于解决形如(A+A^2+...+A^K)%m的矩阵求和问题。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题意：给定一个n*n的矩阵A，要求(A+A^2+....A^K)%m后的矩阵

解析：对于任意的A^x，我们都能够利用矩阵快速幂求出，但是我们现在要求的是和。

仔细观察整个式子，那么我们可以对原式进行变形
如果k为偶数，那么(A+A^2+....A^K) = (A+...+A^K/2)+A^K/2*(A+...+A^K/2)
如果k为奇数，那么(A+A^2+....A^K) = (A+...+A^K/2)+A^K/2*(A+...+A^K/2)+A^k

那么对于上面的式子的变形，就是二分的思想，那么我们可以利用二分来求和，然后对于单个的矩阵的x次方我们利用快速幂

#include <functional>
#include <algorithm>
#include <iostream>
#include <fstream>
#include <sstream>
#include <iomanip>
#include <numeric>
#include <cstring>
#include <cassert>
#include <cstdio>
#include <string>
#include <vector>
#include <bitset>
#include <queue>
#include <stack>
#include <cmath>
#include <ctime>
#include <list>
#include <set>
#include <map>
using namespace std;
const int N=35;

int n,m,k;
struct node
{
    int a[N][N];
    node  operator*(const node &b) const
    {
       node temp;
       for(int i=0;i<n;i++)
       {
          for(int j=0;j<n;j++)
          {
              temp.a[i][j]=0;
              for(int k=0;k<n;k++)
                 temp.a[i][j]+=a[i][k]*b.a[k][j]%m;
              temp.a[i][j]%=m;
          }
       }
       return temp;
    }
     node  operator+(const node &b) const
    {
       node temp;
       for(int i=0;i<n;i++)
          for(int j=0;j<n;j++)
             temp.a[i][j]=(a[i][j]+b.a[i][j])%m;
       return temp;
    }
};

node poww(node x,int t)
{
   node r;
   memset(r.a,0,sizeof(r.a));
   for(int i=0;i<n;i++) r.a[i][i]=1;
   while(t)
   {
       if(t&1) r=r*x;
       x=x*x;
       t>>=1;
   }
   return r;
}

node solvee(node x,int k)
{
    node temp,t;
     if(k == 1)
        return x;
    memset(t.a,0,sizeof(t.a));
    for(int i=0;i<n;i++) t.a[i][i]=1;
    temp=(poww(x,k/2)+t)*solvee(x,k/2);
    if(k&1) temp=temp+poww(x,k);
    return temp;
}
int main()
{
    int i,j;
    node x,res;
    while(~scanf("%d%d%d",&n,&k,&m))
    {
        for(i=0;i<n;i++)
         for(j=0;j<n;j++)
           scanf("%d",&x.a[i][j]);
        res=solvee(x,k);
        for(i=0;i<n;i++)
        {
            cout<<res.a[i][0]%m;
            for(j=1;j<n;j++)
             cout<<" "<<res.a[i][j]%m;
             cout<<endl;
        }
    }
    return 0;
}