算法学习笔记之斜率优化

前言

近几日回归竞赛后，便开始学些新东西了（对于蒟蒻来说）。这几天就连续的更新一下自己对斜率优化的学习过程。

1.思想

在一些常见的DP题中，可能会出现形如\(f[i]=\min/\max(f[j]+(sum[i]-sum[j])^2\)的转移方程式。
这时，我们就可以把后面的二次项展开：

\(f[i]=f[j]+sum[i]^2+sum[j]^2+2*sum[i]*sum[j]\)

\(f[j]+sum[j]^2=2*sum[i]*sum[j]+sum[i]^2+f[i]\)

再将\(2*sum[i]\)看做斜率，\(f[j]+sum[j]\)看做\(y\)，\(sum[j]\)看做\(x\)，则\(f[i]+sum[i]^2\)就成了截距，由于斜率不变，我们就用一个单调队列来维护之前的决策点，再将该斜率的直线从下到上平移，找到斜率最小的点。所以我们既要维护一个下凸包。

2.实现

利用单调队列维护当前的决策点构成的凸包，在遍历到\(i\)时，先将斜率小于当前斜率的决策点全部删去。然后对头的点就是当前的最优决策点，更新当前值后，在队尾加点时保证斜率的单调递增即可。

如果要维护一个上凸包就反着来（emm）

1.去除队头，比较斜率；

2.更新答案；

3.加入队伍，维护单调性质；

3.典例

1.[HNOI2008]玩具装箱 ( 模板斜率优化DP )

先理解一下题意，总而言之就是划分一个序列，满足贡献值最小的情况。

转移方程：\(f(x)=\min_{i=1}^{x-1}(f(i)+(j-i+\sum^j_{k=i}c[k]-L)^2)\)

\[sum[i]=\sum_{k=1}^{k<=i}{c[k]} \]

\[a[i]=sum[i]+i,b[i]=sum[i]+i+L+1 \]

\[dp[i]=dp[j]+(a[i]−b[j])^2 \]

\[dp[i]=dp[j]+a[i]^2+b[j]^2-2⋅a[i]⋅b[j] \]

\[2⋅a[i]⋅b[j]+dp[i]−a[i] ^2 =dp[j]+b[j] ^2\]

把\(b[j]\)看成x,斜率就是\(2⋅a[i]\)，截距就是\(dp[i]-a[i]^2\)。

且斜率是单调递增的。

然后就要寻找最小的截距，就使用单调队列维护下凸包了~~

点击查看代码

while(head<tail&&slope(Q[head],Q[head+1])<2*a(i)) ++head;//维护队首元素
dp[i]=dp[Q[head]]+(a(i)-b(Q[head]))*(a(i)-b(Q[head]));//更新当前
while(head<tail&&slope(i,Q[tail-1])<slope(Q[tail-1],Q[tail])) --tail;//加入队尾
Q[++tail]=i;

完整代码：

点击查看代码

#include<bits/stdc++.h>
using namespace std;
#define db double
#define ll long long
const int maxn=50010;
int n,L;
db sum[maxn],dp[maxn];
int head,tail,Q[maxn];
inline db a(int i){
	return sum[i]+i;
}
inline db b(int i){
	return a(i)+L+1;
}
inline db X(int i){
	return b(i);
}
inline db Y(int i){
	return dp[i]+b(i)*b(i);
}
inline db slope(int i,int j){
	return (Y(i)-Y(j))/(X(i)-X(j));
}
int main(){
	scanf("%d%d",&n,&L);
	for(int i=1;i<=n;i++){
		scanf("%lf",&sum[i]);
		sum[i]+=sum[i-1];
	}	
	head=tail=1;
	for(int i=1;i<=n;i++){
		while(head<tail&&slope(Q[head],Q[head+1])<2*a(i)) ++head;
		dp[i]=dp[Q[head]]+(a(i)-b(Q[head]))*(a(i)-b(Q[head]));
		while(head<tail&&slope(i,Q[tail-1])<slope(Q[tail-1],Q[tail])) --tail;
		Q[++tail]=i;
 	}
	
	printf("%lld\n",(ll)dp[n]);
	return 0;
}

2.[CEOI2004]锯木厂选址 (维护一个上凸包)

状态转移方程有亿点好想：
\(ans=\min(ans,tot-dis[j]*s[j]-dis[i]*(s[i]-s[j]))\)

完整代码：

#include<bits/stdc++.h>
using namespace std;
int n;
const int N=1e5;
long long dp[N],s[N],wdis[N],dis[N];
long long tot;
int q[N];
inline double calc(int j,int k){
	return (double) (dis[j]*s[j]-dis[k]*s[k])/(s[j]-s[k]);
}
int main(){
	cin>>n;
	for(int i=1;i<=n;i++){
		int a,b;
		cin>>s[i]>>dis[i];
		s[i]=s[i-1]+s[i];
		tot+=s[i]*dis[i];
	}
	//cout<<tot<<endl;
	//cout<<s[n]<<endl;
	for(int i=n;i>=1;i--){
		dis[i]=dis[i+1]+dis[i];
	//	cout<<dis[i]<<endl;
	}
	int head=1,tail=1;
	long long ans=1e10;
	for(int i=1;i<=n;i++){
		while(head<tail&&calc(q[head],q[head+1])>=dis[i]) ++head;
		ans=min(ans,tot-dis[q[head]]*s[q[head]]-dis[i]*(s[i]-s[q[head]]));
		while(head<tail&&calc(q[tail],i)>=calc(q[tail-1],q[tail]))--tail;
		q[++tail]=i;
	}
	cout<<ans<<endl;
	return 0;
}

3.[APIO2010]特别行动队 （也是维护一个上凸包）

转移方程有点难：

\[dp[i]=\max(dp[j]+a*(sum[i]-sum[j])^2+b*(sum[i]-sum[j])+c \]

\[2*a*sum[i]*sum[j]+dp[i]-a*sum[i]^2-b*sum[i]-c=a*sum[j]^2+dp[j]-b*sum[j] \]

2.

\[x(i)=sum[i] \]

\[y(i)=a*sum[i]^2+dp[j]-b*sum[j] \]

\[k(i)=2*a*sum[i] \]

\[b(i)=dp[i]-a*sum[i]^2-b*sum[i]-c \]

3.

\[k(i)*x(j)+b(i)=y(i) \]

完整代码：

#include<bits/stdc++.h>
#define y(i) (dp[i]+a*sum[i]*sum[i]-b*sum[i])
#define int long long
#define re register
using namespace std;
int a,b,c;
int n;
const int N=1e6+100;
int sum[N];
int q[N];
int dp[N];
inline double slope(int j,int k){
	return 1.0*(y(j)-y(k))/(sum[j]-sum[k]);
} 
signed main(){
	scanf("%lld",&n);
	scanf("%lld %lld %lld",&a,&b,&c);
	for(re int i=1;i<=n;i++){
		scanf("%lld",&sum[i]);
		sum[i]+=sum[i-1];
	}
	re int head=1,tail=1;
	for(re int i=1;i<=n;i++){
		while(head<tail&&slope(q[head],q[head+1])>2*a*sum[i]) head++;
		dp[i]=(1ll*(sum[i]-sum[q[head]])*(sum[i]-sum[q[head]])*a+1ll*b*(sum[i]-sum[q[head]])+c+dp[q[head]]);
		while(head<tail&&slope(q[tail],i)>=slope(q[tail-1],q[tail])) --tail;
		q[++tail]=i;
	}
	printf("%lld\n",dp[n]);
	
	return 0;
} 

4.P2365 任务安排

一本通上经典的例题，最主要的思想是费用的提前计算，每次开机时把后面增加的费用提前计算好。

\(f[i]=\min _{0\le j \lt i}(f[j]+S*c(j,n)+t_i*c(j,i))\)

\(c(a,b)=\sum _{i\le b}^{i=a+1} c_i\)

这样前缀和优化即可。

\(f_i=f_j+S*(sc_n-sc_j)+st_i*(sc_i-sc_j)\)

\(f_i=f_j+S*sc_n-S*sc_j+st_i*st_i-st_i*sc_j\)

\(f_j=st_i*st_j+S*sc_j-S*sc_n-st_i*sc_i+f_i\)

\(f_j=(S+st_i)*st_j-S*sc_n-st_i*sc_i+f_i\)

\(Y=f_j\)

\(X=st_j\)

\(K=(S+st_i)\)

\(B=f_i-S*sc_n-st_i*sc_i\)

所以维护一个下图宝即可。

#include<bits/stdc++.h>
using namespace std;
long long f[300010],sumt[300010],sumc[300010];
int q[300010],n,s;


int main(){
	cin>>n>>s;
	for(int i=1;i<=n;i++){
		int t,c;
		cin>>t>>c;
		sumt[i]=sumt[i-1]+t;
		sumc[i]=sumc[i-1]+c;
	}
	memset(f,0x3f,sizeof f);
	f[0]=0;
	int l=1;
	int r=1;
	q[1]=0;
	for(int i=1;i<=n;i++){
		while(l<r&&(f[q[l+1]]-f[q[l]])<=(s+sumt[i])*(sumc[q[l+1]]-sumc[q[l]])) l++;
		f[i]=f[q[l]]-(s+sumt[i])*sumc[q[l]]+sumt[i]*sumc[i]+s*sumc[n];
		while(l<r&&(f[q[r]]-f[q[r-1]])*(sumc[i]-sumc[q[r]])>=(f[i]-f[q[r]])*(sumc[q[r]]-sumc[q[r-1]])) r--;
		q[++r]=i;
	
	} 
	cout<<f[n]<<endl;
return 0;
}

5.P5785 [SDOI2012]任务安排

思路：由于T_i的值有正负所以无法保证决策点的y值单调递增，所以会出现这样的图像：

就要二分查找决策点，取代单调队列的队头。

代码：

#include<bits/stdc++.h>
#define int long long
using namespace std;
long long f[300010],sumt[300010],sumc[300010];
int q[300010],n,s;
int Y(int p){
	return f[p];
}
int X(int p){
	return sumc[p];
}
int K(int p){
	return s+sumt[p];
}
int Search(int l,int r,int S)
{
	int mid=0,res=r;
	while(l<r){
		mid=l+r>>1;
		if(Y(q[mid+1])-Y(q[mid])>S*(X(q[mid+1])-X(q[mid]))){
			r=mid;
			res=mid;
		}else{
			l=mid+1;
		}
	}
	return q[res];
}
signed main(){
	cin>>n>>s;
	for(int i=1;i<=n;i++){
		int t,c;
		cin>>t>>c;
		sumt[i]=sumt[i-1]+t;
		sumc[i]=sumc[i-1]+c;
	}
	memset(f,0x3f,sizeof f);
	f[0]=0;
	int l=1,r=1;
	q[1]=0;
	for(int i=1;i<=n;i++){
		int p=Search(l,r,K(i));
		f[i]=f[p]-(s+sumt[i])*sumc[p]+sumt[i]*sumc[i]+s*sumc[n];
		while(l<r&&(f[q[r]]-f[q[r-1]])*(sumc[i]-sumc[q[r]])>=(f[i]-f[q[r]])*(sumc[q[r]]-sumc[q[r-1]])) r--;
		q[++r]=i;
	
	} 
	cout<<f[n]<<endl;
return 0;
}