Epetra Package Browser (Single Doxygen Collection)  Development
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
57 #include "Epetra_MpiDistributor.h"
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
66 template<typename int_type>
67 void 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
99  NumPermuteIDs_ = 0;
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) {
115  PermuteToLIDs_ = new int[NumPermuteIDs_];
116  PermuteFromLIDs_ = new int[NumPermuteIDs_];
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];
130  RemoteLIDs_[NumRemoteIDs_++] = 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  }
161 
162  //Get rid of IDs that don't exist in SourceMap
163  if(NumRemoteIDs_>0) {
164  int cnt = 0;
165  for( i = 0; i < NumRemoteIDs_; ++i )
166  if( RemotePIDs[i] == -1 ) ++cnt;
167  if( cnt ) {
168  if( NumRemoteIDs_-cnt ) {
169  int_type * NewRemoteGIDs = new int_type[NumRemoteIDs_-cnt];
170  int * NewRemotePIDs = new int[NumRemoteIDs_-cnt];
171  int * NewRemoteLIDs = new int[NumRemoteIDs_-cnt];
172  cnt = 0;
173  for( i = 0; i < NumRemoteIDs_; ++i )
174  if( RemotePIDs[i] != -1 ) {
175  NewRemoteGIDs[cnt] = RemoteGIDs[i];
176  NewRemotePIDs[cnt] = RemotePIDs[i];
177  NewRemoteLIDs[cnt] = targetMap.LID(RemoteGIDs[i]);
178  ++cnt;
179  }
180  NumRemoteIDs_ = cnt;
181  delete [] RemoteGIDs;
182  delete [] RemotePIDs;
183  delete [] RemoteLIDs_;
184  RemoteGIDs = NewRemoteGIDs;
185  RemotePIDs = NewRemotePIDs;
186  RemoteLIDs_ = NewRemoteLIDs;
187  ReportError("Warning in Epetra_Import: Target IDs not found in Source Map (Do you want to import to subset of Target Map?)", 1);
188  }
189  else { //valid RemoteIDs empty
190  NumRemoteIDs_ = 0;
191  delete [] RemoteGIDs;
192  RemoteGIDs = 0;
193  delete [] RemotePIDs;
194  RemotePIDs = 0;
195  }
196  }
197  }
198 
199  //Sort Remote IDs by processor so DoReverses will work
200  Epetra_Util util;
201 
202  if(targetMap.GlobalIndicesLongLong())
203  {
204  // FIXME (mfh 11 Jul 2013) This breaks ANSI aliasing rules, if
205  // int_type != long long. On some compilers, it results in
206  // warnings such as this: "dereferencing type-punned pointer
207  // will break strict-aliasing rules".
208  util.Sort(true,NumRemoteIDs_,RemotePIDs,0,0, 1,&RemoteLIDs_, 1,(long long**)&RemoteGIDs);
209  }
210  else if(targetMap.GlobalIndicesInt())
211  {
212  int* ptrs[2] = {RemoteLIDs_, (int*)RemoteGIDs};
213  util.Sort(true,NumRemoteIDs_,RemotePIDs,0,0,2,&ptrs[0], 0, 0);
214  }
215  else
216  {
217  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
218  }
219 
220  // Build distributor & Export lists
221  Distor_ = sourceMap.Comm().CreateDistributor();
222 
223  NumExportIDs_=UserNumExportIDs;
224  ExportLIDs_ = new int[NumExportIDs_];
225  ExportPIDs_ = new int[NumExportIDs_];
226  for(i=0; i<NumExportIDs_; i++) {
227  ExportPIDs_[i] = UserExportPIDs[i];
228  ExportLIDs_[i] = UserExportLIDs[i];
229  }
230 
231  // Send Counts
232  for (i=0; i< NumExportIDs_; i++) {
233  NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
234  }
235 
236 #ifdef HAVE_MPI
237  Epetra_MpiDistributor* MpiDistor = dynamic_cast< Epetra_MpiDistributor*>(Distor_);
238  bool Deterministic = true;
239  if(MpiDistor)
241  NumRemoteIDs_, RemoteGIDs, RemotePIDs,Deterministic);
242  else ierr=-10;
243 #else
244  ierr=-20;
245 #endif
246 
247  if (ierr!=0) throw ReportError("Error in Epetra_Distributor.CreateFromRecvs()", ierr);
248  }
249 
250  if( NumRemoteIDs_>0 ) delete [] RemoteGIDs;
251  if( NumRemoteIDs_>0 ) delete [] RemotePIDs;
252 
253  if (NumTargetIDs>0) delete [] TargetGIDs;
254  if (NumSourceIDs>0) delete [] SourceGIDs;
255 
256 
257 #ifdef EPETRA_ENABLE_DEBUG
258 // Sanity check to make sure we got the import right
259  Epetra_IntVector Source(sourceMap);
260  Epetra_IntVector Target(targetMap);
261 
262  for(i=0; i<Source.MyLength(); i++)
263  Source[i] = (int) (Source.Map().GID(i) % INT_MAX);
264  Target.PutValue(-1);
265 
266  Target.Import(Source,*this,Insert);
267 
268  bool test_passed=true;
269  for(i=0; i<Target.MyLength(); i++){
270  if(Target[i] != Target.Map().GID(i) % INT_MAX) test_passed=false;
271  }
272 
273  if(!test_passed) {
274  printf("[%d] PROCESSOR has a mismatch... prepearing to crash or hang!\n",sourceMap.Comm().MyPID());
275  fflush(stdout);
276  sourceMap.Comm().Barrier();
277  sourceMap.Comm().Barrier();
278  sourceMap.Comm().Barrier();
279  throw ReportError("Epetra_Import: ERROR. User provided IDs do not match what an import generates.");
280  }
281 #endif
282 
283  return;
284 }
285 
286 //==============================================================================
287 // Epetra_Import constructor function for a Epetra_BlockMap object
288 template<typename int_type>
289 void Epetra_Import::Construct( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap, int NumRemotePIDs, const int * UserRemotePIDs)
290 {
291  int i,ierr;
292  // Build three ID lists:
293  // NumSameIDs - Number of IDs in TargetMap and SourceMap that are identical, up to the first
294  // nonidentical ID.
295  // NumPermuteIDs - Number of IDs in SourceMap that must be indirectly loaded but are on this processor.
296  // NumRemoteIDs - Number of IDs that are in SourceMap but not in TargetMap, and thus must be imported.
297 
298  int NumSourceIDs = sourceMap.NumMyElements();
299  int NumTargetIDs = targetMap.NumMyElements();
300 
301  int_type *TargetGIDs = 0;
302  if (NumTargetIDs>0) {
303  TargetGIDs = new int_type[NumTargetIDs];
304  targetMap.MyGlobalElements(TargetGIDs);
305  }
306 
307  int_type * SourceGIDs = 0;
308  if (NumSourceIDs>0) {
309  SourceGIDs = new int_type[NumSourceIDs];
310  sourceMap.MyGlobalElements(SourceGIDs);
311  }
312 
313  int MinIDs = EPETRA_MIN(NumSourceIDs, NumTargetIDs);
314 
315 
316  NumSameIDs_ = 0;
317  for (i=0; i< MinIDs; i++) if (TargetGIDs[i]==SourceGIDs[i]) NumSameIDs_++; else break;
318 
319  // Find count of Target IDs that are truly remote and those that are local but permuted
320  NumPermuteIDs_ = 0;
321  NumRemoteIDs_ = 0;
322  for (i=NumSameIDs_; i< NumTargetIDs; i++)
323  if (sourceMap.MyGID(TargetGIDs[i])) NumPermuteIDs_++; // Check if Target GID is a local Source GID
324  else NumRemoteIDs_++; // If not, then it is remote
325 
326 
327 
328  // Define remote and permutation lists
329  int_type * RemoteGIDs=0;
330  RemoteLIDs_ = 0;
331  if (NumRemoteIDs_>0) {
332  RemoteLIDs_ = new int[NumRemoteIDs_];
333  RemoteGIDs = new int_type[NumRemoteIDs_];
334  }
335  if (NumPermuteIDs_>0) {
336  PermuteToLIDs_ = new int[NumPermuteIDs_];
337  PermuteFromLIDs_ = new int[NumPermuteIDs_];
338  }
339 
340  NumPermuteIDs_ = 0;
341  NumRemoteIDs_ = 0;
342  for (i=NumSameIDs_; i< NumTargetIDs; i++) {
343  if (sourceMap.MyGID(TargetGIDs[i])) {
345  PermuteFromLIDs_[NumPermuteIDs_++] = sourceMap.LID(TargetGIDs[i]);
346  }
347  else {
348  //NumRecv_ +=TargetMap.ElementSize(i); // Count total number of entries to receive
349  NumRecv_ +=targetMap.MaxElementSize(); // Count total number of entries to receive (currently need max)
350  RemoteGIDs[NumRemoteIDs_] = TargetGIDs[i];
351  RemoteLIDs_[NumRemoteIDs_++] = i;
352  }
353  }
354 
355  if( NumRemoteIDs_>0 && !sourceMap.DistributedGlobal() )
356  ReportError("Warning in Epetra_Import: Serial Import has remote IDs. (Importing to Subset of Target Map)", 1);
357 
358  // Test for distributed cases
359  int * RemotePIDs = 0;
360 
361  if (sourceMap.DistributedGlobal()) {
362  if (NumRemoteIDs_>0) RemotePIDs = new int[NumRemoteIDs_];
363 
364 #ifdef EPETRA_ENABLE_DEBUG
365  if(NumRemotePIDs!=-1){
366  int myeq = (NumRemotePIDs==NumRemoteIDs_);
367  int globaleq=0;
368  sourceMap.Comm().MinAll(&myeq,&globaleq,1);
369  if(globaleq!=1) {
370  printf("[%d] UserRemotePIDs count wrong %d != %d\n",sourceMap.Comm().MyPID(),NumRemotePIDs,NumRemoteIDs_);
371  fflush(stdout);
372  sourceMap.Comm().Barrier();
373  sourceMap.Comm().Barrier();
374  sourceMap.Comm().Barrier();
375  throw ReportError("Epetra_Import: UserRemotePIDs count wrong",-1);
376  }
377  }
378 #endif
379 
380  if(NumRemotePIDs!=-1 && NumRemotePIDs==NumRemoteIDs_){
381  // Since I need to sort these, I'll copy them
382  for(i=0; i<NumRemoteIDs_; i++) RemotePIDs[i] = UserRemotePIDs[i];
383  }
384  else{
385  ierr=sourceMap.RemoteIDList(NumRemoteIDs_, RemoteGIDs, RemotePIDs, 0); // Get remote PIDs
386  if (ierr) throw ReportError("Error in sourceMap.RemoteIDList call", ierr);
387  }
388 
389  //Get rid of IDs that don't exist in SourceMap
390  if(NumRemoteIDs_>0) {
391  int cnt = 0;
392  for( i = 0; i < NumRemoteIDs_; ++i )
393  if( RemotePIDs[i] == -1 ) ++cnt;
394  if( cnt ) {
395  if( NumRemoteIDs_-cnt ) {
396  int_type * NewRemoteGIDs = new int_type[NumRemoteIDs_-cnt];
397  int * NewRemotePIDs = new int[NumRemoteIDs_-cnt];
398  int * NewRemoteLIDs = new int[NumRemoteIDs_-cnt];
399  cnt = 0;
400  for( i = 0; i < NumRemoteIDs_; ++i )
401  if( RemotePIDs[i] != -1 ) {
402  NewRemoteGIDs[cnt] = RemoteGIDs[i];
403  NewRemotePIDs[cnt] = RemotePIDs[i];
404  NewRemoteLIDs[cnt] = targetMap.LID(RemoteGIDs[i]);
405  ++cnt;
406  }
407  NumRemoteIDs_ = cnt;
408  delete [] RemoteGIDs;
409  delete [] RemotePIDs;
410  delete [] RemoteLIDs_;
411  RemoteGIDs = NewRemoteGIDs;
412  RemotePIDs = NewRemotePIDs;
413  RemoteLIDs_ = NewRemoteLIDs;
414  ReportError("Warning in Epetra_Import: Target IDs not found in Source Map (Do you want to import to subset of Target Map?)", 1);
415  }
416  else { //valid RemoteIDs empty
417  NumRemoteIDs_ = 0;
418  delete [] RemoteGIDs;
419  RemoteGIDs = 0;
420  delete [] RemotePIDs;
421  RemotePIDs = 0;
422  }
423  }
424  }
425 
426  //Sort Remote IDs by processor so DoReverses will work
427  if(targetMap.GlobalIndicesLongLong())
428  {
429  // FIXME (mfh 11 Jul 2013) This breaks ANSI aliasing rules, if
430  // int_type != long long. On some compilers, it results in
431  // warnings such as this: "dereferencing type-punned pointer
432  // will break strict-aliasing rules".
433  Epetra_Util::Sort(true,NumRemoteIDs_,RemotePIDs,0,0, 1,&RemoteLIDs_, 1,(long long**)&RemoteGIDs);
434  }
435  else if(targetMap.GlobalIndicesInt())
436  {
437  int* ptrs[2] = {RemoteLIDs_, (int*)RemoteGIDs};
438  Epetra_Util::Sort(true,NumRemoteIDs_,RemotePIDs,0,0,2,&ptrs[0], 0, 0);
439  }
440  else
441  {
442  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
443  }
444  Distor_ = sourceMap.Comm().CreateDistributor();
445 
446  // Construct list of exports that calling processor needs to send as a result
447  // of everyone asking for what it needs to receive.
448 
449  bool Deterministic = true;
450  int_type* tmp_ExportLIDs; //Export IDs come in as GIDs
451  ierr = Distor_->CreateFromRecvs( NumRemoteIDs_, RemoteGIDs, RemotePIDs,
452  Deterministic, NumExportIDs_, tmp_ExportLIDs, ExportPIDs_ );
453  if (ierr!=0) throw ReportError("Error in Epetra_Distributor.CreateFromRecvs()", ierr);
454 
455 
456  // Export IDs come in as GIDs, convert to LIDs
457  if(targetMap.GlobalIndicesLongLong())
458  {
459  ExportLIDs_ = new int[NumExportIDs_];
460 
461  for (i=0; i< NumExportIDs_; i++) {
462  if (ExportPIDs_[i] < 0) throw ReportError("targetMap requested a GID that is not in the sourceMap.", -1);
463  ExportLIDs_[i] = sourceMap.LID(tmp_ExportLIDs[i]);
464  NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
465  }
466 
467  delete[] tmp_ExportLIDs;
468  }
469  else if(targetMap.GlobalIndicesInt())
470  {
471  for (i=0; i< NumExportIDs_; i++) {
472  if (ExportPIDs_[i] < 0) throw ReportError("targetMap requested a GID that is not in the sourceMap.", -1);
473  tmp_ExportLIDs[i] = sourceMap.LID(tmp_ExportLIDs[i]);
474  NumSend_ += sourceMap.MaxElementSize(); // Count total number of entries to send (currently need max)
475  }
476 
477  ExportLIDs_ = reinterpret_cast<int *>(tmp_ExportLIDs); // Can't reach here if tmp_ExportLIDs is long long.
478  }
479  else
480  {
481  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndices Internal Error", -1);
482  }
483  }
484 
485  if( NumRemoteIDs_>0 ) delete [] RemoteGIDs;
486  if( NumRemoteIDs_>0 ) delete [] RemotePIDs;
487 
488  if (NumTargetIDs>0) delete [] TargetGIDs;
489  if (NumSourceIDs>0) delete [] SourceGIDs;
490 
491  return;
492 }
493 
494 //==============================================================================
495 Epetra_Import::Epetra_Import( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap,int NumRemotePIDs,const int * RemotePIDs)
496  : Epetra_Object("Epetra::Import"),
497  TargetMap_(targetMap),
498  SourceMap_(sourceMap),
499  NumSameIDs_(0),
500  NumPermuteIDs_(0),
501  PermuteToLIDs_(0),
502  PermuteFromLIDs_(0),
503  NumRemoteIDs_(0),
504  RemoteLIDs_(0),
505  NumExportIDs_(0),
506  ExportLIDs_(0),
507  ExportPIDs_(0),
508  NumSend_(0),
509  NumRecv_(0),
510  Distor_(0)
511 {
512  if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
513  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
514 
515  if(targetMap.GlobalIndicesInt())
516 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
517  Construct<int>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs);
518 #else
519  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
520 #endif
521  else if(targetMap.GlobalIndicesLongLong())
522 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
523  Construct<long long>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs);
524 #else
525  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
526 #endif
527  else
528  throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
529 }
530 
531 
532 
533 //==============================================================================
534 Epetra_Import::Epetra_Import( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap,int NumRemotePIDs,const int * RemotePIDs,
535  const int & numExportIDs, const int * theExportLIDs, const int * theExportPIDs)
536  : Epetra_Object("Epetra::Import"),
537  TargetMap_(targetMap),
538  SourceMap_(sourceMap),
539  NumSameIDs_(0),
540  NumPermuteIDs_(0),
541  PermuteToLIDs_(0),
542  PermuteFromLIDs_(0),
543  NumRemoteIDs_(0),
544  RemoteLIDs_(0),
545  NumExportIDs_(0),
546  ExportLIDs_(0),
547  ExportPIDs_(0),
548  NumSend_(0),
549  NumRecv_(0),
550  Distor_(0)
551 {
552  if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
553  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
554 
555  if(targetMap.GlobalIndicesInt())
556 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
557  Construct_Expert<int>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs,numExportIDs,theExportLIDs,theExportPIDs);
558 #else
559  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
560 #endif
561  else if(targetMap.GlobalIndicesLongLong())
562 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
563  Construct_Expert<long long>(targetMap, sourceMap,NumRemotePIDs,RemotePIDs,numExportIDs,theExportLIDs,theExportPIDs);
564 #else
565  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
566 #endif
567  else
568  throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
569 }
570 
571 
572 //==============================================================================
573 Epetra_Import::Epetra_Import( const Epetra_BlockMap & targetMap, const Epetra_BlockMap & sourceMap)
574  : Epetra_Object("Epetra::Import"),
575  TargetMap_(targetMap),
576  SourceMap_(sourceMap),
577  NumSameIDs_(0),
578  NumPermuteIDs_(0),
579  PermuteToLIDs_(0),
580  PermuteFromLIDs_(0),
581  NumRemoteIDs_(0),
582  RemoteLIDs_(0),
583  NumExportIDs_(0),
584  ExportLIDs_(0),
585  ExportPIDs_(0),
586  NumSend_(0),
587  NumRecv_(0),
588  Distor_(0)
589 {
590  if(!targetMap.GlobalIndicesTypeMatch(sourceMap))
591  throw ReportError("Epetra_Import::Epetra_Import: GlobalIndicesTypeMatch failed", -1);
592 
593  if(targetMap.GlobalIndicesInt())
594 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
595  Construct<int>(targetMap, sourceMap);
596 #else
597  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesInt but no API for it.",-1);
598 #endif
599  else if(targetMap.GlobalIndicesLongLong())
600 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
601  Construct<long long>(targetMap, sourceMap);
602 #else
603  throw ReportError("Epetra_Import::Epetra_Import: ERROR, GlobalIndicesLongLong but no API for it.",-1);
604 #endif
605  else
606  throw ReportError("Epetra_Import::Epetra_Import: Bad global indices type", -1);
607 }
608 
609 //==============================================================================
610 // Epetra_Import copy constructor
612  : Epetra_Object(Importer),
613  TargetMap_(Importer.TargetMap_),
614  SourceMap_(Importer.SourceMap_),
615  NumSameIDs_(Importer.NumSameIDs_),
616  NumPermuteIDs_(Importer.NumPermuteIDs_),
617  PermuteToLIDs_(0),
618  PermuteFromLIDs_(0),
619  NumRemoteIDs_(Importer.NumRemoteIDs_),
620  RemoteLIDs_(0),
621  NumExportIDs_(Importer.NumExportIDs_),
622  ExportLIDs_(0),
623  ExportPIDs_(0),
624  NumSend_(Importer.NumSend_),
625  NumRecv_(Importer.NumRecv_),
626  Distor_(0)
627 {
628  int i;
629  if (NumPermuteIDs_>0) {
630  PermuteToLIDs_ = new int[NumPermuteIDs_];
631  PermuteFromLIDs_ = new int[NumPermuteIDs_];
632  for (i=0; i< NumPermuteIDs_; i++) {
633  PermuteToLIDs_[i] = Importer.PermuteToLIDs_[i];
634  PermuteFromLIDs_[i] = Importer.PermuteFromLIDs_[i];
635  }
636  }
637 
638  if (NumRemoteIDs_>0) {
639  RemoteLIDs_ = new int[NumRemoteIDs_];
640  for (i=0; i< NumRemoteIDs_; i++) RemoteLIDs_[i] = Importer.RemoteLIDs_[i];
641  }
642 
643  if (NumExportIDs_>0) {
644  ExportLIDs_ = new int[NumExportIDs_];
645  ExportPIDs_ = new int[NumExportIDs_];
646  for (i=0; i< NumExportIDs_; i++) {
647  ExportLIDs_[i] = Importer.ExportLIDs_[i];
648  ExportPIDs_[i] = Importer.ExportPIDs_[i];
649  }
650  }
651 
652  if (Importer.Distor_!=0) Distor_ = Importer.Distor_->Clone();
653 
654 }
655 
656 //==============================================================================
657 // Epetra_Import destructor
659 {
660  if( Distor_ != 0 ) delete Distor_;
661  if (RemoteLIDs_ != 0) delete [] RemoteLIDs_;
662  if (PermuteToLIDs_ != 0) delete [] PermuteToLIDs_;
663  if (PermuteFromLIDs_ != 0) delete [] PermuteFromLIDs_;
664 
665  if( ExportPIDs_ != 0 ) delete [] ExportPIDs_; // These were created by Distor_
666  if( ExportLIDs_ != 0 ) delete [] ExportLIDs_;
667 }
668 
669 //==============================================================================
670 // Epetra_Import pseudo-copy constructor.
672  TargetMap_(Exporter.SourceMap_), //reverse
673  SourceMap_(Exporter.TargetMap_),//reverse
674  NumSameIDs_(Exporter.NumSameIDs_),
675  NumPermuteIDs_(Exporter.NumPermuteIDs_),
676  PermuteToLIDs_(0),
677  PermuteFromLIDs_(0),
678  NumRemoteIDs_(Exporter.NumExportIDs_),//reverse
679  RemoteLIDs_(0),
680  NumExportIDs_(Exporter.NumRemoteIDs_),//reverse
681  ExportLIDs_(0),
682  ExportPIDs_(0),
683  NumSend_(Exporter.NumRecv_),//reverse
684  NumRecv_(Exporter.NumSend_),//revsese
685  Distor_(0)
686 {
687  // Reverse the permutes
688  if (NumPermuteIDs_ > 0) {
689  PermuteToLIDs_ = new int[NumPermuteIDs_];
690  PermuteFromLIDs_ = new int[NumPermuteIDs_];
691  for (int i = 0; i < NumPermuteIDs_; ++i) {
692  PermuteFromLIDs_[i] = Exporter.PermuteToLIDs_[i];
693  PermuteToLIDs_[i] = Exporter.PermuteFromLIDs_[i];
694  }
695  }
696 
697  // Copy the exports to the remotes
698  if (NumRemoteIDs_>0) {
699  RemoteLIDs_ = new int[NumRemoteIDs_];
700  for (int i = 0; i < NumRemoteIDs_; ++i) RemoteLIDs_[i] = Exporter.ExportLIDs_[i];
701  }
702 
703  // Copy the remotes to the exports
704  if (NumExportIDs_>0) {
705  ExportLIDs_ = new int[NumExportIDs_];
706  ExportPIDs_ = new int[NumExportIDs_];
707  for (int i = 0; i < NumExportIDs_; ++i) ExportLIDs_[i] = Exporter.RemoteLIDs_[i];
708 
709  // Extract the RemotePIDs from the Distributor
710 #ifdef HAVE_MPI
711  Epetra_MpiDistributor *D=dynamic_cast<Epetra_MpiDistributor*>(&Exporter.Distributor());
712  if(!D) throw ReportError("Epetra_Import: Can't have ExportPIDs w/o an Epetra::MpiDistributor.",-1);
713 
714  // Get the distributor's data
715  const int NumReceives = D->NumReceives();
716  const int *ProcsFrom = D->ProcsFrom();
717  const int *LengthsFrom = D->LengthsFrom();
718 
719  // Now, for each remote ID, record who actually owns it. This loop follows the operation order in the
720  // MpiDistributor so it ought to duplicate that effect.
721  int i = 0, j = 0;
722  for (i = 0, j = 0; i < NumReceives; ++i) {
723  const int pid = ProcsFrom[i];
724  for (int k = 0; k < LengthsFrom[i]; ++k) {
725  ExportPIDs_[j] = pid;
726  j++;
727  }
728  }
729 #else
730  throw ReportError("Epetra_Import: Can't have ExportPIDs w/o an Epetra::MpiDistributor.",-2);
731 #endif
732  }//end NumExportIDs>0
733 
734  if (Exporter.Distor_!=0) Distor_ = Exporter.Distor_->ReverseClone();
735 }
736 
737 //=============================================================================
738 void Epetra_Import::Print(std::ostream & os) const
739 {
740  // mfh 14 Dec 2011: The implementation of Print() I found here
741  // previously didn't print much at all, and it included a message
742  // saying that it wasn't finished ("Epetra_Import::Print needs
743  // attention!!!"). What you see below is a port of
744  // Tpetra::Import::print, which does have a full implementation.
745  // This should allow a side-by-side comparison of Epetra_Import with
746  // Tpetra::Import.
747 
748  // If true, then copy the array data and sort it before printing.
749  // Otherwise, leave the data in its original order.
750  //
751  // NOTE: Do NOT sort the arrays in place! Only sort in the copy.
752  // Epetra depends on the order being preserved, and some arrays'
753  // orders are coupled.
754  const bool sortIDs = false;
755 
756  const Epetra_Comm& comm = SourceMap_.Comm();
757  const int myRank = comm.MyPID();
758  const int numProcs = comm.NumProc();
759 
760  if (myRank == 0) {
761  os << "Import Data Members:" << std::endl;
762  }
763  // We don't need a barrier before this for loop, because Proc 0 is
764  // the first one to do anything in the for loop anyway.
765  for (int p = 0; p < numProcs; ++p) {
766  if (myRank == p) {
767  os << "Image ID : " << myRank << std::endl;
768 
769  os << "permuteFromLIDs:";
770  if (PermuteFromLIDs_ == NULL) {
771  os << " NULL";
772  } else {
773  std::vector<int> permuteFromLIDs (NumPermuteIDs_);
775  permuteFromLIDs.begin());
776  if (sortIDs) {
777  std::sort (permuteFromLIDs.begin(), permuteFromLIDs.end());
778  }
779  os << " {";
780  for (int i = 0; i < NumPermuteIDs_; ++i) {
781  os << permuteFromLIDs[i];
782  if (i < NumPermuteIDs_ - 1) {
783  os << ", ";
784  }
785  }
786  os << "}";
787  }
788  os << std::endl;
789 
790  os << "permuteToLIDs :";
791  if (PermuteToLIDs_ == NULL) {
792  os << " NULL";
793  } else {
794  std::vector<int> permuteToLIDs (NumPermuteIDs_);
796  permuteToLIDs.begin());
797  if (sortIDs) {
798  std::sort (permuteToLIDs.begin(), permuteToLIDs.end());
799  }
800  os << " {";
801  for (int i = 0; i < NumPermuteIDs_; ++i) {
802  os << permuteToLIDs[i];
803  if (i < NumPermuteIDs_ - 1) {
804  os << ", ";
805  }
806  }
807  os << "}";
808  }
809  os << std::endl;
810 
811  os << "remoteLIDs :";
812  if (RemoteLIDs_ == NULL) {
813  os << " NULL";
814  } else {
815  std::vector<int> remoteLIDs (NumRemoteIDs_);
816  std::copy (RemoteLIDs_, RemoteLIDs_ + NumRemoteIDs_,
817  remoteLIDs.begin());
818  if (sortIDs) {
819  std::sort (remoteLIDs.begin(), remoteLIDs.end());
820  }
821  os << " {";
822  for (int i = 0; i < NumRemoteIDs_; ++i) {
823  os << remoteLIDs[i];
824  if (i < NumRemoteIDs_ - 1) {
825  os << ", ";
826  }
827  }
828  os << "}";
829  }
830  os << std::endl;
831 
832  // If sorting for output, the export LIDs and export PIDs have
833  // to be sorted together. We can use Epetra_Util::Sort, using
834  // the PIDs as the keys to match Tpetra::Import.
835  std::vector<int> exportLIDs (NumExportIDs_);
836  std::vector<int> exportPIDs (NumExportIDs_);
837  if (ExportLIDs_ != NULL) {
838  std::copy (ExportLIDs_, ExportLIDs_ + NumExportIDs_, exportLIDs.begin());
839  std::copy (ExportPIDs_, ExportPIDs_ + NumExportIDs_, exportPIDs.begin());
840 
841  if (sortIDs && NumExportIDs_ > 0) {
842  int* intCompanions[1]; // Input for Epetra_Util::Sort().
843  intCompanions[0] = &exportLIDs[0];
844  Epetra_Util::Sort (true, NumExportIDs_, &exportPIDs[0],
845  0, (double**) NULL, 1, intCompanions, 0, 0);
846  }
847  }
848 
849  os << "exportLIDs :";
850  if (ExportLIDs_ == NULL) {
851  os << " NULL";
852  } else {
853  os << " {";
854  for (int i = 0; i < NumExportIDs_; ++i) {
855  os << exportLIDs[i];
856  if (i < NumExportIDs_ - 1) {
857  os << ", ";
858  }
859  }
860  os << "}";
861  }
862  os << std::endl;
863 
864  os << "exportImageIDs :";
865  if (ExportPIDs_ == NULL) {
866  os << " NULL";
867  } else {
868  os << " {";
869  for (int i = 0; i < NumExportIDs_; ++i) {
870  os << exportPIDs[i];
871  if (i < NumExportIDs_ - 1) {
872  os << ", ";
873  }
874  }
875  os << "}";
876  }
877  os << std::endl;
878 
879  os << "numSameIDs : " << NumSameIDs_ << std::endl;
880  os << "numPermuteIDs : " << NumPermuteIDs_ << std::endl;
881  os << "numRemoteIDs : " << NumRemoteIDs_ << std::endl;
882  os << "numExportIDs : " << NumExportIDs_ << std::endl;
883 
884  // Epetra keeps NumSend_ and NumRecv_, whereas in Tpetra, these
885  // are stored in the Distributor object. This is why we print
886  // them here.
887  os << "Number of sends: " << NumSend_ << std::endl;
888  os << "Number of recvs: " << NumRecv_ << std::endl;
889  } // if my rank is p
890 
891  // A few global barriers give I/O a chance to complete.
892  comm.Barrier();
893  comm.Barrier();
894  comm.Barrier();
895  } // for each rank p
896 
897  const bool printMaps = false;
898  if (printMaps) {
899  // The original implementation printed the Maps first. We moved
900  // printing the Maps to the end, for easy comparison with the
901  // output of Tpetra::Import::print().
902  if (myRank == 0) {
903  os << std::endl << std::endl << "Source Map:" << std::endl << std::flush;
904  }
905  comm.Barrier();
906  SourceMap_.Print(os);
907  comm.Barrier();
908 
909  if (myRank == 0) {
910  os << std::endl << std::endl << "Target Map:" << std::endl << std::flush;
911  }
912  comm.Barrier();
913  TargetMap_.Print(os);
914  comm.Barrier();
915  }
916 
917  if (myRank == 0) {
918  os << std::endl << std::endl << "Distributor:" << std::endl << std::flush;
919  }
920  comm.Barrier();
921  if (Distor_ == NULL) {
922  if (myRank == 0) {
923  os << " is NULL." << std::endl;
924  }
925  } else {
926  Distor_->Print(os); // Printing the Distributor is itself distributed.
927  }
928  comm.Barrier();
929 }
930 
const Epetra_BlockMap & Map() const
Returns the address of the Epetra_BlockMap for this multi-vector.
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
virtual ~Epetra_Import(void)
Epetra_Import destructor.
int * PermuteToLIDs_
bool DistributedGlobal() const
Returns true if map is defined across more than one processor.
Epetra_Distributor & Distributor() const
Epetra_IntVector: A class for constructing and using dense integer vectors on a parallel computer...
Epetra_BlockMap TargetMap_
virtual void Print(std::ostream &os) const =0
int MyLength() const
Returns the local vector length on the calling processor of vectors in the multi-vector.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements...
Definition: Epetra_Export.h:62
int PutValue(int Value)
Set all elements of the vector to Value.
int * PermuteFromLIDs_
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.
#define EPETRA_MIN(x, y)
MPI implementation of Epetra_Distributor.
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
virtual int MinAll(double *PartialMins, double *GlobalMins, int Count) const =0
Epetra_Comm Global Min function.
Epetra_Import: This class builds an import object for efficient importing of off-processor elements...
Definition: Epetra_Import.h:63
virtual Epetra_Distributor * ReverseClone()=0
Create and extract the reverse version of the distributor.
virtual void Barrier() const =0
Epetra_Comm Barrier function.
virtual int MyPID() const =0
Return my process ID.
Epetra_Util: The Epetra Util Wrapper Class.
Definition: Epetra_Util.h:79
bool MyGID(int GID_in) const
Returns true if the GID passed in belongs to the calling processor in this map, otherwise returns fal...
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
void Construct(const Epetra_BlockMap &targetMap, const Epetra_BlockMap &sourceMap, int NumRemotePIDs=-1, const int *UserRemotePIDs=0)
Epetra_Distributor * Distor_
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:57
virtual void Print(std::ostream &os) const
Print object to an output stream.
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
int NumMyElements() const
Number of elements on the calling processor.
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
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 * PermuteToLIDs_
int LID(int GID) const
Returns local ID of global ID, return -1 if not found on this processor.
Epetra_Distributor * Distor_
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 int NumProc() const =0
Returns total number of processes.
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)
bool GlobalIndicesTypeMatch(const Epetra_BlockMap &other) const
virtual Epetra_Distributor * CreateDistributor() const =0
Create a distributor object.
int * PermuteFromLIDs_
int MaxElementSize() const
Maximum element size across all processors.
virtual void Print(std::ostream &os) const
Print object to an output stream Print method.
static void 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)
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
virtual Epetra_Distributor * Clone()=0
Epetra_Distributor clone constructor.
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.
int GID(int LID) const
Returns global ID of local ID, return IndexBase-1 if not found on this processor. ...
Epetra_BlockMap SourceMap_
Epetra_Import(const Epetra_BlockMap &TargetMap, const Epetra_BlockMap &SourceMap)
Constructs a Epetra_Import object from the source and target maps.