EpetraExt  Development
EpetraExt_MMHelpers.h
Go to the documentation of this file.
1 //@HEADER
2 // ***********************************************************************
3 //
4 // EpetraExt: Epetra Extended - Linear Algebra Services Package
5 // Copyright (2011) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 
42 #ifndef EPETRAEXT_MMHELPERS_H
43 #define EPETRAEXT_MMHELPERS_H
44 
45 #include "EpetraExt_ConfigDefs.h"
46 #include "Epetra_ConfigDefs.h"
47 #include "Epetra_DistObject.h"
48 #include "Epetra_Map.h"
49 #include "Teuchos_RCP.hpp"
50 
51 #include <vector>
52 #include <set>
53 #include <map>
54 
55 
56 class Epetra_CrsMatrix;
57 class Epetra_Import;
58 class Epetra_Distributor;
59 
60 namespace EpetraExt {
61 class LightweightCrsMatrix;
62 
63 //#define HAVE_EPETRAEXT_DEBUG // for extra sanity checks
64 
65 // ==============================================================
66 //struct that holds views of the contents of a CrsMatrix. These
67 //contents may be a mixture of local and remote rows of the
68 //actual matrix.
70 public:
72 
73  virtual ~CrsMatrixStruct();
74 
75  void deleteContents();
76 
77  int numRows;
78 
79  // The following class members get used in the transpose modes of the MMM
80  // but not in the A*B mode.
82  int** indices;
83  double** values;
84  bool* remote;
85  int numRemote;
86  const Epetra_BlockMap* importColMap;
87 
88  // Maps and matrices (all modes)
89  const Epetra_Map* origRowMap;
90  const Epetra_Map* rowMap;
91  const Epetra_Map* colMap;
92  const Epetra_Map* domainMap;
94  const Epetra_CrsMatrix *origMatrix;
95 
96  // The following class members are only used for A*B mode
97  std::vector<int> targetMapToOrigRow;
98  std::vector<int> targetMapToImportRow;
99 };
100 
102 
103 // ==============================================================
104 class CrsWrapper {
105  public:
106  virtual ~CrsWrapper(){}
107 
108  virtual const Epetra_Map& RowMap() const = 0;
109 
110  virtual bool Filled() = 0;
111 
112 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
113  virtual int InsertGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices) = 0;
114 
115  virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices) = 0;
116 #endif
117 
118 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
119  virtual int InsertGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices) = 0;
120 
121  virtual int SumIntoGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices) = 0;
122 #endif
123 };
124 
125 // ==============================================================
127  public:
128  CrsWrapper_Epetra_CrsMatrix(Epetra_CrsMatrix& epetracrsmatrix);
129  virtual ~CrsWrapper_Epetra_CrsMatrix();
130 
131  const Epetra_Map& RowMap() const;
132 
133  bool Filled();
134 
135 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
136  int InsertGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices);
137  int SumIntoGlobalValues(int GlobalRow, int NumEntries, double* Values, int* Indices);
138 #endif
139 
140 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
141  int InsertGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices);
142  int SumIntoGlobalValues(long long GlobalRow, int NumEntries, double* Values, long long* Indices);
143 #endif
144 
145  private:
146  Epetra_CrsMatrix& ecrsmat_;
147 };
148 
149 // ==============================================================
150 template<typename int_type>
152  public:
153  CrsWrapper_GraphBuilder(const Epetra_Map& emap);
154  virtual ~CrsWrapper_GraphBuilder();
155 
156  const Epetra_Map& RowMap() const {return rowmap_; }
157 
158  bool Filled();
159 
160  int InsertGlobalValues(int_type GlobalRow, int NumEntries, double* Values, int_type* Indices);
161  int SumIntoGlobalValues(int_type GlobalRow, int NumEntries, double* Values, int_type* Indices);
162 
163  std::map<int_type,std::set<int_type>*>& get_graph();
164 
165  int get_max_row_length() { return max_row_length_; }
166 
167  private:
168  std::map<int_type,std::set<int_type>*> graph_;
169  const Epetra_Map& rowmap_;
170  int max_row_length_;
171 };
172 
173 // ==============================================================
174 template<typename int_type>
176  Epetra_CrsMatrix& C);
177 
178 template<typename int_type>
179 void Tpack_outgoing_rows(const Epetra_CrsMatrix& mtx,
180  const std::vector<int_type>& proc_col_ranges,
181  std::vector<int_type>& send_rows,
182  std::vector<int>& rows_per_send_proc);
183 
184 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
185 void pack_outgoing_rows(const Epetra_CrsMatrix& mtx,
186  const std::vector<int>& proc_col_ranges,
187  std::vector<int>& send_rows,
188  std::vector<int>& rows_per_send_proc);
189 #endif
190 
191 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
192 void pack_outgoing_rows(const Epetra_CrsMatrix& mtx,
193  const std::vector<long long>& proc_col_ranges,
194  std::vector<long long>& send_rows,
195  std::vector<int>& rows_per_send_proc);
196 #endif
197 
198 template<typename int_type>
199 std::pair<int_type,int_type> get_col_range(const Epetra_Map& emap);
200 
201 template<typename int_type>
202 std::pair<int_type,int_type> get_col_range(const Epetra_CrsMatrix& mtx);
203 
204 // ==============================================================
205 class LightweightMapData : Epetra_Data {
206  friend class LightweightMap;
207  public:
210  long long IndexBase_;
211 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
212  std::vector<int> MyGlobalElements_int_;
213  Epetra_HashTable<int> * LIDHash_int_;
214 #endif
215 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
216  std::vector<long long> MyGlobalElements_LL_;
217  Epetra_HashTable<long long> * LIDHash_LL_;
218 #endif
219 
220  // For "copy" constructor only...
221  Epetra_Map * CopyMap_;
222 
223  };
224 
226  public:
227  LightweightMap();
228 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
229  LightweightMap(int NumGlobalElements,int NumMyElements, const int * MyGlobalElements, int IndexBase, bool GenerateHash=true);
230 #endif
231 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
232  LightweightMap(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, int IndexBase, bool GenerateHash=true);
233  LightweightMap(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
234 #endif
235  LightweightMap(const Epetra_Map & Map);
236  LightweightMap(const LightweightMap & Map);
237  ~LightweightMap();
238 
239  LightweightMap & operator=(const LightweightMap & map);
240 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
241  int LID(int GID) const;
242  int GID(int LID) const;
243 #endif
244 
245 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
246  int LID(long long GID) const;
247 #endif
248  long long GID64(int LID) const;
249  int NumMyElements() const;
250 
251 #if defined(EPETRA_NO_32BIT_GLOBAL_INDICES) && defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
252  // default implementation so that no compiler/linker error in case neither 32 nor 64
253  // bit indices present.
254  int LID(long long GID) const { return -1; }
255 #endif
256 
257 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
258  int* MyGlobalElements() const;
259  int IndexBase() const {
260  if(IndexBase64() == (long long) static_cast<int>(IndexBase64()))
261  return (int) IndexBase64();
262  throw "EpetraExt::LightweightMap::IndexBase: IndexBase cannot fit an int.";
263  }
264  void MyGlobalElementsPtr(int *& MyGlobalElementList) const;
265 #endif
266 
267 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
268  long long* MyGlobalElements64() const;
269  void MyGlobalElementsPtr(long long *& MyGlobalElementList) const;
270 #endif
271  long long IndexBase64() const {return Data_->IndexBase_;}
272 
273  int MinLID() const;
274  int MaxLID() const;
275 
276  bool GlobalIndicesInt() const { return IsInt; }
277  bool GlobalIndicesLongLong() const { return IsLongLong; }
278  private:
279  void CleanupData();
280  LightweightMapData *Data_;
281  bool IsLongLong;
282  bool IsInt;
283  //Epetra_BlockMapData* Data_;
284 
285 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
286  void Construct_int(int NumGlobalElements,int NumMyElements, const int * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
287 #endif
288 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
289  void Construct_LL(long long NumGlobalElements,int NumMyElements, const long long * MyGlobalElements, long long IndexBase, bool GenerateHash=true);
290 #endif
291 };
292 
293 
294 // ==============================================================
296  public:
297  RemoteOnlyImport(const Epetra_Import & Importer, LightweightMap & RemoteOnlyTargetMap);
298  ~RemoteOnlyImport();
299 
300  int NumSameIDs() {return 0;}
301 
302  int NumPermuteIDs() {return 0;}
303 
304  int NumRemoteIDs() {return NumRemoteIDs_;}
305 
306  int NumExportIDs() {return NumExportIDs_;}
307 
308  int* ExportLIDs() {return ExportLIDs_;}
309 
310  int* ExportPIDs() {return ExportPIDs_;}
311 
312  int* RemoteLIDs() {return RemoteLIDs_;}
313 
314  int* PermuteToLIDs() {return 0;}
315 
316  int* PermuteFromLIDs() {return 0;}
317 
318  int NumSend() {return NumSend_;}
319 
320  Epetra_Distributor & Distributor() {return *Distor_;}
321 
322  const Epetra_BlockMap & SourceMap() const {return *SourceMap_;}
323  const LightweightMap & TargetMap() const {return *TargetMap_;}
324 
325  private:
326  int NumSend_;
327  int NumRemoteIDs_;
328  int NumExportIDs_;
329  int * ExportLIDs_;
330  int * ExportPIDs_;
331  int * RemoteLIDs_;
332  Epetra_Distributor* Distor_;
333  const Epetra_BlockMap* SourceMap_;
334  const LightweightMap *TargetMap_;
335 };
336 
337 // ==============================================================
339  public:
340  LightweightCrsMatrix(const Epetra_CrsMatrix & A, RemoteOnlyImport & RowImporter);
341  LightweightCrsMatrix(const Epetra_CrsMatrix & A, Epetra_Import & RowImporter);
343 
344  // Standard crs data structures
345  std::vector<int> rowptr_;
346  std::vector<int> colind_;
347  std::vector<double> vals_;
348 
349 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
350  // Colind in LL-GID space (if needed)
351  std::vector<long long> colind_LL_;
352 #endif
353 
354  // Epetra Maps
355  bool use_lw;
357  Epetra_BlockMap *RowMapEP_;
359  Epetra_Map DomainMap_;
360 
361 
362  // List of owning PIDs (from the DomainMap) as ordered by entries in the column map.
363  std::vector<int> ColMapOwningPIDs_;
364 
365  // For building the final importer for C
366  std::vector<int> ExportLIDs_;
367  std::vector<int> ExportPIDs_;
368 
369  private:
370 
371  template <typename ImportType, typename int_type>
372  void Construct(const Epetra_CrsMatrix & A, ImportType & RowImporter);
373 
374  // Templated versions of MakeColMapAndReindex (to prevent code duplication)
375  template <class GO>
376  int MakeColMapAndReindex(std::vector<int> owningPIDs,std::vector<GO> Gcolind);
377 
378  template<typename int_type>
379  std::vector<int_type>& getcolind();
380 
381  template<typename ImportType, typename int_type>
382  int PackAndPrepareReverseComm(const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
383  std::vector<int> &ReverseSendSizes, std::vector<int_type> &ReverseSendBuffer);
384 
385  template<typename ImportType, typename int_type>
386  int MakeExportLists(const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
387  std::vector<int> &ReverseRecvSizes, const int_type *ReverseRecvBuffer,
388  std::vector<int> & ExportPIDs, std::vector<int> & ExportLIDs);
389 
390 };
391 
392 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
393 template<> inline std::vector<int>& LightweightCrsMatrix::getcolind() { return colind_; }
394 #endif
395 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
396 template<> inline std::vector<long long>& LightweightCrsMatrix::getcolind() { return colind_LL_; }
397 #endif
398 
399 }//namespace EpetraExt
400 
401 #endif
402 
LightweightCrsMatrix * importMatrix
Epetra_HashTable< long long > * LIDHash_LL_
std::vector< int > targetMapToOrigRow
void pack_outgoing_rows(const Epetra_CrsMatrix &mtx, const std::vector< int > &proc_col_ranges, std::vector< int > &send_rows, std::vector< int > &rows_per_send_proc)
const Epetra_CrsMatrix * origMatrix
const Epetra_Map & RowMap() const
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.
std::vector< int > targetMapToImportRow
std::vector< int > MyGlobalElements_int_
const Epetra_BlockMap & SourceMap() const
long long IndexBase64() const
const LightweightMap & TargetMap() const
void Tpack_outgoing_rows(const Epetra_CrsMatrix &mtx, const std::vector< int_type > &proc_col_ranges, std::vector< int_type > &send_rows, std::vector< int > &rows_per_send_proc)
Epetra_HashTable< int > * LIDHash_int_
std::vector< long long > MyGlobalElements_LL_
Epetra_Distributor & Distributor()
const Epetra_BlockMap * importColMap
void insert_matrix_locations(CrsWrapper_GraphBuilder< int_type > &graphbuilder, Epetra_CrsMatrix &C)
int dumpCrsMatrixStruct(const CrsMatrixStruct &M)
std::vector< long long > colind_LL_
std::pair< int_type, int_type > get_col_range(const Epetra_Map &emap)