Ifpack2 Templated Preconditioning Package Version 1.0
Ifpack2_RILUK_decl.hpp
Go to the documentation of this file.
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 
00045 
00046 #ifndef IFPACK2_CRSRILUK_DECL_HPP
00047 #define IFPACK2_CRSRILUK_DECL_HPP
00048 
00049 #include <Ifpack2_ConfigDefs.hpp>
00050 #include <Ifpack2_Preconditioner.hpp>
00051 #include <Ifpack2_Details_CanChangeMatrix.hpp>
00052 #include <Tpetra_CrsMatrix.hpp>
00053 
00054 #include "Ifpack2_ScalingType.hpp"
00055 #include "Ifpack2_IlukGraph.hpp"
00056 
00057 namespace Teuchos {
00058   class ParameterList; // forward declaration
00059 }
00060 
00061 namespace Ifpack2 {
00062 
00221 template<class MatrixType>
00222 class RILUK:
00223     virtual public Ifpack2::Preconditioner<typename MatrixType::scalar_type,
00224                                            typename MatrixType::local_ordinal_type,
00225                                            typename MatrixType::global_ordinal_type,
00226                                            typename MatrixType::node_type>,
00227     virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type,
00228                                                                        typename MatrixType::local_ordinal_type,
00229                                                                        typename MatrixType::global_ordinal_type,
00230                                                                        typename MatrixType::node_type> >
00231 {
00232  public:
00234   typedef typename MatrixType::scalar_type scalar_type;
00235 
00237   TEUCHOS_DEPRECATED typedef typename MatrixType::scalar_type Scalar;
00238 
00239 
00241   typedef typename MatrixType::local_ordinal_type local_ordinal_type;
00242 
00244   TEUCHOS_DEPRECATED typedef typename MatrixType::local_ordinal_type LocalOrdinal;
00245 
00246 
00248   typedef typename MatrixType::global_ordinal_type global_ordinal_type;
00249 
00251   TEUCHOS_DEPRECATED typedef typename MatrixType::global_ordinal_type GlobalOrdinal;
00252 
00253 
00255   typedef typename MatrixType::node_type node_type;
00256 
00258   TEUCHOS_DEPRECATED typedef typename MatrixType::node_type Node;
00259 
00261   typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitude_type;
00262 
00264   TEUCHOS_DEPRECATED typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitudeType;
00265 
00267   typedef Tpetra::RowMatrix<scalar_type,
00268                             local_ordinal_type,
00269                             global_ordinal_type,
00270                             node_type> row_matrix_type;
00271 
00273   typedef Tpetra::CrsMatrix<scalar_type,
00274                             local_ordinal_type,
00275                             global_ordinal_type,
00276                             node_type> crs_matrix_type;
00277 
00278   template <class NewMatrixType> friend class RILUK;
00279 
00283   RILUK (const Teuchos::RCP<const row_matrix_type>& A_in);
00284 
00292   RILUK (const Teuchos::RCP<const crs_matrix_type>& A_in);
00293 
00294  private:
00297   RILUK (const RILUK<MatrixType> & src);
00298 
00299  public:
00305   template <typename NewMatrixType>
00306   Teuchos::RCP< RILUK< NewMatrixType > >
00307   clone (const Teuchos::RCP<const NewMatrixType>& A_newnode) const;
00308 
00310   virtual ~RILUK ();
00311 
00316   void TEUCHOS_DEPRECATED SetRelaxValue (const magnitude_type RelaxValue) {
00317     RelaxValue_ = RelaxValue;
00318   }
00323   void TEUCHOS_DEPRECATED SetAbsoluteThreshold (const magnitude_type Athresh) {
00324     Athresh_ = Athresh;
00325   }
00330   void TEUCHOS_DEPRECATED SetRelativeThreshold (const magnitude_type Rthresh) {
00331     Rthresh_ = Rthresh;
00332   }
00337   void TEUCHOS_DEPRECATED SetOverlapMode (const Tpetra::CombineMode OverlapMode) {
00338     TEUCHOS_TEST_FOR_EXCEPTION(
00339       true, std::logic_error, "Ifpack2::RILUK::SetOverlapMode: "
00340       "RILUK no longer implements overlap on its own.  "
00341       "Use RILUK with AdditiveSchwarz if you want overlap.");
00342   }
00343 
00351   void setParameters (const Teuchos::ParameterList& params);
00352 
00354   void initialize ();
00355 
00364   void compute ();
00365 
00367   bool isInitialized () const {
00368     return isInitialized_;
00369   }
00371   bool isComputed () const {
00372     return isComputed_;
00373   }
00374 
00376   int getNumInitialize () const {
00377     return numInitialize_;
00378   }
00380   int getNumCompute () const {
00381     return numCompute_;
00382   }
00384   int getNumApply () const {
00385     return numApply_;
00386   }
00387 
00389   double getInitializeTime () const {
00390     return initializeTime_;
00391   }
00393   double getComputeTime () const {
00394     return computeTime_;
00395   }
00397   double getApplyTime () const {
00398     return applyTime_;
00399   }
00400 
00402 
00403 
00426   virtual void
00427   setMatrix (const Teuchos::RCP<const row_matrix_type>& A);
00428 
00430 
00431 
00432 
00434   std::string description () const;
00435 
00437 
00438 
00439 
00441   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00442   getDomainMap () const;
00443 
00445   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00446   getRangeMap () const;
00447 
00477   void
00478   apply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00479          Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00480          Teuchos::ETransp mode = Teuchos::NO_TRANS,
00481          scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one (),
00482          scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero ()) const;
00484 
00485 private:
00507   void
00508   multiply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00509             Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00510             const Teuchos::ETransp mode = Teuchos::NO_TRANS) const;
00511 public:
00512 
00522   magnitude_type TEUCHOS_DEPRECATED computeCondEst (Teuchos::ETransp mode) const;
00523 
00533   virtual magnitude_type TEUCHOS_DEPRECATED
00534   computeCondEst (CondestType CT = Ifpack2::Cheap,
00535                   local_ordinal_type MaxIters = 1550,
00536                   magnitude_type Tol = 1e-9,
00537                   const Teuchos::Ptr<const row_matrix_type>& Matrix = Teuchos::null);
00538 
00539   magnitude_type getCondEst () const { return Condest_; }
00540 
00542   Teuchos::RCP<const row_matrix_type> getMatrix () const;
00543 
00544   // Attribute access functions
00545 
00547   magnitude_type getRelaxValue () const { return RelaxValue_; }
00548 
00550   magnitude_type getAbsoluteThreshold () const { return Athresh_; }
00551 
00553   magnitude_type getRelativeThreshold () const {return Rthresh_;}
00554 
00556   int getLevelOfFill () const { return LevelOfFill_; }
00557 
00559   Tpetra::CombineMode getOverlapMode () {
00560     TEUCHOS_TEST_FOR_EXCEPTION(
00561       true, std::logic_error, "Ifpack2::RILUK::SetOverlapMode: "
00562       "RILUK no longer implements overlap on its own.  "
00563       "Use RILUK with AdditiveSchwarz if you want overlap.");
00564   }
00565 
00567   Tpetra::global_size_t getGlobalNumEntries () const {
00568     return getL ().getGlobalNumEntries () + getU ().getGlobalNumEntries ();
00569   }
00570 
00572   Teuchos::RCP<Ifpack2::IlukGraph<Tpetra::CrsGraph<local_ordinal_type,global_ordinal_type,node_type> > > getGraph () const {
00573     return Graph_;
00574   }
00575 
00577   const crs_matrix_type& getL () const;
00578 
00580   const Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>&
00581   getD () const;
00582 
00584   const crs_matrix_type& getU () const;
00585 
00587   Teuchos::RCP<const crs_matrix_type> getCrsMatrix () const;
00588 
00589 private:
00590   typedef Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> MV;
00591   typedef Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> vec_type;
00592   typedef Teuchos::ScalarTraits<scalar_type> STS;
00593   typedef Teuchos::ScalarTraits<magnitude_type> STM;
00594 
00595   void allocate_L_and_U();
00596   void initAllValues (const row_matrix_type& A);
00597 
00603   static Teuchos::RCP<const row_matrix_type>
00604   makeLocalFilter (const Teuchos::RCP<const row_matrix_type>& A);
00605 
00607   Teuchos::RCP<const row_matrix_type> A_;
00608 
00610   Teuchos::RCP<Ifpack2::IlukGraph<Tpetra::CrsGraph<local_ordinal_type,
00611                                                    global_ordinal_type,
00612                                                    node_type> > > Graph_;
00621   Teuchos::RCP<const crs_matrix_type> A_local_crs_;
00622 
00624   Teuchos::RCP<crs_matrix_type> L_;
00626   Teuchos::RCP<crs_matrix_type> U_;
00628   Teuchos::RCP<vec_type> D_;
00629 
00630   int LevelOfFill_;
00631 
00632   bool isAllocated_;
00633   bool isInitialized_;
00634   bool isComputed_;
00635 
00636   int numInitialize_;
00637   int numCompute_;
00638   mutable int numApply_;
00639 
00640   double initializeTime_;
00641   double computeTime_;
00642   mutable double applyTime_;
00643 
00644   magnitude_type RelaxValue_;
00645   magnitude_type Athresh_;
00646   magnitude_type Rthresh_;
00647 
00648   mutable magnitude_type Condest_;
00649 };
00650 
00651 //Set necessary local solve parameters when using ThrustGPU node
00652 namespace detail {
00653   template<class MatrixType, class NodeType, class MatSolveType>
00654   struct setLocalSolveParams{
00655     static Teuchos::RCP<Teuchos::ParameterList>
00656     setParams (const Teuchos::RCP<Teuchos::ParameterList>& param) {
00657       return param;
00658     }
00659   };
00660 #if defined(HAVE_KOKKOSCLASSIC_THRUST) && defined(HAVE_KOKKOSCLASSIC_CUSPARSE)
00661   template<class MatrixType, class Scalar>
00662   struct setLocalSolveParams<MatrixType,
00663                              KokkosClassic::ThrustGPUNode,
00664                              KokkosClassic::CUSPARSEOps<Scalar, KokkosClassic::ThrustGPUNode> >
00665   {
00666     static Teuchos::RCP<Teuchos::ParameterList>
00667     setParams (const Teuchos::RCP<Teuchos::ParameterList>& param) {
00668       param->sublist ("fillComplete").sublist ("Local Sparse Ops").set ("Prepare Solve", true);
00669       return param;
00670     }
00671   };
00672 #endif
00673 } //end namespace detail
00674 
00675 template <class MatrixType>
00676 template <typename NewMatrixType>
00677 Teuchos::RCP<RILUK<NewMatrixType> >
00678 RILUK<MatrixType>::
00679 clone (const Teuchos::RCP<const NewMatrixType>& A_newnode) const
00680 {
00681   using Teuchos::ParameterList;
00682   using Teuchos::RCP;
00683   using Teuchos::rcp;
00684   typedef typename NewMatrixType::node_type new_node_type;
00685   typedef typename NewMatrixType::mat_solve_type mat_solve_type;
00686   typedef RILUK<NewMatrixType> new_riluk_type;
00687 
00688   RCP<new_riluk_type> new_riluk = rcp (new new_riluk_type (A_newnode));
00689 
00690   RCP<ParameterList> plClone = Teuchos::parameterList ();
00691   plClone = detail::setLocalSolveParams<NewMatrixType,
00692     new_node_type, mat_solve_type>::setParams (plClone);
00693 
00694   RCP<new_node_type> new_node = A_newnode->getNode ();
00695   new_riluk->L_ = L_->clone (new_node, plClone);
00696   new_riluk->U_ = U_->clone (new_node, plClone);
00697   new_riluk->D_ = D_->clone (new_node);
00698 
00699   new_riluk->LevelOfFill_ = LevelOfFill_;
00700 
00701   new_riluk->isAllocated_ = isAllocated_;
00702   new_riluk->isInitialized_ = isInitialized_;
00703   new_riluk->isComputed_ = isComputed_;
00704 
00705   new_riluk->numInitialize_ = numInitialize_;
00706   new_riluk->numCompute_ = numCompute_;
00707   new_riluk->numApply_ =  numApply_;
00708 
00709   new_riluk->RelaxValue_ = RelaxValue_;
00710   new_riluk->Athresh_ = Athresh_;
00711   new_riluk->Rthresh_ = Rthresh_;
00712   new_riluk->Condest_ = Condest_;
00713 
00714   return new_riluk;
00715 }
00716 
00717 } // namespace Ifpack2
00718 
00719 #endif /* IFPACK2_CRSRILUK_DECL_HPP */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends