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
intergrid_map.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1999 - 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
18
21
24
25#include <deal.II/fe/fe.h>
26
28#include <deal.II/grid/tria.h>
31
32
34
35
36template <class MeshType>
39 : source_grid(nullptr, typeid(*this).name())
40 , destination_grid(nullptr, typeid(*this).name())
41{}
42
43
44
45template <class MeshType>
47void InterGridMap<MeshType>::make_mapping(const MeshType &source_grid,
48 const MeshType &destination_grid)
49{
50 // first delete all contents
51 clear();
52
53 // next store pointers to grids
54 this->source_grid = &source_grid;
55 this->destination_grid = &destination_grid;
56
57 // then set up the meshes from
58 // scratch and fill them with end-iterators
59 const unsigned int n_levels = source_grid.get_triangulation().n_levels();
60 mapping.resize(n_levels);
61 for (unsigned int level = 0; level < n_levels; ++level)
62 {
63 // first find out about the highest
64 // index used on this level. We could
65 // in principle ask the triangulation
66 // about this, but we would have to
67 // know the underlying data structure
68 // for this and we would like to
69 // avoid such knowledge here
70 unsigned int n_cells = 0;
71 cell_iterator cell = source_grid.begin(level),
72 endc = source_grid.end(level);
73 for (; cell != endc; ++cell)
74 if (static_cast<unsigned int>(cell->index()) > n_cells)
75 n_cells = cell->index();
76
77 // note: n_cells is now the largest
78 // zero-based index, but we need the
79 // number of cells, which is one larger
80 mapping[level].resize(n_cells + 1, destination_grid.end());
81 }
82
83 // now make up the mapping
84 // loop over all cells and set the user
85 // pointers as well as the contents of
86 // the two arrays. note that the function
87 // takes a *reference* to the int and
88 // this may change it
89 cell_iterator src_cell = source_grid.begin(0),
90 dst_cell = destination_grid.begin(0), endc = source_grid.end(0);
91 for (; src_cell != endc; ++src_cell, ++dst_cell)
92 set_mapping(src_cell, dst_cell);
93
94 // little assertion that the two grids
95 // are indeed related:
96 Assert(dst_cell == destination_grid.end(0), ExcIncompatibleGrids());
97}
98
99
100
101template <class MeshType>
103void InterGridMap<MeshType>::set_mapping(const cell_iterator &src_cell,
104 const cell_iterator &dst_cell)
105{
106 // first set the map for this cell
107 mapping[src_cell->level()][src_cell->index()] = dst_cell;
108
109 // if both cells have children, we may
110 // recurse further into the hierarchy
111 if (src_cell->has_children() && dst_cell->has_children())
112 {
113 Assert(src_cell->n_children() ==
116 Assert(dst_cell->n_children() ==
119 Assert(src_cell->refinement_case() == dst_cell->refinement_case(),
121 for (unsigned int c = 0;
122 c < GeometryInfo<MeshType::dimension>::max_children_per_cell;
123 ++c)
124 set_mapping(src_cell->child(c), dst_cell->child(c));
125 }
126 else if (src_cell->has_children() && !dst_cell->has_children())
127 // src grid is more refined here.
128 // set entries for all children
129 // of this cell to the one
130 // dst_cell
131 for (unsigned int c = 0; c < src_cell->n_children(); ++c)
132 set_entries_to_cell(src_cell->child(c), dst_cell);
133 // else (no cell is refined or
134 // dst_cell is refined): no pointers
135 // to be set
136}
137
138
139
140template <class MeshType>
142void InterGridMap<MeshType>::set_entries_to_cell(const cell_iterator &src_cell,
143 const cell_iterator &dst_cell)
144{
145 // first set the map for this cell
146 mapping[src_cell->level()][src_cell->index()] = dst_cell;
147
148 // then do so for the children as well
149 // if there are any
150 if (src_cell->has_children())
151 for (unsigned int c = 0; c < src_cell->n_children(); ++c)
152 set_entries_to_cell(src_cell->child(c), dst_cell);
153}
154
155
156template <class MeshType>
160{
161 Assert(source_cell.state() == IteratorState::valid,
162 ExcInvalidKey(source_cell));
163 Assert(source_cell->level() <= static_cast<int>(mapping.size()),
164 ExcInvalidKey(source_cell));
165 Assert(source_cell->index() <=
166 static_cast<int>(mapping[source_cell->level()].size()),
167 ExcInvalidKey(source_cell));
168
169 return mapping[source_cell->level()][source_cell->index()];
170}
171
172
173
174template <class MeshType>
176void InterGridMap<MeshType>::clear()
177{
178 mapping.clear();
179 source_grid = nullptr;
180 destination_grid = nullptr;
181}
182
183
184
185template <class MeshType>
187const MeshType &InterGridMap<MeshType>::get_source_grid() const
188{
189 return *source_grid;
190}
191
192
193
194template <class MeshType>
196const MeshType &InterGridMap<MeshType>::get_destination_grid() const
197{
198 return *destination_grid;
199}
200
201
202
203template <class MeshType>
205std::size_t InterGridMap<MeshType>::memory_consumption() const
206{
209 MemoryConsumption::memory_consumption(destination_grid));
210}
211
212
213
214// explicit instantiations
215#include "intergrid_map.inst"
216
typename MeshType::cell_iterator cell_iterator
cell_iterator operator[](const cell_iterator &source_cell) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
@ valid
Iterator points to a valid object.
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
STL namespace.