Using libigl for cell cutting in immersed boundary methods

[manavbhatia]

Hi,

I work in scientific computations and one of my projects deals with level-set based immersed boundary finite element analysis. This requires the following:

• definition of a surface implicitly defined by a level set function equal to zero.
• intersection of this surface with a geometric object (triangle, quadrilateral in 2D and tetrahedron or hexahedron in 3D).
• creation of elements in the sub-domains of the geometric object on either side of the surface.

In continuation to my previous question #899 I am trying to identify if this should be done with Cgal or libigl. It appears that libigl uses Cgal for some functionality (not sure which ones), and I don't have the clarity to decide between the two. Would the finite-precision vs arbitrary precision computation become an issue for robustness?

What is the best way to accomplish the tasks identified above using libigl.

Thanks again!
-Manav

[alecjacobson]

libigl doesn't have exactly this, but has some of the components to do this. I know that others have used libigl's mesh booleans to do this (for clipping Voronoi cells in a power diagram to use the "power particles" fluid simulation method).

As far as I know, cgal doesn't have exactly what you want either.

So in either case you'll have to code up something from raw pieces. libigl uses very basic data-structures that make it very easy to move between geometric operations and the linear algebra and assembly code for FEM. CGAL tends to want to "own" your mesh data and augmenting its representation with your specific data can be cumbersome.

In either case, you'll make a decision about wither arbitrary precision is necessary (if you want to stitch together a mesh after cutting then probably yes, though FEM numerics might still freak out). libigl uses CGAL's rational arithmetic types under the hood.

[manavbhatia]

I am guessing that representing a quad or hex in libigl will be relatively straight forward.

Since my zero level-set surface is implicitly defined, what would be the best way to defined this in libigl? Is there a method that allows the creation of implicitly defined surfaces? Or do they need to be created explicitly.

If it is the latter, then should I find a set of points on the surface in my code and provide them to libigl to create a surface?

Would you recommend any other approach to define this surface?

[jdumas]

Most of libigl functions deal with triangle or tetreahedral meshes, which are represented using simple matrices. You could decide to use a similar approach for quad and hex, although you'd have to decide of a numbering convention for hexes probably. libigl also doesn't do anything special about implicit functions: in iso-surface extraction methods (such as marching cubes, or marching tets -- see method slice_tets I think), the user provides a vector where the function is evaluated on the vertices/points of interest.

[manavbhatia]

Do the tris and tets need to be with 3 nodes and 4 nodes, respectively? Or is there also support for high-order tris with 6 nodes and likewise with tets?

[jdumas]

Only linear triangles/tetrahedra.

[manavbhatia]

Thanks.
Do you know if this is also true for Cgal?

Seems like I will need to represent my surface as a set of tris before I can request the slice operations?

I am wondering if it is typical for geometric operations to be done on tris and tets. For example, in isogeometric finite element analysis it is typical for domains to be represented with NURBS. Is there a library out there (may be Cgal?) that can directly work with such high-order representations?

[jdumas]

Well you can always try to use libs like http://mfem.org/ that offer NURBS FEM. And if you want high order elements you should look into using a more specialized FEM library. libigl and CGAL are geometry processing libraries, so they operate mostly on linear triangle/tet meshes (for libigl), and while CGAL also supports polyhedral meshes and more complicated data structures. I don't think CGAL does quadratic geometric elements, it's already complex enough to deal with linear geometry as it is.

[alecjacobson]

https://github.com/mtao/mandoline