PhD proposal / 2.5D Frame Fields for Hexahedral Meshing
PhD thesis proposal: 2.5D Frame Fields for Hexahedral Meshing
PhD supervisors:
Dmitry Sokolov 
University of Lorraine 
Dmitry.Sokolov@univlorraine.fr 
Nicolas Ray 
INRIA 
Nicolas.Ray@inria.fr 

Proposal’s context
Numerical simulations solve physical equations to better understand complex phenomena. It is a key aspect in both industry and science: the idea is to replace costly physical experiments (e.g. wind tunnel in fluid dynamics) with computer simulation. This dramatically reduces the overall cost of product development, because it allows to prune hypotheses before undertaking expensive real actions such as the realization of a mechanical part, on site investments such as drilling of oil wells, and then submitting it to tests. One of the central points in numerical simulations is the ability to represent the functions (temperature, pressure, speed, etc.) on the studied objects. The most versatile way to do so is to discretize the domain of the interest into the cells of a mesh. But mesh generation is a notoriously difficult task. NASA indicates mesh generation as one of the major challenges for 2030, and estimates that it costs 80% of time and effort in numerical simulation [Slotnick2014]. The main difficulty is due to the need for constructing supports that match both the geometry and the physics of the system to be modeled. More specifically, hexahedral meshes are known to be very hard to produce, if only because of the “rigid” nature of hexahedra limiting the potential combination with their neighbors, but they exhibit unique properties that are mandatory for certain physics simulations (deformation mechanics, fluid dynamics, etc.).
Today, a majority of pure or highly dominant hexahedral meshes is still handmade, because no method has proven yet to be versatile and largely applicable enough with an acceptable accuracy, thus inducing large costs and time overheads.
A new and original approach to generate hexahedral meshes was proposed recently by the team [Sokolov2016, Ray2018]. It comes from the observation that good quality hexahedral meshes look like a deformed grid almost everywhere. The idea is to define a deformation of the object such that if the final hexahedral mesh (the result) undergoes this deformation, it will match a unit, axis aligned grid. The direct application of this idea computes this deformation, applies its inverse to the unit grid, and obtains a hexahedral mesh. In practice, we introduce more degrees of freedom by considering global parameterizations instead of deformation. In this case, the deformation functions have some discontinuities that make it possible to represent a much larger family of hexahedral meshes: the deformed grid can be cut and glued to itself in a nontrivial way.
Global parameterizations are constructed in two steps: the Frame Field (FF) step defines the orientation of the grid at each point of the domain, and the Cube Covering (CC) step generates a global parameterization that will align the final mesh with the orientation (FF) defined in the first step. When both steps succeed, the result is a very regular full hexahedral mesh. Unfortunately, the FF step may generate an orientation field with a topological structure (its singularity graph) that is not valid for the second step, and the CC step often generates invalid mappings with (locally) negative Jacobians.

The goal and the organisation of the project
Frame Fields are actually generated by a heuristic that minimizes the field rotation; it seems fair to minimize the bending of the deformed grid. However, it can produce singularity graphs from which it is impossible to derive a valid parameterization. In fact, valid singularity graphs are extremely constrained: the only valid singularity curves are following one direction of the frame field, and the two other directions define a 2D frame field that is singular on the singularity curve. This constraint can be enforced a posteriori by combinatorial optimizations [Li2012], but this strategy is only local. We want to study an alternative approach that will incorporate the integrability constraints directly during the field generation.
The frame field optimization basically propagates inside the domain the constraints that are defined on the boundary. These constraints enforce one axis of the frame field to be equal to the normal of the surface of the domain. The singularity graph of the frame comes from locations where it must deal with contradicting constraints, but most of them have only two of the three axis that have constraints: the last direction remains well defined. Our idea is to define this “stable” vector field in the domain, and constraint the other two directions to avoid rotating when being advected by this direction. This approach will mimic the sweeping operation done in hexahedral mesh modeling. To reach our objective, we need to define a stable direction field inside the domain, and to smooth the frame field such that it does not rotate when moving in the stable direction.

Our current frame field optimization algorithm relies on spherical harmonics to merge constraints inside the domain, and it will be easier to derive a locally stable vector field from this representation. The challenge will be to handle complex objects where this stable direction is not represented by the same axis everywhere in the domain.

Simply fixing the stable direction and minimizing the frame field curvature will not guarantee that singularities will align to the stable direction; we need to constraint the field to have minimal rotation when moving in the stable direction. We will try to enforce it by two approaches: in the solver by a strong anisotropy term or Lagrange multipliers, or by manipulating a new set of variables that naturally enforce the constraint.
We hope that frame fields generated by this task will exhibit less configurations that cannot be integrated in the CC step.
The project will be split into three main parts:

Months 112: preprocessing the model. Here, we will address the preprocessing of the model to help the GP. We will focus on conflicting configurations of the feature edges. The starting point will be to observe failure cases, define the regions that will be impacted, and remove them from the domain to be remeshed (the postprocessing will deal with it). Then, we will better manage these regions by replacing them by local hexahedral mesh patterns, and we will try to get an exhaustive classification of conflicting cases.
We will manage this problem at two levels: by proposing new solutions to local conflicting configurations, and by determining mesh configurations along feature edges that respect global constraints [Kowalski2012a]. The first aspect builds on top of the results of our results on polycubes, with more degrees of freedom to solve them. Indeed, inner singularities can prevent local conflicting cases to occur, and local meshing patterns can generate more types of constraints on the surrounding frame field.

Months 1224: Locally integrable Frame Fields. The singularities of a 2D Frame Field are simple to characterize: a point is singular if the integral of the field rotation along the edge of an (infinitely small) neighborhood of this point is not null. For quadrangular remeshing (small deformed squares), this implies the presence of an irregular vertex i.e. with a valence different from 4. This characterization unfortunately does not extend to the volume case because the integration of the 3D rotations is poorly defined. However, we know that the only singularity curves corresponding to a hexahedral remeshing are 2D singularities “extruded” in the third direction: they can become edges shared by a number of hexahedrons different from 4.
When you only find the minimum curvature Frame Field, other types of singularity appear (Figure 3right). The idea of this project is to decompose the computation of the Frame Field into two steps: the first one will define a stable direction, and the second one will calculate the two other directions. The stable direction can be imposed by the model, as it is the case in geology with the normal to the geological strata, but it can also be defined automatically by a heuristic. For example, we will try to take the direction of least curvature in a Frame Field (vertical in Figure 3). The two other directions will have to be computed using 2D Frame Field generation algorithms, but by imposing that the field remains constant along the flux lines of the stable direction.

Months 2436: Extract hexahedral mesh. The basic idea of using global parameterization is to define the hexahedral mesh as the preimage of the unit grid by a parameterization function. This process assumes that the parameterization function is locally one to one, with consistent orientation enforced by the positive Jacobian condition. In practice, this is very difficult to enforce, but robust extraction [Lyon2016] accepts some local constraint violation. It is very interesting to consider that this algorithm works well thanks to the variables that are constrained to be integers; this is the same set of variables that defines the combinatorial structure of the hexahedral mesh. We would like to improve the robustness of the process by completely ignoring the geometry of parameterization. In theory, we could extract the hexahedral mesh combinatory using only variables that are constrained to be integers. The geometry can always be optimized directly on the final mesh, that is likely to be a better discretization for this purpose.
If we consider the case of polycubes, the combinatorial structure of the hexahedral mesh is quite easy to extract; the length of the edges of the polycube are directly given by the difference of its vertices coordinates, that are constrained to be integers. It determines directly the combinatorial structure of the final mesh and let us experiment to directly optimize the geometry on the hexahedral mesh. From this point, we can start introducing a singularity graph and extend the combinatorial reconstruction with it. We will experiment it by introducing singularity graphs of increasing complexity until we can generalize the approach to any singularity graph. The risk is that the problem may be too difficult to tackle, but we will at least have a solution for simple cases, plus the geometric optimization of the hexahedral mesh geometry that could help to improve our other results.

References related to the project
[Kowalski2012a] N. Kowalski, F. Ledoux, M.L. Staten, S.J. Owen, Fun sheet matching: towards automatic block decomposition for hexahedral meshes, Engineering with Computers, 2012.
[Lyon2016] M. Lyon, D. Bommes, L. Kobbelt, HexEx: Robust Hexahedral Mesh Extraction, ACM Transactions on Graphics, 2016
[Ray2018] N. Ray, D. Sokolov, M. Reberol, F. Ledoux, B. Lévy, Hexdominant meshing: Mind the gap!, ComputerAided Design, Volume 102, 2018.
[Slotnick2014] J. Slotnick, A. Khodadoust, J. Alonso, D. Darmofal, G. William, L. Elizabeth and D. Mavriplis,CFD Vision 2030 Study: A Path to Revolutionary Computational Aerosciences, NASA/CR2014218178, NF1676L18332.
[Sokolov2016] D. Sokolov, N. Ray, L. Untereiner and B. Lévy. 2016. HexahedralDominant Meshing. ACM Transactions on Graphics 35.