【HDU4741】空间解析几何

1.题目链接。题目大意：给出两条空间中不平行的直线，求出这两条直线的距离和对应的点。

2.分析：在空间中我们知道，直线有三种关系：相交，平行，异面。但是题目中已经说了，是不相交的直线，所以只可能有两种关系，平行或者异面。在空间中，直线方程并不是很好用（行吧，或者说我不会用）。虽然说这是一个计算几何模板题，但是还是希望可以明白其中道理。

3.我们从向量的角度来考虑这个问题那么我们就不用考虑很多的位置关系。首先是求距离，求距离很简单了，我们可以利用叉积求出这两个向量的公共法向量，然后再在这两条直线上各取一点，求这两点之间的距离在法向量上的投影即可。交点该怎么求呢？我们可以通过求线与面的交点求得。根据法向量，我们可以确定每条直线各自所在的一个平面，并且这个平面是包含法向量的，这样的平面根据几何公理，有且仅有一个。（其实也就是空间解析几何中常用的：两点和垂直于这两点连线的直线可以共同确定一个平面）其实最短距离对应的两个点，就是法向量和这两个平面的交点。至于如何求平面和直线的交点，其实好求，就是求这个直线在法向量方向的偏移。

#include <iostream>
#include <string.h>
#include <stdio.h>
#include <math.h>
using namespace std;
const double eps = 1e-8;
//三维空间点
struct Point
{
	double x, y, z;
	Point(double x = 0, double y = 0, double z = 0) : x(x), y(y), z(z) {}
	Point(const Point& a)
	{
		x = a.x;
		y = a.y;
		z = a.z;
		return;
	}
	void Print()
	{
		printf("%lf %lf %lf\n", x, y, z);
	}
	Point operator + (Point &t)
	{
		return Point(x + t.x, y + t.y, z + t.z);
	}
};
//空间直线
struct Line
{
	Point a, b;
};
//空间平面
struct Plane
{
	Point a, b, c;
	Plane() {}
	Plane(Point a, Point b, Point c) :a(a), b(b), c(c) {}
	void showPlane()
	{
		a.Print();
		b.Print();
		c.Print();
		return;
	}
};
double dcmp(double t)
{
	if (fabs(t) < eps) return 0;
	return t < 0 ? -1 : 1;
}
//三维叉积
Point cross(Point u, Point v)
{
	Point ret;
	ret.x = u.y * v.z - v.y * u.z;
	ret.y = u.z * v.x - u.x * v.z;
	ret.z = u.x * v.y - u.y * v.x;
	return ret;
}
//三维点积
double multi(Point u, Point v)
{
	return u.x * v.x + u.y * v.y + u.z * v.z;
}
//矢量差
Point sub(Point u, Point v)
{
	Point ret;
	ret.x = u.x - v.x;
	ret.y = u.y - v.y;
	ret.z = u.z - v.z;
	return ret;
}
//两点距离
double dist(Point p1, Point p2)
{
	return 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));
}
//向量的模
double VectorLength(Point p)
{
	return sqrt(p.x*p.x + p.y*p.y + p.z*p.z);
}
//空间直线距离
double LineToLine(Line u, Line v, Point &tmp)
{
	tmp = cross(sub(u.a, u.b), sub(v.a, v.b));
	return fabs(multi(sub(u.a, v.a), tmp)) / VectorLength(tmp);
}
//取平面法向量
Point normalVector(Plane s)
{
	return cross(sub(s.a, s.b), sub(s.b, s.c));
}
//空间平面与直线的交点
Point Intersection(Line l, Plane s)
{
	Point ret = normalVector(s);
	double t = (ret.x*(s.a.x - l.a.x) + ret.y*(s.a.y - l.a.y) + ret.z*(s.a.z - l.a.z)) / (ret.x*(l.b.x - l.a.x) + ret.y*(l.b.y - l.a.y) + ret.z*(l.b.z - l.a.z));
	ret.x = l.a.x + (l.b.x - l.a.x) * t;
	ret.y = l.a.y + (l.b.y - l.a.y) * t;
	ret.z = l.a.z + (l.b.z - l.a.z) * t;
	return ret;
}
/************以上为模板*************/
void work(Line A, Line B)
{
	Point normal;
	double d = LineToLine(A, B, normal);
	printf("%.6lf\n", d);
	Plane alpha = Plane(A.a, A.b, A.a + normal);
	Plane beta = Plane(B.a, B.b, B.a + normal);
	Point u = Intersection(B, alpha);
	Point v = Intersection(A, beta);
	printf("%.6lf %.6lf %.6lf %.6lf %.6lf %.6lf\n", v.x, v.y, v.z, u.x, u.y, u.z);
}
int main()
{
	int T;
	scanf("%d", &T);
	while (T--)
	{
		Line A, B;
		scanf("%lf%lf%lf", &A.a.x, &A.a.y, &A.a.z);
		scanf("%lf%lf%lf", &A.b.x, &A.b.y, &A.b.z);
		scanf("%lf%lf%lf", &B.a.x, &B.a.y, &B.a.z);
		scanf("%lf%lf%lf", &B.b.x, &B.b.y, &B.b.z);
		work(A, B);
	}
	return 0;
}