From aa2ba0647954c362c469e6a2340934dc05dccc26 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 20 Aug 2004 20:49:00 +0000 Subject: [PATCH] Add umfpack testers. git-svn-id: https://svn.dealii.org/trunk@9568 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/bits/Makefile | 3 +- tests/bits/umfpack_01.cc | 114 ++++++++++++++++++++++++++++++++++ tests/bits/umfpack_02.cc | 120 +++++++++++++++++++++++++++++++++++ tests/bits/umfpack_03.cc | 131 +++++++++++++++++++++++++++++++++++++++ 4 files changed, 367 insertions(+), 1 deletion(-) create mode 100644 tests/bits/umfpack_01.cc create mode 100644 tests/bits/umfpack_02.cc create mode 100644 tests/bits/umfpack_03.cc diff --git a/tests/bits/Makefile b/tests/bits/Makefile index 37f3067f49..d57444bf95 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -62,7 +62,8 @@ tests_x = anna_? \ cylinder_shell_* \ grid_in_* \ compressed_sparsity_pattern_* \ - point_difference_* + point_difference_* \ + umfpack_* # tests for the hp branch: # fe_collection_* diff --git a/tests/bits/umfpack_01.cc b/tests/bits/umfpack_01.cc new file mode 100644 index 0000000000..5c5aba1f94 --- /dev/null +++ b/tests/bits/umfpack_01.cc @@ -0,0 +1,114 @@ +//---------------------------- umfpack_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes +// +// 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. +// +//---------------------------- umfpack_01.cc --------------------------- + +// test the umfpack sparse direct solver on a mass matrix + +#include "../tests.h" +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include +#include +#include + +#include +#include + +#include +#include + + + +template +void test () +{ + deallog << dim << 'd' << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube (tria,0,1); + tria.refine_global (1); + + // destroy the uniformity of the matrix by + // refining one cell + tria.begin_active()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + tria.refine_global(8-2*dim); + + FE_Q fe (1); + DoFHandler dof_handler (tria); + dof_handler.distribute_dofs (fe); + + deallog << "Number of dofs = " << dof_handler.n_dofs() << std::endl; + + SparsityPattern sparsity_pattern; + sparsity_pattern.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern); + sparsity_pattern.compress (); + + SparseMatrix B; + B.reinit (sparsity_pattern); + + QGauss qr (2); + MatrixTools::create_mass_matrix (dof_handler, qr, B); + + // for a number of different solution + // vectors, make up a matching rhs vector + // and check what the UMFPACK solver finds + for (unsigned int i=0; i<3; ++i) + { + Vector solution (dof_handler.n_dofs()); + Vector x (dof_handler.n_dofs()); + Vector b (dof_handler.n_dofs()); + + for (unsigned int j=0; j (); + test<2> (); + test<3> (); +} diff --git a/tests/bits/umfpack_02.cc b/tests/bits/umfpack_02.cc new file mode 100644 index 0000000000..f591be902c --- /dev/null +++ b/tests/bits/umfpack_02.cc @@ -0,0 +1,120 @@ +//---------------------------- umfpack_02.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes +// +// 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. +// +//---------------------------- umfpack_02.cc --------------------------- + +// test the umfpack sparse direct solver on a mass matrix; test that the +// functions that allow for repeated solves yield the correct result even +// though the matrix is decomposed only once + +#include "../tests.h" +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include +#include +#include + +#include +#include + +#include +#include + + + +template +void test () +{ + deallog << dim << 'd' << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube (tria,0,1); + tria.refine_global (1); + + // destroy the uniformity of the matrix by + // refining one cell + tria.begin_active()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + tria.refine_global(8-2*dim); + + FE_Q fe (1); + DoFHandler dof_handler (tria); + dof_handler.distribute_dofs (fe); + + deallog << "Number of dofs = " << dof_handler.n_dofs() << std::endl; + + SparsityPattern sparsity_pattern; + sparsity_pattern.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern); + sparsity_pattern.compress (); + + SparseMatrix B; + B.reinit (sparsity_pattern); + + QGauss qr (2); + MatrixTools::create_mass_matrix (dof_handler, qr, B); + + // compute a decomposition of the matrix + SparseDirectUMFPACK Binv; + Binv.factorize (B); + + // for a number of different solution + // vectors, make up a matching rhs vector + // and check what the UMFPACK solver finds + for (unsigned int i=0; i<3; ++i) + { + Vector solution (dof_handler.n_dofs()); + Vector x (dof_handler.n_dofs()); + Vector b (dof_handler.n_dofs()); + + for (unsigned int j=0; j (); + test<2> (); + test<3> (); +} diff --git a/tests/bits/umfpack_03.cc b/tests/bits/umfpack_03.cc new file mode 100644 index 0000000000..5e8b229262 --- /dev/null +++ b/tests/bits/umfpack_03.cc @@ -0,0 +1,131 @@ +//---------------------------- umfpack_03.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2002, 2003, 2004 by the deal.II authors and Brian Carnes +// +// 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. +// +//---------------------------- umfpack_03.cc --------------------------- + +// test the umfpack sparse direct solver on a mass matrix that is slightly modified to make it nonsymmetric + +#include "../tests.h" +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include +#include +#include + +#include +#include + +#include +#include + + + +template +void test () +{ + deallog << dim << 'd' << std::endl; + + Triangulation tria; + GridGenerator::hyper_cube (tria,0,1); + tria.refine_global (1); + + // destroy the uniformity of the matrix by + // refining one cell + tria.begin_active()->set_refine_flag (); + tria.execute_coarsening_and_refinement (); + tria.refine_global(8-2*dim); + + FE_Q fe (1); + DoFHandler dof_handler (tria); + dof_handler.distribute_dofs (fe); + + deallog << "Number of dofs = " << dof_handler.n_dofs() << std::endl; + + SparsityPattern sparsity_pattern; + sparsity_pattern.reinit (dof_handler.n_dofs(), + dof_handler.n_dofs(), + dof_handler.max_couplings_between_dofs()); + DoFTools::make_sparsity_pattern (dof_handler, sparsity_pattern); + sparsity_pattern.compress (); + + SparseMatrix B; + B.reinit (sparsity_pattern); + + QGauss qr (2); + MatrixTools::create_mass_matrix (dof_handler, qr, B); + + // scale lower left part of the matrix by + // 1/2 and upper right part by 2 to make + // matrix nonsymmetric + for (SparseMatrix::iterator p=B.begin(); + p!=B.end(); ++p) + if (p->column() < p->row()) + p->value() = p->value()/2; + else if (p->column() > p->row()) + p->value() = p->value() * 2; + + // check that we've done it right + for (SparseMatrix::iterator p=B.begin(); + p!=B.end(); ++p) + if (p->column() != p->row()) + Assert (B(p->row(),p->column()) != B(p->column(),p->row()), + ExcInternalError()); + + // for a number of different solution + // vectors, make up a matching rhs vector + // and check what the UMFPACK solver finds + for (unsigned int i=0; i<3; ++i) + { + Vector solution (dof_handler.n_dofs()); + Vector x (dof_handler.n_dofs()); + Vector b (dof_handler.n_dofs()); + + for (unsigned int j=0; j (); + test<2> (); + test<3> (); +} -- 2.39.5