Reference documentation for deal.II version GIT edc7d5c3ce 2023-09-25 07:10:03+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\}}\)
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 
26 namespace 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()),
91  ExcMessage(
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()),
117  ExcMessage(
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 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,
143  ExcMessage(
144  "parallel::CellWeights requires a parallel::TriangulationBase object."));
145 
146  return [&dof_handler, tria, weighting_function](
147  const typename ::Triangulation<dim, spacedim>::cell_iterator
148  &cell,
149  const CellStatus status) -> unsigned int {
151  status,
152  std::cref(
153  dof_handler),
154  std::cref(*tria),
155  weighting_function);
156  };
157  }
158 
159 
160 
161  template <int dim, int spacedim>
162  unsigned int
164  const typename ::Triangulation<dim, spacedim>::cell_iterator &cell_,
165  const CellStatus status,
166  const DoFHandler<dim, spacedim> &dof_handler,
169  &weighting_function)
170  {
171  // Check if we are still working with the correct combination of
172  // Triangulation and DoFHandler.
173  Assert(&triangulation == &(dof_handler.get_triangulation()),
174  ExcMessage(
175  "Triangulation associated with the DoFHandler has changed!"));
176  (void)triangulation;
177 
178  // Skip if the DoFHandler has not been initialized yet.
179  if (dof_handler.get_fe_collection().size() == 0)
180  return 0;
181 
182  // Convert cell type from Triangulation to DoFHandler to be able
183  // to access the information about the degrees of freedom.
184  const typename DoFHandler<dim, spacedim>::cell_iterator cell(*cell_,
185  &dof_handler);
186 
187  // Determine which FiniteElement object will be present on this cell after
188  // refinement and will thus specify the number of degrees of freedom.
190  switch (status)
191  {
195  fe_index = cell->future_fe_index();
196  break;
197 
199 #ifdef DEBUG
200  for (const auto &child : cell->child_iterators())
201  Assert(child->is_active() && child->coarsen_flag_set(),
202  typename ::Triangulation<
203  dim>::ExcInconsistentCoarseningFlags());
204 #endif
205 
206  fe_index = ::internal::hp::DoFHandlerImplementation::
207  dominated_future_fe_on_children<dim, spacedim>(cell);
208  break;
209 
210  default:
211  Assert(false, ExcInternalError());
212  break;
213  }
214 
215  // Return the cell weight determined by the function of choice.
216  return weighting_function(cell, dof_handler.get_fe(fe_index));
217  }
218 } // namespace parallel
219 
220 
221 // explicit instantiations
222 #include "cell_weights.inst"
223 
CellStatus
Definition: cell_status.h:32
@ cell_will_be_refined
@ children_will_be_coarsened
const Triangulation< dim, spacedim > & get_triangulation() const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
Signals signals
Definition: tria.h:2524
unsigned int size() const
Definition: collection.h:265
static unsigned int weighting_callback(const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status, const DoFHandler< dim, spacedim > &dof_handler, const parallel::TriangulationBase< dim, spacedim > &triangulation, const WeightingFunction &weighting_function)
static WeightingFunction constant_weighting(const unsigned int factor=1)
Definition: cell_weights.cc:66
std::function< unsigned int(const typename DoFHandler< dim, spacedim >::cell_iterator &, const FiniteElement< dim, spacedim > &)> WeightingFunction
Definition: cell_weights.h:115
static WeightingFunction ndofs_weighting(const std::pair< float, float > &coefficients)
Definition: cell_weights.cc:78
void reinit(const DoFHandler< dim, spacedim > &dof_handler, const WeightingFunction &weighting_function)
Definition: cell_weights.cc:48
static std::function< unsigned int(const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status)> make_weighting_callback(const DoFHandler< dim, spacedim > &dof_handler, const WeightingFunction &weighting_function)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1616
static ::ExceptionBase & ExcMessage(std::string arg1)
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
static const unsigned int invalid_unsigned_int
Definition: types.h:213
unsigned short int fe_index
Definition: types.h:60
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
const ::Triangulation< dim, spacedim > & tria