题解 P4027 【[NOI2007]货币兑换】
litble
2018-01-03 20:30:34
[蒟蒻的博客](http://blog.csdn.net/litble/article/details/78965163 )
# 斜率优化
首先，由于如果在i天买在j天卖有利可图，那么最优方法就是在i天花完钱在j天卖完。我们令$f_i$为第i天可以得到的最多钱数，然后可以先列方程求出花完钱在第j天得到的两种金券数$x_j=\frac{f_jR_j}{a_jR_j+b_j}$和$y_j=\frac{f_j}{a_jR_j+b_j}$，然后得到状态转移方程：$f_i=x_ja_i+y_jb_i$
将方程稍微变形成直线斜截式方程：$y_j=-\frac{a_ix_j}{b_i}+\frac{f_i}{b_i}$，可以知道，对于平面上众多的点$(x_j,y_j)$，我们每次用一条斜率为$-\frac{a_i}{b_i}$的直线去切它们中的一个，可以求得最大截距的点即为最优决策。
如何求得最大截距呢？以下我们将斜率为$-\frac{a_i}{b_i}$的直线称为当前直线。
现在我们维护点之间的一个凸包，如果对于点j，其左边的线斜率小于当前直线，那么显然把当前直线移到左边的点截距更大（建议自己画图理解）。如果其右边的线斜率大于当前直线，那么将其右移更优。
于是就有两种维护方式：
# splay维护
我们需要优化，那么就需要排除不可能状态。可以发现，在点集中，如果一个点在凸包内部，那么这个点就是不可能是最优决策的点，因为其上一定有一个点，那么把当前直线上移可以获得更大截距。所以我们可以排除这些决策．
维护lk:点x凸包左边那条线的斜率，rk:右边那条线的。
splay按照x排序，寻找最优决策的方法见上，而添加新点的方法：
首先将新点x旋转到根。以寻找其左边最后一个可以与其构成凸包的点为例。对于当前点t，如果t左边斜率大于了直线tx的斜率，那么如果t的右边还有不在凸包内的节点，用tx的连线就不能构成凸包，所以应该继续往右找点。否则，继续往左找点。
然后，删掉在凸包里的点。
找右边第一个可以与其构成凸包的点同理。
最后，如果发现这个点本来就在旧的凸包里面，直接将其删除。（即$lk_x$大于$rk_x$）
```cpp
#include<bits/stdc++.h>
using namespace std;
#define db double
#define eps 1e-9//一定要注意精度
#define inf 1e9
const int N=100005;
int n,rt,sc,f[N],son[N][2];
db dp[N],A[N],B[N],R[N],lk[N],rk[N],X[N],Y[N];//lk：凸包点x左线斜率，rk：右线斜率
int is(int x) {return son[f[x]][1]==x;}
void spin(int x,int &mb) {
int fa=f[x],g=f[fa],t=is(x);
if(fa==mb) mb=x;
else son[g][is(fa)]=x;
f[fa]=x,f[x]=g,f[son[x][t^1]]=fa;
son[fa][t]=son[x][t^1],son[x][t^1]=fa;
}
void splay(int x,int &mb) {
while(x!=mb) {
if(f[x]!=mb) {
if(is(x)^is(f[x])) spin(x,mb);
else spin(f[x],mb);
}
spin(x,mb);
}
}
int find(int x,db num) {//寻找最优解
if(!x) return 0;
if(lk[x]+eps>=num&&rk[x]<=num+eps) return x;
else if(lk[x]<num+eps) return find(son[x][0],num);
else return find(son[x][1],num);
}
db getk(int a,int b) {//获得斜率
if(X[a]-X[b]<eps&&X[a]-X[b]>-eps) return -inf;
return (Y[b]-Y[a])/(X[b]-X[a]);
}
int pre(int x) {//寻找左边最后一个与x可以构成凸包的点
int y=son[x][0],re=y;
while(y) {
if(lk[y]+eps>=getk(y,x)) re=y,y=son[y][1];
else y=son[y][0];
}
return re;
}
int nxt(int x) {//寻找右边第一个与x可以构成凸包的点
int y=son[x][1],re=y;
while(y) {
if(rk[y]<=getk(x,y)+eps) re=y,y=son[y][0];
else y=son[y][1];
}
return re;
}
void newjd(int x) {
splay(x,rt);
if(son[x][0]) {
int kl=pre(x);
splay(kl,son[x][0]),son[kl][1]=0;
lk[x]=rk[kl]=getk(kl,x);
}
else lk[x]=inf;//请勿往左
if(son[x][1]) {
int kl=nxt(x);
splay(kl,son[x][1]),son[kl][0]=0;
rk[x]=lk[kl]=getk(x,kl);
}
else rk[x]=-inf;//请勿往右
if(lk[x]<=rk[x]+eps) {//在原凸包内部,直接删除该点
rt=son[x][0],son[rt][1]=son[x][1],f[son[x][1]]=rt,f[rt]=0;
lk[rt]=rk[son[rt][1]]=getk(rt,son[rt][1]);
}
}
void ins(int &x,int las,int bh) {
if(!x) {x=bh,f[x]=las;return;}
if(X[bh]<=X[x]+eps) ins(son[x][0],x,bh);
else ins(son[x][1],x,bh);
}
int main()
{
scanf("%d%lf",&n,&dp[0]);
for(int i=1;i<=n;++i) {
scanf("%lf%lf%lf",&A[i],&B[i],&R[i]);
int j=find(rt,-A[i]/B[i]);
dp[i]=max(dp[i-1],X[j]*A[i]+Y[j]*B[i]);
Y[i]=dp[i]/(A[i]*R[i]+B[i]),X[i]=Y[i]*R[i];
ins(rt,0,i),newjd(i);
}
printf("%.3lf\n",dp[n]);
return 0;
}
```
# CDQ分治维护
[关于CDQ分治](http://blog.csdn.net/litble/article/details/79016999)
我们可以对决策的时间进行二分，对于左半边区间，首先求出它们的dp值，右半边区间维持$k_i=-\frac{a_i}{b_i}$的有序，左半边区间维持$x_i$的有序。
我们可以开一个栈来维护左半边区间的斜率单调递减的凸包。
对于右半边区间，由于k值排了序，所以可以O(n)查询。查询方法就是关于斜率优化的分析的倒数第二段的内容。
```cpp
#include<bits/stdc++.h>
using namespace std;
#define db double
#define inf 1e9
#define eps 1e-9
const int N=100005;
int n,s[N];db dp[N];
struct node{db k,x,y,a,b,r;int id;}Q[N],kl[N];
db getk(int i,int j) {
if(fabs(Q[i].x-Q[j].x)<=eps) return inf;
return (Q[j].y-Q[i].y)/(Q[j].x-Q[i].x);
}
void merge(int l,int r,int mid) {//归并排序
int t1=l,t2=mid+1;
for(int i=l;i<=r;++i)
if(t1<=mid&&(t2>r||Q[t1].x<Q[t2].x+eps)) kl[i]=Q[t1],++t1;
else kl[i]=Q[t2],++t2;
for(int i=l;i<=r;++i) Q[i]=kl[i];
}
void cdq(int l,int r) {
if(l==r) {//那么在l之前的所有询问都已经处理完毕，可以更新l的答案了
dp[l]=max(dp[l],dp[l-1]);
Q[l].y=dp[l]/(Q[l].a*Q[l].r+Q[l].b),Q[l].x=Q[l].y*Q[l].r;
return;
}
int mid=(l+r)>>1,t1=l-1,t2=mid,top=0;
for(int i=l;i<=r;++i)//把前mid个询问放在左边，后mid个放在右边
if(Q[i].id<=mid) kl[++t1]=Q[i];
else kl[++t2]=Q[i];
for(int i=l;i<=r;++i) Q[i]=kl[i];
cdq(l,mid);//递归处理左边
for(int i=l;i<=mid;++i) {//维护斜率递减的凸包
while(top>=2&&getk(s[top],i)+eps>getk(s[top-1],s[top])) --top;
s[++top]=i;
}
for(int i=mid+1;i<=r;++i) {//处理右边的询问
while(top>=2&&getk(s[top-1],s[top])<=Q[i].k+eps) --top;
int j=s[top];
dp[Q[i].id]=max(dp[Q[i].id],Q[j].x*Q[i].a+Q[j].y*Q[i].b);
}
cdq(mid+1,r),merge(l,r,mid);//递归处理右边后，按照x值为关键字归并排序
}
int cmp1(node t1,node t2) {return t1.k<t2.k;}
int main()
{
scanf("%d%lf",&n,&dp[0]);
for(int i=1;i<=n;++i) {
scanf("%lf%lf%lf",&Q[i].a,&Q[i].b,&Q[i].r);
Q[i].k=-Q[i].a/Q[i].b,Q[i].id=i;
}
sort(Q+1,Q+1+n,cmp1),cdq(1,n);
printf("%.3lf\n",dp[n]);
return 0;
}
```