From 0728ffec67b3908358960206857a96269364e25c Mon Sep 17 00:00:00 2001 From: Marc Fehling Date: Thu, 12 Dec 2019 14:55:46 +0100 Subject: [PATCH] Introduced hp::Refinement::p_adaptivity_fixed_number(). Decide how many cells will be p-adapted from all cells flagged for adaptation. --- include/deal.II/hp/refinement.h | 91 ++++++-- source/hp/refinement.cc | 215 +++++++++++++++++- source/hp/refinement.inst.in | 11 + tests/mpi/p_adaptivity_fixed_number.cc | 152 +++++++++++++ ...xed_number.with_p4est=true.mpirun=1.output | 11 + ...xed_number.with_p4est=true.mpirun=2.output | 23 ++ ...xed_number.with_p4est=true.mpirun=8.output | 95 ++++++++ 7 files changed, 571 insertions(+), 27 deletions(-) create mode 100644 tests/mpi/p_adaptivity_fixed_number.cc create mode 100644 tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=1.output create mode 100644 tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=2.output create mode 100644 tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=8.output diff --git a/include/deal.II/hp/refinement.h b/include/deal.II/hp/refinement.h index d6cc090c0b..cb269b75e7 100644 --- a/include/deal.II/hp/refinement.h +++ b/include/deal.II/hp/refinement.h @@ -150,9 +150,9 @@ namespace hp * Each cell flagged for h-refinement will also be flagged for p-refinement. * The same applies to coarsening. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -167,9 +167,9 @@ namespace hp * Each entry of the parameter @p p_flags needs to correspond to an active * cell. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -196,9 +196,9 @@ namespace hp * Each entry of the parameter @p criteria needs to correspond to an active * cell. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -238,9 +238,9 @@ namespace hp * cell. Parameters @p p_refine_fraction and @p p_coarsen_fraction need to be * in the interval $[0,1]$. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -254,6 +254,49 @@ namespace hp const ComparisonFunction::type> &compare_coarsen = std::less_equal()); + /** + * Adapt which finite element to use on a given fraction of cells. + * + * Out of all cells flagged for a certain type of adaptation, be it + * refinement or coarsening, we will determine a fixed number of cells among + * this subset that will be flagged for the corresponding p-adaptive + * variant. + * + * For each of both refinement and coarsening subsets, we will determine a + * threshold based on the provided parameter @p criteria containing + * indicators for every active cell. In the default case for refinement, all + * cells with an indicator larger than or equal to the corresponding + * threshold will be considered for p-refinement, while for coarsening all + * cells with an indicator less than or equal to the matching threshold are + * taken into account. However, different compare function objects can be + * supplied via the parameters @p compare_refine and @p compare_coarsen to + * impose different decision strategies. + * + * For refinement, the threshold will be associated with the cell that has + * the @p p_refine_fraction times Triangulation::n_active_cells() largest + * indicator, while it is the cell with the @p p_refine_coarsen times + * Triangulation::n_active_cells() lowest indicator for coarsening. + * + * Each entry of the parameter @p criteria needs to correspond to an active + * cell. Parameters @p p_refine_fraction and @p p_coarsen_fraction need to be + * in the interval $[0,1]$. + * + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. + */ + template + void + p_adaptivity_fixed_number( + const hp::DoFHandler &dof_handler, + const Vector & criteria, + const double p_refine_fraction = 0.5, + const double p_coarsen_fraction = 0.5, + const ComparisonFunction::type> + &compare_refine = std::greater_equal(), + const ComparisonFunction::type> + &compare_coarsen = std::less_equal()); + /** * Adapt which finite element to use on cells based on the regularity of the * (unknown) analytical solution. @@ -289,9 +332,9 @@ namespace hp * } * @endcode * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -314,9 +357,9 @@ namespace hp * Each entry of the parameters @p criteria and @p references needs to * correspond to an active cell. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -515,9 +558,9 @@ namespace hp * Removes all refine and coarsen flags on cells that have a * @p future_fe_index assigned. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void @@ -559,9 +602,9 @@ namespace hp * the decision that Triangulation::prepare_coarsening_and_refinement() * would have made later on. * - * @note Preceeding calls of Triangulation::prepare_for_coarsening_and_refinement() - * may change refine and coarsen flags, which will ultimately change the - * results of this function. + * @note Triangulation::prepare_for_coarsening_and_refinement() may change + * refine and coarsen flags. Avoid calling it before this particular + * function. */ template void diff --git a/source/hp/refinement.cc b/source/hp/refinement.cc index 491206ec34..6004b53780 100644 --- a/source/hp/refinement.cc +++ b/source/hp/refinement.cc @@ -18,6 +18,8 @@ #include +#include +#include #include #include @@ -29,6 +31,25 @@ DEAL_II_NAMESPACE_OPEN +namespace +{ + /** + * ComparisonFunction returning 'true' or 'false' for any set of parameters. + * + * These will be used to overwrite user-provided comparison functions whenever + * no actual comparison is required in the decision process, i.e. when no or + * all cells will be refined or coarsened. + */ + template + hp::Refinement::ComparisonFunction compare_false = + [](const Number &, const Number &) { return false; }; + template + hp::Refinement::ComparisonFunction compare_true = + [](const Number &, const Number &) { return true; }; +} // namespace + + + namespace hp { namespace Refinement @@ -168,9 +189,12 @@ namespace hp } } - if (const parallel::TriangulationBase *parallel_tria = - dynamic_cast *>( - &dof_handler.get_triangulation())) + const parallel::TriangulationBase *parallel_tria = + dynamic_cast *>( + &dof_handler.get_triangulation()); + if (parallel_tria != nullptr && + dynamic_cast *>( + &dof_handler.get_triangulation()) == nullptr) { max_criterion_refine = Utilities::MPI::max(max_criterion_refine, @@ -207,6 +231,191 @@ namespace hp + template + void + p_adaptivity_fixed_number( + const hp::DoFHandler &dof_handler, + const Vector & criteria, + const double p_refine_fraction, + const double p_coarsen_fraction, + const ComparisonFunction::type> &compare_refine, + const ComparisonFunction::type> + &compare_coarsen) + { + AssertDimension(dof_handler.get_triangulation().n_active_cells(), + criteria.size()); + Assert((p_refine_fraction >= 0) && (p_refine_fraction <= 1), + dealii::GridRefinement::ExcInvalidParameterValue()); + Assert((p_coarsen_fraction >= 0) && (p_coarsen_fraction <= 1), + dealii::GridRefinement::ExcInvalidParameterValue()); + + // 1.) First extract from the vector of indicators the ones that + // correspond to cells that we locally own. + unsigned int n_flags_refinement = 0; + unsigned int n_flags_coarsening = 0; + Vector indicators_refinement( + dof_handler.get_triangulation().n_active_cells()); + Vector indicators_coarsening( + dof_handler.get_triangulation().n_active_cells()); + for (const auto &cell : + dof_handler.get_triangulation().active_cell_iterators()) + if (!cell->is_artificial() && cell->is_locally_owned()) + { + if (cell->refine_flag_set()) + indicators_refinement(n_flags_refinement++) = + criteria(cell->active_cell_index()); + else if (cell->coarsen_flag_set()) + indicators_coarsening(n_flags_coarsening++) = + criteria(cell->active_cell_index()); + } + indicators_refinement.grow_or_shrink(n_flags_refinement); + indicators_coarsening.grow_or_shrink(n_flags_coarsening); + + // 2.) Determine the number of cells for p-refinement and p-coarsening on + // basis of the flagged cells. + // + // 3.) Find thresholds for p-refinment and p-coarsening on only those + // cells flagged for adaptation. + // + // For cases in which no or all cells flagged for refinement and/or + // coarsening are subject to p-adaptation, we usually pick thresholds + // that apply to all or none of the cells at once. However here, we + // do not know which threshold would suffice for this task because the + // user could provide any comparison function. Thus if necessary, we + // overwrite the user's choice with suitable functions simplying + // returning 'true' and 'false' for any cell with reference wrappers. + // Thus, no function object copies are stored. + // + // 4.) Perform p-adaptation with absolute thresholds. + Number threshold_refinement = 0.; + Number threshold_coarsening = 0.; + auto reference_compare_refine = std::cref(compare_refine); + auto reference_compare_coarsen = std::cref(compare_coarsen); + + const parallel::TriangulationBase *parallel_tria = + dynamic_cast *>( + &dof_handler.get_triangulation()); + if (parallel_tria != nullptr && + dynamic_cast *>( + &dof_handler.get_triangulation()) == nullptr) + { +#ifndef DEAL_II_WITH_P4EST + Assert(false, ExcInternalError()); +#else + // + // parallel implementation with distributed memory + // + + MPI_Comm mpi_communicator = parallel_tria->get_communicator(); + + // 2.) Communicate the number of cells scheduled for p-adaptation + // globally. + const unsigned int n_global_flags_refinement = + Utilities::MPI::sum(n_flags_refinement, mpi_communicator); + const unsigned int n_global_flags_coarsening = + Utilities::MPI::sum(n_flags_coarsening, mpi_communicator); + + const unsigned int target_index_refinement = + static_cast( + std::floor(p_refine_fraction * n_global_flags_refinement)); + const unsigned int target_index_coarsening = + static_cast( + std::ceil((1 - p_coarsen_fraction) * n_global_flags_coarsening)); + + // 3.) Figure out the global max and min of the criteria. We don't + // need it here, but it's a collective communication call. + const std::pair global_min_max_refinement = + dealii::internal::parallel::distributed::GridRefinement:: + compute_global_min_and_max_at_root(indicators_refinement, + mpi_communicator); + + const std::pair global_min_max_coarsening = + dealii::internal::parallel::distributed::GridRefinement:: + compute_global_min_and_max_at_root(indicators_coarsening, + mpi_communicator); + + // 3.) Compute thresholds if necessary. + if (target_index_refinement == 0) + reference_compare_refine = std::cref(compare_false); + else if (target_index_refinement == n_global_flags_refinement) + reference_compare_refine = std::cref(compare_true); + else + threshold_refinement = dealii::internal::parallel::distributed:: + GridRefinement::RefineAndCoarsenFixedNumber::compute_threshold( + indicators_refinement, + global_min_max_refinement, + target_index_refinement, + mpi_communicator); + + if (target_index_coarsening == n_global_flags_coarsening) + reference_compare_coarsen = std::cref(compare_false); + else if (target_index_coarsening == 0) + reference_compare_coarsen = std::cref(compare_true); + else + threshold_coarsening = dealii::internal::parallel::distributed:: + GridRefinement::RefineAndCoarsenFixedNumber::compute_threshold( + indicators_coarsening, + global_min_max_coarsening, + target_index_coarsening, + mpi_communicator); +#endif + } + else + { + // + // serial implementation (and parallel::shared implementation) + // + + // 2.) Determine the number of cells scheduled for p-adaptation. + const unsigned int n_p_refine_cells = static_cast( + std::floor(p_refine_fraction * n_flags_refinement)); + const unsigned int n_p_coarsen_cells = static_cast( + std::floor(p_coarsen_fraction * n_flags_coarsening)); + + // 3.) Compute thresholds if necessary. + if (n_p_refine_cells == 0) + reference_compare_refine = std::cref(compare_false); + else if (n_p_refine_cells == n_flags_refinement) + reference_compare_refine = std::cref(compare_true); + else + { + std::nth_element(indicators_refinement.begin(), + indicators_refinement.begin() + + n_p_refine_cells - 1, + indicators_refinement.end(), + std::greater()); + threshold_refinement = + *(indicators_refinement.begin() + n_p_refine_cells - 1); + } + + if (n_p_coarsen_cells == 0) + reference_compare_coarsen = std::cref(compare_false); + else if (n_p_coarsen_cells == n_flags_coarsening) + reference_compare_coarsen = std::cref(compare_true); + else + { + std::nth_element(indicators_coarsening.begin(), + indicators_coarsening.begin() + + n_p_coarsen_cells - 1, + indicators_coarsening.end(), + std::less()); + threshold_coarsening = + *(indicators_coarsening.begin() + n_p_coarsen_cells - 1); + } + } + + // 4.) Finally perform adaptation. + p_adaptivity_from_absolute_threshold(dof_handler, + criteria, + threshold_refinement, + threshold_coarsening, + std::cref(reference_compare_refine), + std::cref( + reference_compare_coarsen)); + } + + + template void p_adaptivity_from_regularity( diff --git a/source/hp/refinement.inst.in b/source/hp/refinement.inst.in index fdbf312fe8..031966e669 100644 --- a/source/hp/refinement.inst.in +++ b/source/hp/refinement.inst.in @@ -73,6 +73,17 @@ for (deal_II_dimension : DIMENSIONS; S : REAL_SCALARS; const ComparisonFunction::type> &, const ComparisonFunction::type> &); + template void + p_adaptivity_fixed_number( + const hp::DoFHandler &, + const Vector &, + const double, + const double, + const ComparisonFunction::type> &, + const ComparisonFunction::type> &); + template void p_adaptivity_from_regularity + +#include + +#include + +#include + +#include +#include +#include + +#include + +#include "../tests.h" + + + +template +void +validate(const hp::DoFHandler &dh) +{ + deallog << " (cellid,feidx):"; + for (const auto &cell : dh.active_cell_iterators()) + if (!cell->is_artificial() && cell->is_locally_owned()) + { + const std::string cellid = cell->id().to_string(); + const unsigned int coarse_cellid = + std::stoul(cellid.substr(0, cellid.find("_"))); + + deallog << " (" << coarse_cellid << "," << cell->future_fe_index() + << ")"; + } + deallog << std::endl; +} + + + +template +void +setup(Triangulation & tria, + hp::DoFHandler & dh, + const hp::FECollection &fes) +{ + // Initialize triangulation and dofhandler. + GridGenerator::subdivided_hyper_cube(tria, 4); + Assert(tria.n_cells(0) == tria.n_global_active_cells(), ExcInternalError()); + + dh.initialize(tria, fes); + + // Set all active fe indices to 1. + // Flag first half of cells for refinement, and the other half for coarsening. + for (const auto &cell : dh.active_cell_iterators()) + if (!cell->is_artificial() && cell->is_locally_owned()) + { + cell->set_active_fe_index(1); + + const std::string cellid = cell->id().to_string(); + const unsigned int coarse_cellid = + std::stoul(cellid.substr(0, cellid.find('_'))); + + if (coarse_cellid < 0.5 * tria.n_global_active_cells()) + cell->set_refine_flag(); + else + cell->set_coarsen_flag(); + } +} + + + +template +void +test() +{ + hp::FECollection fes; + for (unsigned int d = 1; d <= 3; ++d) + fes.push_back(FE_Q(d)); + + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + hp::DoFHandler dh; + setup(tria, dh, fes); + + deallog << "starting situation" << std::endl; + validate(dh); + + + // We flag the first half of all cells to be refined and the last half of all + // cells to be coarsened for p adapativity. Ultimately, the first quarter of + // all cells will be flagged for p refinement, and the last quarter for p + // coarsening. + + Vector indicators(tria.n_active_cells()); + for (const auto &cell : tria.active_cell_iterators()) + if (!cell->is_artificial() && cell->is_locally_owned()) + { + const std::string cellid = cell->id().to_string(); + const unsigned int coarse_cellid = + std::stoul(cellid.substr(0, cellid.find('_'))); + + if (coarse_cellid < .25 * tria.n_global_active_cells()) + indicators[cell->active_cell_index()] = 2.; + else if (coarse_cellid < .75 * tria.n_global_active_cells()) + indicators[cell->active_cell_index()] = 1.; + else + indicators[cell->active_cell_index()] = 0.; + } + + hp::Refinement::p_adaptivity_fixed_number(dh, indicators); + + deallog << "p-adaptivity fixed number" << std::endl; + validate(dh); + + + deallog << "OK" << std::endl; +} + + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + deallog.push("2d"); + test<2>(); + deallog.pop(); + deallog.push("3d"); + test<3>(); + deallog.pop(); +} diff --git a/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=1.output b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=1.output new file mode 100644 index 0000000000..54a72f0b66 --- /dev/null +++ b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=1.output @@ -0,0 +1,11 @@ + +DEAL:0:2d::starting situation +DEAL:0:2d:: (cellid,feidx): (0,1) (1,1) (4,1) (5,1) (2,1) (3,1) (6,1) (7,1) (8,1) (9,1) (12,1) (13,1) (10,1) (11,1) (14,1) (15,1) +DEAL:0:2d::p-adaptivity fixed number +DEAL:0:2d:: (cellid,feidx): (0,2) (1,2) (4,1) (5,1) (2,2) (3,2) (6,1) (7,1) (8,1) (9,1) (12,0) (13,0) (10,1) (11,1) (14,0) (15,0) +DEAL:0:2d::OK +DEAL:0:3d::starting situation +DEAL:0:3d:: (cellid,feidx): (0,1) (1,1) (8,1) (9,1) (2,1) (3,1) (10,1) (11,1) (16,1) (17,1) (24,1) (25,1) (18,1) (19,1) (26,1) (27,1) (4,1) (5,1) (12,1) (13,1) (6,1) (7,1) (14,1) (15,1) (20,1) (21,1) (28,1) (29,1) (22,1) (23,1) (30,1) (31,1) (32,1) (33,1) (40,1) (41,1) (34,1) (35,1) (42,1) (43,1) (48,1) (49,1) (56,1) (57,1) (50,1) (51,1) (58,1) (59,1) (36,1) (37,1) (44,1) (45,1) (38,1) (39,1) (46,1) (47,1) (52,1) (53,1) (60,1) (61,1) (54,1) (55,1) (62,1) (63,1) +DEAL:0:3d::p-adaptivity fixed number +DEAL:0:3d:: (cellid,feidx): (0,2) (1,2) (8,2) (9,2) (2,2) (3,2) (10,2) (11,2) (16,1) (17,1) (24,1) (25,1) (18,1) (19,1) (26,1) (27,1) (4,2) (5,2) (12,2) (13,2) (6,2) (7,2) (14,2) (15,2) (20,1) (21,1) (28,1) (29,1) (22,1) (23,1) (30,1) (31,1) (32,1) (33,1) (40,1) (41,1) (34,1) (35,1) (42,1) (43,1) (48,0) (49,0) (56,0) (57,0) (50,0) (51,0) (58,0) (59,0) (36,1) (37,1) (44,1) (45,1) (38,1) (39,1) (46,1) (47,1) (52,0) (53,0) (60,0) (61,0) (54,0) (55,0) (62,0) (63,0) +DEAL:0:3d::OK diff --git a/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=2.output b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..b355ad91b2 --- /dev/null +++ b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=2.output @@ -0,0 +1,23 @@ + +DEAL:0:2d::starting situation +DEAL:0:2d:: (cellid,feidx): (0,1) (1,1) (4,1) (5,1) (2,1) (3,1) (6,1) (7,1) +DEAL:0:2d::p-adaptivity fixed number +DEAL:0:2d:: (cellid,feidx): (0,2) (1,2) (4,1) (5,1) (2,2) (3,2) (6,1) (7,1) +DEAL:0:2d::OK +DEAL:0:3d::starting situation +DEAL:0:3d:: (cellid,feidx): (0,1) (1,1) (8,1) (9,1) (2,1) (3,1) (10,1) (11,1) (16,1) (17,1) (24,1) (25,1) (18,1) (19,1) (26,1) (27,1) (4,1) (5,1) (12,1) (13,1) (6,1) (7,1) (14,1) (15,1) (20,1) (21,1) (28,1) (29,1) (22,1) (23,1) (30,1) (31,1) +DEAL:0:3d::p-adaptivity fixed number +DEAL:0:3d:: (cellid,feidx): (0,2) (1,2) (8,2) (9,2) (2,2) (3,2) (10,2) (11,2) (16,1) (17,1) (24,1) (25,1) (18,1) (19,1) (26,1) (27,1) (4,2) (5,2) (12,2) (13,2) (6,2) (7,2) (14,2) (15,2) (20,1) (21,1) (28,1) (29,1) (22,1) (23,1) (30,1) (31,1) +DEAL:0:3d::OK + +DEAL:1:2d::starting situation +DEAL:1:2d:: (cellid,feidx): (8,1) (9,1) (12,1) (13,1) (10,1) (11,1) (14,1) (15,1) +DEAL:1:2d::p-adaptivity fixed number +DEAL:1:2d:: (cellid,feidx): (8,1) (9,1) (12,0) (13,0) (10,1) (11,1) (14,0) (15,0) +DEAL:1:2d::OK +DEAL:1:3d::starting situation +DEAL:1:3d:: (cellid,feidx): (32,1) (33,1) (40,1) (41,1) (34,1) (35,1) (42,1) (43,1) (48,1) (49,1) (56,1) (57,1) (50,1) (51,1) (58,1) (59,1) (36,1) (37,1) (44,1) (45,1) (38,1) (39,1) (46,1) (47,1) (52,1) (53,1) (60,1) (61,1) (54,1) (55,1) (62,1) (63,1) +DEAL:1:3d::p-adaptivity fixed number +DEAL:1:3d:: (cellid,feidx): (32,1) (33,1) (40,1) (41,1) (34,1) (35,1) (42,1) (43,1) (48,0) (49,0) (56,0) (57,0) (50,0) (51,0) (58,0) (59,0) (36,1) (37,1) (44,1) (45,1) (38,1) (39,1) (46,1) (47,1) (52,0) (53,0) (60,0) (61,0) (54,0) (55,0) (62,0) (63,0) +DEAL:1:3d::OK + diff --git a/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=8.output b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=8.output new file mode 100644 index 0000000000..bcf60d347e --- /dev/null +++ b/tests/mpi/p_adaptivity_fixed_number.with_p4est=true.mpirun=8.output @@ -0,0 +1,95 @@ + +DEAL:0:2d::starting situation +DEAL:0:2d:: (cellid,feidx): (0,1) (1,1) +DEAL:0:2d::p-adaptivity fixed number +DEAL:0:2d:: (cellid,feidx): (0,2) (1,2) +DEAL:0:2d::OK +DEAL:0:3d::starting situation +DEAL:0:3d:: (cellid,feidx): (0,1) (1,1) (2,1) (3,1) (4,1) (5,1) (6,1) (7,1) +DEAL:0:3d::p-adaptivity fixed number +DEAL:0:3d:: (cellid,feidx): (0,2) (1,2) (2,2) (3,2) (4,2) (5,2) (6,2) (7,2) +DEAL:0:3d::OK + +DEAL:1:2d::starting situation +DEAL:1:2d:: (cellid,feidx): (2,1) (3,1) +DEAL:1:2d::p-adaptivity fixed number +DEAL:1:2d:: (cellid,feidx): (2,2) (3,2) +DEAL:1:2d::OK +DEAL:1:3d::starting situation +DEAL:1:3d:: (cellid,feidx): (8,1) (9,1) (10,1) (11,1) (12,1) (13,1) (14,1) (15,1) +DEAL:1:3d::p-adaptivity fixed number +DEAL:1:3d:: (cellid,feidx): (8,2) (9,2) (10,2) (11,2) (12,2) (13,2) (14,2) (15,2) +DEAL:1:3d::OK + + +DEAL:2:2d::starting situation +DEAL:2:2d:: (cellid,feidx): (4,1) (5,1) +DEAL:2:2d::p-adaptivity fixed number +DEAL:2:2d:: (cellid,feidx): (4,1) (5,1) +DEAL:2:2d::OK +DEAL:2:3d::starting situation +DEAL:2:3d:: (cellid,feidx): (16,1) (17,1) (18,1) (19,1) (20,1) (21,1) (22,1) (23,1) +DEAL:2:3d::p-adaptivity fixed number +DEAL:2:3d:: (cellid,feidx): (16,1) (17,1) (18,1) (19,1) (20,1) (21,1) (22,1) (23,1) +DEAL:2:3d::OK + + +DEAL:3:2d::starting situation +DEAL:3:2d:: (cellid,feidx): (6,1) (7,1) +DEAL:3:2d::p-adaptivity fixed number +DEAL:3:2d:: (cellid,feidx): (6,1) (7,1) +DEAL:3:2d::OK +DEAL:3:3d::starting situation +DEAL:3:3d:: (cellid,feidx): (24,1) (25,1) (26,1) (27,1) (28,1) (29,1) (30,1) (31,1) +DEAL:3:3d::p-adaptivity fixed number +DEAL:3:3d:: (cellid,feidx): (24,1) (25,1) (26,1) (27,1) (28,1) (29,1) (30,1) (31,1) +DEAL:3:3d::OK + + +DEAL:4:2d::starting situation +DEAL:4:2d:: (cellid,feidx): (8,1) (9,1) +DEAL:4:2d::p-adaptivity fixed number +DEAL:4:2d:: (cellid,feidx): (8,1) (9,1) +DEAL:4:2d::OK +DEAL:4:3d::starting situation +DEAL:4:3d:: (cellid,feidx): (32,1) (33,1) (34,1) (35,1) (36,1) (37,1) (38,1) (39,1) +DEAL:4:3d::p-adaptivity fixed number +DEAL:4:3d:: (cellid,feidx): (32,1) (33,1) (34,1) (35,1) (36,1) (37,1) (38,1) (39,1) +DEAL:4:3d::OK + + +DEAL:5:2d::starting situation +DEAL:5:2d:: (cellid,feidx): (10,1) (11,1) +DEAL:5:2d::p-adaptivity fixed number +DEAL:5:2d:: (cellid,feidx): (10,1) (11,1) +DEAL:5:2d::OK +DEAL:5:3d::starting situation +DEAL:5:3d:: (cellid,feidx): (40,1) (41,1) (42,1) (43,1) (44,1) (45,1) (46,1) (47,1) +DEAL:5:3d::p-adaptivity fixed number +DEAL:5:3d:: (cellid,feidx): (40,1) (41,1) (42,1) (43,1) (44,1) (45,1) (46,1) (47,1) +DEAL:5:3d::OK + + +DEAL:6:2d::starting situation +DEAL:6:2d:: (cellid,feidx): (12,1) (13,1) +DEAL:6:2d::p-adaptivity fixed number +DEAL:6:2d:: (cellid,feidx): (12,0) (13,0) +DEAL:6:2d::OK +DEAL:6:3d::starting situation +DEAL:6:3d:: (cellid,feidx): (48,1) (49,1) (50,1) (51,1) (52,1) (53,1) (54,1) (55,1) +DEAL:6:3d::p-adaptivity fixed number +DEAL:6:3d:: (cellid,feidx): (48,0) (49,0) (50,0) (51,0) (52,0) (53,0) (54,0) (55,0) +DEAL:6:3d::OK + + +DEAL:7:2d::starting situation +DEAL:7:2d:: (cellid,feidx): (14,1) (15,1) +DEAL:7:2d::p-adaptivity fixed number +DEAL:7:2d:: (cellid,feidx): (14,0) (15,0) +DEAL:7:2d::OK +DEAL:7:3d::starting situation +DEAL:7:3d:: (cellid,feidx): (56,1) (57,1) (58,1) (59,1) (60,1) (61,1) (62,1) (63,1) +DEAL:7:3d::p-adaptivity fixed number +DEAL:7:3d:: (cellid,feidx): (56,0) (57,0) (58,0) (59,0) (60,0) (61,0) (62,0) (63,0) +DEAL:7:3d::OK + -- 2.39.5