计算几何基础
(以下内容均默认在二维平面下讨论的)
精度误差
程序设计中,浮点数 double
存在精度误差,所以在进行大小比较时,通常允许一定的误差,即对于两个数x, y 如果 |x-y|< eps,则认为 x=y 。eps是根据题目要求设置的极小值,一般为10-8。
(d为x, y的差值,相等时函数返回0,x>y时返回1,x<y时返回-1)
int sgn(double d){ //精度控制 if(fabs(d)<eps) return 0; if(d>0) return 1; return -1; }
点和矢量
点的定义:可以直接用横纵坐标表示
struct Point{ //二维平面点 double x,y; Point(double x=0, double y=0 ):x(x),y(y) { } };
矢量的定义:矢量是有大小和方向的量,其表示方式依旧可以用横纵坐标表示,与点的表示方式一致
typedef Point Vector; //二维平面向量;
重载相关的运算符:
//向量加法 Vector operator + (Vector A, Vector B){ return Vector(A.x+B.x, A.y+B.y); } //向量减法 Vector operator - (Vector A, Vector B){ return Vector(A.x-B.x, A.y-B.y); } //向量与常数的乘积 Vector operator * (Vector A, double p){ return Vector(A.x*p, A.y*p); } //向量与常数的除法 Vector operator / (Vector A, double p){ return Vector(A.x/p, A.y/p); } //将点集按照x坐标升序排序 bool operator < (const Point &a, const Point &b){ if(a.x==b.x) return a.y<b.y; return a.x<b.x; } //点坐标是否相等判断 bool operator == (const Point &a, const Point &b){ if(dcmp(a.x,b.x)==0&&dcmp(a.y,b.y)==0) return true; else return false; }
叉积和点积
点积
几何意义:点积的结果是一个标量,等于向量大小与夹角的cos值的乘积。
a•b = |a||b|cosθ
在坐标系中:设 A⃗ =(x1, y1), B⃗ =(x2, y2),则点积的计算如下
//向量的点积 double Dot(Vector A,Vector B){ //两向量夹角为锐角时为正,钝角为负,直角为0 return A.x*B.x + A.y*B.y; }
叉积
叉积的几何意义:叉积的结果是一个矢量, |c|=|a×b|=|a| |b|sinα (α为a,b向量之间的夹角)
叉积的应用有很多,后续会提出
//向量的叉积 double Cross(Vector A, Vector B){ return A.x*B.y - A.y*B.x; }
一些矢量的常见操作:
取模
//取模(长度) double Length(Vector A){ return sqrt(Dot(A,A)); //sqrt(A.x*A.x+A.y*A.y); }
向量夹角
弧度制与角度制的换算 :
弧度制2π = 角度制360°
所以 弧度制1rad=360°/(2π)=(180/π)°
角度制1°=2π/360=π/180
//计算量向量的夹角 double Angle(Vector A, Vector B){
return acos(Dot(A,B)/Length(A)/Length(B)); //由点积 A•B = |A||B|cosθ,可以求出cosθ
}
两向量构成的平行四边的有向面积(叉积的应用)
double Area2(Point A, Point B, Point C){ //也是三角形面积的两倍 return Cross(B-A,C-A); }
计算向量逆时针旋转后的向量
Vector Rotate(Vector A, double rad){ return Vector(A.x*cos(rad)-A.y*sin(rad), A.x*sin(rad)+A.y*cos(rad)); }
计算向量逆时针旋转90°的单位法向量
Vector Normal(Vector A){ double L=Length(A); return Vector(-A.y/L, A.x/L); //由垂直可知叉积为1 }
判断两向量的转向
(叉积的应用)
//判断折线bc是不是ab的逆时针方向转向 bool ToLeftTest(Point a, Point b, Point c){ return Cross(b-a, c-b)>0; //如果bc在ab的逆时针方向则叉积大于0,顺时针方向小于,共线等于0 }
点和线
线的定义:两点可以确定一条直线或者线段
struct Line{ //直线的点向式定义 Point v, p; Line (Point v, Point p):v(v),p(p) {} Point point(double t){ //返回点p=v+(p-v)*t return v + (p-v)*t; } };
线的相关操作:
点到直线的距离
double DistanceToLine(Point P, Point A, Point B){ //点P到直线AB的距离 Vector v1=B-A, v2=P-A; return fabs(Cross(v1,v2)/Length(v1)); //不去绝对值得到的是有向距离 }
点到线段的距离
double DistanceToSegment(Point P, Point A, Point B){ //点p到线段AB的距离 if(A==B) return Length(P-A); //A=B,返回两点间距离 Vector v1=B-A, v2=P-A, v3=P-B; if(sgn(Dot(v1, v2)) < 0) //AB与AP为钝角,则A与P的距离最小,返回|AP| return Length(v2); if(sgn(Dot(v1, v3)) > 0) //BA与BP与钝角,则B与P的距离最下,返回|BP| return Length(v3); return DistanceToLine(P, A, B); //否则直接计算P到直线的距离 }
求点在直线上的投影点
Point GetLineProjection(Point P, Point A, Point B){ //点P在直线AB上的投影点 Vector v=B-A; return A + v*(Dot(v, P-A)/Dot(v, v)); }
判断点是否在线段上
bool OnSegment(Point p, Point a1, Point a2){ //判断点P是否在线段a1a2上 return sgn(Cross(a1-p, a2-p))==0&&sgn(Dot(a1-p,a2-p))<0; }
判断两线段是否相交
bool SegmentProperIntersection(Point a1, Point a2, Point b1, Point b2){ double c1 = Cross(a2-a1, b1-a1), c2 = Cross(a2-a1, b2-a1); double c3 = Cross(b2-b1, a1-b1), c4 = Cross(b2-b1, a2-b1); //if语句判断控制是否允许线段在端点处相交,如果允许则添加,不允许删去即可 if(!sgn(c1) || !sgn(c2) || !sgn(c3) || !sgn(c4)){ bool f1 = OnSegment(b1, a1, a2); bool f2 = OnSegment(b2, a1, a2); bool f3 = OnSegment(a1, b1, b2); bool f4 = OnSegment(a2, b1, b2); bool f = (f1|f2|f3|f4); return f; } return (sgn(c1)*sgn(c2) < 0 && sgn(c3)*sgn(c4) < 0); }
两直线的交点
(有交点情况下)
//计算两直线的交点(保证有交点下) Point GetLineIntersection(Point P, Vector v, Point Q, Vector w){ //v,w为直线的方向向量 Vector u = P-Q; double t = Cross(w,u)/Cross(v,w); return P+v*t; }
多边形
(以下讨论的是简单多边形)
多边形的表示:一般用顺时针或者逆时针排列的点集表示
多边形的面积:
(叉积的应用)
double PolygonArea(Point *p, int n){ //p为端点集合(一般按逆时针存储),n为端点个数 double s=0; for(int i=1;i<n-1;i++) s+=Cross(p[i]-p[0], p[i+1]-p[0]); //固定点p[0],将多边形分为n-2个三角形,利用叉积求其面积 return s; }
这里返回的是多边形的面积的两倍,根据叉积的性质,如果返回面积的负数则说明点是按顺时针排列,若为正数则为逆时针排列。
判断点与多边形的关系:
//判断点是否在多边形内,若点在多边形内返回1,在多边形外部返回0,在多边形上返回-1 int isPointInPolygon(Point p, Point *poly, int size){ int wn = 0; for(int i = 0; i < size; ++i){ if(OnSegment(p, poly[i], poly[(i+1)%size])) return -1; int k = sgn(Cross(poly[(i+1)%size] - poly[i], p - poly[i])); int d1 = sgn(poly[i].y - p.y); int d2 = sgn(poly[(i+1)%size].y - p.y); if(k > 0 && d1 <= 0 && d2 > 0) wn++; if(k < 0 && d2 <= 0 && d1 > 0) wn--; } if(wn != 0) return 1; return 0; }
凸多边形的判断方法:
凸多边形的概念:凸多边形就是多边形所有的点都在任意两个相邻顶点所连的直线的左边(按逆时针连接)。
bool is_covbag(){ int dir=0; p[n]=p[0],p[n+1]=p[1]; for(int i=0;i<=n-1;++i){ int tem=sgn(Cross(p[i+1]-p[i]),(p[i+2]-p[i+1]))); if(!dir) dir=tem; if(dir*tem<0) return 0; } return 1; }
圆
圆的表示方法:
圆一般用圆心和半径表示即可
struct Circle{ Point c; double r; Circle(Point c, double r):c(c), r(r) {}; //通过圆心角求坐标 Point point(double a){ return Point(c.x+cos(a)*r, c.y+sin(a)*r); } };
两圆的关系:
1.相交
两圆的圆心距离内之和小于两圆的半径之和。
2.相切
外切:两圆的圆心距离之和等于两圆的半径之和。
内容切:两圆的圆心距离之和等于两圆的半径之差。
3.相离
外离:两圆的圆心距离之和大于两圆的半径之和。
内含:两圆的圆心距离之和小于两圆的半径之差
圆与直线的交点:
//求圆与直线的交点 int getLineCircleIntersection(Line L, Circle C, double& t1, double& t2, vector<Point>& sol){ double a = L.v.x, b = L.p.x - C.c.x, c = L.v.y, d = L.p.y - C.c.y; double e = a*a + c*c, f = 2*(a*b + c*d), g = b*b + d*d - C.r*C.r; double delta = f*f - 4*e*g;//判别式 if(sgn(delta) < 0)//相离 return 0; if(sgn(delta) == 0){//相切 t1 = -f /(2*e); t2 = -f /(2*e); sol.push_back(L.point(t1));//sol存放交点本身 return 1; } //相交 t1 = (-f - sqrt(delta))/(2*e); sol.push_back(L.point(t1)); t2 = (-f + sqrt(delta))/(2*e); sol.push_back(L.point(t2)); return 2; }
结束语
基础篇的就写到这里了,本篇基本没有对结论的证明,更倾向于总结,给大家提供一个模板,同时对算法竞赛中的计算几何有一个大概的模型,不会的大家可以自己度娘 。