Tpetra parallel linear algebra Version of the Day
Loading...
Searching...
No Matches
Tpetra_Import_Util2.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Tpetra: Templated Linear Algebra Services Package
5// Copyright (2008) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ************************************************************************
40// @HEADER
41
42#ifndef TPETRA_IMPORT_UTIL2_HPP
43#define TPETRA_IMPORT_UTIL2_HPP
44
49
50#include "Tpetra_ConfigDefs.hpp"
51#include "Tpetra_Import.hpp"
52#include "Tpetra_HashTable.hpp"
53#include "Tpetra_Map.hpp"
54#include "Tpetra_Util.hpp"
55#include "Tpetra_Distributor.hpp"
58#include "Tpetra_Vector.hpp"
59#include "Kokkos_DualView.hpp"
60#include <Teuchos_Array.hpp>
61#include <utility>
62#include <set>
63
65
66namespace Tpetra {
67namespace Import_Util {
68
71template<typename Scalar, typename Ordinal>
72void
73sortCrsEntries (const Teuchos::ArrayView<size_t>& CRS_rowptr,
74 const Teuchos::ArrayView<Ordinal>& CRS_colind,
75 const Teuchos::ArrayView<Scalar>&CRS_vals);
76
77template<typename Ordinal>
78void
79sortCrsEntries (const Teuchos::ArrayView<size_t>& CRS_rowptr,
80 const Teuchos::ArrayView<Ordinal>& CRS_colind);
81
82template<typename rowptr_array_type, typename colind_array_type, typename vals_array_type>
83void
84sortCrsEntries (const rowptr_array_type& CRS_rowptr,
85 const colind_array_type& CRS_colind,
86 const vals_array_type& CRS_vals);
87
88template<typename rowptr_array_type, typename colind_array_type>
89void
90sortCrsEntries (const rowptr_array_type& CRS_rowptr,
91 const colind_array_type& CRS_colind);
92
97template<typename Scalar, typename Ordinal>
98void
99sortAndMergeCrsEntries (const Teuchos::ArrayView<size_t>& CRS_rowptr,
100 const Teuchos::ArrayView<Ordinal>& CRS_colind,
101 const Teuchos::ArrayView<Scalar>& CRS_vals);
102
103template<typename Ordinal>
104void
105sortAndMergeCrsEntries (const Teuchos::ArrayView<size_t>& CRS_rowptr,
106 const Teuchos::ArrayView<Ordinal>& CRS_colind);
107
123template <typename LocalOrdinal, typename GlobalOrdinal, typename Node>
124void
125lowCommunicationMakeColMapAndReindex (const Teuchos::ArrayView<const size_t> &rowPointers,
126 const Teuchos::ArrayView<LocalOrdinal> &columnIndices_LID,
127 const Teuchos::ArrayView<GlobalOrdinal> &columnIndices_GID,
128 const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & domainMap,
129 const Teuchos::ArrayView<const int> &owningPids,
130 Teuchos::Array<int> &remotePids,
131 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & colMap);
132
133
134
135
149 template <typename LocalOrdinal, typename GlobalOrdinal, typename Node>
150 void getTwoTransferOwnershipVector(const ::Tpetra::Details::Transfer<LocalOrdinal, GlobalOrdinal, Node>& transferThatDefinesOwnership,
151 bool useReverseModeForOwnership,
152 const ::Tpetra::Details::Transfer<LocalOrdinal, GlobalOrdinal, Node>& transferForMigratingData,
153 bool useReverseModeForMigration,
155
156} // namespace Import_Util
157} // namespace Tpetra
158
159
160//
161// Implementations
162//
163
164namespace Tpetra {
165namespace Import_Util {
166
167
168template<typename PID, typename GlobalOrdinal>
169bool sort_PID_then_GID(const std::pair<PID,GlobalOrdinal> &a,
170 const std::pair<PID,GlobalOrdinal> &b)
171{
172 if(a.first!=b.first)
173 return (a.first < b.first);
174 return (a.second < b.second);
175}
176
177template<typename PID,
178 typename GlobalOrdinal,
179 typename LocalOrdinal>
180bool sort_PID_then_pair_GID_LID(const std::pair<PID, std::pair< GlobalOrdinal, LocalOrdinal > > &a,
181 const std::pair<PID, std::pair< GlobalOrdinal, LocalOrdinal > > &b)
182{
183 if(a.first!=b.first)
184 return a.first < b.first;
185 else
186 return (a.second.first < b.second.first);
187}
188
189template<typename Scalar,
190 typename LocalOrdinal,
191 typename GlobalOrdinal,
192 typename Node>
193void
194reverseNeighborDiscovery(const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & SourceMatrix,
195 const typename CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::row_ptrs_host_view_type & rowptr,
196 const typename CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type & colind,
198 Teuchos::RCP<const Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> > MyImporter,
199 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > MyDomainMap,
200 Teuchos::ArrayRCP<int>& type3PIDs,
201 Teuchos::ArrayRCP<LocalOrdinal>& type3LIDs,
202 Teuchos::RCP<const Teuchos::Comm<int> >& rcomm)
203{
204#ifdef HAVE_TPETRACORE_MPI
205 using Teuchos::TimeMonitor;
206 using ::Tpetra::Details::Behavior;
207 using Kokkos::AllowPadding;
208 using Kokkos::view_alloc;
209 using Kokkos::WithoutInitializing;
210 typedef LocalOrdinal LO;
211 typedef GlobalOrdinal GO;
212 typedef std::pair<GO,GO> pidgidpair_t;
213 using Teuchos::RCP;
214 const std::string prefix {" Import_Util2::ReverseND:: "};
215 const std::string label ("IU2::Neighbor");
216
217 // There can be no neighbor discovery if you don't have an importer
218 if(MyImporter.is_null()) return;
219
220 std::ostringstream errstr;
221 bool error = false;
222 auto const comm = MyDomainMap->getComm();
223
224 MPI_Comm rawComm = getRawMpiComm(*comm);
225 const int MyPID = rcomm->getRank ();
226
227 // Things related to messages I am sending in forward mode (RowTransfer)
228 // *** Note: this will be incorrect for transferAndFillComplete if it is in reverse mode. FIXME cbl.
229 auto ExportPIDs = RowTransfer.getExportPIDs();
230 auto ExportLIDs = RowTransfer.getExportLIDs();
231 auto NumExportLIDs = RowTransfer.getNumExportIDs();
232
233 Distributor & Distor = MyImporter->getDistributor();
234 const size_t NumRecvs = Distor.getNumReceives();
235 const size_t NumSends = Distor.getNumSends();
236 auto RemoteLIDs = MyImporter->getRemoteLIDs();
237 auto const ProcsFrom = Distor.getProcsFrom();
238 auto const ProcsTo = Distor.getProcsTo();
239
240 auto LengthsFrom = Distor.getLengthsFrom();
241 auto MyColMap = SourceMatrix.getColMap();
242 const size_t numCols = MyColMap->getLocalNumElements ();
243 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > target = MyImporter->getTargetMap();
244
245 // Get the owning pids in a special way,
246 // s.t. ProcsFrom[RemotePIDs[i]] is the proc that owns RemoteLIDs[j]....
247 Teuchos::Array<int> RemotePIDOrder(numCols,-1);
248
249 // For each remote ID, record index into ProcsFrom, who owns it.
250 for (size_t i = 0, j = 0; i < NumRecvs; ++i) {
251 for (size_t k = 0; k < LengthsFrom[i]; ++k) {
252 const int pid = ProcsFrom[i];
253 if (pid != MyPID) {
254 RemotePIDOrder[RemoteLIDs[j]] = i;
255 }
256 j++;
257 }
258 }
259
260 // Step One: Start tacking the (GID,PID) pairs on the std sets
261 //
262 // For each index in ProcsFrom, we will insert into a set of (PID,
263 // GID) pairs, in order to build a list of such pairs for each of
264 // those processes. Since this is building a reverse, we will send
265 // to these processes.
266 Teuchos::Array<int> ReverseSendSizes(NumRecvs,0);
267 // do this as C array to avoid Teuchos::Array value initialization of all reserved memory
268 Teuchos::Array< Teuchos::ArrayRCP<pidgidpair_t > > RSB(NumRecvs);
269
270 {
271#ifdef HAVE_TPETRA_MMM_TIMINGS
272 TimeMonitor set_all(*TimeMonitor::getNewTimer(prefix + std::string("isMMallSetRSB")));
273#endif
274
275 // 25 Jul 2018: CBL
276 // todo:std::unordered_set (hash table),
277 // with an adequate prereservation ("bucket count").
278 // An onordered_set has to have a custom hasher for pid/gid pair
279 // However, when pidsets is copied to RSB, it will be in key
280 // order _not_ in pid,gid order. (unlike std::set).
281 // Impliment this with a reserve, and time BOTH building pidsets
282 // _and_ the sort after the receive. Even if unordered_set saves
283 // time, if it causes the sort to be longer, it's not a win.
284
285 Teuchos::Array<std::set<pidgidpair_t>> pidsets(NumRecvs);
286 {
287#ifdef HAVE_TPETRA_MMM_TIMINGS
288 TimeMonitor set_insert(*TimeMonitor::getNewTimer(prefix + std::string("isMMallSetRSBinsert")));
289#endif
290 for(size_t i=0; i < NumExportLIDs; i++) {
291 LO lid = ExportLIDs[i];
292 GO exp_pid = ExportPIDs[i];
293 for(auto j=rowptr[lid]; j<rowptr[lid+1]; j++){
294 int pid_order = RemotePIDOrder[colind[j]];
295 if(pid_order!=-1) {
296 GO gid = MyColMap->getGlobalElement(colind[j]); //Epetra SM.GCID46 =>sm->graph-> {colmap(colind)}
297 auto tpair = pidgidpair_t(exp_pid,gid);
298 pidsets[pid_order].insert(pidsets[pid_order].end(),tpair);
299 }
300 }
301 }
302 }
303
304 {
305#ifdef HAVE_TPETRA_MMM_TIMINGS
306 TimeMonitor set_cpy(*TimeMonitor::getNewTimer(prefix + std::string("isMMallSetRSBcpy")));
307#endif
308 int jj = 0;
309 for(auto &&ps : pidsets) {
310 auto s = ps.size();
311 RSB[jj] = Teuchos::arcp(new pidgidpair_t[ s ],0, s ,true);
312 std::copy(ps.begin(),ps.end(),RSB[jj]);
313 ReverseSendSizes[jj]=s;
314 ++jj;
315 }
316 }
317 } // end of set based packing.
318
319 Teuchos::Array<int> ReverseRecvSizes(NumSends,-1);
320 Teuchos::Array<MPI_Request> rawBreq(ProcsFrom.size()+ProcsTo.size(), MPI_REQUEST_NULL);
321 // 25 Jul 2018: MPI_TAG_UB is the largest tag value; could be < 32768.
322 const int mpi_tag_base_ = 3;
323
324 int mpireq_idx=0;
325 for(int i=0;i<ProcsTo.size();++i) {
326 int Rec_Tag = mpi_tag_base_ + ProcsTo[i];
327 int * thisrecv = (int *) (&ReverseRecvSizes[i]);
328 MPI_Request rawRequest = MPI_REQUEST_NULL;
329 MPI_Irecv(const_cast<int*>(thisrecv),
330 1,
331 MPI_INT,
332 ProcsTo[i],
333 Rec_Tag,
334 rawComm,
335 &rawRequest);
336 rawBreq[mpireq_idx++]=rawRequest;
337 }
338 for(int i=0;i<ProcsFrom.size();++i) {
339 int Send_Tag = mpi_tag_base_ + MyPID;
340 int * mysend = ( int *)(&ReverseSendSizes[i]);
341 MPI_Request rawRequest = MPI_REQUEST_NULL;
342 MPI_Isend(mysend,
343 1,
344 MPI_INT,
345 ProcsFrom[i],
346 Send_Tag,
347 rawComm,
348 &rawRequest);
349 rawBreq[mpireq_idx++]=rawRequest;
350 }
351 Teuchos::Array<MPI_Status> rawBstatus(rawBreq.size());
352#ifdef HAVE_TPETRA_DEBUG
353 const int err1 =
354#endif
355 MPI_Waitall (rawBreq.size(), rawBreq.getRawPtr(),
356 rawBstatus.getRawPtr());
357
358
359#ifdef HAVE_TPETRA_DEBUG
360 if(err1) {
361 errstr <<MyPID<< "sE1 reverseNeighborDiscovery Mpi_Waitall error on send ";
362 error=true;
363 std::cerr<<errstr.str()<<std::flush;
364 }
365#endif
366
367 int totalexportpairrecsize = 0;
368 for (size_t i = 0; i < NumSends; ++i) {
369 totalexportpairrecsize += ReverseRecvSizes[i];
370#ifdef HAVE_TPETRA_DEBUG
371 if(ReverseRecvSizes[i]<0) {
372 errstr << MyPID << "E4 reverseNeighborDiscovery: got < 0 for receive size "<<ReverseRecvSizes[i]<<std::endl;
373 error=true;
374 }
375#endif
376 }
377 Teuchos::ArrayRCP<pidgidpair_t >AllReverseRecv= Teuchos::arcp(new pidgidpair_t[totalexportpairrecsize],0,totalexportpairrecsize,true);
378 int offset = 0;
379 mpireq_idx=0;
380 for(int i=0;i<ProcsTo.size();++i) {
381 int recv_data_size = ReverseRecvSizes[i]*2;
382 int recvData_MPI_Tag = mpi_tag_base_*2 + ProcsTo[i];
383 MPI_Request rawRequest = MPI_REQUEST_NULL;
384 GO * rec_bptr = (GO*) (&AllReverseRecv[offset]);
385 offset+=ReverseRecvSizes[i];
386 MPI_Irecv(rec_bptr,
387 recv_data_size,
388 ::Tpetra::Details::MpiTypeTraits<GO>::getType(rec_bptr[0]),
389 ProcsTo[i],
390 recvData_MPI_Tag,
391 rawComm,
392 &rawRequest);
393 rawBreq[mpireq_idx++]=rawRequest;
394 }
395 for(int ii=0;ii<ProcsFrom.size();++ii) {
396 GO * send_bptr = (GO*) (RSB[ii].getRawPtr());
397 MPI_Request rawSequest = MPI_REQUEST_NULL;
398 int send_data_size = ReverseSendSizes[ii]*2; // 2 == count of pair
399 int sendData_MPI_Tag = mpi_tag_base_*2+MyPID;
400 MPI_Isend(send_bptr,
401 send_data_size,
402 ::Tpetra::Details::MpiTypeTraits<GO>::getType(send_bptr[0]),
403 ProcsFrom[ii],
404 sendData_MPI_Tag,
405 rawComm,
406 &rawSequest);
407
408 rawBreq[mpireq_idx++]=rawSequest;
409 }
410#ifdef HAVE_TPETRA_DEBUG
411 const int err =
412#endif
413 MPI_Waitall (rawBreq.size(),
414 rawBreq.getRawPtr(),
415 rawBstatus.getRawPtr());
416#ifdef HAVE_TPETRA_DEBUG
417 if(err) {
418 errstr <<MyPID<< "E3.r reverseNeighborDiscovery Mpi_Waitall error on receive ";
419 error=true;
420 std::cerr<<errstr.str()<<std::flush;
421 }
422#endif
423 std::sort(AllReverseRecv.begin(), AllReverseRecv.end(), Tpetra::Import_Util::sort_PID_then_GID<GlobalOrdinal, GlobalOrdinal>);
424
425 auto newEndOfPairs = std::unique(AllReverseRecv.begin(), AllReverseRecv.end());
426// don't resize to remove non-unique, just use the end-of-unique iterator
427 if(AllReverseRecv.begin() == newEndOfPairs) return;
428 int ARRsize = std::distance(AllReverseRecv.begin(),newEndOfPairs);
429 auto rPIDs = Teuchos::arcp(new int[ARRsize],0,ARRsize,true);
430 auto rLIDs = Teuchos::arcp(new LocalOrdinal[ARRsize],0,ARRsize,true);
431
432 int tsize=0;
433 for(auto itr = AllReverseRecv.begin(); itr!=newEndOfPairs; ++itr ) {
434 if((int)(itr->first) != MyPID) {
435 rPIDs[tsize]=(int)itr->first;
436 LocalOrdinal lid = MyDomainMap->getLocalElement(itr->second);
437 rLIDs[tsize]=lid;
438 tsize++;
439 }
440 }
441
442 type3PIDs=rPIDs.persistingView(0,tsize);
443 type3LIDs=rLIDs.persistingView(0,tsize);
444
445 if(error){
446 std::cerr<<errstr.str()<<std::flush;
447 comm->barrier();
448 comm->barrier();
449 comm->barrier();
450 MPI_Abort (MPI_COMM_WORLD, -1);
451 }
452#endif
453}
454
455// Note: This should get merged with the other Tpetra sort routines eventually.
456template<typename Scalar, typename Ordinal>
457void
458sortCrsEntries (const Teuchos::ArrayView<size_t> &CRS_rowptr,
459 const Teuchos::ArrayView<Ordinal> & CRS_colind,
460 const Teuchos::ArrayView<Scalar> &CRS_vals)
461{
462 // For each row, sort column entries from smallest to largest.
463 // Use shell sort. Stable sort so it is fast if indices are already sorted.
464 // Code copied from Epetra_CrsMatrix::SortEntries()
465 size_t NumRows = CRS_rowptr.size()-1;
466 size_t nnz = CRS_colind.size();
467
468 const bool permute_values_array = CRS_vals.size() > 0;
469
470 for(size_t i = 0; i < NumRows; i++){
471 size_t start=CRS_rowptr[i];
472 if(start >= nnz) continue;
473
474 size_t NumEntries = CRS_rowptr[i+1] - start;
475 Teuchos::ArrayRCP<Scalar> locValues;
476 if (permute_values_array)
477 locValues = Teuchos::arcp<Scalar>(&CRS_vals[start], 0, NumEntries, false);
478 Teuchos::ArrayRCP<Ordinal> locIndices(&CRS_colind[start], 0, NumEntries, false);
479
480 Ordinal n = NumEntries;
481 Ordinal m = 1;
482 while (m<n) m = m*3+1;
483 m /= 3;
484
485 while(m > 0) {
486 Ordinal max = n - m;
487 for(Ordinal j = 0; j < max; j++) {
488 for(Ordinal k = j; k >= 0; k-=m) {
489 if(locIndices[k+m] >= locIndices[k])
490 break;
491 if (permute_values_array) {
492 Scalar dtemp = locValues[k+m];
493 locValues[k+m] = locValues[k];
494 locValues[k] = dtemp;
495 }
496 Ordinal itemp = locIndices[k+m];
497 locIndices[k+m] = locIndices[k];
498 locIndices[k] = itemp;
499 }
500 }
501 m = m/3;
502 }
503 }
504}
505
506template<typename Ordinal>
507void
508sortCrsEntries (const Teuchos::ArrayView<size_t> &CRS_rowptr,
509 const Teuchos::ArrayView<Ordinal> & CRS_colind)
510{
511 // Generate dummy values array
512 Teuchos::ArrayView<Tpetra::Details::DefaultTypes::scalar_type> CRS_vals;
513 sortCrsEntries (CRS_rowptr, CRS_colind, CRS_vals);
514}
515
516namespace Impl {
517
518template<class RowOffsetsType, class ColumnIndicesType, class ValuesType>
519class SortCrsEntries {
520private:
521 typedef typename ColumnIndicesType::non_const_value_type ordinal_type;
522 typedef typename ValuesType::non_const_value_type scalar_type;
523
524public:
525 SortCrsEntries (const RowOffsetsType& ptr,
526 const ColumnIndicesType& ind,
527 const ValuesType& val) :
528 ptr_ (ptr),
529 ind_ (ind),
530 val_ (val)
531 {
532 static_assert (std::is_signed<ordinal_type>::value, "The type of each "
533 "column index -- that is, the type of each entry of ind "
534 "-- must be signed in order for this functor to work.");
535 }
536
537 KOKKOS_FUNCTION void operator() (const size_t i) const
538 {
539 const size_t nnz = ind_.extent (0);
540 const size_t start = ptr_(i);
541 const bool permute_values_array = val_.extent(0) > 0;
542
543 if (start < nnz) {
544 const size_t NumEntries = ptr_(i+1) - start;
545
546 const ordinal_type n = static_cast<ordinal_type> (NumEntries);
547 ordinal_type m = 1;
548 while (m<n) m = m*3+1;
549 m /= 3;
550
551 while (m > 0) {
552 ordinal_type max = n - m;
553 for (ordinal_type j = 0; j < max; j++) {
554 for (ordinal_type k = j; k >= 0; k -= m) {
555 const size_t sk = start+k;
556 if (ind_(sk+m) >= ind_(sk)) {
557 break;
558 }
559 if (permute_values_array) {
560 const scalar_type dtemp = val_(sk+m);
561 val_(sk+m) = val_(sk);
562 val_(sk) = dtemp;
563 }
564 const ordinal_type itemp = ind_(sk+m);
565 ind_(sk+m) = ind_(sk);
566 ind_(sk) = itemp;
567 }
568 }
569 m = m/3;
570 }
571 }
572 }
573
574 static void
575 sortCrsEntries (const RowOffsetsType& ptr,
576 const ColumnIndicesType& ind,
577 const ValuesType& val)
578 {
579 // For each row, sort column entries from smallest to largest.
580 // Use shell sort. Stable sort so it is fast if indices are already sorted.
581 // Code copied from Epetra_CrsMatrix::SortEntries()
582 // NOTE: This should not be taken as a particularly efficient way to sort
583 // rows of matrices in parallel. But it is correct, so that's something.
584 if (ptr.extent (0) == 0) {
585 return; // no rows, so nothing to sort
586 }
587 const size_t NumRows = ptr.extent (0) - 1;
588
589 typedef size_t index_type; // what this function was using; not my choice
590 typedef typename ValuesType::execution_space execution_space;
591 // Specify RangePolicy explicitly, in order to use correct execution
592 // space. See #1345.
593 typedef Kokkos::RangePolicy<execution_space, index_type> range_type;
594 Kokkos::parallel_for ("sortCrsEntries", range_type (0, NumRows),
595 SortCrsEntries (ptr, ind, val));
596 }
597
598private:
599 RowOffsetsType ptr_;
600 ColumnIndicesType ind_;
601 ValuesType val_;
602};
603
604} // namespace Impl
605
606template<typename rowptr_array_type, typename colind_array_type, typename vals_array_type>
607void
608sortCrsEntries (const rowptr_array_type& CRS_rowptr,
609 const colind_array_type& CRS_colind,
610 const vals_array_type& CRS_vals)
611{
612 Impl::SortCrsEntries<rowptr_array_type, colind_array_type,
613 vals_array_type>::sortCrsEntries (CRS_rowptr, CRS_colind, CRS_vals);
614}
615
616template<typename rowptr_array_type, typename colind_array_type>
617void
618sortCrsEntries (const rowptr_array_type& CRS_rowptr,
619 const colind_array_type& CRS_colind)
620{
621 // Generate dummy values array
622 typedef typename colind_array_type::execution_space execution_space;
623 typedef Tpetra::Details::DefaultTypes::scalar_type scalar_type;
624 typedef typename Kokkos::View<scalar_type*, execution_space> scalar_view_type;
625 scalar_view_type CRS_vals;
626 sortCrsEntries<rowptr_array_type, colind_array_type,
627 scalar_view_type>(CRS_rowptr, CRS_colind, CRS_vals);
628}
629
630// Note: This should get merged with the other Tpetra sort routines eventually.
631template<typename Scalar, typename Ordinal>
632void
633sortAndMergeCrsEntries (const Teuchos::ArrayView<size_t> &CRS_rowptr,
634 const Teuchos::ArrayView<Ordinal> & CRS_colind,
635 const Teuchos::ArrayView<Scalar> &CRS_vals)
636{
637 // For each row, sort column entries from smallest to largest,
638 // merging column ids that are identify by adding values. Use shell
639 // sort. Stable sort so it is fast if indices are already sorted.
640 // Code copied from Epetra_CrsMatrix::SortEntries()
641
642 if (CRS_rowptr.size () == 0) {
643 return; // no rows, so nothing to sort
644 }
645 const size_t NumRows = CRS_rowptr.size () - 1;
646 const size_t nnz = CRS_colind.size ();
647 size_t new_curr = CRS_rowptr[0];
648 size_t old_curr = CRS_rowptr[0];
649
650 const bool permute_values_array = CRS_vals.size() > 0;
651
652 for(size_t i = 0; i < NumRows; i++){
653 const size_t old_rowptr_i=CRS_rowptr[i];
654 CRS_rowptr[i] = old_curr;
655 if(old_rowptr_i >= nnz) continue;
656
657 size_t NumEntries = CRS_rowptr[i+1] - old_rowptr_i;
658 Teuchos::ArrayRCP<Scalar> locValues;
659 if (permute_values_array)
660 locValues = Teuchos::arcp<Scalar>(&CRS_vals[old_rowptr_i], 0, NumEntries, false);
661 Teuchos::ArrayRCP<Ordinal> locIndices(&CRS_colind[old_rowptr_i], 0, NumEntries, false);
662
663 // Sort phase
664 Ordinal n = NumEntries;
665 Ordinal m = n/2;
666
667 while(m > 0) {
668 Ordinal max = n - m;
669 for(Ordinal j = 0; j < max; j++) {
670 for(Ordinal k = j; k >= 0; k-=m) {
671 if(locIndices[k+m] >= locIndices[k])
672 break;
673 if (permute_values_array) {
674 Scalar dtemp = locValues[k+m];
675 locValues[k+m] = locValues[k];
676 locValues[k] = dtemp;
677 }
678 Ordinal itemp = locIndices[k+m];
679 locIndices[k+m] = locIndices[k];
680 locIndices[k] = itemp;
681 }
682 }
683 m = m/2;
684 }
685
686 // Merge & shrink
687 for(size_t j=old_rowptr_i; j < CRS_rowptr[i+1]; j++) {
688 if(j > old_rowptr_i && CRS_colind[j]==CRS_colind[new_curr-1]) {
689 if (permute_values_array) CRS_vals[new_curr-1] += CRS_vals[j];
690 }
691 else if(new_curr==j) {
692 new_curr++;
693 }
694 else {
695 CRS_colind[new_curr] = CRS_colind[j];
696 if (permute_values_array) CRS_vals[new_curr] = CRS_vals[j];
697 new_curr++;
698 }
699 }
700 old_curr=new_curr;
701 }
702
703 CRS_rowptr[NumRows] = new_curr;
704}
705
706template<typename Ordinal>
707void
708sortAndMergeCrsEntries (const Teuchos::ArrayView<size_t> &CRS_rowptr,
709 const Teuchos::ArrayView<Ordinal> & CRS_colind)
710{
711 Teuchos::ArrayView<Tpetra::Details::DefaultTypes::scalar_type> CRS_vals;
712 return sortAndMergeCrsEntries(CRS_rowptr, CRS_colind, CRS_vals);
713}
714
715
716template <typename LocalOrdinal, typename GlobalOrdinal, typename Node>
717void
718lowCommunicationMakeColMapAndReindex (const Teuchos::ArrayView<const size_t> &rowptr,
719 const Teuchos::ArrayView<LocalOrdinal> &colind_LID,
720 const Teuchos::ArrayView<GlobalOrdinal> &colind_GID,
721 const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMapRCP,
722 const Teuchos::ArrayView<const int> &owningPIDs,
723 Teuchos::Array<int> &remotePIDs,
724 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > & colMap)
725{
726 using Teuchos::rcp;
727 typedef LocalOrdinal LO;
728 typedef GlobalOrdinal GO;
729 typedef Tpetra::global_size_t GST;
730 typedef Tpetra::Map<LO, GO, Node> map_type;
731 const char prefix[] = "lowCommunicationMakeColMapAndReindex: ";
732
733 // The domainMap is an RCP because there is a shortcut for a
734 // (common) special case to return the columnMap = domainMap.
735 const map_type& domainMap = *domainMapRCP;
736
737 // Scan all column indices and sort into two groups:
738 // Local: those whose GID matches a GID of the domain map on this processor and
739 // Remote: All others.
740 const size_t numDomainElements = domainMap.getLocalNumElements ();
741 Teuchos::Array<bool> LocalGIDs;
742 if (numDomainElements > 0) {
743 LocalGIDs.resize (numDomainElements, false); // Assume domain GIDs are not local
744 }
745
746 // In principle it is good to have RemoteGIDs and RemotGIDList be as
747 // long as the number of remote GIDs on this processor, but this
748 // would require two passes through the column IDs, so we make it
749 // the max of 100 and the number of block rows.
750 //
751 // FIXME (mfh 11 Feb 2015) Tpetra::Details::HashTable can hold at
752 // most INT_MAX entries, but it's possible to have more rows than
753 // that (if size_t is 64 bits and int is 32 bits).
754 const size_t numMyRows = rowptr.size () - 1;
755 const int hashsize = std::max (static_cast<int> (numMyRows), 100);
756
757 Tpetra::Details::HashTable<GO, LO> RemoteGIDs (hashsize);
758 Teuchos::Array<GO> RemoteGIDList;
759 RemoteGIDList.reserve (hashsize);
760 Teuchos::Array<int> PIDList;
761 PIDList.reserve (hashsize);
762
763 // Here we start using the *LocalOrdinal* colind_LID array. This is
764 // safe even if both columnIndices arrays are actually the same
765 // (because LocalOrdinal==GO). For *local* GID's set
766 // colind_LID with with their LID in the domainMap. For *remote*
767 // GIDs, we set colind_LID with (numDomainElements+NumRemoteColGIDs)
768 // before the increment of the remote count. These numberings will
769 // be separate because no local LID is greater than
770 // numDomainElements.
771
772 size_t NumLocalColGIDs = 0;
773 LO NumRemoteColGIDs = 0;
774 for (size_t i = 0; i < numMyRows; ++i) {
775 for(size_t j = rowptr[i]; j < rowptr[i+1]; ++j) {
776 const GO GID = colind_GID[j];
777 // Check if GID matches a row GID
778 const LO LID = domainMap.getLocalElement (GID);
779 if(LID != -1) {
780 const bool alreadyFound = LocalGIDs[LID];
781 if (! alreadyFound) {
782 LocalGIDs[LID] = true; // There is a column in the graph associated with this domain map GID
783 NumLocalColGIDs++;
784 }
785 colind_LID[j] = LID;
786 }
787 else {
788 const LO hash_value = RemoteGIDs.get (GID);
789 if (hash_value == -1) { // This means its a new remote GID
790 const int PID = owningPIDs[j];
791 TEUCHOS_TEST_FOR_EXCEPTION(
792 PID == -1, std::invalid_argument, prefix << "Cannot figure out if "
793 "PID is owned.");
794 colind_LID[j] = static_cast<LO> (numDomainElements + NumRemoteColGIDs);
795 RemoteGIDs.add (GID, NumRemoteColGIDs);
796 RemoteGIDList.push_back (GID);
797 PIDList.push_back (PID);
798 NumRemoteColGIDs++;
799 }
800 else {
801 colind_LID[j] = static_cast<LO> (numDomainElements + hash_value);
802 }
803 }
804 }
805 }
806
807 // Possible short-circuit: If all domain map GIDs are present as
808 // column indices, then set ColMap=domainMap and quit.
809 if (domainMap.getComm ()->getSize () == 1) {
810 // Sanity check: When there is only one process, there can be no
811 // remoteGIDs.
812 TEUCHOS_TEST_FOR_EXCEPTION(
813 NumRemoteColGIDs != 0, std::runtime_error, prefix << "There is only one "
814 "process in the domain Map's communicator, which means that there are no "
815 "\"remote\" indices. Nevertheless, some column indices are not in the "
816 "domain Map.");
817 if (static_cast<size_t> (NumLocalColGIDs) == numDomainElements) {
818 // In this case, we just use the domainMap's indices, which is,
819 // not coincidently, what we clobbered colind with up above
820 // anyway. No further reindexing is needed.
821 colMap = domainMapRCP;
822 return;
823 }
824 }
825
826 // Now build the array containing column GIDs
827 // Build back end, containing remote GIDs, first
828 const LO numMyCols = NumLocalColGIDs + NumRemoteColGIDs;
829 Teuchos::Array<GO> ColIndices;
830 GO* RemoteColIndices = NULL;
831 if (numMyCols > 0) {
832 ColIndices.resize (numMyCols);
833 if (NumLocalColGIDs != static_cast<size_t> (numMyCols)) {
834 RemoteColIndices = &ColIndices[NumLocalColGIDs]; // Points to back half of ColIndices
835 }
836 }
837
838 for (LO i = 0; i < NumRemoteColGIDs; ++i) {
839 RemoteColIndices[i] = RemoteGIDList[i];
840 }
841
842 // Build permute array for *remote* reindexing.
843 Teuchos::Array<LO> RemotePermuteIDs (NumRemoteColGIDs);
844 for (LO i = 0; i < NumRemoteColGIDs; ++i) {
845 RemotePermuteIDs[i]=i;
846 }
847
848 // Sort External column indices so that all columns coming from a
849 // given remote processor are contiguous. This is a sort with two
850 // auxilary arrays: RemoteColIndices and RemotePermuteIDs.
851 Tpetra::sort3 (PIDList.begin (), PIDList.end (),
852 ColIndices.begin () + NumLocalColGIDs,
853 RemotePermuteIDs.begin ());
854
855 // Stash the RemotePIDs.
856 //
857 // Note: If Teuchos::Array had a shrink_to_fit like std::vector,
858 // we'd call it here.
859 remotePIDs = PIDList;
860
861 // Sort external column indices so that columns from a given remote
862 // processor are not only contiguous but also in ascending
863 // order. NOTE: I don't know if the number of externals associated
864 // with a given remote processor is known at this point ... so I
865 // count them here.
866
867 // NTS: Only sort the RemoteColIndices this time...
868 LO StartCurrent = 0, StartNext = 1;
869 while (StartNext < NumRemoteColGIDs) {
870 if (PIDList[StartNext]==PIDList[StartNext-1]) {
871 StartNext++;
872 }
873 else {
874 Tpetra::sort2 (ColIndices.begin () + NumLocalColGIDs + StartCurrent,
875 ColIndices.begin () + NumLocalColGIDs + StartNext,
876 RemotePermuteIDs.begin () + StartCurrent);
877 StartCurrent = StartNext;
878 StartNext++;
879 }
880 }
881 Tpetra::sort2 (ColIndices.begin () + NumLocalColGIDs + StartCurrent,
882 ColIndices.begin () + NumLocalColGIDs + StartNext,
883 RemotePermuteIDs.begin () + StartCurrent);
884
885 // Reverse the permutation to get the information we actually care about
886 Teuchos::Array<LO> ReverseRemotePermuteIDs (NumRemoteColGIDs);
887 for (LO i = 0; i < NumRemoteColGIDs; ++i) {
888 ReverseRemotePermuteIDs[RemotePermuteIDs[i]] = i;
889 }
890
891 // Build permute array for *local* reindexing.
892 bool use_local_permute = false;
893 Teuchos::Array<LO> LocalPermuteIDs (numDomainElements);
894
895 // Now fill front end. Two cases:
896 //
897 // (1) If the number of Local column GIDs is the same as the number
898 // of Local domain GIDs, we can simply read the domain GIDs into
899 // the front part of ColIndices, otherwise
900 //
901 // (2) We step through the GIDs of the domainMap, checking to see if
902 // each domain GID is a column GID. we want to do this to
903 // maintain a consistent ordering of GIDs between the columns
904 // and the domain.
905 Teuchos::ArrayView<const GO> domainGlobalElements = domainMap.getLocalElementList();
906 if (static_cast<size_t> (NumLocalColGIDs) == numDomainElements) {
907 if (NumLocalColGIDs > 0) {
908 // Load Global Indices into first numMyCols elements column GID list
909 std::copy (domainGlobalElements.begin (), domainGlobalElements.end (),
910 ColIndices.begin ());
911 }
912 }
913 else {
914 LO NumLocalAgain = 0;
915 use_local_permute = true;
916 for (size_t i = 0; i < numDomainElements; ++i) {
917 if (LocalGIDs[i]) {
918 LocalPermuteIDs[i] = NumLocalAgain;
919 ColIndices[NumLocalAgain++] = domainGlobalElements[i];
920 }
921 }
922 TEUCHOS_TEST_FOR_EXCEPTION(
923 static_cast<size_t> (NumLocalAgain) != NumLocalColGIDs,
924 std::runtime_error, prefix << "Local ID count test failed.");
925 }
926
927 // Make column Map
928 const GST minus_one = Teuchos::OrdinalTraits<GST>::invalid ();
929 colMap = rcp (new map_type (minus_one, ColIndices, domainMap.getIndexBase (),
930 domainMap.getComm ()));
931
932 // Low-cost reindex of the matrix
933 for (size_t i = 0; i < numMyRows; ++i) {
934 for (size_t j = rowptr[i]; j < rowptr[i+1]; ++j) {
935 const LO ID = colind_LID[j];
936 if (static_cast<size_t> (ID) < numDomainElements) {
937 if (use_local_permute) {
938 colind_LID[j] = LocalPermuteIDs[colind_LID[j]];
939 }
940 // In the case where use_local_permute==false, we just copy
941 // the DomainMap's ordering, which it so happens is what we
942 // put in colind_LID to begin with.
943 }
944 else {
945 colind_LID[j] = NumLocalColGIDs + ReverseRemotePermuteIDs[colind_LID[j]-numDomainElements];
946 }
947 }
948 }
949}
950
951
952
953
954// Generates an list of owning PIDs based on two transfer (aka import/export objects)
955// Let:
956// OwningMap = useReverseModeForOwnership ? transferThatDefinesOwnership.getTargetMap() : transferThatDefinesOwnership.getSourceMap();
957// MapAo = useReverseModeForOwnership ? transferThatDefinesOwnership.getSourceMap() : transferThatDefinesOwnership.getTargetMap();
958// MapAm = useReverseModeForMigration ? transferThatDefinesMigration.getTargetMap() : transferThatDefinesMigration.getSourceMap();
959// VectorMap = useReverseModeForMigration ? transferThatDefinesMigration.getSourceMap() : transferThatDefinesMigration.getTargetMap();
960// Precondition:
961// 1) MapAo.isSameAs(*MapAm) - map compatibility between transfers
962// 2) VectorMap->isSameAs(*owningPIDs->getMap()) - map compabibility between transfer & vector
963// 3) OwningMap->isOneToOne() - owning map is 1-to-1
964// --- Precondition 3 is only checked in DEBUG mode ---
965// Postcondition:
966// owningPIDs[VectorMap->getLocalElement(GID i)] = j iff (OwningMap->isLocalElement(GID i) on rank j)
967template <typename LocalOrdinal, typename GlobalOrdinal, typename Node>
968void getTwoTransferOwnershipVector(const ::Tpetra::Details::Transfer<LocalOrdinal, GlobalOrdinal, Node>& transferThatDefinesOwnership,
969 bool useReverseModeForOwnership,
970 const ::Tpetra::Details::Transfer<LocalOrdinal, GlobalOrdinal, Node>& transferThatDefinesMigration,
971 bool useReverseModeForMigration,
975
976 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > OwningMap = useReverseModeForOwnership ?
977 transferThatDefinesOwnership.getTargetMap() :
978 transferThatDefinesOwnership.getSourceMap();
979 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > MapAo = useReverseModeForOwnership ?
980 transferThatDefinesOwnership.getSourceMap() :
981 transferThatDefinesOwnership.getTargetMap();
982 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > MapAm = useReverseModeForMigration ?
983 transferThatDefinesMigration.getTargetMap() :
984 transferThatDefinesMigration.getSourceMap();
985 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > VectorMap = useReverseModeForMigration ?
986 transferThatDefinesMigration.getSourceMap() :
987 transferThatDefinesMigration.getTargetMap();
988
989 TEUCHOS_TEST_FOR_EXCEPTION(!MapAo->isSameAs(*MapAm),std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector map mismatch between transfers");
990 TEUCHOS_TEST_FOR_EXCEPTION(!VectorMap->isSameAs(*owningPIDs.getMap()),std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector map mismatch transfer and vector");
991#ifdef HAVE_TPETRA_DEBUG
992 TEUCHOS_TEST_FOR_EXCEPTION(!OwningMap->isOneToOne(),std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector owner must be 1-to-1");
993#endif
994
995 int rank = OwningMap->getComm()->getRank();
996 // Generate "A" vector and fill it with owning information. We can read this from transferThatDefinesOwnership w/o communication
997 // Note: Due to the 1-to-1 requirement, several of these options throw
999 const import_type* ownAsImport = dynamic_cast<const import_type*> (&transferThatDefinesOwnership);
1000 const export_type* ownAsExport = dynamic_cast<const export_type*> (&transferThatDefinesOwnership);
1001
1002 Teuchos::ArrayRCP<int> pids = temp.getDataNonConst();
1003 Teuchos::ArrayView<int> v_pids = pids();
1004 if(ownAsImport && useReverseModeForOwnership) {TEUCHOS_TEST_FOR_EXCEPTION(1,std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector owner must be 1-to-1");}
1005 else if(ownAsImport && !useReverseModeForOwnership) getPids(*ownAsImport,v_pids,false);
1006 else if(ownAsExport && useReverseModeForMigration) {TEUCHOS_TEST_FOR_EXCEPTION(1,std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector this option not yet implemented");}
1007 else {TEUCHOS_TEST_FOR_EXCEPTION(1,std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector owner must be 1-to-1");}
1008
1009 const import_type* xferAsImport = dynamic_cast<const import_type*> (&transferThatDefinesMigration);
1010 const export_type* xferAsExport = dynamic_cast<const export_type*> (&transferThatDefinesMigration);
1011 TEUCHOS_TEST_FOR_EXCEPTION(!xferAsImport && !xferAsExport,std::runtime_error,"Tpetra::Import_Util::getTwoTransferOwnershipVector transfer undefined");
1012
1013 // Migrate from "A" vector to output vector
1014 owningPIDs.putScalar(rank);
1015 if(xferAsImport && useReverseModeForMigration) owningPIDs.doExport(temp,*xferAsImport,Tpetra::REPLACE);
1016 else if(xferAsImport && !useReverseModeForMigration) owningPIDs.doImport(temp,*xferAsImport,Tpetra::REPLACE);
1017 else if(xferAsExport && useReverseModeForMigration) owningPIDs.doImport(temp,*xferAsExport,Tpetra::REPLACE);
1018 else owningPIDs.doExport(temp,*xferAsExport,Tpetra::REPLACE);
1019
1020}
1021
1022
1023
1024} // namespace Import_Util
1025} // namespace Tpetra
1026
1027#endif // TPETRA_IMPORT_UTIL_HPP
Declaration of the Tpetra::CrsMatrix class.
Add specializations of Teuchos::Details::MpiTypeTraits for Kokkos::complex<float> and Kokkos::complex...
Declaration and definition of Tpetra::Details::reallocDualViewIfNeeded, an implementation detail of T...
void lowCommunicationMakeColMapAndReindex(const Teuchos::ArrayView< const size_t > &rowPointers, const Teuchos::ArrayView< LocalOrdinal > &columnIndices_LID, const Teuchos::ArrayView< GlobalOrdinal > &columnIndices_GID, const Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::ArrayView< const int > &owningPids, Teuchos::Array< int > &remotePids, Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap)
lowCommunicationMakeColMapAndReindex
void getTwoTransferOwnershipVector(const ::Tpetra::Details::Transfer< LocalOrdinal, GlobalOrdinal, Node > &transferThatDefinesOwnership, bool useReverseModeForOwnership, const ::Tpetra::Details::Transfer< LocalOrdinal, GlobalOrdinal, Node > &transferForMigratingData, bool useReverseModeForMigration, Tpetra::Vector< int, LocalOrdinal, GlobalOrdinal, Node > &owningPIDs)
Generates an list of owning PIDs based on two transfer (aka import/export objects) Let: OwningMap = u...
void sortAndMergeCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< Ordinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals)
Sort and merge the entries of the (raw CSR) matrix by column index within each row.
void sortCrsEntries(const Teuchos::ArrayView< size_t > &CRS_rowptr, const Teuchos::ArrayView< Ordinal > &CRS_colind, const Teuchos::ArrayView< Scalar > &CRS_vals)
Sort the entries of the (raw CSR) matrix by column index within each row.
void getPids(const Tpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &Importer, Teuchos::Array< int > &pids, bool use_minus_one_for_local)
Like getPidGidPairs, but just gets the PIDs, ordered by the column Map.
Stand-alone utility functions and macros.
ValueType get(const KeyType key)
Get the value corresponding to the given key.
void add(const KeyType key, const ValueType value)
Add a key and its value to the hash table.
Common base class of Import and Export.
Teuchos::RCP< const map_type > getTargetMap() const
The target Map used to construct this Export or Import.
Teuchos::ArrayView< const LO > getExportLIDs() const
List of entries in the source Map that will be sent to other processes.
size_t getNumExportIDs() const
Number of entries that must be sent by the calling process to other processes.
Teuchos::ArrayView< const int > getExportPIDs() const
List of processes to which entries will be sent.
void doImport(const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false)
Import data into this object using an Import object ("forward mode").
void doExport(const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false)
Export data into this object using an Export object ("forward mode").
virtual Teuchos::RCP< const map_type > getMap() const
The Map describing the parallel distribution of this object.
Communication plan for data redistribution from a (possibly) multiply-owned to a uniquely-owned distr...
Communication plan for data redistribution from a uniquely-owned to a (possibly) multiply-owned distr...
A parallel distribution of indices over processes.
size_t getLocalNumElements() const
The number of elements belonging to the calling process.
void putScalar(const Scalar &value)
Set all values in the multivector with the given value.
A distributed dense vector.
Teuchos::ArrayRCP< Scalar > getDataNonConst()
View of the local values of this vector.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void sort3(const IT1 &first1, const IT1 &last1, const IT2 &first2, const IT3 &first3)
Sort the first array, and apply the same permutation to the second and third arrays.
void sort2(const IT1 &first1, const IT1 &last1, const IT2 &first2)
Sort the first array, and apply the resulting permutation to the second array.
size_t global_size_t
Global size_t object.
@ REPLACE
Replace existing values with new values.