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
tria_description.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2019 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_grid_construction_utilities_h
17#define dealii_grid_construction_utilities_h
18
19#include <deal.II/base/config.h>
20
23
26#include <deal.II/grid/tria.h>
27
28
30#ifndef DOXYGEN
31namespace LinearAlgebra
32{
33 namespace distributed
34 {
35 template <typename Number, typename MemorySpace>
36 class Vector;
37 }
38} // namespace LinearAlgebra
39#endif
40
41/*------------------------------------------------------------------------*/
42
78template <int structdim>
80{
100 std::vector<unsigned int> vertices;
101
109 union
110 {
120
131 };
132
142
160 CellData(const unsigned int n_vertices =
161 ReferenceCells::get_hypercube<structdim>().n_vertices());
162
166 bool
167 operator==(const CellData<structdim> &other) const;
168
174 template <class Archive>
175 void
176 serialize(Archive &ar, const unsigned int version);
177
178 static_assert(structdim > 0,
179 "The class CellData can only be used for structdim>0.");
180};
181
182
183
239{
247 std::vector<CellData<1>> boundary_lines;
248
263 std::vector<CellData<2>> boundary_quads;
264
273 bool
274 check_consistency(const unsigned int dim) const;
275};
276
277
278template <int structdim>
279template <class Archive>
280void
281CellData<structdim>::serialize(Archive &ar, const unsigned int /*version*/)
282{
283 ar &vertices;
284 ar &material_id;
285 ar &boundary_id;
286 ar &manifold_id;
287}
288
294{
300 {
311 };
312
325 template <int dim>
326 struct CellData
327 {
333 template <class Archive>
334 void
335 serialize(Archive &ar, const unsigned int /*version*/);
336
340 bool
341 operator==(const CellData<dim> &other) const;
342
347
352
357
362
368 std::array<types::manifold_id, GeometryInfo<dim>::lines_per_cell>
370
376 std::array<types::manifold_id,
377 dim == 1 ? 1 : GeometryInfo<3>::quads_per_cell>
379
384 std::vector<std::pair<unsigned int, types::boundary_id>> boundary_ids;
385 };
386
390 template <int dim, int spacedim>
392 {
398 template <class Archive>
399 void
400 serialize(Archive &ar, const unsigned int /*version*/);
401
405 bool
406 operator==(const Description<dim, spacedim> &other) const;
407
411 std::vector<::CellData<dim>> coarse_cells;
412
416 std::vector<Point<spacedim>> coarse_cell_vertices;
417
423 std::vector<types::coarse_cell_id> coarse_cell_index_to_coarse_cell_id;
424
430 std::vector<std::vector<CellData<dim>>> cell_infos;
431
442
447
452 };
453
454
460 namespace Utilities
461 {
527 template <int dim, int spacedim = dim>
530 const ::Triangulation<dim, spacedim> &tria,
531 const MPI_Comm comm,
534 const unsigned int my_rank_in = numbers::invalid_unsigned_int);
535
556 template <int dim, int spacedim>
561 & partition,
564
569 template <int dim, int spacedim>
574 &partition,
575 const std::vector<
577 & mg_partitions,
580
638 template <int dim, int spacedim = dim>
641 const std::function<void(::Triangulation<dim, spacedim> &)>
642 & serial_grid_generator,
643 const std::function<void(::Triangulation<dim, spacedim> &,
644 const MPI_Comm,
645 const unsigned int)> &serial_grid_partitioner,
646 const MPI_Comm comm,
647 const int group_size = 1,
648 const typename Triangulation<dim, spacedim>::MeshSmoothing smoothing =
652
653 } // namespace Utilities
654
655
656
657 template <int dim>
658 template <class Archive>
659 void
660 CellData<dim>::serialize(Archive &ar, const unsigned int /*version*/)
661 {
662 ar &id;
663 ar &subdomain_id;
664 ar &level_subdomain_id;
665 ar &manifold_id;
666 if (dim >= 2)
667 ar &manifold_line_ids;
668 if (dim >= 3)
669 ar &manifold_quad_ids;
670 ar &boundary_ids;
671 }
672
673
674 template <int dim, int spacedim>
675 template <class Archive>
676 void
678 const unsigned int /*version*/)
679 {
680 ar &coarse_cells;
683 ar &cell_infos;
684 ar &settings;
685 ar &smoothing;
686 }
687
688
689
690 template <int dim>
691 bool
693 {
694 if (this->id != other.id)
695 return false;
696 if (this->subdomain_id != other.subdomain_id)
697 return false;
698 if (this->level_subdomain_id != other.level_subdomain_id)
699 return false;
700 if (this->manifold_id != other.manifold_id)
701 return false;
702 if (dim >= 2 && this->manifold_line_ids != other.manifold_line_ids)
703 return false;
704 if (dim >= 3 && this->manifold_quad_ids != other.manifold_quad_ids)
705 return false;
706 if (this->boundary_ids != other.boundary_ids)
707 return false;
708
709 return true;
710 }
711
712
713
714 template <int dim, int spacedim>
715 bool
717 const Description<dim, spacedim> &other) const
718 {
719 if (this->coarse_cells != other.coarse_cells)
720 return false;
722 return false;
725 return false;
726 if (this->cell_infos != other.cell_infos)
727 return false;
728 if (this->settings != other.settings)
729 return false;
730 if (this->smoothing != other.smoothing)
731 return false;
732
733 return true;
734 }
735} // namespace TriangulationDescription
736
737
739
740#endif
std::array< unsigned int, 4 > binary_type
Definition cell_id.h:81
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
Description< dim, spacedim > create_description_from_triangulation(const ::Triangulation< dim, spacedim > &tria, const MPI_Comm comm, const TriangulationDescription::Settings settings=TriangulationDescription::Settings::default_setting, const unsigned int my_rank_in=numbers::invalid_unsigned_int)
Description< dim, spacedim > create_description_from_triangulation_in_groups(const std::function< void(::Triangulation< dim, spacedim > &)> &serial_grid_generator, const std::function< void(::Triangulation< dim, spacedim > &, const MPI_Comm, const unsigned int)> &serial_grid_partitioner, const MPI_Comm comm, const int group_size=1, const typename Triangulation< dim, spacedim >::MeshSmoothing smoothing=::Triangulation< dim, spacedim >::none, const TriangulationDescription::Settings setting=TriangulationDescription::Settings::default_setting)
static const unsigned int invalid_unsigned_int
Definition types.h:213
unsigned int manifold_id
Definition types.h:153
std::vector< unsigned int > vertices
bool operator==(const CellData< structdim > &other) const
types::manifold_id manifold_id
types::material_id material_id
types::boundary_id boundary_id
void serialize(Archive &ar, const unsigned int version)
std::vector< CellData< 2 > > boundary_quads
bool check_consistency(const unsigned int dim) const
std::vector< CellData< 1 > > boundary_lines
std::array< types::manifold_id, dim==1 ? 1 :GeometryInfo< 3 >::quads_per_cell > manifold_quad_ids
std::array< types::manifold_id, GeometryInfo< dim >::lines_per_cell > manifold_line_ids
void serialize(Archive &ar, const unsigned int)
std::vector< std::pair< unsigned int, types::boundary_id > > boundary_ids
bool operator==(const CellData< dim > &other) const
std::vector< std::vector< CellData< dim > > > cell_infos
Triangulation< dim, spacedim >::MeshSmoothing smoothing
void serialize(Archive &ar, const unsigned int)
std::vector<::CellData< dim > > coarse_cells
bool operator==(const Description< dim, spacedim > &other) const
std::vector< Point< spacedim > > coarse_cell_vertices
std::vector< types::coarse_cell_id > coarse_cell_index_to_coarse_cell_id
std::vector< std::pair< unsigned int, Point< spacedim > > > coarse_cell_vertices
const TriangulationDescription::Settings settings
std::vector< types::coarse_cell_id > coarse_cell_index_to_coarse_cell_id
const MPI_Comm comm
std::vector< std::vector< CellData< dim > > > cell_infos
const ::Triangulation< dim, spacedim > & tria
std::vector<::CellData< dim > > coarse_cells