]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
SparseMatrixEZ added to binary libs
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Feb 2003 10:03:36 +0000 (10:03 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Feb 2003 10:03:36 +0000 (10:03 +0000)
git-svn-id: https://svn.dealii.org/trunk@7120 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/lac/include/lac/sparse_matrix_ez.2.templates [new file with mode: 0644]
deal.II/lac/include/lac/sparse_matrix_ez.h
deal.II/lac/include/lac/sparse_matrix_ez.templates.h
deal.II/lac/source/sparse_matrix_ez.double.cc [new file with mode: 0644]
deal.II/lac/source/sparse_matrix_ez.float.cc [new file with mode: 0644]
tests/lac/sparse_matrices.cc
tests/lac/testmatrix.cc

diff --git a/deal.II/lac/include/lac/sparse_matrix_ez.2.templates b/deal.II/lac/include/lac/sparse_matrix_ez.2.templates
new file mode 100644 (file)
index 0000000..9423287
--- /dev/null
@@ -0,0 +1,65 @@
+//---------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//---------------------------------------------------------------------------
+
+
+// File copied from sparse_matrix.2.templates and modified
+
+// Driver file for SparseMatrixEZ functions with two types.
+
+// TYPEMAT and TYPE2 are defined in sparsematrix?.cc
+
+//template SparseMatrixEZ<TYPEMAT> &
+//SparseMatrixEZ<TYPEMAT>::copy_from<TYPE2> (const SparseMatrixEZ<TYPE2> &);
+
+//template 
+//void SparseMatrixEZ<TYPEMAT>::copy_from<TYPE2> (const FullMatrix<TYPE2> &);
+
+//template void SparseMatrixEZ<TYPEMAT>::add_scaled<TYPE2> (const TYPEMAT,
+//                                                       const SparseMatrixEZ<TYPE2> &);
+
+template void SparseMatrixEZ<TYPEMAT>::vmult<TYPE2> (Vector<TYPE2> &,
+                                                    const Vector<TYPE2> &) const;
+template void SparseMatrixEZ<TYPEMAT>::Tvmult<TYPE2> (Vector<TYPE2> &,
+                                                     const Vector<TYPE2> &) const;
+template void SparseMatrixEZ<TYPEMAT>::vmult_add<TYPE2> (Vector<TYPE2> &,
+                                                        const Vector<TYPE2> &) const;
+template void SparseMatrixEZ<TYPEMAT>::Tvmult_add<TYPE2> (Vector<TYPE2> &,
+                                                         const Vector<TYPE2> &) const;
+
+//template TYPE2
+//SparseMatrixEZ<TYPEMAT>::matrix_norm_square<TYPE2> (const Vector<TYPE2> &) const;
+
+//template TYPE2
+//SparseMatrixEZ<TYPEMAT>::matrix_scalar_product<TYPE2> (const Vector<TYPE2> &,
+//                                                    const Vector<TYPE2> &) const;
+
+//template TYPE2 SparseMatrixEZ<TYPEMAT>::residual<TYPE2> (Vector<TYPE2> &,
+//                                                      const Vector<TYPE2> &,
+//                                                      const Vector<TYPE2> &) const;
+
+template void SparseMatrixEZ<TYPEMAT>::precondition_SSOR<TYPE2> (Vector<TYPE2> &,
+                                                                const Vector<TYPE2> &,
+                                                                const TYPEMAT) const;
+
+template void SparseMatrixEZ<TYPEMAT>::precondition_SOR<TYPE2> (Vector<TYPE2> &,
+                                                               const Vector<TYPE2> &,
+                                                               const TYPEMAT) const;
+
+template void SparseMatrixEZ<TYPEMAT>::precondition_TSOR<TYPE2> (Vector<TYPE2> &,
+                                                                const Vector<TYPE2> &,
+                                                                const TYPEMAT) const;
+
+template void SparseMatrixEZ<TYPEMAT>::precondition_Jacobi<TYPE2> (Vector<TYPE2> &,
+                                                                  const Vector<TYPE2> &,
+                                                                  const TYPEMAT) const;
+
index 72e163fc583c97c32ce8b7c7c7a7b9893c806892..517625d1ec932cf4983ff0eaa54482fc7cbe1b2f 100644 (file)
@@ -72,8 +72,12 @@ template<typename number> class FullMatrix;
  * be achieved. @p{default_reserve} should then be an estimate for the
  * number of hanging nodes times @p{default_increment}.
  *
+ * Letting @p{default_increment} zero causes an exception whenever a
+ * row overflows.
+ *
  * If the rows are expected to be filled more or less from first to
- * last, using a @p{default_row_length} may not be such a bad idea.
+ * last, using a @p{default_row_length} of zero may not be such a bad
+ * idea.
  *
  * The name of this matrix is in reverence to a publication of the
  * Internal Revenue Service of the United States of America. I hope
@@ -838,6 +842,18 @@ class SparseMatrixEZ : public Subscriptor
                                      */
     unsigned int memory_consumption () const;
 
+                                    /**
+                                     * Print statistics. If @p{full}
+                                     * is @p{true}, prints a
+                                     * histogram of all existing row
+                                     * lengths and allocated row
+                                     * lengths. Otherwise, just the
+                                     * relation of allocated and used
+                                     * entries is shown.
+                                     */
+    template <class STREAM>
+    void print_statistics (STREAM& s, bool full = false);
+    
                                     /**
                                      * Exception for applying
                                      * inverse-type operators to
@@ -1269,7 +1285,7 @@ inline
 number SparseMatrixEZ<number>::el (const unsigned int i,
                                   const unsigned int j) const
 {
-  Entry* entry = locate(i,j);
+  const Entry* entry = locate(i,j);
   if (entry)
     return entry->value;
   return 0.;
@@ -1400,5 +1416,49 @@ SparseMatrixEZ<number>::conjugate_add (const MATRIXA& A,
     }
 }
 
+
+template <typename number>
+template <class STREAM>
+void
+SparseMatrixEZ<number>::print_statistics(STREAM& out, bool full)
+{
+  typename std::vector<RowInfo>::const_iterator row = row_info.begin();
+  const typename std::vector<RowInfo>::const_iterator endrow = row_info.end();
+
+                                  // Add up entries actually used
+  unsigned int entries_used = 0;
+  unsigned int max_length = 0;
+  for (; row != endrow ; ++ row)
+    {
+      entries_used += row->length;
+      if (max_length < row->length)
+       max_length = row->length;
+    }
+  
+                                  // Number of entries allocated is
+                                  // position of last entry used
+  --row;
+  unsigned int entries_alloc = row->start + row->length;
+
+  out << "SparseMatrixEZ:used     entries:" << entries_used << std::endl
+      << "SparseMatrixEZ:alloc    entries:" << entries_alloc << std::endl
+      << "SparseMatrixEZ:reserved entries:" << data.capacity() << std::endl;
+  
+  if (full)
+    {
+      std::vector<unsigned int> length_used (max_length+1);
+      
+      for (row = row_info.begin() ; row != endrow; ++row)
+       {
+         ++length_used[row->length];
+       }
+      for (unsigned int i=0; i< length_used.size();++i)
+       if (length_used[i] != 0)
+         out << "SparseMatrixEZ:entries\t" << i
+             << "\trows\t" << length_used[i]
+             << std::endl;
+    }
+}
+
 #endif
 /*----------------------------   sparse_matrix.h     ---------------------------*/
index 87ab3794b947e04957fbd6878f1858534c77c3e2..945e18324aae51457ea8ee9d6205a589100acf0f 100644 (file)
@@ -1,4 +1,5 @@
 #include <lac/sparse_matrix_ez.h>
+#include <lac/vector.h>
 
 #include <algorithm>
 
@@ -82,7 +83,7 @@ template <typename number>
 bool
 SparseMatrixEZ<number>::empty() const
 {
-  return ((n_columns == 0) && (row_start.size()==0));
+  return ((n_columns == 0) && (row_info.size()==0));
 }
 
 
@@ -155,7 +156,7 @@ SparseMatrixEZ<number>::vmult_add (Vector<somenumber>& dst,
   const unsigned int end_row = row_info.size();
   for (unsigned int row = 0; row < end_row; ++row)
     {
-      const RowInfo& ri = row_info[i];
+      const RowInfo& ri = row_info[row];
       typename std::vector<Entry>::const_iterator
        entry = data.begin() + ri.start;
       double s = 0.;
@@ -333,9 +334,8 @@ template <typename number>
 unsigned int
 SparseMatrixEZ<number>::memory_consumption() const
 {
-  unsigned int result =
+  return
     sizeof (*this)
-    + sizeof(unsigned int) * row_start.capacity();
-    + sizeof(SparseMatrixEZ<number>::Entry) * data.capacity();
-  return result;
+    + sizeof(unsigned int) * row_info.capacity()
+    + sizeof(typename SparseMatrixEZ<number>::Entry) * data.capacity();
 }
diff --git a/deal.II/lac/source/sparse_matrix_ez.double.cc b/deal.II/lac/source/sparse_matrix_ez.double.cc
new file mode 100644 (file)
index 0000000..90ad8d5
--- /dev/null
@@ -0,0 +1,38 @@
+//----------------------------  sparse_matrix.double.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  sparse_matrix.double.cc  ---------------------------
+
+
+#include <lac/sparse_matrix_ez.templates.h>
+
+#define TYPEMAT double
+
+template class SparseMatrixEZ<TYPEMAT>;
+
+#define TYPE2 float
+
+#include <lac/sparse_matrix_ez.2.templates>
+
+#undef TYPE2
+#define TYPE2 double
+
+#include <lac/sparse_matrix_ez.2.templates>
+
+                                // a prerelease of gcc3.0 fails to
+                                // compile this due to long double
+//  #undef TYPE2
+//  #define TYPE2 long double
+
+//  #include <lac/sparse_matrix.2.templates>
+
+#undef TYPE2
+#undef TYPEMAT
diff --git a/deal.II/lac/source/sparse_matrix_ez.float.cc b/deal.II/lac/source/sparse_matrix_ez.float.cc
new file mode 100644 (file)
index 0000000..88f5d32
--- /dev/null
@@ -0,0 +1,38 @@
+//----------------------------  sparse_matrix.double.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  sparse_matrix.double.cc  ---------------------------
+
+
+#include <lac/sparse_matrix_ez.templates.h>
+
+#define TYPEMAT float
+
+template class SparseMatrixEZ<TYPEMAT>;
+
+#define TYPE2 float
+
+#include <lac/sparse_matrix_ez.2.templates>
+
+#undef TYPE2
+#define TYPE2 double
+
+#include <lac/sparse_matrix_ez.2.templates>
+
+                                // a prerelease of gcc3.0 fails to
+                                // compile this due to long double
+//  #undef TYPE2
+//  #define TYPE2 long double
+
+//  #include <lac/sparse_matrix.2.templates>
+
+#undef TYPE2
+#undef TYPEMAT
index 060c4006ad08773099b3b579e764b443321fd3d2..cb73649a57923ae31cd8a5129a1d0cc8a6587fa9 100644 (file)
@@ -159,6 +159,7 @@ int main()
 #ifdef DEBUG
   check_iterator(E);
 #endif
+  E.print_statistics(deallog, true);
   E.add_scaled(-1., A);
   if (E.l2_norm() < 1.e-14)
     deallog << "Matrices are equal" << std::endl;
@@ -180,6 +181,7 @@ int main()
   deallog << "Assemble" << std::endl;
   testproblem.nine_point(E);
   check_vmult_quadratic(E_res, E, "9-SparseMatrixEZ<double>");
+  E.print_statistics(deallog, true);
 
   for (unsigned int i=0;i<E_res.size();++i)
     if (std::fabs(A_res[i] - E_res[i]) > 1.e-14)
index 8bde56b5c70367dbf0537837894921fe0381c439..dce736b5385fa83c4b2c30edaeb979cc3327109c 100644 (file)
@@ -17,8 +17,6 @@
 #include <lac/sparse_matrix_ez.h>
 #include <lac/vector.h>
 
-#include <lac/sparse_matrix_ez.templates.h>
-
 FDMatrix::FDMatrix(unsigned int nx, unsigned int ny)
                :
                nx(nx), ny(ny)

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.