Dune::Amg::KAMG< M, X, S, PI, K, A > Class Template Reference
[Parallel Algebraic Multigrid]

an algebraic multigrid method using a Krylov-cycle. More...

#include <dune/istl/paamg/kamg.hh>

Inheritance diagram for Dune::Amg::KAMG< M, X, S, PI, K, A >:
Inheritance graph

List of all members.

Public Types

enum  { category = Amg::category }
typedef AMG< M, X, S, PI, AAmg
 The type of the underlying AMG.
typedef K KrylovSolver
 The type of the Krylov solver for the cycle.
typedef Amg::OperatorHierarchy OperatorHierarchy
 The type of the hierarchy of operators.
typedef Amg::CoarseSolver CoarseSolver
 The type of the coarse solver.
typedef Amg::ParallelInformation ParallelInformation
 the type of the parallelinformation to use.
typedef Amg::SmootherArgs SmootherArgs
 The type of the arguments for construction of the smoothers.
typedef Amg::Operator Operator
 the type of the lineatr operator.
typedef Amg::Domain Domain
 the type of the domain.
typedef Amg::Range Range
 The type of the range.
typedef
Amg::ParallelInformationHierarchy 
ParallelInformationHierarchy
 The type of the hierarchy of parallel information.
typedef Amg::ScalarProduct ScalarProduct
 The type of the scalar product.
typedef X domain_type
 The domain type of the preconditioner.
typedef X range_type
 The range type of the preconditioner.
typedef X::field_type field_type
 The field type of the preconditioner.

Public Member Functions

 KAMG (const OperatorHierarchy &matrices, CoarseSolver &coarseSolver, const SmootherArgs &smootherArgs, std::size_t gamma, std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1, std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1)
 Construct a new amg with a specific coarse solver.
template<class C >
 KAMG (const Operator &fineOperator, const C &criterion, const SmootherArgs &smootherArgs, std::size_t gamma=1, std::size_t preSmoothingSteps=1, std::size_t postSmoothingSteps=1, std::size_t maxLevelKrylovSteps=3, double minDefectReduction=1e-1, const ParallelInformation &pinfo=ParallelInformation())
 Construct an AMG with an inexact coarse solver based on the smoother.
void pre (Domain &x, Range &b)
void post (Domain &x)
void apply (Domain &x, const Range &b)
std::size_t maxlevels ()
virtual void pre (X &x, X &b)=0
 Prepare the preconditioner.
virtual void apply (X &v, const X &d)=0
 Apply one step of the preconditioner to the system A(v)=d.
virtual void post (X &x)=0
 Clean up.

Detailed Description

template<class M, class X, class S, class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
class Dune::Amg::KAMG< M, X, S, PI, K, A >

an algebraic multigrid method using a Krylov-cycle.

Template Parameters:
M The type of the linear operator.
X The type of the range and domain.
PI The parallel information object. Use SequentialInformation (default) for a sequential AMG, OwnerOverlapCopyCommunication for the parallel case.
K The type of the Krylov method to use for the cycle.
A The type of the allocator to use.

Member Typedef Documentation

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef AMG<M,X,S,PI,A> Dune::Amg::KAMG< M, X, S, PI, K, A >::Amg

The type of the underlying AMG.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::CoarseSolver Dune::Amg::KAMG< M, X, S, PI, K, A >::CoarseSolver

The type of the coarse solver.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::Domain Dune::Amg::KAMG< M, X, S, PI, K, A >::Domain

the type of the domain.

typedef X Dune::Preconditioner< X, X >::domain_type [inherited]

The domain type of the preconditioner.

Reimplemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.

typedef X::field_type Dune::Preconditioner< X, X >::field_type [inherited]

The field type of the preconditioner.

Reimplemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef K Dune::Amg::KAMG< M, X, S, PI, K, A >::KrylovSolver

The type of the Krylov solver for the cycle.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::Operator Dune::Amg::KAMG< M, X, S, PI, K, A >::Operator

the type of the lineatr operator.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::OperatorHierarchy Dune::Amg::KAMG< M, X, S, PI, K, A >::OperatorHierarchy

The type of the hierarchy of operators.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::ParallelInformation Dune::Amg::KAMG< M, X, S, PI, K, A >::ParallelInformation

the type of the parallelinformation to use.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::ParallelInformationHierarchy Dune::Amg::KAMG< M, X, S, PI, K, A >::ParallelInformationHierarchy

The type of the hierarchy of parallel information.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::Range Dune::Amg::KAMG< M, X, S, PI, K, A >::Range

The type of the range.

typedef X Dune::Preconditioner< X, X >::range_type [inherited]

The range type of the preconditioner.

Reimplemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::ScalarProduct Dune::Amg::KAMG< M, X, S, PI, K, A >::ScalarProduct

The type of the scalar product.

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
typedef Amg::SmootherArgs Dune::Amg::KAMG< M, X, S, PI, K, A >::SmootherArgs

The type of the arguments for construction of the smoothers.


Member Enumeration Documentation

template<class M , class X , class S , class PI = SequentialInformation, class K = BiCGSTABSolver<X>, class A = std::allocator<X>>
anonymous enum
Enumerator:
category 

The solver category.


Member Function Documentation

virtual void Dune::Preconditioner< X, X >::apply ( X &  v,
const X &  d 
) [pure virtual, inherited]

Apply one step of the preconditioner to the system A(v)=d.

On entry v=0 and d=b-A(x) (although this might not be computed in that way. On exit v contains the update, i.e one step computes $ v = M^{-1} d $ where $ M $ is the approximate inverse of the operator $ A $ characterizing the preconditioner.

Parameters:
[out] v The update to be computed
d The current defect.

Implemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.

virtual void Dune::Preconditioner< X, X >::post ( X &  x  )  [pure virtual, inherited]

Clean up.

This method is called after the last apply call for the linear system to be solved. Memory may be deallocated safely here. x is the solution of the linear equation.

Parameters:
x The right hand side of the equation.

Implemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.

virtual void Dune::Preconditioner< X, X >::pre ( X &  x,
X &  b 
) [pure virtual, inherited]

Prepare the preconditioner.

A solver solves a linear operator equation A(x)=b by applying one or several steps of the preconditioner. The method pre() is called before the first apply operation. b and x are right hand side and solution vector of the linear system respectively. It may. e.g., scale the system, allocate memory or compute a (I)LU decomposition. Note: The ILU decomposition could also be computed in the constructor or with a separate method of the derived method if several linear systems with the same matrix are to be solved.

Parameters:
x The left hand side of the equation.
b The right hand side of the equation.

Implemented in Dune::SeqOverlappingSchwarz< M, X, TM, TD, TA >, and Dune::SeqOverlappingSchwarz< M, X, TM, TS, TA >.


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

Generated on Fri Apr 29 2011 with Doxygen (ver 1.7.1) [doxygen-log,error-log].