高级算法指北——浅谈网格图优化 dp
0.简述
对于一类二维或更高维度的 dp,我们把转移的过程看作在二维平面上的一种行动,于是可以直观地通过某些观察,找到加速 dp 转移的方式。
通常在状态过多导致复杂度不正确的情况下,我们会考虑使用这一类 dp 方式进行优化。
下面通过一些例题来介绍此类 dp 的具体优化方式。
1.利用网格图游走的组合意义完成计数
e.g.01. 卡特兰数通项公式的推导
卡特兰数的一个经典意义是“长度为 \(2n\) 的合法括号序列数量”,这个东西显然也可以用 dp 表示:\(f_{i,j}\) 表示填了 \(i\) 个括号,设左括号为 \(1\),右括号为 \(-1\) 时,前缀和为 \(j\) 的序列数量。一个合法的括号序列需要前缀和时刻 \(\geq 0\)。这个 dp 的转移是 \(f_{i,j}\gets f_{i-1,j-1}\) 和 \(f_{i,j}\gets f_{i-1,j+1}\)。
接下来,考虑把 dp 的转移放上平面直角坐标系的第一象限进行考虑。通常,我们喜欢把转移通过换元变成向 \(x\) 或 \(y\) 的正半轴走一步的情形,因此我们旋转坐标轴 \(45^{\circ}\),令 \((i,j)\to (\frac{i-j}{2},\frac{i+j}{2})\),则加入一个左括号可以看作是向上走一步;加入一个右括号可以看作是向右走一步。把图像画出来,不难发现前缀和非负的限制被转化为了 \(i\leq j\) 时刻成立,也即不能碰到 \(y=x-1\) 这条直线。
此时,问题变为求从 \((0,0)\) 走到 \((n,n)\),全过程中不能碰到直线 \(y=x-1\) 的方案数。注意到任意一条碰到 \(y=x-1\) 的路径,都可以从首次碰到的点开始,让接下来的路径关于这条直线对称,得到一个到 \((n+1,n-1)\) 的路径。也即一条到 \((n+1,n-1)\) 的路径与到 \((n,n)\) 的不合法路径一一对应。做一个简单容斥,我们知道从 \((0,0)\) 走到 \((n,m)\),每次向右或上走一步,不设其他限制的方案数为 \(\binom {n+m}{n}\),因此有答案为 \(\binom{2n}{n}-\binom{2n}{n+1}=\frac{\binom{2n}{n}}{n+1}\)。
此即为卡特兰数 \(C_n\) 的通项公式。
e.g.02. P10868 [HBCPC2024] Points on the Number Axis B
简要题意:给定数轴上的 \(n\) 个点,每次随机选择一对相邻的点删去,加入坐标为两个点中点的一个新点。不断执行这个操作直到只剩下一个点为止,问最后点的坐标期望。\(n\leq 10^6\)。
用线性性拆期望,那么每个点的坐标对最终答案的贡献都会带一个系数 \(k_i\)。不难发现一个点的贡献系数仅与其左右两侧点的数量有关。
考虑 dp 求系数。设 \(dp_{i,j}\) 表示对于某个点,若其左边有 \(i\) 个点,右边有 \(j\) 个点时,它的坐标的贡献系数。首先有 \(dp_{0,0}=1\),转移时枚举下一步删掉的是哪两个点,有:
注意到在网格图游走的时候, \(i+j\) 不是一个很好看的形式;又由于每一种 \(i+j\) 在转移的时候都会被计算一次,因此我们考虑把它整体提出来:
不妨设 \(f_{i,j}=(i+j)!dp_{i,j}\),有 \(f_{i,j}=(i-\frac12)f_{i-1,j}+(j-\frac12)f_{i,j-1}\)。
这是一个经典的网格图优化形式,等价于一个带权路径和问题:一条路径从 \((i-1,j)\) 走到 \((i,j)\),会乘上一个权值 \(i-\frac12\),从 \((i,j-1)\) 走到 \((i,j)\),会乘上一个权值 \(j-\frac12\),所有路径初始权值均为 \(1\),\(f_{i,j}\) 就表示从 \((0,0)\) 走到 \((i,j)\) 的所有路径的权值和。
不难发现从 \((0,0)\) 到 \((i,j)\) 的所有路径权值均为 \(\prod_{x=1}^i(x-\frac12)\prod_{y=1}^j(y-\frac12)\),总合法路径数量为 \(\binom{i+j}{i}\)。于是预处理 \(\prod_{x=1}^i(x-\frac12)\) 之后,\(f\) 和 \(dp\) 都可以 \(O(1)\) 求了。复杂度可以做到 \(O(n)\)。
int n,a[N];
int jc[N],qj[N],val[N];
int quick_power(int base,int x){
int res=1;
while(x){
if(x&1)res*=base,res%=mo;
base*=base,base%=mo;
x>>=1;
}
return res;
}
int C(int x,int y){
if(x<y)return 0;
return jc[x]*qj[y]%mo*qj[x-y]%mo;
}
int inv2;
signed main(){
read(n),inv2=quick_power(2,mo-2);
rep(i,1,n)
read(a[i]);
jc[0]=qj[0]=val[0]=1;
rep(i,1,n){
val[i]=val[i-1]*(i-inv2+mo)%mo;
jc[i]=jc[i-1]*i%mo;
qj[i]=quick_power(jc[i],mo-2);
}
int ans=0;
rep(i,1,n){
int le=i-1,ri=n-i;
int k=val[le]*val[ri]%mo*C(le+ri,le)%mo*qj[le+ri]%mo;
ans+=k*a[i]%mo,ans%=mo;
}
printf("%lld\n",ans);
return 0;
}
e.g.03.P7386 「EZEC-6」0-1 Trie
网格图优化 dp。
设 \(f_{n,m}\) 表示剩余 \(n\) 个 \(1\),\(m\) 个 \(0\) 时的答案。我们知道根节点一定是 \(0\),考虑下一位若是 \(0\),则递归到 \(f_{n,m-1}\) 的子问题;若下一位是 \(1\),则再下一位一定是 \(0\),递归到 \(f_{n-1,m-1}\) 的子问题。于是可以写出 dp 式子:\(f_{n,m}=f_{n,m-1}+f_{n-1,m-1}+2\),其中加的 \(2\) 是根节点和根节点的 \(1\) 儿子这两个点。边界是 \(f_{1,m}=m+1\)。
首先可以换元消去 \(+2\) 的项,令 \(g_{n,m}=f_{n,m}+2\),代入可以发现 \(g_{n,m}=g_{n,m-1}+g_{n-1,m-1}\),边界为 \(g_{1,m}=m+3\),\(g_{i,i-1}=2\)。注意加了 \(2\) 后,原来为 \(0\) 的位置会变成 \(2\),因此我们新增了 \(g_{i,i-1}\) 这个边界。
将 dp 转移放在平面直角坐标系的第一象限,等价于从某个边界出发,每次向上或向右上走一步,计算所有路径的权值和。下面分别考虑从两种不同边界出发的路径:
-
从边界 \((i,i-1)(i\ge 2)\) 出发的路径,我们要求第一步必须向上走到 \((i,i)\) 处。这样的路径权值均为 \(2\),总贡献为 \(\sum_{j=2}^n\binom{m-j}{n-j}\times 2=2\binom{m-1}{m-n+1}\)。
-
从边界 \((1,i)(i\ge 1)\) 出发的路径,我们要求第一步必须向右上走到 \((2,i+1)\) 处。这样的路径权值为 \(i+3\)。总贡献为 \(\sum_{j=1}^m\binom{m-j-1}{n-2}\times (j+3)\)。
注意到这里的求和由于带了一个乘积,是不好算的。因此我们考虑把每个乘积表示为另一种路线的种数。你发现若我们将起点设为 \((0,-3)\),每次要求往上或往右走一步,走到 \((1,i)\) 的方案数恰好为 \(i+3\)。因此,我们考虑将这部分计算为从 \((0,-3)\) 出发,走到 \((n,m)\) 的方案数,每条路径在首次离开第 \(1\) 列时作分割,与原来的路径即可形成一一对应关系。
需要注意的是,我们在走的过程中经过了 \((2,1),(2,0),(2,-1)\) 这三个点的路径是不合法的,需要减去。为了避免算重,我们只计数从右下方走到这三个点的路径,分别减去。最后可得出贡献为 \(\binom{n+3}{n}-\binom{m+1}{n-2}-2\times \binom{m}{n-2}-3\times \binom{m-1}{n-2}\)。
故 \(g(n,m)\) 为 \(2\binom{m-1}{m-n+1}+\binom{n+3}{n}-\binom{m+1}{n-2}-2\times \binom{m}{n-2}-3\times \binom{m-1}{n-2}\),\(f\) 还要再减 \(2\)。预处理模数以下的阶乘并使用 Lucas 定理即可。
int quick_power(int base,int x){
int res=1;
while(x){
if(x&1)res*=base,res%=mo;
base*=base,base%=mo;
x>>=1;
}
return res;
}
int T,jc[N],qj[N];
int c(int x,int y){
if(x<y)return 0;
return jc[x]*qj[y]%mo*qj[x-y]%mo;
}
int lucas(int x,int y){
if(x<mo&&y<mo)return c(x,y);
return c(x%mo,y%mo)*lucas(x/mo,y/mo)%mo;
}
int solve(){
int n,m;
read(n),read(m);
if(n==1)return m+1;
if(n>m)return 0;
int res=0;
res=2*lucas(m-1,m-n+1)+lucas(m+3,n)-lucas(m+1,n-2)-2*lucas(m,n-2)-3*lucas(m-1,n-2)-2;
res=(res%mo+mo)%mo;
return res;
}
signed main(){
jc[0]=qj[0]=1;
rep(i,1,mo-1)
jc[i]=jc[i-1]*i%mo;
qj[mo-1]=quick_power(jc[mo-1],mo-2);
repp(i,mo-2,1)
qj[i]=qj[i+1]*(i+1)%mo;
read(T);
int ans=0;
while(T--)
ans^=solve();
printf("%lld\n",ans);
return 0;
}
2.转为网格图上选点的问题
e.g.01.gym102059K Interesting Drug
简要题意:数轴上有 \(n\) 个关键点,你需要从某个关键点出发,在数轴上任意游走。记录 \(P_i\) 表示第 \(i\) 个首次经过的点的编号,对于每一个 \(i\in[1,n]\),若 \(P_{c_i}=i\),则你可以获得 \(d_i\) 的收益。对于每个关键点作为起点的情况,分别求可能的收益最大值,\(n\leq 3\times 10^5\)。
显然在任意时刻,走到过的点是一个区间。考虑倒着 dp:\(f_{i,j}\) 表示 \([i,j]\) 区间内的点被走过了,此后的最大收益,转移有 \(f_{i,j}\gets f_{i-1,j}+[c_i=j-i+2]\times d_i\) 和 \(f_{i,j}\gets f_{i,j+1}+[c_j=j-i+2]\)。初始状态是 \(f_{1,n}=0\)。
注意到有值的转移只有 \(O(n)\) 个,对应在网格图上是 \(O(n)\) 条边,每条边有一个边权。于是问题转化为从 \((1,n)\) 开始走到某个 \((i,i-1)\) 结束,每次可以向右或向下走一步,求经过边的权值和的最大值。直接设 \(g_i\) 表示从起点走到 \(i\) 这条边的最大值,将边排序后使用 BIT 优化即可做到 \(O(n\log n)\)。
int n,d[N],c[N],cntl,cntn;
pii lsh[N];
struct lines{
pii st,ed;
int sid,eid;
int v;
friend bool operator<(lines x,lines y){
return x.sid<y.sid;
}
}l[N];
bool cmpp(pii x,pii y){
if(x.sec!=y.sec)return x.sec>y.sec;
return x.fir<y.fir;
}
struct BIT{
int t[N];
void add(int x,int v){
if(!x)return;
for(int i=x;i<=n;i+=lowbit(i))
t[i]=max(t[i],v);
}
int query(int x){
int res=0;
for(int i=x;i;i-=lowbit(i))
res=max(res,t[i]);
return res;
}
}T;
int f[N],ans[N];
vector<int>stpos[N],endpos[N];
signed main(){
read(n);
rep(i,1,n)
read(c[i]);
rep(i,1,n)
read(d[i]);
rep(i,1,n){
if(i+1<=n&&i+c[i]-1<=n)l[++cntl]=(lines){mp(i,i+c[i]-1),mp(i+1,i+c[i]-1),0,0,d[i]};
if(i-c[i]+1>=1&&i-1>=1)l[++cntl]=(lines){mp(i-c[i]+1,i),mp(i-c[i]+1,i-1),0,0,d[i]};
l[++cntl]=(lines){mp(i,i),mp(i,i-1),0,0,0};
}
rep(i,1,cntl)
lsh[++cntn]=l[i].st,lsh[++cntn]=l[i].ed;
sort(lsh+1,lsh+cntn+1,cmpp),cntn=unique(lsh+1,lsh+cntn+1)-lsh-1;
rep(i,1,cntl){
l[i].sid=lower_bound(lsh+1,lsh+cntn+1,l[i].st,cmpp)-lsh;
l[i].eid=lower_bound(lsh+1,lsh+cntn+1,l[i].ed,cmpp)-lsh;
stpos[l[i].sid].push_back(i),endpos[l[i].eid].push_back(i);
}
rep(i,1,cntn){
for(auto j:endpos[i])
T.add(lsh[i].fir,f[j]);
int res=T.query(lsh[i].fir);
for(auto j:stpos[i])
f[j]=res+l[j].v;
}
rep(i,1,cntl)
if(l[i].st.fir==l[i].st.sec)ans[l[i].st.fir]=f[i]+(c[l[i].st.fir]==1)*d[l[i].st.fir];
rep(i,1,n)
printf("%lld ",ans[i]);
return 0;
}
e.g.02.P10395 青蛙寻青。
简要题意:给定长度分别为 \(n,m\) 的两个序列 \(a,b\),值域在 \([1,k]\),保证 \(1\sim k\) 都在 \(a\) 中出现过。考虑构造一个二分图 \((N,M)\),两侧点的编号均按升序排列,对于任意 \((i,j)\) 满足 \(i\in[1,n],j\in[1,m]\) 且 \(a_i=b_j\),我们在 \((i,j)\) 之间连一条边。现在,你需要修改若干个 \(b\) 中的值,使得 \(1\sim k\) 都在 \(b\) 中出现过,且生成的二分图的边不在非顶点处相交。求修改次数的最小值。\(k\leq n,m\leq 2\times 10^6\)。
首先要会判断无解。当 \(a\) 中存在下标 \(i,j,k\) 使得 \(i<j<k\) 且 \(a_i=a_k\neq a_j\),则引出的连线一定会相交,无解。否则一定有解。在有解的情况下,不妨可以对 \(a\) 序列去重,使得一个数最多在 \(a\) 中出现 \(1\) 次。下文出现的 \(n\) 均表示去重后序列 \(a\) 的长度。
朴素的做法是设 \(dp_{i,j}\) 表示考虑到 \(b\) 序列的第 \(i\) 个数,已经匹配到 \(a\) 序列的第 \(j\) 个数的最小修改次数。转移是 \(dp_{i,j}\gets \min(dp_{i-1,j},dp_{i-1,j-1})+[a_j\neq b_i]\)。整个 dp 的过程显然可以被看作是从 \((0,0)\) 开始,每次往右或右上走一步,最终走到 \((m,n)\) 的一条路径,要求最小化路径上格子权值和。
不难发现,从起始点 \((0,0)\) 开始走到 \((m,n)\) 结束的路径,其长度必定为 \(m\)。同时,由于去重后的序列 \(a\) 中所有元素均不相同,不满足 \(a_y\neq b_x\) 的格子 \((x,y)\) 仅有 \(m\) 个。那么,一条路径的权值和可以表示为 \(m\) 减去途径的满足 \(a_y=b_x\) 的格子 \((x,y)\) 的个数所得到的值。因此,问题转为最大化经过的满足 \(a_y=b_x\) 的格子 \((x,y)\) 的数量。
考虑把所有合法的格子坐标按 \(x\) 从小到大排序,设 \(f_i\) 表示以第 \(i\) 个格子结尾,能走到的格子数量的最大值。有转移 \(f_i\gets f_j+1\)。此时转移的限制比较复杂,需要所在斜率为 \(1\) 的直线的截距不降。于是我们考虑把向右上方走一步转换成向上走一步,把原本坐标为 \((x,y)\) 的格子的坐标重新设为 \((x-y,y)\),这样向右走一步仍然是 \(x\) 坐标 \(+1\),向右上方走一步就转化为了 \(y\) 坐标 \(+1\) 而 \(x\) 坐标不变。此时直接利用 BIT 做最长上升子序列类似物,即可在 \(O(\log n)\) 的复杂度内完成一步转移。最终答案为 \(m-\max_i f_i\)。
时间复杂度 \(O(m\log n)\)。
int n,m,a[N],b[N],dp[N],pos[N],cntp;
bool vis[N];
struct BIT{
int t[N];
void add(int x,int v){
for(int i=x;i<=n;i+=lowbit(i))
t[i]=max(t[i],v);
}
int query(int x){
int res=0;
for(int i=x;i;i-=lowbit(i))
res=max(res,t[i]);
return res;
}
}T;
pii tp[N];
int main(){
read(n),read(m);
rep(i,1,n)
read(a[i]);
rep(i,1,m)
read(b[i]);
rep(i,1,n){
if(!vis[a[i]]||a[i]==a[i-1]){
vis[a[i]]=1;
continue;
}
puts("-1");
return 0;
}
n=unique(a+1,a+n+1)-a-1;
rep(i,1,n)
pos[a[i]]=i;
rep(i,1,m)
if(i>=pos[b[i]])tp[++cntp]=mp(i-pos[b[i]],pos[b[i]]);
if(n>m){
puts("-1");
return 0;
}
tp[++cntp]=mp(m-n,n);
sort(tp+1,tp+cntp+1);
int ans=0;
rep(i,1,cntp){
dp[i]=T.query(tp[i].sec)+1;
T.add(tp[i].sec,dp[i]);
if(tp[i].fir==m-n&&tp[i].sec==n){
ans=m-dp[i]+1;
break;
}
}
printf("%d\n",ans);
return 0;
}
3.借用网格图的其他性质进行优化
e.g.P9870 [NOIP2023] 双序列拓展
不妨只考虑 \(f_i<g_i\) 的情况,另一种情况交换 \(n,m\) 以及 \(X,Y\) 即可类似做。
朴素的 dp 是 \(f_{i,j}\) 表示当前两个拓展的末尾分别对应 \(X_i\) 和 \(Y_j\) 是否合法。则 \(f_{i,j}=1\) 的条件首先是 \(X_i<Y_j\),若满足则有 \(f_{i,j}=f_{i,j-1}\operatorname{or}f_{i-1,j}\operatorname{or}f_{i-1,j-1}\)。对应在网格图上,等价于找一条从 \((1,1)\) 走到 \((n,m)\) 的路径,每次可以往上,往右或往右上走一步,要求走到的点 \((i,j)\) 均满足 \(X_i<Y_j\),判断是否可行。
注意到若 \(x_{min}\geq y_{min}\) 或 \(x_{max}\geq y_{max}\),则此种情况下不合法。
否则我们找到 \(x_{min}\) 所在的列和 \(y_{max}\) 所在的行,这些位置上的点一定合法,并且一条合法路径一定经过这些点中的至少一个。于是可以划分成两个子问题:从 \((1,1)\) 是否存在一条路径到达这个合法边框,以及从 \((n,m)\) 是否存在一条合法路径到达这个合法边框。以 \((1,1)\) 出发为例,此时若想要穿过所有的 \(x\) 到达右边框,则至少需要 \(x_{max}<y_{max}\),此时 \(y_{max}\) 所在行一定合法且与上一层的合法边框连通,于是可以递归到一个子问题;若想穿过所有的 \(y\) 到达上边框,则至少需要 \(y_{min}>x_{min}\),此时 \(x_{min}\) 所在列一定合法且与上一层的合法边框连通,同样可以递归到一个子问题。若两者均不满足则显然不存在合法方案。
预处理 \(X\) 和 \(Y\) 的前后缀 min 和 max,每层的复杂度可以做到 \(O(1)\)。由于矩形大小都在不断缩减,故递归层数为 \(O(n)\),总复杂度也是 \(O(n)\)。
int cid,n,m,q;
int a[N],b[N],aa[N],bb[N];
vector<int>ans;
bool swp=0;
int mina[N],mnpa[N],minb[N],mnpb[N];
int maxa[N],mxpa[N],maxb[N],mxpb[N];
int check1(int x,int y){
if(x==1||y==1)return 1;
if(mina[x-1]<minb[y-1])return check1(mnpa[x-1],y);
else if(maxa[x-1]<maxb[y-1])return check1(x,mxpb[y-1]);
else return 0;
}
int check2(int x,int y){//不断递归
if(x==n||y==m)return 1;
if(mina[x+1]<minb[y+1])return check2(mnpa[x+1],y);
else if(maxa[x+1]<maxb[y+1])return check2(x,mxpb[y+1]);
else return 0;
}
int solve(){//找到初始时候的两个子问题
if(a[1]>=b[1])return 0;
rep(i,1,n){
if(i==1||a[i]<mina[i-1])mina[i]=a[i],mnpa[i]=i;
else mina[i]=mina[i-1],mnpa[i]=mnpa[i-1];
if(i==1||a[i]>maxa[i-1])maxa[i]=a[i],mxpa[i]=i;
else maxa[i]=maxa[i-1],mxpa[i]=mxpa[i-1];
}
rep(i,1,m){
if(i==1||b[i]<minb[i-1])minb[i]=b[i],mnpb[i]=i;
else minb[i]=minb[i-1],mnpb[i]=mnpb[i-1];
if(i==1||b[i]>maxb[i-1])maxb[i]=b[i],mxpb[i]=i;
else maxb[i]=maxb[i-1],mxpb[i]=mxpb[i-1];
}
if(maxa[n]>=maxb[m]||mina[n]>=minb[m])return 0;
if(!check1(mnpa[n],mxpb[m]))return 0;
repp(i,n,1){
if(i==n||a[i]<mina[i+1])mina[i]=a[i],mnpa[i]=i;
else mina[i]=mina[i+1],mnpa[i]=mnpa[i+1];
if(i==n||a[i]>maxa[i+1])maxa[i]=a[i],mxpa[i]=i;
else maxa[i]=maxa[i+1],mxpa[i]=mxpa[i+1];
}
repp(i,m,1){
if(i==m||b[i]<minb[i+1])minb[i]=b[i],mnpb[i]=i;
else minb[i]=minb[i+1],mnpb[i]=mnpb[i+1];
if(i==m||b[i]>maxb[i+1])maxb[i]=b[i],mxpb[i]=i;
else maxb[i]=maxb[i+1],mxpb[i]=mxpb[i+1];
}
return check2(mnpa[1],mxpb[1]);
}
int main(){
read(cid),read(n),read(m),read(q);
rep(i,1,n)
read(a[i]),aa[i]=a[i];
rep(i,1,m)
read(b[i]),bb[i]=b[i];
if(a[1]>b[1]){
swp=1;
rep(i,1,max(n,m))
swap(a[i],b[i]);
swap(n,m);
}
ans.push_back(solve());
if(swp)swap(n,m),swp=0;
while(q--){
rep(i,1,n)
a[i]=aa[i];
rep(i,1,m)
b[i]=bb[i];
int num1,num2;
read(num1),read(num2);
rep(i,1,num1){
int x,v;
read(x),read(v),a[x]=v;
}
rep(i,1,num2){
int x,v;
read(x),read(v),b[x]=v;
}
if(a[1]>b[1]){
swp=1;
rep(i,1,max(n,m))
swap(a[i],b[i]);
swap(n,m);
}
ans.push_back(solve());
if(swp)swap(n,m),swp=0;
}
rep(i,0,(int)ans.size()-1)
printf("%d",ans[i]);
return 0;
}