From 02906ea23d02a418c8e3d22ea7483f488f929e5b Mon Sep 17 00:00:00 2001 From: guido Date: Mon, 17 Feb 2003 10:03:36 +0000 Subject: [PATCH] SparseMatrixEZ added to binary libs git-svn-id: https://svn.dealii.org/trunk@7120 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/lac/sparse_matrix_ez.2.templates | 65 +++++++++++++++++++ deal.II/lac/include/lac/sparse_matrix_ez.h | 64 +++++++++++++++++- .../include/lac/sparse_matrix_ez.templates.h | 12 ++-- deal.II/lac/source/sparse_matrix_ez.double.cc | 38 +++++++++++ deal.II/lac/source/sparse_matrix_ez.float.cc | 38 +++++++++++ tests/lac/sparse_matrices.cc | 2 + tests/lac/testmatrix.cc | 2 - 7 files changed, 211 insertions(+), 10 deletions(-) create mode 100644 deal.II/lac/include/lac/sparse_matrix_ez.2.templates create mode 100644 deal.II/lac/source/sparse_matrix_ez.double.cc create mode 100644 deal.II/lac/source/sparse_matrix_ez.float.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 index 0000000000..9423287001 --- /dev/null +++ b/deal.II/lac/include/lac/sparse_matrix_ez.2.templates @@ -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 & +//SparseMatrixEZ::copy_from (const SparseMatrixEZ &); + +//template +//void SparseMatrixEZ::copy_from (const FullMatrix &); + +//template void SparseMatrixEZ::add_scaled (const TYPEMAT, +// const SparseMatrixEZ &); + +template void SparseMatrixEZ::vmult (Vector &, + const Vector &) const; +template void SparseMatrixEZ::Tvmult (Vector &, + const Vector &) const; +template void SparseMatrixEZ::vmult_add (Vector &, + const Vector &) const; +template void SparseMatrixEZ::Tvmult_add (Vector &, + const Vector &) const; + +//template TYPE2 +//SparseMatrixEZ::matrix_norm_square (const Vector &) const; + +//template TYPE2 +//SparseMatrixEZ::matrix_scalar_product (const Vector &, +// const Vector &) const; + +//template TYPE2 SparseMatrixEZ::residual (Vector &, +// const Vector &, +// const Vector &) const; + +template void SparseMatrixEZ::precondition_SSOR (Vector &, + const Vector &, + const TYPEMAT) const; + +template void SparseMatrixEZ::precondition_SOR (Vector &, + const Vector &, + const TYPEMAT) const; + +template void SparseMatrixEZ::precondition_TSOR (Vector &, + const Vector &, + const TYPEMAT) const; + +template void SparseMatrixEZ::precondition_Jacobi (Vector &, + const Vector &, + const TYPEMAT) const; + diff --git a/deal.II/lac/include/lac/sparse_matrix_ez.h b/deal.II/lac/include/lac/sparse_matrix_ez.h index 72e163fc58..517625d1ec 100644 --- a/deal.II/lac/include/lac/sparse_matrix_ez.h +++ b/deal.II/lac/include/lac/sparse_matrix_ez.h @@ -72,8 +72,12 @@ template 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 + void print_statistics (STREAM& s, bool full = false); + /** * Exception for applying * inverse-type operators to @@ -1269,7 +1285,7 @@ inline number SparseMatrixEZ::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::conjugate_add (const MATRIXA& A, } } + +template +template +void +SparseMatrixEZ::print_statistics(STREAM& out, bool full) +{ + typename std::vector::const_iterator row = row_info.begin(); + const typename std::vector::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 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 ---------------------------*/ diff --git a/deal.II/lac/include/lac/sparse_matrix_ez.templates.h b/deal.II/lac/include/lac/sparse_matrix_ez.templates.h index 87ab3794b9..945e18324a 100644 --- a/deal.II/lac/include/lac/sparse_matrix_ez.templates.h +++ b/deal.II/lac/include/lac/sparse_matrix_ez.templates.h @@ -1,4 +1,5 @@ #include +#include #include @@ -82,7 +83,7 @@ template bool SparseMatrixEZ::empty() const { - return ((n_columns == 0) && (row_start.size()==0)); + return ((n_columns == 0) && (row_info.size()==0)); } @@ -155,7 +156,7 @@ SparseMatrixEZ::vmult_add (Vector& 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::const_iterator entry = data.begin() + ri.start; double s = 0.; @@ -333,9 +334,8 @@ template unsigned int SparseMatrixEZ::memory_consumption() const { - unsigned int result = + return sizeof (*this) - + sizeof(unsigned int) * row_start.capacity(); - + sizeof(SparseMatrixEZ::Entry) * data.capacity(); - return result; + + sizeof(unsigned int) * row_info.capacity() + + sizeof(typename SparseMatrixEZ::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 index 0000000000..90ad8d5076 --- /dev/null +++ b/deal.II/lac/source/sparse_matrix_ez.double.cc @@ -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 + +#define TYPEMAT double + +template class SparseMatrixEZ; + +#define TYPE2 float + +#include + +#undef TYPE2 +#define TYPE2 double + +#include + + // a prerelease of gcc3.0 fails to + // compile this due to long double +// #undef TYPE2 +// #define TYPE2 long double + +// #include + +#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 index 0000000000..88f5d323cd --- /dev/null +++ b/deal.II/lac/source/sparse_matrix_ez.float.cc @@ -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 + +#define TYPEMAT float + +template class SparseMatrixEZ; + +#define TYPE2 float + +#include + +#undef TYPE2 +#define TYPE2 double + +#include + + // a prerelease of gcc3.0 fails to + // compile this due to long double +// #undef TYPE2 +// #define TYPE2 long double + +// #include + +#undef TYPE2 +#undef TYPEMAT diff --git a/tests/lac/sparse_matrices.cc b/tests/lac/sparse_matrices.cc index 060c4006ad..cb73649a57 100644 --- a/tests/lac/sparse_matrices.cc +++ b/tests/lac/sparse_matrices.cc @@ -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"); + E.print_statistics(deallog, true); for (unsigned int i=0;i 1.e-14) diff --git a/tests/lac/testmatrix.cc b/tests/lac/testmatrix.cc index 8bde56b5c7..dce736b538 100644 --- a/tests/lac/testmatrix.cc +++ b/tests/lac/testmatrix.cc @@ -17,8 +17,6 @@ #include #include -#include - FDMatrix::FDMatrix(unsigned int nx, unsigned int ny) : nx(nx), ny(ny) -- 2.39.5