EpetraExt Development
EpetraExt_MMHelpers.h
Go to the documentation of this file.
00001 //@HEADER
00002 // ***********************************************************************
00003 //
00004 //     EpetraExt: Epetra Extended - Linear Algebra Services Package
00005 //                 Copyright (2011) Sandia Corporation
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ***********************************************************************
00040 //@HEADER
00041 
00042 #ifndef EPETRAEXT_MMHELPERS_H
00043 #define EPETRAEXT_MMHELPERS_H
00044 
00045 #include "EpetraExt_ConfigDefs.h"
00046 #include "Epetra_ConfigDefs.h"
00047 #include "Epetra_DistObject.h"
00048 #include "Epetra_Map.h"
00049 #include "Teuchos_RCP.hpp"
00050 
00051 #include <vector>
00052 #include <set>
00053 #include <map>
00054 
00055 
00056 class Epetra_CrsMatrix;
00057 class Epetra_Import;
00058 class Epetra_Distributor;
00059 
00060 namespace EpetraExt {
00061 class LightweightCrsMatrix;
00062 
00063 // Turn on to enable the special MMM timers
00064 //#define ENABLE_MMM_TIMINGS
00065 
00066 //#define HAVE_EPETRAEXT_DEBUG // for extra sanity checks
00067 
00068 // ==============================================================
00069 //struct that holds views of the contents of a CrsMatrix. These
00070 //contents may be a mixture of local and remote rows of the
00071 //actual matrix.
00072 class CrsMatrixStruct {
00073 public:
00074   CrsMatrixStruct();
00075 
00076   virtual ~CrsMatrixStruct();
00077 
00078   void deleteContents();
00079 
00080   int numRows;
00081 
00082   // The following class members get used in the transpose modes of the MMM
00083   // but not in the A*B mode.  
00084   int* numEntriesPerRow;
00085   int** indices;
00086   double** values;
00087   bool* remote;
00088   int numRemote;
00089   const Epetra_BlockMap* importColMap;
00090 
00091   // Maps and matrices (all modes)
00092   const Epetra_Map* origRowMap;
00093   const Epetra_Map* rowMap;
00094   const Epetra_Map* colMap;
00095   const Epetra_Map* domainMap;
00096   LightweightCrsMatrix* importMatrix;
00097   const Epetra_CrsMatrix *origMatrix;
00098 
00099   // The following class members are only used for A*B mode
00100   std::vector<int> targetMapToOrigRow;
00101   std::vector<int> targetMapToImportRow;
00102 };
00103 
00104 int dumpCrsMatrixStruct(const CrsMatrixStruct& M);
00105 
00106 // ==============================================================
00107 class CrsWrapper {
00108  public:
00109   virtual ~CrsWrapper(){}
00110 
00111   virtual const Epetra_Map& RowMap() const = 0;
00112 
00113   virtual bool Filled() = 0;
00114 
00115 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00116   virtual int InsertGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices) = 0;
00117 
00118   virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices) = 0;
00119 #endif
00120 
00121 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00122   virtual int InsertGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices) = 0;
00123 
00124   virtual int SumIntoGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices) = 0;
00125 #endif
00126 };
00127 
00128 // ==============================================================
00129 class CrsWrapper_Epetra_CrsMatrix : public CrsWrapper {
00130  public:
00131   CrsWrapper_Epetra_CrsMatrix(Epetra_CrsMatrix& epetracrsmatrix);
00132   virtual ~CrsWrapper_Epetra_CrsMatrix();
00133 
00134   const Epetra_Map& RowMap() const;
00135 
00136   bool Filled();
00137 
00138 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00139   int InsertGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices);
00140   int SumIntoGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices);
00141 #endif
00142 
00143 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00144   int InsertGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices);
00145   int SumIntoGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices);
00146 #endif
00147 
00148  private:
00149   Epetra_CrsMatrix& ecrsmat_;
00150 };
00151 
00152 // ==============================================================
00153 template<typename int_type>
00154 class CrsWrapper_GraphBuilder : public CrsWrapper {
00155  public:
00156   CrsWrapper_GraphBuilder(const Epetra_Map& emap);
00157   virtual ~CrsWrapper_GraphBuilder();
00158 
00159   const Epetra_Map& RowMap() const {return rowmap_; }
00160 
00161   bool Filled();
00162 
00163   int InsertGlobalValues(int_type GlobalRow, int NumEntries, double* Values, int_type* Indices);
00164   int SumIntoGlobalValues(int_type GlobalRow, int NumEntries, double* Values, int_type* Indices);
00165 
00166   std::map<int_type,std::set<int_type>*>& get_graph();
00167 
00168   int get_max_row_length() { return max_row_length_; }
00169 
00170  private:
00171   std::map<int_type,std::set<int_type>*> graph_;
00172   const Epetra_Map& rowmap_;
00173   int max_row_length_;
00174 };
00175 
00176 // ==============================================================
00177 template<typename int_type>
00178 void insert_matrix_locations(CrsWrapper_GraphBuilder<int_type>& graphbuilder,
00179                               Epetra_CrsMatrix& C);
00180 
00181 template<typename int_type>
00182 void Tpack_outgoing_rows(const Epetra_CrsMatrix& mtx,
00183                         const std::vector<int_type>& proc_col_ranges,
00184                         std::vector<int_type>& send_rows,
00185                         std::vector<int>& rows_per_send_proc);
00186 
00187 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00188 void pack_outgoing_rows(const Epetra_CrsMatrix& mtx,
00189                         const std::vector<int>& proc_col_ranges,
00190                         std::vector<int>& send_rows,
00191                         std::vector<int>& rows_per_send_proc);
00192 #endif
00193 
00194 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00195 void pack_outgoing_rows(const Epetra_CrsMatrix& mtx,
00196                         const std::vector<long long>& proc_col_ranges,
00197                         std::vector<long long>& send_rows,
00198                         std::vector<int>& rows_per_send_proc);
00199 #endif
00200 
00201 template<typename int_type>
00202 std::pair<int_type,int_type> get_col_range(const Epetra_Map& emap);
00203 
00204 template<typename int_type>
00205 std::pair<int_type,int_type> get_col_range(const Epetra_CrsMatrix& mtx);
00206 
00207 // ==============================================================
00208 class LightweightMapData : Epetra_Data {
00209   friend class LightweightMap; 
00210  public:
00211   LightweightMapData();
00212   ~LightweightMapData();
00213   long long IndexBase_;
00214 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00215   std::vector<int> MyGlobalElements_int_; 
00216   Epetra_HashTable<int> * LIDHash_int_;
00217 #endif
00218 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00219   std::vector<long long> MyGlobalElements_LL_; 
00220   Epetra_HashTable<long long> * LIDHash_LL_;
00221 #endif
00222 
00223   // For "copy" constructor only...
00224   Epetra_Map * CopyMap_;
00225   
00226  };
00227 
00228 class LightweightMap {
00229  public:
00230   LightweightMap();
00231 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00232   LightweightMap(int NumGlobalElements,int NumMyElements, const int * MyGlobalElements, int IndexBase, bool GenerateHash=true);
00233 #endif
00234 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00235   LightweightMap(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, int IndexBase, bool GenerateHash=true);
00236   LightweightMap(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
00237 #endif
00238   LightweightMap(const Epetra_Map & Map);
00239   LightweightMap(const LightweightMap & Map);
00240   ~LightweightMap();
00241 
00242   LightweightMap & operator=(const LightweightMap & map);
00243 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00244   int LID(int GID) const;
00245   int GID(int LID) const;      
00246 #endif
00247 
00248 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00249   int  LID(long long GID) const;
00250 #endif
00251   long long GID64(int LID) const;
00252   int NumMyElements() const;
00253 
00254 #if defined(EPETRA_NO_32BIT_GLOBAL_INDICES) && defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
00255   // default implementation so that no compiler/linker error in case neither 32 nor 64
00256   // bit indices present.
00257   int  LID(long long GID) const { return -1; }
00258 #endif
00259 
00260 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00261   int* MyGlobalElements() const; 
00262   int IndexBase() const {
00263     if(IndexBase64() == (long long) static_cast<int>(IndexBase64()))
00264       return (int) IndexBase64();
00265     throw "EpetraExt::LightweightMap::IndexBase: IndexBase cannot fit an int.";
00266   }
00267   void MyGlobalElementsPtr(int *& MyGlobalElementList) const;
00268 #endif
00269 
00270 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00271   long long* MyGlobalElements64() const;
00272   void MyGlobalElementsPtr(long long *& MyGlobalElementList) const;
00273 #endif
00274   long long IndexBase64() const {return Data_->IndexBase_;}
00275 
00276   int MinLID() const;
00277   int MaxLID() const;
00278 
00279   bool GlobalIndicesInt() const { return IsInt; }
00280   bool GlobalIndicesLongLong() const { return IsLongLong; }
00281  private:
00282   void CleanupData();
00283   LightweightMapData *Data_;
00284   bool IsLongLong;
00285   bool IsInt;
00286   //Epetra_BlockMapData* Data_;
00287 
00288 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00289   void Construct_int(int NumGlobalElements,int NumMyElements, const int * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
00290 #endif
00291 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00292   void Construct_LL(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
00293 #endif
00294 };
00295 
00296 
00297 // ==============================================================
00298 class RemoteOnlyImport {
00299  public:
00300   RemoteOnlyImport(const Epetra_Import & Importer, LightweightMap & RemoteOnlyTargetMap);
00301   ~RemoteOnlyImport();
00302 
00303   int NumSameIDs() {return 0;}
00304 
00305   int NumPermuteIDs() {return 0;} 
00306 
00307   int NumRemoteIDs() {return NumRemoteIDs_;}
00308 
00309   int NumExportIDs() {return NumExportIDs_;}
00310 
00311   int* ExportLIDs() {return ExportLIDs_;}
00312 
00313   int* ExportPIDs() {return ExportPIDs_;}
00314 
00315   int* RemoteLIDs() {return RemoteLIDs_;}  
00316 
00317   int* PermuteToLIDs() {return 0;}
00318 
00319   int* PermuteFromLIDs() {return 0;}
00320 
00321   int NumSend() {return NumSend_;}
00322   
00323   Epetra_Distributor & Distributor() {return *Distor_;}  
00324 
00325   const Epetra_BlockMap & SourceMap() const {return *SourceMap_;}
00326   const LightweightMap & TargetMap() const {return *TargetMap_;}
00327 
00328  private:
00329   int NumSend_;
00330   int NumRemoteIDs_;
00331   int NumExportIDs_;
00332   int * ExportLIDs_;
00333   int * ExportPIDs_;
00334   int * RemoteLIDs_;
00335   Epetra_Distributor* Distor_;
00336   const Epetra_BlockMap* SourceMap_;
00337   const LightweightMap  *TargetMap_;
00338 };
00339    
00340 // ==============================================================
00341 class LightweightCrsMatrix {
00342  public:
00343   LightweightCrsMatrix(const Epetra_CrsMatrix & A, RemoteOnlyImport & RowImporter);
00344   LightweightCrsMatrix(const Epetra_CrsMatrix & A, Epetra_Import & RowImporter);
00345   ~LightweightCrsMatrix();
00346 
00347   // Standard crs data structures
00348   std::vector<int>    rowptr_;
00349   std::vector<int>    colind_;
00350   std::vector<double> vals_;
00351 
00352 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00353   // Colind in LL-GID space (if needed)
00354   std::vector<long long>   colind_LL_;
00355 #endif
00356 
00357   // Epetra Maps
00358   bool                     use_lw;
00359   LightweightMap           *RowMapLW_;
00360   Epetra_BlockMap          *RowMapEP_;
00361   LightweightMap           ColMap_;
00362   Epetra_Map               DomainMap_;
00363 
00364 
00365   // List of owning PIDs (from the DomainMap) as ordered by entries in the column map.
00366   std::vector<int>    ColMapOwningPIDs_;
00367 
00368   // For building the final importer for C
00369   std::vector<int>    ExportLIDs_;
00370   std::vector<int>    ExportPIDs_;
00371 
00372  private: 
00373 
00374   template <typename ImportType, typename int_type>
00375   void Construct(const Epetra_CrsMatrix & A, ImportType & RowImporter);
00376 
00377   // Templated versions of MakeColMapAndReindex (to prevent code duplication)
00378   template <class GO>
00379   int MakeColMapAndReindex(std::vector<int> owningPIDs,std::vector<GO> Gcolind);
00380 
00381   template<typename int_type>
00382   std::vector<int_type>& getcolind();
00383 
00384   template<typename ImportType, typename int_type>
00385   int PackAndPrepareReverseComm(const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
00386         std::vector<int> &ReverseSendSizes, std::vector<int_type> &ReverseSendBuffer);
00387 
00388   template<typename ImportType, typename int_type>
00389   int MakeExportLists(const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
00390           std::vector<int> &ReverseRecvSizes, const int_type *ReverseRecvBuffer,
00391           std::vector<int> & ExportPIDs, std::vector<int> & ExportLIDs);
00392 
00393 };
00394 
00395 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
00396 template<> inline std::vector<int>& LightweightCrsMatrix::getcolind() { return colind_; }
00397 #endif
00398 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
00399 template<> inline std::vector<long long>& LightweightCrsMatrix::getcolind() { return colind_LL_; }
00400 #endif
00401 
00402 }//namespace EpetraExt
00403 
00404 #endif
00405 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines