From: Daniel Arndt Date: Tue, 19 Jun 2018 10:26:13 +0000 (+0200) Subject: Fix back_interpolate overload resolution X-Git-Tag: v9.1.0-rc1~994^2~2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=25dbd7e1201fd68066634f4c509b0f80b952e6a8;p=dealii.git Fix back_interpolate overload resolution --- diff --git a/include/deal.II/fe/fe_tools_interpolate.templates.h b/include/deal.II/fe/fe_tools_interpolate.templates.h index 27ed8f0329..2bc7c462e7 100644 --- a/include/deal.II/fe/fe_tools_interpolate.templates.h +++ b/include/deal.II/fe/fe_tools_interpolate.templates.h @@ -368,7 +368,7 @@ namespace FETools namespace internal { template - void + typename std::enable_if::value>::type back_interpolate( const DoFHandler & dof1, const AffineConstraints &constraints1, @@ -412,6 +412,19 @@ namespace FETools u2 = u2_out; interpolate(dof2, u2, dof1, constraints1, u1_interpolated); } + + template + void + back_interpolate( + const DoFHandler &, + const AffineConstraints &, + const PETScWrappers::MPI::BlockVector &, + const DoFHandler &, + const AffineConstraints &, + PETScWrappers::MPI::BlockVector &) + { + Assert(false, ExcNotImplemented()); + } #endif // special version for Trilinos @@ -444,6 +457,50 @@ namespace FETools u2 = u2_out; interpolate(dof2, u2, dof1, constraints1, u1_interpolated); } + + template + void + back_interpolate( + const DoFHandler &dof1, + const AffineConstraints< + typename TrilinosWrappers::MPI::BlockVector::value_type> &constraints1, + const TrilinosWrappers::MPI::BlockVector & u1, + const DoFHandler & dof2, + const AffineConstraints< + typename TrilinosWrappers::MPI::BlockVector::value_type> &constraints2, + TrilinosWrappers::MPI::BlockVector &u1_interpolated) + { + if (u1.n_blocks() == 0) + return; + const MPI_Comm &mpi_communicator = u1.block(0).get_mpi_communicator(); + const IndexSet &dof2_locally_owned_dofs = dof2.locally_owned_dofs(); + IndexSet dof2_locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof2, dof2_locally_relevant_dofs); + + TrilinosWrappers::MPI::Vector u2_out(dof2_locally_owned_dofs, + mpi_communicator); + interpolate(dof1, u1, dof2, constraints2, u2_out); + TrilinosWrappers::MPI::Vector u2(dof2_locally_owned_dofs, + dof2_locally_relevant_dofs, + mpi_communicator); + u2 = u2_out; + interpolate(dof2, u2, dof1, constraints1, u1_interpolated); + } + + template + void + back_interpolate( + const DoFHandler &, + const AffineConstraints< + typename LinearAlgebra::EpetraWrappers::Vector::value_type> &, + const LinearAlgebra::EpetraWrappers::Vector &, + const DoFHandler &, + const AffineConstraints< + typename LinearAlgebra::EpetraWrappers::Vector::value_type> &, + LinearAlgebra::EpetraWrappers::Vector &) + { + AssertThrow(false, ExcNotImplemented()); + } #endif // special version for LinearAlgebra::distributed::Vector @@ -469,6 +526,19 @@ namespace FETools u2.update_ghost_values(); interpolate(dof2, u2, dof1, constraints1, u1_interpolated); } + + // special version for LinearAlgebra::distributed::BlockVector + template + void + back_interpolate(const DoFHandler &, + const AffineConstraints &, + const LinearAlgebra::distributed::BlockVector &, + const DoFHandler &, + const AffineConstraints &, + LinearAlgebra::distributed::BlockVector &) + { + AssertThrow(false, ExcNotImplemented()); + } } // namespace internal diff --git a/tests/mpi/parallel_vector_back_interpolate_02.cc b/tests/mpi/parallel_vector_back_interpolate_02.cc new file mode 100644 index 0000000000..185c0ebba2 --- /dev/null +++ b/tests/mpi/parallel_vector_back_interpolate_02.cc @@ -0,0 +1,164 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2018 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. +// +// --------------------------------------------------------------------- + + +// check FETools::back_interpolate on Trilinos vectors + +#include +#include +#include + +#include + +#include +#include +#include + +#include +#include +#include + +#include + +#include +#include +#include + +#include + +#include +#include + +#include "../tests.h" + + +void +test() +{ + const unsigned int dim = 2; + + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tria); + tria.refine_global(3); + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + + FESystem fe1(FE_Q(1), 1, FE_Q(1), 1); + FESystem fe2(FE_Q(2), 1, FE_Q(2), 1); + DoFHandler dof1(tria), dof2(tria); + dof1.distribute_dofs(fe1); + dof2.distribute_dofs(fe2); + + DoFRenumbering::component_wise(dof1); + DoFRenumbering::component_wise(dof2); + + IndexSet locally_owned_dofs1 = dof1.locally_owned_dofs(); + IndexSet locally_relevant_dofs1; + DoFTools::extract_locally_relevant_dofs(dof1, locally_relevant_dofs1); + + std::vector owned_partitioning1; + std::vector relevant_partitioning1; + { + std::vector dofs_per_block(2); + DoFTools::count_dofs_per_block(dof1, dofs_per_block); + const unsigned int n1 = dofs_per_block[0], n2 = dofs_per_block[1]; + + owned_partitioning1.push_back(locally_owned_dofs1.get_view(0, n1)); + owned_partitioning1.push_back(locally_owned_dofs1.get_view(n1, n1 + n2)); + + relevant_partitioning1.push_back(locally_relevant_dofs1.get_view(0, n1)); + relevant_partitioning1.push_back( + locally_relevant_dofs1.get_view(n1, n1 + n2)); + } + + AffineConstraints c1, c2; + DoFTools::make_hanging_node_constraints(dof1, c1); + c1.close(); + DoFTools::make_hanging_node_constraints(dof2, c2); + c2.close(); + + { + TrilinosWrappers::MPI::Vector v1_distributed(locally_owned_dofs1, + MPI_COMM_WORLD); + TrilinosWrappers::MPI::Vector v1(locally_owned_dofs1, + locally_relevant_dofs1, + MPI_COMM_WORLD); + TrilinosWrappers::MPI::Vector v1_interpolated(v1_distributed); + + for (const auto &el : locally_owned_dofs1) + v1_distributed(el) = random_value(); + c1.distribute(v1_distributed); + v1 = v1_distributed; + + FETools::back_interpolate(dof1, c1, v1, dof2, c2, v1_interpolated); + for (const auto &el : locally_owned_dofs1) + { + if (std::abs(v1_interpolated(el) - v1(el)) > 1.e-10) + { + std::cout << v1_interpolated(el) << " should be " << v1(el) + << std::endl; + AssertThrow(false, ExcInternalError()); + } + } + deallog << "TrilinosWrappers::MPI::Vector: OK" << std::endl; + } + { + TrilinosWrappers::MPI::BlockVector v1_distributed(owned_partitioning1, + MPI_COMM_WORLD); + TrilinosWrappers::MPI::BlockVector v1(owned_partitioning1, + relevant_partitioning1, + MPI_COMM_WORLD); + TrilinosWrappers::MPI::BlockVector v1_interpolated(v1_distributed); + + for (const auto &el : locally_owned_dofs1) + v1_distributed(el) = random_value(); + c1.distribute(v1_distributed); + v1 = v1_distributed; + + FETools::back_interpolate(dof1, c1, v1, dof2, c2, v1_interpolated); + for (const auto &el : locally_owned_dofs1) + { + if (std::abs(v1_interpolated(el) - v1(el)) > 1.e-10) + { + std::cout << v1_interpolated(el) << " should be " << v1(el) + << std::endl; + AssertThrow(false, ExcInternalError()); + } + } + deallog << "TrilinosWrappers::MPI::BlockVector: OK" << std::endl; + } +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + deallog.push(Utilities::int_to_string(myid)); + + if (myid == 0) + { + initlog(); + deallog << std::setprecision(4); + + test(); + } + else + test(); +} diff --git a/tests/mpi/parallel_vector_back_interpolate_02.with_trilinos=true.mpirun=4.output b/tests/mpi/parallel_vector_back_interpolate_02.with_trilinos=true.mpirun=4.output new file mode 100644 index 0000000000..0becb58586 --- /dev/null +++ b/tests/mpi/parallel_vector_back_interpolate_02.with_trilinos=true.mpirun=4.output @@ -0,0 +1,3 @@ + +DEAL:0::TrilinosWrappers::MPI::Vector: OK +DEAL:0::TrilinosWrappers::MPI::BlockVector: OK