From 38d6d947be4704eb737f93485b2cc5fe3eb1ab4f Mon Sep 17 00:00:00 2001 From: heister Date: Thu, 25 Aug 2011 14:25:14 +0000 Subject: [PATCH] fix DoFRenumbering::hierarchical for a distributed Triangulation git-svn-id: https://svn.dealii.org/trunk@24194 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/distributed/tria.h | 15 ++++++++ deal.II/source/distributed/tria.cc | 16 +++++++++ deal.II/source/dofs/dof_renumbering.cc | 42 +++++++++++++++++++--- 3 files changed, 68 insertions(+), 5 deletions(-) diff --git a/deal.II/include/deal.II/distributed/tria.h b/deal.II/include/deal.II/distributed/tria.h index 888afd1d62..13c877498c 100644 --- a/deal.II/include/deal.II/distributed/tria.h +++ b/deal.II/include/deal.II/distributed/tria.h @@ -495,6 +495,16 @@ namespace parallel const CellStatus, const void*)> & unpack_callback); + /** + * Returns a permutation vector for the order the coarse + * cells are handed of to p4est. For example the first + * element i in this vector denotes that the first cell + * in hierarchical ordering is the ith deal cell starting + * from begin(0). + */ + const std::vector & + get_p4est_tree_to_coarse_cell_permutation() const; + private: /** * MPI communicator to be @@ -699,6 +709,11 @@ namespace parallel */ MPI_Comm get_communicator () const; + /** + * */ + const std::vector & + get_p4est_tree_to_coarse_cell_permutation() const; + /** * Return the subdomain id of * those cells that are owned diff --git a/deal.II/source/distributed/tria.cc b/deal.II/source/distributed/tria.cc index 47191ee625..a548df84de 100644 --- a/deal.II/source/distributed/tria.cc +++ b/deal.II/source/distributed/tria.cc @@ -2803,6 +2803,15 @@ namespace parallel + template + const std::vector & + Triangulation::get_p4est_tree_to_coarse_cell_permutation() const + { + return p4est_tree_to_coarse_cell_permutation; + } + + + template MPI_Comm Triangulation::get_communicator () const @@ -2943,6 +2952,13 @@ namespace parallel return MPI_COMM_WORLD; } + template <> + const std::vector & + Triangulation<1,1>::get_p4est_tree_to_coarse_cell_permutation() const + { + static std::vector a; + return a; + } template <> diff --git a/deal.II/source/dofs/dof_renumbering.cc b/deal.II/source/dofs/dof_renumbering.cc index 7af3b620e2..9f306b1a53 100644 --- a/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/source/dofs/dof_renumbering.cc @@ -961,11 +961,42 @@ namespace DoFRenumbering unsigned int next_free = 0; const IndexSet locally_owned = dof_handler.locally_owned_dofs(); - for(cell = dof_handler.begin(0); cell != dof_handler.end(0); ++cell) - next_free = compute_hierarchical_recursive (next_free, - renumbering, - cell, - locally_owned); + const parallel::distributed::Triangulation * tria + = dynamic_cast*> + (&dof_handler.get_tria()); + + if (tria) + { +#ifdef DEAL_II_USE_P4EST + //this is a distributed Triangulation. We need to traverse the coarse + //cells in the order p4est does + for (unsigned int c = 0; c < tria->n_cells (0); ++c) + { + unsigned int coarse_cell_index = + tria->get_p4est_tree_to_coarse_cell_permutation() [c]; + + const typename DoFHandler::cell_iterator + cell (tria, 0, coarse_cell_index, &dof_handler); + + next_free = compute_hierarchical_recursive (next_free, + renumbering, + cell, + locally_owned); + } +#else + Assert (false, ExcNotImplemented()); +#endif + } + else + { + //this is not a distributed Triangulation. Traverse coarse cells in the + //normal order + for (cell = dof_handler.begin (0); cell != dof_handler.end (0); ++cell) + next_free = compute_hierarchical_recursive (next_free, + renumbering, + cell, + locally_owned); + } // verify that the last numbered // degree of freedom is either @@ -986,6 +1017,7 @@ namespace DoFRenumbering numbers::invalid_unsigned_int) == renumbering.end(), ExcInternalError()); + dof_handler.renumber_dofs(renumbering); } -- 2.39.5