CgalMesh for dtk is the name of an Action de Développement Technologique, aka ADT, that aims at simplifying the use of the Cgal Mesh features so that Inria researchers and engineers can handle them in a quite fast and easy way.
After six months, it appears that using Cgal Mesh features in a runtime fashion is not so straightforward. In order to make it feasible, we have been developping a kind of wrapper on top of Cgal Mesh named cgalMeshBridge. This article details the problem that we tackle and the manner we use to circumvent it.
CgalMesh brief review
Cgal has been developped in the manner of the Standard Template Library (STL) and makes the use of generic programming. Generic programming is based on the notion of concept that tells the methods that a class must implement so as to satisfy the concept. Such classes are called models of the concept.
Cgal Mesh package features several concepts for generation of isotropic simplicial meshes discretizing 3D domains. Among these concepts, the main ones, from a end-user point of view, are:
- Mesh Domain that describes the knowledge required on the object to be discretized
- Mesh Criteria that encapsulates the criteria to tune the meshing process
- Mesh Complex In Triangulation that contains the result of the meshing process
Let us consider the following code that enables to generate the volume mesh of an ellipsoid defined by an implicit function.
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Mesh_triangulation_3.h>
#include <CGAL/Mesh_complex_3_in_triangulation_3.h>
#include <CGAL/Mesh_criteria_3.h>
#include <CGAL/make_mesh_3.h>
#include <CGAL/perturb_mesh_3.h>
#include <CGAL/exude_mesh_3.h>
#include <CGAL/Implicit_mesh_domain_3.h>
// /////////////////////////////////////////////////////////////////
//Kernel or K
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef K::FT FT;
typedef K::Point_3 Point;
// Mesh domain or MD
typedef FT (Function)(const Point&);
typedef CGAL::Implicit_mesh_domain_3<Function, K> MD;
// Complex Triangulation or C3T3
typedef CGAL::Mesh_triangulation_3<MD>::type Tr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3T3;
// Mesh Criteriaor or MC
typedef CGAL::Mesh_criteria_3<Tr> MC;
// /////////////////////////////////////////////////////////////////
FT ellipsoid_function (const Point& p)
{
const FT x2 = p.x() * p.x();
const FT y2 = p.y() * p.y();
const FT z2 = p.z() * p.z();
return (x2 + 2 * y2 + 4 * z2 - 1);
}
// /////////////////////////////////////////////////////////////////
int main (int argc, char** argv)
{
MD domain(ellipsoid_function, K::Sphere_3(CGAL::ORIGIN, 2.));
MC criteria(facet_angle=25, facet_size=0.15, facet_distance=0.008,
cell_radius_edge_ratio=3);
C3T3 c3t3 = CGAL::make_mesh_3<C3T3>(domain, criteria);
std::ofstream medit_file("implicit_mesh_domain.mesh");
c3t3.output_to_medit(medit_file);
medit_file.close();
return 1;
}
// /////////////////////////////////////////////////////////////////The above code involves the three concepts we talked about. Firstly, it creates an implicit mesh domain based on a function that describes an ellipsoid:
MD domain(ellipsoid_function, K::Sphere_3(CGAL::ORIGIN, 2.));Then, it defines the Mesh Criteria:
MC criteria(facet_angle=25, facet_size=0.15, facet_distance=0.008,
cell_radius_edge_ratio=3);Eventually, the mesh is generated and stored into a Mesh Complex:
C3T3 c3t3 = CGAL::make_mesh_3<C3T3>(domain, criteria);It is thus quite simple to create a mesh using only these three concepts. Nevertheless, one can notice the significant amount of lines of code that are needed to write the three lines that do the real job. These additional lines are not straighforward to write for people having few software engineering skills. Moreover, this code is not fully reusable to mesh an object that is not described by an implicit function. In order to make CGal Mesh easier to use by Inria researchers and engineers, it is necessary to provide new ways to handle it. In the current ADT, one of these ways is to resort to the visual programming framework of dtk.
CgalMesh and dtkComposer: impossible union ?
Polymorphic requirements vs generic programming
The dtk visual programming framework, also called dtkComposer, enables:
- to draw workflows in a runtime application
- to execute these workflows on the local machine or on remote cluster
- to save a workflow and insert it into a more complex one
- etc.
In the case of CgalMesh, we would like to draw a pipeline, or a composition, that chains a meshing process followed by a mesh refiner then a mesh optimizer as shown in the image below:
The mesher node takes a mesh domain as input data, then transfers a mesh complex to the refiner node that sends this mesh complex to the optimizer. The connexions between the nodes are carried out by edges connecting output ports to input ports. Each of these ports are typed according to the type of the data they handle. In the manner of typed programming languages, connecting two ports that handle different types is not allowed. For instance, an output port sending a string cannot be connected to an input port waiting for a mesh complex. And this where the problem arises since it is therefore unlikely to connect a port handling a mesh complex templated by a polyhedral mesh domain to a port dealing with a mesh complex templated by an implicit mesh domain because they do not handle the same type ! We have to deal with the same issue as far as the input of the mesher node is concerned: we cannot have a unique input port to handle different kinds of mesh domain. Having a port for each kind of mesh domain is not relevant, it would imply to modify the interface of the node every time a new mesh domain is introduced. This issue makes the visual programming inefficient to handle CgalMesh functionalities.
When using object oriented programming, such a difficulty is solved by resorting to a common parent class that is used as unique type for the ports. Unfortunately, when using generic programming paradigm, the notion of concept defines the API that a class must fullfill to be a model of this concept. There is no need of a common parent class even if it is not forbidden. In practice, this is how the three main concepts of CgalMesh were designed so there are no common parent class.
The problem of the chain of template parameters
In the previous section, we pointed out that dynamic polymorphism is not part of the three main concepts of CgalMesh. But another issue lays in the imbrication of the template parameters between the concepts. Let us take a look on the typedef that are written in the first example.
typedef FT (Function)(const Point&);
typedef CGAL::Implicit_mesh_domain_3<Function, K> MD;This typedef defines the kind of mesh domain and we notice that this latter has two template parameters, the type of kernel and the type of the implicit function which is itself templated by the kernel.
typedef CGAL::Mesh_triangulation_3<MD>::type Tr;The above typedef defines the kind of mesh triangulation. This latter relies on the type of the mesh domain as template parameter.
typedef CGAL::Mesh_complex_3_in_triangulation_3<Tr> C3T3;
typedef CGAL::Mesh_criteria_3<Tr> MC;The mesh complex and the mesh criteria are both templated by the mesh triangulation. The chain is thus highlighted and can be simplified as follows:
Input<K>, MD< Input<K> >, TR< MD< Input<K> > >, C3T3< TR< MD< Input<K> > > >where Input is the type governing the mesh domain, e.g. an implicit function, a polyhedron, a nurbs, etc.
Such a chain of nested template parameters has two main drawbacks:
-
changing the type of the input modifies the type of all the concepts in the chain.
- in a runtime environment, such as the visual programming framework, user cannot make generic pipelines since the type of the nodes and the type of the data that will travel through them would rely on the type of the input and there is no way to forecast all the cases.
From the software engineering point of view, the reasons for such a chain remain unclear since alternate design would have avoid it. Anyway, we had to come up with a solution in order to fullfill the objectives of the work. Hopefully, a feature of the next C++17 standard dealing with polymorphic allocators in STL containers helped us to find a remedy.
Bridge pattern to resist to templation
Polymorphic allocators
In the paper N3816 proposed to the C++ standard committees, Pablo Halpern makes the following observations about the inability to use allocators in non-generic contexts:
In large software systems, most of the application program consists of non-generic procedural or object-oriented code that is compiled once and linked many times. Allocators in C++, however, have historically relied solely on compile-time polymorphism, and therefore have not been suitable for use in vocabulary types, which are passed through interfaces between separately-compiled modules, because the allocator type necessarily affects the type of the object that uses it.
The last part should remind you some situation we encountered previously, shouldn’t it. To better understand the commonality, let us say that you conceived a class that provides a STL vector of double as return type of an accessor.
#include <vector>
class dummyObject
{
public:
const std::vector<double>& realValues(void) const;
};When no additional template argument is given, the STL vector uses the default allocator as follows:
template < typename T, typename Allocator = std::allocator<T> > class vector;Now let us say that you want to use a pool management of the memory in your dummyObject. You may try to replace the default allocator by a dummyPoolAllocator. The dummyObject then reads:
#include <vector>
template < typename T > class dummyPoolAllocator;
class dummyObject
{
public:
const std::vector<double, dummyPoolAllocator<double> >& realValues(void) const;
};Patatra ! You broke the binary compatibility between your library and the codes that use it. Yet, you just tried to improve your library, but to modify its implementation you had to modify its interface just because the allocator is a template parameter of the vector. This is why nobody used allocators and this fits perfectly with what Pablo Halpern says in his paper.
The proposition of Halpern is to replace the default allocator by a wrapper, or a bridge, around a pure abstract class defining the interface for all the derived allocators. Derived classes contain the machinery for actually allocating and deallocating memory according to user’s demand. Here is a very simplified version of this design:
template < typename T > class abstractAllocatorData
{
public:
virtual ~abstractAllocatorData(void) = default;
virtual T *allocate(std::size_t n) = 0;
virtual void deallocate(T *ptr) = 0;
};
template < typename T > class polymorphicAllocator
{
public:
polymorphicAllocator(abstractAllocatorData *data) : m_data(data) {}
T *allocate(std::size_t n) { return m_data->allocate(n); }
void deallocate(T *ptr) { m_data->deallocate(ptr); }
private:
abstractAllocatorData<T> *m_data;
};
template < typename T > class dummyAllocatorData : public abstractAllocatorData
{
public:
T *allocate(std::size_t n) { return my_dummy_allocation_process(); }
void deallocate(T *ptr) { my_dummy_deallocation_process(); }
}Now, it is likely to have a unique type of vector when setting polymorphicAllocator as the second template argument.
template < typename T > using generic_vector = vector< T, polymorphicAllocator<T> >;It follows that changing the way the memory is allocated does not modify the interface of the dummyObject and the binary compatibility is preserved. The goal is clearly reached.
Polymorphic mesh domain
We can learn two things from the bridge pattern of the polymorphic allocator:
-
firstly, the bridge enables to discouple the interface from the implementation.
- secondly, the use of the bridge as template parameter enables to keep the same type of the templated class whatever the implementation that is chosen by the user.
Breaking the chain
We can adapt the same strategy in the context of CgalMesh so as to circumvent the issue due to the chain of nested template parameters. Indeed, using the bridge pattern to design a polymorphic mesh domain, we can obtain a single chain of template parameters that do not vary whatever the mesh domain’s implementation.
As the input of the mesh domain (e.g. implicit function, polyhedron, nurbs, etc) is related to the implementation, the new polymorphic mesh domain can only rely on the kernel. We thus have:
template < typename K > class cgalPolymorphicMeshDomain
{
public:
... API to model the concept of Mesh Domain ...
};Of course, the cgalPolymorphicMeshDomain is a bridge and wraps a pure abstract class that serves as base class for all the implementations. But we will detail it further. For now, let us focus on the consequence of such a design on the rest of the chain.
typedef cgalPolymorphicMeshDomain<K> PMD;
typedef CGAL::Mesh_triangulation_3<PMD>::type PTr;
typedef CGAL::Mesh_complex_3_in_triangulation_3<PTr> PC3T3;
typedef CGAL::Mesh_criteria_3<PTr> PMC;In the manner of the polymorphic allocator, the polymorphic Mesh Domain used as template parameter enables to keep the same type for all the concepts that rely on it. We can now design separately-compiled modules that can exchange data that always have the same type.
Hence, provided we suceed in implementing it, the polymorphic Mesh Domain fullfills the requirements that the visual programming framework imposes to use Cgal Mesh.
The Bridge
The proposed solution aims at discoupling the implementation of any mesh domain classes (governed both by the input data and the way this data is treated) from the interface defined by the concept of Mesh Domain. In order to achieve this, we resort to the Bridge design pattern with some additional refinements.
The bridge has an API consistent with the concept of Mesh Domain.
template < typename K > class cgalPolymorphicMeshDomain
{
public:
typedef typename K::FT FT;
typedef K R;
typedef CGAL::Bbox_3 Bbox_3;
public:
typedef int Subdomain_index;
typedef boost::optional<Subdomain_index> Subdomain;
typedef std::pair<Subdomain_index, Subdomain_index> Surface_patch_index;
typedef boost::optional<Surface_patch_index> Surface_patch;
typedef boost::variant<Subdomain_index, Surface_patch_index> Index;
typedef CGAL::cpp11::tuple<Point_3, Index, int> Intersection;
public:
cgalPolymorphicMeshDomain() = default;
cgalPolymorphicMeshDomain(cgalAbstractMeshDomainData<K> *data);
public:
void setData(cgalAbstractMeshDomainData<K> *data);
cgalAbstractMeshDomainData<K> *data(void);
public:
struct Construct_initial_points;
Construct_initial_points construct_initial_points_object(void) const;
public:
struct Construct_intersection;
Construct_intersection construct_intersection_object(void) const;
public:
struct Do_intersect_surface;
Do_intersect_surface do_intersect_surface_object(void) const;
public:
struct Is_in_domain;
Is_in_domain is_in_domain_object(void) const;
public:
Bbox_3 bbox(void) const;
Index index_from_surface_patch_index(const Surface_patch_index& index) const;
Index index_from_subdomain_index(const Subdomain_index& index) const;
Surface_patch_index surface_patch_index(const Index& index) const;
Subdomain_index subdomain_index(const Index& index) const;
private:
cgalAbstractMeshDomainData<K> *d;
};As far as the implementation part is concerned, the bridge wraps a pure abstract class cgalAbstractMeshDomainData which serves as the base class for all the likely implementations of the mesh domain concept.
template < typename K > class cgalAbstractMeshDomainData
{
public:
typedef typename K::Point_3 Point_3;
typedef typename K::Segment_3 Segment_3;
typedef typename K::Ray_3 Ray_3;
typedef typename K::Line_3 Line_3;
public:
typedef typename K::FT FT;
typedef K R;
typedef CGAL::Bbox_3 Bbox_3;
public:
typedef int Subdomain_index;
typedef boost::optional<Subdomain_index> Subdomain;
typedef std::pair<Subdomain_index, Subdomain_index> Surface_patch_index;
typedef boost::optional<Surface_patch_index> Surface_patch;
typedef boost::variant<Subdomain_index, Surface_patch_index> Index;
typedef CGAL::cpp11::tuple<Point_3, Index, int> Intersection;
public:
typedef typename std::vector<std::pair<Point_3, Index> >::iterator Initial_points_iterator;
typedef std::back_insert_iterator<std::vector<std::pair<Point_3, Index> > > Initial_points_back_insert_iterator;
public:
virtual ~cgalAbstractMeshDomainData(void) {}
public:
virtual Initial_points_iterator eval_construct_initial_points(Initial_points_iterator pts, const int n) const = 0;
virtual Initial_points_iterator eval_construct_initial_points(Initial_points_iterator pts) const = 0;
virtual Initial_points_back_insert_iterator eval_construct_initial_points(Initial_points_back_insert_iterator pts, const int n) const = 0;
virtual Initial_points_back_insert_iterator eval_construct_initial_points(Initial_points_back_insert_iterator pts) const = 0;
public:
virtual Intersection eval_construct_intersection(Segment_3 s) const = 0;
virtual Intersection eval_construct_intersection(Ray_3 r) const = 0;
virtual Intersection eval_construct_intersection(Line_3 l) const = 0;
public:
virtual Surface_patch eval_do_intersect_surface(Segment_3 s) const = 0;
virtual Surface_patch eval_do_intersect_surface(Ray_3 r) const = 0;
virtual Surface_patch eval_do_intersect_surface(Line_3 l) const = 0;
public:
virtual Subdomain eval_is_in_domain(Point_3 p) const = 0;
public:
virtual Bbox_3 boundingBox(void) const = 0;
};This abstraction has several pure virtual methods to make the connection between the derived class and the bridge. Implementing all these methods in derived classes is a tedious task. That’s why, an intermediate class based on the CRTP pattern is used to simplify users’ own implementation of the mesh domain. The implementation therefore derived from the CRTP class and not from the abstraction directly.
template < typename Derive, typename K > class cgalMeshDomainDataCRTP : public cgalAbstractMeshDomainData<K>
{
public:
typedef typename K::Point_3 Point_3;
typedef typename K::Segment_3 Segment_3;
typedef typename K::Ray_3 Ray_3;
typedef typename K::Line_3 Line_3;
public:
typedef typename K::FT FT;
typedef K R;
typedef CGAL::Bbox_3 Bbox_3;
public:
typedef int Subdomain_index;
typedef boost::optional<Subdomain_index> Subdomain;
typedef std::pair<Subdomain_index, Subdomain_index> Surface_patch_index;
typedef boost::optional<Surface_patch_index> Surface_patch;
typedef boost::variant<Subdomain_index, Surface_patch_index> Index;
typedef CGAL::cpp11::tuple<Point_3, Index, int> Intersection;
public:
typedef typename std::vector<std::pair<Point_3, Index> >::iterator Initial_points_iterator;
typedef std::back_insert_iterator<std::vector<std::pair<Point_3, Index> > > Initial_points_back_insert_iterator;
public:
cgalMeshDomainDataCRTP(void);
~cgalMeshDomainDataCRTP(void);
public:
Initial_points_iterator eval_construct_initial_points(Initial_points_iterator pts, const int n) const final;
Initial_points_iterator eval_construct_initial_points(Initial_points_iterator pts) const final;
Initial_points_back_insert_iterator eval_construct_initial_points(Initial_points_back_insert_iterator pts, const int n) const final;
Initial_points_back_insert_iterator eval_construct_initial_points(Initial_points_back_insert_iterator pts) const final;
public:
Intersection eval_construct_intersection(Segment_3 s) const final;
Intersection eval_construct_intersection(Ray_3 r) const final;
Intersection eval_construct_intersection(Line_3 l) const final;
public:
Surface_patch eval_do_intersect_surface(Segment_3 s) const final;
Surface_patch eval_do_intersect_surface(Ray_3 r) const final;
Surface_patch eval_do_intersect_surface(Line_3 l) const final;
public:
Subdomain eval_is_in_domain(Point_3 p) const final;
public:
Bbox_3 boundingBox(void) const final;
};Of course, it remains to deal with all the mesh domain classes that have been developped by so many people for a long time. The solution that consists in modifying theses classes by adding the CRTP as parent class is clearly unfeasible. Any modification of theses classes cannnot be allowed since it will break the compatibilty with all the existant codes using Cgal. So, one has to resort to composition to circumvent this issu. A wrapper has thus been written, it is templated by the mesh domain it wraps and it inherits from the CRTP class.
template < typename MD > class cgalMeshDomainDataWrapper : public cgalMeshDomainDataCRTP< cgalMeshDomainDataWrapper<MD>, typename CGAL::Kernel_traits<typename MD::R::Point_3>::type >
{
MD *m_md;
public:
typedef typename MD::R K;
public:
typedef typename K::FT FT;
typedef K R;
typedef CGAL::Bbox_3 Bbox_3;
public:
typedef typename MD::Subdomain_index Subdomain_index;
typedef typename MD::Subdomain Subdomain;
typedef typename MD::Surface_patch_index Surface_patch_index;
typedef typename MD::Surface_patch Surface_patch;
typedef typename MD::Index Index;
typedef typename MD::Intersection Intersection;
public:
cgalMeshDomainDataWrapper(MD *md);
~cgalMeshDomainDataWrapper(void);
public:
struct Construct_initial_points;
Construct_initial_points construct_initial_points_object(void) const;
public:
struct Construct_intersection;
Construct_intersection construct_intersection_object(void) const;
public:
struct Do_intersect_surface;
Do_intersect_surface do_intersect_surface_object(void) const;
public:
struct Is_in_domain;
Is_in_domain is_in_domain_object(void) const;
public:
Bbox_3 bbox(void) const;
};Findings
The reasons that drove these developments could seem to be some quite artificial and we agree with that partially. CgalMesh is used by a very large number of people so why would we like to change the way they use it? As said previously, the work we do for the ADT does not target people who are already convinced by CgalMesh but it aims at making it easier to use by phD and post-docs who do not master neither programming nor discrete geometry. Visual programming is a way to make CgalMesh easiest to use. The bridge pattern mixed with CRTP is the solution we propose. We have been using it for 3 months and the bridge was found to be interesting for other problems that remained unfixed until today. Upcoming posts will detail these unexpected additional possibilities.
Cheers !
Thibaud and Côme.
