From: guido Date: Thu, 3 Oct 2002 08:30:57 +0000 (+0000) Subject: begin new BlockSparseMatrixEZ X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=4f7900f7713e344dc21da4e9ffa82910181cf5e0;p=dealii-svn.git begin new BlockSparseMatrixEZ git-svn-id: https://svn.dealii.org/trunk@6598 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/block_sparse_matrix_ez.h b/deal.II/lac/include/lac/block_sparse_matrix_ez.h new file mode 100644 index 0000000000..5caa44322a --- /dev/null +++ b/deal.II/lac/include/lac/block_sparse_matrix_ez.h @@ -0,0 +1,485 @@ +//--------------------------------------------------------------------------- +// $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. +// +//--------------------------------------------------------------------------- +#ifndef __deal2__block_sparse_matrix_ez_h +#define __deal2__block_sparse_matrix_ez_h + + +#include +#include +#include +#include +#include + +template class BlockVector; + +/** + * A block matrix consisting of blocks of type @p{SparseMatrixEZ}. + * + * Like the other Block-objects, this matrix can be used like a + * @p{SparseMatrixEZ}, when it comes to access to entries. Then, there + * are functions for the multiplication with @p{BlockVector} and + * access to the individual blocks. + * + * @author Guido Kanschat, 2002 + */ +template +class BlockSparseMatrixEZ : public Subscriptor +{ + public: + /** + * Default constructor. The + * result is an empty object with + * zero dimensions. + */ + BlockSparseMatrixEZ (); + + /** + * Constructor setting up an + * object with given unmber of + * block rows and columns. The + * blocks themselves still have + * zero dimension. + */ + BlockSparseMatrixEZ (const unsigned int block_rows, + const unsigned int block_cols); + + /** + * Copy constructor. This is + * needed for some container + * classes. It creates an object + * of the same number of block + * rows and columns. Since it + * calls the copy constructor of + * @ref{SparseMatrixEZ}, the + * block s must be empty. + */ + BlockSparseMatrixEZ (const BlockSparseMatrixEZ&); + + /** + * Copy operator. Like the copy + * constructor, this may be + * called for objects with empty + * blocks only. + */ + BlockSparseMatrixEZ operator = (const BlockSparseMatrixEZ&); + + /** + * Initialize to given block + * numbers. After this + * operation, the matrix will + * have the block dimensions + * provided. Each block will have + * zero dimensions and must be + * initialized + * subsequently. After setting + * the sizes of the blocks, + * @ref{collect_sizes} must be + * called to update internal data + * structures. + */ + void reinit (const unsigned int n_block_rows, + const unsigned int n_block_cols); + /** + * This function collects the + * sizes of the sub-objects and + * stores them in internal + * arrays, in order to be able to + * relay global indices into the + * matrix to indices into the + * subobjects. You *must* call + * this function each time after + * you have changed the size of + * the sub-objects. + */ + void collect_sizes (); + + /** + * Access the block with the + * given coordinates. + */ + SparseMatrixEZ& + block (const unsigned int row, + const unsigned int column); + + + /** + * Access the block with the + * given coordinates. Version for + * constant objects. + */ + const SparseMatrixEZ& + block (const unsigned int row, + const unsigned int column) const; + + /** + * Return the number of blocks in a + * column. + */ + unsigned int n_block_rows () const; + + /** + * Return the number of blocks in a + * row. + */ + unsigned int n_block_cols () const; + + /** + * Return whether the object is + * empty. It is empty if no + * memory is allocated, which is + * the same as that both + * dimensions are zero. This + * function is just the + * concatenation of the + * respective call to all + * sub-matrices. + */ + bool empty () const; + + /** + * Return number of rows of this + * matrix, which equals the + * dimension of the image + * space. It is the sum of rows + * of the rows of sub-matrices. + */ + unsigned int n_rows () const; + + /** + * Return number of columns of + * this matrix, which equals the + * dimension of the range + * space. It is the sum of + * columns of the columns of + * sub-matrices. + */ + unsigned int n_cols () const; + + /** + * Return the dimension of the + * image space. To remember: the + * matrix is of dimension + * $m \times n$. + */ + unsigned int m () const; + + /** + * Return the dimension of the + * range space. To remember: the + * matrix is of dimension + * $m \times n$. + */ + unsigned int n () const; + + /** + * Set the element @p{(i,j)} to + * @p{value}. Throws an error if + * the entry does not + * exist. Still, it is allowed to + * store zero values in + * non-existent fields. + */ + void set (const unsigned int i, + const unsigned int j, + const number value); + + /** + * Add @p{value} to the element + * @p{(i,j)}. Throws an error if + * the entry does not + * exist. Still, it is allowed to + * store zero values in + * non-existent fields. + */ + void add (const unsigned int i, const unsigned int j, + const number value); + + + /** + * Matrix-vector multiplication: + * let $dst = M*src$ with $M$ + * being this matrix. + */ + template + void vmult (BlockVector &dst, + const BlockVector &src) const; + + /** + * Matrix-vector multiplication: + * let $dst = M^T*src$ with $M$ + * being this matrix. This + * function does the same as + * @p{vmult} but takes the + * transposed matrix. + */ + template + void Tvmult (BlockVector &dst, + const BlockVector &src) const; + + /** + * Adding Matrix-vector + * multiplication. Add $M*src$ on + * $dst$ with $M$ being this + * matrix. + */ + template + void vmult_add (BlockVector &dst, + const BlockVector &src) const; + + /** + * Adding Matrix-vector + * multiplication. Add $M^T*src$ + * to $dst$ with $M$ being this + * matrix. This function does the + * same as @p{vmult_add} but takes + * the transposed matrix. + */ + template + void Tvmult_add (BlockVector &dst, + const BlockVector &src) const; + + + private: + /** + * Object storing and managing + * the transformation of row + * indices to indices of the + * sub-objects. + */ + BlockIndices row_indices; + + /** + * Object storing and managing + * the transformation of column + * indices to indices of the + * sub-objects. + */ + BlockIndices column_indices; + + /** + * The actual matrices + */ + Table<2, SparseMatrixEZ > blocks; +}; + + +/*----------------------------------------------------------------------*/ + + +template +inline +unsigned int +BlockSparseMatrixEZ::n_block_rows () const +{ + return row_indices.size() +} + + + +template +inline +unsigned int +BlockSparseMatrixEZ::n_rows () const +{ + return row_indices.total_size() +} + + + +template +inline +unsigned int +BlockSparseMatrixEZ::n_block_cols () const +{ + return column_indices.size() +} + + + +template +inline +unsigned int +BlockSparseMatrixEZ::n_cols () const +{ + return column_indices.total_size() +} + + + +template +inline +SparseMatrixEZ & +BlockSparseMatrixEZ::block (const unsigned int row, + const unsigned int column) +{ + Assert (row +inline +const SparseMatrixEZ & +BlockSparseMatrixEZ::block (const unsigned int row, + const unsigned int column) const +{ + Assert (row +inline +unsigned int +BlockSparseMatrixEZ::m () const +{ + return n_rows(); +}; + + + +template +inline +unsigned int +BlockSparseMatrixEZ::n () const +{ + return n_cols(); +}; + + + +template +inline +void +BlockSparseMatrixEZ::set (const unsigned int i, + const unsigned int j, + const number value) +{ + const std::pair + row_index = row_indices.global_to_local (i), + col_index = column_indices.global_to_local (j); + block(row_index.first,col_index.first).set (row_index.second, + col_index.second, + value); +}; + + + +template +inline +void +BlockSparseMatrixEZ::add (const unsigned int i, + const unsigned int j, + const number value) +{ + const std::pair + row_index = row_indices.global_to_local (i), + col_index = column_indices.global_to_local (j); + block(row_index.first,col_index.first).add (row_index.second, + col_index.second, + value); +}; + + +template +template +void +BlockSparseMatrixEZ::vmult (BlockVector &dst, + const BlockVector &src) const +{ + Assert (dst.n_blocks() == n_block_rows(), + ExcDimensionMismatch(dst.n_blocks(), n_block_rows())); + Assert (src.n_blocks() == n_block_cols(), + ExcDimensionMismatch(src.n_blocks(), n_block_cols())); + + dst = 0.; + + for (unsigned int row=0; row +template +void +BlockSparseMatrix::vmult_add (BlockVector &dst, + const BlockVector &src) const +{ + Assert (dst.n_blocks() == n_block_rows(), + ExcDimensionMismatch(dst.n_blocks(), n_block_rows())); + Assert (src.n_blocks() == n_block_cols(), + ExcDimensionMismatch(src.n_blocks(), n_block_cols())); + + for (unsigned int row=0; row +template +void +BlockSparseMatrix::Tvmult (BlockVector &dst, + const BlockVector& src) const +{ + Assert (dst.n_blocks() == n_block_cols(), + ExcDimensionMismatch(dst.n_blocks(), n_block_cols())); + Assert (src.n_blocks() == n_block_rows(), + ExcDimensionMismatch(src.n_blocks(), n_block_rows())); + + dst = 0.; + + for (unsigned int row=0; row +template +void +BlockSparseMatrix::Tvmult_add (BlockVector &dst, + const BlockVector &src) const +{ + Assert (dst.n_blocks() == n_block_cols(), + ExcDimensionMismatch(dst.n_blocks(), n_block_cols())); + Assert (src.n_blocks() == n_block_rows(), + ExcDimensionMismatch(src.n_blocks(), n_block_rows())); + + for (unsigned int row=0; row +#include +#include + + + +template +BlockSparseMatrixEZ::BlockSparseMatrixEZ () : + row_indices (0), + column_indices (0) +{}; + + + +template +BlockSparseMatrixEZ:: +BlockSparseMatrixEZ (unsigned int rows, unsigned int cols) : + row_indices (rows), + column_indices (cols) +{}; + + + +// template +// BlockSparseMatrixEZ::~BlockSparseMatrixEZ () +// { +// // delete previous content of +// // the subobjects array +// for (unsigned int r=0; r *p = sub_objects[r][c]; +// sub_objects[r][c] = 0; +// delete p; +// }; +// }; + + + +template +BlockSparseMatrixEZ & +BlockSparseMatrixEZ:: +operator = (const BlockSparseMatrixEZ &m) +{ + Assert (n_block_rows() == m.n_block_rows(), ExcIncompatibleObjects()); + Assert (n_block_cols() == m.n_block_cols(), ExcIncompatibleObjects()); + // this operator does not do + // anything except than checking + // whether the base objects want to + // do something + for (unsigned int r=0; r +BlockSparseMatrixEZ::BlockSparseMatrixEZ ( + const BlockSparseMatrixEZ &m) + : + row_indices(m.row_indices), + column_indices(m.column_indices), + blocks(m.blocks) +{}; + + + +template +void +BlockSparseMatrixEZ::reinit (unsigned int rows, + unsigned int cols) +{ + row_indices.reinit(rows, 0); + column_indices.reinit(cols, 0); + blocks.reinit(rows, cols); +}; + + + +template +void +BlockSparseMatrixEZ::clear () +{ + row_indices.reinit(0, 0); + column_indices.reinit(0, 0); + blocks.reinit(0, 0); +}; + + + +template +bool +BlockSparseMatrixEZ::empty () const +{ + for (unsigned int r=0; r +void +BlockSparseMatrixEZ::collect_sizes () +{ + const unsigned int rows = n_block_rows(); + const unsigned int columns = n_block_cols(); + std::vector row_sizes (rows); + std::vector col_sizes (columns); + + // first find out the row sizes + // from the first block column + for (unsigned int r=0; r +// unsigned int +// BlockSparseMatrixEZ::n_nonzero_elements () const +// { +// return sparsity_pattern->n_nonzero_elements (); +// }; + + + +// template +// unsigned int +// BlockSparseMatrixEZ::n_actually_nonzero_elements () const +// { +// unsigned int count = 0; +// for (unsigned int i=0; in_actually_nonzero_elements (); +// return count; +// }; + + + +// template +// unsigned int +// BlockSparseMatrixEZ::memory_consumption () const +// { +// unsigned int mem = sizeof(*this); +// mem += MemoryConsumption::memory_consumption (sub_objects); +// for (unsigned int r=0; r