Tesselation of 3D space into truncated octahedron voxels

[rlavendomme]

Current status

Currently MoloVol analyzes the 3D space by dividing it into cubic voxels.

Issue

1. Due to the shape and arrangement of cubes (of side length r), each voxel has 26 neighbors:

• 6 for faces (separated by 1 r)
• 12 for edges (separated by Sqrt(2) r = 1.414 r)
• 8 for vertices (separated by Sqrt(3) r = 1.732 r)
  This caused some issue regarding evaluation of voxel types within a range search. We fixed this issue by slightly increasing the range of the search algorithm to compensate for neighbors that are further than 1 r. This issue and the solution are detailed in the scientific publication for MoloVol.

2. The flood-fill algorithm has to consider all 26 neighbors for each voxel which is time-consuming.

3. The anisotropy of the cubes (far from spherical shape) lead to a structure orientation dependency in the calculated values. Calculation at higher resolution (smaller grid step) reduce the impact of the starting orientation.

Proposed change

Other than cubes, there are other polyhedra that can fill a 3D space. One interesting polyhedron is the truncated octahedron that can be conveniently described by a 3D grid of voxels and symmetric along X, Y and Z axes, similarly to cubic voxels.
Each center of a truncated octahedron voxel would be sitting on two uniform cartesian grids: (x,y,z) and (0.5+x,0.5+y,0.5+z).

From a crystallographic point of view, cubic voxels correspond to a primitive cubic lattice in which each spherical atom is inflated until filling the space as cubes, whereas truncated octahedron voxels correspond to a body-centered cubic lattice (see cubic crystal systems).

Some adjustments would be necessary for crystal unit cell analysis in MoloVol as the truncated octahedron voxels don't define a planar border around the unit cell. Yet it should be simple to fix if we keep only half of the volume of voxels going through faces X=0, Y=0, Z=0, and the necessary ratio for faces X=a, Y=b, Z=c (similarly to cubes that go above each face a, b, c in the current version of MoloVol).

For surface area calculations, two possibilities could be explored:

• We keep the current algorithm focusing only on voxels in the cubic grid (x,y,z) which should give the same results as in the current version.
• We adapt the calculations with an algorithm devoted on truncated octahedron voxels. It is a problem that was studied by Erik Miller who published a paper based on its freely accessible Master thesis. He showed that truncated octahedron lead to better surface area calculations compared to cubes for different types of curves. Note that rhombic dodecahedra perform even better but are not suitable for MoloVol since they are not symmetric along X, Y and Z axes.

Improvement

1. Each truncated octahedron voxel now has 14 neighbors:

• 6 for square faces (separated by 1 r)
• 8 for hexagonal faces (separated by Sqrt(3)/2 r = 0.866 r)
  These distances are much more uniform compared to the previous case (from 1 to 1.732 r) which may allow to completely remove the search range fix currently used.

2. There are fewer neighbors to consider for the flood-fill algorithm.

3. For a same grid step r, the new voxel system would contain twice the number of voxels thus increasing the precision of results but likely increasing calculation time (even if the flood-fill algorithm would be faster). Consequently, larger grid step could be used to maintain the same precision for results with, perhaps, an improved calculation time.

4. The voxel space would be more "isotropic" (each voxel is closer to a sphere), which means that the calculation results would be less dependent on the initial orientation of the chemical structure.