Tpetra parallel linear algebra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Tpetra_Map_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 // ************************************************************************
38 // @HEADER
39 
44 
45 #ifndef TPETRA_MAP_DEF_HPP
46 #define TPETRA_MAP_DEF_HPP
47 
48 #include "Tpetra_Directory.hpp" // must include for implicit instantiation to work
50 #include "Tpetra_Details_FixedHashTable.hpp"
53 #include "Tpetra_Core.hpp"
54 #include "Tpetra_Util.hpp"
55 #include "Teuchos_as.hpp"
56 #include "Teuchos_TypeNameTraits.hpp"
57 #include "Teuchos_CommHelpers.hpp"
58 #include "Tpetra_Details_mpiIsInitialized.hpp"
59 #include "Tpetra_Details_extractMpiCommFromTeuchos.hpp" // teuchosCommIsAnMpiComm
61 #include <memory>
62 #include <sstream>
63 #include <stdexcept>
64 #include <typeinfo>
65 
66 namespace { // (anonymous)
67 
68  template<class ExecutionSpace>
69  void
70  checkMapInputArray (const char ctorName[],
71  const void* indexList,
72  const size_t indexListSize,
73  const ExecutionSpace& execSpace,
74  const Teuchos::Comm<int>* const comm)
75  {
77 
78  const bool debug = Behavior::debug("Map");
79  if (debug) {
80  using Teuchos::outArg;
81  using Teuchos::REDUCE_MIN;
82  using Teuchos::reduceAll;
83  using std::endl;
84 
85  const int myRank = comm == nullptr ? 0 : comm->getRank ();
86  const bool verbose = Behavior::verbose("Map");
87  std::ostringstream lclErrStrm;
88  int lclSuccess = 1;
89 
90  if (indexListSize != 0 && indexList == nullptr) {
91  lclSuccess = 0;
92  if (verbose) {
93  lclErrStrm << "Proc " << myRank << ": indexList is null, "
94  "but indexListSize=" << indexListSize << " != 0." << endl;
95  }
96  }
97  int gblSuccess = 0; // output argument
98  reduceAll (*comm, REDUCE_MIN, lclSuccess, outArg (gblSuccess));
99  if (gblSuccess != 1) {
100  std::ostringstream gblErrStrm;
101  gblErrStrm << "Tpetra::Map constructor " << ctorName <<
102  " detected a problem with the input array "
103  "(raw array, Teuchos::ArrayView, or Kokkos::View) "
104  "of global indices." << endl;
105  if (verbose) {
106  using ::Tpetra::Details::gathervPrint;
107  gathervPrint (gblErrStrm, lclErrStrm.str (), *comm);
108  }
109  TEUCHOS_TEST_FOR_EXCEPTION
110  (true, std::invalid_argument, gblErrStrm.str ());
111  }
112  }
113  }
114 } // namespace (anonymous)
115 
116 namespace Tpetra {
117 
118  template <class LocalOrdinal, class GlobalOrdinal, class Node>
120  Map () :
121  comm_ (new Teuchos::SerialComm<int> ()),
122  indexBase_ (0),
123  numGlobalElements_ (0),
124  numLocalElements_ (0),
125  minMyGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
126  maxMyGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
127  minAllGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
128  maxAllGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
129  firstContiguousGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
130  lastContiguousGID_ (Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ()),
131  uniform_ (false), // trivially
132  contiguous_ (false),
133  distributed_ (false), // no communicator yet
134  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
135  {
137  }
138 
139 
140  template <class LocalOrdinal, class GlobalOrdinal, class Node>
142  Map (const global_size_t numGlobalElements,
143  const global_ordinal_type indexBase,
144  const Teuchos::RCP<const Teuchos::Comm<int> > &comm,
145  const LocalGlobal lOrG) :
146  comm_ (comm),
147  uniform_ (true),
148  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
149  {
150  using Teuchos::as;
151  using Teuchos::broadcast;
152  using Teuchos::outArg;
153  using Teuchos::reduceAll;
154  using Teuchos::REDUCE_MIN;
155  using Teuchos::REDUCE_MAX;
156  using Teuchos::typeName;
157  using std::endl;
158  using GO = global_ordinal_type;
159  using GST = global_size_t;
160  const GST GSTI = Tpetra::Details::OrdinalTraits<GST>::invalid ();
161  const char funcName[] = "Map(gblNumInds,indexBase,comm,LG)";
162  const char exPfx[] =
163  "Tpetra::Map::Map(gblNumInds,indexBase,comm,LG): ";
164 
165  const bool debug = Details::Behavior::debug("Map");
166  const bool verbose = Details::Behavior::verbose("Map");
167  std::unique_ptr<std::string> prefix;
168  if (verbose) {
169  prefix = Details::createPrefix(
170  comm_.getRawPtr(), "Map", funcName);
171  std::ostringstream os;
172  os << *prefix << "Start" << endl;
173  std::cerr << os.str();
174  }
176 
177  // In debug mode only, check whether numGlobalElements and
178  // indexBase are the same over all processes in the communicator.
179  if (debug) {
180  GST proc0NumGlobalElements = numGlobalElements;
181  broadcast(*comm_, 0, outArg(proc0NumGlobalElements));
182  GST minNumGlobalElements = numGlobalElements;
183  GST maxNumGlobalElements = numGlobalElements;
184  reduceAll(*comm, REDUCE_MIN, numGlobalElements,
185  outArg(minNumGlobalElements));
186  reduceAll(*comm, REDUCE_MAX, numGlobalElements,
187  outArg(maxNumGlobalElements));
188  TEUCHOS_TEST_FOR_EXCEPTION
189  (minNumGlobalElements != maxNumGlobalElements ||
190  numGlobalElements != minNumGlobalElements,
191  std::invalid_argument, exPfx << "All processes must "
192  "provide the same number of global elements. Process 0 set "
193  "numGlobalElements="<< proc0NumGlobalElements << ". The "
194  "calling process " << comm->getRank() << " set "
195  "numGlobalElements=" << numGlobalElements << ". The min "
196  "and max values over all processes are "
197  << minNumGlobalElements << " resp. " << maxNumGlobalElements
198  << ".");
199 
200  GO proc0IndexBase = indexBase;
201  broadcast<int, GO> (*comm_, 0, outArg (proc0IndexBase));
202  GO minIndexBase = indexBase;
203  GO maxIndexBase = indexBase;
204  reduceAll(*comm, REDUCE_MIN, indexBase, outArg(minIndexBase));
205  reduceAll(*comm, REDUCE_MAX, indexBase, outArg(maxIndexBase));
206  TEUCHOS_TEST_FOR_EXCEPTION
207  (minIndexBase != maxIndexBase || indexBase != minIndexBase,
208  std::invalid_argument, exPfx << "All processes must "
209  "provide the same indexBase argument. Process 0 set "
210  "indexBase=" << proc0IndexBase << ". The calling process "
211  << comm->getRank() << " set indexBase=" << indexBase
212  << ". The min and max values over all processes are "
213  << minIndexBase << " resp. " << maxIndexBase << ".");
214  }
215 
216  // Distribute the elements across the processes in the given
217  // communicator so that global IDs (GIDs) are
218  //
219  // - Nonoverlapping (only one process owns each GID)
220  // - Contiguous (the sequence of GIDs is nondecreasing, and no two
221  // adjacent GIDs differ by more than one)
222  // - As evenly distributed as possible (the numbers of GIDs on two
223  // different processes do not differ by more than one)
224 
225  // All processes have the same numGlobalElements, but we still
226  // need to check that it is valid. numGlobalElements must be
227  // positive and not the "invalid" value (GSTI).
228  //
229  // This comparison looks funny, but it avoids compiler warnings
230  // for comparing unsigned integers (numGlobalElements_in is a
231  // GST, which is unsigned) while still working if we
232  // later decide to make GST signed.
233  TEUCHOS_TEST_FOR_EXCEPTION(
234  (numGlobalElements < 1 && numGlobalElements != 0),
235  std::invalid_argument, exPfx << "numGlobalElements (= "
236  << numGlobalElements << ") must be nonnegative.");
237 
238  TEUCHOS_TEST_FOR_EXCEPTION
239  (numGlobalElements == GSTI, std::invalid_argument, exPfx <<
240  "You provided numGlobalElements = Teuchos::OrdinalTraits<"
241  "Tpetra::global_size_t>::invalid(). This version of the "
242  "constructor requires a valid value of numGlobalElements. "
243  "You probably mistook this constructor for the \"contiguous "
244  "nonuniform\" constructor, which can compute the global "
245  "number of elements for you if you set numGlobalElements to "
246  "Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid().");
247 
248  size_t numLocalElements = 0; // will set below
249  if (lOrG == GloballyDistributed) {
250  // Compute numLocalElements:
251  //
252  // If numGlobalElements == numProcs * B + remainder,
253  // then Proc r gets B+1 elements if r < remainder,
254  // and B elements if r >= remainder.
255  //
256  // This strategy is valid for any value of numGlobalElements and
257  // numProcs, including the following border cases:
258  // - numProcs == 1
259  // - numLocalElements < numProcs
260  //
261  // In the former case, remainder == 0 && numGlobalElements ==
262  // numLocalElements. In the latter case, remainder ==
263  // numGlobalElements && numLocalElements is either 0 or 1.
264  const GST numProcs = static_cast<GST> (comm_->getSize ());
265  const GST myRank = static_cast<GST> (comm_->getRank ());
266  const GST quotient = numGlobalElements / numProcs;
267  const GST remainder = numGlobalElements - quotient * numProcs;
268 
269  GO startIndex;
270  if (myRank < remainder) {
271  numLocalElements = static_cast<size_t> (1) + static_cast<size_t> (quotient);
272  // myRank was originally an int, so it should never overflow
273  // reasonable GO types.
274  startIndex = as<GO> (myRank) * as<GO> (numLocalElements);
275  } else {
276  numLocalElements = as<size_t> (quotient);
277  startIndex = as<GO> (myRank) * as<GO> (numLocalElements) +
278  as<GO> (remainder);
279  }
280 
281  minMyGID_ = indexBase + startIndex;
282  maxMyGID_ = indexBase + startIndex + numLocalElements - 1;
283  minAllGID_ = indexBase;
284  maxAllGID_ = indexBase + numGlobalElements - 1;
285  distributed_ = (numProcs > 1);
286  }
287  else { // lOrG == LocallyReplicated
288  numLocalElements = as<size_t> (numGlobalElements);
289  minMyGID_ = indexBase;
290  maxMyGID_ = indexBase + numGlobalElements - 1;
291  distributed_ = false;
292  }
293 
294  minAllGID_ = indexBase;
295  maxAllGID_ = indexBase + numGlobalElements - 1;
296  indexBase_ = indexBase;
297  numGlobalElements_ = numGlobalElements;
298  numLocalElements_ = numLocalElements;
299  firstContiguousGID_ = minMyGID_;
300  lastContiguousGID_ = maxMyGID_;
301  contiguous_ = true;
302 
303  // Create the Directory on demand in getRemoteIndexList().
304  //setupDirectory ();
305 
306  if (verbose) {
307  std::ostringstream os;
308  os << *prefix << "Done" << endl;
309  std::cerr << os.str();
310  }
311  }
312 
313 
314  template <class LocalOrdinal, class GlobalOrdinal, class Node>
316  Map (const global_size_t numGlobalElements,
317  const size_t numLocalElements,
318  const global_ordinal_type indexBase,
319  const Teuchos::RCP<const Teuchos::Comm<int> > &comm) :
320  comm_ (comm),
321  uniform_ (false),
322  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
323  {
324  using Teuchos::as;
325  using Teuchos::broadcast;
326  using Teuchos::outArg;
327  using Teuchos::reduceAll;
328  using Teuchos::REDUCE_MIN;
329  using Teuchos::REDUCE_MAX;
330  using Teuchos::REDUCE_SUM;
331  using Teuchos::scan;
332  using std::endl;
333  using GO = global_ordinal_type;
334  using GST = global_size_t;
335  const GST GSTI = Tpetra::Details::OrdinalTraits<GST>::invalid ();
336  const char funcName[] =
337  "Map(gblNumInds,lclNumInds,indexBase,comm)";
338  const char exPfx[] =
339  "Tpetra::Map::Map(gblNumInds,lclNumInds,indexBase,comm): ";
340  const char suffix[] =
341  ". Please report this bug to the Tpetra developers.";
342 
343  const bool debug = Details::Behavior::debug("Map");
344  const bool verbose = Details::Behavior::verbose("Map");
345  std::unique_ptr<std::string> prefix;
346  if (verbose) {
347  prefix = Details::createPrefix(
348  comm_.getRawPtr(), "Map", funcName);
349  std::ostringstream os;
350  os << *prefix << "Start" << endl;
351  std::cerr << os.str();
352  }
354 
355  // Global sum of numLocalElements over all processes.
356  // Keep this for later debug checks.
357  GST debugGlobalSum {};
358  if (debug) {
359  debugGlobalSum = initialNonuniformDebugCheck(exPfx,
360  numGlobalElements, numLocalElements, indexBase, comm);
361  }
362 
363  // Distribute the elements across the nodes so that they are
364  // - non-overlapping
365  // - contiguous
366 
367  // This differs from the first Map constructor (that only takes a
368  // global number of elements) in that the user has specified the
369  // number of local elements, so that the elements are not
370  // (necessarily) evenly distributed over the processes.
371 
372  // Compute my local offset. This is an inclusive scan, so to get
373  // the final offset, we subtract off the input.
374  GO scanResult = 0;
375  scan<int, GO> (*comm, REDUCE_SUM, numLocalElements, outArg (scanResult));
376  const GO myOffset = scanResult - numLocalElements;
377 
378  if (numGlobalElements != GSTI) {
379  numGlobalElements_ = numGlobalElements; // Use the user's value.
380  }
381  else {
382  // Inclusive scan means that the last process has the final sum.
383  // Rather than doing a reduceAll to get the sum of
384  // numLocalElements, we can just have the last process broadcast
385  // its result. That saves us a round of log(numProcs) messages.
386  const int numProcs = comm->getSize ();
387  GST globalSum = scanResult;
388  if (numProcs > 1) {
389  broadcast (*comm, numProcs - 1, outArg (globalSum));
390  }
391  numGlobalElements_ = globalSum;
392 
393  if (debug) {
394  // No need for an all-reduce here; both come from collectives.
395  TEUCHOS_TEST_FOR_EXCEPTION
396  (globalSum != debugGlobalSum, std::logic_error, exPfx <<
397  "globalSum = " << globalSum << " != debugGlobalSum = " <<
398  debugGlobalSum << suffix);
399  }
400  }
401  numLocalElements_ = numLocalElements;
402  indexBase_ = indexBase;
403  minAllGID_ = (numGlobalElements_ == 0) ?
404  std::numeric_limits<GO>::max () :
405  indexBase;
406  maxAllGID_ = (numGlobalElements_ == 0) ?
407  std::numeric_limits<GO>::lowest () :
408  indexBase + GO(numGlobalElements_) - GO(1);
409  minMyGID_ = (numLocalElements_ == 0) ?
410  std::numeric_limits<GO>::max () :
411  indexBase + GO(myOffset);
412  maxMyGID_ = (numLocalElements_ == 0) ?
413  std::numeric_limits<GO>::lowest () :
414  indexBase + myOffset + GO(numLocalElements) - GO(1);
415  firstContiguousGID_ = minMyGID_;
416  lastContiguousGID_ = maxMyGID_;
417  contiguous_ = true;
418  distributed_ = checkIsDist ();
419 
420  // Create the Directory on demand in getRemoteIndexList().
421  //setupDirectory ();
422 
423  if (verbose) {
424  std::ostringstream os;
425  os << *prefix << "Done" << endl;
426  std::cerr << os.str();
427  }
428  }
429 
430  template <class LocalOrdinal, class GlobalOrdinal, class Node>
434  const char errorMessagePrefix[],
435  const global_size_t numGlobalElements,
436  const size_t numLocalElements,
437  const global_ordinal_type indexBase,
438  const Teuchos::RCP<const Teuchos::Comm<int>>& comm) const
439  {
440  const bool debug = Details::Behavior::debug("Map");
441  if (! debug) {
442  return global_size_t(0);
443  }
444 
445  using Teuchos::broadcast;
446  using Teuchos::outArg;
447  using Teuchos::ptr;
448  using Teuchos::REDUCE_MAX;
449  using Teuchos::REDUCE_MIN;
450  using Teuchos::REDUCE_SUM;
451  using Teuchos::reduceAll;
452  using GO = global_ordinal_type;
453  using GST = global_size_t;
454  const GST GSTI = Tpetra::Details::OrdinalTraits<GST>::invalid ();
455 
456  // The user has specified the distribution of indices over the
457  // processes. The distribution is not necessarily contiguous or
458  // equally shared over the processes.
459  //
460  // We assume that the number of local elements can be stored in a
461  // size_t. The instance member numLocalElements_ is a size_t, so
462  // this variable and that should have the same type.
463 
464  GST debugGlobalSum = 0; // Will be global sum of numLocalElements
465  reduceAll<int, GST> (*comm, REDUCE_SUM, static_cast<GST> (numLocalElements),
466  outArg (debugGlobalSum));
467  // In debug mode only, check whether numGlobalElements and
468  // indexBase are the same over all processes in the communicator.
469  {
470  GST proc0NumGlobalElements = numGlobalElements;
471  broadcast<int, GST> (*comm_, 0, outArg (proc0NumGlobalElements));
472  GST minNumGlobalElements = numGlobalElements;
473  GST maxNumGlobalElements = numGlobalElements;
474  reduceAll<int, GST> (*comm, REDUCE_MIN, numGlobalElements,
475  outArg (minNumGlobalElements));
476  reduceAll<int, GST> (*comm, REDUCE_MAX, numGlobalElements,
477  outArg (maxNumGlobalElements));
478  TEUCHOS_TEST_FOR_EXCEPTION
479  (minNumGlobalElements != maxNumGlobalElements ||
480  numGlobalElements != minNumGlobalElements,
481  std::invalid_argument, errorMessagePrefix << "All processes "
482  "must provide the same number of global elements, even if "
483  "that argument is "
484  "Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid() "
485  "(which signals that the Map should compute the global "
486  "number of elements). Process 0 set numGlobalElements"
487  "=" << proc0NumGlobalElements << ". The calling process "
488  << comm->getRank() << " set numGlobalElements=" <<
489  numGlobalElements << ". The min and max values over all "
490  "processes are " << minNumGlobalElements << " resp. " <<
491  maxNumGlobalElements << ".");
492 
493  GO proc0IndexBase = indexBase;
494  broadcast<int, GO> (*comm_, 0, outArg (proc0IndexBase));
495  GO minIndexBase = indexBase;
496  GO maxIndexBase = indexBase;
497  reduceAll<int, GO> (*comm, REDUCE_MIN, indexBase, outArg (minIndexBase));
498  reduceAll<int, GO> (*comm, REDUCE_MAX, indexBase, outArg (maxIndexBase));
499  TEUCHOS_TEST_FOR_EXCEPTION
500  (minIndexBase != maxIndexBase || indexBase != minIndexBase,
501  std::invalid_argument, errorMessagePrefix <<
502  "All processes must provide the same indexBase argument. "
503  "Process 0 set indexBase = " << proc0IndexBase << ". The "
504  "calling process " << comm->getRank() << " set indexBase="
505  << indexBase << ". The min and max values over all "
506  "processes are " << minIndexBase << " resp. " << maxIndexBase
507  << ".");
508 
509  // Make sure that the sum of numLocalElements over all processes
510  // equals numGlobalElements.
511  TEUCHOS_TEST_FOR_EXCEPTION
512  (numGlobalElements != GSTI &&
513  debugGlobalSum != numGlobalElements, std::invalid_argument,
514  errorMessagePrefix << "The sum of each process' number of "
515  "indices over all processes, " << debugGlobalSum << ", != "
516  << "numGlobalElements=" << numGlobalElements << ". If you "
517  "would like this constructor to compute numGlobalElements "
518  "for you, you may set numGlobalElements="
519  "Teuchos::OrdinalTraits<Tpetra::global_size_t>::invalid() "
520  "on input. Please note that this is NOT necessarily -1.");
521  }
522  return debugGlobalSum;
523  }
524 
525  template <class LocalOrdinal, class GlobalOrdinal, class Node>
526  void
527  Map<LocalOrdinal,GlobalOrdinal,Node>::
528  initWithNonownedHostIndexList(
529  const char errorMessagePrefix[],
530  const global_size_t numGlobalElements,
531  const Kokkos::View<const global_ordinal_type*,
532  Kokkos::LayoutLeft,
533  Kokkos::HostSpace,
534  Kokkos::MemoryUnmanaged>& entryList_host,
535  const global_ordinal_type indexBase,
536  const Teuchos::RCP<const Teuchos::Comm<int>>& comm)
537  {
538  using Kokkos::LayoutLeft;
539  using Kokkos::subview;
540  using Kokkos::View;
541  using Kokkos::view_alloc;
542  using Kokkos::WithoutInitializing;
543  using Teuchos::as;
544  using Teuchos::broadcast;
545  using Teuchos::outArg;
546  using Teuchos::ptr;
547  using Teuchos::REDUCE_MAX;
548  using Teuchos::REDUCE_MIN;
549  using Teuchos::REDUCE_SUM;
550  using Teuchos::reduceAll;
551  using LO = local_ordinal_type;
552  using GO = global_ordinal_type;
553  using GST = global_size_t;
554  const GST GSTI = Tpetra::Details::OrdinalTraits<GST>::invalid ();
555 
556  // Make sure that Kokkos has been initialized (Github Issue #513).
557  TEUCHOS_TEST_FOR_EXCEPTION
558  (! Kokkos::is_initialized (), std::runtime_error,
559  errorMessagePrefix << "The Kokkos execution space "
560  << Teuchos::TypeNameTraits<execution_space>::name()
561  << " has not been initialized. "
562  "Please initialize it before creating a Map.")
563 
564  // The user has specified the distribution of indices over the
565  // processes, via the input array of global indices on each
566  // process. The distribution is not necessarily contiguous or
567  // equally shared over the processes.
568 
569  // The length of the input array on this process is the number of
570  // local indices to associate with this process, even though the
571  // input array contains global indices. We assume that the number
572  // of local indices on a process can be stored in a size_t;
573  // numLocalElements_ is a size_t, so this variable and that should
574  // have the same type.
575  const size_t numLocalElements(entryList_host.size());
576 
577  initialNonuniformDebugCheck(errorMessagePrefix, numGlobalElements,
578  numLocalElements, indexBase, comm);
579 
580  // NOTE (mfh 20 Feb 2013, 10 Oct 2016) In some sense, this global
581  // reduction is redundant, since the directory Map will have to do
582  // the same thing. Thus, we could do the scan and broadcast for
583  // the directory Map here, and give the computed offsets to the
584  // directory Map's constructor. However, a reduction costs less
585  // than a scan and broadcast, so this still saves time if users of
586  // this Map don't ever need the Directory (i.e., if they never
587  // call getRemoteIndexList on this Map).
588  if (numGlobalElements != GSTI) {
589  numGlobalElements_ = numGlobalElements; // Use the user's value.
590  }
591  else { // The user wants us to compute the sum.
592  reduceAll(*comm, REDUCE_SUM,
593  static_cast<GST>(numLocalElements),
594  outArg(numGlobalElements_));
595  }
596 
597  // mfh 20 Feb 2013: We've never quite done the right thing for
598  // duplicate GIDs here. Duplicate GIDs have always been counted
599  // distinctly in numLocalElements_, and thus should get a
600  // different LID. However, we've always used std::map or a hash
601  // table for the GID -> LID lookup table, so distinct GIDs always
602  // map to the same LID. Furthermore, the order of the input GID
603  // list matters, so it's not desirable to sort for determining
604  // uniqueness.
605  //
606  // I've chosen for now to write this code as if the input GID list
607  // contains no duplicates. If this is not desired, we could use
608  // the lookup table itself to determine uniqueness: If we haven't
609  // seen the GID before, it gets a new LID and it's added to the
610  // LID -> GID and GID -> LID tables. If we have seen the GID
611  // before, it doesn't get added to either table. I would
612  // implement this, but it would cost more to do the double lookups
613  // in the table (one to check, and one to insert).
614  //
615  // More importantly, since we build the GID -> LID table in (a
616  // thread-) parallel (way), the order in which duplicate GIDs may
617  // get inserted is not defined. This would make the assignment of
618  // LID to GID nondeterministic.
619 
620  numLocalElements_ = numLocalElements;
621  indexBase_ = indexBase;
622 
623  minMyGID_ = indexBase_;
624  maxMyGID_ = indexBase_;
625 
626  // NOTE (mfh 27 May 2015): While finding the initial contiguous
627  // GID range requires looking at all the GIDs in the range,
628  // dismissing an interval of GIDs only requires looking at the
629  // first and last GIDs. Thus, we could do binary search backwards
630  // from the end in order to catch the common case of a contiguous
631  // interval followed by noncontiguous entries. On the other hand,
632  // we could just expose this case explicitly as yet another Map
633  // constructor, and avoid the trouble of detecting it.
634  if (numLocalElements_ > 0) {
635  // Find contiguous GID range, with the restriction that the
636  // beginning of the range starts with the first entry. While
637  // doing so, fill in the LID -> GID table.
638  typename decltype (lgMap_)::non_const_type lgMap
639  (view_alloc ("lgMap", WithoutInitializing), numLocalElements_);
640  auto lgMap_host =
641  Kokkos::create_mirror_view (Kokkos::HostSpace (), lgMap);
642 
643  // The input array entryList_host is already on host, so we
644  // don't need to take a host view of it.
645  // auto entryList_host =
646  // Kokkos::create_mirror_view (Kokkos::HostSpace (), entryList);
647  // Kokkos::deep_copy (entryList_host, entryList);
648 
649  firstContiguousGID_ = entryList_host[0];
650  lastContiguousGID_ = firstContiguousGID_+1;
651 
652  // FIXME (mfh 23 Sep 2015) We need to copy the input GIDs
653  // anyway, so we have to look at them all. The logical way to
654  // find the first noncontiguous entry would thus be to "reduce,"
655  // where the local reduction result is whether entryList[i] + 1
656  // == entryList[i+1].
657 
658  lgMap_host[0] = firstContiguousGID_;
659  size_t i = 1;
660  for ( ; i < numLocalElements_; ++i) {
661  const GO curGid = entryList_host[i];
662  const LO curLid = as<LO> (i);
663 
664  if (lastContiguousGID_ != curGid) break;
665 
666  // Add the entry to the LID->GID table only after we know that
667  // the current GID is in the initial contiguous sequence, so
668  // that we don't repeat adding it in the first iteration of
669  // the loop below over the remaining noncontiguous GIDs.
670  lgMap_host[curLid] = curGid;
671  ++lastContiguousGID_;
672  }
673  --lastContiguousGID_;
674 
675  // [firstContiguousGID_, lastContigousGID_] is the initial
676  // sequence of contiguous GIDs. We can start the min and max
677  // GID using this range.
678  minMyGID_ = firstContiguousGID_;
679  maxMyGID_ = lastContiguousGID_;
680 
681  // Compute the GID -> LID lookup table, _not_ including the
682  // initial sequence of contiguous GIDs.
683  {
684  const std::pair<size_t, size_t> ncRange (i, entryList_host.extent (0));
685  auto nonContigGids_host = subview (entryList_host, ncRange);
686  TEUCHOS_TEST_FOR_EXCEPTION
687  (static_cast<size_t> (nonContigGids_host.extent (0)) !=
688  static_cast<size_t> (entryList_host.extent (0) - i),
689  std::logic_error, "Tpetra::Map noncontiguous constructor: "
690  "nonContigGids_host.extent(0) = "
691  << nonContigGids_host.extent (0)
692  << " != entryList_host.extent(0) - i = "
693  << (entryList_host.extent (0) - i) << " = "
694  << entryList_host.extent (0) << " - " << i
695  << ". Please report this bug to the Tpetra developers.");
696 
697  // FixedHashTable's constructor expects an owned device View,
698  // so we must deep-copy the subview of the input indices.
699  View<GO*, LayoutLeft, device_type>
700  nonContigGids (view_alloc ("nonContigGids", WithoutInitializing),
701  nonContigGids_host.size ());
702  // DEEP_COPY REVIEW - HOST-TO-DEVICE
703  Kokkos::deep_copy (execution_space(), nonContigGids, nonContigGids_host);
704  Kokkos::fence(); // for UVM issues below - which will be refatored soon so FixedHashTable can build as pure CudaSpace - then I think remove this fence
705 
706  glMap_ = global_to_local_table_type(nonContigGids,
707  firstContiguousGID_,
708  lastContiguousGID_,
709  static_cast<LO> (i));
710  // Make host version - when memory spaces match these just do trivial assignment
711  glMapHost_ = global_to_local_table_host_type(glMap_);
712  }
713 
714  // FIXME (mfh 10 Oct 2016) When we construct the global-to-local
715  // table above, we have to look at all the (noncontiguous) input
716  // indices anyway. Thus, why not have the constructor compute
717  // and return the min and max?
718 
719  for ( ; i < numLocalElements_; ++i) {
720  const GO curGid = entryList_host[i];
721  const LO curLid = static_cast<LO> (i);
722  lgMap_host[curLid] = curGid; // LID -> GID table
723 
724  // While iterating through entryList, we compute its
725  // (process-local) min and max elements.
726  if (curGid < minMyGID_) {
727  minMyGID_ = curGid;
728  }
729  if (curGid > maxMyGID_) {
730  maxMyGID_ = curGid;
731  }
732  }
733 
734  // We filled lgMap on host above; now sync back to device.
735  // DEEP_COPY REVIEW - HOST-TO-DEVICE
736  Kokkos::deep_copy (execution_space(), lgMap, lgMap_host);
737 
738  // "Commit" the local-to-global lookup table we filled in above.
739  lgMap_ = lgMap;
740  // We've already created this, so use it.
741  lgMapHost_ = lgMap_host;
742  }
743  else {
744  minMyGID_ = std::numeric_limits<GlobalOrdinal>::max();
745  maxMyGID_ = std::numeric_limits<GlobalOrdinal>::lowest();
746  // This insures tests for GIDs in the range
747  // [firstContiguousGID_, lastContiguousGID_] fail for processes
748  // with no local elements.
749  firstContiguousGID_ = indexBase_+1;
750  lastContiguousGID_ = indexBase_;
751  // glMap_ was default constructed, so it's already empty.
752  }
753 
754  // Compute the min and max of all processes' GIDs. If
755  // numLocalElements_ == 0 on this process, minMyGID_ and maxMyGID_
756  // are both indexBase_. This is wrong, but fixing it would
757  // require either a fancy sparse all-reduce, or a custom reduction
758  // operator that ignores invalid values ("invalid" means
759  // Tpetra::Details::OrdinalTraits<GO>::invalid()).
760  //
761  // Also, while we're at it, use the same all-reduce to figure out
762  // if the Map is distributed. "Distributed" means that there is
763  // at least one process with a number of local elements less than
764  // the number of global elements.
765  //
766  // We're computing the min and max of all processes' GIDs using a
767  // single MAX all-reduce, because min(x,y) = -max(-x,-y) (when x
768  // and y are signed). (This lets us combine the min and max into
769  // a single all-reduce.) If each process sets localDist=1 if its
770  // number of local elements is strictly less than the number of
771  // global elements, and localDist=0 otherwise, then a MAX
772  // all-reduce on localDist tells us if the Map is distributed (1
773  // if yes, 0 if no). Thus, we can append localDist onto the end
774  // of the data and get the global result from the all-reduce.
775  if (std::numeric_limits<GO>::is_signed) {
776  // Does my process NOT own all the elements?
777  const GO localDist =
778  (as<GST> (numLocalElements_) < numGlobalElements_) ? 1 : 0;
779 
780  GO minMaxInput[3];
781  minMaxInput[0] = -minMyGID_;
782  minMaxInput[1] = maxMyGID_;
783  minMaxInput[2] = localDist;
784 
785  GO minMaxOutput[3];
786  minMaxOutput[0] = 0;
787  minMaxOutput[1] = 0;
788  minMaxOutput[2] = 0;
789  reduceAll<int, GO> (*comm, REDUCE_MAX, 3, minMaxInput, minMaxOutput);
790  minAllGID_ = -minMaxOutput[0];
791  maxAllGID_ = minMaxOutput[1];
792  const GO globalDist = minMaxOutput[2];
793  distributed_ = (comm_->getSize () > 1 && globalDist == 1);
794  }
795  else { // unsigned; use two reductions
796  // This is always correct, no matter the signedness of GO.
797  reduceAll<int, GO> (*comm_, REDUCE_MIN, minMyGID_, outArg (minAllGID_));
798  reduceAll<int, GO> (*comm_, REDUCE_MAX, maxMyGID_, outArg (maxAllGID_));
799  distributed_ = checkIsDist ();
800  }
801 
802  contiguous_ = false; // "Contiguous" is conservative.
803 
804  TEUCHOS_TEST_FOR_EXCEPTION(
805  minAllGID_ < indexBase_,
806  std::invalid_argument,
807  "Tpetra::Map constructor (noncontiguous): "
808  "Minimum global ID = " << minAllGID_ << " over all process(es) is "
809  "less than the given indexBase = " << indexBase_ << ".");
810 
811  // Create the Directory on demand in getRemoteIndexList().
812  //setupDirectory ();
813  }
814 
815  template <class LocalOrdinal, class GlobalOrdinal, class Node>
817  Map (const global_size_t numGlobalElements,
818  const GlobalOrdinal indexList[],
819  const LocalOrdinal indexListSize,
820  const GlobalOrdinal indexBase,
821  const Teuchos::RCP<const Teuchos::Comm<int> >& comm) :
822  comm_ (comm),
823  uniform_ (false),
824  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
825  {
826  using std::endl;
827  const char funcName[] =
828  "Map(gblNumInds,indexList,indexListSize,indexBase,comm)";
829 
830  const bool verbose = Details::Behavior::verbose("Map");
831  std::unique_ptr<std::string> prefix;
832  if (verbose) {
833  prefix = Details::createPrefix(
834  comm_.getRawPtr(), "Map", funcName);
835  std::ostringstream os;
836  os << *prefix << "Start" << endl;
837  std::cerr << os.str();
838  }
840  checkMapInputArray ("(GST, const GO[], LO, GO, comm)",
841  indexList, static_cast<size_t> (indexListSize),
842  Kokkos::DefaultHostExecutionSpace (),
843  comm.getRawPtr ());
844  // Not quite sure if I trust all code to behave correctly if the
845  // pointer is nonnull but the array length is nonzero, so I'll
846  // make sure the raw pointer is null if the length is zero.
847  const GlobalOrdinal* const indsRaw = indexListSize == 0 ? NULL : indexList;
848  Kokkos::View<const GlobalOrdinal*,
849  Kokkos::LayoutLeft,
850  Kokkos::HostSpace,
851  Kokkos::MemoryUnmanaged> inds (indsRaw, indexListSize);
852  initWithNonownedHostIndexList(funcName, numGlobalElements, inds,
853  indexBase, comm);
854  if (verbose) {
855  std::ostringstream os;
856  os << *prefix << "Done" << endl;
857  std::cerr << os.str();
858  }
859  }
860 
861  template <class LocalOrdinal, class GlobalOrdinal, class Node>
863  Map (const global_size_t numGlobalElements,
864  const Teuchos::ArrayView<const GlobalOrdinal>& entryList,
865  const GlobalOrdinal indexBase,
866  const Teuchos::RCP<const Teuchos::Comm<int> >& comm) :
867  comm_ (comm),
868  uniform_ (false),
869  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
870  {
871  using std::endl;
872  const char funcName[] =
873  "Map(gblNumInds,entryList(Teuchos::ArrayView),indexBase,comm)";
874 
875  const bool verbose = Details::Behavior::verbose("Map");
876  std::unique_ptr<std::string> prefix;
877  if (verbose) {
878  prefix = Details::createPrefix(
879  comm_.getRawPtr(), "Map", funcName);
880  std::ostringstream os;
881  os << *prefix << "Start" << endl;
882  std::cerr << os.str();
883  }
885  const size_t numLclInds = static_cast<size_t> (entryList.size ());
886  checkMapInputArray ("(GST, ArrayView, GO, comm)",
887  entryList.getRawPtr (), numLclInds,
888  Kokkos::DefaultHostExecutionSpace (),
889  comm.getRawPtr ());
890  // Not quite sure if I trust both ArrayView and View to behave
891  // correctly if the pointer is nonnull but the array length is
892  // nonzero, so I'll make sure it's null if the length is zero.
893  const GlobalOrdinal* const indsRaw =
894  numLclInds == 0 ? NULL : entryList.getRawPtr ();
895  Kokkos::View<const GlobalOrdinal*,
896  Kokkos::LayoutLeft,
897  Kokkos::HostSpace,
898  Kokkos::MemoryUnmanaged> inds (indsRaw, numLclInds);
899  initWithNonownedHostIndexList(funcName, numGlobalElements, inds,
900  indexBase, comm);
901  if (verbose) {
902  std::ostringstream os;
903  os << *prefix << "Done" << endl;
904  std::cerr << os.str();
905  }
906  }
907 
908  template <class LocalOrdinal, class GlobalOrdinal, class Node>
910  Map (const global_size_t numGlobalElements,
911  const Kokkos::View<const GlobalOrdinal*, device_type>& entryList,
912  const GlobalOrdinal indexBase,
913  const Teuchos::RCP<const Teuchos::Comm<int> >& comm) :
914  comm_ (comm),
915  uniform_ (false),
916  directory_ (new Directory<LocalOrdinal, GlobalOrdinal, Node> ())
917  {
918  using Kokkos::LayoutLeft;
919  using Kokkos::subview;
920  using Kokkos::View;
921  using Kokkos::view_alloc;
922  using Kokkos::WithoutInitializing;
923  using Teuchos::arcp;
924  using Teuchos::ArrayView;
925  using Teuchos::as;
926  using Teuchos::broadcast;
927  using Teuchos::outArg;
928  using Teuchos::ptr;
929  using Teuchos::REDUCE_MAX;
930  using Teuchos::REDUCE_MIN;
931  using Teuchos::REDUCE_SUM;
932  using Teuchos::reduceAll;
933  using Teuchos::typeName;
934  using std::endl;
935  using LO = local_ordinal_type;
936  using GO = global_ordinal_type;
937  using GST = global_size_t;
938  const GST GSTI = Tpetra::Details::OrdinalTraits<GST>::invalid ();
939  const char funcName[] =
940  "Map(gblNumInds,entryList(Kokkos::View),indexBase,comm)";
941 
942  const bool verbose = Details::Behavior::verbose("Map");
943  std::unique_ptr<std::string> prefix;
944  if (verbose) {
945  prefix = Details::createPrefix(
946  comm_.getRawPtr(), "Map", funcName);
947  std::ostringstream os;
948  os << *prefix << "Start" << endl;
949  std::cerr << os.str();
950  }
952  checkMapInputArray ("(GST, Kokkos::View, GO, comm)",
953  entryList.data (),
954  static_cast<size_t> (entryList.extent (0)),
955  execution_space (), comm.getRawPtr ());
956 
957  // The user has specified the distribution of indices over the
958  // processes, via the input array of global indices on each
959  // process. The distribution is not necessarily contiguous or
960  // equally shared over the processes.
961 
962  // The length of the input array on this process is the number of
963  // local indices to associate with this process, even though the
964  // input array contains global indices. We assume that the number
965  // of local indices on a process can be stored in a size_t;
966  // numLocalElements_ is a size_t, so this variable and that should
967  // have the same type.
968  const size_t numLocalElements(entryList.size());
969 
970  initialNonuniformDebugCheck(funcName, numGlobalElements,
971  numLocalElements, indexBase, comm);
972 
973  // NOTE (mfh 20 Feb 2013, 10 Oct 2016) In some sense, this global
974  // reduction is redundant, since the directory Map will have to do
975  // the same thing. Thus, we could do the scan and broadcast for
976  // the directory Map here, and give the computed offsets to the
977  // directory Map's constructor. However, a reduction costs less
978  // than a scan and broadcast, so this still saves time if users of
979  // this Map don't ever need the Directory (i.e., if they never
980  // call getRemoteIndexList on this Map).
981  if (numGlobalElements != GSTI) {
982  numGlobalElements_ = numGlobalElements; // Use the user's value.
983  }
984  else { // The user wants us to compute the sum.
985  reduceAll(*comm, REDUCE_SUM,
986  static_cast<GST>(numLocalElements),
987  outArg(numGlobalElements_));
988  }
989 
990  // mfh 20 Feb 2013: We've never quite done the right thing for
991  // duplicate GIDs here. Duplicate GIDs have always been counted
992  // distinctly in numLocalElements_, and thus should get a
993  // different LID. However, we've always used std::map or a hash
994  // table for the GID -> LID lookup table, so distinct GIDs always
995  // map to the same LID. Furthermore, the order of the input GID
996  // list matters, so it's not desirable to sort for determining
997  // uniqueness.
998  //
999  // I've chosen for now to write this code as if the input GID list
1000  // contains no duplicates. If this is not desired, we could use
1001  // the lookup table itself to determine uniqueness: If we haven't
1002  // seen the GID before, it gets a new LID and it's added to the
1003  // LID -> GID and GID -> LID tables. If we have seen the GID
1004  // before, it doesn't get added to either table. I would
1005  // implement this, but it would cost more to do the double lookups
1006  // in the table (one to check, and one to insert).
1007  //
1008  // More importantly, since we build the GID -> LID table in (a
1009  // thread-) parallel (way), the order in which duplicate GIDs may
1010  // get inserted is not defined. This would make the assignment of
1011  // LID to GID nondeterministic.
1012 
1013  numLocalElements_ = numLocalElements;
1014  indexBase_ = indexBase;
1015 
1016  minMyGID_ = indexBase_;
1017  maxMyGID_ = indexBase_;
1018 
1019  // NOTE (mfh 27 May 2015): While finding the initial contiguous
1020  // GID range requires looking at all the GIDs in the range,
1021  // dismissing an interval of GIDs only requires looking at the
1022  // first and last GIDs. Thus, we could do binary search backwards
1023  // from the end in order to catch the common case of a contiguous
1024  // interval followed by noncontiguous entries. On the other hand,
1025  // we could just expose this case explicitly as yet another Map
1026  // constructor, and avoid the trouble of detecting it.
1027  if (numLocalElements_ > 0) {
1028  // Find contiguous GID range, with the restriction that the
1029  // beginning of the range starts with the first entry. While
1030  // doing so, fill in the LID -> GID table.
1031  typename decltype (lgMap_)::non_const_type lgMap
1032  (view_alloc ("lgMap", WithoutInitializing), numLocalElements_);
1033  auto lgMap_host =
1034  Kokkos::create_mirror_view (Kokkos::HostSpace (), lgMap);
1035 
1036  using array_layout =
1037  typename View<const GO*, device_type>::array_layout;
1038  View<GO*, array_layout, Kokkos::HostSpace> entryList_host
1039  (view_alloc ("entryList_host", WithoutInitializing),
1040  entryList.extent(0));
1041  // DEEP_COPY REVIEW - DEVICE-TO-HOST
1042  Kokkos::deep_copy (execution_space(), entryList_host, entryList);
1043  Kokkos::fence(); // UVM follows
1044  firstContiguousGID_ = entryList_host[0];
1045  lastContiguousGID_ = firstContiguousGID_+1;
1046 
1047  // FIXME (mfh 23 Sep 2015) We need to copy the input GIDs
1048  // anyway, so we have to look at them all. The logical way to
1049  // find the first noncontiguous entry would thus be to "reduce,"
1050  // where the local reduction result is whether entryList[i] + 1
1051  // == entryList[i+1].
1052 
1053  lgMap_host[0] = firstContiguousGID_;
1054  size_t i = 1;
1055  for ( ; i < numLocalElements_; ++i) {
1056  const GO curGid = entryList_host[i];
1057  const LO curLid = as<LO> (i);
1058 
1059  if (lastContiguousGID_ != curGid) break;
1060 
1061  // Add the entry to the LID->GID table only after we know that
1062  // the current GID is in the initial contiguous sequence, so
1063  // that we don't repeat adding it in the first iteration of
1064  // the loop below over the remaining noncontiguous GIDs.
1065  lgMap_host[curLid] = curGid;
1066  ++lastContiguousGID_;
1067  }
1068  --lastContiguousGID_;
1069 
1070  // [firstContiguousGID_, lastContigousGID_] is the initial
1071  // sequence of contiguous GIDs. We can start the min and max
1072  // GID using this range.
1073  minMyGID_ = firstContiguousGID_;
1074  maxMyGID_ = lastContiguousGID_;
1075 
1076  // Compute the GID -> LID lookup table, _not_ including the
1077  // initial sequence of contiguous GIDs.
1078  {
1079  const std::pair<size_t, size_t> ncRange (i, entryList.extent (0));
1080  auto nonContigGids = subview (entryList, ncRange);
1081  TEUCHOS_TEST_FOR_EXCEPTION
1082  (static_cast<size_t> (nonContigGids.extent (0)) !=
1083  static_cast<size_t> (entryList.extent (0) - i),
1084  std::logic_error, "Tpetra::Map noncontiguous constructor: "
1085  "nonContigGids.extent(0) = "
1086  << nonContigGids.extent (0)
1087  << " != entryList.extent(0) - i = "
1088  << (entryList.extent (0) - i) << " = "
1089  << entryList.extent (0) << " - " << i
1090  << ". Please report this bug to the Tpetra developers.");
1091 
1092  glMap_ = global_to_local_table_type(nonContigGids,
1093  firstContiguousGID_,
1094  lastContiguousGID_,
1095  static_cast<LO> (i));
1096  // Make host version - when memory spaces match these just do trivial assignment
1097  glMapHost_ = global_to_local_table_host_type(glMap_);
1098  }
1099 
1100  // FIXME (mfh 10 Oct 2016) When we construct the global-to-local
1101  // table above, we have to look at all the (noncontiguous) input
1102  // indices anyway. Thus, why not have the constructor compute
1103  // and return the min and max?
1104 
1105  for ( ; i < numLocalElements_; ++i) {
1106  const GO curGid = entryList_host[i];
1107  const LO curLid = static_cast<LO> (i);
1108  lgMap_host[curLid] = curGid; // LID -> GID table
1109 
1110  // While iterating through entryList, we compute its
1111  // (process-local) min and max elements.
1112  if (curGid < minMyGID_) {
1113  minMyGID_ = curGid;
1114  }
1115  if (curGid > maxMyGID_) {
1116  maxMyGID_ = curGid;
1117  }
1118  }
1119 
1120  // We filled lgMap on host above; now sync back to device.
1121  // DEEP_COPY REVIEW - HOST-TO-DEVICE
1122  Kokkos::deep_copy (execution_space(), lgMap, lgMap_host);
1123 
1124  // "Commit" the local-to-global lookup table we filled in above.
1125  lgMap_ = lgMap;
1126  // We've already created this, so use it.
1127  lgMapHost_ = lgMap_host;
1128  }
1129  else {
1130  minMyGID_ = std::numeric_limits<GlobalOrdinal>::max();
1131  maxMyGID_ = std::numeric_limits<GlobalOrdinal>::lowest();
1132  // This insures tests for GIDs in the range
1133  // [firstContiguousGID_, lastContiguousGID_] fail for processes
1134  // with no local elements.
1135  firstContiguousGID_ = indexBase_+1;
1136  lastContiguousGID_ = indexBase_;
1137  // glMap_ was default constructed, so it's already empty.
1138  }
1139 
1140  // Compute the min and max of all processes' GIDs. If
1141  // numLocalElements_ == 0 on this process, minMyGID_ and maxMyGID_
1142  // are both indexBase_. This is wrong, but fixing it would
1143  // require either a fancy sparse all-reduce, or a custom reduction
1144  // operator that ignores invalid values ("invalid" means
1145  // Tpetra::Details::OrdinalTraits<GO>::invalid()).
1146  //
1147  // Also, while we're at it, use the same all-reduce to figure out
1148  // if the Map is distributed. "Distributed" means that there is
1149  // at least one process with a number of local elements less than
1150  // the number of global elements.
1151  //
1152  // We're computing the min and max of all processes' GIDs using a
1153  // single MAX all-reduce, because min(x,y) = -max(-x,-y) (when x
1154  // and y are signed). (This lets us combine the min and max into
1155  // a single all-reduce.) If each process sets localDist=1 if its
1156  // number of local elements is strictly less than the number of
1157  // global elements, and localDist=0 otherwise, then a MAX
1158  // all-reduce on localDist tells us if the Map is distributed (1
1159  // if yes, 0 if no). Thus, we can append localDist onto the end
1160  // of the data and get the global result from the all-reduce.
1161  if (std::numeric_limits<GO>::is_signed) {
1162  // Does my process NOT own all the elements?
1163  const GO localDist =
1164  (as<GST> (numLocalElements_) < numGlobalElements_) ? 1 : 0;
1165 
1166  GO minMaxInput[3];
1167  minMaxInput[0] = -minMyGID_;
1168  minMaxInput[1] = maxMyGID_;
1169  minMaxInput[2] = localDist;
1170 
1171  GO minMaxOutput[3];
1172  minMaxOutput[0] = 0;
1173  minMaxOutput[1] = 0;
1174  minMaxOutput[2] = 0;
1175  reduceAll<int, GO> (*comm, REDUCE_MAX, 3, minMaxInput, minMaxOutput);
1176  minAllGID_ = -minMaxOutput[0];
1177  maxAllGID_ = minMaxOutput[1];
1178  const GO globalDist = minMaxOutput[2];
1179  distributed_ = (comm_->getSize () > 1 && globalDist == 1);
1180  }
1181  else { // unsigned; use two reductions
1182  // This is always correct, no matter the signedness of GO.
1183  reduceAll<int, GO> (*comm_, REDUCE_MIN, minMyGID_, outArg (minAllGID_));
1184  reduceAll<int, GO> (*comm_, REDUCE_MAX, maxMyGID_, outArg (maxAllGID_));
1185  distributed_ = checkIsDist ();
1186  }
1187 
1188  contiguous_ = false; // "Contiguous" is conservative.
1189 
1190  TEUCHOS_TEST_FOR_EXCEPTION(
1191  minAllGID_ < indexBase_,
1192  std::invalid_argument,
1193  "Tpetra::Map constructor (noncontiguous): "
1194  "Minimum global ID = " << minAllGID_ << " over all process(es) is "
1195  "less than the given indexBase = " << indexBase_ << ".");
1196 
1197  // Create the Directory on demand in getRemoteIndexList().
1198  //setupDirectory ();
1199 
1200  if (verbose) {
1201  std::ostringstream os;
1202  os << *prefix << "Done" << endl;
1203  std::cerr << os.str();
1204  }
1205  }
1206 
1207 
1208  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1210  {
1211  if (! Kokkos::is_initialized ()) {
1212  std::ostringstream os;
1213  os << "WARNING: Tpetra::Map destructor (~Map()) is being called after "
1214  "Kokkos::finalize() has been called. This is user error! There are "
1215  "two likely causes: " << std::endl <<
1216  " 1. You have a static Tpetra::Map (or RCP or shared_ptr of a Map)"
1217  << std::endl <<
1218  " 2. You declare and construct a Tpetra::Map (or RCP or shared_ptr "
1219  "of a Tpetra::Map) at the same scope in main() as Kokkos::finalize() "
1220  "or Tpetra::finalize()." << std::endl << std::endl <<
1221  "Don't do either of these! Please refer to GitHib Issue #2372."
1222  << std::endl;
1223  ::Tpetra::Details::printOnce (std::cerr, os.str (),
1224  this->getComm ().getRawPtr ());
1225  }
1226  else {
1227  using ::Tpetra::Details::mpiIsInitialized;
1228  using ::Tpetra::Details::mpiIsFinalized;
1229  using ::Tpetra::Details::teuchosCommIsAnMpiComm;
1230 
1231  Teuchos::RCP<const Teuchos::Comm<int> > comm = this->getComm ();
1232  if (! comm.is_null () && teuchosCommIsAnMpiComm (*comm) &&
1233  mpiIsInitialized () && mpiIsFinalized ()) {
1234  // Tpetra itself does not require MPI, even if building with
1235  // MPI. It is legal to create Tpetra objects that do not use
1236  // MPI, even in an MPI program. However, calling Tpetra stuff
1237  // after MPI_Finalize() has been called is a bad idea, since
1238  // some Tpetra defaults may use MPI if available.
1239  std::ostringstream os;
1240  os << "WARNING: Tpetra::Map destructor (~Map()) is being called after "
1241  "MPI_Finalize() has been called. This is user error! There are "
1242  "two likely causes: " << std::endl <<
1243  " 1. You have a static Tpetra::Map (or RCP or shared_ptr of a Map)"
1244  << std::endl <<
1245  " 2. You declare and construct a Tpetra::Map (or RCP or shared_ptr "
1246  "of a Tpetra::Map) at the same scope in main() as MPI_finalize() or "
1247  "Tpetra::finalize()." << std::endl << std::endl <<
1248  "Don't do either of these! Please refer to GitHib Issue #2372."
1249  << std::endl;
1250  ::Tpetra::Details::printOnce (std::cerr, os.str (), comm.getRawPtr ());
1251  }
1252  }
1253  // mfh 20 Mar 2018: We can't check Tpetra::isInitialized() yet,
1254  // because Tpetra does not yet require Tpetra::initialize /
1255  // Tpetra::finalize.
1256  }
1257 
1258  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1259  bool
1261  {
1262  TEUCHOS_TEST_FOR_EXCEPTION(
1263  getComm ().is_null (), std::logic_error, "Tpetra::Map::isOneToOne: "
1264  "getComm() returns null. Please report this bug to the Tpetra "
1265  "developers.");
1266 
1267  // This is a collective operation, if it hasn't been called before.
1268  setupDirectory ();
1269  return directory_->isOneToOne (*this);
1270  }
1271 
1272  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1273  LocalOrdinal
1275  getLocalElement (GlobalOrdinal globalIndex) const
1276  {
1277  if (isContiguous ()) {
1278  if (globalIndex < getMinGlobalIndex () ||
1279  globalIndex > getMaxGlobalIndex ()) {
1280  return Tpetra::Details::OrdinalTraits<LocalOrdinal>::invalid ();
1281  }
1282  return static_cast<LocalOrdinal> (globalIndex - getMinGlobalIndex ());
1283  }
1284  else if (globalIndex >= firstContiguousGID_ &&
1285  globalIndex <= lastContiguousGID_) {
1286  return static_cast<LocalOrdinal> (globalIndex - firstContiguousGID_);
1287  }
1288  else {
1289  // If the given global index is not in the table, this returns
1290  // the same value as OrdinalTraits<LocalOrdinal>::invalid().
1291  // glMapHost_ is Host and does not assume UVM
1292  return glMapHost_.get (globalIndex);
1293  }
1294  }
1295 
1296  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1297  GlobalOrdinal
1299  getGlobalElement (LocalOrdinal localIndex) const
1300  {
1301  if (localIndex < getMinLocalIndex () || localIndex > getMaxLocalIndex ()) {
1302  return Tpetra::Details::OrdinalTraits<GlobalOrdinal>::invalid ();
1303  }
1304  if (isContiguous ()) {
1305  return getMinGlobalIndex () + localIndex;
1306  }
1307  else {
1308  // This is a host Kokkos::View access, with no RCP or ArrayRCP
1309  // involvement. As a result, it is thread safe.
1310  //
1311  // lgMapHost_ is a host pointer; this does NOT assume UVM.
1312  return lgMapHost_[localIndex];
1313  }
1314  }
1315 
1316  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1317  bool
1319  isNodeLocalElement (LocalOrdinal localIndex) const
1320  {
1321  if (localIndex < getMinLocalIndex () || localIndex > getMaxLocalIndex ()) {
1322  return false;
1323  } else {
1324  return true;
1325  }
1326  }
1327 
1328  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1329  bool
1331  isNodeGlobalElement (GlobalOrdinal globalIndex) const {
1332  return this->getLocalElement (globalIndex) !=
1333  Tpetra::Details::OrdinalTraits<LocalOrdinal>::invalid ();
1334  }
1335 
1336  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1338  return uniform_;
1339  }
1340 
1341  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1343  return contiguous_;
1344  }
1345 
1346 
1347  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1351  {
1352  return local_map_type (glMap_, lgMap_, getIndexBase (),
1353  getMinGlobalIndex (), getMaxGlobalIndex (),
1354  firstContiguousGID_, lastContiguousGID_,
1355  getLocalNumElements (), isContiguous ());
1356  }
1357 
1358  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1359  bool
1362  {
1363  using Teuchos::outArg;
1364  using Teuchos::REDUCE_MIN;
1365  using Teuchos::reduceAll;
1366  //
1367  // Tests that avoid the Boolean all-reduce below by using
1368  // globally consistent quantities.
1369  //
1370  if (this == &map) {
1371  // Pointer equality on one process always implies pointer
1372  // equality on all processes, since Map is immutable.
1373  return true;
1374  }
1375  else if (getComm ()->getSize () != map.getComm ()->getSize ()) {
1376  // The two communicators have different numbers of processes.
1377  // It's not correct to call isCompatible() in that case. This
1378  // may result in the all-reduce hanging below.
1379  return false;
1380  }
1381  else if (getGlobalNumElements () != map.getGlobalNumElements ()) {
1382  // Two Maps are definitely NOT compatible if they have different
1383  // global numbers of indices.
1384  return false;
1385  }
1386  else if (isContiguous () && isUniform () &&
1387  map.isContiguous () && map.isUniform ()) {
1388  // Contiguous uniform Maps with the same number of processes in
1389  // their communicators, and with the same global numbers of
1390  // indices, are always compatible.
1391  return true;
1392  }
1393  else if (! isContiguous () && ! map.isContiguous () &&
1394  lgMap_.extent (0) != 0 && map.lgMap_.extent (0) != 0 &&
1395  lgMap_.data () == map.lgMap_.data ()) {
1396  // Noncontiguous Maps whose global index lists are nonempty and
1397  // have the same pointer must be the same (and therefore
1398  // contiguous).
1399  //
1400  // Nonempty is important. For example, consider a communicator
1401  // with two processes, and two Maps that share this
1402  // communicator, with zero global indices on the first process,
1403  // and different nonzero numbers of global indices on the second
1404  // process. In that case, on the first process, the pointers
1405  // would both be NULL.
1406  return true;
1407  }
1408 
1409  TEUCHOS_TEST_FOR_EXCEPTION(
1410  getGlobalNumElements () != map.getGlobalNumElements (), std::logic_error,
1411  "Tpetra::Map::isCompatible: There's a bug in this method. We've already "
1412  "checked that this condition is true above, but it's false here. "
1413  "Please report this bug to the Tpetra developers.");
1414 
1415  // Do both Maps have the same number of indices on each process?
1416  const int locallyCompat =
1417  (getLocalNumElements () == map.getLocalNumElements ()) ? 1 : 0;
1418 
1419  int globallyCompat = 0;
1420  reduceAll<int, int> (*comm_, REDUCE_MIN, locallyCompat, outArg (globallyCompat));
1421  return (globallyCompat == 1);
1422  }
1423 
1424  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1425  bool
1428  {
1429  using Teuchos::ArrayView;
1430  using GO = global_ordinal_type;
1431  using size_type = typename ArrayView<const GO>::size_type;
1432 
1433  // If both Maps are contiguous, we can compare their GID ranges
1434  // easily by looking at the min and max GID on this process.
1435  // Otherwise, we'll compare their GID lists. If only one Map is
1436  // contiguous, then we only have to call getLocalElementList() on
1437  // the noncontiguous Map. (It's best to avoid calling it on a
1438  // contiguous Map, since it results in unnecessary storage that
1439  // persists for the lifetime of the Map.)
1440 
1441  if (this == &map) {
1442  // Pointer equality on one process always implies pointer
1443  // equality on all processes, since Map is immutable.
1444  return true;
1445  }
1446  else if (getLocalNumElements () != map.getLocalNumElements ()) {
1447  return false;
1448  }
1449  else if (getMinGlobalIndex () != map.getMinGlobalIndex () ||
1450  getMaxGlobalIndex () != map.getMaxGlobalIndex ()) {
1451  return false;
1452  }
1453  else {
1454  if (isContiguous ()) {
1455  if (map.isContiguous ()) {
1456  return true; // min and max match, so the ranges match.
1457  }
1458  else { // *this is contiguous, but map is not contiguous
1459  TEUCHOS_TEST_FOR_EXCEPTION(
1460  ! this->isContiguous () || map.isContiguous (), std::logic_error,
1461  "Tpetra::Map::locallySameAs: BUG");
1462  ArrayView<const GO> rhsElts = map.getLocalElementList ();
1463  const GO minLhsGid = this->getMinGlobalIndex ();
1464  const size_type numRhsElts = rhsElts.size ();
1465  for (size_type k = 0; k < numRhsElts; ++k) {
1466  const GO curLhsGid = minLhsGid + static_cast<GO> (k);
1467  if (curLhsGid != rhsElts[k]) {
1468  return false; // stop on first mismatch
1469  }
1470  }
1471  return true;
1472  }
1473  }
1474  else if (map.isContiguous ()) { // *this is not contiguous, but map is
1475  TEUCHOS_TEST_FOR_EXCEPTION(
1476  this->isContiguous () || ! map.isContiguous (), std::logic_error,
1477  "Tpetra::Map::locallySameAs: BUG");
1478  ArrayView<const GO> lhsElts = this->getLocalElementList ();
1479  const GO minRhsGid = map.getMinGlobalIndex ();
1480  const size_type numLhsElts = lhsElts.size ();
1481  for (size_type k = 0; k < numLhsElts; ++k) {
1482  const GO curRhsGid = minRhsGid + static_cast<GO> (k);
1483  if (curRhsGid != lhsElts[k]) {
1484  return false; // stop on first mismatch
1485  }
1486  }
1487  return true;
1488  }
1489  else if (this->lgMap_.data () == map.lgMap_.data ()) {
1490  // Pointers to LID->GID "map" (actually just an array) are the
1491  // same, and the number of GIDs are the same.
1492  return this->getLocalNumElements () == map.getLocalNumElements ();
1493  }
1494  else { // we actually have to compare the GIDs
1495  if (this->getLocalNumElements () != map.getLocalNumElements ()) {
1496  return false; // We already checked above, but check just in case
1497  }
1498  else {
1499  ArrayView<const GO> lhsElts = getLocalElementList ();
1500  ArrayView<const GO> rhsElts = map.getLocalElementList ();
1501 
1502  // std::equal requires that the latter range is as large as
1503  // the former. We know the ranges have equal length, because
1504  // they have the same number of local entries.
1505  return std::equal (lhsElts.begin (), lhsElts.end (), rhsElts.begin ());
1506  }
1507  }
1508  }
1509  }
1510 
1511  template <class LocalOrdinal,class GlobalOrdinal, class Node>
1512  bool
1515  {
1516  if (this == &map)
1517  return true;
1518 
1519  // We are going to check if lmap1 is fitted into lmap2:
1520  // Is lmap1 (map) a subset of lmap2 (this)?
1521  // And do the first lmap1.getLocalNumElements() global elements
1522  // of lmap1,lmap2 owned on each process exactly match?
1523  auto lmap1 = map.getLocalMap();
1524  auto lmap2 = this->getLocalMap();
1525 
1526  auto numLocalElements1 = lmap1.getLocalNumElements();
1527  auto numLocalElements2 = lmap2.getLocalNumElements();
1528 
1529  if (numLocalElements1 > numLocalElements2) {
1530  // There are more indices in the first map on this process than in second map.
1531  return false;
1532  }
1533 
1534  if (lmap1.isContiguous () && lmap2.isContiguous ()) {
1535  // When both Maps are contiguous, just check the interval inclusion.
1536  return ((lmap1.getMinGlobalIndex () == lmap2.getMinGlobalIndex ()) &&
1537  (lmap1.getMaxGlobalIndex () <= lmap2.getMaxGlobalIndex ()));
1538  }
1539 
1540  if (lmap1.getMinGlobalIndex () < lmap2.getMinGlobalIndex () ||
1541  lmap1.getMaxGlobalIndex () > lmap2.getMaxGlobalIndex ()) {
1542  // The second map does not include the first map bounds, and thus some of
1543  // the first map global indices are not in the second map.
1544  return false;
1545  }
1546 
1547  using LO = local_ordinal_type;
1548  using range_type =
1549  Kokkos::RangePolicy<LO, typename node_type::execution_space>;
1550 
1551  // Check all elements.
1552  LO numDiff = 0;
1553  Kokkos::parallel_reduce(
1554  "isLocallyFitted",
1555  range_type(0, numLocalElements1),
1556  KOKKOS_LAMBDA (const LO i, LO& diff) {
1557  diff += (lmap1.getGlobalElement(i) != lmap2.getGlobalElement(i));
1558  }, numDiff);
1559 
1560  return (numDiff == 0);
1561  }
1562 
1563  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1564  bool
1567  {
1568  using Teuchos::outArg;
1569  using Teuchos::REDUCE_MIN;
1570  using Teuchos::reduceAll;
1571  //
1572  // Tests that avoid the Boolean all-reduce below by using
1573  // globally consistent quantities.
1574  //
1575  if (this == &map) {
1576  // Pointer equality on one process always implies pointer
1577  // equality on all processes, since Map is immutable.
1578  return true;
1579  }
1580  else if (getComm ()->getSize () != map.getComm ()->getSize ()) {
1581  // The two communicators have different numbers of processes.
1582  // It's not correct to call isSameAs() in that case. This
1583  // may result in the all-reduce hanging below.
1584  return false;
1585  }
1586  else if (getGlobalNumElements () != map.getGlobalNumElements ()) {
1587  // Two Maps are definitely NOT the same if they have different
1588  // global numbers of indices.
1589  return false;
1590  }
1591  else if (getMinAllGlobalIndex () != map.getMinAllGlobalIndex () ||
1592  getMaxAllGlobalIndex () != map.getMaxAllGlobalIndex () ||
1593  getIndexBase () != map.getIndexBase ()) {
1594  // If the global min or max global index doesn't match, or if
1595  // the index base doesn't match, then the Maps aren't the same.
1596  return false;
1597  }
1598  else if (isDistributed () != map.isDistributed ()) {
1599  // One Map is distributed and the other is not, which means that
1600  // the Maps aren't the same.
1601  return false;
1602  }
1603  else if (isContiguous () && isUniform () &&
1604  map.isContiguous () && map.isUniform ()) {
1605  // Contiguous uniform Maps with the same number of processes in
1606  // their communicators, with the same global numbers of indices,
1607  // and with matching index bases and ranges, must be the same.
1608  return true;
1609  }
1610 
1611  // The two communicators must have the same number of processes,
1612  // with process ranks occurring in the same order. This uses
1613  // MPI_COMM_COMPARE. The MPI 3.1 standard (Section 6.4) says:
1614  // "Operations that access communicators are local and their
1615  // execution does not require interprocess communication."
1616  // However, just to be sure, I'll put this call after the above
1617  // tests that don't communicate.
1618  if (! ::Tpetra::Details::congruent (*comm_, * (map.getComm ()))) {
1619  return false;
1620  }
1621 
1622  // If we get this far, we need to check local properties and then
1623  // communicate local sameness across all processes.
1624  const int isSame_lcl = locallySameAs (map) ? 1 : 0;
1625 
1626  // Return true if and only if all processes report local sameness.
1627  int isSame_gbl = 0;
1628  reduceAll<int, int> (*comm_, REDUCE_MIN, isSame_lcl, outArg (isSame_gbl));
1629  return isSame_gbl == 1;
1630  }
1631 
1632  namespace { // (anonymous)
1633  template <class LO, class GO, class DT>
1634  class FillLgMap {
1635  public:
1636  FillLgMap (const Kokkos::View<GO*, DT>& lgMap,
1637  const GO startGid) :
1638  lgMap_ (lgMap), startGid_ (startGid)
1639  {
1640  Kokkos::RangePolicy<LO, typename DT::execution_space>
1641  range (static_cast<LO> (0), static_cast<LO> (lgMap.size ()));
1642  Kokkos::parallel_for (range, *this);
1643  }
1644 
1645  KOKKOS_INLINE_FUNCTION void operator () (const LO& lid) const {
1646  lgMap_(lid) = startGid_ + static_cast<GO> (lid);
1647  }
1648 
1649  private:
1650  const Kokkos::View<GO*, DT> lgMap_;
1651  const GO startGid_;
1652  };
1653 
1654  } // namespace (anonymous)
1655 
1656 
1657  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1658  typename Map<LocalOrdinal,GlobalOrdinal,Node>::global_indices_array_type
1660  {
1661  using std::endl;
1662  using LO = local_ordinal_type;
1663  using GO = global_ordinal_type;
1664  using const_lg_view_type = decltype(lgMap_);
1665  using lg_view_type = typename const_lg_view_type::non_const_type;
1666  const bool debug = Details::Behavior::debug("Map");
1667  const bool verbose = Details::Behavior::verbose("Map");
1668 
1669  std::unique_ptr<std::string> prefix;
1670  if (verbose) {
1671  prefix = Details::createPrefix(
1672  comm_.getRawPtr(), "Map", "getMyGlobalIndices");
1673  std::ostringstream os;
1674  os << *prefix << "Start" << endl;
1675  std::cerr << os.str();
1676  }
1677 
1678  // If the local-to-global mapping doesn't exist yet, and if we
1679  // have local entries, then create and fill the local-to-global
1680  // mapping.
1681  const bool needToCreateLocalToGlobalMapping =
1682  lgMap_.extent (0) == 0 && numLocalElements_ > 0;
1683 
1684  if (needToCreateLocalToGlobalMapping) {
1685  if (verbose) {
1686  std::ostringstream os;
1687  os << *prefix << "Need to create lgMap" << endl;
1688  std::cerr << os.str();
1689  }
1690  if (debug) {
1691  // The local-to-global mapping should have been set up already
1692  // for a noncontiguous map.
1693  TEUCHOS_TEST_FOR_EXCEPTION
1694  (! isContiguous(), std::logic_error,
1695  "Tpetra::Map::getMyGlobalIndices: The local-to-global "
1696  "mapping (lgMap_) should have been set up already for a "
1697  "noncontiguous Map. Please report this bug to the Tpetra "
1698  "developers.");
1699  }
1700  const LO numElts = static_cast<LO> (getLocalNumElements ());
1701 
1702  using Kokkos::view_alloc;
1703  using Kokkos::WithoutInitializing;
1704  lg_view_type lgMap ("lgMap", numElts);
1705  if (verbose) {
1706  std::ostringstream os;
1707  os << *prefix << "Fill lgMap" << endl;
1708  std::cerr << os.str();
1709  }
1710  FillLgMap<LO, GO, no_uvm_device_type> fillIt (lgMap, minMyGID_);
1711 
1712  if (verbose) {
1713  std::ostringstream os;
1714  os << *prefix << "Copy lgMap to lgMapHost" << endl;
1715  std::cerr << os.str();
1716  }
1717 
1718  auto lgMapHost =
1719  Kokkos::create_mirror_view (Kokkos::HostSpace (), lgMap);
1720  // DEEP_COPY REVIEW - DEVICE-TO-HOST
1721  Kokkos::deep_copy (execution_space(), lgMapHost, lgMap);
1722 
1723  // "Commit" the local-to-global lookup table we filled in above.
1724  lgMap_ = lgMap;
1725  lgMapHost_ = lgMapHost;
1726  }
1727 
1728  if (verbose) {
1729  std::ostringstream os;
1730  os << *prefix << "Done" << endl;
1731  std::cerr << os.str();
1732  }
1733  return lgMapHost_;
1734  }
1735 
1736  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1737  Teuchos::ArrayView<const GlobalOrdinal>
1739  {
1740  using GO = global_ordinal_type;
1741 
1742  // If the local-to-global mapping doesn't exist yet, and if we
1743  // have local entries, then create and fill the local-to-global
1744  // mapping.
1745  (void) this->getMyGlobalIndices ();
1746 
1747  // This does NOT assume UVM; lgMapHost_ is a host pointer.
1748  const GO* lgMapHostRawPtr = lgMapHost_.data ();
1749  // The third argument forces ArrayView not to try to track memory
1750  // in a debug build. We have to use it because the memory does
1751  // not belong to a Teuchos memory management class.
1752  return Teuchos::ArrayView<const GO>(
1753  lgMapHostRawPtr,
1754  lgMapHost_.extent (0),
1755  Teuchos::RCP_DISABLE_NODE_LOOKUP);
1756  }
1757 
1758  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1760  return distributed_;
1761  }
1762 
1763  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1765  using Teuchos::TypeNameTraits;
1766  std::ostringstream os;
1767 
1768  os << "Tpetra::Map: {"
1769  << "LocalOrdinalType: " << TypeNameTraits<LocalOrdinal>::name ()
1770  << ", GlobalOrdinalType: " << TypeNameTraits<GlobalOrdinal>::name ()
1771  << ", NodeType: " << TypeNameTraits<Node>::name ();
1772  if (this->getObjectLabel () != "") {
1773  os << ", Label: \"" << this->getObjectLabel () << "\"";
1774  }
1775  os << ", Global number of entries: " << getGlobalNumElements ()
1776  << ", Number of processes: " << getComm ()->getSize ()
1777  << ", Uniform: " << (isUniform () ? "true" : "false")
1778  << ", Contiguous: " << (isContiguous () ? "true" : "false")
1779  << ", Distributed: " << (isDistributed () ? "true" : "false")
1780  << "}";
1781  return os.str ();
1782  }
1783 
1788  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1789  std::string
1791  localDescribeToString (const Teuchos::EVerbosityLevel vl) const
1792  {
1793  using LO = local_ordinal_type;
1794  using std::endl;
1795 
1796  // This preserves current behavior of Map.
1797  if (vl < Teuchos::VERB_HIGH) {
1798  return std::string ();
1799  }
1800  auto outStringP = Teuchos::rcp (new std::ostringstream ());
1801  Teuchos::RCP<Teuchos::FancyOStream> outp =
1802  Teuchos::getFancyOStream (outStringP);
1803  Teuchos::FancyOStream& out = *outp;
1804 
1805  auto comm = this->getComm ();
1806  const int myRank = comm->getRank ();
1807  const int numProcs = comm->getSize ();
1808  out << "Process " << myRank << " of " << numProcs << ":" << endl;
1809  Teuchos::OSTab tab1 (out);
1810 
1811  const LO numEnt = static_cast<LO> (this->getLocalNumElements ());
1812  out << "My number of entries: " << numEnt << endl
1813  << "My minimum global index: " << this->getMinGlobalIndex () << endl
1814  << "My maximum global index: " << this->getMaxGlobalIndex () << endl;
1815 
1816  if (vl == Teuchos::VERB_EXTREME) {
1817  out << "My global indices: [";
1818  const LO minLclInd = this->getMinLocalIndex ();
1819  for (LO k = 0; k < numEnt; ++k) {
1820  out << minLclInd + this->getGlobalElement (k);
1821  if (k + 1 < numEnt) {
1822  out << ", ";
1823  }
1824  }
1825  out << "]" << endl;
1826  }
1827 
1828  out.flush (); // make sure the ostringstream got everything
1829  return outStringP->str ();
1830  }
1831 
1832  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1833  void
1835  describe (Teuchos::FancyOStream &out,
1836  const Teuchos::EVerbosityLevel verbLevel) const
1837  {
1838  using Teuchos::TypeNameTraits;
1839  using Teuchos::VERB_DEFAULT;
1840  using Teuchos::VERB_NONE;
1841  using Teuchos::VERB_LOW;
1842  using Teuchos::VERB_HIGH;
1843  using std::endl;
1844  using LO = local_ordinal_type;
1845  using GO = global_ordinal_type;
1846  const Teuchos::EVerbosityLevel vl =
1847  (verbLevel == VERB_DEFAULT) ? VERB_LOW : verbLevel;
1848 
1849  if (vl == VERB_NONE) {
1850  return; // don't print anything
1851  }
1852  // If this Map's Comm is null, then the Map does not participate
1853  // in collective operations with the other processes. In that
1854  // case, it is not even legal to call this method. The reasonable
1855  // thing to do in that case is nothing.
1856  auto comm = this->getComm ();
1857  if (comm.is_null ()) {
1858  return;
1859  }
1860  const int myRank = comm->getRank ();
1861  const int numProcs = comm->getSize ();
1862 
1863  // Only Process 0 should touch the output stream, but this method
1864  // in general may need to do communication. Thus, we may need to
1865  // preserve the current tab level across multiple "if (myRank ==
1866  // 0) { ... }" inner scopes. This is why we sometimes create
1867  // OSTab instances by pointer, instead of by value. We only need
1868  // to create them by pointer if the tab level must persist through
1869  // multiple inner scopes.
1870  Teuchos::RCP<Teuchos::OSTab> tab0, tab1;
1871 
1872  if (myRank == 0) {
1873  // At every verbosity level but VERB_NONE, Process 0 prints.
1874  // By convention, describe() always begins with a tab before
1875  // printing.
1876  tab0 = Teuchos::rcp (new Teuchos::OSTab (out));
1877  out << "\"Tpetra::Map\":" << endl;
1878  tab1 = Teuchos::rcp (new Teuchos::OSTab (out));
1879  {
1880  out << "Template parameters:" << endl;
1881  Teuchos::OSTab tab2 (out);
1882  out << "LocalOrdinal: " << TypeNameTraits<LO>::name () << endl
1883  << "GlobalOrdinal: " << TypeNameTraits<GO>::name () << endl
1884  << "Node: " << TypeNameTraits<Node>::name () << endl;
1885  }
1886  const std::string label = this->getObjectLabel ();
1887  if (label != "") {
1888  out << "Label: \"" << label << "\"" << endl;
1889  }
1890  out << "Global number of entries: " << getGlobalNumElements () << endl
1891  << "Minimum global index: " << getMinAllGlobalIndex () << endl
1892  << "Maximum global index: " << getMaxAllGlobalIndex () << endl
1893  << "Index base: " << getIndexBase () << endl
1894  << "Number of processes: " << numProcs << endl
1895  << "Uniform: " << (isUniform () ? "true" : "false") << endl
1896  << "Contiguous: " << (isContiguous () ? "true" : "false") << endl
1897  << "Distributed: " << (isDistributed () ? "true" : "false") << endl;
1898  }
1899 
1900  // This is collective over the Map's communicator.
1901  if (vl >= VERB_HIGH) { // VERB_HIGH or VERB_EXTREME
1902  const std::string lclStr = this->localDescribeToString (vl);
1903  Tpetra::Details::gathervPrint (out, lclStr, *comm);
1904  }
1905  }
1906 
1907  template <class LocalOrdinal, class GlobalOrdinal, class Node>
1908  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
1910  replaceCommWithSubset (const Teuchos::RCP<const Teuchos::Comm<int> >& newComm) const
1911  {
1912  using Teuchos::RCP;
1913  using Teuchos::rcp;
1914  using GST = global_size_t;
1915  using LO = local_ordinal_type;
1916  using GO = global_ordinal_type;
1917  using map_type = Map<LO, GO, Node>;
1918 
1919  // mfh 26 Mar 2013: The lazy way to do this is simply to recreate
1920  // the Map by calling its ordinary public constructor, using the
1921  // original Map's data. This only involves O(1) all-reduces over
1922  // the new communicator, which in the common case only includes a
1923  // small number of processes.
1924 
1925  // Create the Map to return.
1926  if (newComm.is_null () || newComm->getSize () < 1) {
1927  return Teuchos::null; // my process does not participate in the new Map
1928  }
1929  else if (newComm->getSize () == 1) {
1930  // The case where the new communicator has only one process is
1931  // easy. We don't have to communicate to get all the
1932  // information we need. Use the default comm to create the new
1933  // Map, then fill in all the fields directly.
1934  RCP<map_type> newMap (new map_type ());
1935 
1936  newMap->comm_ = newComm;
1937  // mfh 07 Oct 2016: Preserve original behavior, even though the
1938  // original index base may no longer be the globally min global
1939  // index. See #616 for why this doesn't matter so much anymore.
1940  newMap->indexBase_ = this->indexBase_;
1941  newMap->numGlobalElements_ = this->numLocalElements_;
1942  newMap->numLocalElements_ = this->numLocalElements_;
1943  newMap->minMyGID_ = this->minMyGID_;
1944  newMap->maxMyGID_ = this->maxMyGID_;
1945  newMap->minAllGID_ = this->minMyGID_;
1946  newMap->maxAllGID_ = this->maxMyGID_;
1947  newMap->firstContiguousGID_ = this->firstContiguousGID_;
1948  newMap->lastContiguousGID_ = this->lastContiguousGID_;
1949  // Since the new communicator has only one process, neither
1950  // uniformity nor contiguity have changed.
1951  newMap->uniform_ = this->uniform_;
1952  newMap->contiguous_ = this->contiguous_;
1953  // The new communicator only has one process, so the new Map is
1954  // not distributed.
1955  newMap->distributed_ = false;
1956  newMap->lgMap_ = this->lgMap_;
1957  newMap->lgMapHost_ = this->lgMapHost_;
1958  newMap->glMap_ = this->glMap_;
1959  newMap->glMapHost_ = this->glMapHost_;
1960  // It's OK not to initialize the new Map's Directory.
1961  // This is initialized lazily, on first call to getRemoteIndexList.
1962 
1963  return newMap;
1964  }
1965  else { // newComm->getSize() != 1
1966  // Even if the original Map is contiguous, the new Map might not
1967  // be, especially if the excluded processes have ranks != 0 or
1968  // newComm->getSize()-1. The common case for this method is to
1969  // exclude many (possibly even all but one) processes, so it
1970  // likely doesn't pay to do the global communication (over the
1971  // original communicator) to figure out whether we can optimize
1972  // the result Map. Thus, we just set up the result Map as
1973  // noncontiguous.
1974  //
1975  // TODO (mfh 07 Oct 2016) We don't actually need to reconstruct
1976  // the global-to-local table, etc. Optimize this code path to
1977  // avoid unnecessary local work.
1978 
1979  // Make Map (re)compute the global number of elements.
1980  const GST RECOMPUTE = Tpetra::Details::OrdinalTraits<GST>::invalid ();
1981  // TODO (mfh 07 Oct 2016) If we use any Map constructor, we have
1982  // to use the noncontiguous Map constructor, since the new Map
1983  // might not be contiguous. Even if the old Map was contiguous,
1984  // some process in the "middle" might have been excluded. If we
1985  // want to avoid local work, we either have to do the setup by
1986  // hand, or write a new Map constructor.
1987 #if 1
1988  // The disabled code here throws the following exception in
1989  // Map's replaceCommWithSubset test:
1990  //
1991  // Throw test that evaluated to true: static_cast<unsigned long long> (numKeys) > static_cast<unsigned long long> (::Kokkos::Details::ArithTraits<ValueType>::max ())
1992  // 10:
1993  // 10: Tpetra::Details::FixedHashTable: The number of keys -3 is greater than the maximum representable ValueType value 2147483647. This means that it is not possible to use this constructor.
1994  // 10: Process 3: origComm->replaceCommWithSubset(subsetComm) threw an exception: /scratch/prj/Trilinos/Trilinos/packages/tpetra/core/src/Tpetra_Details_FixedHashTable_def.hpp:1044:
1995 
1996  auto lgMap = this->getMyGlobalIndices ();
1997  using size_type =
1998  typename std::decay<decltype (lgMap.extent (0)) >::type;
1999  const size_type lclNumInds =
2000  static_cast<size_type> (this->getLocalNumElements ());
2001  using Teuchos::TypeNameTraits;
2002  TEUCHOS_TEST_FOR_EXCEPTION
2003  (lgMap.extent (0) != lclNumInds, std::logic_error,
2004  "Tpetra::Map::replaceCommWithSubset: Result of getMyGlobalIndices() "
2005  "has length " << lgMap.extent (0) << " (of type " <<
2006  TypeNameTraits<size_type>::name () << ") != this->getLocalNumElements()"
2007  " = " << this->getLocalNumElements () << ". The latter, upon being "
2008  "cast to size_type = " << TypeNameTraits<size_type>::name () << ", "
2009  "becomes " << lclNumInds << ". Please report this bug to the Tpetra "
2010  "developers.");
2011 #else
2012  Teuchos::ArrayView<const GO> lgMap = this->getLocalElementList ();
2013 #endif // 1
2014 
2015  const GO indexBase = this->getIndexBase ();
2016  // map stores HostSpace of CudaSpace but constructor is still CudaUVMSpace
2017  auto lgMap_device = Kokkos::create_mirror_view_and_copy(device_type(), lgMap);
2018  return rcp (new map_type (RECOMPUTE, lgMap_device, indexBase, newComm));
2019  }
2020  }
2021 
2022  template <class LocalOrdinal, class GlobalOrdinal, class Node>
2023  Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
2026  {
2027  using Teuchos::Comm;
2028  using Teuchos::null;
2029  using Teuchos::outArg;
2030  using Teuchos::RCP;
2031  using Teuchos::rcp;
2032  using Teuchos::REDUCE_MIN;
2033  using Teuchos::reduceAll;
2034 
2035  // Create the new communicator. split() returns a valid
2036  // communicator on all processes. On processes where color == 0,
2037  // ignore the result. Passing key == 0 tells MPI to order the
2038  // processes in the new communicator by their rank in the old
2039  // communicator.
2040  const int color = (numLocalElements_ == 0) ? 0 : 1;
2041  // MPI_Comm_split must be called collectively over the original
2042  // communicator. We can't just call it on processes with color
2043  // one, even though we will ignore its result on processes with
2044  // color zero.
2045  RCP<const Comm<int> > newComm = comm_->split (color, 0);
2046  if (color == 0) {
2047  newComm = null;
2048  }
2049 
2050  // Create the Map to return.
2051  if (newComm.is_null ()) {
2052  return null; // my process does not participate in the new Map
2053  } else {
2054  RCP<Map> map = rcp (new Map ());
2055 
2056  map->comm_ = newComm;
2057  map->indexBase_ = indexBase_;
2058  map->numGlobalElements_ = numGlobalElements_;
2059  map->numLocalElements_ = numLocalElements_;
2060  map->minMyGID_ = minMyGID_;
2061  map->maxMyGID_ = maxMyGID_;
2062  map->minAllGID_ = minAllGID_;
2063  map->maxAllGID_ = maxAllGID_;
2064  map->firstContiguousGID_= firstContiguousGID_;
2065  map->lastContiguousGID_ = lastContiguousGID_;
2066 
2067  // Uniformity and contiguity have not changed. The directory
2068  // has changed, but we've taken care of that above.
2069  map->uniform_ = uniform_;
2070  map->contiguous_ = contiguous_;
2071 
2072  // If the original Map was NOT distributed, then the new Map
2073  // cannot be distributed.
2074  //
2075  // If the number of processes in the new communicator is 1, then
2076  // the new Map is not distributed.
2077  //
2078  // Otherwise, we have to check the new Map using an all-reduce
2079  // (over the new communicator). For example, the original Map
2080  // may have had some processes with zero elements, and all other
2081  // processes with the same number of elements as in the whole
2082  // Map. That Map is technically distributed, because of the
2083  // processes with zero elements. Removing those processes would
2084  // make the new Map locally replicated.
2085  if (! distributed_ || newComm->getSize () == 1) {
2086  map->distributed_ = false;
2087  } else {
2088  const int iOwnAllGids = (numLocalElements_ == numGlobalElements_) ? 1 : 0;
2089  int allProcsOwnAllGids = 0;
2090  reduceAll<int, int> (*newComm, REDUCE_MIN, iOwnAllGids, outArg (allProcsOwnAllGids));
2091  map->distributed_ = (allProcsOwnAllGids == 1) ? false : true;
2092  }
2093 
2094  map->lgMap_ = lgMap_;
2095  map->lgMapHost_ = lgMapHost_;
2096  map->glMap_ = glMap_;
2097  map->glMapHost_ = glMapHost_;
2098 
2099  // Map's default constructor creates an uninitialized Directory.
2100  // The Directory will be initialized on demand in
2101  // getRemoteIndexList().
2102  //
2103  // FIXME (mfh 26 Mar 2013) It should be possible to "filter" the
2104  // directory more efficiently than just recreating it. If
2105  // directory recreation proves a bottleneck, we can always
2106  // revisit this. On the other hand, Directory creation is only
2107  // collective over the new, presumably much smaller
2108  // communicator, so it may not be worth the effort to optimize.
2109 
2110  return map;
2111  }
2112  }
2113 
2114  template <class LocalOrdinal, class GlobalOrdinal, class Node>
2115  void
2117  {
2118  TEUCHOS_TEST_FOR_EXCEPTION(
2119  directory_.is_null (), std::logic_error, "Tpetra::Map::setupDirectory: "
2120  "The Directory is null. "
2121  "Please report this bug to the Tpetra developers.");
2122 
2123  // Only create the Directory if it hasn't been created yet.
2124  // This is a collective operation.
2125  if (! directory_->initialized ()) {
2126  directory_->initialize (*this);
2127  }
2128  }
2129 
2130  template <class LocalOrdinal, class GlobalOrdinal, class Node>
2131  LookupStatus
2133  getRemoteIndexList (const Teuchos::ArrayView<const GlobalOrdinal>& GIDs,
2134  const Teuchos::ArrayView<int>& PIDs,
2135  const Teuchos::ArrayView<LocalOrdinal>& LIDs) const
2136  {
2137  using Tpetra::Details::OrdinalTraits;
2139  using std::endl;
2140  using size_type = Teuchos::ArrayView<int>::size_type;
2141 
2142  const bool verbose = Details::Behavior::verbose("Map");
2143  const size_t maxNumToPrint = verbose ?
2145  std::unique_ptr<std::string> prefix;
2146  if (verbose) {
2147  prefix = Details::createPrefix(comm_.getRawPtr(),
2148  "Map", "getRemoteIndexList(GIDs,PIDs,LIDs)");
2149  std::ostringstream os;
2150  os << *prefix << "Start: ";
2151  verbosePrintArray(os, GIDs, "GIDs", maxNumToPrint);
2152  os << endl;
2153  std::cerr << os.str();
2154  }
2155 
2156  // Empty Maps (i.e., containing no indices on any processes in the
2157  // Map's communicator) are perfectly valid. In that case, if the
2158  // input GID list is nonempty, we fill the output arrays with
2159  // invalid values, and return IDNotPresent to notify the caller.
2160  // It's perfectly valid to give getRemoteIndexList GIDs that the
2161  // Map doesn't own. SubmapImport test 2 needs this functionality.
2162  if (getGlobalNumElements () == 0) {
2163  if (GIDs.size () == 0) {
2164  if (verbose) {
2165  std::ostringstream os;
2166  os << *prefix << "Done; both Map & input are empty" << endl;
2167  std::cerr << os.str();
2168  }
2169  return AllIDsPresent; // trivially
2170  }
2171  else {
2172  if (verbose) {
2173  std::ostringstream os;
2174  os << *prefix << "Done: Map is empty on all processes, "
2175  "so all output PIDs & LIDs are invalid (-1)." << endl;
2176  std::cerr << os.str();
2177  }
2178  for (size_type k = 0; k < PIDs.size (); ++k) {
2179  PIDs[k] = OrdinalTraits<int>::invalid ();
2180  }
2181  for (size_type k = 0; k < LIDs.size (); ++k) {
2182  LIDs[k] = OrdinalTraits<LocalOrdinal>::invalid ();
2183  }
2184  return IDNotPresent;
2185  }
2186  }
2187 
2188  // getRemoteIndexList must be called collectively, and Directory
2189  // initialization is collective too, so it's OK to initialize the
2190  // Directory on demand.
2191 
2192  if (verbose) {
2193  std::ostringstream os;
2194  os << *prefix << "Call setupDirectory" << endl;
2195  std::cerr << os.str();
2196  }
2197  setupDirectory();
2198  if (verbose) {
2199  std::ostringstream os;
2200  os << *prefix << "Call directory_->getDirectoryEntries" << endl;
2201  std::cerr << os.str();
2202  }
2203  const Tpetra::LookupStatus retVal =
2204  directory_->getDirectoryEntries (*this, GIDs, PIDs, LIDs);
2205  if (verbose) {
2206  std::ostringstream os;
2207  os << *prefix << "Done; getDirectoryEntries returned "
2208  << (retVal == IDNotPresent ? "IDNotPresent" : "AllIDsPresent")
2209  << "; ";
2210  verbosePrintArray(os, PIDs, "PIDs", maxNumToPrint);
2211  os << ", ";
2212  verbosePrintArray(os, LIDs, "LIDs", maxNumToPrint);
2213  os << endl;
2214  std::cerr << os.str();
2215  }
2216  return retVal;
2217  }
2218 
2219  template <class LocalOrdinal, class GlobalOrdinal, class Node>
2220  LookupStatus
2222  getRemoteIndexList (const Teuchos::ArrayView<const GlobalOrdinal> & GIDs,
2223  const Teuchos::ArrayView<int> & PIDs) const
2224  {
2226  using std::endl;
2227 
2228  const bool verbose = Details::Behavior::verbose("Map");
2229  const size_t maxNumToPrint = verbose ?
2231  std::unique_ptr<std::string> prefix;
2232  if (verbose) {
2233  prefix = Details::createPrefix(comm_.getRawPtr(),
2234  "Map", "getRemoteIndexList(GIDs,PIDs)");
2235  std::ostringstream os;
2236  os << *prefix << "Start: ";
2237  verbosePrintArray(os, GIDs, "GIDs", maxNumToPrint);
2238  os << endl;
2239  std::cerr << os.str();
2240  }
2241 
2242  if (getGlobalNumElements () == 0) {
2243  if (GIDs.size () == 0) {
2244  if (verbose) {
2245  std::ostringstream os;
2246  os << *prefix << "Done; both Map & input are empty" << endl;
2247  std::cerr << os.str();
2248  }
2249  return AllIDsPresent; // trivially
2250  }
2251  else {
2252  if (verbose) {
2253  std::ostringstream os;
2254  os << *prefix << "Done: Map is empty on all processes, "
2255  "so all output PIDs are invalid (-1)." << endl;
2256  std::cerr << os.str();
2257  }
2258  for (Teuchos::ArrayView<int>::size_type k = 0; k < PIDs.size (); ++k) {
2259  PIDs[k] = Tpetra::Details::OrdinalTraits<int>::invalid ();
2260  }
2261  return IDNotPresent;
2262  }
2263  }
2264 
2265  // getRemoteIndexList must be called collectively, and Directory
2266  // initialization is collective too, so it's OK to initialize the
2267  // Directory on demand.
2268 
2269  if (verbose) {
2270  std::ostringstream os;
2271  os << *prefix << "Call setupDirectory" << endl;
2272  std::cerr << os.str();
2273  }
2274  setupDirectory();
2275  if (verbose) {
2276  std::ostringstream os;
2277  os << *prefix << "Call directory_->getDirectoryEntries" << endl;
2278  std::cerr << os.str();
2279  }
2280  const Tpetra::LookupStatus retVal =
2281  directory_->getDirectoryEntries(*this, GIDs, PIDs);
2282  if (verbose) {
2283  std::ostringstream os;
2284  os << *prefix << "Done; getDirectoryEntries returned "
2285  << (retVal == IDNotPresent ? "IDNotPresent" : "AllIDsPresent")
2286  << "; ";
2287  verbosePrintArray(os, PIDs, "PIDs", maxNumToPrint);
2288  os << endl;
2289  std::cerr << os.str();
2290  }
2291  return retVal;
2292  }
2293 
2294  template <class LocalOrdinal, class GlobalOrdinal, class Node>
2295  Teuchos::RCP<const Teuchos::Comm<int> >
2297  return comm_;
2298  }
2299 
2300  template <class LocalOrdinal,class GlobalOrdinal, class Node>
2301  bool
2303  checkIsDist() const
2304  {
2305  using Teuchos::as;
2306  using Teuchos::outArg;
2307  using Teuchos::REDUCE_MIN;
2308  using Teuchos::reduceAll;
2309  using std::endl;
2310 
2311  const bool verbose = Details::Behavior::verbose("Map");
2312  std::unique_ptr<std::string> prefix;
2313  if (verbose) {
2314  prefix = Details::createPrefix(
2315  comm_.getRawPtr(), "Map", "checkIsDist");
2316  std::ostringstream os;
2317  os << *prefix << "Start" << endl;
2318  std::cerr << os.str();
2319  }
2320 
2321  bool global = false;
2322  if (comm_->getSize () > 1) {
2323  // The communicator has more than one process, but that doesn't
2324  // necessarily mean the Map is distributed.
2325  int localRep = 0;
2326  if (numGlobalElements_ == as<global_size_t> (numLocalElements_)) {
2327  // The number of local elements on this process equals the
2328  // number of global elements.
2329  //
2330  // NOTE (mfh 22 Nov 2011) Does this still work if there were
2331  // duplicates in the global ID list on input (the third Map
2332  // constructor), so that the number of local elements (which
2333  // are not duplicated) on this process could be less than the
2334  // number of global elements, even if this process owns all
2335  // the elements?
2336  localRep = 1;
2337  }
2338  int allLocalRep;
2339  reduceAll<int, int> (*comm_, REDUCE_MIN, localRep, outArg (allLocalRep));
2340  if (allLocalRep != 1) {
2341  // At least one process does not own all the elements.
2342  // This makes the Map a distributed Map.
2343  global = true;
2344  }
2345  }
2346  // If the communicator has only one process, then the Map is not
2347  // distributed.
2348 
2349  if (verbose) {
2350  std::ostringstream os;
2351  os << *prefix << "Done; global=" << (global ? "true" : "false")
2352  << endl;
2353  std::cerr << os.str();
2354  }
2355  return global;
2356  }
2357 
2358 } // namespace Tpetra
2359 
2360 template <class LocalOrdinal, class GlobalOrdinal>
2361 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal> >
2362 Tpetra::createLocalMap (const size_t numElements,
2363  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
2364 {
2365  typedef LocalOrdinal LO;
2366  typedef GlobalOrdinal GO;
2367  using NT = typename ::Tpetra::Map<LO, GO>::node_type;
2368  return createLocalMapWithNode<LO, GO, NT> (numElements, comm);
2369 }
2370 
2371 template <class LocalOrdinal, class GlobalOrdinal>
2372 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal> >
2373 Tpetra::createUniformContigMap (const global_size_t numElements,
2374  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
2375 {
2376  typedef LocalOrdinal LO;
2377  typedef GlobalOrdinal GO;
2378  using NT = typename ::Tpetra::Map<LO, GO>::node_type;
2379  return createUniformContigMapWithNode<LO, GO, NT> (numElements, comm);
2380 }
2381 
2382 template <class LocalOrdinal, class GlobalOrdinal, class Node>
2383 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >
2384 Tpetra::createUniformContigMapWithNode (const global_size_t numElements,
2385  const Teuchos::RCP<const Teuchos::Comm<int> >& comm
2386 )
2387 {
2388  using Teuchos::rcp;
2390  const GlobalOrdinal indexBase = static_cast<GlobalOrdinal> (0);
2391 
2392  return rcp (new map_type (numElements, indexBase, comm, GloballyDistributed));
2393 }
2394 
2395 template <class LocalOrdinal, class GlobalOrdinal, class Node>
2396 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >
2397 Tpetra::createLocalMapWithNode (const size_t numElements,
2398  const Teuchos::RCP<const Teuchos::Comm<int> >& comm
2399 )
2400 {
2401  using Tpetra::global_size_t;
2402  using Teuchos::rcp;
2404  const GlobalOrdinal indexBase = 0;
2405  const global_size_t globalNumElts = static_cast<global_size_t> (numElements);
2406 
2407  return rcp (new map_type (globalNumElts, indexBase, comm, LocallyReplicated));
2408 }
2409 
2410 template <class LocalOrdinal, class GlobalOrdinal, class Node>
2411 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >
2413  const size_t localNumElements,
2414  const Teuchos::RCP<const Teuchos::Comm<int> >& comm
2415 )
2416 {
2417  using Teuchos::rcp;
2419  const GlobalOrdinal indexBase = 0;
2420 
2421  return rcp (new map_type (numElements, localNumElements, indexBase, comm));
2422 }
2423 
2424 template <class LocalOrdinal, class GlobalOrdinal>
2425 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal> >
2427  const size_t localNumElements,
2428  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
2429 {
2430  typedef LocalOrdinal LO;
2431  typedef GlobalOrdinal GO;
2432  using NT = typename Tpetra::Map<LO, GO>::node_type;
2433 
2434  return Tpetra::createContigMapWithNode<LO, GO, NT> (numElements, localNumElements, comm);
2435 }
2436 
2437 template <class LocalOrdinal, class GlobalOrdinal>
2438 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal> >
2439 Tpetra::createNonContigMap(const Teuchos::ArrayView<const GlobalOrdinal>& elementList,
2440  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
2441 {
2442  typedef LocalOrdinal LO;
2443  typedef GlobalOrdinal GO;
2444  using NT = typename Tpetra::Map<LO, GO>::node_type;
2445 
2446  return Tpetra::createNonContigMapWithNode<LO, GO, NT> (elementList, comm);
2447 }
2448 
2449 template <class LocalOrdinal, class GlobalOrdinal, class Node>
2450 Teuchos::RCP< const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >
2451 Tpetra::createNonContigMapWithNode (const Teuchos::ArrayView<const GlobalOrdinal>& elementList,
2452  const Teuchos::RCP<const Teuchos::Comm<int> >& comm
2453 )
2454 {
2455  using Teuchos::rcp;
2457  using GST = Tpetra::global_size_t;
2458  const GST INV = Tpetra::Details::OrdinalTraits<GST>::invalid ();
2459  // FIXME (mfh 22 Jul 2016) This is what I found here, but maybe this
2460  // shouldn't be zero, given that the index base is supposed to equal
2461  // the globally min global index?
2462  const GlobalOrdinal indexBase = 0;
2463 
2464  return rcp (new map_type (INV, elementList, indexBase, comm));
2465 }
2466 
2467 template<class LO, class GO, class NT>
2468 Teuchos::RCP<const Tpetra::Map<LO, GO, NT> >
2469 Tpetra::createOneToOne (const Teuchos::RCP<const Tpetra::Map<LO, GO, NT> >& M)
2470 {
2472  using Teuchos::Array;
2473  using Teuchos::ArrayView;
2474  using Teuchos::as;
2475  using Teuchos::rcp;
2476  using std::cerr;
2477  using std::endl;
2478  using map_type = Tpetra::Map<LO, GO, NT>;
2479  using GST = global_size_t;
2480 
2481  const bool verbose = Details::Behavior::verbose("Map");
2482  std::unique_ptr<std::string> prefix;
2483  if (verbose) {
2484  auto comm = M.is_null() ? Teuchos::null : M->getComm();
2485  prefix = Details::createPrefix(
2486  comm.getRawPtr(), "createOneToOne(Map)");
2487  std::ostringstream os;
2488  os << *prefix << "Start" << endl;
2489  cerr << os.str();
2490  }
2491  const size_t maxNumToPrint = verbose ?
2492  Details::Behavior::verbosePrintCountThreshold() : size_t(0);
2493  const GST GINV = Tpetra::Details::OrdinalTraits<GST>::invalid ();
2494  const int myRank = M->getComm ()->getRank ();
2495 
2496  // Bypasses for special cases where either M is known to be
2497  // one-to-one, or the one-to-one version of M is easy to compute.
2498  // This is why we take M as an RCP, not as a const reference -- so
2499  // that we can return M itself if it is 1-to-1.
2500  if (! M->isDistributed ()) {
2501  // For a locally replicated Map, we assume that users want to push
2502  // all the GIDs to Process 0.
2503 
2504  // mfh 05 Nov 2013: getGlobalNumElements() does indeed return what
2505  // you think it should return, in this special case of a locally
2506  // replicated contiguous Map.
2507  const GST numGlobalEntries = M->getGlobalNumElements ();
2508  if (M->isContiguous()) {
2509  const size_t numLocalEntries =
2510  (myRank == 0) ? as<size_t>(numGlobalEntries) : size_t(0);
2511  if (verbose) {
2512  std::ostringstream os;
2513  os << *prefix << "Input is locally replicated & contiguous; "
2514  "numLocalEntries=" << numLocalEntries << endl;
2515  cerr << os.str ();
2516  }
2517  auto retMap =
2518  rcp(new map_type(numGlobalEntries, numLocalEntries,
2519  M->getIndexBase(), M->getComm()));
2520  if (verbose) {
2521  std::ostringstream os;
2522  os << *prefix << "Done" << endl;
2523  cerr << os.str ();
2524  }
2525  return retMap;
2526  }
2527  else {
2528  if (verbose) {
2529  std::ostringstream os;
2530  os << *prefix << "Input is locally replicated & noncontiguous"
2531  << endl;
2532  cerr << os.str ();
2533  }
2534  ArrayView<const GO> myGids =
2535  (myRank == 0) ? M->getLocalElementList() : Teuchos::null;
2536  auto retMap =
2537  rcp(new map_type(GINV, myGids(), M->getIndexBase(),
2538  M->getComm()));
2539  if (verbose) {
2540  std::ostringstream os;
2541  os << *prefix << "Done" << endl;
2542  cerr << os.str ();
2543  }
2544  return retMap;
2545  }
2546  }
2547  else if (M->isContiguous ()) {
2548  if (verbose) {
2549  std::ostringstream os;
2550  os << *prefix << "Input is distributed & contiguous" << endl;
2551  cerr << os.str ();
2552  }
2553  // Contiguous, distributed Maps are one-to-one by construction.
2554  // (Locally replicated Maps can be contiguous.)
2555  return M;
2556  }
2557  else {
2558  if (verbose) {
2559  std::ostringstream os;
2560  os << *prefix << "Input is distributed & noncontiguous" << endl;
2561  cerr << os.str ();
2562  }
2564  const size_t numMyElems = M->getLocalNumElements ();
2565  ArrayView<const GO> myElems = M->getLocalElementList ();
2566  Array<int> owner_procs_vec (numMyElems);
2567 
2568  if (verbose) {
2569  std::ostringstream os;
2570  os << *prefix << "Call Directory::getDirectoryEntries: ";
2571  verbosePrintArray(os, myElems, "GIDs", maxNumToPrint);
2572  os << endl;
2573  cerr << os.str();
2574  }
2575  directory.getDirectoryEntries (*M, myElems, owner_procs_vec ());
2576  if (verbose) {
2577  std::ostringstream os;
2578  os << *prefix << "getDirectoryEntries result: ";
2579  verbosePrintArray(os, owner_procs_vec, "PIDs", maxNumToPrint);
2580  os << endl;
2581  cerr << os.str();
2582  }
2583 
2584  Array<GO> myOwned_vec (numMyElems);
2585  size_t numMyOwnedElems = 0;
2586  for (size_t i = 0; i < numMyElems; ++i) {
2587  const GO GID = myElems[i];
2588  const int owner = owner_procs_vec[i];
2589 
2590  if (myRank == owner) {
2591  myOwned_vec[numMyOwnedElems++] = GID;
2592  }
2593  }
2594  myOwned_vec.resize (numMyOwnedElems);
2595 
2596  if (verbose) {
2597  std::ostringstream os;
2598  os << *prefix << "Create Map: ";
2599  verbosePrintArray(os, myOwned_vec, "GIDs", maxNumToPrint);
2600  os << endl;
2601  cerr << os.str();
2602  }
2603  auto retMap = rcp(new map_type(GINV, myOwned_vec(),
2604  M->getIndexBase(), M->getComm()));
2605  if (verbose) {
2606  std::ostringstream os;
2607  os << *prefix << "Done" << endl;
2608  cerr << os.str();
2609  }
2610  return retMap;
2611  }
2612 }
2613 
2614 template<class LocalOrdinal, class GlobalOrdinal, class Node>
2615 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
2618 {
2619  using Details::Behavior;
2621  using Teuchos::Array;
2622  using Teuchos::ArrayView;
2623  using Teuchos::RCP;
2624  using Teuchos::rcp;
2625  using Teuchos::toString;
2626  using std::cerr;
2627  using std::endl;
2628  using LO = LocalOrdinal;
2629  using GO = GlobalOrdinal;
2630  using map_type = Tpetra::Map<LO, GO, Node>;
2631 
2632  const bool verbose = Behavior::verbose("Map");
2633  std::unique_ptr<std::string> prefix;
2634  if (verbose) {
2635  auto comm = M.is_null() ? Teuchos::null : M->getComm();
2636  prefix = Details::createPrefix(
2637  comm.getRawPtr(), "createOneToOne(Map,TieBreak)");
2638  std::ostringstream os;
2639  os << *prefix << "Start" << endl;
2640  cerr << os.str();
2641  }
2642  const size_t maxNumToPrint = verbose ?
2643  Behavior::verbosePrintCountThreshold() : size_t(0);
2644 
2645  // FIXME (mfh 20 Feb 2013) We should have a bypass for contiguous
2646  // Maps (which are 1-to-1 by construction).
2647 
2649  if (verbose) {
2650  std::ostringstream os;
2651  os << *prefix << "Initialize Directory" << endl;
2652  cerr << os.str ();
2653  }
2654  directory.initialize (*M, tie_break);
2655  if (verbose) {
2656  std::ostringstream os;
2657  os << *prefix << "Done initializing Directory" << endl;
2658  cerr << os.str ();
2659  }
2660  size_t numMyElems = M->getLocalNumElements ();
2661  ArrayView<const GO> myElems = M->getLocalElementList ();
2662  Array<int> owner_procs_vec (numMyElems);
2663  if (verbose) {
2664  std::ostringstream os;
2665  os << *prefix << "Call Directory::getDirectoryEntries: ";
2666  verbosePrintArray(os, myElems, "GIDs", maxNumToPrint);
2667  os << endl;
2668  cerr << os.str();
2669  }
2670  directory.getDirectoryEntries (*M, myElems, owner_procs_vec ());
2671  if (verbose) {
2672  std::ostringstream os;
2673  os << *prefix << "getDirectoryEntries result: ";
2674  verbosePrintArray(os, owner_procs_vec, "PIDs", maxNumToPrint);
2675  os << endl;
2676  cerr << os.str();
2677  }
2678 
2679  const int myRank = M->getComm()->getRank();
2680  Array<GO> myOwned_vec (numMyElems);
2681  size_t numMyOwnedElems = 0;
2682  for (size_t i = 0; i < numMyElems; ++i) {
2683  const GO GID = myElems[i];
2684  const int owner = owner_procs_vec[i];
2685  if (myRank == owner) {
2686  myOwned_vec[numMyOwnedElems++] = GID;
2687  }
2688  }
2689  myOwned_vec.resize (numMyOwnedElems);
2690 
2691  // FIXME (mfh 08 May 2014) The above Directory should be perfectly
2692  // valid for the new Map. Why can't we reuse it?
2693  const global_size_t GINV =
2694  Tpetra::Details::OrdinalTraits<global_size_t>::invalid ();
2695  if (verbose) {
2696  std::ostringstream os;
2697  os << *prefix << "Create Map: ";
2698  verbosePrintArray(os, myOwned_vec, "GIDs", maxNumToPrint);
2699  os << endl;
2700  cerr << os.str();
2701  }
2702  RCP<const map_type> retMap
2703  (new map_type (GINV, myOwned_vec (), M->getIndexBase (),
2704  M->getComm ()));
2705  if (verbose) {
2706  std::ostringstream os;
2707  os << *prefix << "Done" << endl;
2708  cerr << os.str();
2709  }
2710  return retMap;
2711 }
2712 
2713 //
2714 // Explicit instantiation macro
2715 //
2716 // Must be expanded from within the Tpetra namespace!
2717 //
2718 
2720 
2721 #define TPETRA_MAP_INSTANT(LO,GO,NODE) \
2722  \
2723  template class Map< LO , GO , NODE >; \
2724  \
2725  template Teuchos::RCP< const Map<LO,GO,NODE> > \
2726  createLocalMapWithNode<LO,GO,NODE> (const size_t numElements, \
2727  const Teuchos::RCP< const Teuchos::Comm< int > >& comm); \
2728  \
2729  template Teuchos::RCP< const Map<LO,GO,NODE> > \
2730  createContigMapWithNode<LO,GO,NODE> (const global_size_t numElements, \
2731  const size_t localNumElements, \
2732  const Teuchos::RCP< const Teuchos::Comm< int > >& comm); \
2733  \
2734  template Teuchos::RCP< const Map<LO,GO,NODE> > \
2735  createNonContigMapWithNode(const Teuchos::ArrayView<const GO> &elementList, \
2736  const Teuchos::RCP<const Teuchos::Comm<int> > &comm); \
2737  \
2738  template Teuchos::RCP< const Map<LO,GO,NODE> > \
2739  createUniformContigMapWithNode<LO,GO,NODE> (const global_size_t numElements, \
2740  const Teuchos::RCP< const Teuchos::Comm< int > >& comm); \
2741  \
2742  template Teuchos::RCP<const Map<LO,GO,NODE> > \
2743  createOneToOne (const Teuchos::RCP<const Map<LO,GO,NODE> >& M); \
2744  \
2745  template Teuchos::RCP<const Map<LO,GO,NODE> > \
2746  createOneToOne (const Teuchos::RCP<const Map<LO,GO,NODE> >& M, \
2747  const Tpetra::Details::TieBreak<LO,GO>& tie_break); \
2748 
2749 
2751 #define TPETRA_MAP_INSTANT_DEFAULTNODE(LO,GO) \
2752  template Teuchos::RCP< const Map<LO,GO> > \
2753  createLocalMap<LO,GO>( const size_t, const Teuchos::RCP< const Teuchos::Comm< int > > &); \
2754  \
2755  template Teuchos::RCP< const Map<LO,GO> > \
2756  createContigMap<LO,GO>( global_size_t, size_t, \
2757  const Teuchos::RCP< const Teuchos::Comm< int > > &); \
2758  \
2759  template Teuchos::RCP< const Map<LO,GO> > \
2760  createNonContigMap(const Teuchos::ArrayView<const GO> &, \
2761  const Teuchos::RCP<const Teuchos::Comm<int> > &); \
2762  \
2763  template Teuchos::RCP< const Map<LO,GO> > \
2764  createUniformContigMap<LO,GO>( const global_size_t, \
2765  const Teuchos::RCP< const Teuchos::Comm< int > > &); \
2766 
2767 #endif // TPETRA_MAP_DEF_HPP
Interface for breaking ties in ownership.
bool congruent(const Teuchos::Comm< int > &comm1, const Teuchos::Comm< int > &comm2)
Whether the two communicators are congruent.
Definition: Tpetra_Util.cpp:64
global_ordinal_type getMaxGlobalIndex() const
The maximum global index owned by the calling process.
LocalOrdinal local_ordinal_type
The type of local indices.
Declaration of Tpetra::Details::extractMpiCommFromTeuchos.
virtual ~Map()
Destructor (virtual for memory safety of derived classes).
GlobalOrdinal global_ordinal_type
The type of global indices.
bool isSameAs(const Map< local_ordinal_type, global_ordinal_type, Node > &map) const
True if and only if map is identical to this Map.
Declaration of Tpetra::Details::printOnce.
global_indices_array_type getMyGlobalIndices() const
Return a view of the global indices owned by this process.
bool isContiguous() const
True if this Map is distributed contiguously, else false.
LookupStatus
Return status of Map remote index lookup (getRemoteIndexList()).
Declaration of a function that prints strings from each process.
size_t getLocalNumElements() const
The number of elements belonging to the calling process.
bool isLocallyFitted(const Map< local_ordinal_type, global_ordinal_type, Node > &map) const
True if and only if map is locally fitted to this Map.
bool isUniform() const
Whether the range of global indices is uniform.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createContigMapWithNode(const global_size_t numElements, const size_t localNumElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Nonmember constructor for a (potentially) nonuniformly distributed, contiguous Map for a user-specifi...
bool isDistributed() const
Whether this Map is globally distributed or locally replicated.
Teuchos::RCP< const Map< local_ordinal_type, global_ordinal_type, Node > > removeEmptyProcesses() const
Advanced methods.
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Accessors for the Teuchos::Comm and Kokkos Node objects.
static bool debug()
Whether Tpetra is in debug mode.
void verbosePrintArray(std::ostream &out, const ArrayType &x, const char name[], const size_t maxNumToPrint)
Print min(x.size(), maxNumToPrint) entries of x.
int local_ordinal_type
Default value of Scalar template parameter.
&quot;Local&quot; part of Map suitable for Kokkos kernels.
Declaration of Tpetra::Details::initializeKokkos.
void gathervPrint(std::ostream &out, const std::string &s, const Teuchos::Comm< int > &comm)
On Process 0 in the given communicator, print strings from each process in that communicator, in rank order.
Teuchos::ArrayView< const global_ordinal_type > getLocalElementList() const
Return a NONOWNING view of the global indices owned by this process.
void initializeKokkos()
Initialize Kokkos, using command-line arguments (if any) given to Teuchos::GlobalMPISession.
size_t global_size_t
Global size_t object.
void deep_copy(MultiVector< DS, DL, DG, DN > &dst, const MultiVector< SS, SL, SG, SN > &src)
Copy the contents of the MultiVector src into dst.
bool isCompatible(const Map< local_ordinal_type, global_ordinal_type, Node > &map) const
True if and only if map is compatible with this Map.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createOneToOne(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &M)
Nonmember constructor for a contiguous Map with user-defined weights and a user-specified, possibly nondefault Kokkos Node type.
global_ordinal_type getIndexBase() const
The index base for this Map.
global_ordinal_type getGlobalElement(local_ordinal_type localIndex) const
The global index corresponding to the given local index.
bool isNodeLocalElement(local_ordinal_type localIndex) const
Whether the given local index is valid for this Map on the calling process.
Teuchos::RCP< const Map< local_ordinal_type, global_ordinal_type, Node > > replaceCommWithSubset(const Teuchos::RCP< const Teuchos::Comm< int > > &newComm) const
Replace this Map&#39;s communicator with a subset communicator.
Functions for initializing and finalizing Tpetra.
bool isNodeGlobalElement(global_ordinal_type globalIndex) const
Whether the given global index is owned by this Map on the calling process.
static bool verbose()
Whether Tpetra is in verbose mode.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal > > createUniformContigMap(const global_size_t numElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Non-member constructor for a uniformly distributed, contiguous Map with the default Kokkos Node...
LookupStatus getDirectoryEntries(const map_type &map, const Teuchos::ArrayView< const GlobalOrdinal > &globalIDs, const Teuchos::ArrayView< int > &nodeIDs) const
Given a global ID list, return the list of their owning process IDs.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal > > createContigMap(const global_size_t numElements, const size_t localNumElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Non-member constructor for a (potentially) non-uniformly distributed, contiguous Map using the defaul...
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal > > createLocalMap(const size_t numElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Nonmember constructor for a locally replicated Map with the default Kokkos Node.
local_map_type getLocalMap() const
Get the local Map for Kokkos kernels.
global_ordinal_type getMinGlobalIndex() const
The minimum global index owned by the calling process.
typename device_type::execution_space execution_space
The Kokkos execution space.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createUniformContigMapWithNode(const global_size_t numElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Non-member constructor for a uniformly distributed, contiguous Map with a user-specified Kokkos Node...
A parallel distribution of indices over processes.
global_ordinal_type getMaxAllGlobalIndex() const
The maximum global index over all processes in the communicator.
local_ordinal_type getLocalElement(global_ordinal_type globalIndex) const
The local index corresponding to the given global index.
Implement mapping from global ID to process ID and local ID.
Stand-alone utility functions and macros.
bool isOneToOne() const
Whether the Map is one to one.
void initialize(const map_type &map)
Initialize the Directory with its Map.
bool locallySameAs(const Map< local_ordinal_type, global_ordinal_type, node_type > &map) const
Is this Map locally the same as the input Map?
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createLocalMapWithNode(const size_t numElements, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Nonmember constructor for a locally replicated Map with a specified Kokkos Node.
static size_t verbosePrintCountThreshold()
Number of entries below which arrays, lists, etc. will be printed in debug mode.
typename Node::device_type device_type
This class&#39; Kokkos::Device specialization.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal > > createNonContigMap(const Teuchos::ArrayView< const GlobalOrdinal > &elementList, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Nonmember constructor for a non-contiguous Map using the default Kokkos::Device type.
LocalGlobal
Enum for local versus global allocation of Map entries.
LookupStatus getRemoteIndexList(const Teuchos::ArrayView< const global_ordinal_type > &GIDList, const Teuchos::ArrayView< int > &nodeIDList, const Teuchos::ArrayView< local_ordinal_type > &LIDList) const
Return the process ranks and corresponding local indices for the given global indices.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Describe this object in a human-readable way to the given output stream.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > createNonContigMapWithNode(const Teuchos::ArrayView< const GlobalOrdinal > &elementList, const Teuchos::RCP< const Teuchos::Comm< int > > &comm)
Nonmember constructor for a noncontiguous Map with a user-specified, possibly nondefault Kokkos Node ...
Node node_type
Legacy typedef that will go away at some point.
global_ordinal_type getMinAllGlobalIndex() const
The minimum global index over all processes in the communicator.
std::unique_ptr< std::string > createPrefix(const int myRank, const char prefix[])
Create string prefix for each line of verbose output.
global_size_t getGlobalNumElements() const
The number of elements in this Map.
Description of Tpetra&#39;s behavior.
std::string description() const
Implementation of Teuchos::Describable.
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra&#39;s behavior.
Map()
Default constructor (that does nothing).