IFPACK Development
Loading...
Searching...
No Matches
Ifpack_IKLU.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_IKLU_H
44#define IFPACK_IKLU_H
45
46#include "Ifpack_ConfigDefs.h"
47#include "Ifpack_CondestType.h"
48#include "Ifpack_ScalingType.h"
49#include "Ifpack_Preconditioner.h"
50#include "Ifpack_IKLU_Utils.h"
51#include "Epetra_Vector.h"
52#include "Epetra_CrsMatrix.h"
53#include "Epetra_Time.h"
54#include "Teuchos_RefCountPtr.hpp"
55
58class Epetra_Comm;
59class Epetra_Map;
61
62namespace Teuchos {
63 class ParameterList;
64}
65
67
80
81public:
82 // @{ Constructors and Destructors
85
87 virtual ~Ifpack_IKLU();
88
89 // @}
90 // @{ Construction methods
92 /* This method is only available if the Teuchos package is enabled.
93 This method recognizes five parameter names: level_fill, drop_tolerance,
94 absolute_threshold, relative_threshold and overlap_mode. These names are
95 case insensitive. For level_fill the ParameterEntry must have type int, the
96 threshold entries must have type double and overlap_mode must have type
97 Epetra_CombineMode.
98 */
99 int SetParameters(Teuchos::ParameterList& parameterlis);
100
102
108 int Initialize();
109
111 bool IsInitialized() const
112 {
113 return(IsInitialized_);
114 }
115
117
125 int Compute();
126
128 bool IsComputed() const {return(IsComputed_);};
129
130 // Mathematical functions.
131
133
141 int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
142
143 int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
144
146 double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
147 const int MaxIters = 1550,
148 const double Tol = 1e-9,
149 Epetra_RowMatrix* Matrix_in = 0);
150
152 double Condest() const
153 {
154 return(Condest_);
155 }
156
158
167 int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
168
170 double NormInf() const {return(0.0);};
171
173 bool HasNormInf() const {return(false);};
174
176 bool UseTranspose() const {return(UseTranspose_);};
177
179 const Epetra_Map & OperatorDomainMap() const {return(A_.OperatorDomainMap());};
180
182 const Epetra_Map & OperatorRangeMap() const{return(A_.OperatorRangeMap());};
183
185 const Epetra_Comm & Comm() const{return(Comm_);};
186
189 {
190 return(A_);
191 }
192
194 const Epetra_CrsMatrix & L() const {return(*L_);};
195
197 const Epetra_CrsMatrix & U() const {return(*U_);};
198
200 const char* Label() const
201 {
202 return(Label_.c_str());
203 }
204
206 int SetLabel(const char* Label_in)
207 {
208 Label_ = Label_in;
209 return(0);
210 }
211
213 virtual std::ostream& Print(std::ostream& os) const;
214
216 virtual int NumInitialize() const
217 {
218 return(NumInitialize_);
219 }
220
222 virtual int NumCompute() const
223 {
224 return(NumCompute_);
225 }
226
228 virtual int NumApplyInverse() const
229 {
230 return(NumApplyInverse_);
231 }
232
234 virtual double InitializeTime() const
235 {
236 return(InitializeTime_);
237 }
238
240 virtual double ComputeTime() const
241 {
242 return(ComputeTime_);
243 }
244
246 virtual double ApplyInverseTime() const
247 {
248 return(ApplyInverseTime_);
249 }
250
252 virtual double InitializeFlops() const
253 {
254 return(0.0);
255 }
256
257 virtual double ComputeFlops() const
258 {
259 return(ComputeFlops_);
260 }
261
262 virtual double ApplyInverseFlops() const
263 {
264 return(ApplyInverseFlops_);
265 }
266
267 inline double LevelOfFill() const {
268 return(LevelOfFill_);
269 }
270
272 inline double RelaxValue() const {
273 return(Relax_);
274 }
275
277 inline double AbsoluteThreshold() const
278 {
279 return(Athresh_);
280 }
281
283 inline double RelativeThreshold() const
284 {
285 return(Rthresh_);
286 }
287
289 inline double DropTolerance() const
290 {
291 return(DropTolerance_);
292 }
293
295#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
296 int NumGlobalNonzeros() const {
297 // FIXME: diagonal of L_ should not be stored
298 return(L().NumGlobalNonzeros() + U().NumGlobalNonzeros() - L().NumGlobalRows());
299 }
300#endif
301 long long NumGlobalNonzeros64() const {
302 // FIXME: diagonal of L_ should not be stored
303 return(L().NumGlobalNonzeros64() + U().NumGlobalNonzeros64() - L().NumGlobalRows64());
304 }
305
307 int NumMyNonzeros() const {
308 return(L().NumMyNonzeros() + U().NumMyNonzeros());
309 }
310
311private:
312
313 // @}
314 // @{ Internal methods
315
317 Ifpack_IKLU(const Ifpack_IKLU& RHS) :
318 A_(RHS.Matrix()),
319 Comm_(RHS.Comm()),
320 Time_(Comm())
321 {};
322
324 Ifpack_IKLU& operator=(const Ifpack_IKLU& /* RHS */)
325 {
326 return(*this);
327 }
328
330 void Destroy();
331
332 // @}
333 // @{ Internal data
334
336 const Epetra_RowMatrix& A_;
338 const Epetra_Comm& Comm_;
340 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
342 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
344 double Condest_;
346 double Relax_;
348 double Athresh_;
350 double Rthresh_;
352 double LevelOfFill_;
354 double DropTolerance_;
356 std::string Label_;
358 bool IsInitialized_;
360 bool IsComputed_;
362 bool UseTranspose_;
364 int NumMyRows_;
366 int NumMyNonzeros_;
368 int NumInitialize_;
370 int NumCompute_;
372 mutable int NumApplyInverse_;
374 double InitializeTime_;
376 double ComputeTime_;
378 mutable double ApplyInverseTime_;
380 double ComputeFlops_;
382 mutable double ApplyInverseFlops_;
384 mutable Epetra_Time Time_;
386 long long GlobalNonzeros_;
387 Teuchos::RefCountPtr<Epetra_SerialComm> SerialComm_;
388 Teuchos::RefCountPtr<Epetra_Map> SerialMap_;
389
391 csr* csrA_;
392 css* cssS_;
394 csrn* csrnN_;
395
396}; // Ifpack_IKLU
397
398#endif /* IFPACK_IKLU_H */
Ifpack_ScalingType enumerable type.
virtual const Epetra_Map & OperatorDomainMap() const=0
virtual const Epetra_Map & OperatorRangeMap() const=0
Ifpack_IKLU: A class for constructing and using an incomplete LU factorization of a given Epetra_RowM...
Definition Ifpack_IKLU.h:79
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
virtual int NumCompute() const
Returns the number of calls to Compute().
double DropTolerance() const
Gets the dropping tolerance.
virtual ~Ifpack_IKLU()
Ifpack_IKLU Destructor.
int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied.
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
int Compute()
Compute IC factor U using the specified graph, diagonal perturbation thresholds and relaxation parame...
int SetParameters(Teuchos::ParameterList &parameterlis)
Set parameters using a Teuchos::ParameterList object.
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
const char * Label() const
Returns the label of this object.
int NumGlobalNonzeros() const
Returns the number of nonzero entries in the global graph.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
double AbsoluteThreshold() const
Get absolute threshold value.
int NumMyNonzeros() const
Returns the number of nonzero entries in the local graph.
Ifpack_IKLU(const Epetra_RowMatrix *A)
Ifpack_IKLU constuctor with variable number of indices per row.
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
double RelativeThreshold() const
Get relative threshold value.
virtual double InitializeTime() const
Returns the time spent in Initialize().
bool UseTranspose() const
Returns the current UseTranspose setting.
double Condest() const
Returns the computed estimated condition number, or -1.0 if no computed.
virtual double ComputeTime() const
Returns the time spent in Compute().
int Initialize()
Initialize L and U with values from user matrix A.
const Epetra_CrsMatrix & L() const
Returns a reference to the L factor.
const Epetra_CrsMatrix & U() const
Returns a reference to the U factor.
double RelaxValue() const
Set relative threshold value.
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Ifpack_IKLU forward/back solve on a Epetra_MultiVector X in Y.
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
virtual std::ostream & Print(std::ostream &os) const
Prints basic information on iostream. This function is used by operator<<.
int SetLabel(const char *Label_in)
Sets the label for this object.
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.