计算几何模板.md · 金铉杰/acm代码模板

# 计算几何代码模板

## 简单的知识板子

## 点与线：

```c++
const double eps = 1e−8;
const double DINF = 1e20;
const double PI = acos(−1.0);
const int maxp = 1010;
//Compares a double to zero
int sgn(double x){
	if(fabs(x) < eps)return 0;
 	if(x < 0)return -1;
 	else return 1;
}
 //square of a double
inline double sqr(double x){return x*x;}
```

### 点的模板：

```c++
struct Point{
 double x,y;
 Point(){}
 Point(double _x,double _y){
 x = _x;
 y = _y;
 }
 void input(){
 scanf("%lf%lf",&x,&y);
 }
 void output(){
 printf("%.2f %.2f\n",x,y);
 }
 bool operator == (Point b)const{
 return sgn(x-b.x) == 0 && sgn(y-b.y) == 0;
 }
 bool operator < (Point b)const{
 return sgn(x-b.x)== 0?sgn(y-b.y)<0:x<b.x;
 }
 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;
 }
 //返回长度
 double len(){
 return hypot(x,y);//库函数
 }
 //返回长度的平方
 double len2(){
 return x*x + y*y;
 }
 //返回两点的距离
 double distance(Point p){
 return hypot(x-p.x,y-p.y);
 }
 Point operator +(const Point &b)const{
 return Point(x+b.x,y+b.y);
 }
 Point operator *(const double &k)const{
 return Point(x*k,y*k);
 }
 Point operator /(const double &k)const{
 return Point(x/k,y/k);
 }
 //`计算pa 和 pb 的夹角`
 //`就是求这个点看a,b 所成的夹角`
 //`测试 LightOJ1203`
 double rad(Point a,Point b){
 Point p = *this;
 return fabs(atan2( fabs((a-p)^(b-p)),(a-p)*(b-p) ));
 }
 //`化为长度为r的向量`
 Point trunc(double r){
 double l = len();
 if(!sgn(l))return *this;
 r /= l;
 return Point(x*r,y*r);
 }
 //`逆时针旋转90度`
 Point rotleft(){
 return Point(-y,x);
 }
 //`顺时针旋转90度`
 Point rotright(){
 return Point(y,-x);
 }
 //`绕着p点逆时针旋转angle`
 Point rotate(Point p,double angle){
 Point v = (*this) - p;
 double c = cos(angle), s = sin(angle);
 return Point(p.x + v.x*c - v.y*s,p.y + v.x*s + v.y*c);
 }
};
```

### 线

```c++
struct Line {
	Point s,e;
	Line() {}
	Line(Point _s,Point _e) {
		s=_s;
		e=_e;
	}
	bool operator ==(Line v) {
		return (s==v.s) && (e==v.e);
	}
	//根据一个点和倾斜角 angle 确定直线,0<=angle<π
	Line(Point p, double angle) {
		s=p;
		if (sgn(angle-PI/2)==0) {
			e=(s+Point(0,1));
		} 
		else {
			e=(s+Point(1,tan(angle)));
		}
	}
	//ax+by+c=0
	Line(double a,double b,double c) {
		if(sgn(a)==0) {
			s=Point(0,-c/b);
			e=Point(1,-c/b);
		} 
		else if(sgn(b)==0) {
			s=Point(-c/a,0);
			e=Point(-c/a,1);
		} 
		else {
			s=Point(0,-c/b);
			e=Point(1,(-c-a)/b);
		}
	}
	void input() {
		s.input();
		e.input();
	}
	void adjust() {
		if(e<s) swap(s,e);	
	}
	//求线段长度
	double length() {
		return s.distance(e);
	}
	//返回直线倾斜角 0<=angle<π 
	double angle() {
		double k=atan2(e.y-s.y,e.x-s.x);
		if(sgn(k)<0) k+=PI;
		if(sgn(k-PI)==0) k-= PI;
		return k;
	}
	//点和直线关系
	// 1 在左侧
	// 2 在右侧
	// 3 在直线上 
	int relation(Point p) {
		int c=sgn((p-s)^(e-s));
		if(c<0) return 1; 
		else if(c>0) return 2; 
		else return 3;	 
	}
	//点在线段上的判断
	bool pointonseg(Point p) {
		return sgn((p-s)^(e-s))==0 && sgn((p-s)*(p-e))<=0;
	}
	//两向量平行 (对应直线平行或重合)
	bool parallel(Line v) {
		return sgn((e-s)^(v.e-v.s))==0;
	}
	//两线段相交判断
	//2 规范相交
	//1 非规范相交
	//0 不相交
	int segcrossseg(Line v) {
		int d1=sgn((e-s)^(v.s-s));
		int d2=sgn((e-s)^(v.e-s));
		int d3=sgn((v.e-v.s)^(s-v.s));
		int d4=sgn((v.e-v.s)^(e-v.s));
		if((d1^d2)==-2&&(d3^d4)==-2)return 2;
		return (d1==0 && sgn((v.s-s)*(v.s-e))<=0) ||
		 (d2==0 && sgn((v.e-s)*(v.e-e))<=0) ||
		 (d3==0 && sgn((s-v.s)*(s-v.e))<=0) ||
		 (d4==0 && sgn((e-v.s)*(e-v.e))<=0);
	}
	//直线和线段相交判断
	//-*this line -v seg
	//2 规范相交
	//1 非规范相交
	//0 不相交
	int linecrossseg(Line v) {
		int d1=sgn((e-s)^(v.s-s));
		int d2=sgn((e-s)^(v.e-s));
		if((d1^d2)==-2) return 2;
		return (d1==0||d2==0);
	}
	//两直线关系
	//0 平行
	//1 重合
	//2 相交
	int linecrossline(Line v) {
		if((*this).parallel(v))
			return v.relation(s)==3;
		return 2;
	}
	//求两直线的交点
	//要保证两直线不平行或重合
	Point crosspoint(Line v) {
		double a1 = (v.e-v.s)^(s-v.s);
		double a2 = (v.e-v.s)^(e-v.s);
		return Point((s.x*a2-e.x*a1)/(a2-a1),(s.y*a2-e.y*a1)/(a2-a1
		 ));
	}
	//点到直线的距离
	double dispointtoline(Point p) {
		return fabs((p-s)^(e-s))/length();
	}
	//点到线段的距离
	double dispointtoseg(Point p) {
		if(sgn((p-s)*(e-s))<0 || sgn((p-e)*(s-e))<0)
			return min(p.distance(s),p.distance(e));
		return dispointtoline(p);
	}
	//返回线段到线段的距离
	//前提是两线段不相交，相交距离就是 0 了
	double dissegtoseg(Line v) {
		return min(min(dispointtoseg(v.s),dispointtoseg(v.e)),min(v
		 .dispointtoseg(s),v.dispointtoseg(e)));
	}
	//返回点 p 在直线上的投影
	Point lineprog(Point p) {
		return s + ( ((e-s)*((e-s)*(p-s)))/((e-s).len2()) );
	}
	//返回点 p 关于直线的对称点
	Point symmetrypoint(Point p) {
		Point q = lineprog(p);
		return Point(2*q.x-p.x,2*q.y-p.y);
	}
};

```

## 多边形

```c++
struct polygon 
{
 int n;
 Point p[N];
 Line l[N];
 void input(int _n) {
 n=_n;
 for (int i=0; i<n; i++)
 p[i].input();
 }
 void add(Point q) {
 p[n++]=q;
 }
 void getline() {
 for(int i=0; i<n; i++) {
 l[i]=Line(p[i],p[(i+1)%n]);
 }
 }
 struct cmp {
 Point p;
 cmp(const Point &p0) {
 p=p0;
 }
 bool operator()(const Point &aa,const Point &bb) {
 Point a=aa, b = bb;
 int d=sgn((a-p)^(b-p));
 if(d==0) {
 return sgn(a.distance(p)-b.distance(p))<0;
 }
 return d>0;
 }
 };
 //进行极角排序
 //首先需要找到最左下角的点
 //需要重载号好 Point 的 < 操作符 (min 函数要用)
 void norm() {
 Point mi=p[0];
 for (int i=1; i<n; i++) 
 mi=min(mi,p[i]);
 sort(p,p+n,cmp(mi));
 }
 //得到凸包
 //得到的凸包里面的点编号是 0～n-1 的
 //两种凸包的方法
 //注意如果有影响，要特判下所有点共点，或者共线的特殊情况
 void getconvex(polygon &convex) {
 sort(p,p+n);
 convex.n=n;
 for (int i=0; i<min(n,2); i++) {
 convex.p[i]=p[i];
 }
 if (convex.n==2 && (convex.p[0]==convex.p[1]))
 convex.n--;//特判
 if (n <= 2) return;
 int &top=convex.n;
 top=1;
 for (int i=2; i<n; i++) {
 while (top && sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<=0)
 top--;
 convex.p[++top]=p[i];
 }
 int temp=top;
 convex.p[++top]=p[n-2];
 for (int i=n-3; i>=0; i--) {
 while (top!=temp && sgn((convex.p[top]-p[i])^(convex.p[top-1]-p[i]))<=0)
 top--;
 convex.p[++top]=p[i];
 }
 if (convex.n==2 && (convex.p[0]==convex.p[1]))
 convex.n--;//特判
 convex.norm();//原来得到的是顺时针的点，排序后逆时针
 }
 //得到凸包的另外一种方法
 void Graham(polygon &convex) {
 norm();
 int &top=convex.n;
 top=0;
 if (n==1) {
 top=1;
 convex.p[0]=p[0];
 return;
 }
 if (n==2) {
 top=2;
 convex.p[0]=p[0];
 convex.p[1]=p[1];
 if(convex.p[0]==convex.p[1]) top--;
 return;
 }
 convex.p[0]=p[0];
 convex.p[1]=p[1];
 top=2;
 for (int i=2; i<n; i++) {
 while ( top>1 && sgn((convex.p[top-1]-convex.p[top-2])^(p[i]-convex.p[top-2]))<=0)
 top--;
 convex.p[top++]=p[i];
 }
 if (convex.n==2 && (convex.p[0]==convex.p[1]))
 convex.n--;//特 判
 }
 //判断是不是凸的
 bool isconvex() {
 bool s[2];
 memset(s,false,sizeof(s));
 for (int i=0; i<n; i++) {
 int j=(i+1)%n;
 int k=(j+1)%n;
 s[sgn((p[j]-p[i])^(p[k]-p[i]))+1]=true;
 if (s[0] && s[2]) return false;
 }
 return true;
 }
 //判断点和任意多边形的关系
 // 3 点上
 // 2 边上
 // 1 内部
 // 0 外部
 int relationpoint(Point q) {
 for (int i=0; i<n; i++) {
 if (p[i]==q) return 3;
 }
 getline();
 for (int i=0; i<n; i++) {
 if (l[i].pointonseg(q)) return 2;
 }
 int cnt=0;
 for (int i=0; i<n; i++) {
 int j=(i+1)%n;
 int k=sgn((q-p[j])^(p[i]-p[j]));
 int u=sgn(p[i].y-q.y);
 int v=sgn(p[j].y-q.y);
 if (k>0 && u<0 && v>=0) cnt++;
 if (k<0 && v<0 && u>=0) cnt--;
 }
 return cnt!=0;
 }
 //直线 u 切割凸多边形左侧
 //注意直线方向
 void convexcut(Line u,polygon &po) {
 int &top=po.n;//注意引用
 top=0;
 for (int i=0; i<n; i++) {
 int d1=sgn((u.e-u.s)^(p[i]-u.s));
 int d2=sgn((u.e-u.s)^(p[(i+1)%n]-u.s));
 if (d1>=0) po.p[top++]=p[i];
 if (d1*d2<0) 
 po.p[top++] = u.crosspoint(Line(p[i],p[(i+1)%n]));
 }
 }
 //得到周长
 double getcircumference() {
 double sum=0;
 for (int i=0; i<n; i++) {
 sum+=p[i].distance(p[(i+1)%n]);
 }
 return sum;
 }
 //得到面积
 double getarea() {
 double sum=0;
 for (int i=0; i<n; i++) {
 sum+=(p[i]^p[(i+1)%n]);
 }
 return fabs(sum)/2;
 }
 //得到方向
 // 1 表示逆时针，0 表示顺时针
 bool getdir() {
 double sum=0;
 for (int i=0; i<n; i++)
 sum+=(p[i]^p[(i+1)%n]);
 if (sgn(sum)>0) return 1;
 return 0;
 }
 //得到重心
 Point getbarycentre() {
 Point ret(0,0);
 double area=0;
 for (int i=1; i<n-1; i++) {
 double tmp=(p[i]-p[0])^(p[i+1]-p[0]);
 if (sgn(tmp)==0)continue;
 area+=tmp;
 ret.x+=(p[0].x+p[i].x+p[i+1].x)/3*tmp;
 ret.y+=(p[0].y+p[i].y+p[i+1].y)/3*tmp;
 }
 if (sgn(area)) ret=ret/area;
 return ret;
 }
};
```

## 半平面交

```c++
struct halfplane:public Line {
	double angle;
	halfplane() {}
	//表示向量 s->e 逆时针 (左侧) 的半平面
	halfplane(Point _s,Point _e) {
		s=_s;
		e=_e;
	}
	halfplane(Line v) {
		s=v.s;
		e=v.e;
	}
	void calcangle() {
		angle=atan2(e.y-s.y,e.x-s.x);
	}
	bool operator <(const halfplane &b)const {
		return angle<b.angle;
	}
};
struct halfplanes {
	int n;
	halfplane hp[maxp];
	Point p[maxp];
	int que[maxp];
	int st,ed;
	void push(halfplane tmp) {
		hp[n++]=tmp;
	}
	//去重
	void unique() {
		int m=1;
		for (int i=1; i<n; i++) {
			if (sgn(hp[i].angle-hp[i-1].angle)!=0)
				hp[m++]=hp[i];
			else if(sgn((hp[m-1].e-hp[m-1].s)^(hp[i].s-hp[m-1].s))>0)
				hp[m-1]=hp[i];
		}
		n=m;
	}
	bool halfplaneinsert() {
		for (int i=0; i<n; i++) hp[i].calcangle();
		sort(hp,hp+n);
		unique();
		que[st=0]=0;
		que[ed=1]=1;
		p[1]=hp[0].crosspoint(hp[1]);
		for (int i=2; i<n; i++) {
			while (st<ed && sgn((hp[i].e-hp[i].s)^(p[ed]-hp[i].s))<0)
				ed--;
			while (st<ed && sgn((hp[i].e-hp[i].s)^(p[st+1]-hp[i].s))<0)
				st++;
			que[++ed]=i;
			if (hp[i].parallel(hp[que[ed-1]])) return false;
			p[ed]=hp[i].crosspoint(hp[que[ed-1]]);
		}
		while (st<ed && sgn((hp[que[st]].e-hp[que[st]].s)^(p[ed]-hp[que[st]].s))<0)
			ed--;
		while (st<ed && sgn((hp[que[ed]].e-hp[que[ed]].s)^(p[st+1]-hp[que[ed]].s))<0)
			st++;
		if (st+1>=ed) return false;
		return true;
	}
	//得到最后半平面交得到的凸多边形
	//需要先调用 halfplaneinsert() 且返回 true
	void getconvex(polygon &con) {
		p[st]=hp[que[st]].crosspoint(hp[que[ed]]);
		con.n=ed-st+1;
		for(int j=st,i=0; j<=ed; i++,j++)
			con.p[i]=p[j];
	}
};

```

# 自己的板子

```c++
template<typename T> struct Point
{
 T x, y;
 int id;
 
 void input(){cin >> x >> y;}
 void output(string str){cout << str << ": " << x << "==" << y << endl;}
 bool operator==(const Point &a) const {return (fabs(x-a.x)<=eps && fabs(y-a.y)<=eps);}
 Point operator+(const Point &a) const {return {x+a.x,y+a.y};}
 Point operator-(const Point &a) const {return {x-a.x,y-a.y};}
 Point operator-() const {return {-x,-y};}
 Point operator*(const T k) const {return {k*x,k*y};}
 Point operator/(const T k) const {return {x/k,y/k};}
 T operator*(const Point &a) const {return x*a.x+y*a.y;}
 T operator^(const Point &a) const {return x*a.y-y*a.x;}
 bool is_ver(const Point &a) const {return fabs((*this)*a)<=eps;}
 bool is_par(const Point &a) const {return fabs((*this)^a)<=eps;}
 bool operator<(const Point &a) const {if (fabs(x-a.x)<=eps) return y<a.y-eps; return x<a.x-eps;}
 int toleft(const Point &a) const {auto t=(*this)^a; return (t>eps)-(t<-eps);}
 T len2() const {return (*this)*(*this);}
 double len() const {return sqrt(len2());}
 T dis2(const Point &a) const {return (a-(*this)).len2();}
 double dis(const Point &a) const {return (a-(*this)).len();}
 double ang(const Point &a) const {return acos(((*this)*a)/(this->len()*a.len()));}
 Point rot(const double rad) const {return {x*cos(rad)-y*sin(rad),x*sin(rad)+y*cos(rad)};}
 Point rot_val(const double s) const {double c = sqrt(1 - s * s);return {x*c-y*s,x*s+y*c};}
 int quad() const { return sign(y) == 1 || (sign(y) == 0 && sign(x) >= 0); }
 
 // 左边x轴和下半平面一个值，上平面和右x轴算一个值
 /* sort(p, p + n, [&](P a, P b) {
 int qa = a.quad(), qb = b.quad();
 if (qa != qb) return qa < qb;
 else return sign(a ^ b) > 0; // 若还要对第二关键字排序，可以在这里修改
 });*/
};
using P = Point<ty>;

template<typename T> struct Line
{
 Point<T> p, t, v; //p+kv
 
 bool operator==(const Line &a) const {return (v.is_par(a.v) && v.is_par(p-a.p));}
 bool is_par(const Line &a) const {return (v.is_par(a.v) && !v.is_par(p-a.p));}
 bool is_ver(const Line &a) const {return v.is_ver(a.v);}
 bool is_on_line(const Point<T> &a) const {return v.is_par(a-p);}
 bool is_on_seg(const Point<T> &a){return is_on_line(a) && (sign(((a-p)*(a-t))<=eps));}
 int toleft(const Point<T> &a) const {return v.toleft(a-p);}
 Point<T> inter(const Line &a) const {return p+v*((a.v^(p-a.p))/(v^a.v));}
 double dis(const Point<T> &a) const {return fabs(v^(a-p))/v.len();}
 double dis_seg(const P &q){if ((q-p)*(t-p)>=-eps && (q-t)*(p-t)>=-eps) return dis(q);return min(q.dis(p),q.dis(t));}
 Point<T> proj(const Point<T> &a) const {return p+v*((v*(a-p))/(v*v));}
 bool by_pass(Line<T> &a){return a.toleft(p)*a.toleft(t)<=0;} // 自己是线段，a是直线，不能4点共线
 bool by_pass_strict(Line<T> &a){return a.toleft(p)*a.toleft(t)<0;} 
 bool collinear(Line<T> &a){return (!a.toleft(p) && !a.toleft(t));} // 四点共线
 double len2(){return t.dis2(p);}
 double len(){return t.dis(p);}
 bool isSS(Line<T> &a){return intersect(p.x, t.x, a.p.x, a.t.x) && intersect(p.y, t.y, a.p.y, a.t.y) && a.pass(*(this)) && pass(a);}
 bool isSS_strict(Line<T> &a){return a.pass_strict(*(this)) && pass_strict(a);}
};
using L = Line<ty>;

template<typename T> struct Polygon
{
 bool convex_hull = false;
 vector<Point<T>> p;
 vector<Line<T>> l;
 vector<T> sum;

void input(int n){p.resize(n), l.resize(n);for (int i=0;i<n;i++)p[i].input();}
 void add(Point<T> pp){p.push_back(pp);}
 inline size_t nxt(const size_t i) const {return i==p.size()-1?0:i+1;} 
 inline size_t pre(const size_t i) const {return i==0?p.size()-1:i-1;} 
 double get_circ(){double res = 0;for (int i = 0; i < l.size(); i ++ ) res += l[i].len();return res;}
 void get_sum(){const auto &p=this->p;vector<T> a(p.size());for (size_t i=0;i<p.size();i++) a[i]=p[this->pre(i)]^p[i];sum.resize(p.size());partial_sum(a.begin(),a.end(),sum.begin());}
 T query_sum(const size_t l,const size_t r) const{const auto &p=this->p;if (l<=r) return sum[r]-sum[l]+(p[r]^p[l]);return sum[p.size()-1]-sum[l]+sum[r]+(p[r]^p[l]);}
 T query_sum() const {return sum.back();}
 
 // 闵可夫斯基和
 Polygon<T> operator+(Polygon<T> &c) 
 {
 auto &p=this->p;
 vector<L> e1(p.size()),e2(c.p.size()),edge(p.size()+c.p.size());
 vector res; res.reserve(p.size()+c.p.size());
 auto cmp=[](L &u,L &v) {return argcmp()(u.t-u.p,v.p-v.t);};
 for (size_t i=0;i<p.size();i++) e1[i]={p[i],p[this->nxt(i)]};
 for (size_t i=0;i<c.p.size();i++) e2[i]={c.p[i],c.p[c.nxt(i)]};
 rotate(e1.begin(),min_element(e1.begin(),e1.end(),cmp),e1.end());
 rotate(e2.begin(),min_element(e2.begin(),e2.end(),cmp),e2.end());
 merge(e1.begin(),e1.end(),e2.begin(),e2.end(),edge.begin(),cmp);
 const auto check=[](const vector &res,const P &u)
 {
 const auto back1=res.back(),back2=*prev(res.end(),2);
 return (back1-back2).toleft(u-back1)==0 && (back1-back2)*(u-back1)>=-eps;
 };
 auto u=e1[0].p+e2[0].p;
 for (const auto &v:edge)
 {
 while (res.size()>1 && check(res,u)) res.pop_back();
 res.push_back(u);
 u=u+v.t-v.p;
 }
 if (res.size()>1 && check(res,res[0])) res.pop_back();
 Polygon<T> tt;
 tt.p = res;
 return tt;
 }

// 回转数判断点与多边形关系
 pair<bool,int> winding(const Point<T> &a) const
 {
 int cnt=0;
 for (size_t i=0;i<p.size();i++)
 {
 Point<T> u=p[i],v=p[nxt(i)];
 if (abs((a-u)^(a-v))<=eps && (a-u)*(a-v)<=eps) return {true,0};
 if (abs(u.y-v.y)<=eps) continue;
 Line<T> uv={u,v-u};
 if (u.y<v.y-eps && uv.toleft(a)<=0) continue;
 if (u.y>v.y+eps && uv.toleft(a)>=0) continue;
 if (u.y<a.y-eps && v.y>=a.y-eps) cnt++;
 if (u.y>=a.y-eps && v.y<a.y-eps) cnt--;
 }
 return {false,cnt};
 }

void get_line()
 {
 l.resize(p.size());
 for (int i = 0; i < p.size(); i ++ )
 {
 Point<T> u = p[i], v = p[nxt(i)];
 l[i] = {u, v, v - u};
 }
 }
	
 // 严格凸包
 void convexHull() 
 {
 auto ps = p;
 int n = ps.size(); if(n <= 1) return;
 sort(ps.begin(), ps.end());
 vector<Point<T>> qs(n * 2); int k = 0;
 for (int i = 0; i < n; qs[k++] = ps[i++]) 
 while (k > 1 && crossOp(qs[k - 2], qs[k - 1], ps[i]) <= 0) --k;
 for (int i = n - 2, t = k; i >= 0; qs[k++] = ps[i--])
 while (k > t && crossOp(qs[k - 2], qs[k - 1], ps[i]) <= 0) --k;
 qs.resize(k - 1);
 p = qs;
 get_line();
 convex_hull = true;
 }

void convexHullNonStrict() 
 {
 auto ps = p;
 int n = ps.size(); if(n <= 1) return;
 sort(ps.begin(), ps.end());
 vector<Point<T>> qs(n * 2); int k = 0;
 for (int i = 0; i < n; qs[k++] = ps[i++]) 
 while (k > 1 && crossOp(qs[k - 2], qs[k - 1], ps[i]) < 0) --k;
 for (int i = n - 2, t = k; i >= 0; qs[k++] = ps[i--])
 while (k > t && crossOp(qs[k - 2], qs[k - 1], ps[i]) < 0) --k;
 qs.resize(k - 1);
 p = qs;
 get_line();
 convex_hull = true;
 }
	
 // 旋转卡壳
 // auto func = [&](P &a, P &b, P &c, int &ans) {ans = max({ans, c.dis2(a), c.dis2(b)});}; 
 template<typename F> void rotcaliper(F &func, int &ans)
 {
 if (!convex_hull) convexHull();

auto &p=this->p;

if (p.size() <= 2) 
 {
 if (p.size() == 2) ans = p[0].dis2(p[1]);
 return;
 }

for (size_t i=0,j=1;i<p.size();i++)
 {
 auto nxti=this->nxt(i);
 func(p[i],p[nxti],p[j], ans);
 while (fabs(cross(p[nxti], p[this->nxt(j)],p[i]))>=fabs(cross(p[nxti], p[j],p[i])))
 {
 j=this->nxt(j);
 func(p[i],p[nxti],p[j], ans);
 }
 }
 }
 
	// 平面最大最小三角形
 pair<T,T> minmax_triangle()
 {
 if (p.size()<=2) return {0,0};
 vector<pair<int,int>> evt;
 evt.reserve(p.size()*p.size());
 const T t=abs((p[0]-p[1])^(p[0]-p[2]));
 T maxans=t,minans=t;
 for (size_t i=0;i<p.size();i++)
 {
 for (size_t j=0;j<p.size();j++)
 {
 if (i==j || p[i]==p[j]) continue;
 evt.push_back({i,j});
 }
 }

sort(evt.begin(),evt.end(),[&](pair<int,int> &u, pair<int,int> &v)
 {
 Point<T> a=p[u.second]-p[u.first],b=p[v.second]-p[v.first];
 int qa = (Point<T>){-a.y, a.x}.quad(), qb = (Point<T>){-b.y, b.x}.quad();
 if (qa != qb) return qa < qb;
 else return sign(a ^ b) > 0; // 若还要对第二关键字排序，可以在这里修改
 });

vector<size_t> vx(p.size()),pos(p.size());
 for (size_t i=0;i<p.size();i++) vx[i]=i;
 sort(vx.begin(),vx.end(),[&](int x,int y){return p[x]<p[y];});
 for (size_t i=0;i<vx.size();i++) pos[vx[i]]=i;
 for (auto [u,v]:evt)
 {
 size_t i=pos[u],j=pos[v];
 if (i>j) swap(u,v),swap(i,j);
 Point<T> vecu=p[u],vecv=p[v];
 if (i>0) minans=min(minans,abs((p[vx[i-1]]-vecu)^(p[vx[i-1]]-vecv)));
 if (j<vx.size()-1) minans=min(minans,abs((p[vx[j+1]]-vecu)^(p[vx[j+1]]-vecv)));
 maxans=max({maxans,abs((p[vx[0]]-vecu)^(p[vx[0]]-vecv)),abs((p[vx.back()]-vecu)^(p[vx.back()]-vecv))});
 swap(vx[i],vx[j]);
 pos[u]=j; pos[v]=i;
 }
 return {minans,maxans};
 }
 
 // 最小矩形覆盖
 pair<double, vector> minRectangleCoverage()
 {
 if (!convex_hull) convexHull();

int n = p.size();

p.resize(p.size() * 2);
 for (int i = 0; i < n; i ++ ) p[i + n] = p[i];

int l, r, u;
        l = r = u = -1;

P d0 = p[1] - p[0];

for (int i = 0; i < n; i ++ )
 {
 if (l == -1 || (d0 * p[i]) <= (d0 * p[l])) l = i;
 if (r == -1 || (d0 * p[i]) >= (d0 * p[r])) r = i;
 if (u == -1 || (d0 ^ p[i]) >= (d0 ^ p[u])) u = i;
 }

double ans = DINF;
 vector ansP;

for (size_t i=0;i<n;i++)
 {
 int nxti = i + 1;
 P d = p[i + 1] - p[i];

auto check1 = [&](int l, P d) {return (d * p[l + 1]) <= (d * p[l]);};
 auto check2 = [&](int r, P d) {return (d * p[r + 1]) >= (d * p[r]);};
 auto check3 = [&](int u, P d) {return (d ^ p[u + 1]) >= (d ^ p[u]);};

while (l + 1 < n * 2 && check1(l, d)) l ++ ;
 while (r + 1 < n * 2 && check2(r, d)) r ++ ;
 while (u + 1 < n * 2 && check3(u, d)) u ++ ;

double s1 = d * (p[r] - p[l]) / d.len();
            double s2 = (d ^ (p[u] - p[i])) / d.len();

if (s1 * s2 < ans)
 {
 L ll = {p[l], p[l] + d.rot(PI / 2), d.rot(PI / 2)}, rr = {p[r], p[r] + d.rot(-PI / 2), d.rot(-PI / 2)}, uu = {p[u], p[u] - d, -d}, dd = {p[i], p[i + 1], d};
 P p1 = ll.inter(dd), p2 = dd.inter(rr), p3 = rr.inter(uu), p4 = uu.inter(ll);

ans = s1 * s2;
 ansP = vector{p1, p2, p3, p4};
 }
 }

rotate(ansP.begin(), min_element(ansP.begin(), ansP.end(), [&](P a, P b){if (!cmp(a.y, b.y)) return cmp(a.x, b.x) < 0;return cmp(a.y, b.y) < 0;}), 
 ansP.end());

return {ans, ansP};
    }
    // 2个多边形都要弄一遍 凸多边形间的最小距离
    void two_convex_mindis(Polygon &p2, double &ans)
    {
        auto func = [&](P a, P b, P c, double &ans){ans = min(ans, (L){a, b, b - a}.dis_seg(c));};

int j = 0;
 L t0 = {p[0], p[1], p[1] - p[0]};
 for (int i = 0; i < p2.p.size(); i ++ ) if ((cross(p[0], p2.p[i],p[1]))<(cross(p[0], p2.p[j],p[1]))) j = i;
 
 func(p[1],p[0],p2.p[j], ans);
 for (size_t i=0;i<p.size();i++)
 {
 auto nxti=this->nxt(i);
 func(p[i],p[nxti],p2.p[j], ans);
 while ((cross(p[i], p2.p[p2.nxt(j)],p[nxti]))<=(cross(p[i], p2.p[j],p[nxti])))
 {
 j=p2.nxt(j);
 func(p[i],p[nxti],p2.p[j], ans);
 }
 }
 }
 
 // -1 点在多边形边上 | 0 点在多边形外 | 1 点在多边形内
 int is_in(P &a)
 {
 if (!convex_hull) convexHull();
 auto &bp=this->p;
 if (bp.size()==1) return a==bp[0]?-1:0;
 if (bp.size()==2) return (L){bp[0],bp[1], bp[1] - bp[0]}.is_on_seg(a)?-1:0; 
 if (a==bp[0]) return -1;
 if ((bp[1]-bp[0]).toleft(a-bp[0])==-1 || (bp.back()-bp[0]).toleft(a-bp[0])==1) return 0;
 auto cmp=[&](P u,P v){return (u-bp[0]).toleft(v-bp[0])==1;};

size_t i=lower_bound(bp.begin()+1,bp.end(),a,cmp)-bp.begin();
 if (i==1) return (L){bp[0],bp[i], bp[i] - bp[0]}.is_on_seg(a)?-1:0;
 if (i==bp.size()-1 && (L){bp[0],bp[i], bp[i] - bp[0]}.is_on_seg(a)) return -1;
 if ((L){bp[i-1],bp[i], bp[i] - bp[i - 1]}.is_on_seg(a)) return -1;
 return (bp[i]-bp[i-1]).toleft(a-bp[i-1])>0;
 }
 template<typename F> size_t extreme(const F &dir)
 {
 const auto &p=this->p;
 const auto check=[&](const size_t i){return dir(p[i]).toleft(p[this->nxt(i)]-p[i])>=0;};
 const auto dir0=dir(p[0]); const auto check0=check(0);
 if (!check0 && check(p.size()-1)) return 0;
 const auto cmp=[&](const P &v)
 {
 const size_t vi=&v-p.data();
 if (vi==0) return 1;
 const auto checkv=check(vi);
 const auto t=dir0.toleft(v-p[0]);
 if (vi==1 && checkv==check0 && t==0) return 1;
 return checkv^(checkv==check0 && t<=0);
 };
 return partition_point(p.begin(),p.end(),cmp)-p.begin();
 }

// 过凸多边形外一点求凸多边形的切线，返回切点下标
 // 复杂度 O(logn)
 // 必须保证点在多边形外
 pair<size_t,size_t> tangent(const P &a) 
 {
 if (!convex_hull) convexHull();
 const size_t i=extreme([&](const P &u){return u-a;});
 const size_t j=extreme([&](const P &u){return a-u;});
 return {i,j};
 }

// 求平行于给定直线的凸多边形的切线，返回切点下标
 // 复杂度 O(logn)
 pair<size_t,size_t> tangent(const L &a)
 {
 if (!convex_hull) convexHull();
 const size_t i=extreme([&](...){return a.v;});
 const size_t j=extreme([&](...){return -a.v;});
 return {i,j};
 }
};
using Poly = Polygon<ty>;

long long MergeArea( Poly &poly, vector &vec)
{//把vec里的点插入凸包后，新凸包的面积
 Poly tmp;
 if (poly.p.size()<=10)
 {
 tmp.p=vec;
 tmp.p.insert(tmp.p.end(),poly.p.begin(),poly.p.end());
 tmp.convexHull(); 
 tmp.get_sum();
 return tmp.query_sum();
 }
 // vector tmp;
 for (auto u:vec)
 {
 if (poly.is_in(u)) continue;
 const auto [x,y]=poly.tangent(u);
 P px=poly.p[x],py=poly.p[y];
 u.id = -1;
 px.id=x; py.id=y;
 tmp.p.push_back(u);
 tmp.p.push_back(px);
 tmp.p.push_back(py);
 }
 if (tmp.p.empty()) return poly.query_sum();
 tmp.convexHull(); tmp.get_sum();
 long long ans=tmp.query_sum();
 for (size_t i=0;i<tmp.p.size();i++)
 {
 const auto j=tmp.nxt(i);
 if (tmp.p[i].id!=-1 && tmp.p[j].id!=-1)
 {
 ans+=poly.query_sum(tmp.p[i].id,tmp.p[j].id);
 }
 }
 return ans;
}

// 建立convex后 convexHull(), get_sum不能忘

struct argcmp
{
 bool operator () (const P & a, const P & b) const
 {
 const auto quad = [] (const P & a)
 {
 if(a.y < -eps) return 1;
 if(a.y > eps) return 4;
 if(a.x < -eps) return 5;
 if(a.x > eps) return 3;
 return 2;
 };
 const int qa = quad(a), qb = quad(b);
 if(qa != qb) return qa < qb;
 const auto t = a ^ b;
 if(abs(t) <= eps) return a * a eps;
 }
};
struct dynamic_convex
{
 set<P,argcmp> p; //坐标扩大三倍，便于整数运算
 P o; //凸包内一点
 double cir;
 
 decltype(p.begin()) nxt (decltype(p.begin()) it) const { it ++ ; return it == p.end() ? p.begin() : it; }
 decltype(p.begin()) pre (decltype(p.begin()) it) const { if(it == p.begin()) it = p.end(); return -- it; }
 size_t size () { return p.size(); }
 double circ () { return cir / 3; }
 
 // 判断是否在凸包内或凸包边界上
 bool is_in (const P & a) const
 {
 if(p.size() == 0) return 0;
 if(p.size() == 1) return *p.begin() == a * 3;
 if(p.size() == 2) return L{*p.begin(), *p.rbegin(), -(*p.begin()) + (*p.rbegin())}.is_on_seg(a * 3) ? 1 : 0;
 auto it = p.lower_bound(a * 3 - o);
 if(it == p.end()) it = p.begin();
 return ( (*it - *pre(it)) ^ ( (a * 3 - o) - *pre(it) ) ) >= -eps;
 }
 
 void insert (P a)
 {
 if(p.size() <= 1) { p.insert(a * 3); return ; }
 if(p.size() == 2) {
 P u = *p.begin(), v = *p.rbegin();
 if ((L){u,v, v - u}.is_on_seg(a*3)) return;
 if ((L){u,a *3, a*3 - u}.is_on_seg(v))
 {
 p.erase(v); p.insert(a*3);
 return;
 }
 if ((L){v,a*3, a*3 - v}.is_on_seg(u))
 {
 p.erase(u); p.insert(a*3);
 return;
 }
 o = (u + v + a * 3) / 3;
 p = {u - o, v - o, a * 3 - o};
 cir = (u - o).dis(v - o) + (u - o).dis(a * 3 - o) + (v - o).dis(a * 3 - o);
 return ;
 }
 if(is_in(a)) return ;
 a = a * 3 - o;
 auto it = p.insert(a).first;

auto nxtp = nxt(it);
 cir -= (*nxtp).dis(*pre(it));
 while (p.size() > 3 && ( (*nxtp - a) ^ (*nxt(nxtp) - *nxtp) ) <= eps){
 cir -= (*nxtp).dis(*nxt(nxtp));
 p.erase(nxtp);
 nxtp = nxt(it);
 }
 auto prep = pre(it);
 while (p.size() > 3 && ( (*prep - *pre(prep)) ^ (a - *prep) ) <= eps){
 cir -= (*prep).dis(*pre(prep));
 p.erase(prep);
 prep = pre(it);
 }
 cir += (*it).dis(*pre(it)) + (*it).dis(*nxt(it));
 }
};
 
dynamic_convex dc;

```

金铉杰/acm代码模板

简介

发行版

贡献者

近期动态

金铉杰/acm代码模板 .gitee-modal { width: 500px !important; }

简介

发行版

贡献者

近期动态

搜索帮助

金铉杰/acm代码模板