IFPACK  Development
 All Classes Namespaces Files Functions Variables Enumerations Friends Pages
Ifpack_PointRelaxation.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_POINTRELAXATION_H
44 #define IFPACK_POINTRELAXATION_H
45 
46 #include "Ifpack_ConfigDefs.h"
47 #include "Ifpack_Preconditioner.h"
48 
49 #include "Epetra_Vector.h"
50 #include "Epetra_Time.h"
51 #include "Epetra_RowMatrix.h"
52 #include "Epetra_Import.h"
53 
54 #include "Teuchos_RefCountPtr.hpp"
55 
56 namespace Teuchos {
57  class ParameterList;
58 }
59 class Epetra_MultiVector;
60 class Epetra_Vector;
61 class Epetra_Map;
62 class Epetra_Comm;
63 class Epetra_CrsMatrix;
64 
66 
131 
132 public:
133 
135 
141 
144 
146 
153  virtual inline int SetUseTranspose(bool UseTranspose_in)
154  {
155  UseTranspose_ = UseTranspose_in;
156  return(0);
157  }
158 
160 
162 
164 
172  virtual inline int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const
173  {
174  if (IsComputed() == false)
175  IFPACK_CHK_ERR(-3);
176 
177  if (X.NumVectors() != Y.NumVectors())
178  IFPACK_CHK_ERR(-2);
179 
180  IFPACK_CHK_ERR(Matrix_->Multiply(UseTranspose(),X,Y));
181  return(0);
182  }
183 
185 
195  virtual int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
196 
198  virtual double NormInf() const
199  {
200  return(-1.0);
201  }
203 
205 
206  virtual const char * Label() const
207  {
208  return(Label_.c_str());
209  }
210 
212  virtual bool UseTranspose() const
213  {
214  return(UseTranspose_);
215  }
216 
218  virtual bool HasNormInf() const
219  {
220  return(false);
221  }
222 
224  virtual const Epetra_Comm & Comm() const;
225 
227  virtual const Epetra_Map & OperatorDomainMap() const;
228 
230  virtual const Epetra_Map & OperatorRangeMap() const;
231 
232  virtual int Initialize();
233 
234  virtual bool IsInitialized() const
235  {
236  return(IsInitialized_);
237  }
238 
240  virtual inline bool IsComputed() const
241  {
242  return(IsComputed_);
243  }
244 
246  virtual int Compute();
247 
249 
251 
252  virtual const Epetra_RowMatrix& Matrix() const
253  {
254  return(*Matrix_);
255  }
256 
258  virtual double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
259  const int MaxIters = 1550,
260  const double Tol = 1e-9,
261  Epetra_RowMatrix* Matrix = 0);
262 
264  virtual double Condest() const
265  {
266  return(Condest_);
267  }
268 
270  virtual int SetParameters(Teuchos::ParameterList& List);
271 
273  virtual std::ostream& Print(std::ostream & os) const;
274 
276 
278 
280  virtual int NumInitialize() const
281  {
282  return(NumInitialize_);
283  }
284 
286  virtual int NumCompute() const
287  {
288  return(NumCompute_);
289  }
290 
292  virtual int NumApplyInverse() const
293  {
294  return(NumApplyInverse_);
295  }
296 
298  virtual double InitializeTime() const
299  {
300  return(InitializeTime_);
301  }
302 
304  virtual double ComputeTime() const
305  {
306  return(ComputeTime_);
307  }
308 
310  virtual double ApplyInverseTime() const
311  {
312  return(ApplyInverseTime_);
313  }
314 
316  virtual double InitializeFlops() const
317  {
318  return(0.0);
319  }
320 
322  virtual double ComputeFlops() const
323  {
324  return(ComputeFlops_);
325  }
326 
328  virtual double ApplyInverseFlops() const
329  {
330  return(ApplyInverseFlops_);
331  }
332 
333  // @}
334 
335 private:
336 
337  // @{ Application of the preconditioner
338 
340  virtual int ApplyInverseJacobi(const Epetra_MultiVector& X,
341  Epetra_MultiVector& Y) const;
342 
344  virtual int ApplyInverseGS(const Epetra_MultiVector& X,
345  Epetra_MultiVector& Y) const;
346 
347  virtual int ApplyInverseGS_RowMatrix(const Epetra_MultiVector& X,
348  Epetra_MultiVector& Y) const;
349 
350  virtual int ApplyInverseGS_CrsMatrix(const Epetra_CrsMatrix* A,
351  const Epetra_MultiVector& X,
352  Epetra_MultiVector& Y) const;
353 
354  virtual int ApplyInverseGS_FastCrsMatrix(const Epetra_CrsMatrix* A,
355  const Epetra_MultiVector& X,
356  Epetra_MultiVector& Y) const;
357  virtual int ApplyInverseGS_LocalFastCrsMatrix(const Epetra_CrsMatrix* A,
358  const Epetra_MultiVector& X,
359  Epetra_MultiVector& Y) const;
360 
362  virtual int ApplyInverseSGS(const Epetra_MultiVector& X,
363  Epetra_MultiVector& Y) const;
364 
365  virtual int ApplyInverseSGS_RowMatrix(const Epetra_MultiVector& X,
366  Epetra_MultiVector& Y) const;
367 
368  virtual int ApplyInverseSGS_CrsMatrix(const Epetra_CrsMatrix* A,
369  const Epetra_MultiVector& X,
370  Epetra_MultiVector& Y) const;
371 
372  virtual int ApplyInverseSGS_FastCrsMatrix(const Epetra_CrsMatrix* A,
373  const Epetra_MultiVector& X,
374  Epetra_MultiVector& Y) const;
375 
376  virtual int ApplyInverseSGS_LocalFastCrsMatrix(const Epetra_CrsMatrix* A,
377  const Epetra_MultiVector& X,
378  Epetra_MultiVector& Y) const;
379 
380 
382 
383 private:
384 
386  virtual void SetLabel();
387 
390  {}
391 
393  Ifpack_PointRelaxation& operator=(const Ifpack_PointRelaxation& /* rhs */)
394  {
395  return(*this);
396  }
397 
398  // @{ Initializations, timing and flops
400  bool IsInitialized_;
402  bool IsComputed_;
404  int NumInitialize_;
406  int NumCompute_;
408  mutable int NumApplyInverse_;
410  double InitializeTime_;
412  double ComputeTime_;
414  mutable double ApplyInverseTime_;
416  double ComputeFlops_;
418  mutable double ApplyInverseFlops_;
419  // @}
420 
421  // @{ Settings
423  int NumSweeps_;
425  double DampingFactor_;
427  bool UseTranspose_;
429  double Condest_;
430 #if 0
431  // Unused; commented out to avoid build warnings
432 
434  bool ComputeCondest_;
435 #endif // 0
436  std::string Label_;
438  int PrecType_;
439  double MinDiagonalValue_;
440  // @}
441 
442  // @{ Other data
444  int NumMyRows_;
446  int NumMyNonzeros_;
448  long long NumGlobalRows_;
450  long long NumGlobalNonzeros_;
452  Teuchos::RefCountPtr<const Epetra_RowMatrix> Matrix_;
454  Teuchos::RefCountPtr<Epetra_Import> Importer_;
456  mutable Teuchos::RefCountPtr<Epetra_Vector> Diagonal_;
458  Teuchos::RefCountPtr<Epetra_Time> Time_;
460  bool IsParallel_;
462  bool ZeroStartingSolution_;
464  bool DoBackwardGS_;
466  bool DoL1Method_;
468  double L1Eta_;
469 
471  int NumLocalSmoothingIndices_;
473  int * LocalSmoothingIndices_;
474 
475  // @}
476 
477 
478 
479 };
480 
481 #endif // IFPACK_POINTRELAXATION_H
virtual std::ostream & Print(std::ostream &os) const
Prints object to an output stream.
virtual int Apply(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the matrix to an Epetra_MultiVector.
virtual const Epetra_Comm & Comm() const
Returns a pointer to the Epetra_Comm communicator associated with this operator.
virtual int NumApplyInverse() const
Returns the number of calls to ApplyInverse().
virtual const Epetra_Map & OperatorDomainMap() const
Returns the Epetra_Map object associated with the domain of this operator.
virtual double ComputeFlops() const
Returns the number of flops in the computation phase.
virtual int SetUseTranspose(bool UseTranspose_in)
virtual double ApplyInverseFlops() const
Returns the number of flops for the application of the preconditioner.
Ifpack_PointRelaxation(const Epetra_RowMatrix *Matrix)
Ifpack_PointRelaxation constructor with given Epetra_RowMatrix.
virtual int NumCompute() const
Returns the number of calls to Compute().
virtual double Condest() const
Returns the condition number estimate, or -1.0 if not computed.
virtual int SetParameters(Teuchos::ParameterList &List)
Sets all the parameters for the preconditioner.
virtual double ComputeTime() const
Returns the time spent in Compute().
virtual double NormInf() const
Returns the infinity norm of the global matrix (not implemented)
virtual bool UseTranspose() const
Returns the current UseTranspose setting.
Ifpack_Preconditioner: basic class for preconditioning in Ifpack.
virtual bool IsInitialized() const
Returns true if the preconditioner has been successfully initialized, false otherwise.
Ifpack_PointRelaxation: a class to define point relaxation preconditioners of for Epetra_RowMatrix&#39;s...
virtual int Compute()
Computes the preconditioners.
virtual int Initialize()
Computes all it is necessary to initialize the preconditioner.
virtual const Epetra_RowMatrix & Matrix() const
Returns a pointer to the matrix to be preconditioned.
virtual bool IsComputed() const
Returns true if the preconditioner has been successfully computed.
virtual int NumInitialize() const
Returns the number of calls to Initialize().
virtual double InitializeFlops() const
Returns the number of flops in the initialization phase.
virtual bool HasNormInf() const
Returns true if the this object can provide an approximate Inf-norm, false otherwise.
virtual const Epetra_Map & OperatorRangeMap() const
Returns the Epetra_Map object associated with the range of this operator.
virtual double InitializeTime() const
Returns the time spent in Initialize().
virtual double ApplyInverseTime() const
Returns the time spent in ApplyInverse().
virtual int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Applies the preconditioner to X, returns the result in Y.
virtual ~Ifpack_PointRelaxation()
Destructor.