计算几何总结

几何总结:
int sgn(double x)
{
    if(fabs(x)<eps) return 0;
    if(x<0) return -1;
    else return 1;
}
struct Point
{
    double x,y;
    Point(){}
    Point(double _x,double _y){x=_x,y=_y;}
    Point operator-(const Point &b)const
    {
        return Point(x-b.x,y-b.y);
    }
    double operator^(const Point &b)const//叉乘
    {
        return x*b.y-y*b.x;
    }
    double operator*(const Point &b)const//点积
    {
        return x*b.x+y*b.y;
    }
    void input()
    {
        scanf("%lf%lf",&x,&y);
    }
};

struct Line
{
    Point s,e;
    Line(){}
    Line(Point _s,Point _e){s=_s,e=_e;}
};

int cross(Point a,Point b,Point c)
{
    return (b-a)^(c-b);
}
int dis(Point a,Point b)
{
    return (a-b)*(a-b);
}
double gcd(double x,double y)
{
    while(fabs(x)>eps && fabs(y)>eps)
    {
        if(x>y) x-=floor(x/y)*y;
        else y-=floor(y/x)*x;
    }
    return x+y;
}
首先在二维坐标下介绍一些定义：
点：A（x1,y1）,B（x2,y2）
向量：向量AB=（x2-x1,y2-y1）= (x,y);
向量的模|AB|=sqrt(x*x+y*y);

向量的点积：x1*x2 + y1*y2。
点积的结果是一个数值。
应用：可以根据集合意义求两向量的夹角
cosθ=（ 向量a*向量b ）/（|a|*|b|）= (x1*x2+y1*y2)/（|a|*|b|）

向量的叉积：x1*y2 - x2*y1
叉积的结果也是一个向量，是垂直于向量a，b所形成的平面，如果看成三维坐标的话是在 z 轴上，上面结果是它的模。
叉积的集合意义：
|向量a x 向量b|=|向量a|*|向量b|*sinθ
1：其结果是a和b为相邻边形成平行四边形的面积。
2：结果有正有负，有sinθ可知和其夹角有关，夹角大于180°为负值。
3：叉积不满足交换律
应用：
1：通过结果的正负判断两矢量之间的顺逆时针关系
若 a x b > 0表示a在b的顺时针方向上
若 a x b < 0表示a在b的逆时针方向上
若 a x b == 0表示a在b共线，但不确定方向是否相同
2：判断折线拐向，可转化为判断第三点在前两的形成直线的顺逆时针方向，然后判断拐向。
3：判断一个点在一条直线的那一侧，同样上面的方法。
4：判断点是否在线段上，可利用叉乘首先判断是否共线，然后在判断是否在其上。
5：判断两条直线是否相交（快速排斥实验和跨立实验）
根据判断点在直线那一侧我们可以判断一个线段的上的两点分别在另一个线段的两侧，当然这是不够的，因为我们画图发现这样只能够让直线想交，而不是线段，所以我们还要对另一条线段也进行相同的判断就ok。

//HDU-5563 
题目大意：给5个点，判断能不能组成一个正五角星。
思路：其实就是判正五边形。在正五边形中任意两点之间的距离只有两种情况
      所以求出所有点之间的距离排序，前后比较记录有多少个不同距离的边长。

//POJ-1269///
题意大意：判断直线ab与直线cd的关系，如果相交就输出交点坐标
void solve(Point a,Point b,Point c,Point d)
{
    if(fabs(cross(a,b,c))<=eps&&fabs(cross(a,b,d))<=eps)puts("LINE");//共线
    else if((b.x-a.x)*(d.y-c.y)==(d.x-c.x)*(b.y-a.y))puts("NONE");//平行
    else{
        double a1=a.y-b.y;
        double b1=b.x-a.x;
        double c1=a.x*b.y-b.x*a.y;

        double a2=c.y-d.y;
        double b2=d.x-c.x;
        double c2=c.x*d.y-d.x*c.y;

        double x=(b1*c2-b2*c1)/(a1*b2-a2*b1);
        double y=(a2*c1-a1*c2)/(a1*b2-a2*b1);
        printf("POINT %.2f %.2f\n",x,y);//相交及交点坐标
    }
}

//POJ-3304
题意：求是否存在一条直线，使所有线段到这条直线的投影至少有一个交点。
思路：把问题转化为是否存在一条直线与每条线段都有交点
int sgn(double x);
int n;
Point p[maxn];
Line line[maxn];
bool Seg_inter_line(Point a,Point b,Line l2) //判断直线ab和线段l2是否相交 相交返回true否则返回false
{
    return sgn((l2.a-b)^(a-b))*sgn((l2.b-b)^(a-b)) <= 0;
}
/*
bool Seg_inter_line(Line l1,Line l2) //判断直线l1和线段l2是否相交
{
    return sgn((l2.s-l1.e)^(l1.s-l1.e))*sgn((l2.e-l1.e)^(l1.s-l1.e)) <= 0;
}
*/
bool check(Point a,Point b)
{
    if(fabs(a.x-b.x)<=eps&&fabs(a.y-b.y)<=eps) return false;
    for(int i=0;i<n;i++)
        //if(cross(line[i].a,a,b)*cross(line[i].b,a,b)>eps) return false;
        if(!Seg_inter_line(a,b,line[i])) return false;
    return true;
}
int main(void)
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d",&n);
        for(int i=0;i<n;i++)
            scanf("%lf%lf%lf%lf",&line[i].a.x,&line[i].a.y,&line[i].b.x,&line[i].b.y);
        int flag=0;
        if(n<3) flag=1;
        for(int i=0;i<n;i++)
        {
            if(check(line[i].a,line[i].b)) flag=1;
            for(int j=i+1;j<n;j++)
            {
                if(check(line[i].a,line[j].a)||
                   check(line[i].a,line[j].b)||
                   check(line[i].b,line[j].a)||
                   check(line[i].b,line[j].b))
                   {flag=1;break;}
            }
            if(flag) break;
        }
        if(flag) printf("Yes!\n");
        else printf("No!\n");
    }
    return 0;
}

//HDU - 5533 
题意：给你几个点，问能否形成正多边形。
思路：①找出所有点之间的最短距离，即边长。然后再在所有点中找出与边长相等的，只要数量为n即说明能够成正多边形
      ②要求输入的为整数，所以只有正四边形才能成立，判断下即可

//H-Intersection POJ-1410/
题意:给了一个线段和矩形。
如果线段和矩形的边相交，或者线段在矩形内。输出'T',否则输出'F'

bool inter(Line l1,Line l2)//判断线段相交
{
    return
    max(l1.s.x,l1.e.x) >= min(l2.s.x,l2.e.x) &&
    max(l2.s.x,l2.e.x) >= min(l1.s.x,l1.e.x) &&
    max(l1.s.y,l1.e.y) >= min(l2.s.y,l2.e.y) &&
    max(l2.s.y,l2.e.y) >= min(l1.s.y,l1.e.y) &&
    sgn((l2.s-l1.e)^(l1.s-l1.e))*sgn((l2.e-l1.e)^(l1.s-l1.e)) <= 0 &&
    sgn((l1.s-l2.e)^(l2.s-l2.e))*sgn((l1.e-l2.e)^(l2.s-l2.e)) <= 0;
}

bool OnSeg(Point P,Line L)//判断点在线段上
{
    return
    sgn((L.s-P)^(L.e-P)) == 0 &&
    sgn((P.x - L.s.x) * (P.x - L.e.x)) <= 0 &&
    sgn((P.y - L.s.y) * (P.y - L.e.y)) <= 0;
}
//判断点在凸多边形内
//点形成一个凸包，而且按逆时针排序（如果是顺时针把里面的<0改为>0）
//点的编号:0~n-1
//返回值：
//-1:点在凸多边形外
//0:点在凸多边形边界上
//1:点在凸多边形内
int inConvexPoly(Point a,Point p[],int n)
{
    for(int i = 0;i < n;i++)
    {
        if(sgn((p[i]-a)^(p[(i+1)%n]-a)) < 0)return -1;
        else if(OnSeg(a,Line(p[i],p[(i+1)%n])))return 0;
    }
    return 1;
}
//判断点在任意多边形内
//射线法，poly[]的顶点数要大于等于3,点的编号0~n-1
//返回值
//-1:点在凸多边形外
//0:点在凸多边形边界上
//1:点在凸多边形内
int inPoly(Point p,Point poly[],int n)
{
    int cnt;
    Line ray,side;
    cnt = 0;
    ray.s = p;
    ray.e.y = p.y;
    ray.e.x = -100000000000.0;//-INF,注意取值防止越界
    for(int i = 0;i < n;i++)
    {
        side.s = poly[i];
        side.e = poly[(i+1)%n];

        if(OnSeg(p,side)) return 0;
        //如果平行轴则不考虑
        if(sgn(side.s.y - side.e.y) == 0) continue;
        if(OnSeg(side.s,ray))
        {
            if(sgn(side.s.y - side.e.y) > 0) cnt++;
        }
        else if(OnSeg(side.e,ray))
        {
            if(sgn(side.e.y - side.s.y) > 0) cnt++;
        }
        else if(inter(ray,side)) cnt++;
    }
    if(cnt % 2 == 1)return 1;
    else return -1;
}
int main()
{
    int T;
    double x1,y1,x2,y2;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lf%lf%lf%lf",&x1,&y1,&x2,&y2);
        Line line = Line(Point(x1,y1),Point(x2,y2));
        scanf("%lf%lf%lf%lf",&x1,&y1,&x2,&y2);
        if(x1 > x2)swap(x1,x2);
        if(y1 > y2)swap(y1,y2);
        Point p[10];
        p[0] = Point(x1,y1);
        p[1] = Point(x2,y1);
        p[2] = Point(x2,y2);
        p[3] = Point(x1,y2);
        if(inter(line,Line(p[0],p[1]))||inter(line,Line(p[1],p[2]))||
           inter(line,Line(p[2],p[3]))||inter(line,Line(p[3],p[0])))
            printf("T\n");
        else if(inConvexPoly(line.s,p,4) >= 0 || inConvexPoly(line.e,p,4) >= 0)
            printf("T\n");
        else printf("F\n");
    }
    return 0;
}

// HDU5120 Intersection 【求圆的面积交】
// 题意：求两个圆环的相交面积
// 思路：S = A大B大 - A大B小 - A小B大 + A小B小。（A表示A环，大表示大圆，B同）。
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
double work(Point c1,double r1,Point c2,double r2)//圆心c1半径r1于圆心c2半径r2的相交面积
{
    double d=dist(c1,c2);
    if(r1+r2<d+eps)return 0;
    if(d<fabs(r1-r2)+eps)
    {
        double r=min(r1,r2);
        return PI*r*r;
    }
    double x=(d*d+r1*r1-r2*r2)/(2*d);
    double t1=acos(x/r1);
    double t2=acos((d-x)/r2);
    return r1*r1*t1+r2*r2*t2-d*r1*sin(t1);
}
int main(void)
{
    int T;
    scanf("%d",&T);
    for(int cas=1;cas<=T;cas++)
    {
        double r,R,x1,y1,x2,y2;
        scanf("%lf%lf%lf%lf%lf%lf", &r, &R, &x1, &y1, &x2, &y2);
        Point c1=Point(x1,y1);
        Point c2=Point(x2,y2);
        double ans=work(c1,R,c2,R)-work(c1,R,c2,r)
                    -work(c1,r,c2,R)+work(c1,r,c2,r);
        printf("Case #%d: %.6lf\n",cas,ans);
    }
    return 0;
}

//hdu1411
//题意：求四面体的体积
//已知任意四面体D-ABC，记DA=a，DB=b，DC=c，cos角ADB=x，cos角BDC=y，cos角CDA=z，则：
//V=1/6*abc*sqrt(1+2xyz-x^2-y^2-z^2)
int main(void)
{
    double ab,ac,ad,bc,bd,cd;
    while(scanf("%lf%lf%lf%lf%lf%lf",&ab,&ac,&ad,&bc,&bd,&cd)!=EOF)//以d为顶点
    {
        double x=(ad*ad+bd*bd-ab*ab)/(2*ad*bd);
        double y=(bd*bd+cd*cd-bc*bc)/(2*bd*cd);
        double z=(ad*ad+cd*cd-ac*ac)/(2*ad*cd);
        double V=ad*bd*cd*sqrt(1+2*x*y*z-x*x-y*y-z*z)/6;
        printf("%.4f\n",V);
    }
    return 0;
}

题意：给出两条线段的两个端点。求一个圆，使其与两条线段都有且只有一个交点；
分析：找找出两个线段的两个端点的最小距离 ，这两个端点的中点为圆心，这个最小距离加上一个很小的值即为半径
//p[0].x,p[0].y为第一条线段的端点坐标 p[1].x,p[1].y为第一条线段的另一个端点坐标。
//p[2].x,p[2].y为第二条线段的端点坐标 p[3].x,p[3].y为第二条线段的另一个端点坐标。
double mi=INF;
int dian1=0,dian2=0;
for(int i=0;i<2;i++)
{
    for(int j=2;j<4;j++)
    {
        double temp=dist(p[i].x,p[i].y,p[j].x,p[j].y);
        if(mi>temp)
        {
            mi=temp;
            dian1=i;
            dian2=j;
        }
    }
}
Point o;
double r=mi/2.0+0.001;
o.x=(p[dian1].x+p[dian2].x)/2.0;
o.y=(p[dian1].y+p[dian2].y)/2.0;
printf("%lf %lf %lf\n",o.x,o.y,r);
///
题意：给你一个无限深矩形的洞和一个立方体，问你这个立方体能否放入这个洞里面
题解：注意可以斜着放，那就枚举角度往里面塞就行了
const double pi = acos(-1.0);
double a[3],d,e;
bool check(int x,int y)
{
    for(double i=0;i<=90.0;i+=0.001)
    {
        double t=i*pi/180.0;
        double l=a[x]*cos(t)+a[y]*sin(t);
        double k=a[x]*sin(t)+a[y]*cos(t);
        if(l<=d&&k<=e) return true;
    }
    return false;
}
int main(void)
{
    for(int i=0;i<3;i++)
        scanf("%lf",&a[i]);
    scanf("%lf%lf",&d,&e);
    int flag=0;
    for(int i=0;i<3;i++)
    {
        for(int j=0;j<3;j++)
        {
            if(i==j) continue;
            if(check(i,j)) flag=1;
        }
        if(flag) break;
    }
    if(flag) puts("YES");
    else puts("NO");
    return 0;
}
///
题意：给你一个等n角形，然后求面积（内接圆的半径r、外接圆的半径R）
const double pi = acos(-1.0);
int main(void)
{
    double n,R,r;
    scanf("%lf%lf%lf",&n,&R,&r);
    printf("%.10f\n",sin(pi/n)*r*R*n);
    return 0;
}
///
题意：给出两条线段的两个端点。求一个圆，使其与两条线段都有且只有一个交点；
分析：找找出两个线段的两个端点的最小距离 ，这两个端点的中点为圆心，这个最小距离加上一个很小的值即为半径
double mi=INF;
int dian1=0,dian2=0;
for(int i=0;i<2;i++)
{
    for(int j=2;j<4;j++)
    {
        double temp=dist(p[i].x,p[i].y,p[j].x,p[j].y);
        if(mi>temp)
        {
            mi=temp;
            dian1=i;
            dian2=j;
        }
    }
}
Point o;
double r=mi/2.0+0.001;
o.x=(p[dian1].x+p[dian2].x)/2.0;
o.y=(p[dian1].y+p[dian2].y)/2.0;
printf("%lf %lf %lf\n",o.x,o.y,r);
/
题意：一个三维空间，给n个点，问这n个点是否在同一个平面上
题解：首先先找到不在同一条直线上的三个点，做法向量，
      然后我们再枚举任意两个点，如果这两个点是和法向量垂直的话，就说明在一个平面上
typedef long long ll;
const double eps = 1e-7;
const double pi = acos(-1.0);
const int maxn = 1e5+1000;
struct Point
{ll x,y,z;
};
Point A[maxn],vi;
bool check(Point a,Point b)//判断两个向量是否共线(三维) 向量成比例
{
    return a.x*b.y==b.x*a.y&&a.y*b.z==b.y*a.z&&a.x*b.z==a.z*b.x;
}
int main(void)
{
    int n;
    scanf("%d",&n);
    for(int i=0;i<n;i++)
        scanf("%I64d%I64d%I64d",&A[i].x,&A[i].y,&A[i].z);
    if(n<=3) printf("TAK\n");
    else{
        int flag=1;
        for(int i=0;i<n;i++)
        {
            int a=i;
            int b=(i+1)%n;
            int c=(i+2)%n;
            Point vi1,vi2;
            vi1.x = A[b].x - A[a].x;
            vi1.y = A[b].y - A[a].y;
            vi1.z = A[b].z - A[a].z;
            vi2.x = A[c].x - A[b].x;
            vi2.y = A[c].y - A[b].y;
            vi2.z = A[c].z - A[b].z;
            if(check(vi1,vi2)) continue;
            else//做向量vi1、vi2的法向量vi。
            {
                vi.x = vi1.y*vi2.z - vi2.y*vi1.z;
                vi.y = vi2.x*vi1.z - vi1.x*vi2.z;
                vi.z = vi1.x*vi2.y - vi1.y*vi2.x;
                flag=0;
                break;
            }
        }
        if(flag) printf("TAK\n");
        else{
            flag=1;
            for(int i=0;i<n;i++)
            {
                int a=i;
                int b=(i+1)%n;
                Point tx;
                tx.x = A[b].x - A[a].x;
                tx.y = A[b].y - A[a].y;
                tx.z = A[b].z - A[a].z;
                if(vi.x*tx.x+vi.y*tx.y+vi.z*tx.z!=0)
                {
                    flag=0;
                    break;
                }
            }
            if(flag) printf("TAK\n");
            else printf("NIE\n");
        }
    }
    return 0;
}
/
题意：求一个不规则简单多边形的重心。
double cross(Point a,Point b)
{
    return a^b;
}
double CalcConvexArea(Point* p,int n){ //凸包面积
    double area = 0.0;
    for(int i=1;i<n-1;i++)
        area += cross(p[i]-p[0],p[i+1]-p[0]);
    return area*0.5;
}
Point p[110];
int main()
{
    int n,cas=1;
    while(scanf("%d",&n)!=EOF,n)
    {
        for(int i=0;i<n;i++) p[i].input();
        double S = CalcConvexArea(p,n);
        double X = 0.0, Y = 0.0;
        for(int i=2;i<n;i++){
            double area = 0.5*cross(p[i-1]-p[0],p[i]-p[0]);
            X += area*(p[0].x+p[i-1].x+p[i].x)/3.0;
            Y += area*(p[0].y+p[i-1].y+p[i].y)/3.0;
        }
        printf("Stage #%d: %.6f %.6f\n",cas++,X/S,Y/S);
    }
    return 0;
}
//
题意:给你n个在x轴的圆，半径都是r，问你总共的面积是多少
题解:用总的面积减去两两相交的面积就好了
const double eps=1e-6;
const double PI = acos(-1.0);
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
double work(Point c1,double r1,Point c2,double r2)//两圆相交的面积
{
    double d=dist(c1,c2);
    if(r1+r2<d+eps)return 0;
    if(d<fabs(r1-r2)+eps)
    {
        double r=min(r1,r2);
        return PI*r*r;
    }
    double x=(d*d+r1*r1-r2*r2)/(2*d);
    double t1=acos(x/r1);
    double t2=acos((d-x)/r2);
    return r1*r1*t1+r2*r2*t2-d*r1*sin(t1);
}
bool cmp(const Point a,const Point b)
{
    return a.x<b.x;
}
Point p[5050];
int main(void)
{
    int n;
    double r;
    scanf("%d%lf",&n,&r);
    for(int i=1;i<=n;i++)
    {
        scanf("%lf",&p[i].x);
        p[i].y=0;
    }
    sort(p+1,p+1+n,cmp);
    double ans=n*PI*r*r;
    for(int i=1;i<n;i++)
    {
        if(2*r>fabs(p[i].x-p[i+1].x))
            ans-=work(p[i],r,p[i+1],r);
    }
    printf("%.12f\n",ans);
    return 0;
}

题意:给你两个点A,B.然后给你n条直线，这n条直线会把平面切成几块，然后问你从A走到B至少跨过多少条直线。
题解：对于每一条直线，判断这两个点是否在同侧就好了。但是乘法会爆 long long。所以就直接判断符号吧
int main(void)
{
    ll x1,x2,y1,y2;
    cin>>x1>>y1>>x2>>y2;
    int n;
    cin>>n;
    int ans=0;
    for(int i=0;i<n;i++)
    {
        ll a,b,c;
        cin>>a>>b>>c;
        ll temp1=a*x1+b*y1+c;
        ll temp2=a*x2+b*y2+c;
        if(temp1>0 != temp2>0) ans++;
    }
    cout<<ans<<endl;
    return 0;
}
///
题意:给你一个多边形,然后给你一个原点,问你这个多边形绕着这个原点旋转一圈所经过的面积是多少,原点一定不在多边形内
题解:很显然是一个大圆减去一个小圆的面积，所以我们直接找原点到该多边形最远的距离和最近的距离就好了
     最远的距离一定是原点到某个顶点的距离（取最大）
     最近的距离一定是原点到每条边（线段）的距离（取最小）
const double eps=1e-6;
const double PI = acos(-1.0);
const double INF = 0x3f3f3f3f;
const int maxn=1e5+100;
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
double point_to_seg(Point a,Point b,Point c)//点a到线段bc的距离
{
    Point v1=c-b,v2=a-b,v3=a-c;
    if((v1*v2)<0) return sqrt(v2*v2);
    else if((v1*v3)>0) return sqrt(v3*v3);
    else return fabs((v1^v2)/sqrt(v1*v1));
}
Point p[maxn],o;
int main(void)
{
    int n;
    scanf("%d",&n);
    o.input();
    for(int i=0;i<n;i++)
        p[i].input();
    p[n]=p[0];
    double ma=0,mi=INF;
    for(int i=0;i<n;i++)
    {
        ma=max(ma,dist(o,p[i]));
        mi=min(mi,point_to_seg(o,p[i],p[i+1]));
    }
    printf("%.10f\n",ma*ma*PI-mi*mi*PI);
    return 0;
}
/
题意：给出三个点，求出以这三个点为顶点的最小正多边形的面积。
(正多边形内一定有一点,这一点到正多边形各个角的顶点距离相等.
    这一点到正多边形各个角的距离我们称之为正多边形的半径)
思路:给你三个点，很容易算出多边形的半径:R=abc/4S，abc是三边边长，S是三角形面积。
     然后能够构成的最小多边形就是n=pi/gcd(A,B,C)。
     gcd是浮点数。然后强行算一波面积就好了。
const double eps=1e-5;
const double PI = acos(-1.0);
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
double gcd(double x,double y)
{
    while(fabs(x)>eps && fabs(y)>eps)
    {
        if(x>y) x-=floor(x/y)*y;
        else y-=floor(y/x)*x;
    }
    return x+y;
}
Point a,b,c;
int main(void)
{
    a.input();b.input();c.input();
    double A,B,C;
    A=dist(b,c);B=dist(a,c);C=dist(a,b);
    double p=(A+B+C)/2.0;
    double S=sqrt(p*(p-A)*(p-B)*(p-C));
    double AA=acos((B*B+C*C-A*A)/(2.0*B*C));
    double BB=acos((C*C+A*A-B*B)/(2.0*A*C));
    double CC=acos((A*A+B*B-C*C)/(2.0*A*B));
    double R=(A*B*C)/(4.0*S);
    double n=(PI/(gcd(AA,gcd(BB,CC))));
    printf("%.12f\n",n*R*R*sin(2.0*PI/n)/2.0);
    return 0;
}
/
题意：求n个圆覆盖的总面积。时间复杂度O(n^2*log(n))
思路：用自适应Simpson积分来求圆并
#include <cstdio>
#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <string>
#include <cstdlib>
#include <stack>
#include <map>
#include <queue>
using namespace std;
const int maxn=1e5+100;
const double eps=1e-6;
int sgn(double x)
{
    if(fabs(x)<eps) return 0;
    if(x<0) return -1;
    else return 1;
}
struct Point
{
    double x,y;
    Point(){}
    Point(double _x,double _y){x=_x,y=_y;}
    Point operator-(const Point &b)const
    {
        return Point(x-b.x,y-b.y);
    }
    double operator^(const Point &b)const//叉乘
    {
        return x*b.y-y*b.x;
    }
    double operator*(const Point &b)const//点积
    {
        return x*b.x+y*b.y;
    }
    bool operator<(const Point &b)const
    {
        return x==b.x?y<b.y:x<b.x;
    }
    void input()
    {
        scanf("%lf%lf",&x,&y);
    }
};
struct Circle
{
    Point o;
    double r;
    Circle(){}
    Circle(Point _o,double _r){o=_o,r=_r;}
    void input()
    {
        scanf("%lf%lf%lf",&o.x,&o.y,&r);
    }
};
bool cmp1(Circle a,Circle b)//按圆最左端升序排序
{
    if(fabs(a.o.x-a.r-b.o.x+b.r)<=eps) return a.o.x+a.r<b.o.x+b.r;
    return a.o.x-a.r<b.o.x-b.r;
}
bool cmp2(Circle a,Circle b)//按半径升序排序
{
    return a.r>b.r;
}
int n;
Circle circle[maxn],temp[maxn];
pair<double,double>intervals[maxn];
int tot=0;
int st,ed;
bool check(Circle a,Circle b)//检查圆b是否被套在圆a里面
{
    return (a.o.x-b.o.x)*(a.o.x-b.o.x)+(a.o.y-b.o.y)*(a.o.y-b.o.y)<=(a.r-b.r)*(a.r-b.r);
}
void prework()
{
    sort(circle+1,circle+1+n,cmp2);//按半径升序排序
    int cnt=0;
    for(int i=1;i<=n;i++)
    {
        bool flag=true;
        for(int j=1;j<=cnt;j++)
        {
            if(check(temp[j],circle[i]))
            {
                flag=false;
                break;
            }
        }
        if(flag) temp[++cnt]=circle[i];
    }
    n=cnt;
    for(int i=1;i<=n;i++)
        circle[i]=temp[i];
}
void getCircleIntersec(Circle a,double x)//求x=x这条平行于y轴的直线穿过圆a的长度(和当前圆的两交点y坐标)
{
    intervals[++tot]=make_pair(a.o.y-sqrt(a.r*a.r-(x-a.o.x)*(x-a.o.x)),a.o.y+sqrt(a.r*a.r-(x-a.o.x)*(x-a.o.x)));
}
double f(double x)//求扫描线x被圆覆盖部分的长度
{
    tot=0;
    for(int i=st;i<=ed;i++)
        if(x<circle[i].o.x+circle[i].r&&x>circle[i].o.x-circle[i].r)
            getCircleIntersec(circle[i],x);
    sort(intervals+1,intervals+tot+1);
    double ans=0,start=-1e12,end=-1e12; //!!!!!!!!!!!
    for(int i=1;i<=tot;i++)
    {
        if(intervals[i].first>=end)
        {
            ans+=end-start;
            start=intervals[i].first;
            end=intervals[i].second;
        }
        else end=max(end,intervals[i].second);
    }
    ans+=end-start;
    return ans;
}
//积分公式
double calc(double len,double fl,double fmid,double fr)//求长度为len的[L,R]区间，中点为M的Simpson近似面积
{
    return (fl+4*fmid+fr)*len/6;
}
//Simpson积分求区间[L,R]的面积并,f(L)=L,f(R)=R,f(M)=M,把[L,R]当成整体来拟合得到的面积是sqr
double Simpson(double l,double mid,double r,double fl,double fm,double fr,double sqr)
{
    double m1=(l+mid)/2,m2=(mid+r)/2;
    double fm1=f(m1),fm2=f(m2);
    double g1=calc(mid-l,fl,fm1,fm),g2=calc(r-mid,fm,fm2,fr);
    if(fabs(sqr-g1-g2)<=eps)//把当前区间分成2半再拟合得到的答案差别很小，就不再递归下去了
        return g1+g2;
    return Simpson(l,m1,mid,fl,fm1,fm,g1)+Simpson(mid,m2,r,fm,fm2,fr,g2);
}
int main(void)
{
    double ans=0;
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
        circle[i].input();
    prework();//要做一个预处理，删掉所有被其他大圆完全覆盖的小圆
    sort(circle+1,circle+n+1,cmp1);
    for(int i=1,j;i<=n;i++)
    {
        double l=circle[i].o.x-circle[i].r,r=circle[i].o.x+circle[i].r;
        for(j=i+1;j<=n;j++)
        {
            if(circle[j].o.x-circle[j].r>r) break;
            else r=max(r,circle[j].o.x+circle[j].r);//!!!!!!!!!!!!!
        }
        double mid=(l+r)/2;
        st=i,ed=j-1;
        i=j-1;
        double fl=f(l),fm=f(mid),fr=f(r);
        ans+=Simpson(l,mid,r,fl,fm,fr,calc(r-l,fl,fm,fr));
    }
    printf("%.3f\n",ans);
    return 0;
}

题意：给出椭圆,l和r求在这之间的椭圆的面积
const double eps=1e-6;
double a,b,l,r;
//simpson公式用到的函数
double F(double x)
{
    return b*sqrt(1-(x*x)/(a*a));
}
//三点simpson法。这里要求F是一个全局函数
double simpson(double a,double b)
{
    double c=a+(b-a)/2;
    return (F(a)+4*F(c)+F(b))*(b-a)/6;
}
// 自适应Simpson公式（递归过程）。已知整个区间[a,b]上的三点simpson值A
double asr(double a,double b,double eps,double A)
{
    double c=a+(b-a)/2;
    double L=simpson(a,c),R=simpson(c,b);
    if(fabs(L+R-A)<=15*eps)return L+R+(L+R-A)/15.0;
    return asr(a,c,eps/2,L) + asr(c,b,eps/2,R);
}
// 自适应Simpson公式（主过程）
double asr(double a,double b,double eps)
{
    return asr(a,b,eps,simpson(a,b));
}
int main(void)
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lf%lf%lf%lf",&a,&b,&l,&r);
        double ans=asr(l,r,eps);
        ans*=2;
        printf("%.3f\n",ans);
    }
    return 0;
}
//
[UVA-10652]
题意：在二维平面给出一些矩形的重心坐标、长于宽，以及顺时针旋转角度，求矩形面积与凸包面积的比例。
题解：计算凸包->凸包面积 顺时针一定要是负.还有给的ang要转化为rad
#include <cstdio>
#include <cstring>
#include <cmath>
#include <iostream>
#include <algorithm>
#include <string>
#include <cstdlib>
#include <stack>
#include <map>
#include <queue>
using namespace std;
const int maxn=1e5+100;
const double PI=acos(-1.0);
const double eps=1e-6;
int sgn(double x)
{
    if(fabs(x)<eps) return 0;
    if(x<0) return -1;
    else return 1;
}
struct Point
{
    double x,y;
    Point(){}
    Point(double _x,double _y){x=_x,y=_y;}
    Point operator-(const Point &b)const
    {
        return Point(x-b.x,y-b.y);
    }
    Point operator+(const Point &b)const
    {
        return Point(x+b.x,y+b.y);
    }
    double operator^(const Point &b)const//叉乘
    {
        return x*b.y-y*b.x;
    }
    double operator*(const Point &b)const//点积
    {
        return x*b.x+y*b.y;
    }
    friend bool operator <(const Point a,const Point b)
    {
        return a.x<b.x||(a.x==b.x&&a.y>b.y);
    }
    void input()
    {
        scanf("%lf%lf",&x,&y);
    }
};
double torad(double deg)
{
    return deg/180*PI;
}
double cross(Point a,Point b)
{
    return a^b;
}
Point Rotate(Point a,double rad)
{
    return Point(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));
}
//多边形有向面积
double CalcConvexArea(Point* p,int n){ //凸包面积
    double area = 0.0;
    for(int i=1;i<n-1;i++)
        area += cross(p[i]-p[0],p[i+1]-p[0]);
    return area*0.5;
}
//计算凸包，输入点数组p，个数为n，输出点数组ch，函数返回凸包顶点数m
//输入不能有重复点，函数执行完后输入的点的顺序被破坏
//如果不希望在凸包的边上有输入点，把两个<=改成<
//在精度要求高时建议用dcmp比较
int ConvexHull(Point* p , int n, Point* ch)
{
    sort(p,p+n);        //先比较x坐标，再比较y坐标
    int m=0;
    for (int i = 0; i < n; i++) {
        while (m > 1 && cross(ch[m-1]-ch[m-2],p[i]-ch[m-2]) <= 0) m--;
        ch[m++] = p[i];
    }
    int k=m;
    for (int i = n - 2; i >= 0; i--) {
        while (m > k && cross(ch[m-1]-ch[m-2],p[i]-ch[m-1]) <= 0) m--;
        ch[m++] = p[i];
    }
    if(n>1) m--;
    return m;
}
Point p[2500],ch[2500];
int main(void)
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        int n,cnt=0;
        double area1=0,area2=0;
        scanf("%d",&n);
        for(int i=0;i<n;i++)
        {
            double x,y,w,h,j,ang;
            scanf("%lf%lf%lf%lf%lf",&x,&y,&w,&h,&j);
            Point o=Point(x,y);
            ang=-torad(j);//顺时针旋转
            p[cnt++]=o+Rotate(Point(-w/2,-h/2),ang);
            p[cnt++]=o+Rotate(Point(w/2,-h/2),ang);
            p[cnt++]=o+Rotate(Point(-w/2,h/2),ang);
            p[cnt++]=o+Rotate(Point(w/2,h/2),ang);
            area1+=w*h;
        }
        int m=ConvexHull(p,cnt,ch);
        area2=CalcConvexArea(ch,m);
        printf("%.1f %%\n",area1*100/area2);
    }
    return 0;
}

题意：给出等腰三角形ABC，AB=AC，M为BC中点。P点为三角形内使min{∠MPB+∠APC,∠MPC+∠APB}
      最大的点。求P点轨迹长度。
题解:|AM|+一段圆弧(点B、C、以及ABC内心这三点外接圆在三角形ABC内的那段圆弧)
double dist(Point a,Point b)
{
    return sqrt((a-b)*(a-b));
}
Point a,b,c;
int main(void)
{
    int T;
    scanf("%d",&T);
    for(int cas=1;cas<=T;cas++)
    {
        a.input();b.input();c.input();
        Point m=Point((b.x+c.x)/2,(b.y+c.y)/2);
        long double ans=dist(a,m);
        long double ab=dist(a,b),ac=dist(a,c),bc=dist(b,c);
        long double bac=acos((ab*ab+ac*ac-bc*bc)/(2*ab*ac));
        long double bo=tan(bac/2.0)*ab;
        ans+=bo*(pi-bac);
        printf("Case #%d: %.4f\n",cas,(double)ans);
    }
    return 0;
}
/
题意:求一个给定的圆(x^2+y^2=r^2)，在圆周上有多少个点的坐标是整数。
typedef long long ll;
ll gcd(ll x,ll y)
{
    return x%y==0? y:gcd(y,x%y);
}
bool check(ll y,double x)
{
    if(x==floor(x))
    {
        ll x1=(ll)floor(x);
        if(gcd(x1*x1,y*y)==1&&x1*x1!=y*y)
            return true;
    }
    return false;
}
int main()
{
    ll r,ans=0;
    scanf("%lld",&r);
    for(ll d=1; d<=sqrt(2*r); d++)
    {
        if((2*r)%d==0)
        {
            for(ll a=1; a<=(ll)sqrt(2*r/(2*d)); a++)
            {
                double b=sqrt(((2*r)/d)-a*a);
                if(check(a,b))ans++;
            }
            if(d!=(2*r)/d)
            {
                for(ll a=1; a<=(ll)sqrt(d/2); a++)
                {
                    double b=sqrt(d-a*a);
                    if(check(a,b))
                        ans++;
                }
            }
        }
    }
    printf("%lld",ans*4+4);
    return 0;
}
/
bool Seg_inter_line(Point a,Point b,Line l2) //判断直线ab和线段l2是否相交 相交返回true否则返回false
{
    return (sgn((l2.a-b)^(a-b))*sgn((l2.b-b)^(a-b))) <= 0;
}
double Intersect(Point a,Point b,Point c,Point d)//返回直线ab和线段cd相交点的横坐标
{
    double area1=cross(a,b,c);
    double area2=cross(a,b,d);
    int c1=sgn(area1);
    int c2=sgn(area2);
    if(c1*c2<0)//如果小于零,规范相交
        return (area2*c.x-area1*d.x)/(area2-area1);//计算交点的公式
    else if(c1*c2==0)//和端点相交,判断和哪一个端点相交
    {
        if(c1==0) return c.x;
        else if(c2==0) return d.x;
    }
    return -INF;
}
/