NOI 2007 货币兑换Cash （bzoj 1492）

Description

小Y最近在一家金券交易所工作。该金券交易所只发行交易两种金券：A纪念券（以下简称A券）和 B纪念券（以下

简称B券）。每个持有金券的顾客都有一个自己的帐户。金券的数目可以是一个实数。每天随着市场的起伏波动，

两种金券都有自己当时的价值，即每一单位金券当天可以兑换的人民币数目。我们记录第 K 天中 A券 和 B券 的

价值分别为 AK 和 BK（元/单位金券）。为了方便顾客，金券交易所提供了一种非常方便的交易方式：比例交易法

。比例交易法分为两个方面：（a）卖出金券：顾客提供一个 [0,100] 内的实数 OP 作为卖出比例，其意义为：将

 OP% 的 A券和 OP% 的 B券 以当时的价值兑换为人民币；（b）买入金券：顾客支付 IP 元人民币，交易所将会兑

换给用户总价值为 IP 的金券，并且，满足提供给顾客的A券和B券的比例在第 K 天恰好为 RateK；例如，假定接

下来 3 天内的 Ak、Bk、RateK 的变化分别为：

假定在第一天时，用户手中有 100元 人民币但是没有任何金券。用户可以执行以下的操作：

注意到，同一天内可以进行多次操作。小Y是一个很有经济头脑的员工，通过较长时间的运作和行情测算，他已经

知道了未来N天内的A券和B券的价值以及Rate。他还希望能够计算出来，如果开始时拥有S元钱，那么N天后最多能

够获得多少元钱。

Input

输入第一行两个正整数N、S，分别表示小Y能预知的天数以及初始时拥有的钱数。接下来N行，第K行三个实数AK、B

K、RateK，意义如题目中所述。对于100%的测试数据，满足：0<AK≤10；0<BK≤10；0<RateK≤100；MaxProfit≤1

0^9。

【提示】

1.输入文件可能很大，请采用快速的读入方式。

2.必然存在一种最优的买卖方案满足：

每次买进操作使用完所有的人民币；

每次卖出操作卖出所有的金券。

 

Output

只有一个实数MaxProfit，表示第N天的操作结束时能够获得的最大的金钱数目。答案保留3位小数。

Sample Input

3 100
1 1 1
1 2 2
2 2 3

Sample Output

225.000

HINT

---

题目大意 每天可以根据一定比例购买或者卖出纪念劵。问最多在$n$天后可以得到多少钱。

 　　由题意易得，一定存在一种最优的方案满足，能买时就全买，能卖时就全卖。因为如果要赚就要尽可能地多赚，会亏就一点都不去碰。

 　　设$f[i]$表示第$i$天后把手上的纪念劵都卖掉后可以得到最多的钱，于是轻松地列出了dp方程：

$f[i] = maxleft { frac{f[j]left ( r_{j}A_{i} + B_{i} ight )}{r_{j}A_{j} + B_{j}}, f[i - 1] ight }$

　　暂时先不管$max$，当成等式，两边同除以$B_{i}$。

$frac{f[i] }{B_{i}}= frac{f[j]left ( r_{j}frac{A_{i}}{B_{i}} + 1 ight )}{r_{j}A_{j} + B_{j}}$

　　然后移移项什么的。。

$frac{f[i] }{B_{i}}= frac{f[j]r_{j}}{r_{j}A_{j} + B_{j}}cdot frac{A_{i}}{B_{i}}+frac{ f[j]}{r_{j}A_{j} + B_{j}}$

$frac{ f[j]}{r_{j}A_{j} + B_{j}}=- frac{f[j]r_{j}}{r_{j}A_{j} + B_{j}}cdot frac{A_{i}}{B_{i}}+frac{f[i] }{B_{i}}$

　　（第二步是把和$j$相关的扔到等号左边去，当做$y$，把形如$Mleft(i ight)Nleft(j ight )$，扔到左边，把其中的$Mleft(i ight)$看作常数项，把$Nleft(j ight )$看作$x$）

　　所以有：

$x_{i} = frac{r_{i}f[i]}{r_{i}A_{i} + B{i}},y_{i}=frac{f[i]}{r_{i}A_{i} + B{i}}$

　　然后可得：

$y_{j} = -frac{A_{i}}{B_{i}}x_{j} + frac{f[i]}{B_{i}}$

　　因为要最大化$f[i]$，所以应该最大化截距。所以维护上凸壳。

　　同时方程也可以写成：

$f[i] =maxleft { A_{i}x_{j}+B_{i}y_{j}, f[i - 1] ight }$

　　考虑如何来求最大的截距。

　　因为这里插入点的$x$坐标不单调，询问的斜率也不单调，所以不能开单调队列暴力移指针了。

Solution 1 平衡树维护动态凸壳

　　（这是什么啊？可以吃吗？）

　　最后维护的凸壳是要长这个样子:

　　由于边不是很稳定，因为插入或者删除一个点，维护边的信息很麻烦。。。所以考虑用平衡树维护点，如果能够找到前驱后继，那么就可以轻松地找到边的信息。

　　由于要能找到边的信息，所以将点按照横坐标进行排序。

　　为了更好地偷懒，所以假设第0个点连接第1个点的直线的斜率为$inf$，最后一条点的后一条直线斜率为$-inf$

　　考虑插入一个点。

　　向左向右分别找到第一个能与它组成凸壳的点（即满足斜率递减，对于它前面的点来说就是，满足连接点$pre$的前驱和$pre$的直线的斜率大于连接点$pre$和当前点的直线的斜率）

　　显然这个东西可以二分。其实是不用的。。。（时间复杂均摊有保证，每个点被访后要么成为了要的点，结束了操作，要么被删了，所以时间复杂度$Oleft(n ight)$）

　　当然还需要判断一下，这种情况

　　这种的话，只需要看看前面的点和后面的点以及当前点是否合法，不合法就说明当前点不该加入，否则就把它加入凸包。

　　至于如何维护前驱后继？

　　其实很简单，因为只在插入和删除的时候会发生改变，就在那个时候维护一下就好了。

　　现在考虑要在凸壳上找一个点，使得一条斜率为$k$的直线经过它，并且截距最大。

　　（这玩意儿显然可以三分）

　　

　　显然，当直线于凸壳相切（或者和凸壳的某一条边重合）的时候截距最大。不如一个不相切的情况：

　　显然可以通过平移使得它更优。

　　考虑什么时候相切呢？当经过的点的前一条直线的斜率大于等于$k$，以及它的后一条直线的斜率小于等于$k$时。

　　这个显然可以二分。于是剩下的就是代码的问题。。

　　（注意精度问题。。）

Code

/**
 * bzoj
 * Problem#1492
 * Accepted
 * Time: 1160ms
 * Memory: 11468k 
 */ 
#include <bits/stdc++.h>
using namespace std;
typedef bool boolean;
#define pnn pair<TreapNode*, TreapNode*>
#define fi first
#define sc second

const double eps = 1e-7;

int dcmp(double x) {
    if(fabs(x) < eps)    return 0;
    return (x > 0) ? (1) : (-1);
}

typedef class TreapNode {
    public:
        double x;
        double y;
        int rd;
        TreapNode *l, *r;
        TreapNode *pre, *suf;
}TreapNode;

#define Limit 200000

TreapNode pool[Limit];
TreapNode *top = pool;

TreapNode* newnode(double x, double y) {
    top->x = x, top->y = y;
    top->l = top->r = NULL;
    top->rd = rand();
    return top++;
}

typedef class DynamicConvexHull {
    public:
        TreapNode* root;
        
        DynamicConvexHull():root(NULL) {        } 
        
        pnn split(TreapNode* p, double x) {
            if(!p)    return pnn(NULL, NULL);
            pnn rt;
            if(p->x - eps > x) {
                rt = split(p->l, x);
                p->l = rt.sc, rt.sc = p;
            } else {
                rt = split(p->r, x);
                p->r = rt.fi, rt.fi = p;
            }
            return rt;
        }
        
        TreapNode* merge(TreapNode* a, TreapNode* b) {
            if(a == NULL)    return b;
            if(b == NULL)    return a;
            if(a->rd > b->rd) {
                a->r = merge(a->r, b);
                return a;
            }
            b->l = merge(a, b->l);
            return b;
        }
        
        double slope(TreapNode* a, TreapNode* b) {
            if(a == NULL)    return 1e100;
            if(b == NULL)    return -1e100;
            if(!dcmp(a->x - b->x))
                return (dcmp(a->y - b->y) == -1) ? (1e100) : (-1e100);
            return (a->y - b->y) / (a->x - b->x);
        }

        TreapNode* findPre(TreapNode* p, TreapNode* np) {
            TreapNode* rt = NULL;
            while(p) {
                if(slope(p->pre, p) - eps > slope(p, np))
                    rt = p, p = p->r;
                else
                    p = p->l;
            }
            return rt;
        }

        TreapNode* findSuf(TreapNode* p, TreapNode* np) {
            TreapNode* rt = NULL;
            while(p) {
                if(slope(np, p) - eps > slope(p, p->suf))
                    rt = p, p = p->l;
                else
                    p = p->r;    
            }
            return rt;
        }

        void insert(double x, double y) {
            TreapNode* pn = newnode(x, y);
            pnn pr = split(root, x);
            TreapNode* pre = findPre(pr.fi, pn);
            TreapNode* suf = findSuf(pr.sc, pn);
            pn->pre = pre, pn->suf = suf;
            if(slope(pre, pn) - eps >= slope(pn, suf)) {
                pr.fi = (pre) ? (split(pr.fi, pre->x + eps + eps).fi) : (NULL);
                pr.sc = (suf) ? (split(pr.sc, suf->x - eps - eps).sc) : (NULL);
                if(pre)    pre->suf = pn;
                if(suf)    suf->pre = pn;
                pr.fi = merge(pr.fi, pn);
            }
            root = merge(pr.fi, pr.sc);
        }
        
        TreapNode* query(double k) {
            TreapNode* p = root, *rt = NULL;
            double cur, cmp;
            while(p) {
                cur = slope(p->pre, p);
                if(cur - eps >= k) {
                    if(!rt || cur + eps < cmp)
                        rt = p, cmp = cur;
                    p = p->r;
                } else
                    p = p->l;
            }
            return rt;
        }
        
//        void debugOut(TreapNode* p) {
//            if(!p)    return;
//            debugOut(p->l);
//            cerr << "(" << p->x << ", " << p->y << ")" << endl;
//            debugOut(p->r);
//        }
}DynamicConvexHull;

int n;
double *A, *B, *rate;
double *f;
DynamicConvexHull dch;

inline void init() {
    scanf("%d", &n);
    A = new double[(n + 1)];
    B = new double[(n + 1)];
    f = new double[(n + 1)];
    rate = new double[(n + 1)];
    fill(f, f + n + 1, 0.0);
    scanf("%lf", f);
    for(int i = 1; i <= n; i++)
        scanf("%lf%lf%lf", A + i, B + i, rate + i);
}

inline void solve() {
    f[1] = f[0];
    double y = f[1] / (rate[1] * A[1] + B[1]);
    double x = rate[1] * y, k;
    dch.insert(x, y);
    for(int i = 2; i <= n; i++) {
        k = -A[i] / B[i];
        TreapNode* rt = dch.query(k);
        f[i] = A[i] * rt->x + B[i] * rt->y;
        f[i] = max(f[i], f[i - 1]);
        y = f[i] / (rate[i] * A[i] + B[i]);
        x = y * rate[i];
        dch.insert(x, y);
    }
    printf("%.3lf
", f[n]);
}

int main() {
    init();
    solve();
    return 0;
}    

Solution 2 CDQ分治

　　考虑当前要求出$[l, r]$中的dp值。

　　根据CDQ分治的常用套路，考虑左区间对右区间的贡献。

　　假设现在已经成功计算出左区间中的dp值，并将这些状态按照横坐标排序。

　　那么就可以用单调队列维护静态凸壳把左区间的凸壳建出来。

　　将右区间按照询问的斜率从大到小排序。

　　于是，这就变成了最智障的斜率优化问题了。。

　　但是$Oleft ( nlog^{2}n ight )$会不会T掉？

　　考虑计算右区间的时候并不需要按照横坐标排序，而是按照询问的斜率排序。

　　所以，在分治前按照询问的斜率排序，然后在回溯的过程中按照横坐标进行归并。

　　于是成功去掉一个$log$，总时间复杂度$Oleft ( nlog n ight )$

　　但是因为自带大常数，比别人的Splay慢好多，sad....

Code

/**
 * bzoj
 * Problem#1492
 * Accepted
 * Time: 1208ms
 * Memory: 8732k
 */ 
#include <bits/stdc++.h>
using namespace std;
typedef bool boolean;

const double eps = 1e-7;

int dcmp(double x) {
    if(fabs(x) < eps)    return 0;
    return (x > 0) ? (1) : (-1);
}

typedef class Query {
    public:
        double k;
        int id;
        
        boolean operator < (Query b) const {
            return k > b.k;
        }
}Query;

int n;
double *A, *B, *rate;
double *xs, *ys;
double *f;
Query *qs, *qbuf;
int* sta;

inline void init() {
    scanf("%d", &n);
    A = new double[(n + 1)];
    B = new double[(n + 1)];
    rate = new double[(n + 1)];
    xs = new double[(n + 1)];
    ys = new double[(n + 1)];
    f = new double[(n + 1)];
    qs = new Query[(n + 1)];
    qbuf = new Query[(n + 1)];
    sta = new int[(n + 1)];
    scanf("%lf", f);
    for(int i = 1; i <= n; i++) {
        scanf("%lf%lf%lf", A + i, B + i, rate + i);
        qs[i].k = -A[i] / B[i], qs[i].id = i;
    }
}

double slope(int s, int t) {
    if(dcmp(xs[s] - xs[t]) == 0)    return (1e100);
    return (ys[t] - ys[s]) / (xs[t] - xs[s]); 
}

boolean cmpPoint(int a, int b) {
    int d = dcmp(xs[a] - xs[b]);
    return (d == -1 || (d == 0 && dcmp(ys[a] - ys[b]) == -1));
}

void CDQDividing(int l, int r, int L, int R) {
    if(l == r) {
        f[l] = max(f[l], f[l - 1]);
        xs[l] = rate[l] * f[l] / (rate[l] * A[l] + B[l]);
        ys[l] = f[l] / (rate[l] * A[l] + B[l]);
        return; 
    }
    
    int mid = (l + r) >> 1, qL = L - 1, qR = mid;
    
    for(int i = L; i <= R; i++)
        if(qs[i].id <= mid)
            qbuf[++qL] = qs[i];
        else
            qbuf[++qR] = qs[i];
    for(int i = L; i <= qR; i++)
        qs[i] = qbuf[i];
    CDQDividing(l, mid, L, qL);
    
    int pl = 1, pr = 0, t = L;
    for(int i = L; i <= qL; i++) {
        while(pr - pl > 0 && dcmp(slope(sta[pr - 1], sta[pr]) - slope(sta[pr], qs[i].id)) != 1)    pr--;
        sta[++pr] = qs[i].id;
    }
    
    for(int i = mid + 1, id; i <= R; i++) {
        id = qs[i].id;
        while(pr - pl > 0 && dcmp(qs[i].k - slope(sta[pl], sta[pl + 1])) == -1)    pl++;
        f[id] = max(f[id], A[id] * xs[sta[pl]] + B[id] * ys[sta[pl]]);
    }
    
    CDQDividing(mid + 1, r, mid + 1, R);
    
    pl = L, pr = mid + 1;
    while(pl <= qL || pr <= R) {
        if((pr > R) || (pl <= qL && cmpPoint(qs[pl].id, qs[pr].id)))
            qbuf[t++] = qs[pl++];
        else
            qbuf[t++] = qs[pr++];
    }
    for(int i = L; i <= R; i++)
        qs[i] = qbuf[i];
}

inline void solve() {
    sort(qs + 1, qs + n + 1);
    fill(f + 1, f + n + 1, 0);
    CDQDividing(1, n, 1, n);
    printf("%.3lf
", f[n]);
}

int main() {
    init();
    solve();
    return 0;
}