• 懒人的福利?教你用set维护斜率优化凸包


    斜率优化题目大家肯定都做得不少了,有一些题目查询插入点的x坐标和查询斜率都不单调,这样就需要维护动态凸包并二分斜率。(例如bzoj1492)

    常规的做法是cdq分治或手写平衡树维护凸包,然而如果我不愿意写分治,也懒得打平衡树,怎么办呢?

    没关系,今天我告诉你怎么用一个set维护这种凸包。

    首先orzLH,没什么特殊意义,只是单纯的orz。

    我们定义f[i]表示在第i天能拥有的金券组数,按照第i天的比例。

    那么,我们要把前面的金券在今天卖出获得最多的钱,并在今天进行买入。

    所以,f[i]=max((f[j]*a[i]+f[j]/rate[j]*b[i])/(a[i]+rate[i]*b[i]))。

    除下去的东西是一个常数,扔掉。

    然后我们就有:

    t=max(a[i]*(f[j])+b[i]*(f[j]/rate[j]))。

    如果我们把f[j]看做x,f[j]/rate[j]看做y,我们有:

    t=a*x+b*y,两边同时除以b,得到:

    t/b=(a/b)x+y

    y=(t/b)-(a/b)*x

    好的,现在我们有一条斜率为-(a/b)的直线,要找一个点使之截距最大。

    这样我们维护一个右上1/4凸壳即可。

    怎么维护?

    我们考虑不用斜率优化,单纯水平序维护凸包,那么点是按照x坐标单增在平衡树上排列的。

    现在我们在每个点维护他与后面点连线斜率,我们会发现:这个斜率是单降的。

    所以,我们可以通过适当地转换cmp函数,来通过一个set完成两种比较。

    我们定义:

     1 int cmp; // 0 compare x , 1 compare slope .
     2 struct Point {
     3     double x,y,slop;
     4     friend bool operator < (const Point &a,const Point &b) {
     5         if( !cmp ) return a.x != b.x ? a.x < b.x : a.y < b.y;
     6         return a.slop > b.slop;
     7     }
     8     friend Point operator - (const Point &a,const Point &b) {
     9         return (Point){a.x-b.x,a.y-b.y};
    10     }
    11     friend double operator * (const Point &a,const Point &b) {
    12         return a.x * b.y - b.x * a.y;
    13     }
    14     inline double calc(double a,double b) const {
    15         return a * x + b * y;
    16     }
    17 };
    18 set<Point> st;

    插入就是正常凸包插入,最后再维护一下斜率就行了。注意弹出左边后迭代器会失效,所以需要重新lower_bound一下。(可能原来你的迭代器是原来的end,结果弹出左边后end改变了,两个end不同,然后你去弹出右边,访问无效迭代器,就直接RE了)

     1 inline void Pop_right(set<Point>::iterator nxt,const Point &p) {
     2     set<Point>::iterator lst;
     3     while(1) {
     4         lst = nxt , ++nxt;
     5         if( nxt == st.end() ) return;
     6         if( (*lst-p) * (*nxt-*lst) <= 0 ) return;
     7         st.erase(lst);
     8     }
     9 }
    10 inline void Pop_left(set<Point>::iterator prv,const Point &p) {
    11     set<Point>::iterator lst;
    12     while(prv!=st.begin()) {
    13         lst = prv , --prv;
    14         if( (*lst-*prv) * (p-*lst) <= 0 ) break;
    15         st.erase(lst);
    16     }
    17 }
    18 inline void insert(const Point &p) {
    19     cmp = 0;
    20     set<Point>::iterator prv,nxt,lst=st.lower_bound(p);
    21     if(lst!=st.begin()) Pop_left(--lst,p);
    22     lst=st.lower_bound(p);
    23     if(lst!=st.end()) Pop_right(lst,p);
    24     st.insert(p) , lst = st.find(p);
    25     if(lst!=st.begin()) {
    26         prv = lst , --prv;
    27         Point x = *prv;
    28         x.slop = ( p.y - x.y ) / ( p.x - x.x );
    29         st.erase(prv) , st.insert(x);
    30     }
    31     nxt = lst , ++nxt;
    32     if(nxt!=st.end()) {
    33         Point x = p , n = *nxt;
    34         x.slop = ( n.y - x.y ) / ( n.x - x.x );
    35         st.erase(lst) , st.insert(x);
    36     } else {
    37         Point x = p;
    38         x.slop = -1e18;
    39         st.erase(lst) , st.insert(x);
    40     }
    41 }

    查询的话就更改一下比较函数,然后特判一下边界防止RE就好。

     1 inline double query(const int id) {
     2     cmp = 1;
     3     const double k = -a[id] / b[id];
     4     set<Point>::iterator it = st.lower_bound((Point){0,0,k}); // it can't be st.end() if st isn't empty .
     5     if( it==st.end() ) return 0;
     6     double ret = it->calc(a[id],b[id]);
     7     if( it != st.begin() ) {
     8         --it;
     9         ret = max( ret , it->calc(a[id],b[id]) );
    10     }
    11     return ret;
    12 }

    所以整体代码:

    Bzoj1492:

     1 #include<cstdio>
     2 #include<algorithm>
     3 #include<set>
     4 #include<cmath>
     5 using namespace std;
     6 const int maxn=1e5+1e2;
     7  
     8 int cmp; // 0 compare x , 1 compare slope .
     9 struct Point {
    10     double x,y,slop;
    11     friend bool operator < (const Point &a,const Point &b) {
    12         if( !cmp ) return a.x != b.x ? a.x < b.x : a.y < b.y;
    13         return a.slop > b.slop;
    14     }
    15     friend Point operator - (const Point &a,const Point &b) {
    16         return (Point){a.x-b.x,a.y-b.y};
    17     }
    18     friend double operator * (const Point &a,const Point &b) {
    19         return a.x * b.y - b.x * a.y;
    20     }
    21     inline double calc(double a,double b) const {
    22         return a * x + b * y;
    23     }
    24 };
    25 set<Point> st;
    26 double a[maxn],b[maxn],rate[maxn],f[maxn],ans;
    27  
    28 inline void Pop_right(set<Point>::iterator nxt,const Point &p) {
    29     set<Point>::iterator lst;
    30     while(1) {
    31         lst = nxt , ++nxt;
    32         if( nxt == st.end() ) return;
    33         if( (*lst-p) * (*nxt-*lst) <= 0 ) return;
    34         st.erase(lst);
    35     }
    36 }
    37 inline void Pop_left(set<Point>::iterator prv,const Point &p) {
    38     set<Point>::iterator lst;
    39     while(prv!=st.begin()) {
    40         lst = prv , --prv;
    41         if( (*lst-*prv) * (p-*lst) <= 0 ) break;
    42         st.erase(lst);
    43     }
    44 }
    45 inline void insert(const Point &p) {
    46     cmp = 0;
    47     set<Point>::iterator prv,nxt,lst=st.lower_bound(p);
    48     if(lst!=st.begin()) Pop_left(--lst,p);
    49     lst=st.lower_bound(p);
    50     if(lst!=st.end()) Pop_right(lst,p);
    51     st.insert(p) , lst = st.find(p);
    52     if(lst!=st.begin()) {
    53         prv = lst , --prv;
    54         Point x = *prv;
    55         x.slop = ( p.y - x.y ) / ( p.x - x.x );
    56         st.erase(prv) , st.insert(x);
    57     }
    58     nxt = lst , ++nxt;
    59     if(nxt!=st.end()) {
    60         Point x = p , n = *nxt;
    61         x.slop = ( n.y - x.y ) / ( n.x - x.x );
    62         st.erase(lst) , st.insert(x);
    63     } else {
    64         Point x = p;
    65         x.slop = -1e18;
    66         st.erase(lst) , st.insert(x);
    67     }
    68 }
    69 inline double query(const int id) {
    70     cmp = 1;
    71     const double k = -a[id] / b[id];
    72     set<Point>::iterator it = st.lower_bound((Point){0,0,k}); // it can't be st.end() if st isn't empty .
    73     if( it==st.end() ) return 0;
    74     double ret = it->calc(a[id],b[id]);
    75     if( it != st.begin() ) {
    76         --it;
    77         ret = max( ret , it->calc(a[id],b[id]) );
    78     }
    79     return ret;
    80 }
    81  
    82 int main() {
    83     static int n;
    84     scanf("%d%lf",&n,&ans);
    85     for(int i=1;i<=n;i++) scanf("%lf%lf%lf",a+i,b+i,rate+i);
    86     for(int i=1;i<=n;i++) {
    87         ans = max( ans , query(i) );
    88         f[i] = ans * rate[i] / ( a[i] * rate[i] + b[i] );
    89         insert((Point){f[i],f[i]/rate[i],0});
    90     }
    91     printf("%0.3lf
    ",ans);
    92     return 0;
    93 }
    View Code

    不得不说STL的set跑的还是挺快的。

    这里是回档后的世界,无论你做什么,你都一定会这样做。

    而我努力改变命运,只是为了防止那一切重现。

    (来自某中二病晚期患者(不))

  • 相关阅读:
    遍历路径下的所有文件
    房间安排(题目168)
    创建BitMap
    字母统计(241)
    DataTable的Select方法
    ArcEngine中Feature对象的Shape属性和ShapeCopy属性
    C# 轻松获取路径中文件名、目录、扩展名等
    TreeList获取节点中的值
    【算法】LeetCode算法题-Two Sum
    JSP(一):初识JSP
  • 原文地址:https://www.cnblogs.com/Cmd2001/p/8492492.html
Copyright © 2020-2023  润新知