poj 3233 Matrix Power Series 矩阵快速幂

缘起

日常浪费生命~ poj 3233 Matrix Power Series

分析

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
给定n*n的矩阵A和正整数k和m, 求矩阵A的幂次的和
S=A+A^2+...+A^k
输出S的各元素对M取余后的结果

【输入】
n,k,m(三个正整数),然后n行每行n个[0,32768]范围内的正整数.
n ≤ 30,k ≤ 10^9,m < 10^4

【输出】
答案

【样例输入】
2 2 4
0 1
1 1

【样例输出】
1 2
2 3

【限制】
3s

首先就闻到了矩阵快速幂的味道~ 但是等等~ 能直接用矩阵快速幂吗? 显然不能. 因为计算A^k的话, 复杂度是O(n^3*logk)的. 这个复杂度在本题数据规模下是可以承受的. 但是如果再做矩阵加法的话, 复杂度就变成O(n^3*k)的了. k高达10亿. 这就是不能承受的了. 所以得另想办法.

用高代里面的分块矩阵即可.

分块矩阵为

1
2
3
4
{
{A,0}
{I,I}
}

则A^{k}的左下角的n*n的块就是I+A+…+A^{k-1},最后再乘以

1
2
3
4
{
{A, 0}
{0, 0}
}

即可,答案在左下角分块. 复杂度依旧是 O(n^3logk)的

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
//#include "stdafx.h"

#include <stdio.h>
#include <string.h>
//#define LOCAL

int n,k,m, a[65][65],t[65][65]; // t用于保存原a,就是 {{A,0},{0,0}}
void m_mul(int a[65][65], int b[65][65])
{
int ans[65][65];
memset(ans, 0, sizeof(ans));
for (int i = 0; i<2*n; i++)
{
for (int j = 0; j<2*n; j++)
{
for (int k = 0; k<2*n; k++)
{
ans[i][j]+=a[i][k]*b[k][j];
ans[i][j]%=m;
}
ans[i][j]%=m;
}
}
memcpy(a, ans, sizeof(ans)); // 切记这里不能写sizeof(a), a只是行指针, 要写sizeof(ans)
}

void ksm(int a[65][65], int b) // a^b快速幂
{
int ans[65][65];
for (int i = 0;i<2*n; i++)
{
for (int j = 0; j<2*n;j++)
{
ans[i][j]=i==j;
}
}
while(b)
{
if (b&1)
{
m_mul(ans, a);
}
if (b>1)
{
m_mul(a,a);
}
b>>=1;
}
memcpy(a,ans, sizeof(ans));
}

int main()
{
#ifdef LOCAL
freopen("d:\\data.in", "r", stdin);
#endif
scanf("%d%d%d", &n,&k,&m);
for (int i = n;i<2*n; i++)
{
for (int j = 0; j<n; j++)
{
a[i][j] = i-n==j;
}
}
for (int i = n; i<2*n; i++)
{
for (int j = n; j<2*n;j++)
{
a[i][j]=i==j;
}
}
for (int i = 0; i<n; i++)
{
for (int j = 0; j<n; j++)
{
scanf("%d", a[i]+j);
}
} // 初始化分块矩阵
memset(t, 0, sizeof(t));
for (int i = 0; i<n;i++)
{
for (int j = 0; j<n;j++)
{
t[i][j]=a[i][j];
}
}
ksm(a,k);
m_mul(a,t);
for (int i = n; i<2*n; i++)
{
for (int j = 0; j<n; j++)
{
if (j)
{
putchar(' ');
}
printf("%d", a[i][j]);
}
puts("");
} // 输出左下角分块
return 0;
}

ac情况

Status Accepted
Time 625ms
Memory 176kB
Length 1438
Lang C++
Submitted 2019-09-20 09:26:46
Shared
RemoteRunId 20877592