Epetra Package Browser (Single Doxygen Collection) Development
Epetra_Map.cpp
Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 //
00005 //               Epetra: Linear Algebra Services Package
00006 //                 Copyright 2011 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00039 //
00040 // ************************************************************************
00041 //@HEADER
00042 
00043 #include "Epetra_ConfigDefs.h"
00044 #include "Epetra_Map.h"
00045 #include "Epetra_Comm.h"
00046 #ifdef HAVE_MPI
00047 #include "Epetra_MpiComm.h"
00048 #endif
00049 #include "Epetra_HashTable.h"
00050 
00051 //==============================================================================
00052 // Epetra_Map constructor for a Epetra-defined uniform linear distribution of elements.
00053 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00054 Epetra_Map::Epetra_Map(int numGlobalElements, int indexBase, const Epetra_Comm& comm)
00055   : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00056 {
00057   SetLabel("Epetra::Map");
00058 }
00059 #endif
00060 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00061 Epetra_Map::Epetra_Map(long long numGlobalElements, int indexBase, const Epetra_Comm& comm)
00062   : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00063 {
00064   SetLabel("Epetra::Map");
00065 }
00066 
00067 Epetra_Map::Epetra_Map(long long numGlobalElements, long long indexBase, const Epetra_Comm& comm)
00068   : Epetra_BlockMap(numGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00069 {
00070   SetLabel("Epetra::Map");
00071 }
00072 #endif
00073 //==============================================================================
00074 // Epetra_Map constructor for a user-defined linear distribution of constant block size elements.
00075 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00076 Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements, int indexBase, const Epetra_Comm& comm)
00077   : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00078 {
00079   SetLabel("Epetra::Map");
00080 }
00081 #endif
00082 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00083 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements, int indexBase, const Epetra_Comm& comm)
00084   : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00085 {
00086   SetLabel("Epetra::Map");
00087 }
00088 
00089 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements, long long indexBase, const Epetra_Comm& comm)
00090   : Epetra_BlockMap(numGlobalElements, numMyElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00091 {
00092   SetLabel("Epetra::Map");
00093 }
00094 #endif
00095 //==============================================================================
00096 // Epetra_Map constructor for a user-defined arbitrary distribution of constant block size elements.
00097 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00098 Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements,
00099                        const int * myGlobalElements,
00100                        int indexBase, const Epetra_Comm& comm)
00101   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00102 {
00103   SetLabel("Epetra::Map");
00104 }
00105 #endif
00106 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00107 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
00108                        const long long * myGlobalElements,
00109                        int indexBase, const Epetra_Comm& comm)
00110   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00111 {
00112   SetLabel("Epetra::Map");
00113 }
00114 
00115 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
00116                        const long long * myGlobalElements,
00117                        long long indexBase, const Epetra_Comm& comm)
00118   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm) // Map is just a special case of BlockMap
00119 {
00120   SetLabel("Epetra::Map");
00121 }
00122 #endif
00123 
00124 //==============================================================================
00125 // Epetra_Map constructor for a user-defined arbitrary distribution of constant block size elements w/ user provided globals.
00126 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00127 Epetra_Map::Epetra_Map(int numGlobalElements, int numMyElements,
00128                        const int * myGlobalElements,
00129                        int indexBase, const Epetra_Comm& comm,
00130            bool UserIsDistributedGlobal,
00131            int UserMinAllGID, int UserMaxAllGID)
00132   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
00133 {
00134   SetLabel("Epetra::Map");
00135 }
00136 #endif
00137 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00138 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
00139                        const long long * myGlobalElements,
00140                        int indexBase, const Epetra_Comm& comm,
00141            bool UserIsDistributedGlobal,
00142            long long UserMinAllGID, long long UserMaxAllGID)
00143   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
00144 {
00145   SetLabel("Epetra::Map");
00146 }
00147 Epetra_Map::Epetra_Map(long long numGlobalElements, int numMyElements,
00148                        const long long * myGlobalElements,
00149                        long long indexBase, const Epetra_Comm& comm,
00150            bool UserIsDistributedGlobal,
00151            long long UserMinAllGID, long long UserMaxAllGID)
00152   : Epetra_BlockMap(numGlobalElements, numMyElements, myGlobalElements, 1, indexBase, comm, UserIsDistributedGlobal, UserMinAllGID, UserMaxAllGID) // Map is just a special case of BlockMap
00153 {
00154   SetLabel("Epetra::Map");
00155 }
00156 #endif
00157 
00158 //==============================================================================
00159 Epetra_Map::Epetra_Map(const Epetra_Map& map)
00160   : Epetra_BlockMap(map) // Map is just a special case of BlockMap
00161 {
00162 }
00163 
00164 //==============================================================================
00165 Epetra_Map::~Epetra_Map(void)
00166 {
00167 }
00168 
00169 //=============================================================================
00170 Epetra_Map & Epetra_Map::operator= (const Epetra_Map& map)
00171 {
00172   if(this != &map) {
00173     Epetra_BlockMap::operator=(map); // call this->Epetra_BlockMap::operator=
00174   }
00175   return(*this);
00176 }
00177 
00178 //=============================================================================
00179 Epetra_Map * Epetra_Map::RemoveEmptyProcesses() const
00180 {
00181 #ifdef HAVE_MPI
00182   const Epetra_MpiComm * MpiComm = dynamic_cast<const Epetra_MpiComm*>(&Comm());
00183 
00184   // If the Comm isn't MPI, just treat this as a copy constructor
00185   if(!MpiComm) return new Epetra_Map(*this);
00186 
00187   MPI_Comm NewComm,MyMPIComm = MpiComm->Comm();
00188 
00189   // Create the new communicator.  MPI_Comm_split returns a valid
00190   // communicator on all processes.  On processes where color == MPI_UNDEFINED,
00191   // ignore the result.  Passing key == 0 tells MPI to order the
00192   // processes in the new communicator by their rank in the old
00193   // communicator.
00194   const int color = (NumMyElements() == 0) ? MPI_UNDEFINED : 1;
00195 
00196   // MPI_Comm_split must be called collectively over the original
00197   // communicator.  We can't just call it on processes with color
00198   // one, even though we will ignore its result on processes with
00199   // color zero.
00200   int rv = MPI_Comm_split(MyMPIComm,color,0,&NewComm);
00201   if(rv!=MPI_SUCCESS) throw ReportError("Epetra_Map::RemoveEmptyProcesses: MPI_Comm_split failed.",-1);
00202 
00203   if(color == MPI_UNDEFINED)
00204     return 0; // We're not in the new map
00205   else {
00206     Epetra_MpiComm * NewEpetraComm = new Epetra_MpiComm(NewComm);
00207 
00208     // Use the copy constructor for a new map, but basically because it does nothing useful
00209     Epetra_Map * NewMap = new Epetra_Map(*this);
00210 
00211     // Get rid of the old BlockMapData, now make a new one from scratch...
00212     NewMap->CleanupData();
00213     if(GlobalIndicesInt()) {
00214 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00215       NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements(),0,IndexBase(),*NewEpetraComm,false);
00216 #endif
00217     }
00218     else {
00219 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00220       NewMap->BlockMapData_ = new Epetra_BlockMapData(NumGlobalElements64(),0,IndexBase64(),*NewEpetraComm,true);
00221 #endif
00222     }
00223 
00224     // Now copy all of the relevent bits of BlockMapData...
00225     //    NewMap->BlockMapData_->Comm_                    = NewEpetraComm;
00226     NewMap->BlockMapData_->LID_                     = BlockMapData_->LID_;
00227 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00228     NewMap->BlockMapData_->MyGlobalElements_int_    = BlockMapData_->MyGlobalElements_int_;
00229 #endif
00230 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00231     NewMap->BlockMapData_->MyGlobalElements_LL_     = BlockMapData_->MyGlobalElements_LL_;
00232 #endif
00233     NewMap->BlockMapData_->FirstPointInElementList_ = BlockMapData_->FirstPointInElementList_;
00234     NewMap->BlockMapData_->ElementSizeList_         = BlockMapData_->ElementSizeList_;
00235     NewMap->BlockMapData_->PointToElementList_      = BlockMapData_->PointToElementList_;
00236 
00237     NewMap->BlockMapData_->NumGlobalElements_       = BlockMapData_->NumGlobalElements_;
00238     NewMap->BlockMapData_->NumMyElements_           = BlockMapData_->NumMyElements_;
00239     NewMap->BlockMapData_->IndexBase_               = BlockMapData_->IndexBase_;
00240     NewMap->BlockMapData_->ElementSize_             = BlockMapData_->ElementSize_;
00241     NewMap->BlockMapData_->MinMyElementSize_        = BlockMapData_->MinMyElementSize_;
00242     NewMap->BlockMapData_->MaxMyElementSize_        = BlockMapData_->MaxMyElementSize_;
00243     NewMap->BlockMapData_->MinElementSize_          = BlockMapData_->MinElementSize_;
00244     NewMap->BlockMapData_->MaxElementSize_          = BlockMapData_->MaxElementSize_;
00245     NewMap->BlockMapData_->MinAllGID_               = BlockMapData_->MinAllGID_;
00246     NewMap->BlockMapData_->MaxAllGID_               = BlockMapData_->MaxAllGID_;
00247     NewMap->BlockMapData_->MinMyGID_                = BlockMapData_->MinMyGID_;
00248     NewMap->BlockMapData_->MaxMyGID_                = BlockMapData_->MaxMyGID_;
00249     NewMap->BlockMapData_->MinLID_                  = BlockMapData_->MinLID_;
00250     NewMap->BlockMapData_->MaxLID_                  = BlockMapData_->MaxLID_;
00251     NewMap->BlockMapData_->NumGlobalPoints_         = BlockMapData_->NumGlobalPoints_;
00252     NewMap->BlockMapData_->NumMyPoints_             = BlockMapData_->NumMyPoints_;
00253     NewMap->BlockMapData_->ConstantElementSize_     = BlockMapData_->ConstantElementSize_;
00254     NewMap->BlockMapData_->LinearMap_               = BlockMapData_->LinearMap_;
00255     NewMap->BlockMapData_->DistributedGlobal_       = NewEpetraComm->NumProc()==1 ? false : BlockMapData_->DistributedGlobal_;
00256     NewMap->BlockMapData_->OneToOneIsDetermined_    = BlockMapData_->OneToOneIsDetermined_;
00257     NewMap->BlockMapData_->OneToOne_                = BlockMapData_->OneToOne_;
00258     NewMap->BlockMapData_->GlobalIndicesInt_        = BlockMapData_->GlobalIndicesInt_;
00259     NewMap->BlockMapData_->GlobalIndicesLongLong_   = BlockMapData_->GlobalIndicesLongLong_;
00260     NewMap->BlockMapData_->LastContiguousGID_       = BlockMapData_->LastContiguousGID_;
00261     NewMap->BlockMapData_->LastContiguousGIDLoc_    = BlockMapData_->LastContiguousGIDLoc_;
00262     NewMap->BlockMapData_->LIDHash_                 = BlockMapData_->LIDHash_ ? new Epetra_HashTable<int>(*BlockMapData_->LIDHash_) : 0;
00263 
00264     // Delay directory construction
00265     NewMap->BlockMapData_->Directory_               = 0;
00266 
00267     // Cleanup
00268     delete NewEpetraComm;
00269 
00270     return NewMap;
00271   }
00272 #else
00273     // MPI isn't compiled, so just treat this as a copy constructor
00274     return new Epetra_Map(*this);
00275 #endif
00276 }
00277 
00278 //=============================================================================
00279 Epetra_Map* Epetra_Map::ReplaceCommWithSubset(const Epetra_Comm * Comm) const
00280 {
00281   // mfh 26 Mar 2013: The lazy way to do this is simply to recreate
00282   // the Map by calling its ordinary public constructor, using the
00283   // original Map's data.  This only involves O(1) all-reduces over
00284   // the new communicator, which in the common case only includes a
00285   // small number of processes.
00286   Epetra_Map * NewMap=0;
00287 
00288   // Create the Map to return (unless Comm is zero, in which case we return zero).
00289   if(Comm) {
00290     // Map requires that the index base equal the global min GID.
00291     // Figuring out the global min GID requires a reduction over all
00292     // processes in the new communicator.  It could be that some (or
00293     // even all) of these processes contain zero entries.  (Recall
00294     // that this method, unlike removeEmptyProcesses(), may remove
00295     // an arbitrary subset of processes.)  We deal with this by
00296     // doing a min over the min GID on each process if the process
00297     // has more than zero entries, or the global max GID, if that
00298     // process has zero entries.  If no processes have any entries,
00299     // then the index base doesn't matter anyway.
00300 
00301 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00302     if(GlobalIndicesInt()) {
00303       int MyMin, IndexBase;
00304       MyMin  = NumMyElements() > 0 ? MinMyGID() : MaxAllGID();
00305       Comm->MinAll(&MyMin,&IndexBase,1);
00306       NewMap = new Epetra_Map(-1,NumMyElements(),MyGlobalElements(),IndexBase,*Comm);
00307     }
00308     else
00309 #endif
00310 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00311     if(GlobalIndicesLongLong()) {
00312       long long MyMin, IndexBase;
00313       MyMin = NumMyElements() > 0 ? MinMyGID64() : MaxAllGID64();
00314       Comm->MinAll(&MyMin,&IndexBase,1);
00315       NewMap = new Epetra_Map(-1,NumMyElements(),MyGlobalElements64(),IndexBase,*Comm);
00316     }
00317     else
00318 #endif
00319     throw ReportError("Epetra_Map::ReplaceCommWithSubset ERROR, GlobalIndices type unknown.",-1);
00320   }
00321   return NewMap;
00322 }
00323 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines