Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Epetra_Map.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_Map.h"
45#include "Epetra_Comm.h"
46#ifdef HAVE_MPI
47#include "Epetra_MpiComm.h"
48#endif
49#include "Epetra_HashTable.h"
50
51//==============================================================================
52// Epetra_Map constructor for a Epetra-defined uniform linear distribution of elements.
53#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
54Epetra_Map::Epetra_Map(int numGlobalElements, int indexBase, const Epetra_Comm& comm)
55 : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
56{
57 SetLabel("Epetra::Map");
58}
59#endif
60#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
61Epetra_Map::Epetra_Map(long long numGlobalElements, int indexBase, const Epetra_Comm& comm)
62 : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
63{
64 SetLabel("Epetra::Map");
65}
66
67Epetra_Map::Epetra_Map(long long numGlobalElements, long long indexBase, const Epetra_Comm& comm)
68 : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
69{
70 SetLabel("Epetra::Map");
71}
72#endif
73//==============================================================================
74// Epetra_Map constructor for a user-defined linear distribution of constant block size elements.
75#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
76Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements, int indexBase, const Epetra_Comm& comm)
77 : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
78{
79 SetLabel("Epetra::Map");
80}
81#endif
82#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
83Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements, int indexBase, const Epetra_Comm& comm)
84 : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
85{
86 SetLabel("Epetra::Map");
87}
88
89Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements, long long indexBase, const Epetra_Comm& comm)
90 : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
91{
92 SetLabel("Epetra::Map");
93}
94#endif
95//==============================================================================
96// Epetra_Map constructor for a user-defined arbitrary distribution of constant block size elements.
97#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
98Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements,
99 const int * myGlobalElements,
100 int indexBase, const Epetra_Comm& comm)
101 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
102{
103 SetLabel("Epetra::Map");
104}
105#endif
106#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
107Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
108 const long long * myGlobalElements,
109 int indexBase, const Epetra_Comm& comm)
110 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
111{
112 SetLabel("Epetra::Map");
113}
114
115Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
116 const long long * myGlobalElements,
117 long long indexBase, const Epetra_Comm& comm)
118 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
119{
120 SetLabel("Epetra::Map");
121}
122#endif
123
124//==============================================================================
125// Epetra_Map constructor for a user-defined arbitrary distribution of constant block size elements w/ user provided globals.
126#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
127Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements,
128 const int * myGlobalElements,
129 int indexBase, const Epetra_Comm& comm,
130 bool UserIsDistributedGlobal,
131 int UserMinAllGID, int UserMaxAllGID)
132 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
133{
134 SetLabel("Epetra::Map");
135}
136#endif
137#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
138Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
139 const long long * myGlobalElements,
140 int indexBase, const Epetra_Comm& comm,
141 bool UserIsDistributedGlobal,
142 long long UserMinAllGID, long long UserMaxAllGID)
143 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
144{
145 SetLabel("Epetra::Map");
146}
147Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
148 const long long * myGlobalElements,
149 long long indexBase, const Epetra_Comm& comm,
150 bool UserIsDistributedGlobal,
151 long long UserMinAllGID, long long UserMaxAllGID)
152 : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
153{
154 SetLabel("Epetra::Map");
155}
156#endif
157
158//==============================================================================
160 : Epetra_BlockMap(map) // Map is just a special case of BlockMap
161{
162}
163
164//==============================================================================
166{
167}
168
169//=============================================================================
171{
172 if(this != &map) {
173 Epetra_BlockMap::operator=(map); // call this->Epetra_BlockMap::operator=
174 }
175 return(*this);
176}
177
178//=============================================================================
180{
181#ifdef HAVE_MPI
182 const Epetra_MpiComm * MpiComm = dynamic_cast<const Epetra_MpiComm*>(&Comm());
183
184 // If the Comm isn't MPI, just treat this as a copy constructor
185 if(!MpiComm) return new Epetra_Map(*this);
186
187 MPI_Comm NewComm,MyMPIComm = MpiComm->Comm();
188
189 // Create the new communicator. MPI_Comm_split returns a valid
190 // communicator on all processes. On processes where color == MPI_UNDEFINED,
191 // ignore the result. Passing key == 0 tells MPI to order the
192 // processes in the new communicator by their rank in the old
193 // communicator.
194 const int color = (NumMyElements() == 0) ? MPI_UNDEFINED : 1;
195
196 // MPI_Comm_split must be called collectively over the original
197 // communicator. We can't just call it on processes with color
198 // one, even though we will ignore its result on processes with
199 // color zero.
200 int rv = MPI_Comm_split(MyMPIComm,color,0,&NewComm);
201 if(rv!=MPI_SUCCESS) throw ReportError("Epetra_Map::RemoveEmptyProcesses: MPI_Comm_split failed.",-1);
202
203 if(color == MPI_UNDEFINED)
204 return 0; // We're not in the new map
205 else {
206 Epetra_MpiComm * NewEpetraComm = new Epetra_MpiComm(NewComm);
207
208 // Use the copy constructor for a new map, but basically because it does nothing useful
209 Epetra_Map * NewMap = new Epetra_Map(*this);
210
211 // Get rid of the old BlockMapData, now make a new one from scratch...
212 NewMap->CleanupData();
213 if(GlobalIndicesInt()) {
214#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
215 NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements(),0,IndexBase(),*NewEpetraComm,false);
216#endif
217 }
218 else {
219#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
220 NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements64(),0,IndexBase64(),*NewEpetraComm,true);
221#endif
222 }
223
224 // Now copy all of the relevent bits of BlockMapData...
225 // NewMap->BlockMapData_->Comm_ = NewEpetraComm;
227#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
229#endif
230#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
232#endif
236
255 NewMap->BlockMapData_->DistributedGlobal_ = NewEpetraComm->NumProc()==1 ? false : BlockMapData_->DistributedGlobal_;
263
264 // Delay directory construction
265 NewMap->BlockMapData_->Directory_ = 0;
266
267 // Cleanup
268 delete NewEpetraComm;
269
270 return NewMap;
271 }
272#else
273 // MPI isn't compiled, so just treat this as a copy constructor
274 return new Epetra_Map(*this);
275#endif
276}
277
278//=============================================================================
280{
281 // mfh 26 Mar 2013: The lazy way to do this is simply to recreate
282 // the Map by calling its ordinary public constructor, using the
283 // original Map's data. This only involves O(1) all-reduces over
284 // the new communicator, which in the common case only includes a
285 // small number of processes.
286 Epetra_Map * NewMap=0;
287
288 // Create the Map to return (unless theComm is NULL, in which case
289 // we return zero).
290 if(theComm) {
291 // Note: we replace the communicator with a subset of that communicator (that is
292 // we do not just remove the empty procs). This does not affect the index
293 // base: we keep the existing index base.
294#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
295 if(GlobalIndicesInt()) {
296 NewMap = new Epetra_Map(-1,NumMyElements(),MyGlobalElements(),IndexBase(),*theComm);
297 }
298 else
299#endif
300#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
302 NewMap = new Epetra_Map(-1,NumMyElements(),MyGlobalElements64(),IndexBase64(),*theComm);
303 }
304 else
305#endif
306 throw ReportError("Epetra_Map::ReplaceCommWithSubset ERROR, GlobalIndices type unknown.",-1);
307 }
308 return NewMap;
309}
310
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_
Epetra_LongLongSerialDenseVector MyGlobalElements_LL_
Epetra_IntSerialDenseVector PointToElementList_
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
long long IndexBase64() const
long long * MyGlobalElements64() const
Epetra_BlockMap & operator=(const Epetra_BlockMap &map)
Assignment Operator.
int IndexBase() const
Index base for this map.
Epetra_BlockMapData * BlockMapData_
long long NumGlobalElements64() const
bool GlobalIndicesInt() const
Returns true if map create with int NumGlobalElements.
int NumGlobalElements() const
Number of elements across all processors.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
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 GlobalIndicesLongLong() const
Returns true if map create with long long NumGlobalElements.
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
Epetra_Map & operator=(const Epetra_Map &map)
Assignment Operator.
Definition: Epetra_Map.cpp:170
virtual ~Epetra_Map(void)
Epetra_Map destructor.
Definition: Epetra_Map.cpp:165
Epetra_Map(int NumGlobalElements, int IndexBase, const Epetra_Comm &Comm)
Epetra_Map constructor for a Epetra-defined uniform linear distribution of elements.
Definition: Epetra_Map.cpp:54
Epetra_Map * RemoveEmptyProcesses() const
Return a new BlockMap with processes with zero elements removed.
Definition: Epetra_Map.cpp:179
Epetra_Map * ReplaceCommWithSubset(const Epetra_Comm *Comm) const
Replace this Map's communicator with a subset communicator.
Definition: Epetra_Map.cpp:279
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.
virtual int ReportError(const std::string Message, int ErrorCode) const
Error reporting method.
virtual void SetLabel(const char *const Label)
Epetra_Object Label definition using char *.