From 0780010602667335e0050fa24a87cc09add20774 Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Sun, 1 Mar 2020 14:58:34 -0500 Subject: [PATCH] use global_cell_index in many places --- .../deal.II/base/parsed_convergence_table.h | 2 +- include/deal.II/distributed/grid_refinement.h | 10 +++---- include/deal.II/distributed/tria_base.h | 4 +-- include/deal.II/grid/grid_refinement.h | 8 ++--- include/deal.II/grid/tria.h | 2 +- source/distributed/grid_refinement.cc | 30 +++++++------------ source/distributed/grid_refinement.inst.in | 6 ++-- source/distributed/tria_base.cc | 4 +-- source/grid/grid_refinement.cc | 10 +++---- source/grid/grid_refinement.inst.in | 4 +-- source/grid/tria.cc | 2 +- 11 files changed, 37 insertions(+), 45 deletions(-) diff --git a/include/deal.II/base/parsed_convergence_table.h b/include/deal.II/base/parsed_convergence_table.h index 76659f449a..4103301136 100644 --- a/include/deal.II/base/parsed_convergence_table.h +++ b/include/deal.II/base/parsed_convergence_table.h @@ -566,7 +566,7 @@ ParsedConvergenceTable::error_from_exact( AssertDimension(exact.n_components, n_components); AssertDimension(dh.get_fe().n_components(), n_components); - const unsigned int n_active_cells = + const types::global_cell_index n_active_cells = dh.get_triangulation().n_global_active_cells(); const unsigned int n_dofs = dh.n_dofs(); diff --git a/include/deal.II/distributed/grid_refinement.h b/include/deal.II/distributed/grid_refinement.h index 6b4358ec89..22fcd22865 100644 --- a/include/deal.II/distributed/grid_refinement.h +++ b/include/deal.II/distributed/grid_refinement.h @@ -57,7 +57,7 @@ namespace internal number compute_threshold(const dealii::Vector & criteria, const std::pair &global_min_and_max, - const types::global_dof_index n_target_cells, + const types::global_cell_index n_target_cells, MPI_Comm mpi_communicator); } // namespace RefineAndCoarsenFixedNumber @@ -131,10 +131,10 @@ namespace parallel refine_and_coarsen_fixed_number( parallel::distributed::Triangulation &tria, const dealii::Vector & criteria, - const double top_fraction_of_cells, - const double bottom_fraction_of_cells, - const types::global_dof_index max_n_cells = - std::numeric_limits::max()); + const double top_fraction_of_cells, + const double bottom_fraction_of_cells, + const types::global_cell_index max_n_cells = + std::numeric_limits::max()); /** * Like dealii::GridRefinement::refine_and_coarsen_fixed_fraction, but diff --git a/include/deal.II/distributed/tria_base.h b/include/deal.II/distributed/tria_base.h index 0e6507defc..458c72699c 100644 --- a/include/deal.II/distributed/tria_base.h +++ b/include/deal.II/distributed/tria_base.h @@ -145,7 +145,7 @@ namespace parallel * by each processor. This equals the overall number of active cells in * the triangulation. */ - virtual types::global_dof_index + virtual types::global_cell_index n_global_active_cells() const override; /** @@ -244,7 +244,7 @@ namespace parallel * The total number of active cells (sum of @p * n_locally_owned_active_cells). */ - types::global_dof_index n_global_active_cells; + types::global_cell_index n_global_active_cells; /** * The global number of levels computed as the maximum number of levels * taken over all MPI ranks, so n_levels()<=n_global_levels = diff --git a/include/deal.II/grid/grid_refinement.h b/include/deal.II/grid/grid_refinement.h index 3ec8e03862..af09fb1c50 100644 --- a/include/deal.II/grid/grid_refinement.h +++ b/include/deal.II/grid/grid_refinement.h @@ -87,10 +87,10 @@ namespace GridRefinement template std::pair adjust_refine_and_coarsen_number_fraction( - const types::global_dof_index current_n_cells, - const types::global_dof_index max_n_cells, - const double top_fraction_of_cells, - const double bottom_fraction_of_cells); + const types::global_cell_index current_n_cells, + const types::global_cell_index max_n_cells, + const double top_fraction_of_cells, + const double bottom_fraction_of_cells); /** * This function provides a strategy to mark cells for refinement and diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h index 621dab4caa..b8902f7225 100644 --- a/include/deal.II/grid/tria.h +++ b/include/deal.II/grid/tria.h @@ -3006,7 +3006,7 @@ public: * may return a value greater than the number of active cells reported by * the triangulation object on the current processor. */ - virtual types::global_dof_index + virtual types::global_cell_index n_global_active_cells() const; diff --git a/source/distributed/grid_refinement.cc b/source/distributed/grid_refinement.cc index fc45b9f954..4eb262fab7 100644 --- a/source/distributed/grid_refinement.cc +++ b/source/distributed/grid_refinement.cc @@ -228,7 +228,7 @@ namespace internal number compute_threshold(const dealii::Vector & criteria, const std::pair &global_min_and_max, - const types::global_dof_index n_target_cells, + const types::global_cell_index n_target_cells, MPI_Comm mpi_communicator) { double interesting_range[2] = {global_min_and_max.first, @@ -257,7 +257,7 @@ namespace internal // Count how many of our own elements would be above this // threshold: - types::global_dof_index my_count = + types::global_cell_index my_count = std::count_if(criteria.begin(), criteria.end(), [test_threshold](const double c) { @@ -265,16 +265,8 @@ namespace internal }); // Potentially accumulate in a 64bit int to avoid overflow: - types::global_dof_index total_count = 0; - - ierr = MPI_Reduce(&my_count, - &total_count, - 1, - DEAL_II_DOF_INDEX_MPI_TYPE, - MPI_SUM, - master_mpi_rank, - mpi_communicator); - AssertThrowMPI(ierr); + types::global_cell_index total_count = + Utilities::MPI::sum(my_count, mpi_communicator); // now adjust the range. if we have too many cells, we take the // upper half of the previous range, otherwise the lower half. @@ -434,9 +426,9 @@ namespace parallel refine_and_coarsen_fixed_number( parallel::distributed::Triangulation &tria, const dealii::Vector & criteria, - const double top_fraction_of_cells, - const double bottom_fraction_of_cells, - const types::global_dof_index max_n_cells) + const double top_fraction_of_cells, + const double bottom_fraction_of_cells, + const types::global_cell_index max_n_cells) { Assert(criteria.size() == tria.n_active_cells(), ExcDimensionMismatch(criteria.size(), tria.n_active_cells())); @@ -478,8 +470,8 @@ namespace parallel GridRefinement::RefineAndCoarsenFixedNumber::compute_threshold( locally_owned_indicators, global_min_and_max, - static_cast(adjusted_fractions.first * - tria.n_global_active_cells()), + static_cast(adjusted_fractions.first * + tria.n_global_active_cells()), mpi_communicator); // compute bottom threshold only if necessary. otherwise use a threshold @@ -489,8 +481,8 @@ namespace parallel GridRefinement::RefineAndCoarsenFixedNumber::compute_threshold( locally_owned_indicators, global_min_and_max, - static_cast( - std::ceil((1 - adjusted_fractions.second) * + static_cast( + std::ceil((1. - adjusted_fractions.second) * tria.n_global_active_cells())), mpi_communicator); else diff --git a/source/distributed/grid_refinement.inst.in b/source/distributed/grid_refinement.inst.in index 688294715d..c0fc46f489 100644 --- a/source/distributed/grid_refinement.inst.in +++ b/source/distributed/grid_refinement.inst.in @@ -34,7 +34,7 @@ for (S : REAL_SCALARS) template S compute_threshold(const dealii::Vector &, const std::pair &, - const types::global_dof_index, + const types::global_cell_index, MPI_Comm); \} namespace RefineAndCoarsenFixedFraction @@ -69,7 +69,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS) const dealii::Vector &, const double, const double, - const types::global_dof_index); + const types::global_cell_index); template void refine_and_coarsen_fixed_fraction &, const double, const double, - const types::global_dof_index); + const types::global_cell_index); template void refine_and_coarsen_fixed_fraction - types::global_dof_index + types::global_cell_index TriangulationBase::n_global_active_cells() const { return number_cache.n_global_active_cells; @@ -219,7 +219,7 @@ namespace parallel // Potentially cast to a 64 bit type before accumulating to avoid overflow: number_cache.n_global_active_cells = - Utilities::MPI::sum(static_cast( + Utilities::MPI::sum(static_cast( number_cache.n_locally_owned_active_cells), this->mpi_communicator); diff --git a/source/grid/grid_refinement.cc b/source/grid/grid_refinement.cc index 29dc3bdf70..93dc2350eb 100644 --- a/source/grid/grid_refinement.cc +++ b/source/grid/grid_refinement.cc @@ -110,10 +110,10 @@ GridRefinement::coarsen(Triangulation &tria, template std::pair GridRefinement::adjust_refine_and_coarsen_number_fraction( - const types::global_dof_index current_n_cells, - const types::global_dof_index max_n_cells, - const double top_fraction, - const double bottom_fraction) + const types::global_cell_index current_n_cells, + const types::global_cell_index max_n_cells, + const double top_fraction, + const double bottom_fraction) { Assert(top_fraction >= 0, ExcInvalidParameterValue()); Assert(top_fraction <= 1, ExcInvalidParameterValue()); @@ -166,7 +166,7 @@ GridRefinement::adjust_refine_and_coarsen_number_fraction( // again, this is true for isotropically // refined cells. we take this as an // approximation of a mixed refinement. - else if (static_cast( + else if (static_cast( current_n_cells + refine_cells * cell_increase_on_refine - coarsen_cells * cell_decrease_on_coarsen) > max_n_cells) { diff --git a/source/grid/grid_refinement.inst.in b/source/grid/grid_refinement.inst.in index 8bda2e72e1..14452345cd 100644 --- a/source/grid/grid_refinement.inst.in +++ b/source/grid/grid_refinement.inst.in @@ -99,8 +99,8 @@ for (deal_II_dimension : DIMENSIONS) { template std::pair GridRefinement::adjust_refine_and_coarsen_number_fraction< - deal_II_dimension>(const types::global_dof_index, - const types::global_dof_index, + deal_II_dimension>(const types::global_cell_index, + const types::global_cell_index, const double, const double); } diff --git a/source/grid/tria.cc b/source/grid/tria.cc index 354fb50985..87dfd8653c 100644 --- a/source/grid/tria.cc +++ b/source/grid/tria.cc @@ -12695,7 +12695,7 @@ Triangulation::n_active_cells() const } template -types::global_dof_index +types::global_cell_index Triangulation::n_global_active_cells() const { return n_active_cells(); -- 2.39.5