From 9c0f3957b0b60a23cdb2e69ee9057aca9552035f Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Mon, 19 Sep 2016 20:48:05 +0200 Subject: [PATCH] Fix bug in threaded matrixfree --- .../deal.II/matrix_free/dof_info.templates.h | 1 + .../parallel_multigrid_adaptive_mf.cc | 27 ++-- ...h_mpi=true.with_p4est=true.mpirun=7.output | 144 ++++++++++++------ 3 files changed, 115 insertions(+), 57 deletions(-) diff --git a/include/deal.II/matrix_free/dof_info.templates.h b/include/deal.II/matrix_free/dof_info.templates.h index 33789b3ccf..1ee634d7a4 100644 --- a/include/deal.II/matrix_free/dof_info.templates.h +++ b/include/deal.II/matrix_free/dof_info.templates.h @@ -1218,6 +1218,7 @@ no_constraint: // adjust end of boundary cells to the remainder size_info.boundary_cells_end += (remainder+vectorization_length-1)/vectorization_length; + start_nonboundary = 0; } else { diff --git a/tests/matrix_free/parallel_multigrid_adaptive_mf.cc b/tests/matrix_free/parallel_multigrid_adaptive_mf.cc index 746302e207..a7493f8165 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_mf.cc +++ b/tests/matrix_free/parallel_multigrid_adaptive_mf.cc @@ -57,11 +57,15 @@ public: const DoFHandler &dof_handler, const MGConstrainedDoFs &mg_constrained_dofs, const typename FunctionMap::type &dirichlet_boundary, - const unsigned int level = numbers::invalid_unsigned_int) + const unsigned int level, + const bool threaded) { const QGauss<1> quad (n_q_points_1d); typename MatrixFree::AdditionalData addit_data; - addit_data.tasks_parallel_scheme = MatrixFree::AdditionalData::none; + if (threaded) + addit_data.tasks_parallel_scheme = MatrixFree::AdditionalData::partition_partition; + else + addit_data.tasks_parallel_scheme = MatrixFree::AdditionalData::none; addit_data.tasks_block_size = 3; addit_data.level_mg_handler = level; addit_data.mpi_communicator = MPI_COMM_WORLD; @@ -431,7 +435,7 @@ public: template -void do_test (const DoFHandler &dof) +void do_test (const DoFHandler &dof, const bool threaded) { deallog << "Testing " << dof.get_fe().get_name(); deallog << std::endl; @@ -453,7 +457,7 @@ void do_test (const DoFHandler &dof) MappingQ mapping(fe_degree+1); LaplaceOperator fine_matrix; fine_matrix.initialize(mapping, dof, mg_constrained_dofs, dirichlet_boundary, - numbers::invalid_unsigned_int); + numbers::invalid_unsigned_int, threaded); LinearAlgebra::distributed::Vector in, sol; fine_matrix.initialize_dof_vector(in); @@ -472,7 +476,7 @@ void do_test (const DoFHandler &dof) for (unsigned int level = 0; level > mg_interface_matrices; mg_interface_matrices.resize(0, dof.get_triangulation().n_global_levels()-1); @@ -524,7 +528,7 @@ void do_test (const DoFHandler &dof) { // avoid output from inner (coarse-level) solver - deallog.depth_file(2); + deallog.depth_file(3); ReductionControl control(30, 1e-20, 1e-7); SolverCG > solver(control); solver.solve(fine_matrix, sol, in, preconditioner); @@ -547,8 +551,8 @@ void test () for (typename Triangulation::active_cell_iterator cell=tria.begin_active(); cell != tria.end(); ++cell) if (cell->is_locally_owned() && - (cell->center().norm() < 0.5 && (cell->level() < 5 || - cell->center().norm() > 0.45) + ((cell->center().norm() < 0.5 && (cell->level() < 5 || + cell->center().norm() > 0.45)) || (dim == 2 && cell->center().norm() > 1.2))) cell->set_refine_flag(); @@ -558,7 +562,12 @@ void test () dof.distribute_dofs(fe); dof.distribute_mg_dofs(fe); - do_test (dof); + deallog.push("nothread"); + do_test (dof, false); + deallog.pop(); + deallog.push("threaded"); + do_test (dof, true); + deallog.pop(); } } diff --git a/tests/matrix_free/parallel_multigrid_adaptive_mf.with_mpi=true.with_p4est=true.mpirun=7.output b/tests/matrix_free/parallel_multigrid_adaptive_mf.with_mpi=true.with_p4est=true.mpirun=7.output index b7f4d1b5bf..b984bc41d2 100644 --- a/tests/matrix_free/parallel_multigrid_adaptive_mf.with_mpi=true.with_p4est=true.mpirun=7.output +++ b/tests/matrix_free/parallel_multigrid_adaptive_mf.with_mpi=true.with_p4est=true.mpirun=7.output @@ -1,49 +1,97 @@ -DEAL::Testing FE_Q<2>(1) -DEAL::Number of degrees of freedom: 507 -DEAL:cg::Starting value 21.9317 -DEAL:cg::Convergence step 5 value 7.81205e-07 -DEAL::Testing FE_Q<2>(1) -DEAL::Number of degrees of freedom: 881 -DEAL:cg::Starting value 27.7669 -DEAL:cg::Convergence step 6 value 9.24482e-07 -DEAL::Testing FE_Q<2>(1) -DEAL::Number of degrees of freedom: 2147 -DEAL:cg::Starting value 43.2551 -DEAL:cg::Convergence step 6 value 1.14520e-06 -DEAL::Testing FE_Q<2>(1) -DEAL::Number of degrees of freedom: 6517 -DEAL:cg::Starting value 76.9220 -DEAL:cg::Convergence step 6 value 1.55925e-06 -DEAL::Testing FE_Q<2>(3) -DEAL::Number of degrees of freedom: 4259 -DEAL:cg::Starting value 64.2573 -DEAL:cg::Convergence step 5 value 1.63117e-06 -DEAL::Testing FE_Q<2>(3) -DEAL::Number of degrees of freedom: 7457 -DEAL:cg::Starting value 83.1084 -DEAL:cg::Convergence step 5 value 7.69301e-06 -DEAL::Testing FE_Q<2>(3) -DEAL::Number of degrees of freedom: 18535 -DEAL:cg::Starting value 130.977 -DEAL:cg::Convergence step 5 value 6.04333e-06 -DEAL::Testing FE_Q<3>(1) -DEAL::Number of degrees of freedom: 186 -DEAL:cg::Starting value 12.3693 -DEAL:cg::Convergence step 3 value 8.40297e-08 -DEAL::Testing FE_Q<3>(1) -DEAL::Number of degrees of freedom: 648 -DEAL:cg::Starting value 21.1424 -DEAL:cg::Convergence step 6 value 1.05183e-07 -DEAL::Testing FE_Q<3>(1) -DEAL::Number of degrees of freedom: 2930 -DEAL:cg::Starting value 47.1699 -DEAL:cg::Convergence step 7 value 1.64836e-06 -DEAL::Testing FE_Q<3>(2) -DEAL::Number of degrees of freedom: 1106 -DEAL:cg::Starting value 30.8707 -DEAL:cg::Convergence step 5 value 5.33421e-08 -DEAL::Testing FE_Q<3>(2) -DEAL::Number of degrees of freedom: 4268 -DEAL:cg::Starting value 57.4891 -DEAL:cg::Convergence step 7 value 1.03521e-06 +DEAL:nothread::Testing FE_Q<2>(1) +DEAL:nothread::Number of degrees of freedom: 507 +DEAL:nothread:cg::Starting value 21.9317 +DEAL:nothread:cg::Convergence step 5 value 7.81205e-07 +DEAL:threaded::Testing FE_Q<2>(1) +DEAL:threaded::Number of degrees of freedom: 507 +DEAL:threaded:cg::Starting value 21.9317 +DEAL:threaded:cg::Convergence step 5 value 7.81205e-07 +DEAL:nothread::Testing FE_Q<2>(1) +DEAL:nothread::Number of degrees of freedom: 881 +DEAL:nothread:cg::Starting value 27.7669 +DEAL:nothread:cg::Convergence step 6 value 9.24482e-07 +DEAL:threaded::Testing FE_Q<2>(1) +DEAL:threaded::Number of degrees of freedom: 881 +DEAL:threaded:cg::Starting value 27.7669 +DEAL:threaded:cg::Convergence step 6 value 9.24482e-07 +DEAL:nothread::Testing FE_Q<2>(1) +DEAL:nothread::Number of degrees of freedom: 2147 +DEAL:nothread:cg::Starting value 43.2551 +DEAL:nothread:cg::Convergence step 6 value 1.14520e-06 +DEAL:threaded::Testing FE_Q<2>(1) +DEAL:threaded::Number of degrees of freedom: 2147 +DEAL:threaded:cg::Starting value 43.2551 +DEAL:threaded:cg::Convergence step 6 value 1.14520e-06 +DEAL:nothread::Testing FE_Q<2>(1) +DEAL:nothread::Number of degrees of freedom: 6517 +DEAL:nothread:cg::Starting value 76.9220 +DEAL:nothread:cg::Convergence step 6 value 1.55925e-06 +DEAL:threaded::Testing FE_Q<2>(1) +DEAL:threaded::Number of degrees of freedom: 6517 +DEAL:threaded:cg::Starting value 76.9220 +DEAL:threaded:cg::Convergence step 6 value 1.55925e-06 +DEAL:nothread::Testing FE_Q<2>(3) +DEAL:nothread::Number of degrees of freedom: 4259 +DEAL:nothread:cg::Starting value 64.2573 +DEAL:nothread:cg::Convergence step 5 value 1.63125e-06 +DEAL:threaded::Testing FE_Q<2>(3) +DEAL:threaded::Number of degrees of freedom: 4259 +DEAL:threaded:cg::Starting value 64.2573 +DEAL:threaded:cg::Convergence step 5 value 1.63120e-06 +DEAL:nothread::Testing FE_Q<2>(3) +DEAL:nothread::Number of degrees of freedom: 7457 +DEAL:nothread:cg::Starting value 83.1084 +DEAL:nothread:cg::Convergence step 5 value 7.69258e-06 +DEAL:threaded::Testing FE_Q<2>(3) +DEAL:threaded::Number of degrees of freedom: 7457 +DEAL:threaded:cg::Starting value 83.1084 +DEAL:threaded:cg::Convergence step 5 value 7.69315e-06 +DEAL:nothread::Testing FE_Q<2>(3) +DEAL:nothread::Number of degrees of freedom: 18535 +DEAL:nothread:cg::Starting value 130.977 +DEAL:nothread:cg::Convergence step 5 value 6.04397e-06 +DEAL:threaded::Testing FE_Q<2>(3) +DEAL:threaded::Number of degrees of freedom: 18535 +DEAL:threaded:cg::Starting value 130.977 +DEAL:threaded:cg::Convergence step 5 value 6.04357e-06 +DEAL:nothread::Testing FE_Q<3>(1) +DEAL:nothread::Number of degrees of freedom: 186 +DEAL:nothread:cg::Starting value 12.3693 +DEAL:nothread:cg::Convergence step 3 value 8.40297e-08 +DEAL:threaded::Testing FE_Q<3>(1) +DEAL:threaded::Number of degrees of freedom: 186 +DEAL:threaded:cg::Starting value 12.3693 +DEAL:threaded:cg::Convergence step 3 value 8.40297e-08 +DEAL:nothread::Testing FE_Q<3>(1) +DEAL:nothread::Number of degrees of freedom: 648 +DEAL:nothread:cg::Starting value 21.1424 +DEAL:nothread:cg::Convergence step 6 value 1.05183e-07 +DEAL:threaded::Testing FE_Q<3>(1) +DEAL:threaded::Number of degrees of freedom: 648 +DEAL:threaded:cg::Starting value 21.1424 +DEAL:threaded:cg::Convergence step 6 value 1.05183e-07 +DEAL:nothread::Testing FE_Q<3>(1) +DEAL:nothread::Number of degrees of freedom: 2930 +DEAL:nothread:cg::Starting value 47.1699 +DEAL:nothread:cg::Convergence step 7 value 1.64836e-06 +DEAL:threaded::Testing FE_Q<3>(1) +DEAL:threaded::Number of degrees of freedom: 2930 +DEAL:threaded:cg::Starting value 47.1699 +DEAL:threaded:cg::Convergence step 7 value 1.64836e-06 +DEAL:nothread::Testing FE_Q<3>(2) +DEAL:nothread::Number of degrees of freedom: 1106 +DEAL:nothread:cg::Starting value 30.8707 +DEAL:nothread:cg::Convergence step 5 value 5.33407e-08 +DEAL:threaded::Testing FE_Q<3>(2) +DEAL:threaded::Number of degrees of freedom: 1106 +DEAL:threaded:cg::Starting value 30.8707 +DEAL:threaded:cg::Convergence step 5 value 5.33448e-08 +DEAL:nothread::Testing FE_Q<3>(2) +DEAL:nothread::Number of degrees of freedom: 4268 +DEAL:nothread:cg::Starting value 57.4891 +DEAL:nothread:cg::Convergence step 7 value 1.03521e-06 +DEAL:threaded::Testing FE_Q<3>(2) +DEAL:threaded::Number of degrees of freedom: 4268 +DEAL:threaded:cg::Starting value 57.4891 +DEAL:threaded:cg::Convergence step 7 value 1.03521e-06 -- 2.39.5