1.算法描述
由于矩阵乘法具有结合律,因此A^4 = A * A * A * A = (A*A) * (A*A) = A^2 * A^2。我们可以得到这样的结论:当n为偶数时,A^n = A^(n/2) * A^(n/2);当n为奇数时,A^n = A^(n/2) * A^(n/2) * A (其中n/2取整)。这就告诉我们,计算A^n也可以使用二分快速求幂的方法。例如,为了算出A^25的值,我们只需要递归地计算出A^12、A^6、A^3的值即可[1]。
在运算过程中,为了避免高精度运算,可采用取模(具体原因参看[1, 2])。
2.Referrence
[1] Matrix67, 十个利用矩阵乘法解决的经典题目
[2] Matrix67, 同余运算及其基本性质
[3]shiwei408,矩阵快速幂 poj3070 3233 3735 3150
3.问题
3.1 POJ 3070
求Fibonacci数Fn=Fn− 1+Fn− 2forn≥ 2 对10000取余。
题目已经给出了优化解决方案:对矩阵求n次幂,取n次幂的[0][1]元素,即为Fn 。
源代码:
3070 | Accepted | 116K | 0MS | C | 1178B | 2013-08-25 22:27:23 |
复制代码
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#include "stdio.h" #include "string.h" #define Mod 10000 typedef struct { int arr[2][2]; }Matrix; /*multiply the two matrices*/ Matrix mMulti(Matrix a, Matrix b) { int i, j, k; Matrix d; memset(d.arr,0,sizeof(d.arr)); for(i=0;i<2;i++) for(j=0;j<2;j++) { for(k=0;k<2;k++) d.arr[i][j]+=a.arr[i][k]*b.arr[k][j]; d.arr[i][j]%=Mod; } return d; } /*the k-th power of matrix a*/ Matrix mPow(Matrix a, int k) { Matrix modd,meven; modd.arr[0][0]=1; modd.arr[0][1]=0; //matrix modd is initialized to be a unit matrix modd.arr[1][0]=0; modd.arr[1][1]=1; meven=a; //matrix meven is initialized if(k==0) return modd; else { while(k!=1) { if(k&1) { k--; modd=mMulti(modd,meven); } else { k/=2; meven=mMulti(meven,meven); } } return mMulti(modd,meven); } } int main() { int n; Matrix fibonacci, goal; fibonacci.arr[0][0]=1; fibonacci.arr[0][1]=1; fibonacci.arr[1][0]=1; fibonacci.arr[1][1]=0; while(scanf("%d",&n)&&n!=-1) { goal=mPow(fibonacci,n); printf("%dn",goal.arr[0][1]%Mod); } return 0; }
3.2 POJ 3233
求n×n方阵A+A2+A3+ … +Ak
在Discuss中给出了优化解决思路:
复制代码
1
2
3| A+A^2+A^3+……A^k | |A A| (k-1)次方 | A | | | = | | | | | E | |0 E| | E |
源代码:
3233 | Accepted | 284K | 719MS | C | 1451B | 2013-08-27 23:04:22 |
复制代码
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#include "stdio.h" #include "string.h" #define MAX 60 typedef struct { int arr[MAX][MAX]; }Matrix; int n,m; /*multiplication*/ Matrix mMulti(Matrix a, Matrix b) { int i, j, k; Matrix d; memset(d.arr,0,sizeof(d.arr)); for(i=0;i<2*n;i++) for(j=0;j<2*n;j++) for(k=0;k<2*n;k++) d.arr[i][j]=(d.arr[i][j]+a.arr[i][k]*b.arr[k][j])%m; return d; } /*the k-th power of matrix a*/ Matrix mPow(Matrix a,int k) { int i; Matrix modd,meven; memset(modd.arr,0,sizeof(modd.arr)); //matrix modd is initialized to be a unit matrix for(i=0;i<n;i++) modd.arr[i][i]=1; meven=a; //matrix meven is initialized if(k==0) return modd; else { while(k!=1) { if(k&1) { k--; modd=mMulti(modd,meven); } else { k/=2; meven=mMulti(meven,meven); } } return mMulti(modd,meven); } } int main() { int i,j,k; Matrix A,B,C; scanf("%d%d%d",&n,&k,&m); memset(A.arr,0,sizeof(A.arr)); for(i=0;i<n;i++) for(j=0;j<n;j++) scanf("%d",&A.arr[i][j]); memcpy(B.arr,A.arr,sizeof(A.arr)); for(i=0;i<n;i++) for(j=0;j<n;j++) A.arr[i][j+n]=A.arr[i][j]; for(i=n;i<2*n;i++) { A.arr[i][i]=1; B.arr[i][i-n]=1; } C=mMulti(mPow(A,k-1),B); /*print the sum*/ for(i=0;i<n;i++) { for(j=0;j<n-1;j++) printf("%d ",C.arr[i][j]); printf("%dn",C.arr[i][n-1]); } return 0; }
最后
以上就是幸福枕头最近收集整理的关于【算法】矩阵快速求幂1.算法描述2.Referrence3.问题的全部内容,更多相关【算法】矩阵快速求幂1内容请搜索靠谱客的其他文章。
本图文内容来源于网友提供,作为学习参考使用,或来自网络收集整理,版权属于原作者所有。
发表评论 取消回复