ML_Epetra::EdgeMatrixFreePreconditioner Class Reference

#include <ml_EdgeMatrixFreePreconditioner.h>

Collaboration diagram for ML_Epetra::EdgeMatrixFreePreconditioner:

Collaboration graph
[legend]

List of all members.

Public Member Functions

 EdgeMatrixFreePreconditioner (const Epetra_Operator_With_MatMat &Operator, const Epetra_Vector &Diagonal, const Epetra_CrsMatrix &D0_Matrix, const Epetra_CrsMatrix &D0_Clean_Matrix, const Epetra_CrsMatrix &TMT_Matrix, const int *BCedges, const int numBCedges, const Teuchos::ParameterList &List, const bool ComputePrec=true)
 Constructs an EdgeMatrixFreePreconditioner.
 ~EdgeMatrixFreePreconditioner ()
 Destructor.
int ComputePreconditioner (const bool CheckFiltering=false)
 Computes the multilevel hierarchy.
int ReComputePreconditioner ()
 Recomputes the preconditioner.
int Apply (const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
 Apply the inverse of the preconditioner to an Epetra_MultiVector (NOT AVAILABLE).
int ApplyInverse (const Epetra_MultiVector &B, Epetra_MultiVector &X) const
 Apply the preconditioner to RHS B to get result X (X is also initial guess).
void Print (const char *whichHierarchy="main")
 Print the individual operators in the multigrid hierarchy.
int DestroyPreconditioner ()
 Destroys all structures allocated in ComputePreconditioner() if the preconditioner has been computed.
int SetUseTranspose (bool UseTranspose)
 Sets use transpose (not implemented).
double NormInf () const
 Returns the infinity norm (not implemented).
bool UseTranspose () const
 Returns the current UseTranspose setting.
bool HasNormInf () const
 Returns true if the this object can provide an approximate Inf-norm, false otherwise.
const Epetra_Comm & Comm () const
 Returns a pointer to the Epetra_Comm communicator associated with this operator.
const Epetra_Map & OperatorDomainMap () const
 Returns the Epetra_Map object associated with the domain of this operator.
const Epetra_Map & OperatorRangeMap () const
 Returns the Epetra_Map object associated with the range of this operator.

Private Member Functions

int SetupSmoother ()
 Sets up the Chebyshev smoother.
Epetra_MultiVector * BuildNullspace ()
 Build the edge nullspace.
int BuildProlongator (const Epetra_MultiVector &nullspace)
 Build the edge-to-vector-node prolongator described in Bochev, Hu, Siefert and Tuminaro (2006).
int FormCoarseMatrix ()
 Forms the coarse matrix, given the build prolongator.

Private Attributes

int dim
 Dimension of space.
ML_Comm * ml_comm_
 ML Communicator.
const Epetra_Operator_With_MatMat * Operator_
 Fine-level operator.
const Epetra_CrsMatrix * D0_Matrix_
 D0 or T matrix w/ dirichlet nodes and edges zero'd. Used to generate prolongator.
const Epetra_CrsMatrix * D0_Clean_Matrix_
 D0 or T matrix w/ nothing zero'd. Needed to generate the nullspace.
const Epetra_CrsMatrix * TMT_Matrix_
 TMT_Matrix. Needed for nodal maps.
const int * BCedges_
 Dirichlet edges.
const int numBCedges_
Epetra_CrsMatrix * Prolongator_
 Prolongator.
EpetraExt::CrsMatrix_SolverMap ProlongatorColMapTrans_
Epetra_Vector * InvDiagonal_
 Inverse Diagonal.
Epetra_CrsMatrix * CoarseMatrix
 Coarse Matrix.
ML_OperatorCoarseMat_ML
MultiLevelPreconditionerCoarsePC
 Level 2+ Preconditioner.
Ifpack_Chebyshev * Smoother_
 Ifpack Chebyshev Smoother.
const Epetra_Map * EdgeDomainMap_
 Edge Domain Map.
const Epetra_Map * EdgeRangeMap_
 Edge Range Map.
const Epetra_Comm * Comm_
 Epetra communicator object.
const Epetra_Map * NodeDomainMap_
 Nodal Domain Map.
const Epetra_Map * NodeRangeMap_
 Nodal Range Map.
Epetra_Map * CoarseMap_
 Coarse Domain/Range Map.
int num_cycles
 Number of V-cycles to run.
int MaxLevels
bool verbose_
bool very_verbose_
bool print_hierarchy


Detailed Description

Matrix-Free preconditioning class for edge Maxwell Problems.

Constructor & Destructor Documentation

ML_Epetra::EdgeMatrixFreePreconditioner::EdgeMatrixFreePreconditioner ( const Epetra_Operator_With_MatMat &  Operator,
const Epetra_Vector &  Diagonal,
const Epetra_CrsMatrix &  D0_Matrix,
const Epetra_CrsMatrix &  D0_Clean_Matrix,
const Epetra_CrsMatrix &  TMT_Matrix,
const int *  BCedges,
const int  numBCedges,
const Teuchos::ParameterList &  List,
const bool  ComputePrec = true 
)

ML_Epetra::EdgeMatrixFreePreconditioner::~EdgeMatrixFreePreconditioner (  ) 


Member Function Documentation

int ML_Epetra::EdgeMatrixFreePreconditioner::Apply ( const Epetra_MultiVector &  X,
Epetra_MultiVector &  Y 
) const [inline]

int ML_Epetra::EdgeMatrixFreePreconditioner::ApplyInverse ( const Epetra_MultiVector &  B,
Epetra_MultiVector &  X 
) const

Epetra_MultiVector* ML_Epetra::EdgeMatrixFreePreconditioner::BuildNullspace (  )  [private]

int ML_Epetra::EdgeMatrixFreePreconditioner::BuildProlongator ( const Epetra_MultiVector &  nullspace  )  [private]

const Epetra_Comm& ML_Epetra::EdgeMatrixFreePreconditioner::Comm (  )  const [inline]

int ML_Epetra::EdgeMatrixFreePreconditioner::ComputePreconditioner ( const bool  CheckFiltering = false  ) 

Computes the multilevel hierarchy. This function retrives the user's defines parameters (as specified in the input ParameterList), or takes default values otherwise, and creates the ML objects for aggregation and hierarchy. Allocated data can be freed used DestroyPreconditioner(), or by the destructor,

In a Newton-type procedure, several linear systems have to be solved, Often, these systems are not too different. In this case, it might be convenient to keep the already computed preconditioner (with hierarchy, coarse solver, smoothers), and use it to precondition the next linear system. ML offers a way to determine whether the already available preconditioner is "good enough" for the next linear system. The user should proceed as follows:

  • define "adaptive: enable" == true
  • solve the first linear system. ML tries to estimate the rate of convergence, and record it;
  • change the values of the linear system matrix (but NOT its structure)
  • compute the new preconditioner as ComputePreconditioner(true) It is supposed that the pointer to the Epetra_RowMatrix remains constant. Currently, it is not possible to modify this pointer (other than creating a new preconditioner) Computes the preconditioner

int ML_Epetra::EdgeMatrixFreePreconditioner::DestroyPreconditioner (  ) 

int ML_Epetra::EdgeMatrixFreePreconditioner::FormCoarseMatrix (  )  [private]

bool ML_Epetra::EdgeMatrixFreePreconditioner::HasNormInf (  )  const [inline]

double ML_Epetra::EdgeMatrixFreePreconditioner::NormInf (  )  const [inline]

const Epetra_Map& ML_Epetra::EdgeMatrixFreePreconditioner::OperatorDomainMap (  )  const [inline]

const Epetra_Map& ML_Epetra::EdgeMatrixFreePreconditioner::OperatorRangeMap (  )  const [inline]

void ML_Epetra::EdgeMatrixFreePreconditioner::Print ( const char *  whichHierarchy = "main"  ) 

int ML_Epetra::EdgeMatrixFreePreconditioner::ReComputePreconditioner (  )  [inline]

int ML_Epetra::EdgeMatrixFreePreconditioner::SetupSmoother (  )  [private]

int ML_Epetra::EdgeMatrixFreePreconditioner::SetUseTranspose ( bool  UseTranspose  )  [inline]

bool ML_Epetra::EdgeMatrixFreePreconditioner::UseTranspose (  )  const [inline]


Member Data Documentation

const Epetra_Comm* ML_Epetra::EdgeMatrixFreePreconditioner::Comm_ [private]

const Epetra_CrsMatrix* ML_Epetra::EdgeMatrixFreePreconditioner::D0_Matrix_ [private]

const Epetra_Operator_With_MatMat* ML_Epetra::EdgeMatrixFreePreconditioner::Operator_ [private]

EpetraExt::CrsMatrix_SolverMap ML_Epetra::EdgeMatrixFreePreconditioner::ProlongatorColMapTrans_ [private]

const Epetra_CrsMatrix* ML_Epetra::EdgeMatrixFreePreconditioner::TMT_Matrix_ [private]


The documentation for this class was generated from the following file: