Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_BlockMap.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#include "Epetra_ConfigDefs.h"
44#include "Epetra_BlockMap.h"
45#include "Epetra_Comm.h"
46#include "Epetra_Directory.h"
48#include "Epetra_HashTable.h"
49#include <limits>
50
51#ifdef HAVE_MPI
52#include "Epetra_MpiComm.h"
53#endif
54
55// Use the new LID hash table approach by default
56#define EPETRA_BLOCKMAP_NEW_LID
57
58//==============================================================================
59// Epetra_BlockMap constructor function for a Epetra-defined uniform linear distribution of constant size elements.
60void Epetra_BlockMap::ConstructAutoUniform(long long NumGlobal_Elements, int Element_Size, long long Index_Base, const Epetra_Comm& comm, bool IsLongLong)
61{
62
63 // Each processor gets roughly numGlobalPoints/p points
64 // This routine automatically defines a linear partitioning of a
65 // map with numGlobalPoints across the processors
66 // specified in the given Epetra_Comm
67
68 if (NumGlobal_Elements < 0)
69 throw ReportError("NumGlobal_Elements = " + toString(NumGlobal_Elements) + ". Should be >= 0.", -1);
70 if (Element_Size <= 0)
71 throw ReportError("ElementSize = " + toString(Element_Size) + ". Should be > 0.", -2);
72
73 BlockMapData_ = new Epetra_BlockMapData(NumGlobal_Elements, Element_Size, Index_Base, comm, IsLongLong);
74 int NumProc = comm.NumProc();
77
78 int MyPID = comm.MyPID();
79
80 if(BlockMapData_->NumGlobalElements_ / NumProc > (long long) std::numeric_limits<int>::max())
81 throw ReportError("Epetra_BlockMap::ConstructAutoUniform: Error. Not enough space for elements on each processor", -99);
82
84 int remainder = (int) (BlockMapData_->NumGlobalElements_ % NumProc); // remainder will fit int
85 long long start_index = ((long long)(MyPID)) * (BlockMapData_->NumMyElements_ + 1);
86
87 if (MyPID < remainder)
89 else
90 start_index -= (MyPID - remainder);
91
94
99
105
107}
108
109//==============================================================================
110#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
111Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int Element_Size, int Index_Base, const Epetra_Comm& comm)
112 : Epetra_Object("Epetra::BlockMap"),
113 BlockMapData_(0)
114{
115 const bool IsLongLong = true;
116 ConstructAutoUniform(NumGlobal_Elements, Element_Size, static_cast<long long>(Index_Base), comm, IsLongLong);
117}
118
119Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int Element_Size, long long Index_Base, const Epetra_Comm& comm)
120 : Epetra_Object("Epetra::BlockMap"),
121 BlockMapData_(0)
122{
123 const bool IsLongLong = true;
124 ConstructAutoUniform(NumGlobal_Elements, Element_Size, Index_Base, comm, IsLongLong);
125}
126#endif
127
128#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
129Epetra_BlockMap::Epetra_BlockMap(int NumGlobal_Elements, int Element_Size, int Index_Base, const Epetra_Comm& comm)
130 : Epetra_Object("Epetra::BlockMap"),
131 BlockMapData_(0)
132{
133 const bool IsLongLong = false;
134 ConstructAutoUniform((long long)NumGlobal_Elements, Element_Size, Index_Base, comm, IsLongLong);
135}
136#endif
137//==============================================================================
138
139// Epetra_BlockMap constructor function for a user-defined linear distribution of constant size elements.
141 long long NumGlobal_Elements, int NumMy_Elements,
142 int Element_Size, long long Index_Base, const Epetra_Comm& comm, bool IsLongLong)
143{
144 if (NumGlobal_Elements < -1)
145 throw ReportError("NumGlobal_Elements = " + toString(NumGlobal_Elements) + ". Should be >= -1.", -1);
146 if (NumMy_Elements < 0)
147 throw ReportError("NumMy_Elements = " + toString(NumMy_Elements) + ". Should be >= 0.", -2);
148 if (Element_Size <= 0)
149 throw ReportError("ElementSize = " + toString(Element_Size) + ". Should be > 0.", -3);
150
151 BlockMapData_ = new Epetra_BlockMapData(NumGlobal_Elements, Element_Size, Index_Base, comm, IsLongLong);
152 BlockMapData_->NumMyElements_ = NumMy_Elements;
159
160 // Each processor gets NumMyElements points
161
162
163 // Get processor information
164
165 int NumProc = comm.NumProc();
166
167 BlockMapData_->DistributedGlobal_ = IsDistributedGlobal(NumGlobal_Elements, NumMy_Elements);
168
169 // Local Map and uniprocessor case: Each processor gets a complete copy of all elements
170 if (!BlockMapData_->DistributedGlobal_ || NumProc==1) {
172 CheckValidNGE(NumGlobal_Elements);
173
176
177 BlockMapData_-> MinAllGID_ = BlockMapData_->IndexBase_;
181 }
182 else if (NumProc > 1) {
183 // Sum up all local element counts to get global count
184 long long tmp_NumMyElements = BlockMapData_->NumMyElements_;
185 BlockMapData_->Comm_->SumAll(&tmp_NumMyElements, &BlockMapData_->NumGlobalElements_, 1);
186
187 CheckValidNGE(NumGlobal_Elements);
188
191
194
195 // Use the ScanSum function to compute a prefix sum of the number of points
196 long long tmp2_NumMyElements = BlockMapData_->NumMyElements_;
197 BlockMapData_->Comm_->ScanSum(&tmp2_NumMyElements, &BlockMapData_->MaxMyGID_, 1);
198
199 long long start_index = BlockMapData_->MaxMyGID_ - BlockMapData_->NumMyElements_;
202 }
203 else
204 throw ReportError("Internal Error. Report to Epetra developer", -99);
205
206
208}
209
210//==============================================================================
211
212#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
213Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
214 int Element_Size, int Index_Base, const Epetra_Comm& comm)
215 : Epetra_Object("Epetra::BlockMap"),
216 BlockMapData_(0)
217{
218 const bool IsLongLong = true;
219 ConstructUserLinear(NumGlobal_Elements, NumMy_Elements, Element_Size,static_cast<long long>(Index_Base), comm, IsLongLong);
220}
221
222Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
223 int Element_Size, long long Index_Base, const Epetra_Comm& comm)
224 : Epetra_Object("Epetra::BlockMap"),
225 BlockMapData_(0)
226{
227 const bool IsLongLong = true;
228 ConstructUserLinear(NumGlobal_Elements, NumMy_Elements, Element_Size,Index_Base, comm, IsLongLong);
229}
230#endif
231
232#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
233Epetra_BlockMap::Epetra_BlockMap(int NumGlobal_Elements, int NumMy_Elements,
234 int Element_Size, int Index_Base, const Epetra_Comm& comm)
235 : Epetra_Object("Epetra::BlockMap"),
236 BlockMapData_(0)
237{
238 const bool IsLongLong = false;
239 ConstructUserLinear((long long)NumGlobal_Elements, NumMy_Elements, Element_Size,Index_Base, comm, IsLongLong);
240}
241#endif
242
243// Epetra_BlockMap constructor for a user-defined arbitrary distribution of constant size elements.
244template<typename int_type>
245void Epetra_BlockMap::ConstructUserConstant(int_type NumGlobal_Elements, int NumMy_Elements,
246 const int_type * myGlobalElements,
247 int Element_Size, int_type indexBase,
248 const Epetra_Comm& comm, bool IsLongLong)
249{
250 int i;
251 // Each processor gets NumMyElements points
252
253 if (NumGlobal_Elements < -1)
254 throw ReportError("NumGlobal_Elements = " + toString(NumGlobal_Elements) + ". Should be >= -1.", -1);
255 if (NumMy_Elements < 0)
256 throw ReportError("NumMy_Elements = " + toString(NumMy_Elements) + ". Should be >= 0.", -2);
257 if (Element_Size <= 0)
258 throw ReportError("ElementSize = " + toString(Element_Size) + ". Should be > 0.", -3);
259
260 // Allocate storage for global index list information
261
262 BlockMapData_ = new Epetra_BlockMapData(NumGlobal_Elements, Element_Size, indexBase, comm, IsLongLong);
263 if (NumMy_Elements > 0) {
264 int errorcode = SizeMyGlobalElement<int_type>(NumMy_Elements);
265 if(errorcode != 0)
266 throw ReportError("Error with MyGlobalElements allocation.", -99);
267 }
268
269 BlockMapData_->NumMyElements_ = NumMy_Elements;
275 BlockMapData_->LinearMap_ = false;
276 // Get processor information
277
278 int NumProc = comm.NumProc();
279 if (NumMy_Elements > 0) {
280 // Compute min/max GID on this processor
281 BlockMapData_->MinMyGID_ = myGlobalElements[0];
282 BlockMapData_->MaxMyGID_ = myGlobalElements[0];
283 for (i = 0; i < NumMy_Elements; i++) {
284 MyGlobalElementVal<int_type>(i) = myGlobalElements[i];
285 BlockMapData_->MinMyGID_ = EPETRA_MIN(BlockMapData_->MinMyGID_, (long long) myGlobalElements[i]);
286 BlockMapData_->MaxMyGID_ = EPETRA_MAX(BlockMapData_->MaxMyGID_, (long long) myGlobalElements[i]);
287 }
288 }
289 else {
292 }
293
294 BlockMapData_->DistributedGlobal_ = IsDistributedGlobal(NumGlobal_Elements, NumMy_Elements);
295
296 // Local Map and uniprocessor case: Each processor gets a complete copy of all elements
297 if (!BlockMapData_->DistributedGlobal_ || NumProc==1) {
299 CheckValidNGE(NumGlobal_Elements);
302
305 }
306 else if (NumProc > 1) {
307 // Sum up all local element counts to get global count
308 long long tmp_NumMyElements = BlockMapData_->NumMyElements_;
309 BlockMapData_->Comm_->SumAll(&tmp_NumMyElements, &BlockMapData_->NumGlobalElements_, 1);
310 CheckValidNGE(NumGlobal_Elements);
311
314
315 // Use the Allreduce function to find min/max GID
316 long long *tmp_send = new long long[2];
317 long long *tmp_recv = new long long[2];
318 if (BlockMapData_->NumMyElements_ > 0 ||
319 BlockMapData_->NumGlobalElements_ == 0) // This test restores behavior for empty map
320 tmp_send[0] = - BlockMapData_->MinMyGID_; // Negative sign lets us do one reduction
321 else
322 tmp_send[0] = -(std::numeric_limits<int_type>::max())-1;
323 tmp_send[1] = BlockMapData_->MaxMyGID_;
324 BlockMapData_->Comm_->MaxAll(tmp_send, tmp_recv, 2);
325 BlockMapData_->MinAllGID_ = - tmp_recv[0];
326 BlockMapData_->MaxAllGID_ = tmp_recv[1];
327 delete [] tmp_send;
328 delete [] tmp_recv;
330 throw ReportError("Minimum global element index = " + toString(BlockMapData_->MinAllGID_) +
331 " is less than index base = " + toString(BlockMapData_->IndexBase_) +".", -5);
332 }
333 else
334 throw ReportError("Internal Error. Report to Epetra developer", -99);
335
336
338}
339
340#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
341Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
342 const long long * myGlobalElements,
343 int Element_Size, int indexBase,
344 const Epetra_Comm& comm)
345 : Epetra_Object("Epetra::BlockMap"),
346 BlockMapData_(0)
347{
348 const bool IsLongLong = true;
349 ConstructUserConstant(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
350 Element_Size, static_cast<long long>(indexBase), comm, IsLongLong);
351}
352
353Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
354 const long long * myGlobalElements,
355 int Element_Size, long long indexBase,
356 const Epetra_Comm& comm)
357 : Epetra_Object("Epetra::BlockMap"),
358 BlockMapData_(0)
359{
360 const bool IsLongLong = true;
361 ConstructUserConstant(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
362 Element_Size, indexBase, comm, IsLongLong);
363}
364#endif
365
366#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
367Epetra_BlockMap::Epetra_BlockMap(int NumGlobal_Elements, int NumMy_Elements,
368 const int * myGlobalElements,
369 int Element_Size, int indexBase,
370 const Epetra_Comm& comm)
371 : Epetra_Object("Epetra::BlockMap"),
372 BlockMapData_(0)
373{
374 const bool IsLongLong = false;
375 ConstructUserConstant(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
376 Element_Size, indexBase, comm, IsLongLong);
377}
378#endif
379
380
381//==============================================================================
382// Epetra_BlockMap constructor function for a user-defined arbitrary distribution of variable size elements.
383template<typename int_type>
384void Epetra_BlockMap::ConstructUserVariable(int_type NumGlobal_Elements, int NumMy_Elements,
385 const int_type * myGlobalElements,
386 const int *elementSizeList, int_type indexBase,
387 const Epetra_Comm& comm, bool IsLongLong)
388{
389
390 int i;
391 // Each processor gets NumMyElements points
392
393 if (NumGlobal_Elements < -1)
394 throw ReportError("NumGlobal_Elements = " + toString(NumGlobal_Elements) + ". Should be >= -1.", -1);
395 if (NumMy_Elements < 0)
396 throw ReportError("NumMy_Elements = " + toString(NumMy_Elements) + ". Should be >= 0.", -2);
397 for (i = 0; i < NumMy_Elements; i++)
398 if (elementSizeList[i] <= 0)
399 throw ReportError("elementSizeList["+toString(i)+"] = " + toString(elementSizeList[i]) + ". Should be > 0.", -3);
400
401 BlockMapData_ = new Epetra_BlockMapData(NumGlobal_Elements, 0, indexBase, comm, IsLongLong);
402 BlockMapData_->NumMyElements_ = NumMy_Elements;
404 BlockMapData_->LinearMap_ = false;
405 // Allocate storage for global index list and element size information
406
407 if (NumMy_Elements > 0) {
408 int errorcode = SizeMyGlobalElement<int_type>(NumMy_Elements);
409 if(errorcode != 0)
410 throw ReportError("Error with MyGlobalElements allocation.", -99);
411 errorcode = BlockMapData_->ElementSizeList_.Size(NumMy_Elements);
412 if(errorcode != 0)
413 throw ReportError("Error with ElementSizeList allocation.", -99);
414 }
415 // Get processor information
416
417 int NumProc = comm.NumProc();
418
419 if (NumMy_Elements > 0) {
420 // Compute min/max GID and element size, number of points on this processor
421 BlockMapData_->MinMyGID_ = myGlobalElements[0];
422 BlockMapData_->MaxMyGID_ = myGlobalElements[0];
423 BlockMapData_->MinMyElementSize_ = elementSizeList[0];
424 BlockMapData_->MaxMyElementSize_ = elementSizeList[0];
426 for (i = 0; i < NumMy_Elements; i++) {
427 MyGlobalElementVal<int_type>(i) = myGlobalElements[i];
428 BlockMapData_->ElementSizeList_[i] = elementSizeList[i];
429 BlockMapData_->MinMyGID_ = EPETRA_MIN(BlockMapData_->MinMyGID_,(long long) myGlobalElements[i]);
430 BlockMapData_->MaxMyGID_ = EPETRA_MAX(BlockMapData_->MaxMyGID_,(long long) myGlobalElements[i]);
433 BlockMapData_->NumMyPoints_ += elementSizeList[i];
434 }
435 }
436 else {
442 }
443
444 BlockMapData_->DistributedGlobal_ = IsDistributedGlobal(NumGlobal_Elements, NumMy_Elements);
445
446 // Local Map and uniprocessor case: Each processor gets a complete copy of all elements
447 if (!BlockMapData_->DistributedGlobal_ || NumProc == 1) {
449 CheckValidNGE(NumGlobal_Elements);
451
456 }
457 else if (NumProc > 1) {
458 // Sum up all local element and point counts to get global counts
459 int_type *tmp_send = new int_type[4];
460 int_type *tmp_recv = new int_type[4];
461 tmp_send[0] = BlockMapData_->NumMyElements_;
462 tmp_send[1] = BlockMapData_->NumMyPoints_;
463 BlockMapData_->Comm_->SumAll(tmp_send, tmp_recv, 2);
464 BlockMapData_->NumGlobalElements_ = tmp_recv[0];
465 BlockMapData_->NumGlobalPoints_ = tmp_recv[1];
466
467 CheckValidNGE(NumGlobal_Elements);
468
469 // Use the MaxAll function to find min/max GID
470 if (BlockMapData_->NumMyElements_ > 0 ||
471 BlockMapData_->NumGlobalElements_ == 0) // This test restores behavior for empty map
472 tmp_send[0] = - static_cast<int_type>(BlockMapData_->MinMyGID_); // Negative sign lets us do one reduction
473 else
474 tmp_send[0] = -(std::numeric_limits<int_type>::max())-1;
475 tmp_send[1] = static_cast<int_type>(BlockMapData_->MaxMyGID_);
476 tmp_send[2] = - BlockMapData_->MinMyElementSize_;
478 tmp_send[2] = - static_cast<int_type>(BlockMapData_->NumGlobalPoints_); // This processor has no elements, so should not sizes.
479 tmp_send[3] = BlockMapData_->MaxMyElementSize_;
480
481 BlockMapData_->Comm_->MaxAll(tmp_send, tmp_recv, 4);
482
483 BlockMapData_->MinAllGID_ = - tmp_recv[0];
484 BlockMapData_->MaxAllGID_ = tmp_recv[1];
485 BlockMapData_->MinElementSize_ = - (int) tmp_recv[2];
486 BlockMapData_->MaxElementSize_ = (int) tmp_recv[3];
487
488 delete [] tmp_send;
489 delete [] tmp_recv;
490
491 // Check for constant element size
495 }
496
498 throw ReportError("Minimum global element index = " + toString(BlockMapData_->MinAllGID_) +
499 " is less than index base = " + toString(BlockMapData_->IndexBase_) +".", -5);
500 }
501 else
502 throw ReportError("Internal Error. Report to Epetra developer", -99);
503
504
506}
507
508#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
509Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
510 const long long * myGlobalElements,
511 const int *elementSizeList, int indexBase,
512 const Epetra_Comm& comm)
513 : Epetra_Object("Epetra::BlockMap"),
514 BlockMapData_(0)
515{
516 const bool IsLongLong = true;
517 ConstructUserVariable(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
518 elementSizeList, static_cast<long long>(indexBase), comm, IsLongLong);
519}
520
521Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
522 const long long * myGlobalElements,
523 const int *elementSizeList, long long indexBase,
524 const Epetra_Comm& comm)
525 : Epetra_Object("Epetra::BlockMap"),
526 BlockMapData_(0)
527{
528 const bool IsLongLong = true;
529 ConstructUserVariable(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
530 elementSizeList, indexBase, comm, IsLongLong);
531}
532#endif
533
534#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
535Epetra_BlockMap::Epetra_BlockMap(int NumGlobal_Elements, int NumMy_Elements,
536 const int * myGlobalElements,
537 const int *elementSizeList, int indexBase,
538 const Epetra_Comm& comm)
539 : Epetra_Object("Epetra::BlockMap"),
540 BlockMapData_(0)
541{
542 const bool IsLongLong = false;
543 ConstructUserVariable(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
544 elementSizeList, indexBase, comm, IsLongLong);
545}
546#endif
547
548
549//==============================================================================
550// Epetra_BlockMap constructor function for a user-defined arbitrary distribution of variable size elements,
551// with all the information on globals provided by the user.
552template<typename int_type>
553void Epetra_BlockMap::ConstructUserConstantNoComm(int_type NumGlobal_Elements, int NumMy_Elements,
554 const int_type * myGlobalElements,
555 int Element_Size, int_type indexBase,
556 const Epetra_Comm& comm, bool IsLongLong,
557 bool UserIsDistributedGlobal,
558 int_type UserMinAllGID, int_type UserMaxAllGID)
559{
560
561
562 int i;
563 // Each processor gets NumMyElements points
564
565 if (NumGlobal_Elements < -1)
566 throw ReportError("NumGlobal_Elements = " + toString(NumGlobal_Elements) + ". Should be >= -1.", -1);
567 if (NumMy_Elements < 0)
568 throw ReportError("NumMy_Elements = " + toString(NumMy_Elements) + ". Should be >= 0.", -2);
569 if (Element_Size <= 0)
570 throw ReportError("ElementSize = " + toString(Element_Size) + ". Should be > 0.", -3);
571
572 // Allocate storage for global index list information
573
574 BlockMapData_ = new Epetra_BlockMapData(NumGlobal_Elements, Element_Size, indexBase, comm, IsLongLong);
575 if (NumMy_Elements > 0) {
576 int errorcode = SizeMyGlobalElement<int_type>(NumMy_Elements);
577 if(errorcode != 0)
578 throw ReportError("Error with MyGlobalElements allocation.", -99);
579 }
580
581 BlockMapData_->NumMyElements_ = NumMy_Elements;
587 BlockMapData_->LinearMap_ = false;
588 // Get processor information
589
590 int NumProc = comm.NumProc();
591 if (NumMy_Elements > 0) {
592 // Compute min/max GID on this processor
593 BlockMapData_->MinMyGID_ = myGlobalElements[0];
594 BlockMapData_->MaxMyGID_ = myGlobalElements[0];
595 for (i = 0; i < NumMy_Elements; i++) {
596 MyGlobalElementVal<int_type>(i) = myGlobalElements[i];
597 BlockMapData_->MinMyGID_ = EPETRA_MIN(BlockMapData_->MinMyGID_, (long long) myGlobalElements[i]);
598 BlockMapData_->MaxMyGID_ = EPETRA_MAX(BlockMapData_->MaxMyGID_, (long long) myGlobalElements[i]);
599 }
600 }
601 else {
604 }
605
606 BlockMapData_->DistributedGlobal_ = UserIsDistributedGlobal;
607
608 // Local Map and uniprocessor case: Each processor gets a complete copy of all elements
609 if (!BlockMapData_->DistributedGlobal_ || NumProc==1) {
611 CheckValidNGE(NumGlobal_Elements);
614
617 }
618 else if (NumProc > 1) {
619 if(NumGlobal_Elements==-1) {
620 // Sum up all local element counts to get global count
621 long long tmp_NumMyElements = BlockMapData_->NumMyElements_;
622 BlockMapData_->Comm_->SumAll(&tmp_NumMyElements, &BlockMapData_->NumGlobalElements_, 1);
623 }
624 else {
625 // User provides this information
626 BlockMapData_->NumGlobalElements_ = NumGlobal_Elements;
627 }
629
632
633 BlockMapData_->MinAllGID_ = UserMinAllGID;
634 BlockMapData_->MaxAllGID_ = UserMaxAllGID;
636 throw ReportError("Minimum global element index = " + toString(BlockMapData_->MinAllGID_) +
637 " is less than index base = " + toString(BlockMapData_->IndexBase_) +".", -5);
638 }
639 else
640 throw ReportError("Internal Error. Report to Epetra developer", -99);
641
642
644}
645
646#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
647Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
648 const long long * myGlobalElements,
649 int theElementSize, int indexBase,
650 const Epetra_Comm& comm,
651 bool UserIsDistributedGlobal,
652 long long UserMinAllGID, long long UserMaxAllGID)
653 : Epetra_Object("Epetra::BlockMap"),
654 BlockMapData_(0)
655{
656 const bool IsLongLong = true;
657 ConstructUserConstantNoComm(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
658 theElementSize, (long long) indexBase, comm, IsLongLong,
659 UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID);
660}
661Epetra_BlockMap::Epetra_BlockMap(long long NumGlobal_Elements, int NumMy_Elements,
662 const long long * myGlobalElements,
663 int theElementSize, long long indexBase,
664 const Epetra_Comm& comm,
665 bool UserIsDistributedGlobal,
666 long long UserMinAllGID, long long UserMaxAllGID)
667 : Epetra_Object("Epetra::BlockMap"),
668 BlockMapData_(0)
669{
670 const bool IsLongLong = true;
671 ConstructUserConstantNoComm(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
672 theElementSize, indexBase, comm, IsLongLong,
673 UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID);
674}
675#endif
676
677#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
678Epetra_BlockMap::Epetra_BlockMap(int NumGlobal_Elements, int NumMy_Elements,
679 const int * myGlobalElements,
680 int theElementSize, int indexBase,
681 const Epetra_Comm& comm,
682 bool UserIsDistributedGlobal,
683 int UserMinAllGID, int UserMaxAllGID)
684 : Epetra_Object("Epetra::BlockMap"),
685 BlockMapData_(0)
686{
687 const bool IsLongLong = false;
688 ConstructUserConstantNoComm(NumGlobal_Elements, NumMy_Elements, myGlobalElements,
689 theElementSize, indexBase, comm, IsLongLong,
690 UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID);
691}
692#endif
693
694
695
696
697//==============================================================================
699 : Epetra_Object(map.Label()),
700 BlockMapData_(map.BlockMapData_)
701{
703
704 // This call appears to be unnecessary overhead. Removed 10-Aug-2004 maherou.
705 // GlobalToLocalSetup(); // Setup any information for making global index to local index translation fast.
706}
707
708//==============================================================================
710
711 // Quickest test: See if both maps share an inner data class
712 if (this->BlockMapData_ == Map.BlockMapData_)
713 return(true);
714 return false;
715}
716
717//==============================================================================
719
720 // Quickest test: See if both maps share an inner data class
721 if (this->BlockMapData_ == Map.BlockMapData_)
722 return(true);
723
724 if(!GlobalIndicesTypeMatch(Map))
725 return(false);
726
727 // Next check other global properties that are easy global attributes
728 if (BlockMapData_->MinAllGID_ != Map.MinAllGID64() ||
732 return(false);
733
734 // Last possible global check for constant element sizes
736 return(false);
737
738 // If we get this far, we need to check local properties and then check across
739 // all processors to see if local properties are all true
740
741 int numMyElements = BlockMapData_->NumMyElements_;
742
743 int MySameMap = 1; // Assume not needed
744
745 // First check if number of element is the same in each map
746 if (numMyElements != Map.NumMyElements()) MySameMap = 0;
747
748 // If numMyElements is the same, check to see that list of GIDs is the same
749 if (MySameMap==1) {
750 if (LinearMap() && Map.LinearMap() ) {
751 // For linear maps, just need to check whether lower bound is the same
752 if (MinMyGID64() != Map.MinMyGID64() )
753 MySameMap = 0;
754 }
755 else {
756 for (int i = 0; i < numMyElements; i++) {
757 if (GID64(i) != Map.GID64(i)) {
758 MySameMap = 0;
759 break;
760 }
761 }
762 }
763 }
764// for (int i = 0; i < numMyElements; i++)
765// if (GID64(i) != Map.GID64(i)) MySameMap = 0;
766
767 // If GIDs are the same, check to see element sizes are the same
768 if (MySameMap==1 && !BlockMapData_->ConstantElementSize_) {
769 int * sizeList1 = ElementSizeList();
770 int * sizeList2 = Map.ElementSizeList();
771 for (int i = 0; i < numMyElements; i++) if (sizeList1[i] != sizeList2[i]) MySameMap=0;
772 }
773 // Now get min of MySameMap across all processors
774
775 int GlobalSameMap = 0;
776#ifdef NDEBUG
777 (void) Comm().MinAll(&MySameMap, &GlobalSameMap, 1);
778#else
779 // Don't declare 'err' unless we actually use it (in the assert(),
780 // which gets defined away in a release build).
781 int err = Comm().MinAll(&MySameMap, &GlobalSameMap, 1);
782 assert(err==0);
783#endif // NDEBUG
784 return(GlobalSameMap==1);
785}
786
787//==============================================================================
789{
790 if (this->BlockMapData_ == Map.BlockMapData_)
791 return(true);
792
793 if(!GlobalIndicesTypeMatch(Map))
794 return(false);
795
797 return(false);
798
799 // If we get this far, we need to check local properties and then check across
800 // all processors to see if local properties are all true
801
802 int MySameMap = 1; // Assume not needed
804 MySameMap = 0;
805
806 int GlobalSameMap = 0;
807
808#ifdef NDEBUG
809 (void) Comm().MinAll(&MySameMap, &GlobalSameMap, 1);
810#else
811 // Don't declare 'err' unless we actually use it (in the assert(),
812 // which gets defined away in a release build).
813 int err = Comm().MinAll(&MySameMap, &GlobalSameMap, 1);
814 assert( err == 0 );
815#endif // NDEBUG
816
817 return(GlobalSameMap==1);
818}
819
820//==============================================================================
821#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
822int Epetra_BlockMap::MyGlobalElements(long long * myGlobalElements) const
823{
824 // Although one can populate long long data from int data, we don't
825 // allow it to maintain int/long long symmetry.
827 throw ReportError("Epetra_BlockMap::MyGlobalElements(long long *) ERROR, Can't call for non long long* map.",-1);
828
829 // If the global element list is not create, then do so. This can only happen when
830 // a linear distribution has been specified. Thus we can easily construct the update
831 // list in this case.
832
833 int i;
834 int numMyElements = BlockMapData_->NumMyElements_;
835
837 for (i = 0; i < numMyElements; i++)
838 myGlobalElements[i] = BlockMapData_->MinMyGID_ + i;
839 else
840 for (i = 0; i < numMyElements; i++)
841 myGlobalElements[i] = BlockMapData_->MyGlobalElements_LL_[i];
842 return(0);
843}
844#endif
845
846//==============================================================================
847#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
848int Epetra_BlockMap::MyGlobalElements(int * myGlobalElements) const
849{
851 throw ReportError("Epetra_BlockMap::MyGlobalElements(int *) ERROR, Can't call for non int* map.",-1);
852
853 // If the global element list is not create, then do so. This can only happen when
854 // a linear distribution has been specified. Thus we can easily construct the update
855 // list in this case.
856
857 int i;
858 int numMyElements = BlockMapData_->NumMyElements_;
859
861 for (i = 0; i < numMyElements; i++)
862 myGlobalElements[i] = (int) BlockMapData_->MinMyGID_ + i;
863 else
864 for (i = 0; i < numMyElements; i++)
865 myGlobalElements[i] = (int) BlockMapData_->MyGlobalElements_int_[i];
866 return(0);
867}
868#endif
869
870#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
871int Epetra_BlockMap::MyGlobalElementsPtr(long long *& MyGlobalElementList) const
872{
873 MyGlobalElementList = MyGlobalElements64();
874 return(0);
875}
876#endif
877
878#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
879int Epetra_BlockMap::MyGlobalElementsPtr(int *& MyGlobalElementList) const
880{
881 MyGlobalElementList = MyGlobalElements();
882 return(0);
883}
884#endif
885//==============================================================================
886#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
889 throw ReportError("Epetra_BlockMap::MyGlobalElements() ERROR, Can't call for non int* map.",-1);
890
891 int numMyElements = BlockMapData_->NumMyElements_;
892
893 // If ElementSizeList not built, do so
894 if(BlockMapData_->MyGlobalElements_int_.Length() == 0 && numMyElements > 0) {
895 int errorcode = BlockMapData_->MyGlobalElements_int_.Size(numMyElements + 1);
896 if(errorcode != 0)
897 throw ReportError("Error with MyGlobalElements allocation.", -99);
898
899 // Build the array
900 for (int i = 0; i < numMyElements; i++)
902 }
904}
905#endif
906//==============================================================================
907#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
910 throw ReportError("Epetra_BlockMap::MyGlobalElements64 ERROR, Can't call for non long long* map.",-1);
911
912 int numMyElements = BlockMapData_->NumMyElements_;
913
914 // If ElementSizeList not built, do so
915 if(BlockMapData_->MyGlobalElements_LL_.Length() == 0 && numMyElements > 0) {
916 int errorcode = BlockMapData_->MyGlobalElements_LL_.Size(numMyElements + 1);
917 if(errorcode != 0)
918 throw ReportError("Error with MyGlobalElements allocation.", -99);
919
920 // Build the array
921 for (int i = 0; i < numMyElements; i++)
923 }
925}
926#endif
927//==============================================================================
929{
930 if (!MyLID(lid))
931 EPETRA_CHK_ERR(-1);
932
933 int entry;
934
936 entry = MaxElementSize() * lid; // convert to vector entry
937 else {
938 int * entrylist = FirstPointInElementList(); // get entry list
939 entry = entrylist[lid];
940 }
941 return(entry);
942}
943
944//==============================================================================
945int Epetra_BlockMap::FirstPointInElementList(int * firstPointInElementList) const
946{
947 // If the first element entry list is not create, then do so.
948
949 // Note: This array is of length NumMyElement+1
950
951 int i;
952 int numMyElements = BlockMapData_->NumMyElements_;
953
955 firstPointInElementList[0] = 0; // First element of first entry is always zero
956
958 for (i = 0; i < numMyElements; i++)
959 firstPointInElementList[i+1] = firstPointInElementList[i] + BlockMapData_->ElementSize_;
960 else
961 for (i = 0; i < numMyElements; i++)
962 firstPointInElementList[i+1] = firstPointInElementList[i] + BlockMapData_->ElementSizeList_[i];
963 }
964 else
965 for (i = 0; i <= numMyElements; i++)
966 firstPointInElementList[i] = BlockMapData_->FirstPointInElementList_[i];
967 return(0);
968}
969
970//==============================================================================
972 int numMyElements = BlockMapData_->NumMyElements_;
973
974 // If ElementSizeList not built, do so
975 if ((BlockMapData_->FirstPointInElementList_.Length() == 0) && (numMyElements > 0)) {
977 BlockMapData_->FirstPointInElementList_[0] = 0; // First element of first entry is always zero
979 for (int i = 0; i < numMyElements; i++)
981 else
982 for (int i = 0; i < numMyElements; i++)
984 }
986}
987
988//==============================================================================
989int Epetra_BlockMap::ElementSizeList(int * elementSizeList) const
990{
991 // If the element size list is not create, then do so. This can only happen when
992 // a constant element size has been specified. Thus we can easily construct the element size
993 // list in this case.
994
995 int i;
996 int numMyElements = BlockMapData_->NumMyElements_;
997
999 for (i = 0; i < numMyElements; i++)
1000 elementSizeList[i] = BlockMapData_->ElementSize_;
1001 else
1002 for (i = 0; i < numMyElements; i++)
1003 elementSizeList[i] = BlockMapData_->ElementSizeList_[i];
1004
1005 return(0);
1006}
1007
1008//==============================================================================
1010 int numMyElements = BlockMapData_->NumMyElements_;
1011
1012 // If ElementSizeList not built, do so
1013 if ((BlockMapData_->ElementSizeList_.Length() == 0) && (numMyElements > 0)) {
1014 BlockMapData_->ElementSizeList_.Size(numMyElements);
1015 for (int i = 0; i < numMyElements; i++)
1017 }
1019}
1020
1021//==============================================================================
1022int Epetra_BlockMap::PointToElementList(int * pointToElementList) const {
1023 // Build an array such that the local element ID is stored for each point
1024
1025 int i;
1027 int numMyElements = BlockMapData_->NumMyElements_;
1028 int * ptr = pointToElementList;
1029 for (i = 0; i < numMyElements; i++) {
1030 int Size = ElementSize(i);
1031 for (int j = 0; j < Size; j++)
1032 *ptr++ = i;
1033 }
1034 }
1035 else {
1036 int numMyPoints = BlockMapData_->NumMyPoints_;
1037 for (i = 0; i < numMyPoints; i++)
1038 pointToElementList[i] = BlockMapData_->PointToElementList_[i];
1039 }
1040 return(0);
1041}
1042
1043//==============================================================================
1045
1046 // If PointToElementList not built, do so
1049 int numMyElements = BlockMapData_->NumMyElements_;
1051 for (int i = 0; i < numMyElements; i++) {
1052 int Size = ElementSize(i);
1053 for (int j = 0; j < Size; j++)
1054 *ptr++ = i;
1055 }
1056 }
1058}
1059
1060//==============================================================================
1062
1063 if (ConstantElementSize())
1064 return(BlockMapData_->ElementSize_);
1065 else
1066 return(BlockMapData_->ElementSizeList_[lid]);
1067}
1068
1069//==============================================================================
1074 }
1075 return(BlockMapData_->OneToOne_);
1076}
1077
1078//==============================================================================
1079template<typename int_type>
1081{
1082 int i;
1083 int numMyElements = BlockMapData_->NumMyElements_;
1084
1086 return; // Nothing to do
1087 }
1088
1089 if (LinearMap() || numMyElements == 0) {
1090 return; // Nothing else to do
1091 }
1092
1093 // Build LID_ vector to make look up of local index values fast
1094
1095#ifdef EPETRA_BLOCKMAP_NEW_LID
1096
1097 // Here follows an optimization that checks for an initial block of
1098 // contiguous GIDs, and stores the GID->LID mapping for those in a
1099 // more efficient way. This supports a common use case for
1100 // overlapping Maps, in which owned entries of a vector are ordered
1101 // before halo entries.
1102 //
1103 // Epetra defines EPETRA_BLOCKMAP_NEW_LID by default (see the top of
1104 // this file).
1105
1106 //check for initial contiguous block
1107 int_type val = MyGlobalElementValGet<int_type>(0);
1108 for( i = 0 ; i < numMyElements; ++i ) {
1109 if (val != MyGlobalElementValGet<int_type>(i)) break;
1110 ++val;
1111 }
1114 BlockMapData_->LastContiguousGID_ = MyGlobalElementValGet<int_type>(0);
1115 }
1116 else {
1118 MyGlobalElementValGet<int_type>(BlockMapData_->LastContiguousGIDLoc_);
1119 }
1120
1121 //Hash everything else
1122 if(i < numMyElements) {
1123 if (BlockMapData_->LIDHash_ != NULL) {
1124 delete BlockMapData_->LIDHash_;
1125 }
1126
1127 BlockMapData_->LIDHash_ = new Epetra_HashTable<int>(numMyElements - i + 1 );
1128 for(; i < numMyElements; ++i )
1129 BlockMapData_->LIDHash_->Add( MyGlobalElementValGet<int_type>(i), i );
1130 }
1131
1132#else
1133
1134 int SpanGID = BlockMapData_->MaxMyGID_ - BlockMapData_->MinMyGID_ + 1;
1135 BlockMapData_->LID_.Size(SpanGID);
1136
1137 for (i = 0; i < SpanGID; i++)
1138 BlockMapData_->LID_[i] = -1; // Fill all locations with -1
1139
1140 for (i = 0; i < numMyElements; i++) {
1141 int tmp = MyGlobalElementValGet<int_type>(i) - BlockMapData_->MinMyGID_;
1142 assert(tmp >= 0);
1143 assert(tmp < SpanGID);
1144 BlockMapData_->LID_[MyGlobalElementValGet<int_type>(i) - BlockMapData_->MinMyGID_] = i; // Spread local indices
1145 }
1146
1147#endif
1148
1149}
1150
1152{
1154 {
1155#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1156 TGlobalToLocalSetup<int>();
1157#else
1158 throw ReportError("Epetra_BlockMap::GlobalToLocalSetup ERROR, GlobalIndices int but no API for it.",-1);
1159#endif
1160 }
1162 {
1163#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1164 TGlobalToLocalSetup<long long>();
1165#else
1166 throw ReportError("Epetra_BlockMap::GlobalToLocalSetup ERROR, GlobalIndices long long but no API for it.",-1);
1167#endif
1168 }
1169 else
1170 {
1171 throw ReportError("Epetra_BlockMap::GlobalToLocalSetup ERROR, GlobalIndices type unknown.",-1);
1172 }
1173}
1174
1175//==============================================================================
1176#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1177int Epetra_BlockMap::LID(long long gid) const
1178{
1179 if ((gid < BlockMapData_->MinMyGID_) ||
1180 (gid > BlockMapData_->MaxMyGID_)) {
1181 return(-1); // Out of range
1182 }
1183
1185 return (int) (gid - BlockMapData_->MinMyGID_); // Can compute with an offset
1186 }
1187
1188#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1190 if( (int) gid >= BlockMapData_->MyGlobalElements_int_[0] &&
1191 (int) gid <= BlockMapData_->LastContiguousGID_ ) {
1192 return (int) gid - BlockMapData_->MyGlobalElements_int_[0];
1193 }
1194 }
1195 else
1196#endif
1197#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1199 if( gid >= BlockMapData_->MyGlobalElements_LL_[0] &&
1200 gid <= BlockMapData_->LastContiguousGID_ ) {
1201 return (int) ( gid - BlockMapData_->MyGlobalElements_LL_[0] );
1202 }
1203 }
1204 else
1205#endif
1206 throw ReportError("Epetra_BlockMap::LID ERROR, GlobalIndices type unknown.",-1);
1207
1208#ifdef EPETRA_BLOCKMAP_NEW_LID
1209 return BlockMapData_->LIDHash_->Get( gid );
1210#else
1211 return(BlockMapData_->LID_[gid - BlockMapData_->MinMyGID_]); // Find it in LID array
1212#endif
1213}
1214#endif
1215
1216//==============================================================================
1217#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1218int Epetra_BlockMap::LID(int gid) const
1219{
1220 if ((gid < (int) BlockMapData_->MinMyGID_) ||
1221 (gid > (int) BlockMapData_->MaxMyGID_)) {
1222 return(-1); // Out of range
1223 }
1224
1226 return (int) (gid - (int) BlockMapData_->MinMyGID_); // Can compute with an offset
1227 }
1228
1230 if( gid >= BlockMapData_->MyGlobalElements_int_[0] &&
1231 gid <= (int) BlockMapData_->LastContiguousGID_ ) {
1232 return (int) ( gid - BlockMapData_->MyGlobalElements_int_[0] );
1233 }
1234 }
1236 throw ReportError("Epetra_BlockMap::LID ERROR, int version called for long long map.",-1);
1237 }
1238 else {
1239 throw ReportError("Epetra_BlockMap::LID ERROR, GlobalIndices type unknown.",-1);
1240 }
1241
1242#ifdef EPETRA_BLOCKMAP_NEW_LID
1243 return BlockMapData_->LIDHash_->Get( gid );
1244#else
1245 return(BlockMapData_->LID_[gid - BlockMapData_->MinMyGID_]); // Find it in LID array
1246#endif
1247}
1248#endif
1249
1250//==============================================================================
1251
1252long long Epetra_BlockMap::GID64(int lid) const
1253{
1254 if ((BlockMapData_->NumMyElements_==0) ||
1255 (lid < BlockMapData_->MinLID_) ||
1256 (lid > BlockMapData_->MaxLID_)) {
1257 return(BlockMapData_->IndexBase_ - 1); // Out of range
1258 }
1259
1260 if (LinearMap()) {
1261 return(lid + BlockMapData_->MinMyGID_); // Can compute with an offset
1262 }
1263
1264#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1266 {
1267 return(BlockMapData_->MyGlobalElements_int_[lid]); // Find it in MyGlobalElements array
1268 }
1269 else
1270#endif
1271#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1273 {
1274 return(BlockMapData_->MyGlobalElements_LL_[lid]); // Find it in MyGlobalElements array
1275 }
1276 else
1277#endif
1278 throw ReportError("Epetra_BlockMap::GID64 ERROR, GlobalIndices type unknown.",-1);
1279}
1280
1281#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1282int Epetra_BlockMap::GID(int lid) const
1283{
1285 {
1286 if ((BlockMapData_->NumMyElements_==0) ||
1287 (lid < BlockMapData_->MinLID_) ||
1288 (lid > BlockMapData_->MaxLID_)) {
1289 return((int) BlockMapData_->IndexBase_ - 1); // Out of range
1290 }
1291
1292 if (LinearMap()) {
1293 return(lid + (int) BlockMapData_->MinMyGID_); // Can compute with an offset
1294 }
1295
1296 return(BlockMapData_->MyGlobalElements_int_[lid]); // Find it in MyGlobalElements array
1297 }
1298
1299 throw ReportError("Epetra_BlockMap::GID ERROR, GlobalIndices type unknown or long long.",-1);
1300}
1301#endif
1302
1303//==============================================================================
1304int Epetra_BlockMap::FindLocalElementID(int PointID, int & ElementID, int & ElementOffset) const {
1305
1306 if (PointID >= BlockMapData_->NumMyPoints_)
1307 return(-1); // Point is out of range
1308
1309 if (ConstantElementSize()) {
1310 ElementID = PointID / BlockMapData_->MaxElementSize_;
1311 ElementOffset = PointID % BlockMapData_->MaxElementSize_;
1312 return(0);
1313 }
1314 else {
1315 int * tmpPointToElementList = PointToElementList();
1316 int * tmpFirstPointInElementList = FirstPointInElementList();
1317 ElementID = tmpPointToElementList[PointID];
1318 ElementOffset = PointID - tmpFirstPointInElementList[ElementID];
1319 return(0);
1320 }
1321}
1322
1323//==============================================================================
1324#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1325int Epetra_BlockMap::RemoteIDList(int NumIDs, const int * GIDList,
1326 int * PIDList, int * LIDList,
1327 int * SizeList) const
1328{
1330 throw ReportError("Epetra_BlockMap::RemoteIDList ERROR, Can't call int* version for non int* map.",-1);
1331
1332 if (BlockMapData_->Directory_ == NULL) {
1334 }
1335
1337 if (directory == NULL) {
1338 return(-1);
1339 }
1340
1341 EPETRA_CHK_ERR( directory->GetDirectoryEntries(*this, NumIDs, GIDList,
1342 PIDList, LIDList, SizeList) );
1343
1344 return(0);
1345}
1346#endif
1347//==============================================================================
1348#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1349int Epetra_BlockMap::RemoteIDList(int NumIDs, const long long * GIDList,
1350 int * PIDList, int * LIDList,
1351 int * SizeList) const
1352{
1354 throw ReportError("Epetra_BlockMap::RemoteIDList ERROR, Can't call long long* version for non long long* map.",-1);
1355
1356 if (BlockMapData_->Directory_ == NULL) {
1358 }
1359
1361 if (directory == NULL) {
1362 return(-1);
1363 }
1364
1365 EPETRA_CHK_ERR( directory->GetDirectoryEntries(*this, NumIDs, GIDList,
1366 PIDList, LIDList, SizeList) );
1367
1368 return(0);
1369}
1370#endif
1371
1372//==============================================================================
1374{
1375 if (Comm().NumProc() < 2) {
1376 return(true);
1377 }
1378
1379 if (BlockMapData_->Directory_ == NULL) {
1381 }
1382
1384 if (directory == NULL) {
1385 throw ReportError("Epetra_BlockMap::IsOneToOne ERROR, CreateDirectory failed.",-1);
1386 }
1387
1388 return(directory->GIDsAllUniquelyOwned());
1389}
1390
1391//==============================================================================
1392bool Epetra_BlockMap::IsDistributedGlobal(long long numGlobalElements, int numMyElements) const {
1393
1394 bool isDistributedGlobal = false; // Assume map is not global distributed
1395 if (BlockMapData_->Comm_->NumProc() > 1) {
1396 int LocalReplicated = 0;
1397 int AllLocalReplicated;
1398 if (numGlobalElements == numMyElements)
1399 LocalReplicated=1;
1400 BlockMapData_->Comm_->MinAll(&LocalReplicated, &AllLocalReplicated, 1);
1401
1402 // If any PE has LocalReplicated=0, then map is distributed global
1403 if (AllLocalReplicated != 1)
1404 isDistributedGlobal = true;
1405 }
1406 return(isDistributedGlobal);
1407}
1408
1409//==============================================================================
1410void Epetra_BlockMap::CheckValidNGE(long long numGlobalElements) {
1411 // Check to see if user's value for numGlobalElements is either -1
1412 // (in which case we use our computed value) or matches ours.
1413 if ((numGlobalElements != -1) && (numGlobalElements != BlockMapData_->NumGlobalElements_)) {
1414 long long BmdNumGlobalElements = BlockMapData_->NumGlobalElements_;
1415 CleanupData();
1416 throw ReportError("Invalid NumGlobalElements. NumGlobalElements = " + toString(numGlobalElements) +
1417 ". Should equal " + toString(BmdNumGlobalElements) +
1418 ", or be set to -1 to compute automatically", -4);
1419 }
1420}
1421
1422//==============================================================================
1426
1427 GlobalToLocalSetup(); // Setup any information for making global index to local index translation fast.
1428}
1429
1430//==============================================================================
1431void Epetra_BlockMap::Print(std::ostream & os) const
1432{
1433 int * FirstPointInElementList1 = 0;
1434 int * ElementSizeList1 = 0;
1435 if (!ConstantElementSize()) {
1436 FirstPointInElementList1 = FirstPointInElementList();
1437 ElementSizeList1 = ElementSizeList();
1438 }
1439 int MyPID = Comm().MyPID();
1440 int NumProc = Comm().NumProc();
1441
1442 for (int iproc = 0; iproc < NumProc; iproc++) {
1443 if (MyPID == iproc) {
1444 if (MyPID == 0) {
1445 os << "\nNumber of Global Elements = "; os << NumGlobalElements64(); os << std::endl;
1446 os << "Number of Global Points = "; os << NumGlobalPoints64(); os << std::endl;
1447 os << "Maximum of all GIDs = "; os << MaxAllGID64(); os << std::endl;
1448 os << "Minimum of all GIDs = "; os << MinAllGID64(); os << std::endl;
1449 os << "Index Base = "; os << IndexBase64(); os << std::endl;
1450 if (ConstantElementSize()) {
1451 os << "Constant Element Size = "; os << ElementSize(); os << std::endl;
1452 }
1453 }
1454 os << std::endl;
1455
1456 os << "Number of Local Elements = "; os << NumMyElements(); os << std::endl;
1457 os << "Number of Local Points = "; os << NumMyPoints(); os << std::endl;
1458 os << "Maximum of my GIDs = "; os << MaxMyGID64(); os << std::endl;
1459 os << "Minimum of my GIDs = "; os << MinMyGID64(); os << std::endl;
1460 os << std::endl;
1461
1462 os.width(14);
1463 os << " MyPID"; os << " ";
1464 os.width(14);
1465 os << " Local Index "; os << " ";
1466 os.width(14);
1467 os << " Global Index "; os << " ";
1468 if (!ConstantElementSize()) {
1469 os.width(14);
1470 os <<" FirstPointInElement "; os << " ";
1471 os.width(14);
1472 os <<" ElementSize "; os << " ";
1473 }
1474 os << std::endl;
1475
1476 for (int i = 0; i < NumMyElements(); i++) {
1477 os.width(14);
1478 os << MyPID; os << " ";
1479 os.width(14);
1480 os << i; os << " ";
1481 os.width(14);
1482
1484 {
1485#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1486 long long * MyGlobalElements1 = MyGlobalElements64();
1487 os << MyGlobalElements1[i]; os << " ";
1488#else
1489 throw ReportError("Epetra_BlockMap::Print: ERROR, GlobalIndicesLongLong but no API for it.",-1);
1490#endif
1491 }
1493 {
1494#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1495 int * MyGlobalElements1 = MyGlobalElements();
1496 os << MyGlobalElements1[i]; os << " ";
1497#else
1498 throw ReportError("Epetra_BlockMap::Print: ERROR, no GlobalIndicesLongLong but no API for it.",-1);
1499#endif
1500 }
1501
1502 if (!ConstantElementSize()) {
1503 os.width(14);
1504 os << FirstPointInElementList1[i]; os << " ";
1505 os.width(14);
1506 os << ElementSizeList1[i]; os << " ";
1507 }
1508 os << std::endl;
1509 }
1510
1511 os << std::flush;
1512
1513 }
1514 // Do a few global ops to give I/O a chance to complete
1515 Comm().Barrier();
1516 Comm().Barrier();
1517 Comm().Barrier();
1518 }
1519 return;
1520}
1521
1522//==============================================================================
1524 CleanupData();
1525}
1526
1527//==============================================================================
1529{
1530 if(BlockMapData_ != 0) {
1531
1533 if(BlockMapData_->ReferenceCount() == 0) {
1534 delete BlockMapData_;
1535 BlockMapData_ = 0;
1536 }
1537 }
1538}
1539
1540//=============================================================================
1542{
1543 if((this != &map) && (BlockMapData_ != map.BlockMapData_)) {
1544 CleanupData();
1547 }
1548
1549 return(*this);
1550}
1551
1552//=============================================================================
1554{
1555#ifdef HAVE_MPI
1556 const Epetra_MpiComm * MpiComm = dynamic_cast<const Epetra_MpiComm*>(&Comm());
1557
1558 // If the Comm isn't MPI, just treat this as a copy constructor
1559 if(!MpiComm) return new Epetra_BlockMap(*this);
1560
1561 MPI_Comm NewComm,MyMPIComm = MpiComm->Comm();
1562
1563 // Create the new communicator. MPI_Comm_split returns a valid
1564 // communicator on all processes. On processes where color == MPI_UNDEFINED,
1565 // ignore the result. Passing key == 0 tells MPI to order the
1566 // processes in the new communicator by their rank in the old
1567 // communicator.
1568 const int color = (NumMyElements() == 0) ? MPI_UNDEFINED : 1;
1569
1570 // MPI_Comm_split must be called collectively over the original
1571 // communicator. We can't just call it on processes with color
1572 // one, even though we will ignore its result on processes with
1573 // color zero.
1574 int rv = MPI_Comm_split(MyMPIComm,color,0,&NewComm);
1575 if(rv!=MPI_SUCCESS) throw ReportError("Epetra_BlockMap::RemoveEmptyProcesses: MPI_Comm_split failed.",-1);
1576
1577 if(color == MPI_UNDEFINED)
1578 return 0; // We're not in the new map
1579 else {
1580 Epetra_MpiComm * NewEpetraComm = new Epetra_MpiComm(NewComm);
1581
1582 // Use the copy constructor for a new map, but basically because it does nothing useful
1583 Epetra_BlockMap * NewMap = new Epetra_BlockMap(*this);
1584
1585 // Get rid of the old BlockMapData, now make a new one from scratch...
1586 NewMap->CleanupData();
1587 if(GlobalIndicesInt()) {
1588#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1589 NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements(),0,IndexBase(),*NewEpetraComm,false);
1590#endif
1591 }
1592 else {
1593#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1594 NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements64(),0,IndexBase64(),*NewEpetraComm,true);
1595#endif
1596 }
1597 // Now copy all of the relevent bits of BlockMapData...
1598 // NewMap->BlockMapData_->Comm_ = NewEpetraComm;
1600#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1602#endif
1603#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1605#endif
1609
1628 NewMap->BlockMapData_->DistributedGlobal_ = NewEpetraComm->NumProc()==1 ? false : BlockMapData_->DistributedGlobal_;
1636
1637 // Delay directory construction
1638 NewMap->BlockMapData_->Directory_ = 0;
1639
1640 // Cleanup
1641 delete NewEpetraComm;
1642
1643 return NewMap;
1644 }
1645#else
1646 // MPI isn't compiled, so just treat this as a copy constructor
1647 return new Epetra_BlockMap(*this);
1648#endif
1649}
1650
1651//=============================================================================
1653{
1654 // mfh 26 Mar 2013: The lazy way to do this is simply to recreate
1655 // the Map by calling its ordinary public constructor, using the
1656 // original Map's data. This only involves O(1) all-reduces over
1657 // the new communicator, which in the common case only includes a
1658 // small number of processes.
1659 Epetra_BlockMap * NewMap=0;
1660
1661 // Create the Map to return (unless theComm is zero, in which case we return zero).
1662 if (theComm) {
1663 // Map requires that the index base equal the global min GID.
1664 // Figuring out the global min GID requires a reduction over all
1665 // processes in the new communicator. It could be that some (or
1666 // even all) of these processes contain zero entries. (Recall
1667 // that this method, unlike removeEmptyProcesses(), may remove
1668 // an arbitrary subset of processes.) We deal with this by
1669 // doing a min over the min GID on each process if the process
1670 // has more than zero entries, or the global max GID, if that
1671 // process has zero entries. If no processes have any entries,
1672 // then the index base doesn't matter anyway.
1673
1674#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
1675 if(GlobalIndicesInt()) {
1676 int MyMin, theIndexBase;
1677 MyMin = NumMyElements() > 0 ? MinMyGID() : MaxAllGID();
1678 theComm->MinAll(&MyMin,&theIndexBase,1);
1679 NewMap = new Epetra_BlockMap(-1,NumMyElements(),MyGlobalElements(),ElementSizeList(),theIndexBase,*theComm);
1680 }
1681 else
1682#endif
1683#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
1684 if(GlobalIndicesLongLong()) {
1685 long long MyMin, theIndexBase;
1686 MyMin = NumMyElements() > 0 ? MinMyGID64() : MaxAllGID64();
1687 theComm->MinAll(&MyMin,&theIndexBase,1);
1688 NewMap = new Epetra_BlockMap(-1,NumMyElements(),MyGlobalElements64(),ElementSizeList(),theIndexBase,*theComm);
1689 }
1690 else
1691#endif
1692 throw ReportError("Epetra_BlockMap::ReplaceCommWithSubset ERROR, GlobalIndices type unknown.",-1);
1693 }
1694 return NewMap;
1695}
#define EPETRA_MIN(x, y)
#define EPETRA_MAX(x, y)
#define EPETRA_CHK_ERR(a)
Epetra_BlockMapData: The Epetra BlockMap Data Class.
Epetra_IntSerialDenseVector MyGlobalElements_int_
Epetra_HashTable< int > * LIDHash_
Epetra_IntSerialDenseVector LID_
Epetra_IntSerialDenseVector ElementSizeList_
Epetra_Directory * Directory_
Epetra_IntSerialDenseVector FirstPointInElementList_
const Epetra_Comm * Comm_
Epetra_LongLongSerialDenseVector MyGlobalElements_LL_
Epetra_IntSerialDenseVector PointToElementList_
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
long long IndexBase64() const
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.
void ConstructUserLinear(long long NumGlobal_Elements, int NumMy_Elements, int Element_Size, long long Index_Base, const Epetra_Comm &comm, bool IsLongLong)
int MaxElementSize() const
Maximum element size across all processors.
bool PointSameAs(const Epetra_BlockMap &Map) const
Returns true if this and Map have identical point-wise structure.
long long * MyGlobalElements64() const
bool SameBlockMapDataAs(const Epetra_BlockMap &Map) const
Returns true if maps share same block map data underneath.
long long MaxAllGID64() const
int MinMyGID() const
Returns the minimum global ID owned by this processor.
long long MaxMyGID64() const
bool DetermineIsOneToOne() const
bool GlobalIndicesTypeMatch(const Epetra_BlockMap &other) const
virtual void Print(std::ostream &os) const
Print object to an output stream.
void ConstructAutoUniform(long long NumGlobal_Elements, int Element_Size, long long Index_Base, const Epetra_Comm &comm, bool IsLongLong)
int * ElementSizeList() const
List of the element sizes corresponding to the array MyGlobalElements().
int * PointToElementList() const
For each local point, indicates the local element ID that the point belongs to.
Epetra_BlockMap(int NumGlobalElements, int ElementSize, int IndexBase, const Epetra_Comm &Comm)
Epetra_BlockMap constructor for a Epetra-defined uniform linear distribution of constant size element...
void ConstructUserVariable(int_type NumGlobal_Elements, int NumMy_Elements, const int_type *myGlobalElements, const int *elementSizeList, int_type indexBase, const Epetra_Comm &comm, bool IsLongLong)
int NumMyPoints() const
Number of local points for this map; equals the sum of all element sizes on the calling processor.
Epetra_BlockMap & operator=(const Epetra_BlockMap &map)
Assignment Operator.
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 LinearMap() const
Returns true if the global ID space is contiguously divided (but not necessarily uniformly) across al...
int IndexBase() const
Index base for this map.
Epetra_BlockMapData * BlockMapData_
int MaxAllGID() const
Returns the maximum global ID across the entire map.
long long GID64(int LID) const
int ElementSize() const
Returns the size of elements in the map; only valid if map has constant element size.
int * FirstPointInElementList() const
Pointer to internal array containing a mapping between the local elements and the first local point n...
long long NumGlobalElements64() const
bool SameAs(const Epetra_BlockMap &Map) const
Returns true if this and Map are identical maps.
bool IsDistributedGlobal(long long NumGlobalElements, int NumMyElements) const
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
Epetra_BlockMap * ReplaceCommWithSubset(const Epetra_Comm *Comm) const
Replace this BlockMap's communicator with a subset communicator.
int FindLocalElementID(int PointID, int &ElementID, int &ElementOffset) const
Returns the LID of the element that contains the given local PointID, and the Offset of the point in ...
int NumGlobalElements() const
Number of elements across all processors.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
virtual ~Epetra_BlockMap(void)
Epetra_BlockMap destructor.
bool ConstantElementSize() const
Returns true if map has constant element size.
long long MinAllGID64() const
long long MinMyGID64() const
int * MyGlobalElements() const
Pointer to internal array containing list of global IDs assigned to the calling processor.
int NumMyElements() const
Number of elements on the calling processor.
bool IsOneToOne() const
bool MyLID(int lid) const
Returns true if the LID passed in belongs to the calling processor in this map, otherwise returns fal...
long long NumGlobalPoints64() const
void ConstructUserConstantNoComm(int_type NumGlobal_Elements, int NumMy_Elements, const int_type *myGlobalElements, int ElementSize, int_type indexBase, const Epetra_Comm &comm, bool IsLongLong, bool UserIsDistributedGlobal, int_type UserMinAllGID, int_type UserMaxAllGID)
int FirstPointInElement(int LID) const
Returns the requested entry in the FirstPointInElementList; see FirstPointInElementList() for details...
void ConstructUserConstant(int_type NumGlobal_Elements, int NumMy_Elements, const int_type *myGlobalElements, int Element_Size, int_type indexBase, const Epetra_Comm &comm, bool IsLongLong)
void CheckValidNGE(long long NumGlobalElements)
bool GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
Epetra_BlockMap * RemoveEmptyProcesses() const
Return a new BlockMap with processes with zero elements removed.
int MyGlobalElementsPtr(int *&MyGlobalElementList) const
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
virtual int MaxAll(double *PartialMaxs, double *GlobalMaxs, int Count) const =0
Epetra_Comm Global Max function.
virtual int NumProc() const =0
Returns total number of processes.
virtual Epetra_Directory * CreateDirectory(const Epetra_BlockMap &Map) const =0
Create a directory object for the given Epetra_BlockMap.
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 int ScanSum(double *MyVals, double *ScanSums, int Count) const =0
Epetra_Comm Scan Sum function.
virtual void Barrier() const =0
Epetra_Comm Barrier function.
virtual int SumAll(double *PartialSums, double *GlobalSums, int Count) const =0
Epetra_Comm Global Sum function.
void DecrementReferenceCount()
Decrement reference count.
Definition: Epetra_Data.cpp:66
int ReferenceCount() const
Get reference count.
Definition: Epetra_Data.cpp:71
void IncrementReferenceCount()
Increment reference count.
Definition: Epetra_Data.cpp:61
Epetra_Directory: This class is a pure virtual class whose interface allows Epetra_Map and Epetr_Bloc...
virtual int GetDirectoryEntries(const Epetra_BlockMap &Map, const int NumEntries, const int *GlobalEntries, int *Procs, int *LocalEntries, int *EntrySizes, bool high_rank_sharing_procs=false) const =0
GetDirectoryEntries : Returns proc and local id info for non-local map entries.
virtual bool GIDsAllUniquelyOwned() const =0
GIDsAllUniquelyOwned: returns true if all GIDs appear on just one processor.
void Add(const long long key, const value_type value)
value_type Get(const long long key)
int Size(int Length_in)
Set length of a Epetra_IntSerialDenseVector object; init values to zero.
int Length() const
Returns length of vector.
int * Values()
Returns pointer to the values in vector.
int Length() const
Returns length of vector.
int Size(int Length_in)
Set length of a Epetra_LongLongSerialDenseVector object; init values to zero.
long long * Values()
Returns pointer to the values in vector.
Epetra_MpiComm: The Epetra MPI Communication Class.
MPI_Comm Comm() const
Extract MPI Communicator from a Epetra_MpiComm object.
int NumProc() const
Returns total number of processes.
Epetra_Object: The base Epetra class.
Definition: Epetra_Object.h:57
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
std::string toString(const int &x) const