]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
begin new BlockSparseMatrixEZ
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 3 Oct 2002 08:30:57 +0000 (08:30 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 3 Oct 2002 08:30:57 +0000 (08:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@6598 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/block_sparse_matrix_ez.h [new file with mode: 0644]
deal.II/lac/include/lac/block_sparse_matrix_ez.templates.h [new file with mode: 0644]

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 (file)
index 0000000..5caa443
--- /dev/null
@@ -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 <base/config.h>
+#include <base/exceptions.h>
+#include <base/subscriptor.h>
+#include <base/smartpointer.h>
+#include <lac/block_indices.h>
+
+template <typename Number> 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<typename Number>
+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<Number>&);
+
+                                    /**
+                                     * Copy operator. Like the copy
+                                     * constructor, this may be
+                                     * called for objects with empty
+                                     * blocks only.
+                                     */
+    BlockSparseMatrixEZ operator = (const BlockSparseMatrixEZ<Number>&);
+
+                                    /**
+                                     * 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<Number>&
+    block (const unsigned int row,
+          const unsigned int column);
+    
+    
+                                    /**
+                                     * Access the block with the
+                                     * given coordinates. Version for
+                                     * constant objects.
+                                     */
+    const SparseMatrixEZ<Number>&
+    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 <typename somenumber>
+    void vmult (BlockVector<somenumber>       &dst,
+               const BlockVector<somenumber> &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 <typename somenumber>
+    void Tvmult (BlockVector<somenumber>       &dst,
+                const BlockVector<somenumber> &src) const;
+  
+                                    /**
+                                     * Adding Matrix-vector
+                                     * multiplication. Add $M*src$ on
+                                     * $dst$ with $M$ being this
+                                     * matrix.
+                                     */
+    template <typename somenumber>
+    void vmult_add (BlockVector<somenumber>       &dst,
+                   const BlockVector<somenumber> &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 <typename somenumber>
+    void Tvmult_add (BlockVector<somenumber>       &dst,
+                    const BlockVector<somenumber> &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<Number> > blocks;
+};
+
+  
+/*----------------------------------------------------------------------*/
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::n_block_rows () const
+{
+  return row_indices.size()
+}
+
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::n_rows () const
+{
+  return row_indices.total_size()
+}
+
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::n_block_cols () const
+{
+  return column_indices.size()
+}
+
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::n_cols () const
+{
+  return column_indices.total_size()
+}
+
+
+
+template <typename number>
+inline
+SparseMatrixEZ<number> &
+BlockSparseMatrixEZ<number>::block (const unsigned int row,
+                                 const unsigned int column)
+{
+  Assert (row<n_rows(), ExcIndexRange (row, 0, n_rows()));
+  Assert (column<n_cols(), ExcIndexRange (column, 0, n_cols()));
+  
+  return blocks[row][column];
+};
+
+
+
+template <typename number>
+inline
+const SparseMatrixEZ<number> &
+BlockSparseMatrixEZ<number>::block (const unsigned int row,
+                                 const unsigned int column) const
+{
+  Assert (row<n_rows(), ExcIndexRange (row, 0, n_rows()));
+  Assert (column<n_cols(), ExcIndexRange (column, 0, n_cols()));
+  
+  return blocks[row][column];
+};
+
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::m () const
+{
+  return n_rows();
+};
+
+
+
+template <typename number>
+inline
+unsigned int
+BlockSparseMatrixEZ<number>::n () const
+{
+  return n_cols();
+};
+
+
+
+template <typename number>
+inline
+void
+BlockSparseMatrixEZ<number>::set (const unsigned int i,
+                                 const unsigned int j,
+                                 const number value)
+{
+  const std::pair<unsigned int,unsigned int>
+    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 <typename number>
+inline
+void
+BlockSparseMatrixEZ<number>::add (const unsigned int i,
+                                 const unsigned int j,
+                                 const number value)
+{
+  const std::pair<unsigned int,unsigned int>
+    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 <typename number>
+template <typename somenumber>
+void
+BlockSparseMatrixEZ<number>::vmult (BlockVector<somenumber>       &dst,
+                                   const BlockVector<somenumber> &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<n_block_rows(); ++row)
+    {
+      for (unsigned int col=0; col<n_block_cols(); ++col)
+       block(row,col).vmult_add (dst.block(row),
+                                 src.block(col));
+    };
+};
+
+
+
+template <typename number>
+template <typename somenumber>
+void
+BlockSparseMatrix<number>::vmult_add (BlockVector<somenumber>       &dst,
+                                     const BlockVector<somenumber> &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<n_block_rows(); ++row)
+    {
+      for (unsigned int col=0; col<n_block_cols(); ++col)
+       block(row,col).vmult_add (dst.block(row),
+                                 src.block(col));
+    };
+};
+
+
+
+
+template <typename number>
+template <typename somenumber>
+void
+BlockSparseMatrix<number>::Tvmult (BlockVector<somenumber>   &dst,
+                                  const BlockVector<somenumber>& 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<n_block_rows(); ++row)
+    {
+      for (unsigned int col=0; col<n_block_cols(); ++col)
+       block(row,col).Tvmult_add (dst.block(col),
+                                  src.block(row));
+    };
+};
+
+
+
+template <typename number>
+template <typename somenumber>
+void
+BlockSparseMatrix<number>::Tvmult_add (BlockVector<somenumber>    &dst,
+                                      const BlockVector<somenumber> &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<n_block_rows(); ++row)
+    {
+      for (unsigned int col=0; col<n_block_cols(); ++col)
+       block(row,col).Tvmult_add (dst.block(col),
+                                  src.block(row));
+    };
+};
diff --git a/deal.II/lac/include/lac/block_sparse_matrix_ez.templates.h b/deal.II/lac/include/lac/block_sparse_matrix_ez.templates.h
new file mode 100644 (file)
index 0000000..b2fbcdb
--- /dev/null
@@ -0,0 +1,200 @@
+//-------------------------------------------------------------------------------
+//    $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_templates_h
+#define __deal2__block_sparse_matrix_ez_templates_h
+
+
+#include <base/config.h>
+#include <base/memory_consumption.h>
+#include <lac/block_sparse_matrix_ez.h>
+
+
+
+template <typename number>
+BlockSparseMatrixEZ<number>::BlockSparseMatrixEZ () :
+               row_indices (0),
+               column_indices (0)
+{};
+
+
+
+template <typename number>
+BlockSparseMatrixEZ<number>::
+BlockSparseMatrixEZ (unsigned int rows, unsigned int cols) :
+               row_indices (rows),
+               column_indices (cols)
+{};
+
+
+
+//  template <typename number>
+//  BlockSparseMatrixEZ<number>::~BlockSparseMatrixEZ ()
+//  {
+//                                // delete previous content of
+//                                // the subobjects array
+//    for (unsigned int r=0; r<rows; ++r)
+//      for (unsigned int c=0; c<columns; ++c)
+//        {
+//     SparseMatrixEZ<number> *p = sub_objects[r][c];
+//     sub_objects[r][c] = 0;
+//     delete p;
+//        };
+//  };
+
+
+
+template <typename number>
+BlockSparseMatrixEZ<number> &
+BlockSparseMatrixEZ<number>::
+operator = (const BlockSparseMatrixEZ<number> &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<n_block_rows(); ++r)
+    for (unsigned int c=0; c<n_block_cols(); ++c)
+      block(r,c) = m.block(r,c);
+  return *this;
+};
+
+
+template <typename number>
+BlockSparseMatrixEZ<number>::BlockSparseMatrixEZ (
+  const BlockSparseMatrixEZ<number> &m)
+               :
+               row_indices(m.row_indices),
+               column_indices(m.column_indices),
+               blocks(m.blocks)
+{};
+
+
+template <typename number>
+void
+BlockSparseMatrixEZ<number>::reinit (unsigned int rows,
+                                    unsigned int cols)
+{
+  row_indices.reinit(rows, 0);
+  column_indices.reinit(cols, 0);
+  blocks.reinit(rows, cols);  
+};
+
+
+
+template <typename number>
+void
+BlockSparseMatrixEZ<number>::clear () 
+{
+  row_indices.reinit(0, 0);
+  column_indices.reinit(0, 0);
+  blocks.reinit(0, 0);  
+};
+
+
+
+template <typename number>
+bool
+BlockSparseMatrixEZ<number>::empty () const
+{
+  for (unsigned int r=0; r<n_block_rows(); ++r)
+    for (unsigned int c=0; c<n_block_cols(); ++c)
+      if (block(r,c).empty () == false)
+       return false;
+  return true;
+};
+
+
+
+template <typename number>
+void
+BlockSparseMatrixEZ<number>::collect_sizes ()
+{
+  const unsigned int rows = n_block_rows();
+  const unsigned int columns = n_block_cols();
+  std::vector<unsigned int> row_sizes (rows);
+  std::vector<unsigned int> col_sizes (columns);
+
+                                  // first find out the row sizes
+                                  // from the first block column
+  for (unsigned int r=0; r<rows; ++r)
+    row_sizes[r] = blocks[r][0].n_rows();
+                                  // then check that the following
+                                  // block columns have the same
+                                  // sizes
+  for (unsigned int c=1; c<columns; ++c)
+    for (unsigned int r=0; r<rows; ++r)
+      Assert (row_sizes[r] == blocks[r][c].n_rows(),
+             ExcIncompatibleRowNumbers (r,0,r,c));
+
+                                  // finally initialize the row
+                                  // indices with this array
+  row_indices.reinit (row_sizes);
+  
+  
+                                  // then do the same with the columns
+  for (unsigned int c=0; c<columns; ++c)
+    col_sizes[c] = blocks[0][c].n_cols();
+  for (unsigned int r=1; r<rows; ++r)
+    for (unsigned int c=0; c<columns; ++c)
+      Assert (col_sizes[c] == blocks[r][c].n_cols(),
+             ExcIncompatibleRowNumbers (0,c,r,c));
+
+                                  // finally initialize the row
+                                  // indices with this array
+  column_indices.reinit (col_sizes);
+};
+
+
+
+
+//  template <typename number>
+//  unsigned int
+//  BlockSparseMatrixEZ<number>::n_nonzero_elements () const
+//  {
+//    return sparsity_pattern->n_nonzero_elements ();
+//  };
+
+
+
+//  template <typename number>
+//  unsigned int
+//  BlockSparseMatrixEZ<number>::n_actually_nonzero_elements () const
+//  {
+//    unsigned int count = 0;
+//    for (unsigned int i=0; i<rows; ++i)
+//      for (unsigned int j=0; j<columns; ++j)
+//        count += sub_objects[i][j]->n_actually_nonzero_elements ();
+//    return count;
+//  };
+
+
+
+//  template <typename number>
+//  unsigned int
+//  BlockSparseMatrixEZ<number>::memory_consumption () const
+//  {
+//    unsigned int mem = sizeof(*this);
+//    mem += MemoryConsumption::memory_consumption (sub_objects);
+//    for (unsigned int r=0; r<rows; ++r)
+//      for (unsigned int c=0; c<columns; ++c)
+//        mem += MemoryConsumption::memory_consumption(*sub_objects[r][c]);
+//    return mem;
+//  };
+
+
+
+#endif // ifdef block_sparse_matrix_templates_h

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.