hdu4679 Convex Hull 计算几何

50个点，每个点有一个恒定速度，给定一个时刻T

求，其中 S(t) 为所有点在t时刻组成的凸包面积。

找到所有特殊的时间点： 只有在存在三点共线的时候凸包上的点才会改变。

如果凸包上的点保持不变，其面积随时间变化是一个一元二次函数（考虑相邻两点的位置向量是一个关于时间的一次函数，叉乘后得到一个一元二次的式子）。

所以求出a，b，c后求积分即可。求a，b，c时可以取开始时间点，结束时间点，再取一个中间点。这个过程在网上看的代码，比较巧妙。。

精度问题，除以2 用乘以0.5代替

#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <cassert>
#include <ctime>
#include <algorithm>
#include <iostream>
#include <vector>
#define pb push_back
#define sqr(x) ((x)*(x))
#define cube(x) ((x)*(x)*(x))
using namespace std;
const double eps = 1e-9;
const int maxN = 50;
int n, T;
vector<double> t;

int dcmp(double x){
    if (fabs(x)<eps) return 0;
    return x>0?1:-1;
}
struct Point{
    double x, y;
    int id;
    Point(){}
    Point(double xx, double yy):x(xx), y(yy){}
    void input(){ scanf("%lf%lf", &x, &y); }
    friend Point operator + (const Point& a, const Point& b){
        return Point(a.x+b.x, a.y+b.y);
    }
    friend Point operator - (const Point& a, const Point& b){
        return Point(a.x-b.x, a.y-b.y);
    }
    friend double operator * (const Point& a, const Point& b){
        return a.x*b.y - a.y+b.x;
    }
    friend Point operator * (const Point& a, const double& k){
        return Point(a.x*k, a.y*k);
    }
    friend bool operator < (const Point& a, const Point& b){
        return a.x<b.x || a.x==b.x&&a.y<b.y;
    }
}p[maxN+2], q[maxN+2], v[maxN+2], ch[maxN+2];
int det(const Point& a, const Point& b){
    return dcmp(a.x*b.y - a.y*b.x);
}
int Hull(int n, Point p[], Point ch[]){
    sort(p+1, p+1+n);
    int m=0;
    for (int i=1; i<=n; i++){
        while (m>1 && det(ch[m]-ch[m-1], p[i]-ch[m-1])<=0) m--;
        ch[++m] = p[i];
    }
    int pre=m;
    for (int i=n; i; i--){
        while (m>pre && det(ch[m]-ch[m-1], p[i]-ch[m-1])<=0) m--;
        ch[++m] = p[i];
    }
    return m;
}
void get_abc(double& a, double& b, double& c, int i, int j, int k){
    Point p1=p[j]-p[i], v1=v[j]-v[i];
    Point p2=p[k]-p[i], v2=v[k]-v[i];
    a = v1*v2;
    b = p1*v2 + p2*v1;
    c = p1*p2;
}
void solve(double a, double b, double c){
    double x;
    if (a==0){
        if (dcmp(b)){
            x = -c/b;
            if (x>=0 && x<=T) t.pb(x);
        }
        return;
    }
    double delta = b*b - 4*a*c;
    if (delta<0) return;
    if (!dcmp(delta)){
        x = -b*0.5/a;
        if (x>=0 && x<=T) t.pb(x);
        return;
    }
    double tmp = sqrt(delta);
    x = (-b+tmp)*0.5/a;
    if (x>=0 && x<=T) t.pb(x);
    x = (-b-tmp)*0.5/a;
    if (x>=0 && x<=T) t.pb(x);
}

int main(){
    while (scanf("%d%d", &n, &T) != EOF){
        for (int i=1; i<=n; i++){
            p[i].input();
            v[i].input();
        }
        t.clear();
        t.pb(0); t.pb(T);
        for (int i=1; i<=n; i++)
            for (int j=i+1; j<=n; j++)
                for (int k=j+1; k<=n; k++){
                    double a, b, c;
                    get_abc(a, b, c, i, j, k);
                    solve(a, b, c);
                }
        int tot=t.size();
        sort(t.begin(), t.end());
        double ans=0;
        for (int i=0; i<tot-1; i++){
            double mt = (t[i]+t[i+1])*0.5;
            for (int i=1; i<=n; i++){
                q[i] = p[i] + v[i]*mt;
                q[i].id = i;
            }
            int m = Hull(n, q, ch);
            for (int i=1; i<=m; i++)
                printf("\t%d\n", ch[i].id);
            for (int i=1; i<m; i++){
                double a, b, c;
                get_abc(a, b, c, ch[1].id, ch[i].id, ch[i+1].id);
                ans += a/3*(cube(t[i+1]) - cube(t[i])) + b*0.5*(sqr(t[i+1]) - sqr(t[i])) + c*(t[i+1]-t[i]);
            }
        }
        printf("%.10f\n", ans/T*0.5);
    }
    return 0;
}