/**
* Constructor.
*
- * @param[in] dof The DoFHandler on which all operations will happen.
- * At the time when this constructor is called, the DoFHandler still
- * points to the Triangulation before the refinement in question
+ * @param[in] dof_handler The DoFHandler on which all operations will
+ * happen. At the time when this constructor is called, the DoFHandler
+ * still points to the Triangulation before the refinement in question
* happens.
+ * @param[in] average_values Average the contribututions to the same
+ * DoF coming from different cells. Note: averaging requires an
+ * additional communication step, since the valence of the DoF has to be
+ * determined.
*/
- SolutionTransfer(const DoFHandler<dim, spacedim> &dof);
+ SolutionTransfer(const DoFHandler<dim, spacedim> &dof_handler,
+ const bool average_values = false);
/**
* Destructor.
SolutionTransfer<dim, VectorType, spacedim>>
dof_handler;
+ /**
+ * Flag indicating if averaging should be performed.
+ */
+ const bool average_values;
+
/**
* A vector that stores pointers to all the vectors we are supposed to
* copy over from the old to the new mesh.
const typename Triangulation<dim, spacedim>::CellStatus status,
const boost::iterator_range<std::vector<char>::const_iterator>
& data_range,
- std::vector<VectorType *> &all_out);
+ std::vector<VectorType *> &all_out,
+ VectorType & valence);
/**
DoFHandler<dimension_, space_dimension_>::invalid_fe_index,
const bool perform_check = false) const;
+ /**
+ * Similar to set_dof_values_by_interpolation() with the difference that
+ * values are added into the vector.
+ *
+ * @note In parallel::distributed::SolutionTransfer, this function is used
+ * to accumulate the contributions of all cells to a DoF; with a
+ * subsequent multiplication with the inverse of the valence, finally,
+ * the average value is obtained.
+ */
+ template <class OutputVector, typename number>
+ void
+ distribute_local_to_global_by_interpolation(
+ const Vector<number> &local_values,
+ OutputVector & values,
+ const unsigned int fe_index =
+ DoFHandler<dimension_, space_dimension_>::invalid_fe_index) const;
+
/**
* Distribute a local (cell based) vector to a global one by mapping the
* local numbering of the degrees of freedom to the global one and entering
{
template <int dim, typename VectorType, int spacedim>
SolutionTransfer<dim, VectorType, spacedim>::SolutionTransfer(
- const DoFHandler<dim, spacedim> &dof)
+ const DoFHandler<dim, spacedim> &dof,
+ const bool average_values)
: dof_handler(&dof, typeid(*this).name())
+ , average_values(average_values)
, handle(numbers::invalid_unsigned_int)
{
Assert(
&dof_handler->get_triangulation())));
Assert(tria != nullptr, ExcInternalError());
+ for (const auto vec : all_out)
+ *vec = 0.0;
+
+ VectorType valence;
+
+ // initialize valence vector only if we need to average
+ if (average_values)
+ valence.reinit(*all_out[0]);
+
tria->notify_ready_to_unpack(
handle,
- [this, &all_out](
+ [this, &all_out, &valence](
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) {
- this->unpack_callback(cell_, status, data_range, all_out);
+ this->unpack_callback(cell_, status, data_range, all_out, valence);
});
- for (typename std::vector<VectorType *>::iterator it = all_out.begin();
- it != all_out.end();
- ++it)
- (*it)->compress(::dealii::VectorOperation::insert);
+ if (average_values)
+ {
+ // finalize valence: compress and invert
+ valence.compress(::dealii::VectorOperation::add);
+ for (const auto i : valence.locally_owned_elements())
+ valence[i] = valence[i] == 0.0 ? 0.0 : (1.0 / valence[i]);
+ valence.compress(::dealii::VectorOperation::insert);
+
+ for (const auto vec : all_out)
+ {
+ // comress and weight with valence
+ vec->compress(::dealii::VectorOperation::add);
+ vec->scale(valence);
+ }
+ }
+ else
+ {
+ for (const auto vec : all_out)
+ vec->compress(::dealii::VectorOperation::insert);
+ }
input_vectors.clear();
}
const typename Triangulation<dim, spacedim>::CellStatus status,
const boost::iterator_range<std::vector<char>::const_iterator>
& data_range,
- std::vector<VectorType *> &all_out)
+ std::vector<VectorType *> &all_out,
+ VectorType & valence)
{
typename DoFHandler<dim, spacedim>::cell_iterator cell(*cell_,
dof_handler);
auto it_input = dof_values.cbegin();
auto it_output = all_out.begin();
for (; it_input != dof_values.cend(); ++it_input, ++it_output)
- cell->set_dof_values_by_interpolation(*it_input,
- *(*it_output),
- fe_index,
- true);
+ if (average_values)
+ cell->distribute_local_to_global_by_interpolation(*it_input,
+ *(*it_output),
+ fe_index);
+ else
+ cell->set_dof_values_by_interpolation(*it_input,
+ *(*it_output),
+ fe_index,
+ true);
+
+ if (average_values)
+ {
+ // compute valence vector if averaging should be performed
+ Vector<typename VectorType::value_type> ones(dofs_per_cell);
+ ones = 1.0;
+ cell->distribute_local_to_global_by_interpolation(ones,
+ valence,
+ fe_index);
+ }
}
} // namespace distributed
} // namespace parallel
}
+template <int dim, int spacedim, bool lda>
+template <class OutputVector, typename number>
+void
+DoFCellAccessor<dim, spacedim, lda>::
+ distribute_local_to_global_by_interpolation(
+ const Vector<number> &local_values,
+ OutputVector & values,
+ const unsigned int fe_index_) const
+{
+ internal::process_by_interpolation<dim, spacedim, lda, OutputVector, number>(
+ *this,
+ local_values,
+ values,
+ fe_index_,
+ [](const DoFCellAccessor<dim, spacedim, lda> &cell,
+ const Vector<number> & local_values,
+ OutputVector & values) {
+ std::vector<types::global_dof_index> dof_indices(
+ cell.get_fe().n_dofs_per_cell());
+ cell.get_dof_indices(dof_indices);
+ AffineConstraints<number>().distribute_local_to_global(local_values,
+ dof_indices,
+ values);
+ });
+}
+
+
// --------------------------------------------------------------------------
// explicit instantiations
#include "dof_accessor_set.inst"
const unsigned int fe_index,
const bool) const;
+ template void DoFCellAccessor<deal_II_dimension, deal_II_dimension, lda>::
+ distribute_local_to_global_by_interpolation(
+ const Vector<VEC::value_type> &, VEC &, const unsigned int fe_index)
+ const;
+
#if deal_II_dimension != 3
template void
const unsigned int fe_index,
const bool) const;
+ template void
+ DoFCellAccessor<deal_II_dimension, deal_II_dimension + 1, lda>::
+ distribute_local_to_global_by_interpolation(
+ const Vector<VEC::value_type> &, VEC &, const unsigned int fe_index)
+ const;
+
#endif
#if deal_II_dimension == 3
const unsigned int fe_index,
const bool) const;
+ template void
+ DoFCellAccessor<1, 3, lda>::distribute_local_to_global_by_interpolation(
+ const Vector<VEC::value_type> &, VEC &, const unsigned int fe_index)
+ const;
+
#endif
}