Fast way to compute volume to surface Voronoi projection in R

Question

I have a 3D surface and volume, and would like to find the volume Voronoi cell for each vertex. So within each cell, there will be n voxels belonging to only one vertex and this vertex is the closest vertex to these n voxels among all the other vertices. There is also a distance constrain which ensures we only look for voxels that are not so far from the surface. I made a function below which works but it takes so long (1k sec) to calculate for a volume with 100k voxels and a surface with 10k vertices. Could you please help to find a more efficient way? Thanks a lot.

# voxel_coord is a 100k by 3 matrix for voxel's x, y, z coordinate
# surf_coord is a 10k by 3 matrix for vertex's x, y, z coordinate
# slabsize is the distance constrain, if a voxel is to far from the surface, 
# it will get a index -1, otherwise, the output is vertex index for each voxel. 
voronoi_cell <- apply(voxel_coord, 1, find_voronoi, surf_coord, 5)

find_voronoi <- function(voxel_coord, surf_coord, slabsize){
  dist_temp <- pdist2(voxel_coord, surf_coord)[1,]
  if(min(dist_temp) <= slabsize){
    voronoi_cell <- which(dist_temp == min(dist_temp))
  }else {
    voronoi_cell <- -1
  }
  return(voronoi_cell)
}