Luogu P2333 [SCOI2006]一孔之见 计算几何

#include <cstdio>#include <cmath>#include <algorithm>#define eps 1e-5#define INF 1e20using namespace std;const double pi=acos(-1.0);inline double sq(double x) { return x*x; }struct Point {    double x,y;    Point() {}    Point(double _x,double _y):x(_x),y(_y) {}    void scan() { scanf("%lf%lf",&x,&y); }    Point operator - (const Point& rhs) const {        return Point(x-rhs.x,y-rhs.y);    }    double operator ^ (const Point& rhs) const {        return x*rhs.x+y*rhs.y;    }    double operator * (const Point& rhs) const {        return x*rhs.y-y*rhs.x;    }    double operator | (const Point& rhs) const {        return sqrt(sq(x-rhs.x)+sq(y-rhs.y));    }    double getlen() {        return sqrt(sq(x)+sq(y));    }    void rotate90() {        double _x=x,_y=y;        x=-_y, y=_x;        return ;    }}p[55],o(0,0);typedef Point Vector;int n;double S;double calc_dist(Point a,Point b,Point c) { /// c to a->b    Vector v1=a-c,v2=b-c;    double tmp=fabs(v1*v2);    return tmp/(a|b);}double calc_angle(Vector a,Vector b) {    return acos((a^b)/a.getlen()/b.getlen());}bool check(double r) {    double s=0;    ///calculate current area    for(int i=1;i<=n;++i) {        Point p1=p[i],p2=p[i+1];        double d1=o|p1,d2=o|p2,d=calc_dist(p1,p2,o);        if(d1>d2) swap(d1,d2), swap(p1,p2);        Vector v1=p1-o,v2=p2-o,v=p1-p2;        v.rotate90();        bool side=(v*v1)*(v*v2)<0?false:true;        double ang=calc_angle(v1,v2);        if(d>=r) s+=ang*sq(r)/2;        else if(d1<=r && d2<=r) s+=fabs(v1*v2)/2;        else if(d1>=r && d2>=r) {            if(!side) {                ///calculate triangle                double _d=sqrt(sq(r)-sq(d));                s+=_d*d;                ///calculate circle                double rest_ang=ang-acos(d/r)*2;                s+=rest_ang*sq(r)/2;            }            else s+=ang*sq(r)/2;        }        else {            double ang1=acos(d/d1),ang2=acos(d/r),ang3=ang;            if(!side) ang2+=ang1;            else ang2-=ang1;            ang3-=ang2;            ///calculate triangle            s+=sin(ang2)*r*d1/2;            ///calculate circle            s+=ang3*sq(r)/2;        }    }    return s>=S;}int main() {    scanf("%d%lf",&n,&S);    for(int i=1;i<=n;++i) p[i].scan();    p[n+1]=p[1];    double l=eps,r=INF,mid,ans;    while(r-l>eps) {        mid=(l+r)/2;        if(check(mid)) ans=mid, r=mid-eps;        else l=mid+eps;    }    printf("%.2f\n",ans);    return 0;}