Reference documentation for deal.II version Git d3aed38b93 2021-10-28 13:33:27 +0200
\(\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\}}\)
solution_transfer.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2009 - 2021 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 
17 #include <deal.II/base/config.h>
18 
19 #ifdef DEAL_II_WITH_P4EST
20 
23 
25 # include <deal.II/dofs/dof_tools.h>
26 
29 
30 # include <deal.II/hp/dof_handler.h>
31 
39 # include <deal.II/lac/vector.h>
40 
41 # include <functional>
42 # include <numeric>
43 
44 
46 
47 
48 namespace
49 {
58  template <typename value_type>
59  std::vector<char>
60  pack_dof_values(std::vector<Vector<value_type>> &dof_values,
61  const unsigned int dofs_per_cell)
62  {
63  for (const auto &values : dof_values)
64  {
65  AssertDimension(values.size(), dofs_per_cell);
66  (void)values;
67  }
68 
69  const std::size_t bytes_per_entry = sizeof(value_type) * dofs_per_cell;
70 
71  std::vector<char> buffer(dof_values.size() * bytes_per_entry);
72  for (unsigned int i = 0; i < dof_values.size(); ++i)
73  std::memcpy(&buffer[i * bytes_per_entry],
74  &dof_values[i](0),
75  bytes_per_entry);
76 
77  return buffer;
78  }
79 
80 
81 
85  template <typename value_type>
86  std::vector<Vector<value_type>>
87  unpack_dof_values(
88  const boost::iterator_range<std::vector<char>::const_iterator> &data_range,
89  const unsigned int dofs_per_cell)
90  {
91  const std::size_t bytes_per_entry = sizeof(value_type) * dofs_per_cell;
92  const unsigned int n_elements = data_range.size() / bytes_per_entry;
93 
94  Assert((data_range.size() % bytes_per_entry == 0), ExcInternalError());
95 
96  std::vector<Vector<value_type>> unpacked_data;
97  unpacked_data.reserve(n_elements);
98  for (unsigned int i = 0; i < n_elements; ++i)
99  {
100  Vector<value_type> dof_values(dofs_per_cell);
101  std::memcpy(&dof_values(0),
102  &(*std::next(data_range.begin(), i * bytes_per_entry)),
103  bytes_per_entry);
104  unpacked_data.emplace_back(std::move(dof_values));
105  }
106 
107  return unpacked_data;
108  }
109 } // namespace
110 
111 
112 
113 namespace parallel
114 {
115  namespace distributed
116  {
117  template <int dim, typename VectorType, int spacedim>
119  const DoFHandler<dim, spacedim> &dof)
120  : dof_handler(&dof, typeid(*this).name())
121  , handle(numbers::invalid_unsigned_int)
122  {
123  Assert(
124  (dynamic_cast<
126  &dof_handler->get_triangulation()) != nullptr),
127  ExcMessage(
128  "parallel::distributed::SolutionTransfer requires a parallel::distributed::Triangulation object."));
129  }
130 
131 
132 
133  template <int dim, typename VectorType, int spacedim>
134  void
137  const std::vector<const VectorType *> &all_in)
138  {
139  for (unsigned int i = 0; i < all_in.size(); ++i)
140  Assert(all_in[i]->size() == dof_handler->n_dofs(),
141  ExcDimensionMismatch(all_in[i]->size(), dof_handler->n_dofs()));
142 
143  input_vectors = all_in;
145  }
146 
147 
148 
149  template <int dim, typename VectorType, int spacedim>
150  void
152  {
153  // TODO: casting away constness is bad
156  const_cast<::Triangulation<dim, spacedim> *>(
157  &dof_handler->get_triangulation())));
158  Assert(tria != nullptr, ExcInternalError());
159 
161  [this](
162  const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
163  const typename Triangulation<dim, spacedim>::CellStatus status) {
164  return this->pack_callback(cell_, status);
165  },
166  /*returns_variable_size_data=*/dof_handler->has_hp_capabilities());
167  }
168 
169 
170 
171  template <int dim, typename VectorType, int spacedim>
172  void
175  {
176  std::vector<const VectorType *> all_in(1, &in);
178  }
179 
180 
181 
182  template <int dim, typename VectorType, int spacedim>
183  void
185  const VectorType &in)
186  {
187  std::vector<const VectorType *> all_in(1, &in);
189  }
190 
191 
192 
193  template <int dim, typename VectorType, int spacedim>
194  void
196  const std::vector<const VectorType *> &all_in)
197  {
199  }
200 
201 
202 
203  template <int dim, typename VectorType, int spacedim>
204  void
206  {
207  std::vector<VectorType *> all_in(1, &in);
208  deserialize(all_in);
209  }
210 
211 
212 
213  template <int dim, typename VectorType, int spacedim>
214  void
216  std::vector<VectorType *> &all_in)
217  {
219 
220  // this makes interpolate() happy
221  input_vectors.resize(all_in.size());
222 
223  interpolate(all_in);
224  }
225 
226 
227  template <int dim, typename VectorType, int spacedim>
228  void
230  std::vector<VectorType *> &all_out)
231  {
232  Assert(input_vectors.size() == all_out.size(),
233  ExcDimensionMismatch(input_vectors.size(), all_out.size()));
234  for (unsigned int i = 0; i < all_out.size(); ++i)
235  Assert(all_out[i]->size() == dof_handler->n_dofs(),
236  ExcDimensionMismatch(all_out[i]->size(), dof_handler->n_dofs()));
237 
238  // TODO: casting away constness is bad
241  const_cast<::Triangulation<dim, spacedim> *>(
242  &dof_handler->get_triangulation())));
243  Assert(tria != nullptr, ExcInternalError());
244 
246  handle,
247  [this, &all_out](
248  const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
249  const typename Triangulation<dim, spacedim>::CellStatus status,
250  const boost::iterator_range<std::vector<char>::const_iterator>
251  &data_range) {
252  this->unpack_callback(cell_, status, data_range, all_out);
253  });
254 
255  for (typename std::vector<VectorType *>::iterator it = all_out.begin();
256  it != all_out.end();
257  ++it)
258  (*it)->compress(::VectorOperation::insert);
259 
260  input_vectors.clear();
261  }
262 
263 
264 
265  template <int dim, typename VectorType, int spacedim>
266  void
268  {
269  std::vector<VectorType *> all_out(1, &out);
270  interpolate(all_out);
271  }
272 
273 
274 
275  template <int dim, typename VectorType, int spacedim>
276  std::vector<char>
278  const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
279  const typename Triangulation<dim, spacedim>::CellStatus status)
280  {
281  typename DoFHandler<dim, spacedim>::cell_iterator cell(*cell_,
282  dof_handler);
283 
284  // create buffer for each individual object
285  std::vector<::Vector<typename VectorType::value_type>> dof_values(
286  input_vectors.size());
287 
288  unsigned int fe_index = 0;
289  if (dof_handler->has_hp_capabilities())
290  {
291  switch (status)
292  {
294  spacedim>::CELL_PERSIST:
296  spacedim>::CELL_REFINE:
297  {
298  fe_index = cell->future_fe_index();
299  break;
300  }
301 
303  spacedim>::CELL_COARSEN:
304  {
305  // In case of coarsening, we need to find a suitable FE index
306  // for the parent cell. We choose the 'least dominant fe'
307  // on all children from the associated FECollection.
308 # ifdef DEBUG
309  for (const auto &child : cell->child_iterators())
310  Assert(child->is_active() && child->coarsen_flag_set(),
311  typename ::Triangulation<
312  dim>::ExcInconsistentCoarseningFlags());
313 # endif
314 
315  fe_index = ::internal::hp::DoFHandlerImplementation::
316  dominated_future_fe_on_children<dim, spacedim>(cell);
317  break;
318  }
319 
320  default:
321  Assert(false, ExcInternalError());
322  break;
323  }
324  }
325 
326  const unsigned int dofs_per_cell =
327  dof_handler->get_fe(fe_index).n_dofs_per_cell();
328 
329  if (dofs_per_cell == 0)
330  return std::vector<char>(); // nothing to do for FE_Nothing
331 
332  auto it_input = input_vectors.cbegin();
333  auto it_output = dof_values.begin();
334  for (; it_input != input_vectors.cend(); ++it_input, ++it_output)
335  {
336  it_output->reinit(dofs_per_cell);
337  cell->get_interpolated_dof_values(*(*it_input), *it_output, fe_index);
338  }
339 
340  return pack_dof_values<typename VectorType::value_type>(dof_values,
341  dofs_per_cell);
342  }
343 
344 
345 
346  template <int dim, typename VectorType, int spacedim>
347  void
349  const typename Triangulation<dim, spacedim>::cell_iterator &cell_,
350  const typename Triangulation<dim, spacedim>::CellStatus status,
351  const boost::iterator_range<std::vector<char>::const_iterator>
352  & data_range,
353  std::vector<VectorType *> &all_out)
354  {
355  typename DoFHandler<dim, spacedim>::cell_iterator cell(*cell_,
356  dof_handler);
357 
358  unsigned int fe_index = 0;
359  if (dof_handler->has_hp_capabilities())
360  {
361  switch (status)
362  {
364  spacedim>::CELL_PERSIST:
366  spacedim>::CELL_COARSEN:
367  {
368  fe_index = cell->active_fe_index();
369  break;
370  }
371 
373  spacedim>::CELL_REFINE:
374  {
375  // After refinement, this particular cell is no longer active,
376  // and its children have inherited its FE index. However, to
377  // unpack the data on the old cell, we need to recover its FE
378  // index from one of the children. Just to be sure, we also
379  // check if all children have the same FE index.
380  fe_index = cell->child(0)->active_fe_index();
381  for (unsigned int child_index = 1;
382  child_index < cell->n_children();
383  ++child_index)
384  Assert(cell->child(child_index)->active_fe_index() ==
385  fe_index,
386  ExcInternalError());
387  break;
388  }
389 
390  default:
391  Assert(false, ExcInternalError());
392  break;
393  }
394  }
395 
396  const unsigned int dofs_per_cell =
397  dof_handler->get_fe(fe_index).n_dofs_per_cell();
398 
399  if (dofs_per_cell == 0)
400  return; // nothing to do for FE_Nothing
401 
402  const std::vector<::Vector<typename VectorType::value_type>>
403  dof_values =
404  unpack_dof_values<typename VectorType::value_type>(data_range,
405  dofs_per_cell);
406 
407  // check if sizes match
408  Assert(dof_values.size() == all_out.size(), ExcInternalError());
409 
410  // check if we have enough dofs provided by the FE object
411  // to interpolate the transferred data correctly
412  for (auto it_dof_values = dof_values.begin();
413  it_dof_values != dof_values.end();
414  ++it_dof_values)
415  Assert(
416  dofs_per_cell == it_dof_values->size(),
417  ExcMessage(
418  "The transferred data was packed with a different number of dofs than the "
419  "currently registered FE object assigned to the DoFHandler has."));
420 
421  // distribute data for each registered vector on mesh
422  auto it_input = dof_values.cbegin();
423  auto it_output = all_out.begin();
424  for (; it_input != dof_values.cend(); ++it_input, ++it_output)
425  cell->set_dof_values_by_interpolation(*it_input,
426  *(*it_output),
427  fe_index);
428  }
429  } // namespace distributed
430 } // namespace parallel
431 
432 
433 // explicit instantiations
434 # include "solution_transfer.inst"
435 
437 
438 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
void notify_ready_to_unpack(const unsigned int handle, const std::function< void(const cell_iterator &, const CellStatus, const boost::iterator_range< std::vector< char >::const_iterator > &)> &unpack_callback)
Definition: tria_base.cc:790
const ::Triangulation< dim, spacedim > & tria
void interpolate(std::vector< VectorType *> &all_out)
SolutionTransfer(const DoFHandler< dim, spacedim > &dof)
static ::ExceptionBase & ExcMessage(std::string arg1)
typename ::Triangulation< dim, spacedim >::CellStatus CellStatus
Definition: tria.h:294
#define Assert(cond, exc)
Definition: exceptions.h:1461
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
std::vector< char > pack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::CellStatus status)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
SmartPointer< const DoFHandler< dim, spacedim >, SolutionTransfer< dim, VectorType, spacedim > > dof_handler
unsigned int register_data_attach(const std::function< std::vector< char >(const cell_iterator &, const CellStatus)> &pack_callback, const bool returns_variable_size_data)
Definition: tria_base.cc:761
std::vector< const VectorType * > input_vectors
void prepare_for_serialization(const VectorType &in)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:466
void prepare_for_coarsening_and_refinement(const std::vector< const VectorType *> &all_in)
void unpack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename Triangulation< dim, spacedim >::CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range, std::vector< VectorType *> &all_out)
typename ::Triangulation< dim, spacedim >::cell_iterator cell_iterator
Definition: tria.h:270
static ::ExceptionBase & ExcInternalError()