EpetraExt Package Browser (Single Doxygen Collection) Development
EpetraExt_RestrictedCrsMatrixWrapper.cpp
Go to the documentation of this file.
00001 #include "EpetraExt_ConfigDefs.h"
00002 
00003 
00004 #ifdef HAVE_MPI
00005 
00006 
00007 #include "EpetraExt_RestrictedCrsMatrixWrapper.h"
00008 #include "Epetra_MpiComm.h"
00009 #include "Epetra_Map.h"
00010 #include "Epetra_CrsMatrix.h"
00011 
00012 
00013 namespace EpetraExt{
00014 
00015 RestrictedCrsMatrixWrapper::RestrictedCrsMatrixWrapper()
00016   : proc_is_active(true),
00017     subcomm_is_set(false),
00018     MPI_SubComm_(MPI_COMM_NULL),
00019     RestrictedComm_(0),
00020     ResRowMap_(0),
00021     ResColMap_(0),
00022     input_matrix_(),
00023     restricted_matrix_()  
00024 {
00025   
00026 }
00027 
00028 RestrictedCrsMatrixWrapper::~RestrictedCrsMatrixWrapper() {
00029   delete ResRowMap_;
00030   delete ResColMap_;
00031   delete RestrictedComm_;
00032 }
00033 
00034 
00035 
00036 int RestrictedCrsMatrixWrapper::SetMPISubComm(MPI_Comm MPI_SubComm){
00037   if(!subcomm_is_set){
00038     MPI_SubComm_=MPI_SubComm; delete RestrictedComm_; subcomm_is_set=true;
00039     return 0;
00040   }
00041   else return -1;
00042 }
00043 
00044 
00045 
00046 int RestrictedCrsMatrixWrapper::restrict_comm(Teuchos::RCP<Epetra_CrsMatrix> input_matrix){
00047   /* Pull the Matrix Info */
00048   input_matrix_=input_matrix;
00049   
00050   const Epetra_MpiComm *InComm = dynamic_cast<const Epetra_MpiComm*>(& input_matrix_->Comm());
00051   const Epetra_Map *InRowMap= dynamic_cast<const Epetra_Map* >(& input_matrix_->RowMap());
00052   const Epetra_Map *InColMap= dynamic_cast<const Epetra_Map* >(& input_matrix_->ColMap());
00053 
00054   if(!InComm || !InRowMap || !InColMap) return (-1);
00055   
00056   int Nrows=InRowMap->NumGlobalElements();
00057   int Ncols=InColMap->NumGlobalElements();
00058   
00059   if(!subcomm_is_set){
00060     /* Build the Split Communicators, If Needed */
00061     int color;
00062     if(InRowMap->NumMyElements()) color=1;
00063     else color=MPI_UNDEFINED;
00064     MPI_Comm_split(InComm->Comm(),color,InComm->MyPID(),&MPI_SubComm_);
00065   }
00066   else{
00067     /* Sanity check user-provided subcomm - drop an error if the MPISubComm
00068        does not include a processor with data. */
00069     if (input_matrix->NumMyRows() && MPI_SubComm_ == MPI_COMM_NULL)
00070       return(-2);
00071   }
00072 
00073   /* Mark active processors */
00074   if(MPI_SubComm_ == MPI_COMM_NULL) proc_is_active=false;
00075   else proc_is_active=true;
00076   
00077 
00078   if(proc_is_active){      
00079     RestrictedComm_=new Epetra_MpiComm(MPI_SubComm_);
00080     
00081     /* Build the Restricted Maps */
00082     ResRowMap_ = new Epetra_Map(Nrows,InRowMap->NumMyElements(),InRowMap->MyGlobalElements(),
00083                                 InRowMap->IndexBase(),*RestrictedComm_);
00084     ResColMap_ = new Epetra_Map(Ncols,InColMap->NumMyElements(),InColMap->MyGlobalElements(),
00085                                 InColMap->IndexBase(),*RestrictedComm_);        
00086     
00087     int *rowptr,*colind,Nr;
00088     double *values;
00089     
00090     /* Allocate the Restricted Matrix */
00091     restricted_matrix_= Teuchos::rcp(new Epetra_CrsMatrix(View,*ResRowMap_,*ResColMap_,0));
00092     for(int i=0;i<input_matrix_->NumMyRows();i++) {
00093       input_matrix_->ExtractMyRowView(i,Nr,values,colind);
00094       restricted_matrix_->InsertMyValues(i,Nr,values,colind);
00095     }
00096     restricted_matrix_->FillComplete();      
00097   }
00098 
00099   return 0;
00100 }/*end restrict_comm*/
00101 
00102 
00103 
00104 }  
00105 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines