Notes on implementation of surface input/reconstruction
The current algorithm for input of complex geometries in Gerris requires the user to provide a surface composed of connected triangles. The new algorithm is tolerant of many defects found in real-world CAD models (small gaps, self-intersections etc...) but only within reasonable limits: i.e. if holes (or self-intersections) are larger than the minimum cell size this will cause the algorithm to fail.
Also in some cases, users may not have a triangulated surface at their disposal but only a cloud of points and normals (e.g. from data acquired using a 3D scanner). Surface reconstruction techniques could be used to extract a triangulated surface from this dataset but this is a non-trivial operation which could also result in topologically incorrect surfaces (depending on the algorithm used).
Generally speaking, it is easy to convert a triangulated surface into a cloud of points and normals but the inverse problem is hard. So, the cloud point representation can be seen as a "minimal" representation easily obtained from any other representation.
Surface reconstruction techniques often use space-partitioning techniques (such as octrees). A very nice paper by F. Duguet et al describes this technique (Level of Detail Continuum for Huge Geometric Data).
Together with Stephane Zaleski, we have also been working on a mixed Volume-Of-Fluid, level set formulation to solve flows with interfaces and surface tension. An important element of the algorithm is an efficient way to compute a level-set function from the VOF interface reconstruction. In Duguet et al this is done (from the cloud of points and normals) by solving a vector diffusion equation (Gradient Vector Flow) and a Poisson equation to get the distance field from the vector field. Although this is done using a multigrid algorithm on the octree (inspired by Gerris!), I think our technique (based on direct hierarchical traversal) could be faster and would yield a closer estimate of the true Euclidean distance (Duguet et al technique does not compute the true Euclidean distance).
Once the distance function is computed, it can be used almost directly in the current implementation of the solid fractions computation: the intersections between cell edges and the surface can be computed using the distance function (interpolated at the cell vertices) rather than using the triangular facets.
To sum up, the algorithm would be as follows:
- Compute distance and normal info on each leaf cell of the octree by inserting (serially) points/normals from the input point cloud (as per Duguet et al) ("triangle clouds" could also be used)
- Compute signed distance function using the hierarchical technique of gerris-csf
- Compute solid fractions using the current algorithm
The advantages would be:
- No topological or geometrical constraints on the input surface: any model (triangulated surfaces, NURBS, cloud of points etc...) could be imported directly into Gerris
- The signed distance function could also be useful for other things: collision detection, differential geometry etc... Also, it is sometimes useful to be able to "displace" the boundary perpendicularly to itself (so-called "displacement surface"). This is used for example to approximate the effects of finite-thickness viscous boundary layers (see for example Applications of a Cartesian Mesh Boundary-Layer Approach for Complex Conﬁgurations). The distance function makes this easy.
- Unified representation of solid boundaries and interfaces (sort of)