题意:一个多边形,A点和B点,满足PB <= k * PA的P的范围与多边形的公共面积。
分析:这是个阿波罗尼斯圆。既然是圆,那么设圆的一般方程:(x + D/2) ^ 2 + (y + E/2) ^ 2 = (D ^ 2 + E ^ 2 - 4 * F) / 4,通过PB == PA * k解方程来求解圆心以及半径。然后就是套模板啦,上海交大的红书。
#include <bits/stdc++.h>
using namespace std;
#define lson l, mid, rt << 1
#define rson mid + 1, r, rt << 1 | 1
typedef long long ll;
const int N = 5e2 + 10;
const int INF = 0x3f3f3f3f;
const int MOD = 1e9 + 7;
const double EPS = 1e-10;
const double PI = acos (-1.0);
int dcmp(double x) {
if (fabs (x) < EPS) return 0;
else return x < 0 ? -1 : 1;
}
struct Point {
double x, y;
Point () {}
Point (double x, double y) : x (x), y (y) {}
Point operator + (const Point &r) const { //向量加法
return Point (x + r.x, y + r.y);
}
Point operator - (const Point &r) const {
return Point (x - r.x, y - r.y);
}
Point operator * (double p) const { //向量乘以标量
return Point (x * p, y * p);
}
Point operator / (double p) const { //向量除以标量
return Point (x / p, y / p);
}
bool operator < (const Point &r) const {
return x < r.x || (x == r.x && y < r.y);
}
bool operator == (const Point &r) const {
return dcmp (x - r.x) == 0 && dcmp (y - r.y) == 0;
}
};
typedef Point Vector;
Point read_point(void) {
double x, y; scanf ("%lf%lf", &x, &y);
return Point (x, y);
}
double polar_angle(Vector V) {
return atan2 (V.y, V.x);
}
double dot(Point a, Point b) {
return a.x * b.x + a.y * b.y;
}
double cross(Point a, Point b) {
return a.x * b.y - a.y * b.x;
}
double length(Vector V) {
return sqrt (dot (V, V));
}
double my_sqrt(double x) {
return sqrt (max (0.0, x));
}
Point ps[N];
double r, k;
double x1, _y1, x2, _y2;
int n;
double sqr(double x) {
return x * x;
}
struct Line {
Point p;
Vector v;
double r;
Line () {}
Line (const Point &p, const Vector &v) : p (p), v (v) {
r = polar_angle (v);
}
Point point(double a) {
return p + v * a;
}
};
struct Circle {
Point c;
double r;
Circle () {}
Circle (Point c, double r) : c (c), r (r) {}
Point point(double a) {
return Point (c.x + cos (a) * r, c.y + sin (a) * r);
}
};
int line_cir_inter(Line L, Circle C, double &t1, double &t2, vector<Point> &P) {
double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y;
double e = a * a + c * c, f = 2 * (a * b + c * d), g = b * b + d * d - C.r * C.r;
double delta = f * f - 4 * e * g;
if (dcmp (delta) < 0) return 0;
if (dcmp (delta) == 0) {
t1 = t2 = -f / (2 * e); P.push_back (L.point (t1));
return 1;
}
t1 = (-f - sqrt (delta)) / (2 * e);
t2 = (-f + sqrt (delta)) / (2 * e);
if (t1 > t2) swap (t1, t2);
if (dcmp (t1) > 0 && dcmp (t1 - 1) < 0) P.push_back (L.point (t1));
if (dcmp (t2) > 0 && dcmp (t2 - 1) < 0) P.push_back (L.point (t2));
return (int) P.size ();
}
double sector_area(Point a, Point b) {
double theta = polar_angle (a) - polar_angle (b);
while (dcmp (theta) <= 0) theta += 2 * PI;
while (theta > 2 * PI) theta -= 2 * PI;
theta = min (theta, 2 * PI - theta);
return r * r * theta / 2;
}
double cal(Point a, Point b) {
double t1, t2;
bool ina = dcmp (length (a) - r) < 0;
bool inb = dcmp (length (b) - r) < 0;
if (ina && inb) return fabs (cross (a, b)) / 2.0;
vector<Point> p;
int num = line_cir_inter (Line (a, b - a), Circle (Point (0, 0), r), t1, t2, p);
if (ina) return sector_area (b, p[0]) + fabs (cross (a, p[0])) / 2.0;
if (inb) return sector_area (p[0], a) + fabs (cross (p[0], b)) / 2.0;
if (num == 2) return sector_area (a, p[0]) + sector_area (p[1], b) + fabs (cross (p[0], p[1])) / 2.0;
return sector_area (a, b);
}
double cir_poly_area() {
double ret = 0;
for (int i=0; i<n; ++i) {
int sgn = dcmp (cross (ps[i], ps[i+1]));
if (sgn != 0) {
ret += sgn * cal (ps[i], ps[i+1]);
}
}
return ret;
}
void init(void) {
double D = (2 * x2 - 2 * sqr (k) * x1) / (1 - sqr (k));
double E = (2 * _y2 - 2 * sqr (k) * _y1) / (1 - sqr (k));
double F = (sqr (k*x1)+sqr (k*_y1)-sqr (x2)-sqr(_y2)) / (1 - sqr (k));
double x0 = D / 2, y0 = E / 2;
r = sqrt (F + sqr (D)/4 + sqr (E)/4);
for (int i=0; i<n; ++i) {
ps[i].x -= x0; ps[i].y -= y0;
}
ps[n] = ps[0];
}
int main(void) {
int cas = 0;
while (scanf ("%d%lf", &n, &k) == 2) {
for (int i=0; i<n; ++i) {
ps[i] = read_point ();
}
scanf ("%lf%lf", &x1, &_y1);
scanf ("%lf%lf", &x2, &_y2);
init ();
printf ("Case %d: %.10f\n", ++cas, fabs (cir_poly_area ()));
}
//cout << "Time elapsed: " << 1.0 * clock() / CLOCKS_PER_SEC << " s.\n";
return 0;
}