From: Martin Kronbichler Date: Wed, 29 Apr 2009 16:34:27 +0000 (+0000) Subject: Add a test that tests the feature of inhomogeneous constraints set only with vectors. X-Git-Tag: v8.0.0~7722 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b4a210d3c16bf9340ce08a424c3c1a168c979544;p=dealii.git Add a test that tests the feature of inhomogeneous constraints set only with vectors. git-svn-id: https://svn.dealii.org/trunk@18793 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/tests/deal.II/Makefile b/tests/deal.II/Makefile index b862879553..c5f681d341 100644 --- a/tests/deal.II/Makefile +++ b/tests/deal.II/Makefile @@ -28,6 +28,7 @@ tests_x = block_matrices \ inhomogeneous_constraints \ inhomogeneous_constraints_nonsymmetric \ inhomogeneous_constraints_block \ + inhomogeneous_constraints_vector \ data_out \ derivative_* \ derivatives \ diff --git a/tests/deal.II/inhomogeneous_constraints_vector.cc b/tests/deal.II/inhomogeneous_constraints_vector.cc new file mode 100644 index 0000000000..d66685a411 --- /dev/null +++ b/tests/deal.II/inhomogeneous_constraints_vector.cc @@ -0,0 +1,296 @@ +//------------------ inhomogeneous_constraints.cc ------------------------ +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009 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. +// +//------------------ inhomogeneous_constraints.cc ------------------------ + + +// this function tests the correctness of the implementation of +// inhomogeneous constraints with +// ConstraintMatrix::distribute_local_to_global operating only on a vector, +// based on a modification of the step-5 tutorial program. It assumes +// correctness of the function ConstraintMatrix::distribute_local_to_global +// operating on both the matrix and vector simultaneously. + +#include "../tests.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +std::ofstream logfile("inhomogeneous_constraints_vector/output"); + +using namespace dealii; + +template +class LaplaceProblem +{ + public: + LaplaceProblem (); + void run (); + + private: + void setup_system (); + void assemble_system (); + void solve (); + + Triangulation triangulation; + FE_Q fe; + DoFHandler dof_handler; + + SparsityPattern sparsity_pattern; + SparseMatrix system_matrix; + + Vector solution; + Vector system_rhs; + ConstraintMatrix constraints; +}; + +template +class Coefficient : public Function +{ + public: + Coefficient () : Function() {} + + virtual double value (const Point &p, + const unsigned int component = 0) const; + + virtual void value_list (const std::vector > &points, + std::vector &values, + const unsigned int component = 0) const; +}; + + +template +double Coefficient::value (const Point &p, + const unsigned int /*component*/) const +{ + if (p.square() < 0.5*0.5) + return 20; + else + return 1; +} + + +template +void Coefficient::value_list (const std::vector > &points, + std::vector &values, + const unsigned int component) const +{ + Assert (values.size() == points.size(), + ExcDimensionMismatch (values.size(), points.size())); + Assert (component == 0, + ExcIndexRange (component, 0, 1)); + + const unsigned int n_points = points.size(); + + for (unsigned int i=0; i +LaplaceProblem::LaplaceProblem () : + fe (1), + dof_handler (triangulation) +{} + + + +template +void LaplaceProblem::setup_system () +{ + dof_handler.distribute_dofs (fe); + + deallog << " Number of degrees of freedom: " + << dof_handler.n_dofs() + << std::endl; + + constraints.clear(); + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ConstantFunction(1), + constraints); + constraints.close(); + 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, constraints, false); + sparsity_pattern.compress(); + + system_matrix.reinit (sparsity_pattern); + + solution.reinit (dof_handler.n_dofs()); + system_rhs.reinit (dof_handler.n_dofs()); +} + + + +template +void LaplaceProblem::assemble_system () +{ + QGauss quadrature_formula(2); + Vector test(dof_handler.n_dofs()); + + FEValues fe_values (fe, quadrature_formula, + update_values | update_gradients | + update_quadrature_points | update_JxW_values); + + const unsigned int dofs_per_cell = fe.dofs_per_cell; + const unsigned int n_q_points = quadrature_formula.size(); + + FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); + Vector cell_rhs (dofs_per_cell); + + std::vector local_dof_indices (dofs_per_cell); + + const Coefficient coefficient; + std::vector coefficient_values (n_q_points); + + typename DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + { + cell_matrix = 0; + cell_rhs = 0; + + fe_values.reinit (cell); + + coefficient.value_list (fe_values.get_quadrature_points(), + coefficient_values); + + for (unsigned int q_point=0; q_pointget_dof_indices (local_dof_indices); + + constraints.distribute_local_to_global(cell_matrix, cell_rhs, local_dof_indices, system_matrix, system_rhs); + constraints.distribute_local_to_global(cell_rhs, local_dof_indices, test, cell_matrix); + + } + test -= system_rhs; + Assert (test.l2_norm() <= 1e-12, ExcInternalError()); +} + +template +void LaplaceProblem::solve () +{ + SolverControl solver_control (1000, 1e-12); + SolverBicgstab<> bicgstab (solver_control); + + PreconditionSSOR<> preconditioner; + preconditioner.initialize(system_matrix, 1.2); + + bicgstab.solve (system_matrix, solution, system_rhs, + preconditioner); + + constraints.distribute (solution); + + deallog << " " << solver_control.last_step() + << " Bicgstab iterations needed to obtain convergence." + << std::endl; +} + + + +template +void LaplaceProblem::run () +{ + for (unsigned int cycle=0; cycle<3; ++cycle) + { + deallog << "Cycle " << cycle << ':' << std::endl; + + if (cycle != 0) + triangulation.refine_global (1); + else + { + GridGenerator::hyper_cube (triangulation, -1, 1); + triangulation.refine_global (2); + { + typename DoFHandler::active_cell_iterator cell = dof_handler.begin_active(); + cell->set_refine_flag(); + cell++; + cell->set_refine_flag(); + } + triangulation.execute_coarsening_and_refinement(); + } + + deallog << " Number of active cells: " + << triangulation.n_active_cells() + << std::endl + << " Total number of cells: " + << triangulation.n_cells() + << std::endl; + + setup_system (); + assemble_system (); + solve (); + } +} + + +int main () +{ + deallog << std::setprecision (2); + logfile << std::setprecision (2); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + LaplaceProblem<2> laplace_problem_2d; + laplace_problem_2d.run (); + + LaplaceProblem<3> laplace_problem_3d; + laplace_problem_3d.run (); + + return 0; +} diff --git a/tests/deal.II/inhomogeneous_constraints_vector/cmp/generic b/tests/deal.II/inhomogeneous_constraints_vector/cmp/generic new file mode 100644 index 0000000000..ac83885d8c --- /dev/null +++ b/tests/deal.II/inhomogeneous_constraints_vector/cmp/generic @@ -0,0 +1,43 @@ + +DEAL::Cycle 0: +DEAL:: Number of active cells: 22 +DEAL:: Total number of cells: 29 +DEAL:: Number of degrees of freedom: 34 +DEAL:Bicgstab::Starting value 4.8 +DEAL:Bicgstab::Convergence step 7 value 0 +DEAL:: 7 Bicgstab iterations needed to obtain convergence. +DEAL::Cycle 1: +DEAL:: Number of active cells: 88 +DEAL:: Total number of cells: 117 +DEAL:: Number of degrees of freedom: 111 +DEAL:Bicgstab::Starting value 6.3 +DEAL:Bicgstab::Convergence step 10 value 0 +DEAL:: 10 Bicgstab iterations needed to obtain convergence. +DEAL::Cycle 2: +DEAL:: Number of active cells: 352 +DEAL:: Total number of cells: 469 +DEAL:: Number of degrees of freedom: 397 +DEAL:Bicgstab::Starting value 8.7 +DEAL:Bicgstab::Convergence step 18 value 0 +DEAL:: 18 Bicgstab iterations needed to obtain convergence. +DEAL::Cycle 0: +DEAL:: Number of active cells: 78 +DEAL:: Total number of cells: 89 +DEAL:: Number of degrees of freedom: 158 +DEAL:Bicgstab::Starting value 4.9 +DEAL:Bicgstab::Convergence step 7 value 0 +DEAL:: 7 Bicgstab iterations needed to obtain convergence. +DEAL::Cycle 1: +DEAL:: Number of active cells: 624 +DEAL:: Total number of cells: 713 +DEAL:: Number of degrees of freedom: 909 +DEAL:Bicgstab::Starting value 4.9 +DEAL:Bicgstab::Convergence step 10 value 0 +DEAL:: 10 Bicgstab iterations needed to obtain convergence. +DEAL::Cycle 2: +DEAL:: Number of active cells: 4992 +DEAL:: Total number of cells: 5705 +DEAL:: Number of degrees of freedom: 6065 +DEAL:Bicgstab::Starting value 4.9 +DEAL:Bicgstab::Convergence step 16 value 0 +DEAL:: 16 Bicgstab iterations needed to obtain convergence.