Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_Import.cpp
Go to the documentation of this file.
1/*
2//@HEADER
3// ************************************************************************
4//
5// Epetra: Linear Algebra Services Package
6// Copyright 2011 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
39//
40// ************************************************************************
41//@HEADER
42*/
43
44#include "Epetra_ConfigDefs.h"
45#include "Epetra_Import.h"
46#include "Epetra_BlockMap.h"
47#include "Epetra_Distributor.h"
48#include "Epetra_Comm.h"
49#include "Epetra_Util.h"
50#include "Epetra_Export.h"
51
52#include <algorithm>
53#include <vector>
54
55
56#ifdef HAVE_MPI
58#endif
59
60#ifdef EPETRA_ENABLE_DEBUG
61#include "Epetra_IntVector.h"
62#endif
63
64//==============================================================================
65// Epetra_Import constructor function for a Epetra_BlockMap object
66template<typename int_type>
67void Epetra_Import::Construct_Expert( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap,
68 int NumRemotePIDs, const int * UserRemotePIDs,
69 const int & UserNumExportIDs, const int * UserExportLIDs, const int * UserExportPIDs)
70{
71 int i,ierr;
72 // Build three ID lists:
73 // NumSameIDs - Number of IDs in TargetMap and SourceMap that are identical, up to the first
74 // nonidentical ID.
75 // NumPermuteIDs - Number of IDs in SourceMap that must be indirectly loaded but are on this processor.
76 // NumRemoteIDs - Number of IDs that are in SourceMap but not in TargetMap, and thus must be imported.
77
78 int NumSourceIDs = sourceMap.NumMyElements();
79 int NumTargetIDs = targetMap.NumMyElements();
80
81 int_type *TargetGIDs = 0;
82 if (NumTargetIDs>0) {
83 TargetGIDs = new int_type[NumTargetIDs];
84 targetMap.MyGlobalElements(TargetGIDs);
85 }
86
87 int_type * SourceGIDs = 0;
88 if (NumSourceIDs>0) {
89 SourceGIDs = new int_type[NumSourceIDs];
90 sourceMap.MyGlobalElements(SourceGIDs);
91 }
92
93 int MinIDs = EPETRA_MIN(NumSourceIDs, NumTargetIDs);
94
95 NumSameIDs_ = 0;
96 for (i=0; i< MinIDs; i++) if (TargetGIDs[i]==SourceGIDs[i]) NumSameIDs_++; else break;
97
98 // Find count of Target IDs that are truly remote and those that are local but permuted
100 NumRemoteIDs_ = 0;
101 for (i=NumSameIDs_; i< NumTargetIDs; i++)
102 if (sourceMap.MyGID(TargetGIDs[i])) NumPermuteIDs_++; // Check if Target GID is a local Source GID
103 else NumRemoteIDs_++; // If not, then it is remote
104
105
106
107 // Define remote and permutation lists
108 int_type * RemoteGIDs=0;
109 RemoteLIDs_ = 0;
110 if (NumRemoteIDs_>0) {
111 RemoteLIDs_ = new int[NumRemoteIDs_];
112 RemoteGIDs = new int_type[NumRemoteIDs_];
113 }
114 if (NumPermuteIDs_>0) {
117 }
118
119 NumPermuteIDs_ = 0;
120 NumRemoteIDs_ = 0;
121 for (i=NumSameIDs_; i< NumTargetIDs; i++) {
122 if (sourceMap.MyGID(TargetGIDs[i])) {
124 PermuteFromLIDs_[NumPermuteIDs_++] = sourceMap.LID(TargetGIDs[i]);
125 }
126 else {
127 //NumRecv_ +=TargetMap.ElementSize(i); // Count total number of entries to receive
128 NumRecv_ +=targetMap.MaxElementSize(); // Count total number of entries to receive (currently need max)
129 RemoteGIDs[NumRemoteIDs_] = TargetGIDs[i];
131 }
132 }
133
134 if( NumRemoteIDs_>0 && !sourceMap.DistributedGlobal() )
135 ReportError("Warning in Epetra_Import: Serial Import has remote IDs. (Importing to Subset of Target Map)", 1);
136
137 // Test for distributed cases
138 int * RemotePIDs = 0;
139
140 if (sourceMap.DistributedGlobal()) {
141 if (NumRemoteIDs_>0) RemotePIDs = new int[NumRemoteIDs_];
142
143#ifdef EPETRA_ENABLE_DEBUG
144 int myeq = (NumRemotePIDs==NumRemoteIDs_);
145 int globaleq=0;
146 sourceMap.Comm().MinAll(&myeq,&globaleq,1);
147 if(globaleq!=1) {
148 printf("[%d] UserRemotePIDs count wrong %d != %d\n",sourceMap.Comm().MyPID(),NumRemotePIDs,NumRemoteIDs_);
149 fflush(stdout);
150 sourceMap.Comm().Barrier();
151 sourceMap.Comm().Barrier();
152 sourceMap.Comm().Barrier();
153 ReportError("Epetra_Import: UserRemotePIDs count wrong",-1);
154 }
155#endif
156
157 if(NumRemotePIDs==NumRemoteIDs_){
158 // Since I need to sort these, I'll copy them
159 for(i=0; i<NumRemoteIDs_; i++) RemotePIDs[i] = UserRemotePIDs[i];
160 } else {
161 // This can happen if the UserRemotePIDs comes from a Map that is not ordered properly
162 ReportError("Epetra_Import: UserRemotePIDs count wrong",-1);
163 }
164
165 //Get rid of IDs that don't exist in SourceMap
166 if(NumRemoteIDs_>0) {
167 int cnt = 0;
168 for( i = 0; i < NumRemoteIDs_; ++i )
169 if( RemotePIDs[i] == -1 ) ++cnt;
170 if( cnt ) {
171 if( NumRemoteIDs_-cnt ) {
172 int_type * NewRemoteGIDs = new int_type[NumRemoteIDs_-cnt];
173 int * NewRemotePIDs = new int[NumRemoteIDs_-cnt];
174 int * NewRemoteLIDs = new int[NumRemoteIDs_-cnt];
175 cnt = 0;
176 for( i = 0; i < NumRemoteIDs_; ++i )
177 if( RemotePIDs[i] != -1 ) {
178 NewRemoteGIDs[cnt] = RemoteGIDs[i];
179 NewRemotePIDs[cnt] = RemotePIDs[i];
180 NewRemoteLIDs[cnt] = targetMap.LID(RemoteGIDs[i]);
181 ++cnt;
182 }
183 NumRemoteIDs_ = cnt;
184 delete [] RemoteGIDs;
185 delete [] RemotePIDs;
186 delete [] RemoteLIDs_;
187 RemoteGIDs = NewRemoteGIDs;
188 RemotePIDs = NewRemotePIDs;
189 RemoteLIDs_ = NewRemoteLIDs;
190 ReportError("Warning in Epetra_Import: Target IDs not found in Source Map (Do you want to import to subset of Target Map?)", 1);
191 }
192 else { //valid RemoteIDs empty
193 NumRemoteIDs_ = 0;
194 delete [] RemoteGIDs;
195 RemoteGIDs = 0;
196 delete [] RemotePIDs;
197 RemotePIDs = 0;
198 }
199 }
200 }
201
202 //Sort Remote IDs by processor so DoReverses will work
203 Epetra_Util util;
204
205 if(targetMap.GlobalIndicesLongLong())
206 {
207 // FIXME (mfh 11 Jul 2013) This breaks ANSI aliasing rules, if
208 // int_type != long long. On some compilers, it results in
209 // warnings such as this: "dereferencing type-punned pointer
210 // will break strict-aliasing rules".
211 util.Sort(true,NumRemoteIDs_,RemotePIDs,0,0, 1,&RemoteLIDs_, 1,(long long**)&RemoteGIDs);
212 }
213 else if(targetMap.GlobalIndicesInt())
214 {
215 int* ptrs[2] = {RemoteLIDs_, (int*)RemoteGIDs};
216 util.Sort(true,NumRemoteIDs_,RemotePIDs,0,0,2,&ptrs[0], 0, 0);
217 }
218 else
219 {
220 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
221 }
222
223 // Build distributor & Export lists
224 Distor_ = sourceMap.Comm().CreateDistributor();
225
226 NumExportIDs_=UserNumExportIDs;
227 ExportLIDs_ = new int[NumExportIDs_];
228 ExportPIDs_ = new int[NumExportIDs_];
229 for(i=0; i<NumExportIDs_; i++) {
230 ExportPIDs_[i] = UserExportPIDs[i];
231 ExportLIDs_[i] = UserExportLIDs[i];
232 }
233
234 // Send Counts
235 for (i=0; i< NumExportIDs_; i++) {
236 NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
237 }
238
239#ifdef HAVE_MPI
240 Epetra_MpiDistributor* MpiDistor = dynamic_cast< Epetra_MpiDistributor*>(Distor_);
241 bool Deterministic = true;
242 if(MpiDistor)
244 NumRemoteIDs_, RemoteGIDs, RemotePIDs,Deterministic);
245 else ierr=-10;
246#else
247 ierr=-20;
248#endif
249
250 if (ierr!=0) throw ReportError("Error in Epetra_Distributor.CreateFromRecvs()", ierr);
251 }
252
253 if( NumRemoteIDs_>0 ) delete [] RemoteGIDs;
254 if( NumRemoteIDs_>0 ) delete [] RemotePIDs;
255
256 if (NumTargetIDs>0) delete [] TargetGIDs;
257 if (NumSourceIDs>0) delete [] SourceGIDs;
258
259
260#ifdef EPETRA_ENABLE_DEBUG
261// Sanity check to make sure we got the import right
262 Epetra_IntVector Source(sourceMap);
263 Epetra_IntVector Target(targetMap);
264
265 for(i=0; i<Source.MyLength(); i++)
266 Source[i] = (int) (Source.Map().GID(i) % INT_MAX);
267 Target.PutValue(-1);
268
269 Target.Import(Source,*this,Insert);
270
271 bool test_passed=true;
272 for(i=0; i<Target.MyLength(); i++){
273 if(Target[i] != Target.Map().GID(i) % INT_MAX) test_passed=false;
274 }
275
276 if(!test_passed) {
277 printf("[%d] PROCESSOR has a mismatch... prepearing to crash or hang!\n",sourceMap.Comm().MyPID());
278 fflush(stdout);
279 sourceMap.Comm().Barrier();
280 sourceMap.Comm().Barrier();
281 sourceMap.Comm().Barrier();
282 throw ReportError("Epetra_Import: ERROR. User provided IDs do not match what an import generates.");
283 }
284#endif
285
286 return;
287}
288
289//==============================================================================
290// Epetra_Import constructor function for a Epetra_BlockMap object
291template<typename int_type>
292void Epetra_Import::Construct( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap, int NumRemotePIDs, const int * UserRemotePIDs)
293{
294 int i,ierr;
295 // Build three ID lists:
296 // NumSameIDs - Number of IDs in TargetMap and SourceMap that are identical, up to the first
297 // nonidentical ID.
298 // NumPermuteIDs - Number of IDs in SourceMap that must be indirectly loaded but are on this processor.
299 // NumRemoteIDs - Number of IDs that are in SourceMap but not in TargetMap, and thus must be imported.
300
301 int NumSourceIDs = sourceMap.NumMyElements();
302 int NumTargetIDs = targetMap.NumMyElements();
303
304 int_type *TargetGIDs = 0;
305 if (NumTargetIDs>0) {
306 TargetGIDs = new int_type[NumTargetIDs];
307 targetMap.MyGlobalElements(TargetGIDs);
308 }
309
310 int_type * SourceGIDs = 0;
311 if (NumSourceIDs>0) {
312 SourceGIDs = new int_type[NumSourceIDs];
313 sourceMap.MyGlobalElements(SourceGIDs);
314 }
315
316 int MinIDs = EPETRA_MIN(NumSourceIDs, NumTargetIDs);
317
318
319 NumSameIDs_ = 0;
320 for (i=0; i< MinIDs; i++) if (TargetGIDs[i]==SourceGIDs[i]) NumSameIDs_++; else break;
321
322 // Find count of Target IDs that are truly remote and those that are local but permuted
323 NumPermuteIDs_ = 0;
324 NumRemoteIDs_ = 0;
325 for (i=NumSameIDs_; i< NumTargetIDs; i++)
326 if (sourceMap.MyGID(TargetGIDs[i])) NumPermuteIDs_++; // Check if Target GID is a local Source GID
327 else NumRemoteIDs_++; // If not, then it is remote
328
329
330
331 // Define remote and permutation lists
332 int_type * RemoteGIDs=0;
333 RemoteLIDs_ = 0;
334 if (NumRemoteIDs_>0) {
335 RemoteLIDs_ = new int[NumRemoteIDs_];
336 RemoteGIDs = new int_type[NumRemoteIDs_];
337 }
338 if (NumPermuteIDs_>0) {
341 }
342
343 NumPermuteIDs_ = 0;
344 NumRemoteIDs_ = 0;
345 for (i=NumSameIDs_; i< NumTargetIDs; i++) {
346 if (sourceMap.MyGID(TargetGIDs[i])) {
348 PermuteFromLIDs_[NumPermuteIDs_++] = sourceMap.LID(TargetGIDs[i]);
349 }
350 else {
351 //NumRecv_ +=TargetMap.ElementSize(i); // Count total number of entries to receive
352 NumRecv_ +=targetMap.MaxElementSize(); // Count total number of entries to receive (currently need max)
353 RemoteGIDs[NumRemoteIDs_] = TargetGIDs[i];
355 }
356 }
357
358 if( NumRemoteIDs_>0 && !sourceMap.DistributedGlobal() )
359 ReportError("Warning in Epetra_Import: Serial Import has remote IDs. (Importing to Subset of Target Map)", 1);
360
361 // Test for distributed cases
362 int * RemotePIDs = 0;
363
364 if (sourceMap.DistributedGlobal()) {
365 if (NumRemoteIDs_>0) RemotePIDs = new int[NumRemoteIDs_];
366
367#ifdef EPETRA_ENABLE_DEBUG
368 if(NumRemotePIDs!=-1){
369 int myeq = (NumRemotePIDs==NumRemoteIDs_);
370 int globaleq=0;
371 sourceMap.Comm().MinAll(&myeq,&globaleq,1);
372 if(globaleq!=1) {
373 printf("[%d] UserRemotePIDs count wrong %d != %d\n",sourceMap.Comm().MyPID(),NumRemotePIDs,NumRemoteIDs_);
374 fflush(stdout);
375 sourceMap.Comm().Barrier();
376 sourceMap.Comm().Barrier();
377 sourceMap.Comm().Barrier();
378 throw ReportError("Epetra_Import: UserRemotePIDs count wrong",-1);
379 }
380 }
381#endif
382
383 if(NumRemotePIDs!=-1 && NumRemotePIDs==NumRemoteIDs_){
384 // Since I need to sort these, I'll copy them
385 for(i=0; i<NumRemoteIDs_; i++) RemotePIDs[i] = UserRemotePIDs[i];
386 }
387 else{
388 ierr=sourceMap.RemoteIDList(NumRemoteIDs_, RemoteGIDs, RemotePIDs, 0); // Get remote PIDs
389 if (ierr) throw ReportError("Error in sourceMap.RemoteIDList call", ierr);
390 }
391
392 //Get rid of IDs that don't exist in SourceMap
393 if(NumRemoteIDs_>0) {
394 int cnt = 0;
395 for( i = 0; i < NumRemoteIDs_; ++i )
396 if( RemotePIDs[i] == -1 ) ++cnt;
397 if( cnt ) {
398 if( NumRemoteIDs_-cnt ) {
399 int_type * NewRemoteGIDs = new int_type[NumRemoteIDs_-cnt];
400 int * NewRemotePIDs = new int[NumRemoteIDs_-cnt];
401 int * NewRemoteLIDs = new int[NumRemoteIDs_-cnt];
402 cnt = 0;
403 for( i = 0; i < NumRemoteIDs_; ++i )
404 if( RemotePIDs[i] != -1 ) {
405 NewRemoteGIDs[cnt] = RemoteGIDs[i];
406 NewRemotePIDs[cnt] = RemotePIDs[i];
407 NewRemoteLIDs[cnt] = targetMap.LID(RemoteGIDs[i]);
408 ++cnt;
409 }
410 NumRemoteIDs_ = cnt;
411 delete [] RemoteGIDs;
412 delete [] RemotePIDs;
413 delete [] RemoteLIDs_;
414 RemoteGIDs = NewRemoteGIDs;
415 RemotePIDs = NewRemotePIDs;
416 RemoteLIDs_ = NewRemoteLIDs;
417 ReportError("Warning in Epetra_Import: Target IDs not found in Source Map (Do you want to import to subset of Target Map?)", 1);
418 }
419 else { //valid RemoteIDs empty
420 NumRemoteIDs_ = 0;
421 delete [] RemoteGIDs;
422 RemoteGIDs = 0;
423 delete [] RemotePIDs;
424 RemotePIDs = 0;
425 }
426 }
427 }
428
429 //Sort Remote IDs by processor so DoReverses will work
430 if(targetMap.GlobalIndicesLongLong())
431 {
432 // FIXME (mfh 11 Jul 2013) This breaks ANSI aliasing rules, if
433 // int_type != long long. On some compilers, it results in
434 // warnings such as this: "dereferencing type-punned pointer
435 // will break strict-aliasing rules".
436 Epetra_Util::Sort(true,NumRemoteIDs_,RemotePIDs,0,0, 1,&RemoteLIDs_, 1,(long long**)&RemoteGIDs);
437 }
438 else if(targetMap.GlobalIndicesInt())
439 {
440 int* ptrs[2] = {RemoteLIDs_, (int*)RemoteGIDs};
441 Epetra_Util::Sort(true,NumRemoteIDs_,RemotePIDs,0,0,2,&ptrs[0], 0, 0);
442 }
443 else
444 {
445 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
446 }
447 Distor_ = sourceMap.Comm().CreateDistributor();
448
449 // Construct list of exports that calling processor needs to send as a result
450 // of everyone asking for what it needs to receive.
451
452 bool Deterministic = true;
453 int_type* tmp_ExportLIDs; //Export IDs come in as GIDs
454 ierr = Distor_->CreateFromRecvs( NumRemoteIDs_, RemoteGIDs, RemotePIDs,
455 Deterministic, NumExportIDs_, tmp_ExportLIDs, ExportPIDs_ );
456 if (ierr!=0) throw ReportError("Error in Epetra_Distributor.CreateFromRecvs()", ierr);
457
458
459 // Export IDs come in as GIDs, convert to LIDs
460 if(targetMap.GlobalIndicesLongLong())
461 {
462 ExportLIDs_ = new int[NumExportIDs_];
463
464 for (i=0; i< NumExportIDs_; i++) {
465 if (ExportPIDs_[i] < 0) throw ReportError("targetMap requested a GID that is not in the sourceMap.", -1);
466 ExportLIDs_[i] = sourceMap.LID(tmp_ExportLIDs[i]);
467 NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
468 }
469
470 delete[] tmp_ExportLIDs;
471 }
472 else if(targetMap.GlobalIndicesInt())
473 {
474 for (i=0; i< NumExportIDs_; i++) {
475 if (ExportPIDs_[i] < 0) throw ReportError("targetMap requested a GID that is not in the sourceMap.", -1);
476 tmp_ExportLIDs[i] = sourceMap.LID(tmp_ExportLIDs[i]);
477 NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
478 }
479
480 ExportLIDs_ = reinterpret_cast<int *>(tmp_ExportLIDs); // Can't reach here if tmp_ExportLIDs is long long.
481 }
482 else
483 {
484 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
485 }
486 }
487
488 if( NumRemoteIDs_>0 ) delete [] RemoteGIDs;
489 if( NumRemoteIDs_>0 ) delete [] RemotePIDs;
490
491 if (NumTargetIDs>0) delete [] TargetGIDs;
492 if (NumSourceIDs>0) delete [] SourceGIDs;
493
494 return;
495}
496
497//==============================================================================
498Epetra_Import::Epetra_Import( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap,int NumRemotePIDs,const int * RemotePIDs)
499 : Epetra_Object("Epetra::Import"),
500 TargetMap_(targetMap),
501 SourceMap_(sourceMap),
502 NumSameIDs_(0),
503 NumPermuteIDs_(0),
504 PermuteToLIDs_(0),
505 PermuteFromLIDs_(0),
506 NumRemoteIDs_(0),
507 RemoteLIDs_(0),
508 NumExportIDs_(0),
509 ExportLIDs_(0),
510 ExportPIDs_(0),
511 NumSend_(0),
512 NumRecv_(0),
513 Distor_(0)
514{
515 if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
516 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
517
518 if(targetMap.GlobalIndicesInt())
519#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
520 Construct<int>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs);
521#else
522 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
523#endif
524 else if(targetMap.GlobalIndicesLongLong())
525#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
526 Construct<long long>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs);
527#else
528 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
529#endif
530 else
531 throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
532}
533
534
535
536//==============================================================================
537Epetra_Import::Epetra_Import( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap,int NumRemotePIDs,const int * RemotePIDs,
538 const int & numExportIDs, const int * theExportLIDs, const int * theExportPIDs)
539 : Epetra_Object("Epetra::Import"),
540 TargetMap_(targetMap),
541 SourceMap_(sourceMap),
542 NumSameIDs_(0),
543 NumPermuteIDs_(0),
544 PermuteToLIDs_(0),
545 PermuteFromLIDs_(0),
546 NumRemoteIDs_(0),
547 RemoteLIDs_(0),
548 NumExportIDs_(0),
549 ExportLIDs_(0),
550 ExportPIDs_(0),
551 NumSend_(0),
552 NumRecv_(0),
553 Distor_(0)
554{
555 if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
556 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
557
558 if(targetMap.GlobalIndicesInt())
559#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
560 Construct_Expert<int>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs,numExportIDs,theExportLIDs,theExportPIDs);
561#else
562 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
563#endif
564 else if(targetMap.GlobalIndicesLongLong())
565#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
566 Construct_Expert<long long>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs,numExportIDs,theExportLIDs,theExportPIDs);
567#else
568 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
569#endif
570 else
571 throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
572}
573
574
575//==============================================================================
577 : Epetra_Object("Epetra::Import"),
578 TargetMap_(targetMap),
579 SourceMap_(sourceMap),
580 NumSameIDs_(0),
581 NumPermuteIDs_(0),
582 PermuteToLIDs_(0),
583 PermuteFromLIDs_(0),
584 NumRemoteIDs_(0),
585 RemoteLIDs_(0),
586 NumExportIDs_(0),
587 ExportLIDs_(0),
588 ExportPIDs_(0),
589 NumSend_(0),
590 NumRecv_(0),
591 Distor_(0)
592{
593 if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
594 throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
595
596 if(targetMap.GlobalIndicesInt())
597#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
598 Construct<int>(targetMap, sourceMap);
599#else
600 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
601#endif
602 else if(targetMap.GlobalIndicesLongLong())
603#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
604 Construct<long long>(targetMap, sourceMap);
605#else
606 throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
607#endif
608 else
609 throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
610}
611
612//==============================================================================
613// Epetra_Import copy constructor
615 : Epetra_Object(Importer),
616 TargetMap_(Importer.TargetMap_),
617 SourceMap_(Importer.SourceMap_),
618 NumSameIDs_(Importer.NumSameIDs_),
619 NumPermuteIDs_(Importer.NumPermuteIDs_),
620 PermuteToLIDs_(0),
621 PermuteFromLIDs_(0),
622 NumRemoteIDs_(Importer.NumRemoteIDs_),
623 RemoteLIDs_(0),
624 NumExportIDs_(Importer.NumExportIDs_),
625 ExportLIDs_(0),
626 ExportPIDs_(0),
627 NumSend_(Importer.NumSend_),
628 NumRecv_(Importer.NumRecv_),
629 Distor_(0)
630{
631 int i;
632 if (NumPermuteIDs_>0) {
635 for (i=0; i< NumPermuteIDs_; i++) {
636 PermuteToLIDs_[i] = Importer.PermuteToLIDs_[i];
637 PermuteFromLIDs_[i] = Importer.PermuteFromLIDs_[i];
638 }
639 }
640
641 if (NumRemoteIDs_>0) {
642 RemoteLIDs_ = new int[NumRemoteIDs_];
643 for (i=0; i< NumRemoteIDs_; i++) RemoteLIDs_[i] = Importer.RemoteLIDs_[i];
644 }
645
646 if (NumExportIDs_>0) {
647 ExportLIDs_ = new int[NumExportIDs_];
648 ExportPIDs_ = new int[NumExportIDs_];
649 for (i=0; i< NumExportIDs_; i++) {
650 ExportLIDs_[i] = Importer.ExportLIDs_[i];
651 ExportPIDs_[i] = Importer.ExportPIDs_[i];
652 }
653 }
654
655 if (Importer.Distor_!=0) Distor_ = Importer.Distor_->Clone();
656
657}
658
659//==============================================================================
660// Epetra_Import destructor
662{
663 if( Distor_ != 0 ) delete Distor_;
664 if (RemoteLIDs_ != 0) delete [] RemoteLIDs_;
665 if (PermuteToLIDs_ != 0) delete [] PermuteToLIDs_;
666 if (PermuteFromLIDs_ != 0) delete [] PermuteFromLIDs_;
667
668 if( ExportPIDs_ != 0 ) delete [] ExportPIDs_; // These were created by Distor_
669 if( ExportLIDs_ != 0 ) delete [] ExportLIDs_;
670}
671
672//==============================================================================
673// Epetra_Import pseudo-copy constructor.
675 TargetMap_(Exporter.SourceMap_), //reverse
676 SourceMap_(Exporter.TargetMap_),//reverse
677 NumSameIDs_(Exporter.NumSameIDs_),
678 NumPermuteIDs_(Exporter.NumPermuteIDs_),
679 PermuteToLIDs_(0),
680 PermuteFromLIDs_(0),
681 NumRemoteIDs_(Exporter.NumExportIDs_),//reverse
682 RemoteLIDs_(0),
683 NumExportIDs_(Exporter.NumRemoteIDs_),//reverse
684 ExportLIDs_(0),
685 ExportPIDs_(0),
686 NumSend_(Exporter.NumRecv_),//reverse
687 NumRecv_(Exporter.NumSend_),//revsese
688 Distor_(0)
689{
690 // Reverse the permutes
691 if (NumPermuteIDs_ > 0) {
694 for (int i = 0; i < NumPermuteIDs_; ++i) {
695 PermuteFromLIDs_[i] = Exporter.PermuteToLIDs_[i];
696 PermuteToLIDs_[i] = Exporter.PermuteFromLIDs_[i];
697 }
698 }
699
700 // Copy the exports to the remotes
701 if (NumRemoteIDs_>0) {
702 RemoteLIDs_ = new int[NumRemoteIDs_];
703 for (int i = 0; i < NumRemoteIDs_; ++i) RemoteLIDs_[i] = Exporter.ExportLIDs_[i];
704 }
705
706 // Copy the remotes to the exports
707 if (NumExportIDs_>0) {
708 ExportLIDs_ = new int[NumExportIDs_];
709 ExportPIDs_ = new int[NumExportIDs_];
710 for (int i = 0; i < NumExportIDs_; ++i) ExportLIDs_[i] = Exporter.RemoteLIDs_[i];
711
712 // Extract the RemotePIDs from the Distributor
713#ifdef HAVE_MPI
714 Epetra_MpiDistributor *D=dynamic_cast<Epetra_MpiDistributor*>(&Exporter.Distributor());
715 if(!D) throw ReportError("Epetra_Import: Can't have ExportPIDs w/o an Epetra::MpiDistributor.",-1);
716
717 // Get the distributor's data
718 const int NumReceives = D->NumReceives();
719 const int *ProcsFrom = D->ProcsFrom();
720 const int *LengthsFrom = D->LengthsFrom();
721
722 // Now, for each remote ID, record who actually owns it. This loop follows the operation order in the
723 // MpiDistributor so it ought to duplicate that effect.
724 int i = 0, j = 0;
725 for (i = 0, j = 0; i < NumReceives; ++i) {
726 const int pid = ProcsFrom[i];
727 for (int k = 0; k < LengthsFrom[i]; ++k) {
728 ExportPIDs_[j] = pid;
729 j++;
730 }
731 }
732#else
733 throw ReportError("Epetra_Import: Can't have ExportPIDs w/o an Epetra::MpiDistributor.",-2);
734#endif
735 }//end NumExportIDs>0
736
737 if (Exporter.Distor_!=0) Distor_ = Exporter.Distor_->ReverseClone();
738}
739
740//=============================================================================
741void Epetra_Import::Print(std::ostream & os) const
742{
743 // mfh 14 Dec 2011: The implementation of Print() I found here
744 // previously didn't print much at all, and it included a message
745 // saying that it wasn't finished ("Epetra_Import::Print needs
746 // attention!!!"). What you see below is a port of
747 // Tpetra::Import::print, which does have a full implementation.
748 // This should allow a side-by-side comparison of Epetra_Import with
749 // Tpetra::Import.
750
751 // If true, then copy the array data and sort it before printing.
752 // Otherwise, leave the data in its original order.
753 //
754 // NOTE: Do NOT sort the arrays in place! Only sort in the copy.
755 // Epetra depends on the order being preserved, and some arrays'
756 // orders are coupled.
757 const bool sortIDs = false;
758
759 const Epetra_Comm& comm = SourceMap_.Comm();
760 const int myRank = comm.MyPID();
761 const int numProcs = comm.NumProc();
762
763 if (myRank == 0) {
764 os << "Import Data Members:" << std::endl;
765 }
766 // We don't need a barrier before this for loop, because Proc 0 is
767 // the first one to do anything in the for loop anyway.
768 for (int p = 0; p < numProcs; ++p) {
769 if (myRank == p) {
770 os << "Image ID : " << myRank << std::endl;
771
772 os << "permuteFromLIDs:";
773 if (PermuteFromLIDs_ == NULL) {
774 os << " NULL";
775 } else {
776 std::vector<int> permuteFromLIDs (NumPermuteIDs_);
778 permuteFromLIDs.begin());
779 if (sortIDs) {
780 std::sort (permuteFromLIDs.begin(), permuteFromLIDs.end());
781 }
782 os << " {";
783 for (int i = 0; i < NumPermuteIDs_; ++i) {
784 os << permuteFromLIDs[i];
785 if (i < NumPermuteIDs_ - 1) {
786 os << ", ";
787 }
788 }
789 os << "}";
790 }
791 os << std::endl;
792
793 os << "permuteToLIDs :";
794 if (PermuteToLIDs_ == NULL) {
795 os << " NULL";
796 } else {
797 std::vector<int> permuteToLIDs (NumPermuteIDs_);
799 permuteToLIDs.begin());
800 if (sortIDs) {
801 std::sort (permuteToLIDs.begin(), permuteToLIDs.end());
802 }
803 os << " {";
804 for (int i = 0; i < NumPermuteIDs_; ++i) {
805 os << permuteToLIDs[i];
806 if (i < NumPermuteIDs_ - 1) {
807 os << ", ";
808 }
809 }
810 os << "}";
811 }
812 os << std::endl;
813
814 os << "remoteLIDs :";
815 if (RemoteLIDs_ == NULL) {
816 os << " NULL";
817 } else {
818 std::vector<int> remoteLIDs (NumRemoteIDs_);
820 remoteLIDs.begin());
821 if (sortIDs) {
822 std::sort (remoteLIDs.begin(), remoteLIDs.end());
823 }
824 os << " {";
825 for (int i = 0; i < NumRemoteIDs_; ++i) {
826 os << remoteLIDs[i];
827 if (i < NumRemoteIDs_ - 1) {
828 os << ", ";
829 }
830 }
831 os << "}";
832 }
833 os << std::endl;
834
835 // If sorting for output, the export LIDs and export PIDs have
836 // to be sorted together. We can use Epetra_Util::Sort, using
837 // the PIDs as the keys to match Tpetra::Import.
838 std::vector<int> exportLIDs (NumExportIDs_);
839 std::vector<int> exportPIDs (NumExportIDs_);
840 if (ExportLIDs_ != NULL) {
841 std::copy (ExportLIDs_, ExportLIDs_ + NumExportIDs_, exportLIDs.begin());
842 std::copy (ExportPIDs_, ExportPIDs_ + NumExportIDs_, exportPIDs.begin());
843
844 if (sortIDs && NumExportIDs_ > 0) {
845 int* intCompanions[1]; // Input for Epetra_Util::Sort().
846 intCompanions[0] = &exportLIDs[0];
848 0, (double**) NULL, 1, intCompanions, 0, 0);
849 }
850 }
851
852 os << "exportLIDs :";
853 if (ExportLIDs_ == NULL) {
854 os << " NULL";
855 } else {
856 os << " {";
857 for (int i = 0; i < NumExportIDs_; ++i) {
858 os << exportLIDs[i];
859 if (i < NumExportIDs_ - 1) {
860 os << ", ";
861 }
862 }
863 os << "}";
864 }
865 os << std::endl;
866
867 os << "exportImageIDs :";
868 if (ExportPIDs_ == NULL) {
869 os << " NULL";
870 } else {
871 os << " {";
872 for (int i = 0; i < NumExportIDs_; ++i) {
873 os << exportPIDs[i];
874 if (i < NumExportIDs_ - 1) {
875 os << ", ";
876 }
877 }
878 os << "}";
879 }
880 os << std::endl;
881
882 os << "numSameIDs : " << NumSameIDs_ << std::endl;
883 os << "numPermuteIDs : " << NumPermuteIDs_ << std::endl;
884 os << "numRemoteIDs : " << NumRemoteIDs_ << std::endl;
885 os << "numExportIDs : " << NumExportIDs_ << std::endl;
886
887 // Epetra keeps NumSend_ and NumRecv_, whereas in Tpetra, these
888 // are stored in the Distributor object. This is why we print
889 // them here.
890 os << "Number of sends: " << NumSend_ << std::endl;
891 os << "Number of recvs: " << NumRecv_ << std::endl;
892 } // if my rank is p
893
894 // A few global barriers give I/O a chance to complete.
895 comm.Barrier();
896 comm.Barrier();
897 comm.Barrier();
898 } // for each rank p
899
900 const bool printMaps = false;
901 if (printMaps) {
902 // The original implementation printed the Maps first. We moved
903 // printing the Maps to the end, for easy comparison with the
904 // output of Tpetra::Import::print().
905 if (myRank == 0) {
906 os << std::endl << std::endl << "Source Map:" << std::endl << std::flush;
907 }
908 comm.Barrier();
909 SourceMap_.Print(os);
910 comm.Barrier();
911
912 if (myRank == 0) {
913 os << std::endl << std::endl << "Target Map:" << std::endl << std::flush;
914 }
915 comm.Barrier();
916 TargetMap_.Print(os);
917 comm.Barrier();
918 }
919
920 if (myRank == 0) {
921 os << std::endl << std::endl << "Distributor:" << std::endl << std::flush;
922 }
923 comm.Barrier();
924 if (Distor_ == NULL) {
925 if (myRank == 0) {
926 os << " is NULL." << std::endl;
927 }
928 } else {
929 Distor_->Print(os); // Printing the Distributor is itself distributed.
930 }
931 comm.Barrier();
932}
933
@ Insert
#define EPETRA_MIN(x, y)
T * Epetra_Util_data_ptr(std::vector< T > &vec)
Function that returns either a pointer to the first entry in the vector or, if the vector is empty,...
Definition: Epetra_Util.h:422
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int RemoteIDList(int NumIDs, const int *GIDList, int *PIDList, int *LIDList) const
Returns the processor IDs and corresponding local index value for a given list of global indices.
int MaxElementSize() const
Maximum element size across all processors.
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
bool GlobalIndicesTypeMatch(const Epetra_BlockMap &other) const
virtual void Print(std::ostream &os) const
Print object to an output stream.
bool DistributedGlobal() const
Returns true if map is defined across more than one processor.
int LID(int GID) const
Returns local ID of global ID, return -1 if not found on this processor.
int GID(int LID) const
Returns global ID of local ID, return IndexBase-1 if not found on this processor.
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
int NumMyElements() const
Number of elements on the calling processor.
bool MyGID(int GID_in) const
Returns true if the GID passed in belongs to the calling processor in this map, otherwise returns fal...
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
virtual int NumProc() const =0
Returns total number of processes.
virtual int MinAll(double *PartialMins, double *GlobalMins, int Count) const =0
Epetra_Comm Global Min function.
virtual int MyPID() const =0
Return my process ID.
virtual Epetra_Distributor * CreateDistributor() const =0
Create a distributor object.
virtual void Barrier() const =0
Epetra_Comm Barrier function.
const Epetra_BlockMap & Map() const
Returns the address of the Epetra_BlockMap for this multi-vector.
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Imports an Epetra_DistObject using the Epetra_Import object.
virtual Epetra_Distributor * Clone()=0
Epetra_Distributor clone constructor.
virtual int CreateFromRecvs(const int &NumRemoteIDs, const int *RemoteGIDs, const int *RemotePIDs, bool Deterministic, int &NumExportIDs, int *&ExportGIDs, int *&ExportPIDs)=0
Create Distributor object using list of Remote global IDs and corresponding PIDs.
virtual void Print(std::ostream &os) const =0
virtual Epetra_Distributor * ReverseClone()=0
Create and extract the reverse version of the distributor.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements.
Definition: Epetra_Export.h:62
int * PermuteToLIDs_
Epetra_Distributor & Distributor() const
Epetra_Distributor * Distor_
int * PermuteFromLIDs_
Epetra_Import: This class builds an import object for efficient importing of off-processor elements.
Definition: Epetra_Import.h:63
Epetra_BlockMap SourceMap_
virtual void Print(std::ostream &os) const
Print object to an output stream Print method.
int * PermuteFromLIDs_
Epetra_BlockMap TargetMap_
void Construct(const Epetra_BlockMap &targetMap, const Epetra_BlockMap &sourceMap, int NumRemotePIDs=-1, const int *UserRemotePIDs=0)
void Construct_Expert(const Epetra_BlockMap &TargetMap, const Epetra_BlockMap &SourceMap, int NumRemotePIDs, const int *RemotePIDs, const int &NumExportIDs, const int *ExportLIDs, const int *ExportPIDs)
int * PermuteToLIDs_
Epetra_Import(const Epetra_BlockMap &TargetMap, const Epetra_BlockMap &SourceMap)
Constructs a Epetra_Import object from the source and target maps.
Epetra_Distributor * Distor_
virtual ~Epetra_Import(void)
Epetra_Import destructor.
Epetra_IntVector: A class for constructing and using dense integer vectors on a parallel computer.
int PutValue(int Value)
Set all elements of the vector to Value.
int MyLength() const
Returns the local vector length on the calling processor of vectors in the multi-vector.
MPI implementation of Epetra_Distributor.
int CreateFromSendsAndRecvs(const int &NumExportIDs, const int *ExportPIDs, const int &NumRemoteIDs, const int *RemoteGIDs, const int *RemotePIDs, bool Deterministic)
Create a communication plan from send list and a recv list.
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:57
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
Epetra_Util: The Epetra Util Wrapper Class.
Definition: Epetra_Util.h:79
static void EPETRA_LIB_DLL_EXPORT Sort(bool SortAscending, int NumKeys, T *Keys, int NumDoubleCompanions, double **DoubleCompanions, int NumIntCompanions, int **IntCompanions, int NumLongLongCompanions, long long **LongLongCompanions)
Epetra_Util Sort Routine (Shell sort)