Title of Your Paper

Anisotropic Meshing of Implicit Surfaces

Sergei Azernikov and Anath Fischer
Laboratory for CAD and LCE
Faculty of Mechanical Engineering
Technion - Israel Institute of Technology
Haifa, Israel
meranath@technion.ac.il

Contents

Abstract:

This paper proposes a new grid-based method for adaptive anisotropic meshing of implicit surfaces. Grid-based methods are considered the major technique for implicit surface meshing, mainly due to their efficiency and simplicity. However, these methods suffer from a number of inherent drawbacks, resulting from the fact that the imposed Cartesian grid is generally not well adapted to the iso-surface, either in size or in orientation. To overcome the above obstacles, we propose a new implicit surface meshing method. The main idea of this method is first to construct a geometric field which is induced by the shape of the surface. This geometric field represents the natural directions and grid cell size for each point in R³. Then, the imposed volumetric grid is deformed toward the object’s shape according to the produced geometric field. The iso-surface meshes can be extracted from the resulting adaptive grid by any conventional grid-based contouring technique, since the topology of the grid is not modified during the adaptation process. These meshes provide better approximation of the unknown surface and exhibit anisotropy, which is present in this surface.

1. Introduction

In nature, a surface usually appears implicitly as a boundary between two environments. In synthetic solid modeling, this implicit form is used as well, since CSG operations with implicit models are rather simple. Recently, several researchers have shown that high quality implicit surfaces can be reconstructed from unorganized point clouds [9, 22] and polygonal meshes [28]. For other applications however, i.e. rendering, analysis and manufacturing, an explicit representation may be beneficial. An implicit representation is usually converted into an explicit one by iso-surface extraction or contouring. Therefore, numerous studies on iso-surface extraction have been published during the last two decades [6]. The existing contouring methods can be classified into two main categories: volume-based methods and surface-based methods.

Volume-based methods impose a structured grid of voxels on the domain. These voxels can be cubical, as in the original Marching Cubes algorithm [20], tetrahedral [5] or octahedral [8], uniform [5] or adaptive [25]. Then, intersection points are found between the surface and edges of the grid cells. These points represent the geometry of the surface. Roughly speaking, this geometry can be interpreted as a projection or mapping of the grid onto the surface. This mapping makes the volume-based methods efficient and robust, since the topology of the surface is dictated a priori by the imposed grid. Moreover, all the portions of the iso-surface are detected. There is also an option for parallel computing acceleration, because each voxel can be processed independently.

Unfortunately, though efficient, volume-based methods often produce low quality meshes, since the imposed grid in general is not well adapted to the iso-surface shape. Many attempts have been made to solve this problem. The main research directions are 1) post-processing or remeshing of the low quality meshes [2, 23]; 2) surface-based methods [1]; and 3) grid adaptation [4, 21].

Surface-based methods start from a seed point on the surface and propagate the mesh by region growing. These methods produce high quality meshes because they consider the local differential geometry of the surface during the propagation process [1,12,18]. However, since the topology of the surface is recovered on the fly, dealing with complex high genus objects might by difficult. Consider, for example, the object shown in Figure 1, which is constructed from several disjoint parts. Extensive investigation of this problem can be found in [26].

The grid adaptation approach for implicit surface contouring was introduced by Moore and Warren [21]. In this work, the imposed grid vertices are perturbed in order to reduce the number of triangles produced by the Marching Cubes algorithm. However, this perturbation is indifferent to the surface geometry. Balmelli et al. [4] proposed to warp the grid according to an importance map. This map makes grid vertices attract to regions of interest chosen by the user. As a result, the extracted iso-surface mesh is denser in these regions and sparser elsewhere. The importance map constitutes a scalar metric field, which guides the adaptation of the grid. Since the field is scalar, the adaptation is isotropic. Recently, Tchon et al. [27] showed that employing an anisotropic geometric field for grid adaptation is very natural and effective. These authors proposed an effective and simple method to adapt a uniform volumetric grid to the triangulated surface. The adaptation is guided by the geometric field, which is constructed from the tensor of curvature and thickness of the solid bounded by the surface. The adapted grid tends to follow the shape of the surface and clearly reflects anisotropic properties of this shape. These results suggest that a geometric field is a very natural and general representation of a shape for applications, where anisotropic properties are important. Indeed, these properties must be considered during high quality surface meshing [14]. Therefore, we employ the anisotropic grid adaptation approach for adaptive meshing of implicit surfaces.


(a)	(b)

Fig. 1. Adaptive anisotropic meshing of a complex implicit solid constructed from three disjoint rings: (a) iso-surface extracted from the uniform Cartesian grid 25x25x25 (1776 faces); (b) iso-surface extracted from the same grid after proposed geometric adaptation (2374 faces).

Our adaptation method is based on the one proposed by Tchon et al. [27], with several important differences due to a completely distinct application. Tchon et al. proposed their method in the context of hexahedral meshing of volumes bounded by a given polygonized surface. In contrast, we are interested in meshing implicit surfaces. Therefore, we must deal with a number of specific problems:

· How to construct the geometric tensor field from an implicitly defined surface?

· How to sample the produced tensor field when the region of interest is the surface only?

In the following sections we describe the proposed approach in detail and give answers to the stated questions.

2. The Approach

In the current work we propose to enhance the performance of grid-based implicit surface contouring methods by adapting the imposed volumetric grid with the anisotropic metric field induced by this surface shape. The motivating idea behind this approach is simple. In terms of signal processing, the imposed grid is a sampling of R³. For the homogeneous and isotropic Euclidean space, uniform sampling is optimal. In our case, however, the presence of the implicit surface introduces curvature into the space. As a result, the space is no longer Euclidean but rather Riemannian. Indeed, the uniform grid in Euclidean space is non-uniform in the Riemannian space. Consequently, this grid should be warped in order to provide uniform sampling of the anisotropic space.

The proposed method consists of four main phases, as shown in Figure 2. The process begins by constructing the background mesh. Then, the geometric tensor field is evaluated and propagated on this mesh. As a result, a metric tensor is set for each point in the problem domain. The grid is adapted geometrically by relaxing the vertex position, while edge length is calculated in the produced anisotropic metric. Finally, the iso-surface mesh is extracted using one of the standard grid-based contouring techniques.

Fig. 2. Flow chart of the proposed implicit surface meshing method.

This approach has the following important advantages:

· The geometric adaptation is shape-driven and not uniform or axis-aligned, as with octrees. Therefore, complex geometric features can be recovered with a lower number of voxels.

· The geometric adaptation only perturbs the grid’s geometry, while maintaining the structured topology. That is, the adaptation can be used as a preprocessor for any existing grid-based contouring method.

· The meshes extracted from the adapted grids are of much higher quality than those extracted from the original uniform grids with the same number of voxels. Moreover, dual contouring [3, 17] of the adapted grids produces all-quad, anisotropic meshes.

· Recovery of all parts of the iso-surface, connected or not, is guaranteed. This is not trivial to achieve using surface-based methods [12].

· In analysis applications, the behavior of some physical phenomena often must be incorporated into the meshing process. This is supported naturally by the proposed approach, by applying metric intersection [14].

Since the geometric tensor field definition requires notation from the differential geometry of implicit surfaces, we start with a brief review of this notation.

3. Differential Geometry of Implicit Surfaces

Although the results of this section can be found in [16], we summarize them here for completeness. Consider F : R³ → R a scalar function defined over the domain of support Ω. Let S be a surface defined implicitly as S = {"x(x, y, z)∈ Ω|F (x)=0}. Our goal is to determine local differential properties for "x ∈ S. These properties are compactly encoded in the Weingarten map [13], or more intuitively in the shape operator L (v). This operator shows how the unit surface normal changes at the point x ∈ S in the tangent direction v. The surface S unit normal at x is equal to the normalized gradient of F at x,

(1)

Since the gradient field ÑF is defined over the whole domain W, it is possible to extend the unit normal field : Â³ → Â³ from S to W as follows,

(2)

(3)

By definition, L (v) is a covariant derivative of the normal field in direction v,

(4)

The Hessian matrix of F is defined as,

(5)

It is easy to show that,

(6)

Thus, the matrix form of the shape operator L (v) at point x is given by,

(7)

where (b₁, b₂) is an arbitrary basis of the tangent space of F at x. Having l(x) in hand we can calculate all the required differential properties at x, as summarized in Table 1.

Property	Expression
principal curvatures κ₁_,₂	eigenvalues of l(x)
principal directions	eigenvectors of l(x)
Gaussian curvature K	det[l(x)]
mean curvature H	1/2 trace[l(x)]

Table 1. Differential properties of F at x
that can be extracted from the shape operator matrix form l(x).

These properties set a tensor of curvature ,

(8)

This tensor is a core of the geometric tensor defined in Section 4. Figure 3 presents evaluated tensors of curvature on torus.


(a)	(b)

Fig. 3. Tensors of curvature can be efficiently visualized as ellipsoids: (a) tensor of curvature represented as ellipse; (b) tensors of curvature evaluated for the rings model.

4. Geometric Tensor Field

Using the results of the previous section, we are able to formally define the geometric field as a Riemannian metric induced by the shape of the surface S on the domain Ω. For each point x∈Ω this metric defines a tensor M_3×3 as follows:

(9)

where Â is an orthonormal rotation matrix and Λ is a diagonal scaling matrix. This tensor defines the mapping of the desired hexahedral voxel to a unit axis-aligned cube for "x Î Ω[14].

As in [27], we associate the lower 2×2 minor of M with the tensor of curvature C of the surface. But for the normal direction we give a different interpretation, which is more appropriate in our context. Tchon et al. [27] deduce the target size of an element at x in the normal direction from the local thickness of the solid enclosed by S. In terms of the normal field produced by S, this is twice the distance from x to the closest singular point of . This thickness is approximated from the digital skeleton of the solid bounded by S. In our case, however, S does not necessarily bound a solid, and in any case we are interested only in a narrow band surrounding S and not in the bounded volume. In our implementation, the size of this band is set to twice the surface-voxel intersection evaluation precision e.

The geometric tensor components can be now rewritten as follows:

(10)

(11)

where are the principal directions which span the tangent space of the surface, is the normal vector and α is a coefficient which allows control of the introduced anisotropy level.

5. Background Mesh Construction

For the geometric field to be represented faithfully, it should be appropriately sampled on the background mesh. Sampling it on a dense uniform grid is possible, but such an approach is extremely memory-intensive. A better solution is to sample this field adaptively. Octrees have been used for adaptive sampling of scalar distance fields by Frisken et al. [15], and geometric tensor fields by Tchon et al. [27]. Octrees are very suit