IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
Ifpack_SILU.h
1 /*@HEADER
2 // ***********************************************************************
3 //
4 // Ifpack: Object-Oriented Algebraic Preconditioner Package
5 // Copyright (2002) 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 
43 #ifndef IFPACK_SILU_H
44 #define IFPACK_SILU_H
45 
46 #include "Ifpack_ConfigDefs.h"
47 
48 #ifdef HAVE_IFPACK_SUPERLU
49 #include "Ifpack_Preconditioner.h"
50 #include "Ifpack_Condest.h"
51 #include "Ifpack_ScalingType.h"
52 #include "Ifpack_IlukGraph.h"
53 #include "Epetra_CompObject.h"
54 #include "Epetra_MultiVector.h"
55 #include "Epetra_Vector.h"
56 #include "Epetra_CrsGraph.h"
57 #include "Epetra_CrsMatrix.h"
58 #include "Epetra_BlockMap.h"
59 #include "Epetra_Map.h"
60 #include "Epetra_Object.h"
61 #include "Epetra_Comm.h"
62 #include "Epetra_RowMatrix.h"
63 #include "Epetra_Time.h"
64 #include "Teuchos_RefCountPtr.hpp"
65 
66 namespace Teuchos {
67  class ParameterList;
68 }
69 
70 // SuperLU includes
71 #include "slu_ddefs.h"
72 
73 
84 
85 public:
87 
90 
93  {
94  Destroy();
95  }
96 
98 
100 
102  int Initialize();
103 
105  bool IsInitialized() const
106  {
107  return(IsInitialized_);
108  }
109 
111  int Compute();
112 
114  bool IsComputed() const
115  {
116  return(IsComputed_);
117  }
118 
120  /* This method is only available if the Teuchos package is enabled.
121  This method recognizes four parameter names: relax_value,
122  absolute_threshold, relative_threshold and overlap_mode. These names are
123  case insensitive, and in each case except overlap_mode, the ParameterEntry
124  must have type double. For overlap_mode, the ParameterEntry must have
125  type Epetra_CombineMode.
126  */
127  int SetParameters(Teuchos::ParameterList& parameterlist);
128 
130 
139  int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
141 
143 
144  // Applies the matrix to X, returns the result in Y.
145  int Apply(const Epetra_MultiVector& X,
146  Epetra_MultiVector& Y) const
147  {
148  return(Multiply(false,X,Y));
149  }
150 
151  int Multiply(bool Trans, const Epetra_MultiVector& X,
152  Epetra_MultiVector& Y) const;
153 
155 
168  int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
169 
171  double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
172  const int MaxIters = 1550,
173  const double Tol = 1e-9,
174  Epetra_RowMatrix* Matrix_in = 0);
175 
177  double Condest() const
178  {
179  return(Condest_);
180  }
181 
183 
185 
186 
188  const char* Label() const {return(Label_);}
189 
191  int SetLabel(const char* Label_in)
192  {
193  strcpy(Label_,Label_in);
194  return(0);
195  }
196 
198  double NormInf() const {return(0.0);};
199 
201  bool HasNormInf() const {return(false);};
202 
204  bool UseTranspose() const {return(UseTranspose_);};
205 
207  const Epetra_Map & OperatorDomainMap() const {return(A_->OperatorDomainMap());};
208 
210  const Epetra_Map & OperatorRangeMap() const{return(A_->OperatorRangeMap());};
211 
213  const Epetra_Comm & Comm() const{return(Comm_);};
214 
216  const Epetra_RowMatrix& Matrix() const
217  {
218  return(*A_);
219  }
220 
222  virtual std::ostream& Print(std::ostream& os) const;
223 
225  virtual int NumInitialize() const
226  {
227  return(NumInitialize_);
228  }
229 
231  virtual int NumCompute() const
232  {
233  return(NumCompute_);
234  }
235 
237  virtual int NumApplyInverse() const
238  {
239  return(NumApplyInverse_);
240  }
241 
243  virtual double InitializeTime() const
244  {
245  return(InitializeTime_);
246  }
247 
249  virtual double ComputeTime() const
250  {
251  return(ComputeTime_);
252  }
253 
255  virtual double ApplyInverseTime() const
256  {
257  return(ApplyInverseTime_);
258  }
259 
260  virtual double InitializeFlops() const
261  {
262  return(0.0);
263  }
264 
265  virtual double ComputeFlops() const
266  {
267  return(0.0);
268  }
269 
270  virtual double ApplyInverseFlops() const
271  {
272  return(0.0);
273  }
274 
275 private:
276 
278 
280 
282  Ifpack_SILU(const Ifpack_SILU& RHS):
283  Comm_(RHS.Comm()),
284  Time_(RHS.Comm())
285  {}
286 
288  Ifpack_SILU& operator=(const Ifpack_SILU& RHS)
289  {
290  return(*this);
291  }
292 
294  void Destroy();
295 
297 
307  int Solve(bool Trans, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
308 
309  int InitAllValues(const Epetra_RowMatrix & A, int MaxNumEntries);
310 
312  double DropTol() const {return DropTol_;}
313 
315  double FillTol() const{return FillTol_;}
316 
318  double FillFactor() const{return FillFactor_;}
319 
321  int DropRule() const{return DropRule_;}
322 
323 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
324  int NumGlobalRows() const {return(Graph().NumGlobalRows());};
326 
328  int NumGlobalCols() const {return(Graph().NumGlobalCols());};
329 
331  int NumGlobalNonzeros() const {return(Graph().NumGlobalNonzeros());};
332 
334  virtual int NumGlobalBlockDiagonals() const {return(Graph().NumGlobalBlockDiagonals());};
335 #endif
336 
338  long long NumGlobalRows64() const {return(Graph().NumGlobalRows64());};
339 
341  long long NumGlobalCols64() const {return(Graph().NumGlobalCols64());};
342 
344  long long NumGlobalNonzeros64() const {return(Graph().NumGlobalNonzeros64());};
345 
347  virtual long long NumGlobalBlockDiagonals64() const {return(Graph().NumGlobalBlockDiagonals64());};
348 
350  int NumMyRows() const {return(Graph().NumMyRows());};
351 
353  int NumMyCols() const {return(Graph().NumMyCols());};
354 
356  int NumMyNonzeros() const {return(Graph().NumMyNonzeros());};
357 
359  virtual int NumMyBlockDiagonals() const {return(Graph().NumMyBlockDiagonals());};
360 
362  virtual int NumMyDiagonals() const {return(NumMyDiagonals_);};
363 
364 #ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
365  int IndexBase() const {return(Graph().IndexBase());};
367 #endif
368  long long IndexBase64() const {return(Graph().IndexBase64());};
369 
371  const Epetra_CrsGraph & Graph() const {return(*Graph_);};
372 
375  {
376  return(*A_);
377  }
378 
380 
382 
383 
385  Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
386  Teuchos::RefCountPtr<Epetra_CrsGraph> Graph_;
387  Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
388  Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
389  Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
390  const Epetra_Comm & Comm_;
392  Teuchos::RefCountPtr<Epetra_CrsMatrix> Aover_;
393  bool UseTranspose_;
394 
395  int NumMyDiagonals_;
396  bool Allocated_;
397  bool ValuesInitialized_;
398  bool Factored_;
399 
401  double DropTol_;
403  double FillTol_;
405  double FillFactor_;
407  int DropRule_;
408 
410  double Condest_;
412  bool IsInitialized_;
414  bool IsComputed_;
416  char Label_[160];
418  int NumInitialize_;
420  int NumCompute_;
422  mutable int NumApplyInverse_;
424  double InitializeTime_;
426  double ComputeTime_;
428  mutable double ApplyInverseTime_;
430  mutable Epetra_Time Time_;
432 #ifdef HAVE_IFPACK_SUPERLU5_API
433  mutable GlobalLU_t lu_;
434 #endif
435  mutable SuperLUStat_t stat_;
438  mutable superlu_options_t options_;
440  mutable SuperMatrix SA_,SAc_,SL_,SU_,SY_;
442  int *etree_,*perm_r_,*perm_c_;
444 
445 
446  template<typename int_type>
447  int TInitialize();
448 };
449 
450 #endif /* HAVE_IFPACK_SUPERLU */
451 #endif /* IFPACK_ILU_H */
int SetLabel(const char *Label_in)
Sets label for this object.
Definition: Ifpack_SILU.h:191
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
Definition: Ifpack_SILU.h:255
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
Definition: Ifpack_SILU.h:260
const char * Label() const
Returns a character string describing the operator.
Definition: Ifpack_SILU.h:188
Ifpack_SILU(Epetra_RowMatrix *A)
Constructor.
Definition: Ifpack_SILU.cpp:70
virtual int NumInitialize() const
Returns the number of calls to Initialize().
Definition: Ifpack_SILU.h:225
bool UseTranspose() const
Returns the current UseTranspose setting.
Definition: Ifpack_SILU.h:204
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
Definition: Ifpack_SILU.h:216
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y...
A wrapper to SuperLU 4.0&#39;s supernodal ILUT w/ partial pivoting.
Definition: Ifpack_SILU.h:83
virtual double InitializeTime() const
Returns the time spent in Initialize().
Definition: Ifpack_SILU.h:243
~Ifpack_SILU()
Destructor.
Definition: Ifpack_SILU.h:92
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
Definition: Ifpack_SILU.h:105
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
Definition: Ifpack_SILU.h:201
virtual std::ostream & Print(std::ostream &os) const
Prints on stream basic information about this object.
Ifpack_ScalingType enumerable type.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
int Compute()
Compute ILU factors L and U using the specified graph, diagonal perturbation thresholds and relaxatio...
int Initialize()
Initialize the preconditioner, does not touch matrix values.
int SetParameters(Teuchos::ParameterList &parameterlist)
Set parameters using a Teuchos::ParameterList object.
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
Definition: Ifpack_SILU.h:207
int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied.
Definition: Ifpack_SILU.h:139
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
Definition: Ifpack_SILU.h:114
virtual double ComputeTime() const
Returns the time spent in Compute().
Definition: Ifpack_SILU.h:249
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
Definition: Ifpack_SILU.h:270
virtual int NumCompute() const
Returns the number of calls to Compute().
Definition: Ifpack_SILU.h:231
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
Definition: Ifpack_SILU.h:210
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
Definition: Ifpack_SILU.h:213
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
Definition: Ifpack_SILU.h:237
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
Definition: Ifpack_SILU.h:198
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
Definition: Ifpack_SILU.h:265
double Condest() const
Returns the computed estimated condition number, or -1.0 if not computed.
Definition: Ifpack_SILU.h:177