From c4a1cd03284e87a906a09736598d4b679f909501 Mon Sep 17 00:00:00 2001
From: Ce Qin <qince168@gmail.com>
Date: Sat, 1 Aug 2020 13:44:00 +0800
Subject: [PATCH] Fix step-18.

---
 source/grid/grid_refinement.cc                | 10 ++-
 tests/sharedtria/refine_and_coarsen_01.cc     | 89 +++++++++++++++++++
 .../refine_and_coarsen_01.mpirun=3.output     |  7 ++
 3 files changed, 104 insertions(+), 2 deletions(-)
 create mode 100644 tests/sharedtria/refine_and_coarsen_01.cc
 create mode 100644 tests/sharedtria/refine_and_coarsen_01.mpirun=3.output

diff --git a/source/grid/grid_refinement.cc b/source/grid/grid_refinement.cc
index f24f82709d..8c3afb2217 100644
--- a/source/grid/grid_refinement.cc
+++ b/source/grid/grid_refinement.cc
@@ -16,6 +16,8 @@
 
 #include <deal.II/base/template_constraints.h>
 
+#include <deal.II/distributed/tria_base.h>
+
 #include <deal.II/grid/grid_refinement.h>
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_accessor.h>
@@ -71,7 +73,9 @@ GridRefinement::refine(Triangulation<dim, spacedim> &tria,
 
   unsigned int marked = 0;
   for (const auto &cell : tria.active_cell_iterators())
-    if (cell->is_locally_owned() &&
+    if ((dynamic_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
+           &tria) == nullptr ||
+         cell->is_locally_owned()) &&
         std::fabs(criteria(cell->active_cell_index())) >= new_threshold)
       {
         if (max_to_mark != numbers::invalid_unsigned_int &&
@@ -95,7 +99,9 @@ GridRefinement::coarsen(Triangulation<dim, spacedim> &tria,
   Assert(criteria.is_non_negative(), ExcNegativeCriteria());
 
   for (const auto &cell : tria.active_cell_iterators())
-    if (cell->is_locally_owned() &&
+    if ((dynamic_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
+           &tria) == nullptr ||
+         cell->is_locally_owned()) &&
         std::fabs(criteria(cell->active_cell_index())) <= threshold)
       if (!cell->refine_flag_set())
         cell->set_coarsen_flag();
diff --git a/tests/sharedtria/refine_and_coarsen_01.cc b/tests/sharedtria/refine_and_coarsen_01.cc
new file mode 100644
index 0000000000..770e9116c6
--- /dev/null
+++ b/tests/sharedtria/refine_and_coarsen_01.cc
@@ -0,0 +1,89 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2010 - 2018 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// The purpose of this test is to ensure that the p::s::T is
+// refined/coarsened in the same way on all participating processors.
+
+#include <deal.II/distributed/shared_tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+
+#include "../tests.h"
+
+
+
+template <int dim>
+void
+check()
+{
+  parallel::shared::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(3);
+
+  Vector<float> estimated_error_per_cell;
+  estimated_error_per_cell.reinit(tria.n_active_cells());
+  for (unsigned int i = 0; i < estimated_error_per_cell.size(); ++i)
+    estimated_error_per_cell(i) = i + 1;
+
+  GridRefinement::refine_and_coarsen_fixed_number(tria,
+                                                  estimated_error_per_cell,
+                                                  0.2,
+                                                  0.1);
+
+  std::vector<bool> refined_cells(tria.n_active_cells() * dim);
+  tria.save_refine_flags(refined_cells);
+  int n_refined_cells =
+    std::count(refined_cells.begin(), refined_cells.end(), true);
+
+  std::vector<bool> coarsened_cells(tria.n_active_cells() * dim);
+  tria.save_coarsen_flags(coarsened_cells);
+  int n_coarsened_cells =
+    std::count(coarsened_cells.begin(), coarsened_cells.end(), true);
+
+  tria.execute_coarsening_and_refinement();
+  int n_cells = tria.n_active_cells();
+
+  if (Utilities::MPI::max(n_refined_cells, MPI_COMM_WORLD) ==
+      Utilities::MPI::min(n_refined_cells, MPI_COMM_WORLD))
+    deallog << "OK" << std::endl;
+
+  if (Utilities::MPI::max(n_coarsened_cells, MPI_COMM_WORLD) ==
+      Utilities::MPI::min(n_coarsened_cells, MPI_COMM_WORLD))
+    deallog << "OK" << std::endl;
+
+  if (Utilities::MPI::max(n_cells, MPI_COMM_WORLD) ==
+      Utilities::MPI::min(n_cells, MPI_COMM_WORLD))
+    deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  mpi_initlog();
+
+  deallog.push("2d");
+  check<2>();
+  deallog.pop();
+
+  deallog.push("3d");
+  check<3>();
+  deallog.pop();
+
+  return 0;
+}
diff --git a/tests/sharedtria/refine_and_coarsen_01.mpirun=3.output b/tests/sharedtria/refine_and_coarsen_01.mpirun=3.output
new file mode 100644
index 0000000000..43fb2379a9
--- /dev/null
+++ b/tests/sharedtria/refine_and_coarsen_01.mpirun=3.output
@@ -0,0 +1,7 @@
+
+DEAL:2d::OK
+DEAL:2d::OK
+DEAL:2d::OK
+DEAL:3d::OK
+DEAL:3d::OK
+DEAL:3d::OK
-- 
2.39.5