Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Namespaces | Classes | Typedefs | Enumerations | Functions
CGALWrappers Namespace Reference

Namespaces

namespace  internal
 

Classes

struct  AdditionalData
 
struct  AdditionalData< 2 >
 
struct  FaceInfo2
 

Typedefs

using ConcurrencyTag = CGAL::Sequential_tag
 
using K = CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt
 
using K_exact = CGAL::Exact_predicates_exact_constructions_kernel
 
using CGALPolygon = CGAL::Polygon_2< K >
 
using Polygon_with_holes_2 = CGAL::Polygon_with_holes_2< K >
 
using CGALTriangle2 = K::Triangle_2
 
using CGALTriangle3 = K::Triangle_3
 
using CGALTriangle3_exact = K_exact::Triangle_3
 
using CGALPoint2 = K::Point_2
 
using CGALPoint3 = K::Point_3
 
using CGALPoint3_exact = K_exact::Point_3
 
using CGALSegment2 = K::Segment_2
 
using Surface_mesh = CGAL::Surface_mesh< K_exact::Point_3 >
 
using CGALSegment3 = K::Segment_3
 
using CGALSegment3_exact = K_exact::Segment_3
 
using CGALTetra = K::Tetrahedron_3
 
using CGALTetra_exact = K_exact::Tetrahedron_3
 
using Triangulation2 = CGAL::Triangulation_2< K >
 
using Triangulation3 = CGAL::Triangulation_3< K >
 
using Triangulation3_exact = CGAL::Triangulation_3< K_exact >
 
using Vb = CGAL::Triangulation_vertex_base_2< K >
 
using Fbb = CGAL::Triangulation_face_base_with_info_2< FaceInfo2, K >
 
using CFb = CGAL::Constrained_triangulation_face_base_2< K, Fbb >
 
using Fb = CGAL::Delaunay_mesh_face_base_2< K, CFb >
 
using Tds = CGAL::Triangulation_data_structure_2< Vb, Fb >
 
using Itag = CGAL::Exact_predicates_tag
 
using CDT = CGAL::Constrained_Delaunay_triangulation_2< K, Tds, Itag >
 
using Criteria = CGAL::Delaunay_mesh_size_criteria_2< CDT >
 
using Vertex_handle = CDT::Vertex_handle
 
using Face_handle = CDT::Face_handle
 

Enumerations

enum class  FacetTopology { facet_vertices_on_surface = CGAL::FACET_VERTICES_ON_SURFACE , facet_vertices_on_same_surface_patch , facet_vertices_on_same_surface_patch_with_adjacency_check }
 
enum class  BooleanOperation { compute_corefinement = 1 << 0 , compute_difference = 1 << 1 , compute_intersection = 1 << 2 , compute_union = 1 << 3 }
 

Functions

template<int dim0, int dim1, int spacedim>
std::vector< std::array< Point< spacedim >, dim1+1 > > compute_intersection_of_cells (const typename Triangulation< dim0, spacedim >::cell_iterator &cell0, const typename Triangulation< dim1, spacedim >::cell_iterator &cell1, const Mapping< dim0, spacedim > &mapping0, const Mapping< dim1, spacedim > &mapping1, const double tol=1e-9)
 
template<int dim0, int dim1, int spacedim>
std::vector< std::array< Point< spacedim >, dim1+1 > > compute_intersection_of_cells (const ArrayView< const Point< spacedim > > &vertices0, const ArrayView< const Point< spacedim > > &vertices1, const double tol=1e-9)
 
template<typename CGALPointType , int dim, int spacedim>
void dealii_cell_to_cgal_surface_mesh (const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const ::Mapping< dim, spacedim > &mapping, CGAL::Surface_mesh< CGALPointType > &mesh)
 
template<typename CGALPointType , int dim, int spacedim>
void dealii_tria_to_cgal_surface_mesh (const ::Triangulation< dim, spacedim > &triangulation, CGAL::Surface_mesh< CGALPointType > &mesh)
 
template<typename C3t3 >
void cgal_surface_mesh_to_cgal_triangulation (CGAL::Surface_mesh< typename C3t3::Point::Point > &surface_mesh, C3t3 &triangulation, const AdditionalData< 3 > &data=AdditionalData< 3 >{})
 
template<typename CGALPointType >
void compute_boolean_operation (const CGAL::Surface_mesh< CGALPointType > &surface_mesh_1, const CGAL::Surface_mesh< CGALPointType > &surface_mesh_2, const BooleanOperation &boolean_operation, CGAL::Surface_mesh< CGALPointType > &output_surface_mesh)
 
template<typename CGALTriangulationType >
::Quadrature< CGALTriangulationType::Point::Ambient_dimension::value > compute_quadrature (const CGALTriangulationType &tria, const unsigned int degree)
 
template<int dim0, int dim1, int spacedim>
::Quadrature< spacedim > compute_quadrature_on_boolean_operation (const typename ::Triangulation< dim0, spacedim >::cell_iterator &cell0, const typename ::Triangulation< dim1, spacedim >::cell_iterator &cell1, const unsigned int degree, const BooleanOperation &bool_op, const Mapping< dim0, spacedim > &mapping0=(ReferenceCells::get_hypercube< dim0 >() .template get_default_linear_mapping< dim0, spacedim >()), const Mapping< dim1, spacedim > &mapping1=(ReferenceCells::get_hypercube< dim1 >() .template get_default_linear_mapping< dim1, spacedim >()))
 
template<int dim0, int dim1, int spacedim>
::Quadrature< spacedim > compute_quadrature_on_intersection (const typename ::Triangulation< dim0, spacedim >::cell_iterator &cell0, const typename ::Triangulation< dim1, spacedim >::cell_iterator &cell1, const unsigned int degree, const Mapping< dim0, spacedim > &mapping0, const Mapping< dim1, spacedim > &mapping1)
 
template<typename CGALPointType >
void remesh_surface (CGAL::Surface_mesh< CGALPointType > &surface_mesh, const AdditionalData< 3 > &data=AdditionalData< 3 >{})
 

Detailed Description

Interface to the Computational Geometry Algorithm Library (CGAL).

CGAL is a software project that provides easy access to efficient and reliable geometric algorithms. The library offers data structures and algorithms like triangulations, Voronoi diagrams, Boolean operations on polygons and polyhedra, point set processing, arrangements of curves, surface and volume mesh generation, geometry processing, alpha shapes, convex hull algorithms, shape reconstruction, AABB and KD trees...

You can learn more about the CGAL library at https://www.cgal.org/

Typedef Documentation

◆ ConcurrencyTag

using CGALWrappers::ConcurrencyTag = typedef CGAL::Sequential_tag

Definition at line 83 of file utilities.h.

◆ K

using CGALWrappers::K = typedef CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt

Definition at line 64 of file intersections.cc.

◆ K_exact

using CGALWrappers::K_exact = typedef CGAL::Exact_predicates_exact_constructions_kernel

Definition at line 65 of file intersections.cc.

◆ CGALPolygon

using CGALWrappers::CGALPolygon = typedef CGAL::Polygon_2<K>

Definition at line 66 of file intersections.cc.

◆ Polygon_with_holes_2

using CGALWrappers::Polygon_with_holes_2 = typedef CGAL::Polygon_with_holes_2<K>

Definition at line 67 of file intersections.cc.

◆ CGALTriangle2

using CGALWrappers::CGALTriangle2 = typedef K::Triangle_2

Definition at line 68 of file intersections.cc.

◆ CGALTriangle3

using CGALWrappers::CGALTriangle3 = typedef K::Triangle_3

Definition at line 69 of file intersections.cc.

◆ CGALTriangle3_exact

using CGALWrappers::CGALTriangle3_exact = typedef K_exact::Triangle_3

Definition at line 70 of file intersections.cc.

◆ CGALPoint2

using CGALWrappers::CGALPoint2 = typedef K::Point_2

Definition at line 71 of file intersections.cc.

◆ CGALPoint3

using CGALWrappers::CGALPoint3 = typedef K::Point_3

Definition at line 72 of file intersections.cc.

◆ CGALPoint3_exact

using CGALWrappers::CGALPoint3_exact = typedef K_exact::Point_3

Definition at line 73 of file intersections.cc.

◆ CGALSegment2

using CGALWrappers::CGALSegment2 = typedef K::Segment_2

Definition at line 74 of file intersections.cc.

◆ Surface_mesh

using CGALWrappers::Surface_mesh = typedef CGAL::Surface_mesh<K_exact::Point_3>

Definition at line 75 of file intersections.cc.

◆ CGALSegment3

using CGALWrappers::CGALSegment3 = typedef K::Segment_3

Definition at line 76 of file intersections.cc.

◆ CGALSegment3_exact

using CGALWrappers::CGALSegment3_exact = typedef K_exact::Segment_3

Definition at line 77 of file intersections.cc.

◆ CGALTetra

using CGALWrappers::CGALTetra = typedef K::Tetrahedron_3

Definition at line 78 of file intersections.cc.

◆ CGALTetra_exact

using CGALWrappers::CGALTetra_exact = typedef K_exact::Tetrahedron_3

Definition at line 79 of file intersections.cc.

◆ Triangulation2

using CGALWrappers::Triangulation2 = typedef CGAL::Triangulation_2<K>

Definition at line 80 of file intersections.cc.

◆ Triangulation3

using CGALWrappers::Triangulation3 = typedef CGAL::Triangulation_3<K>

Definition at line 81 of file intersections.cc.

◆ Triangulation3_exact

using CGALWrappers::Triangulation3_exact = typedef CGAL::Triangulation_3<K_exact>

Definition at line 82 of file intersections.cc.

◆ Vb

using CGALWrappers::Vb = typedef CGAL::Triangulation_vertex_base_2<K>

Definition at line 96 of file intersections.cc.

◆ Fbb

using CGALWrappers::Fbb = typedef CGAL::Triangulation_face_base_with_info_2<FaceInfo2, K>

Definition at line 97 of file intersections.cc.

◆ CFb

using CGALWrappers::CFb = typedef CGAL::Constrained_triangulation_face_base_2<K, Fbb>

Definition at line 98 of file intersections.cc.

◆ Fb

using CGALWrappers::Fb = typedef CGAL::Delaunay_mesh_face_base_2<K, CFb>

Definition at line 99 of file intersections.cc.

◆ Tds

using CGALWrappers::Tds = typedef CGAL::Triangulation_data_structure_2<Vb, Fb>

Definition at line 100 of file intersections.cc.

◆ Itag

using CGALWrappers::Itag = typedef CGAL::Exact_predicates_tag

Definition at line 101 of file intersections.cc.

◆ CDT

using CGALWrappers::CDT = typedef CGAL::Constrained_Delaunay_triangulation_2<K, Tds, Itag>

Definition at line 102 of file intersections.cc.

◆ Criteria

using CGALWrappers::Criteria = typedef CGAL::Delaunay_mesh_size_criteria_2<CDT>

Definition at line 103 of file intersections.cc.

◆ Vertex_handle

using CGALWrappers::Vertex_handle = typedef CDT::Vertex_handle

Definition at line 104 of file intersections.cc.

◆ Face_handle

using CGALWrappers::Face_handle = typedef CDT::Face_handle

Definition at line 105 of file intersections.cc.

Enumeration Type Documentation

◆ FacetTopology

enum class CGALWrappers::FacetTopology
strong
Enumerator
facet_vertices_on_surface 

Each vertex of the facet have to be on the surface, on a curve, or on a corner.

facet_vertices_on_same_surface_patch 

The three vertices of a facet belonging to a surface patch s have to be on the same surface patch s, on a curve or on a corner.

facet_vertices_on_same_surface_patch_with_adjacency_check 

The three vertices of a facet belonging to a surface patch s have to be on the same surface patch s, or on a curve incident to the surface patch s or on a corner incident to the surface patch s.

Definition at line 32 of file additional_data.h.

◆ BooleanOperation

enum class CGALWrappers::BooleanOperation
strong

An enum type given to functions that compute boolean operations between geometrical objects, defined by triangulated surface grids.

As an example, we show all supported boolean operations applied to a hedra (the union of two pyramids) and a cube.

Enumerator
compute_corefinement 

Given two triangulated surfaces, refine the first surface until its intersection with the second surface is captured exactly by the refinement

compute_difference 

Given two triangulated surfaces, compute the boolean difference of the first surface minus the second surface

compute_intersection 

Given two triangulated surfaces, compute their intersection

compute_union 

Given two triangulated surfaces, compute their union

Definition at line 95 of file utilities.h.

Function Documentation

◆ compute_intersection_of_cells() [1/2]

template<int dim0, int dim1, int spacedim>
std::vector< std::array< Point< spacedim >, dim1+1 > > CGALWrappers::compute_intersection_of_cells ( const typename Triangulation< dim0, spacedim >::cell_iterator &  cell0,
const typename Triangulation< dim1, spacedim >::cell_iterator &  cell1,
const Mapping< dim0, spacedim > &  mapping0,
const Mapping< dim1, spacedim > &  mapping1,
const double  tol = 1e-9 
)

Given two deal.II affine cells, compute their intersection and return a vector of simplices, each one identified by an array of deal.II Points. All the simplices together are a subdivision of the intersection. If cells are non-affine, a geometrical error is introduced. If the measure of one of the simplices is below a certain threshold which defaults to 1e-9, then it is discarded. In case the two cells are disjoint, an empty array is returned.

Note
If the two cells have different intrinsic dimensions, than the first iterator is assumed to have with higher intrinsic dimension. The same holds for the two Mapping objects.
Parameters
cell0Iterator to the first cell.
cell1Iterator to the second cell.
mapping0Mapping for the first cell.
mapping1Mapping for the second cell.
tolThreshold to decide whether or not a simplex is included.
Returns
Vector of arrays, where each array identify a simplex by its vertices.

Definition at line 828 of file intersections.cc.

◆ compute_intersection_of_cells() [2/2]

template<int dim0, int dim1, int spacedim>
std::vector< std::array< Point< spacedim >, dim1+1 > > CGALWrappers::compute_intersection_of_cells ( const ArrayView< const Point< spacedim > > &  vertices0,
const ArrayView< const Point< spacedim > > &  vertices1,
const double  tol = 1e-9 
)

Same function as above, but working directly with vertices. It's assumed that the first vertices are coming from a Triangulation<dim0,spacedim>, while the other vertices are from a Triangulation<dim1,spacedim>, with dim0 > dim1.

Note
The vertices have to be given in CGAL order.

Definition at line 758 of file intersections.cc.

◆ dealii_cell_to_cgal_surface_mesh()

template<typename CGALPointType , int dim, int spacedim>
void CGALWrappers::dealii_cell_to_cgal_surface_mesh ( const typename ::Triangulation< dim, spacedim >::cell_iterator &  cell,
const ::Mapping< dim, spacedim > &  mapping,
CGAL::Surface_mesh< CGALPointType > &  mesh 
)

Build a CGAL::Surface_mesh from a deal.II cell.

The class Surface_mesh implements a halfedge data structure and can be used to represent polyhedral surfaces. It is an edge-centered data structure capable of maintaining incidence information of vertices, edges, and faces. Each edge is represented by two halfedges with opposite orientation. The orientation of a face is chosen so that the halfedges around a face are oriented counterclockwise.

More information on this class is available at https://doc.cgal.org/latest/Surface_mesh/index.html

The function will throw an exception in dimension one. In dimension two, it generates a surface mesh of the quadrilateral cell or of the triangle cell, while in dimension three it will generate the surface mesh of the cell, i.e., a polyhedral mesh containing the faces of the input cell.

The generated mesh is useful when performing geometric operations using CGAL::Polygon_mesh_processing, i.e., to compute boolean operations on cells, splitting, cutting, slicing, etc.

For examples on how to use the resulting CGAL::Surface_mesh see https://doc.cgal.org/latest/Polygon_mesh_processing/

Parameters
[in]cellThe input deal.II cell iterator
[in]mappingThe mapping used to map the vertices of the cell
[out]meshThe output CGAL::Surface_mesh

◆ dealii_tria_to_cgal_surface_mesh()

template<typename CGALPointType , int dim, int spacedim>
void CGALWrappers::dealii_tria_to_cgal_surface_mesh ( const ::Triangulation< dim, spacedim > &  triangulation,
CGAL::Surface_mesh< CGALPointType > &  mesh 
)

Convert a deal.II triangulation to a CGAL::Surface_mesh. The output depends on the intrinsic dimension of the input deal.II triangulation.

In 2d, i.e. with a Triangulation<2> or a Triangulation<2,3>, the output is the CGAL::Surface_mesh describing the whole triangulation.

In 3d, the boundary the of the deal.II Triangulation is converted to a CGAL::Surface_mesh by looping over all the boundary faces.

Parameters
[in]triangulationThe input deal.II triangulation.
[out]meshThe output CGAL::Surface_mesh.

◆ cgal_surface_mesh_to_cgal_triangulation()

template<typename C3t3 >
void CGALWrappers::cgal_surface_mesh_to_cgal_triangulation ( CGAL::Surface_mesh< typename C3t3::Point::Point > &  surface_mesh,
C3t3 &  triangulation,
const AdditionalData< 3 > &  data = AdditionalData< 3 >{} 
)

Given a closed CGAL::Surface_mesh, this function fills the region bounded by the surface with tets.

Parameters
[in]surface_meshThe (closed) surface mesh bounding the volume that has to be filled.
[out]triangulationThe output triangulation filled with tetrahedra.
[in]dataAdditionalData object to pass to the CGAL::make_mesh_3 function. See the documentation of that struct for a description of those parameters.

◆ compute_boolean_operation()

template<typename CGALPointType >
void CGALWrappers::compute_boolean_operation ( const CGAL::Surface_mesh< CGALPointType > &  surface_mesh_1,
const CGAL::Surface_mesh< CGALPointType > &  surface_mesh_2,
const BooleanOperation boolean_operation,
CGAL::Surface_mesh< CGALPointType > &  output_surface_mesh 
)

Given two triangulated surface meshes that bound two volumes, execute a boolean operation on them, and store the result in a third surface mesh.

Quoting from CGAL documentation (https://doc.cgal.org/latest/Polygon_mesh_processing/index.html#title14):

‍Given a closed triangulated surface mesh, each connected component splits the 3d space into two subspaces. The vertex sequence of each face of a component is seen either clockwise or counterclockwise from these two subspaces. The subspace that sees the sequence clockwise (resp. counterclockwise) is on the negative (resp. positive) side of the component.

‍Given a closed triangulated surface mesh surface with no self-intersections, we say that surface bounds a volume if each subspace lies exclusively on the positive (or negative) side of all the incident connected components of surface. The volume bounded by surface is the union of all subspaces that are on negative sides of their incident connected components of surface.

‍There is no restriction on the topology of the input volumes. However, there are some requirements on the input to guarantee that the operation is possible. First, the input meshes must not self-intersect. Second, the operation is possible only if the output can be bounded by a manifold triangulated surface mesh. In particular this means that the output volume has no part with zero thickness. Mathematically speaking, the intersection with an infinitesimally small ball centered in the output volume is a topological ball. At the surface level this means that no non-manifold vertex or edge is allowed in the output. For example, it is not possible to compute the union of two cubes that are disjoint but sharing an edge.

See BooleanOperation for a list of available operations, with the corresponding examples.

Parameters
[in]surface_mesh_1The first surface mesh.
[in]surface_mesh_2The second surface mesh.
[in]boolean_operationSee BooleanOperation for the list of the allowed operations.
[out]output_surface_meshThe surface mesh with the result of the boolean operation.

◆ compute_quadrature()

template<typename CGALTriangulationType >
::Quadrature< CGALTriangulationType::Point::Ambient_dimension::value > CGALWrappers::compute_quadrature ( const CGALTriangulationType &  tria,
const unsigned int  degree 
)

Given a CGAL Triangulation describing a polyhedral region, create a Quadrature rule to integrate over the polygon by looping through all the vertices and exploiting QGaussSimplex.

Parameters
[in]triaThe CGAL triangulation object describing the polyhedral region.
[in]degreeDesired degree of the Quadrature rule on each element of the polyhedral.
Returns
A global Quadrature rule that allows to integrate over the polyhedron.

◆ compute_quadrature_on_boolean_operation()

template<int dim0, int dim1, int spacedim>
::Quadrature< spacedim > CGALWrappers::compute_quadrature_on_boolean_operation ( const typename ::Triangulation< dim0, spacedim >::cell_iterator &  cell0,
const typename ::Triangulation< dim1, spacedim >::cell_iterator &  cell1,
const unsigned int  degree,
const BooleanOperation bool_op,
const Mapping< dim0, spacedim > &  mapping0 = (ReferenceCells::get_hypercube< dim0 >() .template get_default_linear_mapping< dim0, spacedim >()),
const Mapping< dim1, spacedim > &  mapping1 = (ReferenceCells::get_hypercube< dim1 >() .template get_default_linear_mapping< dim1, spacedim >()) 
)

Compute a Quadrature formula over the polygonal/polyhedral region described by a BooleanOperation between two deal.II cell_iterator objects. The degree of the Quadrature rule is specified by the argument degree. This function splits the polygon/polyhedron into simplices and collects QGaussSimplex quadrature rules.

Parameters
[in]cell0A cell_iterator to the first deal.II cell.
[in]cell1A cell_iterator to the second deal.II cell.
[in]degreeThe degree of accuracy you wish to get for the global quadrature formula.
[in]bool_opThe BooleanOperation to be performed.
[in]mapping0Mapping object for the first cell.
[in]mapping1Mapping object for the first cell.
Returns
The global quadrature rule on the polygon/polyhedron.

◆ compute_quadrature_on_intersection()

template<int dim0, int dim1, int spacedim>
::Quadrature< spacedim > CGALWrappers::compute_quadrature_on_intersection ( const typename ::Triangulation< dim0, spacedim >::cell_iterator &  cell0,
const typename ::Triangulation< dim1, spacedim >::cell_iterator &  cell1,
const unsigned int  degree,
const Mapping< dim0, spacedim > &  mapping0,
const Mapping< dim1, spacedim > &  mapping1 
)

A specialization of the function above when the BooleanOperation is an intersection. The rationale behind this specialization is that deal.II affine cells are convex sets, and as the intersection of convex sets is itself convex, this function internally exploits this to use a cheaper way to mesh the inside.

Parameters
[in]cell0A cell_iterator to the first deal.II cell.
[in]cell1A cell_iterator to the second deal.II cell.
[in]mapping0Mapping object for the first cell.
[in]mapping1Mapping object for the first cell.
[in]degreeThe degree of accuracy you wish to get for the global quadrature formula.
Returns
The global quadrature rule on the polygon/polyhedron.

◆ remesh_surface()

template<typename CGALPointType >
void CGALWrappers::remesh_surface ( CGAL::Surface_mesh< CGALPointType > &  surface_mesh,
const AdditionalData< 3 > &  data = AdditionalData< 3 >{} 
)

Remesh a CGAL::Surface_mesh.

If the domain has 1-dimensional exposed features, the criteria includes a sizing field to guide the sampling of 1-dimensional features with protecting balls centers.

  • CGAL::parameters::edge_size.

This utility should be exploited to improve the quality of a mesh coming from boolean operations. As an example, consider two CGAL::Surface_mesh objects describing two hyper_spheres that intersect non-trivially. After computing their corefinement and union with compute_boolean_operation(), the resulting mesh is the following:

Clearly, elements (triangles) like the ones arising along the intersection of the two spheres should be avoided as they're quite degenerate and would decrease the accuracy of the results. A call to the present function will try to remesh the surface, according to the given named parameters, to improve its quality. In the present case, the result is the following:

Parameters
surface_meshThe input CGAL::Surface_mesh.
[in]dataAdditionalData object to pass to the CGAL::make_mesh_3 function. See the documentation of that struct for a description of those parameters.