关于计算几何

这里是蒟蒻整理的一些计算几何的知识。挖一个坑，不定期填坑。。。

一.叉积的性质
1.a*b=-b*a;//反交换律
2.不满足结合律，但满足 a*(b*c)+b*(a*c)+c*(a*b)=0
3.a*b=0//向量a，b平行
4.拉格朗日公式：a*(b*c)=b*(a*c)-c*(a*b)
5.设向量P(x1,y1),Q(x2,y2) => P*Q=x1*y2-x2*y1，若P*Q>0,P在Q的顺时针方向； 若P*Q<0,P在Q的逆时针方向； 若=0，两者共线。设P*Q=delta,那么|delta|=|p|*|Q|*sin(theta)

int get_cross(node a,node b,node c){  //calculate the cross product of vector ab and vector ac
    return ((b.x-a.x)*(c.y-a.y)-(b.y-a.y)*(c.x-a.x));  
} 

二.凸包
1.构造：首先找基准点(先按照y从大到小排，若y相同，就按照x从大到小排)；再极角排序；然后将左下角的元素加入栈，对于栈顶的元素o和栈顶第二个元素p若当前元素x满足xo*xp>0弹出栈顶元素。相同的，再将右上角的元素加入栈，重复这种操作。
2.极角排序：
a.atan2函数(要注意可能卡精度)：atan2(y,x)求坐标原点与坐标(x,y)所连的射线与x轴的夹角的大小
   或者是：先排极角(根据叉积来排)，若极角相同，就根据距离基准点的远近来排(根据这个方法写一个bool型的cmp函数，再sort即可)

struct node{  
    int x,y;  
}p[maxn],point,stack[maxn];  
  
int get_cross(node a,node b,node c){  
    return ((b.x-a.x)*(c.y-a.y)-(b.y-a.y)*(c.x-a.x));  
}  
  
double get_dis(node a,node b){  
    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));  
}  
  
void find_point(){  
    int col=0; node temp; point=p[0];  
    for(int i=1;i<N;i++){  
        if(p[i].y<point.y||(p[i].y==point.y&&p[i].x<point.x)){  
            point=p[i]; col=i;  
        }  
    }  
    temp=p[0]; p[0]=p[col]; p[col]=temp;  
}  
  
bool cmp(node a,node b){  
    int K=get_cross(point,a,b);  
    if(K<0)return true;  
    if(K>0)return false;  
    double xx=get_dis(point,a),yy=get_dis(point,b);  
    return xx<yy;  
}  
  
void Graham(){  
    int top=2; p[N]=p[0];  
    stack[0]=p[0]; stack[1]=p[1]; stack[2]=p[2];  
    for(int i=3;i<=N;i++){  
        while(top>1&&get_cross(stack[top-1],stack[top],p[i])>=0)top--;  
        stack[++top]=p[i];  
    }  
    double ans=0;  
    node temp=stack[top],sai=temp; top--;  
    while(top){  
        node xx=stack[top]; top--;  
        ans+=get_dis(xx,temp);  
        temp=xx;  
    }  
    ans+=get_dis(sai,temp);  
    ans+=pi*L*2;  
    printf("%.0lf\n",ans);  
} 

三.旋转卡壳

1.凸包的对踵点个数不超过3N/2个
2.求凸包上最远点对：
a.枚举点：求到该点最远的点的距离(先按照一个方向旋转求解，这时据该点的点的距离是先递增再递减的，如果递减，就结束该次操作)，然后旋转继续求解(接着上次枚举断掉的点枚举)，更新最大值。
b.枚举边：求到一条边的距离最远的点的距离。先枚举一条边，在枚举点，对于每个点，计算该点到边的两个端点的距离，更新最大值。对于这个距离也应是先递增再递减的。当距离开始减小时，就停止这一次的操作，开始按照一定一定方向枚举下一条边，枚举点还是接着上一次的断点进行枚举，这样的复杂度就是线性的。由于点到边的距离增大，点与边构成的三角形的面积是增大的，所以计算可以用叉积做。
3.两个凸包间的最小距离：考虑三种情况。
a.点对点 b.点对边 c.边对边 (此时两个平行线的方向相反)
4.凸多边形外接矩形最小面积：将矩形的一条边与凸多边形的一条边重合，然后找最远点。用两组相互垂直的平行线即可求出此时的矩形面积。

 求凸多边形外接矩形最小周长和求凸多边形外接矩形最小面积依赖于一个结论：凸多边形的外接矩形与该凸多边形的一条边重合。这样就可以降低复杂度，用旋转卡壳做。

void rotating_colipers(){  
    int ans=0,x=1; stack[top+1]=p[0];  
    for(int i=0;i<=top;i++){  
        while(get_cross(stack[x+1],stack[i+1],stack[i])>get_cross(stack[x],stack[i+1],stack[i])){  
            x=(x+1)%top;  
        }  
        ans=max(ans,get_dis(stack[x],stack[i]));  
        ans=max(ans,get_dis(stack[x+1],stack[i+1]));  
    }  
    printf("%d",ans);  
}  

四.半平面交
1.半平面交算法(NlogN)的一般过程
1.将所有的半平面进行极角排序
2.用双端队列维护，加入最开始的两个半平面。对于新加入的一个半平面，如果队列顶端的两个半平面的交点在该平面的外面，删除顶端半平面(重复此操作直到在平面内)。如果队列底端的两个平面的交点在该平面的外面，删除底端的半平面。
3.对于最后的，删除两端多余的半平面。如果队列顶端的两个半平面的交点在底端半平面的外部，删除顶端的半平面。对于底端的半平面的交点在队列顶端半平面外部的，删除队列底端的半平面。(重复此操作直到在平面内)。
2.半平面交的精度控制在1e-10。

3.时间复杂度(N log N)，主要是排序的时间复杂度。

下面是一个还没有检查正确性的半平面交的模板= =。

const double eps = (double)1e-10;
const int maxn = 100 + 10;//根据实际情况决定

struct point{//点 
	int x, y;
}p[maxn];

struct segment{//平面
	point s, e;// s - > e
	double angle;
	void get_angle(){
		angle = atan2(e.y - s.y, e.x - s.x);
	}
}seg[maxn], deq[maxn];

double cross(point a, point b, point c){//vector ab and vertor ac。叉积大于0，说明点c在有向直线的左端
	return (b.x - a.x) * (c.y - a.y) - (b.y - a.y) * (c.x - a.x);
}

point get_intersect(segment a, segment b){//求两个直线的交点
	double u = cross(a.s, a.e, b.s), v = cross(a.e, a.s, b.e);
	point t;
	t.x = (b.s.x * v + b.e.x * u) / (u + v); t.y = (b.s.y * v + b.e.y * u) / (u + v);
	return t;
}

bool cmp(segment a, segment b){
	if(a.angle > b.angle)return true;//先根据极角排序
	else if(abs(a.angle - b.angle) < eps && cross(b.s, b.e, a.e) > -eps)return true;//极角相同，在内侧的排在前面
	else return false;
}

void halfplane(int N){
	sort(seg, seg + N, cmp);
	int temp = 1;
	for(int i = 1; i < N; ++i){//对于极角相同的平面，保留在内侧的平面
		if(abs(seg[i].angle - seg[temp - 1].angle) > eps){
			seg[temp ++] = seg[i];
		}
	}
	int M = temp;//进行删减后的实际平面的个数
	
	deq[0] = seg[0]; deq[1] = seg[1];
	int L = 0, R = 1;
	for(int i = 2; i < M; ++i){
		while(L < R && cross(seg[i].s, seg[i].e, get_intersect(deq[L], deq[L + 1])) < - eps)L ++;
		while(L < R && cross(seg[i].s, seg[i].e, get_intersect(deq[R], deq[R - 1])) < - eps)R --;
		deq[++ R] = seg[i];
	}
	
	while(L < R && cross(deq[L].s, deq[L].e, get_intersect(deq[R], deq[R - 1])) < - eps)R --;
	while(L < R && cross(deq[R].s, deq[R].e, get_intersect(deq[L], deq[L + 1])) < - eps)L ++;
	
	if(L == R)return ;
	int num = 0;
	for(int i = L; i < R; ++i){//计算剩下的直线的交点
		p[++ num] = get_intersect(deq[i], deq[i + 1]);
	}
	if(R > L + 1)p[++ num] = get_intersect(deq[L], deq[R]);//对于构成凸包的情况
}

double get_area(int N){//选取p[0]点进行三角剖分，计算面积即为叉积的1/2
	double area = 0;
	for(int i = 1; i < N - 1; i++){//点是从 0 到 N - 1 标号的
		area += cross(p[0], p[i], p[i + 1]);
	}
	return abs(area) / 2.0;
}