42#ifndef EPETRAEXT_MMHELPERS_H
43#define EPETRAEXT_MMHELPERS_H
46#include "Epetra_ConfigDefs.h"
47#include "Epetra_DistObject.h"
48#include "Epetra_Map.h"
49#include "Teuchos_RCP.hpp"
50#include "Epetra_Comm.h"
51#include "Epetra_Import.h"
52#include "Epetra_CrsMatrix.h"
53#include <Teuchos_TimeMonitor.hpp>
63class LightweightCrsMatrix;
117#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
123#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
124 virtual int InsertGlobalValues(
long long GlobalRow,
int NumEntries,
double* Values,
long long* Indices) = 0;
126 virtual int SumIntoGlobalValues(
long long GlobalRow,
int NumEntries,
double* Values,
long long* Indices) = 0;
140#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
145#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
146 int InsertGlobalValues(
long long GlobalRow,
int NumEntries,
double* Values,
long long* Indices);
147 int SumIntoGlobalValues(
long long GlobalRow,
int NumEntries,
double* Values,
long long* Indices);
155template<
typename int_type>
165 int InsertGlobalValues(int_type GlobalRow,
int NumEntries,
double* Values, int_type* Indices);
166 int SumIntoGlobalValues(int_type GlobalRow,
int NumEntries,
double* Values, int_type* Indices);
168 std::map<int_type,std::set<int_type>*>&
get_graph();
173 std::map<int_type,std::set<int_type>*> graph_;
179template<
typename int_type>
183template<
typename int_type>
185 const std::vector<int_type>& proc_col_ranges,
186 std::vector<int_type>& send_rows,
187 std::vector<int>& rows_per_send_proc);
189#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
191 const std::vector<int>& proc_col_ranges,
192 std::vector<int>& send_rows,
193 std::vector<int>& rows_per_send_proc);
196#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
198 const std::vector<long long>& proc_col_ranges,
199 std::vector<long long>& send_rows,
200 std::vector<int>& rows_per_send_proc);
203template<
typename int_type>
206template<
typename int_type>
216#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
220#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
233#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
236#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
245#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
250#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
251 int LID(
long long GID)
const;
256#if defined(EPETRA_NO_32BIT_GLOBAL_INDICES) && defined(EPETRA_NO_64BIT_GLOBAL_INDICES)
259 int LID(
long long GID)
const {
return -1; }
262#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
267 throw "EpetraExt::LightweightMap::IndexBase: IndexBase cannot fit an int.";
272#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
290#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
293#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
354#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
376 template <
typename ImportType,
typename int_type>
377 void Construct(
const Epetra_CrsMatrix & A, ImportType & RowImporter,
bool SortGhosts=
false,
const char * label=0);
381 int MakeColMapAndReindex(std::vector<int> owningPIDs,std::vector<GO> Gcolind,
bool SortGhosts=
false,
const char * label=0);
383 template<
typename int_type>
384 std::vector<int_type>& getcolind();
386 template<
typename ImportType,
typename int_type>
387 int PackAndPrepareReverseComm(
const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
388 std::vector<int> &ReverseSendSizes, std::vector<int_type> &ReverseSendBuffer);
390 template<
typename ImportType,
typename int_type>
391 int MakeExportLists(
const Epetra_CrsMatrix & SourceMatrix, ImportType & RowImporter,
392 std::vector<int> &ReverseRecvSizes,
const int_type *ReverseRecvBuffer,
393 std::vector<int> & ExportPIDs, std::vector<int> & ExportLIDs);
397#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
398template<>
inline std::vector<int>& LightweightCrsMatrix::getcolind() {
return colind_; }
400#ifndef EPETRA_NO_64BIT_GLOBAL_INDICES
401template<>
inline std::vector<long long>& LightweightCrsMatrix::getcolind() {
return colind_LL_; }
405template<
typename int_type>
409 const Epetra_Import * prototypeImporter=0,
bool SortGhosts=
false,
410 const char * label=0)
421#ifdef ENABLE_MMM_TIMINGS
423 if(label) tpref = std::string(label);
424 using Teuchos::TimeMonitor;
425 Teuchos::RCP<Teuchos::TimeMonitor> MM = Teuchos::rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(tpref + std::string(
"EpetraExt: MMM Ionly Setup"))));
428 Mview.deleteContents();
430 Mview.origMatrix = &M;
432 int numProcs = Mrowmap.Comm().NumProc();
433 Mview.numRows = targetMap.NumMyElements();
435 Mview.origRowMap = &(M.
RowMap());
436 Mview.rowMap = &targetMap;
437 Mview.colMap = &(M.
ColMap());
439 Mview.importColMap = NULL;
443 int N = Mview.rowMap->NumMyElements();
445 if(Mrowmap.SameAs(targetMap)) {
447 Mview.targetMapToOrigRow.resize(N);
448 for(i=0;i<N; i++) Mview.targetMapToOrigRow[i]=i;
451 else if(prototypeImporter && prototypeImporter->SourceMap().SameAs(M.
RowMap()) && prototypeImporter->TargetMap().SameAs(targetMap)){
452 numRemote = prototypeImporter->NumRemoteIDs();
454 Mview.targetMapToOrigRow.resize(N); Mview.targetMapToOrigRow.assign(N,-1);
455 Mview.targetMapToImportRow.resize(N); Mview.targetMapToImportRow.assign(N,-1);
457 const int * PermuteToLIDs = prototypeImporter->PermuteToLIDs();
458 const int * PermuteFromLIDs = prototypeImporter->PermuteFromLIDs();
459 const int * RemoteLIDs = prototypeImporter->RemoteLIDs();
461 for(i=0; i<prototypeImporter->NumSameIDs();i++)
462 Mview.targetMapToOrigRow[i] = i;
465 for(i=0; i<prototypeImporter->NumPermuteIDs();i++)
466 Mview.targetMapToOrigRow[PermuteToLIDs[i]] = PermuteFromLIDs[i];
468 for(i=0; i<prototypeImporter->NumRemoteIDs();i++)
469 Mview.targetMapToImportRow[RemoteLIDs[i]] = i;
473 throw std::runtime_error(
"import_only: This routine only works if you either have the right map or no prototypeImporter");
476 if (Mview.numRemote > 0) {
477 std::cerr <<
"EpetraExt::MatrixMatrix::Multiply ERROR, numProcs < 2 but "
478 <<
"attempting to import remote matrix rows."<<std::endl;
486#ifdef ENABLE_MMM_TIMINGS
487 MM = Teuchos::rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(tpref + std::string(
"EpetraExt: MMM Ionly Import-1"))));
489 const int * RemoteLIDs = prototypeImporter->RemoteLIDs();
492 int_type* MremoteRows = numRemote>0 ?
new int_type[prototypeImporter->NumRemoteIDs()] : 0;
493 for(i=0; i<prototypeImporter->NumRemoteIDs(); i++)
494 MremoteRows[i] = (int_type) targetMap.GID64(RemoteLIDs[i]);
496 LightweightMap MremoteRowMap((int_type) -1, numRemote, MremoteRows, (int_type)Mrowmap.IndexBase64());
498#ifdef ENABLE_MMM_TIMINGS
499 MM = Teuchos::rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(tpref + std::string(
"EpetraExt: MMM Ionly Import-2"))));
506#ifdef ENABLE_MMM_TIMINGS
507 MM = Teuchos::rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(tpref + std::string(
"EpetraExt: MMM Ionly Import-3"))));
512#ifdef ENABLE_MMM_TIMINGS
513 MM = Teuchos::rcp(
new TimeMonitor(*TimeMonitor::getNewTimer(tpref + std::string(
"EpetraExt: MMM Ionly Import-4"))));
516#ifdef ENABLE_MMM_STATISTICS
522 delete [] MremoteRows;
virtual ~CrsMatrixStruct()
const Epetra_CrsMatrix * origMatrix
const Epetra_Map * rowMap
std::vector< int > targetMapToImportRow
const Epetra_Map * colMap
const Epetra_Map * domainMap
LightweightCrsMatrix * importMatrix
const Epetra_BlockMap * importColMap
std::vector< int > targetMapToOrigRow
const Epetra_Map * origRowMap
CrsWrapper_Epetra_CrsMatrix(Epetra_CrsMatrix &epetracrsmatrix)
const Epetra_Map & RowMap() const
int SumIntoGlobalValues(int GlobalRow, int NumEntries, double *Values, int *Indices)
virtual ~CrsWrapper_Epetra_CrsMatrix()
int InsertGlobalValues(int GlobalRow, int NumEntries, double *Values, int *Indices)
std::map< int_type, std::set< int_type > * > & get_graph()
const Epetra_Map & RowMap() const
virtual ~CrsWrapper_GraphBuilder()
int InsertGlobalValues(int_type GlobalRow, int NumEntries, double *Values, int_type *Indices)
int SumIntoGlobalValues(int_type GlobalRow, int NumEntries, double *Values, int_type *Indices)
CrsWrapper_GraphBuilder(const Epetra_Map &emap)
virtual const Epetra_Map & RowMap() const =0
virtual int InsertGlobalValues(long long GlobalRow, int NumEntries, double *Values, long long *Indices)=0
virtual int SumIntoGlobalValues(long long GlobalRow, int NumEntries, double *Values, long long *Indices)=0
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, double *Values, int *Indices)=0
virtual int SumIntoGlobalValues(int GlobalRow, int NumEntries, double *Values, int *Indices)=0
Epetra_BlockMap * RowMapEP_
std::vector< double > vals_
std::vector< int > ExportPIDs_
LightweightMap * RowMapLW_
std::vector< long long > colind_LL_
std::vector< int > ColMapOwningPIDs_
LightweightCrsMatrix(const Epetra_CrsMatrix &A, RemoteOnlyImport &RowImporter, bool SortGhosts=false, const char *label=0)
std::vector< int > colind_
std::vector< int > ExportLIDs_
std::vector< int > rowptr_
std::vector< long long > MyGlobalElements_LL_
Epetra_HashTable< int > * LIDHash_int_
Epetra_HashTable< long long > * LIDHash_LL_
std::vector< int > MyGlobalElements_int_
bool GlobalIndicesInt() const
LightweightMap & operator=(const LightweightMap &map)
long long GID64(int LID) const
long long * MyGlobalElements64() const
int * MyGlobalElements() const
long long IndexBase64() const
bool GlobalIndicesLongLong() const
void MyGlobalElementsPtr(int *&MyGlobalElementList) const
int NumMyElements() const
const LightweightMap & TargetMap() const
Epetra_Distributor & Distributor()
RemoteOnlyImport(const Epetra_Import &Importer, LightweightMap &RemoteOnlyTargetMap)
const Epetra_BlockMap & SourceMap() const
const Epetra_Map & RowMap() const
const Epetra_Map & ColMap() const
const Epetra_Map & DomainMap() const
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.
void insert_matrix_locations(CrsWrapper_GraphBuilder< int_type > &graphbuilder, Epetra_CrsMatrix &C)
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)
void PrintMultiplicationStatistics(Epetra_Import *Transfer, const std::string &label)
void printMultiplicationStatistics(Epetra_Import *Transfer, const std::string &label)
int import_only(const Epetra_CrsMatrix &M, const Epetra_Map &targetMap, CrsMatrixStruct &Mview, const Epetra_Import *prototypeImporter)
std::pair< int_type, int_type > get_col_range(const Epetra_Map &emap)
int dumpCrsMatrixStruct(const CrsMatrixStruct &M)
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)