每日算法 — 使用java实现斜率优化DP:单调队列维护凸包加速动态规划

引言

动态规划(DP)是算法竞赛和工程实践中最核心的技术之一。当DP转移方程形如 $dp[i] = \min_{j < i}{dp[j] + a[i] \cdot b[j] + c[i] + d[j]}$ 时,朴素的 $O(n^2)$ 实现往往难以满足大规模数据的效率要求。斜率优化通过维护下凸壳(Lower Convex Hull)配合单调队列,将时间复杂度降至 $O(n)$。本文将用Java从零实现斜率优化DP,详解其几何直观与代码细节。

问题建模:什么样的DP可以斜率优化

斜率优化适用于形如以下的DP转移:

$$dp[i] = \min_{j < i} {dp[j] + a[i] \cdot b[j] + c[i] + d[j]}$$

或更一般的形式:

$$dp[i] = \min_{j < i} {dp[j] + a[i] \cdot b[j]} + c[i]$$

其中 $a[i]$ 和 $b[j]$ 为已知序列,且 $a[i]$ 具有单调性(递增或递减)。关键观察:若将每个决策点 $j$ 视为二维平面上的点 $(b[j], dp[j])$,则转移方程等价于用斜率为 $-a[i]$ 的直线去”切”这些点,寻找最小的截距。

经典例题:划分序列最小代价

问题描述:给定长度为n的序列和正整数m,将其划分为若干连续段。每段 $[l, r]$ 的代价为 $(\sum_{i=l}^{r} a_i)^2 + m$。求最小总代价。

设前缀和为 $S[i]$,状态转移方程为:

$$dp[i] = \min_{0 \leq j < i} {dp[j] + (S[i] – S[j])^2 + m}$$

展开后:

$$dp[i] = \min_{0 \leq j < i} {dp[j] + S[j]^2 – 2 \cdot S[i] \cdot S[j]} + S[i]^2 + m$$

令 $x_j = S[j]$,$y_j = dp[j] + S[j]^2$,则:

$$dp[i] = \min_{j < i} {y_j – 2 \cdot S[i] \cdot x_j} + S[i]^2 + m$$

这正符合斜率优化的标准形式,其中斜率 $k_i = 2 \cdot S[i]$ 随 $S[i]$ 单调递增。

几何直观:下凸壳与最优决策点

将每个决策点 $j$ 表示为点 $P_j = (x_j, y_j)$。对于固定的 $i$,我们需要最小化 $y_j – k_i \cdot x_j$,即寻找一条斜率为 $k_i$ 的直线,使其经过某个 $P_j$ 时截距最小。

核心洞察:最优决策点只可能出现在这些点的下凸壳上。当 $k_i$ 单调递增时,最优决策点在凸壳上单调向右移动,可以用单调队列维护。

import java.util.*;

/**
 * 斜率优化DP核心实现
 * 以"序列划分最小代价"问题为例
 */
public class SlopeOptimizationDP {

    /**
     * 表示凸壳上的一个点(决策点)
     */
    static class Point {
        long x; // 对应 S[j]
        long y; // 对应 dp[j] + S[j]^2
        int index; // 原始下标j

        Point(long x, long y, int index) {
            this.x = x;
            this.y = y;
            this.index = index;
        }

        /**
         * 计算两点之间的斜率
         * 斜率 = (y2 - y1) / (x2 - x1)
         */
        static double slope(Point p1, Point p2) {
            if (p1.x == p2.x) {
                return p1.y < p2.y ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
            }
            return (double)(p2.y - p1.y) / (double)(p2.x - p1.x);
        }

        /**
         * 使用交叉相乘避免浮点精度问题
         * 判断p1->p2的斜率是否小于等于p2->p3的斜率
         * 即 (y2-y1)*(x3-x2) <= (y3-y2)*(x2-x1)
         * 用于维护下凸壳
         */
        static boolean isBad(Point p1, Point p2, Point p3) {
            // 若三点共线或形成上凸,则p2是"坏点"
            // 使用交叉相乘:slope(p1,p2) >= slope(p2,p3) 时p2在凸壳内部或上凸
            return (p2.y - p1.y) * (p3.x - p2.x) >= (p3.y - p2.y) * (p2.x - p1.x);
        }
    }

    /**
     * 使用斜率优化求解序列划分问题
     * @param a 原序列
     * @param m 每段固定代价
     * @return 最小总代价
     */
    public static long solvePartition(long[] a, long m) {
        int n = a.length;

        // 前缀和
        long[] prefix = new long[n + 1];
        for (int i = 1; i <= n; i++) {
            prefix[i] = prefix[i - 1] + a[i - 1];
        }

        long[] dp = new long[n + 1];
        dp[0] = 0; // 空前缀代价为0

        // 单调队列维护下凸壳上的决策点
        ArrayDeque<Point> deque = new ArrayDeque<>();
        // 初始决策点 j=0: x=S[0]=0, y=dp[0]+S[0]^2=0
        deque.addLast(new Point(0, 0, 0));

        for (int i = 1; i <= n; i++) {
            long si = prefix[i];
            long ki = 2 * si; // 当前斜率

            // Step 1: 从队首弹出非最优决策点
            // 若队首两个点中,第二个更优,则弹出队首
            while (deque.size() >= 2) {
                Point first = deque.peekFirst();
                Point second = deque.get(1);

                // 比较在first和second处计算的dp值
                // dp via first: first.y - ki * first.x
                // dp via second: second.y - ki * second.x
                // 若后者更小,first不再最优
                long val1 = first.y - ki * first.x;
                long val2 = second.y - ki * second.x;

                if (val2 <= val1) {
                    deque.pollFirst();
                } else {
                    break;
                }
            }

            // Step 2: 使用队首最优决策计算dp[i]
            Point best = deque.peekFirst();
            dp[i] = best.y - ki * best.x + si * si + m;

            // Step 3: 将当前点加入凸壳,维护下凸性质
            Point newPoint = new Point(si, dp[i] + si * si, i);

            while (deque.size() >= 2) {
                Point last = deque.peekLast();
                Point secondLast = deque.get(deque.size() - 2);

                // 若加入新点后形成上凸或共线,弹出队尾
                if (Point.isBad(secondLast, last, newPoint)) {
                    deque.pollLast();
                } else {
                    break;
                }
            }

            deque.addLast(newPoint);
        }

        return dp[n];
    }

精度优化:完全避免浮点数

上述实现中,val1 和 val2 的比较已使用整数运算,但 Point.isBad 中的乘法可能导致64位溢出。对于极端数据(如 $n=10^6$, $a_i=10^9$),需要使用 BigInteger 或 __int128(C++)级别的精度。Java中可用 java.math.BigInteger:

    /**
     * 高精度版本:使用BigInteger避免溢出
     */
    static boolean isBadBigInteger(Point p1, Point p2, Point p3) {
        java.math.BigInteger left = java.math.BigInteger.valueOf(p2.y - p1.y)
            .multiply(java.math.BigInteger.valueOf(p3.x - p2.x));
        java.math.BigInteger right = java.math.BigInteger.valueOf(p3.y - p2.y)
            .multiply(java.math.BigInteger.valueOf(p2.x - p1.x));
        return left.compareTo(right) >= 0;
    }

斜率优化的推广形式

斜率优化不仅适用于求最小值,也可用于求最大值(维护上凸壳)。核心区别在于不等号方向:

    /**
     * 求最大值版本:维护上凸壳
     * 判断条件取反
     */
    static boolean isBadForMax(Point p1, Point p2, Point p3) {
        // 维护上凸壳:slope(p1,p2) <= slope(p2,p3) 时p2是坏点
        return (p2.y - p1.y) * (p3.x - p2.x) <= (p3.y - p2.y) * (p2.x - p1.x);
    }

    /**
     * 通用斜率优化DP框架
     * @param n 状态数
     * @param kFunc 斜率函数 k(i)
     * @param xFunc x坐标函数 x(j)
     * @param yFunc y坐标函数 y(j) = dp[j] + 某表达式
     * @param extra extra(i) 只与i有关的项
     */
    public static long genericSlopeOptimization(
            int n,
            java.util.function.IntToLongFunction kFunc,
            java.util.function.IntToLongFunction xFunc,
            java.util.function.LongUnaryOperator yFunc, // 需要dp[j]作为输入
            java.util.function.IntToLongFunction extra) {

        long[] dp = new long[n + 1];
        ArrayDeque<Point> deque = new ArrayDeque<>();
        deque.addLast(new Point(xFunc.applyAsLong(0), yFunc.applyAsLong(dp[0]), 0));

        for (int i = 1; i <= n; i++) {
            long ki = kFunc.applyAsLong(i);

            // 弹出队首非最优点
            while (deque.size() >= 2) {
                Point f = deque.peekFirst();
                Point s = deque.get(1);
                if (s.y - ki * s.x >= f.y - ki * f.x) { // 最大值版本
                    deque.pollFirst();
                } else break;
            }

            Point best = deque.peekFirst();
            dp[i] = best.y - ki * best.x + extra.applyAsLong(i);

            Point np = new Point(xFunc.applyAsLong(i), dp[i] + /* depends on problem */ 0, i);
            // ... 维护凸壳
        }
        return dp[n];
    }

完整测试程序

    /**
     * 主程序:测试斜率优化DP
     */
    public static void main(String[] args) {
        // 测试用例1:小规模验证正确性
        long[] a1 = {1, 2, 3, 4, 5};
        long m1 = 1;
        long result1 = solvePartition(a1, m1);
        System.out.println("测试1结果: " + result1);
        // 手工验证:最优划分为 [1,2,3,4,5] 整体一段
        // cost = (15)^2 + 1 = 226

        // 测试用例2:全负数
        long[] a2 = {-5, -4, -3, -2, -1};
        long m2 = 10;
        long result2 = solvePartition(a2, m2);
        System.out.println("测试2结果: " + result2);

        // 测试用例3:大规模随机数据(性能测试)
        int n = 200000;
        long[] a3 = new long[n];
        Random rand = new Random(42);
        for (int i = 0; i < n; i++) {
            a3[i] = rand.nextInt(10000) + 1;
        }
        long start = System.currentTimeMillis();
        long result3 = solvePartition(a3, 1000);
        long end = System.currentTimeMillis();
        System.out.println("测试3结果: " + result3);
        System.out.println("大规模测试耗时: " + (end - start) + "ms");

        // 测试用例4:与O(n^2)暴力解法对比(小规模)
        long[] a4 = {3, 1, 4, 1, 5, 9, 2, 6};
        long m4 = 5;
        long fastResult = solvePartition(a4, m4);
        long bruteResult = bruteForcePartition(a4, m4);
        System.out.println("斜率优化结果: " + fastResult);
        System.out.println("暴力解法结果: " + bruteResult);
        System.out.println("结果一致性: " + (fastResult == bruteResult));
    }

    /**
     * O(n^2)暴力解法,用于验证斜率优化的正确性
     */
    public static long bruteForcePartition(long[] a, long m) {
        int n = a.length;
        long[] prefix = new long[n + 1];
        for (int i = 1; i <= n; i++) {
            prefix[i] = prefix[i - 1] + a[i - 1];
        }

        long[] dp = new long[n + 1];
        Arrays.fill(dp, Long.MAX_VALUE);
        dp[0] = 0;

        for (int i = 1; i <= n; i++) {
            for (int j = 0; j < i; j++) {
                long sum = prefix[i] - prefix[j];
                long cost = dp[j] + sum * sum + m;
                dp[i] = Math.min(dp[i], cost);
            }
        }
        return dp[n];
    }
}

复杂度分析

实现方式 时间复杂度 空间复杂度 适用场景
暴力DP $O(n^2)$ $O(n)$ n <= 5000
斜率优化 $O(n)$ $O(n)$ n <= 10^6,单调斜率
分治优化 $O(n \log n)$ $O(n)$ 决策点不单调
李超树 $O(n \log n)$ $O(n)$ 斜率不单调

常见错误与调试技巧

  1. 凸壳方向错误:最小值维护下凸壳(斜率递增),最大值维护上凸壳(斜率递减)
  2. 精度溢出:long 乘法溢出时改用 BigInteger 或 double 配合容差
  3. 斜率单调性假设失效:若 $k_i$ 不单调,需改用二分查找或李超线段树
  4. 初始点遗漏:务必加入 $j=0$ 的初始决策点

总结

斜率优化DP的本质是将代数问题几何化:把每个决策点看作平面上的点,DP转移转化为用直线”切”凸壳寻找最优截距。单调队列在此扮演了高效维护凸壳顶点的角色,使得每个点最多入队出队一次,整体达到线性复杂度。掌握这一技巧后,可将大量 $O(n^2)$ 的动态规划优化至 $O(n)$,在竞赛和工程实践中具有极高价值。