19 #ifdef DEAL_II_WITH_P4EST
39 # include <functional>
56 template <
typename value_type>
59 const unsigned int dofs_per_cell)
61 for (
const auto &
values : dof_values)
67 const std::size_t bytes_per_entry =
sizeof(value_type) * dofs_per_cell;
69 std::vector<char> buffer(dof_values.size() * bytes_per_entry);
70 for (
unsigned int i = 0; i < dof_values.size(); ++i)
71 std::memcpy(&buffer[i * bytes_per_entry],
83 template <
typename value_type>
84 std::vector<Vector<value_type>>
86 const boost::iterator_range<std::vector<char>::const_iterator> &data_range,
87 const unsigned int dofs_per_cell)
89 const std::size_t bytes_per_entry =
sizeof(value_type) * dofs_per_cell;
90 const unsigned int n_elements = data_range.size() / bytes_per_entry;
94 std::vector<Vector<value_type>> unpacked_data;
95 unpacked_data.reserve(n_elements);
96 for (
unsigned int i = 0; i < n_elements; ++i)
99 std::memcpy(&dof_values(0),
100 &(*std::next(data_range.begin(), i * bytes_per_entry)),
102 unpacked_data.emplace_back(std::move(dof_values));
105 return unpacked_data;
113 namespace distributed
115 template <
int dim,
typename VectorType,
int spacedim>
118 const bool average_values)
119 : dof_handler(&dof, typeid(*this).name())
120 , average_values(average_values)
128 "parallel::distributed::SolutionTransfer requires a parallel::distributed::Triangulation object."));
133 template <
int dim,
typename VectorType,
int spacedim>
137 const std::vector<const VectorType *> &all_in)
139 for (
unsigned int i = 0; i < all_in.size(); ++i)
140 Assert(all_in[i]->size() == dof_handler->n_dofs(),
143 input_vectors = all_in;
144 register_data_attach();
149 template <
int dim,
typename VectorType,
int spacedim>
157 &dof_handler->get_triangulation())));
161 ExcMessage(
"You can only add one solution per "
162 "SolutionTransfer object."));
168 return this->pack_callback(cell_, status);
170 dof_handler->has_hp_capabilities());
175 template <
int dim,
typename VectorType,
int spacedim>
180 std::vector<const VectorType *> all_in(1, &in);
181 prepare_for_coarsening_and_refinement(all_in);
186 template <
int dim,
typename VectorType,
int spacedim>
189 const VectorType &in)
191 std::vector<const VectorType *> all_in(1, &in);
192 prepare_for_serialization(all_in);
197 template <
int dim,
typename VectorType,
int spacedim>
200 const std::vector<const VectorType *> &all_in)
202 prepare_for_coarsening_and_refinement(all_in);
207 template <
int dim,
typename VectorType,
int spacedim>
211 std::vector<VectorType *> all_in(1, &in);
217 template <
int dim,
typename VectorType,
int spacedim>
220 std::vector<VectorType *> &all_in)
222 register_data_attach();
225 input_vectors.resize(all_in.size());
231 template <
int dim,
typename VectorType,
int spacedim>
234 std::vector<VectorType *> &all_out)
236 Assert(input_vectors.size() == all_out.size(),
238 for (
unsigned int i = 0; i < all_out.size(); ++i)
239 Assert(all_out[i]->size() == dof_handler->n_dofs(),
246 &dof_handler->get_triangulation())));
250 for (
auto *
const vec : all_out)
257 valence.reinit(*all_out[0]);
261 [
this, &all_out, &valence](
264 const boost::iterator_range<std::vector<char>::const_iterator>
266 this->unpack_callback(cell_, status, data_range, all_out, valence);
272 using Number =
typename VectorType::value_type;
274 for (
const auto i : valence.locally_owned_elements())
275 valence[i] = (
static_cast<Number
>(valence[i]) == Number() ?
277 (Number(1.0) /
static_cast<Number
>(valence[i])));
280 for (
auto *
const vec : all_out)
289 for (
auto *
const vec : all_out)
293 input_vectors.clear();
299 template <
int dim,
typename VectorType,
int spacedim>
303 std::vector<VectorType *> all_out(1, &out);
309 template <
int dim,
typename VectorType,
int spacedim>
319 std::vector<::Vector<typename VectorType::value_type>> dof_values(
320 input_vectors.size());
323 if (dof_handler->has_hp_capabilities())
340 for (
const auto &child : cell->child_iterators())
341 Assert(child->is_active() && child->coarsen_flag_set(),
342 typename ::Triangulation<
343 dim>::ExcInconsistentCoarseningFlags());
346 fe_index = ::internal::hp::DoFHandlerImplementation::
347 dominated_future_fe_on_children<dim, spacedim>(cell);
357 const unsigned int dofs_per_cell =
358 dof_handler->get_fe(
fe_index).n_dofs_per_cell();
360 if (dofs_per_cell == 0)
361 return std::vector<char>();
363 auto it_input = input_vectors.cbegin();
364 auto it_output = dof_values.begin();
365 for (; it_input != input_vectors.cend(); ++it_input, ++it_output)
367 it_output->reinit(dofs_per_cell);
368 cell->get_interpolated_dof_values(*(*it_input), *it_output,
fe_index);
371 return pack_dof_values<typename VectorType::value_type>(dof_values,
377 template <
int dim,
typename VectorType,
int spacedim>
382 const boost::iterator_range<std::vector<char>::const_iterator>
384 std::vector<VectorType *> &all_out,
391 if (dof_handler->has_hp_capabilities())
409 fe_index = cell->child(0)->active_fe_index();
410 for (
unsigned int child_index = 1;
411 child_index < cell->n_children();
413 Assert(cell->child(child_index)->active_fe_index() ==
425 const unsigned int dofs_per_cell =
426 dof_handler->get_fe(
fe_index).n_dofs_per_cell();
428 if (dofs_per_cell == 0)
431 const std::vector<::Vector<typename VectorType::value_type>>
433 unpack_dof_values<typename VectorType::value_type>(data_range,
441 for (
auto it_dof_values = dof_values.begin();
442 it_dof_values != dof_values.end();
445 dofs_per_cell == it_dof_values->size(),
447 "The transferred data was packed with a different number of dofs than the "
448 "currently registered FE object assigned to the DoFHandler has."));
451 auto it_input = dof_values.cbegin();
452 auto it_output = all_out.begin();
453 for (; it_input != dof_values.cend(); ++it_input, ++it_output)
455 cell->distribute_local_to_global_by_interpolation(*it_input,
459 cell->set_dof_values_by_interpolation(*it_input,
469 cell->distribute_local_to_global_by_interpolation(ones,
479 # include "solution_transfer.inst"
@ children_will_be_coarsened
unsigned int register_data_attach(const std::function< std::vector< char >(const cell_iterator &, const ::CellStatus)> &pack_callback, const bool returns_variable_size_data)
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)
std::vector< char > pack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status)
void unpack_callback(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range, std::vector< VectorType * > &all_out, VectorType &valence)
void register_data_attach()
SmartPointer< const DoFHandler< dim, spacedim >, SolutionTransfer< dim, VectorType, spacedim > > dof_handler
void deserialize(VectorType &in)
void prepare_for_coarsening_and_refinement(const std::vector< const VectorType * > &all_in)
void prepare_for_serialization(const VectorType &in)
void interpolate(std::vector< VectorType * > &all_out)
SolutionTransfer(const DoFHandler< dim, spacedim > &dof_handler, const bool average_values=false)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMessage(std::string arg1)
static const unsigned int invalid_unsigned_int
unsigned short int fe_index
const ::Triangulation< dim, spacedim > & tria