From: Wolfgang Bangerth Date: Thu, 22 Jun 2017 15:17:14 +0000 (-0600) Subject: Let the DoFHandler policy classes store references to the DoFHandler. X-Git-Tag: v9.0.0-rc1~1475^2~6 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=254bdc078deb10510d2bbe267ba634149453cd65;p=dealii.git Let the DoFHandler policy classes store references to the DoFHandler. This is rather than passing this reference to each and every call to the policy functions. The primary goal of this is so that we can use the Policy base class not only for the 'normal' DoFHandler, but also for the hp::DoFHandler by removing the reference to the actual class type from the interface of the base class. --- diff --git a/include/deal.II/dofs/dof_handler_policy.h b/include/deal.II/dofs/dof_handler_policy.h index c73e78b967..633cc7c7c8 100644 --- a/include/deal.II/dofs/dof_handler_policy.h +++ b/include/deal.II/dofs/dof_handler_policy.h @@ -62,25 +62,24 @@ namespace internal virtual ~PolicyBase () = default; /** - * Distribute degrees of freedom on the object given as first - * argument. The reference to the NumberCache of the DoFHandler object - * has to be passed in a second argument. It could then be modified to + * Distribute degrees of freedom on the DoFHandler object associated + * with this policy object. The argument is a reference to the NumberCache + * of the DoFHandler object. The function may modify it to * make DoFHandler related functions work properly when called within * the policies classes. The updated NumberCache is written to that * argument. */ virtual void - distribute_dofs (dealii::DoFHandler &dof_handler, - NumberCache &number_cache) const = 0; + distribute_dofs (NumberCache &number_cache) const = 0; /** - * Distribute the multigrid dofs on each level + * Distribute the multigrid dofs on each level of the DoFHandler + * associated with this policy object. */ virtual void - distribute_mg_dofs (dealii::DoFHandler &dof_handler, - std::vector &number_caches) const = 0; + distribute_mg_dofs (std::vector &number_caches) const = 0; /** * Renumber degrees of freedom as specified by the first argument. The @@ -93,7 +92,6 @@ namespace internal virtual void renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache) const = 0; }; @@ -107,29 +105,33 @@ namespace internal { public: /** - * Distribute degrees of freedom on the object given as last argument. + * Constructor. + * @param dof_handler The DoFHandler object upon which this + * policy class is supposed to work. */ + Sequential (dealii::DoFHandler &dof_handler); + + // documentation is inherited virtual void - distribute_dofs (dealii::DoFHandler &dof_handler, - NumberCache &number_cache) const; + distribute_dofs (NumberCache &number_cache) const; - /** - * Distribute multigrid DoFs. - */ + // documentation is inherited virtual void - distribute_mg_dofs (dealii::DoFHandler &dof_handler, - std::vector &number_caches) const; + distribute_mg_dofs (std::vector &number_caches) const; - /** - * Renumber degrees of freedom as specified by the first argument. - */ + // documentation is inherited virtual void renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache) const; + + protected: + /** + * The DoFHandler object on which this policy object works. + */ + dealii::DoFHandler &dof_handler; }; /** @@ -140,6 +142,12 @@ namespace internal class ParallelShared : public Sequential { public: + /** + * Constructor. + * @param dof_handler The DoFHandler object upon which this + * policy class is supposed to work. + */ + ParallelShared (dealii::DoFHandler &dof_handler); /** * Distribute degrees of freedom on the object given as first @@ -151,16 +159,14 @@ namespace internal */ virtual void - distribute_dofs (dealii::DoFHandler &dof_handler, - NumberCache &number_cache) const; + distribute_dofs (NumberCache &number_cache) const; /** * This function is not yet implemented. */ virtual void - distribute_mg_dofs (dealii::DoFHandler &dof_handler, - std::vector &number_caches) const; + distribute_mg_dofs (std::vector &number_caches) const; /** * Renumber degrees of freedom as specified by the first argument. @@ -174,7 +180,6 @@ namespace internal virtual void renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache) const; }; @@ -188,29 +193,33 @@ namespace internal { public: /** - * Distribute degrees of freedom on the object given as last argument. + * Constructor. + * @param dof_handler The DoFHandler object upon which this + * policy class is supposed to work. */ + ParallelDistributed (dealii::DoFHandler &dof_handler); + + // documentation is inherited virtual void - distribute_dofs (dealii::DoFHandler &dof_handler, - NumberCache &number_cache) const; + distribute_dofs (NumberCache &number_cache) const; - /** - * Distribute multigrid DoFs. - */ + // documentation is inherited virtual void - distribute_mg_dofs (dealii::DoFHandler &dof_handler, - std::vector &number_caches) const; + distribute_mg_dofs (std::vector &number_caches) const; - /** - * Renumber degrees of freedom as specified by the first argument. - */ + // documentation is inherited virtual void renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache) const; + + private: + /** + * The DoFHandler object on which this policy object works. + */ + dealii::DoFHandler &dof_handler; }; } } diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 5013463546..d90bd71e20 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -659,13 +659,13 @@ DoFHandler::DoFHandler (const Triangulation &tria) if (dynamic_cast*> (&tria) != nullptr) - policy.reset (new internal::DoFHandler::Policy::ParallelShared()); + policy.reset (new internal::DoFHandler::Policy::ParallelShared(*this)); else if (dynamic_cast*> (&tria) == nullptr) - policy.reset (new internal::DoFHandler::Policy::Sequential()); + policy.reset (new internal::DoFHandler::Policy::Sequential(*this)); else - policy.reset (new internal::DoFHandler::Policy::ParallelDistributed()); + policy.reset (new internal::DoFHandler::Policy::ParallelDistributed(*this)); } @@ -682,14 +682,19 @@ DoFHandler::~DoFHandler () { // release allocated memory clear (); + + // also release the policy. this needs to happen before the + // current object disappears because the policy objects + // store references to the DoFhandler object they work on + policy.reset (); } template void -DoFHandler::initialize( - const Triangulation &t, - const FiniteElement &fe) +DoFHandler:: +initialize(const Triangulation &t, + const FiniteElement &fe) { tria = &t; faces = nullptr; @@ -701,13 +706,13 @@ DoFHandler::initialize( if (dynamic_cast*> (&t) != nullptr) - policy.reset (new internal::DoFHandler::Policy::ParallelShared()); + policy.reset (new internal::DoFHandler::Policy::ParallelShared(*this)); else if (dynamic_cast*> (&t) == nullptr) - policy.reset (new internal::DoFHandler::Policy::Sequential()); + policy.reset (new internal::DoFHandler::Policy::Sequential(*this)); else - policy.reset (new internal::DoFHandler::Policy::ParallelDistributed()); + policy.reset (new internal::DoFHandler::Policy::ParallelDistributed(*this)); distribute_dofs(fe); } @@ -1093,7 +1098,7 @@ void DoFHandler::distribute_dofs (const FiniteElementdistribute_dofs (*this,number_cache); + policy->distribute_dofs (number_cache); // initialize the block info object // only if this is a sequential @@ -1127,7 +1132,7 @@ void DoFHandler::distribute_mg_dofs (const FiniteElementn_global_levels()); - policy->distribute_mg_dofs (*this, mg_number_cache); + policy->distribute_mg_dofs (mg_number_cache); // initialize the block info object // only if this is a sequential @@ -1216,7 +1221,7 @@ DoFHandler::renumber_dofs (const std::vectorrenumber_dofs (new_numbers, *this,number_cache); + policy->renumber_dofs (new_numbers, number_cache); } diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index bf859fdc18..c18993589a 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -911,11 +911,20 @@ namespace internal /* --------------------- class Sequential ---------------- */ + + template + Sequential:: + Sequential (dealii::DoFHandler &dof_handler) + : + dof_handler (dof_handler) + {} + + + template void Sequential:: - distribute_dofs (DoFHandler &dof_handler, - NumberCache &number_cache_current ) const + distribute_dofs (NumberCache &number_cache_current) const { const types::global_dof_index n_dofs = Implementation::distribute_dofs (0, @@ -949,8 +958,7 @@ namespace internal template void Sequential:: - distribute_mg_dofs (DoFHandler &dof_handler, - std::vector &number_caches) const + distribute_mg_dofs (std::vector &number_caches) const { std::vector user_flags; @@ -978,7 +986,6 @@ namespace internal void Sequential:: renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache_current) const { Implementation::renumber_dofs (new_numbers, IndexSet(0), @@ -1009,22 +1016,33 @@ namespace internal /* --------------------- class ParallelShared ---------------- */ + + template + ParallelShared:: + ParallelShared (dealii::DoFHandler &dof_handler) + : + Sequential (dof_handler) + {} + + + + + template void ParallelShared:: - distribute_dofs (DoFHandler &dof_handler, - NumberCache &number_cache) const + distribute_dofs (NumberCache &number_cache) const { // If the underlying shared::Tria allows artificial cells, we need to do // some tricks here to make Sequential algorithms play nicely. // Namely, we first restore the original partition (without artificial cells) // and then turn artificial cells on at the end of this function. const parallel::shared::Triangulation *tr = - (dynamic_cast*> (&dof_handler.get_triangulation())); + (dynamic_cast*> (&this->dof_handler.get_triangulation())); Assert(tr != nullptr, ExcInternalError()); typename parallel::shared::Triangulation::active_cell_iterator - cell = dof_handler.get_triangulation().begin_active(), - endc = dof_handler.get_triangulation().end(); + cell = this->dof_handler.get_triangulation().begin_active(), + endc = this->dof_handler.get_triangulation().end(); std::vector current_subdomain_ids(tr->n_active_cells()); const std::vector &true_subdomain_ids = tr->get_true_subdomain_ids_of_cells(); if (tr->with_artificial_cells()) @@ -1036,8 +1054,8 @@ namespace internal // let the sequential algorithm do its magic, then sort DoF indices // by subdomain - Sequential::distribute_dofs (dof_handler,number_cache); - DoFRenumbering::subdomain_wise (dof_handler); + this->Sequential::distribute_dofs (number_cache); + DoFRenumbering::subdomain_wise (this->dof_handler); // dofrenumbering will reset subdomains, this is ugly but we need to do it again: cell = tr->begin_active(); @@ -1045,12 +1063,12 @@ namespace internal for (unsigned int index=0; cell != endc; cell++, index++) cell->set_subdomain_id(true_subdomain_ids[index]); - number_cache.locally_owned_dofs_per_processor = DoFTools::locally_owned_dofs_per_subdomain (dof_handler); - number_cache.locally_owned_dofs = number_cache.locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain()]; + number_cache.locally_owned_dofs_per_processor = DoFTools::locally_owned_dofs_per_subdomain (this->dof_handler); + number_cache.locally_owned_dofs = number_cache.locally_owned_dofs_per_processor[this->dof_handler.get_triangulation().locally_owned_subdomain()]; number_cache.n_locally_owned_dofs_per_processor.resize (number_cache.locally_owned_dofs_per_processor.size()); for (unsigned int i = 0; i < number_cache.n_locally_owned_dofs_per_processor.size(); i++) number_cache.n_locally_owned_dofs_per_processor[i] = number_cache.locally_owned_dofs_per_processor[i].n_elements(); - number_cache.n_locally_owned_dofs = number_cache.n_locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain()]; + number_cache.n_locally_owned_dofs = number_cache.n_locally_owned_dofs_per_processor[this->dof_handler.get_triangulation().locally_owned_subdomain()]; // restore current subdomain ids cell = tr->begin_active(); @@ -1064,11 +1082,10 @@ namespace internal template void ParallelShared:: - distribute_mg_dofs (DoFHandler &dof_handler, - std::vector &number_caches) const + distribute_mg_dofs (std::vector &number_caches) const { // first, call the sequential function to distribute dofs - Sequential::distribute_mg_dofs (dof_handler, number_caches); + this->Sequential::distribute_mg_dofs (number_caches); // now we need to update the number cache. This part is not yet implemented. AssertThrow(false,ExcNotImplemented()); } @@ -1079,25 +1096,23 @@ namespace internal void ParallelShared:: renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache) const { #ifndef DEAL_II_WITH_MPI (void)new_numbers; - (void)dof_handler; (void)number_cache; Assert (false, ExcNotImplemented()); #else // Similar to distribute_dofs() we need to have a special treatment in // case artificial cells are present. const parallel::shared::Triangulation *tr = - (dynamic_cast*> (&dof_handler.get_triangulation())); + (dynamic_cast*> (&this->dof_handler.get_triangulation())); Assert(tr != nullptr, ExcInternalError()); typename parallel::shared::Triangulation::active_cell_iterator - cell = dof_handler.get_triangulation().begin_active(), - endc = dof_handler.get_triangulation().end(); + cell = this->dof_handler.get_triangulation().begin_active(), + endc = this->dof_handler.get_triangulation().end(); std::vector current_subdomain_ids(tr->n_active_cells()); const std::vector &true_subdomain_ids = tr->get_true_subdomain_ids_of_cells(); if (tr->with_artificial_cells()) @@ -1107,22 +1122,22 @@ namespace internal cell->set_subdomain_id(true_subdomain_ids[index]); } - std::vector global_gathered_numbers (dof_handler.n_dofs (), 0); + std::vector global_gathered_numbers (this->dof_handler.n_dofs (), 0); // as we call DoFRenumbering::subdomain_wise (dof_handler) from distribute_dofs(), // we need to support sequential-like input. // Distributed-like input from, for example, component_wise renumbering is also supported. - if (new_numbers.size () == dof_handler.n_dofs ()) + if (new_numbers.size () == this->dof_handler.n_dofs ()) { global_gathered_numbers = new_numbers; } else { - Assert(new_numbers.size() == dof_handler.locally_owned_dofs().n_elements(), + Assert(new_numbers.size() == this->dof_handler.locally_owned_dofs().n_elements(), ExcInternalError()); const unsigned int n_cpu = Utilities::MPI::n_mpi_processes (tr->get_communicator ()); - std::vector gathered_new_numbers (dof_handler.n_dofs (), 0); + std::vector gathered_new_numbers (this->dof_handler.n_dofs (), 0); Assert(Utilities::MPI::this_mpi_process (tr->get_communicator ()) == - dof_handler.get_triangulation().locally_owned_subdomain (), + this->dof_handler.get_triangulation().locally_owned_subdomain (), ExcInternalError()) // gather new numbers among processors into one vector @@ -1165,8 +1180,8 @@ namespace internal // flag_1 and flag_2 are // used to control that there is a // one-to-one relation between old and new DoFs. - std::vector flag_1 (dof_handler.n_dofs (), 0); - std::vector flag_2 (dof_handler.n_dofs (), 0); + std::vector flag_1 (this->dof_handler.n_dofs (), 0); + std::vector flag_2 (this->dof_handler.n_dofs (), 0); for (unsigned int i = 0; i < n_cpu; i++) { const IndexSet &iset = @@ -1176,8 +1191,8 @@ namespace internal { const types::global_dof_index target = iset.nth_index_in_set (ind); const types::global_dof_index value = gathered_new_numbers[shift + ind]; - Assert(target < dof_handler.n_dofs(), ExcInternalError()); - Assert(value < dof_handler.n_dofs(), ExcInternalError()); + Assert(target < this->dof_handler.n_dofs(), ExcInternalError()); + Assert(value < this->dof_handler.n_dofs(), ExcInternalError()); global_gathered_numbers[target] = value; flag_1[target]++; flag_2[value]++; @@ -1195,12 +1210,12 @@ namespace internal ExcInternalError()); } - Sequential::renumber_dofs (global_gathered_numbers, dof_handler, number_cache); + this->Sequential::renumber_dofs (global_gathered_numbers, number_cache); // correct number_cache: number_cache.locally_owned_dofs_per_processor = - DoFTools::locally_owned_dofs_per_subdomain (dof_handler); + DoFTools::locally_owned_dofs_per_subdomain (this->dof_handler); number_cache.locally_owned_dofs = - number_cache.locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain ()]; + number_cache.locally_owned_dofs_per_processor[this->dof_handler.get_triangulation().locally_owned_subdomain ()]; // sequential renumbering returns a vector of size 1 here, // correct this: number_cache.n_locally_owned_dofs_per_processor.resize(number_cache.locally_owned_dofs_per_processor.size()); @@ -1209,7 +1224,7 @@ namespace internal number_cache.n_locally_owned_dofs_per_processor[i] = number_cache.locally_owned_dofs_per_processor[i].n_elements (); number_cache.n_locally_owned_dofs = - number_cache.n_locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain ()]; + number_cache.n_locally_owned_dofs_per_processor[this->dof_handler.get_triangulation().locally_owned_subdomain ()]; // restore artificial cells cell = tr->begin_active(); @@ -2059,16 +2074,26 @@ namespace internal + template + ParallelDistributed:: + ParallelDistributed (dealii::DoFHandler &dof_handler) + : + dof_handler (dof_handler) + {} + + + + + + template void ParallelDistributed:: - distribute_dofs (DoFHandler &dof_handler, - NumberCache &number_cache_current) const + distribute_dofs (NumberCache &number_cache_current) const { NumberCache number_cache; #ifndef DEAL_II_WITH_P4EST - (void)dof_handler; Assert (false, ExcNotImplemented()); #else @@ -2271,11 +2296,9 @@ namespace internal template void ParallelDistributed:: - distribute_mg_dofs (DoFHandler &dof_handler, - std::vector &number_caches) const + distribute_mg_dofs (std::vector &number_caches) const { #ifndef DEAL_II_WITH_P4EST - (void)dof_handler; (void)number_caches; Assert (false, ExcNotImplemented()); #else @@ -2524,11 +2547,9 @@ namespace internal void ParallelDistributed:: renumber_dofs (const std::vector &new_numbers, - dealii::DoFHandler &dof_handler, NumberCache &number_cache_current) const { (void)new_numbers; - (void)dof_handler; Assert (new_numbers.size() == dof_handler.locally_owned_dofs().n_elements(), ExcInternalError());