Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
cell_data_transfer.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 - 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 #ifndef dealii_distributed_cell_data_transfer_h
17 #define dealii_distributed_cell_data_transfer_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/base/smartpointer.h>
22 
23 #include <deal.II/distributed/tria.h>
24 
25 #include <boost/range/iterator_range.hpp>
26 
27 #include <functional>
28 
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 namespace parallel
33 {
34  namespace distributed
35  {
125  template <int dim, int spacedim = dim, typename VectorType = Vector<double>>
127  {
128  public:
139  {
144  static typename VectorType::value_type
147  const VectorType &input_vector)
148  {
149  const typename VectorType::value_type value =
150  input_vector[parent->child(0)->active_cell_index()];
151 
152  for (unsigned int child_index = 1; child_index < parent->n_children();
153  ++child_index)
154  Assert(
155  value ==
156  input_vector[parent->child(child_index)->active_cell_index()],
157  ExcMessage(
158  "Values on cells that will be coarsened are not equal!"));
159 
160  return value;
161  }
162 
167  static typename VectorType::value_type
169  cell_iterator & parent,
170  const VectorType &input_vector)
171  {
172  typename VectorType::value_type sum = 0;
173 
174  for (unsigned int child_index = 0; child_index < parent->n_children();
175  ++child_index)
176  sum +=
177  input_vector[parent->child(child_index)->active_cell_index()];
178 
179  return sum;
180  }
181 
186  static typename VectorType::value_type
187  mean(const typename parallel::distributed::
189  const VectorType & input_vector)
190  {
191  return sum(parent, input_vector) / parent->n_children();
192  }
193  };
194 
209  & triangulation,
210  const bool transfer_variable_size_data = false,
211  const std::function<typename VectorType::value_type(
212  const typename parallel::distributed::
214  const VectorType &input_vector)> coarsening_strategy =
216 
228  void
229  prepare_for_coarsening_and_refinement(const VectorType &in);
230 
234  void
236  const std::vector<const VectorType *> &all_in);
237 
247  void
248  prepare_for_serialization(const VectorType &in);
249 
253  void
254  prepare_for_serialization(const std::vector<const VectorType *> &all_in);
255 
260  void
261  unpack(VectorType &out);
262 
266  void
267  unpack(std::vector<VectorType *> &all_out);
268 
273  void
274  deserialize(VectorType &out);
275 
279  void
280  deserialize(std::vector<VectorType *> &all_out);
281 
282  private:
289 
294 
299  const std::function<typename VectorType::value_type(
301  cell_iterator & parent,
302  const VectorType &input_vector)>
304 
309  std::vector<const VectorType *> input_vectors;
310 
315  unsigned int handle;
316 
321  void
323 
329  std::vector<char>
330  pack_callback(const typename parallel::distributed::
332  const typename parallel::distributed::
333  Triangulation<dim, spacedim>::CellStatus status);
334 
340  void
343  cell_iterator &cell,
345  CellStatus status,
346  const boost::iterator_range<std::vector<char>::const_iterator>
347  & data_range,
348  std::vector<VectorType *> &all_out);
349  };
350  } // namespace distributed
351 } // namespace parallel
352 
353 
354 DEAL_II_NAMESPACE_CLOSE
355 
356 #endif
CellDataTransfer(const parallel::distributed::Triangulation< dim, spacedim > &triangulation, const bool transfer_variable_size_data=false, const std::function< typename VectorType::value_type(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &parent, const VectorType &input_vector)> coarsening_strategy=&CoarseningStrategies::check_equality)
std::vector< const VectorType * > input_vectors
std::vector< char > pack_callback(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &cell, const typename parallel::distributed::Triangulation< dim, spacedim >::CellStatus status)
static VectorType::value_type mean(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &parent, const VectorType &input_vector)
SmartPointer< const parallel::distributed::Triangulation< dim, spacedim >, CellDataTransfer< dim, spacedim, VectorType > > triangulation
static VectorType::value_type sum(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &parent, const VectorType &input_vector)
static ::ExceptionBase & ExcMessage(std::string arg1)
static VectorType::value_type check_equality(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &parent, const VectorType &input_vector)
#define Assert(cond, exc)
Definition: exceptions.h:1407
const std::function< typename VectorType::value_type(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &parent, const VectorType &input_vector)> coarsening_strategy
void prepare_for_coarsening_and_refinement(const VectorType &in)
void unpack_callback(const typename parallel::distributed::Triangulation< dim, spacedim >::cell_iterator &cell, const typename parallel::distributed::Triangulation< dim, spacedim >::CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range, std::vector< VectorType *> &all_out)
void prepare_for_serialization(const VectorType &in)
typename ::Triangulation< dim, spacedim >::cell_iterator cell_iterator
Definition: tria.h:269