Здравствуйте, iperesypkina!
Держите программу. Максимально прокомментировал.
Но без знания теории понять не получится...
Надеюсь, с аналитической геометрией Вы на "ты".
[code h=200]#include <stdlib.h>
#include <stdio.h>
#include <time.h>
#include <malloc.h>
#include <math.h>
//структура, описывающая точку в трехмерном пространстве
typedef struct _point
{
double x;
double y;
double z;
}POINT;
//структура, описывающая линию, заданную двумя точками
typedef struct _line
{
POINT p1;
POINT p2;
} LINE;
double e = 0.00001;
//нахождение смешанного произведения векторов
double mixesPrud(POINT *v1, POINT *v2, POINT *v3)
{
return (v1->x*(v2->y*v3->z - v2->z*v3->y) - v1->y*(v2->x*v3->z - v2->z*v3->x) + v1->z*(v2->x*v3->y - v2->y*v3->x));
}
//проверка коллинеарности векторов
//вещественное число сравниваем не c 0, а на меньше числа e=0.00001
bool isCollinear(POINT *v1, POINT *v2)
{
return (((v1->x/v2->x - v1->y/v2->y) < e) && ((v1->x/v2->x - v1->z/v2->z) < e));
}
//проверка прямых на совпадение
//подставляем значение одной координаты одной прямой в уравнение втророй и сравниваем с остальными
//если совпадает, то прямые совпадают
//учитываем, что прямые могут быть перпендикулярны осям координат
//вещественное число сравниваем не c 0, а на меньше числа e=0.00001
bool isEqualLines(LINE *pline1, LINE *pline2)
{
double x, y, z;
if (pline2->p2.x - pline2->p1.x > e) //прямая не перпендикулярна оси Х
{
y = pline2->p1.y + (pline1->p1.x - pline2->p1.x) * (pline2->p2.y - pline2->p1.y) / (pline2->p2.x - pline2->p1.x);
z = pline2->p1.z + (pline1->p1.x - pline2->p1.x) * (pline2->p2.z - pline2->p1.z) / (pline2->p2.x - pline2->p1.x);
return (((y - pline1->p1.y) < e) && ((z - pline1->p1.z) < e));
}
else if (pline2->p2.y - pline2->p1.y > e) //прямая не перпендикулярна оси Y
{
x = pline2->p1.x + (pline1->p1.y - pline2->p1.y) * (pline2->p2.x - pline2->p1.x) / (pline2->p2.y - pline2->p1.y);
z = pline2->p1.z + (pline1->p1.y - pline2->p1.y) * (pline2->p2.z - pline2->p1.z) / (pline2->p2.y - pline2->p1.y);
return (((x - pline1->p1.x) < e) && ((z - pline1->p1.z) < e));
}
else if (pline2->p2.z - pline2->p1.z > e) //прямая не перпендикулярна оси Z
{
x = pline2->p1.x + (pline1->p1.z - pline2->p1.z) * (pline2->p2.x - pline2->p1.x) / (pline2->p2.z - pline2->p1.z);
y = pline2->p1.y + (pline1->p1.z - pline2->p1.z) * (pline2->p2.y - pline2->p1.y) / (pline2->p2.z - pline2->p1.z);
return (((x - pline1->p1.x) < e) && ((y - pline1->p1.y) < e));
}
else
return false; //прямая фактически задана двумя одинаковыми точками
// при желании можно отфильтровать при задании прямых
}
//нахождение расстояния между прямыми
double distance(LINE *pline1, LINE *pline2)
{
POINT v, s, n1, n2;
double mp;
//сначала найдем направляющие вектора n1 и n2 прямых pline1 и pline2
n1.x = pline1->p2.x - pline1->p1.x;
n1.y = pline1->p2.y - pline1->p1.y;
n1.z = pline1->p2.z - pline1->p1.z;
n2.x = pline2->p2.x - pline2->p1.x;
n2.y = pline2->p2.y - pline2->p1.y;
n2.z = pline2->p2.z - pline2->p1.z;
//Найдем вектор, соединяющий точку pline1->p2 с точкой pline2->p1
v.x = pline1->p2.x - pline2->p1.x;
v.y = pline1->p2.y - pline2->p1.y;
v.z = pline1->p2.z - pline2->p1.z;
//прямые скрещиваются, если векторы n1, n2 и v не компланарны
//найдем модуль смешанного произведения и сравним на 0
mp = fabs(mixesPrud(&n1, &n2, &v));
if (mp < e)
{ //компланарны, т.е. в одной плоскости
if (isCollinear(&n1, &n2)) //коллинеарны? (параллельны?)
{ //да - значит, прямые не пересекаются
if (isEqualLines(pline1, pline2)) //прямые совпадают?
return 0.; //да, расстояние = 0
else
{ //прямые параллельны
//найдем векторное произведение вектора v,
// и направляющего вектора, например, прямой pline2
s.x = v.y * n2.z - v.z * n2.y;
s.y = - v.x * n2.z + v.z * n2.x;
s.z = v.x * n2.y - v.y * n2.x;
//получим расстояние, как отношение длины вектора s к длине n2
return (sqrt(s.x*s.x + s.y*s.y+ s.z*s.z) / sqrt(n2.x*n2.x + n2.y*n2.y+ n2.z*n2.z));
}
}
else
return 0.; //прямые пересекаются, расстояние = 0
}
else
{ //некомпланарны - прямые скрещиваются
//найдем векторное произведение направляющих векторов
s.x = n1.y * n2.z - n1.z * n2.y;
s.y = - n1.x * n2.z + n1.z * n2.x;
s.z = n1.x * n2.y - n1.y * n2.x;
//получим расстояние, как отношение модуля смешанного произведения
// и длины векторного произведения направляющих векторов
return (mp / sqrt(s.x*s.x + s.y*s.y + s.z*s.z));
}
}
int main()
{
const int N = 10; //количество прямых
int i, j;
LINE* pLines = (LINE*)malloc(N*sizeof(LINE)); //память под прямые
srand(time(0)); //инициируем генератор псевдослучайных чисел
for(i=0; i<N; i++) //генерируем случайные прямые, величины берем в интервале [-100.0:99.9]
{
pLines[i].p1.x = ((double)(rand()%2000) / 10) - 100;
pLines[i].p1.y = ((double)(rand()%2000) / 10) - 100;
pLines[i].p1.z = ((double)(rand()%2000) / 10) - 100;
pLines[i].p2.x = ((double)(rand()%2000) / 10) - 100;
pLines[i].p2.y = ((double)(rand()%2000) / 10) - 100;
pLines[i].p2.z = ((double)(rand()%2000) / 10) - 100;
}
double d, max = 0; //текущее и максимальное расстояние
int i1, i2; //индексы прямых, при которых максимальное расстояние
//начинаем искать
for(i=0; i<N-1; i++) //от первой прямой до предпоследней
{
for(j=i+1; j<N; j++) //от следующей до конца
{
d = distance(&pLines[i], &pLines[j]); //получаем расстояние
if (d > max)
{
max = d; //находим максимальное
i1 = i;
i2 = j;
}
}
}
//выводим максимальное расстояние и точки прямых, максимально удаленных друг от друга
printf("d = %7.4f\n[%4.1f,%4.1f,%4.1f]:[%4.1f,%4.1f,%4.1f]\n[%4.1f,%4.1f,%4.1f]:[%4.1f,%4.1f,%4.1f]\n", max,
pLines[i1].p1.x, pLines[i1].p1.y, pLines[i1].p1.z,
pLines[i1].p2.x, pLines[i1].p2.y, pLines[i1].p2.z,
pLines[i2].p1.x, pLines[i2].p1.y, pLines[i2].p1.z,
pLines[i2].p2.x, pLines[i2].p2.y, pLines[i2].p2.z);
free(pLines); //освобождаем память
return 0;
}
[/code]