HDU4697 Convex hull

X姐说没那么难写。

但是A完这题，我觉得我可以去死了= =

该死的精度问题。

基本解题思路是这样的，对任意三点共线的情况求出其发生的时间t，将其标记为一个事件，可以证明，只有发生事件的时候才有可能使凸包的形状发生变化。

对所有事件按照时间排序，求出所有相邻事件间的凸包形状（在凸包上的点集），并根据该点集和凸包面积的积分求出答案。

本题的难点在于：

1. 精度问题（1e-9丧心病狂！！！double和int的转换，写template有优势）

2. 解方程问题（学习了套公式求参数）

3. 求积分问题（也是套公式求参数，这个求参数的函数真有用）

#pragma comment(linker,"/STACK:102400000,102400000")
#include <cstdio>
#include <iostream>
#include <cstring>
#include <cmath>
#include <ctime>
#include <vector>
#include <set>
#include <queue>
#include <stack>
#include <map>
#include <algorithm>
using namespace std;
typedef long long ll;

const int maxn = 60;
const double eps = 1e-9;

template <class Type>
struct Point
{
	Type x, y;
	int id;
	Point(Type x=0, Type y=0): x(x), y(y) {}
};

template <class Type>
bool operator < (const Point<Type> &a, const Point<Type> &b) {return a.x < b.x || (a.x == b.x && a.y <b.y);}
int dcmp(double x) {if (fabs(x) < eps) return 0; else return x < 0 ? -1 : 1;}

template <class Type>
Point<Type> operator + (Point<Type> A, Point<Type> B) {return Point<Type>(A.x+B.x, A.y+B.y);}
template <class Type>
Point<Type> operator - (Point<Type> A, Point<Type> B) {return Point<Type>(A.x-B.x, A.y-B.y);}
template <class Type>
Point<Type> operator * (Point<Type> A, Type p) {return Point<Type>(A.x*p, A.y*p);}
template <class Type>
Type Dot(Point<Type> A, Point<Type> B) {return A.x*B.x + A.y*B.y;}
template <class Type>
Type Cross(Point<Type> A, Point<Type> B) {return A.x*B.y - A.y*B.x;}

template <class Type>
int ConvexHull(Point<Type> *A, int n, Point<Type> *ch)
{
	sort(A, A+n);
	int m = 0;
	for (int i = 0; i < n; i++)
	{
		while (m > 1 && dcmp(Cross(ch[m-1]-ch[m-2], A[i]-ch[m-1])) <= 0) m--;
		ch[m++] = A[i];
	}
	int k = m;
	for (int i = n-2; i >= 0; i--)
	{
		while (m > k && dcmp(Cross(ch[m-1]-ch[m-2], A[i]-ch[m-1])) <= 0) m--;
		ch[m++] = A[i];
	}
	return m;
}

int n, T, m;
Point<int> p[maxn], v[maxn];
Point<double> q[maxn], ch[maxn];
double t[250000];
int tot;
double ans;

template <class Type>
void get_coef(Type &a, Type &b, Type &c, Type d, Type k_1, Type b_1, Type k_2, Type b_2) // d * (k_1 x + b_1) * (k_2 x + b_2)
{
    a += d * (k_1 * k_2);
    b += d * (k_1 * b_2 + k_2 * b_1);
    c += d * b_1 * b_2;
}

inline void solve(int a, int b, int c)
{
    if (a == 0)
    {
        if (b != 0) t[tot++] = -1.0 * c / b;
        return;
    }

    int delta = b*b-4*a*c;
    if (delta < 0) return;
    if (delta == 0)
    {
        t[tot++] = 0.5 * (-b) / a;
        return;
    }
    double tmp = sqrt(1.0*delta);
    t[tot++] = 0.5 * (-b - tmp) / a;
    t[tot++] = 0.5 * (-b + tmp) / a;
    return;
}

inline double integral(int i, int j, double t1, double t2)
{
    int a = 0, b = 0, c = 0;
    get_coef(a, b, c, 1, v[i].x, p[i].x, v[j].y, p[j].y);
    get_coef(a, b, c, -1, v[j].x, p[j].x, v[i].y, p[i].y);
    return a * (t2*t2*t2-t1*t1*t1)/3.0 + b * (t2+t1)*(t2-t1)/2.0 + c * (t2-t1);
}

int main()
{
    while (scanf("%d%d", &n, &T) == 2)
    {
        for (int i=0;i<n;i++)
        {
            scanf("%d%d%d%d", &p[i].x, &p[i].y, &v[i].x, &v[i].y);
            p[i].id = i;
        }

        if (n <= 2) {printf("%.10f\n", 0); continue;}

        t[0] = 0; t[1] = T;
        tot = 2;
        for (int i=0;i<n;i++)
            for (int j=i+1;j<n;j++)
                for (int k=j+1;k<n;k++)
                {
                    int a = 0, b = 0, c = 0;
                    get_coef(a, b, c, 1, v[i].x-v[j].x, p[i].x-p[j].x, v[i].y-v[k].y, p[i].y-p[k].y);
                    get_coef(a, b, c, -1, v[i].x-v[k].x, p[i].x-p[k].x, v[i].y-v[j].y, p[i].y-p[j].y);
                    solve(a, b, c);
                }

        sort(t, t+tot);
        ans = 0;
        for (int k=0;k<tot-1;k++)
        {
            double mt = 0.5 * (t[k] + t[k+1]);
            if (mt < 0 || mt > T) continue;
            for (int i=0;i<n;i++)
            {
                q[i] = Point<double>(p[i].x + v[i].x * mt, p[i].y + v[i].y * mt);
                q[i].id = p[i].id;
            }
            m = ConvexHull(q, n, ch);
            for (int i=1;i<m;i++) ans += integral(ch[i-1].id, ch[i].id, t[k], t[k+1]);
        }
        printf("%.10f\n", abs(ans * 0.5 / T));
    }
    return 0;
}