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