From 632a69b3f354b8f33178d46c675b9732e97b7f29 Mon Sep 17 00:00:00 2001 From: Marco Feder Date: Thu, 4 May 2023 15:37:40 +0200 Subject: [PATCH] first part of review --- .../major/20230426FederMunchKronbichler | 3 + .../matrix_free/tensor_product_kernels.h | 6 +- .../multigrid/mg_transfer_global_coarsening.h | 10 +++- .../mg_transfer_global_coarsening.templates.h | 56 +++++++++---------- .../mg_transfer_util.h | 23 -------- ...rid_01.mpirun=1.with_trilinos=true.output} | 0 ...rid_02.mpirun=1.with_trilinos=true.output} | 0 ...rid_03.mpirun=1.with_trilinos=true.output} | 0 8 files changed, 42 insertions(+), 56 deletions(-) create mode 100644 doc/news/changes/major/20230426FederMunchKronbichler rename tests/multigrid-global-coarsening/{non_nested_multigrid_01.output => non_nested_multigrid_01.mpirun=1.with_trilinos=true.output} (100%) rename tests/multigrid-global-coarsening/{non_nested_multigrid_02.output => non_nested_multigrid_02.mpirun=1.with_trilinos=true.output} (100%) rename tests/multigrid-global-coarsening/{non_nested_multigrid_03.output => non_nested_multigrid_03.mpirun=1.with_trilinos=true.output} (100%) diff --git a/doc/news/changes/major/20230426FederMunchKronbichler b/doc/news/changes/major/20230426FederMunchKronbichler new file mode 100644 index 0000000000..753fcf86fd --- /dev/null +++ b/doc/news/changes/major/20230426FederMunchKronbichler @@ -0,0 +1,3 @@ +New: Add a two level transfer operator between non-nested multigrid levels. +
+(Marco Feder, Peter Munch, Martin Kronbichler, 2023/04/26) diff --git a/include/deal.II/matrix_free/tensor_product_kernels.h b/include/deal.II/matrix_free/tensor_product_kernels.h index b2f32a40ba..e02030183a 100644 --- a/include/deal.II/matrix_free/tensor_product_kernels.h +++ b/include/deal.II/matrix_free/tensor_product_kernels.h @@ -3488,10 +3488,10 @@ namespace internal constexpr unsigned int array_size = length > 0 ? length : 1; std::array shape_values_x; std::array shape_derivs_x; - for (unsigned int i = 0; i < array_size; ++i) + for (unsigned int j = 0; j < array_size; ++j) { - shape_values_x[i] = shapes[i][0][0]; - shape_derivs_x[i] = shapes[i][1][0]; + shape_values_x[j] = shapes[j][0][0]; + shape_derivs_x[j] = shapes[j][1][0]; } for (unsigned int i1 = 0; i1 < (dim > 1 ? length : 1); ++i1) { diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.h index 55c64c1f22..c6b2980810 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.h @@ -31,6 +31,8 @@ #include #include +#include + DEAL_II_NAMESPACE_OPEN @@ -771,7 +773,7 @@ private: /** * MappingInfo object needed as Mapping argument by FEPointEvaluation. */ - std::shared_ptr> mapping_info; + std::shared_ptr> mapping_info; /** * Helper class for reading from and writing to global vectors and for @@ -784,6 +786,12 @@ private: * Finite element of the coarse DoFHandler passed to reinit(). */ std::unique_ptr> fe_coarse; + + /** + * DoF indices of the fine cells, expressed in indices local to the MPI + * rank. + */ + std::vector level_dof_indices_fine; }; diff --git a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h index d15aa7a0a1..2b6115c832 100644 --- a/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h +++ b/include/deal.II/multigrid/mg_transfer_global_coarsening.templates.h @@ -3474,8 +3474,8 @@ namespace internal { namespace { - template - std::shared_ptr> + template + std::shared_ptr> fill_mapping_info(const Utilities::MPI::RemotePointEvaluation &rpe) { const auto &cell_data = rpe.get_cell_data(); @@ -3486,10 +3486,10 @@ namespace internal for (unsigned int i = 0; i < cell_data.cells.size(); ++i) { - typename Triangulation::active_cell_iterator cell = { + typename Triangulation::active_cell_iterator cell( &rpe.get_triangulation(), cell_data.cells[i].first, - cell_data.cells[i].second}; + cell_data.cells[i].second); const ArrayView> unit_points( cell_data.reference_point_values.data() + @@ -3497,19 +3497,14 @@ namespace internal cell_data.reference_point_ptrs[i + 1] - cell_data.reference_point_ptrs[i]); - std::vector> unit_points_vector_cell; - - unit_points_vector_cell.insert(unit_points_vector_cell.begin(), - unit_points.begin(), - unit_points.end()); - cell_iterators.emplace_back(cell); - unit_points_vector.emplace_back(unit_points_vector_cell); + unit_points_vector.emplace_back(unit_points.begin(), + unit_points.end()); } auto mapping_info = - std::make_shared>(rpe.get_mapping(), - update_values); + std::make_shared>( + rpe.get_mapping(), update_values); mapping_info->reinit_cells(cell_iterators, unit_points_vector); return mapping_info; @@ -3541,12 +3536,11 @@ MGTwoLevelTransferNonNested>:: // create partitioners and internal vectors { - IndexSet locally_relevant_dofs; - DoFTools::extract_locally_relevant_dofs(dof_handler_coarse, - locally_relevant_dofs); + IndexSet locally_active_dofs = + DoFTools::extract_locally_active_dofs(dof_handler_coarse); this->partitioner_coarse.reset( new Utilities::MPI::Partitioner(dof_handler_coarse.locally_owned_dofs(), - locally_relevant_dofs, + locally_active_dofs, dof_handler_coarse.get_communicator())); this->vec_coarse.reinit(this->partitioner_coarse); @@ -3610,7 +3604,7 @@ MGTwoLevelTransferNonNested>:: rpe.reinit(points, dof_handler_coarse.get_triangulation(), mapping_coarse); // set up MappingInfo for easier data access - mapping_info = internal::fill_mapping_info(rpe); + mapping_info = internal::fill_mapping_info(rpe); // set up constraints const auto &cell_data = rpe.get_cell_data(); @@ -3621,11 +3615,11 @@ MGTwoLevelTransferNonNested>:: for (unsigned int i = 0; i < cell_data.cells.size(); ++i) { - typename DoFHandler::active_cell_iterator cell = { + typename DoFHandler::active_cell_iterator cell( &rpe.get_triangulation(), cell_data.cells[i].first, cell_data.cells[i].second, - &dof_handler_coarse}; + &dof_handler_coarse); constraint_info.read_dof_indices(i, numbers::invalid_unsigned_int, @@ -3697,23 +3691,24 @@ MGTwoLevelTransferNonNested>:: if (rpe.is_map_unique() == false) { const auto evaluation_point_results_temp = evaluation_point_results; - evaluation_point_results.assign(rpe.get_point_ptrs().size() - 1, 0); + evaluation_point_results.clear(); + evaluation_point_results.reserve(rpe.get_point_ptrs().size() - 1); const auto &ptr = rpe.get_point_ptrs(); for (unsigned int i = 0; i < ptr.size() - 1; ++i) { const auto n_entries = ptr[i + 1] - ptr[i]; - if (n_entries == 0) - continue; Number result = 0.0; - for (unsigned int j = 0; j < n_entries; ++j) - result += evaluation_point_results_temp[ptr[i] + j]; - result /= n_entries; - - evaluation_point_results[i] = result; + if (n_entries > 0) + { + for (unsigned int j = 0; j < n_entries; ++j) + result += evaluation_point_results_temp[ptr[i] + j]; + result /= n_entries; + } + evaluation_point_results.push_back(result); } } @@ -3817,7 +3812,10 @@ MGTwoLevelTransferNonNested>:: &external_partitioner_fine) { this->internal_enable_inplace_operations_if_possible( - external_partitioner_coarse, external_partitioner_fine, constraint_info); + external_partitioner_coarse, + external_partitioner_fine, + constraint_info, + this->level_dof_indices_fine); } diff --git a/tests/multigrid-global-coarsening/mg_transfer_util.h b/tests/multigrid-global-coarsening/mg_transfer_util.h index 931981f9ae..1031b7651f 100644 --- a/tests/multigrid-global-coarsening/mg_transfer_util.h +++ b/tests/multigrid-global-coarsening/mg_transfer_util.h @@ -208,29 +208,6 @@ test_non_nested_transfer( // via constraint matrix constraint_fine.distribute(dst); - // Visualize the result - { - DataOut<2> data_out; - data_out.attach_dof_handler(dof_handler_coarse); - data_out.add_data_vector( - src, - "coarse_solution", - DataOut_DoFData::DataVectorType::type_dof_data); - data_out.build_patches(); - std::ofstream output_coarse("coarse_sol.vtk"); - data_out.write_vtk(output_coarse); - data_out.clear(); - std::ofstream output_fine("prolonged_sol.vtk"); - data_out.attach_dof_handler(dof_handler_fine); - data_out.add_data_vector( - dst, - "prolonged_solution", - DataOut_DoFData::DataVectorType::type_dof_data); - - data_out.build_patches(); - data_out.write_vtk(output_fine); - } - // print norms if (true) { diff --git a/tests/multigrid-global-coarsening/non_nested_multigrid_01.output b/tests/multigrid-global-coarsening/non_nested_multigrid_01.mpirun=1.with_trilinos=true.output similarity index 100% rename from tests/multigrid-global-coarsening/non_nested_multigrid_01.output rename to tests/multigrid-global-coarsening/non_nested_multigrid_01.mpirun=1.with_trilinos=true.output diff --git a/tests/multigrid-global-coarsening/non_nested_multigrid_02.output b/tests/multigrid-global-coarsening/non_nested_multigrid_02.mpirun=1.with_trilinos=true.output similarity index 100% rename from tests/multigrid-global-coarsening/non_nested_multigrid_02.output rename to tests/multigrid-global-coarsening/non_nested_multigrid_02.mpirun=1.with_trilinos=true.output diff --git a/tests/multigrid-global-coarsening/non_nested_multigrid_03.output b/tests/multigrid-global-coarsening/non_nested_multigrid_03.mpirun=1.with_trilinos=true.output similarity index 100% rename from tests/multigrid-global-coarsening/non_nested_multigrid_03.output rename to tests/multigrid-global-coarsening/non_nested_multigrid_03.mpirun=1.with_trilinos=true.output -- 2.39.5