diff --git a/include/deal.II/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index b6974fa198a6..8065efa3106b 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -689,15 +689,13 @@ class DoFHandler : public Subscriptor * function set_fe(). */ void - distribute_dofs(const FiniteElement &fe, - const dealii::types::global_dof_index &virtual_dofs = 0); + distribute_dofs(const FiniteElement &fe); /** * Same as above but taking an hp::FECollection object. */ void - distribute_dofs(const hp::FECollection &fe, - const dealii::types::global_dof_index &virtual_dofs = 0); + distribute_dofs(const hp::FECollection &fe); /** * Distribute level degrees of freedom on each level for geometric @@ -707,6 +705,16 @@ class DoFHandler : public Subscriptor void distribute_mg_dofs(); + /** + * FIXME: documentation + * + * Distribute virtual degrees of freedom. [...] + * + * @pre The locally owned index set must be contiguous. + */ + void + distribute_virtual_dofs(const dealii::types::global_dof_index virtual_dofs); + /** * Returns whether this DoFHandler has hp-capabilities. */ @@ -1190,6 +1198,8 @@ class DoFHandler : public Subscriptor locally_owned_dofs() const; /** + * FIXME: documentation + * * Return an IndexSet describing the subset of locally owned virtual DoFs. */ const IndexSet & diff --git a/include/deal.II/dofs/dof_handler_policy.h b/include/deal.II/dofs/dof_handler_policy.h index fc8c3c75b3d6..4d356976b691 100644 --- a/include/deal.II/dofs/dof_handler_policy.h +++ b/include/deal.II/dofs/dof_handler_policy.h @@ -73,7 +73,7 @@ namespace internal * argument. */ virtual NumberCache - distribute_dofs(const types::global_dof_index virtual_dofs) const = 0; + distribute_dofs() const = 0; /** * Distribute the multigrid dofs on each level of the DoFHandler @@ -83,6 +83,13 @@ namespace internal virtual std::vector distribute_mg_dofs() const = 0; + /** + * FIXME: documentation + */ + virtual NumberCache + distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const = 0; + /** * Renumber degrees of freedom as specified by the first argument. * @@ -90,7 +97,8 @@ namespace internal */ virtual NumberCache renumber_dofs( - const std::vector &new_numbers) const = 0; + const std::vector &new_numbers, + const types::global_dof_index &n_locally_owned_dofs) const = 0; /** * Renumber multilevel degrees of freedom on one level of a multigrid @@ -124,17 +132,21 @@ namespace internal // documentation is inherited virtual NumberCache - distribute_dofs( - const types::global_dof_index virtual_dofs) const override; + distribute_dofs() const override; // documentation is inherited virtual std::vector distribute_mg_dofs() const override; + virtual NumberCache + distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const override; + // documentation is inherited virtual NumberCache - renumber_dofs(const std::vector &new_numbers) - const override; + renumber_dofs( + const std::vector &new_numbers, + const types::global_dof_index &n_locally_owned_dofs) const override; // documentation is inherited virtual NumberCache @@ -175,8 +187,7 @@ namespace internal * number_cache.locally_owned_dofs are updated consistently. */ virtual NumberCache - distribute_dofs( - const types::global_dof_index virtual_dofs) const override; + distribute_dofs() const override; /** * This function is not yet implemented. @@ -184,6 +195,13 @@ namespace internal virtual std::vector distribute_mg_dofs() const override; + /** + * This function is not yet implemented. + */ + virtual NumberCache + distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const override; + /** * Renumber degrees of freedom as specified by the first argument. * @@ -194,8 +212,9 @@ namespace internal * parallel::distributed case. */ virtual NumberCache - renumber_dofs(const std::vector &new_numbers) - const override; + renumber_dofs( + const std::vector &new_numbers, + const types::global_dof_index &n_locally_owned_dofs) const override; // documentation is inherited virtual NumberCache @@ -228,8 +247,7 @@ namespace internal // documentation is inherited virtual NumberCache - distribute_dofs( - const types::global_dof_index virtual_dofs) const override; + distribute_dofs() const override; // documentation is inherited virtual std::vector @@ -237,8 +255,14 @@ namespace internal // documentation is inherited virtual NumberCache - renumber_dofs(const std::vector &new_numbers) - const override; + distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const override; + + // documentation is inherited + virtual NumberCache + renumber_dofs( + const std::vector &new_numbers, + const types::global_dof_index &n_locally_owned_dofs) const override; // documentation is inherited virtual NumberCache diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index f181a308044f..a34c9eec5097 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -2133,10 +2133,9 @@ std::size_t DoFHandler::memory_consumption() const template DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim)) void DoFHandler::distribute_dofs( - const FiniteElement &fe, - const dealii::types::global_dof_index &virtual_dofs) + const FiniteElement &fe) { - this->distribute_dofs(hp::FECollection(fe), virtual_dofs); + this->distribute_dofs(hp::FECollection(fe)); } @@ -2144,8 +2143,7 @@ void DoFHandler::distribute_dofs( template DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim)) void DoFHandler::distribute_dofs( - const hp::FECollection &ff, - const dealii::types::global_dof_index &virtual_dofs) + const hp::FECollection &ff) { Assert(this->tria != nullptr, ExcMessage( @@ -2256,7 +2254,7 @@ void DoFHandler::distribute_dofs( } // hand the actual work over to the policy - this->number_cache = this->policy->distribute_dofs(virtual_dofs); + this->number_cache = this->policy->distribute_dofs(); // do some housekeeping: compress indices // if(hp_capability_enabled) @@ -2313,6 +2311,16 @@ void DoFHandler::distribute_mg_dofs() +template +DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim)) +void DoFHandler::distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) +{ + this->number_cache = this->policy->distribute_virtual_dofs(virtual_dofs); +} + + + template DEAL_II_CXX20_REQUIRES((concepts::is_valid_dim_spacedim)) void DoFHandler::initialize_local_block_info() @@ -2410,11 +2418,10 @@ void DoFHandler::renumber_dofs( // [0...n_dofs()) into itself but only globally, not on each processor if (this->n_locally_owned_dofs() == this->n_dofs()) { - std::vector tmp(new_numbers); + auto tmp = new_numbers; std::sort(tmp.begin(), tmp.end()); - std::vector::const_iterator p = tmp.begin(); - types::global_dof_index i = 0; - for (; p != tmp.end(); ++p, ++i) + types::global_dof_index i = 0; + for (auto p = tmp.begin(); p != tmp.end(); ++p, ++i) Assert(*p == i, ExcNewNumbersNotConsecutive(i)); } else diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index eb4657bf7264..7bc8306364d8 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -2649,8 +2649,7 @@ namespace internal template NumberCache - Sequential::distribute_dofs( - const types::global_dof_index n_virtual_dofs) const + Sequential::distribute_dofs() const { const types::global_dof_index n_initial_dofs = Implementation::distribute_dofs(numbers::invalid_subdomain_id, @@ -2661,7 +2660,34 @@ namespace internal n_initial_dofs, /*check_validity=*/true); + return NumberCache(n_dofs); + } + + + + template + NumberCache + Sequential::distribute_virtual_dofs( + const types::global_dof_index n_virtual_dofs) const + { + const auto &old_locally_owned_dofs = + dof_handler->locally_owned_dofs(); + + Assert(old_locally_owned_dofs.is_contiguous(), + ExcMessage( + "Virtual degrees of freedom can only be distributed when the " + "locally owned index range is contiguous.")); + + const auto &old_locally_owned_virtual_dofs = + dof_handler->locally_owned_virtual_dofs(); + + Assert(old_locally_owned_virtual_dofs.n_elements() == 0, + ExcMessage( + "Can distribute virtual degrees of freedom only once after " + "distribute_dofs()")); + // return a sequential, complete index set + const auto n_dofs = old_locally_owned_dofs.size(); NumberCache number_cache(n_dofs + n_virtual_dofs); number_cache.locally_owned_virtual_dofs = @@ -2703,50 +2729,45 @@ namespace internal template NumberCache Sequential::renumber_dofs( - const std::vector &new_numbers) const + const std::vector &new_numbers, + const dealii::types::global_dof_index &n_total_dofs) const { + AssertDimension(new_numbers.size(), dof_handler->n_dofs()); + Implementation::renumber_dofs(new_numbers, IndexSet(0), *dof_handler, /*check_validity=*/true); - // FIXME begin: refactor into own function - IndexSet my_locally_owned_new_virtual_dofs(dof_handler->n_dofs()); - { - const auto &old_virtual_dofs = - dof_handler->locally_owned_virtual_dofs(); - std::vector new_numbers_sorted( - old_virtual_dofs.n_elements()); - - std::transform(std::begin(old_virtual_dofs), - std::end(old_virtual_dofs), - std::begin(new_numbers_sorted), - [&](const auto &it) { return new_numbers[it]; }); - - std::sort(new_numbers_sorted.begin(), new_numbers_sorted.end()); - - my_locally_owned_new_virtual_dofs.add_indices( - new_numbers_sorted.begin(), new_numbers_sorted.end()); - my_locally_owned_new_virtual_dofs.compress(); - } - // FIXME end - - // return a sequential, complete index set. take into account that the - // number of DoF indices may in fact be smaller than there were before - // if some previously separately numbered dofs have been identified. - // this is, for example, what we do when the DoFHandler has hp- - // capabilities enabled: it first enumerates all DoFs on cells - // independently, and then unifies some located at vertices or faces; - // this leaves us with fewer DoFs than there were before, so use the - // largest index as the one to determine the size of the index space if (new_numbers.empty()) return NumberCache(); else { - const auto n_total_dofs = - *std::max_element(new_numbers.begin(), new_numbers.end()) + 1; +#ifdef DEBUG + const auto largest_index = + *std::max_element(new_numbers.begin(), new_numbers.end()); + Assert(largest_index < n_total_dofs, ExcInternalError()); +#endif + + const auto &old_virtual_dofs = + dof_handler->locally_owned_virtual_dofs(); + std::vector new_virtual_sorted( + old_virtual_dofs.n_elements()); + + std::transform(std::begin(old_virtual_dofs), + std::end(old_virtual_dofs), + std::begin(new_virtual_sorted), + [&](const auto &it) { return new_numbers[it]; }); + + std::sort(new_virtual_sorted.begin(), new_virtual_sorted.end()); + + IndexSet my_locally_owned_new_virtual_dofs(n_total_dofs); + my_locally_owned_new_virtual_dofs.add_indices( + new_virtual_sorted.begin(), new_virtual_sorted.end()); + my_locally_owned_new_virtual_dofs.compress(); // return a sequential, complete index set + NumberCache number_cache(n_total_dofs); number_cache.locally_owned_virtual_dofs = my_locally_owned_new_virtual_dofs; @@ -2920,11 +2941,8 @@ namespace internal template NumberCache - ParallelShared::distribute_dofs( - const types::global_dof_index virtual_dofs [[maybe_unused]]) const + ParallelShared::distribute_dofs() const { - Assert(virtual_dofs == 0, ExcNotImplemented()); - const dealii::parallel::shared::Triangulation *tr = (dynamic_cast< const dealii::parallel::shared::Triangulation *>( @@ -3066,6 +3084,18 @@ namespace internal + template + NumberCache + ParallelShared::distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const + { + // FIXME: implement me + AssertThrow(virtual_dofs == 0, ExcNotImplemented()); + return NumberCache(); + } + + + template std::vector ParallelShared::distribute_mg_dofs() const @@ -3644,8 +3674,7 @@ namespace internal template NumberCache - ParallelDistributed::distribute_dofs( - const types::global_dof_index n_locally_virtual [[maybe_unused]]) const + ParallelDistributed::distribute_dofs() const { #ifndef DEAL_II_WITH_MPI DEAL_II_NOT_IMPLEMENTED(); @@ -3741,13 +3770,9 @@ namespace internal Utilities::MPI::partial_and_total_sum( n_locally_owned_dofs, triangulation->get_communicator()); - const auto [my_shift_virtual, n_global_virtual] = - Utilities::MPI::partial_and_total_sum( - n_locally_virtual, triangulation->get_communicator()); - // make dof indices globally consecutive Implementation::enumerate_dof_indices_for_renumbering( - renumbering, all_constrained_indices, my_shift + my_shift_virtual); + renumbering, all_constrained_indices, my_shift); // now re-enumerate all dofs to this shifted and condensed // numbering form. we renumber some dofs as invalid, so @@ -3757,26 +3782,17 @@ namespace internal *dof_handler, /*check_validity=*/false); - NumberCache number_cache; - number_cache.n_global_dofs = n_global_dofs + n_global_virtual; - number_cache.n_locally_owned_dofs = - n_locally_owned_dofs + n_locally_virtual; + number_cache.n_global_dofs = n_global_dofs; + number_cache.n_locally_owned_dofs = n_locally_owned_dofs; - number_cache.locally_owned_dofs = - IndexSet(n_global_dofs + n_global_virtual); + number_cache.locally_owned_dofs = IndexSet(n_global_dofs); number_cache.locally_owned_dofs.add_range( // - my_shift + my_shift_virtual, - my_shift + my_shift_virtual + n_locally_owned_dofs + - n_locally_virtual); + my_shift, + my_shift + n_locally_owned_dofs); number_cache.locally_owned_dofs.compress(); - number_cache.locally_owned_virtual_dofs = - IndexSet(n_global_dofs + n_global_virtual); - number_cache.locally_owned_virtual_dofs.add_range( - my_shift + my_shift_virtual + n_locally_owned_dofs, - my_shift + my_shift_virtual + n_locally_owned_dofs + - n_locally_virtual); + number_cache.locally_owned_virtual_dofs = IndexSet(n_global_dofs); number_cache.locally_owned_virtual_dofs.compress(); @@ -3868,6 +3884,70 @@ namespace internal + template + NumberCache + ParallelDistributed::distribute_virtual_dofs( + const types::global_dof_index virtual_dofs) const + { + const auto &old_locally_owned_dofs = + dof_handler->locally_owned_dofs(); + + Assert(old_locally_owned_dofs.is_contiguous(), + ExcMessage( + "Virtual degrees of freedom can only be distributed when the " + "locally owned index range is contiguous.")); + + const auto &old_locally_owned_virtual_dofs = + dof_handler->locally_owned_virtual_dofs(); + + Assert(old_locally_owned_virtual_dofs.n_elements() == 0, + ExcMessage( + "Can distribute virtual degrees of freedom only once after " + "distribute_dofs()")); + + dealii::parallel::DistributedTriangulationBase + *triangulation = + (dynamic_cast< + dealii::parallel::DistributedTriangulationBase *>( + const_cast *>( + &dof_handler->get_triangulation()))); + Assert(triangulation != nullptr, ExcInternalError()); + + // + // Renumber degrees of freedom by adding a global shift to the + // locally owned index range: + // + + const auto [my_shift, n_global_virtual_dofs] = + Utilities::MPI::partial_and_total_sum( + virtual_dofs, triangulation->get_communicator()); + + std::vector renumbering( + old_locally_owned_dofs.n_elements()); + std::generate(renumbering.begin(), + renumbering.end(), + [my_shift = my_shift, + n = *old_locally_owned_dofs.begin()]() mutable { + return my_shift + n++; + }); + + NumberCache number_cache = this->renumber_dofs(renumbering); + + return number_cache; + +#if 0 + // FIXME resume here + number_cache.locally_owned_virtual_dofs = + number_cache.locally_owned_virtual_dofs.add_range( + my_shift + my_shift_virtual + n_locally_owned_dofs, + my_shift + my_shift_virtual + n_locally_owned_dofs + + n_locally_virtual); + number_cache.locally_owned_virtual_dofs.compress(); +#endif + } + + + template std::vector ParallelDistributed::distribute_mg_dofs() const @@ -4047,10 +4127,9 @@ namespace internal template NumberCache ParallelDistributed::renumber_dofs( - const std::vector &new_numbers) const + const std::vector &new_numbers, + const dealii::types::global_dof_index &n_locally_owned_dofs) const { - (void)new_numbers; - Assert(new_numbers.size() == dof_handler->n_locally_owned_dofs(), ExcInternalError()); @@ -4072,12 +4151,13 @@ namespace internal // ranks changed, in which case we do not need to find a new index // set. const IndexSet &owned_dofs = dof_handler->locally_owned_dofs(); - const bool locally_owned_set_changes = + bool locally_owned_set_changes = std::any_of(new_numbers.cbegin(), new_numbers.cend(), [&owned_dofs](const types::global_dof_index i) { return owned_dofs.is_element(i) == false; - }); + }) || + true; IndexSet my_locally_owned_new_dof_indices = owned_dofs; if (locally_owned_set_changes && owned_dofs.n_elements() > 0)