[原创]POJ 3233 Matrix Power Series 【矩阵快速幂+等比矩阵】

2016-07-17 12:56:31 Tabris_ 阅读数：295

博客爬取于2020-06-14 22:44:11
以下为正文

版权声明：本文为Tabris原创文章，未经博主允许不得私自转载。
https://blog.csdn.net/qq_33184171/article/details/51931853

题目连接 : http://poj.org/problem?id=3233

Matrix Power Series
Time Limit: 3000MS Memory Limit: 131072K
Total Submissions: 20125 Accepted: 8449
Description

Given a n × n matrix A and a positive integer k, find the sum S = A + A2 + A3 + … + Ak.

Input

The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 109) and m (m < 104). 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

题目大意 ： 就是求S %m；

题解 : 其实很容易知道A^k 之后只需加和就行了

但是直接加和还是不行 k的范围是在太大 会超时 所以就构造一个矩阵

因为S可以看成S=A(I+A(I+A(I+…A(I+A)))) (I是单位矩阵)

拿k=3举例S=A(I+A(I+A))

那么我们想，可不可以构造一个矩阵T使得TT（因为是k次幂）这样乘下去每次可以得到A(A+I)

那么肯定T有个两个元素就是A与I

那么假设：T={A I }
I I
那么T=TT={AA+II AI+II}
AI+II II+II
这样存在一个I（A+I）的式子 ，当T再乘以T的时候会出现A（A+I）

这个时候我们可以简化将T={A I}

                                        0   I

这样可以简化很多计算TT={AA AI+II}
0 I

那么容易得到T^(K+1)={A^(K+1) I+A+A^2+A^3+…+A^K}
0 I

这样我们只需要算T的k+1次幂就可以了

上文摘自 http://blog.csdn.net/yihuikang/article/details/7722765

实在懒得写这些 其实就是简单的等比矩阵求和 学过线性代数 就应该能构造出来矩阵的

错将输入 n,k,m 当成的 n,m,k 错了一上午 。。汗2333333333

附本提代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
# include <iostream>
# include <stdio.h>
# include <string.h>
# include <string>
# include <algorithm>
# include <math.h>

using namespace std;

/******************************/
# define LL long long int
# define _LL __int64
/*****************************/

/*
const int M = 100;
const int MOD = 1e4;
*/

# define M   n*2
# define MOD m

LL n,m,k;


struct Matrix
{
   // int r,c;   //C行R列
    LL m[66][66];
};

Matrix operator * (Matrix a,Matrix b)
{
    Matrix c;
    for(int i=0; i<M; i++) //初始化矩阵
        for(int j=0; j<M; j++)
            c.m[i][j]= 0;

    for(int k=0; k<M; k++)
        for(int i=0; i<M; i++) //实现矩阵乘法
        {
            if(a.m[i][k] <= 0)  continue;  //剪枝
            for(int j=0; j<M; j++)
            {
               // if(b.m[k][j] <= 0)    continue;  //剪枝
                c.m[i][j]=(c.m[i][j]+a.m[i][k]*b.m[k][j]+MOD)%MOD;
            }
        }
    return c;
}

Matrix operator ^ (Matrix a,LL b)
{
    Matrix c;
    for(int i=0; i<M; i++) //初始化单位矩阵
        for(int j=0; j<M; j++)
            c.m[i][j]= ( i == j );

    while(b)
    {
        if(b&1) c= c * a ;
        b >>= 1;
        a = a * a ;
    }

    return c;
}

Matrix operator + (Matrix a,Matrix b)
{
    for(int i=0; i<M; i++)
        for(int j=0; j<M; j++)
            a.m[i][j]=(a.m[i][j]+b.m[i][j]+MOD)%MOD;

    return a;
}

int main()
{
    ios::sync_with_stdio(false);

    while( cin>>n>>k>>m )
    {
        Matrix a;
        memset(a.m,0,sizeof(a.m));

        for(int i=0; i<n; i++)
            for(int j=0; j<n; j++)
                cin >> a.m[i][j];

        for(int i=0;i<n;i++)
            a.m[i+n][i+n]=a.m[i][i+n]=1;

        Matrix c;

        c = a ^ ( k + 1 );

        for(int i=0;i<n;i++)  //减掉单位矩阵。。
            c.m[i][i+n] = (c.m[i][i+n]-1+MOD)%MOD;

        for(int i=0; i<n; i++)
        {
            for(int j=0; j<n; j++)
            {
                if(j) cout<<" ";
                cout<<c.m[i][j+n];
            }
            cout<<endl;
        }

    }
    return 0;
}