【动态规划】常见基础DP优化技巧（二）

· · 算法·理论

upd on 26.8.16：修改了一处笔误。

前言

感谢 @jucason_xu 讲解的部分芝士点，让我们一起膜拜他！

有问题可以在评论区或私信提出。

一些前置芝士：单调队列、单调栈、线段树、矩阵。

正文

一般来说，DP 转移的优化一般需要保证：状态合理的设计。

其实就是数组开的下，否则一切对于时间上的优化都是空谈。

I. 线段树优化 DP

有时候我们需要快速进行区间求值来帮助 DP 转移，而线段树就很好的支持了这种操作。

先看例题吧。

1. P1848 [USACO12OPEN] Bookshelf G

考虑朴素 DP。设计 f_i 为划分出 [1,i] 的答案。我们枚举 j 表示上一个区间的右端点，那么转移是简单的：

f_i=\min\{f_j+\max_{k=j+1}^i\{h_k\}\},\sum_{k=j+1}^iw_k\le L

我们从 \max_{k=j+1}^i\{h_k\} 下手优化转移。我们考虑维护一个递减的单调栈，这样就能很好的统计出对于每个最大值，它所覆盖的区间。同时我们用数据结构维护 \min\{val_j=f_j+\max_{k=j+1}^i\{h_k\}\}。

我们把假设某个时刻单调栈内情况如下图：

枚举到 i 时，我们把 i 加入单调栈。首先我们发现 h_{top} 没有 h_i 大（即蓝色部分是没有橙色部分高的），说明栈顶所对应区间的最大值会变成 h_i（即蓝色部分会被橙色部分覆盖），于是把栈顶覆盖区间的 val_j 内加上 h_i-h_{top} ，并且将栈顶弹出……重复进行刚刚的操作直到 h_i<h_{top}，然后 i 入栈。

而 \sum_{k=j+1}^iw_k\le L 的限制，显然 j 具有单调性。我们用二分找到最小的 j 满足这个限制，然后在 [j,i-1] 取最大值即可。

综上，我们需要一个支持区间加，区间查询的数据结构，这不就是线段树吗？时间复杂度 O(n\log n)。

:::success[AC code]

const int N=1e5+8;
int n,m;
int h[N],w[N];
int st[N],top;
struct SGT {
#define mid ((l+r)>>1)
#define ls o<<1
#define rs o<<1|1
    int mn[N<<2],tag[N<<2];
    void pd(int o) {
        tag[ls]+=tag[o];
        tag[rs]+=tag[o];
        mn[ls]+=tag[o];
        mn[rs]+=tag[o];
        tag[o]=0;
    }
    void upd(int o,int l,int r,int L,int R,int v) {
        if(L<=l&&r<=R) {
            tag[o]+=v;
            mn[o]+=v;
            return ;
        }
        pd(o);
        if(L<=mid) upd(ls,l,mid,L,R,v);
        if(mid<R) upd(rs,mid+1,r,L,R,v);
        mn[o]=min(mn[ls],mn[rs]);
    }
    int qry(int o,int l,int r,int L,int R) {
        if(L<=l&&r<=R) return mn[o];
        pd(o);
        int res=INF;
        if(L<=mid) res=min(res,qry(ls,l,mid,L,R));
        if(mid<R) res=min(res,qry(rs,mid+1,r,L,R));
        return res;
    }
}sgt;
int f[N];
signed main()
{
    read(n);
    read(m);
    rep(i,1,n) read(h[i]),read(w[i]),w[i]=w[i-1]+w[i];
    st[++top]=0;
    h[0]=INF;
    int pos=0;
    rep(i,1,n) {
        sgt.upd(1,0,n,i-1,i-1,f[i-1]+h[i]);
        while(top&&h[st[top]]<h[i]) {
            sgt.upd(1,0,n,st[top-1],st[top]-1,h[i]-h[st[top]]);
            top--;
        }
        st[++top]=i;
        while(pos<i&&w[i]-w[pos]>m) pos++;
        f[i]=sgt.qry(1,0,n,pos,i-1);
    }
    print(f[n]);
    return 0;
}

:::

这里的线段树是普通的下标维护，我们再看一道权值线段树优化 DP 的题目。

2. AT_dp_q Flowers

仍先考虑朴素 DP，设 f_{i} 表示以 i 为结尾的上升子序列的权值最大值。枚举上一个结尾 j，于是：

f_i=\max_j\{f_j+a_i\},(h_j\le h_i)

我们发现：与转移相关的东西其实不是下标，它只决定了你 DP 的顺序。真正与转移相关的是值域（即 h_j），于是考虑重新设计状态：设计 f_x 表示：以高度为 x 的花结尾的最大值，转移则变成：

f_{h_i}=\max_{j=1}^{h_i}\{f_j\}+a_i

我们发现 j 是一段区间，于是用权值线段树维护 f_x 的区间最大值。不过，这个 j 也是一段前缀，可以用常熟更优的树状数组维护。

时间复杂度 O(n\log n)。

:::success[AC code]

const int N=2e5+8;
int n;
int h[N],a[N],f[N];
struct BIT {
#define lowbit(x) x&(-x)
    int w[N];
    void upd(int pos,int v) {
        while(pos<=n) {
            w[pos]=max(w[pos],v);
            pos+=lowbit(pos);
        }
    }
    int qry(int pos) {
        int res=0;
        while(pos) {
            res=max(w[pos],res);
            pos-=lowbit(pos);
        }
        return res;
    }
}fenwick;
signed main()
{
    read(n);
    rep(i,1,n) read(h[i]);
    rep(i,1,n) read(a[i]);
    rep(i,1,n) {
        f[h[i]]=fenwick.qry(h[i]-1)+a[i];
        fenwick.upd(h[i],f[h[i]]);
    }
    print(fenwick.qry(n));
    return 0;
}

:::

总结

为啥第一题是以普通下标建立线段树，第二题却是值域呢？因为第一题中，以 h_i 为最大值的区间 [l,r]，在固定 r 时的 l（即下标）是一段区间；但第二题中，真正是区间的是值域。即使用线段树优化 DP 时，确定区间的构成是一个必要的步骤。

II. 斜率优化 DP

感谢 @jucason_xu 对于此部分芝士的讲解。

斜率优化，其实就是单调队列优化的升级版。只要你能推出式子，看出是到斜率优化的题目，敲出板子，那这道题就没了。

概括的讲，如果转移式子形如 f_i=\max_{j=1}^{i-1}\{f_j+a_j\times b_i+c_j\}+d_i 的形式，那么它就可以使用斜率优化成 O(1) 或 O(\log n) 的单次转移。

3. AT_dp_z Frog 3

依旧考虑暴力 DP：设计 f_i 表示到达 i 的答案，枚举上一个点 j，所以：

f_i=\min_{j=1}^{i-1}\{f_j+(h_i-h_j)^2+C\}

考虑优化。我们把转移式中的括号展开：

f_i=\min_{j=1}^{i-1}\{f_j+h_i^2-2\times h_i\times h_j+h_j^2+C\}

我们发现 h_i^2+C 与 j 无关，将其提出来，那么：

f_i=\min_{j=1}^{i-1}\{f_j+h_j^2-2\times h_i\times h_j\}+h_i^2+C

此时有人就能发现了：如果我们设 b_i=f_i+h_i^2,k_i=-2\times h_i，那么转移式子又变成了：

f_i=\min_{j=1}^{i-1}\{k_j\times h_i+b_j\}+h_i^2+C

这不就是在一些一次函数中，查询 x=b_i 的最大值吗？这时就有 DS 大手子说话了：直接无脑李超线段树啊！但我不会啊。我们思考一下怎么去优化。

我们可以画一个图：

紫色部分就是我们要求的答案。容易发现，这个是凸包的一部分。于是：每个 j 最多只会贡献一段区间，或是不贡献。且如果 j 贡献出了 [l_j,r_j]，那么就一定有 r_j<l_k(j<k)。

这就和单调队列优化 DP 很像。前面也提到过：

斜率优化，其实就是单调队列优化的升级版。

那么我们就可以利用相同的方法优化 DP。先考虑什么时候弹出队头：假设某个时候 x 比 y 优秀，那么显然：

k_x\times h_i+b_x\le k_y\times h_i+b_y\\ b_x-b_y\le (k_y-k_x)\times h_i\\ \frac{b_x-b_y}{k_y-k_x}\le h_i\\ \frac{f_x+h_x^2-f_y-h_y^2}{h_x-h_y}\le 2\times h_i\\

由于 \frac{b_x-b_y}{k_y-k_x} 的形式很像斜率，这也便是斜率优化名字的由来。

接下来看如何弹出队尾。首先你得知道为啥要弹出对尾：

如图，橙色为我们维护的凸包。但是绿线并没有经过这个凸包，这就是为啥要弹出队尾：因为有的直线根本就做不了贡献。都说单调队列优化 DP 是：比你小的比你强，你就可以退役了。那么斜率优化就是：比你小的发力比你早，你也可以退役了。（学 OI 太难了 TwT）。

考虑怎么弹出队尾：我们找出 i 和队尾 q_{tl} 的交点，这就是 i 发力的时候；然后是我们看 q_{tl} 与 q_{tl-1} 的交点，这就是 q_{tl} 发力的时候。如果 q_{tl} 发力的时候，i 已经发力了，那么 q_{tl} 就没地方发力了，只能退役。交点判断也是简单的，初中知识不细说。

那么这道题就没了。时空复杂度 O(n)。

:::success[AC code]

const int N=2e5+8;
int n,c;
int h[N];
int f[N];
int hd,tl,q[N];
double spole(int x,int y) {
    return 1.0*(f[x]+h[x]*h[x]-f[y]-h[y]*h[y])/(h[x]-h[y]);
}
signed main()
{
    read(n);read(c);
    rep(i,1,n) read(h[i]);
    q[hd=tl=1]=1;
    rep(i,2,n) {
        while(hd<tl&&spole(q[hd+1],q[hd])<=2*h[i]) hd++;
        f[i]=f[q[hd]]+(h[i]-h[q[hd]])*(h[i]-h[q[hd]])+c;
        while(hd<tl&&spole(q[tl-1],q[tl])>=spole(q[tl],i)) tl--;
        q[++tl]=i;
    }
    print(f[n]);
    return 0;
}

:::

总结

其实上面的斜率优化是 O(1) 转移的是因为 h_i 是单调递增的。但是如果 h_i 不是单调递增，我们就得二分查找第一个使斜率大于等于 2\times h_i 的 j，因此是 O(\log n) 转移的。

斜率优化理解了是很好写的。当然如果你是 DS 大手子你就可以忽略了。

III.决策单调性

什么是决策单调性呢？设 p(i) 表示转移 i 的最优决策点，那么有 p(i)\le p(i+1)。决策单调性一定满足一个性质：四边形不等式：

w(a,c)+w(b,d)\le w(a,d)+w(b,c)

如果说某个题目需要你划分出 k 个区间，并且其满足四边形不等式时，那么它可能是一道决策单调性优化 DP 的题目。怎么看满足不满足四边形不等式呢？@jucason_xu 说：

打个表就可以了。

其实决策单调性是个很大的东西：单调队列，斜率优化也算一种决策单调性优化。不过有的时候，式子并不是很好用以上两种方法维护，那就是决策单调性优化的时候了。

那么怎么用上这个性质呢？考虑分治。

假设我们处理 j\in[l,r]，并且知道 p(j)\in[p(l),p(r)] 中：

显然时间复杂度是 O(nk\log n)，其中 k 表示你一共要重复分治 k 次。

但假如我们区间没有个数限制，即只有一维呢？考虑 CDQ 分治。假设处理区间 [l,r]：

其实关于决策单调性的优化还有二分栈、二分单调队列的方式，不过笔者并没有掌握很熟悉，所以不会提及。

4. CF868F Yet Another Minimization Problem

考虑暴力 DP。设 f_{i,j} 表示第 i 个区间的结尾为 j 的最小值，枚举上一个区间右端点 k，于是：

f_{i,j}=\min\{f_k+w(k+1,j)\}

显然这个东西满足四边形不等式，按照上面分治即可。那么难点就在于如何快速计算 w(l,r)。

我们借用莫队的思想，每一层均摊是 O(n) 的。总时间复杂度还是 O(nk\log n)。

:::success[AC code]

const int N=1e5+8;
int n,k;
int a[N];
int f[N][28];
int res,L=1,R;
int tot[N];
int w(int ql,int qr) {
    while(ql<L) res+=tot[a[--L]]++;
    while(R<qr) res+=tot[a[++R]]++;
    while(ql>L) res-=--tot[a[L++]];
    while(R>qr) res-=--tot[a[R--]];
    return res;
}
void solve(int l,int r,int pl,int pr,int j) {
    if(l>r) return ;
    int mid=(l+r)>>1;
    int pos;
    f[mid][j]=INF;
    rep(i,pl,min(pr,mid)) 
        if(f[mid][j]>f[i-1][j-1]+w(i,mid)) 
            f[mid][j]=f[i-1][j-1]+w(i,mid),pos=i;
    solve(l,mid-1,pl,pos,j);
    solve(mid+1,r,pos,pr,j);
}
signed main()
{
    read(n);read(k);
    rep(i,1,n) read(a[i]);
    rep(i,1,n) f[i][0]=INF;
    f[0][0]=0;
    rep(j,1,k) solve(1,n,1,n,j);
    print(f[n][k]);
    return 0;
}

:::

IV. 根号分治优化 DP

有的时候 DP 复杂度为 O(n\times \frac{n}{a_i})，可以使用根号分治优化。

具体的，对于 a_i\le \sqrt{n}，我们可以使用一些性质或是数组来记录某些值达到优化的目的；对于 a_i>\sqrt{n}，我们直接暴力 DP。

当然，有的时候分治是相反的，怎么搞还得看具体的题目。

5. P10761 [BalticOI 2024] Trains

设计 f_i 表示到达 i 的答案。转移显然：

f_{i+j\times d_i}\Leftarrow f_i,(1\le j\le x_i)

本题直接暴力是 O(n^2) 的，考虑优化。

对于 d_i>\sqrt{n}，我们直接暴力转移；否则我们记录 g_{i,j} 表示通过 d_k=j\le \sqrt{n} 到达 i 的方案数。转移显然：g_{i,j}\Leftarrow g_{i-j,j}。

接下来我们考虑 f 对 g 的贡献。假设没有 x_i 的限制，显然 g_{i+d_i,d_i}\Leftarrow f_i；我们考虑一个差分的思想：g_{i+d_i\times(x_i+1),d_i}\Leftarrow -f_i。而统计 g 对 f 的贡献是简单的，这里不多赘述。

时间复杂度 O(n\sqrt{n})。

:::success[AC code]

const int N=1e5+8;
const int p=1e9+7;
int n,ans;
int a[N],b[N];
int f[N],g[N][408];
signed main()
{
    n=read();
    rep(i,1,n) a[i]=read(),b[i]=read();
    f[1]=1;
    rep(i,1,n) { 
        rep(j,1,sqrt(n)) if(i-j>0) g[i][j]=(g[i][j]+g[i-j][j])%p; 
        rep(j,1,sqrt(n)) f[i]=(f[i]+g[i][j])%p;
        ans=(ans+f[i])%p;
        if(a[i]>sqrt(n)) {
            rep(j,1,b[i]) if(i+j*a[i]<=n) {
                    f[i+j*a[i]]+=f[i];
                    f[i+j*a[i]]%=p;
                } else { break; }
        }
        else if(i+a[i]<=n&&a[i]) {
            g[i+a[i]][a[i]]+=f[i];
            g[i+a[i]][a[i]]%=p;
            if(i+a[i]*(b[i]+1)<=n) {
                g[i+a[i]*(b[i]+1)][a[i]]+=p-f[i];
                g[i+a[i]*(b[i]+1)][a[i]]%=p;
            }
        }
    }
    write(ans);
    return 0;
}

:::

6. P6189 [NOI Online #1 入门组] 跑步

整数拆分问题，考虑两种暴力做法：

  1. 设计 f_{i,j} 表示最大数为 i 得出 j 的方案数，那么有：f_{i,j}=\sum_{k=1}^i f_{k,j-i}，可以用前缀和优化。

  2. 设计 f_{i,j} 表示 i 个数得出 j 的方案数，即完全背包。对于 i,j，我们可以创造 0 为 i，那么 f_{i,j}\Leftarrow f_{i-1,j}；我们也可以选择给 i 的前缀加 1，那么 f_{i,j}\Leftarrow f_{i,j-i}。

两种方法不由得让我们想道根号分治。考虑一个性质：最多只会有 \sqrt{n} 个数大于等于 \sqrt{n}。因此，我们把组成部分以 \sqrt{n} 为分割点，划分成两个段计算贡献。

对于第一个段为最大值小于 \sqrt{n}。我们可以直接用方法 1 求出。

对于第二个段为最小值大于等于 \sqrt{n}。我们发现：方法 2 是以 0 作为初值的，我们可以把这个初值设为 \sqrt{n}，那么对于第一个转移方法就变成了 f_{i,j}\Leftarrow f_{i-1,j-\sqrt{n}}，即创造 \sqrt{n} 为第 i 个数。

最后计算答案，我们枚举第一段得出的数 j，那么第二段就得出了 n-j，将两种方案相乘即可。 :::success[AC code]

const int N=1e5+8;
int n,p,B;
int s[338][N],f[N];
int g[338][N];
signed main()
{
    n=read(),p=read();
    B=sqrt(n);
    f[0]=1;
    rep(i,0,B-1) s[i][0]=1;
    rep(j,1,n) {
        rep(i,1,min(j,B-1)) f[i]=s[i][j-i];
        rep(i,1,B-1) s[i][j]=((ll)(s[i-1][j]+f[i]))%p;
    }
    g[0][0]=1;
    rep(i,1,n/B) rep(j,i,n) {
        if(j>=B) g[i][j]=g[i-1][j-B];
        g[i][j]=((ll)(g[i][j]+g[i][j-i]))%p; 
    } 
    int ans=0;
    rep(j,0,n) {
        int cur=0;
        rep(i,0,n/B) cur=(cur+g[i][n-j])%p;
        ans=(ans+((ll)s[B-1][j]*cur)%p)%p;
    }
    write(ans);
    return 0;
}

:::

总结

其实大多数根号分治优化 DP 是极典的，只要你看出来时间复杂度是和 a_i 值呈反比例的基本都是根号分治优化的。

V. 矩阵优化 DP

差点忘了这个。最好想最难调的一集。

我们知道：无论是普通的矩阵还是 min-plus 矩阵都是满足结合律的。有时 DP 的某一维状态极少，并且与前面相关的状态也很少，那么就可以用矩阵加速从 O(nm) 优化到 O(m^3\log n)。

比如说斐波那契数列：f_{i}=f_{i-1}+f_{i-2}。我们考虑设计矩阵转移：

V\times\begin{bmatrix}f_i \\ f_{i-1}\end{bmatrix}=\begin{bmatrix}f_i+f_{i-1} \\ f_i\end{bmatrix}=\begin{bmatrix}f_{i+1} \\ f_i\end{bmatrix}

我们发现当 V=\begin{bmatrix}1&1 \\ 1&0\end{bmatrix} 时满足上方等式。那么就有：

\begin{bmatrix}1&1 \\ 1&0\end{bmatrix}\times\cdots\times\begin{bmatrix}1&1 \\ 1&0\end{bmatrix}\times\begin{bmatrix}1\\1\end{bmatrix}=\begin{bmatrix}f_{n} \\ f_{n-1}\end{bmatrix}

（一共 n-2 个 \begin{bmatrix}1&1 \\ 1&0\end{bmatrix}）

由于矩阵满足结合律，于是：

\begin{bmatrix}1&1 \\ 1&0\end{bmatrix}^{n-2}\times\begin{bmatrix}1\\1\end{bmatrix}=\begin{bmatrix}f_{n} \\ f_{n-1}\end{bmatrix}

通过矩阵快速幂即可优化到 O(\log n) 求出。不过这个式子太好推了，我们可以看一道不是一眼秒的题目。

7. P15234 「CROI · R3」浣熊的长木桥

首先考虑开始没有砖块的方案数。

由于砖块的样子，不难发现其一定是由若干个 2\times k 的长方形组成的。所以其一定是一列一列填的。

考虑状压 DP。设 f_{i,0/1/2/3} 表示第 i 列的状态为 00/01/10/11 时的方案数。则有转移：

\begin{cases}f_{i,0}=f_{i-1,3} \\ f_{i,1}=f_{i-1,0}+f_{i-1,3} \\ f_{i,2}=f_{i-1,0}+f_{i-1,3} \\ f_{i,3}=2\times f_{i-1,0}+f_{i-1,1}+f_{i-1,2}+f_{i-1,3}\end{cases}

不难发现：这个 f_{i,state} 只和前面的 f_{i-1,state'} 有关，可以用矩阵乘法表示（注意这里的矩阵下标是从 0 开始的）。设现在的状态为：

\begin{bmatrix}f_{i,0} \\ f_{i,1} \\ f_{i,2} \\ f_{i,3}\end{bmatrix}

则有转移为：

\begin{bmatrix}0&0&0&1 \\ 1&0&0&1 \\ 1&0&0&1 \\ 2&1&1&1\end{bmatrix}\times\begin{bmatrix}f_{i-1,0} \\ f_{i-1,1} \\ f_{i-1,2} \\ f_{i-1,3}\end{bmatrix}=\begin{bmatrix}f_{i,0} \\ f_{i,1} \\ f_{i,2} \\ f_{i,3}\end{bmatrix}

这个转移可以用矩阵快速幂优化。

接着考虑开始有砖块的方案数。设当前矩阵为 a。对于第 i 列每一列被提前填好的格子的状态 p，如果 p 和 state 有重复部分，则说明 state 是非法情况；否则，就把新矩阵的 a'_{p+state,0}\Leftarrow a_{state,0}。

最终答案就为 a_{3,0}。时间复杂度是 O(m+\log n) 的，不过有个矩阵自带的 64 的常数。

:::success[AC code]

const int N=1e5+8;
const int p=1e9+7;
struct matrix {
    int r,c;
    int a[4][4];
    matrix() { memset(a,0,sizeof(a)); }
    matrix operator *(const matrix &qwq) const {
        matrix res;
        res.r=r;res.c=qwq.c;
        rep(i,0,res.r-1) rep(j,0,res.c-1)
            rep(k,0,qwq.r-1) res.a[i][j]=(res.a[i][j]+a[i][k]*qwq.a[k][j]%p)%p;
        return res;
    }
}base,f;
matrix qpow(matrix base,int k) {
    matrix res;
    res.r=res.c=base.r;
    for(int i=0;i<res.r;i++)
        res.a[i][i]=1;
    while(k) {
        if(k&1) res=res*base;
        base=base*base;
        k>>=1;
    }
    return res;
} 
int n,m;
struct node {
    int x,y;
    bool operator <(const node &qwq) const {
        if(y==qwq.y) return x<qwq.x;
        return y<qwq.y;
    }
}a[N];
signed main()
{
    read(n);read(m);
    rep(i,1,m) read(a[i].x),read(a[i].y),a[i].x--;
    sort(a+1,a+m+1);
    base.r=base.c=4;
    base.a[0][3]=base.a[1][3]=
    base.a[2][3]=base.a[3][3]=
    base.a[1][0]=base.a[2][0]=
    base.a[3][1]=base.a[3][2]=1;
    base.a[3][0]=2;
    f.r=4;f.c=1;
    f.a[3][0]=1;
    int lst=0;
    rep(i,1,m) {
        int p=1<<a[i].x;
        if(a[i].y==a[i+1].y) p|=1<<a[i+1].x,i++;
        f=qpow(base,a[i].y-lst)*f;
        matrix tmp;
        tmp.r=4;tmp.c=1;
        rep(st,0,3) if((st&p)==0) 
            tmp.a[st|p][0]+=f.a[st][0];
        f=tmp;
        lst=a[i].y;
    }
    f=qpow(base,n-lst)*f;
    print(f.a[3][0]);
    return 0;
}

:::

本题是一般矩阵的优化。接下来我们看一道 min-plus 矩阵优化的题目：

8. P2886 [USACO07NOV] Cow Relays G

定义 f_{i,j} 表示经过 i 条边到达 j 的方案数。显然有转移：

f_{i,j}=\min_{k=1}^n\{f_{i-1,k}+dis_{k,j}\}

其中 dis_{x,y} 表示 x,y 的边的长度。如果没有边则为无穷大。发现边数很少，考虑对出现过的点离散化，设计 min-plus 矩阵：

\begin{bmatrix}f_{i,1} \\ f_{i,2} \\ \vdots \\ f_{i,n}\end{bmatrix}

那么转移就变成了：

V\times\begin{bmatrix}f_{i-1,1} \\ f_{i-1,2} \\ \vdots \\ f_{i-1,n}\end{bmatrix}=\begin{bmatrix}f_{i,1} \\ f_{i,2} \\ \vdots \\ f_{i,n}\end{bmatrix}

如果你的直觉较好，那么就可以发现：其实 V 就是邻接矩阵。因为 f 的转移形式不就是 min-plus 的转移。

那么这道题就没了。时间复杂度 O(m\log m+m^3\log n)。

:::success[AC code]

const int N=108;
int n,m,s,t,l;
int dot[N*10];
int x[N*10],y[N*10],z[N*10];
struct matrix {
    int r,c;
    int a[N][N];
    matrix() {memset(a,0x3f3f3f,sizeof(a));}
    matrix operator *(const matrix &qwq) const {
        matrix res;res.r=r;res.c=qwq.c;
        rep(i,1,res.r) rep(j,1,res.c) rep(k,1,qwq.r)  
            res.a[i][j]=min(res.a[i][j],a[i][k]+qwq.a[k][j]);
        return res;
    }
}base;
matrix qpow(matrix x,int y) {
    matrix res=x;
    y--;
    while(y) {
        if(y&1) res=res*x;
        y>>=1;
        x=x*x;
    }
    return res;
}
signed main()
{
    l=read();m=read();
    s=read();t=read();
    rep(i,1,m) x[i]=read(),dot[++n]=y[i]=read(),dot[++n]=z[i]=read(); 
    sort(dot+1,dot+n+1);
    n=unique(dot+1,dot+n+1)-dot-1;
    s=lower_bound(dot+1,dot+n+1,s)-dot;
    t=lower_bound(dot+1,dot+n+1,t)-dot;
    base.r=base.c=n;
    rep(i,1,m) {
        y[i]=lower_bound(dot+1,dot+n+1,y[i])-dot,
        z[i]=lower_bound(dot+1,dot+n+1,z[i])-dot;
        base.a[y[i]][z[i]]=base.a[z[i]][y[i]]=x[i];
    }
    base=qpow(base,l);
    write(base.a[s][t]);
    return 0;
}

:::

总结

矩阵加速是较板的，毕竟就推个式子、写个快速幂。当然矩阵还可以与线段树、树链剖分组成动态 DP，不过这里并不会过多说明。（其实就是笔者太菜了）

VI. 非常见 DP 优化

其实就是剪枝卡常（？

不过在大多数时候这个是必要的，良好的卡常技巧其实可以起到很大的作用。

9. CF2041C Cube

看到 n\le 12，考虑状态压缩 DP。设计 f_{i,p_1,p_2} 表示选完第 i 层后，行的状态为 p_1，列的状态为 p_2。转移是显然的：枚举行、列可选的的状态中选择 (i,j,k)，那么：

f_{i+1,p_1|2^j,p_2|2^k}=\min(f_{i+1,p_1|2^j,p_2|2^k},f_{i,p_1,p_2}+a_{i+1,j,k})

如果说我们对于每个 i 暴力枚举状态 p_1,p_2，判断是否合法，这样的时间复杂度是 O(2^{2n}n^3) 的；我们考虑一个优化：我们先提前存下每个 i 对应的合法状态，然后 p_1,p_2 从这些合法状态中枚举，时间复杂度是 O(2^{2n}n^2) 的，可以通过。

:::success[AC code]

const int N=12;
int n;
int a[N+1][N][N];
vector<int> p[N+1];
int f[N+1][1<<N][1<<N];
signed main()
{
    read(n);    
    rep(i,1,n) rep(j,0,n-1) rep(k,0,n-1)
        read(a[i][j][k]);
    rep(i,0,(1<<n)-1) p[__builtin_popcount(i)].pb(i);
    rep(i,1,n) for(auto p1:p[i]) for(auto p2:p[i])
        f[i][p1][p2]=INF;
    f[0][0][0]=0;
    rep(i,0,n-1) for(auto p1:p[i]) rep(j,0,n-1) if(((1<<j)&p1)==0)
        for(auto p2:p[i]) rep(k,0,n-1) if(((1<<k)&p2)==0) f[i+1][p1|(1<<j)][p2|(1<<k)]=
            min(f[i+1][p1|(1<<j)][p2|(1<<k)],f[i][p1][p2]+a[i+1][j][k]);
    print(f[n][(1<<n)-1][(1<<n)-1]);
    return 0;
}

:::

总结

状态的合理剪枝可以帮助你在考场上获得更多的分数，并且极为实用。还有记得：unsigned long long 的乘法取模比 long long 快。

最后

感谢 @jucason_xu！！！

有问题可以在评论区、私信说。感谢各位阅读。