[復習]斜率優化
好久沒寫過了，跟忘了沒啥區別了。
然后重新理解一遍這個東西，感覺我原來對於斜率優化的想法有着很大的問題。
所以這些東西舉例子重新推一推吧QwQ。
[HNOI2010]玩具裝箱
首先寫暴力\(O(n^2)\)的轉移，設\(S_i\)是\(C_i\)的前綴和。
然后把式子拆開，和\(j\)無關的直接移出去，只和\(j\)相關的放在一起，同時和\(i,j\)相關的放在一起。
那么分類之后就是這樣的：
- 和\(j\)無關的：\(i^2-2i+1+S_i^2+L^2-2LS_i+2iS_i-2iL-2S_i+2L\)
- 只和\(j\)有關的：\(f[j]+j^2+2j+S_{j}^2+2LS_j+2jS_j+2jL+2S_j\)
- 同時和\(i,j\)相關的：\(-2ij-2S_iS_j-2iS_j-2jS_i=-2(i+S_i)(j+S_j)\)
一共\(22\)項，似乎沒有什么問題。(其實可以直接令\(M_i=i-1+S_i-L\)，但是為了鍛煉拆式子能力就這樣吧......算了，我編不下去了.....)
那么把和\(j\)無關的部分記做\(pre[i]\)，只和\(j\)有關的記做\(y[j]\)，\(j+S_j\)記做\(x[j]\)，\(2(i+S_i)\)記做\(k_i\)。
那么轉移方程可以改寫成：
而\(k_i\)是一個常數，我們把后面這個式子理解為一個一次函數\(y=kx+b\)的形式，得到\(b=y-kx\)。
什么意思呢？平面上有若干個點\((x[j],y[j])\)，你要過這些點畫一條斜率為\(k_i\)的直線，使得其截距最小。
不難發現滿足條件的\(j\)一定在下凸殼上。
這里有個很優秀的性質，也就是\(k_i,x[j]\)都是單增的。
這樣子凸殼可以直接用單調隊列維護，而取最優值只需要每次找到凸殼左側最優位置就好啦。
#include<iostream>
#include<cstdio>
using namespace std;
#define ll long long
#define MAX 50050
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int n,c[MAX];
ll L,S[MAX],f[MAX],pre[MAX],x[MAX],y[MAX];
ll Sqr(ll x){return x*x;}
ll Calc(int i,int j){return f[j]+Sqr(i-j-1+S[i]-S[j]-L);}
int Q[MAX],h,t;
double Slope(int i,int j){return 1.0*((y[i]+f[i])-(y[j]+f[j]))/(x[i]-x[j]);}
int main()
{
n=read();L=read();
for(int i=1;i<=n;++i)c[i]=read(),S[i]=S[i-1]+c[i];
for(int i=1;i<=n;++i)pre[i]=1ll*i*i-2*i+1+S[i]*S[i]+L*L-2*L*S[i]+2*i*S[i]-2*i*L-2*S[i]+2*L;
for(int i=1;i<=n;++i)y[i]=1ll*i*i+2*i+S[i]*S[i]+2*L*S[i]+2*i*S[i]+2*i*L+2*S[i];
for(int i=1;i<=n;++i)x[i]=i+S[i];
f[0]=0;Q[h=t=1]=0;
for(int i=1;i<=n;++i)
{
while(h<t&&Calc(i,Q[h])>=Calc(i,Q[h+1]))++h;
f[i]=Calc(i,Q[h]);
while(h<t&&Slope(Q[t-1],i)<=Slope(Q[t-1],Q[t]))--t;
Q[++t]=i;
}
printf("%lld\n",f[n]);
return 0;
}
[APIO2010]特別行動隊
首先還是可以寫出\(O(n^2)\)的暴力\(dp\)。
還是把式子拆開之后按照前面的三類分開。
- 與\(j\)無關：\(aS_i^2+bS_i+c\)
- 只與\(j\)有關：\(f[j]+aS_j^2-bS_j\)
- 與\(i,j\)有關：\(-2aS_iS_j\)
還是和前面一樣，與\(j\)無關記做\(pre[i]\)，只與\(j\)有關記做\(y[j]\)，同時有關的剔除與\(i\)有關的部分之后就是\(x[j]=S_j\)。
轉移改寫成
似乎又是前面那種給定點和斜率求最大值截距的問題了。
這次是取最大值，且因為\(a<0\)，所以\(-2aS_i\)單增，因此維護的是上凸殼。
並且\(x[j]\)單調，所以可以直接拿單調隊列維護。
#include<iostream>
#include<cstdio>
using namespace std;
#define MAX 1000100
#define ll long long
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int n;ll A,B,C,f[MAX],S[MAX],x[MAX],y[MAX];
double Slope(int a,int b){return 1.0*((y[a]+f[a])-(y[b]+f[b]))/(x[a]-x[b]);}
ll Calc(int i,int j){return f[j]+A*(S[i]-S[j])*(S[i]-S[j])+B*(S[i]-S[j])+C;}
int Q[MAX],h,t;
int main()
{
n=read();A=read();B=read();C=read();
for(int i=1;i<=n;++i)S[i]=S[i-1]+read();
for(int i=1;i<=n;++i)y[i]=A*S[i]*S[i]-B*S[i];
for(int i=1;i<=n;++i)x[i]=S[i];
Q[h=t=1]=0;
for(int i=1;i<=n;++i)
{
while(h<t&&Calc(i,Q[h])<=Calc(i,Q[h+1]))++h;
f[i]=Calc(i,Q[h]);
while(h<t&&Slope(Q[t-1],i)>=Slope(Q[t-1],Q[t]))--t;
Q[++t]=i;
}
printf("%lld\n",f[n]);
return 0;
}
[SDOI2012]任務安排
考慮\(O(N^2)\)暴力\(dp\)，顯然每次把后面所有的費用在當前的時間中的貢獻一起計算。
記\(SF[i]\)為\(F[i]\)的后綴和，\(ST[i]\)為\(T[i]\)的前綴和，得到轉移：
假裝\(SF\)數組被向左平移了一位，這樣子也可以看做只和\(j\)相關的部分。
然后來分類：
- 和\(j\)無關：似乎並沒有
- 只和\(j\)相關：\(f[j]+SF[j+1]*(S-ST[j])\)
- 和\(i,j\)相關：\(SF[j+1]*ST[i]\)
然后就可以很套路的令\(y[j]=f[j]+SF[j+1]*(S-ST[j])\)，然后\(x[j]=SF[j+1]\)，然后每次就可以看做在斜率\(-ST[i]\)上詢問最小截距？
看起來很美好，但是我們發現了一個問題，\(T\)這個東西可以是負數，所以\(ST\)並不是單調的，所以不能直接用單調隊列來維護這個東西。而且\(x[j]\)這玩意是從右往左單減的，所以凸包還要反着維護。
那么凸殼還是可以直接維護的，只是每次詢問的斜率不單調而已，所以我們只需要每次在凸殼上二分這個斜率就好了。
注意一下這里因為要求的是截距的最小值，且斜率為負，所以需要求的是一個下凸殼。
這里為了防止掉精度(比如\(\Delta x=0\)這樣子)所以斜率的比較全部都乘過去變成乘積的比較。
#include<iostream>
#include<cstdio>
using namespace std;
#define ll long long
#define MAX 300300
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int n;
ll T[MAX],F[MAX],S,f[MAX],x[MAX],y[MAX];
ll Calc(int i,int j){return f[j]+F[j+1]*(S+T[i]-T[j]);}
bool chk1(int i,int j,int k){return 1.0*((y[i]+f[i])-(y[j]+f[j]))*(x[i]-x[k])<=1.0*((y[i]+f[i])-(y[k]+f[k]))*(x[i]-x[j]);}
bool chk2(int i,int j,ll k){return 1.0*(y[i]+f[i])-(y[j]+f[j])<=1.0*k*(x[i]-x[j]);}
int Q[MAX],top;
int Bound(ll k)
{
int l=1,r=top-1,ret=top;
while(l<=r)
{
int mid=(l+r)>>1;
if(chk2(Q[mid],Q[mid+1],k))r=mid-1,ret=mid;
else l=mid+1;
}
return Q[ret];
}
int main()
{
n=read();S=read();
for(int i=1;i<=n;++i)T[i]=read(),F[i]=read();
for(int i=n;i>=0;--i)F[i]+=F[i+1];
for(int i=1;i<=n;++i)T[i]+=T[i-1];
for(int i=0;i<=n;++i)y[i]=F[i+1]*(S-T[i]);
for(int i=0;i<=n;++i)x[i]=F[i+1];
Q[top=1]=0;
for(int i=1;i<=n;++i)
{
f[i]=Calc(i,Bound(-T[i]));
while(top>1&&chk1(Q[top-1],Q[top],i))--top;
Q[++top]=i;
}
printf("%lld\n",f[n]);
return 0;
}
[NOI2007]貨幣兌換
首先題目里的提示已經告訴我們每次要么花完所有錢，要么賣出所有金券。
我們設\(f[i]\)表示第\(i\)天結束時能夠得到的最大價值的錢，考慮這個東西怎么轉移。那么我們枚舉上一次買入金券是哪一天，把那一天結束時能夠得到的錢全部換成金券，然后到當前第\(i\)天賣掉，那么我們就可以寫出一個\(O(n^2)\)的暴力了。
令\(VB[j]=f[j]/(SA[j]*Rate[j]+SB[j])\)，得到轉移：
其中\(SA[i]\)表示第\(i\)天金券\(A\)的價值，\(SB[i]\)同理，\(Rate[i]\)同題面的含義。
類似前面分類：
和\(j\)無關：似乎並沒有
之和\(j\)有關：似乎也沒有
同時和\(i,j\)有關：似乎都是的。
看起來似乎並不好寫成斜率的形式了？
直接令\(x[j]=VB[j]*Rate[j],y[j]=VB[j]\)，那么轉移可以寫成：
然后強行把其中一項\(i\)給提出來，變成：
於是問題又變成了斜率的問題，即給定一對點\((x[j],y[j])\)，求用斜率\(-\frac{SA[i]}{SB[i]}\)的直線穿過這些點之后的最大截距。
於是我們就要來維護一個上凸殼。
那么問題來了，\(x[j]\)完全不單調，應該怎么維護呢？
方法有兩種，第一種是用\(CDQ\)分治處理，另外一種是用\(Splay\)動態維護凸殼。當然，你如果願意也可以使用二進制分組來強行把\(CDQ\)給在線化，這里就不討論這個東西。
先說第一種。
首先明白這樣一個事實：我們斜率優化並不是只能從一個凸包上取答案，我們可以把所有的點分組建立凸包，分別在每個凸包上詢問取最優值。那么\(CDQ\)分治就可以很好的處理這個問題。
因為我們在詢問一個點之前，它前面的所有點的值都必須算出，所以\(CDQ\)的步驟是：遞歸處理左側，處理當前整個區間，處理右側，這個順序不要弄錯。
這樣子只需要給左半邊的所有點排序構建一個上凸殼，然后給右側的所有點依次二分做詢問就好了。
注意這里寫的時候不要忘記還有\(f[i]=\max\{f[i],f[i-1],f[0]\}\)這個轉移了。
#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;
#define MAX 100100
int n;double f[MAX];
double SA[MAX],SB[MAX],Rate[MAX];
double Calc(int i,int j)
{
double B=f[j]/(SA[j]*Rate[j]+SB[j]);
double A=Rate[j]*B;
return SA[i]*A+SB[i]*B;
}
struct Node{double x,y;int id;}p[MAX],tmp[MAX];
int Q[MAX],top;
double Slope(int i,int j)
{
if(fabs(p[i].x-p[j].x)<1e-9)return 1e-15;
return (p[i].y-p[j].y)/(p[i].x-p[j].x);
}
int Bound(double K)
{
int l=1,r=top-1,ret=top;
while(l<=r)
{
int mid=(l+r)>>1;
if(Slope(Q[mid],Q[mid+1])<=K)r=mid-1,ret=mid;
else l=mid+1;
}
return Q[ret];
}
void CDQ(int l,int r)
{
if(l==r)
{
f[l]=max(f[l],f[l-1]);p[l].id=l;
p[l].y=f[l]/(SA[l]*Rate[l]+SB[l]);
p[l].x=p[l].y*Rate[l];
return;
}
int mid=(l+r)>>1;CDQ(l,mid);top=0;
for(int i=l;i<=mid;++i)
{
while(top>1&&Slope(Q[top-1],Q[top])<=Slope(Q[top-1],i))--top;
Q[++top]=i;
}
for(int i=mid+1;i<=r;++i)f[i]=max(f[i],Calc(i,p[Bound(-SA[i]/SB[i])].id));
CDQ(mid+1,r);int t1=l,t2=mid+1,t=l;
for(;t1<=mid||t2<=r;)
if(t1<=mid&&(t2>r||(p[t1].x<p[t2].x)||(fabs(p[t1].x-p[t2].x)<1e-9&&p[t1].y<p[t2].y)))tmp[t++]=p[t1++];
else tmp[t++]=p[t2++];
for(int i=l;i<=r;++i)p[i]=tmp[i];
}
int main()
{
scanf("%d%lf",&n,&f[0]);
for(int i=1;i<=n;++i)scanf("%lf%lf%lf",&SA[i],&SB[i],&Rate[i]);
for(int i=1;i<=n;++i)f[i]=f[0];
CDQ(1,n);printf("%.3lf\n",f[n]);
return 0;
}
第二種方法是平衡樹，用平衡樹維護這個凸殼，每個節點按照\(x\)軸排序，每次插入一個新點的時候先檢查這個點是否存在於凸殼內，如果存在，把兩側不合法的點全部刪掉。
這個東西就是動態維護凸殼而已，似乎和單純的斜率優化關系不是很大了。
#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;
#define MAX 100100
const double inf=1e9;
const double eps=1e-7;
inline int read()
{
int x=0;bool t=false;char ch=getchar();
while((ch<'0'||ch>'9')&&ch!='-')ch=getchar();
if(ch=='-')t=true,ch=getchar();
while(ch<='9'&&ch>='0')x=x*10+ch-48,ch=getchar();
return t?-x:x;
}
int n;
double f[MAX],x[MAX],y[MAX];
double SA[MAX],SB[MAX],Rate[MAX];
double Calc(int i,int j)
{
double B=f[j]/(SA[j]*Rate[j]+SB[j]);
double A=Rate[j]*B;
return SA[i]*A+SB[i]*B;
}
double Slope(int i,int j)
{
if(fabs(x[i]-x[j])<eps)return -1e20;
return (y[i]-y[j])/(x[i]-x[j]);
}
struct SplayTree
{
#define ls (t[x].ch[0])
#define rs (t[x].ch[1])
int root,tot;
struct Node{int ch[2],ff;double k1,k2;}t[MAX];
void rotate(int x)
{
int y=t[x].ff,z=t[y].ff;
int k=t[y].ch[1]==x;
if(z)t[z].ch[t[z].ch[1]==y]=x;t[x].ff=z;
t[y].ch[k]=t[x].ch[k^1];if(t[x].ch[k^1])t[t[x].ch[k^1]].ff=y;
t[x].ch[k^1]=y;t[y].ff=x;
}
void Splay(int x,int goal)
{
while(t[x].ff!=goal)
{
int y=t[x].ff,z=t[y].ff;
if(z!=goal)
(t[y].ch[0]==x)^(t[z].ch[0]==y)?rotate(x):rotate(y);
rotate(x);
}
if(!goal)root=x;
}
void Insert(int u)
{
if(!root){root=u;return;}
int nw=root;
while(233)
{
int p=x[nw]+eps<=x[u];
if(!t[nw].ch[p]){t[nw].ch[p]=u;t[u].ff=nw;break;}
else nw=t[nw].ch[p];
}
}
int Pre(int x)
{
Splay(x,0);int u=ls,ret=u;
while(u)
if(Slope(u,x)<=t[u].k1)ret=u,u=t[u].ch[1];
else u=t[u].ch[0];
return ret;
}
int Suf(int x)
{
Splay(x,0);int u=rs,ret=u;
while(u)
if(Slope(u,x)>=t[u].k2)ret=u,u=t[u].ch[0];
else u=t[u].ch[1];
return ret;
}
void Update(int x)
{
Splay(x,0);
if(ls){int p=Pre(x);Splay(p,x);t[p].ch[1]=0;t[p].k2=t[x].k1=Slope(x,p);}
else t[x].k1=inf;
if(rs){int p=Suf(x);Splay(p,x);t[p].ch[0]=0;t[p].k1=t[x].k2=Slope(x,p);}
else t[x].k2=-inf;
if(t[x].k1<=t[x].k2)
t[ls].ch[1]=rs,t[ls].ff=0,t[rs].ff=ls,root=ls,t[ls].k2=t[rs].k1=Slope(ls,rs);
}
int Bound(double K)
{
int x=root;
while(233)
{
if(!x)return 0;
if(t[x].k1>=K&&t[x].k2<=K)return x;
if(t[x].k1<K)x=ls;else x=rs;
}
}
}T;
int main()
{
scanf("%d%lf",&n,&f[0]);
for(int i=1;i<=n;++i)scanf("%lf%lf%lf",&SA[i],&SB[i],&Rate[i]);
for(int i=1;i<=n;++i)
{
int j=T.Bound(-SA[i]/SB[i]);
f[i]=max(f[i-1],Calc(i,j));
y[i]=f[i]/(SA[i]*Rate[i]+SB[i]);
x[i]=y[i]*Rate[i];
T.Insert(i);T.Update(i);
}
printf("%.3lf\n",f[n]);
return 0;
}
Ending
似乎就這么多了QwQ
總結一下，其實推式子都是大同小異的。
主要的區別就在於坐標是否具有單調性。
如果坐標和詢問的斜率都單調，直接使用單調隊列就好了。
如果坐標單調而詢問的斜率不單調，那么就需要維護單調棧，在其中二分斜率。
如果坐標不單調，那么用\(CDQ\)分治或者\(Splay\)維護凸殼(二進制分組多一個\(log\))。
upd：還可以用李超線段樹來維護。
除了這幾道裸的\(dp\)題之外，顯然還是有很多套路的，比如說套上線段樹分治讓你維護凸殼之類的，這里都懶得寫了。
就這樣啦QaQ。