动态规划的斜率优化

发布于 2018-12-01  4,340 次阅读

在某些情况下, 动态规划的转移数目过多导致即使使用动态规划算法也会令复杂度过高. 但往往, 在某些情况下状态转移方程具有某些特别的性质, 对转移方程进行变形就可以转化为某种斜率的形式, 利用斜率的几何意义构造凸集, 通过维护一个单调双端队列使得​\( O(n^2) \)​复杂度的状态转移优化至​\( O(n) \)​.

下面结合几个例子讲一讲。

有关动态规划

        动态规划最经典的应用就是背包问题了, 与背包问题相关的问题可见: 背包九讲

       另外的可参阅 算法导论 的相关章节, 上面有详细深刻的讲解和许多极具启发性的例子.

       其余内容不再赘述。

两个例题

POJ2018

链接: POJ2018

题目大意: 给出一个长为N的整数序列, 求一个区间, 区间的长度至少为F, 并且使得这个区间的均值最大.

        严格地说, 这不是一个动态规划的题目, 但是它很好地体现了斜率优化的思想. 如果用​\( sum[i] \)​表示序列中前i个数的和, 那么区间​\( [i, j] \)​的均值即为 ​\( k(i,j)=\frac{sum[i]-sum[j]}{i-j} \)​, 从这个形式可以看出它具有斜率的形式, 即是说如果区间​\( [a, b] \)​的均值高于区间​\( [c, d] \)​的均值, 那么斜率​\( k(a,b)>k(c,d) \)​, 如果将\( (i,sum[i]) \)画在坐标系里, ​\( (a,sum[a]),(b,sum[b])) \)​两点的连线线斜率大于​\( (c,sum[c]),(d,sum[d]) \)​.

        以区间的右端点对不同的区间进行划分, 对于以 i 为右端点的所有区间中, 要求最大的区间均值, 想象坐标系中一条经过\( (i,sum[i]) \)的直线, 其倾角\( \theta \)从​\( \frac{\pi}{2} \)​逐渐减小到​\( 0 \)​, 与这条直线最先相交的点对应的横坐标, 就是我们要求的区间的左端点. 这个几何直观给了我们一个hint: 如果横坐标分别为\( i,j,k \)的三点在坐标系中形成上凸, 即​\( k(i,j)>k(j,k) \)​, 那么无论区间右端点是那个点, 其对应的最优左端点都不可能是上凸起来的那个点j, 也就是说, 我们可以把所有出现上凸的点都不予考虑, 那么只需要维护一个下凸集, 从1开始向右逐个枚举区间的右端点, 那么将点​\( (i,sum[i]) \)​和 横坐标依次为​\( i-1,i-2,\cdots,1 \)​的点​\( (k,sum[k]) \)​相连, 如果某个处导致\( i,k,k-1 \)对应的几个点出现上凸, 那么删去k点, 直到出现满足下凸的点, 就可以跳出. 借由这样的操作, 我们维护得到了一个上凸集q, 我们只需要在这个上凸集中寻找该区间的左端点就可以了. 在寻找左端点时, 也可以利用考察凸性的办法, 因为旋转过程中第一个相交的点​\( q[k] \)​必然恰好是​\( (i,sum[i]),\ (q[k],sum[q[k]]),\ (q[k-1],sum[q[k-1]]) \)​凸性反转的点.需要注意的是, 在这个题目中区间有最小长度限制, 所以略有不同.

下面是我的AC代码:

#include <iostream>
#include <cstdio>

using namespace std;
const int MAXN = 100005;
int sum[MAXN], q[MAXN];

inline double k(int i, int j) { return (double) (sum[i] - sum[j]) / (double) (i - j); }

int main() {
    int N, F;
    double ans = -1000;
    scanf("%d%d", &N, &F);
    int x;
    for (int i = 1; i <= N; ++i) {
        scanf("%d", &x);
        sum[i] = sum[i - 1] + x;
    }
    int tail = 0, head = 0;
    for (int i = 0; i <= N - F; ++i) {
        for (; tail > head && k(i, q[tail]) < k(q[tail], q[tail - 1]); --tail);
        q[++tail] = i;
        int j;
        for (j = tail; j > head && k(i + F, q[j]) < k(q[j], q[j - 1]); --j);
        head = j;
        ans = max(ans, k(q[j], i + F));
    }
    printf("%d\n", (int) (ans * 1000));
    return 0;
}

CodeForces674C

链接: CodeForces674C

题意简述:

        (略QAQ)

        此题一开始有个小难点, 就是如何求期望. 但是稍微观察一下就会发现, 按照题目所设定的规则, 玩家肯定会从第一个level开始一个一个地通关, 直到通过一个region, 再直到通过所有的region. 由于期望具有线性性, 一个region的通关时间期望等于各个level通关时间期望之和. 当然, 一个level的通关时间期望是与region的划分方法有关的, 本题正是需要找到一个给出最佳划分的快速算法. 下面先给出一个求期望的方法:

        不失一般性, 用​\( S[i] \)​表示总的区间​\( [1,N] \)​的部分和. 现考察region​\( [r,i] \)​. 设该区间的前面若干level都被打败, 现在正在打第​\( j \)​个level(指标都按照整个区间数), 假设所需时间为​\( T_j \)​, 那么

\[ T_j=\frac{t_j}{S_j-S_r}\times 1+\frac{S_j-S_r-t_j}{S_j-S_r}\times(T_j+1) \]

        则通关该region的期望为诸​\( T_j=\frac{S_j-S_r}{t_j} \)​在该区间内求和, 不妨将其记为​\( E(r,i) \)​. 则​\( E(r,i)=SSh[i]-SSh[r]+S[r]harmon[r]-S[r]harmon[i] \)​, 其中 S[] 表示部分和, harmon[] 表示​\( \frac{1}{t_i} \)​的部分和, SSh[] 表示​\( \frac{S_i}{t_i} \)​ 的部分和, 引入这些部分和是为了简化问题.

        接下来构造递推关系. 对于待划分的区间, 考察划分方案中的最后一个区间的左端点, 以​\( dp[i][j] \)​记前 i 个level划分为 j 个region的最小期望时间, 则可得到

\[ dp[i][j]=min_{1\le k\le i-1}\{dp[k][j-1]+E(k+1,i)\} \]

        该递推式需要遍历​\( i,j \)​, 每一次遍历时需要计算​\( E(r,s) \)​, 复杂度​\( O(n) \)​, 故而总复杂度​\( O(n^2k) \)​. 我们尝试通过斜率优化的方法对于对这个动态规划进行优化. 假设​\( 1\le r < s\le i-1 \)​, 且在 i 点处, r 的决策优于 s 点, 即

\[ dp[r][j-1]+SSh[i]-SSh[r]+S[r]harmon[r]-S[r]harmon[i]\le dp[s][j-1]+SSh[i]-SSh[s]+S[s]harmon[s]-S[s]harmon[i] \]

, 代入表达式就可转化为斜率的形式:

\[ \frac{G[r][j-1]-G[s][j-1]}{S[r]-S[s]}>h[i] \]

其中

\[ G[x][j]=dp[x][j]+S[x]harmon[x]-SSh[x] \]

        然后就是很套路的操作了, 通过单调双端队列维护一个上凸集, 每次通过找切线找到当前最优解, 每个点仅入队出队一次, 故而复杂度为​\( O(nk) \)​.

        有一个需要注意的地方emmm, 由于需要预处理出部分和, 题目给的数据范围会导致部分和超过 int 上限, 需要用 long long 储存部分和.

我的AC代码:

#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>

using namespace std;
const int MAXN = 200005;
const int MAXK = 52;
double dp[MAXN][MAXK];
long long _t[MAXN], S[MAXN];
double harmon[MAXN], SSh[MAXN];

inline double G(int r, int j) { return dp[r][j] + S[r] * harmon[r] - SSh[r]; }

inline double K(int r, int s, int j) { return (G(r, j) - G(s, j)) / (S[r] - S[s]); }

inline void init(int n) {
    S[1] = _t[1];
    harmon[1] = 1.0 / _t[1];
    SSh[1] = 1.0;
    for (int i = 2; i <= n; ++i) {
        S[i] = S[i - 1] + _t[i];
        harmon[i] = harmon[i - 1] + 1.0 / _t[i];
        SSh[i] = SSh[i - 1] + (S[i] * (double) 1.0 / (double) _t[i]);
    }
    for (int i = 1; i <= n; ++i) { dp[i][1] = ((double) S[i] / (double) _t[i]); }
    for (int i = 2; i <= n; ++i) { dp[i][1] += dp[i - 1][1]; }
}

inline double get_ans(int n, int k) {
    init(n);
    int q[MAXN];
    int head = 0, tail = 0;
    for (int j = 2; j <= k; ++j) {
        head = tail = 0;
        for (int i = j; i <= n; ++i) {
            for (; tail > head + 1 && K(q[tail - 2], q[tail - 1], j - 1) > K(q[tail - 1], i - 1, j - 1); --tail);
            q[tail++] = i - 1;
            int r;
            for (r = head; tail > r + 1 && K(q[r], q[r + 1], j - 1) <= harmon[i]; ++r);
            head = r;
            dp[i][j] = dp[q[r]][j - 1] + SSh[i] - SSh[q[r]] - S[q[r]] * harmon[i] + S[q[r]] * harmon[q[r]];
        }
    }
    return dp[n][k];
}

int main() {
    int n, k;
    scanf("%d%d", &n, &k);
    for (int i = 1; i <= n; ++i) { cin >> _t[i]; }
    printf("%lf\n", get_ans(n, k));
    return 0;
}

终有一日, 仰望星空