Ifpack Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
Ifpack_ILU.h
Go to the documentation of this file.
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_ILU_H
44#define IFPACK_ILU_H
45
46#include "Ifpack_ConfigDefs.h"
48#include "Ifpack_Condest.h"
49#include "Ifpack_ScalingType.h"
50#include "Ifpack_IlukGraph.h"
51#include "Epetra_CompObject.h"
52#include "Epetra_MultiVector.h"
53#include "Epetra_Vector.h"
54#include "Epetra_CrsGraph.h"
55#include "Epetra_CrsMatrix.h"
56#include "Epetra_BlockMap.h"
57#include "Epetra_Map.h"
58#include "Epetra_Object.h"
59#include "Epetra_Comm.h"
60#include "Epetra_RowMatrix.h"
61#include "Epetra_Time.h"
62#include "Teuchos_RefCountPtr.hpp"
63
64namespace Teuchos {
65 class ParameterList;
66}
67
69
84
85public:
86 // @{ Constructors and destructors.
89
92 {
93 Destroy();
94 }
95
96 // @}
97 // @{ Construction methods
98
100 int Initialize();
101
103 bool IsInitialized() const
104 {
105 return(IsInitialized_);
106 }
107
109
117 int Compute();
118
120 bool IsComputed() const
121 {
122 return(IsComputed_);
123 }
124
126 /* This method is only available if the Teuchos package is enabled.
127 This method recognizes four parameter names: relax_value,
128 absolute_threshold, relative_threshold and overlap_mode. These names are
129 case insensitive, and in each case except overlap_mode, the ParameterEntry
130 must have type double. For overlap_mode, the ParameterEntry must have
131 type Epetra_CombineMode.
132 */
133 int SetParameters(Teuchos::ParameterList& parameterlist);
134
136
144 int SetUseTranspose(bool UseTranspose_in) {UseTranspose_ = UseTranspose_in; return(0);};
145 // @}
146
147 // @{ Mathematical functions.
148 // Applies the matrix to X, returns the result in Y.
150 Epetra_MultiVector& Y) const
151 {
152 return(Multiply(false,X,Y));
153 }
154
155 int Multiply(bool Trans, const Epetra_MultiVector& X,
156 Epetra_MultiVector& Y) const;
157
159
172 int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
173
175 double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
176 const int MaxIters = 1550,
177 const double Tol = 1e-9,
178 Epetra_RowMatrix* Matrix_in = 0);
179
181 double Condest() const
182 {
183 return(Condest_);
184 }
185
186 // @}
187 // @{ Query methods
188
190 const Epetra_CrsMatrix & L() const {return(*L_);};
191
193 const Epetra_Vector & D() const {return(*D_);};
194
196 const Epetra_CrsMatrix & U() const {return(*U_);};
197
199 const char* Label() const {return(Label_);}
200
202 int SetLabel(const char* Label_in)
203 {
204 strcpy(Label_,Label_in);
205 return(0);
206 }
207
209 double NormInf() const {return(0.0);};
210
212 bool HasNormInf() const {return(false);};
213
215 bool UseTranspose() const {return(UseTranspose_);};
216
218 const Epetra_Map & OperatorDomainMap() const {return(U_->OperatorDomainMap());};
219
221 const Epetra_Map & OperatorRangeMap() const{return(L_->OperatorRangeMap());};
222
224 const Epetra_Comm & Comm() const{return(Comm_);};
225
228 {
229 return(*A_);
230 }
231
233 virtual std::ostream& Print(std::ostream& os) const;
234
236 virtual int NumInitialize() const
237 {
238 return(NumInitialize_);
239 }
240
242 virtual int NumCompute() const
243 {
244 return(NumCompute_);
245 }
246
248 virtual int NumApplyInverse() const
249 {
250 return(NumApplyInverse_);
251 }
252
254 virtual double InitializeTime() const
255 {
256 return(InitializeTime_);
257 }
258
260 virtual double ComputeTime() const
261 {
262 return(ComputeTime_);
263 }
264
266 virtual double ApplyInverseTime() const
267 {
268 return(ApplyInverseTime_);
269 }
270
272 virtual double InitializeFlops() const
273 {
274 return(0.0);
275 }
276
277 virtual double ComputeFlops() const
278 {
279 return(ComputeFlops_);
280 }
281
282 virtual double ApplyInverseFlops() const
283 {
284 return(ApplyInverseFlops_);
285 }
286
287private:
288
289 // @}
290 // @{ Private methods
291
294 Comm_(RHS.Comm()),
295 Time_(RHS.Comm())
296 {}
297
299 Ifpack_ILU& operator=(const Ifpack_ILU& /* RHS */)
300 {
301 return(*this);
302 }
303
305 void Destroy();
306
308
318 int Solve(bool Trans, const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
319
320 int ComputeSetup();
321 int InitAllValues(const Epetra_RowMatrix & A, int MaxNumEntries);
322
324 int LevelOfFill() const {return LevelOfFill_;}
325
327 double RelaxValue() const {return RelaxValue_;}
328
330 double AbsoluteThreshold() const {return Athresh_;}
331
333 double RelativeThreshold() const {return Rthresh_;}
334
335#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
337 int NumGlobalRows() const {return(Graph().NumGlobalRows());};
338
340 int NumGlobalCols() const {return(Graph().NumGlobalCols());};
341
344
346 virtual int NumGlobalBlockDiagonals() const {return(Graph().NumGlobalBlockDiagonals());};
347#endif
348
349 long long NumGlobalRows64() const {return(Graph().NumGlobalRows64());};
350
351 long long NumGlobalCols64() const {return(Graph().NumGlobalCols64());};
352
353 long long NumGlobalNonzeros64() const {return(L().NumGlobalNonzeros64()+U().NumGlobalNonzeros64());};
354
355 virtual long long NumGlobalBlockDiagonals64() const {return(Graph().NumGlobalBlockDiagonals64());};
356
358 int NumMyRows() const {return(Graph().NumMyRows());};
359
361 int NumMyCols() const {return(Graph().NumMyCols());};
362
364 int NumMyNonzeros() const {return(L().NumMyNonzeros()+U().NumMyNonzeros());};
365
367 virtual int NumMyBlockDiagonals() const {return(Graph().NumMyBlockDiagonals());};
368
370 virtual int NumMyDiagonals() const {return(NumMyDiagonals_);};
371
373#ifndef EPETRA_NO_32BIT_GLOBAL_INDICES
374 int IndexBase() const {return(Graph().IndexBase());};
375#endif
376 long long IndexBase64() const {return(Graph().IndexBase64());};
377
379 const Ifpack_IlukGraph & Graph() const {return(*Graph_);};
380
383 {
384 return(*A_);
385 }
386
387 // @}
388 // @{ Internal data
389
391 Teuchos::RefCountPtr<Epetra_RowMatrix> A_;
392 Teuchos::RefCountPtr<Ifpack_IlukGraph> Graph_;
393 Teuchos::RefCountPtr<Epetra_CrsGraph> CrsGraph_;
394 Teuchos::RefCountPtr<Epetra_Map> IlukRowMap_;
395 Teuchos::RefCountPtr<Epetra_Map> IlukDomainMap_;
396 Teuchos::RefCountPtr<Epetra_Map> IlukRangeMap_;
401 Teuchos::RefCountPtr<Epetra_CrsMatrix> L_;
403 Teuchos::RefCountPtr<Epetra_CrsMatrix> U_;
404 Teuchos::RefCountPtr<Epetra_CrsGraph> L_Graph_;
405 Teuchos::RefCountPtr<Epetra_CrsGraph> U_Graph_;
407 Teuchos::RefCountPtr<Epetra_Vector> D_;
409
417 double Athresh_;
419 double Rthresh_;
421 double Condest_;
429 char Label_[160];
435 mutable int NumApplyInverse_;
441 mutable double ApplyInverseTime_;
445 mutable double ApplyInverseFlops_;
448
449};
450
451#endif /* IFPACK_ILU_H */
Ifpack_CondestType
Ifpack_CondestType: enum to define the type of condition number estimate.
@ Ifpack_Cheap
cheap estimate
Ifpack_ScalingType enumerable type.
#define RHS(a)
Definition MatGenFD.c:60
Ifpack_ILU: A class for constructing and using an incomplete lower/upper (ILU) factorization of a giv...
Definition Ifpack_ILU.h:83
virtual double InitializeTime() const
Returns the time spent in Initialize().
Definition Ifpack_ILU.h:254
const Epetra_CrsMatrix & U() const
Returns the address of the L factor associated with this factored matrix.
Definition Ifpack_ILU.h:196
bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized.
Definition Ifpack_ILU.h:103
double ComputeFlops_
Contains the number of flops for Compute().
Definition Ifpack_ILU.h:443
int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Definition Ifpack_ILU.h:149
virtual int NumMyDiagonals() const
Returns the number of nonzero diagonal values found in matrix.
Definition Ifpack_ILU.h:370
long long NumGlobalRows64() const
Definition Ifpack_ILU.h:349
virtual double ApplyInverseFlops() const
Returns the number of flops in the application of the preconditioner.
Definition Ifpack_ILU.h:282
const Epetra_Map * L_RangeMap_
Definition Ifpack_ILU.h:398
double Athresh_
absolute threshold
Definition Ifpack_ILU.h:417
int NumGlobalRows() const
Returns the number of global matrix rows.
Definition Ifpack_ILU.h:337
int SetLabel(const char *Label_in)
Sets label for this object.
Definition Ifpack_ILU.h:202
Teuchos::RefCountPtr< Epetra_Map > IlukRangeMap_
Definition Ifpack_ILU.h:396
int NumGlobalNonzeros() const
Returns the number of nonzero entries in the global graph.
Definition Ifpack_ILU.h:343
bool ValuesInitialized_
Definition Ifpack_ILU.h:412
double ComputeTime_
Contains the time for all successful calls to Compute().
Definition Ifpack_ILU.h:439
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
Definition Ifpack_ILU.h:272
bool UseTranspose_
Definition Ifpack_ILU.h:408
Teuchos::RefCountPtr< Epetra_CrsMatrix > U_
Contains the U factors.
Definition Ifpack_ILU.h:403
const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
Definition Ifpack_ILU.h:221
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
Definition Ifpack_ILU.h:277
long long NumGlobalCols64() const
Definition Ifpack_ILU.h:351
virtual double ComputeTime() const
Returns the time spent in Compute().
Definition Ifpack_ILU.h:260
double RelaxValue_
Relaxation value.
Definition Ifpack_ILU.h:415
Teuchos::RefCountPtr< Ifpack_IlukGraph > Graph_
Definition Ifpack_ILU.h:392
int SetParameters(Teuchos::ParameterList &parameterlist)
Set parameters using a Teuchos::ParameterList object.
Ifpack_ILU & operator=(const Ifpack_ILU &)
operator= (should never be used)
Definition Ifpack_ILU.h:299
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
Definition Ifpack_ILU.h:248
const Epetra_RowMatrix & Matrix() const
Returns a reference to the matrix to be preconditioned.
Definition Ifpack_ILU.h:227
const Epetra_Map * U_DomainMap_
Definition Ifpack_ILU.h:397
int Compute()
Compute ILU factors L and U using the specified graph, diagonal perturbation thresholds and relaxatio...
Ifpack_ILU(const Ifpack_ILU &RHS)
Copy constructor (should never be used)
Definition Ifpack_ILU.h:293
int Solve(bool Trans, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Ifpack_ILU forward/back solve on a Epetra_MultiVector X in Y.
const Epetra_Comm & Comm() const
Returns the Epetra_BlockMap object associated with the range of this matrix operator.
Definition Ifpack_ILU.h:224
virtual std::ostream & Print(std::ostream &os) const
Prints on stream basic information about this object.
const Epetra_CrsMatrix & L() const
Returns the address of the L factor associated with this factored matrix.
Definition Ifpack_ILU.h:190
bool Factored_
Definition Ifpack_ILU.h:413
Epetra_RowMatrix & Matrix()
Returns a reference to the matrix.
Definition Ifpack_ILU.h:382
const Epetra_Comm & Comm_
Definition Ifpack_ILU.h:399
Teuchos::RefCountPtr< Epetra_RowMatrix > A_
Pointer to the Epetra_RowMatrix to factorize.
Definition Ifpack_ILU.h:391
int NumApplyInverse_
Contains the number of successful call to ApplyInverse().
Definition Ifpack_ILU.h:435
Teuchos::RefCountPtr< Epetra_CrsMatrix > L_
Contains the L factors.
Definition Ifpack_ILU.h:401
int ComputeSetup()
virtual int NumMyBlockDiagonals() const
Returns the number of diagonal entries found in the local input graph.
Definition Ifpack_ILU.h:367
double AbsoluteThreshold() const
Get absolute threshold value.
Definition Ifpack_ILU.h:330
int IndexBase() const
Returns the index base for row and column indices for this graph.
Definition Ifpack_ILU.h:374
bool IsComputed_
If true, the preconditioner has been successfully computed.
Definition Ifpack_ILU.h:427
Teuchos::RefCountPtr< Epetra_Map > IlukDomainMap_
Definition Ifpack_ILU.h:395
const Epetra_Vector & D() const
Returns the address of the D factor associated with this factored matrix.
Definition Ifpack_ILU.h:193
double RelativeThreshold() const
Get relative threshold value.
Definition Ifpack_ILU.h:333
int NumMyRows() const
Returns the number of local matrix rows.
Definition Ifpack_ILU.h:358
int InitAllValues(const Epetra_RowMatrix &A, int MaxNumEntries)
bool Allocated_
Definition Ifpack_ILU.h:411
int NumMyCols() const
Returns the number of local matrix columns.
Definition Ifpack_ILU.h:361
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
Definition Ifpack_ILU.h:266
bool UseTranspose() const
Returns the current UseTranspose setting.
Definition Ifpack_ILU.h:215
int NumCompute_
Contains the number of successful call to Compute().
Definition Ifpack_ILU.h:433
long long NumGlobalNonzeros64() const
Definition Ifpack_ILU.h:353
double Condest() const
Returns the computed estimated condition number, or -1.0 if not computed.
Definition Ifpack_ILU.h:181
const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
Definition Ifpack_ILU.h:218
long long IndexBase64() const
Definition Ifpack_ILU.h:376
bool HasNormInf() const
Returns false because this class cannot compute an Inf-norm.
Definition Ifpack_ILU.h:212
Teuchos::RefCountPtr< Epetra_Map > IlukRowMap_
Definition Ifpack_ILU.h:394
Teuchos::RefCountPtr< Epetra_CrsGraph > CrsGraph_
Definition Ifpack_ILU.h:393
virtual int NumInitialize() const
Returns the number of calls to Initialize().
Definition Ifpack_ILU.h:236
int LevelOfFill_
Level of fill.
Definition Ifpack_ILU.h:423
double Rthresh_
relative threshold
Definition Ifpack_ILU.h:419
Teuchos::RefCountPtr< Epetra_Vector > D_
Diagonal of factors.
Definition Ifpack_ILU.h:407
void Destroy()
Destroys all internal data.
int SetUseTranspose(bool UseTranspose_in)
If set true, transpose of this operator will be applied.
Definition Ifpack_ILU.h:144
const Ifpack_IlukGraph & Graph() const
Returns the address of the Ifpack_IlukGraph associated with this factored matrix.
Definition Ifpack_ILU.h:379
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.
int Initialize()
Initialize the preconditioner, does not touch matrix values.
double Condest_
condition number estimate
Definition Ifpack_ILU.h:421
virtual int NumCompute() const
Returns the number of calls to Compute().
Definition Ifpack_ILU.h:242
double RelaxValue() const
Get ILU(k) relaxation parameter.
Definition Ifpack_ILU.h:327
int NumGlobalCols() const
Returns the number of global matrix columns.
Definition Ifpack_ILU.h:340
int Multiply(bool Trans, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
char Label_[160]
Label of this object.
Definition Ifpack_ILU.h:429
double NormInf() const
Returns 0.0 because this class cannot compute Inf-norm.
Definition Ifpack_ILU.h:209
Epetra_Time Time_
Used for timing issues.
Definition Ifpack_ILU.h:447
~Ifpack_ILU()
Destructor.
Definition Ifpack_ILU.h:91
virtual long long NumGlobalBlockDiagonals64() const
Definition Ifpack_ILU.h:355
bool IsComputed() const
If factor is completed, this query returns true, otherwise it returns false.
Definition Ifpack_ILU.h:120
double InitializeTime_
Contains the time for all successful calls to Initialize().
Definition Ifpack_ILU.h:437
int LevelOfFill() const
Returns the level of fill.
Definition Ifpack_ILU.h:324
Teuchos::RefCountPtr< Epetra_CrsGraph > U_Graph_
Definition Ifpack_ILU.h:405
Ifpack_ILU(Epetra_RowMatrix *A)
Constructor.
const char * Label() const
Returns a character string describing the operator.
Definition Ifpack_ILU.h:199
virtual int NumGlobalBlockDiagonals() const
Returns the number of diagonal entries found in the global input graph.
Definition Ifpack_ILU.h:346
Teuchos::RefCountPtr< Epetra_CrsGraph > L_Graph_
Definition Ifpack_ILU.h:404
int NumMyNonzeros() const
Returns the number of nonzero entries in the local graph.
Definition Ifpack_ILU.h:364
double ApplyInverseTime_
Contains the time for all successful calls to ApplyInverse().
Definition Ifpack_ILU.h:441
int NumInitialize_
Contains the number of successful calls to Initialize().
Definition Ifpack_ILU.h:431
bool IsInitialized_
If true, the preconditioner has been successfully initialized.
Definition Ifpack_ILU.h:425
int NumMyDiagonals_
Definition Ifpack_ILU.h:410
double ApplyInverseFlops_
Contain sthe number of flops for ApplyInverse().
Definition Ifpack_ILU.h:445
Ifpack_IlukGraph: A class for constructing level filled graphs for use with ILU(k) class precondition...
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.