From: Svenja Schoeder Date: Wed, 20 Jun 2018 13:26:08 +0000 (+0200) Subject: improved the description for the handling of initialized ghost elements X-Git-Tag: v9.1.0-rc1~945^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=48b19b386699c0fd2151ff867c90810a047c3176;p=dealii.git improved the description for the handling of initialized ghost elements extended test to show the behavior in the presence of initialized and set ghost elements added test for ReadWriteVector discussed the case that compress is called with VectorOperation::max two times and added test coverage of this case indention --- diff --git a/include/deal.II/lac/la_parallel_vector.h b/include/deal.II/lac/la_parallel_vector.h index 0a7fbc9dfd..97c345574d 100644 --- a/include/deal.II/lac/la_parallel_vector.h +++ b/include/deal.II/lac/la_parallel_vector.h @@ -417,12 +417,18 @@ namespace LinearAlgebra * ghost element is found, it is compared to the value on the owning * processor and an exception is thrown if these elements do not agree. * If called with VectorOperation::min or VectorOperation::max the - * minimum or maximum on all elements across the processors is set. If - * a processor does not set values for the ghosts, the value 0.0 is - * assumed, which can be the minimal or maximal value across all - * processors. The operations min and max are reasonable if all - * processors set relevant values or set values that do not alter the - * results, e.g. small values for the max operation. + * minimum or maximum on all elements across the processors is set. + * @note This vector class has a fixed set of ghost entries attached to + * the local representation. As a consequence, all ghost entries are + * assumed to be valid and will be unconditionally exchanged according + * to the given VectorOperation. Make sure to initialize all ghost + * entries with the neutral element of the given VectorOperation. The + * neutral element is zero for VectorOperation::add and + * VectorOperation::insert, `+inf` for VectorOperation::min, and `-inf` + * for VectorOperation::max or touch all ghost entries. If all values + * are initialized with values below zero and compress is called with + * VectorOperation::max two times subsequently, the maximal value + * after the second calculation will be zero. */ virtual void compress(::dealii::VectorOperation::values operation) override; diff --git a/tests/mpi/parallel_vector_22.cc b/tests/mpi/parallel_vector_22.cc index ecae1101db..d01784514f 100644 --- a/tests/mpi/parallel_vector_22.cc +++ b/tests/mpi/parallel_vector_22.cc @@ -62,7 +62,7 @@ test() deallog << myid << ":" << "second owned entry: " << v(myid * 2 + 1) << std::endl; - // set ghost dof on not owning processor and maximize + // set ghost dof on owning processor and maximize if (myid) v(1) = 7. * myid; v.compress(VectorOperation::max); @@ -72,7 +72,7 @@ test() // check deallog << myid << ":" - << "ghost entry: " << v(1) << std::endl; + << "ghost entry after max from owner: " << v(1) << std::endl; // ghosts are set to zero v.zero_out_ghosts(); @@ -83,9 +83,9 @@ test() // check deallog << myid << ":" - << "ghost entry: " << v(1) << std::endl; + << "ghost entry after min from zero: " << v(1) << std::endl; - // update of ghost value from owner and minimize + // set ghost dof on non-owning processors and minimize v.zero_out_ghosts(); if (!myid) v(1) = -1.; @@ -94,7 +94,44 @@ test() // check deallog << myid << ":" - << "ghost entry: " << v(1) << std::endl; + << "ghost entry after min from : " << v(1) << std::endl; + + // set vector to 1, zeros in ghosts except on owner where -1. is set + v.zero_out_ghosts(); + v = 1.0; + if (!myid) + v(1) = -1.0; + + // maximize + v.compress(VectorOperation::max); + v.update_ghost_values(); + + // even if only one value is set (-1. on owner), the other values + // contribute a "0" and maximization receives zero and returns it + deallog << myid << ":" + << "ghost entry after max and partly init: " << v(1) << std::endl; + + // however, if the ghost value is set on all processors, the + // maximum is -1: + v.zero_out_ghosts(); + v = 1.0; + v(1) = -1.0; + v.compress(VectorOperation::max); + v.update_ghost_values(); + deallog << myid << ":" + << "ghost entry after max and full init: " << v(1) << std::endl; + + // what happens in case max is called two times and all values were smaller + // than zero + v.zero_out_ghosts(); + v = -1.0; + v(1) = -1.0; + v.compress(VectorOperation::max); + deallog << myid << ":" + << "ghost entry after first max: " << v(1) << std::endl; + v.compress(VectorOperation::max); + deallog << myid << ":" + << "ghost entry after second max: " << v(1) << std::endl; if (myid == 0) deallog << "OK" << std::endl; diff --git a/tests/mpi/parallel_vector_22.mpirun=10.output b/tests/mpi/parallel_vector_22.mpirun=10.output index 8a13c7f82a..d8c6b01b01 100644 --- a/tests/mpi/parallel_vector_22.mpirun=10.output +++ b/tests/mpi/parallel_vector_22.mpirun=10.output @@ -2,70 +2,110 @@ DEAL:0::numproc=10 DEAL:0::0:first owned entry: 0.00000 DEAL:0::0:second owned entry: 2.00000 -DEAL:0::0:ghost entry: 63.0000 -DEAL:0::0:ghost entry: 0.00000 -DEAL:0::0:ghost entry: -1.00000 +DEAL:0::0:ghost entry after max from owner: 63.0000 +DEAL:0::0:ghost entry after min from zero: 0.00000 +DEAL:0::0:ghost entry after min from : -1.00000 +DEAL:0::0:ghost entry after max and partly init: 0.00000 +DEAL:0::0:ghost entry after max and full init: -1.00000 +DEAL:0::0:ghost entry after first max: -1.00000 +DEAL:0::0:ghost entry after second max: 0.00000 DEAL:0::OK DEAL:1::1:first owned entry: 4.00000 DEAL:1::1:second owned entry: 6.00000 -DEAL:1::1:ghost entry: 63.0000 -DEAL:1::1:ghost entry: 0.00000 -DEAL:1::1:ghost entry: -1.00000 +DEAL:1::1:ghost entry after max from owner: 63.0000 +DEAL:1::1:ghost entry after min from zero: 0.00000 +DEAL:1::1:ghost entry after min from : -1.00000 +DEAL:1::1:ghost entry after max and partly init: 0.00000 +DEAL:1::1:ghost entry after max and full init: -1.00000 +DEAL:1::1:ghost entry after first max: 0.00000 +DEAL:1::1:ghost entry after second max: 0.00000 DEAL:2::2:first owned entry: 8.00000 DEAL:2::2:second owned entry: 10.0000 -DEAL:2::2:ghost entry: 63.0000 -DEAL:2::2:ghost entry: 0.00000 -DEAL:2::2:ghost entry: -1.00000 +DEAL:2::2:ghost entry after max from owner: 63.0000 +DEAL:2::2:ghost entry after min from zero: 0.00000 +DEAL:2::2:ghost entry after min from : -1.00000 +DEAL:2::2:ghost entry after max and partly init: 0.00000 +DEAL:2::2:ghost entry after max and full init: -1.00000 +DEAL:2::2:ghost entry after first max: 0.00000 +DEAL:2::2:ghost entry after second max: 0.00000 DEAL:3::3:first owned entry: 12.0000 DEAL:3::3:second owned entry: 14.0000 -DEAL:3::3:ghost entry: 63.0000 -DEAL:3::3:ghost entry: 0.00000 -DEAL:3::3:ghost entry: -1.00000 +DEAL:3::3:ghost entry after max from owner: 63.0000 +DEAL:3::3:ghost entry after min from zero: 0.00000 +DEAL:3::3:ghost entry after min from : -1.00000 +DEAL:3::3:ghost entry after max and partly init: 0.00000 +DEAL:3::3:ghost entry after max and full init: -1.00000 +DEAL:3::3:ghost entry after first max: 0.00000 +DEAL:3::3:ghost entry after second max: 0.00000 DEAL:4::4:first owned entry: 16.0000 DEAL:4::4:second owned entry: 18.0000 -DEAL:4::4:ghost entry: 63.0000 -DEAL:4::4:ghost entry: 0.00000 -DEAL:4::4:ghost entry: -1.00000 +DEAL:4::4:ghost entry after max from owner: 63.0000 +DEAL:4::4:ghost entry after min from zero: 0.00000 +DEAL:4::4:ghost entry after min from : -1.00000 +DEAL:4::4:ghost entry after max and partly init: 0.00000 +DEAL:4::4:ghost entry after max and full init: -1.00000 +DEAL:4::4:ghost entry after first max: 0.00000 +DEAL:4::4:ghost entry after second max: 0.00000 DEAL:5::5:first owned entry: 20.0000 DEAL:5::5:second owned entry: 22.0000 -DEAL:5::5:ghost entry: 63.0000 -DEAL:5::5:ghost entry: 0.00000 -DEAL:5::5:ghost entry: -1.00000 +DEAL:5::5:ghost entry after max from owner: 63.0000 +DEAL:5::5:ghost entry after min from zero: 0.00000 +DEAL:5::5:ghost entry after min from : -1.00000 +DEAL:5::5:ghost entry after max and partly init: 0.00000 +DEAL:5::5:ghost entry after max and full init: -1.00000 +DEAL:5::5:ghost entry after first max: 0.00000 +DEAL:5::5:ghost entry after second max: 0.00000 DEAL:6::6:first owned entry: 24.0000 DEAL:6::6:second owned entry: 26.0000 -DEAL:6::6:ghost entry: 63.0000 -DEAL:6::6:ghost entry: 0.00000 -DEAL:6::6:ghost entry: -1.00000 +DEAL:6::6:ghost entry after max from owner: 63.0000 +DEAL:6::6:ghost entry after min from zero: 0.00000 +DEAL:6::6:ghost entry after min from : -1.00000 +DEAL:6::6:ghost entry after max and partly init: 0.00000 +DEAL:6::6:ghost entry after max and full init: -1.00000 +DEAL:6::6:ghost entry after first max: 0.00000 +DEAL:6::6:ghost entry after second max: 0.00000 DEAL:7::7:first owned entry: 28.0000 DEAL:7::7:second owned entry: 30.0000 -DEAL:7::7:ghost entry: 63.0000 -DEAL:7::7:ghost entry: 0.00000 -DEAL:7::7:ghost entry: -1.00000 +DEAL:7::7:ghost entry after max from owner: 63.0000 +DEAL:7::7:ghost entry after min from zero: 0.00000 +DEAL:7::7:ghost entry after min from : -1.00000 +DEAL:7::7:ghost entry after max and partly init: 0.00000 +DEAL:7::7:ghost entry after max and full init: -1.00000 +DEAL:7::7:ghost entry after first max: 0.00000 +DEAL:7::7:ghost entry after second max: 0.00000 DEAL:8::8:first owned entry: 32.0000 DEAL:8::8:second owned entry: 34.0000 -DEAL:8::8:ghost entry: 63.0000 -DEAL:8::8:ghost entry: 0.00000 -DEAL:8::8:ghost entry: -1.00000 +DEAL:8::8:ghost entry after max from owner: 63.0000 +DEAL:8::8:ghost entry after min from zero: 0.00000 +DEAL:8::8:ghost entry after min from : -1.00000 +DEAL:8::8:ghost entry after max and partly init: 0.00000 +DEAL:8::8:ghost entry after max and full init: -1.00000 +DEAL:8::8:ghost entry after first max: 0.00000 +DEAL:8::8:ghost entry after second max: 0.00000 DEAL:9::9:first owned entry: 36.0000 DEAL:9::9:second owned entry: 38.0000 -DEAL:9::9:ghost entry: 63.0000 -DEAL:9::9:ghost entry: 0.00000 -DEAL:9::9:ghost entry: -1.00000 +DEAL:9::9:ghost entry after max from owner: 63.0000 +DEAL:9::9:ghost entry after min from zero: 0.00000 +DEAL:9::9:ghost entry after min from : -1.00000 +DEAL:9::9:ghost entry after max and partly init: 0.00000 +DEAL:9::9:ghost entry after max and full init: -1.00000 +DEAL:9::9:ghost entry after first max: 0.00000 +DEAL:9::9:ghost entry after second max: 0.00000 diff --git a/tests/mpi/parallel_vector_22.mpirun=4.output b/tests/mpi/parallel_vector_22.mpirun=4.output index 4c36987598..f1940df5e2 100644 --- a/tests/mpi/parallel_vector_22.mpirun=4.output +++ b/tests/mpi/parallel_vector_22.mpirun=4.output @@ -2,28 +2,44 @@ DEAL:0::numproc=4 DEAL:0::0:first owned entry: 0.00000 DEAL:0::0:second owned entry: 2.00000 -DEAL:0::0:ghost entry: 21.0000 -DEAL:0::0:ghost entry: 0.00000 -DEAL:0::0:ghost entry: -1.00000 +DEAL:0::0:ghost entry after max from owner: 21.0000 +DEAL:0::0:ghost entry after min from zero: 0.00000 +DEAL:0::0:ghost entry after min from : -1.00000 +DEAL:0::0:ghost entry after max and partly init: 0.00000 +DEAL:0::0:ghost entry after max and full init: -1.00000 +DEAL:0::0:ghost entry after first max: -1.00000 +DEAL:0::0:ghost entry after second max: 0.00000 DEAL:0::OK DEAL:1::1:first owned entry: 4.00000 DEAL:1::1:second owned entry: 6.00000 -DEAL:1::1:ghost entry: 21.0000 -DEAL:1::1:ghost entry: 0.00000 -DEAL:1::1:ghost entry: -1.00000 +DEAL:1::1:ghost entry after max from owner: 21.0000 +DEAL:1::1:ghost entry after min from zero: 0.00000 +DEAL:1::1:ghost entry after min from : -1.00000 +DEAL:1::1:ghost entry after max and partly init: 0.00000 +DEAL:1::1:ghost entry after max and full init: -1.00000 +DEAL:1::1:ghost entry after first max: 0.00000 +DEAL:1::1:ghost entry after second max: 0.00000 DEAL:2::2:first owned entry: 8.00000 DEAL:2::2:second owned entry: 10.0000 -DEAL:2::2:ghost entry: 21.0000 -DEAL:2::2:ghost entry: 0.00000 -DEAL:2::2:ghost entry: -1.00000 +DEAL:2::2:ghost entry after max from owner: 21.0000 +DEAL:2::2:ghost entry after min from zero: 0.00000 +DEAL:2::2:ghost entry after min from : -1.00000 +DEAL:2::2:ghost entry after max and partly init: 0.00000 +DEAL:2::2:ghost entry after max and full init: -1.00000 +DEAL:2::2:ghost entry after first max: 0.00000 +DEAL:2::2:ghost entry after second max: 0.00000 DEAL:3::3:first owned entry: 12.0000 DEAL:3::3:second owned entry: 14.0000 -DEAL:3::3:ghost entry: 21.0000 -DEAL:3::3:ghost entry: 0.00000 -DEAL:3::3:ghost entry: -1.00000 +DEAL:3::3:ghost entry after max from owner: 21.0000 +DEAL:3::3:ghost entry after min from zero: 0.00000 +DEAL:3::3:ghost entry after min from : -1.00000 +DEAL:3::3:ghost entry after max and partly init: 0.00000 +DEAL:3::3:ghost entry after max and full init: -1.00000 +DEAL:3::3:ghost entry after first max: 0.00000 +DEAL:3::3:ghost entry after second max: 0.00000 diff --git a/tests/mpi/parallel_vector_23.cc b/tests/mpi/parallel_vector_23.cc new file mode 100644 index 0000000000..ae79e2112b --- /dev/null +++ b/tests/mpi/parallel_vector_23.cc @@ -0,0 +1,92 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + +// check LA::Vector::compress(VectorOperation::min/max) from ghosts + +#include +#include + +#include +#include + +#include +#include + +#include "../tests.h" + + +void +test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + unsigned int numproc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + if (myid == 0) + deallog << "numproc=" << numproc << std::endl; + + + // each processor owns 2 indices and all + // are ghosting element 1 (the second) + IndexSet local_owned(numproc * 2); + local_owned.add_range(myid * 2, myid * 2 + 2); + IndexSet local_relevant(numproc * 2); + local_relevant = local_owned; + local_relevant.add_range(1, 2); + + // create vector + LinearAlgebra::distributed::Vector v(local_owned, + local_relevant, + MPI_COMM_WORLD); + + // the read write vector additionally has ghost elements + IndexSet read_write_owned(numproc * 2); + LinearAlgebra::ReadWriteVector read_write_vector(local_relevant); + + read_write_vector.local_element(0) = myid; + read_write_vector(1) = 2. * myid; + + v.import(read_write_vector, VectorOperation::max); + v.update_ghost_values(); + + deallog << myid << ":" + << "ghost entry after max: " << v(1) << std::endl; + + if (!myid) + read_write_vector(1) = -1.0; + + v.import(read_write_vector, VectorOperation::min); + v.update_ghost_values(); + + deallog << myid << ":" + << "ghost entry after min: " << v(1) << std::endl; + + + if (myid == 0) + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + MPILogInitAll log; + + test(); +} diff --git a/tests/mpi/parallel_vector_23.mpirun=4.output b/tests/mpi/parallel_vector_23.mpirun=4.output new file mode 100644 index 0000000000..d1ab3a1164 --- /dev/null +++ b/tests/mpi/parallel_vector_23.mpirun=4.output @@ -0,0 +1,17 @@ + +DEAL:0::numproc=4 +DEAL:0::0:ghost entry after max: 6.00000 +DEAL:0::0:ghost entry after min: -1.00000 +DEAL:0::OK + +DEAL:1::1:ghost entry after max: 6.00000 +DEAL:1::1:ghost entry after min: -1.00000 + + +DEAL:2::2:ghost entry after max: 6.00000 +DEAL:2::2:ghost entry after min: -1.00000 + + +DEAL:3::3:ghost entry after max: 6.00000 +DEAL:3::3:ghost entry after min: -1.00000 +