From: Wolfgang Bangerth Date: Tue, 9 Sep 2003 16:09:15 +0000 (+0000) Subject: Allow some block matrix-usual vector operations. X-Git-Tag: v8.0.0~16242 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b29b622a78c9b1307f0ab5250b7bc7768da86c3a;p=dealii.git Allow some block matrix-usual vector operations. git-svn-id: https://svn.dealii.org/trunk@7968 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 8e3c3cba99..6686c6773e 100644 --- a/deal.II/doc/news/c-4-0.html +++ b/deal.II/doc/news/c-4-0.html @@ -142,6 +142,16 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

lac

    +
  1. New: Some of the member matrix-vector functions of the + BlockSparseMatrix class that + previously could only be used with arguments of type BlockVector can now also be used with the + usual Vector class provided the + block matrix has only one block row or column. +
    + (WB 2002/09/09) +

    +
  2. Fixed: FullMatrix::copy_from didn't compile when copying from a sparse matrix. This is now fixed. diff --git a/tests/lac/sparse_matrices.cc b/tests/lac/sparse_matrices.cc index 8f8ffc379e..1dafad0e11 100644 --- a/tests/lac/sparse_matrices.cc +++ b/tests/lac/sparse_matrices.cc @@ -17,6 +17,7 @@ #include "testmatrix.h" #include #include +#include #include #include #include @@ -93,6 +94,49 @@ check_vmult_quadratic(std::vector& residuals, } + +void +check_vmult_quadratic(std::vector& residuals, + const BlockSparseMatrix & A, + const char* prefix) +{ + deallog.push(prefix); + + Vector u(A.n()); + Vector f(A.m()); + GrowingVectorMemory<> mem; + + SolverControl control(10, 1.e-13, false); + SolverRichardson<> rich(control, mem, .01); + SolverRichardson<> prich(control, mem, 1.); + PreconditionIdentity identity; + PreconditionJacobi > jacobi; + jacobi.initialize(A, .5); + + PreconditionBlock, float>::AdditionalData + data((unsigned int) std::sqrt(A.n()+.3), 1.2); + + PreconditionBlockJacobi, float> block_jacobi; + block_jacobi.initialize(A, data); + + u = 0.; + f = 1.; + + PREC_CHECK(rich, solve, identity); + PREC_CHECK(prich, solve, jacobi); + u = 0.; + PREC_CHECK(prich, solve, block_jacobi); + + u = 0.; + deallog << "Transpose" << std::endl; + PREC_CHECK(rich, Tsolve, identity); + PREC_CHECK(prich, Tsolve, jacobi); + u = 0.; + PREC_CHECK(prich, Tsolve, block_jacobi); + deallog.pop(); +} + + template void check_iterator (const MATRIX& A) @@ -238,43 +282,71 @@ int main() std::vector A_res; std::vector E_res; - - deallog << "Structure" << std::endl; + + // usual sparse matrix SparsityPattern structure(dim, dim, 5); - testproblem.five_point_structure(structure); - structure.compress(); - SparseMatrix A(structure); - deallog << "Assemble" << std::endl; - testproblem.five_point(A, true); - check_vmult_quadratic(A_res, A, "5-SparseMatrix"); + SparseMatrix A; + { + deallog << "Structure" << std::endl; + testproblem.five_point_structure(structure); + structure.compress(); + A.reinit (structure); + deallog << "Assemble" << std::endl; + testproblem.five_point(A, true); + check_vmult_quadratic(A_res, A, "5-SparseMatrix"); +#ifdef DEBUG + check_iterator(A); +#endif + } + + // block sparse matrix with only + // one block + { + deallog << "Structure" << std::endl; + BlockSparsityPattern block_structure(1,1); + block_structure.block(0,0).reinit (dim, dim, 5); + block_structure.collect_sizes (); + testproblem.five_point_structure(block_structure.block(0,0)); + block_structure.compress(); + BlockSparseMatrix block_A(block_structure); + deallog << "Assemble" << std::endl; + testproblem.five_point(block_A, true); + std::vector block_A_res; + check_vmult_quadratic(block_A_res, block_A, "5-BlockSparseMatrix"); #ifdef DEBUG - check_iterator(A); + check_iterator(block_A); #endif + } + // ez sparse matrix SparseMatrixEZ E(dim,dim,row_length,2); - deallog << "Assemble" << std::endl; - testproblem.five_point(E, true); - check_vmult_quadratic(E_res, E, "5-SparseMatrixEZ"); + { + deallog << "Assemble" << std::endl; + testproblem.five_point(E, true); + check_vmult_quadratic(E_res, E, "5-SparseMatrixEZ"); #ifdef DEBUG - check_iterator(E); + 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; - else - deallog << "Matrices differ!!" << std::endl; + E.print_statistics(deallog, true); + E.add_scaled(-1., A); + if (E.l2_norm() < 1.e-14) + deallog << "Matrices are equal" << std::endl; + else + deallog << "Matrices differ!!" << std::endl; + } + + { + structure.reinit(dim, dim, 5); + deallog << "Structure" << std::endl; + structure.reinit(dim, dim, 9); + testproblem.nine_point_structure(structure); + structure.compress(); + A.reinit (structure); + deallog << "Assemble" << std::endl; + testproblem.nine_point(A); + check_vmult_quadratic(A_res, A, "9-SparseMatrix"); + } - A.clear(); - deallog << "Structure" << std::endl; - structure.reinit(dim, dim, 9); - testproblem.nine_point_structure(structure); - structure.compress(); - A.reinit(structure); - deallog << "Assemble" << std::endl; - testproblem.nine_point(A); - check_vmult_quadratic(A_res, A, "9-SparseMatrix"); - E.clear(); E.reinit(dim,dim,row_length,2); deallog << "Assemble" << std::endl; diff --git a/tests/lac/testmatrix.cc b/tests/lac/testmatrix.cc index ad5ff5f51c..d22d4f25f6 100644 --- a/tests/lac/testmatrix.cc +++ b/tests/lac/testmatrix.cc @@ -14,6 +14,7 @@ #include "testmatrix.h" #include +#include #include #include @@ -249,6 +250,8 @@ FDMatrix::gnuplot_print(std::ostream& s, const Vector& V) const template void FDMatrix::five_point(SparseMatrix&, bool) const; template void FDMatrix::five_point(SparseMatrix&, bool) const; +template void FDMatrix::five_point(BlockSparseMatrix&, bool) const; +template void FDMatrix::five_point(BlockSparseMatrix&, bool) const; template void FDMatrix::five_point(SparseMatrixEZ&, bool) const; template void FDMatrix::five_point(SparseMatrix&, bool) const; template void FDMatrix::nine_point(SparseMatrix&, bool) const;