第九届河南理工大学算法程序设计大赛 I. 安全距离 (计算几何)

 题目链接

在一个三维空间中给定一个三角形薄片和一个球形安全区，现在需要求出这个三角形薄片到安全区的最短距离。若三角形薄片上的某点在安全区内，则最短距离为零。

输入

第一行一个整数 ?，（ ?≤103）。
之后的 ? 行，每行 13 个整数 ?1 ?1 ?1 ?2 ?2 ?2 ?3 ?3 ?3 ?? ?? ?? ?。
(?1,?1,?1),(?2,?2,?2),(?3,?3,?3) 表示三角形薄片的三个顶点，(??????),? 表示球的球心和半径。
数据保证所有的数均不超过 103 。

输出

输出 ? 行，每行一个实数为所求答案，请保留6位小数（四舍五入）。

样例

input

3
0 0 0 0 0 1 0 1 0 2 0 0 1
1 0 5 0 0 1 0 1 0 2 0 0 1
0 0 0 0 0 1 0 1 0 0 0 0 1

output

1.000000
1.121320
0.000000

 简化一下题意：

求空间中一个点到一个三角形的最短距离 

大致思路分为两种情况

一： 过这个点向三角形所在平面引引垂线，如果垂足在三角形封闭区域内，那么最短距离即为这个点到垂足的距离，两点距离就很好求了

二：如果垂足不在三角形的封闭区域内，那么这个点到三角形的最短距离即为这个点到三条线段的最短距离中最小的那个

按照这个思路就很好写了，求点到线段的最短距离可以用三分，也可以用计算几何公式

 

公式解法，复杂度比较低

#include <bits/stdc++.h>
using namespace std;
#define eps 1e-8
inline double max(double a,double b)
{
	return a > b ? a : b;
}
inline double min(double a,double b)
{
	return a > b ? b : a;
}
struct point
{
	double x,y,z;
	point(){}
	point(double _x,double _y,double _z)
	{
		x=_x;
		y=_y;
		z=_z;
	}
};
double dist(point p1,point p2)
{
    double ans=sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y)+(p1.z-p2.z)*(p1.z-p2.z));
    return ans;
}
double point_line(point p,point a,point b)
{
	double t=((a.x-p.x)*(a.x-b.x)+(a.y-p.y)*(a.y-b.y)+(a.z-p.z)*(a.z-b.z))/((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y)+(a.z-b.z)*(a.z-b.z));
	point q;
	q.x=a.x+t*(b.x-a.x);
	q.y=a.y+t*(b.y-a.y);
	q.z=a.z+t*(b.z-a.z);
	if(t>=0&&t<=1)
		return dist(p,q);
	else if(t<0)
		return dist(p,a);
	else
		return dist(p,b);
}
void solve(double R,point l,point r,point m,point p)
{
	double A=(r.y-l.y)*(m.z-l.z)-(r.z-l.z)*(m.y-l.y);
	double B=(m.x-l.x)*(r.z-l.z)-(r.x-l.x)*(m.z-l.z);
	double C=(r.x-l.x)*(m.y-l.y)-(m.x-l.x)*(r.y-l.y);
	double t=(A*(l.x-p.x)+B*(l.y-p.y)+C*(l.z-p.z))/(A*A+B*B+C*C);
	point Q;
	Q.x=p.x+t*A;
	Q.y=p.y+t*B;
 	Q.z=p.z+t*C;
	double dis1=Q.x*(r.y-m.y)-Q.y*(r.x-m.x)+(r.x*m.y-r.y*m.x);
	double dis2=l.x*(Q.y-m.y)-l.y*(Q.x-m.x)+(Q.x*m.y-Q.y*m.x);
	double dis3=l.x*(r.y-Q.y)-l.y*(r.x-Q.x)+(r.x*Q.y-Q.x*r.y);
	if((fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(A))<eps)||(fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(B))<eps)||(fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(C))<eps))
	{
		double ans1=dist(Q,p);
        printf("%.6lf\n",max(0,(ans1-R)));
		return ;
	}
	double ans1=(point_line(p,l,r)-R);
	double ans2=(point_line(p,l,m)-R);
	double ans3=(point_line(p,r,m)-R);
	double out=min(ans1,min(ans2,ans3));
	printf("%.6f\n",max(0,out));
}
int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
    	point p,l,r,m;double R;
        scanf("%lf%lf%lf",&l.x,&l.y,&l.z);
        scanf("%lf%lf%lf",&r.x,&r.y,&r.z);
        scanf("%lf%lf%lf",&m.x,&m.y,&m.z);
        scanf("%lf%lf%lf",&p.x,&p.y,&p.z);
        scanf("%lf",&R);
        solve(R, l, r, m, p);
    }
    return 0;
}

三分解法

#include <bits/stdc++.h>
using namespace std;
#define max(a,b) a>b?a:b
#define eps 1e-8
struct point{//点结构 
	double x,y,z;
};
int sgn(double a)
{
	return (a>eps)-(a<-eps);
}
point p,l,r,m,p1,p2;
point Q;// 垂足
double d1,d2;
double R;
double dist(point p1,point p2){//计算两点间距离 
    double ans=sqrt((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y)+(p1.z-p2.z)*(p1.z-p2.z));
    return ans;
}
void solve()
{
	double A=(r.y-l.y)*(m.z-l.z)-(r.z-l.z)*(m.y-l.y);
	double B=(m.x-l.x)*(r.z-l.z)-(r.x-l.x)*(m.z-l.z);
	double C=(r.x-l.x)*(m.y-l.y)-(m.x-l.x)*(r.y-l.y);
	double D=-(A*l.x+B*l.y+C*l.z);
	double d=(fabs(A*p.x+B*p.y+C*p.z+D)/sqrt(A*A+B*B+C*C));
	double t=(A*(l.x-p.x)+B*(l.y-p.y)+C*(l.z-p.z))/(A*A+B*B+C*C);
	Q.x=p.x+t*A;
	Q.y=p.y+t*B;
 	Q.z=p.z+t*C;
	double dis1=Q.x*(r.y-m.y)-Q.y*(r.x-m.x)+(r.x*m.y-r.y*m.x);
	double dis2=l.x*(Q.y-m.y)-l.y*(Q.x-m.x)+(Q.x*m.y-Q.y*m.x);
	double dis3=l.x*(r.y-Q.y)-l.y*(r.x-Q.x)+(r.x*Q.y-Q.x*r.y);
	if((fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(A))<eps)||(fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(B))<eps)||(fabs((fabs(dis1)+fabs(dis2)+fabs(dis3))-fabs(C))<eps))
	{
		double ans1=sqrt((Q.x-p.x)*(Q.x-p.x)+(Q.y-p.y)*(Q.y-p.y)+(Q.z-p.z)*(Q.z-p.z));
        printf("%.6lf\n",max(0,(ans1-R)));
		//if((d-R)<0) printf("0.000000\n");
		//else printf("%.6f\n", (d-R));
		return ;
	}
	point aa,bb,cc;
	point dd,ee,ff;
	point gg,hh,ii;
	aa=l,bb=r,cc=m;
	dd=l,ee=r,ff=m;
	gg=l,hh=r,ii=m;
	for(int i=0;i<100;i++)
            p1.x=(l.x+r.x)/2.0;
            p1.y=(l.y+r.y)/2.0;
            p1.z=(l.z+r.z)/2.0;
 
            p2.x=(p1.x+r.x)/2.0;
            p2.y=(p1.y+r.y)/2.0;
            p2.z=(p1.z+r.z)/2.0;
            d1=dist(p1,p);
	       d2=dist(p2,p);
            if(d1<d2) r=p2;
            else l=p1;
    }
    double vis1=dist(p,l);
    for(int i=0;i<100;i++)
    {
            p1.x=(dd.x+ff.x)/2.0;
            p1.y=(dd.y+ff.y)/2.0;
            p1.z=(dd.z+ff.z)/2.0;
 
            p2.x=(p1.x+ff.x)/2.0;
            p2.y=(p1.y+ff.y)/2.0;
            p2.z=(p1.z+ff.z)/2.0;
            d1=dist(p1,p);
	       d2=dist(p2,p);
            if(d1<d2) ff=p2;
            else dd=p1;
    }
    double vis2=dist(p,dd);
    for(int i=0;i<100;i++)
    {
            p1.x=(hh.x+ii.x)/2.0;
            p1.y=(hh.y+ii.y)/2.0;
            p1.z=(hh.z+ii.z)/2.0;
            p2.x=(p1.x+ii.x)/2.0;
            p2.y=(p1.y+ii.y)/2.0;
            p2.z=(p1.z+ii.z)/2.0;
            d1=dist(p1,p);
	        d2=dist(p2,p);
            if(d1<d2) ii=p2;
            else hh=p1;
    }
    double vis3=dist(p,hh);
    //printf("%.6f\n",vis1-R);
    //printf("%.6f\n",vis2-R);
    //printf("%.6f\n",vis3-R);
    double ans=(min(vis1,min(vis2,vis3))-R);
    printf("%.6f\n",max(0,ans));
}
int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lf%lf%lf",&l.x,&l.y,&l.z);
        scanf("%lf%lf%lf",&r.x,&r.y,&r.z);
        scanf("%lf%lf%lf",&m.x,&m.y,&m.z);
        scanf("%lf%lf%lf",&p.x,&p.y,&p.z);
        scanf("%lf",&R);

        solve();

    }
    return 0;
}