摘要 初识矩阵快速幂已是很久之前的事情了,我的忘性有那么好,现在自然是忘得一干二净啦,简要写一下矩阵快速幂相关,以备后用.
事情要从很久很久以前说起,我记得那是第一节算法课,老师说了这样一个题目:某人下楼梯可以一次走两阶,也可以一次走一阶,现在他在第 k 阶楼梯,问题是他要走到第0阶有多少种不同的走法?显然,这是fibonacci数列,即:
$$ f(n) = \begin{cases} 1& \text{n=1}\\ 1& \text{n=2}\\ f(n-1) + f(n-2)& \text(n \geq 3) \end{cases} $$
这个问题并不复杂,我们甚至能求出它的通项公式:
$$ \frac{1}{\sqrt{5}}[(\frac{1+\sqrt{5}}{2})^{n}- (\frac{1-\sqrt{5}}{2})^{n}] $$
不过,我们知道斐波那契数列的数字总是整数,这个通项公式并没有显式的表明这是一个整数,而且,即使利用这个通项公式求解 $ f(n) $ 依然要计算至少 $ n $次,反倒不如直接递推!
事实上我们能在 log 级的时间复杂度内解决这个问题。
首先,我们构造矩阵:
$$ \begin{equation} \left| \begin{array}{cc} 1 & 1 \\ 1 & 0 \end{array} \right| \end{equation} $$
$$ \begin{equation} \left| \begin{array}{c} f(2) \\ f(1) \end{array} \right| \end{equation} $$
矩阵一和矩阵二相乘会得到一个新的矩阵:
$$ \begin{equation} \left| \begin{array}{c} f(2) + f(1) \\ f(2) \end{array} \right| \end{equation} $$
显然,矩阵三等价于:
$$ \begin{equation} \left| \begin{array}{c} f(3) \\ f(2) \end{array} \right| \end{equation} $$
如果把得到的矩阵三再和矩阵一相乘呢?我们会得到含有f(4)和f(3)的矩阵,依此类推,我们总能得到这样的矩阵:
$$ \begin{equation} \left| \begin{array}{c} f(n) \\ f(n-1) \end{array} \right| \end{equation} $$
记这个矩阵为$ M_{n} $ ,矩阵一为$ M_{0} $,矩阵二为$ M_{2} $则有:
$$ M_{n} = M_{0}*M_{n-1} = M_{0}*M_{0}*M_{n-2} = ... = M_{0}^{n-2}*M_{2} $$
作为一个"伪ACMer"我还是知道快速幂取模这个东西的:
快速幂取模能在log级的时间复杂度内解决 $ A^{n} mod p $的问题, A , p , n 均为正整数
新得到的式子和快速幂取模能解决的问题竟如此相似,事实上是完全一样的,如果你了解快速幂取模的过程的话,唔---还是简要的说一下吧!
代码是这样的:
#define LL long longLL Quick_Mod(LL a,LL k,LL p){ LL ans = 1 ; a %= p ; while (k) { if (k&1) ans = ans*a%p; a = a*a%p ; k = k>>1 ; } return ans ;}
思路是这样的:如果 n 是奇数,那么$ answer = A^{n-1}*A mod p $,如果 n 是偶数,那么$ answer = (A^{2})^{(n/2)} mod p $。这样,就会把时间复杂度降到log级别
咦><整数和矩阵的区别很大么?把快速幂取模代码中的乘号重载成矩阵的乘法不就ok了吗?没错,就是这样!!!现在,已经成功的解决了这个问题!!
顺带说一下昨天 :
对,就是这个,读完题目我就知道要构造矩阵了,瞎写了一个多小时,都没能正确的构造矩阵,也因此才写了这篇随笔。
分析 我们给它的递推式换一种写法:
$$ f_{n} = \begin{cases} f_{2}^{0}-->0& \text{n=1}\\ f_{2}^{1}-->1& \text{n=2}\\ f_{2}*f_{n-1}^{c}*f_{n-2}-->(1+c*g_{n-1}+g_{n-2})& \text(n \geq 3) \end{cases} $$
这样写的意思就是都表示为$ f_{2} $的某次方,最后得到$ g_{n} = 1 + c*g_{n-1} + g_{n-2} $且有了$ g_{1} = 0, g_{2} = 1 $
于是构造矩阵:
$$ \left| \begin{array}{ccc} c & 1 & 1 \\ 1 & 0 & 0 \\ 0 & 0 & 1 \end{array} \right| $$
$$ \left| \begin{array}{c} g_{2} \\ g_{1} \\ 1 \end{array} \right| $$
可以发现,相乘一次就能得到$ g_{3} $。同样,到这里这个问题就解决了,下面是代码。
我的代码在HDoj上提交可以 AC ,但是直觉告诉我应该有点什么问题==希望路过的同学能帮忙指出!
1 #include2 #include 3 #include 4 #include 5 #include 6 using namespace std ; 7 #define LL long long 8 #define rep(i,n) for (int i = 1 ; i <= n ; ++ i) 9 const int n = 3 ;10 struct Matrix11 {12 LL mat[n<<1][n<<1] ;13 }e;14 15 void InitMatrix()16 {17 rep(i,n) rep(j,n) e.mat[i][j] = (i == j) ;18 }19 20 Matrix mul(Matrix a,Matrix b,LL Mod)21 {22 Matrix c ;23 rep(i,n) {24 rep(j,n) {25 c.mat[i][j] = 0 ;26 rep(k,n) {27 c.mat[i][j] += (a.mat[i][k]*b.mat[k][j])%Mod ;28 while (c.mat[i][j] >= Mod) c.mat[i][j] -= Mod ;29 }30 }31 }32 return c ;33 }34 35 Matrix Quick_Mod(Matrix a,LL k,LL Mod)36 {37 Matrix c = e ;38 while (k) {39 if (k&1) c = mul(c,a,Mod) ;40 a = mul(a,a,Mod) ;41 k = k/2 ;42 }43 return c ;44 }45 46 LL Quick_Mod(LL a,LL k,LL Mod)47 {48 LL ans = 1 ;49 a %= Mod ;50 while (k) {51 if (k&1) ans = ans*a%Mod ;52 a = a*a%Mod ;53 k = k>>1 ;54 }55 return ans ;56 }57 58 59 int main()60 {61 // freopen("in.txt","r",stdin) ;62 LL N , A , B , C , Mod ;63 int T ;64 Matrix c ;65 InitMatrix() ;66 scanf("%d",&T) ;67 while (T --) {68 scanf("%I64d%I64d%I64d%I64d%I64d",&N,&A,&B,&C,&Mod) ;69 if (A%Mod == 0) {70 if (N == 1) puts("1") ;71 else puts("0") ;72 continue ;73 }74 if (N < 3) {75 if (N == 1) printf("1\n") ;76 else printf("%I64d\n",Quick_Mod(A,B,Mod)) ;77 continue ;78 }79 memset(c.mat,0,sizeof(c.mat)) ;80 c.mat[1][1] = C ;81 c.mat[1][2] = c.mat[1][3] = c.mat[2][1] = c.mat[3][3] = 1 ;82 c = Quick_Mod(c,N-2,Mod-1) ;83 LL k = c.mat[1][1] + c.mat[1][3] ;84 k %= (Mod-1) ;85 k = k*B%(Mod - 1) ;86 printf("%I64d\n",Quick_Mod(A,k,Mod)) ;87 }88 return 0 ;89 }
待学习链接:
matrix67的
~end~