From: Guido Kanschat Date: Thu, 2 Feb 2006 09:29:11 +0000 (+0000) Subject: new class PoiterMatrixAux X-Git-Tag: v8.0.0~12456 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=784dd932e9a9de5871f03538fda1591dc6784601;p=dealii.git new class PoiterMatrixAux git-svn-id: https://svn.dealii.org/trunk@12223 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index c14604d344..0576cb49bf 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -298,6 +298,23 @@ inconvenience this causes.

lac

    + + + >li>

    Improved: BlockMatrixArray::enter_aux allows using matrices without adding vector + muiltplication by using PointerMatrixAux. +
    + (GK 2006/02/02) +

    + +
  1. New: The class PointerMatrixAux was + introduced for use with matrix classes lacking the adding vector + multiplication functions. It implements these by using its own auxiliary + vector. +
    + (GK 2006/02/02) +

    +
  2. Improved: The FilteredMatrix class was able to filter only the SparseMatrix diff --git a/deal.II/lac/include/lac/block_matrix_array.h b/deal.II/lac/include/lac/block_matrix_array.h index 4a867ecfd5..87320b854b 100644 --- a/deal.II/lac/include/lac/block_matrix_array.h +++ b/deal.II/lac/include/lac/block_matrix_array.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors +// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -185,6 +185,21 @@ class BlockMatrixArray : public Subscriptor const unsigned int col, const double prefix = 1., const bool transpose = false); + + /** + * Add an entry like with enter, + * but use PointerMatrixAux for + * matrices not having functions + * vmult_add() and TVmult_add(). + */ + template + void enter_aux (VectorMemory >& mem, + const MATRIX &matrix, + const unsigned int row, + const unsigned int col, + const double prefix = 1., + const bool transpose = false); + /** * Delete all entries, i.e. reset @@ -284,13 +299,30 @@ class BlockMatrixArray : public Subscriptor { public: /** - * Constructor initializing all - * data fields. + * Constructor initializing + * all data fields. A + * PointerMatrix object is + * generated for + * matrix. */ template Entry (const MATRIX& matrix, unsigned row, unsigned int col, double prefix, bool transpose); + + /** + * Constructor initializing + * all data fields. A + * PointerMatrixAux object is + * generated for + * matrix. + */ + template + Entry (const MATRIX& matrix, + unsigned row, unsigned int col, + double prefix, bool transpose, + VectorMemory >& mem); + /** * Copy constructor * invalidating the old @@ -484,6 +516,21 @@ class BlockTrianglePrecondition const double prefix = 1., const bool transpose = false); + /** + * Enter a block. This calls + * BlockMatrixArray::enter_aux(). Remember + * that the diagonal blocks + * should actually be inverse + * matrices or preconditioners. + */ + template + void enter_aux (VectorMemory >& mem, + const MATRIX &matrix, + const unsigned int row, + const unsigned int col, + const double prefix = 1., + const bool transpose = false); + /** * Preconditioning. */ @@ -578,15 +625,38 @@ class BlockTrianglePrecondition template template inline -BlockMatrixArray::Entry::Entry (const MATRIX& matrix, - unsigned row, unsigned int col, - double prefix, bool transpose) +BlockMatrixArray::Entry::Entry ( + const MATRIX& m, + unsigned int row, + unsigned int col, + double prefix, + bool transpose) : row (row), col (col), prefix (prefix), transpose (transpose), - matrix (new PointerMatrix >(&matrix, typeid(*this).name())) + matrix (new PointerMatrix >(&m, typeid(*this).name())) +{} + + + +template +template +inline +BlockMatrixArray::Entry::Entry ( + const MATRIX& m, + unsigned row, + unsigned int col, + double prefix, + bool transpose, + VectorMemory >& mem) + : + row (row), + col (col), + prefix (prefix), + transpose (transpose), + matrix (new PointerMatrixAux >(&mem, &m, typeid(*this).name())) {} @@ -595,9 +665,12 @@ template template inline void -BlockMatrixArray::enter (const MATRIX& matrix, - unsigned row, unsigned int col, - double prefix, bool transpose) +BlockMatrixArray::enter ( + const MATRIX& matrix, + unsigned int row, + unsigned int col, + double prefix, + bool transpose) { Assert (mem != 0, ExcNotInitialized()); Assert(row::enter (const MATRIX& matrix, } +template +template +inline +void +BlockMatrixArray::enter_aux ( + VectorMemory >& mem, + const MATRIX& matrix, + unsigned int row, + unsigned int col, + double prefix, + bool transpose) +{ + Assert(row struct BlockMatrixArrayPointerMatrixLess { @@ -704,6 +795,23 @@ BlockTrianglePrecondition::enter (const MATRIX& matrix, +template +template +inline +void +BlockTrianglePrecondition::enter_aux ( + VectorMemory >& mem, + const MATRIX& matrix, + unsigned int row, + unsigned int col, + double prefix, + bool transpose) +{ + BlockMatrixArray::enter_aux(mem, matrix, row, col, prefix, transpose); +} + + + #endif // DOXYGEN #endif diff --git a/deal.II/lac/include/lac/pointer_matrix.h b/deal.II/lac/include/lac/pointer_matrix.h index 690da5a9c5..c38ea22c5d 100644 --- a/deal.II/lac/include/lac/pointer_matrix.h +++ b/deal.II/lac/include/lac/pointer_matrix.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors +// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -16,6 +16,8 @@ #include #include +template class VectorMemory; + /*! @addtogroup Matrix2 *@{ */ @@ -233,6 +235,122 @@ class PointerMatrix : public PointerMatrixBase }; +/** + * A pointer to be used as a matrix. This class stores a pointer to a + * matrix and can be used as a matrix itself in iterative methods. + * + * The main purpose for the existence of this class is its base class, + * which only has a vector as template argument. Therefore, this + * interface provides an abstract base class for matrices. + * + * This class differs form PointerMatrix by its additional + * VectorMemory object and by the fact that it implements the + * functions vmult_add() and Tvmult_add() only using vmult() and + * Tvmult() of the MATRIX. + * + * @author Guido Kanschat 2006 + */ +template +class PointerMatrixAux : public PointerMatrixBase +{ + public: + /** + * Constructor. The pointer in the + * argument is stored in this + * class. As usual, the lifetime of + * *M must be longer than the + * one of the PointerMatrixAux. + * + * If M is zero, no + * matrix is stored. + */ + PointerMatrixAux (VectorMemory* mem, + const MATRIX* M=0); + + /** + * Constructor. The name argument + * is used to identify the + * SmartPointer for this object. + */ + PointerMatrixAux(VectorMemory* mem, + const char* name); + + /** + * Constructor. M points + * to a matrix which must live + * longer than the + * PointerMatrixAux. The name + * argument is used to identify + * the SmartPointer for this + * object. + */ + PointerMatrixAux(VectorMemory* mem, + const MATRIX* M, + const char* name); + + // Use doc from base class + virtual void clear(); + + /** + * Return whether the object is + * empty. + */ + bool empty () const; + + /** + * Assign a new matrix + * pointer. Deletes the old pointer + * and releases its matrix. + * @see SmartPointer + */ + const PointerMatrixAux& operator= (const MATRIX* M); + + /** + * Matrix-vector product. + */ + virtual void vmult (VECTOR& dst, + const VECTOR& src) const; + + /** + * Tranposed matrix-vector product. + */ + virtual void Tvmult (VECTOR& dst, + const VECTOR& src) const; + + /** + * Matrix-vector product, adding to + * dst. + */ + virtual void vmult_add (VECTOR& dst, + const VECTOR& src) const; + + /** + * Tranposed matrix-vector product, + * adding to dst. + */ + virtual void Tvmult_add (VECTOR& dst, + const VECTOR& src) const; + + private: + /** + * Return the address of the + * matrix for comparison. + */ + virtual const void* get() const; + + /** + * Object for getting the + * auxiliary vector. + */ + SmartPointer > mem; + + /** + * The pointer to the actual matrix. + */ + SmartPointer m; +}; + + /** * This function helps you creating a PointerMatrixBase object if you * do not want to provide the full template arguments of @@ -418,4 +536,117 @@ PointerMatrix::get () const } +//----------------------------------------------------------------------// + + +template +PointerMatrixAux::PointerMatrixAux ( + VectorMemory* mem, + const MATRIX* M) + : mem(mem, typeid(*this).name()), + m(M, typeid(*this).name()) +{} + + +template +PointerMatrixAux::PointerMatrixAux ( + VectorMemory* mem, + const char* name) + : mem(mem, name), + m(0, name) +{} + + +template +PointerMatrixAux::PointerMatrixAux ( + VectorMemory* mem, + const MATRIX* M, + const char* name) + : mem(mem, name), + m(M, name) +{} + + +template +inline void +PointerMatrixAux::clear () +{ + m = 0; +} + + +template +inline const PointerMatrixAux& +PointerMatrixAux::operator= (const MATRIX* M) +{ + m = M; + return *this; +} + + +template +inline bool +PointerMatrixAux::empty () const +{ + if (m == 0) + return true; + return m->empty(); +} + +template +inline void +PointerMatrixAux::vmult (VECTOR& dst, + const VECTOR& src) const +{ + Assert (m != 0, ExcNotInitialized()); + m->vmult (dst, src); +} + + +template +inline void +PointerMatrixAux::Tvmult (VECTOR& dst, + const VECTOR& src) const +{ + Assert (m != 0, ExcNotInitialized()); + m->Tvmult (dst, src); +} + + +template +inline void +PointerMatrixAux::vmult_add (VECTOR& dst, + const VECTOR& src) const +{ + Assert (m != 0, ExcNotInitialized()); + VECTOR* v = mem->alloc(); + v->reinit(dst); + m->vmult (*v, src); + dst += *v; + mem->free(v); +} + + +template +inline void +PointerMatrixAux::Tvmult_add (VECTOR& dst, + const VECTOR& src) const +{ + Assert (m != 0, ExcNotInitialized()); + VECTOR* v = mem->alloc(); + v->reinit(dst); + m->Tvmult (*v, src); + dst += *v; + mem->free(v); +} + + +template +inline const void* +PointerMatrixAux::get () const +{ + return m; +} + + #endif