关于斜率优化,我就是一个傻子啊,真的一直没弄懂……
状态和方程还是很好出来的啊:
(f[i]=min(f[j]+(s[i]-s[j]+i-j-L-1)^2))其中(s[i])表示前缀和,(f[i])表示前(i)个处理后的最小值。
但是我们发现,这个东西要转移的话是个(O(n^2))而(N<=50000),显然转移不了。
接下来就是斜率优化的天下了!
我们要找的是在(i)之前的一个(j)使得(f[i])最小,考虑怎样可以快速的找到。
令(a[i]=s[i]+i,b[i]=s[i]+i-L-1)。
我们将原方程变形,就可以得到这样的形式:(先不要管为啥好吧)
(2*a[i]*b[j]+f[i]-a[i]^2=f[j]+b[j]^2)。先明确一点,当我们在转移(i)时(a[i])是确定的。
对应一下(kx+b=y),我们令(2*a[i])为(k),(b[j])为(x),(f[j]+b[j]^2)为(y),那么,我们所求为直线在(y)上的最小截距。
当我们在转移(i)时,前面的(1~i-1)可表示为一堆点((P_j(b[j],f[j]+b[j]^2)))。
转移即为,找一条过(P_j)的斜率为(2*a[i])的直线,使得其在(y)轴上的截距最小。
不妨拿出三个点,其构成一个向上的凸包(如右上的图)我们发现,当我们平移直线是,(B)一定不为最优决策。
所以我们可以舍掉(B)点,即维护一个向下的凸包(如下面的图)这是我们发现,此时的(A,B,C)三点都有可能成为最优决策。
此时考虑用一个数据结构维护所有可以构成一个向下的凸包的点。
//回来看下之前我们说的斜率(2*a[i]),显然它具有单调性(如果没有单调性就二分)
在前(1~i-1)中,不妨设有两个状态:(j,k)且(j)比(k)更优,则有:
(f[j]+(s[i]-s[j]+i-j-L-1)^2<f[k]+(s[i]-s[k]+i-k-L-1)^2)
经过一系列毫不人道的化简(风骨傲天很懒所以他没把将过程放上来)可以得到:
(2*a[i]>frac{f[j]-f[k]+b[j]^2-b[k]^2}{b[j]-b[k]})(好丑的式子……)
也就是说只要满足这个式子,就有(j)比(k)更优。
同样对应斜率,就有:(P_j)和(P_k)的斜率小于(2*a[i])时,就有(j)比(k)更优。
又因为我们维护的是一个向下的凸包,所以我们就只用考虑相邻的两个点即可。
那么我们就可以用一个单调队列来维护,队列中相邻点间连线的斜率递增。
这张图应该说很好的反应了转移和维护的过程,上面的两张就是转移,下面的两张是维护。
转移:由图中我们可以发现,我们要求的(j)为斜率第一个大于(2*a[i])的点,因此舍掉(A,B)
维护:因为单桥队列中斜率的单调性,删掉(E),由转移后的(i)得出的点(F(b[i],f[i]+b[i]^2))加入队尾。
而实际上,我们的斜率优化有一定的公式:
当方程形如:(f[i]=min(f[j]+S(i,j))+k)((k)为常数)时,可以使用斜率优化。我们的斜率为在当次转移中的一个不变量,维护的是一堆点。不过通常用上面讲到的“假设两个状态法”来求斜率。
所以斜率优化的思考过程应该是和上面讲到的相反,先考虑斜率再变形方程并得出(x,y)。
啊~终于打完了,累成狗
又是快乐的代码时间:
#include<bits/stdc++.h>
using namespace std;
inline int read()
{
int f=1,w=0;char x=0;
while(x<'0'||x>'9') {if(x=='-') f=-1; x=getchar();}
while(x!=EOF&&x>='0'&&x<='9') {w=(w<<3)+(w<<1)+(x^48);x=getchar();}
return w*f;
}
const int N=50010;
int n,q[N],l,r,L;
double s[N],f[N];
inline double A(int i) {return s[i]+i;}
inline double B(int i) {return s[i]+i+L+1;}
inline double S(double x) {return x*x;}
inline double K(int i,int j) {return (f[i]-f[j]+S(B(i))-S(B(j)))/(B(i)-B(j));}
int main(){
#ifndef ONLINE_JUDGE
freopen("A.in","r",stdin);
#endif
n=read();L=read();l=r=1;
for(int i=1,c;i<=n;i++) c=read(),s[i]=c*1.0+s[i-1];
for(int i=1;i<=n;i++)
{
while(l<r&&2*(s[i]+i)>K(q[l],q[l+1])) l++;
f[i]=f[q[l]]+S(A(i)-B(q[l]));
while(l<r&&K(q[r-1],q[r])>K(i,q[r-1])) r--;
q[++r]=i;
}
printf("%lld",(long long)f[n]);
}
这一题也是一个经典的斜率优化,稍微有点不同。
状态:(f[i]=min(f[j]+a*(s[i]-s[j])^2+b*(s[i]-s[j])+c))
化简:(2*a*s[i]*s[j]+f[i]-a*s[i]^2-b*s[i]-c=f[j]+a*s[j]^2-b*s[j])
因为这一题的(a<0),所以我们维护一个向上的凸包即可。(但是(a<0)时仍有单调性)
直接上代码:
#include <cstdio>
using namespace std;
#define F(x) ((x)*(x))
inline int read()
{
int f=1,w=0;char x=0;
while(x<'0'||x>'9') {if(x=='-') f=-1; x=getchar();}
while(x!=EOF&&x>='0'&&x<='9') {w=(w<<3)+(w<<1)+(x^48);x=getchar();}
return w*f;
}
const int N=1000010;
#define int long long
int n,l,r,f[N],Q[N],a,b,c,s[N];
inline double Work(int x,int y)
{
return 1.*(f[x]-f[y]+(F(s[x])-F(s[y]))*a)/(s[x]-s[y])-b;
}
main(){
n=read(),a=read(),b=read(),c=read();
for(int i=1;i<=n;i++) s[i]=read(),s[i]+=s[i-1];
for(int i=1;i<=n;i++)
{
while(l<r&&Work(Q[l+1],Q[l])>s[i]*2*a) l++;
//用于从队列中选出最值
f[i]=f[Q[l]]+a*F(s[i]-s[Q[l]])+b*(s[i]-s[Q[l]])+c;
while(l<r&&Work(Q[r],Q[r-1])<Work(i,Q[r])) r--;
//维护队列单调性
Q[++r]=i;
}
printf("%lld",f[n]);
}
关于不满足单调性时的特殊情况,下面来一个例题。
(BZOJ2726)(慎重声明,这一题和(Luogu)上的不一样)
数据范围:
([1, 4] 0<N<=1000,0<=S<=2^8,0<=Ti<=2^8,0<=Fi<=2^8)
([5, 12] 0<N<=300000,0<=S<=2^8,0<=Ti<=2^8,0<=Fi<=2^8)
([13, 20] 0<N<=100000,0<=S<=2^8,-(2^8)<=Ti<=2^8,0<=Fi<=2^8)
以下的(F[i]),(T[i])表示相应数组的前缀和
因为我们在转移中需要前面分成的批数,所以有一个极为直接的状态:(f[i][j])表示前(i)个,分为(j)组的答案。
但实际上空间上完全不行(你(100000)怎么开二维……)
考虑一维的状态,实际上这里用上里一个叫“费用提前计算”的思想。(先把方程写出来吧)
(f[i]=min(f[j]+T[i]*(F[i]-F[j])+S*(F[n]-F[j])))
考虑这个(S)会对什么产生影响:显然是(j+1~i)的任务产生影响,因此我们将它提出来计算。
然鹅这个怎么看都不是一个好的方法。(你(100000)怎么跑这个……)
用斜率优化考虑转移,式子可以变形为:
(f[j]=(S+T[i])*F[j]+f[i]-T[i]*F[i]-S*F[n])
设(f[j])为(y),(F[j])为(x),但我们在转移时就傻眼了,因为这一题的特殊性出题人的毒瘤性,(T[i])可能为负
即我们的每次转移时那个固定的斜率不单调,我们不能将队头的点删掉!
但是至少我们的(F[i])有单调性,加入的点有单调性。
这时我们可以不删除队列中的点,利用二分在队列中找最适合转移的点,再进行转移。
代码应该会感觉有些奇怪,主要是智障太懒了,直接在弱化版上魔改了……
#include<bits/stdc++.h>
using namespace std;
#define int long long
inline int read()
{
int f=1,w=0;char x=0;
while(x<'0'||x>'9') {if(x=='-') f=-1; x=getchar();}
while(x!=EOF&&x>='0'&&x<='9') {w=(w<<3)+(w<<1)+(x^48);x=getchar();}
return w*f;
}
const int N=1000010;
int n,s,T[N],F[N],q[N],top=1,f[N];
inline bool check(int j,int i)
{
if(j<top) return f[q[j+1]]-f[q[j]]<=(T[i]+s)*(F[q[j+1]]-F[q[j]]);
else return 0;
}
inline bool Check(int j,int i)
{
int a=(f[i]-f[q[j]])*(F[q[j]]-F[q[j-1]]);
int b=(f[q[j]]-f[q[j-1]])*(F[i]-F[q[j]]);
return a<=b;
}
main(){
#ifndef ONLINE_JUDGE
freopen("Text1.in","r",stdin);
#endif
n=read(),s=read();
for(int i=1;i<=n;i++)
T[i]=read(),F[i]=read(),T[i]+=T[i-1],F[i]+=F[i-1];
for(int i=1;i<=n;i++)
{
int L=1,R=top;
while(L<R)
{
int mid=(L+R)>>1;
if(L==R) break ;
if(check(mid,i)) L=mid+1;
else R=mid;
}
int j=q[L];
f[i]=f[j]+T[i]*F[i]+s*F[n]-F[j]*(s+T[i]);
while(top>=2&&Check(top,i)) q[top--]=0; q[++top]=i;
}
printf("%lld",f[n]);
}
所以我们做个总结,斜率优化(DP)适用于状态转移方程为:
(f[i]=min(f[j]+S(i,j))+k)((k)为常数)且(S(i,j))计算时有((i)*(j))的部分。
同时我们设的(x,k)必定要有单调性,如果(k)没有,就不删点二分,如果都没有,就用(CDQ)(相当于是动态插点,动态查找)当然没人拦你用平衡树……
据说……这道题被许多人嘲讽为水题,可我不这么想啊(果真是因为我太弱了吗……
但考场上就真的要……(WOC)你给我解释一下(O_{(mt)})都能过是什么情况啊!!
你确定你真的不是用脚在造数据?!
正解(你(™)别给老子想什么暴力卡常):斜率优化(DP),推方程完全不难,设(f[i])为最后乘编号为(i)的车的最小烦躁值,转移方程为:
决策点为((q_i,f[j]+A*q_j^2-B*q_j)),当(j)比(k)更优时,满足:
但是,我们转移时要满足(p_i>=q_j,y_j=x_i),所以我们不能像以前一样维护一个凸包,从所有的决策中转移。
仔细思考一下,我们的决策来自部分满足条件的前面已经做出的决策,不妨我们对决策按(i)到达的位置分个组。
即在每个节点处维护一个凸包(这样一定满足单调性,不解释了),凸包中的点为(f[j])((j)为目的地为该节点的列车编号),这样我们转移时就可以满足空间限制了。
再考虑时间限制怎么做,我们可以枚举时间,再开一个等待队列,存在(q[i])时到的列车(i),但他们不能被利用,因为还没枚举到他们的到达时间,然后枚举到时间(t)时,将等待队列中所有到达时间为(t)的列车加进相应的凸包中(并维护凸包的单调性),说明他们可以被利用来转移。
最后在每次转移后将该次转移加入等待队列中,判断是否到达(n),如果到达,就更新答案(记得加上(q_i))
给泥萌看我丑陋的代码:
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define S(x) ((x)*(x))
inline int read()
{
int f=1,w=0;char x=0;
while(x<'0'||x>'9') {if(x=='-') f=-1; x=getchar();}
while(x!=EOF&&x>='0'&&x<='9') {w=(w<<3)+(w<<1)+(x^48);x=getchar();}
return w*f;
}
const int N=200010,M=1001;
queue<int> res[M];
int n,m,A,B,C,MaxT,ans=1e18;
vector<int> Tbg[M],Q[N];
int q[N],p[N],x[N],y[N],head[N],f[N];
inline double K(int j,int k)
{
return (double)(f[j]-f[k]+A*(S(q[j])-S(q[k]))-B*(q[j]-q[k]))/(double)(q[j]-q[k]);
}
main(){
#ifndef ONLINE_JUDGE
//freopen("A.in","r",stdin);//Ans=94;
freopen("B.in","r",stdin);//Ans=34;
#endif
n=read(),m=read(),A=read(),B=read(),C=read();
for(int i=1;i<=m;i++)
{
x[i]=read(),y[i]=read(),p[i]=read();
q[i]=read(),Tbg[p[i]].push_back(i);
MaxT=max(MaxT,q[i]);
}
Q[1].push_back(0);
for(int t=0;t<=MaxT;t++)
{
while(!res[t].empty())
{
int pos=y[res[t].front()];
while(Q[pos].size()-head[pos]>=2)
{
int len=Q[pos].size();
if(K(Q[pos][len-1],Q[pos][len-2])<K(Q[pos][len-2],res[t].front())) break;
Q[pos].pop_back();
}
Q[pos].push_back(res[t].front()),res[t].pop();
}
for(int i=0;i<(int)Tbg[t].size();i++)
if((int)Q[x[Tbg[t][i]]].size()>head[x[Tbg[t][i]]])
{
int id=Tbg[t][i],pos=x[id];
while((int)Q[pos].size()-head[pos]>=2)
{
if(K(Q[pos][head[pos]],Q[pos][head[pos]+1])>2.0*A*p[id]) break ;
head[pos]++;
}
int j=Q[pos][head[pos]];
f[id]=f[j]+A*S(p[id]-q[j])+B*(p[id]-q[j])+C;
res[q[id]].push(id);if(y[id]==n) ans=min(ans,f[id]+q[id]);
}
}
printf("%lld",ans);
}