Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix SolutionTransfer to properly compress() vectors. #15577

Merged
merged 3 commits into from
Jul 2, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
8 changes: 8 additions & 0 deletions doc/news/changes/minor/20230701Bangerth
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
Fixed: The SolutionTransfer class writes into output vectors, but does
not call compress() on them. This is of no consequence for deal.II
vectors for which this class is mostly used (in contrast to the
parallel::distributed::SolutionTransfer class), but leads to awkward
downstream failures with, for example, PETSc vectors. This is now
fixed.
<br>
(Wolfgang Bangerth, 2023/07/01)
5 changes: 5 additions & 0 deletions source/numerics/solution_transfer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -583,6 +583,11 @@ SolutionTransfer<dim, VectorType, spacedim>::interpolate(
Assert(false, ExcInternalError());
}
}

// We have written into the output vectors. If this was a PETSc vector, for
// example, then we need to compress these to make future operations safe:
for (auto &vec : all_out)
vec.compress(VectorOperation::insert);
}


Expand Down
100 changes: 100 additions & 0 deletions tests/petsc/solution_transfer_01.cc
Original file line number Diff line number Diff line change
@@ -0,0 +1,100 @@
// ---------------------------------------------------------------------
//
// Copyright (C) 2023 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------



// SolutionTransfer writes into vectors, but used to forget to
// compress the resulting vectors appropriately. This didn't matter
// most of the time because the (sequential) SolutionTransfer was
// apparently used only with our own vector types, for which
// compress() is a no-op. But it failed with PETSc vectors.
//
// This variation of the function deals with the
// SolutionTransfer::interpolate() variant that takes multiple vectors
// to transfer.


#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/dofs/dof_handler.h>

#include <deal.II/fe/fe_q.h>

#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/tria.h>

#include <deal.II/lac/petsc_vector.h>

#include <deal.II/numerics/solution_transfer.h>

#include <iostream>

#include "../tests.h"


template <int dim>
void
test()
{
Triangulation<dim> triangulation;

GridGenerator::hyper_cube(triangulation, 0, 1);
triangulation.refine_global(2);

FE_Q<dim> fe(1);
DoFHandler<dim> dof_handler(triangulation);
dof_handler.distribute_dofs(fe);


PETScWrappers::MPI::Vector solution;
solution.reinit(dof_handler.locally_owned_dofs(), MPI_COMM_SELF);

// Do a fake transfer
SolutionTransfer<dim, PETScWrappers::MPI::Vector> solution_trans(dof_handler);

PETScWrappers::MPI::Vector previous_solution;
previous_solution = solution;
triangulation.prepare_coarsening_and_refinement();
solution_trans.prepare_for_coarsening_and_refinement(previous_solution);

triangulation.execute_coarsening_and_refinement();

solution_trans.interpolate(previous_solution, solution);

deallog << "OK" << std::endl;
}



int
main(int argc, char **argv)
{
Utilities::MPI::MPI_InitFinalize mpi_initialization(
argc, argv, testing_max_num_threads());

mpi_initlog();
deallog << std::setprecision(2);

{
deallog.push("2d");
test<2>();
deallog.pop();
}

{
deallog.push("3d");
test<3>();
deallog.pop();
}
}
3 changes: 3 additions & 0 deletions tests/petsc/solution_transfer_01.output
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@

DEAL:2d::OK
DEAL:3d::OK
100 changes: 100 additions & 0 deletions tests/petsc/solution_transfer_02.cc
Original file line number Diff line number Diff line change
@@ -0,0 +1,100 @@
// ---------------------------------------------------------------------
//
// Copyright (C) 2023 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------



// SolutionTransfer writes into vectors, but used to forget to
// compress the resulting vectors appropriately. This didn't matter
// most of the time because the (sequential) SolutionTransfer was
// apparently used only with our own vector types, for which
// compress() is a no-op. But it failed with PETSc vectors.
//
// This variation of the function deals with the
// SolutionTransfer::interpolate() variant that takes multiple vectors
// to transfer.


#include <deal.II/dofs/dof_accessor.h>
#include <deal.II/dofs/dof_handler.h>

#include <deal.II/fe/fe_q.h>

#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/tria.h>

#include <deal.II/lac/petsc_vector.h>

#include <deal.II/numerics/solution_transfer.h>

#include <iostream>

#include "../tests.h"


template <int dim>
void
test()
{
Triangulation<dim> triangulation;

GridGenerator::hyper_cube(triangulation, 0, 1);
triangulation.refine_global(2);

FE_Q<dim> fe(1);
DoFHandler<dim> dof_handler(triangulation);
dof_handler.distribute_dofs(fe);


PETScWrappers::MPI::Vector solution;
solution.reinit(dof_handler.locally_owned_dofs(), MPI_COMM_SELF);

// Do a fake transfer
SolutionTransfer<dim, PETScWrappers::MPI::Vector> solution_trans(dof_handler);

PETScWrappers::MPI::Vector previous_solution;
previous_solution = solution;
triangulation.prepare_coarsening_and_refinement();
solution_trans.prepare_for_coarsening_and_refinement(previous_solution);

triangulation.execute_coarsening_and_refinement();

solution_trans.interpolate(previous_solution, solution);

deallog << "OK" << std::endl;
}



int
main(int argc, char **argv)
{
Utilities::MPI::MPI_InitFinalize mpi_initialization(
argc, argv, testing_max_num_threads());

mpi_initlog();
deallog << std::setprecision(2);

{
deallog.push("2d");
test<2>();
deallog.pop();
}

{
deallog.push("3d");
test<3>();
deallog.pop();
}
}
3 changes: 3 additions & 0 deletions tests/petsc/solution_transfer_02.output
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@

DEAL:2d::OK
DEAL:3d::OK