zoukankan html css js c++ java

计算几何基础模板

基础部分
凸包
极角排序
旋转卡壳
半平面交
最小圆覆盖

计算几何拖了这么久，终于拖到省选前了。
emm不求能A题，只求写暴力。

参考：
https://oi.men.ci/geometry-notes/
https://www.cnblogs.com/fly-in-milkyway/p/10569895.html
https://blog.csdn.net/clover_hxy/article/details/53966405
https://www.cnblogs.com/lstoi/p/9791654.html
https://blog.csdn.net/qq_34921856/article/details/80690822

ps: 刚开始写的代码，读入用的p[i]=(Point){read(),read()}，因为初始化了构造函数所以这么读进去横纵坐标可能是反的。。
有些题并没有影响所以以前的就不改了。

基础部分

#include <cmath>
#include <cstdio>
#include <cctype>
#include <vector>
#include <cstring>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=1e5+5;

const double eps=1e-10;

inline int dcmp(double x) {return fabs(x)<eps?0:x<0?-1:1;}

struct Vec
{
	double x,y;
	Vec(double x=0,double y=0):x(x),y(y) {}
	Vec operator +(const Vec &a)const {return Vec(x+a.x, y+a.y);}
	Vec operator -(const Vec &a)const {return Vec(x-a.x, y-a.y);}
	Vec operator *(const double p)const {return Vec(x*p, y*p);}

	double operator *(const Vec &a)const {return x*a.y-y*a.x;}//cross product
	bool operator <(const Vec &a)const {return x<a.x||(x==a.x&&y<a.y);}
	bool operator ==(const Vec &a)const {return !dcmp(x-a.x)&&!dcmp(y-a.y);}

	double Norm() {return x*x+y*y;}//范数 
	double Length() {return sqrt(x*x+y*y);}//模长 
	double Dot(Vec a) {return x*a.x+y*a.y;}//dot product
	double Angle(Vec a) {return acos(Dot(a)/Length()/a.Length());}//两向量夹角 
	Vec Normal() {double t=Length(); return Vec(-y/t,x/t);}//单位法向量 
	Vec Rotate(double rad) {return Vec(x*cos(rad)-y*sin(rad),x*sin(rad)+y*cos(rad));}//顺时针旋转rad度 
};
typedef Vec Point;

struct Line
{
	Point a; Vec v;
	Line(Point a,Vec v):a(a),v(v) {}

	bool OnLine(const Point &p) {return !dcmp((a-p)*v);}//!dcmp((a-p)*(b-p))
	bool OnSegment(const Point &p) {return !dcmp((a-p)*v)&&dcmp((a-p).Dot(a+v-p))<=0;}//PA*PB<=0
	int Relation(const Line &l)//直线之间的关系 0:平行 1:相交 2:重合(无数个交点) 
	{
		return dcmp(v*l.v)?1:dcmp(v*(a-l.a))?0:2;
	}
	Point Intersection(const Line &l)//直线交点 
	{
		return a+v*((l.v*(a-l.a))/(v*l.v));//注意方向（叉积的正负） 
	}
};

inline bool cmp(const Point &a,const Point &b) {return a.x==b.x?a.y<b.y:a.x<b.x;}

struct Polygon
{
	int sk[N];
	std::vector<Point> ps;

	bool Include(const Point &p)//点在多边形内 
	{
		int cnt=0;
		for(int i=0,lim=ps.size(); i<lim; ++i)
		{
			const Point a=ps[i],b=ps[i+1==lim?0:i+1];
			if(Line(a,b-a).OnSegment(p)) return 1;
			double d1=a.y-p.y,d2=b.y-p.y,tmp=(a-p)*(b-p);
			if((tmp<0&&d1<0&&d2>=0)||(tmp>0&&d1>=0&&d2<0)) ++cnt;
		}
		return cnt&1;
	}
	double Area()//多边形有向面积（逆时针为正，顺时针为负） 
	{
		double res=0;
		for(int i=0,lim=ps.size(); i<lim; ++i)
			res+=ps[i]*ps[i+1==lim?0:i+1];
		return res*0.5;
	}
	int Convex()//求凸包 存在sk[]里 
	{
		std::sort(ps.begin(),ps.end(),cmp);
		int top=1,n=ps.size(); sk[1]=0;
		for(int i=1; i<n; ++i)
		{
			while(top>=2 && (ps[sk[top]]-ps[sk[top-1]])*(ps[i]-ps[sk[top-1]])<=0) --top;
			sk[++top]=i;
		}
		int k=top;
		for(int i=n-2; ~i; --i)
		{
			while(top>k && (ps[sk[top]]-ps[sk[top-1]])*(ps[i]-ps[sk[top-1]])<=0) --top;
			sk[++top]=i;
		}
		return top;
	}
};

int main()
{
	return 0;
}

凸包

洛谷.2742.[模板]二维凸包。

(Granham's Scan)：

选出所有点中横坐标最小（如果相同取纵坐标最小的）的点(x)，将所有点按与(x)的极角序排序。
将(x)放到栈里，跑一遍单调栈求凸包。

这个做法似乎精度不太好。

另一种做法是：

将所有点按横坐标为第一关键字，纵坐标为第二关键字排序。
将(1)号点放到栈里，单调栈求一遍下凸壳。
保留原先栈中的元素，单调栈求一遍上凸壳。

复杂度都是(O(nlog n))。

#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
//#define gc() getchar()
#define MAXIN 500000
#define gc() (SS==TT&&(TT=(SS=IN)+fread(IN,1,MAXIN,stdin),SS==TT)?EOF:*SS++)
typedef long long LL;
const int N=1e4+6;

int sk[N];
char IN[MAXIN],*SS=IN,*TT=IN;
struct Vec
{
	double x,y;
	Vec(double x=0,double y=0):x(x),y(y) {}
	Vec operator +(const Vec &a)const {return Vec(x+a.x, y+a.y);}
	Vec operator -(const Vec &a)const {return Vec(x-a.x, y-a.y);}
	Vec operator *(const double a)const {return Vec(x*a, y*a);}
	bool operator <(const Vec &a)const {return x<a.x||(x==a.x&&y<a.y);}
	double operator *(const Vec &a)const {return x*a.y-y*a.x;}
	double Length() {return sqrt(x*x+y*y);}
}p[N];
typedef Vec Point;

inline double read()
{
	double x=0,y=0.1,f=1;register char c=gc();
	for(;!isdigit(c);c=='-'&&(f=-1),c=gc());
	for(;isdigit(c);x=x*10+c-48,c=gc());
	for(c=='.'&&(c=gc());isdigit(c);x+=y*(c-48),y*=0.1,c=gc());
	return x*f;
}
double Convex(int n)
{
	std::sort(p+1,p+1+n);
	int top=1; sk[1]=1;
	for(int i=2; i<=n; ++i)
	{
		while(top>1 && (p[sk[top]]-p[sk[top-1]])*(p[i]-p[sk[top-1]])<=0) --top;
		sk[++top]=i;
	}
	int m=top;
	for(int i=n-1; i; --i)
	{
		while(top>m && (p[sk[top]]-p[sk[top-1]])*(p[i]-p[sk[top-1]])<=0) --top;
		sk[++top]=i;
	}
	double ans=0;
	for(int i=1; i<top; ++i) ans+=(p[sk[i+1]]-p[sk[i]]).Length();
	return ans;
}

int main()
{
	int n=read();
	for(int i=1; i<=n; ++i) p[i]=(Point){read(),read()};
	printf("%.2f
",Convex(n));

	return 0;
}

极角排序

一道例题。

旋转卡壳

BZOJ.1069.[SCOI2007]最大土地面积

求四边形最大面积，枚举一条对角线，求对角线两边最大三角形面积即可。
暴力是(n^3)的，容易发现固定对角线一个端点，另一个端点移动时两边最远点的移动也是单调的。所以用旋转卡壳优化一下就是(O(n^2))啦。

旋转卡壳还有很多用途，见最上面的链接。
经常要找和当前直线平行的直线切在哪里，注意用底相同，平行线之间高也相同，以及叉积有正负的性质，求叉积。
求凸多边形最小面积外接矩形，找三个边界点时（枚举下边界点(d)），上边界点(u)可以如上所说用叉积判（判断((p_u-p_d)*(p_{d+1}-p_d))与((p_{u+1}-p_d)*(p_{d+1}-p_d))的大小关系，比它小则(u+1)）。
点积可以求向量在另一个向量上的映射长度，所以右边界点(r)可以用点积求（在(p_{d+1}-p_d)上映射长度最长的）。
同理，左边界点(l)也用点积求，不过注意刚开始要将(l)赋值为(r)。

[BZOJ1069]的代码：

//884kb	120ms
#include <cstdio>
#include <algorithm>
typedef long long LL;
const int N=2005;

struct Vec
{
	double x,y;
	Vec(double x=0,double y=0):x(x),y(y) {}
	Vec operator +(const Vec &a)const {return Vec(x+a.x, y+a.y);}
	Vec operator -(const Vec &a)const {return Vec(x-a.x, y-a.y);}
	double operator *(const Vec &a)const {return x*a.y-y*a.x;}
	bool operator <(const Vec &a)const {return x<a.x||(x==a.x&&y<a.y);}
}p[N],sk[N];


int main()
{
	int n; scanf("%d",&n);
	for(int i=1; i<=n; ++i) scanf("%lf%lf",&p[i].x,&p[i].y);
	std::sort(p+1,p+1+n);
	int top=1; sk[1]=p[1];
	for(int i=2; i<=n; ++i)
	{
		while(top>1 && (sk[top]-sk[top-1])*(p[i]-sk[top-1])<=0) --top;
		sk[++top]=p[i];
	}
	int k=top;
	for(int i=n-1; i>1; --i)
	{
		while(top>k && (sk[top]-sk[top-1])*(p[i]-sk[top-1])<=0) --top;
		sk[++top]=p[i];
	}
	double ans=0; sk[top+1]=sk[1];
	#define Inc(x) (x+1>top?1:x+1)
	for(int i=1; i+2<=top; ++i)
	{
		int x=i+1,y=i+3>top?1:i+3;
		for(int j=i+2; j<=top; ++j)
		{
			while((sk[x+1]-sk[i])*(sk[j]-sk[i])>(sk[x]-sk[i])*(sk[j]-sk[i])) x=Inc(x);
			while((sk[j]-sk[i])*(sk[y+1]-sk[i])>(sk[j]-sk[i])*(sk[y]-sk[i])) y=Inc(y);
			ans=std::max(ans,(sk[x]-sk[i])*(sk[j]-sk[i])+(sk[j]-sk[i])*(sk[y]-sk[i]));
		}
	}
	printf("%.3f
",ans*0.5);

	return 0;
}

半平面交

咕咕了。

最小圆覆盖

具体见这里。

#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
#define MAXIN 500000
//#define gc() (SS==TT&&(TT=(SS=IN)+fread(IN,1,MAXIN,stdin),SS==TT)?EOF:*SS++)
typedef long long LL;
const int N=1e6+5;

char IN[MAXIN],*SS=IN,*TT=IN;
struct Vec
{
	double x,y;
	Vec(double x=0,double y=0):x(x),y(y) {}
	Vec operator +(const Vec &a)const {return Vec(x+a.x, y+a.y);}
	Vec operator -(const Vec &a)const {return Vec(x-a.x, y-a.y);}
	Vec operator *(const double a)const {return Vec(x*a, y*a);}
	double operator *(const Vec &a)const {return x*a.y-y*a.x;}
	Vec Rotate_90()const {return Vec(y,-x);}
	double len()const {return sqrt(x*x+y*y);}
	double len2()const {return x*x+y*y;}
}p[N];
typedef Vec Point;
struct Line
{
	Point p; Vec v;
	Line(Point p,Vec v):p(p),v(v) {}
	Line PerpendicularBisector()const//垂直平分线=-=
	{
		return Line((p+p+v)*0.5,v.Rotate_90());
	}
	Point Intersection(const Line &l)const
	{
		return p+v*((l.v*(p-l.p))/(v*l.v));
	}
};

inline double read()
{
	double x=0,y=0.1,f=1;register char c=gc();
	for(;!isdigit(c);c=='-'&&(f=-1),c=gc());
	for(;isdigit(c);x=x*10+c-48,c=gc());
	for(c=='.'&&(c=gc());isdigit(c);x+=y*(c-48),y*=0.1,c=gc());
	return x*f;
}
Point CalcCircle(const Point &a,const Point &b,const Point &c)
{
//	Line A=Line(a,b-a).PerpendicularBisector(),B=Line(a,c-a).PerpendicularBisector();
	Line A=Line((a+b)*0.5,(b-a).Rotate_90()),B=Line((a+c)*0.5,(c-a).Rotate_90());
	return A.Intersection(B);
}
void Solve(const int n)
{
	srand(330), std::random_shuffle(p+1,p+1+n);//话说这个srand不够随机啊= = 
	Point O=p[1]; double R=0;
	for(int i=2; i<=n; ++i)
		if((p[i]-O).len2()>R)
		{
			O=p[i], R=0;
			for(int j=1; j<i; ++j)
				if((p[j]-O).len2()>R)
				{
					O=(p[i]+p[j])*0.5, R=(p[i]-O).len2();
					for(int k=1; k<j; ++k)
						if((p[k]-O).len2()>R)
							O=CalcCircle(p[i],p[j],p[k]), R=(p[k]-O).len2();
				}
		}
	printf("%.2f %.2f %.2f
",O.x,O.y,sqrt(R));
}

int main()
{
	int n=read();
	for(int i=1; i<=n; ++i) p[i].x=read(),p[i].y=read();
//	for(int i=1; i<=n; ++i) p[i]=(Point){read(),read()};//声明构造函数之后再这么用，貌似。。= = 不同编译器结果不同。。
	Solve(n);

	return 0;
}

查看全文

相关阅读:
三大主流负载均衡软件对比（LVS+Nginx+HAproxy）
nginx 提示the "ssl" directive is deprecated, use the "listen ... ssl" directive instead
centos安装nginx并配置SSL证书
 hadoop创建目录文件失败
 The server time zone value 'EDT' is unrecognized or represents more than one time zone.
脚本启动SpringBoot（jar)
centos做免密登录
 数据库远程连接配置
 Bash 快捷键
 TCP三次握手四次断开

原文地址：https://www.cnblogs.com/SovietPower/p/10593432.html

计算几何基础 模板