Before we talk about the fast power of a matrix, let's talk about what a matrix is.
About matrix
As the name suggests, a matrix is to fill some numbers into a rectangle in the following form:
So?
So let's talk about matrix multiplication.
Matrix multiplication
Matrix multiplication is to multiply two matrices, but how do we do it? Baidu learned from the book that matrix multiplication has the following definitions:
Let A and B be two matrices, let C=A × B. Then:
-
The number of columns of A must be equal to the number of rows of B;
-
Let A be n × The matrix of R, B is r × Matrix of m; Then the product C of A and B is an n × Matrix of m;
-
\(C_{i,j}=A_{i,1}×B_{1,j}+A_{i,2}×B_{2,j}+A_{i,3}×B_{3,j}+...+A_{i,n}×B_{n,j}=\sum^{n}_{k=1}A_{i,k}×B_{k,j}\)
For example:
So we can write the code:
struct matrix{
int n,m,g[N][N];//N is the number of rows, m is the number of columns, g is the matrix, and N is the size of the matrix
matrix operator*(const matrix&b)const{
matrix c;
memset(c.g,0,sizeof(c.g));//initialization
for(int i=1;i<=n;i++)
for(int j=1;j<=m;j++)
for(int k=1;k<=b.m;k++)
c.g[i][j]=(c.g[i][j]+g[i][k]*b.g[k][j]%mod)%mod;
return c;
}
}
matrix multiply(matrix a,matrix b){
matrix c;
memset(c.g,0,sizeof(c.g));
for(int i=1;i<=a.n;i++)
for(int j=1;j<=a.m;j++)
for(int k=1;k<=b.m;k++)
c.g[i][j]=(c.g[i][j]+a.g[i][k]*b.g[k][j]%mod)%mod;
return c;
}
So does matrix multiplication satisfy the operation law?
Yes, the matrix multiplication satisfies the associative law, which is proved as follows:
Let n-order matrix be:
\(A=(a_{i,j})\)
\(B=(b_{i,j})\)
\(C=(c_{i,j})\)
\(A\times{B}=(d_{i,j})\)
\(B\times{C}=(e_{i,j})\)
\((A\times{B})\times{C}=(f_{i,j})\)
$A\times{(B\times{C})}=(g_{i,j}) $
From the multiplication of the matrix
\(d_{i,j}=a_{i,1}\times{b_{1,j}}+a_{i,2}\times{b_{2,j}}+...+a_{i,n}\times{b_{n,j}}\)
\(e_{i,j}=b_{i,1}\times{c_{1,j}}+b_{i,2}\times{c_{2,j}}+...+b_{i,n}\times{c_{n,j}}\)
\(f_{i,j}=d_{i,1}\times{c_{1,j}}+d_{i,2}\times{c_{2,j}}+...+d_{i,n}\times{c_{n,j}}\)
\(g_{i,j}=a_{i,1}\times{e_{1,j}}+a_{i,2}\times{e_{2,j}}+...+a_{i,n}\times{e_{n,j}}\)
Therefore, for any \ (i,j=1,2,...,n \),
$f_{i,j}=d_{i,1}\times{c_{1,j}}+d_{i,2}\times{c_{2,j}}+...+d_{i,n}\times{c_{n,j}} $
\(=(a_{i,1}\times{b_{1,1}}+a_{i,2}\times{b_{2,1}}+...+a_{i,n}\times{b_{n,1}})\times{c_{1,j}}+(a_{i,1}\times{b_{1,1}}+a_{i,2}\times{b_{2,1}}+...+a_{i,n}\times{b_{n,1}})\times{c_{2,j}}\)
$+...+(a_{i,1}\times{b_{1,n}}+a_{i,2}\times{b_{2,n}+...+a_{i,n}}\times{b_{n,n}})\times{c_{n,j}} $
\(=a_{i,1}\times{(b_{1,1}\times{c_{1,j}}+b_{1,2}\times{c_{2,j}}+...+b_{1,n}\times{c_{n,j}})}+a_{i,2}\times{(b_{2,1}\times{c_{1,j}}+b_{2,2}\times{c_{2,j}}+...+b_{2,n}\times{c_{n,j}})}\)
$+...+a_{i,n}\times{(b_{n,1}\times{c_{1,j}}+b_{n,2}\times{c_{2,j}}+...+b_{n,n}\times{c_{n,j}})} $
$=a_{i,1}\times{e_{1,j}}+a_{i,2}\times{e_{2,j}}+...+a_{i,n}\times{e_{n,j}}=g_{i,j} $
So \ ((A\times{B})\times{C}=A\times{(B\times{C})}. \)
Proof from the Internet: https://www.cnblogs.com/Jakson/articles/4557558.html
What's the use? Let's get to the point.
Matrix fast power
For a \ (n × N \), we call it a square matrix. A square matrix can perform power operation, which is defined as \ (C=A^n \). Obviously, for any matrix, only if it is a square matrix can the power operation be carried out (otherwise N and m are not equal).
Because the matrix multiplication satisfies the associative law, we can use the fast power to calculate the power of the matrix.
The advantage of matrix multiplication is that it can store useful states in a matrix and obtain the value of DP once through one multiplication (DP state transition equation must be linear. Although the state transition equation in some problems is not linear, we can convert it into linear equation through certain transformation)
The general forms are row vector left multiplication matrix and column vector right multiplication matrix, for example
When implemented, it is generally supplemented by 0 in the vector to become a square matrix
Why must we be linear? See the following breakdown:
We know linear state transition equations, such as Fibonacci series
\(f_i=f_{i-1}+f_{i-2}\)
We can get the nth term of Fibonacci sequence by single recurrence
It can also be expressed as:
\(f_i=1\times{f_{i-1}}+1\times{f_{i-2}}\)
In matrix multiplication, we just use the last answer to transform the new state, so we can use matrix multiplication to optimize the linear DP
Combined fast power:
matrix operator^(int k){
matrix ans,x=*this;
memset(ans.g,0,sizeof(ans.g));
for(int i=1;i<=N;i++)ans.g[i][i]=1;
for(;k;k>>=1,x=x*x)if(k&1)ans=ans*x;
return ans;
}
matrix power(matrix a,int b){
matrix ans;
memset(ans.g,0,sizeof(ans.g));
for(int i=1;i<=N;i++)ans.g[i][i]=1;
while(b){
if(b&1)ans=multiply(ans,a);
a=multiply(a,a);
b>>=1;
}
return ans;
}
About identity matrix: identity matrix is equivalent to 1 in multiplication. Before each fast power, the ans matrix needs to be initialized. The code is as follows:
memset(ans.g,0,sizeof(ans.g)); for(int i=1;i<=N;i++)ans.g[i][i]=1;
Because the fast power is used, the time complexity of the fast power of the matrix is basically \ (O(N^3\times{log(n)}) \), where N is the size of the matrix
Here are some examples (the explanation may not be as clear as the problem solution, so emm, I hope to explain more clearly? Click to view the problem solution):
Matrix fast power template problem.
because
\(f_{i}=1\times{f_{i-1}}+1\times{f_{i-2}}\)
\(f_{i-1}=1\times{f_{i-1}}+0\times{f_{i-2}}\)
Therefore, the transfer matrix is obtained:
Transfer to:
Use column vector
Multiply the transfer matrix \ (n \) times to get \ (g {1,1} = f {n + 1}, G {2,1} = f {n} \). Output \ (g {2,1} \).
The code is as follows:
#include<bits/stdc++.h>
#define int long long
#define inf 0x3f3f3f3f
using namespace std;
int read(){
int w=0,h=1;char ch=getchar();
while(ch<'0'||ch>'9'){if(ch=='-')h=-h;ch=getchar();}
while(ch>='0'&&ch<='9'){w=w*10+ch-'0';ch=getchar();}
return w*h;
}
const int mod=1e9+7;
int n,q;
struct matrix{
int g[105][105];
matrix operator*(const matrix&b)const{
matrix c;
for(int i=1;i<=2;i++)
for(int j=1;j<=2;j++){
c.g[i][j]=0;
for(int k=1;k<=2;k++){
c.g[i][j]=(c.g[i][j]+g[i][k]*b.g[k][j]%mod)%mod;
}
}
return c;
}
matrix operator^(int k){
matrix ans,x=*this;memset(ans.g,0,sizeof(ans.g));
for(int i=1;i<=2;i++)ans.g[i][i]=1;
for(;k;k>>=1,x=x*x)if(k&1)ans=ans*x;
return ans;
}
}f;
signed main(){
q=read();
f.g[1][1]=1;f.g[1][2]=1;
f.g[2][1]=1;f.g[2][2]=0;
f=f^q;
printf("%lld\n",f.g[2][1]%mod);
return 0;
}
Application of fast power advanced order of matrix. The state transfer equation is given in the title, and the transfer matrix is listed according to the transfer equation:
Initialize \ (RES \) matrix: \ (res.g {1,1} = x 0, res.g {2,2} = 1 \); Since we want to push \ (X_n \) and the first item is \ (X_0 \), it should not be initialized to the identity matrix
It is found that the modulus m is a little large. Use turtle speed to prevent explosion \ (what, you don't know Turtle speed ride？)
#include<bits/stdc++.h>
#pragma GCC optimize(2)
#define int long long
using namespace std
int read(){
int w=0,h=1;char ch=getchar();
while(ch<'0'||ch>'9'){if(ch=='-')h=-h;ch=getchar();}
while(ch>='0'&&ch<='9'){w=w*10+ch-'0';ch=getchar();}
return w*h;
}
struct node{
int g[5][5];
}f,res,ans;
int m,a,c,x0,n,mod;
void init(){
res.g[1][1]=x0%m;res.g[1][2]=1;
f.g[1][1]=a%m;f.g[1][2]=0;
f.g[2][1]=c%m;f.g[2][2]=1；
}
int cnt(int x,int y){
int s=0;
while(y){
if(y&1)s=(s+x)%m;
x=x*2%m;
y>>=1;
}
return s%m;
}
node multiple(node x,node y){
memset(ans.g,0,sizeof(ans.g))；
for(int i=1;i<=2;i++)
for(int j=1;j<=2;j++)
for(int k=1;k<=2;k++)
ans.g[i][j]=(ans.g[i][j]+cnt(x.g[i][k]，y.g[k][j]))%m;
return ans;
}
void Fast(int n){
while(n){
if(n&1)res=multiple(res,f);
f=multiple(f,f);
n>>=1;
}
}
signed main(){
m=read()；a=read();c=read();x0=read();n=read();mod=read();
init();
Fast(n);
printf("%lld",res.g[1][1]%mod);
return 0;
}
Rogu's cold problem, the problem surface of erger
This is the baked test question given to us by the coach when I first entered the Changjun junior squadron. Of course I won't, so I always put it in the collection. Then one day I clicked on this question and found that it was actually the label of matrix multiplication, so I did it happily. Of course I won't, so I dropped it.
When you see this question, the first reaction is graph theory, of course, but when you see the requirements and data range of k, you immediately know that it is matrix multiplication, but what is the relationship between graph theory and matrix multiplication? \ (QwQ \), so you think about it.
What is the relationship between graph theory and matrix multiplication? Please see the following decomposition:
We know: \ (C {I, J} = \ sum {k = 1} ^ {n} a {I, K} × B_{k,j}\)
In the path counting of graph theory, how do we calculate the number of paths from i to j? Similar to the multi-source shortest path algorithm Floyd, we enumerate the transit points. According to the multiplication principle and addition principle, we get the following code:
for(int i=1;i<=N;i++) for(int j=1;j<=N;j++) for(int k=1;k<=N;k++) c[i][j]=c[i][j]+c[i][k]*c[k][j];
Similar to matrix multiplication whale man!!!
So what does the power of a square matrix mean?
When the exponent of graph G is one, it is obviously the number of schemes from i to j in one step. Let's consider \ (G^2 \).
\(G^2 \), we enumerate a transit point, so the original path is broken into two parts, and this is the number of solutions that take two steps from i to j! At this time, we solve the problem that we just reach the end at time k.
But there is another difficulty in this problem: Ke AI's little monster. But as long as we grasp the key of this problem: data range, we can easily solve it.
It is found that Ke AI's little monster's walking law is no more than 4, while the least common multiple of 1,2,3,4 is only 12, so we consider discussing it by case. We mark the point that the monster will reach in each case diagram, judge whether the monster is at that point and whether the blood volume is enough to pass through the current edge, quickly multiply the 12 cases together before the power to achieve the optimization effect, and pay attention to the moment There is no commutative law in matrix multiplication, so the last multiplication in case 0. Since the question asks for the minimum value, put in two points.
See for a clearer explanation of monsters Swamp fish
Freshly baked Code:
#include<bits/stdc++.h>
#define int long long
using namespace std;
int read(){
int w=0,h=1;char ch=getchar();
while(ch<'0'||ch>'9'){if(ch=='-')h=-1;ch=getchar();}
while(ch>='0'&&ch<='9'){w=w*10+ch-'0';ch=getchar();}
return w*h;
}
const int mod=10000;
int n,m,St,En,K,NaCly_Fish;
struct node{
int g[55][55];
node operator * (node b){
node c;memset(c.g,0,sizeof(c.g));
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
for(int k=1;k<=n;k++){
c.g[i][j]|=g[i][k]&b.g[k][j];//Bisection only needs to judge whether it can be reached, so fast bit operation is adopted
}
}
}
return c;
}
node operator ^ (int x){
node ans,t=*this;memset(ans.g,0,sizeof(ans.g));
for(int i=1;i<=n;i++)ans.g[i][i]=1;
for(int i=x;i;i>>=1,t=t*t)if(i&1)ans=ans*t;
return ans;
}
}a[12],ans,res,check[13],mp;
signed main(){
n=read();m=read();St=read();En=read();K=read();
for(int i=1;i<=m;i++){
int u=read(),v=read(),val=read();
mp.g[u][v]=mp.g[v][u]=val;
}
NaCly_Fish=read();//Do you know why it's nacl_fish? Because I made swamp fish. The number of piranhas in swamp fish is NFish
for(int i=0;i<12;i++)check[i]=mp;
for(int i=1;i<=NaCly_Fish;i++){
int T=read();
for(int j=0;j<T;j++){
int cur=read();
for(int k=j;k<12;k+=T)
for(int l=1;l<=n;l++)
check[k].g[l][cur]=0;//Mark the point that the monster will reach in a cycle (reincarnation)
}
}
int l=0,r=3e9;
while(l+1<r){
int mid=(l+r)>>1;
for(int i=0;i<12;i++)a[i]=check[i]；
for(int i=0;i<12;i++){
for(int j=1;j<=n;j++){
for(int k=1;k<=n;k++){
if(a[i].g[j][k])a[i].g[j][k]=a[i].g[j][k]<=mid;//Judge whether the point reached at the ith unit time in a cycle is strange and whether the blood volume is sufficient
}
}
}
memset(res.g,0,sizeof(res.g));
for(int i=1;i<=n;i++)res.g[i][i]=1;
for(int i=1;i<12;i++)res=res*a[i];res=res*a[0];
int t=K/12,step=K%12;//t is the number of times the 12 cases are multiplied together, and step is the number of remaining cases to be multiplied
memset(ans.g,0,sizeof(ans.g));
for(int i=1;i<=n;i++)ans.g[i][i]=1;
ans=ans*(res^t)；
for(int i=1;i<=step;i++)ans=ans*a[i];
if(ans.g[St][En])r=mid;
else l=mid;
}
if(l>=(1<<31-3))cout<<"'IMP0SSBLE!!";
else cout<<r;
return 0;
}