Ifpack2 Templated Preconditioning Package: Ifpack2_Relaxation

00001 /*@HEADER
00002 // ***********************************************************************
00003 //
00004 //       Ifpack2: Tempated Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2009) 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 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ***********************************************************************
00040 //@HEADER
00041 */
00042 
00043 #ifndef IFPACK2_RELAXATION_DECL_HPP
00044 #define IFPACK2_RELAXATION_DECL_HPP
00045 
00046 #include <Ifpack2_ConfigDefs.hpp>
00047 #include <Ifpack2_Preconditioner.hpp>
00048 #include <Ifpack2_Details_CanChangeMatrix.hpp>
00049 #include <Ifpack2_Condest.hpp>
00050 #include <Ifpack2_Parameters.hpp>
00051 #include <Tpetra_Vector.hpp>
00052 #include <Teuchos_ScalarTraits.hpp>
00053 #include <Tpetra_CrsMatrix_decl.hpp> // Don't need the definition here
00054 #include <Tpetra_Experimental_BlockCrsMatrix_decl.hpp>
00055 
00056 namespace Teuchos {
00057   // forward declarations
00058   class ParameterList;
00059   class Time;
00060 } // namespace Teuchos
00061 
00062 namespace Ifpack2 {
00063 
00227 template<class MatrixType>
00228 class Relaxation :
00229   virtual public Ifpack2::Preconditioner<typename MatrixType::scalar_type,
00230                                          typename MatrixType::local_ordinal_type,
00231                                          typename MatrixType::global_ordinal_type,
00232                                          typename MatrixType::node_type>,
00233   virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type,
00234                                                                      typename MatrixType::local_ordinal_type,
00235                                                                      typename MatrixType::global_ordinal_type,
00236                                                                      typename MatrixType::node_type> >
00237 {
00238 public:
00240 
00241 
00243   typedef typename MatrixType::scalar_type scalar_type;
00244 
00246   TEUCHOS_DEPRECATED typedef typename MatrixType::scalar_type Scalar;
00247 
00248 
00250   typedef typename MatrixType::local_ordinal_type local_ordinal_type;
00251 
00253   TEUCHOS_DEPRECATED typedef typename MatrixType::local_ordinal_type LocalOrdinal;
00254 
00255 
00257   typedef typename MatrixType::global_ordinal_type global_ordinal_type;
00258 
00260   TEUCHOS_DEPRECATED typedef typename MatrixType::global_ordinal_type GlobalOrdinal;
00261 
00262 
00264   typedef typename MatrixType::node_type node_type;
00265 
00267   TEUCHOS_DEPRECATED typedef typename MatrixType::node_type Node;
00268 
00269 
00271   typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitude_type;
00272 
00274   TEUCHOS_DEPRECATED typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitudeType;
00275 
00277   typedef Tpetra::RowMatrix<scalar_type, local_ordinal_type,
00278                             global_ordinal_type, node_type> row_matrix_type;
00280 
00281 
00282 
00316   explicit Relaxation (const Teuchos::RCP<const row_matrix_type>& A);
00317 
00319   virtual ~Relaxation();
00320 
00322 
00323 
00324 
00402   void setParameters (const Teuchos::ParameterList& params);
00403 
00405   Teuchos::RCP<const Teuchos::ParameterList>
00406   getValidParameters () const;
00407 
00416   void initialize ();
00417 
00419   inline bool isInitialized() const {
00420     return isInitialized_;
00421   }
00422 
00430   void compute ();
00431 
00432 
00434   inline bool isComputed() const {
00435     return(IsComputed_);
00436   }
00437 
00439 
00440 
00441 
00464   virtual void
00465   setMatrix (const Teuchos::RCP<const row_matrix_type>& A);
00466 
00468 
00469 
00470 
00486   void
00487   apply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00488          Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00489          Teuchos::ETransp mode = Teuchos::NO_TRANS,
00490          scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
00491          scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero()) const;
00492 
00494   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00495   getDomainMap () const;
00496 
00498   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00499   getRangeMap () const;
00500 
00502   bool hasTransposeApply () const;
00503 
00516   void
00517   applyMat (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00518             Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00519             Teuchos::ETransp mode = Teuchos::NO_TRANS) const;
00520 
00522 
00523 
00524 
00534   virtual magnitude_type TEUCHOS_DEPRECATED
00535   computeCondEst (CondestType CT = Cheap,
00536                   local_ordinal_type MaxIters = 1550,
00537                   magnitude_type Tol = 1e-9,
00538                   const Teuchos::Ptr<const row_matrix_type>& matrix = Teuchos::null);
00539 
00541 
00542 
00543 
00547   virtual magnitude_type TEUCHOS_DEPRECATED getCondEst() const;
00548 
00550   Teuchos::RCP<const Teuchos::Comm<int> > getComm() const;
00551 
00553   Teuchos::RCP<const row_matrix_type> getMatrix () const;
00554 
00556   double getComputeFlops() const;
00557 
00559   double getApplyFlops() const;
00560 
00562   int getNumInitialize() const;
00563 
00565   int getNumCompute() const;
00566 
00568   int getNumApply() const;
00569 
00571   double getInitializeTime() const;
00572 
00574   double getComputeTime() const;
00575 
00577   double getApplyTime() const;
00578 
00580 
00581 
00582 
00589   std::string description () const;
00590 
00613   void
00614   describe (Teuchos::FancyOStream &out,
00615             const Teuchos::EVerbosityLevel verbLevel =
00616             Teuchos::Describable::verbLevel_default) const;
00618 
00619 private:
00621 
00622 
00623   typedef Teuchos::ScalarTraits<scalar_type> STS;
00624   typedef Teuchos::ScalarTraits<magnitude_type> STM;
00625 
00630   typedef Tpetra::CrsMatrix<scalar_type, local_ordinal_type,
00631                             global_ordinal_type, node_type> crs_matrix_type;
00632   typedef Tpetra::Experimental::BlockCrsMatrix<scalar_type, local_ordinal_type,
00633                             global_ordinal_type, node_type> block_crs_matrix_type;
00634   typedef Tpetra::Experimental::BlockMultiVector<scalar_type, local_ordinal_type,
00635                             global_ordinal_type, node_type> block_multivector_type;
00637 
00638 
00639 
00641   Relaxation (const Relaxation<MatrixType>& RHS);
00642 
00644   Relaxation<MatrixType>& operator= (const Relaxation<MatrixType>& RHS);
00645 
00647 
00648 
00649 
00654   void setParametersImpl (Teuchos::ParameterList& params);
00655 
00657   void ApplyInverseJacobi(
00658         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00659               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00660 
00662   void ApplyInverseJacobi_BlockCrsMatrix(
00663         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00664               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00665 
00667   void ApplyInverseGS(
00668         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00669               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00670 
00672   void ApplyInverseGS_RowMatrix(
00673         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00674               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00675 
00677   void
00678   ApplyInverseGS_CrsMatrix (const crs_matrix_type& A,
00679                             const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00680                             Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00681 
00683   void
00684   ApplyInverseGS_BlockCrsMatrix (const block_crs_matrix_type& A,
00685                             const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00686                             Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y);
00687 
00689   void ApplyInverseSGS(
00690         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00691               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00692 
00694   void ApplyInverseSGS_RowMatrix(
00695         const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00696               Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00697 
00699   void
00700   ApplyInverseSGS_CrsMatrix (const crs_matrix_type& A,
00701                              const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00702                              Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y) const;
00703 
00705   void
00706   ApplyInverseSGS_BlockCrsMatrix (const block_crs_matrix_type& A,
00707                              const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00708                              Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y);
00709 
00710   void computeBlockCrs ();
00711 
00712 
00714 
00715 
00716 
00723   mutable Teuchos::RCP<const Teuchos::ParameterList> validParams_;
00724 
00726   Teuchos::RCP<const row_matrix_type> A_;
00728   Teuchos::RCP<Teuchos::Time> Time_;
00730   Teuchos::RCP<const Tpetra::Import<local_ordinal_type,global_ordinal_type,node_type> > Importer_;
00732   mutable Teuchos::RCP<Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> > Diagonal_;
00733 
00734   Teuchos::RCP<block_crs_matrix_type> BlockDiagonal_;
00735   Teuchos::Array<int> blockDiagonalFactorizationPivots;
00736   Teuchos::RCP<block_multivector_type> yBlockColumnPointMap_;
00737 
00738 
00740   int NumSweeps_;
00742   Details::RelaxationType PrecType_;
00744   scalar_type DampingFactor_;
00746   bool IsParallel_;
00748   bool ZeroStartingSolution_;
00750   bool DoBackwardGS_;
00752   bool DoL1Method_;
00754   magnitude_type L1Eta_;
00756   scalar_type MinDiagonalValue_;
00758   bool fixTinyDiagEntries_;
00760   bool checkDiagEntries_;
00761 
00763   magnitude_type Condest_;
00765   bool isInitialized_;
00767   bool IsComputed_;
00769   int NumInitialize_;
00771   int NumCompute_;
00773   mutable int NumApply_;
00775   double InitializeTime_;
00777   double ComputeTime_;
00779   mutable double ApplyTime_;
00781   double ComputeFlops_;
00783   mutable double ApplyFlops_;
00784 
00786   magnitude_type globalMinMagDiagEntryMag_;
00788   magnitude_type globalMaxMagDiagEntryMag_;
00790   size_t globalNumSmallDiagEntries_;
00792   size_t globalNumZeroDiagEntries_;
00794   size_t globalNumNegDiagEntries_;
00799   magnitude_type globalDiagNormDiff_;
00800 
00806   Teuchos::ArrayRCP<size_t> diagOffsets_;
00812   bool savedDiagOffsets_;
00813 
00814   bool hasBlockCrsMatrix_;
00815 
00817   Teuchos::ArrayRCP<local_ordinal_type> localSmoothingIndices_;
00818 
00820 }; //class Relaxation
00821 
00822 }//namespace Ifpack2
00823 
00824 #endif // IFPACK2_RELAXATION_DECL_HPP
00825