From 0014c1a973def7f47057a2abb0af58684e731cc9 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Fri, 17 Jan 2020 17:32:49 +0100 Subject: [PATCH] Add Utilities::MPI::NoncontiguousPartitioner --- doc/news/changes/minor/20200304PeterMunch | 7 + .../base/mpi_noncontiguous_partitioner.h | 171 +++++++++ .../mpi_noncontiguous_partitioner.templates.h | 328 ++++++++++++++++++ include/deal.II/base/mpi_tags.h | 2 + source/base/CMakeLists.txt | 2 + source/base/mpi_noncontiguous_partitioner.cc | 27 ++ .../mpi_noncontiguous_partitioner.inst.in | 52 +++ .../base/mpi_noncontiguous_partitioner_01.cc | 76 ++++ ...ncontiguous_partitioner_01.mpirun=2.output | 7 + .../base/mpi_noncontiguous_partitioner_02.cc | 161 +++++++++ ...ncontiguous_partitioner_02.mpirun=1.output | 13 + ...ncontiguous_partitioner_02.mpirun=2.output | 27 ++ ...ncontiguous_partitioner_02.mpirun=4.output | 55 +++ ...ncontiguous_partitioner_02.mpirun=5.output | 69 ++++ 14 files changed, 997 insertions(+) create mode 100644 doc/news/changes/minor/20200304PeterMunch create mode 100644 include/deal.II/base/mpi_noncontiguous_partitioner.h create mode 100644 include/deal.II/base/mpi_noncontiguous_partitioner.templates.h create mode 100644 source/base/mpi_noncontiguous_partitioner.cc create mode 100644 source/base/mpi_noncontiguous_partitioner.inst.in create mode 100644 tests/base/mpi_noncontiguous_partitioner_01.cc create mode 100644 tests/base/mpi_noncontiguous_partitioner_01.mpirun=2.output create mode 100644 tests/base/mpi_noncontiguous_partitioner_02.cc create mode 100644 tests/base/mpi_noncontiguous_partitioner_02.mpirun=1.output create mode 100644 tests/base/mpi_noncontiguous_partitioner_02.mpirun=2.output create mode 100644 tests/base/mpi_noncontiguous_partitioner_02.mpirun=4.output create mode 100644 tests/base/mpi_noncontiguous_partitioner_02.mpirun=5.output diff --git a/doc/news/changes/minor/20200304PeterMunch b/doc/news/changes/minor/20200304PeterMunch new file mode 100644 index 0000000000..92744807b3 --- /dev/null +++ b/doc/news/changes/minor/20200304PeterMunch @@ -0,0 +1,7 @@ +New: Introduce a new communication-pattern class Utilities::MPI::NoncontiguousPartitioner. +It is similar to Utilities::MPI::Partitioner, however, does not make any +restrictions regarding to the ordering of the underlying index sets. This class enables +efficient repartitioning of vectors, which might be benefitial in interfacing with +external libraries that expect a certain fixed order (like checkerboard partitioning). +
+(Peter Munch, 2020/03/04) diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.h b/include/deal.II/base/mpi_noncontiguous_partitioner.h new file mode 100644 index 0000000000..e04ac51190 --- /dev/null +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.h @@ -0,0 +1,171 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + +#ifndef dealii_mpi_noncontiguous_vector_h +#define dealii_mpi_noncontiguous_vector_h + +#include + +#include +#include +#include + +#include +#include + + +DEAL_II_NAMESPACE_OPEN + +namespace Utilities +{ + namespace MPI + { + /** + * A flexible Partitioner class, which does not impose an restrictions + * regarding the order of the underlying index sets. + * + * @author Peter Munch, 2020 + */ + template + class NoncontiguousPartitioner + : public LinearAlgebra::CommunicationPatternBase + { + public: + /** + * Constructor. Set up point-to-point communication pattern based on the + * IndexSets arguments @p indexset_has and @p indexset_want for the MPI + * communicator @p communicator. + */ + NoncontiguousPartitioner(const IndexSet &indexset_has, + const IndexSet &indexset_want, + const MPI_Comm &communicator); + + /** + * Constructor. Same as above but for vectors of indices @p indices_has + * and @p indices_want. This allows that the indices do not have to be + * sorted and the values are read and written automatically at the right + * position of the vector during update_values(), update_values_start(), + * and update_values_finish(). + */ + NoncontiguousPartitioner( + const std::vector &indices_has, + const std::vector &indices_want, + const MPI_Comm & communicator); + + /** + * Fill the vector @p dst according to the precomputed communication + * pattern with values from @src. + * + * @pre The vectors only have to provide a method begin(), which allows + * to access their raw data. + * + * @note This function calls the methods update_values_start() and + * update_values_finish() in sequence. Users can call these two + * functions separately and hereby overlap communication and + * computation. + */ + template + void + update_values(VectorType &dst, const VectorType &src) const; + + /** + * Start update. Data is packed as well as non-blocking send and receives + * are started. + */ + template + void + update_values_start(const VectorType &src, const unsigned int tag) const; + + /** + * Finish update. The method waits until all data has been sent and + * received. Once data from any process is received it is processed and + * placed at the right position of the vector @p dst. + */ + template + void + update_values_finish(VectorType &dst, const unsigned int tag) const; + + /** + * Returns the number of processes this process sends data to and + * number of processes this process received data from. + */ + std::pair + n_targets(); + + /** + * Return memory consumption in Byte. + */ + types::global_dof_index + memory_consumption(); + + /** + * Return the underlying communicator. + */ + const MPI_Comm & + get_mpi_communicator() const override; + + /** + * Initialize the inner data structures. + */ + void + reinit(const IndexSet &indexset_has, + const IndexSet &indexset_want, + const MPI_Comm &communicator) override; + + private: + /// MPI communicator + MPI_Comm communicator; + + /// CRS and MPI data structures for sending + /// The ranks this process sends data to. + std::vector send_ranks; + + /// Offset of each process within send_buffer. + std::vector send_ptr; + + /// Local index of each entry in send_buffer + /// with in the destination vector. + std::vector send_indices; + + /// Buffer containing the values sorted accoding to the ranks. + mutable std::vector send_buffers; + + /// MPI requests. + mutable std::vector send_requests; + + /// CRS and MPI data structures for receiving + //// The ranks this process receives data from. + std::vector recv_ranks; + + /// Offset of each process within recv_buffer. + std::vector recv_ptr; + + /// Local index of each entry in recv_buffer + /// with in the destination vector. + std::vector recv_indices; + + /// Buffer containing the values sorted accoding to the ranks. + mutable std::vector recv_buffers; + + /// MPI requests. + mutable std::vector recv_requests; + }; + + } // namespace MPI +} // namespace Utilities + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h new file mode 100644 index 0000000000..e8b6ca134f --- /dev/null +++ b/include/deal.II/base/mpi_noncontiguous_partitioner.templates.h @@ -0,0 +1,328 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + +#ifndef dealii_mpi_noncontiguous_vector_template_h +#define dealii_mpi_noncontiguous_vector_template_h + +#include + +#include +#include +#include +#include + +#include +#include + + +DEAL_II_NAMESPACE_OPEN + +namespace Utilities +{ + namespace MPI + { + template + NoncontiguousPartitioner::NoncontiguousPartitioner( + const IndexSet &indexset_has, + const IndexSet &indexset_want, + const MPI_Comm &communicator) + { + this->reinit(indexset_has, indexset_want, communicator); + } + + template + NoncontiguousPartitioner::NoncontiguousPartitioner( + const std::vector &indices_has, + const std::vector &indices_want, + const MPI_Comm & communicator) + { + // step 0) determine "number of degrees of freedom" needed for IndexSet + const types::global_dof_index n_dofs = Utilities::MPI::max( + std::max( + [&indices_has]() { + const auto it = max_element(indices_has.begin(), indices_has.end()); + return it == indices_has.end() ? 0 : (*it + 1); + }(), + [&indices_want]() { + const auto it = + max_element(indices_want.begin(), indices_want.end()); + return it == indices_want.end() ? 0 : (*it + 1); + }()), + communicator); + + // step 1) convert vectors to indexsets (sorted!) + IndexSet index_set_has(n_dofs); + index_set_has.add_indices(indices_has.begin(), indices_has.end()); + + IndexSet index_set_want(n_dofs); + index_set_want.add_indices(indices_want.begin(), indices_want.end()); + + // step 2) setup internal data structures with indexset + this->reinit(index_set_has, index_set_want, communicator); + + // step 3) fix inner data structures so that it is sorted as + // in the original vector + { + std::vector temp_map_send( + index_set_has.n_elements()); + + for (types::global_dof_index i = 0; i < indices_has.size(); i++) + temp_map_send[index_set_has.index_within_set(indices_has[i])] = i; + + for (auto &i : send_indices) + i = temp_map_send[i]; + } + + { + std::vector temp_map_recv( + index_set_want.n_elements()); + + for (types::global_dof_index i = 0; i < indices_want.size(); i++) + temp_map_recv[index_set_want.index_within_set(indices_want[i])] = i; + + for (auto &i : recv_indices) + i = temp_map_recv[i]; + } + } + + template + std::pair + NoncontiguousPartitioner::n_targets() + { + return {send_ranks.size(), recv_ranks.size()}; + } + + template + types::global_dof_index + NoncontiguousPartitioner::memory_consumption() + { + return MemoryConsumption::memory_consumption(send_ranks) + + MemoryConsumption::memory_consumption(send_ptr) + + MemoryConsumption::memory_consumption(send_indices) + + MemoryConsumption::memory_consumption(send_buffers) + + MemoryConsumption::memory_consumption(send_requests) + + MemoryConsumption::memory_consumption(recv_ranks) + + MemoryConsumption::memory_consumption(recv_ptr) + + MemoryConsumption::memory_consumption(recv_indices) + + MemoryConsumption::memory_consumption(recv_buffers) + + MemoryConsumption::memory_consumption(recv_requests); + } + + template + const MPI_Comm & + NoncontiguousPartitioner::get_mpi_communicator() const + { + return communicator; + } + + template + void + NoncontiguousPartitioner::reinit(const IndexSet &indexset_has, + const IndexSet &indexset_want, + const MPI_Comm &communicator) + { + this->communicator = communicator; + + // clean up + send_ranks.clear(); + send_ptr.clear(); + send_indices.clear(); + send_buffers.clear(); + send_requests.clear(); + recv_ranks.clear(); + recv_ptr.clear(); + recv_indices.clear(); + recv_buffers.clear(); + recv_requests.clear(); + + // setup communication pattern + std::vector owning_ranks_of_ghosts( + indexset_want.n_elements()); + + // set up dictionary + Utilities::MPI::internal::ComputeIndexOwner::ConsensusAlgorithmPayload + process(indexset_has, + indexset_want, + communicator, + owning_ranks_of_ghosts, + true); + + Utilities::MPI::ConsensusAlgorithmSelector< + std::pair, + unsigned int> + consensus_algorithm(process, communicator); + consensus_algorithm.run(); + + // setup map of processes from where this rank will receive values + { + std::map> recv_map; + + for (const auto &owner : owning_ranks_of_ghosts) + recv_map[owner] = std::vector(); + + for (types::global_dof_index i = 0; i < owning_ranks_of_ghosts.size(); + i++) + recv_map[owning_ranks_of_ghosts[i]].push_back(i); + + recv_ptr.push_back(recv_indices.size() /*=0*/); + for (const auto &target_with_indexset : recv_map) + { + recv_ranks.push_back(target_with_indexset.first); + + for (const auto &cell_index : target_with_indexset.second) + recv_indices.push_back(cell_index); + + recv_ptr.push_back(recv_indices.size()); + } + + recv_buffers.resize(recv_indices.size()); + recv_requests.resize(recv_map.size()); + } + + { + const auto targets_with_indexset = process.get_requesters(); + + send_ptr.push_back(send_indices.size() /*=0*/); + for (const auto &target_with_indexset : targets_with_indexset) + { + send_ranks.push_back(target_with_indexset.first); + + for (const auto &cell_index : target_with_indexset.second) + send_indices.push_back(indexset_has.index_within_set(cell_index)); + + send_ptr.push_back(send_indices.size()); + } + + send_buffers.resize(send_indices.size()); + send_requests.resize(targets_with_indexset.size()); + } + } + + + + template + template + void + NoncontiguousPartitioner::update_values(VectorType & dst, + const VectorType &src) const + { + const auto tag = internal::Tags::noncontiguous_partitioner_update_values; + + this->update_values_start(src, tag); + this->update_values_finish(dst, tag); + } + + + + template + template + void + NoncontiguousPartitioner::update_values_start( + const VectorType & src, + const unsigned int tag) const + { +#ifndef DEAL_II_WITH_MPI + (void)src; + (void)tag; + Assert(false, ExcNeedsMPI()); +#else + // post recv + for (types::global_dof_index i = 0; i < recv_ranks.size(); i++) + { + const auto ierr = MPI_Irecv(recv_buffers.data() + recv_ptr[i], + recv_ptr[i + 1] - recv_ptr[i], + Utilities::MPI::internal::mpi_type_id( + recv_buffers.data()), + recv_ranks[i], + tag, + communicator, + &recv_requests[i]); + AssertThrowMPI(ierr); + } + + auto src_iterator = src.begin(); + + // post send + for (types::global_dof_index i = 0; i < send_ranks.size(); i++) + { + // collect data to be send + for (types::global_dof_index j = send_ptr[i], c = 0; + j < send_ptr[i + 1]; + j++) + send_buffers[send_ptr[i] + c++] = src_iterator[send_indices[j]]; + + // send data + const auto ierr = MPI_Isend(send_buffers.data() + send_ptr[i], + send_ptr[i + 1] - send_ptr[i], + Utilities::MPI::internal::mpi_type_id( + send_buffers.data()), + send_ranks[i], + tag, + communicator, + &send_requests[i]); + AssertThrowMPI(ierr); + } +#endif + } + + + + template + template + void + NoncontiguousPartitioner::update_values_finish( + VectorType & dst, + const unsigned int tag) const + { + (void)tag; + +#ifndef DEAL_II_WITH_MPI + (void)dst; + Assert(false, ExcNeedsMPI()); +#else + auto dst_iterator = dst.begin(); + + // receive all data packages and copy data from buffers + for (types::global_dof_index proc = 0; proc < recv_requests.size(); + proc++) + { + int i; + MPI_Status status; + const auto ierr = MPI_Waitany(recv_requests.size(), + recv_requests.data(), + &i, + &status); + AssertThrowMPI(ierr); + + for (types::global_dof_index j = recv_ptr[i], c = 0; + j < recv_ptr[i + 1]; + j++) + dst_iterator[recv_indices[j]] = recv_buffers[recv_ptr[i] + c++]; + } + + // wait that all data packages have been sent + const auto ierr = MPI_Waitall(send_requests.size(), + send_requests.data(), + MPI_STATUSES_IGNORE); + AssertThrowMPI(ierr); +#endif + } + + } // namespace MPI +} // namespace Utilities + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/include/deal.II/base/mpi_tags.h b/include/deal.II/base/mpi_tags.h index 538b9e5f6b..75ec944938 100644 --- a/include/deal.II/base/mpi_tags.h +++ b/include/deal.II/base/mpi_tags.h @@ -129,6 +129,8 @@ namespace Utilities partitioner_export_start, partitioner_export_end = partitioner_export_start + 200, + /// NoncontiguousPartitioner::update_values + noncontiguous_partitioner_update_values, }; } // namespace Tags diff --git a/source/base/CMakeLists.txt b/source/base/CMakeLists.txt index ed2f0f819f..ee453a1f6e 100644 --- a/source/base/CMakeLists.txt +++ b/source/base/CMakeLists.txt @@ -46,6 +46,7 @@ SET(_unity_include_src hdf5.cc mpi.cc mpi_consensus_algorithm.cc + mpi_noncontiguous_partitioner.cc mu_parser_internal.cc multithread_info.cc named_selection.cc @@ -125,6 +126,7 @@ SET(_inst geometric_utilities.inst.in hdf5.inst.in mpi.inst.in + mpi_noncontiguous_partitioner.inst.in partitioner.inst.in partitioner.cuda.inst.in polynomials_rannacher_turek.inst.in diff --git a/source/base/mpi_noncontiguous_partitioner.cc b/source/base/mpi_noncontiguous_partitioner.cc new file mode 100644 index 0000000000..a4e043b2ef --- /dev/null +++ b/source/base/mpi_noncontiguous_partitioner.cc @@ -0,0 +1,27 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + +#include +#include +#include + +#include +#include + +DEAL_II_NAMESPACE_OPEN + +#include "mpi_noncontiguous_partitioner.inst" + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/base/mpi_noncontiguous_partitioner.inst.in b/source/base/mpi_noncontiguous_partitioner.inst.in new file mode 100644 index 0000000000..e7719c667c --- /dev/null +++ b/source/base/mpi_noncontiguous_partitioner.inst.in @@ -0,0 +1,52 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 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. +// +// --------------------------------------------------------------------- + + + +for (S : REAL_SCALARS) + { + namespace Utilities + \{ + namespace MPI + \{ + template class NoncontiguousPartitioner; + + template void + NoncontiguousPartitioner::update_values( + std::vector & dst, + const std::vector &src) const; + + template void + NoncontiguousPartitioner::update_values( + AlignedVector & dst, + const AlignedVector &src) const; + + template void + NoncontiguousPartitioner::update_values( + ArrayView & dst, + const ArrayView &src) const; + + template void + NoncontiguousPartitioner::update_values( + LinearAlgebra::Vector & dst, + const LinearAlgebra::Vector &src) const; + + template void + NoncontiguousPartitioner::update_values( + LinearAlgebra::distributed::Vector & dst, + const LinearAlgebra::distributed::Vector &src) const; + \} + \} + } diff --git a/tests/base/mpi_noncontiguous_partitioner_01.cc b/tests/base/mpi_noncontiguous_partitioner_01.cc new file mode 100644 index 0000000000..23ba353d5a --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_01.cc @@ -0,0 +1,76 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 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 Utilities::MPI::NoncontiguousPartitioner for non-contiguous index space. + +#include +#include + +#include "../tests.h" + +using namespace dealii; + +void +test(const MPI_Comm comm) +{ + IndexSet index_set_has(4); + IndexSet index_set_want(4); + + if (Utilities::MPI::this_mpi_process(comm) == 0) + { + index_set_has.add_index(1); + index_set_want.add_index(2); + } + else + { + index_set_has.add_index(2); + index_set_want.add_index(1); + index_set_want.add_index(2); + } + + Utilities::MPI::NoncontiguousPartitioner vector(index_set_has, + index_set_want, + comm); + + AlignedVector src(index_set_has.n_elements()); + AlignedVector dst(index_set_want.n_elements()); + + src[0] = Utilities::MPI::this_mpi_process(comm) * 100 + 1; + + vector.update_values(dst, src); + + for (size_t i = 0; i < src.size(); i++) + deallog << static_cast(src[i]) << " "; + deallog << std::endl; + for (size_t i = 0; i < dst.size(); i++) + deallog << static_cast(dst[i]) << " "; + deallog << std::endl; +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + const MPI_Comm comm = MPI_COMM_WORLD; + + { + deallog.push("all"); + test(comm); + deallog.pop(); + } +} diff --git a/tests/base/mpi_noncontiguous_partitioner_01.mpirun=2.output b/tests/base/mpi_noncontiguous_partitioner_01.mpirun=2.output new file mode 100644 index 0000000000..1ceffc9f3f --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_01.mpirun=2.output @@ -0,0 +1,7 @@ + +DEAL:0:all::1 +DEAL:0:all::101 + +DEAL:1:all::101 +DEAL:1:all::1 101 + diff --git a/tests/base/mpi_noncontiguous_partitioner_02.cc b/tests/base/mpi_noncontiguous_partitioner_02.cc new file mode 100644 index 0000000000..bc36dcadb4 --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_02.cc @@ -0,0 +1,161 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 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 Utilities::MPI::NoncontiguousPartitioner for repartitioning of a +// degrees of freedom (Morton order -> layer partitioning) + +#include +#include + +#include + +#include + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +template +void +test(const MPI_Comm &comm, const bool do_revert, const unsigned int dir) +{ + const unsigned int degree = 2; + const unsigned int n_refinements = 2; + const unsigned int n_cells_1D = Utilities::pow(2, n_refinements); + const unsigned int n_points_1D = (degree + 1) * n_cells_1D; + const double delta = 1.0 / n_cells_1D; + const unsigned int n_points_cell_1D = degree + 1; + const unsigned int n_points_cell = Utilities::pow(degree + 1, dim); + const unsigned int n_points_face = Utilities::pow(n_points_1D, dim - 1); + + const unsigned int n_procs = Utilities::MPI::n_mpi_processes(comm); + const unsigned int my_rank = + do_revert ? (n_procs - 1 - Utilities::MPI::this_mpi_process(comm)) : + Utilities::MPI::this_mpi_process(comm); + + parallel::distributed::Triangulation tria(comm); + + GridGenerator::hyper_cube(tria); + tria.refine_global(n_refinements); + + FE_DGQ fe(degree); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + const unsigned n_dofs = dof_handler.n_dofs(); + + + std::vector indices_has, indices_want; + + auto norm_point_to_lex = [&](const auto c) { + // convert normalized point [0, 1] to lex + if (dim == 2) + return std::floor(c[0]) + n_cells_1D * std::floor(c[1]); + else + return std::floor(c[0]) + n_cells_1D * std::floor(c[1]) + + n_cells_1D * n_cells_1D * std::floor(c[2]); + }; + + // ... has (symm) + for (const auto &cell : dof_handler.active_cell_iterators()) + if (cell->active() && cell->is_locally_owned()) + { + auto c = cell->center(); + for (unsigned int i = 0; i < dim; i++) + c[i] = c[i] / delta; + + const unsigned int lid = norm_point_to_lex(c); + + for (unsigned int i = lid * n_points_cell; + i < (lid + 1) * n_points_cell; + i++) + indices_has.push_back(i); + } + + + const unsigned int div = n_points_1D / n_procs; + const unsigned int rem = n_points_1D % n_procs; + + const unsigned int start = div * (my_rank + 0) + std::min(my_rank + 0, rem); + const unsigned int end = div * (my_rank + 1) + std::min(my_rank + 1, rem); + + if (dim == 2 && dir == 0) + { + for (unsigned int j = 0, c = start * n_points_face; j < end; j++) + for (unsigned int i = start; i < n_points_1D; i++) + indices_want.push_back(c++); + } + else if (dim == 2 && dir == 1) + { + for (unsigned int j = 0; j < n_points_1D; j++) + for (unsigned int i = start; i < end; i++) + indices_want.push_back(j * n_points_face + i); + } + else + Assert(false, StandardExceptions::ExcNotImplemented()); + + if (do_revert) + std::reverse(indices_want.begin(), indices_want.end()); + + Utilities::MPI::NoncontiguousPartitioner vector(indices_has, + indices_want, + comm); + + AlignedVector src(indices_has.size()); + for (unsigned int i = 0; i < indices_has.size(); i++) + src[i] = indices_has[i]; + + + AlignedVector dst(indices_want.size()); + + vector.update_values(dst, src); + + for (size_t i = 0; i < src.size(); i++) + deallog << static_cast(src[i]) << " "; + deallog << std::endl; + for (size_t i = 0; i < dst.size(); i++) + deallog << static_cast(dst[i]) << " "; + deallog << std::endl << std::endl; + + + for (size_t i = 0; i < dst.size(); i++) + AssertDimension(dst[i], indices_want[i]); +} + +template +void +test_dim(const MPI_Comm &comm, const bool do_revert) +{ + for (int dir = 0; dir < dim; dir++) + test(comm, do_revert, dir); +} + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + const MPI_Comm comm = MPI_COMM_WORLD; + + test_dim<2>(comm, /*do_revert=*/false); + test_dim<2>(comm, /*do_revert=*/true); +} diff --git a/tests/base/mpi_noncontiguous_partitioner_02.mpirun=1.output b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=1.output new file mode 100644 index 0000000000..20818dd82b --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=1.output @@ -0,0 +1,13 @@ + +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0::143 142 141 140 139 138 137 136 135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 95 94 93 92 91 90 89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 71 70 69 68 67 66 65 64 63 62 61 60 59 58 57 56 55 54 53 52 51 50 49 48 47 46 45 44 43 42 41 40 39 38 37 36 35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:0::143 142 141 140 139 138 137 136 135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 95 94 93 92 91 90 89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 71 70 69 68 67 66 65 64 63 62 61 60 59 58 57 56 55 54 53 52 51 50 49 48 47 46 45 44 43 42 41 40 39 38 37 36 35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0 +DEAL:0:: diff --git a/tests/base/mpi_noncontiguous_partitioner_02.mpirun=2.output b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=2.output new file mode 100644 index 0000000000..adc0250aad --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=2.output @@ -0,0 +1,27 @@ + +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:0::0 1 2 3 4 5 12 13 14 15 16 17 24 25 26 27 28 29 36 37 38 39 40 41 48 49 50 51 52 53 60 61 62 63 64 65 72 73 74 75 76 77 84 85 86 87 88 89 96 97 98 99 100 101 108 109 110 111 112 113 120 121 122 123 124 125 132 133 134 135 136 137 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:0::143 142 141 140 139 138 137 136 135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 95 94 93 92 91 90 89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:0::143 142 141 140 139 138 131 130 129 128 127 126 119 118 117 116 115 114 107 106 105 104 103 102 95 94 93 92 91 90 83 82 81 80 79 78 71 70 69 68 67 66 59 58 57 56 55 54 47 46 45 44 43 42 35 34 33 32 31 30 23 22 21 20 19 18 11 10 9 8 7 6 +DEAL:0:: + +DEAL:1::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:1::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:1:: +DEAL:1::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:1::6 7 8 9 10 11 18 19 20 21 22 23 30 31 32 33 34 35 42 43 44 45 46 47 54 55 56 57 58 59 66 67 68 69 70 71 78 79 80 81 82 83 90 91 92 93 94 95 102 103 104 105 106 107 114 115 116 117 118 119 126 127 128 129 130 131 138 139 140 141 142 143 +DEAL:1:: +DEAL:1::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:1::71 70 69 68 67 66 65 64 63 62 61 60 59 58 57 56 55 54 53 52 51 50 49 48 47 46 45 44 43 42 41 40 39 38 37 36 35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0 +DEAL:1:: +DEAL:1::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:1::137 136 135 134 133 132 125 124 123 122 121 120 113 112 111 110 109 108 101 100 99 98 97 96 89 88 87 86 85 84 77 76 75 74 73 72 65 64 63 62 61 60 53 52 51 50 49 48 41 40 39 38 37 36 29 28 27 26 25 24 17 16 15 14 13 12 5 4 3 2 1 0 +DEAL:1:: + diff --git a/tests/base/mpi_noncontiguous_partitioner_02.mpirun=4.output b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=4.output new file mode 100644 index 0000000000..4860113984 --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=4.output @@ -0,0 +1,55 @@ + +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::0 1 2 12 13 14 24 25 26 36 37 38 48 49 50 60 61 62 72 73 74 84 85 86 96 97 98 108 109 110 120 121 122 132 133 134 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::143 142 141 140 139 138 137 136 135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::143 142 141 131 130 129 119 118 117 107 106 105 95 94 93 83 82 81 71 70 69 59 58 57 47 46 45 35 34 33 23 22 21 11 10 9 +DEAL:0:: + +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::3 4 5 15 16 17 27 28 29 39 40 41 51 52 53 63 64 65 75 76 77 87 88 89 99 100 101 111 112 113 123 124 125 135 136 137 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 95 94 93 92 91 90 89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::140 139 138 128 127 126 116 115 114 104 103 102 92 91 90 80 79 78 68 67 66 56 55 54 44 43 42 32 31 30 20 19 18 8 7 6 +DEAL:1:: + + +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:2:: +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:2::6 7 8 18 19 20 30 31 32 42 43 44 54 55 56 66 67 68 78 79 80 90 91 92 102 103 104 114 115 116 126 127 128 138 139 140 +DEAL:2:: +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:2::89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 71 70 69 68 67 66 65 64 63 62 61 60 59 58 57 56 55 54 53 52 51 50 49 48 47 46 45 44 43 42 41 40 39 38 37 36 +DEAL:2:: +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:2::137 136 135 125 124 123 113 112 111 101 100 99 89 88 87 77 76 75 65 64 63 53 52 51 41 40 39 29 28 27 17 16 15 5 4 3 +DEAL:2:: + + +DEAL:3::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:3::108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:3:: +DEAL:3::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:3::9 10 11 21 22 23 33 34 35 45 46 47 57 58 59 69 70 71 81 82 83 93 94 95 105 106 107 117 118 119 129 130 131 141 142 143 +DEAL:3:: +DEAL:3::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:3::35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0 +DEAL:3:: +DEAL:3::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:3::134 133 132 122 121 120 110 109 108 98 97 96 86 85 84 74 73 72 62 61 60 50 49 48 38 37 36 26 25 24 14 13 12 2 1 0 +DEAL:3:: + diff --git a/tests/base/mpi_noncontiguous_partitioner_02.mpirun=5.output b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=5.output new file mode 100644 index 0000000000..0cdd4ea23c --- /dev/null +++ b/tests/base/mpi_noncontiguous_partitioner_02.mpirun=5.output @@ -0,0 +1,69 @@ + +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::0 1 2 12 13 14 24 25 26 36 37 38 48 49 50 60 61 62 72 73 74 84 85 86 96 97 98 108 109 110 120 121 122 132 133 134 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::143 142 141 140 139 138 137 136 135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 +DEAL:0:: +DEAL:0::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 +DEAL:0::143 142 131 130 119 118 107 106 95 94 83 82 71 70 59 58 47 46 35 34 23 22 11 10 +DEAL:0:: + +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::3 4 5 15 16 17 27 28 29 39 40 41 51 52 53 63 64 65 75 76 77 87 88 89 99 100 101 111 112 113 123 124 125 135 136 137 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::135 134 133 132 131 130 129 128 127 126 125 124 123 122 121 120 119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 +DEAL:1:: +DEAL:1::18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 +DEAL:1::141 140 129 128 117 116 105 104 93 92 81 80 69 68 57 56 45 44 33 32 21 20 9 8 +DEAL:1:: + + +DEAL:2:: +DEAL:2::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 +DEAL:2:: +DEAL:2:: +DEAL:2::6 7 18 19 30 31 42 43 54 55 66 67 78 79 90 91 102 103 114 115 126 127 138 139 +DEAL:2:: +DEAL:2:: +DEAL:2::119 118 117 116 115 114 113 112 111 110 109 108 107 106 105 104 103 102 101 100 99 98 97 96 95 94 93 92 91 90 89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 +DEAL:2:: +DEAL:2:: +DEAL:2::139 138 127 126 115 114 103 102 91 90 79 78 67 66 55 54 43 42 31 30 19 18 7 6 +DEAL:2:: + + +DEAL:3::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:3::96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 +DEAL:3:: +DEAL:3::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:3::8 9 20 21 32 33 44 45 56 57 68 69 80 81 92 93 104 105 116 117 128 129 140 141 +DEAL:3:: +DEAL:3::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:3::89 88 87 86 85 84 83 82 81 80 79 78 77 76 75 74 73 72 71 70 69 68 67 66 65 64 63 62 61 60 59 58 57 56 55 54 53 52 51 50 49 48 47 46 45 44 43 42 41 40 39 38 37 36 +DEAL:3:: +DEAL:3::72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 +DEAL:3::137 136 135 125 124 123 113 112 111 101 100 99 89 88 87 77 76 75 65 64 63 53 52 51 41 40 39 29 28 27 17 16 15 5 4 3 +DEAL:3:: + + +DEAL:4::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:4::120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:4:: +DEAL:4::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:4::10 11 22 23 34 35 46 47 58 59 70 71 82 83 94 95 106 107 118 119 130 131 142 143 +DEAL:4:: +DEAL:4::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:4::35 34 33 32 31 30 29 28 27 26 25 24 23 22 21 20 19 18 17 16 15 14 13 12 11 10 9 8 7 6 5 4 3 2 1 0 +DEAL:4:: +DEAL:4::90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 +DEAL:4::134 133 132 122 121 120 110 109 108 98 97 96 86 85 84 74 73 72 62 61 60 50 49 48 38 37 36 26 25 24 14 13 12 2 1 0 +DEAL:4:: + -- 2.39.5