DP 1
DP 是 OI 界十分特殊的一样东西。
当我们获得一个题面时,如果我们能够定义出一个状态,并且能够直接或间接导向答案,那么这个 DP 就是可行的。
如果复杂度过高,可用优化状态,优化转移的方式减小复杂度。
矩阵优化
当转移较为固定,但数量较多时使用矩阵优化。
注意矩阵不满足交换律(你就说状态转移有没有交换律吧),使用的时候要注意。
P1357 花园
对于 \(m\le 5\),我们可以直接状压。\(f_{i,j}\) 表示第 \(i\) 个盆,后 \(m-1\) 个花盆的状态为 \(j\)。将其右移一位便可以作为现在的状态,改变最前面的 \(0/1\) 表示第 \(i\) 个位置是否放花。如果 pop_count 合法就可以转移。
\(i\) 这一维可以省掉,\(j\) 与 \(j\) 的转移使用矩阵,毕竟矩阵是固定的。由于这是一个环,所以我们要枚举 \(1\) 到 \(m-1\) 位的状态 \(s\),然后计算 \(f_{n+m-1,s}\)。当然不用每次都做一边快速幂,将转移矩阵预处理出来即可。
#include <bits/stdc++.h>
#define int long long
#define pp __builtin_popcount
using namespace std;
const int N=2.5e5+5,P=1e9+7;
int n,m,k,ans;
struct mat{
int a[16][16];
void init(){
memset(a,0,sizeof(a));
}
}f,g;
mat operator *(mat x,mat y){
mat z;z.init();
for(int i=0;i<(1<<m);++i)
for(int j=0;j<(1<<m);++j)
for(int k=0;k<(1<<m);++k)
(z.a[i][j]+=x.a[i][k]*y.a[k][j])%=P;
return z;
}mat qp(mat x,int y){
mat z=x;--y;
for(;y;y>>=1,x=x*x)
if(y&1)z=z*x;
return z;
}signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>m>>k;--m;
for(int i=0;i<(1<<m);++i)
if(pp(i)<=k)f.a[i][i>>1]=1,f.a[i][(i>>1)|(1<<m-1)]=(pp(i)<k);
f=qp(f,n);
for(int i=0;i<(1<<m);++i)
if(pp(i)<=k){
g.init();g.a[0][i]=1;g=g*f;
(ans+=g.a[0][i])%=P;
}
cout<<ans<<'\n';
return 0;
}
P1707 刷题比赛
一个超大转移矩阵,转移一个 \(a_{i+1},b_{i+1},c_{i+1},a_i,b_i,c_i,k^2,k,1,w^k,z^k\) 就行。
P3193 [HNOI2008] GTY考试
先考虑朴素做法,\(f_{i,j}\) 表示第 \(i\) 位匹配到 \(A_j\) 的方案数,转移很显然。
看到数据范围上肯定上矩阵,对于这个矩阵的构造,我们要知道这个位置填上某个数之后后缀与 \(A\) 的前缀最长的一个,于是显然的 KMP。
#include<bits/stdc++.h>
using namespace std;
const int N=2e6+5,INF=0x3f3f3f3f;
int n,m,P;string s;
struct mat{
int a[23][23];
void init(int x){
memset(a,0,sizeof(a));
for(int i=0;i<m;++i)
a[i][i]=x;
}
}g;
inline mat operator *(mat x,mat y){
mat z;z.init(0);
for(int i=0;i<m;++i)
for(int j=0;j<m;++j)
for(int k=0;k<m;++k)
(z.a[i][j]+=x.a[i][k]*y.a[k][j])%=P;
return z;
}
inline mat qp(mat x,int y){
mat ans;ans.init(1);
for(;y;y>>=1,x=x*x)
if(y&1)ans=ans*x;
return ans;
}int kmp[23];
void Kmp(){
for(int i=2,j=0;i<=m;++i){
while(j&&s[i]!=s[j+1])j=kmp[j];
if(s[i]==s[j+1])++j;
kmp[i]=j;
}
for(int i=0;i<m;++i){
for(int c='0';c<='9';++c){
int j=i;
while(s[j+1]!=c&&j)j=kmp[j];
if(s[j+1]==c)++j;
++g.a[i][j];
}
}
}
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>m>>P>>s,s=" "+s,Kmp();
mat res;res.init(0);res.a[0][0]=1;
res=res*qp(g,n);int ans=0;
for(int i=0;i<m;++i)(ans+=res.a[0][i])%=P;
cout<<ans<<'\n';
return 0;
}
P3746 [六省联考 2017] 组合数问题
\(\infty\) 是假的。
考虑 \(r\) 十分的小,不妨维护一个
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=55;
int n,p,k,r;
int a[N][N],b[N][N],ans[N][N];
signed main(){
cin>>n>>p>>k>>r;
for(int i=1;i<=k;++i)
++a[i][i],++a[i][(i+k)%k+1],ans[i][i]=1;
int base=n*k;
for(;base;base>>=1){
if(base&1){
for(int i=1;i<=k;++i)
for(int j=1;j<=k;++j)
b[i][j]=ans[i][j],ans[i][j]=0;
for(int i=1;i<=k;++i)
for(int j=1;j<=k;++j)
for(int l=1;l<=k;++l)
(ans[i][j]+=a[i][l]*b[l][j])%=p;
}
for(int i=1;i<=k;++i)
for(int j=1;j<=k;++j)
b[i][j]=a[i][j],a[i][j]=0;
for(int i=1;i<=k;++i)
for(int j=1;j<=k;++j)
for(int l=1;l<=k;++l)
(a[i][j]+=b[i][l]*b[l][j])%=p;
}cout<<ans[1][r+1]<<'\n';
return 0;
}
P5772 [JSOI2016] 位运算
设选出 \(x_1,x_2,...,x_n\),不妨令 \(x_0=R>x_1>x_2>...>x_n\),以解决互补相同与不重复。
由于 \(n\) 很小,直接状压。设 \(f_{i,s}\) 位从高到低第 \(i\) 位,\(s\) 为一个长度为 \(n\) 的二进制数,如果 \(x_j=x_{j-1}\),则 \(s\) 的第 \(j\) 位为 \(1\)。转移的时候枚举每个数的状态加上即可。
考虑这个过程是重复的,所以使用一个 \(2^n\times2^n\) 的矩阵加速即可。
#include <bits/stdc++.h>
#define int long long
#define pp __builtin_popcount
using namespace std;
const int N=128,M=55,P=1e9+7;
int n,m,k,dp[M][N];
int b[8];
string str;struct mat{
int a[N][N];
void init(){
memset(a,0,sizeof(a));
}
}f,g;
mat operator *(mat x,mat y){
mat z;z.init();
for(int i=0;i<(1<<n);++i)
for(int j=0;j<(1<<n);++j)
for(int k=0;k<(1<<n);++k)
(z.a[i][j]+=x.a[i][k]*y.a[k][j])%=P;
return z;
}mat qp(mat x,int y){
mat z=x;--y;
for(;y;y>>=1,x=x*x)
if(y&1)z=z*x;
return z;
}
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>k>>str,m=str.size(),str=" "+str;
for(int s=0;s<(1<<n);++s){
memset(dp,0,sizeof(dp));
dp[0][s]=1;
for(int i=1;i<=m;++i){
for(int j=0;j<(1<<n);++j){
if(dp[i-1][j]){
for(int k=0;k<(1<<n);++k){
if(pp(k)&1)continue;
b[0]=str[i]-'0';
for(int l=1;l<=n;++l)
b[l]=((k>>l-1)&1);
int zyq=1,to=0;
for(int l=1;l<=n;++l){
if((j>>l-1)&1){
if(b[l]>b[l-1])zyq=0;
if(b[l]==b[l-1])to|=(1<<l-1);
}
}if(zyq)(dp[i][to]+=dp[i-1][j])%=P;
}
}
}
}for(int t=0;t<(1<<n);++t)
f.a[s][t]=dp[m][t];
}g.a[0][(1<<n)-1]=1;
g=g*qp(f,k);
cout<<g.a[0][0]<<'\n';
return 0;
}
P5371 [SNOI2019] 纸牌
\(f_{i,j,k}\) 表示有 \(j\) 种 \((i-1,i,i+1)\),\(k\) 种 \((i,i+1,i+2)\) 的组合,显然 \(0\le j,k\le 2\)。
我们要将其转移到 \((i+1,i+2,i+3)\) 去,所以要枚举 \(i+3\) 的个数,显然每个顺子都有 \(i+1\),所以保证 \(i+1\) 的个数 \(\le c\) 即可,剩余分配一下系数就行。
通过这种方式我们能是实现没有 \(\ge a_i\) 的限制,并且能够处理出转移矩阵。对于有限制的中途特殊计算,没有限制的用矩阵就行。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1005,P=998244353;
struct mat{
int a[9][9];
void init(){
memset(a,0,sizeof(a));
}
void init(int x){
memset(a,0,sizeof(a));
for(int i=0;i<9;++i)
a[i][i]=x;
}
}base,ans,sum;
mat operator *(mat x,mat y){
mat z;z.init();
for(int i=0;i<9;++i)
for(int j=0;j<9;++j)
for(int k=0;k<9;++k)
(z.a[i][j]+=x.a[i][k]*y.a[k][j])%=P;
return z;
}
mat qp(mat x,int y){
mat res;res.init(1);
for(;y;y>>=1,x=x*x)
if(y&1)res=res*x;
return res;
}
int n,c,m;
int a[N],b[N];
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>c>>m;
for(int i=0;i<3;++i)
for(int j=0;j<3;++j)
for(int k=0;k<3;++k)
if(i+j+k<=c)
base.a[i*3+j][j*3+k]=(c-i-j-k)/3+1;
ans.init(1);
for(int t=1;t<=m;++t){
cin>>a[t]>>b[t];sum.init();
ans=ans*qp(base,a[t]-a[t-1]-1);
for(int i=0;i<3;++i)
for(int j=0;j<3;++j)
for(int k=0;k<3;++k)
if(i+j+k<=c){
int s=i+j+k;
if(s<b[t])s=b[t]+((s-b[t])%3+3)%3;
if(s<=c)sum.a[i*3+j][j*3+k]=(c-s)/3+1;
}
ans=ans*sum;
}ans=ans*qp(base,n-a[m]);
cout<<ans.a[0][0]<<'\n';
return 0;
}
P3176 [HAOI2015] 数字串拆分
\(f\) 数组的矩阵转移是显然的。可以发现 \(f(i)=G^i\),\(G\) 是转移矩阵,其满足结合律。我们可以求出每个区间 \([l,r]\) 的 \(f(\overline{a_la_{l+1}…a_r})\),人后设一个 \(g_i\) 为 \(\overline{a_1a_2…,a_i}\) 的答案,然后通过结合律计算出 \(g_n\)。对于每个区间的计算,考虑从 \([l+1,r]\) 推过来,也就是要乘上 \(G^{s_l10^{r-l}}\),这也是可以预处理的。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=505,P=998244353;
int n,m;string s;
struct mat{
int a[5][5];
void init(int x){
memset(a,0,sizeof(a));
for(int i=0;i<m;++i)
a[i][i]=x;
}
}z,res,f[10][N],g[N][N],ans[N];
inline mat operator *(mat x,mat y){
z.init(0);
for(int i=0;i<m;++i)
for(int j=0;j<m;++j)
for(int k=0;k<m;++k)
(z.a[i][j]+=x.a[i][k]*y.a[k][j])%=P;
return z;
}
inline mat operator +(mat x,mat y){
for(int i=0;i<m;++i)
for(int j=0;j<m;++j)
(x.a[i][j]+=y.a[i][j])%=P;
return x;
}
inline mat qp(mat x,int y){
res.init(1);
for(;y;y>>=1,x=x*x)
if(y&1)res=res*x;
return res;
}
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>s>>m,n=s.size(),s=" "+s;
for(int i=1;i<m;++i)
f[1][0].a[i-1][i]=1;
for(int i=0;i<m;++i)
f[1][0].a[i][0]=1;
f[0][0].init(1);
for(int i=2;i<=9;++i)
f[i][0]=f[1][0]*f[i-1][0];
for(int i=0;i<=9;++i)
for(int j=1;j<=n;++j)
f[i][j]=qp(f[i][j-1],10);
for(int r=n;r>=1;--r){
g[r][r]=f[s[r]-'0'][0];
for(int l=r-1;l>=1;--l)
g[l][r]=g[l+1][r]*f[s[l]-'0'][r-l];
}ans[0].a[0][0]=1;
for(int i=1;i<=n;++i){
for(int j=0;j<i;++j)
ans[i]=ans[i]+ans[j]*g[j+1][i];
}cout<<ans[n].a[0][0]<<'\n';
return 0;
}
P6772 [NOI2020] 美食家
考虑一种广义的矩阵乘法。
可以证明其满足结合律。
\(dp_{i,j}\) 为第 \(i\) 天在 \(j\) 城市的最大收益。
如果没有美食节,我们可以推出转移矩阵,因为 \(n\le 50,w\le 5\),所以可以造出一个 \(250\times 250\) 的矩阵,对于美食节我们可以断开特判,但是复杂度 \(O(kw^3n^3\log T)\),但是我们矩阵转移是一个 \(1\times n\) 的矩阵和 \(n\times n\) 的矩阵相乘,所以乘法我们可以做到 \(n^2\) 水平,预处理出转移矩阵 \(G,G^2,G^4,G^8\),所以复杂度 \(O(kw^2n^2\log T)\)。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=53,INF=1e17;
int n,m,t,k,c[N];
int id[N][6],cnt;
struct day{int t,x,y;}a[205];
bool cmp(day x,day y){return x.t<y.t;}
struct mat{
int a[253][253];
void init(int x){
for(int i=1;i<=n*5;++i)
for(int j=1;j<=n*5;++j)
a[i][j]=-INF;
for(int i=1;i<=n*5;++i)
a[i][i]=x;
}
}z,res,ans,base[31];
mat operator *(mat x,mat y){
z.init(-INF);
for(int i=1;i<=n*5;++i)
for(int j=1;j<=n*5;++j)
for(int k=1;k<=n*5;++k)
z.a[i][j]=max(z.a[i][j],x.a[i][k]+y.a[k][j]);
return z;
}
void mul(int delta){
for(int h=0;h<=30;++h)
if((delta>>h)&1){
res.init(-INF);
for(int j=1;j<=n*5;++j)
for(int k=1;k<=n*5;++k)
res.a[1][j]=max(res.a[1][j],ans.a[1][k]+base[h].a[k][j]);
ans=res;
}
}
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>m>>t>>k;
base[0].init(-INF);
for(int i=1;i<=n;++i){
cin>>c[i];
for(int j=1;j<=5;++j)
id[i][j]=++cnt;
for(int j=1;j<5;++j)
base[0].a[id[i][j]][id[i][j+1]]=0;
}
for(int i=1,u,v,w;i<=m;++i)
cin>>u>>v>>w,base[0].a[id[u][w]][id[v][1]]=c[v];
for(int i=1;i<=30;++i)
base[i]=base[i-1]*base[i-1];
for(int i=1;i<=k;++i)
cin>>a[i].t>>a[i].x>>a[i].y;
sort(a+1,a+k+1,cmp);
ans.init(-INF),ans.a[1][1]=c[1];
for(int i=1;i<=k;++i){
mul(a[i].t-a[i-1].t);
if(ans.a[1][id[a[i].x][1]]!=-INF)
ans.a[1][id[a[i].x][1]]+=a[i].y;
}if(a[k].t!=t)mul(t-a[k].t);
if(ans.a[1][1]>0)
cout<<ans.a[1][1]<<'\n';
else cout<<-1<<'\n';
return 0;
}
CF917C Pollywog
和上面的类似。
这个一看 \(k\le 8\) 直接状压,\(f_{i,s}\) 表示所观察的视野在 \([i,i+k-1]\),蝌蚪所在的集合为 \(s\),然后枚举最前面的蝌蚪移动,得到一个转移矩阵,形如 \(f_{i+1,t} \gets \text{getmin} (f_{i,s}+g_{s,t})\)。对于特殊点停下来特别做一下即可。
首先状态可以优化一下,合法状态最多有 \(H=C^4_8=70\) 个。
使用更广义的矩阵来做,预处理一个 \(g^{2^k}\),即可每次矩阵快速幂做到 \(O(H^2\log n)\),总复杂度 \(O(qH^2\log n)\)
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N=30;
int x,k,n,m,q;
int c[N],mp[1<<8],pm[70];
map<int,int> co;
struct mat{
int a[70][70];
mat(){memset(a,0x3f,sizeof(a));}
}f,g[N];
mat operator *(mat x,mat y){mat z;
for(int i=0;i<m;++i)
for(int j=0;j<m;++j)
for(int k=0;k<m;++k)
z.a[i][j]=min(z.a[i][j],x.a[i][k]+y.a[k][j]);
return z;
}struct node{int p,w;}a[N];
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>x>>k>>n>>q;
for(int i=1;i<=k;++i)
cin>>c[i];
for(int s=0;s<(1<<k);++s)
if(__builtin_popcount(s)==x)
pm[mp[s]=m]=s,++m;
for(int i=1;i<=q;++i)
cin>>a[i].p>>a[i].w,co[a[i].p]=a[i].w;
sort(a+1,a+q+1,[](node x,node y){return x.p<y.p;});
while(q&&a[q].p>n-x)--q;
a[0].p=1,a[++q]=(node){n-x+1,0};
for(int i=0;i<m;++i)
if(pm[i]&1){
for(int j=1;j<=k;++j)
if(!((pm[i]>>j)&1))
g[0].a[i][mp[(1<<j|pm[i])>>1]]=c[j];
}else g[0].a[i][mp[pm[i]>>1]]=0;
for(int l=1;l<=26;++l){
for(int i=0;i<m;++i)
for(int j=0;j<m;++j)
for(int k=0;k<m;++k)
g[l].a[i][j]=min(g[l].a[i][j],g[l-1].a[i][k]+g[l-1].a[k][j]);
}f.a[0][mp[(1<<x)-1]]=0;
for(int i=1;i<=q;++i){
for(int l=0;l<=26;++l){
if((a[i].p-a[i-1].p)>>l&1){
mat h;
for(int j=0;j<m;++j)
for(int k=0;k<m;++k)
h.a[0][j]=min(h.a[0][j],f.a[0][k]+g[l].a[k][j]);
f=h;
}
}
for(int j=0;j<m;++j)
if(pm[j]&1)f.a[0][j]+=a[i].w;
}int ans=f.a[0][mp[(1<<x)-1]];
for(int i=1;i<=x;++i)
ans+=co[n-i+1];
cout<<ans<<'\n';
return 0;
}
动态 DP
既然 DP 过程可以用矩阵维护,那么当 DP 中的参数发生变化时,我们不用重新 DP,而是修改某处的转移矩阵计算总乘积即可。
P4719 【模板】动态 DP
\(f_{i,0/1}\) 表示 \(i\) 子树中包含/不包含节点 \(i\) 的最大独立集,我们有转移:
如果用转移矩阵,我们使用这样一种矩阵:
设节点 \(u\) 的重儿子为 \(w\),\(g_{i,0}=\sum_{ v\in son(u),v\not=w}\max(f_{v,0},f_{v,1}),g_{i,1}=a_i+\sum_{v\in son(u),v\not=w}f_{v,0}\)。
所以方程可以写成
这样就可以用上面的矩阵了。
带修改的话,我们会发现只会更改结点到根的 DP 值,对于一条重链,我们维护其矩阵“乘积”,这样我们一次修改只会更改 \(\log n\) 次。
#include<bits/stdc++.h>
#define ls(p) (p<<1)
#define rs(p) (p<<1|1)
#define mid (l+r>>1)
#define pb push_back
using namespace std;
const int N=1e5+5,INF=0x3f3f3f3f;
int n,m,a[N];
int fa[N],dep[N],dfn[N],to[N],cnt;
int siz[N],son[N],top[N],pot[N];
int f[N][2];
vector<int> e[N];
struct mat{
int a[2][2];
void init(){
a[0][0]=a[0][1]=a[1][0]=a[1][1]=-INF;
}
}g[N];
mat operator *(mat x,mat y){
mat z;z.init();
for(int i=0;i<2;++i)
for(int j=0;j<2;++j)
for(int k=0;k<2;++k)
z.a[i][j]=max(z.a[i][j],x.a[i][k]+y.a[k][j]);
return z;
}struct SEG{
mat s[N<<2];
void build(int p,int l,int r){
if(l==r)return s[p]=g[to[l]],void();
build(ls(p),l,mid);build(rs(p),mid+1,r);
s[p]=s[ls(p)]*s[rs(p)];
}void update(int p,int l,int r,int x){
if(l==r)return s[p]=g[to[x]],void();
if(x<=mid)update(ls(p),l,mid,x);
else update(rs(p),mid+1,r,x);
s[p]=s[ls(p)]*s[rs(p)];
}mat query(int p,int l,int r,int L,int R){
if(L<=l&&r<=R)return s[p];
if(L>mid)return query(rs(p),mid+1,r,L,R);
if(R<=mid)return query(ls(p),l,mid,L,R);
return query(ls(p),l,mid,L,R)*query(rs(p),mid+1,r,L,R);
}
}seg;
void dfs(int u){
for(int v:e[u]){
if(v==fa[u])continue;
fa[v]=u,dep[v]=dep[u]+1;
dfs(v),siz[u]+=siz[v];
if(siz[v]>siz[son[u]])son[u]=v;
}siz[u]=1;
}void dfs(int u,int topf){
to[dfn[u]=++cnt]=u;
pot[top[u]=topf]=cnt;
f[u][1]=a[u];g[u].a[0][0]=g[u].a[0][1]=0;
g[u].a[1][0]=a[u];
if(son[u]){
dfs(son[u],topf);f[u][1]+=f[son[u]][0];
f[u][0]+=max(f[son[u]][0],f[son[u]][1]);
}for(int v:e[u])
if(v!=fa[u]&&v!=son[u])
dfs(v,v),f[u][1]+=f[v][0],
f[u][0]+=max(f[v][0],f[v][1]),
g[u].a[0][0]+=max(f[v][0],f[v][1]),
g[u].a[0][1]+=max(f[v][0],f[v][1]),
g[u].a[1][0]+=f[v][0];
}
signed main(){
ios::sync_with_stdio(0);
cin.tie(0),cout.tie(0);
cin>>n>>m;
for(int i=1;i<=n;++i)
cin>>a[i];
for(int i=1,u,v;i<n;++i)
cin>>u>>v,e[u].pb(v),e[v].pb(u);
dfs(1),dfs(1,1);
seg.build(1,1,n);
while(m--){
int u,y;cin>>u>>y;
g[u].a[1][0]+=y-a[u];
a[u]=y;
while(u){
mat m1=seg.query(1,1,n,dfn[top[u]],pot[top[u]]);
seg.update(1,1,n,dfn[u]);
mat m2=seg.query(1,1,n,dfn[top[u]],pot[top[u]]);
u=fa[top[u]];g[u].a[1][0]+=m2.a[0][0]-m1.a[0][0];
g[u].a[0][0]+=max(m2.a[0][0],m2.a[1][0])-max(m1.a[0][0],m1.a[1][0]);
g[u].a[0][1]+=max(m2.a[0][0],m2.a[1][0])-max(m1.a[0][0],m1.a[1][0]);
}mat gty=seg.query(1,1,n,dfn[top[1]],pot[top[1]]);
cout<<max(gty.a[1][0],gty.a[0][0])<<'\n';
}
return 0;
}

浙公网安备 33010602011771号