先来看一道题(HDU3507):
题意:给出N个单词,每个单词有个非负权值Ci,将k个单词排在一行的费用为(∑Ci)^2+M.求最优方案,使得总费用最小.
我们很容易得到一个O(N^2)的算法:
s[i]表示前i个单词的权值和.
先写个东西在这:所有元素非负的数组的前缀和值随下标增加单调递增.后面会用到.
f[i]表示将前i个单词排版完毕后的最优值,f[i]=min{f[j]+(s[i]-s[j])^2+M}.
但题目中N的范围是500000.这个算法明显不行.考虑如何优化.
我们固定i,考虑它的两个一般决策点j,k(j<k).
记g[pos]=f[pos]+(s[i]-s[pos])^2+M,即i从pos转移的代价.
如果决策点k优于j,那么就有g[k]<g[j].展开来:
f[k]+(s[i]-s[k])^2+M<f[j]+(s[i]-s[j])^2+M,化简得
f[k]-f[j]+s[k]^2-s[j]^2<2*s[i]*(s[k]-s[j])
注意到s[k]>s[j],我们在不等式两边除以(s[k]-s[j]).
不等式化为(f[k]-f[j]+s[k]^2-s[j]^2)/(s[k]-s[j])<2*s[i].
方便起见,我们将左边分式的分子分母同时变号.
(f[j]-f[k]+s[j]^2-s[k]^2)/(s[j]-s[k])<2*s[i].
可以看到不等式左边与i无关,右边只与i有关.(而且左边像一个两点间的斜率式).
记slope[j,k]=(f[j]-f[k]+s[j]^2-s[k]^2)/(s[j]-s[k]).
好了,现在我们有一个结论.
△对于i的两个决策点j,k(j<k),决策k优于决策j就等价于slope[j,k]<2*s[i].
换句话说,如果slope[j,k]<2*s[i],那么决策k优于j,反之决策j不比k差.
其实我们还可以知道,决策点k永远会比决策点j优.因为对于以后的i',s[i']>s[i]>slope[j,k].
因此这里的优劣应该是全局的,而不只限于i.
我们再来考虑三个点j,k,l(j<k<l)之间的优劣关系.
还是通过斜率:
如果slope[j,k]>slope[k,l],我们看看能得到什么.
1.若slope[k,l]<2*s[i].那么由之前的结论(△),l比k优.
2.若slope[k,l]>2*s[i],则slope[j,k]>2*s[i],那么由之前的结论(△),决策j不比k差.
综上,如果slope[j,k]>slope[k,l],k是可以淘汰掉的.
我们又得到一个结论.
△对于三个决策点j,k,l(j<k<l),如果slope[j,k]>slope[k,l],那么k永远不会成为某个点的最优决策.
现在我们有了这两个结论,怎样来优化呢?
我们可以将决策放到一个队列中,利用以上两个结论剔除无用决策点,达到快速转移的目的.
记队列的头指针为h,尾指针为t.
对于队列的头部,如果slope[q[h],q[h+1]]<2*s[i],那么,q[h]一定可以去掉了.h=h+1.
事实上经过这样的调整后,q[h]就是i的最有决策,直接取来更新就是了.
更新出f[i]后,将f[i]从尾部加入队列,并用i去剔除无用决策.
对于队列的尾部,如果slope[q[t-1],q[t]]>slope[q[t],i],那么q[t]可以去掉.t=t-1.
(这里我是按照我自己写程序的习惯写的,先用i去更新队尾,再加入i.还可以有不同的写法)
顺便说一句,这样维护以后,队列中的"点"形成一个上凸包(联想上面说的斜率).
程序大致过程
for i=1 to n do
begin
当队列不为空时,更新队头;
取当前队头更新f[i];
用i去更新队尾;
将i加入队尾.
end;
可以看到外层循环是O(N)的,内层里每个元素进出队列仅一次,所以总效率为O(N).
code:(HDU3507)
const oo=1e100;
maxn=500001;
var s,f:Array[0..maxn] of int64;
q:array[0..maxn] of longint;
n,m,i,h,t,c:longint;
function slope(j,k:longint):extended;
begin
if s[j]=s[k] then
begin
if f[j]>=f[k] then slope:=-oo
else slope:=oo;
exit;
end;
slope:=(f[j]-f[k]+s[j]*s[j]-s[k]*s[k])/(s[j]-s[k]);
end;
begin
while not seekeof do
begin
readln(n,m);
for i:=1 to n do
begin
readln(c);
s[i]:=s[i-1]+c;
end;
h:=0; t:=0;
for i:=1 to n do
begin
while (h<t)and(slope(q[h],q[h+1])<2*s[i]) do inc(h);
f[i]:=f[q[h]]+(s[i]-s[q[h]])*(s[i]-s[q[h]])+M;
while (h<t)and(slope(q[t-1],q[t])>slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
writeln(f[n]);
end;
end.
斜率优化相关题目推荐:
HNOI2008 玩具装箱
APIO2010 特别行动队
CEOI2004 锯木场选址
POJ 1180
POJ 3709
HNOI2008 玩具装箱
const oo=1<<60;
var f,s,q:array[0..50001] of int64;
n,l,i,j,c,h,t:longint;
function min(a,b:int64):int64;
begin
if a>b then exit(b); exit(a);
end;
function slope(j,k:longint):extended;
begin
exit((f[j]-f[k]-sqr(s[k]+k)+sqr(s[j]+j))/(s[j]+j-s[k]-k))
end;
begin
readln(n,l);
for i:=1 to n do
begin
readln(c);
s[i]:=s[i-1]+c;
end;
h:=0; t:=0;
for i:=1 to n do
begin
while (h<t)and(slope(q[h],q[h+1])<(i+s[i]-l-1)*2) do inc(h);
f[i]:=f[q[h]]+sqr(i-q[h]-1+s[i]-s[q[h]]-l);
while (h<t)and(slope(q[t-1],q[t])>slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
writeln(f[n]);
end.
APIO2010 特别行动队
const oo=2000000000;
var q,num:array[0..1000001] of longint;
f,s:array[0..1000001] of extended;
n,i,h,t:longint;
w,a,b,c:extended;
function slope(j,k:longint):extended;
begin
exit((f[j]-f[k]+a*s[j]*s[j]-a*s[k]*s[k]-b*s[j]+b*s[k])/(s[j]-s[k]))
end;
begin
readln(n);
readln(a,b,c);
for i:=1 to n do
begin
read(num[i]);
s[i]:=s[i-1]+num[i];
end;
h:=0; t:=0;
for i:=1 to n do
begin
while (h<t)and(slope(q[h],q[h+1])>=2*a*s[i]) do inc(h);
w:=s[i]-s[q[h]];
f[i]:=f[q[h]]+a*w*w+b*w+c;
while (h<t)and(slope(q[t-1],q[t])<=slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
writeln(f[n]:0:0);
end.
CEOI2004 锯木场选址
const maxn=20001;
var f,w,d,sw,sd,c,q:array[0..maxn] of longint;
n,i,j,tmp,ans,h,t:longint;
function cost(l,r:longint):longint;
begin
cost:=c[r]-c[l-1]-sw[l-1]*(sd[r]-sd[l-1]);
end;
function slope(j,k:longint):extended;
begin
slope:=(sw[k]*sd[k]-sw[j]*sd[j])/(sw[k]-sw[j]);
end;
begin
readln(n);
for i:=1 to n do readln(w[i],d[i]);
for i:=1 to n+1 do
begin
sw[i]:=sw[i-1]+w[i];
sd[i]:=sd[i-1]+d[i-1];
end;
for i:=1 to n+1 do c[i]:=c[i-1]+sw[i-1]*d[i-1];
h:=1; t:=1; q[1]:=1;
for i:=2 to n do
begin
while (h<t)and(slope(q[h],q[h+1])<sd[i]) do inc(h);
f[i]:=c[q[h]]+cost(q[h]+1,i)+cost(i+1,n+1);
while (h<t)and(slope(q[t-1],q[t])>slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
ans:=maxlongint;
for i:=2 to n do
if ans>f[i] then ans:=f[i];
writeln(Ans);
end.
POJ1180
var f,tx,fx,st,sf,q:array[0..100001] of longint;
n,s,i,h,t:longint;
function slope(i,j:longint):extended;
begin
exit((f[i]-f[j])/(st[i]-st[j]));
end;
begin
readln(n);
readln(s);
for i:=1 to n do readln(tx[i],fx[i]);
for i:=n downto 1 do
begin
st[i]:=st[i+1]+tx[i];
sf[i]:=sf[i+1]+fx[i];
end;
h:=0; t:=0;
for i:=n downto 1 do
begin
while (h<t)and(slope(q[h],q[h+1])<sf[i]) do inc(h);
f[i]:=f[q[h]]+(s+st[i]-st[q[h]])*sf[i];
while (h<t)and(slope(q[t-1],q[t])>slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
writeln(f[1]);
end.
POJ3709
var q:array[0..500001] of longint;
f,s,num:array[0..500001] of extended;
dnum,d,i,j,m,n,h,t:longint;
function slope(j,k:longint):extended;
var x,y:extended;
begin
y:=f[j]-f[k]-s[j]+s[k]+num[j+1]*j-num[k+1]*k;
x:=num[j+1]-num[k+1];
if x=0 then
if y>=0 then exit(k+m)
else exit(n+1);
slope:=y/x;
if slope<k+m then exit(k+m);
end;
begin
readln(dnum);
for d:=1 to dnum do
begin
readln(n,m);
fillchar(s,sizeof(s),0);
fillchar(f,sizeof(f),0);
for i:=1 to n do
begin
read(num[i]);
s[i]:=s[i-1]+num[i];
end;
readln;
h:=0; t:=0;
for i:=m to n do
begin
while (h<t)and(slope(q[h],q[h+1])<=i) do inc(h);
f[i]:=f[q[h]]+s[i]-s[q[h]]-num[q[h]+1]*(i-q[h]);
while (h<t)and(slope(q[t-1],q[t])>=slope(q[t],i)) do dec(t);
inc(t); q[t]:=i;
end;
writeln(f[n]:0:0);
end;
end.
更多题目:http://www.notonlysuccess.com/?p=740
参考资料:
2004国家集训队论文 周源 《浅谈数形结合思想在信息学竞赛中的应用》
2009JSOI集训队论文 《用单调性优化动态规划》
《动态规划的斜率优化》
