ug4
ug::AssembledMultiGridCycle< TDomain, TAlgebra > Class Template Reference

geometric multi grid preconditioner More...

#include <mg_solver.h>

+ Inheritance diagram for ug::AssembledMultiGridCycle< TDomain, TAlgebra >:

Classes

struct  LevData
 
struct  LevelIndex
 
struct  SurfLevelMap
 

Public Types

typedef TAlgebra algebra_type
 Algebra type. More...
 
typedef TDomain domain_type
 Domain. More...
 
typedef GridFunction< TDomain, TAlgebra > GF
 Grid Function type. More...
 
typedef algebra_type::matrix_type matrix_type
 Matrix type. More...
 
typedef MGStats< TDomain, TAlgebra > mg_stats_type
 MGStats type. More...
 
typedef algebra_type::vector_type vector_type
 Vector type. More...
 
- Public Types inherited from ug::ILinearIterator< TAlgebra::vector_type >
typedef TAlgebra::vector_type codomain_function_type
 Range space. More...
 
typedef TAlgebra::vector_type domain_function_type
 Domain space. More...
 

Public Member Functions

void add_prolongation_post_process (SmartPtr< ITransferPostProcess< TDomain, TAlgebra > > PP)
 add prolongation post process More...
 
void add_restriction_post_process (SmartPtr< ITransferPostProcess< TDomain, TAlgebra > > PP)
 add restriction post process More...
 
 AssembledMultiGridCycle ()
 constructor without arguments More...
 
 AssembledMultiGridCycle (SmartPtr< ApproximationSpace< TDomain > > approxSpace)
 constructor setting approximation space More...
 
void clear_transfer_post_process ()
 clears all transfer post process More...
 
virtual std::string config_string () const
 returns information about configuration parameters More...
 
virtual bool init (SmartPtr< ILinearOperator< vector_type > > J, const vector_type &u)
 Prepare for Operator J(u) and linearization point u (current solution) More...
 
virtual bool init (SmartPtr< ILinearOperator< vector_type > > L)
 Prepare for Linear Operator L. More...
 
virtual const char * name () const
 name More...
 
void set_approximation_space (SmartPtr< ApproximationSpace< TDomain > > approxSpace)
 sets the approximation space More...
 
void set_base_level (int baseLevel)
 sets the level where exact solving is performed in the mg cycle More...
 
void set_base_solver (SmartPtr< ILinearOperatorInverse< vector_type > > baseSolver)
 sets the base solver that is used More...
 
void set_comm_comp_overlap (bool bOverlap)
 sets if communication and computation should be overlaped More...
 
void set_cycle_type (const std::string &type)
 sets the cycle type (1 = V-cycle, 2 = W-cycle, ...) More...
 
void set_cycle_type (int type)
 sets the cycle type (1 = V-cycle, 2 = W-cycle, ...) More...
 
void set_debug (SmartPtr< GridFunctionDebugWriter< TDomain, TAlgebra > > spDebugWriter)
 set debug output More...
 
void set_discretization (SmartPtr< IAssemble< TAlgebra > > spAss)
 sets the assembling procedure that is used to compute coarse grid matrices More...
 
void set_emulate_full_refined_grid (bool bEmulate)
 sets if copies should be used to emulate a full-refined grid More...
 
void set_gathered_base_solver_if_ambiguous (bool bGathered)
 sets if the base solver is applied in parallel More...
 
void set_mg_stats (SmartPtr< mg_stats_type > mgstats)
 sets MGStats, an object which records statistics for individual iterations More...
 
void set_num_postsmooth (int num)
 sets the number of post-smoothing steps to be performed More...
 
void set_num_presmooth (int num)
 sets the number of pre-smoothing steps to be performed More...
 
void set_postsmoother (SmartPtr< ILinearIterator< vector_type > > smoother)
 sets the post-smoother that is used More...
 
void set_presmoother (SmartPtr< ILinearIterator< vector_type > > smoother)
 sets the pre-smoother that is used More...
 
void set_projection (SmartPtr< ITransferOperator< TDomain, TAlgebra > > P)
 sets the projection operator More...
 
void set_prolongation (SmartPtr< ITransferOperator< TDomain, TAlgebra > > P)
 sets the prolongation operator More...
 
void set_rap (bool bRAP)
 sets if RAP - Product used to build coarse grid matrices More...
 
void set_restriction (SmartPtr< ITransferOperator< TDomain, TAlgebra > > P)
 sets the restriction operator More...
 
void set_smooth_on_surface_rim (bool bSmooth)
 sets if smoothing is performed on surface rim More...
 
void set_smoother (SmartPtr< ILinearIterator< vector_type > > smoother)
 sets the smoother that is used More...
 
void set_surface_level (int topLevel)
 sets the surface level (default is top-surface) More...
 
void set_transfer (SmartPtr< ITransferOperator< TDomain, TAlgebra > > P)
 sets the transfer operator More...
 
virtual bool supports_parallel () const
 returns if parallel solving is supported More...
 
- Public Member Functions inherited from ug::ILinearIterator< TAlgebra::vector_type >
virtual bool apply (TAlgebra::vector_type &c, const TAlgebra::vector_type &d)=0
 compute new correction c = B*d More...
 
virtual bool apply_update_defect (TAlgebra::vector_type &c, TAlgebra::vector_type &d)=0
 compute new correction c = B*d and update defect d := d - A*c More...
 
SmartPtr< IDamping< TAlgebra::vector_type, TAlgebra::vector_type > > damping ()
 returns the scaling More...
 
 ILinearIterator ()
 constructor More...
 
 ILinearIterator (const ILinearIterator< TAlgebra::vector_type, TAlgebra::vector_type > &parent)
 copy constructor More...
 
virtual bool init (SmartPtr< ILinearOperator< TAlgebra::vector_type, TAlgebra::vector_type > > J, const TAlgebra::vector_type &u)=0
 initialize for operator J(u) and linearization point u More...
 
virtual bool init (SmartPtr< ILinearOperator< TAlgebra::vector_type, TAlgebra::vector_type > > L)=0
 initialize for linear operator L More...
 
void set_damp (number factor)
 sets the damping to a constant factor More...
 
void set_damp (SmartPtr< IDamping< TAlgebra::vector_type, TAlgebra::vector_type > > spScaling)
 sets a scaling for the correction More...
 
virtual ~ILinearIterator ()
 virtual destructor More...
 

Protected Member Functions

void copy_ghost_to_noghost (SmartPtr< GF > spVecTo, ConstSmartPtr< GF > spVecFrom, const std::vector< size_t > &vMapPatchToGlobal)
 copies vector to smoothing patch using cached mapping More...
 
void copy_noghost_to_ghost (SmartPtr< GF > spVecTo, ConstSmartPtr< GF > spVecFrom, const std::vector< size_t > &vMapPatchToGlobal)
 copies vector from smoothing patch using cached mapping More...
 
void copy_noghost_to_ghost (SmartPtr< matrix_type > spMatTo, ConstSmartPtr< matrix_type > spMatFrom, const std::vector< size_t > &vMapPatchToGlobal)
 copies matrix from smoothing patch using cached mapping More...
 
void enter_debug_writer_section (GridLevel &orig_gl, const char *sec_name, int lev, int cycleNo=-1, int callNo=-1)
 enters a new debugger section for smoothers, base solver etc More...
 
bool gathered_base_master () const
 returns if gathered base master More...
 
void leave_debug_writer_section (GridLevel &orig_gl)
 leaves the current debugger section More...
 
void log_debug_data (int lvl, int cycleNo, std::string name)
 logs a level-data-struct to the terminal More...
 
void mg_stats_defect (GF &gf, int lvl, typename mg_stats_type::Stage stage)
 Calls MGStats::set_defect (if available) with the given parameters. More...
 
void write_debug (const GF &rGF, std::string name, int cycleNo=-1)
 
void write_debug (ConstSmartPtr< GF > spGF, std::string name, int cycleNo=-1)
 writes debug output for a level vector only on smooth path More...
 
template<typename TElem >
void init_noghost_to_ghost_mapping (std::vector< size_t > &vNoGhostToGhostMap, ConstSmartPtr< DoFDistribution > spNoGhostDD, ConstSmartPtr< DoFDistribution > spGhostDD)
 
void init_noghost_to_ghost_mapping (int lev)
 
void write_debug (const matrix_type &mat, std::string name, const GridLevel &glTo, const GridLevel &glFrom)
 writes debug output for a level matrix only on smooth path More...
 
void write_debug (const matrix_type &mat, std::string name, const GF &rTo, const GF &rFrom)
 writes debug output for a level matrix only on smooth path More...
 

Protected Attributes

bool m_bGatheredBaseIfAmbiguous
 flag, if to solve base problem in parallel when gathered and (!) parallel possible More...
 
bool m_bGatheredBaseUsed
 flag if using parallel base solver More...
 
bool m_bMatrixStructureIsConst
 flag indicating whether the matrix structure is supposed to stay the same in the next init phase More...
 
pcl::InterfaceCommunicator< IndexLayoutm_Com
 communicator More...
 
int m_dbgIterCnt
 counter for debug, to distinguish the iterations More...
 
bool m_ignoreInitForBaseSolver
 flag indicating whether the base-solver shall be initialized during a call to init() More...
 
SmartPtr< mg_stats_typem_mgstats
 MGStats are used for debugging to record statistics on individual iterations. More...
 
GFm_pC
 current surface correction More...
 
SmartPtr< GridFunctionDebugWriter< TDomain, TAlgebra > > m_spDebugWriter
 Debug Writer. More...
 
std::vector< SmartPtr< LevData > > m_vLevData
 storage for all level More...
 
SmartPtr< GFspGatheredBaseCorr
 vector for gathered base solver More...
 
SmartPtr< MatrixOperator< matrix_type, vector_type > > spGatheredBaseMat
 Matrix for gathered base solver. More...
 
- Protected Attributes inherited from ug::ILinearIterator< TAlgebra::vector_type >
SmartPtr< IDamping< TAlgebra::vector_type, TAlgebra::vector_type > > m_spDamping
 the scaling More...
 
ConstSmartPtr< matrix_typem_spSurfaceMat
 operator to invert (surface grid) More...
 
const vector_typem_pSurfaceSol
 Solution on surface grid. More...
 
SmartPtr< IAssemble< TAlgebra > > m_spAss
 assembling routine for coarse grid matrices More...
 
SmartPtr< ApproximationSpace< TDomain > > m_spApproxSpace
 approximation space for level and surface grid More...
 
int m_topLev
 
int m_surfaceLev
 operator to invert (surface grid) More...
 
int m_baseLev
 base level (where exact inverse is computed) More...
 
int m_cycleType
 cylce type (1 = V-cycle, 2 = W-cylcle, ...) More...
 
int m_numPreSmooth
 number of Presmooth steps More...
 
int m_numPostSmooth
 number of Postsmooth steps More...
 
int m_LocalFullRefLevel
 lowest level containing surface geom obj (proc-locally) More...
 
GridLevel::ViewType m_GridLevelType
 grid-view for level vectors More...
 
bool m_bUseRAP
 using RAP-Product (assemble coarse-grid matrices otherwise) More...
 
bool m_bSmoothOnSurfaceRim
 flag if smoothing on surface rim More...
 
bool m_bCommCompOverlap
 flag if overlapping communication and computation More...
 
RevisionCounter m_ApproxSpaceRevision
 approximation space revision of cached values More...
 
SmartPtr< ILinearIterator< vector_type > > m_spPreSmootherPrototype
 prototype for pre-smoother More...
 
SmartPtr< ILinearIterator< vector_type > > m_spPostSmootherPrototype
 prototype for post-smoother More...
 
SmartPtr< ITransferOperator< TDomain, TAlgebra > > m_spProjectionPrototype
 prototype for projection operator More...
 
SmartPtr< ITransferOperator< TDomain, TAlgebra > > m_spProlongationPrototype
 prototype for prolongation operator More...
 
SmartPtr< ITransferOperator< TDomain, TAlgebra > > m_spRestrictionPrototype
 prototype for restriction operator More...
 
std::vector< SmartPtr< ITransferPostProcess< TDomain, TAlgebra > > > m_vspProlongationPostProcess
 prototpe for transfer post process More...
 
std::vector< SmartPtr< ITransferPostProcess< TDomain, TAlgebra > > > m_vspRestrictionPostProcess
 operator to invert (surface grid) More...
 
SmartPtr< ILinearOperatorInverse< vector_type > > m_spBaseSolver
 base solver for the coarse problem More...
 
std::vector< LevelIndexm_vSurfToLevelMap
 operator to invert (surface grid) More...
 
static const int _V_ = 1
 operator to invert (surface grid) More...
 
static const int _W_ = 2
 operator to invert (surface grid) More...
 
static const int _F_ = -1
 operator to invert (surface grid) More...
 
virtual void ignore_init_for_base_solver (bool ignore)
 does not call init on base-solver during initialization More...
 
virtual bool ignore_init_for_base_solver () const
 operator to invert (surface grid) More...
 
void set_matrix_structure_is_const (bool b)
 
void force_reinit ()
 reinit transfer operators More...
 
virtual bool apply (vector_type &c, const vector_type &d)
 Compute new correction c = B*d. More...
 
virtual bool apply_update_defect (vector_type &c, vector_type &d)
 Compute new correction c = B*d and return new defect d := d - A*c. More...
 
SmartPtr< ILinearIterator< vector_type > > clone ()
 Clone. More...
 
 ~AssembledMultiGridCycle ()
 Destructor. More...
 
void lmgc (int lev, int cycleType)
 compute correction on level and update defect More...
 
void presmooth_and_restriction (int lev)
 performs presmoothing on the given level More...
 
void prolongation_and_postsmooth (int lev)
 performs prolongation to the level above More...
 
void base_solve (int lev)
 compute base solver More...
 
void init_level_memory (int baseLev, int topLev)
 allocates the memory More...
 
void init ()
 initializes common part More...
 
void init_smoother ()
 initializes the smoother and base solver More...
 
void assemble_level_operator ()
 initializes the coarse grid matrices More...
 
void init_rap_operator ()
 operator to invert (surface grid) More...
 
void init_base_solver ()
 initializes the smoother and base solver More...
 
void init_transfer ()
 initializes the transfers More...
 
void init_projection ()
 initializes the prolongation More...
 
void assemble_rim_cpl (const vector_type *u)
 
void init_rap_rim_cpl ()
 operator to invert (surface grid) More...
 
template<typename TElem >
void init_index_mappings ()
 operator to invert (surface grid) More...
 
void init_index_mappings ()
 operator to invert (surface grid) More...
 

Detailed Description

template<typename TDomain, typename TAlgebra>
class ug::AssembledMultiGridCycle< TDomain, TAlgebra >

geometric multi grid preconditioner

This class implements one step of the geometric multi grid as a preconditioner for linear iteration schemes such as linear iteration, cg or bicgstab.

The coarse grid spaces are build up according to the Approximation Space that is set from outside. In addition an Assembling routine must be specified that is used to assemble the coarse grid matrices.

Template Parameters
TApproximationSpaceType of Approximation Space
TAlgebraType of Algebra

Member Typedef Documentation

◆ algebra_type

template<typename TDomain , typename TAlgebra >
typedef TAlgebra ug::AssembledMultiGridCycle< TDomain, TAlgebra >::algebra_type

Algebra type.

◆ domain_type

template<typename TDomain , typename TAlgebra >
typedef TDomain ug::AssembledMultiGridCycle< TDomain, TAlgebra >::domain_type

◆ GF

template<typename TDomain , typename TAlgebra >
typedef GridFunction<TDomain, TAlgebra> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::GF

Grid Function type.

◆ matrix_type

template<typename TDomain , typename TAlgebra >
typedef algebra_type::matrix_type ug::AssembledMultiGridCycle< TDomain, TAlgebra >::matrix_type

Matrix type.

◆ mg_stats_type

template<typename TDomain , typename TAlgebra >
typedef MGStats<TDomain, TAlgebra> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::mg_stats_type

MGStats type.

◆ vector_type

template<typename TDomain , typename TAlgebra >
typedef algebra_type::vector_type ug::AssembledMultiGridCycle< TDomain, TAlgebra >::vector_type

Vector type.

Constructor & Destructor Documentation

◆ AssembledMultiGridCycle() [1/2]

template<typename TDomain , typename TAlgebra >
ug::AssembledMultiGridCycle< TDomain, TAlgebra >::AssembledMultiGridCycle

constructor without arguments

◆ AssembledMultiGridCycle() [2/2]

template<typename TDomain , typename TAlgebra >
ug::AssembledMultiGridCycle< TDomain, TAlgebra >::AssembledMultiGridCycle ( SmartPtr< ApproximationSpace< TDomain > >  approxSpace)

constructor setting approximation space

◆ ~AssembledMultiGridCycle()

template<typename TDomain , typename TAlgebra >
ug::AssembledMultiGridCycle< TDomain, TAlgebra >::~AssembledMultiGridCycle

Destructor.

Member Function Documentation

◆ add_prolongation_post_process()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::add_prolongation_post_process ( SmartPtr< ITransferPostProcess< TDomain, TAlgebra > >  PP)
inline

◆ add_restriction_post_process()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::add_restriction_post_process ( SmartPtr< ITransferPostProcess< TDomain, TAlgebra > >  PP)
inline

◆ apply()

◆ apply_update_defect()

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::apply_update_defect ( vector_type c,
vector_type d 
)
virtual

Compute new correction c = B*d and return new defect d := d - A*c.

References GMG_PROFILE_FUNC, and UG_THROW.

◆ assemble_level_operator()

◆ assemble_rim_cpl()

◆ base_solve()

◆ clear_transfer_post_process()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::clear_transfer_post_process ( )
inline

◆ clone()

template<typename TDomain , typename TAlgebra >
SmartPtr< ILinearIterator< typename TAlgebra::vector_type > > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::clone
virtual

◆ config_string()

template<typename TDomain , typename TAlgebra >
std::string ug::AssembledMultiGridCycle< TDomain, TAlgebra >::config_string
virtual

returns information about configuration parameters

Reimplemented from ug::ILinearIterator< TAlgebra::vector_type >.

References ug::ConfigShift().

◆ copy_ghost_to_noghost()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::copy_ghost_to_noghost ( SmartPtr< GF spVecTo,
ConstSmartPtr< GF spVecFrom,
const std::vector< size_t > &  vMapPatchToGlobal 
)
protected

copies vector to smoothing patch using cached mapping

References GMG_PROFILE_FUNC, UG_ASSERT, and UG_THROW.

◆ copy_noghost_to_ghost() [1/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::copy_noghost_to_ghost ( SmartPtr< GF spVecTo,
ConstSmartPtr< GF spVecFrom,
const std::vector< size_t > &  vMapPatchToGlobal 
)
protected

copies vector from smoothing patch using cached mapping

References GMG_PROFILE_FUNC, UG_ASSERT, and UG_THROW.

◆ copy_noghost_to_ghost() [2/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::copy_noghost_to_ghost ( SmartPtr< matrix_type spMatTo,
ConstSmartPtr< matrix_type spMatFrom,
const std::vector< size_t > &  vMapPatchToGlobal 
)
protected

◆ enter_debug_writer_section()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::enter_debug_writer_section ( GridLevel orig_gl,
const char *  sec_name,
int  lev,
int  cycleNo = -1,
int  callNo = -1 
)
inlineprotected

enters a new debugger section for smoothers, base solver etc

References PROFILE_FUNC_GROUP.

◆ force_reinit()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::force_reinit

reinit transfer operators

◆ gathered_base_master()

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::gathered_base_master
protected

returns if gathered base master

References UG_THROW.

◆ ignore_init_for_base_solver() [1/2]

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::ignore_init_for_base_solver
virtual

operator to invert (surface grid)

◆ ignore_init_for_base_solver() [2/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::ignore_init_for_base_solver ( bool  ignore)
virtual

does not call init on base-solver during initialization

Use this method with care. It can be useful e.g. during repeated adaption of a linear problem, where the base-level does not change.

Warning
Currently only works in serial environments
Note
This only affects calls to init(...) and has no effect on direct calls to init_base_solver.

References pcl::NumProcs(), and UG_COND_THROW.

◆ init() [1/3]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init
protected

◆ init() [2/3]

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init ( SmartPtr< ILinearOperator< vector_type > >  J,
const vector_type u 
)
virtual

Prepare for Operator J(u) and linearization point u (current solution)

References GMG_PROFILE_FUNC, init(), ug::LIB_DISC_MULTIGRID, UG_DLOG, and SmartPtr< T, FreePolicy >::valid().

◆ init() [3/3]

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init ( SmartPtr< ILinearOperator< vector_type > >  L)
virtual

Prepare for Linear Operator L.

References GMG_PROFILE_FUNC, init(), ug::LIB_DISC_MULTIGRID, UG_DLOG, and SmartPtr< T, FreePolicy >::valid().

◆ init_base_solver()

◆ init_index_mappings() [1/2]

template<typename TDomain , typename TAlgebra >
template<typename TElem >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init_index_mappings
protected

◆ init_index_mappings() [2/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init_index_mappings
protected

operator to invert (surface grid)

References ug::EDGE, ug::FACE, GMG_PROFILE_FUNC, ug::GridLevel::SURFACE, UG_THROW, ug::VERTEX, and ug::VOLUME.

◆ init_level_memory()

◆ init_noghost_to_ghost_mapping() [1/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init_noghost_to_ghost_mapping ( int  lev)
protected

◆ init_noghost_to_ghost_mapping() [2/2]

template<typename TDomain , typename TAlgebra >
template<typename TElem >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::init_noghost_to_ghost_mapping ( std::vector< size_t > &  vNoGhostToGhostMap,
ConstSmartPtr< DoFDistribution spNoGhostDD,
ConstSmartPtr< DoFDistribution spGhostDD 
)
protected

init mapping from noghost -> w/ ghost

References UG_ASSERT.

◆ init_projection()

◆ init_rap_operator()

◆ init_rap_rim_cpl()

◆ init_smoother()

◆ init_transfer()

◆ leave_debug_writer_section()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::leave_debug_writer_section ( GridLevel orig_gl)
inlineprotected

leaves the current debugger section

◆ lmgc()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::lmgc ( int  lev,
int  cycleType 
)
protected

compute correction on level and update defect

References GMG_PROFILE_FUNC, ug::LIB_DISC_MULTIGRID, UG_CATCH_THROW, UG_DLOG, and UG_THROW.

◆ log_debug_data()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::log_debug_data ( int  lvl,
int  cycleNo,
std::string  name 
)
protected

◆ mg_stats_defect()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::mg_stats_defect ( GF gf,
int  lvl,
typename mg_stats_type::Stage  stage 
)
protected

Calls MGStats::set_defect (if available) with the given parameters.

◆ name()

template<typename TDomain , typename TAlgebra >
virtual const char* ug::AssembledMultiGridCycle< TDomain, TAlgebra >::name ( ) const
inlinevirtual

◆ presmooth_and_restriction()

◆ prolongation_and_postsmooth()

◆ set_approximation_space()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_approximation_space ( SmartPtr< ApproximationSpace< TDomain > >  approxSpace)

sets the approximation space

References make_sp().

◆ set_base_level()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_base_level ( int  baseLevel)
inline

sets the level where exact solving is performed in the mg cycle

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_baseLev.

◆ set_base_solver()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_base_solver ( SmartPtr< ILinearOperatorInverse< vector_type > >  baseSolver)
inline

sets the base solver that is used

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spBaseSolver.

◆ set_comm_comp_overlap()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_comm_comp_overlap ( bool  bOverlap)
inline

sets if communication and computation should be overlaped

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bCommCompOverlap.

◆ set_cycle_type() [1/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type ( const std::string &  type)
inline

◆ set_cycle_type() [2/2]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type ( int  type)
inline

sets the cycle type (1 = V-cycle, 2 = W-cycle, ...)

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_cycleType.

◆ set_debug()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_debug ( SmartPtr< GridFunctionDebugWriter< TDomain, TAlgebra > >  spDebugWriter)
inline

set debug output

If a DebugWriter is passed by this method, the multi grid cycle writes the level/surface vectors and matrices for debug purposes.

Parameters
[in]debugWriterDebug Writer to use

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spDebugWriter.

◆ set_discretization()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_discretization ( SmartPtr< IAssemble< TAlgebra > >  spAss)
inline

sets the assembling procedure that is used to compute coarse grid matrices

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spAss.

◆ set_emulate_full_refined_grid()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_emulate_full_refined_grid ( bool  bEmulate)
inline

sets if copies should be used to emulate a full-refined grid

References ug::GridLevel::LEVEL, ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_GridLevelType, and ug::GridLevel::SURFACE.

◆ set_gathered_base_solver_if_ambiguous()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_gathered_base_solver_if_ambiguous ( bool  bGathered)
inline

sets if the base solver is applied in parallel

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bGatheredBaseIfAmbiguous.

◆ set_matrix_structure_is_const()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_matrix_structure_is_const ( bool  b)
inline

\ If the nonzero entry pattern of the matrix does not change (which is usually the case in time-dependent problems without space adaptation), then keeping the entry structure between re-inits saves time during level matrix assembling.

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bMatrixStructureIsConst.

◆ set_mg_stats()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_mg_stats ( SmartPtr< mg_stats_type mgstats)
inline

sets MGStats, an object which records statistics for individual iterations

Setting MGStats is optional. If set, a runtime overhead is incurred. Only useful for debugging.

See also
MGStats

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_mgstats.

◆ set_num_postsmooth()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_num_postsmooth ( int  num)
inline

sets the number of post-smoothing steps to be performed

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_numPostSmooth.

◆ set_num_presmooth()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_num_presmooth ( int  num)
inline

sets the number of pre-smoothing steps to be performed

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_numPreSmooth.

◆ set_postsmoother()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_postsmoother ( SmartPtr< ILinearIterator< vector_type > >  smoother)
inline

◆ set_presmoother()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_presmoother ( SmartPtr< ILinearIterator< vector_type > >  smoother)
inline

◆ set_projection()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_projection ( SmartPtr< ITransferOperator< TDomain, TAlgebra > >  P)
inline

◆ set_prolongation()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_prolongation ( SmartPtr< ITransferOperator< TDomain, TAlgebra > >  P)
inline

◆ set_rap()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_rap ( bool  bRAP)
inline

sets if RAP - Product used to build coarse grid matrices

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bUseRAP.

◆ set_restriction()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_restriction ( SmartPtr< ITransferOperator< TDomain, TAlgebra > >  P)
inline

◆ set_smooth_on_surface_rim()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_smooth_on_surface_rim ( bool  bSmooth)
inline

sets if smoothing is performed on surface rim

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bSmoothOnSurfaceRim.

◆ set_smoother()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_smoother ( SmartPtr< ILinearIterator< vector_type > >  smoother)
inline

◆ set_surface_level()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_surface_level ( int  topLevel)
inline

sets the surface level (default is top-surface)

References ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_surfaceLev.

◆ set_transfer()

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_transfer ( SmartPtr< ITransferOperator< TDomain, TAlgebra > >  P)
inline

◆ supports_parallel()

template<typename TDomain , typename TAlgebra >
virtual bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::supports_parallel ( ) const
inlinevirtual

◆ write_debug() [1/4]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::write_debug ( const GF rGF,
std::string  name,
int  cycleNo = -1 
)
protected

◆ write_debug() [2/4]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::write_debug ( const matrix_type mat,
std::string  name,
const GF rTo,
const GF rFrom 
)
protected

writes debug output for a level matrix only on smooth path

This method writes the level matrix to a debug file, if a debug writer has been set.

Parameters
[in]matLevel Matrix to write for debug purpose
[in]nameFilename

References ug::GridFunction< TDomain, TAlgebra >::grid_level(), and name.

◆ write_debug() [3/4]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::write_debug ( const matrix_type mat,
std::string  name,
const GridLevel glTo,
const GridLevel glFrom 
)
protected

writes debug output for a level matrix only on smooth path

This method writes the level matrix to a debug file, if a debug writer has been set.

Parameters
[in]matLevel Matrix to write for debug purpose
[in]nameFilename

References ug::GridLevelAppendix(), name, and PROFILE_FUNC_GROUP.

◆ write_debug() [4/4]

template<typename TDomain , typename TAlgebra >
void ug::AssembledMultiGridCycle< TDomain, TAlgebra >::write_debug ( ConstSmartPtr< GF spGF,
std::string  name,
int  cycleNo = -1 
)
inlineprotected

writes debug output for a level vector only on smooth path

This method writes the level vector to a debug file, if a debug writer has been set.

Parameters
[in]spGFLevel Vector to write for debug purpose
[in]nameFilename

References name.

Member Data Documentation

◆ _F_

template<typename TDomain , typename TAlgebra >
const int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::_F_ = -1
staticprotected

operator to invert (surface grid)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type().

◆ _V_

template<typename TDomain , typename TAlgebra >
const int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::_V_ = 1
staticprotected

operator to invert (surface grid)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type().

◆ _W_

template<typename TDomain , typename TAlgebra >
const int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::_W_ = 2
staticprotected

operator to invert (surface grid)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type().

◆ m_ApproxSpaceRevision

template<typename TDomain , typename TAlgebra >
RevisionCounter ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_ApproxSpaceRevision
protected

approximation space revision of cached values

◆ m_baseLev

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_baseLev
protected

base level (where exact inverse is computed)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_base_level().

◆ m_bCommCompOverlap

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bCommCompOverlap
protected

flag if overlapping communication and computation

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_comm_comp_overlap().

◆ m_bGatheredBaseIfAmbiguous

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bGatheredBaseIfAmbiguous
protected

flag, if to solve base problem in parallel when gathered and (!) parallel possible

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_gathered_base_solver_if_ambiguous().

◆ m_bGatheredBaseUsed

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bGatheredBaseUsed
protected

flag if using parallel base solver

◆ m_bMatrixStructureIsConst

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bMatrixStructureIsConst
protected

flag indicating whether the matrix structure is supposed to stay the same in the next init phase

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_matrix_structure_is_const().

◆ m_bSmoothOnSurfaceRim

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bSmoothOnSurfaceRim
protected

◆ m_bUseRAP

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_bUseRAP
protected

using RAP-Product (assemble coarse-grid matrices otherwise)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_rap().

◆ m_Com

template<typename TDomain , typename TAlgebra >
pcl::InterfaceCommunicator<IndexLayout> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_Com
protected

communicator

◆ m_cycleType

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_cycleType
protected

cylce type (1 = V-cycle, 2 = W-cylcle, ...)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_cycle_type().

◆ m_dbgIterCnt

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_dbgIterCnt
protected

counter for debug, to distinguish the iterations

◆ m_GridLevelType

template<typename TDomain , typename TAlgebra >
GridLevel::ViewType ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_GridLevelType
protected

◆ m_ignoreInitForBaseSolver

template<typename TDomain , typename TAlgebra >
bool ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_ignoreInitForBaseSolver
protected

flag indicating whether the base-solver shall be initialized during a call to init()

◆ m_LocalFullRefLevel

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_LocalFullRefLevel
protected

lowest level containing surface geom obj (proc-locally)

◆ m_mgstats

template<typename TDomain , typename TAlgebra >
SmartPtr<mg_stats_type> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_mgstats
protected

MGStats are used for debugging to record statistics on individual iterations.

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_mg_stats().

◆ m_numPostSmooth

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_numPostSmooth
protected

◆ m_numPreSmooth

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_numPreSmooth
protected

◆ m_pC

template<typename TDomain , typename TAlgebra >
GF* ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_pC
protected

current surface correction

◆ m_pSurfaceSol

template<typename TDomain , typename TAlgebra >
const vector_type* ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_pSurfaceSol
protected

Solution on surface grid.

◆ m_spApproxSpace

template<typename TDomain , typename TAlgebra >
SmartPtr<ApproximationSpace<TDomain> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spApproxSpace
protected

approximation space for level and surface grid

◆ m_spAss

template<typename TDomain , typename TAlgebra >
SmartPtr<IAssemble<TAlgebra> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spAss
protected

assembling routine for coarse grid matrices

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_discretization().

◆ m_spBaseSolver

template<typename TDomain , typename TAlgebra >
SmartPtr<ILinearOperatorInverse<vector_type> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spBaseSolver
protected

base solver for the coarse problem

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_base_solver().

◆ m_spDebugWriter

template<typename TDomain , typename TAlgebra >
SmartPtr<GridFunctionDebugWriter<TDomain, TAlgebra> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spDebugWriter
protected

◆ m_spPostSmootherPrototype

template<typename TDomain , typename TAlgebra >
SmartPtr<ILinearIterator<vector_type> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spPostSmootherPrototype
protected

◆ m_spPreSmootherPrototype

template<typename TDomain , typename TAlgebra >
SmartPtr<ILinearIterator<vector_type> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spPreSmootherPrototype
protected

◆ m_spProjectionPrototype

template<typename TDomain , typename TAlgebra >
SmartPtr<ITransferOperator<TDomain, TAlgebra> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spProjectionPrototype
protected

prototype for projection operator

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_projection().

◆ m_spProlongationPrototype

template<typename TDomain , typename TAlgebra >
SmartPtr<ITransferOperator<TDomain, TAlgebra> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spProlongationPrototype
protected

prototype for prolongation operator

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_prolongation().

◆ m_spRestrictionPrototype

template<typename TDomain , typename TAlgebra >
SmartPtr<ITransferOperator<TDomain, TAlgebra> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spRestrictionPrototype
protected

prototype for restriction operator

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_restriction().

◆ m_spSurfaceMat

template<typename TDomain , typename TAlgebra >
ConstSmartPtr<matrix_type> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_spSurfaceMat
protected

operator to invert (surface grid)

◆ m_surfaceLev

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_surfaceLev
protected

operator to invert (surface grid)

Referenced by ug::AssembledMultiGridCycle< TDomain, TAlgebra >::set_surface_level().

◆ m_topLev

template<typename TDomain , typename TAlgebra >
int ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_topLev
protected

top level (i.e. highest level in hierarchy. This is the surface level in case of non-adaptive refinement)

◆ m_vLevData

template<typename TDomain , typename TAlgebra >
std::vector<SmartPtr<LevData> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_vLevData
protected

storage for all level

◆ m_vspProlongationPostProcess

template<typename TDomain , typename TAlgebra >
std::vector<SmartPtr<ITransferPostProcess<TDomain, TAlgebra> > > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_vspProlongationPostProcess
protected

◆ m_vspRestrictionPostProcess

template<typename TDomain , typename TAlgebra >
std::vector<SmartPtr<ITransferPostProcess<TDomain, TAlgebra> > > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_vspRestrictionPostProcess
protected

◆ m_vSurfToLevelMap

template<typename TDomain , typename TAlgebra >
std::vector<LevelIndex> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::m_vSurfToLevelMap
protected

operator to invert (surface grid)

◆ spGatheredBaseCorr

template<typename TDomain , typename TAlgebra >
SmartPtr<GF> ug::AssembledMultiGridCycle< TDomain, TAlgebra >::spGatheredBaseCorr
protected

vector for gathered base solver

◆ spGatheredBaseMat

template<typename TDomain , typename TAlgebra >
SmartPtr<MatrixOperator<matrix_type, vector_type> > ug::AssembledMultiGridCycle< TDomain, TAlgebra >::spGatheredBaseMat
protected

Matrix for gathered base solver.


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