snl_fei_Broker_LinSysCore.hpp

00001 /*--------------------------------------------------------------------*/
00002 /*    Copyright 2005 Sandia Corporation.                              */
00003 /*    Under the terms of Contract DE-AC04-94AL85000, there is a       */
00004 /*    non-exclusive license for use of this work by or on behalf      */
00005 /*    of the U.S. Government.  Export of this program may require     */
00006 /*    a license from the United States Government.                    */
00007 /*--------------------------------------------------------------------*/
00008 
00009 #ifndef _snl_fei_Broker_LinSysCore_hpp_
00010 #define _snl_fei_Broker_LinSysCore_hpp_
00011 
00012 #include <fei_macros.hpp>
00013 #include <fei_mpi.h>
00014 #include <fei_CommUtils.hpp>
00015 #include <snl_fei_Broker.hpp>
00016 #include <fei_LinearSystemCore.hpp>
00017 #include <fei_VectorSpace.hpp>
00018 #include <fei_Lookup_Impl.hpp>
00019 #include <fei_MatrixGraph.hpp>
00020 #include <fei_SparseRowGraph.hpp>
00021 #include <fei_Vector_Impl.hpp>
00022 #include <fei_Matrix_Impl.hpp>
00023 #include <fei_MatrixReducer.hpp>
00024 #include <fei_VectorReducer.hpp>
00025 #include <fei_Reducer.hpp>
00026 #include <snl_fei_LinearSystem_General.hpp>
00027 
00028 #undef fei_file
00029 #define fei_file "snl_fei_Broker_LinSysCore.hpp"
00030 #include <fei_ErrMacros.hpp>
00031 
00032 namespace snl_fei {
00033 
00037   class Broker_LinSysCore : public snl_fei::Broker {
00038   public:
00040     Broker_LinSysCore(fei::SharedPtr<LinearSystemCore> lsc,
00041                       fei::SharedPtr<fei::MatrixGraph> matrixGraph,
00042                       fei::SharedPtr<fei::Reducer> reducer,
00043                       bool blockMatrix);
00044 
00046     virtual ~Broker_LinSysCore();
00047 
00057     fei::SharedPtr<fei::Vector> createVector(bool isSolutionVector=false)
00058       {
00059         fei::SharedPtr<fei::Vector> vptr;
00060         if (matrixGraph_.get() == NULL) return(vptr);
00061         if (setGlobalOffsets() != 0) return(vptr);
00062 
00063         fei::SharedPtr<fei::Vector> tmpvec;
00064         tmpvec.reset(new fei::Vector_Impl<LinearSystemCore >(matrixGraph_->getRowSpace(),
00065                                                        linsyscore_.get(),
00066                                                        numLocalEqns_,
00067                                                        isSolutionVector));
00068 
00069         fei::SharedPtr<fei::Vector> vec;
00070         if (reducer_.get() != NULL) {
00071           vec.reset(new fei::VectorReducer(reducer_, tmpvec, isSolutionVector));
00072         }
00073         else {
00074           vec = tmpvec;
00075         }
00076 
00077         return(vec);
00078       }
00079 
00082     fei::SharedPtr<fei::Matrix> createMatrix()
00083       {
00084         fei::SharedPtr<fei::Matrix> mptr;
00085         if (matrixGraph_.get() == NULL) return(mptr);
00086 
00087         if (setMatrixStructure() != 0) return(mptr);
00088 
00089         fei::SharedPtr<fei::Matrix> tmpmat;
00090         tmpmat.reset(new fei::Matrix_Impl<LinearSystemCore>(linsyscore_,
00091                                                          matrixGraph_,
00092                                                            numLocalEqns_));
00093         fei::SharedPtr<fei::Matrix> matptr;
00094         if (reducer_.get() != NULL) {
00095           matptr.reset(new fei::MatrixReducer(reducer_, tmpmat));
00096         }
00097         else {
00098           matptr = tmpmat;
00099         }
00100 
00101         return(matptr);
00102       }
00103 
00106     fei::SharedPtr<fei::LinearSystem> createLinearSystem()
00107       {
00108         fei::SharedPtr<fei::LinearSystem> lsptr;
00109         if (matrixGraph_.get() == NULL) return(lsptr);
00110 
00111         if (setMatrixStructure() != 0) return(lsptr);
00112 
00113         fei::SharedPtr<fei::LinearSystem>
00114           linSys(new snl_fei::LinearSystem_General(matrixGraph_));
00115         return(linSys);
00116       }
00117 
00119     void setMatrixGraph(fei::SharedPtr<fei::MatrixGraph> matrixGraph)
00120       { matrixGraph_ = matrixGraph; }
00121 
00122   private:
00123     int setGlobalOffsets()
00124       {
00125         //only set the global offsets once.
00126         if (setGlobalOffsets_) return(0);
00127 
00128         if (matrixGraph_.get() == NULL) return(-1);
00129 
00130         MPI_Comm comm = matrixGraph_->getRowSpace()->getCommunicator();
00131         int num_procs = fei::numProcs(comm);
00132         int local_proc = fei::localProc(comm);
00133 
00134         std::vector<int> globalOffsets;
00135         std::vector<int> globalBlkOffsets;
00136 
00137         if (reducer_.get() != NULL) {
00138           int localsize = reducer_->getLocalReducedEqns().size();
00139           numLocalEqns_ = localsize;
00140           std::vector<int> lsizes(num_procs, 0);
00141           std::vector<int> gsizes(num_procs, 0);
00142           lsizes[local_proc] = localsize;
00143           fei::GlobalMax(comm, lsizes, gsizes);
00144           globalOffsets.resize(num_procs+1);
00145           int offset = 0;
00146           for(int p=0; p<num_procs; ++p) {
00147             globalOffsets[p] = offset;
00148             offset += gsizes[p];
00149           }
00150           globalOffsets[num_procs] = offset;
00151           globalBlkOffsets = globalOffsets;
00152         }
00153         else {
00154           fei::SharedPtr<fei::VectorSpace> vecSpace = 
00155             matrixGraph_->getRowSpace();
00156 
00157           vecSpace->getGlobalIndexOffsets(globalOffsets);
00158           vecSpace->getGlobalBlkIndexOffsets(globalBlkOffsets);
00159 
00160           numLocalEqns_ = globalOffsets[local_proc+1]-globalOffsets[local_proc];
00161         }
00162 
00163         CHK_ERR(linsyscore_->setGlobalOffsets(num_procs+1,
00164                                               &globalBlkOffsets[0],
00165                                               &globalOffsets[0],
00166                                               &globalBlkOffsets[0]));
00167 
00168         setGlobalOffsets_ = true;
00169         return(0);
00170       }
00171 
00172     int setMatrixStructure()
00173       {
00174         //only set the matrix matrixGraph once.
00175         if (setMatrixStructure_) return(0);
00176 
00177         if (matrixGraph_.get() == NULL) return(-1);
00178 
00179         lookup_ = new fei::Lookup_Impl(matrixGraph_, 0);
00180 
00181         CHK_ERR( linsyscore_->setLookup(*lookup_) );
00182 
00183         CHK_ERR( setGlobalOffsets() );
00184 
00185         MPI_Comm comm = matrixGraph_->getRowSpace()->getCommunicator();
00186 
00187         fei::SharedPtr<fei::SparseRowGraph> localSRGraph =
00188           matrixGraph_->createGraph(blockMatrix_);
00189 
00190         std::vector<int>& rowNumbers = localSRGraph->rowNumbers;
00191         int numLocalRows = rowNumbers.size();
00192         int* rowOffsets = &(localSRGraph->rowOffsets[0]);
00193         int numLocalNonzeros = localSRGraph->packedColumnIndices.size();
00194         int* nonzeros = &(localSRGraph->packedColumnIndices[0]);
00195 
00196         int numGlobalNonzeros = 0;
00197         fei::GlobalSum(comm, numLocalNonzeros, numGlobalNonzeros);
00198 
00199         std::vector<int*> colPtrs(numLocalRows);
00200         std::vector<int> ptRowsPerBlkRow(numLocalRows, 1);
00201         std::vector<int> rowLengths(numLocalRows);
00202         int* rowLengthsPtr = &rowLengths[0];
00203 
00204         for(int i=0; i<numLocalRows; ++i) {
00205           colPtrs[i] = &(nonzeros[rowOffsets[i]]);
00206           rowLengthsPtr[i] = rowOffsets[i+1]-rowOffsets[i];
00207           if (blockMatrix_ == true) {
00208             ptRowsPerBlkRow[i] = lookup_->getBlkEqnSize(rowNumbers[i]);
00209           }
00210         }
00211 
00212         CHK_ERR( linsyscore_->setMatrixStructure(&colPtrs[0],
00213                                                  rowLengthsPtr,
00214                                                  &colPtrs[0],
00215                                                  rowLengthsPtr,
00216                                                  &ptRowsPerBlkRow[0]));
00217 
00218         setMatrixStructure_ = true;
00219 
00220         return(0);
00221       }
00222 
00223 
00224     fei::SharedPtr<LinearSystemCore> linsyscore_;
00225     fei::SharedPtr<fei::MatrixGraph> matrixGraph_;
00226     fei::SharedPtr<fei::Reducer> reducer_;
00227     Lookup* lookup_;
00228 
00229     bool setGlobalOffsets_;
00230     int numLocalEqns_;
00231     bool setMatrixStructure_;
00232     bool blockMatrix_;
00233   };//class Broker_LinSysCore
00234 }//namespace snl_fei
00235 
00236 #endif // _snl_fei_Broker_LinSysCore_hpp_
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Friends
Generated on Wed Apr 13 10:08:24 2011 for FEI by  doxygen 1.6.3