三维几何之求空间俩线段的公垂线,以及分数类

分析:

首先判断线段俩直线是否平行(或重合),如果是的话直接求。考虑4个端点到另外一条线段的距离,取最小值即可。

如果不平行或重合,说明俩条直线是异面直线,这时最短距离既可能是某端点到另外一条线段的距离,也可能是异面直线的最短距离。

如何求异面直线的最短距离?假设俩条直线分别为l1=(p1,v1)和l2=(p2,v2),那么最短距离会在某个q1=p1+sv1和q2=p2+tv2上取到,其中q1和q2分别在l1和l2上,且q1q2是这俩条异面直线的公垂线。

向量q1q1=q2-q1=p2-p1+tv2-sv1垂直于V1,因此Dot(p2-p1+tv2-sv1,v1)=0,根据分配率,有Dot(p2-p1,v1)+t*Dot(v2,v1)-s*Dot(v1,v1)=0.注意这里的三个点积都是可以直接算出来的,因此实际上得到的是一个关于t和s的一次方程。根据q1q2垂直于v2还可以得到一个一次方程,联立求解即可。下面的代码直接使用了经过复杂化简以后的结果。

注意本题比较特殊,要求以分数方式输出。所以可以考虑定义一个Rat类(代表Rational)来保存和计算有理数,并且重载加法、减法、和乘法,然后把本题用到的俩个关键函数:点到线段距离Distace2Tosegment和异面直线的最小距离LineDistace3D改写成返回有理数的版本。

代码:

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cstdlib>
#include<cctype>
#include<cmath>
#include<string>
#include<map>
#include<set>
#include<vector>
#include<queue>
#include<stack>
#define LL long long
using namespace std;
const double eps=1e-6;
int dcmp(double x)
{
    return fabs(x)<eps?0:(x>0?1:-1);
}
struct Point3
{
    LL x,y,z;
    Point3(LL x,LL y,LL z):x(x),y(y),z(z){}
    Point3() {}
    void in()
    {
        cin>>x>>y>>z;
    }
};
typedef Point3 Vector3;
Vector3 operator +(Vector3 A,Vector3 B)
{
    return Vector3(A.x+B.x,A.y+B.y,A.z+B.z);
}
Vector3 operator -(Vector3 A,Vector3 B)
{
    return Vector3(A.x-B.x,A.y-B.y,A.z-B.z);
}
Vector3 operator * (Vector3 A,int p)
{
    return Vector3(A.x*p,A.y*p,A.z*p);
}
Vector3 operator / (Vector3 A,int p)
{
    return Vector3(A.x/p,A.y/p,A.z/p);
}
bool operator == (Vector3 A,Vector3 B)
{
    return A.x==B.x&&A.y==B.y&&A.z==B.z;
}
int Dot(Vector3 A,Vector3 B)
{
    return A.x*B.x+A.y*B.y+A.z*B.z;
}
double Length(Vector3 A)
{
    return sqrt(Dot(A,A));
}
int Length2(Vector3 A)
{
    return Dot(A,A);
}
Vector3 Cross(Vector3 A,Vector3 B)
{
    return Vector3(A.y*B.z-A.z*B.y,A.z*B.x-A.x*B.z,A.x*B.y-A.y*B.x);
}
LL gcd(LL a,LL b){return b?gcd(b,a%b):a;}
LL lcm(LL a,LL b){return a/gcd(a,b)*b;}
struct Rat
{
    LL a,b;
    Rat(LL a=0):a(a),b(1) {}
    Rat(LL x,LL y):a(x),b(y){
        if(b<0) a=-a,b=-b;
        LL d=gcd(a,b);
        if(d<0) d=-d;
        a/=d;
        b/=d;
    }
};
Rat operator + (const Rat &A,const Rat &B)
{
    LL x=lcm(A.b,B.b);
    return Rat(A.a*(x/A.b)+B.a*(x/B.b),x);
}
Rat operator - (const Rat &A,const Rat &B)
{
    return  A+Rat(-B.a,B.b);
}
Rat  operator *(const Rat &A,const Rat &B)
{
    return Rat(A.a*B.a,A.b*B.b);
}

void updatemin(Rat &A,const Rat &B)
{
    if(A.a*B.b>B.a*A.b) A.a=B.a,A.b=B.b;
}
Rat Rat_Distace2ToSegment(const Point3 &P,const Point3 &A,const Point3 &B)
{
    if(A==B) return Length2(P-A);
    Vector3 v1=B-A,v2=P-A,v3=P-B;
    if(Dot(v1,v2)<0) return Length2(v2);
    else if(Dot(v1,v3)>0) return Length2(v3);
    else return Rat(Length2(Cross(v1,v2)),Length2(v1));
}
bool Rat_LineDistace3D(const Point3 &p1,const Point3 &u,const Point3 &p2,const Vector3 &v, Rat &s)
{
    LL b =(LL)Dot(u,u)*Dot(v,v)-(LL)Dot(u,v)*Dot(u,v);
    if(b==0) return false;
    LL a=(LL)Dot(u,v)*Dot(v,p1-p2)-(LL)Dot(v,v)*Dot(u,p1-p2);
    s=Rat(a,b);
    return true;
}
void Rat_GetPointOnLine(const Point3 &A,const Point3 &B,const Rat &t,Rat &x,Rat &y,Rat &z)
{
    x=Rat(A.x)+Rat(B.x-A.x)*t;
    y=Rat(A.y)+Rat(B.y-A.y)*t;
    z=Rat(A.z)+Rat(B.z-A.z)*t;
}
Rat Rat_Distance2(const Rat &x1,const Rat &y1,const Rat &z1,const Rat &x2,const Rat &y2,const Rat &z2)
{
    return (x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)+(z1-z2)*(z1-z2);
}
int main()
{
    int T;
    cin>>T;
    while(T--)
    {
        Point3 A,B,C,D;
        A.in();
        B.in();
        C.in();
        D.in();
        Rat s,t;
        bool ok=false;
        Rat ans=Rat(1000000000);
        if(Rat_LineDistace3D(A,B-A,C,D-C,s))
            if(s.a>0&&s.a<s.b&&Rat_LineDistace3D(C,D-C,A,B-A,t))
                if(t.a>0&&t.a<t.b)
        {
            ok=true;
            Rat x1,y1,z1,x2,y2,z2;
            Rat_GetPointOnLine(A,B,s,x1,y1,z1);
            Rat_GetPointOnLine(C,D,t,x2,y2,z2);
            ans=Rat_Distance2(x1,y1,z1,x2,y2,z2);
        }
        if(!ok)
        {
            updatemin(ans,Rat_Distace2ToSegment(A,C,D));
            updatemin(ans,Rat_Distace2ToSegment(B,C,D));
            updatemin(ans,Rat_Distace2ToSegment(C,A,B));
            updatemin(ans,Rat_Distace2ToSegment(D,A,B));
        }
        printf("%lld %lld\n",ans.a,ans.b);
    }
    return 0;
}

  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值