Determining if a sphere is enclosed completely by other spheres placed around it

Question

Problem: Given a list of spheres, find all empty spaces that are completely enclosed by spheres.

Detail: This is a problem I am working on in which I try to determine the cavities located in a protein. I am given a list of atoms that make up the protein ((x,y,z) coordinates and radius). I then run my algorithm to find all empty spaces that lie within the bounds of the protein by checking if a probe (of given radius) can be placed at a location without colliding with other spheres. There are two types of empty spaces, void spaces and cavities. Void spaces are spaces that can lead to or on the outside of the protein. Cavities are empty spaces that are completely enclosed by protein atoms. Here is an image of the sample "protein" we are working with.

It can be viewed in three dimensions here.

There is a cavity located near the center of the protein, the tunnel you see going through the protein would be considered a void space because it is not fully enclosed by atoms.

Example: Given a list of 26 atoms, these atoms are evenly spaced from (0,0,0) to (1,1,1) in a 3-dimensional grid. Each atom has a radius of 0.25 and is placed on either 0, 0.5, or 1 on any axis. There is no atom at the point (0.5, 0.5, 0.5). If we were to draw a 3D figure of these atoms, it would be a cube like shape with the center missing. A cavity would be designated at (0.5,0.5,0.5) with a radius of 0.25. It can be assumed that this cavity is surrounded by proteins on all sides.

Example image:

Note that the above is only a 2D representation of the cube and protein. It is actually 3D.

How would one go about determining void spaces vs. cavities for a much larger and irregularly shaped group of atoms?

I was thinking about implementing a recursive algorithm that checks every direction to see if it can reach the maximum and minimum bounds of the graph but I am not sure if this is the correct way to go about doing it.

Extra: Is there a different algorithm that would say the cavity in the example is actually a void space because there are very small "paths" to reach the outside of the protein? A cavity would have to be completely enclosed by atoms to exist. Any void spaces that have a path (in any direction, not necessarily straight) to the outside of the protein would not be considered cavities.

Answer 1

Cool question. Here's an algorithm that should do the trick:

Notation:

• Let's call our moveable sphere S.
• Write diam(X) for the diameter of a sphere X
• Write dist(X,Y) for the distance from X to Y; this is the same as the distance from the center of X to the center of Y minus the sum of the radii.

The algorithm:

1. For any two unmoveable spheres A and B, check whether S can pass directly between the centers of A and B (i.e. is diam(S) <= dist(A,B)?).
2. If so, for each other sphere C, check whether S could simultaneously touch all three spheres A, B, and C, if there were no other spheres present. If S could simultaneously touch all 3, draw a triangle between the centers of A, B, and C.
   • This can be checked in several ways. One fairly easy way: the possible positions of the center of S while touching both A and B form a circle. You want to know whether this circle has a point on it which is less than diam(S) + diam(C) away from the center of C. This is easy geometry.
3. The problem now reduces to the question: do the triangles separate the initial position of the center of S from infinity? You can answer this one connected component at a time. In fact, you can even answer this one "edge-connected" component at a time, where a component is edge-connected if any two non-vertex points can be linked by a path that doesn't pass through any vertices. You can calculate these components through a simple graph search.
4. For a given edge-connected component, you need to decide whether the component separates the center of S from infinity. There are a few ways you might do this:
   • Calculate the 2-homology of the component, choose effective generators, and for each, ask whether your point and infinity are on the same side of the cycle or not, which can be checked using the orientation class.
   • Or, just start painting the component:
     • Start with a triangle that you can reach from S, and paint every face that can be reached from there. This is slightly subtle, but the algorithm is just "start anywhere, queue up the edges, cross each edge onto the face forming the smallest angle with that edge, and stop when there are no edges left." Keep in mind that the opposite side of the same triangle might be the face forming the smallest angle.
     • Do the same from infinity. Did you cross any painted triangles? If yes, your sphere can escape. If no, it can't.

Why it works

Step 3 is true because if you don't hit any sphere C while "rolling around" the edge between A and B, then you can reach any side of that edge. Put another way, any position that stops you from going to infinity must involve S touching at least 3 spheres.

Note that there are some subtleties that arise from "exceptional" situations, like when S touches 4 spheres at once. You can avoid these subtleties by re-triangulating your shape before performing steps 3 and 4.