Epetra_CrsSingletonFilter.h

Go to the documentation of this file.
00001 
00002 //@HEADER
00003 /*
00004 ************************************************************************
00005 
00006               Epetra: Linear Algebra Services Package 
00007                 Copyright (2001) Sandia Corporation
00008 
00009 Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00010 license for use of this work by or on behalf of the U.S. Government.
00011 
00012 This library is free software; you can redistribute it and/or modify
00013 it under the terms of the GNU Lesser General Public License as
00014 published by the Free Software Foundation; either version 2.1 of the
00015 License, or (at your option) any later version.
00016  
00017 This library is distributed in the hope that it will be useful, but
00018 WITHOUT ANY WARRANTY; without even the implied warranty of
00019 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00020 Lesser General Public License for more details.
00021  
00022 You should have received a copy of the GNU Lesser General Public
00023 License along with this library; if not, write to the Free Software
00024 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00025 USA
00026 Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00027 
00028 ************************************************************************
00029 */
00030 //@HEADER
00031 
00032 #ifndef EPETRA_CRSSINGLETONFILTER_H
00033 #define EPETRA_CRSSINGLETONFILTER_H
00034 
00035 #include "Epetra_Object.h"
00036 #include "Epetra_CrsMatrix.h"
00037 #include "Epetra_MapColoring.h"
00038 #include "Epetra_SerialDenseVector.h"
00039 class Epetra_LinearProblem;
00040 class Epetra_Map;
00041 class Epetra_MultiVector;
00042 class Epetra_Import;
00043 class Epetra_Export;
00044 class Epetra_IntVector;
00045 
00047 
00103 class Epetra_CrsSingletonFilter {
00104       
00105  public:
00106 
00108 
00109 
00110   Epetra_CrsSingletonFilter();
00111 
00113   virtual ~Epetra_CrsSingletonFilter();
00115 
00116 
00117 
00118 
00127   int Analyze(Epetra_RowMatrix * FullMatrix);
00128 
00130   bool SingletonsDetected() const {if (!AnalysisDone_) return(false); else return(RowMapColors_->MaxNumColors()>1);};
00132 
00134 
00135 
00136 
00141   int ConstructReducedProblem(Epetra_LinearProblem * Problem);
00142 
00144 
00149   int UpdateReducedProblem(Epetra_LinearProblem * Problem);
00150 
00152 
00153 
00154 
00155 
00160   int ComputeFullSolution();
00162 
00163 
00164 
00165   int NumRowSingletons() const {return(NumGlobalRowSingletons_);};
00166 
00168   int NumColSingletons() const {return(NumGlobalColSingletons_);};
00169 
00171 
00176   int NumSingletons() const {return(NumColSingletons()+NumRowSingletons());};
00177 
00179   double RatioOfDimensions() const {return(RatioOfDimensions_);};
00180 
00182   double RatioOfNonzeros() const {return(RatioOfNonzeros_);};
00183 
00185 
00186 
00187 
00189   Epetra_LinearProblem * FullProblem() const {return(FullProblem_);};
00190 
00192   Epetra_LinearProblem * ReducedProblem() const {return(ReducedProblem_);};
00193 
00195   Epetra_RowMatrix * FullMatrix() const {return(FullMatrix_);};
00196 
00198   Epetra_CrsMatrix * ReducedMatrix() const {return(ReducedMatrix_);};
00199 
00201   Epetra_MapColoring * RowMapColors() const {return(RowMapColors_);};
00202 
00204   Epetra_MapColoring * ColMapColors() const {return(ColMapColors_);};
00205 
00207   Epetra_Map * ReducedMatrixRowMap() const {return(ReducedMatrixRowMap_);};
00208 
00210   Epetra_Map * ReducedMatrixColMap() const {return(ReducedMatrixColMap_);};
00211 
00213   Epetra_Map * ReducedMatrixDomainMap() const {return(ReducedMatrixDomainMap_);};
00214 
00216   Epetra_Map * ReducedMatrixRangeMap() const {return(ReducedMatrixRangeMap_);};
00218 
00219  protected:
00220 
00221  
00222 
00223   //  This pointer will be zero if full matrix is not a CrsMatrix.
00224   Epetra_CrsMatrix * FullCrsMatrix() const {return(FullCrsMatrix_);};
00225 
00226   const Epetra_Map & FullMatrixRowMap() const {return(FullMatrix()->RowMatrixRowMap());};
00227   const Epetra_Map & FullMatrixColMap() const {return(FullMatrix()->RowMatrixColMap());};
00228   const Epetra_Map & FullMatrixDomainMap() const {return((FullMatrix()->OperatorDomainMap()));};
00229   const Epetra_Map & FullMatrixRangeMap() const {return((FullMatrix()->OperatorRangeMap()));};
00230   void InitializeDefaults();
00231   int ComputeEliminateMaps();
00232   int Setup(Epetra_LinearProblem * Problem);
00233   int InitFullMatrixAccess();
00234   int GetRow(int Row, int & NumIndices, int * & Indices);
00235   int GetRowGCIDs(int Row, int & NumIndices, double * & Values, int * & GlobalIndices);
00236   int GetRow(int Row, int & NumIndices, double * & Values, int * & Indices);
00237   int CreatePostSolveArrays(const Epetra_IntVector & RowIDs,
00238           const Epetra_MapColoring & RowMapColors,
00239           const Epetra_IntVector & ColProfiles,
00240           const Epetra_IntVector & NewColProfiles,
00241           const Epetra_IntVector & ColHasRowWithSingleton);
00242   
00243   int ConstructRedistributeExporter(Epetra_Map * SourceMap, Epetra_Map * TargetMap,
00244             Epetra_Export * & RedistributeExporter,
00245             Epetra_Map * & RedistributeMap);
00246   
00247   Epetra_LinearProblem * FullProblem_;
00248   Epetra_LinearProblem * ReducedProblem_;
00249   Epetra_RowMatrix * FullMatrix_;
00250   Epetra_CrsMatrix * FullCrsMatrix_;
00251   Epetra_CrsMatrix * ReducedMatrix_;
00252   Epetra_MultiVector * ReducedRHS_;
00253   Epetra_MultiVector * ReducedLHS_;
00254   
00255   Epetra_Map * ReducedMatrixRowMap_;
00256   Epetra_Map * ReducedMatrixColMap_;
00257   Epetra_Map * ReducedMatrixDomainMap_;
00258   Epetra_Map * ReducedMatrixRangeMap_;
00259   Epetra_Map * OrigReducedMatrixDomainMap_;
00260   Epetra_Import * Full2ReducedRHSImporter_;
00261   Epetra_Import * Full2ReducedLHSImporter_;
00262   Epetra_Export * RedistributeDomainExporter_;
00263   
00264   int * ColSingletonRowLIDs_;
00265   int * ColSingletonColLIDs_;
00266   int * ColSingletonPivotLIDs_;
00267   double * ColSingletonPivots_;
00268   
00269   
00270   int AbsoluteThreshold_;
00271   double RelativeThreshold_;
00272 
00273   int NumMyRowSingletons_;
00274   int NumMyColSingletons_;
00275   int NumGlobalRowSingletons_;
00276   int NumGlobalColSingletons_;
00277   double RatioOfDimensions_;
00278   double RatioOfNonzeros_;
00279   
00280   bool HaveReducedProblem_;
00281   bool UserDefinedEliminateMaps_;
00282   bool AnalysisDone_;
00283   bool SymmetricElimination_;
00284   
00285   Epetra_MultiVector * tempExportX_;
00286   Epetra_MultiVector * tempX_;
00287   Epetra_MultiVector * tempB_;
00288   Epetra_MultiVector * RedistributeReducedLHS_;
00289   int * Indices_;
00290   Epetra_SerialDenseVector Values_;
00291   
00292   Epetra_MapColoring * RowMapColors_;
00293   Epetra_MapColoring * ColMapColors_;
00294   bool FullMatrixIsCrsMatrix_;
00295   int MaxNumMyEntries_;
00296   
00297   
00298  private:
00300   Epetra_CrsSingletonFilter(const Epetra_CrsSingletonFilter & Problem);
00301   Epetra_CrsSingletonFilter & operator=(const Epetra_CrsSingletonFilter & Problem);
00302 };
00303 #endif /* EPETRA_CRSSINGLETONFILTER_H */

Generated on Thu Sep 18 12:37:57 2008 for Epetra Package Browser (Single Doxygen Collection) by doxygen 1.3.9.1