This paper introduces a new algorithm for surface reconstruction from unorganized sample points in \(\mathbb{R}^3\). The algorithm is the first to provide provable guarantees on the topological correctness and convergence to the original surface as the sampling density increases. The key idea is to use the three-dimensional Voronoi diagram and Delaunay triangulation to generate a set of triangles called the "crust" of the sample points. The output mesh interpolates the input points, ensuring that the reconstruction is accurate and detailed in feature-rich areas while being sparse in featureless regions. The algorithm's effectiveness is demonstrated through practical examples and theoretical analysis. It is shown that the crust algorithm can handle noisy data and undersampled regions, although it struggles with sharp edges and boundaries. The authors propose heuristics to improve the reconstruction in these challenging cases and discuss future research directions, including handling noise, sharp edges, and boundaries, and exploring applications in mesh compression. The implementation of the algorithm is efficient, capable of handling large datasets in minutes, and robust to numerical issues due to the use of exact arithmetic in Delaunay triangulation.This paper introduces a new algorithm for surface reconstruction from unorganized sample points in \(\mathbb{R}^3\). The algorithm is the first to provide provable guarantees on the topological correctness and convergence to the original surface as the sampling density increases. The key idea is to use the three-dimensional Voronoi diagram and Delaunay triangulation to generate a set of triangles called the "crust" of the sample points. The output mesh interpolates the input points, ensuring that the reconstruction is accurate and detailed in feature-rich areas while being sparse in featureless regions. The algorithm's effectiveness is demonstrated through practical examples and theoretical analysis. It is shown that the crust algorithm can handle noisy data and undersampled regions, although it struggles with sharp edges and boundaries. The authors propose heuristics to improve the reconstruction in these challenging cases and discuss future research directions, including handling noise, sharp edges, and boundaries, and exploring applications in mesh compression. The implementation of the algorithm is efficient, capable of handling large datasets in minutes, and robust to numerical issues due to the use of exact arithmetic in Delaunay triangulation.