From: Wolfgang Bangerth Date: Sat, 10 Mar 2018 05:38:06 +0000 (-0700) Subject: Break the main part of a loop out of the loop and put it into a lambda. X-Git-Tag: v9.0.0-rc1~324^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0bac8300bd32f77318b0c80386c98d870df274ea;p=dealii.git Break the main part of a loop out of the loop and put it into a lambda. --- diff --git a/include/deal.II/fe/fe_tools.templates.h b/include/deal.II/fe/fe_tools.templates.h index 0a5a3c8787..a931c07b46 100644 --- a/include/deal.II/fe/fe_tools.templates.h +++ b/include/deal.II/fe/fe_tools.templates.h @@ -2157,106 +2157,109 @@ namespace FETools mass.gauss_jordan(); } - // loop over all possible - // refinement cases - unsigned int ref_case = (isotropic_only) - ? RefinementCase::isotropic_refinement - : RefinementCase::cut_x; - for (; ref_case <= RefinementCase::isotropic_refinement; ++ref_case) - { - const unsigned int - nc = GeometryInfo::n_children(RefinementCase(ref_case)); - for (unsigned int i=0; i &inverse_mass_matrix, + std::vector> &matrices) + { + const unsigned int + nc = GeometryInfo::n_children(RefinementCase(ref_case)); - // create a respective refinement on the - // triangulation - Triangulation tr; - GridGenerator::hyper_cube (tr, 0, 1); - tr.begin_active()->set_refine_flag(RefinementCase(ref_case)); - tr.execute_coarsening_and_refinement(); + for (unsigned int i=0; i fine (StaticMappingQ1::mapping, fe, q_fine, - update_quadrature_points | update_JxW_values | - update_values); + // create a respective refinement on the triangulation + Triangulation tr; + GridGenerator::hyper_cube (tr, 0, 1); + tr.begin_active()->set_refine_flag(RefinementCase(ref_case)); + tr.execute_coarsening_and_refinement(); - typename Triangulation::cell_iterator coarse_cell - = tr.begin(0); + FEValues fine (StaticMappingQ1::mapping, fe, q_fine, + update_quadrature_points | update_JxW_values | + update_values); - Vector v_coarse(n); - Vector v_fine(n); + typename Triangulation::cell_iterator coarse_cell + = tr.begin(0); - for (unsigned int cell_number=0; cell_number &this_matrix = matrices[ref_case-1][cell_number]; - - // Compute right hand side, - // which is a fine level basis - // function tested with the - // coarse level functions. - fine.reinit(coarse_cell->child(cell_number)); - const std::vector > &q_points_fine = fine.get_quadrature_points(); - std::vector > q_points_coarse(q_points_fine.size()); - for (unsigned int q=0; q q_coarse (q_points_coarse, - fine.get_JxW_values()); - FEValues coarse (StaticMappingQ1::mapping, fe, q_coarse, update_values); - coarse.reinit(coarse_cell); - - // Build RHS - - const std::vector &JxW = fine.get_JxW_values(); - - // Outer loop over all fine - // grid shape functions phi_j - for (unsigned int j=0; j v_coarse(n); + Vector v_fine(n); + + for (unsigned int cell_number=0; cell_number &this_matrix = matrices[cell_number]; + + // Compute right hand side, which is a fine level basis + // function tested with the coarse level functions. + fine.reinit(coarse_cell->child(cell_number)); + const std::vector > &q_points_fine = fine.get_quadrature_points(); + std::vector > q_points_coarse(q_points_fine.size()); + for (unsigned int q=0; q q_coarse (q_points_coarse, + fine.get_JxW_values()); + FEValues coarse (StaticMappingQ1::mapping, fe, q_coarse, + update_values); + coarse.reinit(coarse_cell); - double update = 0; + // Build RHS + + const std::vector &JxW = fine.get_JxW_values(); + + // Outer loop over all fine grid shape functions phi_j + for (unsigned int j=0; j::isotropic_refinement + : RefinementCase::cut_x; + for (; ref_case <= RefinementCase::isotropic_refinement; ++ref_case) + { + compute_one_case (ref_case, mass, matrices[ref_case-1]); } }