Circumcenter of Tetrahedron (in 4D)

760 Views Asked by At

I am trying to calculate the circumcenter of a tetrahedron in 4 dimensional space. Basically what I am looking for is the center of the smallest sphere which passes through all 4 vertices of the tetrahedron. I have searched online but can't seem to find any specific formula for this. My overall aim is to find the circumcenter and check if any other point within a given data set lies within the sphere constructed around the vertices of the tetrahedron. Similar to how the Delaunay triangulation works. Note that the tetrahedron can be a regular tetrahedron and also an irregular tetrahedron.

Currently I am using a custom optimization function which uses a GA to locate a point which is equidistant from all 4 vertices. However, this doesn't always find the smallest enclosing sphere. I was hoping for some concrete mathematical formula which can make this calculation more accurate.

1

There are 1 best solutions below

6
On BEST ANSWER

I don't know an explicit formula (and if one exists, it's probably not easy to digest), but it's easy to get the center as a result of solving a small linear system -- without optimization algorithms. Given points P1, P2, P3, P4, apply translation so that P1 becomes the origin. The center of any sphere containing P1,...,P4 will have the property of being equidistant: from P1 and P2; from P1 and P3; from P1 and P4. Each of these is a linear equation. There are 3 equations for 4 unknowns, so the system is underdetermined. The solution you want is the one of smallest norm (the least squares solution). This is not the solution that backslash operator gives, but one can use pinv to get it, which for small matrices isn't costly.

p1 = rand(1,4); p2 = rand(1,4); p3 = rand(1,4); p4 = rand(1,4);  % test input
v = [p2-p1; p3-p1; p4-p1];       % matrix of linear system
b = 0.5*[v(1,:)*v(1,:)' ; v(2,:)*v(2,:)' ; v(3,:)*v(3,:)'];  % RHS of the system
x = p1 + (pinv(v)*b)'    % least squares solution, translated back

You can check that norm(x-p1)...norm(x-p4) are equal.