文章目錄

• 概述
• 快速冪
• • 解析
  • 代碼
• 矩陣運算定義
• • 加法
  • 乘法
  • 單位矩陣
• 一、斐波拉契（基礎模板）
• • 題目描述
  • 解析
  • 代碼
• 二、行為方案（實際應用）
• • 題目描述
  • 解析
  • 代碼
• 三、矩陣求和（子矩陣作為矩陣元素）
• • 題目描述
  • 解析
  • 代碼
• 四、最短路徑（對矩陣乘法靈活的定義方式）
• • 題目描述
  • 解析
  • 代碼：
• thanks for reading！

概述

矩陣快速冪，顧名思義，就是矩陣乘法與快速冪的結合,可以將時間復雜度優化到log級別。
當數據范圍中圖的尺寸較小，而時間、變換次數、步數等的范圍較大（通常在10^9左右）時，可以考慮矩陣快速冪問題
矩陣快速冪的關鍵在于構造轉移矩陣

快速冪

解析

由于

x^a * x^b = x^(a+b)

逆向來看,對于一個較大的指數k，我們可以進行拆分：

x^k = x^a1 * x^a2 *…(a1+a2+…=k)

由于二進制拆分的唯一存在性，我們就可以把k拆成若干個2的冪數作為a序列，將其對應的乘冪求出，并乘在一起就能得到x^k了
時間復雜度log k

代碼

ll ksm(ll x, int w) {//求x的w次冪ll ans = 1;while (w) {if (w & 1)ans = (ll)(ans * x) % mod;//mod是取模用的x = (x * x) % mod;w >>= 1;}return ans; }

矩陣運算定義

加法

對于兩個大小完全相同的矩陣，定義二者相加為每個對應元素相加（比較簡單，不細說了）

乘法

對于兩個矩陣AB，若A的行數等于B的列數，則二者可以相乘
換言之，ab的矩陣A可以和bc的矩陣B相乘，結果是一個a*c的矩陣C，單次相乘時間復雜度為abc 的乘積
對于每一個c中的元素：

c[i][j]=∑a[i][p]+b[p][j](1<=p<=b)

簡單說就是頂針式的乘積累加
容易證明，矩陣乘法滿足結合律，但是不滿足交換律
對于一個行列相等的方陣，可以不斷累乘自身

單位矩陣

如果一個矩陣的結構為：
1 0 0 0 0
0 1 0 0 0
0 0 1 0 0
0 0 0 1 0
0 0 0 0 1
則稱它為一個單位矩陣
對于一個a*b的矩陣A，它與一個邊長為b的矩陣相乘，結果還是A本身
單位矩陣在矩陣運算中起著單位1的作用

矩陣簡單的運算就是這些，接下來我們通過一些例題由淺入深，看看矩陣快速冪如何應用吧

一、斐波拉契（基礎模板）

萬物之始

題目描述

斐波拉契的定義為：

f[1]=f[2]=1; f[n]=f[n-1]+f[n-2](n>=3)

（不會真的有人不知道吧）
現在要你求出第k項取模后的結果
(k<=2^63)

解析

我們發現：
對于一個行矩陣：

f[n] f[n-1]

嘗試構造出一個轉移矩陣M：
（轉移矩陣是矩陣題目的靈魂，多為方陣）

1 1
1 0

那么兩者相乘后就會得到：

f[n]+f[n-1] f[n]

根據斐波拉契的定義，上面的結果就是：

f[n+1] f[n]

所以，要求第k項的值，就只需要把矩陣f[2] f[1]乘上k-2次轉移矩陣即可
這樣就可以用快速冪提速了
（由于本題原矩陣過于簡單，其實代碼實現時我們都不需要讓那個行矩陣出現）

代碼

#include<bits/stdc++.h> using namespace std; #define ll long long const int mod=1e9+7; ll a,b,c,k; ll ans[4][4],res[4][4],trans[4][4]; void mul1(){//res*resmemset(trans,0,sizeof(trans));for(int i=1;i<=2;i++){for(int j=1;j<=2;j++){for(int p=1;p<=2;p++){trans[i][j]+=res[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=1;i<=2;i++){for(int j=1;j<=2;j++){res[i][j]=trans[i][j];}} } void mul2(){//res*ansmemset(trans,0,sizeof(trans));for(int i=1;i<=2;i++){for(int j=1;j<=2;j++){for(int p=1;p<=2;p++){trans[i][j]+=ans[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=1;i<=2;i++){for(int j=1;j<=2;j++){ans[i][j]=trans[i][j];}} } void ksm(ll w){while(w){if(w&1) mul2();mul1();w>>=1;}return; } void print(){for(int i=1;i<=2;i++){for(int j=1;j<=2;j++){printf("%d ",ans[i][j]);}printf("\n");} } int main() {ans[1][1]=ans[2][2]=1;res[1][1]=res[1][2]=res[2][1]=1;scanf("%lld",&a);if(a<=2) {printf("1");return 0;}ksm(a-2); // print();printf("%lld",(ans[1][1]+ans[2][1])%mod);return 0; }

二、行為方案（實際應用）

題目描述

解析

定義dp[i][j][t]為經過t時間從i走到j的方案數
不難想到轉移

dp[i][j][t]=∑dp[i][p][t-1]*dp[p][j][1]

我們發現出現了矩陣乘法標志性的頂針結構
也就是說，想要從t-1轉移到t，只需要把t-1狀態的矩陣與狀態1的矩陣相乘就行了
換言之，詢問t狀態的矩陣，就是求狀態1矩陣的t次冪
這樣在提速的同時，也就把第三維空間優化掉了
這樣就轉化為快速冪問題了
接下來就是原始矩陣的構造問題

首先肯定要把存在邊的點連上：

dp[from][to]=1

還可以停留在原點，也就相當于一個自環：

dp[i][i]=1;(1<=i<=n)

對于自爆，我們可以把0結點當成自爆的狀態，在任何一個點均可自爆，所以：

dp[i][0]=1(1<=i<=n)

由于自爆后不能向任何狀態轉移了，所以0結點沒有任何出邊
最后的答案就是∑dp[1][i] (0<=i<=n）
這樣本題就結束了

代碼

#include<bits/stdc++.h> using namespace std; #define ll long long const int mod=2017; const int N=150; ll a,b,c,k; int n,m; ll ans[N][N],res[N][N],trans[N][N]; void mul1(){//res*resmemset(trans,0,sizeof(trans));for(int i=0;i<=n;i++){for(int j=0;j<=n;j++){for(int p=0;p<=n;p++){trans[i][j]+=res[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=0;i<=n;i++){for(int j=0;j<=n;j++){res[i][j]=trans[i][j];}} } void mul2(){//res*ansmemset(trans,0,sizeof(trans));for(int i=0;i<=n;i++){for(int j=0;j<=n;j++){for(int p=0;p<=n;p++){trans[i][j]+=ans[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=0;i<=n;i++){for(int j=0;j<=n;j++){ans[i][j]=trans[i][j];}} } void ksm(ll w){while(w){if(w&1) mul2();mul1();w>>=1;}return; } void print(){for(int i=0;i<=n;i++){for(int j=0;j<=n;j++){printf("%d ",ans[i][j]);}printf("\n");} } int main() {scanf("%d%d",&n,&m);for(int i=0;i<=n;i++) ans[i][i]=1;for(int i=1;i<=m;i++){int a,b;scanf("%d%d",&a,&b);res[a][b]=res[b][a]=1;}for(int i=0;i<=n;i++) res[i][i]=1;for(int i=1;i<=n;i++) res[i][0]=1;int t;scanf("%d",&t);ksm(t);int tot=0;for(int i=0;i<=n;i++){tot+=ans[1][i];} // print();printf("%d",tot%mod);return 0; }

三、矩陣求和（子矩陣作為矩陣元素）

題目描述

解析

對于矩陣A，我們構造轉移矩陣：

AI

0	I

（I表示與A邊長相等的單位矩陣，矩陣套矩陣表示把小矩陣值套進大矩陣里）
比如，當A為：

2 3
1 5

時，轉移矩陣為：

2 3 1 0
1 5 0 1
0 0 1 0
0 0 0 1

使它與自身相乘,得到：

A^2I+A

0	I

再乘一次：

A^3I+A+A^2

0	I

這樣我們就可以得到規律了
只需要把這個轉移矩陣構造出來后乘上k次冪在把左上角減去單位矩陣輸出即可
注意：取模減法需要判斷一下會不會減成負的！
（連ybt評測的答案似乎都沒有注意。。。）

代碼

#include<bits/stdc++.h> using namespace std; #define ll long long int mod=2017; const int N=150; ll a,b,c,k; int n,m; ll ans[N][N],res[N][N],trans[N][N]; void mul1(){//res*resmemset(trans,0,sizeof(trans));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){for(int p=1;p<=n;p++){trans[i][j]+=res[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){res[i][j]=trans[i][j];}} } void mul2(){//res*ansmemset(trans,0,sizeof(trans));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){for(int p=1;p<=n;p++){trans[i][j]+=ans[i][p]*res[p][j];trans[i][j]%=mod;}}}for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){ans[i][j]=trans[i][j];}} } void ksm(ll w){while(w){if(w&1) mul2();mul1();w>>=1;}return; } void print(){for(int i=1;i<=n;i++){for(int j=1+n;j<=2*n;j++){if(i+n==j) ans[i][j]--; // if(ans[i][j]<0) ans[i][j]+=mod;printf("%lld ",ans[i][j]);}printf("\n");} } int main() {scanf("%d%d%d",&n,&m,&mod);for(int i=1;i<=n;i++){for(int j=1;j<=n;j++) scanf("%lld",&res[i][j]);}for(int j=n+1;j<=2*n;j++){res[j][j]=res[j-n][j]=1;}n<<=1;for(int i=1;i<=n;i++) ans[i][i]=1;ksm(m+1);n>>=1;print();return 0; }

四、最短路徑（對矩陣乘法靈活的定義方式）

題目描述

解析

定義dp[i][j][k]：從i到j，經過k條邊的最短路徑
易得轉移方程：

dp[i][j][k]=min(dp[i][j][k],dp[i][p][k-1]+dp[p][j][1];

我們發現出現了頂針結構，但似乎運算從乘積累加變成了求min
那么把運算的定義改一下不就好了！
（這樣修改之后不再有單位矩陣的概念，好在本題也不需要）
還有一些細節問題：

1.在原始矩陣中dp[i][i]應該為正無窮而不是0，換言之，經過一條邊走到原處的方案在沒有自環的情況下應該是不可能的
2.由于邊數<=100，所以點不會超過200個，但編號會達到1000。如果直接用原編號，會超時。所以本題需要對點進行離散化處理

代碼：

#include<bits/stdc++.h> using namespace std; #define ll long long int mod=2017; const int N=1500; int s,e,n,m,k; int mx; int a,b,c,d; ll ans[N][N],res[N][N],trans[N][N]; int id[N*100],tot=0; void mul1(){//res*resmemset(trans,0x3f,sizeof(trans));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){for(int p=1;p<=n;p++){trans[i][j]=min(trans[i][j],res[i][p]+res[p][j]); // trans[i][j]%=mod;}}}for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){res[i][j]=trans[i][j];}} } void mul2(){//res*ansmemset(trans,0x3f,sizeof(trans));for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){for(int p=1;p<=n;p++){trans[i][j]=min(trans[i][j],ans[i][p]+res[p][j]); // trans[i][j]%=mod;}}}for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){ans[i][j]=trans[i][j];}} } void ksm(ll w){while(w){if(w&1) mul2();mul1();w>>=1;}return; } void print(){for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){ // if(ans[i][j]<0) ans[i][j]+=mod;if(ans[i][j]>2e9) printf("-1 ");else printf("%lld ",ans[i][j]);}printf("\n");}printf("\n"); } void print_res(){for(int i=1;i<=n;i++){for(int j=1;j<=n;j++){ // if(ans[i][j]<0) ans[i][j]+=mod;if(res[i][j]>2e9) printf("-1 ");else printf("%lld ",res[i][j]);}printf("\n");} } int main() {scanf("%d%d%d%d",&k,&n,&s,&e);memset(res,0x3f,sizeof(res));memset(ans,0x3f,sizeof(ans));for(int i=1;i<=n;i++){scanf("%d%d%d",&c,&a,&b);if(!id[a]) id[a]=++tot;if(!id[b]) id[b]=++tot;a=id[a],b=id[b];res[a][b]=res[b][a]=c; // mx=max(mx,max(a,b)); // printf("a=%d b=%d v=%d\n",a,b,c);}n=tot; // print_res();for(int i=1;i<=n;i++) ans[i][i]=0;ksm(k); // mul2(); // print(); // mul2(); // print();printf("%lld",ans[id[s]][id[e]]);return 0; } /* 1 6 6 4 11 4 6 4 4 8 8 4 9 6 6 8 2 6 9 3 8 9 */

thanks for reading！

總結

以上是生活随笔為你收集整理的矩阵快速幂一篇通的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。