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
cell_weights.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2018 - 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
18
20
21#include <limits>
22
24
25
26namespace parallel
27{
28 template <int dim, int spacedim>
30 const DoFHandler<dim, spacedim> &dof_handler,
31 const WeightingFunction & weighting_function)
32 {
33 reinit(dof_handler, weighting_function);
34 }
35
36
37
38 template <int dim, int spacedim>
40 {
41 connection.disconnect();
42 }
43
44
45
46 template <int dim, int spacedim>
47 void
49 const DoFHandler<dim, spacedim> &dof_handler,
51 &weighting_function)
52 {
53 connection.disconnect();
54 connection = dof_handler.get_triangulation().signals.weight.connect(
55 make_weighting_callback(dof_handler, weighting_function));
56 }
57
58
59
60 // ---------- static member functions ----------
61
62 // ---------- selection of weighting functions ----------
63
64 template <int dim, int spacedim>
67 {
68 return [factor](const typename DoFHandler<dim, spacedim>::cell_iterator &,
69 const FiniteElement<dim, spacedim> &) -> unsigned int {
70 return factor;
71 };
72 }
73
74
75
76 template <int dim, int spacedim>
79 const std::pair<float, float> &coefficients)
80 {
81 return [coefficients](
83 const FiniteElement<dim, spacedim> &future_fe) -> unsigned int {
84 const float result =
85 std::trunc(coefficients.first *
86 std::pow(future_fe.n_dofs_per_cell(), coefficients.second));
87
88 Assert(result >= 0. &&
89 result <=
90 static_cast<float>(std::numeric_limits<unsigned int>::max()),
92 "Cannot cast determined weight for this cell to unsigned int!"));
93
94 return static_cast<unsigned int>(result);
95 };
96 }
97
98
99
100 template <int dim, int spacedim>
103 const std::vector<std::pair<float, float>> &coefficients)
104 {
105 return [coefficients](
107 const FiniteElement<dim, spacedim> &future_fe) -> unsigned int {
108 float result = 0;
109 for (const auto &pair : coefficients)
110 result +=
111 pair.first * std::pow(future_fe.n_dofs_per_cell(), pair.second);
112 result = std::trunc(result);
113
114 Assert(result >= 0. &&
115 result <=
116 static_cast<float>(std::numeric_limits<unsigned int>::max()),
118 "Cannot cast determined weight for this cell to unsigned int!"));
119
120 return static_cast<unsigned int>(result);
121 };
122 }
123
124
125
126 // ---------- handling callback functions ----------
127
128 template <int dim, int spacedim>
129 std::function<unsigned int(
130 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell,
131 const typename ::Triangulation<dim, spacedim>::CellStatus status)>
133 const DoFHandler<dim, spacedim> &dof_handler,
135 &weighting_function)
136 {
138 dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
139 &(dof_handler.get_triangulation()));
140
141 Assert(
142 tria != nullptr,
144 "parallel::CellWeights requires a parallel::TriangulationBase object."));
145
146 return
147 [&dof_handler, tria, weighting_function](
148 const typename ::Triangulation<dim, spacedim>::cell_iterator
149 & cell,
150 const typename ::Triangulation<dim, spacedim>::CellStatus status)
151 -> unsigned int {
153 cell,
154 status,
155 std::cref(dof_handler),
156 std::cref(*tria),
157 weighting_function);
158 };
159 }
160
161
162
163 template <int dim, int spacedim>
164 unsigned int
166 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell_,
167 const typename ::Triangulation<dim, spacedim>::CellStatus status,
168 const DoFHandler<dim, spacedim> & dof_handler,
171 &weighting_function)
172 {
173 // Check if we are still working with the correct combination of
174 // Triangulation and DoFHandler.
175 Assert(&triangulation == &(dof_handler.get_triangulation()),
177 "Triangulation associated with the DoFHandler has changed!"));
178 (void)triangulation;
179
180 // Skip if the DoFHandler has not been initialized yet.
181 if (dof_handler.get_fe_collection().size() == 0)
182 return 0;
183
184 // Convert cell type from Triangulation to DoFHandler to be able
185 // to access the information about the degrees of freedom.
186 const typename DoFHandler<dim, spacedim>::cell_iterator cell(*cell_,
187 &dof_handler);
188
189 // Determine which FiniteElement object will be present on this cell after
190 // refinement and will thus specify the number of degrees of freedom.
191 unsigned int fe_index = numbers::invalid_unsigned_int;
192 switch (status)
193 {
197 fe_index = cell->future_fe_index();
198 break;
199
201#ifdef DEBUG
202 for (const auto &child : cell->child_iterators())
203 Assert(child->is_active() && child->coarsen_flag_set(),
204 typename ::Triangulation<
205 dim>::ExcInconsistentCoarseningFlags());
206#endif
207
208 fe_index = ::internal::hp::DoFHandlerImplementation::
209 dominated_future_fe_on_children<dim, spacedim>(cell);
210 break;
211
212 default:
213 Assert(false, ExcInternalError());
214 break;
215 }
216
217 // Return the cell weight determined by the function of choice.
218 return weighting_function(cell, dof_handler.get_fe(fe_index));
219 }
220} // namespace parallel
221
222
223// explicit instantiations
224#include "cell_weights.inst"
225
const hp::FECollection< dim, spacedim > & get_fe_collection() const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const Triangulation< dim, spacedim > & get_triangulation() const
Signals signals
Definition tria.h:2479
unsigned int size() const
Definition collection.h:265
static unsigned int weighting_callback(const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const typename ::Triangulation< dim, spacedim >::CellStatus status, const DoFHandler< dim, spacedim > &dof_handler, const parallel::TriangulationBase< dim, spacedim > &triangulation, const WeightingFunction &weighting_function)
static std::function< unsigned int(const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const typename ::Triangulation< dim, spacedim >::CellStatus status)> make_weighting_callback(const DoFHandler< dim, spacedim > &dof_handler, const WeightingFunction &weighting_function)
static WeightingFunction constant_weighting(const unsigned int factor=1)
std::function< unsigned int(const typename DoFHandler< dim, spacedim >::cell_iterator &, const FiniteElement< dim, spacedim > &)> WeightingFunction
static WeightingFunction ndofs_weighting(const std::pair< float, float > &coefficients)
void reinit(const DoFHandler< dim, spacedim > &dof_handler, const WeightingFunction &weighting_function)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ActiveSelector::cell_iterator cell_iterator
static const unsigned int invalid_unsigned_int
Definition types.h:213
::VectorizedArray< Number, width > pow(const ::VectorizedArray< Number, width > &, const Number p)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
const ::Triangulation< dim, spacedim > & tria