Skip to content

Commit

Permalink
MGTransferMF: check compatibility of DoFHandlers
Browse files Browse the repository at this point in the history
  • Loading branch information
@peterrum
peterrum committed Sep 29, 2023
1 parent 877a30a commit e3cfd05
Show file tree
Hide file tree
Showing 2 changed files with 192 additions and 25 deletions.
87 changes: 62 additions & 25 deletions include/deal.II/multigrid/mg_transfer_global_coarsening.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,9 @@ namespace RepartitioningPolicyTools
template <int dim, int spacedim>
class Base;
}

template <int dim, typename Number>
class MGTransferMF;
#endif


Expand Down Expand Up @@ -699,7 +702,19 @@ class MGTwoLevelTransfer<dim, LinearAlgebra::distributed::Vector<Number>>
*/
unsigned int n_components;

/**
* Pointer to the DoFHandler object used during initialization.
*/
SmartPointer<const DoFHandler<dim>> dof_handler_fine;

/**
* Muligird level used during initialization.
*/
unsigned int mg_level_fine;

friend class internal::MGTwoLevelTransferImplementation;

friend class MGTransferMF<dim, Number>;
};


Expand Down Expand Up @@ -912,6 +927,16 @@ class MGTwoLevelTransferNonNested<dim,
LinearAlgebra::distributed::Vector<Number> &dst,
const LinearAlgebra::distributed::Vector<Number> &src) const;

/**
* Pointer to the DoFHandler object used during initialization.
*/
SmartPointer<const DoFHandler<dim>> dof_handler_fine;

/**
* Multigrid level used during initialization.
*/
unsigned int mg_level_fine;

/**
* Object to evaluate shape functions on one mesh on visited support points of
* the other mesh.
Expand Down Expand Up @@ -947,6 +972,8 @@ class MGTwoLevelTransferNonNested<dim,
* point.
*/
std::vector<unsigned int> level_dof_indices_fine_ptrs;

friend class MGTransferMF<dim, Number>;
};


Expand Down Expand Up @@ -1156,19 +1183,19 @@ class MGTransferMF : public dealii::MGLevelGlobalTransfer<
*/
template <class InVector>
void
interpolate_to_mg(MGLevelObject<VectorType> &dst, const InVector &src) const;
interpolate_to_mg(const DoFHandler<dim> &dof_handler,
MGLevelObject<VectorType> &dst,
const InVector &src) const;

/**
* Like the above function but with a user-provided DoFHandler as
* additional argument. However, this DoFHandler is not used internally, but
* is required to be able to use MGTransferMF and
* MGTransferMatrixFree as template argument.
* Interpolate fine-mesh field @p src to each multigrid level in
* @p dof_handler and store the result in @p dst.
*/
template <class InVector>
void
interpolate_to_mg(const DoFHandler<dim> &dof_handler,
MGLevelObject<VectorType> &dst,
const InVector &src) const;
DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
"Use the other version that takes a DoFHandler!")
void interpolate_to_mg(MGLevelObject<VectorType> &dst,
const InVector &src) const;

/** @} */

Expand Down Expand Up @@ -1236,6 +1263,14 @@ class MGTransferMF : public dealii::MGLevelGlobalTransfer<
const InVector &vector_reference,
const bool omit_zeroing_entries) const;

/**
* Check that the internal DoFHandler is compatible with the external one
* provided by copy_to_mg(), copy_from_mg() and interpolate_to_mg()
* used, e.g., by PreconditionMG.
*/
void
assert_dof_handler(const DoFHandler<dim> &dof_handler_out) const;

/**
* Internal transfer operator.
*
Expand Down Expand Up @@ -1458,7 +1493,7 @@ MGTransferMF<dim, Number>::copy_to_mg(const DoFHandler<dim> &dof_handler,
MGLevelObject<VectorType> &dst,
const InVector &src) const
{
(void)dof_handler;
assert_dof_handler(dof_handler);

for (unsigned int level = dst.min_level(); level <= dst.max_level(); ++level)
{
Expand Down Expand Up @@ -1517,7 +1552,7 @@ MGTransferMF<dim, Number>::copy_from_mg(
OutVector &dst,
const MGLevelObject<VectorType> &src) const
{
(void)dof_handler;
assert_dof_handler(dof_handler);

if (this->perform_plain_copy)
{
Expand Down Expand Up @@ -1569,6 +1604,22 @@ void
MGTransferMF<dim, Number>::interpolate_to_mg(MGLevelObject<VectorType> &dst,
const InVector &src) const
{
DoFHandler<dim> dof_handler_dummy;

this->interpolate_to_mg(dof_handler_dummy, dst, src);
}



template <int dim, typename Number>
template <class InVector>
void
MGTransferMF<dim, Number>::interpolate_to_mg(const DoFHandler<dim> &dof_handler,
MGLevelObject<VectorType> &dst,
const InVector &src) const
{
assert_dof_handler(dof_handler);

const unsigned int min_level = transfer.min_level();
const unsigned int max_level = transfer.max_level();

Expand Down Expand Up @@ -1629,20 +1680,6 @@ MGTransferMF<dim, Number>::interpolate_to_mg(MGLevelObject<VectorType> &dst,
}
}



template <int dim, typename Number>
template <class InVector>
void
MGTransferMF<dim, Number>::interpolate_to_mg(const DoFHandler<dim> &dof_handler,
MGLevelObject<VectorType> &dst,
const InVector &src) const
{
(void)dof_handler;

this->interpolate_to_mg(dst, src);
}

#endif

DEAL_II_NAMESPACE_CLOSE
Expand Down
130 changes: 130 additions & 0 deletions include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1783,6 +1783,9 @@ namespace internal
(mg_level_coarse + 1 == mg_level_fine),
ExcNotImplemented());

transfer.dof_handler_fine = &dof_handler_fine;
transfer.mg_level_fine = mg_level_fine;

std::unique_ptr<FineDoFHandlerViewBase<dim>> dof_handler_fine_view;

if (internal::h_transfer_uses_first_child_policy(dof_handler_fine,
Expand Down Expand Up @@ -2323,6 +2326,9 @@ namespace internal
"(numbers::invalid_unsigned_int) or on refinement levels without "
"hanging nodes."));

transfer.dof_handler_fine = &dof_handler_fine;
transfer.mg_level_fine = mg_level_fine;

std::unique_ptr<FineDoFHandlerViewBase<dim>> dof_handler_fine_view;

if (internal::p_transfer_involves_repartitioning(dof_handler_fine,
Expand Down Expand Up @@ -4237,6 +4243,127 @@ MGTransferMF<dim, Number>::restrict_and_add(const unsigned int from_level,



template <int dim, typename Number>
void
MGTransferMF<dim, Number>::assert_dof_handler(
const DoFHandler<dim> &dof_handler_out) const
{
#ifndef DEBUG
(void)dof_handler_out;
#else

if (dof_handler_out.n_dofs() == 0)
return; // provided DoFHandler object is dummy

if (this->transfer.n_levels() <= 1)
return; // single level: no check possible

const DoFHandler<dim> *dof_handler_in = nullptr;
unsigned int level_in = 0;

if (const auto t = dynamic_cast<
const MGTwoLevelTransfer<dim,
LinearAlgebra::distributed::Vector<Number>> *>(
this->transfer[this->transfer.max_level()].get()))
{
dof_handler_in = t->dof_handler_fine;
level_in = t->mg_level_fine;
}
else if (const auto t = dynamic_cast<const MGTwoLevelTransferNonNested<
dim,
LinearAlgebra::distributed::Vector<Number>> *>(
this->transfer[this->transfer.max_level()].get()))
{
dof_handler_in = t->dof_handler_fine;
level_in = t->mg_level_fine;
}
else
{
Assert(false, ExcNotImplemented());
}

if (this->perform_plain_copy)
{
// global-coarsening path: compare indices of cells

if (dof_handler_in == &dof_handler_out)
return; // nothing to do

std::vector<types::global_dof_index> dof_indices_in;
std::vector<types::global_dof_index> dof_indices_out;

internal::loop_over_active_or_level_cells(
dof_handler_in->get_triangulation(), level_in, [&](const auto &cell) {
const auto cell_id = cell->id();

Assert(
dof_handler_out.get_triangulation().contains_cell(cell_id),
ExcMessage(
"DoFHandler instanes used for set up of MGTransferMF and copy_to_mg(), "
"copy_from_mg(), or interpolate_to_mg() are not compatible."));

if (level_in == numbers::invalid_unsigned_int)
{
const auto cell_in =
cell->as_dof_handler_iterator(*dof_handler_in);
const auto cell_out =
dof_handler_out.get_triangulation()
.create_cell_iterator(cell_id)
->as_dof_handler_iterator(dof_handler_out);

AssertDimension(cell_in->get_fe().n_dofs_per_cell(),
cell_out->get_fe().n_dofs_per_cell());

dof_indices_in.resize(cell_in->get_fe().n_dofs_per_cell());
dof_indices_out.resize(cell_out->get_fe().n_dofs_per_cell());

cell_in->get_dof_indices(dof_indices_in);
cell_out->get_dof_indices(dof_indices_out);
}
else
{
const auto cell_in =
cell->as_dof_handler_level_iterator(*dof_handler_in);
const auto cell_out =
dof_handler_out.get_triangulation()
.create_cell_iterator(cell_id)
->as_dof_handler_level_iterator(dof_handler_out);

AssertDimension(cell_in->get_fe().n_dofs_per_cell(),
cell_out->get_fe().n_dofs_per_cell());

dof_indices_in.resize(cell_in->get_fe().n_dofs_per_cell());
dof_indices_out.resize(cell_out->get_fe().n_dofs_per_cell());

cell_in->get_mg_dof_indices(dof_indices_in);
cell_out->get_mg_dof_indices(dof_indices_out);
}

Assert(
dof_indices_in == dof_indices_out,
ExcMessage(
"DoFHandler instanes used for set up of MGTransferMF and copy_to_mg(), "
"copy_from_mg(), or interpolate_to_mg() are not compatible."));
});
}
else
{
// local-smoothing path: compare DoFHandlers

Assert(
dof_handler_in == &dof_handler_out,
ExcMessage(
"DoFHandler instances used during MGTransferMF::build() and copy_to_mg(), "
"copy_from_mg(), or interpolate_to_mg() are not the same."));

// TODO: criterion could be refined by checking the DoFs on all level that
// have active cells
}
#endif
}



template <int dim, typename Number>
std::size_t
MGTransferMF<dim, Number>::memory_consumption() const
Expand Down Expand Up @@ -4795,6 +4922,9 @@ MGTwoLevelTransferNonNested<dim, LinearAlgebra::distributed::Vector<Number>>::
dof_handler_fine.get_fe().n_components() > 0,
ExcNotImplemented());

this->dof_handler_fine = &dof_handler_fine;
this->mg_level_fine = numbers::invalid_unsigned_int;

this->fine_element_is_continuous =
dof_handler_fine.get_fe().n_dofs_per_vertex() > 0;

Expand Down

0 comments on commit e3cfd05

Please sign in to comment.