From: Daniel Arndt Date: Fri, 29 Jun 2018 19:10:51 +0000 (+0200) Subject: examples/step-13: Update indenting and modernize X-Git-Tag: v9.1.0-rc1~965^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a98428a0f95037864dbfc0e3de9bb3bf24ffad08;p=dealii.git examples/step-13: Update indenting and modernize --- diff --git a/examples/step-13/doc/results.dox b/examples/step-13/doc/results.dox index 1d95429259..66f1dae6fe 100644 --- a/examples/step-13/doc/results.dox +++ b/examples/step-13/doc/results.dox @@ -164,22 +164,20 @@ of the adjacent cells have hanging nodes; this destroys some superapproximation effects of which the globally refined mesh can profit. Answer 2: this quick hack @code - bool refinement_indicated = false; - typename Triangulation::active_cell_iterator cell; - for (cell=triangulation->begin_active(); - cell!=triangulation->end(); ++cell) - for (unsigned int v=0; v::vertices_per_cell; ++v) - if (cell->vertex(v) == Point(.5,.5)) - { - cell->clear_coarsen_flag(); - refinement_indicated |= cell->refine_flag_set(); - } - if (refinement_indicated) - for (cell=triangulation->begin_active(); - cell!=triangulation->end(); ++cell) - for (unsigned int v=0; v::vertices_per_cell; ++v) + bool refinement_indicated = false; + typename Triangulation::active_cell_iterator cell; + for (const auto &cell : triangulation.active_cell_iterators()) + for (unsigned int v=0; v::vertices_per_cell; ++v) if (cell->vertex(v) == Point(.5,.5)) - cell->set_refine_flag (); + { + cell->clear_coarsen_flag(); + refinement_indicated |= cell->refine_flag_set(); + } + if (refinement_indicated) + for (const auto &cell : triangulation.active_cell_iterators()) + for (unsigned int v=0; v::vertices_per_cell; ++v) + if (cell->vertex(v) == Point(.5,.5)) + cell->set_refine_flag (); @endcode in the refinement function of the Kelly refinement class right before executing refinement would improve the results (exercise: what does diff --git a/examples/step-13/step-13.cc b/examples/step-13/step-13.cc index d65b252ece..3e4fa29104 100644 --- a/examples/step-13/step-13.cc +++ b/examples/step-13/step-13.cc @@ -119,7 +119,7 @@ namespace Step13 class EvaluationBase { public: - virtual ~EvaluationBase(); + virtual ~EvaluationBase() = default; void set_refinement_cycle(const unsigned int refinement_cycle); @@ -131,13 +131,6 @@ namespace Step13 }; - // After the declaration has been discussed above, the implementation is - // rather straightforward: - template - EvaluationBase::~EvaluationBase() - {} - - template void EvaluationBase::set_refinement_cycle(const unsigned int step) @@ -225,54 +218,52 @@ namespace Step13 // vertex matches the evaluation point. If this is the case, then // extract the point value, set a flag that we have found the point of // interest, and exit the loop. - typename DoFHandler::active_cell_iterator cell = dof_handler - .begin_active(), - endc = dof_handler.end(); - bool evaluation_point_found = false; - for (; (cell != endc) && !evaluation_point_found; ++cell) - for (unsigned int vertex = 0; - vertex < GeometryInfo::vertices_per_cell; - ++vertex) - if (cell->vertex(vertex) == evaluation_point) - { - // In order to extract the point value from the global solution - // vector, pick that component that belongs to the vertex of - // interest, and, in case the solution is vector-valued, take - // the first component of it: - point_value = solution(cell->vertex_dof_index(vertex, 0)); - // Note that by this we have made an assumption that is not - // valid always and should be documented in the class - // declaration if this were code for a real application rather - // than a tutorial program: we assume that the finite element - // used for the solution we try to evaluate actually has degrees - // of freedom associated with vertices. This, for example, does - // not hold for discontinuous elements, were the support points - // for the shape functions happen to be located at the vertices, - // but are not associated with the vertices but rather with the - // cell interior, since association with vertices would imply - // continuity there. It would also not hold for edge oriented - // elements, and the like. - // - // Ideally, we would check this at the beginning of the - // function, for example by a statement like Assert - // (dof_handler.get_fe().dofs_per_vertex @> 0, - // ExcNotImplemented()), which should make it quite clear - // what is going wrong when the exception is triggered. In this - // case, we omit it (which is indeed bad style), but knowing - // that that does not hurt here, since the statement - // cell-@>vertex_dof_index(vertex,0) would fail if - // we asked it to give us the DoF index of a vertex if there - // were none. - // - // We stress again that this restriction on the allowed finite - // elements should be stated in the class documentation. - - // Since we found the right point, we now set the respective - // flag and exit the innermost loop. The outer loop will the - // also be terminated due to the set flag. - evaluation_point_found = true; - break; - }; + bool evaluation_point_found = false; + for (const auto &cell : dof_handler.active_cell_iterators()) + if (!evaluation_point_found) + for (unsigned int vertex = 0; + vertex < GeometryInfo::vertices_per_cell; + ++vertex) + if (cell->vertex(vertex) == evaluation_point) + { + // In order to extract the point value from the global solution + // vector, pick that component that belongs to the vertex of + // interest, and, in case the solution is vector-valued, take + // the first component of it: + point_value = solution(cell->vertex_dof_index(vertex, 0)); + // Note that by this we have made an assumption that is not + // valid always and should be documented in the class + // declaration if this were code for a real application rather + // than a tutorial program: we assume that the finite element + // used for the solution we try to evaluate actually has degrees + // of freedom associated with vertices. This, for example, does + // not hold for discontinuous elements, were the support points + // for the shape functions happen to be located at the vertices, + // but are not associated with the vertices but rather with the + // cell interior, since association with vertices would imply + // continuity there. It would also not hold for edge oriented + // elements, and the like. + // + // Ideally, we would check this at the beginning of the + // function, for example by a statement like Assert + // (dof_handler.get_fe().dofs_per_vertex @> 0, + // ExcNotImplemented()), which should make it quite clear + // what is going wrong when the exception is triggered. In this + // case, we omit it (which is indeed bad style), but knowing + // that that does not hurt here, since the statement + // cell-@>vertex_dof_index(vertex,0) would fail if + // we asked it to give us the DoF index of a vertex if there + // were none. + // + // We stress again that this restriction on the allowed finite + // elements should be stated in the class documentation. + + // Since we found the right point, we now set the respective + // flag and exit the innermost loop. The outer loop will also be + // terminated due to the set flag. + evaluation_point_found = true; + break; + }; // Finally, we'd like to make sure that we have indeed found the // evaluation point, since if that were not so we could not give a @@ -507,7 +498,7 @@ namespace Step13 { public: Base(Triangulation &coarse_grid); - virtual ~Base(); + virtual ~Base() = default; virtual void solve_problem() = 0; virtual void postprocess( @@ -528,11 +519,6 @@ namespace Step13 {} - template - Base::~Base() - {} - - // @sect4{A general solver class} // Following now the main class that implements assembling the matrix of @@ -618,10 +604,10 @@ namespace Step13 void solve(Vector &solution) const; - ConstraintMatrix hanging_node_constraints; - SparsityPattern sparsity_pattern; - SparseMatrix matrix; - Vector rhs; + AffineConstraints hanging_node_constraints; + SparsityPattern sparsity_pattern; + SparseMatrix matrix; + Vector rhs; }; @@ -997,7 +983,7 @@ namespace Step13 { hanging_node_constraints.clear(); - void (*mhnc_p)(const DoFHandler &, ConstraintMatrix &) = + void (*mhnc_p)(const DoFHandler &, AffineConstraints &) = &DoFTools::make_hanging_node_constraints; // Start a side task then continue on the main thread @@ -1111,11 +1097,7 @@ namespace Step13 std::vector rhs_values(n_q_points); std::vector local_dof_indices(dofs_per_cell); - typename DoFHandler::active_cell_iterator cell = this->dof_handler - .begin_active(), - endc = - this->dof_handler.end(); - for (; cell != endc; ++cell) + for (const auto &cell : this->dof_handler.active_cell_iterators()) { cell_rhs = 0; fe_values.reinit(cell); @@ -1124,8 +1106,9 @@ namespace Step13 for (unsigned int q_point = 0; q_point < n_q_points; ++q_point) for (unsigned int i = 0; i < dofs_per_cell; ++i) - cell_rhs(i) += (fe_values.shape_value(i, q_point) * - rhs_values[q_point] * fe_values.JxW(q_point)); + cell_rhs(i) += fe_values.shape_value(i, q_point) * // + rhs_values[q_point] * // + fe_values.JxW(q_point); cell->get_dof_indices(local_dof_indices); for (unsigned int i = 0; i < dofs_per_cell; ++i) @@ -1290,9 +1273,11 @@ namespace Step13 template - double Solution::value(const Point &p, - const unsigned int /*component*/) const + double Solution::value(const Point & p, + const unsigned int component) const { + (void)component; + AssertIndexRange(component, 1); double q = p(0); for (unsigned int i = 1; i < dim; ++i) q += std::sin(10 * p(i) + 5 * p(0) * p(0)); @@ -1316,9 +1301,11 @@ namespace Step13 template - double RightHandSide::value(const Point &p, - const unsigned int /*component*/) const + double RightHandSide::value(const Point & p, + const unsigned int component) const { + (void)component; + AssertIndexRange(component, 1); double q = p(0); for (unsigned int i = 1; i < dim; ++i) q += std::sin(10 * p(i) + 5 * p(0) * p(0)); @@ -1376,14 +1363,10 @@ namespace Step13 // annoying, but could be shortened by an alias, if so desired. solver.solve_problem(); - for (typename std::list< - Evaluation::EvaluationBase *>::const_iterator i = - postprocessor_list.begin(); - i != postprocessor_list.end(); - ++i) + for (const auto &postprocessor : postprocessor_list) { - (*i)->set_refinement_cycle(step); - solver.postprocess(**i); + postprocessor->set_refinement_cycle(step); + solver.postprocess(*postprocessor); }; @@ -1433,12 +1416,18 @@ namespace Step13 // Create a solver object of the kind indicated by the argument to this // function. If the name is not recognized, throw an exception! - LaplaceSolver::Base *solver = nullptr; + // The respective solver object is stored in a std::unique_ptr to avoid + // having to delete the pointer after use. For initializing, we want to use + // the C++14 function std::make_unique. Since deal.II only requires C++11 up + // to now, we define this function in a separate namespace called + // `std_cxx14`. In case the compiler supports C++14, this just calls + // std::make_unique. + std::unique_ptr> solver; if (solver_name == "global") - solver = new LaplaceSolver::RefinementGlobal( + solver = std_cxx14::make_unique>( triangulation, fe, quadrature, rhs_function, boundary_values); else if (solver_name == "kelly") - solver = new LaplaceSolver::RefinementKelly( + solver = std_cxx14::make_unique>( triangulation, fe, quadrature, rhs_function, boundary_values); else AssertThrow(false, ExcNotImplemented()); @@ -1464,10 +1453,8 @@ namespace Step13 // simulation on successively refined grids: run_simulation(*solver, postprocessor_list); - // When this all is done, write out the results of the point evaluations, - // and finally delete the solver object: + // When this all is done, write out the results of the point evaluations: results_table.write_text(std::cout); - delete solver; // And one blank line after all results: std::cout << std::endl;