From 2210755b36f12c96c55a99f6214ef2308c38ad14 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 14 Aug 2014 07:52:59 -0500 Subject: [PATCH] Add a test that uses range-based fors. --- tests/deal.II/range_based_for_step-6.cc | 358 ++++++++++++++++++ ...6.with_cxx11=true.with_umfpack=true.output | 16 + 2 files changed, 374 insertions(+) create mode 100644 tests/deal.II/range_based_for_step-6.cc create mode 100644 tests/deal.II/range_based_for_step-6.with_cxx11=true.with_umfpack=true.output diff --git a/tests/deal.II/range_based_for_step-6.cc b/tests/deal.II/range_based_for_step-6.cc new file mode 100644 index 0000000000..cec1eee680 --- /dev/null +++ b/tests/deal.II/range_based_for_step-6.cc @@ -0,0 +1,358 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2014 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + +// A variant of step-6 that replaces loops over cells with range-based for +// loops a la C++11 + +#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 + +#include +#include + + +#include + +#include + +#include + +using namespace dealii; + + + +template +class Step6 +{ +public: + Step6 (); + ~Step6 (); + + void run (); + +private: + void setup_system (); + void assemble_system (); + void solve (); + void refine_grid (); + void output_results (const unsigned int cycle) const; + + Triangulation triangulation; + + DoFHandler dof_handler; + FE_Q fe; + + ConstraintMatrix constraints; + + SparsityPattern sparsity_pattern; + SparseMatrix system_matrix; + + Vector solution; + Vector system_rhs; +}; + + + + +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) 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 +{ + const unsigned int n_points = points.size(); + + Assert (values.size() == n_points, + ExcDimensionMismatch (values.size(), n_points)); + + Assert (component == 0, + ExcIndexRange (component, 0, 1)); + + for (unsigned int i=0; i +Step6::Step6 () + : + dof_handler (triangulation), + fe (2) +{} + + + +template +Step6::~Step6 () +{ + dof_handler.clear (); +} + + + +template +void Step6::setup_system () +{ + dof_handler.distribute_dofs (fe); + + solution.reinit (dof_handler.n_dofs()); + system_rhs.reinit (dof_handler.n_dofs()); + + + constraints.clear (); + DoFTools::make_hanging_node_constraints (dof_handler, + constraints); + + + VectorTools::interpolate_boundary_values (dof_handler, + 0, + ZeroFunction(), + constraints); + + + constraints.close (); + + CompressedSparsityPattern c_sparsity(dof_handler.n_dofs()); + DoFTools::make_sparsity_pattern(dof_handler, + c_sparsity, + constraints, + /*keep_constrained_dofs = */ false); + + sparsity_pattern.copy_from(c_sparsity); + + system_matrix.reinit (sparsity_pattern); +} + + + +template +void Step6::assemble_system () +{ + const QGauss quadrature_formula(3); + + 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); + + for (auto cell : dof_handler.active_cell_iterators()) + { + cell_matrix = 0; + cell_rhs = 0; + + fe_values.reinit (cell); + + coefficient.value_list (fe_values.get_quadrature_points(), + coefficient_values); + + for (unsigned int q_index=0; q_indexget_dof_indices (local_dof_indices); + constraints.distribute_local_to_global (cell_matrix, + cell_rhs, + local_dof_indices, + system_matrix, + system_rhs); + } +} + + + + +template +void Step6::solve () +{ + SolverControl solver_control (1000, 1e-12); + SolverCG<> solver (solver_control); + + PreconditionSSOR<> preconditioner; + preconditioner.initialize(system_matrix, 1.2); + + solver.solve (system_matrix, solution, system_rhs, + preconditioner); + + constraints.distribute (solution); +} + + + +template +void Step6::refine_grid () +{ + Vector estimated_error_per_cell (triangulation.n_active_cells()); + + KellyErrorEstimator::estimate (dof_handler, + QGauss(3), + typename FunctionMap::type(), + solution, + estimated_error_per_cell); + + GridRefinement::refine_and_coarsen_fixed_number (triangulation, + estimated_error_per_cell, + 0.3, 0.03); + + triangulation.execute_coarsening_and_refinement (); +} + + + +template +void Step6::output_results (const unsigned int cycle) const +{ + Assert (cycle < 10, ExcNotImplemented()); + + std::string filename = "grid-"; + filename += ('0' + cycle); + filename += ".eps"; + + std::ofstream output (filename.c_str()); + + GridOut grid_out; + grid_out.write_eps (triangulation, output); +} + + + +template +void Step6::run () +{ + for (unsigned int cycle=0; cycle<3; ++cycle) + { + deallog << "Cycle " << cycle << ':' << std::endl; + + if (cycle == 0) + { + GridGenerator::hyper_ball (triangulation); + + static const HyperBallBoundary boundary; + triangulation.set_boundary (0, boundary); + + triangulation.refine_global (1); + } + else + refine_grid (); + + + deallog << " Number of active cells: " + << triangulation.n_active_cells() + << std::endl; + + setup_system (); + + deallog << " Number of degrees of freedom: " + << dof_handler.n_dofs() + << std::endl; + + assemble_system (); + solve (); + } +} + + + +int main () +{ + std::ofstream logfile ("output"); + deallog << std::setprecision (3); + + deallog.attach(logfile); + deallog.depth_console (0); + deallog.threshold_double(1.e-7); + + Step6<2> laplace_problem_2d; + laplace_problem_2d.run (); +} diff --git a/tests/deal.II/range_based_for_step-6.with_cxx11=true.with_umfpack=true.output b/tests/deal.II/range_based_for_step-6.with_cxx11=true.with_umfpack=true.output new file mode 100644 index 0000000000..1add5861bd --- /dev/null +++ b/tests/deal.II/range_based_for_step-6.with_cxx11=true.with_umfpack=true.output @@ -0,0 +1,16 @@ + +DEAL::Cycle 0: +DEAL:: Number of active cells: 20 +DEAL:: Number of degrees of freedom: 89 +DEAL:cg::Starting value 0.352 +DEAL:cg::Convergence step 22 value 0 +DEAL::Cycle 1: +DEAL:: Number of active cells: 44 +DEAL:: Number of degrees of freedom: 209 +DEAL:cg::Starting value 0.271 +DEAL:cg::Convergence step 30 value 0 +DEAL::Cycle 2: +DEAL:: Number of active cells: 92 +DEAL:: Number of degrees of freedom: 449 +DEAL:cg::Starting value 0.204 +DEAL:cg::Convergence step 51 value 0 -- 2.39.5