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 
00239 template<class MatrixType>
00240 class RILUK:
00241     virtual public Ifpack2::Preconditioner<typename MatrixType::scalar_type,
00242                                            typename MatrixType::local_ordinal_type,
00243                                            typename MatrixType::global_ordinal_type,
00244                                            typename MatrixType::node_type>,
00245     virtual public Ifpack2::Details::CanChangeMatrix<Tpetra::RowMatrix<typename MatrixType::scalar_type,
00246                                                                        typename MatrixType::local_ordinal_type,
00247                                                                        typename MatrixType::global_ordinal_type,
00248                                                                        typename MatrixType::node_type> >
00249 {
00250  public:
00252   typedef typename MatrixType::scalar_type scalar_type;
00253 
00255   TEUCHOS_DEPRECATED typedef typename MatrixType::scalar_type Scalar;
00256 
00257 
00259   typedef typename MatrixType::local_ordinal_type local_ordinal_type;
00260 
00262   TEUCHOS_DEPRECATED typedef typename MatrixType::local_ordinal_type LocalOrdinal;
00263 
00264 
00266   typedef typename MatrixType::global_ordinal_type global_ordinal_type;
00267 
00269   TEUCHOS_DEPRECATED typedef typename MatrixType::global_ordinal_type GlobalOrdinal;
00270 
00271 
00273   typedef typename MatrixType::node_type node_type;
00274 
00276   TEUCHOS_DEPRECATED typedef typename MatrixType::node_type Node;
00277 
00279   typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitude_type;
00280 
00282   TEUCHOS_DEPRECATED typedef typename Teuchos::ScalarTraits<scalar_type>::magnitudeType magnitudeType;
00283 
00285   typedef Tpetra::RowMatrix<scalar_type,
00286                             local_ordinal_type,
00287                             global_ordinal_type,
00288                             node_type> row_matrix_type;
00289 
00291   typedef Tpetra::CrsMatrix<scalar_type,
00292                             local_ordinal_type,
00293                             global_ordinal_type,
00294                             node_type> crs_matrix_type;
00295 
00296   template <class NewMatrixType> friend class RILUK;
00297 
00301   RILUK (const Teuchos::RCP<const row_matrix_type>& A_in);
00302 
00310   RILUK (const Teuchos::RCP<const crs_matrix_type>& A_in);
00311 
00312  private:
00315   RILUK (const RILUK<MatrixType> & src);
00316 
00317  public:
00323   template <typename NewMatrixType>
00324   Teuchos::RCP< RILUK< NewMatrixType > >
00325   clone (const Teuchos::RCP<const NewMatrixType>& A_newnode) const;
00326 
00328   virtual ~RILUK ();
00329 
00334   void TEUCHOS_DEPRECATED SetRelaxValue (const magnitude_type RelaxValue) {
00335     RelaxValue_ = RelaxValue;
00336   }
00341   void TEUCHOS_DEPRECATED SetAbsoluteThreshold (const magnitude_type Athresh) {
00342     Athresh_ = Athresh;
00343   }
00348   void TEUCHOS_DEPRECATED SetRelativeThreshold (const magnitude_type Rthresh) {
00349     Rthresh_ = Rthresh;
00350   }
00355   void TEUCHOS_DEPRECATED SetOverlapMode (const Tpetra::CombineMode OverlapMode) {
00356     TEUCHOS_TEST_FOR_EXCEPTION(
00357       true, std::logic_error, "Ifpack2::RILUK::SetOverlapMode: "
00358       "RILUK no longer implements overlap on its own.  "
00359       "Use RILUK with AdditiveSchwarz if you want overlap.");
00360   }
00361 
00369   void setParameters (const Teuchos::ParameterList& params);
00370 
00372   void initialize ();
00373 
00382   void compute ();
00383 
00385   bool isInitialized () const {
00386     return isInitialized_;
00387   }
00389   bool isComputed () const {
00390     return isComputed_;
00391   }
00392 
00394   int getNumInitialize () const {
00395     return numInitialize_;
00396   }
00398   int getNumCompute () const {
00399     return numCompute_;
00400   }
00402   int getNumApply () const {
00403     return numApply_;
00404   }
00405 
00407   double getInitializeTime () const {
00408     return initializeTime_;
00409   }
00411   double getComputeTime () const {
00412     return computeTime_;
00413   }
00415   double getApplyTime () const {
00416     return applyTime_;
00417   }
00418 
00420 
00421 
00444   virtual void
00445   setMatrix (const Teuchos::RCP<const row_matrix_type>& A);
00446 
00448 
00449 
00450 
00452   std::string description () const;
00453 
00455 
00456 
00457 
00459   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00460   getDomainMap () const;
00461 
00463   Teuchos::RCP<const Tpetra::Map<local_ordinal_type,global_ordinal_type,node_type> >
00464   getRangeMap () const;
00465 
00495   void
00496   apply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00497          Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00498          Teuchos::ETransp mode = Teuchos::NO_TRANS,
00499          scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one (),
00500          scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero ()) const;
00502 
00503 private:
00525   void
00526   multiply (const Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& X,
00527             Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>& Y,
00528             const Teuchos::ETransp mode = Teuchos::NO_TRANS) const;
00529 public:
00530 
00540   magnitude_type TEUCHOS_DEPRECATED computeCondEst (Teuchos::ETransp mode) const;
00541 
00551   virtual magnitude_type TEUCHOS_DEPRECATED
00552   computeCondEst (CondestType CT = Ifpack2::Cheap,
00553                   local_ordinal_type MaxIters = 1550,
00554                   magnitude_type Tol = 1e-9,
00555                   const Teuchos::Ptr<const row_matrix_type>& Matrix = Teuchos::null);
00556 
00557   magnitude_type getCondEst () const { return Condest_; }
00558 
00560   Teuchos::RCP<const row_matrix_type> getMatrix () const;
00561 
00562   // Attribute access functions
00563 
00565   magnitude_type getRelaxValue () const { return RelaxValue_; }
00566 
00568   magnitude_type getAbsoluteThreshold () const { return Athresh_; }
00569 
00571   magnitude_type getRelativeThreshold () const {return Rthresh_;}
00572 
00574   int getLevelOfFill () const { return LevelOfFill_; }
00575 
00577   Tpetra::CombineMode getOverlapMode () {
00578     TEUCHOS_TEST_FOR_EXCEPTION(
00579       true, std::logic_error, "Ifpack2::RILUK::SetOverlapMode: "
00580       "RILUK no longer implements overlap on its own.  "
00581       "Use RILUK with AdditiveSchwarz if you want overlap.");
00582   }
00583 
00585   Tpetra::global_size_t getGlobalNumEntries () const {
00586     return getL ().getGlobalNumEntries () + getU ().getGlobalNumEntries ();
00587   }
00588 
00590   Teuchos::RCP<Ifpack2::IlukGraph<Tpetra::CrsGraph<local_ordinal_type,global_ordinal_type,node_type> > > getGraph () const {
00591     return Graph_;
00592   }
00593 
00595   const crs_matrix_type& getL () const;
00596 
00598   const Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type>&
00599   getD () const;
00600 
00602   const crs_matrix_type& getU () const;
00603 
00605   Teuchos::RCP<const crs_matrix_type> getCrsMatrix () const;
00606 
00607 private:
00608   typedef Tpetra::MultiVector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> MV;
00609   typedef Tpetra::Vector<scalar_type,local_ordinal_type,global_ordinal_type,node_type> vec_type;
00610   typedef Teuchos::ScalarTraits<scalar_type> STS;
00611   typedef Teuchos::ScalarTraits<magnitude_type> STM;
00612 
00613   void allocate_L_and_U();
00614   void initAllValues (const row_matrix_type& A);
00615 
00621   static Teuchos::RCP<const row_matrix_type>
00622   makeLocalFilter (const Teuchos::RCP<const row_matrix_type>& A);
00623 
00625   Teuchos::RCP<const row_matrix_type> A_;
00626 
00628   Teuchos::RCP<Ifpack2::IlukGraph<Tpetra::CrsGraph<local_ordinal_type,
00629                                                    global_ordinal_type,
00630                                                    node_type> > > Graph_;
00639   Teuchos::RCP<const crs_matrix_type> A_local_crs_;
00640 
00642   Teuchos::RCP<crs_matrix_type> L_;
00644   Teuchos::RCP<crs_matrix_type> U_;
00646   Teuchos::RCP<vec_type> D_;
00647 
00648   int LevelOfFill_;
00649 
00650   bool isAllocated_;
00651   bool isInitialized_;
00652   bool isComputed_;
00653 
00654   int numInitialize_;
00655   int numCompute_;
00656   mutable int numApply_;
00657 
00658   double initializeTime_;
00659   double computeTime_;
00660   mutable double applyTime_;
00661 
00662   magnitude_type RelaxValue_;
00663   magnitude_type Athresh_;
00664   magnitude_type Rthresh_;
00665 
00666   mutable magnitude_type Condest_;
00667 };
00668 
00669 //Set necessary local solve parameters when using ThrustGPU node
00670 namespace detail {
00671   template<class MatrixType, class NodeType, class MatSolveType>
00672   struct setLocalSolveParams{
00673     static Teuchos::RCP<Teuchos::ParameterList>
00674     setParams (const Teuchos::RCP<Teuchos::ParameterList>& param) {
00675       return param;
00676     }
00677   };
00678 #if defined(HAVE_KOKKOSCLASSIC_THRUST) && defined(HAVE_KOKKOSCLASSIC_CUSPARSE)
00679   template<class MatrixType, class Scalar>
00680   struct setLocalSolveParams<MatrixType,
00681                              KokkosClassic::ThrustGPUNode,
00682                              KokkosClassic::CUSPARSEOps<Scalar, KokkosClassic::ThrustGPUNode> >
00683   {
00684     static Teuchos::RCP<Teuchos::ParameterList>
00685     setParams (const Teuchos::RCP<Teuchos::ParameterList>& param) {
00686       param->sublist ("fillComplete").sublist ("Local Sparse Ops").set ("Prepare Solve", true);
00687       return param;
00688     }
00689   };
00690 #endif
00691 } //end namespace detail
00692 
00693 template <class MatrixType>
00694 template <typename NewMatrixType>
00695 Teuchos::RCP<RILUK<NewMatrixType> >
00696 RILUK<MatrixType>::
00697 clone (const Teuchos::RCP<const NewMatrixType>& A_newnode) const
00698 {
00699   using Teuchos::ParameterList;
00700   using Teuchos::RCP;
00701   using Teuchos::rcp;
00702   typedef typename NewMatrixType::node_type new_node_type;
00703   typedef typename NewMatrixType::mat_solve_type mat_solve_type;
00704   typedef RILUK<NewMatrixType> new_riluk_type;
00705 
00706   RCP<new_riluk_type> new_riluk = rcp (new new_riluk_type (A_newnode));
00707 
00708   RCP<ParameterList> plClone = Teuchos::parameterList ();
00709   plClone = detail::setLocalSolveParams<NewMatrixType,
00710     new_node_type, mat_solve_type>::setParams (plClone);
00711 
00712   RCP<new_node_type> new_node = A_newnode->getNode ();
00713   new_riluk->L_ = L_->clone (new_node, plClone);
00714   new_riluk->U_ = U_->clone (new_node, plClone);
00715   new_riluk->D_ = D_->clone (new_node);
00716 
00717   new_riluk->LevelOfFill_ = LevelOfFill_;
00718 
00719   new_riluk->isAllocated_ = isAllocated_;
00720   new_riluk->isInitialized_ = isInitialized_;
00721   new_riluk->isComputed_ = isComputed_;
00722 
00723   new_riluk->numInitialize_ = numInitialize_;
00724   new_riluk->numCompute_ = numCompute_;
00725   new_riluk->numApply_ =  numApply_;
00726 
00727   new_riluk->RelaxValue_ = RelaxValue_;
00728   new_riluk->Athresh_ = Athresh_;
00729   new_riluk->Rthresh_ = Rthresh_;
00730   new_riluk->Condest_ = Condest_;
00731 
00732   return new_riluk;
00733 }
00734 
00735 } // namespace Ifpack2
00736 
00737 #endif /* IFPACK2_CRSRILUK_DECL_HPP */
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends