deal.II version GIT relicensing-2645-g9d42f38b55 2025-02-15 04:20:00+00:00
\(\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// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2020 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_grid_tria_description_h
16#define dealii_grid_tria_description_h
17
18#include <deal.II/base/config.h>
19
22
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
43
59{
77
94 template <int dim>
95 struct CellData
96 {
100 CellData();
101
107 template <class Archive>
108 void
109 serialize(Archive &ar, const unsigned int /*version*/);
110
114 bool
115 operator==(const CellData<dim> &other) const;
116
121
126
131
136
142 std::array<types::manifold_id, GeometryInfo<dim>::lines_per_cell>
144
150 std::array<types::manifold_id,
151 dim == 1 ? 1 : GeometryInfo<3>::quads_per_cell>
153
158 std::vector<std::pair<unsigned int, types::boundary_id>> boundary_ids;
159 };
160
161
168 template <int dim, int spacedim = dim>
170 {
174 Description();
175
181 template <class Archive>
182 void
183 serialize(Archive &ar, const unsigned int /*version*/);
184
188 bool
189 operator==(const Description<dim, spacedim> &other) const;
190
194 std::vector<::CellData<dim>> coarse_cells;
195
199 std::vector<Point<spacedim>> coarse_cell_vertices;
200
206 std::vector<types::coarse_cell_id> coarse_cell_index_to_coarse_cell_id;
207
213 std::vector<std::vector<CellData<dim>>> cell_infos;
214
225
230
235 };
236
237
241 namespace Utilities
242 {
329 template <int dim, int spacedim = dim>
332 const ::Triangulation<dim, spacedim> &tria,
333 const MPI_Comm comm,
336 const unsigned int my_rank_in = numbers::invalid_unsigned_int);
337
367 template <int dim, int spacedim>
372 &partition,
375
383 template <int dim, int spacedim>
388 &partition,
389 const std::vector<
391 &mg_partitions,
394
452 template <int dim, int spacedim = dim>
455 const std::function<void(::Triangulation<dim, spacedim> &)>
456 &serial_grid_generator,
457 const std::function<void(::Triangulation<dim, spacedim> &,
458 const MPI_Comm,
459 const unsigned int)> &serial_grid_partitioner,
460 const MPI_Comm comm,
461 const int group_size = 1,
462 const typename Triangulation<dim, spacedim>::MeshSmoothing smoothing =
466
467 } // namespace Utilities
468
469
470
471 template <int dim>
473 : subdomain_id(numbers::invalid_subdomain_id)
474 , level_subdomain_id(numbers::invalid_subdomain_id)
475 , manifold_id(numbers::flat_manifold_id)
476 {
477 std::fill(id.begin(), id.end(), numbers::invalid_unsigned_int);
478 std::fill(manifold_line_ids.begin(),
479 manifold_line_ids.end(),
481 std::fill(manifold_quad_ids.begin(),
482 manifold_quad_ids.end(),
484 }
485
486
487
488 template <int dim>
489 template <class Archive>
490 void
491 CellData<dim>::serialize(Archive &ar, const unsigned int /*version*/)
492 {
493 ar &id;
494 ar &subdomain_id;
495 ar &level_subdomain_id;
496 ar &manifold_id;
497 if (dim >= 2)
498 ar &manifold_line_ids;
499 if (dim >= 3)
500 ar &manifold_quad_ids;
501 ar &boundary_ids;
502 }
503
504
505
506 template <int dim, int spacedim>
508 : comm(MPI_COMM_NULL)
509 , settings(Settings::default_setting)
510 , smoothing(Triangulation<dim, spacedim>::MeshSmoothing::none)
511 {}
512
513
514
515 template <int dim, int spacedim>
516 template <class Archive>
517 void
519 const unsigned int /*version*/)
520 {
521 ar &coarse_cells;
524 ar &cell_infos;
525 ar &settings;
526 ar &smoothing;
527 }
528
529
530
531 template <int dim>
532 bool
534 {
535 if (this->id != other.id)
536 return false;
537 if (this->subdomain_id != other.subdomain_id)
538 return false;
539 if (this->level_subdomain_id != other.level_subdomain_id)
540 return false;
541 if (this->manifold_id != other.manifold_id)
542 return false;
543 if (dim >= 2 && this->manifold_line_ids != other.manifold_line_ids)
544 return false;
545 if (dim >= 3 && this->manifold_quad_ids != other.manifold_quad_ids)
546 return false;
547 if (this->boundary_ids != other.boundary_ids)
548 return false;
549
550 return true;
551 }
552
553
554
555 template <int dim, int spacedim>
556 bool
558 const Description<dim, spacedim> &other) const
559 {
560 if (this->coarse_cells != other.coarse_cells)
561 return false;
563 return false;
566 return false;
567 if (this->cell_infos != other.cell_infos)
568 return false;
569 if (this->settings != other.settings)
570 return false;
571 if (this->smoothing != other.smoothing)
572 return false;
573
574 return true;
575 }
576} // namespace TriangulationDescription
577
578
580
581#endif
std::array< unsigned int, 4 > binary_type
Definition cell_id.h:80
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:518
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:519
const MPI_Comm comm
Definition mpi.cc:918
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)
constexpr unsigned int invalid_unsigned_int
Definition types.h:232
constexpr types::manifold_id flat_manifold_id
Definition types.h:336
unsigned int manifold_id
Definition types.h:165
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
std::vector< types::coarse_cell_id > coarse_cell_index_to_coarse_cell_id
std::vector< std::vector< CellData< dim > > > cell_infos
std::vector<::CellData< dim > > coarse_cells