Matrix Power Series 矩阵快速幂 + 因式分解（分治）

最新推荐文章于 2024-08-30 20:59:12 发布

13vickey

最新推荐文章于 2024-08-30 20:59:12 发布

阅读量739

点赞数

CC 4.0 BY-SA版权

分类专栏：数学

本文链接：https://blog.youkuaiyun.com/vickey1018/article/details/27377597

数学专栏收录该内容

28 篇文章

订阅专栏

本文介绍了一种使用矩阵快速幂算法解决特定矩阵求和问题的方法。对于给定的n×n矩阵A和正整数k，算法高效计算从A到Ak的所有矩阵相加的和，并输出结果模m后的矩阵。

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

Description

Given a n × n matrix A and a positive integer k, find the sum S = A + A² + A³ + … + A^k.

Input

The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 10⁹) and m (m < 10⁴). Then follow n lines each containing n nonnegative integers below 32,768, giving A’s elements in row-major order.

Output

Output the elements of S modulo m in the same way as A is given.

思路 A1+A2+A3+A4+A5+A6= A1+A2+A3 + A3(A1+A2+A3)

= (A1+A2+A3) *（ A3+I) // I 是单位矩阵；

=（ A1*（A1+I）+A3) * （A3+ I ）;

采用分治的思想，不断二分……

#include<stdio.h>
#include<string.h>
#include<iostream>
#include<algorithm>
#include<set>
#include<queue>
#include<map>
#include<vector>
using namespace std;
int k,n,mod;
struct node{
  int w[40][40];
  void init(){
      for(int i=0;i<n;i++)
        for(int j=0;j<n;j++)
          scanf("%d",&w[i][j]);
  }
  void prf(){
    for(int i=0;i<n;i++){
      for(int j=0;j<n-1;j++) printf("%d ",w[i][j]);
      printf("%d\n",w[i][n-1]);
    }
  }
  node operator * (const node &b) const{
     node cur;
     memset(cur.w,0,sizeof(cur.w));
     for(int i=0;i<n;i++){
       for(int j=0;j<n;j++){
         for(int k=0;k<n;k++){
            cur.w[i][j]=(cur.w[i][j]+(w[i][k]*b.w[k][j])%mod)%mod;
         }
       }
     }
     return cur;
  }
  node operator + (const node &b) const{
    node cur ;
    for(int i=0;i<n;i++)
      for(int j=0;j<n;j++) cur.w[i][j]=(w[i][j]+b.w[i][j])%mod;
    return cur;
  }
};
node a,ans,I,tp;
node check(int x){
    node cur=I,p=a;
    while(x){
       if(x&1) cur=cur*p;
       p=p*p;
       x>>=1;
    }
    return cur;
}
void dfs(int num){
    if(num==1)  return;
    dfs(num/2);
    tp=check(num/2);
    ans=ans*(tp+I);
    if(num&1){
      tp=check(num);
      ans=ans+tp;
    }
}
int  main()
{
//    freopen("in.in","r",stdin);
    memset(I.w,0,sizeof(I.w));
    while(~scanf("%d%d%d",&n,&k,&mod)){
       a.init();
       ans=a;
       for(int i=0;i<n;i++) I.w[i][i]=1;
       dfs(k);
       ans.prf();
    }
    return 0;
}