给出 $n$ 个三角形，求它们并的面积。
第一行包含一个正整数 $n$ ($n \leq 100$)，表示三角形的个数。
接下来的 $n$ 行，每行六个实数 $x_1, y_1, x_2, y_2, x_3, y_3$ ($|x_i|, |y_i| \leq 10^6$)，表示三角形的顶点坐标，输入保留 $1$ 位小数。
输出一行一个实数 $U$，表示三角形并的面积，保留两位小数。
说到面积，很多人都可以想到 Simpson 积分，它是曲边积分的一种近似。但是这道题中所有的边都是直线 (线段)，因此存在精确的算法，其实用积分计算面积也是可以做到精确的。
注意到，对于一个很小的区间 $U = \left[ x - \delta, x + \delta \right]$，用 (平行于 $y$ 轴的) 扫描线去扫 (横坐标在) 区间 $U$ 上的区域，所截得的线段长度大多数都是一个线性函数。只有在 $U$ 中包含某两条线的交点的时候，才会出现分段。
因此，我们可以将 $n$ 个三角形的 $3n$ 条线两两相交的交点处理出来并按 $x$ 轴排序，设为 $x_1, x_2, \cdots, x_k$。则对于每个区间 $U_i = \left[ x_i, x_{i+1} \right]$，截得的线段长是一个线性函数，或者说它是若干个梯形的不交并。
由于梯形的面积等于中位线 × 高，因此，我们只需用 $x = \dfrac {x_i + x_{i+1}} 2$ 的扫描线去扫描整个图形，用所截得的线段长乘以 $x_{i+1} - x_i$ (高) 即可得到这部分的面积。
将所有部分的面积相加即可得到答案。由于线段求并是 $O \left( n \log n \right)$ 的，因此总复杂度为 $O \left( n^3 \log n \right)$。
#include <bits/stdc++.h>
const int N = 1000;
const double eps = 1e-8;
#define lt(x, y) ((x) < (y) - eps)
#define gt(x, y) ((x) > (y) + eps)
#define le(x, y) ((x) <= (y) + eps)
#define ge(x, y) ((x) >= (y) - eps)
#define eq(x, y) (le(x, y) && ge(x, y))
#define dot(x, y, z) (((y) - (x)) * ((z) - (x)))
#define cross(x, y, z) (((y) - (x)) ^ ((z) - (x)))
#define def(g) inline bool f##g(const double &x, const double &y) {return g(x, y);}
def(eq)
struct vec2 {
double x, y;
vec2 (double x0 = 0.0, double y0 = 0.0) : x(x0), y(y0) {}
vec2 * read() {scanf("%lf%lf", &x, &y); return this;}
inline vec2 operator - () const {return vec2(-x, -y);}
inline vec2 operator + (const vec2 &B) const {return vec2(x + B.x, y + B.y);}
inline vec2 operator - (const vec2 &B) const {return vec2(x - B.x, y - B.y);}
inline vec2 operator * (double k) const {return vec2(x * k, y * k);}
inline vec2 operator / (double k) const {return *this * (1.0 / k);}
inline double operator * (const vec2 &B) const {return x * B.x + y * B.y;}
inline double operator ^ (const vec2 &B) const {return x * B.y - y * B.x;}
inline double norm2() const {return x * x + y * y;}
inline double norm() const {return sqrt(x * x + y * y);}
inline bool operator < (const vec2 &B) const {return lt(x, B.x) || le(x, B.x) && lt(y, B.y);}
inline bool operator == (const vec2 &B) const {return eq(x, B.x) && eq(y, B.y);}
inline bool operator << (const vec2 &B) const {return lt(y, 0) ^ lt(B.y, 0) ? lt(B.y, 0) : gt(*this ^ B, 0) || ge(*this ^ B, 0) && ge(x, 0) && lt(B.x, 0);}
inline vec2 trans(double a11, double a12, double a21, double a22) const {return vec2(x * a11 + y * a12, x * a21 + y * a22);}
};
struct line {
double A, B, C;
line (double A0 = 0.0, double B0 = 0.0, double C0 = 0.0) : A(A0), B(B0), C(C0) {}
line (const vec2 &u, const vec2 &v) : A(u.y - v.y), B(v.x - u.x), C(u ^ v) {}
inline vec2 normVec() const {return vec2(A, B);}
inline double norm2() const {return A * A + B * B;}
inline double operator () (const vec2 &P) const {return A * P.x + B * P.y + C;}
} l[N];
inline vec2 intersection(const line u, const line v) {
double Det = 1.0 / (u.A * v.B - u.B * v.A);
return vec2(u.B * v.C - u.C * v.B, u.C * v.A - u.A * v.C) * Det;
}
// ---------------- templates end ---------------- //
struct triangle {
vec2 A[3];
triangle * read() {return A[0].read(), A[1].read(), A[2].read(), std::sort(A, A + 3), this;}
} t[N];
typedef std::pair <double, bool> pdb;
int n;
double x[N * N];
pdb p[N];
inline double proj(double x, const vec2 &A, const vec2 &B) {return ((A ^ B) + x * (A.y - B.y)) / (A.x - B.x);}
bool intersect(double x, const triangle &t, pdb *r1, pdb *r2) {
if (le(x, t.A[0].x) || ge(x, t.A[2].x)) return false;
r1->first = proj(x, t.A[0], t.A[2]);
r2->first = proj(x, t.A[le(x, t.A[1].x) ? 0 : 2], t.A[1]);
if (r1->first > r2->first) std::swap(r1->first, r2->first);
return r2->second = false, r1->second = true;
}
double compute(double x) {
int i, cur = 0, tot = 0; double ret = 0.0; pdb g1, g2;
for (i = 0; i < n; ++i)
if (intersect(x, t[i], &g1, &g2)) p[tot++] = g1, p[tot++] = g2;
std::sort(p, p + tot);
for (i = 0; i < tot; p[i++].second ? ++cur : --cur)
ret += cur ? p[i].first - p[i - 1].first : 0.0;
return ret;
}
int main() {
int i, j, cnt = 0, tot = 0; double ans = 0.0;
scanf("%d", &n);
for (i = 0; i < n; ++i)
t[i].read(), l[tot++] = line(t[i].A[0], t[i].A[1]), l[tot++] = line(t[i].A[1], t[i].A[2]), l[tot++] = line(t[i].A[2], t[i].A[0]);
for (i = 0; i < tot; ++i)
for (j = i + 1; j < tot; ++j)
x[cnt++] = intersection(l[i], l[j]).x, std::isfinite(x[cnt - 1]) ? 0 : --cnt;
std::sort(x, x + cnt), cnt = std::unique(x, x + cnt, feq) - x;
for (i = 1; i < cnt; ++i) ans += compute((x[i - 1] + x[i]) * 0.5) * (x[i] - x[i - 1]);
printf("%.2lf\n", ans);
return 0;
}
坑1：注意 lydsy 中需要减 eps 后再输出，不知道是为什么。(可能是它那边是五舍六入？)
坑2：用扫描线的时候注意要对一开始的横坐标去重，去重时注意 eps。