计算几何说起来都是一套一套的,写起来却让人二楞二楞的qwq
接下来直接粘代码,稍微高级的算法讲解详见超链接
【基本的定义&算法】
#include<cmath>
const double eps = 1e-7;//误差允许范围,广泛用于大小比较
inline double dcmp(double x)//强烈建议浮点数的大小比较都用这个函数
{
if(fabs(x) < eps) return 0;
return x < 0 ? -1 : 1;
}
struct point{//由于点与向量在表示形式方面一样,所以就不额外定义一个结构体表示向量了。
double x, y, ang;//x坐标,y坐标,极角
point(){}
point(double a, double b) : x(a), y(b){ang = atan2(y, x);}//算极角
friend double dots(point a, point b) {return a.x * b.x + a.y * b.y;}//点积
friend double cross(point a, point b){return a.x * b.y - a.y * b.x;}//叉积
friend double dist(point a, point b) {return sqrt(dots(a - b, a - b));} //模长
point operator +(const point a) const {return point(x + a.x, y + a.y);}//向量加法
point operator -(const point a) const {return point(x - a.x, y - a.y);}//两点间的向量、向量的减法
point operator *(const double a)const {return point(a * x, a * y);} //向量的数乘
bool operator <(const point a) const {return dcmp(ang - a.ang) == 0 ? dcmp(dots(point(x, y), point(x, y)) - dots(a, a)) < 0 : dcmp(ang - a.ang) < 0;}
bool operator ==(const point a) const {return dcmp(ang - a.ang) == 0 && dcmp(dots(point(x, y), point(x, y)) - dots(a, a)) == 0;}
/*关系运算符因题而异(此处的重载用于极角排序)*/
};
struct line{
point p, v;//一个定点,方向向量
double ang;//极角
line(){}
line(point a, point b) : p(a), v(b - a){ang = atan2(v.y, v.x);}//算极角
friend point get_cross(line a, line b) {return a.p + a.v * (cross(b.v, a.p - b.p) / cross(a.v, b.v));} //求交点
friend bool onleft(line a, point b) {return dcmp(cross(a.v, b - a.p)) > 0;}//判方向
bool operator < (const line a) const {return ang < a.ang;}
/*关系运算符因题而异(此处的重载用于极角排序)*/
};
inline double get_S(int n)//多边形面积
{
double ret = 0;
for(int i = 2; i < n; i++)
ret += cross(p[i] - p[1], p[i + 1] - p[1]);
return ret / 2.0;
}
int main()
{
//这里写东西
}
求凸包
//由于极角排序的算法导致点的储存不稳定(即无论如何处理,都可能会有一些边上点无法完全保留或完全剔除,故采用以x为第一关键字,y为第二关键字的方法)
//本模板将所有凸包上的点都保留了下来
point stk[mn], p[mn]//为保证好读,此处的point已根据需要重载运算符“<”和“==”;
int top;
inline void make_polygon(int n)
{
sort(p + 1, p + 1 + n), top = 0;//排序
stk[++top] = p[1], stk[++top] = p[2];//放入前两个点
for(i = 3; i <= n; i++)
{
while(top > 1 && cross(stk[top] - stk[top - 1], p[i] - stk[top]) > 0)//判断凹凸性
--top;
stk[++top] = p[i];
}
/*上下凸壳各做一次,然后合并得凸包*/
for(i = n; i; i--)
{
while(top > 1 && cross(stk[top] - stk[top - 1], p[i] - stk[top]) > 0)
--top;
stk[++top] = p[i];
}
int siz = unique(stk + 1, stk + 1 + top) - stk - 1;//去重
}
半平面交
point p[mn];
int cnt, h, t;
line q[mn], l[mn];
inline double merge_halfplane(int cnt)
{
sort(l + 1, l + 1 + cnt), q[h = t = 1] = l[1];//初始化
for(i = 2; i <= cnt; i++)
{
while(h < t && !onleft(l[i], p[t - 1])) t--;
while(h < t && !onleft(l[i], p[h])) h++;//剔除无用直线
q[++t] = l[i]; //塞入直线
if(dcmp(cross(q[t - 1].v, q[t].v)) == 0) //剔除先进入的同斜率直线
{
t--;
if(onleft(q[t], l[i].p)) //此时可能不会用上这条直线
q[t] = l[i];
}
if(h < t)
p[t - 1] = get_cross(q[t - 1], q[t]); //获取交点
}
while(h < t && !onleft(q[h], p[t - 1])) t--;//判断首尾连接处
p[t] = get_cross(q[t], q[h]);
return get_S(h, t); //算面积
}