From 08c70c40ea2ebc01ddadff1399f3f667643aba1e Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Wed, 11 Jul 2018 18:54:19 +0200 Subject: [PATCH] Allow iterating over TrilinosSparsityPattern in parallel and after compress --- .../deal.II/lac/trilinos_sparsity_pattern.h | 23 +++-- source/lac/trilinos_sparsity_pattern.cc | 20 ++-- tests/trilinos/sparsity_pattern_06.cc | 94 +++++++++++++++++++ ...rilinos=true.with_mpi=true.mpirun=3.output | 48 ++++++++++ 4 files changed, 168 insertions(+), 17 deletions(-) create mode 100644 tests/trilinos/sparsity_pattern_06.cc create mode 100644 tests/trilinos/sparsity_pattern_06.with_trilinos=true.with_mpi=true.mpirun=3.output diff --git a/include/deal.II/lac/trilinos_sparsity_pattern.h b/include/deal.II/lac/trilinos_sparsity_pattern.h index a61b358c21..e2bf25e04d 100644 --- a/include/deal.II/lac/trilinos_sparsity_pattern.h +++ b/include/deal.II/lac/trilinos_sparsity_pattern.h @@ -1288,6 +1288,7 @@ namespace TrilinosWrappers {} + inline Iterator::Iterator(const Iterator &) = default; @@ -1300,18 +1301,23 @@ namespace TrilinosWrappers ++accessor.a_index; - // If at end of line: do one - // step, then cycle until we - // find a row with a nonzero - // number of entries. + // If at end of line: do one step, then cycle until we find a row with a + // nonzero number of entries that is stored locally. if (accessor.a_index >= accessor.colnum_cache->size()) { accessor.a_index = 0; ++accessor.a_row; - while ((accessor.a_row < accessor.sparsity_pattern->n_rows()) && - (accessor.sparsity_pattern->row_length(accessor.a_row) == 0)) - ++accessor.a_row; + while (accessor.a_row < accessor.sparsity_pattern->n_rows()) + { + const auto row_length = + accessor.sparsity_pattern->row_length(accessor.a_row); + if (row_length == 0 || + row_length == static_cast(-1)) + ++accessor.a_row; + else + break; + } accessor.visit_present_row(); } @@ -1376,7 +1382,8 @@ namespace TrilinosWrappers inline SparsityPattern::const_iterator SparsityPattern::begin() const { - return const_iterator(this, 0, 0); + const size_type first_valid_row = this->local_range().first; + return const_iterator(this, first_valid_row, 0); } diff --git a/source/lac/trilinos_sparsity_pattern.cc b/source/lac/trilinos_sparsity_pattern.cc index 6960fc94d7..5a57541144 100644 --- a/source/lac/trilinos_sparsity_pattern.cc +++ b/source/lac/trilinos_sparsity_pattern.cc @@ -41,12 +41,12 @@ namespace TrilinosWrappers if (this->a_row == sparsity_pattern->n_rows()) { colnum_cache.reset(); - return; } - // otherwise first flush Trilinos caches - sparsity_pattern->compress(); + // otherwise first flush Trilinos caches if necessary + if (!sparsity_pattern->is_compressed()) + sparsity_pattern->compress(); colnum_cache = std::make_shared>( sparsity_pattern->row_length(this->a_row)); @@ -815,14 +815,16 @@ namespace TrilinosWrappers ierr = graph->FillComplete(*column_space_map, static_cast( graph->RangeMap())); + AssertThrow(ierr == 0, ExcTrilinosError(ierr)); } else - ierr = graph->GlobalAssemble(*column_space_map, - static_cast( - graph->RangeMap()), - true); - - AssertThrow(ierr == 0, ExcTrilinosError(ierr)); + { + ierr = graph->GlobalAssemble(*column_space_map, + static_cast( + graph->RangeMap()), + true); + AssertThrow(ierr == 0, ExcTrilinosError(ierr)); + } ierr = graph->OptimizeStorage(); AssertThrow(ierr == 0, ExcTrilinosError(ierr)); diff --git a/tests/trilinos/sparsity_pattern_06.cc b/tests/trilinos/sparsity_pattern_06.cc new file mode 100644 index 0000000000..d4200b25ce --- /dev/null +++ b/tests/trilinos/sparsity_pattern_06.cc @@ -0,0 +1,94 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// Test accessing the entries of TrilinosWrappers::SparsityPattern +// after compress has been called. +// The sparsity pattern used is the same as in sparse_matrix_add_03. + +#include +#include +#include + +#include + +#include "../tests.h" + + +void +test() +{ + const unsigned int MyPID = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + const unsigned int NumProc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + if (MyPID == 0) + deallog << "NumProc=" << NumProc << std::endl; + + // create non-contiguous index set for NumProc > 3 + dealii::IndexSet parallel_partitioning(NumProc * 2); + + // non-contiguous + parallel_partitioning.add_index(MyPID); + parallel_partitioning.add_index(NumProc + MyPID); + + // create sparsity pattern from parallel_partitioning + + // The sparsity pattern corresponds to a [FE_DGQ<1>(p=0)]^2 FESystem, + // on a periodic triangulation in which each MPI process owns 2 cells, + // with reordered dofs by its components, such that the rows in the + // final matrix are locally not in a contiguous set. + + dealii::TrilinosWrappers::SparsityPattern sp(parallel_partitioning, + MPI_COMM_WORLD, + 2); + + sp.add(MyPID, (NumProc + MyPID - 1) % NumProc); + sp.add(MyPID, MyPID); + sp.add(MyPID, (MyPID + 1) % NumProc); + sp.add(MyPID + NumProc, (NumProc + MyPID - 1) % NumProc + NumProc); + sp.add(MyPID + NumProc, MyPID + NumProc); + sp.add(MyPID + NumProc, (MyPID + 1) % NumProc + NumProc); + + deallog << "before compress:" << std::endl; + + for (const auto &el : sp) + { + deallog << "index: " << el.index() << " = " + << "(" << el.row() << " , " << el.column() << ")" << std::endl; + } + + sp.compress(); + + deallog << "after compress:" << std::endl; + for (const auto &el : sp) + { + deallog << "index: " << el.index() << " = " + << "(" << el.row() << " , " << el.column() << ")" << std::endl; + } +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + MPILogInitAll mpi_log; + + test(); +} diff --git a/tests/trilinos/sparsity_pattern_06.with_trilinos=true.with_mpi=true.mpirun=3.output b/tests/trilinos/sparsity_pattern_06.with_trilinos=true.with_mpi=true.mpirun=3.output new file mode 100644 index 0000000000..6c15a08999 --- /dev/null +++ b/tests/trilinos/sparsity_pattern_06.with_trilinos=true.with_mpi=true.mpirun=3.output @@ -0,0 +1,48 @@ + +DEAL:0::NumProc=3 +DEAL:0::before compress: +DEAL:0::index: 0 = (0 , 0) +DEAL:0::index: 1 = (0 , 1) +DEAL:0::index: 2 = (0 , 2) +DEAL:0::index: 0 = (3 , 3) +DEAL:0::index: 1 = (3 , 4) +DEAL:0::index: 2 = (3 , 5) +DEAL:0::after compress: +DEAL:0::index: 0 = (0 , 0) +DEAL:0::index: 1 = (0 , 1) +DEAL:0::index: 2 = (0 , 2) +DEAL:0::index: 0 = (3 , 3) +DEAL:0::index: 1 = (3 , 4) +DEAL:0::index: 2 = (3 , 5) + +DEAL:1::before compress: +DEAL:1::index: 0 = (1 , 1) +DEAL:1::index: 1 = (1 , 0) +DEAL:1::index: 2 = (1 , 2) +DEAL:1::index: 0 = (4 , 4) +DEAL:1::index: 1 = (4 , 3) +DEAL:1::index: 2 = (4 , 5) +DEAL:1::after compress: +DEAL:1::index: 0 = (1 , 1) +DEAL:1::index: 1 = (1 , 0) +DEAL:1::index: 2 = (1 , 2) +DEAL:1::index: 0 = (4 , 4) +DEAL:1::index: 1 = (4 , 3) +DEAL:1::index: 2 = (4 , 5) + + +DEAL:2::before compress: +DEAL:2::index: 0 = (2 , 2) +DEAL:2::index: 1 = (2 , 0) +DEAL:2::index: 2 = (2 , 1) +DEAL:2::index: 0 = (5 , 5) +DEAL:2::index: 1 = (5 , 3) +DEAL:2::index: 2 = (5 , 4) +DEAL:2::after compress: +DEAL:2::index: 0 = (2 , 2) +DEAL:2::index: 1 = (2 , 0) +DEAL:2::index: 2 = (2 , 1) +DEAL:2::index: 0 = (5 , 5) +DEAL:2::index: 1 = (5 , 3) +DEAL:2::index: 2 = (5 , 4) + -- 2.39.5