Matrix Power Series

最新推荐文章于 2020-01-21 15:40:01 发布

原创最新推荐文章于 2020-01-21 15:40:01 发布 · 376 阅读

0 ·

CC 4.0 BY-SA版权

文章标签：

#分治 #矩阵乘法

Problems 专栏收录该内容

10 篇文章

订阅专栏

本文介绍了一种利用矩阵快速幂解决特定数学问题的方法：给定一个 n×n 的矩阵 A 和正整数 k，计算从 A 到 Ak 的所有矩阵幂次之和 S，并输出各元素对 m 取模后的结果。文章提供了两种解决方案，一种是通过分治思想进行优化，另一种则是构造一个元素为矩阵的矩阵来加速计算。

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

Description

Time Limit: 3000MS Memory Limit: 131072K
Given a $n*n$ matrix $A$ and a positive integer $k$ , find the sum $S = A + A^2 + A^3 + … + 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^9)$ and $m (m < 10^4)$ . 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.

Sample Input

2 2 4
0 1
1 1

Sample Output

1 2
2 3

Solution & Code

法一： $k$ 的范围是 $10^9$ ，考虑分治。
$S_i=A+A^2+...+A^i$
$S_{2i}=A^{i+1}+A^{i+2}+...+A^{2i}$
$S_{2i}=S_i*A^i$

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;

typedef long long ll;

const int maxr = 35;

int n, k;
ll m;

struct matrix{ll val[maxr][maxr];};
matrix U, V, Z, I;

matrix multiply(matrix A, matrix B){

    matrix C = Z;
    for(int i = 1; i <= n; ++i)
        for(int j = 1; j <= n; ++j)
            for(int k = 1; k <= n; ++k){
                C.val[i][j] += A.val[i][k] * B.val[k][j];
                C.val[i][j] %= m;
            }
    return C;
}

matrix add(matrix A, matrix B){

    matrix C = Z;
    for(int i = 1; i <= n; ++i)
        for(int j = 1; j <= n; ++j){
            C.val[i][j] = A.val[i][j] + B.val[i][j];
            C.val[i][j] %= m;
        }
    return C;
}

matrix powmod(matrix A, int x){

    matrix B = I;
    while(x){
        if(x & 1) B = multiply(B, A);
        x = x >> 1;
        A = multiply(A, A);
    }
    return B;
}

matrix calc(int x){

    if(x == 1) return U;
    int y = x >> 1;
    matrix L = calc(y);
    matrix R = multiply(L, powmod(U, y));
    if(y * 2 != x) return add(add(L, R), powmod(U, x));
        else return add(L, R);
}

int main(){

    cin >> n >> k >> m;
    for(int i = 1; i <= n; ++i){
        for(int j = 1; j <= n; ++j) cin >> U.val[i][j];
        I.val[i][i] = 1;
    }
    V = calc(k);
    for(int i = 1; i <= n; ++i){
        for(int j = 1; j <= n; ++j) printf("%lld ", V.val[i][j]);
        puts("");
    }

    return 0;
}