题解:P4406 [CQOI2005] 三角形面积并
你说得对,但让我们先抛开这题不谈,来建立一个平面直角坐标系,让我们已知两个点的坐标,求出这两个点构成直线的横截式(如果你已经会了,就略过这里)。
设两个点的坐标分别是 \(\left ( x_1,y_1 \right )\)、\(\left ( x_2,y_2 \right )\),因为我们要求的是直线的横截式,所以它一定能被表示为 \(y = ax + b\)。我们可以根据点的坐标列出方程组:
你的消元法确实很高级,但是我更想要用克莱姆法则,可得:
综上,横截式就是 \(y=\frac{ y_1-y_2}{ x_1-x_2}x+\frac{x_1y_2-x_2y_1}{ x_1-x_2}\)。
但是,等一下。克莱姆法则说明了:当 \(D\) 等于 \(0\) 时,这个方程无解或无数解。那时候怎么作横截式呢?
我们思考一下就可以知道:此时做不出横截式。好吧就到这里,没有就没有,不影响我们计算。
现在回归这题。这题要求我们求三角形面积并。想到面积,我们不由得想到一个东西——积分。
……直观地说,对于一个给定的正实值函数,在一个实数区间上的定积分可以理解为在坐标平面上,由曲线、直线以及轴围成的曲边梯形的面积值(一种确定的实数值)(摘自百度百科)
在这题,我们可以使用定积分求解。说到求解积分,你想到了什么?对!自适应辛普森法!不会了话建议左转进入 P4526 学习包教包会。
那么我们要求的就是 \(\int_{-\infty}^{\infty} f(x)\mathrm{d}x\),因为本题的数据范围,这个积分又变成了 \(\int_{-10^6}^{10^6} f(x)\mathrm{d}x\)。
我们是不是可以聪明一点呢?
我们将所有点的的 \(x\) 坐标存进数组 point 中,让它从小到大排序。
然后呢?使用聪明的大脑我们可以发现,答案就是 \(\sum_{i=1}^{3n-1} \int_{point_{i-1}}^{point_{i}} f(x)\mathrm{d}x\)。
现在我们开考虑 \(f(x)\) 的计算,这个确实挺恶心的。
我们分别考虑每一个三角形:

那我问你,遇到这种普通情况,怎么求出 \(y_1\) 和 \(y_2\)?看看上面。这不是随随便便就可以用横截式一下子算出来。
看看特殊情况:

那我再问你,这个怎么算。有的人蒙了,但是我们通过 AB 和 AC 的横截式算出 B 和 C 的 \(y\) 坐标后,还考虑什么 BC 的横截式呢。
所以,我们可以实现代码了(这个代码的实现真的很折磨人,请在手打的时候做好心理准备)。
AC 代码:
//
#include <bits/stdc++.h>
#define int long long
using namespace std;
const int N = 100 + 5;
const double eps = 1e-9;
int n;
struct node {
double x, y;
}a[N], b[N], c[N];
struct equ {
bool flag = false;
double a, b;
// y = ax + b;
}ab[N], bc[N], ac[N];
struct nod {
double fi;
int se;
};
inline double f(register double x) {
vector <nod> xp;
int sum = 0;
for (int i = 1; i <= n; i++) {
int tot = 0;
// AB
vector <double > tmp;
if (min(a[i].x, b[i].x) <= x && x <= max(a[i].x, b[i].x) && ab[i].flag) {
double y_AB = ab[i].a * x + ab[i].b;
tmp.push_back(y_AB);tot++;
}
// AC
if (min(a[i].x, c[i].x) <= x && x <= max(a[i].x, c[i].x) && ac[i].flag) {
double y_AC = ac[i].a * x + ac[i].b;
tmp.push_back(y_AC);tot++;
}
// BC
if (min(b[i].x, c[i].x) <= x && x <= max(b[i].x, c[i].x) && bc[i].flag) {
double y_BC = bc[i].a * x + bc[i].b;
tmp.push_back(y_BC);tot++;
}
if (tot == 2) {
xp.push_back({min(tmp[0], tmp[1]), 1});
xp.push_back({max(tmp[0], tmp[1]), -1});
}
}
sort(xp.begin(), xp.end(), [](nod a, nod b) {
if (a.fi != b.fi) return a.fi < b.fi;
return a.se > b.se;
});
register double ans = 0;
register int kl = 0;
register double la = xp[0].fi;
for (auto& k : xp) {
register double p = k.fi;
register int type = k.se;
if (kl > 0) {
ans += p - la;
}
kl += type;
la = p;
}
return ans;
}
inline double ss(register double l, register double r) {
register double mid = (l + r) / 2;
return (r - l) * (f(l) + 4 * f(mid) + f(r)) / 6;
}
inline double g(register double l, register double r, register double eps, register double ans, register int step) {
register double mid = (l + r) / 2;
register double tl = ss(l, mid), tr = ss(mid, r);
if (abs(tl + tr - ans) <= 15 * eps && step < 0) {
return tl + tr + (tl + tr - ans) / 15;
}
return g(l, mid, eps / 2, tl, step - 1) + g(mid, r, eps / 2, tr, step - 1);
}
signed main() {
cin >> n;
vector <double> point;
for (int i = 1; i <= n; i++) {
cin >> a[i].x >> a[i].y;
cin >> b[i].x >> b[i].y;
cin >> c[i].x >> c[i].y;
if(a[i].x != b[i].x) {
ab[i].flag = true;
ab[i].a = (a[i].y - b[i].y) / (a[i].x - b[i].x);
ab[i].b = (a[i].x * b[i].y - a[i].y * b[i].x) / (a[i].x - b[i].x);
}
if(a[i].x != c[i].x) {
ac[i].flag = true;
ac[i].a = (a[i].y - c[i].y) / (a[i].x - c[i].x);
ac[i].b = (a[i].x * c[i].y - a[i].y * c[i].x) / (a[i].x - c[i].x);
}
if(b[i].x != c[i].x) {
bc[i].flag = true;
bc[i].a = (b[i].y - c[i].y) / (b[i].x - c[i].x);
bc[i].b = (b[i].x * c[i].y - b[i].y * c[i].x) / (b[i].x - c[i].x);
}
point.push_back(min({a[i].x, b[i].x, c[i].x}));
point.push_back(max({a[i].x, b[i].x, c[i].x}));
}
sort(point.begin(), point.end());
register double ans = 0;
register double last = point[0];
for (int i = 1; i < 2 * n; i++) {
ans += g(last, point[i], eps, ss(last, point[i]), 5);
last = point[i];
}
cout << fixed << setprecision(2) << ans;
return 0;
}

浙公网安备 33010602011771号