From 62c045b9e7135d8e56b1a8130a1d0dc7fb0472e9 Mon Sep 17 00:00:00 2001 From: jthano Date: Sun, 3 Jun 2018 22:56:00 -0700 Subject: [PATCH] Modified reinit function to check if parallel_partitioner is_ascending_and_one_to_one and then pass this result to make_trilinos_map. The reason for this check is that if the check returns true, then make_trilinos_map may be able to create a linear epetra map. added draft test Added test. Test passes Modified reinit function to check if parallel_partitioner is_ascending_and_one_to_one and then pass this result to make_trilinos_map. The reason for this check is that if the check returns true, then make_trilinos_map may be able to create a linear epetra map. Also added test and test output. Fixed issues related to test and added change log Renamed test files so that mpi test will happen. Also made minor edits to test file --- doc/news/changes/minor/20180616JohsuaHanophy | 5 ++ source/lac/trilinos_vector.cc | 5 +- tests/trilinos/vector_reinit_to_linear_map.cc | 80 +++++++++++++++++++ ...ector_reinit_to_linear_map.mpirun=2.output | 2 + 4 files changed, 91 insertions(+), 1 deletion(-) create mode 100644 doc/news/changes/minor/20180616JohsuaHanophy create mode 100644 tests/trilinos/vector_reinit_to_linear_map.cc create mode 100644 tests/trilinos/vector_reinit_to_linear_map.mpirun=2.output diff --git a/doc/news/changes/minor/20180616JohsuaHanophy b/doc/news/changes/minor/20180616JohsuaHanophy new file mode 100644 index 0000000000..373297774a --- /dev/null +++ b/doc/news/changes/minor/20180616JohsuaHanophy @@ -0,0 +1,5 @@ +Fixed: TrilinosWrappers::MPI::Vector::reinit check now checks if +parallel_partitioner is_ascending_and_one_to_one before calling +make_trilinos_map because make_trilinos_map may be able to make +a linear map if is_ascending_and_one_to_one is true. +(Joshua Hanophy, 2018/06/16) diff --git a/source/lac/trilinos_vector.cc b/source/lac/trilinos_vector.cc index 3ad0b5472f..f2acfcec4b 100644 --- a/source/lac/trilinos_vector.cc +++ b/source/lac/trilinos_vector.cc @@ -157,8 +157,11 @@ namespace TrilinosWrappers { nonlocal_vector.reset(); + const bool overlapping = + !parallel_partitioner.is_ascending_and_one_to_one(communicator); + Epetra_Map map = - parallel_partitioner.make_trilinos_map(communicator, true); + parallel_partitioner.make_trilinos_map(communicator, overlapping); vector = std_cxx14::make_unique(map); diff --git a/tests/trilinos/vector_reinit_to_linear_map.cc b/tests/trilinos/vector_reinit_to_linear_map.cc new file mode 100644 index 0000000000..b1964d455f --- /dev/null +++ b/tests/trilinos/vector_reinit_to_linear_map.cc @@ -0,0 +1,80 @@ +// --------------------------------------------------------------------- +// +// 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// run TrilinosWrappers::MPI::Vector::reinit() in MPI environment +// and with an IndexSet that is_ascending_and_one_to_one to check +// that make_trilinos_map will construct the object with a Map that +// is linear + +#include + +#include +#include + +#include +#include + +#include + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +void +test() +{ + MPI_Comm mpi_communicator(MPI_COMM_WORLD); + + parallel::distributed::Triangulation<2> tria(mpi_communicator); + + GridGenerator::hyper_cube(tria, -1, 0); + tria.refine_global(2); + + const unsigned int poly_degree = 1; + FE_Q<2> fe(poly_degree); + + DoFHandler<2> dof_handler(tria); + dof_handler.distribute_dofs(fe); + + IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs(); + + TrilinosWrappers::MPI::Vector vector_linear; + vector_linear.reinit(locally_owned_dofs, mpi_communicator); + + const int is_linear_index_set = + locally_owned_dofs.is_ascending_and_one_to_one(mpi_communicator); + + const int is_linear_map = vector_linear.trilinos_vector().Map().LinearMap(); + + if (is_linear_index_set == 1 && is_linear_map == 1) + { + deallog << "OK" << std::endl; + } +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + mpi_initlog(); + + test(); +} diff --git a/tests/trilinos/vector_reinit_to_linear_map.mpirun=2.output b/tests/trilinos/vector_reinit_to_linear_map.mpirun=2.output new file mode 100644 index 0000000000..0fd8fc12f0 --- /dev/null +++ b/tests/trilinos/vector_reinit_to_linear_map.mpirun=2.output @@ -0,0 +1,2 @@ + +DEAL::OK -- 2.39.5