From: Timo Heister <timo.heister@gmail.com> Date: Tue, 4 Jul 2017 07:53:40 +0000 (+0200) Subject: Fix TrilinosWrappers::PreconditionBlock* with no local rows X-Git-Tag: v9.0.0-rc1~1453^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0232b200263873f6c422592a78128563730c084b;p=dealii.git Fix TrilinosWrappers::PreconditionBlock* with no local rows work around the ifpack error by pretending to use a point smoother on processors without any local rows. --- diff --git a/doc/news/changes/minor/20170704Heister b/doc/news/changes/minor/20170704Heister new file mode 100644 index 0000000000..7d12d2a4fc --- /dev/null +++ b/doc/news/changes/minor/20170704Heister @@ -0,0 +1,3 @@ +Fixed: TrilinosWrappers::PreconditionBlock* classes now work correctly if one of the processors does not own any rows. +<br> +(Timo Heister, 2017/07/04) diff --git a/source/lac/trilinos_precondition.cc b/source/lac/trilinos_precondition.cc index f0dfb4f38e..5172593a14 100644 --- a/source/lac/trilinos_precondition.cc +++ b/source/lac/trilinos_precondition.cc @@ -283,10 +283,15 @@ namespace TrilinosWrappers { // release memory before reallocation preconditioner.reset (); + + // Block relaxation setup fails if we have no locally owned rows. As a work-around + // we just pretend to use point relaxation on those processors: preconditioner.reset (Ifpack().Create - ("block relaxation", - const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), - 0)); + ( + (matrix.trilinos_matrix().NumMyRows()==0) ? + "point relaxation" : "block relaxation", + const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), + 0)); Ifpack_Preconditioner *ifpack = static_cast<Ifpack_Preconditioner *> (preconditioner.get()); @@ -343,10 +348,15 @@ namespace TrilinosWrappers const AdditionalData &additional_data) { preconditioner.reset (); + + // Block relaxation setup fails if we have no locally owned rows. As a work-around + // we just pretend to use point relaxation on those processors: preconditioner.reset (Ifpack().Create - ("block relaxation", - const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), - additional_data.overlap)); + ( + (matrix.trilinos_matrix().NumMyRows()==0) ? + "point relaxation" : "block relaxation", + const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), + additional_data.overlap)); Ifpack_Preconditioner *ifpack = static_cast<Ifpack_Preconditioner *> (preconditioner.get()); @@ -404,10 +414,15 @@ namespace TrilinosWrappers const AdditionalData &additional_data) { preconditioner.reset (); + + // Block relaxation setup fails if we have no locally owned rows. As a work-around + // we just pretend to use point relaxation on those processors: preconditioner.reset (Ifpack().Create - ("block relaxation", - const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), - additional_data.overlap)); + ( + (matrix.trilinos_matrix().NumMyRows()==0) ? + "point relaxation" : "block relaxation", + const_cast<Epetra_CrsMatrix *>(&matrix.trilinos_matrix()), + additional_data.overlap)); Ifpack_Preconditioner *ifpack = static_cast<Ifpack_Preconditioner *> (preconditioner.get()); diff --git a/tests/trilinos/block_relax_no_rows.cc b/tests/trilinos/block_relax_no_rows.cc new file mode 100644 index 0000000000..91d4f2c400 --- /dev/null +++ b/tests/trilinos/block_relax_no_rows.cc @@ -0,0 +1,77 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2013 - 2016 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. +// +// --------------------------------------------------------------------- + + +// Trilinos PreconditionBlock* used to fail if one processor has no locally +// owned rows + +#include "../tests.h" +#include <deal.II/lac/trilinos_sparse_matrix.h> +#include <deal.II/lac/trilinos_precondition.h> + + +template <class Prec> +void test() +{ + const unsigned int n_procs = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + + // only proc 0 and 1 own rows + IndexSet rows(4*n_procs); + if (myid==0) + rows.add_range(0, 2*n_procs); + else if (myid==1) + rows.add_range(2*n_procs, 4*n_procs); + rows.compress(); + + TrilinosWrappers::MPI::Vector src (rows, MPI_COMM_WORLD), + dst (rows, MPI_COMM_WORLD); + + TrilinosWrappers::SparseMatrix mat(rows, rows, MPI_COMM_WORLD); + for (const auto &row : rows) + { + const unsigned int i = row; + mat.set(i, i, 100.); + for (unsigned int j=0; j<mat.n(); ++j) + if (i!=j) + mat.set (i,j, i*j*.5+.5); + } + + for (unsigned int i=0; i<src.size(); ++i) + src(i) = i; + + mat.compress (VectorOperation::insert); + src.compress (VectorOperation::insert); + + Prec preconditioner; + typename Prec::AdditionalData data; + data.block_size = 4; + + preconditioner.initialize(mat, data); + preconditioner.vmult(dst, src); + + deallog << "dst: " << dst.l2_norm() << std::endl; +} + + +int main (int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + MPILogInitAll log; + + test<TrilinosWrappers::PreconditionBlockJacobi>(); + test<TrilinosWrappers::PreconditionBlockSOR>(); + test<TrilinosWrappers::PreconditionBlockSSOR>(); +} diff --git a/tests/trilinos/block_relax_no_rows.mpirun=3.output b/tests/trilinos/block_relax_no_rows.mpirun=3.output new file mode 100644 index 0000000000..beaca9e669 --- /dev/null +++ b/tests/trilinos/block_relax_no_rows.mpirun=3.output @@ -0,0 +1,14 @@ + +DEAL:0::dst: 0.137153 +DEAL:0::dst: 0.112530 +DEAL:0::dst: 0.0911212 + +DEAL:1::dst: 0.137153 +DEAL:1::dst: 0.112530 +DEAL:1::dst: 0.0911212 + + +DEAL:2::dst: 0.137153 +DEAL:2::dst: 0.112530 +DEAL:2::dst: 0.0911212 +