Ifpack2 Templated Preconditioning Package  Version 1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
Ifpack2_Details_FastILU_Base_decl.hpp
Go to the documentation of this file.
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2009) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ***********************************************************************
40 //@HEADER
41 */
42 
47 
48 #ifndef __IFPACK2_FASTILU_BASE_DECL_HPP__
49 #define __IFPACK2_FASTILU_BASE_DECL_HPP__
50 
51 #include <Tpetra_RowMatrix.hpp>
52 #include <Tpetra_CrsMatrix.hpp>
53 #include <Tpetra_KokkosCompat_DefaultNode.hpp>
54 #include <KokkosSparse_CrsMatrix.hpp>
57 #include <shylu_fastutil.hpp>
58 
59 #ifdef HAVE_IFPACK2_METIS
60 #include "metis.h"
61 #endif
62 
63 namespace Ifpack2
64 {
65 namespace Details
66 {
67 
70 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
71  class FastILU_Base : public Ifpack2::Preconditioner<Scalar, LocalOrdinal, GlobalOrdinal, Node>,
73  Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
74 {
75  public:
77  typedef typename Node::device_type device_type;
79  typedef typename device_type::execution_space execution_space;
81  typedef typename Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>::impl_scalar_type ImplScalar;
83  typedef Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> TRowMatrix;
85  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> TCrsMatrix;
87  typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> TMultiVec;
89  typedef KokkosSparse::CrsMatrix<Scalar, LocalOrdinal, execution_space> KCrsMatrix;
91  typedef Kokkos::View<LocalOrdinal *, execution_space> OrdinalArray;
93  typedef typename Kokkos::View<LocalOrdinal *, execution_space>::HostMirror OrdinalArrayHost;
95  typedef Kokkos::View< ImplScalar *, execution_space> ImplScalarArray;
96  typedef Kokkos::View< Scalar *, execution_space> ScalarArray;
97  typedef Kokkos::View<const Scalar *, execution_space> ConstScalarArray;
98  #ifdef HAVE_IFPACK2_METIS
99  typedef Kokkos::View<idx_t*, Kokkos::HostSpace> MetisArrayHost;
100  #endif
101 
104 
107  getDomainMap () const;
108 
111  getRangeMap () const;
112 
114  void
115  apply (const TMultiVec& X,
116  TMultiVec& Y,
118  Scalar alpha = Teuchos::ScalarTraits<Scalar>::one(),
119  Scalar beta = Teuchos::ScalarTraits<Scalar>::zero()) const;
121 
132  void setParameters (const Teuchos::ParameterList& List);
133 
134  bool isBlockCrs() const;
135 
137  void initialize();
138 
140  bool isInitialized() const;
141 
143  void compute();
144 
146  bool isComputed() const;
147 
150 
152  int getNumInitialize() const;
153 
155  int getNumCompute() const;
156 
158  int getNumApply() const;
159 
161  double getInitializeTime() const;
162 
164  double getComputeTime() const;
165 
167  double getApplyTime() const;
168 
170  double getCopyTime() const;
171 
173  virtual int getSweeps() const = 0;
174 
176  virtual std::string getSpTrsvType() const = 0;
177 
179  virtual int getNTrisol() const = 0;
180 
182  virtual void checkLocalILU() const;
183 
185  virtual void checkLocalIC() const;
186 
188  std::string description() const;
189 
194 
195  protected:
197  bool initFlag_;
198  bool computedFlag_;
199  int nInit_;
200  int nComputed_;
201  mutable int nApply_;
202  //store the local CRS components
203  ImplScalarArray localValues_; //set at beginning of compute()
204  OrdinalArray localRowPtrs_; //set in initialize()
205  OrdinalArray localColInds_; //set in initialize()
206  OrdinalArrayHost localRowPtrsHost_; //set in initialize() and used to get localValues_ in compute()
207 
208  double initTime_;
209  double computeTime_;
210  mutable double applyTime_;
211  double crsCopyTime_; //total time spent deep copying values, rowptrs, colinds out of mat
212 
213  //Store validated parameter values (instead of keeping a ParameterList around)
214  struct Params
215  {
216  Params() {}
217  Params(const Teuchos::ParameterList& pL, std::string precType);
218  bool use_metis;
219  FastILU::SpTRSV sptrsv_algo;
220  int nFact;
221  int nTrisol;
222  int level;
223  int blkSize;
224  double omega;
225  double shift;
226  bool guessFlag;
227  int blockSizeILU;
228  int blockSize;
229  bool blockCrs;
230  int blockCrsSize;
231  bool fillBlocks;
232  static Params getDefaults();
233  };
234 
235  Params params_;
236 
237  #ifdef HAVE_IFPACK2_METIS
238  MetisArrayHost metis_perm_;
239  MetisArrayHost metis_iperm_;
240  #endif
241 
243  // \pre !mat_.is_null()
244  virtual void initLocalPrec() = 0;
246  virtual void computeLocalPrec() = 0;
248  virtual void applyLocalPrec(ImplScalarArray x, ImplScalarArray y) const = 0;
250  virtual std::string getName() const = 0;
251 };
252 
253 } //namespace Details
254 } //namespace Ifpack2
255 
256 #endif
257 
virtual std::string getName() const =0
Get the name of the underlying preconditioner (&quot;Filu&quot;, &quot;Fildl&quot; or &quot;Fic&quot;)
Mix-in interface for preconditioners that can change their matrix after construction.
Definition: Ifpack2_Details_CanChangeMatrix.hpp:93
int getNumCompute() const
Get the number of times compute() was called.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:338
virtual void checkLocalIC() const
Verify and print debug information about the underlying IC preconditioner.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:388
double getComputeTime() const
Get the time spent in the last compute() call.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:359
Kokkos::View< LocalOrdinal *, execution_space >::HostMirror OrdinalArrayHost
Array of LocalOrdinal on host.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:93
virtual int getNTrisol() const =0
Get the &quot;triangular solve iterations&quot; parameter.
virtual void applyLocalPrec(ImplScalarArray x, ImplScalarArray y) const =0
Apply the local preconditioner with 1-D views of the local parts of X and Y (one vector only) ...
device_type::execution_space execution_space
Kokkos execution space.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:79
double getCopyTime() const
Get the time spent deep copying local 3-array CRS out of the matrix.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:373
virtual int getSweeps() const =0
Get the &quot;sweeps&quot; parameter.
double getInitializeTime() const
Get the time spent in the last initialize() call.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:352
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Get the range map of the matrix.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:88
double getApplyTime() const
Get the time spent in the last apply() call.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:366
Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TCrsMatrix
Tpetra CRS matrix.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:85
Teuchos::RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Get the domain map of the matrix.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:80
Teuchos::RCP< const TRowMatrix > getMatrix() const
Get the current matrix.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:324
virtual void checkLocalILU() const
Verify and print debug information about the underlying ILU preconditioner (only supported if this is...
Definition: Ifpack2_Details_FastILU_Base_def.hpp:380
Kokkos::View< LocalOrdinal *, execution_space > OrdinalArray
Array of LocalOrdinal on device.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:91
int getNumApply() const
Get the number of times apply() was called.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:345
The base class of the Ifpack2 FastILU wrappers (Filu, Fildl and Fic)
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:71
void setParameters(const Teuchos::ParameterList &List)
Validate parameters, and set defaults when parameters are not provided.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:153
std::string description() const
Return a brief description of the preconditioner, in YAML format.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:395
void compute()
Compute the preconditioner.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:291
Node::device_type device_type
Kokkos device type.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:77
Interface for all Ifpack2 preconditioners.
Definition: Ifpack2_Preconditioner.hpp:107
Kokkos::View< ImplScalar *, execution_space > ImplScalarArray
Array of Scalar on device.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:95
FastILU_Base(Teuchos::RCP< const TRowMatrix > mat_)
Constructor.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:64
bool isComputed() const
Whether compute() has been called since the last time the matrix&#39;s values or structure were changed...
Definition: Ifpack2_Details_FastILU_Base_def.hpp:316
Declaration of interface for preconditioners that can change their matrix after construction.
bool isInitialized() const
Whether initialize() has been called since the last time the matrix&#39;s structure was changed...
Definition: Ifpack2_Details_FastILU_Base_def.hpp:284
KokkosSparse::CrsMatrix< Scalar, LocalOrdinal, execution_space > KCrsMatrix
Kokkos CRS matrix.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:89
void initialize()
Initialize the preconditioner.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:168
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TMultiVec
Tpetra multivector.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:87
Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >::impl_scalar_type ImplScalar
Kokkos scalar type.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:81
Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TRowMatrix
Tpetra row matrix.
Definition: Ifpack2_Details_FastILU_Base_decl.hpp:83
void setMatrix(const Teuchos::RCP< const TRowMatrix > &A)
Definition: Ifpack2_Details_FastILU_Base_def.hpp:421
virtual void initLocalPrec()=0
Construct the underlying preconditioner (localPrec_) using given params and then call localPrec_-&gt;ini...
void apply(const TMultiVec &X, TMultiVec &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
Apply the preconditioner.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:95
int getNumInitialize() const
Get the number of times initialize() was called.
Definition: Ifpack2_Details_FastILU_Base_def.hpp:331
virtual void computeLocalPrec()=0
Get values array from the matrix and then call compute() on the underlying preconditioner.
virtual std::string getSpTrsvType() const =0
Get the name of triangular solve algorithm.