Skip to content

Commit

Permalink
Use AdditionalData and fix headers inclusion
Browse files Browse the repository at this point in the history
  • Loading branch information
fdrmrc committed Jun 7, 2022
1 parent f82669c commit be72338
Show file tree
Hide file tree
Showing 7 changed files with 226 additions and 156 deletions.
105 changes: 105 additions & 0 deletions include/deal.II/cgal/point_conversion.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,105 @@
// ---------------------------------------------------------------------
//
// Copyright (C) 2022 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------

#ifndef dealii_cgal_point_conversion_h
#define dealii_cgal_point_conversion_h

#include <deal.II/base/config.h>

#include <deal.II/base/point.h>

#ifdef DEAL_II_WITH_CGAL

# include <CGAL/Cartesian.h>
# include <CGAL/Exact_predicates_exact_constructions_kernel.h>
# include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
# include <CGAL/Simple_cartesian.h>


DEAL_II_NAMESPACE_OPEN
namespace CGALWrappers
{
/**
* Convert from a deal.II Point to any compatible CGAL point.
*
* @tparam CGALPointType Any of the CGAL point types
* @tparam dim Dimension of the point
* @param [in] p An input deal.II Point<dim>
* @return CGALPointType A CGAL point
*/
template <typename CGALPointType, int dim>
inline CGALPointType
dealii_point_to_cgal_point(const dealii::Point<dim> &p);

/**
* Convert from various CGAL point types to deal.II Point.
*
* @tparam dim Dimension of the point
* @tparam CGALPointType Any of the CGAL point types
* @param p An input CGAL point type
* @return dealii::Point<dim> The corresponding deal.II point.
*/
template <int dim, typename CGALPointType>
inline dealii::Point<dim>
cgal_point_to_dealii_point(const CGALPointType &p);


# ifndef DOXYGEN
// Template implementations

template <typename CGALPointType, int dim>
inline CGALPointType
dealii_point_to_cgal_point(const dealii::Point<dim> &p)
{
constexpr int cdim = CGALPointType::Ambient_dimension::value;
static_assert(dim <= cdim, "Only dim <= cdim supported");
if constexpr (cdim == 1)
return CGALPointType(p[0]);
else if constexpr (cdim == 2)
return CGALPointType(p[0], dim > 1 ? p[1] : 0);
else if constexpr (cdim == 3)
return CGALPointType(p[0], dim > 1 ? p[1] : 0, dim > 2 ? p[2] : 0);
else
Assert(false, dealii::ExcNotImplemented());
return CGALPointType();
}



template <int dim, typename CGALPointType>
inline dealii::Point<dim>
cgal_point_to_dealii_point(const CGALPointType &p)
{
constexpr int cdim = CGALPointType::Ambient_dimension::value;
if constexpr (dim == 1)
return dealii::Point<dim>(CGAL::to_double(p.x()));
else if constexpr (dim == 2)
return dealii::Point<dim>(CGAL::to_double(p.x()),
cdim > 1 ? CGAL::to_double(p.y()) : 0);
else if constexpr (dim == 3)
return dealii::Point<dim>(CGAL::to_double(p.x()),
cdim > 1 ? CGAL::to_double(p.y()) : 0,
cdim > 2 ? CGAL::to_double(p.z()) : 0);
else
Assert(false, dealii::ExcNotImplemented());
}
} // namespace CGALWrappers

# endif

DEAL_II_NAMESPACE_CLOSE

#endif
#endif
4 changes: 1 addition & 3 deletions include/deal.II/cgal/surface_mesh.h
Original file line number Diff line number Diff line change
Expand Up @@ -23,12 +23,10 @@

#include <deal.II/grid/tria.h>

#include <deal.II/cgal/utilities.h>

#ifdef DEAL_II_WITH_CGAL
# include <CGAL/Polygon_mesh_processing/stitch_borders.h>
# include <CGAL/Surface_mesh.h>

# include <deal.II/cgal/point_conversion.h>

DEAL_II_NAMESPACE_OPEN

Expand Down
4 changes: 2 additions & 2 deletions include/deal.II/cgal/triangulation.h
Original file line number Diff line number Diff line change
Expand Up @@ -257,8 +257,8 @@ namespace CGALWrappers
// A face is identified by a cell and by the index within the cell
// of the opposite vertex. Loop over vertices, and retain only those
// that belong to this face.
unsigned int j = 0;
for (unsigned int i = 0; i < 4; ++i)
int j = 0;
for (int i = 0; i < 4; ++i)
if (i != cgal_vertex_face_index)
dealii_face.vertices[j++] =
cgal_to_dealii_vertex_map[cgal_cell->vertex(i)];
Expand Down
106 changes: 43 additions & 63 deletions include/deal.II/cgal/utilities.h
Original file line number Diff line number Diff line change
Expand Up @@ -209,30 +209,6 @@ namespace CGALWrappers



/**
* Convert from a deal.II Point to any compatible CGAL point.
*
* @tparam CGALPointType Any of the CGAL point types
* @tparam dim Dimension of the point
* @param [in] p An input deal.II Point<dim>
* @return CGALPointType A CGAL point
*/
template <typename CGALPointType, int dim>
inline CGALPointType
dealii_point_to_cgal_point(const dealii::Point<dim> &p);

/**
* Convert from various CGAL point types to deal.II Point.
*
* @tparam dim Dimension of the point
* @tparam CGALPointType Any of the CGAL point types
* @param p An input CGAL point type
* @return dealii::Point<dim> The corresponding deal.II point.
*/
template <int dim, typename CGALPointType>
inline dealii::Point<dim>
cgal_point_to_dealii_point(const CGALPointType &p);

/**
* Given a closed CGAL::Surface_mesh, this function fills the
* region bounded by the surface with tets.
Expand Down Expand Up @@ -375,45 +351,6 @@ namespace CGALWrappers
// Template implementations
namespace CGALWrappers
{
template <typename CGALPointType, int dim>
inline CGALPointType
dealii_point_to_cgal_point(const dealii::Point<dim> &p)
{
constexpr int cdim = CGALPointType::Ambient_dimension::value;
static_assert(dim <= cdim, "Only dim <= cdim supported");
if constexpr (cdim == 1)
return CGALPointType(p[0]);
else if constexpr (cdim == 2)
return CGALPointType(p[0], dim > 1 ? p[1] : 0);
else if constexpr (cdim == 3)
return CGALPointType(p[0], dim > 1 ? p[1] : 0, dim > 2 ? p[2] : 0);
else
Assert(false, dealii::ExcNotImplemented());
return CGALPointType();
}



template <int dim, typename CGALPointType>
inline dealii::Point<dim>
cgal_point_to_dealii_point(const CGALPointType &p)
{
constexpr int cdim = CGALPointType::Ambient_dimension::value;
if constexpr (dim == 1)
return dealii::Point<dim>(CGAL::to_double(p.x()));
else if constexpr (dim == 2)
return dealii::Point<dim>(CGAL::to_double(p.x()),
cdim > 1 ? CGAL::to_double(p.y()) : 0);
else if constexpr (dim == 3)
return dealii::Point<dim>(CGAL::to_double(p.x()),
cdim > 1 ? CGAL::to_double(p.y()) : 0,
cdim > 2 ? CGAL::to_double(p.z()) : 0);
else
Assert(false, dealii::ExcNotImplemented());
}



template <typename C3t3>
void
cgal_surface_mesh_to_cgal_triangulation(
Expand Down Expand Up @@ -686,5 +623,48 @@ namespace CGALWrappers

DEAL_II_NAMESPACE_CLOSE

#else

DEAL_II_NAMESPACE_OPEN
namespace CGALWrappers
{
template <int dim>
struct AdditionalData
{};

template <int dim>
void
implicit_function(Triangulation<dim, 3> & tria,
const Function<3> & implicit_function,
const AdditionalData<dim> &data = AdditionalData<dim>{},
const Point<3> & interior_point = Point<3>(),
const double & outer_ball_radius = 1.0)
{
(void)tria;
(void)implicit_function;
(void)data;
(void)interior_point;
(void)outer_ball_radius;
AssertThrow(false, ExcMessage("This function needs CGAL to be installed."));
}



void
surface_mesh_to_volumetric_mesh(
const Triangulation<2, 3> &surface_tria,
Triangulation<3> & vol_tria,
const AdditionalData<3> & data = AdditionalData<3>{})
{
(void)surface_tria;
(void)vol_tria;
(void)data;
AssertThrow(false, ExcMessage("This function needs CGAL to be installed."));
}

} // namespace CGALWrappers

DEAL_II_NAMESPACE_CLOSE

#endif
#endif
90 changes: 10 additions & 80 deletions include/deal.II/grid/grid_generator.h
Original file line number Diff line number Diff line change
Expand Up @@ -30,25 +30,18 @@

#ifdef DEAL_II_WITH_CGAL
// All functions needed by the CGAL mesh generation utilities
<<<<<<< HEAD
# include <CGAL/Complex_2_in_triangulation_3.h>
# include <CGAL/IO/facets_in_complex_2_to_triangle_mesh.h>
# include <CGAL/Implicit_surface_3.h>
=======
>>>>>>> 4b1e1f45a5 (Surface deal.II Tria to volumetric deal.II Tria)
# include <CGAL/Labeled_mesh_domain_3.h>
# include <CGAL/Mesh_complex_3_in_triangulation_3.h>
# include <CGAL/Mesh_criteria_3.h>
# include <CGAL/Mesh_triangulation_3.h>
<<<<<<< HEAD
# include <CGAL/Surface_mesh.h>
# include <CGAL/Surface_mesh_default_triangulation_3.h>
# include <CGAL/make_mesh_3.h>
# include <CGAL/make_surface_mesh.h>
=======
# include <CGAL/make_mesh_3.h>
# include <deal.II/cgal/surface_mesh.h>
>>>>>>> 4b1e1f45a5 (Surface deal.II Tria to volumetric deal.II Tria)
# include <deal.II/cgal/triangulation.h>
# include <deal.II/cgal/utilities.h>
#endif
Expand Down Expand Up @@ -1831,38 +1824,26 @@ namespace GridGenerator
CGALWrappers::AdditionalData<dim>{},
const Point<3> &interior_point = Point<3>(),
const double & outer_ball_radius = 1.0);
/**

/**
* Create a deal.II Triangulation<3> out of a deal.II Triangulation<2,3>
* by filling it with tetrahedrons.
* by filling it with tetrahedra.
*
* The optional variable arguments @p cgal_args can be used to pass additional
* The last optional argument @p data can be used to pass additional
* arguments to the CGAL::Mesh_criteria_3 class (see
* https://doc.cgal.org/latest/Mesh_3/index.html) for more information.
*
* The arguments allow fine control on the size, quality, and distribution of
* the cells of the final triangulation. CGAL uses named parameters for these
* arguments, i.e., they must be specified with the syntax
* `CGAL::parameters::parameter_name=parameter_value`. Accepted parameters
* are:
*
* - CGAL::parameters::facet_angle,
* - CGAL::parameters::facet_size,
* - CGAL::parameters::facet_distance,
* - CGAL::parameters::cell_radius_edge_ratio,
* - CGAL::parameters::cell_size;
* https://doc.cgal.org/latest/Mesh_3/index.html for more information).
*
*
* @param [in] surface_tria The input deal.II Triangulation<2,3>.
* @param [out] vol_tria The output deal.II Triangulation<3>.
* @param [in] cgal_args Additional parameters to pass to the CGAL::make_mesh_3
* function.
* @param[in] data Additional parameters to pass to the CGAL::make_mesh_3
* function and to the CGAL::make_surface_mesh functions.
*/
template <typename... Args>
void
surface_mesh_to_volumetric_mesh(const Triangulation<2, 3> &surface_tria,
Triangulation<3> &vol_tria,
Args... cgal_args);
Triangulation<3> & vol_tria,
const CGALWrappers::AdditionalData<3> &data =
CGALWrappers::AdditionalData<3>{});
#endif
///@}

Expand Down Expand Up @@ -2876,57 +2857,6 @@ namespace GridGenerator
Assert(false, ExcImpossibleInDim(dim));
}
}



template <typename... Args>
void
surface_mesh_to_volumetric_mesh(const Triangulation<2, 3> &surface_tria,
Triangulation<3> &vol_tria,
Args... cgal_args)
{
Assert(
surface_tria.n_cells() > 0,
ExcMessage(
"The input triangulation cannot be empty when calling this function."));
Assert(
vol_tria.n_cells() == 0,
ExcMessage(
"The output triangulation must be empty when calling this function."));
using K = CGAL::Exact_predicates_inexact_constructions_kernel;
using Point_3 = K::Point_3;

using Mesh_domain =
CGAL::Polyhedral_mesh_domain_with_features_3<K,
CGAL::Surface_mesh<Point_3>>;
using Tr = CGAL::Mesh_triangulation_3<Mesh_domain,
CGAL::Default,
CGALWrappers::ConcurrencyTag>::type;
using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
using C3t3 =
CGAL::Mesh_complex_3_in_triangulation_3<Tr,
Mesh_domain::Corner_index,
Mesh_domain::Curve_index>;

CGAL::Surface_mesh<Point_3> mesh;
// This function "fills" the missing arrow of the following diagram.
// Tria<2,3> Tria<3>
// | ^
// | |
// | |
// | |
// V |
// CGAL::Surface_mesh -----------> CGAL::C3t3
CGALWrappers::dealii_tria_to_cgal_surface_mesh(surface_tria, mesh);
CGAL::Polygon_mesh_processing::triangulate_faces(mesh);
CGAL::Polygon_mesh_processing::stitch_borders(mesh);
Mesh_domain domain(mesh);
domain.detect_features();
Mesh_criteria criteria(cgal_args...);
const auto cgal_triangulation = CGAL::make_mesh_3<C3t3>(domain, criteria);
CGALWrappers::cgal_triangulation_to_dealii_triangulation(cgal_triangulation,
vol_tria);
}
# endif

#endif
Expand Down

0 comments on commit be72338

Please sign in to comment.