From: wolf Date: Fri, 21 Sep 2001 06:07:35 +0000 (+0000) Subject: FullMatrix::symmetrize X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=3877a575813299ed1428535581bfb211b25ef316;p=dealii-svn.git FullMatrix::symmetrize git-svn-id: https://svn.dealii.org/trunk@5030 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/2001/c-3-1.html b/deal.II/doc/news/2001/c-3-1.html index a1c704260c..2bae941a21 100644 --- a/deal.II/doc/news/2001/c-3-1.html +++ b/deal.II/doc/news/2001/c-3-1.html @@ -356,7 +356,12 @@ documentation, etc.

  1. + New: Function FullMatrix::symmetrize(). +
    + (WB 2001/09/20) +

    +

  2. Improved: the stopping criterion of SolverBicgstab without computing the exact residual is implemented diff --git a/deal.II/lac/include/lac/full_matrix.h b/deal.II/lac/include/lac/full_matrix.h index 40a62d9379..de5ac4995b 100644 --- a/deal.II/lac/include/lac/full_matrix.h +++ b/deal.II/lac/include/lac/full_matrix.h @@ -318,6 +318,17 @@ class FullMatrix : public vector2d number2 matrix_scalar_product (const Vector &u, const Vector &v) const; + /** + * Symmetrize the matrix by + * forming the mean value between + * the existing matrix and its + * transpose, $A = \frac 12(A+A^T)$. + * + * Obviously the matrix must be + * square for this operation. + */ + void symmetrize (); + /** * Return the $l_1$-norm of the matrix, i.e. * $|M|_1=max_{all columns j}\sum_{all @@ -374,13 +385,13 @@ class FullMatrix : public vector2d /** * A=Inverse(A). Inversion of - * this matrix by - * Gauss-Jordan algorithm with - * partial pivoting. This - * process is well-behaved for - * positive definite matrices, - * but be aware of round-off - * errors in the indefinite case. + * this matrix by Gauss-Jordan + * algorithm with partial + * pivoting. This process is + * well-behaved for positive + * definite matrices, but be + * aware of round-off errors in + * the indefinite case. * * The numerical effort to invert * an @p{n x n} matrix is of the @@ -389,12 +400,14 @@ class FullMatrix : public vector2d void gauss_jordan (); /** - * Computes the determinant of a matrix. - * This is only implemented for one two and - * three dimensions, since for higher - * dimensions the numerical work explodes. - * Obviously, the matrix needs to be square - * for this function. + * Computes the determinant of a + * matrix. This is only + * implemented for one, two, and + * three dimensions, since for + * higher dimensions the + * numerical work explodes. + * Obviously, the matrix needs to + * be square for this function. */ double determinant () const; @@ -482,15 +495,17 @@ class FullMatrix : public vector2d const Vector& b) const; /** - * Forward elimination of lower triangle. - * Inverts the lower triangle of a - * quadratic matrix. + * Forward elimination of lower + * triangle. Inverts the lower + * triangle of a quadratic + * matrix. * - * If the matrix has more columns than rows, - * this function only operates on the left - * square submatrix. If there are more rows, - * the upper square part of the matrix - * is considered + * If the matrix has more columns + * than rows, this function only + * operates on the left square + * submatrix. If there are more + * rows, the upper square part of + * the matrix is considered */ template void forward (Vector &dst, @@ -506,12 +521,15 @@ class FullMatrix : public vector2d /** * QR-factorization of a matrix. - * The orthogonal transformation Q is - * applied to the vector y and this matrix. + * The orthogonal transformation + * Q is applied to the vector y + * and this matrix. * - * After execution of householder, the upper - * triangle contains the resulting matrix R, - * the lower the incomplete factorization + * After execution of + * householder, the upper + * triangle contains the + * resulting matrix R, the lower + * the incomplete factorization * matrices. */ template diff --git a/deal.II/lac/include/lac/full_matrix.templates.h b/deal.II/lac/include/lac/full_matrix.templates.h index 0a0494c8f1..0fd628dad9 100644 --- a/deal.II/lac/include/lac/full_matrix.templates.h +++ b/deal.II/lac/include/lac/full_matrix.templates.h @@ -292,15 +292,16 @@ double FullMatrix::residual (Vector& dst, } + template template -void FullMatrix::forward (Vector& dst, - const Vector& src) const +void FullMatrix::forward (Vector &dst, + const Vector &src) const { Assert (data() != 0, ExcEmptyMatrix()); - Assert(dst.size() == m(), ExcDimensionMismatch(dst.size(), m())); - Assert(src.size() == n(), ExcDimensionMismatch(src.size(), n())); + Assert (dst.size() == m(), ExcDimensionMismatch(dst.size(), m())); + Assert (src.size() == n(), ExcDimensionMismatch(src.size(), n())); unsigned int i,j; unsigned int nu = ( (m()::forward (Vector& dst, } + template template -void FullMatrix::backward (Vector& dst, - const Vector& src) const +void FullMatrix::backward (Vector &dst, + const Vector &src) const { Assert (data() != 0, ExcEmptyMatrix()); @@ -333,20 +335,12 @@ void FullMatrix::backward (Vector& dst, } -/* template */ -/* template */ -/* FullMatrix& */ -/* FullMatrix::operator = (const FullMatrix& m) */ -/* { */ - -/* } */ - template template -void FullMatrix::fill (const FullMatrix& src, - const unsigned int i, - const unsigned int j) +void FullMatrix::fill (const FullMatrix &src, + const unsigned int i, + const unsigned int j) { Assert (n() >= src.n() + j, ExcInvalidDestination(n(), src.n(), j)); Assert (m() >= src.m() + i, ExcInvalidDestination(m(), src.m(), i)); @@ -598,6 +592,23 @@ number2 FullMatrix::matrix_scalar_product (const Vector &u, }; + +template +void +FullMatrix::symmetrize () +{ + Assert (m() == n(), ExcNotQuadratic()); + + const unsigned int N = m(); + for (unsigned int i=0; i number FullMatrix::l1_norm () const {