From a126e8178d042b074c62b309c5a4b365e1f9525e Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Tue, 17 Apr 2018 11:34:58 -0400 Subject: [PATCH] add ConstraintMatrix::is_consistent_in_parallel --- include/deal.II/lac/constraint_matrix.h | 25 + source/lac/constraint_matrix.cc | 95 ++++ tests/mpi/constraints_consistent_01.cc | 149 ++++++ .../constraints_consistent_01.mpirun=1.output | 275 ++++++++++ .../constraints_consistent_01.mpirun=3.output | 500 ++++++++++++++++++ 5 files changed, 1044 insertions(+) create mode 100644 tests/mpi/constraints_consistent_01.cc create mode 100644 tests/mpi/constraints_consistent_01.mpirun=1.output create mode 100644 tests/mpi/constraints_consistent_01.mpirun=3.output diff --git a/include/deal.II/lac/constraint_matrix.h b/include/deal.II/lac/constraint_matrix.h index 7c12c0d0c5..a7dd469771 100644 --- a/include/deal.II/lac/constraint_matrix.h +++ b/include/deal.II/lac/constraint_matrix.h @@ -1203,6 +1203,31 @@ public: const LineRange get_lines() const; + /** + * Check if the current object is consistent on all processors + * in a distributed computation. + * + * This method checks if all processors agree on the constraints for their + * local lines as returned by get_local_lines(). This method is a collective + * operation and will return @p true only if all processors are consistent. + * + * Please supply the owned DoFs per processor as return by + * DoFHandler::locally_owned_dofs_per_processor() as @p locally_owned_dofs. + * + * If @p verbose is set to @p true, additional debug information is written + * to std::cout. + * + * @note This method exchanges all constraint information of locally relevant + * lines and is as such slow for large computations and should probably + * only be used in debug mode. + * + * @return Whether all ConstraintMatrix objects are consistent. + */ + bool is_consistent_in_parallel(const std::vector &locally_owned_dofs, + const MPI_Comm mpi_communicator, + const bool verbose=false) const; + + /** * Exception * diff --git a/source/lac/constraint_matrix.cc b/source/lac/constraint_matrix.cc index 5cd0e16dd6..9c1d98e3d9 100644 --- a/source/lac/constraint_matrix.cc +++ b/source/lac/constraint_matrix.cc @@ -42,6 +42,7 @@ #include #include #include +#include DEAL_II_NAMESPACE_OPEN @@ -103,6 +104,100 @@ const ConstraintMatrix::LineRange ConstraintMatrix::get_lines() const } +namespace +{ + struct line + { + types::global_dof_index row; + double inhom; + std::vector< std::pair > entries; + + template + void serialize(Archive &ar, const unsigned int) + { + ar &row &entries &inhom; + } + + line() {} + + line(const ConstraintMatrix &cm, const types::global_dof_index row) + : row (row) + { + inhom = cm.get_inhomogeneity(row); + if (auto v=cm.get_constraint_entries(row)) + entries = *v; + } + }; +} + +bool ConstraintMatrix::is_consistent_in_parallel(const std::vector &locally_owned_dofs, + const MPI_Comm mpi_communicator, + const bool verbose) const +{ + // identify non-owned rows and send to owner: + std::map< unsigned int, std::vector > to_send; + + const unsigned int myid = dealii::Utilities::MPI::this_mpi_process(mpi_communicator); + IndexSet non_owned = this->get_local_lines(); + non_owned.subtract_set(locally_owned_dofs[myid]); + for (auto it = non_owned.begin(); it != non_owned.end(); ++it) + { + line l(*this, *it); + unsigned int owner; + // this is inefficient, but I guess we don't care + for (owner=0; owner > received = Utilities::MPI::some_to_some (MPI_COMM_WORLD, to_send); + + unsigned int inconsistent = 0; + + + // from each processor: + for (auto &kv : received) + { + // for each incoming line: + for (auto &lineit : kv.second) + { + line ref(*this, lineit.row); + + if (lineit.inhom!=ref.inhom) + { + ++inconsistent; + + if (verbose) + std::cout << "Proc " << myid + << " got line " << lineit.row + << " from " << kv.first + << " inhomogeneity " << lineit.inhom << " != " << ref.inhom << std::endl; + } + else if (ref.entries != lineit.entries) + { + ++inconsistent; + if (verbose) + std::cout << "Proc " << myid + << " got line " << lineit.row + << " from " << kv.first + << " wrong values!" + << std::endl; + } + } + } + + const unsigned int total = Utilities::MPI::sum(inconsistent, mpi_communicator); + if (verbose && total>0 && myid==0) + std::cout << total << " inconsistent lines discovered!" << std::endl; + return total==0; +} + + void ConstraintMatrix::add_lines (const std::set &lines) diff --git a/tests/mpi/constraints_consistent_01.cc b/tests/mpi/constraints_consistent_01.cc new file mode 100644 index 0000000000..bef1b8d747 --- /dev/null +++ b/tests/mpi/constraints_consistent_01.cc @@ -0,0 +1,149 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2009 - 2017 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. +// +// --------------------------------------------------------------------- + + + +// test ConstraintMatrix::is_consistent_in_parallel + +#include "../tests.h" + + +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include +#include +#include +#include + +template +void check(parallel::distributed::Triangulation &tria) +{ + MPILogInitAll all; + + DoFHandler<2> dof_handler(tria); + FESystem fe(FE_Q(2),dim); + + dof_handler.distribute_dofs (fe); + + IndexSet locally_owned_dofs = dof_handler.locally_owned_dofs (); + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs (dof_handler, + locally_relevant_dofs); + + ConstraintMatrix constraints; + + constraints.clear (); + constraints.reinit (locally_relevant_dofs); + DoFTools::make_hanging_node_constraints (dof_handler, constraints); + + for (unsigned int id=0; id < 1; ++id) + dealii::VectorTools::interpolate_boundary_values(dof_handler, + id, + ConstantFunction(static_cast(id)+42.0,dim), + constraints); + + VectorTools::compute_no_normal_flux_constraints(dof_handler, + 0, /*first component*/ + std::set {2}, + constraints); + + constraints.close (); + deallog << "LocallyOwned = " << std::flush; + locally_owned_dofs.print(deallog.get_file_stream()); + deallog << "constraints:" << std::endl; + constraints.print(deallog.get_file_stream()); + deallog << + "consistent? " + << + constraints.is_consistent_in_parallel(dof_handler.locally_owned_dofs_per_processor(), + MPI_COMM_WORLD, + true) + << std::endl; + + +} + + + + +template +void test() +{ + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tr,0.0, 1.0, false); + tr.refine_global (2); + check(tr); + } + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tr,0.0, 1.0, true); + tr.refine_global (2); + check(tr); + } + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::subdivided_hyper_cube(tr, 2); + tr.refine_global(1); + for (auto cell = tr.begin_active(); cell != tr.end(); ++cell) + { + if (cell->id().to_string()=="0_1:0") + cell->set_refine_flag(); + else if (cell->parent()->id().to_string()=="3_0:") + cell->set_coarsen_flag(); + } + + tr.execute_coarsening_and_refinement(); + check(tr); + } + { + deallog << "this should be inconsistent with more than one rank:" << std::endl; + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tr,0.0, 1.0, false); + tr.refine_global (2); + + for (auto cell = tr.begin_active(); cell != tr.end(); ++cell) + if (cell->is_locally_owned()) + { + for (unsigned int f(0); f::faces_per_cell; ++f) + { + if (cell->face(f)->at_boundary() && cell->face(f)->center()[0]<1e-10) + cell->face(f)->set_all_boundary_ids(1); + } + } + + check(tr); + } + + deallog << "OK" << std::endl; +} + + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization (argc, argv, 1); + + test<2>(); +} diff --git a/tests/mpi/constraints_consistent_01.mpirun=1.output b/tests/mpi/constraints_consistent_01.mpirun=1.output new file mode 100644 index 0000000000..c859d7a554 --- /dev/null +++ b/tests/mpi/constraints_consistent_01.mpirun=1.output @@ -0,0 +1,275 @@ + +DEAL:0::LocallyOwned = {[0,161]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 4 = 42.0000 + 5 = 42.0000 + 8 = 42.0000 + 9 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 62 = 42.0000 + 63 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 66 = 42.0000 + 67 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 114 = 42.0000 + 115 = 42.0000 + 118 = 42.0000 + 119 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 + 154 = 42.0000 + 155 = 42.0000 + 156 = 42.0000 + 157 = 42.0000 + 158 = 42.0000 + 159 = 42.0000 +DEAL:0::consistent? 1 +DEAL:0::LocallyOwned = {[0,161]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 3 = 0 + 4 = 42.0000 + 5 = 42.0000 + 8 = 42.0000 + 9 = 42.0000 + 13 = 0 + 19 = 0 + 25 = 0 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 51 = 0 + 57 = 0 + 63 = 0 + 69 = 0 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 114 = 42.0000 + 115 = 42.0000 +DEAL:0::consistent? 1 +DEAL:0::LocallyOwned = {[0,177]} +constraints: + 2 = 42.0000 + 3 = 42.0000 + 4 = 42.0000 + 5 = 42.0000 + 6 = 42.0000 + 7 = 42.0000 + 10 = 42.0000 + 11 = 42.0000 + 14 = 42.0000 + 15 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 20 = 42.0000 + 21 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 36 = 42.0000 + 37 = 42.0000 + 38 = 42.0000 + 39 = 42.0000 + 42 = 42.0000 + 43 = 42.0000 + 58 = 42.0000 + 59 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 70 = 42.0000 + 71 = 42.0000 + 72 = 42.0000 + 73 = 42.0000 + 74 = 42.0000 + 75 = 42.0000 + 76 = 42.0000 + 77 = 42.0000 + 82 12: 1.00000 + 83 13: 1.00000 + 86 0: 0.375000 + 86 12: 0.750000 + 86: -5.25000 + 87 1: 0.375000 + 87 13: 0.750000 + 87: -5.25000 + 90 = 42.0000 + 91 = 42.0000 + 92 0: -0.125000 + 92 12: 0.750000 + 92: 15.7500 + 93 1: -0.125000 + 93 13: 0.750000 + 93: 15.7500 + 96 = 42.0000 + 97 = 42.0000 + 100 = 42.0000 + 101 = 42.0000 + 108 8: 1.00000 + 109 9: 1.00000 + 110 0: 0.375000 + 110 8: 0.750000 + 110: -5.25000 + 111 1: 0.375000 + 111 9: 0.750000 + 111: -5.25000 + 116 = 42.0000 + 117 = 42.0000 + 118 = 42.0000 + 119 = 42.0000 + 120 = 42.0000 + 121 = 42.0000 + 124 = 42.0000 + 125 = 42.0000 + 128 0: -0.125000 + 128 8: 0.750000 + 128: 15.7500 + 129 1: -0.125000 + 129 9: 0.750000 + 129: 15.7500 + 130 = 42.0000 + 131 = 42.0000 + 134 = 42.0000 + 135 = 42.0000 + 136 = 42.0000 + 137 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 142 = 42.0000 + 143 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 152 26: 1.00000 + 153 27: 1.00000 + 154 22: -0.125000 + 154 26: 0.750000 + 154: 15.7500 + 155 23: -0.125000 + 155 27: 0.750000 + 155: 15.7500 + 156 = 42.0000 + 157 = 42.0000 + 162 46: 1.00000 + 163 47: 1.00000 + 164 = 42.0000 + 165 = 42.0000 + 168 22: -0.125000 + 168 46: 0.750000 + 168: 15.7500 + 169 23: -0.125000 + 169 47: 0.750000 + 169: 15.7500 + 172 22: 0.375000 + 172 26: 0.750000 + 172: -5.25000 + 173 23: 0.375000 + 173 27: 0.750000 + 173: -5.25000 + 174 22: 0.375000 + 174 46: 0.750000 + 174: -5.25000 + 175 23: 0.375000 + 175 47: 0.750000 + 175: -5.25000 +DEAL:0::consistent? 1 +DEAL::this should be inconsistent with more than one rank: +DEAL:0::LocallyOwned = {[0,161]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 62 = 42.0000 + 63 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 66 = 42.0000 + 67 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 118 = 42.0000 + 119 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 + 154 = 42.0000 + 155 = 42.0000 + 156 = 42.0000 + 157 = 42.0000 + 158 = 42.0000 + 159 = 42.0000 +DEAL:0::consistent? 1 +DEAL::OK diff --git a/tests/mpi/constraints_consistent_01.mpirun=3.output b/tests/mpi/constraints_consistent_01.mpirun=3.output new file mode 100644 index 0000000000..8151cee2de --- /dev/null +++ b/tests/mpi/constraints_consistent_01.mpirun=3.output @@ -0,0 +1,500 @@ + +DEAL:0::LocallyOwned = {[0,49]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 4 = 42.0000 + 5 = 42.0000 + 8 = 42.0000 + 9 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 +DEAL:0::consistent? 1 + +DEAL:1::LocallyOwned = {[50,129]} +constraints: + 2 = 42.0000 + 3 = 42.0000 + 4 = 42.0000 + 5 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 62 = 42.0000 + 63 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 66 = 42.0000 + 67 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 114 = 42.0000 + 115 = 42.0000 + 118 = 42.0000 + 119 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 +DEAL:1::consistent? 1 + + +DEAL:2::LocallyOwned = {[130,161]} +constraints: + 64 = 42.0000 + 65 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 + 154 = 42.0000 + 155 = 42.0000 + 156 = 42.0000 + 157 = 42.0000 + 158 = 42.0000 + 159 = 42.0000 +DEAL:2::consistent? 1 + +DEAL:0::LocallyOwned = {[0,49]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 3 = 0 + 4 = 42.0000 + 5 = 42.0000 + 8 = 42.0000 + 9 = 42.0000 + 13 = 0 + 19 = 0 + 25 = 0 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 51 = 0 + 57 = 0 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 +DEAL:0::consistent? 1 + +DEAL:1::LocallyOwned = {[50,129]} +constraints: + 3 = 0 + 4 = 42.0000 + 5 = 42.0000 + 19 = 0 + 25 = 0 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 51 = 0 + 57 = 0 + 63 = 0 + 69 = 0 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 114 = 42.0000 + 115 = 42.0000 +DEAL:1::consistent? 1 + + +DEAL:2::LocallyOwned = {[130,161]} +constraints: +DEAL:2::consistent? 1 + +DEAL:0::LocallyOwned = {[0,63]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 20 = 42.0000 + 21 = 42.0000 + 22 = 42.0000 + 23 = 42.0000 + 26 = 42.0000 + 27 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 36 8: 1.00000 + 37 9: 1.00000 + 38 4: -0.125000 + 38 8: 0.750000 + 38: 15.7500 + 39 5: -0.125000 + 39 9: 0.750000 + 39: 15.7500 + 40 = 42.0000 + 41 = 42.0000 + 46 = 42.0000 + 47 = 42.0000 + 48 72: 1.00000 + 49 73: 1.00000 + 50 = 42.0000 + 51 = 42.0000 + 54 4: -0.125000 + 54 72: 0.750000 + 54: 15.7500 + 55 5: -0.125000 + 55 73: 0.750000 + 55: 15.7500 + 58 4: 0.375000 + 58 8: 0.750000 + 58: -5.25000 + 59 5: 0.375000 + 59 9: 0.750000 + 59: -5.25000 + 60 4: 0.375000 + 60 72: 0.750000 + 60: -5.25000 + 61 5: 0.375000 + 61 73: 0.750000 + 61: -5.25000 + 64 = 42.0000 + 65 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 86 = 42.0000 + 87 = 42.0000 + 92 = 42.0000 + 93 = 42.0000 +DEAL:0::consistent? 1 + +DEAL:1::LocallyOwned = {[64,125]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 22 = 42.0000 + 23 = 42.0000 + 36 8: 1.00000 + 37 9: 1.00000 + 46 = 42.0000 + 47 = 42.0000 + 48 72: 1.00000 + 49 73: 1.00000 + 50 = 42.0000 + 51 = 42.0000 + 54 4: -0.125000 + 54 72: 0.750000 + 54: 15.7500 + 55 5: -0.125000 + 55 73: 0.750000 + 55: 15.7500 + 58 4: 0.375000 + 58 8: 0.750000 + 58: -5.25000 + 59 5: 0.375000 + 59 9: 0.750000 + 59: -5.25000 + 60 4: 0.375000 + 60 72: 0.750000 + 60: -5.25000 + 61 5: 0.375000 + 61 73: 0.750000 + 61: -5.25000 + 64 = 42.0000 + 65 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 86 = 42.0000 + 87 = 42.0000 + 92 = 42.0000 + 93 = 42.0000 + 98 = 42.0000 + 99 = 42.0000 + 100 = 42.0000 + 101 = 42.0000 + 102 = 42.0000 + 103 = 42.0000 + 104 = 42.0000 + 105 = 42.0000 + 110 134: 1.00000 + 111 135: 1.00000 + 114 78: 0.375000 + 114 134: 0.750000 + 114: -5.25000 + 115 79: 0.375000 + 115 135: 0.750000 + 115: -5.25000 + 118 = 42.0000 + 119 = 42.0000 + 120 = 42.0000 + 121 = 42.0000 + 122 78: -0.125000 + 122 134: 0.750000 + 122: 15.7500 + 123 79: -0.125000 + 123 135: 0.750000 + 123: 15.7500 + 126 = 42.0000 + 127 = 42.0000 + 128 = 42.0000 + 129 = 42.0000 + 132 = 42.0000 + 133 = 42.0000 + 136 = 42.0000 + 137 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 144 = 42.0000 + 145 = 42.0000 + 152 130: 1.00000 + 153 131: 1.00000 + 154 78: 0.375000 + 154 130: 0.750000 + 154: -5.25000 + 155 79: 0.375000 + 155 131: 0.750000 + 155: -5.25000 +DEAL:1::consistent? 1 + + +DEAL:2::LocallyOwned = {[126,177]} +constraints: + 46 = 42.0000 + 47 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 100 = 42.0000 + 101 = 42.0000 + 110 134: 1.00000 + 111 135: 1.00000 + 114 78: 0.375000 + 114 134: 0.750000 + 114: -5.25000 + 115 79: 0.375000 + 115 135: 0.750000 + 115: -5.25000 + 118 = 42.0000 + 119 = 42.0000 + 120 = 42.0000 + 121 = 42.0000 + 122 78: -0.125000 + 122 134: 0.750000 + 122: 15.7500 + 123 79: -0.125000 + 123 135: 0.750000 + 123: 15.7500 + 126 = 42.0000 + 127 = 42.0000 + 128 = 42.0000 + 129 = 42.0000 + 132 = 42.0000 + 133 = 42.0000 + 136 = 42.0000 + 137 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 144 = 42.0000 + 145 = 42.0000 + 152 130: 1.00000 + 153 131: 1.00000 + 154 78: 0.375000 + 154 130: 0.750000 + 154: -5.25000 + 155 79: 0.375000 + 155 131: 0.750000 + 155: -5.25000 + 160 = 42.0000 + 161 = 42.0000 + 162 = 42.0000 + 163 = 42.0000 + 164 = 42.0000 + 165 = 42.0000 + 168 = 42.0000 + 169 = 42.0000 + 172 78: -0.125000 + 172 130: 0.750000 + 172: 15.7500 + 173 79: -0.125000 + 173 131: 0.750000 + 173: 15.7500 + 174 = 42.0000 + 175 = 42.0000 +DEAL:2::consistent? 1 + +DEAL::this should be inconsistent with more than one rank: +DEAL:0::LocallyOwned = {[0,49]} +constraints: + 0 = 42.0000 + 1 = 42.0000 + 2 = 42.0000 + 3 = 42.0000 + 12 = 42.0000 + 13 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 90 = 42.0000 + 91 = 42.0000 + 94 = 42.0000 + 95 = 42.0000 +DEAL:0::consistent? 0 + +DEAL:1::LocallyOwned = {[50,129]} +constraints: + 2 = 42.0000 + 3 = 42.0000 + 4 = 42.0000 + 5 = 42.0000 + 18 = 42.0000 + 19 = 42.0000 + 24 = 42.0000 + 25 = 42.0000 + 30 = 42.0000 + 31 = 42.0000 + 34 = 42.0000 + 35 = 42.0000 + 50 = 42.0000 + 51 = 42.0000 + 56 = 42.0000 + 57 = 42.0000 + 62 = 42.0000 + 63 = 42.0000 + 64 = 42.0000 + 65 = 42.0000 + 66 = 42.0000 + 67 = 42.0000 + 68 = 42.0000 + 69 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 110 = 42.0000 + 111 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 118 = 42.0000 + 119 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 +DEAL:1::consistent? 0 + + +DEAL:2::LocallyOwned = {[130,161]} +constraints: + 64 = 42.0000 + 65 = 42.0000 + 82 = 42.0000 + 83 = 42.0000 + 84 = 42.0000 + 85 = 42.0000 + 112 = 42.0000 + 113 = 42.0000 + 122 = 42.0000 + 123 = 42.0000 + 126 = 42.0000 + 127 = 42.0000 + 138 = 42.0000 + 139 = 42.0000 + 140 = 42.0000 + 141 = 42.0000 + 146 = 42.0000 + 147 = 42.0000 + 150 = 42.0000 + 151 = 42.0000 + 154 = 42.0000 + 155 = 42.0000 + 156 = 42.0000 + 157 = 42.0000 + 158 = 42.0000 + 159 = 42.0000 +DEAL:2::consistent? 0 + +DEAL::OK -- 2.39.5