Ifpack_OverlapGraph.cpp

Go to the documentation of this file.
00001 /*@HEADER
00002 // ***********************************************************************
00003 // 
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2002) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 //@HEADER
00028 */
00029 
00030 #include "Ifpack_OverlapGraph.h"
00031 #include "Epetra_CrsGraph.h"
00032 #include "Epetra_RowMatrix.h"
00033 #include "Epetra_BlockMap.h"
00034 #include "Epetra_Map.h"
00035 #include "Epetra_Import.h"
00036 
00037 #ifdef HAVE_IFPACK_TEUCHOS
00038 #include <Teuchos_ParameterList.hpp>
00039 #include <ifp_parameters.h>
00040 #endif
00041 
00042 //==============================================================================
00043 Ifpack_OverlapGraph::Ifpack_OverlapGraph(const Epetra_CrsGraph * UserMatrixGraph, int OverlapLevel)
00044   : OverlapGraph_(0),
00045     UserMatrixGraph_(UserMatrixGraph),
00046     UserMatrix_(0),
00047     OverlapRowMap_(0),
00048     OverlapLevel_(OverlapLevel)
00049 {
00050   // Test for non-trivial overlap here so we can use it later.
00051   IsOverlapped_ = (OverlapLevel>0 && UserMatrixGraph->DomainMap().DistributedGlobal());
00052 
00053   ConstructOverlapGraph(UserMatrixGraph);
00054 
00055 }
00056 //==============================================================================
00057 Ifpack_OverlapGraph::Ifpack_OverlapGraph(const Epetra_RowMatrix * UserMatrix, int OverlapLevel)
00058   : OverlapGraph_(0),
00059     UserMatrixGraph_(0),
00060     UserMatrix_(UserMatrix),
00061     OverlapRowMap_(0),
00062     OverlapLevel_(OverlapLevel)
00063 {
00064   // Test for non-trivial overlap here so we can use it later.
00065   IsOverlapped_ = (OverlapLevel>0 && UserMatrix->OperatorDomainMap().DistributedGlobal());
00066   
00067   throw ReportError("This constructor is not implemented yet.  Need to add Epetra_SrcObject support to Epetra_Import/Export", -1);
00068 }
00069 //==============================================================================
00070 Ifpack_OverlapGraph::Ifpack_OverlapGraph(const Ifpack_OverlapGraph & Source)
00071   : OverlapGraph_(Source.OverlapGraph_),
00072     UserMatrixGraph_(Source.UserMatrixGraph_),
00073     UserMatrix_(Source.UserMatrix_),
00074     OverlapRowMap_(Source.OverlapRowMap_),
00075     OverlapLevel_(Source.OverlapLevel_),
00076     IsOverlapped_(Source.IsOverlapped_)
00077 {
00078   if (IsOverlapped_) {
00079     if (OverlapGraph_!=0) OverlapGraph_ = new Epetra_CrsGraph(*OverlapGraph_);
00080     if (OverlapRowMap_!=0) OverlapRowMap_ = new Epetra_BlockMap(*OverlapRowMap_);
00081   }
00082   
00083 }
00084 //==============================================================================
00085 Ifpack_OverlapGraph::~Ifpack_OverlapGraph() {
00086   
00087   if (IsOverlapped_) {
00088     if (OverlapGraph_!=0) delete OverlapGraph_;
00089     if (OverlapRowMap_!=0) delete OverlapRowMap_;
00090   }
00091 }
00092 
00093 #ifdef HAVE_IFPACK_TEUCHOS
00094 //==========================================================================
00095 int Ifpack_OverlapGraph::SetParameters(const Teuchos::ParameterList& parameterlist,
00096                bool cerr_warning_if_unused)
00097 {
00098   Ifpack::param_struct params;
00099   params.int_params[Ifpack::level_overlap-FIRST_INT_PARAM] = OverlapLevel_;
00100 
00101   Ifpack::set_parameters(parameterlist, params, cerr_warning_if_unused);
00102 
00103   OverlapLevel_ = params.int_params[Ifpack::level_overlap-FIRST_INT_PARAM];
00104   return(0);
00105 }
00106 #endif
00107 
00108 //==============================================================================
00109 int Ifpack_OverlapGraph::ConstructOverlapGraph(const Epetra_CrsGraph * UserMatrixGraph) {
00110 
00111   OverlapGraph_ = (Epetra_CrsGraph *) UserMatrixGraph;
00112   OverlapRowMap_ = (Epetra_BlockMap *) &UserMatrixGraph->RowMap();
00113 
00114   if (!IsOverlapped_) return(0); // All done
00115 
00116   Epetra_CrsGraph * OldGraph;
00117   Epetra_BlockMap * OldRowMap;
00118   Epetra_BlockMap * DomainMap = (Epetra_BlockMap *) &UserMatrixGraph->DomainMap();
00119   Epetra_BlockMap * RangeMap = (Epetra_BlockMap *) &UserMatrixGraph->RangeMap();
00120 
00121   for (int level=1; level <= OverlapLevel_; level++) {
00122     OldGraph = OverlapGraph_;
00123     OldRowMap = OverlapRowMap_;
00124 
00125     OverlapImporter_ = (Epetra_Import *) OldGraph->Importer();
00126     OverlapRowMap_ = new Epetra_BlockMap(OverlapImporter_->TargetMap());
00127 
00128     if (level<OverlapLevel_)
00129       OverlapGraph_ = new Epetra_CrsGraph(Copy, *OverlapRowMap_, 0);
00130     else
00131       // On last iteration, we want to filter out all columns except those that correspond
00132       // to rows in the graph.  This assures that our matrix is square
00133       OverlapGraph_ = new Epetra_CrsGraph(Copy, *OverlapRowMap_, *OverlapRowMap_, 0);
00134 
00135     EPETRA_CHK_ERR(OverlapGraph_->Import( *UserMatrixGraph, *OverlapImporter_, Insert));
00136     if (level<OverlapLevel_) {
00137       EPETRA_CHK_ERR(OverlapGraph_->FillComplete(*DomainMap, *RangeMap));
00138     }
00139     else {
00140       // Copy last OverlapImporter because we will use it later
00141       OverlapImporter_ = new Epetra_Import(*OverlapRowMap_, *DomainMap);
00142       EPETRA_CHK_ERR(OverlapGraph_->FillComplete(*DomainMap, *RangeMap));
00143     }
00144 
00145     if (level>1) {
00146       delete OldGraph;
00147       delete OldRowMap;
00148     }
00149   }
00150 
00151   return(0);
00152 }
00153 

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