Более быстрое пересечение тетраэдр-тетраэдр

Для одного моего проекта мне требуется надежное обнаружение пересечения между двумя тетраэдрами в трехмерном пространстве. Мне не нужны точки / линии / грани, чтобы просто знать, есть ли пересечение или нет. Касание тоже считается пересечением, но общая грань треугольника не считается пересечением. После целой борьбы, чтобы добиться этого как можно быстрее, мое решение закипало до такого ужаса:

1. пусть есть тетраэдры v0,v1

   каждый тетраэдр имеет 4 треугольника t[4] где каждый треугольник имеет 3 точки p0,p1,p2 и нормальный вектор n,

2. вычислить плоскости всех четырех сторон обоих тетраэдров

   так что любая точка p на плоскости задается уравнением

   dot(p,n) + d = 0
   

   где n это нормально плоскости. Как это известно, это сводится к вычислительной d

   D0[i] = -dot(v0.t[i].p0,v0.t[i].n)
   D1[i] = -dot(v1.t[i].p0,v1.t[i].n)
   i = { 0,1,2,3 }
   

   для каждого треугольника каждого тетраэдра

3. проверить любую комбинацию треугольника против пересечения треугольника между v0, v1

   так что просто цикл между всеми 16 комбинациями и использовать треугольник против пересечения треугольника.

Треугольник v0.t[i] против треугольника v1.t[j] пересечение сводится к этому:

1. вычислить пересечение между плоскостями

   это, очевидно, луч (для непараллельных плоскостей), поэтому простое перекрестное произведение нормалей плоскости даст направление луча

   dir = cross(v0.t[i].n,v1.t[j].n)
   

   теперь нужно просто найти точку пересечения, принадлежащую обеим плоскостям. Используя вычисление детерминанта непосредственно из перекрестного произведения нормалей, вычисление луча выглядит следующим образом:

   // determinants
   det=vector_len2(dir);
   d0=D0[i]/det;
   d1=D1[j]/det;
   // position
   pos = d0*cross(dir,v1.t[j].n) + d1*cross(v0.t[i].n,dir);
   

   для получения дополнительной информации см .:

   • SO / SE: линия пересечения между двумя плоскостями
   • Вольфрам: пересечение плоскость-плоскость

2. составить знаковые интервалы расстояния пересечения треугольных лучей для каждого треугольника

   так просто вычислите пересечение между лучом и каждой линией ребра треугольника, помня минимальное и максимальное расстояние от pos, Нам не нужна фактическая точка, только скалярное расстояние от pos который является параметром, возвращаемым пересечением линии / луча.

3. проверьте, перекрываются ли диапазоны обоих треугольников или нет

   если перекрывается чем v0,v1 пересекаются … если для всех 16 тестов перекрытия не происходит, чем v0,v1 не пересекается

Как вы можете видеть, это очень много вещей для вычисления. Мои знания по линейной алгебре и векторной математике очень ограничены вещами, которые я использую, поэтому есть большая вероятность, что для этого может быть гораздо лучший подход. Я много чего пытался упростить, но безуспешно (например, с помощью bbox, bsphere, с помощью более простого теста, использующего грани лучей и треугольников в одной плоскости и т. Д.), Но результат был медленнее или даже неправильным (не считая крайние случаи правильно).

Вот мой актуальный C ++ реализация:

//---------------------------------------------------------------------------
bool tetrahedrons::intersect_lin_ray(double *a0,double *a1,double *b0,double *db,double &tb)
{
int i0,i1;
double da[3],ta,q;
vector_sub(da,a1,a0); ta=0.0; tb=0.0; i0=0; i1=1;
if (fabs(da[i0])+fabs(db[i0])<=_zero) i0=2;
else if (fabs(da[i1])+fabs(db[i1])<=_zero) i1=2;
q=(da[i0]*db[i1])-(db[i0]*da[i1]);
if (fabs(q)<=_zero) return 0;       // no intersection
// intersection ta,tb parameters
ta=divide(db[i0]*(a0[i1]-b0[i1])+db[i1]*(b0[i0]-a0[i0]),q);
tb=divide(da[i0]*(a0[i1]-b0[i1])+da[i1]*(b0[i0]-a0[i0]),q);
if ((ta<0.0)||(ta>1.0)) return 0;   // inside line check
return 1;
}
//---------------------------------------------------------------------------
bool tetrahedrons::intersect_vol_vol(_vol4 &v0,_vol4 &v1)   // tetrahedron v0 intersect tetrahedron v1 ?
{
int i,j,_ti,_tj;
_fac3 *f0,*f1;
double pos[3],dir[3],p[3],det,D0[4],D1[4],d0,d1,t,ti0,ti1,tj0,tj1;
// planes offset:  dot(p,v0.t[i].n)+D0[i] = 0  ,  dot(p,v1.t[j].n)+D1[j] = 0
for (i=0;i<4;i++)
{
D0[i]=-vector_mul(pnt.pnt.dat+fac.dat[v0.t[i]].p0,fac.dat[v0.t[i]].n);
D1[i]=-vector_mul(pnt.pnt.dat+fac.dat[v1.t[i]].p0,fac.dat[v1.t[i]].n);
}
// plane plane intersection -> ray
for (i=0;i<4;i++)
for (j=0;j<4;j++)
{
f0=fac.dat+v0.t[i];
f1=fac.dat+v1.t[j];
// no common vertex
if ((f0->p0==f1->p0)||(f0->p0==f1->p1)||(f0->p0==f1->p2)) continue;
if ((f0->p1==f1->p0)||(f0->p1==f1->p1)||(f0->p1==f1->p2)) continue;
if ((f0->p2==f1->p0)||(f0->p2==f1->p1)||(f0->p2==f1->p2)) continue;
// direction
vector_mul(dir,f0->n,f1->n);
det=vector_len2(dir);
if (fabs(det)<=_zero) continue; // parallel planes?
d0=D0[i]/det;
d1=D1[j]/det;
// position
vector_mul(p,dir,f1->n); vector_mul(pos,p,d0);
vector_mul(p,f0->n,dir); vector_mul(p,p,d1);
vector_add(pos,pos,p);
// compute intersection edge points
_ti=1; _tj=1;
if (intersect_lin_ray(pnt.pnt.dat+f0->p0,pnt.pnt.dat+f0->p1,pos,dir,t)){ if (_ti) { _ti=0; ti0=t; ti1=t; } if (ti0>t) ti0=t; if (ti1<t) ti1=t; }
if (intersect_lin_ray(pnt.pnt.dat+f0->p1,pnt.pnt.dat+f0->p2,pos,dir,t)){ if (_ti) { _ti=0; ti0=t; ti1=t; } if (ti0>t) ti0=t; if (ti1<t) ti1=t; }
if (intersect_lin_ray(pnt.pnt.dat+f0->p2,pnt.pnt.dat+f0->p0,pos,dir,t)){ if (_ti) { _ti=0; ti0=t; ti1=t; } if (ti0>t) ti0=t; if (ti1<t) ti1=t; }
if (intersect_lin_ray(pnt.pnt.dat+f1->p0,pnt.pnt.dat+f1->p1,pos,dir,t)){ if (_tj) { _tj=0; tj0=t; tj1=t; } if (tj0>t) tj0=t; if (tj1<t) tj1=t; }
if (intersect_lin_ray(pnt.pnt.dat+f1->p1,pnt.pnt.dat+f1->p2,pos,dir,t)){ if (_tj) { _tj=0; tj0=t; tj1=t; } if (tj0>t) tj0=t; if (tj1<t) tj1=t; }
if (intersect_lin_ray(pnt.pnt.dat+f1->p2,pnt.pnt.dat+f1->p0,pos,dir,t)){ if (_tj) { _tj=0; tj0=t; tj1=t; } if (tj0>t) tj0=t; if (tj1<t) tj1=t; }
if ((_ti)||(_tj)) continue;
if ((ti0>=tj0)&&(ti0<=tj1)) return 1;
if ((ti1>=tj0)&&(ti1<=tj1)) return 1;
if ((tj0>=ti0)&&(tj0<=ti1)) return 1;
if ((tj1>=ti0)&&(tj1<=ti1)) return 1;
}
return 0;
};
//---------------------------------------------------------------------------

Это часть гораздо большей программы. _zero это просто порог для нуля на основе минимального размера детали. _fac3 это треугольник и _vol4 это тетраэдр. И точки, и треугольники индексируются из pnt.pnt.dat[] а также fac.dat[] динамические списки. Я знаю, что это странно, но за этим многое происходит (например, пространственное разделение на сегменты и многое другое для ускорения процессов, для которых это используется).

vector_mul(a,b,c) является a=cross(b,c) а также a=dot(b,c) продукт (который зависит от c если это вектор или нет).

Я бы предпочел избегать любых предварительно вычисленных значений для каждого треугольника / тетраэдра, так как даже сейчас классы уже содержат довольно много информации (например, родительский корабль, счет использования и т. Д.). И поскольку я связан с Win32, память ограничена только вокруг 1.2 GB поэтому любой дополнительный материал ограничит максимальный размер используемой сетки.

Так что я ищу что-нибудь из этого:

1. некоторый математический или кодовый трюк для ускорения текущего подхода, если это возможно
2. другой быстрый подход к этому

Я обязан BDS2006 Win32 C ++ и лучше избегать использования сторонних библиотек.

[Edit1] пример данных

Вот тетраэдрическое облако точек в качестве примера данных для тестирования:

double pnt[192]=    // pnt.pnt.dat[pnt.n*3] = { x,y,z, ... }
{
-0.227,0.108,-0.386,
-0.227,0.153,-0.386,
0.227,0.108,-0.386,
0.227,0.153,-0.386,
0.227,0.108,-0.431,
0.227,0.153,-0.431,
-0.227,0.108,-0.431,
-0.227,0.153,-0.431,
-0.227,0.108,0.429,
-0.227,0.153,0.429,
0.227,0.108,0.429,
0.227,0.153,0.429,
0.227,0.108,0.384,
0.227,0.153,0.384,
-0.227,0.108,0.384,
-0.227,0.153,0.384,
-0.023,0.108,0.409,
-0.023,0.153,0.409,
0.023,0.108,0.409,
0.023,0.153,0.409,
0.023,0.108,-0.409,
0.023,0.153,-0.409,
-0.023,0.108,-0.409,
-0.023,0.153,-0.409,
-0.318,0.210,0.500,
-0.318,0.233,0.500,
0.318,0.210,0.500,
0.318,0.233,0.500,
0.318,0.210,-0.500,
0.318,0.233,-0.500,
-0.318,0.210,-0.500,
-0.318,0.233,-0.500,
-0.273,-0.233,0.432,
-0.273,0.222,0.432,
-0.227,-0.233,0.432,
-0.227,0.222,0.432,
-0.227,-0.233,0.386,
-0.227,0.222,0.386,
-0.273,-0.233,0.386,
-0.273,0.222,0.386,
0.227,-0.233,0.432,
0.227,0.222,0.432,
0.273,-0.233,0.432,
0.273,0.222,0.432,
0.273,-0.233,0.386,
0.273,0.222,0.386,
0.227,-0.233,0.386,
0.227,0.222,0.386,
-0.273,-0.233,-0.386,
-0.273,0.222,-0.386,
-0.227,-0.233,-0.386,
-0.227,0.222,-0.386,
-0.227,-0.233,-0.432,
-0.227,0.222,-0.432,
-0.273,-0.233,-0.432,
-0.273,0.222,-0.432,
0.227,-0.233,-0.386,
0.227,0.222,-0.386,
0.273,-0.233,-0.386,
0.273,0.222,-0.386,
0.273,-0.233,-0.432,
0.273,0.222,-0.432,
0.227,-0.233,-0.432,
0.227,0.222,-0.432,
};

struct _fac3 { int p0,p1,p2; double n[3]; };

_fac3 fac[140]= // fac.dat[fac.num] = { p0,p1,p2,n(x,y,z), ... }
{
78, 84, 96, 0.600,-0.800,-0.000,
72, 84, 96, -0.844,-0.003,-0.537,
72, 78, 84, -0.000,1.000,-0.000,
72, 78, 96, -0.000,-0.152,0.988,
6, 84, 96, -0.859,0.336,-0.385,
6, 78, 96, 0.597,-0.801,0.031,
6, 78, 84, 0.746,-0.666,0.000,
6, 72, 96, -0.852,-0.006,-0.523,
6, 72, 84, -0.834,0.151,-0.530,
78, 84,147, 0.020,1.000,-0.000,
72, 84,147, -0.023,-1.000,-0.015,
72, 78,147, -0.000,1.000,0.014,
78, 96,186, 0.546,-0.776,0.316,
6, 96,186, -0.864,0.067,-0.500,
6, 78,186, 0.995,0.014,-0.104,
78, 84,186, 0.980,-0.201,0.000,
6, 84,186, -0.812,0.078,-0.578,
72, 96,186, -0.865,-0.011,-0.501,
6, 72,186, -0.846,0.071,-0.529,
6, 84,147, -0.153,-0.672,-0.724,
6, 72,147, -0.222,-0.975,-0.024,
84,135,147, 0.018,1.000,-0.013,
78,135,147, -0.311,0.924,0.220,
78, 84,135, 0.258,0.966,-0.000,
72,135,147, -0.018,1.000,0.013,
72, 78,135, -0.000,0.995,0.105,
96,132,186, -0.000,-1.000,-0.000,
78,132,186, 0.995,-0.087,-0.056,
78, 96,132, 0.081,-0.256,0.963,
84,132,186, 0.976,-0.209,-0.055,
78, 84,132, 0.995,-0.101,0.000,
84,147,186, -0.190,-0.111,-0.975,
6,147,186, -0.030,-0.134,0.991,
0, 96,186, -0.587,-0.735,-0.339,
0, 72,186, 0.598,0.801,-0.031,
0, 72, 96, -0.992,-0.087,-0.092,
72,147,186, -0.675,-0.737,-0.044,
135,147,189, 0.000,1.000,-0.000,
84,147,189, -0.018,0.980,-0.197,
84,135,189, 0.126,0.992,-0.007,
81, 84,135, -0.183,0.983,-0.023,
78, 81,135, -0.930,-0.000,0.367,
78, 81, 84, 1.000,-0.000,0.000,
105,135,147, -0.000,1.000,0.000,
72,105,147, -0.126,0.992,0.007,
72,105,135, 0.018,0.980,0.197,
72, 81,135, -0.036,0.996,-0.082,
72, 78, 81, -0.000,-0.000,1.000,
96,120,132, -0.000,-1.000,-0.000,
78,120,132, 0.685,-0.246,0.685,
78, 96,120, -0.000,-0.152,0.988,
132,180,186, -0.000,-1.000,0.000,
84,180,186, 0.000,-0.152,-0.988,
84,132,180, 0.995,-0.101,-0.000,
147,150,186, 0.101,0.010,0.995,
84,150,186, -0.100,-0.131,-0.986,
84,147,150, -0.190,-0.019,-0.982,
96,114,186, 0.000,-1.000,0.000,
0,114,186, -0.584,-0.729,-0.357,
0, 96,114, -0.991,0.134,0.000,
0,147,186, -0.144,-0.058,-0.988,
0, 72,147, -0.926,-0.374,-0.052,
72, 96,114, -0.995,-0.101,0.000,
0, 72,114, -0.993,-0.077,-0.093,
75,147,189, -0.001,1.000,-0.012,
75,135,189, 0.018,1.000,-0.001,
75,135,147, -0.016,-1.000,0.012,
147,159,189, -0.000,1.000,-0.000,
84,159,189, -0.000,0.985,-0.174,
84,147,159, -0.025,-0.999,-0.025,
81,135,189, -0.274,0.962,0.015,
81, 84,189, 0.114,0.993,-0.023,
75,105,147, -0.115,-0.993,0.006,
75,105,135, 0.017,-0.983,0.181,
72, 75,147, -0.999,-0.000,-0.051,
72, 75,105, 0.599,-0.000,0.801,
81,105,135, -0.009,0.996,-0.093,
72, 81,105, -0.036,0.991,0.127,
120,126,132, -0.000,-1.000,-0.000,
78,126,132, 0.995,-0.101,-0.000,
78,120,126, -0.000,-0.152,0.988,
0,150,186, 0.101,-0.000,0.995,
0,147,150, -0.000,-0.000,1.000,
144,150,186, 0.000,-1.000,0.000,
84,144,186, -0.091,-0.133,-0.987,
84,144,150, -0.000,0.249,0.968,
147,150,159, -0.705,-0.071,-0.705,
84,150,159, -0.125,-0.100,-0.987,
114,150,186, 0.000,-1.000,0.000,
0,114,150, -0.998,-0.000,-0.059,
72,114,147, -0.995,-0.088,-0.052,
0,114,147, -0.906,-0.365,-0.215,
93,147,189, -0.009,-0.996,-0.093,
75, 93,189, 0.020,1.000,0.000,
75, 93,147, -0.237,-0.971,-0.000,
75, 81,189, -0.000,1.000,-0.012,
75, 81,135, -0.000,-0.995,0.096,
93,159,189, -0.000,-0.987,-0.160,
93,147,159, -0.069,-0.995,-0.069,
84, 93,189, 0.036,0.991,-0.127,
84, 93,159, -0.036,-0.993,-0.113,
84, 87,189, -0.599,-0.000,-0.801,
81, 87,189, -0.120,0.993,-0.000,
81, 84, 87, 1.000,0.000,0.000,
75, 81,105, -0.000,-0.987,0.160,
72, 93,147, -0.183,-0.983,-0.023,
72, 75, 93, -1.000,0.000,-0.000,
72, 75, 81, 0.000,-0.000,1.000,
114,147,150, -0.993,-0.100,-0.059,
144,162,186, 0.000,-1.000,0.000,
84,162,186, -0.000,-0.152,-0.988,
84,144,162, -0.600,0.800,0.000,
144,150,159, 0.000,0.101,0.995,
84,144,159, -0.125,-0.087,-0.988,
144,147,159, -0.707,0.000,-0.707,
144,147,150, -0.000,0.000,1.000,
93,114,147, 0.732,-0.587,-0.346,
72, 93,114, -0.995,-0.100,-0.002,
81, 93,189, 0.022,1.000,-0.014,
75, 81, 93, -0.000,1.000,0.000,
93,144,159, 0.582,-0.140,-0.801,
93,144,147, -0.930,0.000,0.367,
87, 93,189, -0.000,0.987,0.160,
84, 87, 93, -0.000,0.000,-1.000,
84, 93,144, -0.009,-0.238,-0.971,
81, 87, 93, -0.000,1.000,0.000,
114,144,150, -0.000,-1.000,-0.000,
114,144,147, -1.000,0.000,-0.000,
93,144,162, -0.995,-0.096,0.000,
84, 93,162, -0.005,-0.145,-0.989,
93,114,144, -0.995,-0.096,0.000,
72,114,144, -0.995,-0.101,-0.000,
72, 93,144, -0.995,-0.097,-0.002,
90,144,162, -0.995,-0.101,0.000,
90, 93,162, 0.834,0.000,-0.552,
90, 93,144, -0.930,0.000,0.367,
84, 90,162, 0.000,-0.152,-0.988,
84, 90, 93, 0.000,0.000,-1.000,
72, 90,144, -0.995,-0.101,-0.000,
72, 90, 93, -1.000,0.000,-0.000,
};

struct _vol4 { int p0,p1,p2,p3,t[4]; double s[4]; };

_vol4 vol[62]=  // vol.dat[vol.num] = { p0,p1,p2,p3,t[0],t[1],t[2],t[3],s[0],s[1],s[2],s[3], ... }
{
72, 78, 96, 84,   0,  1,  2,  3, 1,1,1,1,
78, 84, 96,  6,   4,  5,  6,  0, 1,1,1,-1,
72, 84, 96,  6,   4,  7,  8,  1, -1,1,1,-1,
72, 78, 84,147,   9, 10, 11,  2, 1,1,1,-1,
6, 78, 96,186,  12, 13, 14,  5, 1,1,1,-1,
6, 78, 84,186,  15, 16, 14,  6, 1,1,-1,-1,
6, 72, 96,186,  17, 13, 18,  7, 1,-1,1,-1,
6, 72, 84,147,  10, 19, 20,  8, -1,1,1,-1,
78, 84,147,135,  21, 22, 23,  9, 1,1,1,-1,
72, 78,147,135,  22, 24, 25, 11, -1,1,1,-1,
78, 96,186,132,  26, 27, 28, 12, 1,1,1,-1,
78, 84,186,132,  29, 27, 30, 15, 1,-1,1,-1,
6, 84,186,147,  31, 32, 19, 16, 1,1,-1,-1,
72, 96,186,  0,  33, 34, 35, 17, 1,1,1,-1,
6, 72,186,147,  36, 32, 20, 18, 1,-1,-1,-1,
84,135,147,189,  37, 38, 39, 21, 1,1,1,-1,
78, 84,135, 81,  40, 41, 42, 23, 1,1,1,-1,
72,135,147,105,  43, 44, 45, 24, 1,1,1,-1,
72, 78,135, 81,  41, 46, 47, 25, -1,1,1,-1,
78, 96,132,120,  48, 49, 50, 28, 1,1,1,-1,
84,132,186,180,  51, 52, 53, 29, 1,1,1,-1,
84,147,186,150,  54, 55, 56, 31, 1,1,1,-1,
0, 96,186,114,  57, 58, 59, 33, 1,1,1,-1,
0, 72,186,147,  36, 60, 61, 34, -1,1,1,-1,
0, 72, 96,114,  62, 59, 63, 35, 1,-1,1,-1,
135,147,189, 75,  64, 65, 66, 37, 1,1,1,-1,
84,147,189,159,  67, 68, 69, 38, 1,1,1,-1,
84,135,189, 81,  70, 71, 40, 39, 1,1,-1,-1,
105,135,147, 75,  66, 72, 73, 43, -1,1,1,-1,
72,105,147, 75,  72, 74, 75, 44, -1,1,1,-1,
72,105,135, 81,  76, 46, 77, 45, 1,-1,1,-1,
78,120,132,126,  78, 79, 80, 49, 1,1,1,-1,
147,150,186,  0,  81, 60, 82, 54, 1,-1,1,-1,
84,150,186,144,  83, 84, 85, 55, 1,1,1,-1,
84,147,150,159,  86, 87, 69, 56, 1,1,-1,-1,
0,114,186,150,  88, 81, 89, 58, 1,-1,1,-1,
0, 72,147,114,  90, 91, 63, 61, 1,1,-1,-1,
75,147,189, 93,  92, 93, 94, 64, 1,1,1,-1,
75,135,189, 81,  70, 95, 96, 65, -1,1,1,-1,
147,159,189, 93,  97, 92, 98, 67, 1,-1,1,-1,
84,159,189, 93,  97, 99,100, 68, -1,1,1,-1,
81, 84,189, 87, 101,102,103, 71, 1,1,1,-1,
75,105,135, 81,  76, 96,104, 73, -1,-1,1,-1,
72, 75,147, 93,  94,105,106, 74, -1,1,1,-1,
72, 75,105, 81, 104, 77,107, 75, -1,-1,1,-1,
0,147,150,114, 108, 89, 91, 82, 1,-1,-1,-1,
84,144,186,162, 109,110,111, 84, 1,1,1,-1,
84,144,150,159, 112, 87,113, 85, 1,-1,1,-1,
147,150,159,144, 112,114,115, 86, -1,1,1,-1,
72,114,147, 93, 116,105,117, 90, 1,-1,1,-1,
75, 93,189, 81, 118, 95,119, 93, 1,-1,1,-1,
93,147,159,144, 114,120,121, 98, -1,1,1,-1,
84, 93,189, 87, 122,101,123, 99, 1,-1,1,-1,
84, 93,159,144, 120,113,124,100, -1,-1,1,-1,
81, 87,189, 93, 122,118,125,102, -1,-1,1,-1,
114,147,150,144, 115,126,127,108, -1,1,1,-1,
84,144,162, 93, 128,129,124,111, 1,1,-1,-1,
93,114,147,144, 127,121,130,116, -1,-1,1,-1,
72, 93,114,144, 130,131,132,117, -1,1,1,-1,
93,144,162, 90, 133,134,135,128, 1,1,1,-1,
84, 93,162, 90, 134,136,137,129, -1,1,1,-1,
72, 93,144, 90, 135,138,139,132, -1,1,1,-1,
};

p? точечные индексы 0,3,6,9... от pnt n нормально s является признаком нормали (если треугольник является общим, так что нормали указывают одинаково) и t[4] являются индексами треугольников 0,1,2,3,... от fac,
А вот образец теста:

bool tetrahedrons::vols_intersect()                     // test if vol[] intersects each other
{
int i,j;
for (i=0;i<vol.num;i++)
for (j=i+1;j<vol.num;j++,dbg_cnt++)
if (intersect_vol_vol(vol.dat[i],vol.dat[j]))
{
linc=0x800000FF;
if (intersect_vol_vol(vol.dat[j],vol.dat[i])) linc=0x8000FFFF;
lin_add_vol(vol.dat[i]);
lin_add_vol(vol.dat[j]);
return 1;
}
return 0;
}

где dbg_cnt счетчик перекрестных испытаний. Для этой сетки я получил следующие результаты:

tests | time
------+-------------
18910 | 190-215 [ms]

Я позвонил vols_intersect проверить 10 раз, чтобы сделать измерения достаточно долго. Конечно, ни один из размещенных тетраэдров в этом наборе данных не пересекается (что приводит к максимальному времени). В реальном процессе (слишком большом, чтобы публиковать), который приводит к этой сетке, счетчик выглядит так:

intersecting        5
non intersecting 1766
all tests        1771