From: Guido Kanschat Date: Mon, 25 Apr 2005 20:34:34 +0000 (+0000) Subject: a block smoother class X-Git-Tag: v8.0.0~13978 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=aeaa22f2a9853b917420df00d098b4f4e09d6451;p=dealii.git a block smoother class git-svn-id: https://svn.dealii.org/trunk@10577 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_block_smoother.h b/deal.II/deal.II/include/multigrid/mg_block_smoother.h new file mode 100644 index 0000000000..902beb99d3 --- /dev/null +++ b/deal.II/deal.II/include/multigrid/mg_block_smoother.h @@ -0,0 +1,330 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2005 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- +#ifndef __deal2__mg_block_smoother_h +#define __deal2__mg_block_smoother_h + + +#include +#include +#include +#include +#include +#include +#include +#include + +template class MGDoFHandler; + + +/* + * MGSmootherBase is defined in mg_base.h + */ + +/*!@addtogroup mg */ +/*@{*/ + +/** + * General smoother class for block vectors. This class gives complete + * freedom to the choice of a block smoother by being initialized with + * a matrix and a smoother object. Therefore, the smoother object for + * each level must be constructed by hand. + * + * @author Guido Kanschat, 2005 + */ +template +class MGSmootherBlock + : public MGSmootherBase > +{ + public: + /** + * Constructor. Sets memory and + * smoothing parameters. + */ + MGSmootherBlock(VectorMemory >& mem, + const unsigned int steps = 1, + const bool variable = false, + const bool symmetric = false, + const bool transpose = false, + const bool reverse = false); + + /** + * Initialize for matrices. The + * parameter matrices can be + * any object having functions + * get_minlevel() and + * get_maxlevel() as well as + * an operator[] returning a + * reference to @p MATRIX. + * + * The same convention is used + * for the parameter + * smoothers, such that + * operator[] returns + * the object doing the + * block-smoothing on a single + * level. + * + * This function stores pointers + * to the level matrices and + * smoothing operator for each + * level. + */ + template + void initialize (const MGMATRIX& matrices, + const MGRELAX& smoothers); + + /** + * Empty all vectors. + */ + void clear (); + + /** + * Modify the number of smoothing + * steps on finest level. + */ + void set_steps (const unsigned int); + + /** + * Switch on/off variable + * smoothing. + */ + void set_variable (const bool); + + /** + * Switch on/off symmetric + * smoothing. + */ + void set_symmetric (const bool); + + /** + * Switch on/off transposed. This + * is mutually exclusive with + * reverse(). + */ + void set_transpose (const bool); + + /** + * Switch on/off reversed. This + * is mutually exclusive with + * transpose(). + */ + void set_reverse (const bool); + + /** + * Implementation of the + * interface for @p Multigrid. + * This function does nothing, + * which by comparison with the + * definition of this function + * means that the the smoothing + * operator equals the null + * operator. + */ + virtual void smooth (const unsigned int level, + BlockVector& u, + const BlockVector& rhs) const; + private: + /** + * Pointer to the matrices. + */ + MGLevelObject > > matrices; + + /** + * Pointer to the matrices. + */ + MGLevelObject > > smoothers; + + /** + * Number of smoothing steps. + */ + unsigned int steps; + + /** + * Variable smoothing? + */ + bool variable; + + /** + * Symmetric smoothing? + */ + bool symmetric; + + /* + * Transposed? + */ + bool transpose; + + /** + * Reverse? + */ + bool reverse; + + /** + * Memory for auxiliary vectors. + */ + VectorMemory >& mem; + +}; + + +//--------------------------------------------------------------------------- + +template +inline +MGSmootherBlock::MGSmootherBlock( + VectorMemory >& mem, + const unsigned int steps, + const bool variable, + const bool symmetric, + const bool transpose, + const bool reverse) + : + steps(steps), + variable(variable), + symmetric(symmetric), + transpose(transpose), + reverse(reverse), + mem(mem) +{} + + +template +inline void +MGSmootherBlock::clear () +{ + unsigned int i=matrices.get_minlevel(), + max_level=matrices.get_maxlevel(); + for (; i<=max_level; ++i) + { + smoothers[i] = 0; + matrices[i] = 0; + } +} + + +template +template +inline void +MGSmootherBlock::initialize ( + const MGMATRIX& m, + const MGRELAX& s) +{ + const unsigned int min = m.get_minlevel(); + const unsigned int max = m.get_maxlevel(); + + matrices.resize(min, max); + smoothers.resize(min, max); + + for (unsigned int i=min;i<=max;++i) + { + matrices[i] = &m[i]; + smoothers[i] = &s[i]; + } +} + +template +inline void +MGSmootherBlock:: +set_steps (const unsigned int s) +{ + steps = s; +} + + +template +inline void +MGSmootherBlock:: +set_variable (const bool flag) +{ + variable = flag; +} + + +template +inline void +MGSmootherBlock:: +set_symmetric (const bool flag) +{ + symmetric = flag; +} + + +template +inline void +MGSmootherBlock:: +set_transpose (const bool flag) +{ + transpose = flag; +} + + +template +inline void +MGSmootherBlock:: +set_reverse (const bool flag) +{ + reverse = flag; +} + + +template +inline void +MGSmootherBlock::smooth( + const unsigned int level, + BlockVector& u, + const BlockVector& rhs) const +{ + unsigned int maxlevel = matrices.get_maxlevel(); + unsigned int steps2 = steps; + + if (variable) + steps2 *= (1<<(maxlevel-level)); + + BlockVector* r = mem.alloc(); + BlockVector* d = mem.alloc(); + r->reinit(u); + d->reinit(u); + + bool T = transpose; + if (symmetric && (steps2 % 2 == 0)) + T = false; +// cerr << 'S' << level; +// cerr << '(' << matrices[level]->m() << ',' << matrices[level]->n() << ')'; + + for (unsigned int i=0; isadd(-1.,1.,rhs); + smoothers[level].Tvmult(*d, *r); + } else { +// cerr << 'N'; + matrices[level].vmult(*r,u); + r->sadd(-1.,1.,rhs); + smoothers[level].vmult(*d, *r); + } +// cerr << '{' << r->l2_norm() << '}'; + u += *d; + if (symmetric) + T = !T; + } + + mem.free(r); + mem.free(d); +} + + + +#endif