From: guido Date: Mon, 9 Feb 2004 10:13:21 +0000 (+0000) Subject: new add function X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a155ab35dcc9d0262ff08bcbff4cfe7994b4ed5d;p=dealii-svn.git new add function git-svn-id: https://svn.dealii.org/trunk@8425 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/doc/news/c-4-0.html b/deal.II/doc/news/c-4-0.html index 90d6625323..bc6dc5cfdb 100644 --- a/deal.II/doc/news/c-4-0.html +++ b/deal.II/doc/news/c-4-0.html @@ -276,6 +276,12 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
    +
  1. New:FullMatrix has a new function add, + allowing to add to a selected block of the matrix. +
    + (GK 2004/02/09) +

    +
  2. Improved: Initialization routines of class SparseMatrix have an additional parameter controlling the storage of diagonal entries. diff --git a/deal.II/lac/include/lac/full_matrix.h b/deal.II/lac/include/lac/full_matrix.h index 71cb8d0295..416efca4d9 100644 --- a/deal.II/lac/include/lac/full_matrix.h +++ b/deal.II/lac/include/lac/full_matrix.h @@ -52,7 +52,7 @@ template class Vector; * * @ref Instantiations: some (@ @) * - * @author Guido Kanschat, Franz-Theo Suttmeier, Wolfgang Bangerth, 1993-2001 + * @author Guido Kanschat, Franz-Theo Suttmeier, Wolfgang Bangerth, 1993-2004 */ template class FullMatrix : public Table<2,number> @@ -396,6 +396,32 @@ class FullMatrix : public Table<2,number> void add (const number s, const FullMatrix &B); + /** + * Add rectangular block. + * + * A rectangular block of the + * matrix src is added to + * this. The upper left + * corner of the block being + * copied is + * (src_offset_i,src_offset_j). + * The upper left corner of the + * copied block is + * (dst_offset_i,dst_offset_j). + * The size of the rectangular + * block being copied is the + * maximum size possible, + * determined either by the size + * of this or src. + */ + template + void add (const FullMatrix &src, + const double factor, + const unsigned int dst_offset_i = 0, + const unsigned int dst_offset_j = 0, + const unsigned int src_offset_i = 0, + const unsigned int src_offset_j = 0); + /** * Weighted addition of the * transpose of B to this. diff --git a/deal.II/lac/include/lac/full_matrix.templates.h b/deal.II/lac/include/lac/full_matrix.templates.h index 1df2423dc2..01d7ac2262 100644 --- a/deal.II/lac/include/lac/full_matrix.templates.h +++ b/deal.II/lac/include/lac/full_matrix.templates.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -865,6 +865,30 @@ FullMatrix::add (const number s, } } +template +template +void FullMatrix::add (const FullMatrix &src, + const double factor, + const unsigned int dst_offset_i, + const unsigned int dst_offset_j, + const unsigned int src_offset_i, + const unsigned int src_offset_j) +{ + // Compute maximal size of copied block + const unsigned int rows = (m() - dst_offset_i >= src.m() - src_offset_i) + ? src.m() + : m(); + const unsigned int cols = (n() - dst_offset_j >= src.n() - src_offset_j) + ? src.n() + : n(); + + for (unsigned int i=0; iel(dst_offset_i+i,dst_offset_j+j) + += factor * src.el(src_offset_i+i,src_offset_j+j); +} + + template template diff --git a/deal.II/lac/source/full_matrix.double.cc b/deal.II/lac/source/full_matrix.double.cc index 0923aa128f..b74cc61cb2 100644 --- a/deal.II/lac/source/full_matrix.double.cc +++ b/deal.II/lac/source/full_matrix.double.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -19,83 +19,79 @@ template class FullMatrix; -template FullMatrix & FullMatrix::operator =(const FullMatrix&); +template FullMatrix& FullMatrix::operator =( + const FullMatrix&); #define TYPEMAT2 double template void FullMatrix::fill ( - const FullMatrix&, const unsigned, const unsigned, const unsigned, const unsigned); + const FullMatrix&, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::add (const TYPEMAT, const FullMatrix&); +template void FullMatrix::add ( + const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); template void FullMatrix::mmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::Tmmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::add_diag (const TYPEMAT, const FullMatrix&); template void FullMatrix::invert (const FullMatrix&); - #define TYPEVEC double #define TYPERES double -template void FullMatrix::fill_permutation (const FullMatrix&, - const std::vector&, - const std::vector&); -template void FullMatrix::vmult(Vector&, - const Vector&, - const bool) const; -template void FullMatrix::Tvmult(Vector&, - const Vector&, - - const bool) const; -template double FullMatrix::residual(Vector&, const Vector&, - const Vector&) const; -template TYPEVEC FullMatrix::matrix_norm_square (const Vector &) const; -template TYPEVEC FullMatrix::matrix_scalar_product(const Vector&, - const Vector&) const; -template void FullMatrix::forward(Vector&, - const Vector&) const; -template void FullMatrix::backward(Vector&, - const Vector&) const; +template void FullMatrix::fill_permutation ( + const FullMatrix&, + const std::vector&, + const std::vector&); +template void FullMatrix::vmult( + Vector&, const Vector&, bool) const; +template void FullMatrix::Tvmult( + Vector&, const Vector&, bool) const; +template double FullMatrix::residual( + Vector&, const Vector&, const Vector&) const; +template TYPEVEC FullMatrix::matrix_norm_square ( + const Vector &) const; +template TYPEVEC FullMatrix::matrix_scalar_product( + const Vector&, const Vector&) const; +template void FullMatrix::forward( + Vector&, const Vector&) const; +template void FullMatrix::backward( + Vector&, const Vector&) const; template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares(Vector&, - Vector&); +template double FullMatrix::least_squares( + Vector&, Vector&); template -void FullMatrix::precondition_Jacobi (Vector &, - const Vector &, - const TYPEMAT) const; - +void FullMatrix::precondition_Jacobi ( + Vector &, const Vector &, const TYPEMAT) const; #undef TYPEVEC #define TYPEVEC float -template void FullMatrix::fill_permutation (const FullMatrix&, - const std::vector&, - const std::vector&); -template void FullMatrix::vmult(Vector&, - const Vector&, - const bool) const; -template void FullMatrix::Tvmult(Vector&, - const Vector&, - const bool) const; -template double FullMatrix::residual(Vector&, - const Vector&, - const Vector&) const; -template TYPEVEC FullMatrix::matrix_norm_square (const Vector &) const; -template TYPEVEC FullMatrix::matrix_scalar_product(const Vector&, - const Vector&) const; -template void FullMatrix::forward(Vector&, - const Vector&) const; -template void FullMatrix::backward(Vector&, - const Vector&) const; +template void FullMatrix::fill_permutation ( + const FullMatrix&, + const std::vector&, + const std::vector&); +template void FullMatrix::vmult( + Vector&, const Vector&, bool) const; +template void FullMatrix::Tvmult( + Vector&, const Vector&, bool) const; +template double FullMatrix::residual( + Vector&, const Vector&, const Vector&) const; +template TYPEVEC FullMatrix::matrix_norm_square ( + const Vector &) const; +template TYPEVEC FullMatrix::matrix_scalar_product( + const Vector&, const Vector&) const; +template void FullMatrix::forward( + Vector&, const Vector&) const; +template void FullMatrix::backward( + Vector&, const Vector&) const; template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares(Vector&, - Vector&); - +template double FullMatrix::least_squares( + Vector&, Vector&); template -void -FullMatrix::precondition_Jacobi (Vector &, - const Vector &, - const TYPEMAT) const; +void FullMatrix::precondition_Jacobi ( + Vector &, const Vector &, const TYPEMAT) const; + #undef TYPERES #define TYPERES float diff --git a/deal.II/lac/source/full_matrix.float.cc b/deal.II/lac/source/full_matrix.float.cc index fddd4eb756..c130f26d3b 100644 --- a/deal.II/lac/source/full_matrix.float.cc +++ b/deal.II/lac/source/full_matrix.float.cc @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors +// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -18,14 +18,17 @@ template class FullMatrix; -template FullMatrix & FullMatrix::operator =(const FullMatrix&); +template FullMatrix& FullMatrix::operator =( + const FullMatrix&); #define TYPEMAT2 float //template FullMatrix& FullMatrix::operator =<>(const FullMatrix&); template void FullMatrix::fill ( - const FullMatrix&, const unsigned, const unsigned, const unsigned, const unsigned); + const FullMatrix&, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::add (const TYPEMAT, const FullMatrix&); +template void FullMatrix::add ( + const FullMatrix&, double, unsigned, unsigned, unsigned, unsigned); template void FullMatrix::Tadd (const TYPEMAT, const FullMatrix&); template void FullMatrix::mmult (FullMatrix&, const FullMatrix&, const bool) const; template void FullMatrix::Tmmult (FullMatrix&, const FullMatrix&, const bool) const; @@ -35,45 +38,63 @@ template void FullMatrix::invert (const FullMatrix& #define TYPEVEC double #define TYPERES double -template void FullMatrix::fill_permutation (const FullMatrix&, - const std::vector&, - const std::vector&); -template void FullMatrix::vmult(Vector&, const Vector&, const bool) const; -template void FullMatrix::Tvmult(Vector&, const Vector&, const bool) const; -template double FullMatrix::residual(Vector&, const Vector&, const Vector&) const; -template TYPEVEC FullMatrix::matrix_norm_square (const Vector &) const; -template TYPEVEC FullMatrix::matrix_scalar_product(const Vector&, const Vector&) const; -template void FullMatrix::forward(Vector&, const Vector&) const; -template void FullMatrix::backward(Vector&, const Vector&) const; +template void FullMatrix::fill_permutation ( + const FullMatrix&, + const std::vector&, + const std::vector&); +template void FullMatrix::vmult( + Vector&, const Vector&, bool) const; +template void FullMatrix::Tvmult( + Vector&, const Vector&, bool) const; +template double FullMatrix::residual( + Vector&, const Vector&, const Vector&) const; +template TYPEVEC FullMatrix::matrix_norm_square ( + const Vector &) const; +template TYPEVEC FullMatrix::matrix_scalar_product( + const Vector&, const Vector&) const; +template void FullMatrix::forward( + Vector&, const Vector&) const; +template void FullMatrix::backward( + Vector&, const Vector&) const; template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares(Vector&, Vector&); +template double FullMatrix::least_squares( + Vector&, Vector&); + template -void FullMatrix::precondition_Jacobi (Vector &, - const Vector &, - const TYPEMAT) const; +void FullMatrix::precondition_Jacobi ( + Vector &, const Vector &, const TYPEMAT) const; #undef TYPEVEC #define TYPEVEC float -template void FullMatrix::fill_permutation (const FullMatrix&, - const std::vector&, - const std::vector&); -template void FullMatrix::vmult(Vector&, const Vector&, const bool) const; -template void FullMatrix::Tvmult(Vector&, const Vector&, const bool) const; -template double FullMatrix::residual(Vector&, const Vector&, const Vector&) const; -template TYPEVEC FullMatrix::matrix_norm_square (const Vector &) const; -template TYPEVEC FullMatrix::matrix_scalar_product(const Vector&, const Vector&) const; -template void FullMatrix::forward(Vector&, const Vector&) const; -template void FullMatrix::backward(Vector&, const Vector&) const; +template void FullMatrix::fill_permutation ( + const FullMatrix&, + const std::vector&, + const std::vector&); +template void FullMatrix::vmult( + Vector&, const Vector&, bool) const; +template void FullMatrix::Tvmult( + Vector&, const Vector&, bool) const; +template double FullMatrix::residual( + Vector&, const Vector&, const Vector&) const; +template TYPEVEC FullMatrix::matrix_norm_square ( + const Vector &) const; +template TYPEVEC FullMatrix::matrix_scalar_product( + const Vector&, const Vector&) const; +template void FullMatrix::forward( + Vector&, const Vector&) const; +template void FullMatrix::backward( + Vector&, const Vector&) const; template void FullMatrix::householder(Vector&); -template double FullMatrix::least_squares(Vector&, Vector&); +template double FullMatrix::least_squares( + Vector&, Vector&); template -void FullMatrix::precondition_Jacobi (Vector &, - const Vector &, - const TYPEMAT) const; +void FullMatrix::precondition_Jacobi ( + Vector &, const Vector &, const TYPEMAT) const; #undef TYPERES #define TYPERES float -template double FullMatrix::residual(Vector&, const Vector&, const Vector&) const; +template double FullMatrix::residual( + Vector&, const Vector&, const Vector&) const;