15#ifndef dealii_grid_tools_h
16#define dealii_grid_tools_h
48#include <boost/archive/binary_iarchive.hpp>
49#include <boost/archive/binary_oarchive.hpp>
50#include <boost/random/mersenne_twister.hpp>
51#include <boost/serialization/array.hpp>
52#include <boost/serialization/vector.hpp>
54#ifdef DEAL_II_WITH_ZLIB
55# include <boost/iostreams/device/back_inserter.hpp>
56# include <boost/iostreams/filter/gzip.hpp>
57# include <boost/iostreams/filtering_stream.hpp>
58# include <boost/iostreams/stream.hpp>
64#ifdef DEAL_II_HAVE_CXX20
77 template <
int dim,
int spacedim>
86 class MappingCollection;
93 template <
int dim,
int spacedim>
100 template <
int dim,
int spacedim,
typename MeshType>
106 using type =
typename MeshType::active_cell_iterator;
113 template <
int dim,
int spacedim>
219 template <
int dim,
typename Transformation,
int spacedim>
222 std::assignable_from<
545 std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
546 std::vector<std::vector<Point<dim>>>,
547 std::vector<std::vector<unsigned int>>>
591 template <
int dim,
int spacedim>
594 std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
595 std::vector<std::vector<Point<dim>>>,
596 std::vector<std::vector<unsigned int>>,
597 std::vector<unsigned int>>
687 template <
int dim,
int spacedim>
690 std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>,
691 std::vector<std::vector<Point<dim>>>,
692 std::vector<std::vector<unsigned int>>,
693 std::vector<std::vector<Point<spacedim>>>,
694 std::vector<std::vector<unsigned int>>>
702 const double tolerance = 1e-10,
703 const std::vector<bool> &marked_vertices = {},
704 const bool enforce_unique_mapping =
true);
722 template <
int dim,
int spacedim>
747 std::vector<std::tuple<std::pair<int, int>,
777 std::vector<std::tuple<unsigned int, unsigned int, unsigned int>>
802 template <
int dim,
int spacedim>
808 const std::vector<bool> &marked_vertices,
809 const double tolerance,
811 const bool enforce_unique_mapping =
false);
821 template <
int structdim,
int spacedim>
828 std::array<::Point<spacedim>, structdim + 1>;
838 std::vector<std::tuple<std::pair<int, int>,
852 std::vector<std::tuple<unsigned int, unsigned int, IntersectionType>>
899 std::map<unsigned int, std::vector<unsigned int>>
901 const std::vector<std::tuple<unsigned int, unsigned int, unsigned int>>
913 template <
int structdim,
int dim,
int spacedim>
919 const std::vector<bool> &marked_vertices,
920 const double tolerance);
933 template <
int spacedim>
964 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
970 const std::vector<
bool> &marked_vertices = {});
998 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1005 const std::vector<
bool> &marked_vertices = {});
1030 template <
typename MeshType>
1033 std::vector<typename MeshType::active_cell_iterator>
1037 MeshType::space_dimension,
1041 const unsigned int vertex_index);
1108 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1112 std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
Point<dim>>
1114 std::pair<typename ::internal::
1115 ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
1121 const std::vector<bool> &marked_vertices = {},
1122 const double tolerance = 1.e-10);
1134 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1140 typename ::internal::
1141 ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type
1145 const std::vector<bool> &marked_vertices = {},
1146 const double tolerance = 1.e-10);
1154 template <
int dim,
int spacedim>
1155 std::pair<typename DoFHandler<dim, spacedim>::active_cell_iterator,
1161 const double tolerance = 1.e-10);
1214 template <
int dim,
int spacedim>
1215 std::pair<typename Triangulation<dim, spacedim>::active_cell_iterator,
1222 const std::vector<bool> &marked_vertices = {},
1223 const double tolerance = 1.e-10);
1241 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1245 std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
Point<dim>>
1247 std::pair<typename ::internal::
1248 ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
1259 &vertex_to_cell_centers,
1262 const std::vector<bool> &marked_vertices = {},
1265 const double tolerance = 1.e-10,
1267 std::pair<BoundingBox<spacedim>,
1300 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1304 std::vector<std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
1307 std::vector<std::pair<
1308 typename ::internal::
1309 ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
1316 const double tolerance,
1321 *vertex_to_cells =
nullptr);
1332 template <
int dim,
template <
int,
int>
class MeshType,
int spacedim>
1336 std::vector<std::pair<typename MeshType<dim, spacedim>::active_cell_iterator,
1339 std::vector<std::pair<
1340 typename ::internal::
1341 ActiveCellIterator<dim, spacedim, MeshType<dim, spacedim>>::type,
1348 const double tolerance = 1e-10,
1349 const std::vector<bool> &marked_vertices = {});
1374 template <
typename MeshType>
1712 std::tuple<std::vector<std::vector<unsigned int>>,
1713 std::map<unsigned int, unsigned int>,
1714 std::map<unsigned int, std::vector<unsigned int>>>
1757 template <
int spacedim>
1759 std::tuple<std::map<unsigned int, std::vector<unsigned int>>,
1760 std::map<unsigned int, unsigned int>,
1761 std::map<unsigned int, std::vector<unsigned int>>>
1781 template <
int dim,
int spacedim>
1782 std::vector<std::vector<Tensor<1, spacedim>>>
1797 template <
int dim,
int spacedim>
1803 (ReferenceCells::get_hypercube<dim>()
1807 .ReferenceCell::get_default_linear_mapping<dim, spacedim>()
1822 template <
int dim,
int spacedim>
1823 std::map<unsigned int, types::global_vertex_index>
1853 template <
int dim,
int spacedim>
1870 template <
int dim,
int spacedim>
1923 template <
int dim,
int spacedim>
1941 template <
int dim,
int spacedim>
1964 template <
int dim,
int spacedim>
1981 template <
int dim,
int spacedim>
1992 template <
int dim,
int spacedim>
1993 std::vector<types::subdomain_id>
1995 const std::vector<CellId> &cell_ids);
2007 template <
int dim,
int spacedim>
2010 std::vector<types::subdomain_id> &
subdomain);
2026 template <
int dim,
int spacedim>
2061 template <
int dim,
int spacedim>
2091 template <
int dim,
int spacedim>
2149 template <
typename MeshType>
2152 const typename MeshType::active_cell_iterator &cell);
2176 template <
class Container>
2177 std::vector<typename Container::cell_iterator>
2179 const std::vector<typename Container::active_cell_iterator> &
patch_cells);
2247 template <
class Container>
2250 const std::vector<typename Container::active_cell_iterator> &
patch,
2255 Container::space_dimension>::active_cell_iterator,
2289 template <
int dim,
int spacedim>
2292 std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>>
2310 template <
typename CellIterator>
2374 template <
typename FaceIterator>
2375 std::optional<types::geometric_orientation>
2379 const unsigned int direction,
2442 template <
typename MeshType>
2448 const unsigned int direction,
2476 template <
typename MeshType>
2481 const unsigned int direction,
2577 template <
int dim,
int spacedim>
2606 template <
int dim,
int spacedim>
2612 [](
const std::set<types::manifold_id> &manifold_ids) {
2613 if (manifold_ids.size() == 1)
2614 return *manifold_ids.begin();
2707 template <
typename DataType,
typename MeshType>
2714 const DataType &)> &unpack,
2731 template <
typename DataType,
typename MeshType>
2738 const DataType &)> &unpack,
2754 template <
int spacedim>
2755 std::vector<std::vector<BoundingBox<spacedim>>>
2792 template <
int spacedim>
2815 template <
int dim,
int spacedim>
2841 template <
int dim,
int spacedim>
2842 std::map<unsigned int, std::set<::types::subdomain_id>>
2866 template <
int dim,
typename VectorType>
2961 const std::vector<unsigned int> &,
2978 const std::vector<
Point<2>> &points,
2979 const std::vector<unsigned int> &mask,
2990 const std::vector<
Point<3>> &points,
2991 const std::vector<unsigned int> &mask,
3028 <<
"The number of partitions you gave is " <<
arg1
3029 <<
", but must be greater than zero.");
3035 <<
"The subdomain id " <<
arg1
3036 <<
" has no cells associated with it.");
3047 <<
"The scaling factor must be positive, but it is " <<
arg1
3055 <<
"The given vertex with index " <<
arg1
3056 <<
" is not used in the given triangulation.");
3095 template <
int dim,
typename Transformation,
int spacedim>
3098 std::assignable_from<
3115 for (; cell !=
endc; ++cell)
3132 for (; cell !=
endc; ++cell)
3133 for (
const unsigned int face : cell->face_indices())
3134 if (cell->face(face)->has_children() &&
3135 !cell->face(face)->at_boundary())
3137 Assert(cell->reference_cell() ==
3138 ReferenceCells::get_hypercube<dim>(),
3142 cell->face(face)->child(0)->vertex(1) =
3143 (cell->face(face)->vertex(0) + cell->face(face)->vertex(1)) /
3152 for (; cell !=
endc; ++cell)
3153 for (
const unsigned int face : cell->face_indices())
3154 if (cell->face(face)->has_children() &&
3155 !cell->face(face)->at_boundary())
3157 if (
static_cast<std::uint8_t
>(
3158 cell->face(face)->refinement_case()) ==
3161 Assert(cell->reference_cell() ==
3162 ReferenceCells::get_hypercube<dim>(),
3166 cell->face(face)->child(0)->vertex(1) =
3167 (cell->face(face)->vertex(0) +
3168 cell->face(face)->vertex(1)) /
3170 cell->face(face)->child(0)->vertex(2) =
3171 (cell->face(face)->vertex(0) +
3172 cell->face(face)->vertex(2)) /
3174 cell->face(face)->child(1)->vertex(3) =
3175 (cell->face(face)->vertex(1) +
3176 cell->face(face)->vertex(3)) /
3178 cell->face(face)->child(2)->vertex(3) =
3179 (cell->face(face)->vertex(2) +
3180 cell->face(face)->vertex(3)) /
3184 cell->face(face)->child(0)->vertex(3) =
3185 (cell->face(face)->vertex(0) +
3186 cell->face(face)->vertex(1) +
3187 cell->face(face)->vertex(2) +
3188 cell->face(face)->vertex(3)) /
3194 for (
unsigned int line = 0;
3197 if (cell->face(face)->line(line)->has_children())
3198 cell->face(face)->line(line)->child(0)->vertex(1) =
3199 (cell->face(face)->line(line)->vertex(0) +
3200 cell->face(face)->line(line)->vertex(1)) /
3212 template <
typename MeshType>
3215 const typename MeshType::cell_iterator &cell)
3217 std::vector<typename MeshType::active_cell_iterator>
child_cells;
3219 if (cell->has_children())
3221 for (
unsigned int child = 0; child < cell->n_children(); ++child)
3222 if (cell->child(child)->has_children())
3224 const std::vector<typename MeshType::active_cell_iterator>
3239 template <
typename MeshType>
3242 const typename MeshType::active_cell_iterator &cell,
3246 for (
const unsigned int n : cell->face_indices())
3247 if (!cell->at_boundary(n))
3270 if (cell->face(n)->has_children())
3274 for (
unsigned int c = 0;
3275 c < cell->face(n)->n_active_descendants();
3278 cell->neighbor_child_on_subface(n, c));
3292 template <
typename CellIterator>
3296 return sizeof(*this) +
matrix.memory_consumption();
3303 template <
typename DataType,
3314 const std::function<
void(
3317 const std::function<std::set<types::subdomain_id>(
3319 MeshType::space_dimension> &)>
3322# ifndef DEAL_II_WITH_MPI
3332 constexpr int spacedim = MeshType::space_dimension;
3335 &
mesh.get_triangulation());
3339 "The function exchange_cell_data_to_ghosts() only works with parallel triangulations."));
3341 if (
const auto tria =
dynamic_cast<
3343 &
mesh.get_triangulation()))
3346 tria->with_artificial_cells(),
3348 "The functions GridTools::exchange_cell_data_to_ghosts() and "
3349 "GridTools::exchange_cell_data_to_level_ghosts() can only "
3350 "operate on a single layer of ghost cells. However, you have "
3351 "given a Triangulation object of type "
3352 "parallel::shared::Triangulation without artificial cells "
3353 "resulting in an arbitrary number of ghost layers. "
3354 "To use this function for a Triangulation object of type "
3355 "parallel::shared::Triangulation, make sure to create the "
3356 "Triangulation object with allow_artificial_cells set to true. "
3357 "This results in a parallel::shared::Triangulation with only "
3358 "a single layer of ghost cells."));
3362 std::set<::types::subdomain_id> ghost_owners =
3365 std::vector<typename CellId::binary_type>>
3374 constexpr int spacedim = MeshType::space_dimension;
3388 mutex, tria->get_mpi_communicator());
3396 std::vector<MPI_Request> requests(ghost_owners.size());
3398 unsigned int idx = 0;
3403 it.second.size() *
sizeof(
it.second[0]),
3407 tria->get_mpi_communicator(),
3416 std::vector<std::vector<char>>
sendbuffers(ghost_owners.size());
3418 for (
unsigned int idx = 0; idx < ghost_owners.size(); ++idx)
3423 tria->get_mpi_communicator(),
3446 tria->get_mpi_communicator(),
3484 sizeof(std::size_t));
3499 tria->get_mpi_communicator(),
3511 tria->get_mpi_communicator(),
3526 tria->get_mpi_communicator(),
3539 auto received_data = Utilities::unpack<std::vector<DataType>>(
3548 const std::vector<typename CellId::binary_type> &
this_cell_list =
3580 if (requests.size() > 0)
3600 template <
typename DataType,
typename MeshType>
3604 const std::function<std::optional<DataType>(
3605 const typename MeshType::active_cell_iterator &)> &pack,
3606 const std::function<
void(
const typename MeshType::active_cell_iterator &,
3607 const DataType &)> &unpack,
3608 const std::function<
bool(
const typename MeshType::active_cell_iterator &)>
3611# ifndef DEAL_II_WITH_MPI
3618 internal::exchange_cell_data<DataType,
3620 typename MeshType::active_cell_iterator>(
3625 [&](
const auto &process) {
3626 for (
const auto &cell :
mesh.active_cell_iterators())
3627 if (cell->is_ghost())
3630 [](
const auto &tria) {
return tria.ghost_owners(); });
3636 template <
typename DataType,
typename MeshType>
3640 const std::function<std::optional<DataType>(
3641 const typename MeshType::level_cell_iterator &)> &pack,
3642 const std::function<
void(
const typename MeshType::level_cell_iterator &,
3643 const DataType &)> &unpack,
3644 const std::function<
bool(
const typename MeshType::level_cell_iterator &)>
3647# ifndef DEAL_II_WITH_MPI
3654 internal::exchange_cell_data<DataType,
3656 typename MeshType::level_cell_iterator>(
3661 [&](
const auto &process) {
3662 for (
const auto &cell :
mesh.cell_iterators())
3663 if (cell->is_ghost_on_level())
3664 process(cell, cell->level_subdomain_id());
3666 [](
const auto &tria) {
return tria.level_ghost_owners(); });
std::array< unsigned int, 4 > binary_type
Abstract base class for mapping classes.
static constexpr unsigned int dimension
typename MeshType::active_cell_iterator type
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_CXX20_REQUIRES(condition)
#define DEAL_II_NAMESPACE_CLOSE
#define DeclException0(Exception0)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcNeedsMPI()
static ::ExceptionBase & ExcNonExistentSubdomain(int arg1)
#define Assert(cond, exc)
static ::ExceptionBase & ExcVertexNotUsed(unsigned int arg1)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcTriangulationHasBeenRefined()
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcInvalidNumberOfPartitions(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcScalingFactorNotPositive(double arg1)
#define AssertThrow(cond, exc)
typename ActiveSelector::active_cell_iterator active_cell_iterator
void copy_boundary_to_manifold_id(Triangulation< dim, spacedim > &tria, const bool reset_boundary_ids=false)
void copy_material_to_manifold_id(Triangulation< dim, spacedim > &tria, const bool compute_face_ids=false)
void map_boundary_to_manifold_ids(const std::vector< types::boundary_id > &src_boundary_ids, const std::vector< types::manifold_id > &dst_manifold_ids, Triangulation< dim, spacedim > &tria, const std::vector< types::boundary_id > &reset_boundary_ids={})
void assign_co_dimensional_manifold_indicators(Triangulation< dim, spacedim > &tria, const std::function< types::manifold_id(const std::set< types::manifold_id > &)> &disambiguation_function=[](const std::set< types::manifold_id > &manifold_ids) { if(manifold_ids.size()==1) return *manifold_ids.begin();else return numbers::flat_manifold_id;}, bool overwrite_only_flat_manifold_ids=true)
std::vector< index_type > data
@ matrix
Contents is actually a matrix.
std::size_t pack(const T &object, std::vector< char > &dest_buffer, const bool allow_compression=true)
T unpack(const std::vector< char > &buffer, const bool allow_compression=true)
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
constexpr types::manifold_id flat_manifold_id
unsigned int global_dof_index
unsigned int subdomain_id
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
boost::geometry::index::rtree< LeafType, IndexType, IndexableGetter > RTree