MueLu  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
MueLu_RepartitionHeuristicFactory_def.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // MueLu: A package for multigrid based preconditioning
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 // Tobias Wiesner (tawiesn@sandia.gov)
43 //
44 // ***********************************************************************
45 //
46 // @HEADER
47 #ifndef PACKAGES_MUELU_SRC_REBALANCING_MUELU_REPARTITIONHEURISTICFACTORY_DEF_HPP_
48 #define PACKAGES_MUELU_SRC_REBALANCING_MUELU_REPARTITIONHEURISTICFACTORY_DEF_HPP_
49 
50 #include <algorithm>
51 #include <iostream>
52 #include <sstream>
53 
54 #ifdef HAVE_MPI
56 #include <Teuchos_CommHelpers.hpp>
57 
58 //#include <Xpetra_Map.hpp>
59 #include <Xpetra_Matrix.hpp>
60 
61 #include "MueLu_Utilities.hpp"
62 
63 #include "MueLu_RAPFactory.hpp"
64 #include "MueLu_BlockedRAPFactory.hpp"
65 #include "MueLu_SubBlockAFactory.hpp"
66 #include "MueLu_Level.hpp"
67 #include "MueLu_MasterList.hpp"
68 #include "MueLu_Monitor.hpp"
69 
70 #include "MueLu_RepartitionHeuristicFactory.hpp"
71 
72 namespace MueLu {
73 
74  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
76  RCP<ParameterList> validParamList = rcp(new ParameterList());
77 
78 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
79  SET_VALID_ENTRY("repartition: start level");
80  SET_VALID_ENTRY("repartition: node repartition level");
81  SET_VALID_ENTRY("repartition: min rows per proc");
82  SET_VALID_ENTRY("repartition: target rows per proc");
83  SET_VALID_ENTRY("repartition: min rows per thread");
84  SET_VALID_ENTRY("repartition: target rows per thread");
85  SET_VALID_ENTRY("repartition: max imbalance");
86 #undef SET_VALID_ENTRY
87 
88  validParamList->set< RCP<const FactoryBase> >("A", Teuchos::null, "Factory of the matrix A");
89  validParamList->set< RCP<const FactoryBase> >("Node Comm", Teuchos::null, "Generating factory of the node level communicator");
90 
91  return validParamList;
92  }
93 
94  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
96  Input(currentLevel, "A");
97  const Teuchos::ParameterList & pL = GetParameterList();
98  if(pL.isParameter("repartition: node repartition level")) {
99  const int nodeRepartLevel = pL.get<int>("repartition: node repartition level");
100  if(currentLevel.GetLevelID() == nodeRepartLevel) {
101  Input(currentLevel,"Node Comm");
102  }
103  }
104  }
105 
106  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
108  FactoryMonitor m(*this, "Build", currentLevel);
109 
110  const Teuchos::ParameterList & pL = GetParameterList();
111  // Access parameters here to make sure that we set the parameter entry flag to "used" even in case of short-circuit evaluation.
112  // TODO (JG): I don't really know if we want to do this.
113  const int startLevel = pL.get<int> ("repartition: start level");
114  const int nodeRepartLevel = pL.get<int> ("repartition: node repartition level");
115  LO minRowsPerProcess = pL.get<LO> ("repartition: min rows per proc");
116  LO targetRowsPerProcess = pL.get<LO> ("repartition: target rows per proc");
117  LO minRowsPerThread = pL.get<LO> ("repartition: min rows per thread");
118  LO targetRowsPerThread = pL.get<LO> ("repartition: target rows per thread");
119  const double nonzeroImbalance = pL.get<double>("repartition: max imbalance");
120 
121  int thread_per_mpi_rank = 1;
122 #if defined(HAVE_MUELU_KOKKOSCORE) && defined(KOKKOS_ENABLE_OPENMP)
123  using execution_space = typename Node::device_type::execution_space;
124  if (std::is_same<execution_space, Kokkos::OpenMP>::value)
125  thread_per_mpi_rank = execution_space::concurrency();
126 #endif
127 
128  if (minRowsPerThread > 0)
129  // We ignore the value given by minRowsPerProcess and repartition based on threads instead
130  minRowsPerProcess = minRowsPerThread*thread_per_mpi_rank;
131 
132  if (targetRowsPerThread == 0)
133  targetRowsPerThread = minRowsPerThread;
134 
135  if (targetRowsPerThread > 0)
136  // We ignore the value given by targetRowsPerProcess and repartition based on threads instead
137  targetRowsPerProcess = targetRowsPerThread*thread_per_mpi_rank;
138 
139  if (targetRowsPerProcess == 0)
140  targetRowsPerProcess = minRowsPerProcess;
141 
142  // Stick this on the level so Zoltan2Interface can use this later
143  Set<LO>(currentLevel,"repartition: heuristic target rows per process",targetRowsPerProcess);
144 
145  // Check for validity of the node repartition option
146  TEUCHOS_TEST_FOR_EXCEPTION(nodeRepartLevel >= startLevel, Exceptions::RuntimeError, "MueLu::RepartitionHeuristicFactory::Build(): If 'repartition: node repartition level' is set, it must be less than or equal to 'repartition: start level'");
147 
148 
149  RCP<const FactoryBase> Afact = GetFactory("A");
150  if(!Afact.is_null() && Teuchos::rcp_dynamic_cast<const RAPFactory>(Afact) == Teuchos::null &&
151  Teuchos::rcp_dynamic_cast<const BlockedRAPFactory>(Afact) == Teuchos::null &&
152  Teuchos::rcp_dynamic_cast<const SubBlockAFactory>(Afact) == Teuchos::null) {
153  GetOStream(Warnings) <<
154  "MueLu::RepartitionHeuristicFactory::Build: The generation factory for A must " \
155  "be a RAPFactory or a SubBlockAFactory providing the non-rebalanced matrix information! " \
156  "It specifically must not be of type Rebalance(Blocked)AcFactory or similar. " \
157  "Please check the input. Make also sure that \"number of partitions\" is provided to " \
158  "the Interface class and the RepartitionFactory instance. Instead, we have a "<<Afact->description() << std::endl;
159  }
160  // TODO: We only need a CrsGraph. This class does not have to be templated on Scalar types.
161  RCP<Matrix> A = Get< RCP<Matrix> >(currentLevel, "A");
162 
163  // ======================================================================================================
164  // Determine whether partitioning is needed
165  // ======================================================================================================
166  // NOTE: most tests include some global communication, which is why we currently only do tests until we make
167  // a decision on whether to repartition. However, there is value in knowing how "close" we are to having to
168  // rebalance an operator. So, it would probably be beneficial to do and report *all* tests.
169 
170  // Test0: Should we do node repartitioning?
171  if (currentLevel.GetLevelID() == nodeRepartLevel && A->getMap()->getComm()->getSize() > 1) {
172  RCP<const Teuchos::Comm<int> > NodeComm = Get< RCP<const Teuchos::Comm<int> > >(currentLevel, "Node Comm");
173  TEUCHOS_TEST_FOR_EXCEPTION(NodeComm.is_null(), Exceptions::RuntimeError, "MueLu::RepartitionHeuristicFactory::Build(): NodeComm is null.");
174 
175  // If we only have one node, then we don't want to pop down to one rank
176  if(NodeComm()->getSize() != A->getMap()->getComm()->getSize()) {
177  GetOStream(Statistics1) << "Repartitioning? YES: \n Within node only"<<std::endl;
178  int nodeRank = NodeComm->getRank();
179 
180  // Do a reduction to get the total number of nodes
181  int isZero = (nodeRank == 0);
182  int numNodes=0;
183  Teuchos::reduceAll(*A->getMap()->getComm(), Teuchos::REDUCE_SUM, isZero, Teuchos::outArg(numNodes));
184  Set(currentLevel, "number of partitions", numNodes);
185  return;
186  }
187  }
188 
189  // Test1: skip repartitioning if current level is less than the specified minimum level for repartitioning
190  if (currentLevel.GetLevelID() < startLevel) {
191  GetOStream(Statistics1) << "Repartitioning? NO:" <<
192  "\n current level = " << Teuchos::toString(currentLevel.GetLevelID()) <<
193  ", first level where repartitioning can happen is " + Teuchos::toString(startLevel) << std::endl;
194 
195  // a negative number of processors means: no repartitioning
196  Set(currentLevel, "number of partitions", -1);
197 
198  return;
199  }
200 
201  RCP<const Map> rowMap = A->getRowMap();
202 
203  RCP<const Teuchos::Comm<int> > origComm = rowMap->getComm();
204  RCP<const Teuchos::Comm<int> > comm = origComm;
205 
206  // Test 2: check whether A is actually distributed, i.e. more than one processor owns part of A
207  // TODO: this global communication can be avoided if we store the information with the matrix (it is known when matrix is created)
208  // TODO: further improvements could be achieved when we use subcommunicator for the active set. Then we only need to check its size
209 
210  // TODO: The block transfer factories do not check correctly whether or not repartitioning actually took place.
211  {
212  if (comm->getSize() == 1 && Teuchos::rcp_dynamic_cast<const RAPFactory>(Afact) != Teuchos::null) {
213  GetOStream(Statistics1) << "Repartitioning? NO:" <<
214  "\n comm size = 1" << std::endl;
215 
216  Set(currentLevel, "number of partitions", -1);
217  return;
218  }
219 
220  int numActiveProcesses = 0;
221  MueLu_sumAll(comm, Teuchos::as<int>((A->getNodeNumRows() > 0) ? 1 : 0), numActiveProcesses);
222 
223  if (numActiveProcesses == 1) {
224  GetOStream(Statistics1) << "Repartitioning? NO:" <<
225  "\n # processes with rows = " << Teuchos::toString(numActiveProcesses) << std::endl;
226 
227  Set(currentLevel, "number of partitions", 1);
228  return;
229  }
230  }
231 
232  bool test3 = false, test4 = false;
233  std::string msg3, msg4;
234 
235  // Test3: check whether number of rows on any processor satisfies the minimum number of rows requirement
236  // NOTE: Test2 ensures that repartitionning is not done when there is only one processor (it may or may not satisfy Test3)
237  if (minRowsPerProcess > 0) {
238  LO numMyRows = Teuchos::as<LO>(A->getNodeNumRows()), minNumRows, LOMAX = Teuchos::OrdinalTraits<LO>::max();
239  LO haveFewRows = (numMyRows < minRowsPerProcess ? 1 : 0), numWithFewRows = 0;
240  MueLu_sumAll(comm, haveFewRows, numWithFewRows);
241  MueLu_minAll(comm, (numMyRows > 0 ? numMyRows : LOMAX), minNumRows);
242 
243  // TODO: we could change it to repartition only if the number of processors with numRows < minNumRows is larger than some
244  // percentage of the total number. This way, we won't repartition if 2 out of 1000 processors don't have enough elements.
245  // I'm thinking maybe 20% threshold. To implement, simply add " && numWithFewRows < .2*numProcs" to the if statement.
246  if (numWithFewRows > 0)
247  test3 = true;
248 
249  msg3 = "\n min # rows per proc = " + Teuchos::toString(minNumRows) + ", min allowable = " + Teuchos::toString(minRowsPerProcess);
250  }
251 
252  // Test4: check whether the balance in the number of nonzeros per processor is greater than threshold
253  if (!test3) {
254  GO minNnz, maxNnz, numMyNnz = Teuchos::as<GO>(A->getNodeNumEntries());
255  MueLu_maxAll(comm, numMyNnz, maxNnz);
256  MueLu_minAll(comm, (numMyNnz > 0 ? numMyNnz : maxNnz), minNnz); // min nnz over all active processors
257  double imbalance = Teuchos::as<double>(maxNnz)/minNnz;
258 
259  if (imbalance > nonzeroImbalance)
260  test4 = true;
261 
262  msg4 = "\n nonzero imbalance = " + Teuchos::toString(imbalance) + ", max allowable = " + Teuchos::toString(nonzeroImbalance);
263  }
264 
265  if (!test3 && !test4) {
266  GetOStream(Statistics1) << "Repartitioning? NO:" << msg3 + msg4 << std::endl;
267 
268  // A negative number of partitions means: no repartitioning
269  Set(currentLevel, "number of partitions", -1);
270  return;
271  }
272 
273  GetOStream(Statistics1) << "Repartitioning? YES:" << msg3 + msg4 << std::endl;
274 
275  // ======================================================================================================
276  // Calculate number of partitions
277  // ======================================================================================================
278  // FIXME Quick way to figure out how many partitions there should be (same algorithm as ML)
279  // FIXME Should take into account nnz? Perhaps only when user is using min #nnz per row threshold.
280 
281  // The number of partitions is calculated by the RepartitionFactory and stored in "number of partitions" variable on
282  // the current level. If this variable is already set (e.g., by another instance of RepartitionFactory) then this number
283  // is used. The "number of partitions" variable serves as basic communication between the RepartitionFactory (which
284  // requests a certain number of partitions) and the *Interface classes which call the underlying partitioning algorithms
285  // and produce the "Partition" array with the requested number of partitions.
286  const auto globalNumRows = Teuchos::as<GO>(A->getGlobalNumRows());
287  int numPartitions = 1;
288  if (globalNumRows >= targetRowsPerProcess) {
289  // Make sure that each CPU thread has approximately targetRowsPerProcess
290  numPartitions = std::max(Teuchos::as<int>(globalNumRows / targetRowsPerProcess), 1);
291  }
292  numPartitions = std::min(numPartitions, comm->getSize());
293 
294  Set(currentLevel, "number of partitions", numPartitions);
295 
296  GetOStream(Statistics1) << "Number of partitions to use = " << numPartitions << std::endl;
297  } // Build
298 } // namespace MueLu
299 
300 #endif //ifdef HAVE_MPI
301 #endif /* PACKAGES_MUELU_SRC_REBALANCING_MUELU_REPARTITIONHEURISTICFACTORY_DEF_HPP_ */
#define MueLu_sumAll(rcpComm, in, out)
#define MueLu_maxAll(rcpComm, in, out)
void Build(Level &currentLevel) const
Build an object with this factory.
GlobalOrdinal GO
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Print more statistics.
LocalOrdinal LO
#define SET_VALID_ENTRY(name)
#define MueLu_minAll(rcpComm, in, out)
bool isParameter(const std::string &name) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
Definition: MueLu_Level.hpp:99
Factory for building a thresholded operator.
TEUCHOS_DEPRECATED void reduceAll(const Comm< Ordinal > &comm, const EReductionType reductType, const Packet &send, Packet *globalReduct)
void DeclareInput(Level &currentLevel) const
Determines the data that RepartitionHeuristicFactory needs, and the factories that generate that data...
int GetLevelID() const
Return level number.
Definition: MueLu_Level.cpp:76
Exception throws to report errors in the internal logical of the program.
Print all warning messages.
Factory for building coarse matrices.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
virtual std::string description() const
Return a simple one-line description of this object.
std::string toString(const T &t)
bool is_null() const