HDU4752 Polygon

点击打开链接

题目大意是给定一个多边形，一条抛物线，求抛物线在多边形里面部分的长度的总和。

基本思想：对抛物线和多边形求交点，求出来交点之后判断那几段在多边形内，然后对那几段进行曲线积分，然后加和即可。

我的最初思想是求出来所有的交点之后，然后判断两个交点的中点是不是在多边形内，然后确定区间之后直接曲线积分，但是如果这么做的话，如果多边形的边在抛物线上面交错相交，时间复杂度为O(nn)，稳稳地不可以，于是乎这道题被卡住了……

然后几何大神适牛大大放出了这道题的解题报告，我便去拜读了一下，发现不止姿势优美，而且非常好写，尤其是判断抛物线恰好落在多边形顶点上面的时候，那个交点是否选择，他的判断方法真的是神思想……详情参见适牛大大的解题报告

至于code……几乎是仿写了适牛的写法……TAT

（同时学会了自适应辛普森积分法……）

#include <cstdio>#include <cstring>#include <algorithm>#include <cmath>#include <cstdlib>#include <vector>using namespace std;const double eps = 1e-8;int dblcmp(double x){if (fabs(x) < eps) return 0;return x > eps ? 1 : -1;}const double pi = acos(-1.0);inline double sqr(double x){return x * x;}struct Point2D{double x, y;Point2D(){}Point2D(double a, double b) : x(a), y(b){}void input(){scanf("%lf%lf", &x, &y);}void output(){printf("%.2lf %.2lf", x, y);}friend Point2D operator + (const Point2D &a, const Point2D &b){return Point2D(a.x + b.x, a.y + b.y);}friend Point2D operator - (const Point2D &a, const Point2D &b){return Point2D(a.x - b.x, a.y - b.y);}friend bool operator == (const Point2D &a, const Point2D &b){return dblcmp(a.x - b.x) == 0 && dblcmp(a.y - b.y) == 0;}friend Point2D operator * (const double &a, const Point2D &b){return Point2D(a * b.x, a * b.y);}friend Point2D operator / (const Point2D &a, const double &b){return Point2D(a.x / b, a.y / b);}double norm(){return hypot(x, y);}Point2D trunc(double r){double l = norm();if (!dblcmp(l)) return *this;r /= l;return Point2D(x * r, y * r);}};double det(const Point2D &a, const Point2D &b){return a.x * b.y - a.y * b.x;}double dot(const Point2D &a, const Point2D &b){return a.x * b.x + a.y * b.y;}double dist(const Point2D &a, const Point2D &b){return (a - b).norm();}Point2D rotate_point(const Point2D &p, double A){double tx = p.x, ty = p.y;return Point2D(tx * cos(A) - ty * sin(A), tx * sin(A) + ty * cos(A));}double a, b, c, l, r;double f0(double x){return sqrt(1 + sqr(2 * a * x + b));}double adaptive_simpsons_aux(double (*f)(double), double a, double b, double eps, double s, double fa, double fb, double fc, int depth){double c = (a + b) / 2, h = b - a;double d = (a + c) / 2, e = (c + b) / 2;double fd = f(d), fe = f(e);double sl = (fa + 4 * fd + fc) * h / 12;double sr = (fc + 4 * fe + fb) * h / 12;double s2 = sl + sr;if (depth <= 0 || fabs(s2 - s) <= 15 * eps)return s2 + (s2 - s) / 15;return adaptive_simpsons_aux(f, a, c, eps / 2, sl, fa, fc, fd, depth - 1) +    adaptive_simpsons_aux(f, c, b, eps / 2, sr, fc, fb, fe, depth - 1);}double adaptive_simpsons(double (*f)(double), double a, double b, double eps, int depth){double c = (a + b) / 2, h = b - a;double fa = f(a), fb = f(b), fc = f(c);double s = (fa + 4 * fc + fb) * h / 6;return adaptive_simpsons_aux(f, a, b, eps, s, fa, fb, fc, depth);}double f(double x){return a * x * x + b * x + c;}int n;const int maxn = 2e4 + 10;Point2D p[maxn];vector<double> h;bool inSeg(Point2D l, Point2D r, Point2D mid){return dblcmp(mid.x - fmin(l.x, r.x)) >= 0&& dblcmp(fmax(l.x, r.x) - mid.x) >= 0&& dblcmp(mid.y - fmin(l.y, r.y)) >= 0&& dblcmp(fmax(l.y, r.y) - mid.y) >= 0;}bool check(Point2D pre, Point2D now, Point2D nxt){pre = pre - now;nxt = nxt - now;pre = pre.trunc(0.1);nxt = nxt.trunc(0.1);pre = pre + now;nxt = nxt + now;if (dblcmp(pre.y - f(pre.x)) * dblcmp(nxt.y - f(nxt.x)) < 0) return true;return false;}void solve(){scanf("%lf%lf%lf%lf%lf", &a, &b, &c, &l, &r);for (int i = 0; i < n; i++)p[i].input();h.resize(0);if (n <= 2){puts("0.00");return;}for (int i = 0; i < n; ++i){Point2D now, pre, nxt, Inter, Inter2;pre = p[(i + n - 1) % n];now = p[i];nxt = p[(i + 1) % n];if (dblcmp(now.x - nxt.x) == 0){Inter = Point2D(now.x, f(now.x));if (inSeg(now, nxt, Inter)){h.push_back(Inter.x);if (Inter == now)if (!check(pre, now, nxt))h.push_back(Inter.x);}}else{double kk = (nxt.y - now.y) / (nxt.x - now.x);double bb = nxt.y - nxt.x * kk;double A = sqr(b - kk) - 4.0 * (c - bb) * a;if (dblcmp(A) > 0){double x1 = (-(b - kk) + sqrt(A)) / 2.0 / a;double x2 = (-(b - kk) - sqrt(A)) / 2.0 / a;Inter = Point2D(x1, f(x1));Inter2 = Point2D(x2, f(x2));if (inSeg(now, nxt, Inter)){if (!(Inter == pre) && !(Inter == nxt)) h.push_back(Inter.x);if (Inter == now && !check(pre, now, nxt)) h.push_back(Inter.x);}if (inSeg(now, nxt, Inter2)){if (!(Inter2 == pre) && !(Inter2 == nxt)) h.push_back(Inter2.x);if (Inter2 == now && !check(pre, now, nxt)) h.push_back(Inter2.x);}}}}sort(h.begin(), h.end());double ans = 0;for (int i = 0; i + 1 < h.size(); i += 2){double L = fmax(l, h[i]), R = fmin(r, h[i + 1]);if (dblcmp(R - L) > 0)ans += adaptive_simpsons(f0, L, R, eps, 10);if (dblcmp(R - L) < 0) break;}double out = ans;printf("%.2lf\n", out + eps);}int main(){#ifndef ONLINE_JUDGEfreopen("in.txt", "rt", stdin);#endifwhile(scanf("%d", &n) == 1)solve();return 0;}