Tpetra parallel linear algebra  Version of the Day
Tpetra_DirectoryImpl_def.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef __Tpetra_DirectoryImpl_def_hpp
43 #define __Tpetra_DirectoryImpl_def_hpp
44 
47 
49 #include <Tpetra_Distributor.hpp>
50 #include <Tpetra_Map.hpp>
51 #include <Tpetra_TieBreak.hpp>
52 
53 #include <Tpetra_Details_FixedHashTable.hpp>
54 #include <Tpetra_HashTable.hpp>
55 
56 
57 // FIXME (mfh 16 Apr 2013) GIANT HACK BELOW
58 #ifdef HAVE_MPI
59 # include "mpi.h"
60 #endif // HAVE_MPI
61 // FIXME (mfh 16 Apr 2013) GIANT HACK ABOVE
62 
63 
64 namespace Tpetra {
65  namespace Details {
66  template<class LO, class GO, class NT>
69 
70  template<class LO, class GO, class NT>
73  getEntries (const map_type& map,
74  const Teuchos::ArrayView<const GO> &globalIDs,
75  const Teuchos::ArrayView<int> &nodeIDs,
76  const Teuchos::ArrayView<LO> &localIDs,
77  const bool computeLIDs) const
78  {
79  // Ensure that globalIDs, nodeIDs, and localIDs (if applicable)
80  // all have the same size, before modifying any output arguments.
81  TEUCHOS_TEST_FOR_EXCEPTION(nodeIDs.size() != globalIDs.size(),
82  std::invalid_argument, Teuchos::typeName(*this) << "::getEntries(): "
83  "Output arrays do not have the right sizes. nodeIDs.size() = "
84  << nodeIDs.size() << " != globalIDs.size() = " << globalIDs.size()
85  << ".");
86  TEUCHOS_TEST_FOR_EXCEPTION(
87  computeLIDs && localIDs.size() != globalIDs.size(),
88  std::invalid_argument, Teuchos::typeName(*this) << "::getEntries(): "
89  "Output array do not have the right sizes. localIDs.size() = "
90  << localIDs.size() << " != globalIDs.size() = " << globalIDs.size()
91  << ".");
92 
93  // Initially, fill nodeIDs and localIDs (if applicable) with
94  // invalid values. The "invalid" process ID is -1 (this means
95  // the same thing as MPI_ANY_SOURCE to Teuchos, so it's an
96  // "invalid" process ID); the invalid local ID comes from
97  // OrdinalTraits.
98  std::fill (nodeIDs.begin(), nodeIDs.end(), -1);
99  if (computeLIDs) {
100  std::fill (localIDs.begin(), localIDs.end(),
101  Teuchos::OrdinalTraits<LO>::invalid ());
102  }
103  // Actually do the work.
104  return this->getEntriesImpl (map, globalIDs, nodeIDs, localIDs, computeLIDs);
105  }
106 
107 
108  template<class LO, class GO, class NT>
111  numProcs_ (map.getComm ()->getSize ())
112  {}
113 
114 
115  template<class LO, class GO, class NT>
118  numProcs_ (0) // to be set later
119  {}
120 
121 
122  template<class LO, class GO, class NT>
123  bool
125  isOneToOne (const Teuchos::Comm<int>& comm) const
126  {
127  // A locally replicated Map is one-to-one only if there is no
128  // replication, that is, only if the Map's communicator only has
129  // one process.
130  return (numProcs_ == 1);
131  }
132 
133 
134  template<class LO, class GO, class NT>
135  std::string
137  {
138  std::ostringstream os;
139  os << "ReplicatedDirectory"
140  << "<" << Teuchos::TypeNameTraits<LO>::name ()
141  << ", " << Teuchos::TypeNameTraits<GO>::name ()
142  << ", " << Teuchos::TypeNameTraits<NT>::name () << ">";
143  return os.str ();
144  }
145 
146 
147  template<class LO, class GO, class NT>
150  {
151  TEUCHOS_TEST_FOR_EXCEPTION(! map.isContiguous (), std::invalid_argument,
152  Teuchos::typeName (*this) << " constructor: Map is not contiguous.");
153  TEUCHOS_TEST_FOR_EXCEPTION(! map.isUniform (), std::invalid_argument,
154  Teuchos::typeName (*this) << " constructor: Map is not uniform.");
155  }
156 
157 
158  template<class LO, class GO, class NT>
159  std::string
161  {
162  std::ostringstream os;
163  os << "ContiguousUniformDirectory"
164  << "<" << Teuchos::TypeNameTraits<LO>::name ()
165  << ", " << Teuchos::TypeNameTraits<GO>::name ()
166  << ", " << Teuchos::TypeNameTraits<NT>::name () << ">";
167  return os.str ();
168  }
169 
170 
171  template<class LO, class GO, class NT>
175  const Teuchos::ArrayView<const GO> &globalIDs,
176  const Teuchos::ArrayView<int> &nodeIDs,
177  const Teuchos::ArrayView<LO> &localIDs,
178  const bool computeLIDs) const
179  {
180  using Teuchos::Comm;
181  using Teuchos::RCP;
182  typedef typename Teuchos::ArrayView<const GO>::size_type size_type;
183  const LO invalidLid = Teuchos::OrdinalTraits<LO>::invalid ();
185 
186  RCP<const Comm<int> > comm = map.getComm ();
187  const GO g_min = map.getMinAllGlobalIndex ();
188 
189  // Let N_G be the global number of elements in the Map,
190  // and P be the number of processes in its communicator.
191  // Then, N_G = P * N_L + R = R*(N_L + 1) + (P - R)*N_L.
192  //
193  // The first R processes own N_L+1 elements.
194  // The remaining P-R processes own N_L elements.
195  //
196  // Let g be the current GID, g_min be the global minimum GID,
197  // and g_0 = g - g_min. If g is a valid GID in this Map, then
198  // g_0 is in [0, N_G - 1].
199  //
200  // If g is a valid GID in this Map and g_0 < R*(N_L + 1), then
201  // the rank of the process that owns g is floor(g_0 / (N_L +
202  // 1)), and its corresponding local index on that process is g_0
203  // mod (N_L + 1).
204  //
205  // Let g_R = g_0 - R*(N_L + 1). If g is a valid GID in this Map
206  // and g_0 >= R*(N_L + 1), then the rank of the process that
207  // owns g is then R + floor(g_R / N_L), and its corresponding
208  // local index on that process is g_R mod N_L.
209 
210  const size_type N_G =
211  static_cast<size_type> (map.getGlobalNumElements ());
212  const size_type P = static_cast<size_type> (comm->getSize ());
213  const size_type N_L = N_G / P;
214  const size_type R = N_G - N_L * P; // N_G mod P
215  const size_type N_R = R * (N_L + static_cast<size_type> (1));
216 
217 #ifdef HAVE_TPETRA_DEBUG
218  TEUCHOS_TEST_FOR_EXCEPTION(
219  N_G != P*N_L + R, std::logic_error,
220  "Tpetra::ContiguousUniformDirectory::getEntriesImpl: "
221  "N_G = " << N_G << " != P*N_L + R = " << P << "*" << N_L << " + " << R
222  << " = " << P*N_L + R << ". "
223  "Please report this bug to the Tpetra developers.");
224 #endif // HAVE_TPETRA_DEBUG
225 
226  const size_type numGids = globalIDs.size (); // for const loop bound
227  // Avoid signed/unsigned comparisons below, in case GO is
228  // unsigned. (Integer literals are generally signed.)
229  const GO ONE = static_cast<GO> (1);
230 
231  if (computeLIDs) {
232  for (size_type k = 0; k < numGids; ++k) {
233  const GO g_0 = globalIDs[k] - g_min;
234 
235  // The first test is a little strange just in case GO is
236  // unsigned. Compilers raise a warning on tests like "x <
237  // 0" if x is unsigned, but don't usually raise a warning if
238  // the expression is a bit more complicated than that.
239  if (g_0 + ONE < ONE || g_0 >= static_cast<GO> (N_G)) {
240  nodeIDs[k] = -1;
241  localIDs[k] = invalidLid;
242  res = IDNotPresent;
243  }
244  else if (g_0 < static_cast<GO> (N_R)) {
245  // The GID comes from the initial sequence of R processes.
246  nodeIDs[k] = static_cast<int> (g_0 / static_cast<GO> (N_L + 1));
247  localIDs[k] = static_cast<LO> (g_0 % static_cast<GO> (N_L + 1));
248  }
249  else if (g_0 >= static_cast<GO> (N_R)) {
250  // The GID comes from the remaining P-R processes.
251  const GO g_R = g_0 - static_cast<GO> (N_R);
252  nodeIDs[k] = static_cast<int> (R + g_R / N_L);
253  localIDs[k] = static_cast<int> (g_R % N_L);
254  }
255 #ifdef HAVE_TPETRA_DEBUG
256  else {
257  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
258  "Tpetra::ContiguousUniformDirectory::getEntriesImpl: "
259  "should never get here. "
260  "Please report this bug to the Tpetra developers.");
261  }
262 #endif // HAVE_TPETRA_DEBUG
263  }
264  }
265  else { // don't compute local indices
266  for (size_type k = 0; k < numGids; ++k) {
267  const GO g_0 = globalIDs[k] - g_min;
268  // The first test is a little strange just in case GO is
269  // unsigned. Compilers raise a warning on tests like "x <
270  // 0" if x is unsigned, but don't usually raise a warning if
271  // the expression is a bit more complicated than that.
272  if (g_0 + ONE < ONE || g_0 >= static_cast<GO> (N_G)) {
273  nodeIDs[k] = -1;
274  res = IDNotPresent;
275  }
276  else if (g_0 < static_cast<GO> (N_R)) {
277  // The GID comes from the initial sequence of R processes.
278  nodeIDs[k] = static_cast<int> (g_0 / static_cast<GO> (N_L + 1));
279  }
280  else if (g_0 >= static_cast<GO> (N_R)) {
281  // The GID comes from the remaining P-R processes.
282  const GO g_R = g_0 - static_cast<GO> (N_R);
283  nodeIDs[k] = static_cast<int> (R + g_R / N_L);
284  }
285 #ifdef HAVE_TPETRA_DEBUG
286  else {
287  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
288  "Tpetra::ContiguousUniformDirectory::getEntriesImpl: "
289  "should never get here. "
290  "Please report this bug to the Tpetra developers.");
291  }
292 #endif // HAVE_TPETRA_DEBUG
293  }
294  }
295  return res;
296  }
297 
298  template<class LO, class GO, class NT>
301  {
302  using Teuchos::gatherAll;
303  using Teuchos::RCP;
304 
305  RCP<const Teuchos::Comm<int> > comm = map.getComm ();
306 
307  TEUCHOS_TEST_FOR_EXCEPTION(! map.isDistributed (), std::invalid_argument,
308  Teuchos::typeName (*this) << " constructor: Map is not distributed.");
309  TEUCHOS_TEST_FOR_EXCEPTION(! map.isContiguous (), std::invalid_argument,
310  Teuchos::typeName (*this) << " constructor: Map is not contiguous.");
311 
312  const int numProcs = comm->getSize ();
313 
314  // Make room for the min global ID on each process, plus one
315  // entry at the end for the "max cap."
316  allMinGIDs_ = arcp<GO> (numProcs + 1);
317  // Get my process' min global ID.
318  GO minMyGID = map.getMinGlobalIndex ();
319  // Gather all of the min global IDs into the first numProcs
320  // entries of allMinGIDs_.
321 
322  // FIXME (mfh 16 Apr 2013) GIANT HACK BELOW
323  //
324  // The purpose of this giant hack is that gatherAll appears to
325  // interpret the "receive count" argument differently than
326  // MPI_Allgather does. Matt Bettencourt reports Valgrind issues
327  // (memcpy with overlapping data) with MpiComm<int>::gatherAll,
328  // which could relate either to this, or to OpenMPI.
329 #ifdef HAVE_MPI
330  MPI_Datatype rawMpiType = MPI_INT;
331  bool useRawMpi = true;
332  if (typeid (GO) == typeid (int)) {
333  rawMpiType = MPI_INT;
334  } else if (typeid (GO) == typeid (long)) {
335  rawMpiType = MPI_LONG;
336  } else {
337  useRawMpi = false;
338  }
339  if (useRawMpi) {
340  using Teuchos::rcp_dynamic_cast;
341  using Teuchos::MpiComm;
342  RCP<const MpiComm<int> > mpiComm =
343  rcp_dynamic_cast<const MpiComm<int> > (comm);
344  // It could be a SerialComm instead, even in an MPI build, so
345  // be sure to check.
346  if (! comm.is_null ()) {
347  MPI_Comm rawMpiComm = * (mpiComm->getRawMpiComm ());
348  const int err =
349  MPI_Allgather (&minMyGID, 1, rawMpiType,
350  allMinGIDs_.getRawPtr (), 1, rawMpiType,
351  rawMpiComm);
352  TEUCHOS_TEST_FOR_EXCEPTION(err != MPI_SUCCESS, std::runtime_error,
353  "Tpetra::DistributedContiguousDirectory: MPI_Allgather failed");
354  } else {
355  gatherAll<int, GO> (*comm, 1, &minMyGID, numProcs, allMinGIDs_.getRawPtr ());
356  }
357  } else {
358  gatherAll<int, GO> (*comm, 1, &minMyGID, numProcs, allMinGIDs_.getRawPtr ());
359  }
360 #else // NOT HAVE_MPI
361  gatherAll<int, GO> (*comm, 1, &minMyGID, numProcs, allMinGIDs_.getRawPtr ());
362 #endif // HAVE_MPI
363  // FIXME (mfh 16 Apr 2013) GIANT HACK ABOVE
364 
365  //gatherAll<int, GO> (*comm, 1, &minMyGID, numProcs, allMinGIDs_.getRawPtr ());
366 
367  // Put the max cap at the end. Adding one lets us write loops
368  // over the global IDs with the usual strict less-than bound.
369  allMinGIDs_[numProcs] = map.getMaxAllGlobalIndex ()
370  + Teuchos::OrdinalTraits<GO>::one ();
371  }
372 
373  template<class LO, class GO, class NT>
374  std::string
376  {
377  std::ostringstream os;
378  os << "DistributedContiguousDirectory"
379  << "<" << Teuchos::TypeNameTraits<LO>::name ()
380  << ", " << Teuchos::TypeNameTraits<GO>::name ()
381  << ", " << Teuchos::TypeNameTraits<NT>::name () << ">";
382  return os.str ();
383  }
384 
385  template<class LO, class GO, class NT>
389  const Teuchos::ArrayView<const GO> &globalIDs,
390  const Teuchos::ArrayView<int> &nodeIDs,
391  const Teuchos::ArrayView<LO> &localIDs,
392  const bool computeLIDs) const
393  {
394  using Teuchos::Array;
395  using Teuchos::ArrayRCP;
396  using Teuchos::ArrayView;
397  using Teuchos::as;
398  using Teuchos::Comm;
399  using Teuchos::RCP;
400 
402  RCP<const Teuchos::Comm<int> > comm = map.getComm ();
403  const int myRank = comm->getRank ();
404 
405  // Map is on one process or is locally replicated.
406  typename ArrayView<int>::iterator procIter = nodeIDs.begin();
407  typename ArrayView<LO>::iterator lidIter = localIDs.begin();
408  typename ArrayView<const GO>::iterator gidIter;
409  for (gidIter = globalIDs.begin(); gidIter != globalIDs.end(); ++gidIter) {
410  if (map.isNodeGlobalElement (*gidIter)) {
411  *procIter++ = myRank;
412  if (computeLIDs) {
413  *lidIter++ = map.getLocalElement (*gidIter);
414  }
415  }
416  else {
417  // Advance the pointers, leaving these values set to invalid
418  procIter++;
419  if (computeLIDs) {
420  lidIter++;
421  }
422  res = IDNotPresent;
423  }
424  }
425  return res;
426  }
427 
428  template<class LO, class GO, class NT>
432  const Teuchos::ArrayView<const GO> &globalIDs,
433  const Teuchos::ArrayView<int> &nodeIDs,
434  const Teuchos::ArrayView<LO> &localIDs,
435  const bool computeLIDs) const
436  {
437  using Teuchos::Array;
438  using Teuchos::ArrayRCP;
439  using Teuchos::ArrayView;
440  using Teuchos::as;
441  using Teuchos::Comm;
442  using Teuchos::RCP;
443 
444  RCP<const Teuchos::Comm<int> > comm = map.getComm ();
445  const int numProcs = comm->getSize ();
446  const global_size_t nOverP = map.getGlobalNumElements () / numProcs;
447  const LO LINVALID = Teuchos::OrdinalTraits<LO>::invalid();
449 
450  // Map is distributed but contiguous.
451  typename ArrayView<int>::iterator procIter = nodeIDs.begin();
452  typename ArrayView<LO>::iterator lidIter = localIDs.begin();
453  typename ArrayView<const GO>::iterator gidIter;
454  for (gidIter = globalIDs.begin(); gidIter != globalIDs.end(); ++gidIter) {
455  LO LID = LINVALID; // Assume not found until proven otherwise
456  int image = -1;
457  GO GID = *gidIter;
458  // Guess uniform distribution and start a little above it
459  // TODO: replace by a binary search
460  int curRank;
461  { // We go through all this trouble to avoid overflow and
462  // signed / unsigned casting mistakes (that were made in
463  // previous versions of this code).
464  const GO one = as<GO> (1);
465  const GO two = as<GO> (2);
466  const GO nOverP_GID = as<GO> (nOverP);
467  const GO lowerBound = GID / std::max(nOverP_GID, one) + two;
468  curRank = as<int>(std::min(lowerBound, as<GO>(numProcs - 1)));
469  }
470  bool found = false;
471  while (curRank >= 0 && curRank < numProcs) {
472  if (allMinGIDs_[curRank] <= GID) {
473  if (GID < allMinGIDs_[curRank + 1]) {
474  found = true;
475  break;
476  }
477  else {
478  curRank++;
479  }
480  }
481  else {
482  curRank--;
483  }
484  }
485  if (found) {
486  image = curRank;
487  LID = as<LO> (GID - allMinGIDs_[image]);
488  }
489  else {
490  res = IDNotPresent;
491  }
492  *procIter++ = image;
493  if (computeLIDs) {
494  *lidIter++ = LID;
495  }
496  }
497  return res;
498  }
499 
500  template<class LO, class GO, class NT>
503  oneToOneResult_ (ONE_TO_ONE_NOT_CALLED_YET), // to be revised below
504  locallyOneToOne_ (true), // to be revised below
505  useHashTables_ (false) // to be revised below
506  {
507  initialize (map, Teuchos::null);
508  }
509 
510  template<class LO, class GO, class NT>
513  const tie_break_type& tie_break) :
514  oneToOneResult_ (ONE_TO_ONE_NOT_CALLED_YET), // to be revised below
515  locallyOneToOne_ (true), // to be revised below
516  useHashTables_ (false) // to be revised below
517  {
518  initialize (map, Teuchos::ptrFromRef (tie_break));
519  }
520 
521  template<class LO, class GO, class NT>
522  void
524  initialize (const map_type& map,
525  Teuchos::Ptr<const tie_break_type> tie_break)
526  {
527  using Teuchos::Array;
528  using Teuchos::ArrayView;
529  using Teuchos::as;
530  using Teuchos::RCP;
531  using Teuchos::rcp;
532  using Teuchos::typeName;
533  using Teuchos::TypeNameTraits;
534  using std::cerr;
535  using std::endl;
536  typedef Array<int>::size_type size_type;
537 
538  // This class' implementation of getEntriesImpl() currently
539  // encodes the following assumptions:
540  //
541  // 1. global_size_t >= GO
542  // 2. global_size_t >= int
543  // 3. global_size_t >= LO
544  //
545  // We check these assumptions here.
546  TEUCHOS_TEST_FOR_EXCEPTION(sizeof(global_size_t) < sizeof(GO),
547  std::logic_error, typeName (*this) << ": sizeof(Tpetra::"
548  "global_size_t) = " << sizeof(global_size_t) << " < sizeof(Global"
549  "Ordinal = " << TypeNameTraits<LO>::name () << ") = " << sizeof(GO)
550  << ".");
551  TEUCHOS_TEST_FOR_EXCEPTION(sizeof(global_size_t) < sizeof(int),
552  std::logic_error, typeName (*this) << ": sizeof(Tpetra::"
553  "global_size_t) = " << sizeof(global_size_t) << " < sizeof(int) = "
554  << sizeof(int) << ".");
555  TEUCHOS_TEST_FOR_EXCEPTION(sizeof(global_size_t) < sizeof(LO),
556  std::logic_error, typeName (*this) << ": sizeof(Tpetra::"
557  "global_size_t) = " << sizeof(global_size_t) << " < sizeof(Local"
558  "Ordinal = " << TypeNameTraits<LO>::name () << ") = " << sizeof(LO)
559  << ".");
560 
561  RCP<const Teuchos::Comm<int> > comm = map.getComm ();
562  const LO LINVALID = Teuchos::OrdinalTraits<LO>::invalid ();
563  const GO minAllGID = map.getMinAllGlobalIndex ();
564  const GO maxAllGID = map.getMaxAllGlobalIndex ();
565 
566  // The "Directory Map" (see below) will have a range of elements
567  // from the minimum to the maximum GID of the user Map, and a
568  // minimum GID of minAllGID from the user Map. It doesn't
569  // actually have to store all those entries, though do beware of
570  // calling getNodeElementList on it (see Bug 5822).
571  const global_size_t numGlobalEntries = maxAllGID - minAllGID + 1;
572 
573  // We can't afford to replicate the whole directory on each
574  // process, so create the "Directory Map", a uniform contiguous
575  // Map that describes how we will distribute the directory over
576  // processes.
577  //
578  // FIXME (mfh 08 May 2012) Here we're setting minAllGID to be
579  // the index base. The index base should be separate from the
580  // minimum GID.
581  directoryMap_ = rcp (new map_type (numGlobalEntries, minAllGID, comm,
582  GloballyDistributed, map.getNode ()));
583  // The number of Directory elements that my process owns.
584  const size_t dir_numMyEntries = directoryMap_->getNodeNumElements ();
585 
586  // Fix for Bug 5822: If the input Map is "sparse," that is if
587  // the difference between the global min and global max GID is
588  // much larger than the global number of elements in the input
589  // Map, then it's possible that the Directory Map might have
590  // many more entries than the input Map on this process. This
591  // can cause memory scalability issues. In that case, we switch
592  // from the array-based implementation of Directory storage to
593  // the hash table - based implementation. We don't use hash
594  // tables all the time, because they are slower in the common
595  // case of a nonsparse Map.
596  //
597  // NOTE: This is a per-process decision. Some processes may use
598  // array-based storage, whereas others may use hash table -
599  // based storage.
600 
601  // A hash table takes a constant factor more space, more or
602  // less, than an array. Thus, it's not worthwhile, even in
603  // terms of memory usage, always to use a hash table.
604  // Furthermore, array lookups are faster than hash table
605  // lookups, so it may be worthwhile to use an array even if it
606  // takes more space. The "sparsity threshold" governs when to
607  // switch to a hash table - based implementation.
608  const size_t inverseSparsityThreshold = 10;
609  useHashTables_ =
610  dir_numMyEntries >= inverseSparsityThreshold * map.getNodeNumElements ();
611 
612  // Get list of process IDs that own the directory entries for the
613  // Map GIDs. These will be the targets of the sends that the
614  // Distributor will do.
615  const int myRank = comm->getRank ();
616  const size_t numMyEntries = map.getNodeNumElements ();
617  Array<int> sendImageIDs (numMyEntries);
618  ArrayView<const GO> myGlobalEntries = map.getNodeElementList ();
619  // An ID not present in this lookup indicates that it lies outside
620  // of the range [minAllGID,maxAllGID] (from map_). this means
621  // something is wrong with map_, our fault.
622  const LookupStatus lookupStatus =
623  directoryMap_->getRemoteIndexList (myGlobalEntries, sendImageIDs);
624  TEUCHOS_TEST_FOR_EXCEPTION(
625  lookupStatus == IDNotPresent, std::logic_error, Teuchos::typeName(*this)
626  << " constructor: the Directory Map could not find out where one or "
627  "more of my Map's indices should go. The input to getRemoteIndexList "
628  "is " << Teuchos::toString (myGlobalEntries) << ", and the output is "
629  << Teuchos::toString (sendImageIDs ()) << ". The input Map itself has "
630  "the following entries on the calling process " <<
631  map.getComm ()->getRank () << ": " <<
632  Teuchos::toString (map.getNodeElementList ()) << ", and has "
633  << map.getGlobalNumElements () << " total global indices in ["
634  << map.getMinAllGlobalIndex () << "," << map.getMaxAllGlobalIndex ()
635  << "]. The Directory Map has "
636  << directoryMap_->getGlobalNumElements () << " total global indices in "
637  "[" << directoryMap_->getMinAllGlobalIndex () << "," <<
638  directoryMap_->getMaxAllGlobalIndex () << "], and the calling process "
639  "has GIDs [" << directoryMap_->getMinGlobalIndex () << "," <<
640  directoryMap_->getMaxGlobalIndex () << "]. "
641  "This probably means there is a bug in Map or Directory. "
642  "Please report this bug to the Tpetra developers.");
643 
644  // Initialize the distributor using the list of process IDs to
645  // which to send. We'll use the distributor to send out triples
646  // of (GID, process ID, LID). We're sending the entries to the
647  // processes that the Directory Map says should own them, which is
648  // why we called directoryMap_->getRemoteIndexList() above.
649  Distributor distor (comm);
650  const size_t numReceives = distor.createFromSends (sendImageIDs);
651 
652  // NOTE (mfh 21 Mar 2012) The following code assumes that
653  // sizeof(GO) >= sizeof(int) and sizeof(GO) >= sizeof(LO).
654  //
655  // Create and fill buffer of (GID, PID, LID) triples to send
656  // out. We pack the (GID, PID, LID) triples into a single Array
657  // of GO, casting the PID from int to GO and the LID from LO to
658  // GO as we do so.
659  //
660  // FIXME (mfh 23 Mar 2014) This assumes that sizeof(LO) <=
661  // sizeof(GO) and sizeof(int) <= sizeof(GO). The former is
662  // required, and the latter is generally the case, but we should
663  // still check for this.
664  const int packetSize = 3; // We're sending triples, so packet size is 3.
665  Array<GO> exportEntries (packetSize * numMyEntries); // data to send out
666  {
667  size_type exportIndex = 0;
668  for (size_type i = 0; i < static_cast<size_type> (numMyEntries); ++i) {
669  exportEntries[exportIndex++] = myGlobalEntries[i];
670  exportEntries[exportIndex++] = as<GO> (myRank);
671  exportEntries[exportIndex++] = as<GO> (i);
672  }
673  }
674  // Buffer of data to receive. The Distributor figured out for
675  // us how many packets we're receiving, when we called its
676  // createFromSends() method to set up the distribution plan.
677  Array<GO> importElements (packetSize * distor.getTotalReceiveLength ());
678 
679  // Distribute the triples of (GID, process ID, LID).
680  distor.doPostsAndWaits (exportEntries ().getConst (), packetSize, importElements ());
681 
682  // Unpack the redistributed data. Both implementations of
683  // Directory storage map from an LID in the Directory Map (which
684  // is the LID of the GID to store) to either a PID or an LID in
685  // the input Map. Each "packet" (contiguous chunk of
686  // importElements) contains a triple: (GID, PID, LID).
687  if (useHashTables_) {
688  // Create the hash tables. We know exactly how many elements
689  // to expect in each hash table. FixedHashTable's constructor
690  // currently requires all the keys and values at once, so we
691  // have to extract them in temporary arrays. It may be
692  // possible to rewrite FixedHashTable to use a "start fill" /
693  // "end fill" approach that avoids the temporary arrays, but
694  // we won't try that for now.
695 
696  // The constructors of Array and ArrayRCP that take a number
697  // of elements all initialize the arrays. Instead, allocate
698  // raw arrays, then hand them off to ArrayRCP, to avoid the
699  // initial unnecessary initialization without losing the
700  // benefit of exception safety (and bounds checking, in a
701  // debug build).
702  LO* tableKeysRaw = NULL;
703  LO* tableLidsRaw = NULL;
704  int* tablePidsRaw = NULL;
705  try {
706  tableKeysRaw = new LO [numReceives];
707  tableLidsRaw = new LO [numReceives];
708  tablePidsRaw = new int [numReceives];
709  } catch (...) {
710  if (tableKeysRaw != NULL) {
711  delete [] tableKeysRaw;
712  }
713  if (tableLidsRaw != NULL) {
714  delete [] tableLidsRaw;
715  }
716  if (tablePidsRaw != NULL) {
717  delete [] tablePidsRaw;
718  }
719  throw;
720  }
721  ArrayRCP<LO> tableKeys (tableKeysRaw, 0, numReceives, true);
722  ArrayRCP<LO> tableLids (tableLidsRaw, 0, numReceives, true);
723  ArrayRCP<int> tablePids (tablePidsRaw, 0, numReceives, true);
724 
725  if (tie_break.is_null ()) {
726  // Fill the temporary arrays of keys and values.
727  size_type importIndex = 0;
728  for (size_type i = 0; i < static_cast<size_type> (numReceives); ++i) {
729  const GO curGID = importElements[importIndex++];
730  const LO curLID = directoryMap_->getLocalElement (curGID);
731  TEUCHOS_TEST_FOR_EXCEPTION(
732  curLID == LINVALID, std::logic_error,
733  Teuchos::typeName(*this) << " constructor: Incoming global index "
734  << curGID << " does not have a corresponding local index in the "
735  "Directory Map. Please report this bug to the Tpetra developers.");
736  tableKeys[i] = curLID;
737  tablePids[i] = importElements[importIndex++];
738  tableLids[i] = importElements[importIndex++];
739  }
740  // Set up the hash tables. The hash tables' constructor
741  // detects whether there are duplicates, so that we can set
742  // locallyOneToOne_.
743  typedef Kokkos::Device<typename NT::execution_space,
744  typename NT::memory_space> DT;
745  lidToPidTable_ =
746  rcp (new Details::FixedHashTable<LO, int, DT> (tableKeys (),
747  tablePids ()));
748  locallyOneToOne_ = ! (lidToPidTable_->hasDuplicateKeys ());
749  lidToLidTable_ =
750  rcp (new Details::FixedHashTable<LO, LO, DT> (tableKeys (),
751  tableLids ()));
752  }
753  else { // tie_break is NOT null
754 
755  // For each directory Map LID received, collect all the
756  // corresponding (PID,LID) pairs. If the input Map is not
757  // one-to-one, corresponding directory Map LIDs will have
758  // more than one pair. In that case, we will use the
759  // TieBreak object to pick exactly one pair.
760  typedef std::map<LO, std::vector<std::pair<int, LO> > > pair_table_type;
761  pair_table_type ownedPidLidPairs;
762 
763  // For each directory Map LID received, collect the zero or
764  // more input Map (PID,LID) pairs into ownedPidLidPairs.
765  size_type importIndex = 0;
766  for (size_type i = 0; i < static_cast<size_type> (numReceives); ++i) {
767  const GO curGID = importElements[importIndex++];
768  const LO dirMapLid = directoryMap_->getLocalElement (curGID);
769  TEUCHOS_TEST_FOR_EXCEPTION(
770  dirMapLid == LINVALID, std::logic_error,
771  Teuchos::typeName(*this) << " constructor: Incoming global index "
772  << curGID << " does not have a corresponding local index in the "
773  "Directory Map. Please report this bug to the Tpetra developers.");
774  tableKeys[i] = dirMapLid;
775  const int PID = importElements[importIndex++];
776  const int LID = importElements[importIndex++];
777 
778  // These may change below. We fill them in just to ensure
779  // that they won't have invalid values.
780  tablePids[i] = PID;
781  tableLids[i] = LID;
782 
783  // For every directory Map LID, we have to remember all
784  // (PID, LID) pairs. The TieBreak object will arbitrate
785  // between them in the loop below.
786  ownedPidLidPairs[dirMapLid].push_back (std::make_pair (PID, LID));
787  }
788 
789  // Use TieBreak to arbitrate between (PID,LID) pairs
790  // corresponding to each directory Map LID.
791  //
792  // FIXME (mfh 23 Mar 2014) How do I know that i is the same
793  // as the directory Map LID?
794  const size_type numPairs =
795  static_cast<size_type> (ownedPidLidPairs.size ());
796  for (size_type i = 0; i < numPairs; ++i) {
797  const LO dirMapLid = static_cast<LO> (i);
798  const GO dirMapGid = directoryMap_->getGlobalElement (dirMapLid);
799  const std::vector<std::pair<int, LO> >& pidLidList =
800  ownedPidLidPairs[i];
801  const size_t listLen = pidLidList.size ();
802  if (listLen > 0) {
803  if (listLen > 1) {
804  locallyOneToOne_ = false;
805  }
806  // If there is some (PID,LID) pair for the current input
807  // Map LID, then it makes sense to invoke the TieBreak
808  // object to arbitrate between the options. Even if
809  // there is only one (PID,LID) pair, we still want to
810  // give the TieBreak object a chance to do whatever it
811  // likes to do, in terms of side effects (e.g., track
812  // (PID,LID) pairs).
813  const size_type index =
814  static_cast<size_type> (tie_break->selectedIndex (dirMapGid,
815  pidLidList));
816  tablePids[i] = pidLidList[index].first;
817  tableLids[i] = pidLidList[index].second;
818  }
819  }
820 
821  // Set up the hash tables.
822  typedef Kokkos::Device<typename NT::execution_space,
823  typename NT::memory_space> DT;
824  lidToPidTable_ =
825  rcp (new Details::FixedHashTable<LO, int, DT> (tableKeys (),
826  tablePids ()));
827  lidToLidTable_ =
828  rcp (new Details::FixedHashTable<LO, LO, DT> (tableKeys (),
829  tableLids ()));
830  }
831  }
832  else {
833  if (tie_break.is_null ()) {
834  // Use array-based implementation of Directory storage.
835  // Allocate these arrays and fill them with invalid values,
836  // in case the input Map's GID list is sparse (i.e., does
837  // not populate all GIDs from minAllGID to maxAllGID).
838  PIDs_ = arcp<int> (dir_numMyEntries);
839  std::fill (PIDs_.begin (), PIDs_.end (), -1);
840  LIDs_ = arcp<LO> (dir_numMyEntries);
841  std::fill (LIDs_.begin (), LIDs_.end (), LINVALID);
842  // Fill in the arrays with PIDs resp. LIDs.
843  size_type importIndex = 0;
844  for (size_type i = 0; i < static_cast<size_type> (numReceives); ++i) {
845  const GO curGID = importElements[importIndex++];
846  const LO curLID = directoryMap_->getLocalElement (curGID);
847  TEUCHOS_TEST_FOR_EXCEPTION(curLID == LINVALID, std::logic_error,
848  Teuchos::typeName(*this) << " constructor: Incoming global index "
849  << curGID << " does not have a corresponding local index in the "
850  "Directory Map. Please report this bug to the Tpetra developers.");
851 
852  // If PIDs_[curLID] is not -1, then curGID is a duplicate
853  // on the calling process, so the Directory is not locally
854  // one-to-one.
855  if (PIDs_[curLID] != -1) {
856  locallyOneToOne_ = false;
857  }
858  PIDs_[curLID] = importElements[importIndex++];
859  LIDs_[curLID] = importElements[importIndex++];
860  }
861  }
862  else {
863  PIDs_ = arcp<int> (dir_numMyEntries);
864  LIDs_ = arcp<LO> (dir_numMyEntries);
865  std::fill (PIDs_.begin (), PIDs_.end (), -1);
866 
867  // All received (PID, LID) pairs go into ownedPidLidPairs.
868  // This is a map from the directory Map's LID to the (PID,
869  // LID) pair (where the latter LID comes from the input Map,
870  // not the directory Map). If the input Map is not
871  // one-to-one, corresponding LIDs will have
872  // ownedPidLidPairs[curLID].size() > 1. In that case, we
873  // will use the TieBreak object to pick exactly one pair.
874  Array<std::vector<std::pair<int, LO> > > ownedPidLidPairs (dir_numMyEntries);
875  size_type importIndex = 0;
876  for (size_type i = 0; i < static_cast<size_type> (numReceives); ++i) {
877  const GO GID = importElements[importIndex++];
878  const int PID = importElements[importIndex++];
879  const LO LID = importElements[importIndex++];
880 
881  const LO dirMapLid = directoryMap_->getLocalElement (GID);
882  TEUCHOS_TEST_FOR_EXCEPTION(
883  dirMapLid == LINVALID, std::logic_error,
884  Teuchos::typeName(*this) << " constructor: Incoming global index "
885  << GID << " does not have a corresponding local index in the "
886  "Directory Map. Please report this bug to the Tpetra developers.");
887  ownedPidLidPairs[dirMapLid].push_back (std::make_pair (PID, LID));
888  }
889 
890  // Use TieBreak to arbitrate between (PID,LID) pairs
891  // corresponding to each directory Map LID.
892  //
893  // FIXME (mfh 23 Mar 2014) How do I know that i is the same
894  // as the directory Map LID?
895  const size_type numPairs =
896  static_cast<size_type> (ownedPidLidPairs.size ());
897  for (size_type i = 0; i < numPairs; ++i) {
898  const LO dirMapLid = static_cast<LO> (i);
899  const GO dirMapGid = directoryMap_->getGlobalElement (dirMapLid);
900  const std::vector<std::pair<int, LO> >& pidLidList =
901  ownedPidLidPairs[i];
902  const size_t listLen = pidLidList.size ();
903  if (listLen > 0) {
904  if (listLen > 1) {
905  locallyOneToOne_ = false;
906  }
907  // If there is some (PID,LID) pair for the current input
908  // Map LID, then it makes sense to invoke the TieBreak
909  // object to arbitrate between the options. Even if
910  // there is only one (PID,LID) pair, we still want to
911  // give the TieBreak object a chance to do whatever it
912  // likes to do, in terms of side effects (e.g., track
913  // (PID,LID) pairs).
914  const size_type index =
915  static_cast<size_type> (tie_break->selectedIndex (dirMapGid,
916  pidLidList));
917  PIDs_[i] = pidLidList[index].first;
918  LIDs_[i] = pidLidList[index].second;
919  }
920  // else no GID specified by source map
921  }
922  }
923  }
924  }
925 
926  template<class LO, class GO, class NT>
927  std::string
929  {
930  std::ostringstream os;
931  os << "DistributedNoncontiguousDirectory"
932  << "<" << Teuchos::TypeNameTraits<LO>::name ()
933  << ", " << Teuchos::TypeNameTraits<GO>::name ()
934  << ", " << Teuchos::TypeNameTraits<NT>::name () << ">";
935  return os.str ();
936  }
937 
938  template<class LO, class GO, class NT>
942  const Teuchos::ArrayView<const GO> &globalIDs,
943  const Teuchos::ArrayView<int> &nodeIDs,
944  const Teuchos::ArrayView<LO> &localIDs,
945  const bool computeLIDs) const
946  {
947  using Teuchos::Array;
948  using Teuchos::ArrayRCP;
949  using Teuchos::ArrayView;
950  using Teuchos::as;
951  using Teuchos::RCP;
952  using std::cerr;
953  using std::endl;
954  typedef typename Array<GO>::size_type size_type;
955 
956  RCP<const Teuchos::Comm<int> > comm = map.getComm ();
957  const size_t numEntries = globalIDs.size ();
958  const LO LINVALID = Teuchos::OrdinalTraits<LO>::invalid();
960 
961  //
962  // Set up directory structure.
963  //
964 
965  // If we're computing LIDs, we also have to include them in each
966  // packet, along with the GID and process ID.
967  const int packetSize = computeLIDs ? 3 : 2;
968 
969  // For data distribution, we use: Surprise! A Distributor!
970  Distributor distor (comm);
971 
972  // Get directory locations for the requested list of entries.
973  Array<int> dirImages (numEntries);
974  res = directoryMap_->getRemoteIndexList (globalIDs, dirImages ());
975  // Check for unfound globalIDs and set corresponding nodeIDs to -1
976  size_t numMissing = 0;
977  if (res == IDNotPresent) {
978  for (size_t i=0; i < numEntries; ++i) {
979  if (dirImages[i] == -1) {
980  nodeIDs[i] = -1;
981  if (computeLIDs) {
982  localIDs[i] = LINVALID;
983  }
984  numMissing++;
985  }
986  }
987  }
988 
989  Array<GO> sendGIDs;
990  Array<int> sendImages;
991  distor.createFromRecvs (globalIDs, dirImages (), sendGIDs, sendImages);
992  const size_type numSends = sendGIDs.size ();
993 
994  //
995  // mfh 13 Nov 2012:
996  //
997  // The code below temporarily stores LO, GO, and int values in
998  // an array of global_size_t. If one of the signed types (LO
999  // and GO should both be signed) happened to be -1 (or some
1000  // negative number, but -1 is the one that came up today), then
1001  // conversion to global_size_t will result in a huge
1002  // global_size_t value, and thus conversion back may overflow.
1003  // (Teuchos::as doesn't know that we meant it to be an LO or GO
1004  // all along.)
1005  //
1006  // The overflow normally would not be a problem, since it would
1007  // just go back to -1 again. However, Teuchos::as does range
1008  // checking on conversions in a debug build, so it throws an
1009  // exception (std::range_error) in this case. Range checking is
1010  // generally useful in debug mode, so we don't want to disable
1011  // this behavior globally.
1012  //
1013  // We solve this problem by forgoing use of Teuchos::as for the
1014  // conversions below from LO, GO, or int to global_size_t, and
1015  // the later conversions back from global_size_t to LO, GO, or
1016  // int.
1017  //
1018  // I've recorded this discussion as Bug 5760.
1019  //
1020 
1021  // global_size_t >= GO
1022  // global_size_t >= size_t >= int
1023  // global_size_t >= size_t >= LO
1024  // Therefore, we can safely store all of these in a global_size_t
1025  Array<global_size_t> exports (packetSize * numSends);
1026  {
1027  // Packet format:
1028  // - If computing LIDs: (GID, PID, LID)
1029  // - Otherwise: (GID, PID)
1030  //
1031  // "PID" means "process ID" (a.k.a. "node ID," a.k.a. "rank").
1032 
1033  // Current position to which to write in exports array. If
1034  // sending pairs, we pack the (GID, PID) pair for gid =
1035  // sendGIDs[k] in exports[2*k], exports[2*k+1]. If sending
1036  // triples, we pack the (GID, PID, LID) pair for gid =
1037  // sendGIDs[k] in exports[3*k, 3*k+1, 3*k+2].
1038  size_type exportsIndex = 0;
1039 
1040  if (useHashTables_) {
1041  for (size_type gidIndex = 0; gidIndex < numSends; ++gidIndex) {
1042  const GO curGID = sendGIDs[gidIndex];
1043  // Don't use as() here (see above note).
1044  exports[exportsIndex++] = static_cast<global_size_t> (curGID);
1045  const LO curLID = directoryMap_->getLocalElement (curGID);
1046  TEUCHOS_TEST_FOR_EXCEPTION(curLID == LINVALID, std::logic_error,
1047  Teuchos::typeName (*this) << "::getEntriesImpl(): The Directory "
1048  "Map's global index " << curGID << " does not have a corresponding "
1049  "local index. Please report this bug to the Tpetra developers.");
1050  // Don't use as() here (see above note).
1051  exports[exportsIndex++] = static_cast<global_size_t> (lidToPidTable_->get (curLID));
1052  if (computeLIDs) {
1053  // Don't use as() here (see above note).
1054  exports[exportsIndex++] = static_cast<global_size_t> (lidToLidTable_->get (curLID));
1055  }
1056  }
1057  } else {
1058  for (size_type gidIndex = 0; gidIndex < numSends; ++gidIndex) {
1059  const GO curGID = sendGIDs[gidIndex];
1060  // Don't use as() here (see above note).
1061  exports[exportsIndex++] = static_cast<global_size_t> (curGID);
1062  const LO curLID = directoryMap_->getLocalElement (curGID);
1063  TEUCHOS_TEST_FOR_EXCEPTION(curLID == LINVALID, std::logic_error,
1064  Teuchos::typeName (*this) << "::getEntriesImpl(): The Directory "
1065  "Map's global index " << curGID << " does not have a corresponding "
1066  "local index. Please report this bug to the Tpetra developers.");
1067  // Don't use as() here (see above note).
1068  exports[exportsIndex++] = static_cast<global_size_t> (PIDs_[curLID]);
1069  if (computeLIDs) {
1070  // Don't use as() here (see above note).
1071  exports[exportsIndex++] = static_cast<global_size_t> (LIDs_[curLID]);
1072  }
1073  }
1074  }
1075 
1076  TEUCHOS_TEST_FOR_EXCEPTION(
1077  exportsIndex > exports.size (), std::logic_error,
1078  Teuchos::typeName (*this) << "::getEntriesImpl(): On Process " <<
1079  comm->getRank () << ", exportsIndex = " << exportsIndex <<
1080  " > exports.size() = " << exports.size () <<
1081  ". Please report this bug to the Tpetra developers.");
1082  }
1083 
1084  TEUCHOS_TEST_FOR_EXCEPTION(
1085  numEntries < numMissing, std::logic_error,
1086  Teuchos::typeName (*this) << "::getEntriesImpl(): On Process "
1087  << comm->getRank () << ", numEntries = " << numEntries
1088  << " < numMissing = " << numMissing
1089  << ". Please report this bug to the Tpetra developers.");
1090 
1091  //
1092  // mfh 13 Nov 2012: See note above on conversions between
1093  // global_size_t and LO, GO, or int.
1094  //
1095  const size_t numRecv = numEntries - numMissing;
1096 
1097  {
1098  const size_t importLen = packetSize * distor.getTotalReceiveLength ();
1099  const size_t requiredImportLen = numRecv * packetSize;
1100  const int myRank = comm->getRank ();
1101  TEUCHOS_TEST_FOR_EXCEPTION(
1102  importLen < requiredImportLen, std::logic_error,
1103  "Tpetra::Details::DistributedNoncontiguousDirectory::getEntriesImpl: "
1104  "On Process " << myRank << ": The 'imports' array must have length "
1105  "at least " << requiredImportLen << ", but its actual length is " <<
1106  importLen << ". numRecv: " << numRecv << ", packetSize: " <<
1107  packetSize << ", numEntries (# GIDs): " << numEntries <<
1108  ", numMissing: " << numMissing << ": distor.getTotalReceiveLength(): "
1109  << distor.getTotalReceiveLength () << ". " << std::endl <<
1110  "Distributor description: " << distor.description () << ". "
1111  << std::endl <<
1112  "Please report this bug to the Tpetra developers.");
1113  }
1114 
1115  Array<global_size_t> imports (packetSize * distor.getTotalReceiveLength ());
1116  // FIXME (mfh 20 Mar 2014) One could overlap the sort2() below
1117  // with communication, by splitting this call into doPosts and
1118  // doWaits. The code is still correct in this form, however.
1119  distor.doPostsAndWaits (exports ().getConst (), packetSize, imports ());
1120 
1121  Array<GO> sortedIDs (globalIDs); // deep copy (for later sorting)
1122  Array<GO> offset (numEntries); // permutation array (sort2 output)
1123  for (GO ii = 0; ii < static_cast<GO> (numEntries); ++ii) {
1124  offset[ii] = ii;
1125  }
1126  sort2 (sortedIDs.begin(), sortedIDs.begin() + numEntries, offset.begin());
1127 
1128  size_t importsIndex = 0;
1129  //typename Array<global_size_t>::iterator ptr = imports.begin();
1130  typedef typename Array<GO>::iterator IT;
1131 
1132  // we know these conversions are in range, because we loaded this data
1133  for (size_t i = 0; i < numRecv; ++i) {
1134  // Don't use as() here (see above note).
1135  const GO curGID = static_cast<GO> (imports[importsIndex++]);
1136  std::pair<IT, IT> p1 = std::equal_range (sortedIDs.begin(), sortedIDs.end(), curGID);
1137  if (p1.first != p1.second) {
1138  const size_t j = p1.first - sortedIDs.begin();
1139  // Don't use as() here (see above note).
1140  nodeIDs[offset[j]] = static_cast<int> (imports[importsIndex++]);
1141  if (computeLIDs) {
1142  // Don't use as() here (see above note).
1143  localIDs[offset[j]] = static_cast<LO> (imports[importsIndex++]);
1144  }
1145  if (nodeIDs[offset[j]] == -1) {
1146  res = IDNotPresent;
1147  }
1148  }
1149  }
1150 
1151  TEUCHOS_TEST_FOR_EXCEPTION(
1152  static_cast<size_t> (importsIndex) > static_cast<size_t> (imports.size ()),
1153  std::logic_error,
1154  "Tpetra::Details::DistributedNoncontiguousDirectory::getEntriesImpl: "
1155  "On Process " << comm->getRank () << ": importsIndex = " <<
1156  importsIndex << " > imports.size() = " << imports.size () << ". "
1157  "numRecv: " << numRecv << ", packetSize: " << packetSize << ", "
1158  "numEntries (# GIDs): " << numEntries << ", numMissing: " << numMissing
1159  << ": distor.getTotalReceiveLength(): "
1160  << distor.getTotalReceiveLength () << ". Please report this bug to "
1161  "the Tpetra developers.");
1162 
1163  return res;
1164  }
1165 
1166 
1167  template<class LO, class GO, class NT>
1168  bool
1170  isOneToOne (const Teuchos::Comm<int>& comm) const
1171  {
1172  if (oneToOneResult_ == ONE_TO_ONE_NOT_CALLED_YET) {
1173  const int lcl121 = isLocallyOneToOne () ? 1 : 0;
1174  int gbl121 = 0;
1175  Teuchos::reduceAll<int, int> (comm, Teuchos::REDUCE_MIN, lcl121,
1176  Teuchos::outArg (gbl121));
1177  oneToOneResult_ = (gbl121 == 1) ? ONE_TO_ONE_TRUE : ONE_TO_ONE_FALSE;
1178  }
1179  return (oneToOneResult_ == ONE_TO_ONE_TRUE);
1180  }
1181  } // namespace Details
1182 } // namespace Tpetra
1183 
1184 //
1185 // Explicit instantiation macro
1186 //
1187 // Must be expanded from within the Tpetra::Details namespace!
1188 //
1189 #define TPETRA_DIRECTORY_IMPL_INSTANT(LO,GO,NODE) \
1190  template class Directory< LO , GO , NODE >; \
1191  template class ReplicatedDirectory< LO , GO , NODE >; \
1192  template class ContiguousUniformDirectory< LO, GO, NODE >; \
1193  template class DistributedContiguousDirectory< LO , GO , NODE >; \
1194  template class DistributedNoncontiguousDirectory< LO , GO , NODE >; \
1195 
1196 #endif // __Tpetra_DirectoryImpl_def_hpp
Interface for breaking ties in ownership.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
GlobalOrdinal getMinGlobalIndex() const
The minimum global index owned by the calling process.
Interface for breaking ties in ownership.
std::string description() const
A one-line human-readable description of this object.
std::string description() const
A simple one-line description of this object.
bool isNodeGlobalElement(GlobalOrdinal globalIndex) const
Whether the given global index is owned by this Map on the calling process.
bool isContiguous() const
True if this Map is distributed contiguously, else false.
LookupStatus
Return status of Map remote index lookup (getRemoteIndexList()).
LookupStatus getEntriesImpl(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs, const Teuchos::ArrayView< LocalOrdinal > &localIDs, const bool computeLIDs) const
Find process IDs and (optionally) local IDs for the given global IDs.
bool isUniform() const
Whether the range of global indices is uniform.
Implementation of Directory for a distributed noncontiguous Map.
bool isDistributed() const
Whether this Map is globally distributed or locally replicated.
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Accessors for the Teuchos::Comm and Kokkos Node objects.
Implementation of Directory for a distributed contiguous Map.
GlobalOrdinal getMinAllGlobalIndex() const
The minimum global index over all processes in the communicator.
virtual bool isOneToOne(const Teuchos::Comm< int > &comm) const
Whether the Directory&#39;s input Map is (globally) one to one.
ReplicatedDirectory()
Constructor (that takes no arguments).
Implementation details of Tpetra.
GlobalOrdinal getMaxAllGlobalIndex() const
The maximum global index over all processes in the communicator.
virtual bool isOneToOne(const Teuchos::Comm< int > &comm) const
Whether the Directory&#39;s input Map is (globally) one to one.
size_t global_size_t
Global size_t object.
LookupStatus getEntriesImpl(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs, const Teuchos::ArrayView< LocalOrdinal > &localIDs, const bool computeLIDs) const
Find process IDs and (optionally) local IDs for the given global IDs.
void doPostsAndWaits(const ArrayView< const Packet > &exports, size_t numPackets, const ArrayView< Packet > &imports)
Execute the (forward) communication plan.
LookupStatus getEntries(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs, const Teuchos::ArrayView< LocalOrdinal > &localIDs, const bool computeLIDs) const
Sets up and executes a communication plan for a Tpetra DistObject.
size_t getTotalReceiveLength() const
Total number of values this process will receive from other processes.
Teuchos::RCP< Node > getNode() const
Get this Map&#39;s Node object.
LookupStatus getEntriesImpl(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs, const Teuchos::ArrayView< LocalOrdinal > &localIDs, const bool computeLIDs) const
Find process IDs and (optionally) local IDs for the given global IDs.
std::string description() const
A one-line human-readable description of this object.
std::string description() const
A one-line human-readable description of this object.
LocalOrdinal getLocalElement(GlobalOrdinal globalIndex) const
The local index corresponding to the given global index.
Teuchos::ArrayView< const GlobalOrdinal > getNodeElementList() const
Return a view of the global indices owned by this process.
void sort2(const IT1 &first1, const IT1 &last1, const IT2 &first2)
Sort the first array, and apply the resulting permutation to the second array.
size_t getNodeNumElements() const
The number of elements belonging to the calling process.
Describes a parallel distribution of objects over processes.
std::string description() const
A one-line human-readable description of this object.
Declaration of implementation details of Tpetra::Directory.
Implementation of Directory for a contiguous, uniformly distributed Map.
LookupStatus getEntriesImpl(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs, const Teuchos::ArrayView< LocalOrdinal > &localIDs, const bool computeLIDs) const
Find process IDs and (optionally) local IDs for the given global IDs.
void createFromRecvs(const ArrayView< const Ordinal > &remoteIDs, const ArrayView< const int > &remoteNodeIDs, Array< Ordinal > &exportIDs, Array< int > &exportNodeIDs)
Set up Distributor using list of process ranks from which to receive.
global_size_t getGlobalNumElements() const
The number of elements in this Map.