45 #include "Epetra_ConfigDefs.h" 52 #include "Epetra_MpiComm.h" 53 #include "Epetra_Map.h" 54 #include "Epetra_CrsMatrix.h" 59 RestrictedCrsMatrixWrapper::RestrictedCrsMatrixWrapper()
60 : proc_is_active(true),
61 subcomm_is_set(false),
62 MPI_SubComm_(MPI_COMM_NULL),
72 RestrictedCrsMatrixWrapper::~RestrictedCrsMatrixWrapper() {
75 delete RestrictedComm_;
80 int RestrictedCrsMatrixWrapper::SetMPISubComm(MPI_Comm MPI_SubComm){
82 MPI_SubComm_=MPI_SubComm;
delete RestrictedComm_; subcomm_is_set=
true;
89 template<
typename int_type>
90 int RestrictedCrsMatrixWrapper::Trestrict_comm(Teuchos::RCP<Epetra_CrsMatrix> input_matrix){
92 input_matrix_=input_matrix;
94 const Epetra_MpiComm *InComm =
dynamic_cast<const Epetra_MpiComm*
>(& input_matrix_->Comm());
95 const Epetra_Map *InRowMap=
dynamic_cast<const Epetra_Map*
>(& input_matrix_->RowMap());
96 const Epetra_Map *InColMap=
dynamic_cast<const Epetra_Map*
>(& input_matrix_->ColMap());
98 if(!InComm || !InRowMap || !InColMap)
return (-1);
100 int_type Nrows = (int_type) InRowMap->NumGlobalElements64();
101 int_type Ncols = (int_type) InColMap->NumGlobalElements64();
106 if(InRowMap->NumMyElements()) color=1;
107 else color=MPI_UNDEFINED;
108 MPI_Comm_split(InComm->Comm(),color,InComm->MyPID(),&MPI_SubComm_);
113 if (input_matrix->NumMyRows() && MPI_SubComm_ == MPI_COMM_NULL)
118 if(MPI_SubComm_ == MPI_COMM_NULL) proc_is_active=
false;
119 else proc_is_active=
true;
123 RestrictedComm_=
new Epetra_MpiComm(MPI_SubComm_);
125 int_type* RowMapGlobalElements = 0;
126 InRowMap->MyGlobalElementsPtr(RowMapGlobalElements);
127 int_type* ColMapGlobalElements = 0;
128 InColMap->MyGlobalElementsPtr(ColMapGlobalElements);
131 ResRowMap_ =
new Epetra_Map(Nrows,InRowMap->NumMyElements(),RowMapGlobalElements,
132 (int_type) InRowMap->IndexBase64(),*RestrictedComm_);
133 ResColMap_ =
new Epetra_Map(Ncols,InColMap->NumMyElements(),ColMapGlobalElements,
134 (int_type) InColMap->IndexBase64(),*RestrictedComm_);
140 restricted_matrix_= Teuchos::rcp(
new Epetra_CrsMatrix(View,*ResRowMap_,*ResColMap_,0));
141 for(
int i=0;i<input_matrix_->NumMyRows();i++) {
142 input_matrix_->ExtractMyRowView(i,Nr,values,colind);
143 restricted_matrix_->InsertMyValues(i,Nr,values,colind);
145 restricted_matrix_->FillComplete();
151 int RestrictedCrsMatrixWrapper::restrict_comm(Teuchos::RCP<Epetra_CrsMatrix> input_matrix)
153 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES 154 if(input_matrix->RowMap().GlobalIndicesInt()) {
155 return Trestrict_comm<int>(input_matrix);
159 #ifndef EPETRA_NO_64BIT_GLOBAL_INDICES 160 if(input_matrix->RowMap().GlobalIndicesLongLong()) {
161 return Trestrict_comm<long long>(input_matrix);
165 throw "EpetraExt::Trestrict_comm: ERROR, GlobalIndices type unknown.";
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.