牛客网暑期ACM多校训练营（第三场）J Distance to Work【二分+圆交】

https://www.nowcoder.com/acm/contest/141/J

给出一个多边形，多边形内区域为城市。
给出几个工作地点，一个工作地点给出一个比例p/q，表示相比于所有的可能，我这个工作范围恰好站总生活面积的p/q。

我们可以二分圆的半径，使公共部分面积=$p/q*$城市面积。

圆与多边形交可以套用板子。

#include<bits/stdc++.h>
using namespace std;
#define M 500
#define eps 1e-8
const double PI = acos(-1.0);
const double INF = 1e4+7;
class pnt_type
{
public:
    double x,y;
};
class state_type
{
public:
    double angle;
    double CoverArea;
};
pnt_type pnt[M];
pnt_type center;
int n;int m;
double R;
bool read_data()
{
    int i;
    if (cin >> pnt[1].x >> pnt[1].y)
    {
        for (i=2; i<=n; i++) cin >> pnt[i].x >> pnt[i].y;
        return true;
    }
    return false;
}
inline double Area2(pnt_type &a,pnt_type &b,pnt_type &c)
{
    return (b.x - a.x) * (c.y - a.y) - (c.x - a.x) * (b.y - a.y);
}
inline double dot(pnt_type &a,pnt_type &b,pnt_type &c)
{
    return (b.x - a.x) * (c.x - a.x) + (b.y - a.y) * (c.y - a.y);
}
inline double dist(pnt_type &a,pnt_type &b)
{
    return sqrt((b.x - a.x) * (b.x - a.x) + (b.y - a.y) * (b.y - a.y));
}

void init()
{
    int i;
    double temp,sum;
    for (i=2; i<n; i++)
    {
        temp = Area2(pnt[1],pnt[i],pnt[i + 1]);
        sum += temp;
    }
    if (sum < 0) reverse(pnt + 1,pnt + n + 1);
    pnt[n + 1] = pnt[1];
}

inline bool inCircle(pnt_type &s)
{
    return dist(center,s) <= R;
}

bool SameSide(pnt_type a,pnt_type b)
{
    if (dist(a,center) > dist(b,center)) swap(a,b);
    return dot(a,b,center) < eps;
}

double ShadomOnCircle(pnt_type a,pnt_type b)
{
    double flag = Area2(center,a,b),res = 0;
    if (fabs(flag) < eps) return 0;

    bool ina = inCircle(a),inb = inCircle(b);
    if (ina && inb)
    {
        res = fabs(Area2(center,a,b)) / 2;
    }
    else if (!ina && !inb)
    {
        if (SameSide(a,b))
        {
            double theta = acos(dot(center,a,b) / dist(center,a) / dist(center,b));
            res = R * R * theta / 2;
        }
        else
        {
            double height = fabs(Area2(center,a,b)) / dist(a,b);
            double theta = acos(dot(center,a,b) / dist(center,a) / dist(center,b));
            if (height >= R)
            {
                res = R * R * theta / 2;
            }
            else
            {
                double _theta = 2 * acos(height / R);
                res = R * R * (theta - _theta) / 2 + R * R * sin(_theta) / 2;
            }
        }
    }
    else
    {
        if (!ina && inb) swap(a,b);
        double height = fabs(Area2(center,a,b)) / dist(a,b);
        double temp = dot(a,center,b);
        double theta = acos(dot(center,a,b) / dist(center,a) / dist(center,b)),theta1,theta2;
        if (fabs(temp) < eps)
        {
            double _theta = acos(height / R);
            res += R * height / 2 * sin(_theta);
            res += R * R / 2 * (theta - _theta);
        }
        else
        {
            theta1 = asin(height / R);
            theta2 = asin(height / dist(a,center));
            if (temp > 0)
            {
                res += dist(center,a) * R / 2 * sin(PI - theta1 - theta2);
                res += R * R / 2 * (theta + theta1 + theta2 - PI);
            }
            else
            {
                res += dist(center,a) * R / 2 * sin(theta2 - theta1);
                res += R * R / 2 * (theta - theta2 + theta1);
            }
        }
    }
    if (flag < 0) return -res;
    else return res;
}
double Cover()
{
    int i;
    double res = 0;
    for (i=1; i<=n; i++)
        res += ShadomOnCircle(pnt[i],pnt[i + 1]);
    return res;
}

int main()
{
    ios::sync_with_stdio(false);
    double ans;double A=0.0;
    pnt[0].x=pnt[0].y=0;
    cin>>n;
    init();
    read_data();
    for(int i=1;i<=n;i++){
        A+=(pnt[i-1].x*pnt[i].y-pnt[i].x*pnt[i-1].y);
    }
    A+=pnt[n].x*pnt[1].y-pnt[1].x*pnt[n].y;
    A=A/2;
    A=fabs(A);
    //cout<<A<<endl;
    cin>>m;
    double p,q;
    for(int i=1;i<=m;i++){
        cin >> center.x >> center.y;
        cin>>p>>q;
        double rate=p/q;
        double L1=0.0,R1=INF;
        int cs=0;
        while(cs<=100){
            cs++;
            R=(L1+R1)/2;
            ans=Cover();
            ans=fabs(ans);
            //double tmp=R*R*PI;
            if((A-ans)*q<=p*A){
                R1=R;
            }
            else{
                L1=R;
            }
        }
        printf("%.7f\n",L1);
    }

    //printf("%.2lf\n",ans);

    return 0;
}