矩陣快速冪優化的動態規劃
阿新 • • 發佈:2018-12-24
因為最近寫矩陣快速冪總是搞反相乘的順序所以來寫一發部落格
不過突然寫這麼簡單的東西似乎會被鄙視?
矩陣和矩陣乘法
(沒學過矩乘的同學最好不要通過此部落格學習,它的目的不在於此)
矩陣是一個由陣列成的矩形,特殊地,行數為1的矩陣被稱作行向量,列數為1的矩陣被稱作列向量。
矩陣A(m*n)和B(n*p)可以進行乘法得到矩陣C(m*p),也就是說當前者的列數等於後者行數的時候。乘法的定義是
cij=∑mk=1Aik∗Bkj
按照定義,可以寫出矩陣和矩乘的程式碼,其中矩乘的效率是O(n*m*p)
struct matrix { int m[maxm][maxm]; int x,y; inline int *operator [](int x) { return m[x]; } matrix(int x,int y):x(x),y(y) { memset(m,0,sizeof(m)); } matrix operator *(matrix &b) { matrix ans(x,b.y); for(int i=0;i<x;i++) for(int j=0;j<b.y;j++) for(int k=0;k<y;k++) ans[i][j]+=m[i][k]*b[k][j];//一般的題目都會在此處要求取模,請注意相乘時int溢位 return ans; } };
矩陣乘法是滿足結合律的,也就是說,(A*B)*C=A*(B*C)
矩陣乘法與DP
對於滿足DP[i][j]需要從若干個DP[i-1][k]處取值乘上常數求和的,且每一輪從i-1到i轉移方式固定的動態規劃(如果上面的解釋很難懂可以看公式) DP[i][j]=∑mk=1DP[i−1][k]∗A[j][k]可以看做是每次把初始的DP[i-1]向量乘上一個行數和列數相等的轉移矩陣A,變換成一個長度相等的DP[i]向量。若使得DP陣列對應向量是行向量,那麼A[j][k]表示從DP[i-1][j]轉移到DP[i][k]時需要乘的數,需要將行向量乘以轉移矩陣(順序很重要!不滿足交換律)。若DP陣列對應向量是列向量,那麼A[j][k]表示從DP[i-1][k]轉移到DP[i][j]次需要乘以的數,需要將轉移矩陣乘以向量
矩陣快速冪及其優化
前面提到矩乘滿足結合律,那麼設DP[0]對應的向量是B,轉移矩陣為A,B*A*A*...*A=B*(A*A*...*A)=B*(A^n) 其中A^n可以在O(m^3*log n)時間內完成,其中m^3是矩乘的時間 程式碼如下但是上文中的乘法,傳了一整個matrix結構體作為引數,又返回了一個結構體,這使得賦值佔用了大量的時間,可能會被卡常數。那麼能不能通過傳指標的方法避免這部分常數呢?答案是肯定的matrix operator ^(matrix b,long t) { matrix ans;//此處ans應初始化為單位矩陣 while(t) { if(t&1) ans=ans*b; b=b*b; t>>=1; } return ans; }
void quickmul(matrix *a,matrix *b,matrix *ans)
{
memset(ans->m,0,sizeof(ans->m));
for(int i=0;i<a->x;i++)
for(int j=0;j<b->y;j++)
for(int k=0;k<a->y;k++)
ans->m[i][j]+=a->m[i][k]*b->m[k][j];
}
matrix operator ^(matrix in,int t)
{
matrix tmpbase(in.x,in.y);
matrix tmp1(in.x),tmp2(in.x);
matrix *ans=&tmp1,*swapans=&tmp2;
matrix *b=&in,*swapb=&tmpbase;
while(t)
{
if(t&1)
{
quickmul(ans,b,swapans);
swap(ans,swapans);//交換指標地址,不改變其中的值
}
quickmul(b,b,swapb);
swap(b,swapb);
t>>=1;
}
return *ans;
}