斜率優化DP總結

liuchanglc發表於2020-10-22

前言

(本文中的圖片都由\(WPS\)出品)
\(DP\)\(OI\) 中重要的一部分
一般來說,因為 \(DP\) 會把之前的結果儲存下來,所以時間複雜度還是比較優秀的
但是在某些情況下,時間複雜度仍然超出了題目的限制,這是我們就要考慮對其進行優化
\(DP\) 的優化一般從狀態、決策、轉移三個方面去考慮
而斜率優化則是對決策進行優化的一種方法
它適用於類似 \(f[i]=min/max(a[i] \times b[j]+c[i]+d[j])\) 的方程

例題(洛谷 P4072 [SDOI2016]征途 )

題目描述

\(Pine\) 開始了從 \(S\) 地到 \(T\) 地的征途。

\(S\)地到\(T\)地的路可以劃分成 \(n\) 段,相鄰兩段路的分界點設有休息站。

\(Pine\)計劃用\(m\)天到達\(T\)地。除第\(m\)天外,每一天晚上\(Pine\)都必須在休息站過夜。所以,一段路必須在同一天中走完。

\(Pine\)希望每一天走的路長度儘可能相近,所以他希望每一天走的路的長度的方差儘可能小。

幫助\(Pine\)求出最小方差是多少。

設方差是\(v\),可以證明,\(v\times m^2\)是一個整數。為了避免精度誤差,輸出結果時輸出\(v\times m^2\)

輸入格式

第一行兩個數 \(n\)\(m\)

第二行 \(n\) 個數,表示 \(n\) 段路的長度

輸出格式

一個數,最小方差乘以 \(m^2\) 後的值

輸入輸出樣例

輸入 #1

5 2
1 2 5 8 6

輸出 #1

36

說明/提示

對於 \(30\%\) 的資料,\(1 \le n \le 10\)

對於 \(60\%\) 的資料,\(1 \le n \le 100\)

對於 \(100\%\) 的資料,\(1 \le n \le 3000\)

保證從 \(S\)\(T\) 的總路程不超過 \(30000\)

分析

要對一個狀態轉移方程進行優化,首先要把最樸素的方程寫出來

在本題中,稍加推導即可寫出時間複雜度為 \(O(m \times n^2)\)的狀態轉移方程

\[s^2 \times m^2=\frac{(v_1-\overline{v})^2+(v_2-\overline{v})^2+...+(v_m-\overline{v})^2}m \times m^2 \]

\[s^2 \times m^2=((v_1-\overline{v})^2+(v_2-\overline{v})^2+...+(v_m-\overline{v})^2) \times m \]

\[s^2 \times m^2=m\times \sum_{i=1}^mv_i^2+m^2 \times \overline {v}^2-2\times m \times \overline {v} \times sum[n] \]

又因為 \(\overline{v}=\frac{sum[n]}{m}\)

所以

\[s^2\times m^2=m\times \sum_{i=1}^mv_i^2-sum[n]\times sum[n] \]

後面的值是固定的,所以我們只需要讓前面的值最小化即可

我們設\(f[i][j]\)為前\(i\)天分成\(j\)段所得到的最小值

那麼就有

\[f[i][k]=\min(f[j][k-1]+(sum[i]-sum[j])^2) \]

展開就有

\[f[i][k]=f[j][k-1]+sum[i]^2+sum[j]^2-2\times sum[i] \times sum[j] \]

接下來就是本文的重點:如何用斜率優化這類方程

首先,你需要掌握一次函式 \(y=kx+b\) 的影像和性質

這應該問題不大

下面我們就要對方程進行移項,使其變成易於優化的形式

\[f[j][k-1]+sum[j]^2=2\times sum[i] \times sum[j]+f[i][k]-sum[i]^2 \]

我們發現,這和一次函式的解析式完全吻合

我們把\(f[j][k-1]+sum[j]^2\)看成\(y\)

\(2 \times sum[i]\)看成\(k\)

\(sum[j]\)看成\(x\)

\(f[i][k]-sum[i]^2\)看成\(b\)

這樣,對於每一個\(i\)來說,直線的\(k\)是確定的

我們要使\(f[i][k]\)最小,也就是要使\(b\)最小

我們可以把所有的\(j\)想象成空間座標為\((sum[j],f[j][k-1]+sum[j]^2)\)中的點

知道了斜率,知道了直線上的點,那麼這條直線就確定了

那麼我們考慮什麼樣的點使直線的\(b\)最小

直線\(l\)是我們要移動的直線,平面中的點是可以轉移的\(j\)

我們會發現噹噹前點和後一個點形成的直線的斜率恰好大於直線\(l\)的斜率時,由當前點轉移決策是最優的

在這裡要特別強調一下:本題中斜率\(k\) 和 橫座標 \(x\) 均為單調遞增的,對於 \(k\)\(x\)不單調遞增的情況,處理方式不同

這就是程式碼裡面的

while(head<tail && xl(q[head],q[head+1])<2*sum[j]) head++;

我們再去考慮什麼樣的點肯定不會對結果產生貢獻

上面的圖中\(2\)號節點是無論如何也不會更新其它節點的

因為\(1\)號節點或\(3\)號節點總會比它更優

這就是程式碼裡的

while(head<tail && xl(j,q[tail-1])>=xl(j,q[tail])) tail--;

整個過程就相當於維護了一個下凸包

在求斜率的函式中,我們要判掉 \(x\) 相等的情況,在某些時候,還要判掉 \(y\) 相等的情況

double X(int id){
	return (double)sum[id];
}
double Y(int id){
	return (double)(g[id]+sum[id]*sum[id]);
}
double xl(int i,int j){
	if(std::fabs(X(i)-X(j))<eps){
		if(std::fabs(Y(i)-Y(j))<eps) return 0;
		else if(Y(i)>Y(j)) return 1e18;
		else return -1e18;
	}
	return (Y(i)-Y(j))/(X(i)-X(j));
}

擴充一:斜率不單調但x單調

如果斜率不是單調遞增,我們就不能從前面清空佇列直接轉移

比如上面這幅圖如果在遇到直線 \(m\) 時一直從前清空佇列的話那麼就會把\(3\)號決策點彈出佇列

但是如果之後遇到一個斜率比較小的直線\(l\)那麼就不能轉移到最優解

典型的例題是 洛谷P5785 [SDOI2012]任務安排

樸素的狀態轉移方程為

$ f[i] = f[j] + (sumc[i] - sumc[j]) * sumt[i] + s * (sumc[n] - sumc[j]);$

在這一道題中,作為斜率的 \(sumt\) 不再單調

但是 \(x\) 之仍然是單調的

所以我們可以用維護一個斜率單調的佇列

每次在佇列中二分答案

值得一提的是,出題人精心準備了卡精度的資料,所以我們要把二分時的除法改為乘法

程式碼

#include <cstdio>
#define rg register
inline int read() {
    rg int x = 0, fh = 1;
    rg char ch = getchar();
    while (ch < '0' || ch > '9') {
        if (ch == '-')
            fh = -1;
        ch = getchar();
    }
    while (ch >= '0' && ch <= '9') {
        x = (x << 1) + (x << 3) + (ch ^ 48);
        ch = getchar();
    }
    return x * fh;
}
typedef long long ll;
const int maxn = 1e6 + 5;
int t[maxn], c[maxn], n, s, q[maxn], head, tail;
ll f[maxn], sumt[maxn], sumc[maxn];
double Y(int i) { return (double)(f[i] - s * sumc[i]); }
double X(int i) { return (double)sumc[i]; }
double xl(int i, int j) {
    if (X(i) == X(j)) {
        if (Y(i) > Y(j))
            return 1e18;
        else
            return -1e18;
    }
    return (double)(Y(i) - Y(j)) / (X(i) - X(j));
}
int ef(double now) {
    int l = head, r = tail, mids;
    while (l < r) {
        mids = (l + r) >> 1;
        if ((X(q[mids]) > X(q[mids + 1]) &&
             Y(q[mids]) - Y(q[mids + 1]) < now * (X(q[mids]) - X(q[mids + 1]))) ||
            (X(q[mids]) < X(q[mids + 1]) &&
             (Y(q[mids + 1]) - Y(q[mids]) < now * (X(q[mids + 1]) - X(q[mids])))))
            l = mids + 1;
        else
            r = mids;
    }
    return q[l];
}
int main() {
    n = read(), s = read();
    for (rg int i = 1; i <= n; i++) {
        t[i] = read();
        c[i] = read();
        sumt[i] = sumt[i - 1] + t[i];
        sumc[i] = sumc[i - 1] + c[i];
    }
    head = tail = 1;
    for (rg int i = 1; i <= n; i++) {
        rg int wz = ef((double)(sumt[i]));
        f[i] = f[wz] + (sumc[i] - sumc[wz]) * sumt[i] + s * (sumc[n] - sumc[wz]);
        while (head < tail && xl(i, q[tail - 1]) >= xl(i, q[tail])) tail--;
        q[++tail] = i;
    }
    printf("%lld\n", f[n]);
    return 0;
}

擴充套件二、x不單調

沒有找到 \(x\) 不單調但是 \(k\) 單調的題

但是卻有一道 \(x\) 不單調 \(k\) 也不單調的題洛谷P4655 [CEOI2017]Building Bridges

對於這道題,我們同樣可以寫出最樸素的方程

\(f[i]f[j]+(h[i]-h[j])*(h[i]-h[j])+sum[i-1]-sum[j]\)

神奇的一點是 \(h\) 陣列既作為直線的斜率又作為 \(x\)

而且 \(h\) 並不單調

這是我們就不能再用單調佇列去維護,因為凸包的形狀在不斷改變

本題可以用李超線段樹或者平衡樹動態維護凸包解決

但是還有一種 \(CDQ\) 分治離線處理的方法

直接人為地排出單調性,像普通單調佇列那樣維護就可以了

如果左側斜率遞增,並且左側編號小於右側,那麼可以通過單調佇列維護左側的凸包來更新右側答案

並且這樣一定能夠遍歷出每個節點的所有決策點

程式碼

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cstring>
#include<cmath>
#define rg register
inline int read(){
	rg int x=0,fh=1;
	rg char ch=getchar();
	while(ch<'0' || ch>'9'){
		if(ch=='-') fh=-1;
		ch=getchar();
	}
	while(ch>='0' && ch<='9'){
		x=(x<<1)+(x<<3)+(ch^48);
		ch=getchar();
	}
	return x*fh;
}
const int maxn=1e5+5;
const double eps=1e-6;
typedef long long ll;
int n,h[maxn],w[maxn];
ll f[maxn],sum[maxn];
double X(int id){
	return (double)h[id];
}
double Y(int id){
	return (double)(f[id]+1LL*h[id]*h[id]-sum[id]);
}
double xl(int i,int j){
	if(std::fabs(X(i)-X(j))<eps){
		if(std::fabs(Y(i)-Y(j))<eps) return 0;
		else if(Y(i)>Y(j)) return 1e18;
		else return -1e18;
	} else {
		return (Y(i)-Y(j))/(X(i)-X(j));
	}
}
bool cmp(int aa,int bb){
	return h[aa]<h[bb];
}
int tmp[maxn],p[maxn],q[maxn];
void solve(int l,int r){
	if(l==r) return;
	rg int mids=(l+r)>>1,head=l-1,tail=mids;
	for(rg int i=l;i<=r;i++){
		if(p[i]<=mids) tmp[++head]=p[i];
		else tmp[++tail]=p[i];
	}
	for(rg int i=l;i<=r;i++){
		p[i]=tmp[i];
	}
	solve(l,mids);
	head=1,tail=0;
	for(rg int i=l;i<=mids;i++){
		while(head<tail && xl(p[i],q[tail-1])>=xl(p[i],q[tail])) tail--;
		q[++tail]=p[i];
	}
	for(rg int i=mids+1;i<=r;i++){
		while(head<tail && xl(q[head+1],q[head])<=2.0*h[p[i]]) head++;
		f[p[i]]=std::min(f[p[i]],f[q[head]]+1LL*(h[p[i]]-h[q[head]])*(h[p[i]]-h[q[head]])+sum[p[i]-1]-sum[q[head]]);
	}
	solve(mids+1,r);
	head=l,tail=mids+1;
	for(rg int i=l;i<=r;i++){
		if(tail>r || (head<=mids && h[p[head]]<h[p[tail]])) tmp[i]=p[head++];
		else tmp[i]=p[tail++];
	}
	for(rg int i=l;i<=r;i++){
		p[i]=tmp[i];
	}
}
int main(){
	memset(f,0x3f,sizeof(f));
	n=read();
	for(rg int i=1;i<=n;i++){
		h[i]=read();
		p[i]=i;
	}
	for(rg int i=1;i<=n;i++){
		w[i]=read();
		sum[i]=sum[i-1]+w[i];
	}
	std::sort(p+1,p+1+n,cmp);
	f[1]=0;
	solve(1,n);
	printf("%lld\n",f[n]);
	return 0;
}

總結

雖然斜率優化看起來很難,但是熟能生巧,多打幾遍就能掌握

相關文章