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 
142 
145 
147 
154  virtual inline int SetUseTranspose(bool UseTranspose_in)
155  {
156  UseTranspose_ = UseTranspose_in;
157  return(0);
158  }
159 
161 
163 
165 
173  virtual inline int Apply(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const
174  {
175  if (IsComputed() == false)
176  IFPACK_CHK_ERR(-3);
177 
178  if (X.NumVectors() != Y.NumVectors())
179  IFPACK_CHK_ERR(-2);
180 
181  IFPACK_CHK_ERR(Matrix_->Multiply(UseTranspose(),X,Y));
182  return(0);
183  }
184 
186 
196  virtual int ApplyInverse(const Epetra_MultiVector& X, Epetra_MultiVector& Y) const;
197 
199  virtual double NormInf() const
200  {
201  return(-1.0);
202  }
204 
206 
207  virtual const char * Label() const
208  {
209  return(Label_.c_str());
210  }
211 
213  virtual bool UseTranspose() const
214  {
215  return(UseTranspose_);
216  }
217 
219  virtual bool HasNormInf() const
220  {
221  return(false);
222  }
223 
225  virtual const Epetra_Comm & Comm() const;
226 
228  virtual const Epetra_Map & OperatorDomainMap() const;
229 
231  virtual const Epetra_Map & OperatorRangeMap() const;
232 
233  virtual int Initialize();
234 
235  virtual bool IsInitialized() const
236  {
237  return(IsInitialized_);
238  }
239 
241  virtual inline bool IsComputed() const
242  {
243  return(IsComputed_);
244  }
245 
247  virtual int Compute();
248 
250 
252 
253  virtual const Epetra_RowMatrix& Matrix() const
254  {
255  return(*Matrix_);
256  }
257 
259  virtual double Condest(const Ifpack_CondestType CT = Ifpack_Cheap,
260  const int MaxIters = 1550,
261  const double Tol = 1e-9,
262  Epetra_RowMatrix* Matrix = 0);
263 
265  virtual double Condest() const
266  {
267  return(Condest_);
268  }
269 
271  virtual int SetParameters(Teuchos::ParameterList& List);
272 
274  virtual std::ostream& Print(std::ostream & os) const;
275 
277 
279 
281  virtual int NumInitialize() const
282  {
283  return(NumInitialize_);
284  }
285 
287  virtual int NumCompute() const
288  {
289  return(NumCompute_);
290  }
291 
293  virtual int NumApplyInverse() const
294  {
295  return(NumApplyInverse_);
296  }
297 
299  virtual double InitializeTime() const
300  {
301  return(InitializeTime_);
302  }
303 
305  virtual double ComputeTime() const
306  {
307  return(ComputeTime_);
308  }
309 
311  virtual double ApplyInverseTime() const
312  {
313  return(ApplyInverseTime_);
314  }
315 
317  virtual double InitializeFlops() const
318  {
319  return(0.0);
320  }
321 
323  virtual double ComputeFlops() const
324  {
325  return(ComputeFlops_);
326  }
327 
329  virtual double ApplyInverseFlops() const
330  {
331  return(ApplyInverseFlops_);
332  }
333 
334  // @}
335 
336 private:
337 
338  // @{ Application of the preconditioner
339 
341  virtual int ApplyInverseJacobi(const Epetra_MultiVector& X,
342  Epetra_MultiVector& Y) const;
343 
345  virtual int ApplyInverseGS(const Epetra_MultiVector& X,
346  Epetra_MultiVector& Y) const;
347 
348  virtual int ApplyInverseGS_RowMatrix(const Epetra_MultiVector& X,
349  Epetra_MultiVector& Y) const;
350 
351  virtual int ApplyInverseGS_CrsMatrix(const Epetra_CrsMatrix* A,
352  const Epetra_MultiVector& X,
353  Epetra_MultiVector& Y) const;
354 
355  virtual int ApplyInverseGS_FastCrsMatrix(const Epetra_CrsMatrix* A,
356  const Epetra_MultiVector& X,
357  Epetra_MultiVector& Y) const;
358  virtual int ApplyInverseGS_LocalFastCrsMatrix(const Epetra_CrsMatrix* A,
359  const Epetra_MultiVector& X,
360  Epetra_MultiVector& Y) const;
361 
363  virtual int ApplyInverseSGS(const Epetra_MultiVector& X,
364  Epetra_MultiVector& Y) const;
365 
366  virtual int ApplyInverseSGS_RowMatrix(const Epetra_MultiVector& X,
367  Epetra_MultiVector& Y) const;
368 
369  virtual int ApplyInverseSGS_CrsMatrix(const Epetra_CrsMatrix* A,
370  const Epetra_MultiVector& X,
371  Epetra_MultiVector& Y) const;
372 
373  virtual int ApplyInverseSGS_FastCrsMatrix(const Epetra_CrsMatrix* A,
374  const Epetra_MultiVector& X,
375  Epetra_MultiVector& Y) const;
376 
377  virtual int ApplyInverseSGS_LocalFastCrsMatrix(const Epetra_CrsMatrix* A,
378  const Epetra_MultiVector& X,
379  Epetra_MultiVector& Y) const;
380 
381 
383 
384 private:
385 
387  virtual void SetLabel();
388 
391  {}
392 
394  Ifpack_PointRelaxation& operator=(const Ifpack_PointRelaxation& /* rhs */)
395  {
396  return(*this);
397  }
398 
399  // @{ Initializations, timing and flops
401  bool IsInitialized_;
403  bool IsComputed_;
405  int NumInitialize_;
407  int NumCompute_;
409  mutable int NumApplyInverse_;
411  double InitializeTime_;
413  double ComputeTime_;
415  mutable double ApplyInverseTime_;
417  double ComputeFlops_;
419  mutable double ApplyInverseFlops_;
420  // @}
421 
422  // @{ Settings
424  int NumSweeps_;
426  double DampingFactor_;
428  bool UseTranspose_;
430  double Condest_;
431 #if 0
432  // Unused; commented out to avoid build warnings
433 
435  bool ComputeCondest_;
436 #endif // 0
437  std::string Label_;
439  int PrecType_;
440  double MinDiagonalValue_;
441  // @}
442 
443  // @{ Other data
445  int NumMyRows_;
447  int NumMyNonzeros_;
449  long long NumGlobalRows_;
451  long long NumGlobalNonzeros_;
453  Teuchos::RefCountPtr<const Epetra_RowMatrix> Matrix_;
455  Teuchos::RefCountPtr<Epetra_Import> Importer_;
457  mutable Teuchos::RefCountPtr<Epetra_Vector> Diagonal_;
459  Teuchos::RefCountPtr<Epetra_Time> Time_;
461  bool IsParallel_;
463  bool ZeroStartingSolution_;
465  bool DoBackwardGS_;
467  bool DoL1Method_;
469  double L1Eta_;
470 
472  int NumLocalSmoothingIndices_;
474  int * LocalSmoothingIndices_;
475 
476  // @}
477 
478 
479 
480 };
481 
482 #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.