计算几何模板

本文深入探讨了几何算法中的关键问题,包括凸包构建、判断点是否位于多边形内部及半平面交集等,并提供了详细的代码实现。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >


 

poj1113(凸包)

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
const double eps=1e-8;
const double pi=acos(-1.0);

inline int dcmp(double x) 
{
	if (fabs(x)<eps) return 0;
	return x<0 ? -1:1;
}
struct vec
{
	double x,y;
	vec (double X=0,double Y=0)
	{
		x=X,y=Y;
	}
};
typedef vec point;
vec operator +(vec a,vec b) {return vec(a.x+b.x,a.y+b.y);}
vec operator -(vec a,vec b) {return vec(a.x-b.x,a.y-b.y);}
inline double cross(vec a,vec b) {return a.x*b.y-a.y*b.x;}
inline double dot(vec a,vec b) {return a.x*b.x+a.y*b.y;} 
inline double len(vec a) {return sqrt(dot(a,a));}

int n,num;
point p[1005],poly[1005];double L;
bool cmp(point a,point b)
{
	if (a.x==b.x) return a.y<b.y;
	return a.x<b.x;
}

void tubao()
{
	sort(p+1,p+n+1,cmp);
	for (int i=1;i<=n;i++)
	{
		while (num>1&&dcmp(cross(p[i]-poly[num],poly[num]-poly[num-1]))<=0) num--;
		poly[++num]=p[i];
	}
	int k=num;
	for (int i=n-1;i>=1;i--) 
	{
		while (num>k&&dcmp(cross(p[i]-poly[num],poly[num]-poly[num-1]))<=0) num--;
		poly[++num]=p[i];
	}
	if (num>1) num--;
}

int main()
{
	scanf("%d%lf",&n,&L);
	for (int i=1;i<=n;i++) scanf("%lf%lf",&p[i].x,&p[i].y);	
	tubao();
	double ans=2.0*pi*L;
	for (int i=1;i<=num;i++) ans+=len(poly[i%num+1]-poly[i]);
	printf("%.0lf",ans);
	return 0;
}


 

poj1584(点(圆)在多边形内)

#include<cstdio>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
const double eps=1e-8;

inline int dcmp(double x) 
{
	if (fabs(x)<eps) return 0;
	return x<0 ? -1:1;
}
struct vec
{
	double x,y;
	vec (double X=0,double Y=0)
	{
		x=X,y=Y;
	}
};
typedef vec point;
vec operator +(vec a,vec b) {return vec(a.x+b.x,a.y+b.y);}
vec operator -(vec a,vec b) {return vec(a.x-b.x,a.y-b.y);}
inline double cross(vec a,vec b) {return a.x*b.y-a.y*b.x;}
inline double dot(vec a,vec b) {return a.x*b.x+a.y*b.y;} 
inline double len(vec a) {return sqrt(dot(a,a));}











int n;
point poly[100005],O;double R;

bool pan() 
{
	point p1,p2,p3;int ty=0,tmp;
	for (int i=1;i<=n;i++)
	{
		p1=poly[i],p2=poly[i%n+1],p3=poly[(i%n+1)%n+1];
		tmp=dcmp(cross(p2-p1,p3-p2));
		if (tmp==0) continue;
		
		if (ty==0) ty=tmp;
		else if (ty!=tmp) return false;
	}
	return true;
}
double dis(point p,point p1,point p2)
{
	vec v1=p-p1,v2=p2-p1;
	return fabs(cross(v2,v1)/len(v2));
}
bool inpoly()
{
	int cnt=0,y1,y2,k;point p1,p2;
	for (int i=1;i<=n;i++)
	{
		p1=poly[i],p2=poly[i%n+1];
		if (dcmp(dis(O,p1,p2)-R)<0) return false;
		
		y1=dcmp(p1.y-O.y),y2=dcmp(p2.y-O.y);
		k=dcmp(cross(p2-p1,O-p1));
		
		if (k>0&&y2>=0&&y1<0) cnt++;
		if (k<0&&y1>=0&&y2<0) cnt--;
	}
	return cnt;
}

int main()
{
	while (scanf("%d",&n)&&n>=3) 
	{
		scanf("%lf%lf%lf",&R,&O.x,&O.y);
		for (int i=1;i<=n;i++) scanf("%lf%lf",&poly[i].x,&poly[i].y);
		
		if (!pan()) printf("HOLE IS ILL-FORMED\n");
		else
		if (inpoly()) printf("PEG WILL FIT\n");else printf("PEG WILL NOT FIT\n");
	}
	return 0;
}


 n^2半平面交,就是一个一个插入半平面,同时维护当前的凸包。。。

本题维护的是右半平面。。

poj3384

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<cstdlib>
using namespace std;
const double eps=1e-8;
const double inf=1e9;
const double pi=acos(-1.0);

inline int dcmp(double x) 
{
	if (fabs(x)<eps) return 0;
	return x<0? -1:1;
}
struct vec
{
	double x,y;
	vec(double X=0,double Y=0)
	{
		x=X,y=Y;
	}
};
typedef vec point;
vec operator +(vec a,vec b) {return vec(a.x+b.x,a.y+b.y);}
vec operator -(vec a,vec b) {return vec(a.x-b.x,a.y-b.y);}
vec operator *(vec a,double t) {return vec(a.x*t,a.y*t);}
vec operator /(vec a,double t) {return vec(a.x/t,a.y/t);}

inline double dot(vec a,vec b) {return a.x*b.x+a.y*b.y;}
inline double cross(vec a,vec b) {return a.x*b.y-a.y*b.x;}
inline double len(vec a) {return sqrt(dot(a,a));}

point GIL(point p1,vec v1,point p2,vec v2)
{
	vec w=p1-p2;
	double t=cross(v2,w)/cross(v1,v2);
	return p1+v1*t;
}
bool jiaodian(point A,point B,point C,point D)
{
	vec v1=B-A,v2=C-A,v3=D-A;
	return dcmp(cross(v1,v2))!=dcmp(cross(v1,v3));
}
vec rotate(vec v,double rad)
{
	return vec(v.x*cos(rad)-v.y*sin(rad),v.y*cos(rad)+v.x*sin(rad));
}
int n;double r;

point poly[105],polyn[105],p[105];int cnt,cntn;
void init()
{
	poly[++cnt]=point(inf,inf);
	poly[++cnt]=point(inf,-inf);
	poly[++cnt]=point(-inf,-inf);
	poly[++cnt]=point(-inf,inf);
}
void half(point A,point B)
{
	cntn=0;
	for (int i=1;i<=cnt;i++)
	{
		point C=poly[i],D=poly[i%cnt+1];
		if (dcmp(cross(B-A,C-A))<=0) 
				polyn[++cntn]=C;
		if (jiaodian(A,B,C,D)) 
				polyn[++cntn]=GIL(A,B-A,C,D-C);
	}
	cnt=cntn;
	for (int i=1;i<=cnt;i++) poly[i]=polyn[i];
}

int main()
{
	scanf("%d%lf",&n,&r);
	for (int i=1;i<=n;i++) scanf("%lf%lf",&p[i].x,&p[i].y);
	init();
	for (int i=1;i<=n;i++)
	{
		point A=p[i],B=p[i%n+1];
		vec v=(B-A)/len(B-A)*r;
		v=rotate(v,-pi/2.0);
		A=A+v,B=B+v;
		half(A,B);
	}
	
	double ans=0;int ansi=1,ansj=1;
	for (int i=1;i<=cnt;i++)
		for (int j=i+1;j<=cnt;j++) 
		{
			double dis=len(poly[i]-poly[j]);
			if (dcmp(dis-ans)>0) 
			{
				ans=dis;
				ansi=i,ansj=j;
			}
		}
	printf("%.4lf %.4lf %.4lf %.4lf",poly[ansi].x,poly[ansi].y,poly[ansj].x,poly[ansj].y);
	return 0;
}


 bzoj2732

二分,抽象为一次函数不等式问题,转化为半平面求交。

这里用的n log n的半平面求交,极角排序,然后依次加入,双端队列维护

#include<cstdio>
#include<cstring>
#include<cmath>
#include<cstdlib>
#include<algorithm>
using namespace std;
const int N=200015;
const long double inf=1e11;
const long double eps=1e-12;

inline int dcmp(long double x) 
{
	if (fabs(x)<eps) return 0;
	return x<0 ? -1:1;
}
struct vec
{
	long double x,y;
	vec(long double X=0,long double Y=0)
	{
		x=X,y=Y;
	} 
};
typedef vec point;
struct line
{
	point p;
	vec v;
	long double ang;
	line(point P=point(0,0),vec V=vec(0,0))
	{
		p=P,v=V;
		ang=atan2(v.y,v.x);
	}
	bool operator <(const line &b) const 
	{
		return ang<b.ang;
	}
};
vec operator +(vec a,vec b) {return vec(a.x+b.x,a.y+b.y);}
vec operator -(vec a,vec b) {return vec(a.x-b.x,a.y-b.y);}
vec operator *(vec a,long double t) {return vec(a.x*t,a.y*t);}

inline long double dot(vec a,vec b) {return a.x*b.x+a.y*b.y;}
inline long double cross(vec a,vec b) {return a.x*b.y-a.y*b.x;}

inline bool onleft(line l,point p)
{
	return dcmp(cross(l.v,p-l.p))>=0;
}
point GIL(line a,line b) 
{
	vec v=a.p-b.p;
	long double t=cross(b.v,v)/cross(a.v,b.v);
	return a.p+a.v*t;
}

int n;
line s[N],L[N],tmp[N];int num;
point p[N];

bool pan(int mid)
{
	for (int i=1;i<=mid;i++) L[i]=s[i];
	sort(L+1,L+mid+1);
	int l,r;
	tmp[l=r=1]=L[1];
	for (int i=2;i<=mid;i++)
	{
		while (l<r&&!onleft(L[i],p[r-1])) r--;
		while (l<r&&!onleft(L[i],p[l])) l++;
		tmp[++r]=L[i];
		if (dcmp(cross(L[i].v,tmp[r-1].v))==0)
		{
			r--;
			if (onleft(tmp[r],L[i].p)) tmp[r]=L[i];
		}
		if (l<r) p[r-1]=GIL(tmp[r],tmp[r-1]);
	}
	while (l<r&&!onleft(tmp[l],p[r-1])) r--;
	if (r-l<=1) return false;
	return true; 
}

int main()
{
	scanf("%d",&n);
	s[++num]=line(point(0,0),vec(0,inf));
	s[++num]=line(point(0,inf),vec(-inf,0));
	s[++num]=line(point(-inf,inf),vec(0,-inf));
	s[++num]=line(point(-inf,0),vec(inf,0));
	
	long double tmp1=-1.0,tmp2=4.0;
	double x,y1,y2;point P,Q;
	for (int i=1;i<=n;i++)
	{
		scanf("%lf%lf%lf",&x,&y1,&y2);
		if (x<=0) break;
		
		P=point(tmp1,(y2-tmp1*x*x)/x);
		Q=point(tmp2,(y2-tmp2*x*x)/x);
		s[++num]=line(Q,P-Q);
		
		P=point(tmp1,(y1-tmp1*x*x)/x);
		Q=point(tmp2,(y1-tmp2*x*x)/x);
		s[++num]=line(P,Q-P);
	}
	
	int mid,l=5,r=num,ans=4;
	while (l<=r)
	{
		mid=(l+r)>>1;
		if (pan(mid)) l=mid+1,ans=mid;else r=mid-1;
	}
	printf("%d",(ans-4)/2);
	return 0;
}


 

 

 

 

 

 

 

 

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值