MultilevelLinearOp< T > Class Template Reference

#include <MultilevelLinearOp.H>

Inheritance diagram for MultilevelLinearOp< T >:

Inheritance graph
[legend]
Collaboration diagram for MultilevelLinearOp< T >:

Collaboration graph
[legend]

List of all members.


Detailed Description

template<class T>
class MultilevelLinearOp< T >

Operator class for Linear solver for solving L(phi) = rhs for the case where <T> is a multilevel (AMR) data structure; as such, T is actually a Vector<T>. The MultilevelLinearOp is also defined with a Vector<AMRLevelOp> which applies the operator in an AMR way. For now, this is hardwired for T = LevelData<FArrayBox>, but we can hopefully get back redo this once I get the kinks worked out...

Public Member Functions

 MultilevelLinearOp ()
 default constructor (sets defaults for preconditioner)
void define (const Vector< DisjointBoxLayout > &a_vectGrids, const Vector< int > &a_refRatios, const Vector< ProblemDomain > &a_domains, const Vector< RealVect > &a_vectDx, RefCountedPtr< AMRLevelOpFactory< LevelData< T > > > &a_opFactory, int a_lBase)
 define function -- opFactory should be able to define all required levels
virtual ~MultilevelLinearOp ()
virtual void residual (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_phi, const Vector< LevelData< T > * > &a_rhs, bool a_homogeneous=false)
 compute residual using AMRLevelOps AMRResidual functionality
virtual void preCond (Vector< LevelData< T > * > &a_cor, const Vector< LevelData< T > * > &a_residual)
 Apply preconditioner.
virtual void applyOp (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_phi, bool a_homogeneous=false)
virtual void create (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
virtual void clear (Vector< LevelData< T > * > &a_lhs)
virtual void assign (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_rhs)
virtual Real dotProduct (const Vector< LevelData< T > * > &a_1, const Vector< LevelData< T > * > &a_2)
virtual void incr (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_x, Real a_scale)
virtual void axby (Vector< LevelData< T > * > &a_lhs, const Vector< LevelData< T > * > &a_x, const Vector< LevelData< T > * > &a_y, Real a_a, Real a_b)
virtual void scale (Vector< LevelData< T > * > &a_lhs, const Real &a_scale)
virtual Real norm (const Vector< LevelData< T > * > &a_rhs, int a_ord)
virtual void setToZero (Vector< LevelData< T > * > &a_lhs)
virtual void write (const Vector< LevelData< T > * > *a_data, const char *a_filename)

Public Attributes

bool m_use_multigrid_preconditioner
 if true, preCond(...) function uses multigrid v-cycle(s)
int m_num_mg_iterations
 number of mg v-cycles for each preCond call (if using mg preconditioner)
int m_num_mg_smooth
 number of smoothings in the multigrid preconditioner
int m_preCondSolverDepth
 max depth for multigrid preconditioner (if using mg preconditioner)
Vector< DisjointBoxLayoutm_vectGrids
Vector< int > m_refRatios
Vector< ProblemDomainm_domains
Vector< RealVectm_vectDx
Vector< RefCountedPtr
< AMRLevelOp< LevelData< T > > > > 
m_vectOperators
 vector of level operators
int m_lBase
AMRMultiGrid< LevelData< T > > m_preCondSolver
 Solver which is used by the preCond function.
LinearSolver< LevelData< T > > * m_precondBottomSolverPtr
 bottom solver for m_preCondSolver
bool m_isPrecondSolverInitialized
 has the preconditioner solver been initialized?

Constructor & Destructor Documentation

template<class T>
MultilevelLinearOp< T >::MultilevelLinearOp (  )  [inline]

template<class T>
MultilevelLinearOp< T >::~MultilevelLinearOp (  )  [inline, virtual]


Member Function Documentation

template<class T>
void MultilevelLinearOp< T >::define ( const Vector< DisjointBoxLayout > &  a_vectGrids,
const Vector< int > &  a_refRatios,
const Vector< ProblemDomain > &  a_domains,
const Vector< RealVect > &  a_vectDx,
RefCountedPtr< AMRLevelOpFactory< LevelData< T > > > &  a_opFactory,
int  a_lBase 
) [inline]

template<class T>
void MultilevelLinearOp< T >::residual ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_phi,
const Vector< LevelData< T > * > &  a_rhs,
bool  a_homogeneous = false 
) [inline, virtual]

compute residual using AMRLevelOps AMRResidual functionality

a_lhs = L(a_phi) - a_rhs. If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Say you are solving L(phi) = rhs. Make a_lhs = L(a_phi) - a_rhs. If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, MultilevelLinearOp< T >::m_vectOperators, and Vector< T >::size().

Referenced by MultilevelLinearOp< T >::preCond().

template<class T>
void MultilevelLinearOp< T >::preCond ( Vector< LevelData< T > * > &  a_cor,
const Vector< LevelData< T > * > &  a_residual 
) [inline, virtual]

template<class T>
void MultilevelLinearOp< T >::applyOp ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_phi,
bool  a_homogeneous = false 
) [inline, virtual]

In the context of solving L(phi) = rhs, set a_lhs = L(a_phi). If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

In the context of solving L(phi) = rhs, set a_lhs = L(a_phi). If a_homogeneous is true, evaluate the operator using homogeneous boundary conditions.

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, and MultilevelLinearOp< T >::m_vectOperators.

template<class T>
void MultilevelLinearOp< T >::create ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_rhs 
) [inline, virtual]

Create data holder a_lhs that mirrors a_rhs. You do not need to copy the data of a_rhs, just make a holder the same size.

Create data holder a_lhs that mirrors a_rhs. You do not need to copy the data of a_rhs, just make a holder the same size.

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, MultilevelLinearOp< T >::m_vectOperators, and Max().

Referenced by MultilevelLinearOp< T >::dotProduct(), and MultilevelLinearOp< T >::preCond().

template<class T>
void MultilevelLinearOp< T >::clear ( Vector< LevelData< T > * > &  a_lhs  )  [inline, virtual]

Clean up data holder before it goes out of scope. This is necessary because create calls new.

Clean up data holder before it goes out of scope. This is necessary because create calls new.

Reimplemented from LinearOp< Vector< LevelData< T > * > >.

References CH_TIME.

Referenced by MultilevelLinearOp< T >::dotProduct(), and MultilevelLinearOp< T >::preCond().

template<class T>
void MultilevelLinearOp< T >::assign ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_rhs 
) [inline, virtual]

template<class T>
Real MultilevelLinearOp< T >::dotProduct ( const Vector< LevelData< T > * > &  a_1,
const Vector< LevelData< T > * > &  a_2 
) [inline, virtual]

template<class T>
void MultilevelLinearOp< T >::incr ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_x,
Real  a_scale 
) [inline, virtual]

Increment by scaled amount (a_lhs += a_scale*a_x).

Increment by scaled amount (a_lhs += a_scale*a_x).

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, and MultilevelLinearOp< T >::m_vectOperators.

Referenced by MultilevelLinearOp< T >::preCond().

template<class T>
void MultilevelLinearOp< T >::axby ( Vector< LevelData< T > * > &  a_lhs,
const Vector< LevelData< T > * > &  a_x,
const Vector< LevelData< T > * > &  a_y,
Real  a_a,
Real  a_b 
) [inline, virtual]

Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).

Set input to a scaled sum (a_lhs = a_a*a_x + a_b*a_y).

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, and MultilevelLinearOp< T >::m_vectOperators.

template<class T>
void MultilevelLinearOp< T >::scale ( Vector< LevelData< T > * > &  a_lhs,
const Real a_scale 
) [inline, virtual]

Multiply the input by a given scale (a_lhs *= a_scale).

Multiply the input by a given scale (a_lhs *= a_scale).

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MultilevelLinearOp< T >::m_lBase, and MultilevelLinearOp< T >::m_vectOperators.

Referenced by MultilevelLinearOp< T >::dotProduct().

template<class T>
Real MultilevelLinearOp< T >::norm ( const Vector< LevelData< T > * > &  a_rhs,
int  a_ord 
) [inline, virtual]

Return the norm of a_rhs. a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.

Return the norm of a_rhs. a_ord == 0 max norm, a_ord == 1 sum(abs(a_rhs)), else, L(a_ord) norm.

Implements LinearOp< Vector< LevelData< T > * > >.

References CH_TIME, MayDay::Error(), MultilevelLinearOp< T >::m_lBase, MultilevelLinearOp< T >::m_refRatios, and MultilevelLinearOp< T >::m_vectOperators.

template<class T>
void MultilevelLinearOp< T >::setToZero ( Vector< LevelData< T > * > &  a_lhs  )  [inline, virtual]

template<class T>
void MultilevelLinearOp< T >::write ( const Vector< LevelData< T > * > *  a,
const char *  filename 
) [inline, virtual]

Debugging aid for solvers. Print out a "T" to a file named "filename" default implementation is to print out a message saying "LinearOp::write not implemented"

Reimplemented from LinearOp< Vector< LevelData< T > * > >.

References MultilevelLinearOp< T >::m_refRatios, MayDay::Warning(), and writeVectorLevelName().


Member Data Documentation

template<class T>
bool MultilevelLinearOp< T >::m_use_multigrid_preconditioner

if true, preCond(...) function uses multigrid v-cycle(s)

Referenced by MultilevelLinearOp< T >::define(), MultilevelLinearOp< T >::MultilevelLinearOp(), and MultilevelLinearOp< T >::preCond().

template<class T>
int MultilevelLinearOp< T >::m_num_mg_iterations

number of mg v-cycles for each preCond call (if using mg preconditioner)

Referenced by MultilevelLinearOp< T >::MultilevelLinearOp(), and MultilevelLinearOp< T >::preCond().

template<class T>
int MultilevelLinearOp< T >::m_num_mg_smooth

number of smoothings in the multigrid preconditioner

Referenced by MultilevelLinearOp< T >::define(), and MultilevelLinearOp< T >::MultilevelLinearOp().

template<class T>
int MultilevelLinearOp< T >::m_preCondSolverDepth

max depth for multigrid preconditioner (if using mg preconditioner)

Referenced by MultilevelLinearOp< T >::define(), and MultilevelLinearOp< T >::MultilevelLinearOp().

template<class T>
Vector<int> MultilevelLinearOp< T >::m_refRatios

template<class T>
Vector<RealVect> MultilevelLinearOp< T >::m_vectDx

template<class T>
int MultilevelLinearOp< T >::m_lBase

Solver which is used by the preCond function.

if m_use_multigrid_preconditioner is true, this solver is used to apply a AMR multigrid v-cycle as the preconditioner

Referenced by MultilevelLinearOp< T >::define(), and MultilevelLinearOp< T >::preCond().

template<class T>
bool MultilevelLinearOp< T >::m_isPrecondSolverInitialized

has the preconditioner solver been initialized?

Referenced by MultilevelLinearOp< T >::define(), and MultilevelLinearOp< T >::preCond().


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

Generated on Tue Apr 14 14:23:45 2009 for Chombo + EB by  doxygen 1.5.5