hdu 5027 Help!（计算几何 三分求极值）

hdu 5027 Help!

根据题意，在陆地上的路线解分为两部分：①从起点直线到达；②无法直线到达。救到人后返回陆地的时间可以先求出来。

可以求出起点到凸多边形的相切点，在两点间的边可直线到达。找出使总时间最小的入水点，可以用三分搜索求极值（列出公式求导后可发现是一个单峰函数）。

考虑入水点无法直线到达的情况，原理与上类似

计算几何的实现还是相当蛋疼的。有几个比较烦人的地方：

1、找切点。对于每个点 p 求出 sz[p] = ( vector(p, p+1) ^vector(p, point(Xo, Yo)） )>=0，若sz[p] != sz[p+1]，则说明点 p+1 是切点。而且，若起点就在多边形的边上，则切点为其所处的边的两端点。

2、水点无法直线到达的情况下，应先到达切点，再沿着边寻找入水点

3、对于求出的切点 a，b。最好改变其位置使 a->b 之间的点属于①，b->a间的点属于②。

#include <cstdio>#include <cstring>#include <iostream>#include <algorithm>#include <cmath>using namespace std;const int MAXN = 50010;const double eps = 1e-4, inf = 1e100;double Ts, Vr, Vs, rres, aans, len;int n, sz[MAXN];inline int cmp(double a, double b){    return ((a-b)>eps) - ((a-b)<-eps);}struct point{    double x, y;    point() {}    point(double x, double y)    {        this->x = x;        this->y = y;    }    double operator ^ (const point & a)    {        return x*a.y-y*a.x;    }    double operator * (const point & a)    {        return x*a.x + y*a.y;    }    double dis2()    {        return (x*x + y*y);    }    double dis()    {        return sqrt(x*x + y*y);    }    point operator - (const point &a )    {        return point(x-a.x, y-a.y);    }    point operator + (const point &a )    {        return point(x+a.x, y+a.y);    }    point operator / (const double &a )    {        return point(x/a, y/a);    }    point operator * (const double a )    {        return point(x*a, y*a);    }    bool operator == (const point &a) const    {        return cmp(x,a.x)==0 && cmp(y,a.y) == 0;    }    point trunc(double r)    {        double l=dis();        if (!cmp(l,0))return *this;        return point(x*r/l,y*r/l);    }} da[MAXN], OO, PP, frm;double persToedge(point o, point p1, point p2){    if (p1 == p2)    {        return (o-p1).dis();    }    point l1 = o-p1, l2 = o-p2, l3 = p1-p2;    if (cmp(l1*l3, 0)==1)    {        return l1.dis();    }    if (cmp(l2*l3, 0)==-1)    {        return l2.dis();    }    return fabs(l1^l2)/(l3.dis());}double dis_Point_to_Line(point o, point a, point b){    return (fabs((a-o)^(b-o)))/((a-b).dis());}void circle_X_line(point o, double r, point a, point b, int &nc, point *res){    int k ;    double dis = dis_Point_to_Line(o, a, b);    k = cmp(r, dis);    if (k==-1)    {        nc = 0;        return;    }    point xx = a+((b-a)*(((b-a)*(o-a))/ ((b-a).dis2()) ));    if (k == 0)    {        nc = 1;        res[0]=xx;        return;    }    dis = r*r - dis*dis;    nc = 2;    res[0] = xx-((b-a).trunc(dis));    res[1] = xx+((b-a).trunc(dis));    if (res[0].x > res[1].x) swap(res[0], res[1]);}inline bool isPointOnSeg(point o, point a, point b){    if (a.x > b.x) swap(a, b);    return cmp(b.x, o.x)>=0 && cmp(o.x, a.x)>=0;}double getNum(point x, int cs = -1){    double tp = (PP-x).dis()/Vs;    if (cs > 0)        return tp += ((len + (x-frm).dis()) / Vr);    return (OO-x).dis()/Vr + tp;}double tripleSearch(point bg, point ed, int cs=-1){    point ml, mr;    while ((bg-ed).dis() > eps)    {        ml = (bg*2+ed)/3.0;        mr = (bg+ed*2)/3.0;        if (getNum(ml, cs) < getNum(mr, cs)) ed = mr;        else bg = ml;    }    return getNum(bg, cs);}void cal(int bg, int ed, int cs){    int nc, i = bg;    point ap[2], p1, p2;    while (true)    {        nc = 0;        p1 = da[(i+n)%n];        p2 = da[(i+1+n)%n];        if (cs == 1) frm = p1;        else frm = p2;        circle_X_line(PP, (Ts-aans)*Vs, p1, p2, nc, ap);        if (nc == 1)        {            if (isPointOnSeg(ap[0], p1, p2))                rres = min(rres, (len + (ap[0]-p1).dis())/Vr + (PP-ap[0]).dis()/Vs);        }        else if (nc == 2)        {            point e1, e2;            if (cmp(ap[0].x, max(p1.x, p2.x))>0 ||                    cmp(ap[1].x, min(p1.x, p2.x))<0) continue;            if (isPointOnSeg(p1, ap[0], ap[1])) e1 = p1;            else e1 = ap[0];            if (isPointOnSeg(p2, ap[0], ap[1])) e2 = p2;            else e2 = ap[1];            rres = min(rres, tripleSearch(e1, e2, 1));        }        len += (p1-p2).dis();        if (i == ed) break;        i = (i+cs+n)%n;    }}void solve(){    aans = inf;    rres = inf;    int bg, ed, nc, ssa[2], hc = 0;    point ap[2], p1, p2, ssp[2];    if (cmp( (da[1]-da[0]) ^ (da[2]-da[1]), 0) == -1)        reverse(da, da+n);    for (int i = 0; i< n; ++i)        aans = min(aans, persToedge(PP, da[i], da[(i+1)%n])/Vs*2);    if (aans > Ts)    {        puts("-1");        return;    }    for (int i = 0; i< n; ++i)        sz[i] = cmp( (da[(i+1)%n]-da[i]) ^ (OO-da[i]), 0)>=0;    for (int i = 0, j = 0; i< n; ++i)        if (sz[i] != sz[(i+1)%n])            ssa[j] = (i+1)%n, ssp[j] = da[ssa[j]], ++j, hc = 1;    if (!hc)    //    {        // rres = 0;        for (int i = 0; i< n; ++i)            if (cmp(dis_Point_to_Line(OO, da[i], da[(i+1)%n]), 0) == 0)            {                ssa[0] = i;                ssa[1] = (i+1)%n;                ssp[0] = da[i];                ssp[1] = da[(i+1)%n];                break;            }    }    else    {        if (sz[ssa[0]] != 0) swap(ssa[0], ssa[1]), swap(ssp[0], ssp[1]);        bg = ssa[0], ed = ssa[1] < ssa[0]?ssa[1]+n:ssa[1];        for (int i = bg; i<ed; ++i)        {            nc = 0;            p1 = da[i%n];            p2 = da[(i+1)%n];            if (p1.x > p2.x) swap(p1, p2);            circle_X_line(PP, (Ts-aans)*Vs, p1, p2, nc, ap);            if (nc == 0) continue;            else if (nc == 1)            {                if (isPointOnSeg(ap[0], p1, p2))                    rres = min(rres, (ap[0]-OO).dis()/Vr+(ap[0]-PP).dis()/Vs);            }            else            {                point e1, e2;                if (cmp(ap[0].x, p2.x)>0 || cmp(ap[1].x, p1.x) < 0) continue;                if (isPointOnSeg(p1, ap[0], ap[1])) e1 = p1;                else e1 = ap[0];                if (isPointOnSeg(p2, ap[0], ap[1])) e2 = p2;                else e2 = ap[1];                rres = min(rres, tripleSearch(e1, e2));            }        }    }    len = (ssp[1]-OO).dis();    cal(ssa[1], (ssa[0]-1+n)%n, 1);    len = (ssp[0]-OO).dis();    cal((ssa[0]-1+n)%n, ssa[1], -1);    if (rres == inf) puts("-1");    else printf("%.2lf\n", rres + aans);}int main(){#ifndef ONLINE_JUDGE    freopen("in.txt", "r", stdin);    //freopen("out.txt", "w", stdout);#endif // ONLINE_JUDGE    int t;    scanf("%d", &t);    while (t--)    {        scanf("%lf%lf%lf",  &Ts, &Vr, &Vs);        scanf("%lf%lf%lf%lf", &OO.x, &OO.y, &PP.x, &PP.y);        scanf("%d", &n);        for (int i = 0; i< n; ++i)            scanf("%lf%lf", &da[i].x, &da[i].y);        solve();    }    return 0;}