Circumsphere of a tetrahedron undefined?

Question

I am trying to find 3D alpha shapes from my data-set. In doing so, I am keeping only those tetrahedra that have circumradius below a certain threshold. However, while finding the circumradius of the tetrahedra making up the delaunay triangulation of the data, I am finding that for some cases, it becomes imaginary. I am using the standard formulae, given for example, in http://mathworld.wolfram.com/Circumsphere.html. Are these formualae always valid for any (non-degenerate) tetrahedron ? An example of the vertices of a tetrahedra for which it fails:

1.0e+03 *

-0.882361572000000 1.832846680000000 8.039920898000000

-0.871205933000000 2.190948975000000 7.713502440999999

-0.874571533000000 1.637495972000000 7.953884766000000

-0.945120239000000 1.753712891000000 8.093748535000000

I am using Matlab here, along with its delaunayn function. Essentially in the formula $r = \frac{ \sqrt{D_x^2 + D_y^2 + D_z^2 - 4 a c} }{2 |a| }$, the discriminant under the square root becomes negative. My computed values are:

[$D_x \, D_y \, D_z \, a \, c$] =

1.0e+14 *

0.000152327124454

0.000160507051388

-0.000890797019744

-0.000000057932589

-3.705507881061755

Answer 1

If this was StackOverflow they would say "show me your code." It looks like you've at least done something with the minus sign in calculating Dy. Here's my code where the four vertices of the tetrahedron are stored in the 4-by-3 matrix xyz:

w = ones(4,1);
s = sum(xyz.^2,2);
Dx = det([s xyz(:,2:3) w])
Dy = -det([s xyz(:,[1 3]) w])
Dz = det([s xyz(:,1:2) w])
a = det([xyz w])
c = det([s xyz])
r = 0.5*sqrt(Dx^2+Dy^2+Dz^2-4*a*c)/abs(c)

The resulting radius is a positive real value: about 8.72e-6. Yes, the negative Dy goes away when it gets squared (and could be omitted in your case) – did you maybe put the negative inside of the square root?

I can't answer your your more mathematical question, but I'd hazard to guess that this formula is fine. I'd think that numerical issues would potentially be a bigger issue. There's potential for overflow in your case if the numbers get bigger.