From fc28a5d61313534564219b509c817eecd21a1470 Mon Sep 17 00:00:00 2001 From: guido Date: Tue, 13 Aug 2002 16:01:43 +0000 Subject: [PATCH] matrix library started git-svn-id: https://svn.dealii.org/trunk@6323 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/matrix_lib.h | 245 ++++++++++++++++++ .../lac/include/lac/matrix_lib.templates.h | 91 +++++++ deal.II/lac/include/lac/pointer_matrix.h | 1 + deal.II/lac/source/matrix_lib.cc | 39 +++ 4 files changed, 376 insertions(+) create mode 100644 deal.II/lac/include/lac/matrix_lib.h create mode 100644 deal.II/lac/include/lac/matrix_lib.templates.h create mode 100644 deal.II/lac/source/matrix_lib.cc diff --git a/deal.II/lac/include/lac/matrix_lib.h b/deal.II/lac/include/lac/matrix_lib.h new file mode 100644 index 0000000000..4aa5d48c10 --- /dev/null +++ b/deal.II/lac/include/lac/matrix_lib.h @@ -0,0 +1,245 @@ +// $Id$ + +#ifndef __deal2__matrix_lib_h +#define __deal2__matrix_lib_h + +#include +#include +#include + +template class Vector; +template class BlockVector; + +/** + * Poor man's matrix product. Stores two matrices $m_1$ and $m_2$ and + * implements matrix-vector multiplications for the product $m_1 m_2$ + * by performing multiplication with both factors consecutively. + * + * @author Guido Kanschat, 2000, 2001, 2002 + */ +template +class ProductMatrix : public PointerMatrixBase +{ + public: + /** + * Constructor. Additionally to + * the two constituting matrices, a + * memory pool for the auxiliary + * vector must be provided. + */ + template + ProductMatrix(const MATRIX1& m1, + const MATRIX2& m2, + VectorMemory& mem); + + /** + * Destructor. + */ + ~ProductMatrix(); + + /** + * 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 + * @p{dst}. + */ + virtual void vmult_add (VECTOR& dst, + const VECTOR& src) const; + + /** + * Tranposed matrix-vector product, + * adding to @p{dst}. + */ + virtual void Tvmult_add (VECTOR& dst, + const VECTOR& src) const; + + private: + /** + * The left matrix of the product. + */ + PointerMatrixBase* m1; + + /** + * The right matrix of the product. + */ + PointerMatrixBase* m2; + + /** + * Memory for auxiliary vector. + */ + SmartPointer > mem; +}; + + +/** + * Mean value filter. The @p{vmult} functions of this matrix filter + * out mean values of the vector. If the vector is of type + * @p{BlockVector}, then an additional parameter selects a single + * component for this operation. + * + * @author Guido Kanschat, 2002 + */ +class MeanValueFilter : public Subscriptor +{ + public: + /** + * Constructor, optionally + * selecting a component. + */ + MeanValueFilter(unsigned int component = static_cast(-1)); + + /** + * Return the source vector with + * subtracted mean value. + */ + template + void vmult (Vector& dst, + const Vector& src) const; + + /** + * Add source vector with + * subtracted mean value to dest. + */ + template + void vmult_add (Vector& dst, + const Vector& src) const; + + /** + * Return the source vector with + * subtracted mean value in + * selected component. + */ + template + void vmult (BlockVector& dst, + const BlockVector& src) const; + + /** + * Add a soruce to dest, where + * the mean value in the selected + * component is subtracted. + */ + template + void vmult_add (BlockVector& dst, + const BlockVector& src) const; + + + /** + * Not implemented. + */ + template + void Tvmult(VECTOR&, const VECTOR&) const; + + /** + * Not implemented. + */ + template + void Tvmult_add(VECTOR&, const VECTOR&) const; + + private: + /** + * Component for filtering block vectors. + */ + unsigned int component; +}; + + +//----------------------------------------------------------------------// + +template +template +ProductMatrix::ProductMatrix ( + const MATRIX1& mat1, + const MATRIX2& mat2, + VectorMemory& m) + : mem(&m) +{ + m1 = new PointerMatrix(&mat1); + m2 = new PointerMatrix(&mat2); +} + + +template +ProductMatrix::~ProductMatrix () +{ + delete m1; + delete m2; +} + + +template +void +ProductMatrix::vmult (VECTOR& dst, const VECTOR& src) const +{ + VECTOR* v = mem->alloc(); + v->reinit(dst); + m2->vmult (*v, src); + m1->vmult (dst, *v); + mem->free(v); +} + + +template +void +ProductMatrix::vmult_add (VECTOR& dst, const VECTOR& src) const +{ + VECTOR* v = mem->alloc(); + v->reinit(dst); + m2->vmult (*v, src); + m1->vmult_add (dst, *v); + mem->free(v); +} + + +template +void +ProductMatrix::Tvmult (VECTOR& dst, const VECTOR& src) const +{ + VECTOR* v = mem->alloc(); + v->reinit(dst); + m1->Tvmult (*v, src); + m2->Tvmult (dst, *v); + mem->free(v); +} + + +template +void +ProductMatrix::Tvmult_add (VECTOR& dst, const VECTOR& src) const +{ + VECTOR* v = mem->alloc(); + v->reinit(dst); + m1->Tvmult (*v, src); + m2->Tvmult_add (dst, *v); + mem->free(v); +} + + +//----------------------------------------------------------------------// + +template +inline void +MeanValueFilter::Tvmult(VECTOR&, const VECTOR&) const +{ + Assert(false, ExcNotImplemented()); +} + + +template +inline void +MeanValueFilter::Tvmult_add(VECTOR&, const VECTOR&) const +{ + Assert(false, ExcNotImplemented()); +} + + +#endif diff --git a/deal.II/lac/include/lac/matrix_lib.templates.h b/deal.II/lac/include/lac/matrix_lib.templates.h new file mode 100644 index 0000000000..9dd11e20d2 --- /dev/null +++ b/deal.II/lac/include/lac/matrix_lib.templates.h @@ -0,0 +1,91 @@ +//---------------------------- vector.templates.h --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 1998, 1999, 2000, 2001, 2002 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. +// +//---------------------------- vector.templates.h --------------------------- +#ifndef __deal2__matrix_lib_templates_h +#define __deal2__matrix_lib_templates_h + +#include +#include +#include + + +template +void +MeanValueFilter::vmult(Vector& dst, + const Vector& src) const +{ + Assert (dst.size() == src.size(), + ExcDimensionMismatch(dst.size(), src.size())); + + number mean = src.mean_value(); + + for (unsigned int i=0;i +void +MeanValueFilter::vmult_add(Vector& dst, + const Vector& src) const +{ + Assert (dst.size() == src.size(), + ExcDimensionMismatch(dst.size(), src.size())); + + number mean = src.mean_value(); + + for (unsigned int i=0;i +void +MeanValueFilter::vmult(BlockVector& dst, + const BlockVector& src) const +{ + Assert (component != static_cast(-1), + ExcNotInitialized()); + + Assert (dst.n_blocks() == src.n_blocks(), + ExcDimensionMismatch(dst.n_blocks(), src.n_blocks())); + + for (unsigned int i=0;i +void +MeanValueFilter::vmult_add(BlockVector& dst, + const BlockVector& src) const +{ + Assert (component != static_cast(-1), + ExcNotInitialized()); + + Assert (dst.n_blocks() == src.n_blocks(), + ExcDimensionMismatch(dst.n_blocks(), src.n_blocks())); + + for (unsigned int i=0;i +#include #include /** diff --git a/deal.II/lac/source/matrix_lib.cc b/deal.II/lac/source/matrix_lib.cc new file mode 100644 index 0000000000..0422d0dd05 --- /dev/null +++ b/deal.II/lac/source/matrix_lib.cc @@ -0,0 +1,39 @@ +//---------------------------- vector.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 1998, 1999, 2000, 2001, 2002 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. +// +//---------------------------- vector.cc --------------------------- + + +#include + +MeanValueFilter::MeanValueFilter(unsigned int component) + : + component(component) +{} + +template void MeanValueFilter::vmult(Vector&, + const Vector&) const; +template void MeanValueFilter::vmult(Vector&, + const Vector&) const; +template void MeanValueFilter::vmult(BlockVector&, + const BlockVector&) const; +template void MeanValueFilter::vmult(BlockVector&, + const BlockVector&) const; + +template void MeanValueFilter::vmult_add(Vector&, + const Vector&) const; +template void MeanValueFilter::vmult_add(Vector&, + const Vector&) const; +template void MeanValueFilter::vmult_add(BlockVector&, + const BlockVector&) const; +template void MeanValueFilter::vmult_add(BlockVector&, + const BlockVector&) const; + -- 2.39.5