From 98619f091ec3164d8ddca383f75d4cf46fadf395 Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Mon, 4 Sep 2023 22:23:42 +0200 Subject: [PATCH] MatrixFree: do not update ghost values of a ghosted src vector --- include/deal.II/matrix_free/matrix_free.h | 31 ++++++++++++++++--- tests/matrix_free/step-48.cc | 1 + tests/matrix_free/step-48c.cc | 1 + tests/matrix_free/stokes_computation.cc | 3 -- .../stokes_computation.cc | 3 -- 5 files changed, 29 insertions(+), 10 deletions(-) diff --git a/include/deal.II/matrix_free/matrix_free.h b/include/deal.II/matrix_free/matrix_free.h index 78410ba48d..f109651ec6 100644 --- a/include/deal.II/matrix_free/matrix_free.h +++ b/include/deal.II/matrix_free/matrix_free.h @@ -3413,7 +3413,10 @@ namespace internal ExcNotImplemented()); if (ghosts_set) - ghosts_were_set = true; + { + ghosts_were_set = true; + return; + } vec.update_ghost_values(); } @@ -3441,7 +3444,10 @@ namespace internal ExcNotImplemented()); if (ghosts_set) - ghosts_were_set = true; + { + ghosts_were_set = true; + return; + } vec.update_ghost_values_start(component_in_block_vector + channel_shift); } @@ -3472,7 +3478,10 @@ namespace internal ExcNotImplemented()); if (ghosts_set) - ghosts_were_set = true; + { + ghosts_were_set = true; + return; + } if (vec.size() != 0) { @@ -3536,6 +3545,10 @@ namespace internal const VectorType &vec) { (void)component_in_block_vector; + + if (ghosts_were_set) + return; + vec.update_ghost_values_finish(); } @@ -3559,6 +3572,9 @@ namespace internal "Type mismatch between VectorType and VectorDataExchange"); (void)component_in_block_vector; + if (ghosts_were_set) + return; + if (vec.size() != 0) { # ifdef DEAL_II_WITH_MPI @@ -4092,7 +4108,10 @@ namespace internal ExcNotImplemented()); if (ghosts_set) - exchanger.ghosts_were_set = true; + { + exchanger.ghosts_were_set = true; + return; + } vec.update_ghost_values(); } @@ -4463,6 +4482,10 @@ namespace internal const VectorStruct &vec, VectorDataExchange &exchanger) { + // return immediately if there is nothing to do. + if (exchanger.ghosts_were_set == true) + return; + exchanger.reset_ghost_values(vec); } diff --git a/tests/matrix_free/step-48.cc b/tests/matrix_free/step-48.cc index 8c62ae8258..165e0811e9 100644 --- a/tests/matrix_free/step-48.cc +++ b/tests/matrix_free/step-48.cc @@ -329,6 +329,7 @@ namespace Step48 norm_per_cell, QGauss(fe_degree + 1), VectorTools::L2_norm); + solution.zero_out_ghost_values(); const double solution_norm = std::sqrt(Utilities::MPI::sum(norm_per_cell.norm_sqr(), MPI_COMM_WORLD)); diff --git a/tests/matrix_free/step-48c.cc b/tests/matrix_free/step-48c.cc index 330eb55f84..a8d7479946 100644 --- a/tests/matrix_free/step-48c.cc +++ b/tests/matrix_free/step-48c.cc @@ -299,6 +299,7 @@ namespace Step48 norm_per_cell, QGauss(fe_degree + 1), VectorTools::L2_norm); + solution.zero_out_ghost_values(); const double solution_norm = std::sqrt(Utilities::MPI::sum(norm_per_cell.norm_sqr(), MPI_COMM_WORLD)); diff --git a/tests/matrix_free/stokes_computation.cc b/tests/matrix_free/stokes_computation.cc index c685e92c4d..2661b52c6c 100644 --- a/tests/matrix_free/stokes_computation.cc +++ b/tests/matrix_free/stokes_computation.cc @@ -1148,10 +1148,7 @@ namespace StokesClass stokes_matrix.initialize_dof_vector(solution); stokes_matrix.initialize_dof_vector(system_rhs); - solution.update_ghost_values(); solution.collect_sizes(); - - system_rhs.update_ghost_values(); system_rhs.collect_sizes(); } diff --git a/tests/multigrid-global-coarsening/stokes_computation.cc b/tests/multigrid-global-coarsening/stokes_computation.cc index 9810667cbc..c0ab2be5bb 100644 --- a/tests/multigrid-global-coarsening/stokes_computation.cc +++ b/tests/multigrid-global-coarsening/stokes_computation.cc @@ -1148,10 +1148,7 @@ namespace StokesClass stokes_matrix.initialize_dof_vector(solution); stokes_matrix.initialize_dof_vector(system_rhs); - solution.update_ghost_values(); solution.collect_sizes(); - - system_rhs.update_ghost_values(); system_rhs.collect_sizes(); } -- 2.39.5