<ol>
+ <li> <p> New:<code class="class">FullMatrix</code> has a new function add,
+ allowing to add to a selected block of the matrix.
+ <br>
+ (GK 2004/02/09)
+ </p>
+
<li> <p> Improved: Initialization routines of class <code
class="class">SparseMatrix</code> have an additional parameter
controlling the storage of diagonal entries.
*
* @ref Instantiations: some (<tt>@<float@> @<double@></tt>)
*
- * @author Guido Kanschat, Franz-Theo Suttmeier, Wolfgang Bangerth, 1993-2001
+ * @author Guido Kanschat, Franz-Theo Suttmeier, Wolfgang Bangerth, 1993-2004
*/
template<typename number>
class FullMatrix : public Table<2,number>
void add (const number s,
const FullMatrix<number2> &B);
+ /**
+ * Add rectangular block.
+ *
+ * A rectangular block of the
+ * matrix <tt>src</tt> is added to
+ * <tt>this</tt>. The upper left
+ * corner of the block being
+ * copied is
+ * <tt>(src_offset_i,src_offset_j)</tt>.
+ * The upper left corner of the
+ * copied block is
+ * <tt>(dst_offset_i,dst_offset_j)</tt>.
+ * The size of the rectangular
+ * block being copied is the
+ * maximum size possible,
+ * determined either by the size
+ * of <tt>this</tt> or <tt>src</tt>.
+ */
+ template<typename number2>
+ void add (const FullMatrix<number2> &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 <tt>B</tt> to <tt>this</tt>.
// $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
}
}
+template <typename number>
+template <typename number2>
+void FullMatrix<number>::add (const FullMatrix<number2> &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; i<rows ; ++i)
+ for (unsigned int j=0; j<cols ; ++j)
+ this->el(dst_offset_i+i,dst_offset_j+j)
+ += factor * src.el(src_offset_i+i,src_offset_j+j);
+}
+
+
template <typename number>
template <typename number2>
// $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
template class FullMatrix<TYPEMAT>;
-template FullMatrix<double> & FullMatrix<double>::operator =(const FullMatrix<float>&);
+template FullMatrix<TYPEMAT>& FullMatrix<TYPEMAT>::operator =(
+ const FullMatrix<float>&);
#define TYPEMAT2 double
template void FullMatrix<TYPEMAT>::fill<TYPEMAT2> (
- const FullMatrix<TYPEMAT2>&, const unsigned, const unsigned, const unsigned, const unsigned);
+ const FullMatrix<TYPEMAT2>&, unsigned, unsigned, unsigned, unsigned);
template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (
+ const FullMatrix<TYPEMAT2>&, double, unsigned, unsigned, unsigned, unsigned);
template void FullMatrix<TYPEMAT>::Tadd<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
template void FullMatrix<TYPEMAT>::mmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
template void FullMatrix<TYPEMAT>::Tmmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
template void FullMatrix<TYPEMAT>::add_diag<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
template void FullMatrix<TYPEMAT>::invert<TYPEMAT2> (const FullMatrix<TYPEMAT2>&);
-
#define TYPEVEC double
#define TYPERES double
-template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (const FullMatrix<TYPEVEC>&,
- const std::vector<unsigned int>&,
- const std::vector<unsigned int>&);
-template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&,
- const bool) const;
-template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&,
-
- const bool) const;
-template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&,
- const Vector<TYPERES>&) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (const Vector<TYPEVEC> &) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(const Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+ const FullMatrix<TYPEVEC>&,
+ const std::vector<unsigned int>&,
+ const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+ const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+ const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
template void FullMatrix<TYPEMAT>::householder<TYPEVEC>(Vector<TYPEVEC>&);
-template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(Vector<TYPEVEC>&,
- Vector<TYPEVEC>&);
+template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(
+ Vector<TYPEVEC>&, Vector<TYPEVEC>&);
template
-void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (Vector<TYPEVEC> &,
- const Vector<TYPEVEC> &,
- const TYPEMAT) const;
-
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+ Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
#undef TYPEVEC
#define TYPEVEC float
-template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (const FullMatrix<TYPEVEC>&,
- const std::vector<unsigned int>&,
- const std::vector<unsigned int>&);
-template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&,
- const bool) const;
-template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&,
- const bool) const;
-template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&,
- const Vector<TYPERES>&) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (const Vector<TYPEVEC> &) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(const Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(Vector<TYPEVEC>&,
- const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+ const FullMatrix<TYPEVEC>&,
+ const std::vector<unsigned int>&,
+ const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+ const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+ const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
template void FullMatrix<TYPEMAT>::householder<TYPEVEC>(Vector<TYPEVEC>&);
-template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(Vector<TYPEVEC>&,
- Vector<TYPEVEC>&);
-
+template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(
+ Vector<TYPEVEC>&, Vector<TYPEVEC>&);
template
-void
-FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (Vector<TYPEVEC> &,
- const Vector<TYPEVEC> &,
- const TYPEMAT) const;
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+ Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
+
#undef TYPERES
#define TYPERES float
// $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
template class FullMatrix<TYPEMAT>;
-template FullMatrix<float> & FullMatrix<float>::operator =(const FullMatrix<double>&);
+template FullMatrix<TYPEMAT>& FullMatrix<TYPEMAT>::operator =(
+ const FullMatrix<double>&);
#define TYPEMAT2 float
//template FullMatrix<TYPEMAT>& FullMatrix<TYPEMAT>::operator =<>(const FullMatrix<TYPEMAT2>&);
template void FullMatrix<TYPEMAT>::fill<TYPEMAT2> (
- const FullMatrix<TYPEMAT2>&, const unsigned, const unsigned, const unsigned, const unsigned);
+ const FullMatrix<TYPEMAT2>&, unsigned, unsigned, unsigned, unsigned);
template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
+template void FullMatrix<TYPEMAT>::add<TYPEMAT2> (
+ const FullMatrix<TYPEMAT2>&, double, unsigned, unsigned, unsigned, unsigned);
template void FullMatrix<TYPEMAT>::Tadd<TYPEMAT2> (const TYPEMAT, const FullMatrix<TYPEMAT2>&);
template void FullMatrix<TYPEMAT>::mmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
template void FullMatrix<TYPEMAT>::Tmmult<TYPEMAT2> (FullMatrix<TYPEMAT2>&, const FullMatrix<TYPEMAT2>&, const bool) const;
#define TYPEVEC double
#define TYPERES double
-template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (const FullMatrix<TYPEVEC>&,
- const std::vector<unsigned int>&,
- const std::vector<unsigned int>&);
-template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const bool) const;
-template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const bool) const;
-template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (const Vector<TYPEVEC> &) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+ const FullMatrix<TYPEVEC>&,
+ const std::vector<unsigned int>&,
+ const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+ const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+ const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
template void FullMatrix<TYPEMAT>::householder<TYPEVEC>(Vector<TYPEVEC>&);
-template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(Vector<TYPEVEC>&, Vector<TYPEVEC>&);
+template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(
+ Vector<TYPEVEC>&, Vector<TYPEVEC>&);
+
template
-void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (Vector<TYPEVEC> &,
- const Vector<TYPEVEC> &,
- const TYPEMAT) const;
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+ Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
#undef TYPEVEC
#define TYPEVEC float
-template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (const FullMatrix<TYPEVEC>&,
- const std::vector<unsigned int>&,
- const std::vector<unsigned int>&);
-template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const bool) const;
-template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const bool) const;
-template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (const Vector<TYPEVEC> &) const;
-template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
-template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::fill_permutation<TYPEVEC> (
+ const FullMatrix<TYPEVEC>&,
+ const std::vector<unsigned int>&,
+ const std::vector<unsigned int>&);
+template void FullMatrix<TYPEMAT>::vmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template void FullMatrix<TYPEMAT>::Tvmult<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, bool) const;
+template double FullMatrix<TYPEMAT>::residual<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_norm_square<TYPEVEC> (
+ const Vector<TYPEVEC> &) const;
+template TYPEVEC FullMatrix<TYPEMAT>::matrix_scalar_product<TYPEVEC>(
+ const Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::forward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
+template void FullMatrix<TYPEMAT>::backward<TYPEVEC>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&) const;
template void FullMatrix<TYPEMAT>::householder<TYPEVEC>(Vector<TYPEVEC>&);
-template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(Vector<TYPEVEC>&, Vector<TYPEVEC>&);
+template double FullMatrix<TYPEMAT>::least_squares<TYPEVEC>(
+ Vector<TYPEVEC>&, Vector<TYPEVEC>&);
template
-void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (Vector<TYPEVEC> &,
- const Vector<TYPEVEC> &,
- const TYPEMAT) const;
+void FullMatrix<TYPEMAT>::precondition_Jacobi<TYPEVEC> (
+ Vector<TYPEVEC> &, const Vector<TYPEVEC> &, const TYPEMAT) const;
#undef TYPERES
#define TYPERES float
-template double FullMatrix<TYPEMAT>::residual<TYPEVEC,TYPERES>(Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;
+template double FullMatrix<TYPEMAT>::residual<TYPEVEC,TYPERES>(
+ Vector<TYPEVEC>&, const Vector<TYPEVEC>&, const Vector<TYPERES>&) const;