From: wolf Date: Sat, 2 Dec 2000 21:30:13 +0000 (+0000) Subject: Add functions symmetrize to SparsityPattern and SparseMatrix. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=68f032a1f667c7d2a76a1c7d74b0c06a8f4bc273;p=dealii-svn.git Add functions symmetrize to SparsityPattern and SparseMatrix. git-svn-id: https://svn.dealii.org/trunk@3520 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/2000/c-3-0.html b/deal.II/doc/news/2000/c-3-0.html index 218cac5a10..b7885f62c8 100644 --- a/deal.II/doc/news/2000/c-3-0.html +++ b/deal.II/doc/news/2000/c-3-0.html @@ -250,6 +250,15 @@ documentation, etc.

lac

    +

  1. + New: There are now functions + SparsityPattern::symmetrize and + SparseMatrix::symmetrize that + generate a symmetric matrix from a non-symmetric one. +
    + (WB 2000/12/02) +

    +

  2. New: almost all classes that store data now have a function memory_consumption that returns an diff --git a/deal.II/lac/include/lac/sparse_matrix.h b/deal.II/lac/include/lac/sparse_matrix.h index ffd0d5f5ca..09d9ad660a 100644 --- a/deal.II/lac/include/lac/sparse_matrix.h +++ b/deal.II/lac/include/lac/sparse_matrix.h @@ -254,6 +254,34 @@ class SparseMatrix : public Subscriptor void add (const unsigned int i, const unsigned int j, const number value); + /** + * Symmetrize the matrix by + * forming the mean value between + * the existing matrix and its + * transpose, $A = \frac 12(A+A^T)$. + * + * This operation assumes that + * the underlying sparsity + * pattern represents a symmetric + * object. If this is not the + * case, then the result of this + * operation will not be a + * symmetric matrix, since it + * only explicitely symmetrizes + * by looping over the lower left + * triangular part for efficiency + * reasons; if there are entries + * in the upper right triangle, + * then these elements are missed + * in the + * symmetrization. Symmetrization + * of the sparsity pattern can be + * obtain by the + * @ref{SparsityPattern}@p{::symmetrize} + * function. + */ + void symmetrize (); + /** * Copy the given matrix to this * one. The operation throws an diff --git a/deal.II/lac/include/lac/sparse_matrix.templates.h b/deal.II/lac/include/lac/sparse_matrix.templates.h index c1db80dc55..4bf1714dde 100644 --- a/deal.II/lac/include/lac/sparse_matrix.templates.h +++ b/deal.II/lac/include/lac/sparse_matrix.templates.h @@ -40,6 +40,7 @@ SparseMatrix::SparseMatrix () : {}; + template SparseMatrix::SparseMatrix (const SparseMatrix &m) : Subscriptor (m), @@ -53,6 +54,7 @@ SparseMatrix::SparseMatrix (const SparseMatrix &m) : }; + template SparseMatrix& SparseMatrix::operator = (const SparseMatrix &m) @@ -65,6 +67,7 @@ SparseMatrix::operator = (const SparseMatrix &m) }; + template SparseMatrix::SparseMatrix (const SparsityPattern &c) : cols(&c), @@ -75,6 +78,7 @@ SparseMatrix::SparseMatrix (const SparsityPattern &c) : }; + template SparseMatrix::~SparseMatrix () { @@ -85,6 +89,7 @@ SparseMatrix::~SparseMatrix () }; + template void SparseMatrix::reinit () @@ -113,6 +118,7 @@ SparseMatrix::reinit () } + template void SparseMatrix::reinit (const SparsityPattern &sparsity) @@ -122,6 +128,7 @@ SparseMatrix::reinit (const SparsityPattern &sparsity) }; + template void SparseMatrix::clear () @@ -133,6 +140,7 @@ SparseMatrix::clear () }; + template bool SparseMatrix::empty () const @@ -144,6 +152,7 @@ SparseMatrix::empty () const }; + template unsigned int SparseMatrix::n_nonzero_elements () const @@ -153,6 +162,7 @@ SparseMatrix::n_nonzero_elements () const }; + template unsigned int SparseMatrix::n_actually_nonzero_elements () const @@ -163,6 +173,44 @@ SparseMatrix::n_actually_nonzero_elements () const }; + +template +void +SparseMatrix::symmetrize () +{ + Assert (cols != 0, ExcMatrixNotInitialized()); + Assert (cols->rows == cols->cols, ExcMatrixNotSquare()); + + const unsigned int n_rows = m(); + for (unsigned int row=0; rowrowstart[row]+1]; + const unsigned int *colnum_ptr = &cols->colnums[cols->rowstart[row]+1]; + const number *const val_end_of_row = &val[cols->rowstart[row+1]]; + + // treat lower left triangle + while ((val_ptr != val_end_of_row) && (*colnum_ptr template SparseMatrix & diff --git a/deal.II/lac/include/lac/sparsity_pattern.h b/deal.II/lac/include/lac/sparsity_pattern.h index fbf44987ce..c29fe1110b 100644 --- a/deal.II/lac/include/lac/sparsity_pattern.h +++ b/deal.II/lac/include/lac/sparsity_pattern.h @@ -351,7 +351,8 @@ class SparsityPattern : public Subscriptor * because the column numbers are * sorted. */ - unsigned int operator() (const unsigned int i, const unsigned int j) const; + unsigned int operator() (const unsigned int i, + const unsigned int j) const; /** * Add a nonzero entry to the matrix. @@ -361,7 +362,21 @@ class SparsityPattern : public Subscriptor * If the entry already exists, nothing * bad happens. */ - void add (const unsigned int i, const unsigned int j); + void add (const unsigned int i, + const unsigned int j); + + /** + * Make the sparsity pattern + * symmetric by adding the + * sparsity pattern of the + * transpose object. + * + * This function throws an + * exception if the sparsity + * pattern does not represent a + * square matrix. + */ + void symmetrize (); /** * Print the sparsity of the matrix diff --git a/deal.II/lac/source/sparsity_pattern.cc b/deal.II/lac/source/sparsity_pattern.cc index 81bb2999ee..d875283768 100644 --- a/deal.II/lac/source/sparsity_pattern.cc +++ b/deal.II/lac/source/sparsity_pattern.cc @@ -567,7 +567,42 @@ SparsityPattern::add (const unsigned int i, // wrong: there was not enough space // in this line Assert (false, ExcNotEnoughSpace(i, rowstart[i+1]-rowstart[i])); -} +}; + + + +void +SparsityPattern::symmetrize () +{ + Assert ((rowstart!=0) && (colnums!=0), ExcEmptyObject()); + Assert (compressed==false, ExcMatrixIsCompressed()); + Assert (rows==cols, ExcNotSquare()); + + // loop over all elements presently + // in the sparsity pattern and add + // the transpose element. note: + // + // 1. that the sparsity pattern + // changes which we work on + // + // 2. that the @p{add} function can + // be called on elements that + // already exist without any harm + for (unsigned int row=0; row