]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new add function
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Feb 2004 10:13:21 +0000 (10:13 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Feb 2004 10:13:21 +0000 (10:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@8425 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/c-4-0.html
deal.II/lac/include/lac/full_matrix.h
deal.II/lac/include/lac/full_matrix.templates.h
deal.II/lac/source/full_matrix.double.cc
deal.II/lac/source/full_matrix.float.cc

index 90d662532383fe6357ac3f6b5dbcb00caea49e08..bc6dc5cfdb630249b6a1596bbc87be16d41d60ea 100644 (file)
@@ -276,6 +276,12 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 
 <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.
index 71cb8d02958cae189c27fbfc9655220fd0be0bec..416efca4d91f7b983760b2feadd26247ea2dc546 100644 (file)
@@ -52,7 +52,7 @@ template<typename number> class Vector;
  *
  * @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>
@@ -396,6 +396,32 @@ 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>.
index 1df2423dc2ab43285ce9e20e12859fc244723f82..01d7ac2262d1e5b020421f7c54353cfed79cee3a 100644 (file)
@@ -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<number>::add (const number               s,
   }
 }
 
+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>
index 0923aa128ff1fe3d9bf358bb86b977ab6ed86004..b74cc61cb243bd4b6cfe9696f124ce60143767c4 100644 (file)
@@ -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
 
 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
index fddd4eb756faeb293b6c2ae1921191614e954e5c..c130f26d3b3de94be871d6789e1c7c7529b64984 100644 (file)
@@ -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
 
 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;
@@ -35,45 +38,63 @@ 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
 
-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;

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.