Ifpack_AdditiveSchwarz.h

00001 #ifndef IFPACK_ADDITIVESCHWARZ_H
00002 #define IFPACK_ADDITIVESCHWARZ_H
00003 
00004 #include "Ifpack_ConfigDefs.h"
00005 #include "Ifpack_Preconditioner.h"
00006 #include "Ifpack_ConfigDefs.h"
00007 #include "Ifpack_Preconditioner.h"
00008 #include "Ifpack_Reordering.h"
00009 #include "Ifpack_RCMReordering.h"
00010 #include "Ifpack_METISReordering.h"
00011 #include "Ifpack_LocalFilter.h"
00012 #include "Ifpack_SingletonFilter.h"
00013 #include "Ifpack_ReorderFilter.h"
00014 #include "Ifpack_Utils.h"
00015 #include "Ifpack_OverlappingRowMatrix.h"
00016 #include "Epetra_CombineMode.h"
00017 #include "Epetra_MultiVector.h"
00018 #include "Epetra_Map.h"
00019 #include "Epetra_Comm.h"
00020 #include "Epetra_Time.h"
00021 #include "Epetra_LinearProblem.h"
00022 #include "Epetra_RowMatrix.h"
00023 #include "Epetra_CrsMatrix.h"
00024 #include "Teuchos_ParameterList.hpp"
00025 #include "Teuchos_RefCountPtr.hpp"
00026 
00028 
00081 template<typename T>
00082 class Ifpack_AdditiveSchwarz : public virtual Ifpack_Preconditioner {
00083       
00084 public:
00085 
00087 
00088 
00099   Ifpack_AdditiveSchwarz(Epetra_RowMatrix* Matrix_in,
00100              int OverlapLevel_in = 0);
00101   
00103   virtual ~Ifpack_AdditiveSchwarz() {};
00105 
00107 
00109 
00118     virtual int SetUseTranspose(bool UseTranspose_in);
00120   
00122 
00124 
00134     virtual int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00135 
00137 
00148     virtual int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
00149 
00151     virtual double NormInf() const;
00153   
00155 
00157     virtual const char * Label() const;
00158 
00160     virtual bool UseTranspose() const;
00161 
00163     virtual bool HasNormInf() const;
00164 
00166     virtual const Epetra_Comm & Comm() const;
00167 
00169     virtual const Epetra_Map & OperatorDomainMap() const;
00170 
00172     virtual const Epetra_Map & OperatorRangeMap() const;
00174 
00176   virtual bool IsInitialized() const
00177   {
00178     return(IsInitialized_);
00179   }
00180 
00182   virtual bool IsComputed() const
00183   {
00184     return(IsComputed_);
00185   }
00186 
00188 
00207   virtual int SetParameters(Teuchos::ParameterList& List);
00208 
00209   // @}
00210 
00211   // @{ Query methods
00212   
00214   virtual int Initialize();
00215 
00217   virtual int Compute();
00218 
00220   virtual double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
00221                          const int MaxIters = 1550,
00222                          const double Tol = 1e-9,
00223              Epetra_RowMatrix* Matrix_in = 0);
00224 
00226   virtual double Condest() const
00227   {
00228     return(Condest_);
00229   }
00230 
00232   virtual const Epetra_RowMatrix& Matrix() const
00233   {
00234     return(*Matrix_);
00235   }
00236 
00238   virtual bool IsOverlapping() const
00239   {
00240     return(IsOverlapping_);
00241   }
00242 
00244   virtual std::ostream& Print(std::ostream&) const;
00245   
00246   virtual const T* Inverse() const
00247   {
00248     return(&*Inverse_);
00249   }
00250 
00252   virtual int NumInitialize() const
00253   {
00254     return(NumInitialize_);
00255   }
00256 
00258   virtual int NumCompute() const
00259   {
00260     return(NumCompute_);
00261   }
00262 
00264   virtual int NumApplyInverse() const
00265   {
00266     return(NumApplyInverse_);
00267   }
00268 
00270   virtual double InitializeTime() const
00271   {
00272     return(InitializeTime_);
00273   }
00274 
00276   virtual double ComputeTime() const
00277   {
00278     return(ComputeTime_);
00279   }
00280 
00282   virtual double ApplyInverseTime() const
00283   {
00284     return(ApplyInverseTime_);
00285   }
00286 
00288   virtual double InitializeFlops() const
00289   {
00290     return(InitializeFlops_);
00291   }
00292 
00293   virtual double ComputeFlops() const
00294   {
00295     return(ComputeFlops_);
00296   }
00297 
00298   virtual double ApplyInverseFlops() const
00299   {
00300     return(ApplyInverseFlops_);
00301   }
00302 
00304   virtual int OverlapLevel() const
00305   {
00306     return(OverlapLevel_);
00307   }
00308 
00310   virtual const Teuchos::ParameterList& List() const
00311   {
00312     return(List_);
00313   }
00314 
00315 protected:
00316 
00317   // @}
00318 
00319   // @{ Internal merhods.
00320   
00322   Ifpack_AdditiveSchwarz(const Ifpack_AdditiveSchwarz& RHS)
00323   { }
00324 
00326   int Setup();
00327   
00328   // @}
00329 
00330   // @{ Internal data.
00331   
00333   Teuchos::RefCountPtr<const Epetra_RowMatrix> Matrix_;
00335   Teuchos::RefCountPtr<Ifpack_OverlappingRowMatrix> OverlappingMatrix_;
00337   Teuchos::RefCountPtr<Ifpack_LocalFilter> LocalizedMatrix_;
00339   string Label_;
00341   bool IsInitialized_;
00343   bool IsComputed_;
00345   bool UseTranspose_;
00347   bool IsOverlapping_;
00349   int OverlapLevel_;
00351   Teuchos::ParameterList List_;
00353   Epetra_CombineMode CombineMode_;
00355   double Condest_;
00357   bool ComputeCondest_;
00359   bool UseReordering_;
00361   string ReorderingType_;
00363   Teuchos::RefCountPtr<Ifpack_Reordering> Reordering_;
00365   Teuchos::RefCountPtr<Ifpack_ReorderFilter> ReorderedLocalizedMatrix_;
00367   bool FilterSingletons_;
00369   Teuchos::RefCountPtr<Ifpack_SingletonFilter> SingletonFilter_;
00371   int NumInitialize_;
00373   int NumCompute_;
00375   mutable int NumApplyInverse_;
00377   double InitializeTime_;
00379   double ComputeTime_;
00381   mutable double ApplyInverseTime_;
00383   double InitializeFlops_;
00385   double ComputeFlops_;
00387   mutable double ApplyInverseFlops_;
00389   Teuchos::RefCountPtr<Epetra_Time> Time_;
00391   Teuchos::RefCountPtr<T> Inverse_;
00392 
00393 }; // class Ifpack_AdditiveSchwarz<T>
00394 
00395 //==============================================================================
00396 template<typename T>
00397 Ifpack_AdditiveSchwarz<T>::
00398 Ifpack_AdditiveSchwarz(Epetra_RowMatrix* Matrix_in,
00399                int OverlapLevel_in) :
00400   IsInitialized_(false),
00401   IsComputed_(false),
00402   UseTranspose_(false),
00403   IsOverlapping_(false),
00404   OverlapLevel_(OverlapLevel_in),
00405   CombineMode_(Zero),
00406   Condest_(-1.0),
00407   ComputeCondest_(true),
00408   UseReordering_(false),
00409   ReorderingType_("none"),
00410   FilterSingletons_(false),
00411   NumInitialize_(0),
00412   NumCompute_(0),
00413   NumApplyInverse_(0),
00414   InitializeTime_(0.0),
00415   ComputeTime_(0.0),
00416   ApplyInverseTime_(0.0),
00417   InitializeFlops_(0.0),
00418   ComputeFlops_(0.0),
00419   ApplyInverseFlops_(0.0)
00420 {
00421   // Construct a reference-counted pointer with the input matrix, don't manage the memory.
00422   Matrix_ = Teuchos::rcp( Matrix_in, false );
00423 
00424   if (Matrix_->Comm().NumProc() == 1)
00425     OverlapLevel_ = 0;
00426 
00427   if ((OverlapLevel_ != 0) && (Matrix_->Comm().NumProc() > 1))
00428     IsOverlapping_ = true;
00429   // Sets parameters to default values
00430   Teuchos::ParameterList List_in;
00431   SetParameters(List_in);
00432 }
00433 
00434 //==============================================================================
00435 template<typename T>
00436 int Ifpack_AdditiveSchwarz<T>::Setup()
00437 {
00438 
00439   Epetra_RowMatrix* MatrixPtr;
00440 
00441   if (OverlappingMatrix_ != Teuchos::null)
00442     LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(OverlappingMatrix_) );
00443   else
00444     LocalizedMatrix_ = Teuchos::rcp( new Ifpack_LocalFilter(Matrix_) );
00445 
00446   if (LocalizedMatrix_ == Teuchos::null)
00447     IFPACK_CHK_ERR(-5);
00448 
00449   // users may want to skip singleton check
00450   if (FilterSingletons_) {
00451     SingletonFilter_ = Teuchos::rcp( new Ifpack_SingletonFilter(LocalizedMatrix_) );
00452     MatrixPtr = &*SingletonFilter_;
00453   }
00454   else
00455     MatrixPtr = &*LocalizedMatrix_;
00456 
00457   if (UseReordering_) {
00458 
00459     // create reordering and compute it
00460     if (ReorderingType_ == "rcm")
00461       Reordering_ = Teuchos::rcp( new Ifpack_RCMReordering() );
00462     else if (ReorderingType_ == "metis")
00463       Reordering_ = Teuchos::rcp( new Ifpack_METISReordering() );
00464     else {
00465       cerr << "reordering type not correct (" << ReorderingType_ << ")" << endl;
00466       exit(EXIT_FAILURE);
00467     }
00468     if (Reordering_ == Teuchos::null) IFPACK_CHK_ERR(-5);
00469 
00470     IFPACK_CHK_ERR(Reordering_->SetParameters(List_));
00471     IFPACK_CHK_ERR(Reordering_->Compute(*MatrixPtr));
00472 
00473     // now create reordered localized matrix
00474     ReorderedLocalizedMatrix_ = 
00475       Teuchos::rcp( new Ifpack_ReorderFilter(Teuchos::rcp( MatrixPtr, false ), Reordering_) );
00476 
00477     if (ReorderedLocalizedMatrix_ == Teuchos::null) IFPACK_CHK_ERR(-5);
00478 
00479     MatrixPtr = &*ReorderedLocalizedMatrix_;
00480   }
00481 
00482   Inverse_ = Teuchos::rcp( new T(MatrixPtr) );
00483 
00484   if (Inverse_ == Teuchos::null)
00485     IFPACK_CHK_ERR(-5);
00486 
00487   return(0);
00488 }
00489 
00490 //==============================================================================
00491 template<typename T>
00492 int Ifpack_AdditiveSchwarz<T>::SetParameters(Teuchos::ParameterList& List_in)
00493 {
00494  
00495   // compute the condition number each time Compute() is invoked.
00496   ComputeCondest_ = List_in.get("schwarz: compute condest", ComputeCondest_);
00497   // combine mode
00498   if( Teuchos::ParameterEntry *combineModeEntry = List_in.getEntryPtr("schwarz: combine mode") )
00499   {
00500     if( typeid(std::string) == combineModeEntry->getAny().type() )
00501     {
00502       std::string mode = List_in.get("schwarz: combine mode", "Add");
00503       if (mode == "Add")
00504         CombineMode_ = Add;
00505       else if (mode == "Zero")
00506         CombineMode_ = Zero;
00507       else if (mode == "Insert")
00508         CombineMode_ = Insert;
00509       else if (mode == "InsertAdd")
00510         CombineMode_ = InsertAdd;
00511       else if (mode == "Average")
00512         CombineMode_ = Average;
00513       else if (mode == "AbsMax")
00514         CombineMode_ = AbsMax;
00515       else
00516       {
00517         TEST_FOR_EXCEPTION(
00518           true,std::logic_error
00519           ,"Error, The (Epetra) combine mode of \""<<mode<<"\" is not valid!  Only the values"
00520           " \"Add\", \"Zero\", \"Insert\", \"InsertAdd\", \"Average\", and \"AbsMax\" are accepted!"
00521           );
00522       }
00523     }
00524     else if ( typeid(Epetra_CombineMode) == combineModeEntry->getAny().type() )
00525     {
00526       CombineMode_ = Teuchos::any_cast<Epetra_CombineMode>(combineModeEntry->getAny());
00527     }
00528     else
00529     {
00530       // Throw exception with good error message!
00531       Teuchos::getParameter<std::string>(List_in,"schwarz: combine mode");
00532     }
00533   }
00534   else
00535   {
00536     // Make the default be a string to be consistent with the valid parameters!
00537     List_in.get("schwarz: combine mode","Zero");
00538   }
00539   // type of reordering
00540   ReorderingType_ = List_in.get("schwarz: reordering type", ReorderingType_);
00541   if (ReorderingType_ == "none")
00542     UseReordering_ = false;
00543   else 
00544     UseReordering_ = true;
00545   // if true, filter singletons. NOTE: the filtered matrix can still have
00546   // singletons! A simple example: upper triangular matrix, if I remove
00547   // the lower node, I still get a matrix with a singleton! However, filter
00548   // singletons should help for PDE problems with Dirichlet BCs.
00549   FilterSingletons_ = List_in.get("schwarz: filter singletons", FilterSingletons_);
00550 
00551   // This copy may be needed by Amesos or other preconditioners.
00552   List_ = List_in;
00553 
00554   return(0);
00555 }
00556 
00557 //==============================================================================
00558 template<typename T>
00559 int Ifpack_AdditiveSchwarz<T>::Initialize()
00560 {
00561   IsInitialized_ = false;
00562   IsComputed_ = false; // values required
00563   Condest_ = -1.0; // zero-out condest
00564 
00565   if (Time_ == Teuchos::null)
00566     Time_ = Teuchos::rcp( new Epetra_Time(Comm()) );
00567 
00568   Time_->ResetStartTime();
00569 
00570   // compute the overlapping matrix if necessary
00571   if (IsOverlapping_) {
00572     OverlappingMatrix_ = 
00573       Teuchos::rcp( new Ifpack_OverlappingRowMatrix(Matrix_, OverlapLevel_) );
00574     if (OverlappingMatrix_ == Teuchos::null)
00575       IFPACK_CHK_ERR(-5);
00576   }
00577 
00578   IFPACK_CHK_ERR(Setup());
00579 
00580   if (Inverse_ == Teuchos::null)
00581     IFPACK_CHK_ERR(-5);
00582 
00583   if (LocalizedMatrix_ == Teuchos::null)
00584     IFPACK_CHK_ERR(-5);
00585 
00586   IFPACK_CHK_ERR(Inverse_->SetUseTranspose(UseTranspose()));
00587   IFPACK_CHK_ERR(Inverse_->SetParameters(List_));
00588   IFPACK_CHK_ERR(Inverse_->Initialize());
00589 
00590   // Label is for Aztec-like solvers
00591   Label_ = "Ifpack_AdditiveSchwarz, ";
00592   if (UseTranspose())
00593     Label_ += ", transp";
00594   Label_ += ", ov = " + Ifpack_toString(OverlapLevel_)
00595     + ", local solver = \n\t\t***** `" + string(Inverse_->Label()) + "'";
00596 
00597   IsInitialized_ = true;
00598   ++NumInitialize_;
00599   InitializeTime_ += Time_->ElapsedTime();
00600 
00601   // count flops by summing up all the InitializeFlops() in each
00602   // Inverse. Each Inverse() can only give its flops -- it acts on one
00603   // process only
00604   double partial = Inverse_->InitializeFlops();
00605   double total;
00606   Comm().SumAll(&partial, &total, 1);
00607   InitializeFlops_ += total;
00608 
00609   return(0);
00610 }
00611 
00612 //==============================================================================
00613 template<typename T>
00614 int Ifpack_AdditiveSchwarz<T>::Compute()
00615 {
00616 
00617   if (IsInitialized() == false)
00618     IFPACK_CHK_ERR(Initialize());
00619 
00620   Time_->ResetStartTime();
00621   IsComputed_ = false;
00622   Condest_ = -1.0;
00623   
00624   IFPACK_CHK_ERR(Inverse_->Compute());
00625 
00626   IsComputed_ = true; // need this here for Condest(Ifpack_Cheap)
00627   ++NumCompute_;
00628   ComputeTime_ += Time_->ElapsedTime();
00629 
00630   // sum up flops
00631   double partial = Inverse_->ComputeFlops();
00632    double total;
00633   Comm().SumAll(&partial, &total, 1);
00634   ComputeFlops_ += total;
00635 
00636   // reset the Label
00637   string R = "";
00638   if (UseReordering_)
00639     R = ReorderingType_ + " reord, ";
00640 
00641   if (ComputeCondest_)
00642     Condest(Ifpack_Cheap);
00643   
00644   // add Condest() to label
00645   Label_ = "Ifpack_AdditiveSchwarz, ov = " + Ifpack_toString(OverlapLevel_)
00646     + ", local solver = \n\t\t***** `" + string(Inverse_->Label()) + "'"
00647     + "\n\t\t***** " + R + "Condition number estimate = "
00648     + Ifpack_toString(Condest());
00649 
00650   return(0);
00651 }
00652 
00653 //==============================================================================
00654 template<typename T>
00655 int Ifpack_AdditiveSchwarz<T>::SetUseTranspose(bool UseTranspose_in)
00656 {
00657   // store the flag -- it will be set in Initialize() if Inverse_ does not
00658   // exist.
00659   UseTranspose_ = UseTranspose_in;
00660 
00661   // If Inverse_ exists, pass it right now.
00662   if (Inverse_!=Teuchos::null)
00663     IFPACK_CHK_ERR(Inverse_->SetUseTranspose(UseTranspose_in));
00664   return(0);
00665 }
00666 
00667 //==============================================================================
00668 template<typename T>
00669 int Ifpack_AdditiveSchwarz<T>::
00670 Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const
00671 {
00672   IFPACK_CHK_ERR(Matrix_->Apply(X,Y));
00673   return(0);
00674 }
00675 
00676 //==============================================================================
00677 template<typename T>
00678 double Ifpack_AdditiveSchwarz<T>::NormInf() const
00679 {
00680   return(-1.0);
00681 }
00682 
00683 //==============================================================================
00684 template<typename T>
00685 const char * Ifpack_AdditiveSchwarz<T>::Label() const
00686 {
00687   return(Label_.c_str());
00688 }
00689 
00690 //==============================================================================
00691 template<typename T>
00692 bool Ifpack_AdditiveSchwarz<T>::UseTranspose() const
00693 {
00694   return(UseTranspose_);
00695 }
00696 
00697 //==============================================================================
00698 template<typename T>
00699 bool Ifpack_AdditiveSchwarz<T>::HasNormInf() const
00700 {
00701   return(false);
00702 }
00703 
00704 //==============================================================================
00705 template<typename T>
00706 const Epetra_Comm & Ifpack_AdditiveSchwarz<T>::Comm() const
00707 {
00708   return(Matrix_->Comm());
00709 }
00710 
00711 //==============================================================================
00712 template<typename T>
00713 const Epetra_Map & Ifpack_AdditiveSchwarz<T>::OperatorDomainMap() const
00714 {
00715   return(Matrix_->OperatorDomainMap());
00716 }
00717 
00718 //==============================================================================
00719 template<typename T>
00720 const Epetra_Map & Ifpack_AdditiveSchwarz<T>::OperatorRangeMap() const
00721 {
00722   return(Matrix_->OperatorRangeMap());
00723 }
00724 
00725 //==============================================================================
00726 template<typename T>
00727 int Ifpack_AdditiveSchwarz<T>::
00728 ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const
00729 {
00730   // compute the preconditioner is not done by the user
00731   if (!IsComputed())
00732     IFPACK_CHK_ERR(-3);
00733   
00734   int NumVectors = X.NumVectors();
00735 
00736   if (NumVectors != Y.NumVectors())
00737     IFPACK_CHK_ERR(-2); // wrong input
00738 
00739   Time_->ResetStartTime();
00740 
00741   Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingX;
00742   Teuchos::RefCountPtr<Epetra_MultiVector> OverlappingY;
00743   Teuchos::RefCountPtr<Epetra_MultiVector> Xtmp;
00744 
00745   // for flop count, see bottom of this function
00746   double pre_partial = Inverse_->ApplyInverseFlops();
00747   double pre_total;
00748   Comm().SumAll(&pre_partial, &pre_total, 1);
00749 
00750   // process overlap, may need to create vectors and import data
00751   if (IsOverlapping()) {
00752     OverlappingX = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
00753                             X.NumVectors()) );
00754     OverlappingY = Teuchos::rcp( new Epetra_MultiVector(OverlappingMatrix_->RowMatrixRowMap(),
00755                             Y.NumVectors()) );
00756     if (OverlappingY == Teuchos::null) IFPACK_CHK_ERR(-5);
00757 
00758     OverlappingY->PutScalar(0.0);
00759     OverlappingX->PutScalar(0.0);
00760     IFPACK_CHK_ERR(OverlappingMatrix_->ImportMultiVector(X,*OverlappingX,Insert));
00761     // FIXME: this will not work with overlapping and non-zero starting
00762     // solutions. The same for other cases below.
00763     // IFPACK_CHK_ERR(OverlappingMatrix_->ImportMultiVector(Y,*OverlappingY,Insert));
00764   }
00765   else {
00766     Xtmp = Teuchos::rcp( new Epetra_MultiVector(X) );
00767     OverlappingX = Xtmp;
00768     OverlappingY = Teuchos::rcp( &Y, false );
00769   }
00770 
00771   if (FilterSingletons_) {
00772     // process singleton filter
00773     Epetra_MultiVector ReducedX(SingletonFilter_->Map(),NumVectors);
00774     Epetra_MultiVector ReducedY(SingletonFilter_->Map(),NumVectors);
00775     IFPACK_CHK_ERR(SingletonFilter_->SolveSingletons(*OverlappingX,*OverlappingY));
00776     IFPACK_CHK_ERR(SingletonFilter_->CreateReducedRHS(*OverlappingY,*OverlappingX,ReducedX));
00777 
00778     // process reordering
00779     if (!UseReordering_) {
00780       IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReducedX,ReducedY));
00781     }
00782     else {
00783       Epetra_MultiVector ReorderedX(ReducedX);
00784       Epetra_MultiVector ReorderedY(ReducedY);
00785       IFPACK_CHK_ERR(Reordering_->P(ReducedX,ReorderedX));
00786       IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReorderedX,ReorderedY));
00787       IFPACK_CHK_ERR(Reordering_->Pinv(ReorderedY,ReducedY));
00788     }
00789 
00790     // finish up with singletons
00791     IFPACK_CHK_ERR(SingletonFilter_->UpdateLHS(ReducedY,*OverlappingY));
00792   }
00793   else {
00794     // process reordering
00795     if (!UseReordering_) {
00796       IFPACK_CHK_ERR(Inverse_->ApplyInverse(*OverlappingX,*OverlappingY));
00797     }
00798     else {
00799       Epetra_MultiVector ReorderedX(*OverlappingX);
00800       Epetra_MultiVector ReorderedY(*OverlappingY);
00801       IFPACK_CHK_ERR(Reordering_->P(*OverlappingX,ReorderedX));
00802       IFPACK_CHK_ERR(Inverse_->ApplyInverse(ReorderedX,ReorderedY));
00803       IFPACK_CHK_ERR(Reordering_->Pinv(ReorderedY,*OverlappingY));
00804     }
00805   }
00806 
00807   if (IsOverlapping()) {
00808     IFPACK_CHK_ERR(OverlappingMatrix_->ExportMultiVector(*OverlappingY,Y,
00809                              CombineMode_));
00810   }
00811 
00812   // add flops. Note the we only have to add the newly counted
00813   // flops -- and each Inverse returns the cumulative sum
00814   double partial = Inverse_->ApplyInverseFlops();
00815   double total;
00816   Comm().SumAll(&partial, &total, 1);
00817   ApplyInverseFlops_ += total - pre_total;
00818 
00819   // FIXME: right now I am skipping the overlap and singletons
00820   ++NumApplyInverse_;
00821   ApplyInverseTime_ += Time_->ElapsedTime();
00822 
00823   return(0);
00824  
00825 }
00826 
00827 //==============================================================================
00828 template<typename T>
00829 std::ostream& Ifpack_AdditiveSchwarz<T>::
00830 Print(std::ostream& os) const
00831 {
00832   double IF = InitializeFlops();
00833   double CF = ComputeFlops();
00834   double AF = ApplyInverseFlops();
00835 
00836   double IFT = 0.0, CFT = 0.0, AFT = 0.0;
00837   if (InitializeTime() != 0.0)
00838     IFT = IF / InitializeTime();
00839   if (ComputeTime() != 0.0)
00840     CFT = CF / ComputeTime();
00841   if (ApplyInverseTime() != 0.0)
00842     AFT = AF / ApplyInverseTime();
00843 
00844   if (Matrix().Comm().MyPID())
00845     return(os);
00846 
00847   os << endl;
00848   os << "================================================================================" << endl;
00849   os << "Ifpack_AdditiveSchwarz, overlap level = " << OverlapLevel_ << endl;
00850   if (CombineMode_ == Insert)
00851     os << "Combine mode                          = Insert" << endl;
00852   else if (CombineMode_ == Add)
00853     os << "Combine mode                          = Add" << endl;
00854   else if (CombineMode_ == Zero)
00855     os << "Combine mode                          = Zero" << endl;
00856   else if (CombineMode_ == Average)
00857     os << "Combine mode                          = Average" << endl;
00858   else if (CombineMode_ == AbsMax)
00859     os << "Combine mode                          = AbsMax" << endl;
00860 
00861   os << "Condition number estimate             = " << Condest_ << endl;
00862   os << "Global number of rows                 = " << Matrix_->NumGlobalRows() << endl;
00863 
00864   os << endl;
00865   os << "Phase           # calls   Total Time (s)       Total MFlops     MFlops/s" << endl;
00866   os << "-----           -------   --------------       ------------     --------" << endl;
00867   os << "Initialize()    "   << std::setw(5) << NumInitialize()
00868      << "  " << std::setw(15) << InitializeTime() 
00869      << "  " << std::setw(15) << 1.0e-6 * IF 
00870      << "  " << std::setw(15) << 1.0e-6 * IFT << endl;
00871   os << "Compute()       "   << std::setw(5) << NumCompute() 
00872      << "  " << std::setw(15) << ComputeTime()
00873      << "  " << std::setw(15) << 1.0e-6 * CF
00874      << "  " << std::setw(15) << 1.0e-6 * CFT << endl;
00875   os << "ApplyInverse()  "   << std::setw(5) << NumApplyInverse() 
00876      << "  " << std::setw(15) << ApplyInverseTime()
00877      << "  " << std::setw(15) << 1.0e-6 * AF
00878      << "  " << std::setw(15) << 1.0e-6 * AFT << endl;
00879   os << "================================================================================" << endl;
00880   os << endl;
00881 
00882   return(os);
00883 }
00884 
00885 #include "Ifpack_Condest.h"
00886 //==============================================================================
00887 template<typename T>
00888 double Ifpack_AdditiveSchwarz<T>::
00889 Condest(const Ifpack_CondestType CT, const int MaxIters, 
00890         const double Tol, Epetra_RowMatrix* Matrix_in)
00891 {
00892   if (!IsComputed()) // cannot compute right now
00893     return(-1.0);
00894 
00895   Condest_ = Ifpack_Condest(*this, CT, MaxIters, Tol, Matrix_in);
00896 
00897   return(Condest_);
00898 }
00899 
00900 #endif // IFPACK_ADDITIVESCHWARZ_H

Generated on Wed May 12 21:46:03 2010 for IFPACK by  doxygen 1.4.7