计算几何 2017.4.10

1、POJ 1569 Myacm Triangles

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <queue>
#include <vector>
#include <stack>
#include <map>
#include <set>
#include <cmath>
#include <cctype>
#include <ctime>
#include <cassert>

using namespace std;

#define REP(i, n) for (int i = 0; i < (n); ++i)

struct Point {
  double x, y;
  Point(double x = 0.0, double y = 0.0):x(x),y(y) {}
};
typedef Point Vector;
struct Seg { Point a, b; };

const double INF = 1e20;
const double eps = 1e-10;
const double PI = acos(-1.0);
char str[2];
Point point[20];
int n;

Point operator + (const Point& a, const Vector& b) { return Point(a.x + b.x, a.y + b.y); }
Vector operator - (const Point& a, const Point& b) { return Vector(a.x - b.x, a.y - b.y); }
Vector operator * (const Vector& a, double k) { return Vector(a.x * k, a.y * k); }
int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; }
Point read_point() { double x, y; scanf("%lf %lf", &x, &y); return Point{x, y}; }
double cross(const Vector& a, const Vector& b) { return a.x * b.y - a.y * b.x; }
bool is_parallel(const Seg& p, const Seg& q) {
    return (abs((p.a.y - p.b.y) * (q.a.x - q.b.x) - (p.a.x - p.b.x) * (q.a.y - q.b.y)) < eps) ? true : false;
}
double area2(const Point& a, const Point& b, const Point& c) { return cross(b - a, c - a); }

int main() {
#ifdef __AiR_H
    freopen("in.txt", "r", stdin);
//    freopen("out.txt", "w", stdout);
#endif // __AiR_H
    while (scanf("%d", &n) && n != 0) {
        int ansA, ansB, ansC; double Max = 0;
        REP(i, n) { scanf("%s", str); point[i] = read_point(); }
        bool flag; double area_sum, area_1, area_2, area_3;
        for (int i = 0; i < n - 2; ++i) {
            for (int j = i + 1; j < n - 1; ++j) {
                for (int k = j + 1; k < n; ++k) {
                    if (is_parallel(Seg{point[i], point[j]}, Seg{point[j], point[k]})) { continue; }
                    flag = true; area_sum = fabs(area2(point[i], point[j], point[k]));
                    REP(u, n) if (u != i && u != j && u != k) {
                        area_1 = fabs(area2(point[u], point[j], point[k]));
                        area_2 = fabs(area2(point[i], point[u], point[k]));
                        area_3 = fabs(area2(point[i], point[j], point[u]));
                        if (dcmp(area_sum - area_1 - area_2 - area_3) == 0) { flag = false; break; }
                    }
                    if (!flag) { continue; }
                    if (area_sum > Max) { Max = area_sum; ansA = i; ansB = j; ansC = k; }
                }
            }
        }
        printf("%c%c%c\n", ansA + 'A', ansB + 'A', ansC + 'A');
    }
#ifdef __AiR_H
    printf("Time used = %.2fs\n", (double)clock() / CLOCKS_PER_SEC);
#endif // __AiR_H
    return 0;
}

2、POJ 3304 Segments

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <queue>
#include <vector>
#include <stack>
#include <map>
#include <set>
#include <cmath>
#include <cctype>
#include <ctime>
#include <cassert>

using namespace std;

#define REP(i, n) for (int i = 0; i < (n); ++i)

struct Point {
  double x, y;
  Point(double x = 0.0, double y = 0.0):x(x),y(y) {}
};
struct Seg { Point a, b; };
typedef Point Vector;
typedef Seg Line;

const double INF = 1e20;
const double eps = 1e-10;
const double PI = acos(-1.0);
char str[2];
Point point[20];
int T, n;
Seg seg[110];
vector<Point> v;

Point operator + (const Point& a, const Vector& b) { return Point(a.x + b.x, a.y + b.y); }
Vector operator - (const Point& a, const Point& b) { return Vector(a.x - b.x, a.y - b.y); }
Vector operator * (const Vector& a, double k) { return Vector(a.x * k, a.y * k); }
int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; }
Point read_point() { double x, y; scanf("%lf %lf", &x, &y);  v.push_back(Point{x, y}); return Point{x, y}; }
Seg read_seg() { return Seg{read_point(), read_point()}; }
double dot(const Vector& a, const Vector& b) { return a.x * b.x + a.y * b.y; }
double cross(const Vector& a, const Vector& b) { return a.x * b.y - a.y * b.x; }
bool is_parallel(const Seg& p, const Seg& q) {
    return (abs((p.a.y - p.b.y) * (q.a.x - q.b.x) - (p.a.x - p.b.x) * (q.a.y - q.b.y)) < eps) ? true : false;
}
bool on_line(const Point& a, const Line& p) { return is_parallel(Seg{a, p.a}, Seg{a, p.b}); }
bool is_intersect(const Line& p, const Seg& q) {
    double c1 = cross(p.b - p.a, q.a - p.a), c2 = cross(p.b - p.a, q.b - p.a);
    if (dcmp(c1) * dcmp(c2) < 0) { return true; }
    return on_line(q.a, p) || on_line(q.b, p);
}

int main() {
#ifdef __AiR_H
    freopen("in.txt", "r", stdin);
//    freopen("out.txt", "w", stdout);
#endif // __AiR_H
    scanf("%d", &T);
    while (T--) {
        scanf("%d", &n);
        v.clear();
        bool flag = false, all_its = true; Point its;
        REP(i, n) { seg[i] = read_seg(); }
        if (n <= 2) { printf("Yes!\n"); continue; }
        for (int i = 0; i < n * 2 - 1; ++i) {
            for (int j = i + 1; j < n * 2; ++j) {
                if (dcmp(v[i].x - v[j].x) == 0 && dcmp(v[i].y - v[j].y) == 0) { continue; }
                all_its = true;
                for (int k = 0; k < n; ++k) {
                    if (is_parallel(seg[k], Seg{v[i], v[j]})) {
                        if (!is_parallel(seg[k], Seg{seg[k].a, v[j]})) { all_its = false; break; }
                        else { continue; }
                    }
                    if (!is_intersect(Line{v[i], v[j]}, seg[k])) { all_its = false; break; }
                }
                if (all_its) { flag = true; break; }
            }
            if (flag) { break; }
        }
        if (flag) { printf("Yes!\n"); }
        else { printf("No!\n"); }
    }
#ifdef __AiR_H
    printf("Time used = %.2fs\n", (double)clock() / CLOCKS_PER_SEC);
#endif // __AiR_H
    return 0;
}

3、POJ 1039 Pipe

参考：《北京大学暑期课《 ACM/ICPC竞赛训练》》

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <queue>
#include <vector>
#include <stack>
#include <map>
#include <set>
#include <cmath>
#include <cctype>
#include <ctime>
#include <cassert>

using namespace std;

#define REP(i, n) for (int i = 0; i < (n); ++i)

struct Point {
  double x, y;
  Point(double x = 0.0, double y = 0.0):x(x),y(y) {}
};
struct Seg { Point a, b; };
typedef Point Vector;
typedef Seg Line;

const double INF = 1e20;
const double eps = 1e-10;
const double PI = acos(-1.0);
int n;
Point point1[30], point2[30];

Point operator + (const Point& a, const Vector& b) { return Point(a.x + b.x, a.y + b.y); }
Vector operator - (const Point& a, const Point& b) { return Vector(a.x - b.x, a.y - b.y); }
Vector operator * (const Vector& a, double k) { return Vector(a.x * k, a.y * k); }
int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; }
Point read_point() { double x, y; scanf("%lf %lf", &x, &y); return Point{x, y}; }
Seg read_seg() { return Seg{read_point(), read_point()}; }
double dot(const Vector& a, const Vector& b) { return a.x * b.x + a.y * b.y; }
double cross(const Vector& a, const Vector& b) { return a.x * b.y - a.y * b.x; }
bool is_parallel(const Seg& p, const Seg& q) {
    return (abs((p.a.y - p.b.y) * (q.a.x - q.b.x) - (p.a.x - p.b.x) * (q.a.y - q.b.y)) < eps) ? true : false;
}
bool on_line(const Point& a, const Line& p) { return is_parallel(Seg{a, p.a}, Seg{a, p.b}); }
bool is_intersect(const Line& p, const Seg& q) {
    double c1 = cross(p.b - p.a, q.a - p.a), c2 = cross(p.b - p.a, q.b - p.a);
    if (dcmp(c1) * dcmp(c2) < 0) { return true; }
    return on_line(q.a, p) || on_line(q.b, p);
}
Point line_intersection(const Point& p, const Vector& v, const Point& q, const Vector& w) {
    Vector u = p - q; double t = cross(w, u) / cross(v, w); return p + v * t;
}

int main() {
#ifdef __AiR_H
    freopen("in.txt", "r", stdin);
//    freopen("out.txt", "w", stdout);
#endif // __AiR_H
    while (scanf("%d", &n) && n != 0) {
        bool flag = false, all_its = true;
        REP(i, n) { point1[i] = read_point(); point2[i] = Point{point1[i].x, point1[i].y - 1.0}; }
        double Max = point1[0].x; int cur = 1; Point its;
        REP(i, n) REP(j, n) if (i != j) {
            all_its = true;
            REP(k, n) {
                if (!is_intersect(Line{point2[i], point1[j]}, Seg{point1[k], point2[k]})) {
                    if (k >= cur) {
                        cur = k;
                        its = line_intersection(point2[i], point1[j] - point2[i], point2[k], point1[k] - point2[k]);
                        if (its.y > point1[k].y) {
                            its = line_intersection(point2[i], point1[j] - point2[i], point1[k - 1], point1[k] - point1[k - 1]);
                        } else {
                            its = line_intersection(point2[i], point1[j] - point2[i], point2[k - 1], point2[k] - point2[k - 1]);
                        }
                        Max = max(Max, its.x);
                    }
                    all_its = false; break;
                }
            }
            if (all_its) { flag = true; }
        }
        if (!flag) { printf("%.2f\n", Max); }
        else { printf("Through all the pipe.\n"); }
    }
#ifdef __AiR_H
    printf("Time used = %.2fs\n", (double)clock() / CLOCKS_PER_SEC);
#endif // __AiR_H
    return 0;
}

4、POJ 1066 Treasure Hunt

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <queue>
#include <vector>
#include <stack>
#include <map>
#include <set>
#include <cmath>
#include <cctype>
#include <ctime>
#include <cassert>

using namespace std;

#define REP(i, n) for (int i = 0; i < (n); ++i)

struct Point {
  double x, y;
  Point(double x = 0.0, double y = 0.0):x(x),y(y) {}
};
struct Seg { Point a, b; };
typedef Point Vector;
typedef Seg Line;

const double INF = 1e20;
const double eps = 1e-10;
const double PI = acos(-1.0);
int n;
Seg seg[40];
Point E;
vector<Point> v;

Point operator + (const Point& a, const Vector& b) { return Point(a.x + b.x, a.y + b.y); }
Vector operator - (const Point& a, const Point& b) { return Vector(a.x - b.x, a.y - b.y); }
Vector operator * (const Vector& a, double k) { return Vector(a.x * k, a.y * k); }
int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; }
Point read_point() { double x, y; scanf("%lf %lf", &x, &y); v.push_back(Point{x, y}); return Point{x, y}; }
Seg read_seg() { return Seg{read_point(), read_point()}; }
double dot(const Vector& a, const Vector& b) { return a.x * b.x + a.y * b.y; }
double cross(const Vector& a, const Vector& b) { return a.x * b.y - a.y * b.x; }
bool on_seg(const Point& a, const Seg& p) {
    return dcmp(cross(p.a - a, p.b - a)) == 0 && dcmp(dot(p.a - a, p.b - a)) <= 0;
}
bool is_seg_intersect(const Seg& p, const Seg& q) {
    double c1 = cross(p.b - p.a, q.a - p.a), c2 = cross(p.b - p.a, q.b - p.a);
    double c3 = cross(q.b - q.a, p.a - q.a), c4 = cross(q.b - q.a, p.b - q.a);
    if (dcmp(c1) * dcmp(c2) < 0 && dcmp(c3) * dcmp(c4) < 0) { return true; }
    return on_seg(p.a, q) || on_seg(p.b, q) || on_seg(q.a, p) || on_seg(q.b, p);
}

int main() {
#ifdef __AiR_H
    freopen("in.txt", "r", stdin);
//    freopen("out.txt", "w", stdout);
#endif // __AiR_H
    scanf("%d", &n);
    REP(i, n) { seg[i] = read_seg(); } E = read_point();
    if (n == 0) { printf("Number of doors = 1\n"); return 0; }
    int ans = 100, cnt;
    REP(i, 2 * n) {
        cnt = 0;
        REP(j, n) {
            if (is_seg_intersect(seg[j], Seg{v[i], E})) { ++cnt; }
        }
        ans = min(ans, cnt);
    }
    printf("Number of doors = %d\n", ans);
#ifdef __AiR_H
    printf("Time used = %.2fs\n", (double)clock() / CLOCKS_PER_SEC);
#endif // __AiR_H
    return 0;
}

5、UVa 12304 2D Geometry 110 in 1!

参考：《算法竞赛入门经典-训练指南》

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <queue>
#include <vector>
#include <stack>
#include <map>
#include <set>
#include <cmath>
#include <cctype>
#include <ctime>
#include <cassert>

using namespace std;

#define REP(i, n) for (int i = 0; i < (n); ++i)

struct Point {
  double x, y;
  Point(double x = 0.0, double y = 0.0):x(x),y(y) {}
};
struct Seg { Point a, b; };
typedef Point Vector;
typedef Seg Line;
struct Circle {
    Point c; double r;
    Point point(double a) {
        return Point{c.x + cos(a) * r, c.y + sin(a) * r};
    }
};

const double INF = 1e20;
const double eps = 1e-10;
const double PI = acos(-1.0);
char str[1000];

Point operator + (const Point& a, const Vector& b) { return Point(a.x + b.x, a.y + b.y); }
Vector operator - (const Point& a, const Point& b) { return Vector(a.x - b.x, a.y - b.y); }
Vector operator * (const Vector& a, double k) { return Vector(a.x * k, a.y * k); }
Vector operator / (const Vector& a, double k) { return Vector{a.x / k, a.y / k}; }
int dcmp(double x) { if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1; }
bool operator == (const Point& a, const Point& b) {
    return dcmp(a.x - b.x) == 0 && dcmp(a.y - b.y) == 0;
}
bool operator < (const Point& a, const Point& b) {
    return a.x < b.x || (a.x == b.x && a.y < b.y);
}
Point read_point() { double x, y; scanf("%lf %lf", &x, &y); return Point{x, y}; }
Seg read_seg() { return Seg{read_point(), read_point()}; }
Line read_line() { return Line{read_point(), read_point()}; }
double dot(const Vector& a, const Vector& b) { return a.x * b.x + a.y * b.y; }
double cross(const Vector& a, const Vector& b) { return a.x * b.y - a.y * b.x; }
double length(const Vector& a) { return sqrt(dot(a, a)); }
Vector rotate(const Vector& a, double rad) {
    return Vector(a.x * cos(rad) - a.y * sin(rad), a.x * sin(rad) + a.y * cos(rad));
}
Point line_intersection(const Point& p, const Vector& v, const Point& q, const Vector& w) {
    Vector u = p - q; double t = cross(w, u) / cross(v, w); return p + v * t;
}
double distance(const Point& p, const Line& l) {
    Vector v1 = l.b - l.a, v2 = p - l.a;
    return fabs(cross(v1, v2) / length(v1));
}
int get_line_circle_its(const Line& l, const Circle& C, double& t1, double& t2, vector<Point>& sol) {
    double a = l.b.x - l.a.x, b = l.a.x - C.c.x, c = l.b.y - l.a.y, d = l.a.y - C.c.y;
    double e = a * a + c * c, f = 2.0 * (a * b + c * d), g = b * b + d * d - C.r * C.r;
    double delta = f * f - 4.0 * e * g; double dist = distance(C.c, l);
    if (dcmp(dist - C.r) > 0) { return 0; }
    if (dcmp(dist - C.r) == 0) {
        t1 = -f / (2.0 * e); t2 = t1; sol.push_back(l.a + (l.b - l.a) * t1); return 1;
    }
    t1 = (-f - sqrt(delta)) / (2.0 * e); sol.push_back(l.a + (l.b - l.a) * t1);
    t2 = (-f + sqrt(delta)) / (2.0 * e); sol.push_back(l.a + (l.b - l.a) * t2);
    return 2;
}
double angle(const Vector& v) { return atan2(v.y, v.x); }
int get_circle_circle_its(Circle c1, Circle c2, vector<Point>& sol) {
    double d = length(c1.c - c2.c);
    if (dcmp(d) == 0) {
        if (dcmp(c1.r - c2.r) == 0) { return -1; } return 0;
    }
    if (dcmp(c1.r + c2.r - d) < 0) { return 0; }
    if (dcmp(fabs(c1.r - c2.r) - d) > 0) { return 0; }
    double a = angle(c2.c - c1.c);
    double da = acos((c1.r * c1.r + d * d - c2.r * c2.r) / (2.0 * c1.r * d));
    Point p1 = c1.point(a - da), p2 = c1.point(a + da);
    sol.push_back(p1); if (p1 == p2) { return 1; }
    sol.push_back(p2); return 2;
}
Circle circum_circle(const Point& p1, const Point& p2, const Point& p3) {
    double Bx = p2.x - p1.x, By = p2.y - p1.y;
    double Cx = p3.x - p1.x, Cy = p3.y - p1.y;
    double D = 2.0 * (Bx * Cy - By * Cx);
    double cx = (Cy * (Bx * Bx + By * By) - By * (Cx * Cx + Cy * Cy)) / D + p1.x;
    double cy = (Bx * (Cx * Cx + Cy * Cy) - Cx * (Bx * Bx + By * By)) / D + p1.y;
    Point p = Point(cx, cy);
    return Circle{p, length(p1 - p)};
}
Circle inscribed_circle(const Point& p1, const Point& p2, const Point& p3) {
    double a = length(p2 - p3), b = length(p3 - p1), c = length(p1 - p2);
    Point p = (p1 * a + p2 * b + p3 * c) / (a + b + c);
    return Circle{p, distance(p, Line{p1, p2})};
}
int get_tangents(Point p, Circle C, Vector* v) {
    Vector u = C.c - p; double dist = length(u);
    if (dist < C.r) { return 0; }
    if (dcmp(dist - C.r) == 0) { v[0] = rotate(u, PI / 2); return 1; }
    double ang = asin(C.r / dist);
    v[0] = rotate(u, ang); v[1] = rotate(u, -ang); return 2;
}
double line_angle_degree(Vector v) {
    double ang = angle(v) * 180.0 / PI;
    while (dcmp(ang) < 0) { ang += 360.0; } while (dcmp(ang - 180.0) >= 0) { ang -= 180.0; }
    return ang;
}
Line move(Line l, double d) {
    Vector v = l.b - l.a; double len = length(v); l.a = l.a + Vector{-v.y / len, v.x / len} * d;
    return Line{l.a, l.a + v};
}

int main() {
#ifdef __AiR_H
    freopen("in.txt", "r", stdin);
//    freopen("out.txt", "w", stdout);
#endif // __AiR_H
    double xc, yc, xp, yp, r1, r2, r; Point a, b, c; Circle ans_c;
    while (scanf("%s", str) != EOF) {
        if (strcmp(str, "CircumscribedCircle") == 0) {
            a = read_point(); b = read_point(); c = read_point(); ans_c = circum_circle(a, b, c);
            printf("(%.6f,%.6f,%.6f)\n", ans_c.c.x, ans_c.c.y, ans_c.r);
        } else if (strcmp(str, "InscribedCircle") == 0) {
            a = read_point(); b = read_point(); c = read_point(); ans_c = inscribed_circle(a, b, c);
            printf("(%.6f,%.6f,%.6f)\n", ans_c.c.x, ans_c.c.y, ans_c.r);
        } else if (strcmp(str, "TangentLineThroughPoint") == 0) {
            scanf("%lf %lf %lf %lf %lf", &xc, &yc, &r, &xp, &yp);
            Vector v[2]; int cnt = get_tangents(Point{xp, yp}, Circle{Point{xc, yc}, r}, v);
            if (cnt == 0) { printf("[]\n"); continue; }
            if (cnt == 1) { printf("[%.6f]\n", line_angle_degree(v[0])); continue; }
            double ans0 = line_angle_degree(v[0]), ans1 = line_angle_degree(v[1]);
            if (ans0 > ans1) { swap(ans0, ans1); } printf("[%.6f,%.6f]\n", ans0, ans1);
        } else if (strcmp(str, "CircleThroughAPointAndTangentToALineWithRadius") == 0) {
            scanf("%lf %lf", &xp, &yp); Line l = read_line(); scanf("%lf", &r);
            vector<Point> ans; double t1, t2;
            get_line_circle_its(move(l, -r), Circle{Point{xp, yp}, r}, t1, t2, ans);
            get_line_circle_its(move(l, r), Circle{Point{xp, yp}, r}, t1, t2, ans);
            sort(ans.begin(), ans.end()); printf("[");
            if ((int)ans.size() > 0) {
                printf("(%.6f,%.6f)", ans[0].x, ans[0].y);
                for (int i = 1; i < (int)ans.size(); ++i) { printf(",(%.6f,%.6f)", ans[i].x, ans[i].y);  }
            }
            printf("]\n");
        } else if (strcmp(str, "CircleTangentToTwoLinesWithRadius") == 0) {
            Line l1 = read_line(), l2 = read_line(); scanf("%lf", &r);
            vector<Point> ans;
            Line l1_1 = move(l1, -r), l1_2 = move(l1, r), l2_1 = move(l2, -r), l2_2 = move(l2, r);
            ans.push_back(line_intersection(l1_1.a, l1_1.b - l1_1.a, l2_1.a, l2_1.b - l2_1.a));
            ans.push_back(line_intersection(l1_1.a, l1_1.b - l1_1.a, l2_2.a, l2_2.b - l2_2.a));
            ans.push_back(line_intersection(l1_2.a, l1_2.b - l1_2.a, l2_2.a, l2_2.b - l2_2.a));
            ans.push_back(line_intersection(l1_2.a, l1_2.b - l1_2.a, l2_1.a, l2_1.b - l2_1.a));
            sort(ans.begin(), ans.end()); printf("[(%.6f,%.6f)", ans[0].x, ans[0].y);
            for (int i = 1; i < 4; ++i) { printf(",(%.6f,%.6f)", ans[i].x, ans[i].y); }
            printf("]\n");
        } else {
            a = read_point(); scanf("%lf", &r1); b = read_point(); scanf("%lf", &r2);
            scanf("%lf", &r); vector<Point> ans;
            int cnt = get_circle_circle_its(Circle{a, r1 + r}, Circle{b, r2 + r}, ans);
            if (cnt == 0) { printf("[]\n"); continue; }
            if (cnt == 1) { printf("[(%.6f,%.6f)]\n", ans[0].x, ans[0].y); continue; }
            sort(ans.begin(), ans.end());
            printf("[(%.6f,%.6f),(%.6f,%.6f)]\n", ans[0].x, ans[0].y, ans[1].x, ans[1].y);
        }
    }
#ifdef __AiR_H
    printf("Time used = %.2fs\n", (double)clock() / CLOCKS_PER_SEC);
#endif // __AiR_H
    return 0;
}