From: Daniel Arndt Date: Wed, 18 Apr 2018 19:56:35 +0000 (+0200) Subject: Store a copy of the FECollection object in hp::DoFHandler X-Git-Tag: v9.0.0-rc1~137^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=641bb0fcdd31b5bfe2eca37c1f4bc8fba367a043;p=dealii.git Store a copy of the FECollection object in hp::DoFHandler --- diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index b460a1f793..376a04bd18 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -1164,13 +1164,13 @@ namespace internal Assert ( (fe_index != dealii::hp::DoFHandler::default_fe_index), ExcMessage ("You need to specify a FE index when working " "with hp DoFHandlers")); - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection != nullptr, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); Assert (local_index < dof_handler.get_fe(fe_index).dofs_per_vertex, ExcIndexRange(local_index, 0, dof_handler.get_fe(fe_index).dofs_per_vertex)); - Assert (fe_index < dof_handler.finite_elements->size(), + Assert (fe_index < dof_handler.fe_collection->size(), ExcInternalError()); Assert (dof_handler.vertex_dof_offsets[vertex_index] != numbers::invalid_unsigned_int, @@ -1195,7 +1195,7 @@ namespace internal Assert (this_fe_index != numbers::invalid_dof_index, ExcInternalError()); - Assert (this_fe_index < dof_handler.finite_elements->size(), + Assert (this_fe_index < dof_handler.fe_collection->size(), ExcInternalError()); if (this_fe_index == fe_index) @@ -1249,7 +1249,7 @@ namespace internal Assert ( (fe_index != dealii::hp::DoFHandler::default_fe_index), ExcMessage ("You need to specify a FE index when working " "with hp DoFHandlers")); - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); Assert (local_index < dof_handler.get_fe(fe_index).dofs_per_vertex, @@ -1280,7 +1280,7 @@ namespace internal Assert (this_fe_index != numbers::invalid_dof_index, ExcInternalError()); - Assert (this_fe_index < dof_handler.finite_elements->size(), + Assert (this_fe_index < dof_handler.fe_collection->size(), ExcInternalError()); if (this_fe_index == fe_index) @@ -1302,7 +1302,7 @@ namespace internal n_active_vertex_fe_indices (const dealii::hp::DoFHandler &dof_handler, const unsigned int vertex_index) { - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection != nullptr, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); @@ -1350,7 +1350,7 @@ namespace internal const unsigned int vertex_index, const unsigned int n) { - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection != nullptr, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); Assert (n < n_active_vertex_fe_indices(dof_handler, vertex_index), @@ -1379,7 +1379,7 @@ namespace internal Assert((*pointer)::max(), ExcInternalError()); const types::global_dof_index this_fe_index = *pointer; - Assert (this_fe_index < dof_handler.finite_elements->size(), + Assert (this_fe_index < dof_handler.fe_collection->size(), ExcInternalError()); if (counter == n) @@ -1410,10 +1410,10 @@ namespace internal Assert ( (fe_index != dealii::hp::DoFHandler::default_fe_index), ExcMessage ("You need to specify a FE index when working " "with hp DoFHandlers")); - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection != nullptr, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - Assert (fe_index < dof_handler.finite_elements->size(), + Assert (fe_index < dof_handler.fe_collection->size(), ExcInternalError()); // make sure we don't ask on @@ -1438,7 +1438,7 @@ namespace internal Assert((*pointer)::max(), ExcInternalError()); const types::global_dof_index this_fe_index = *pointer; - Assert (this_fe_index < dof_handler.finite_elements->size(), + Assert (this_fe_index < dof_handler.fe_collection->size(), ExcInternalError()); if (this_fe_index == numbers::invalid_dof_index) diff --git a/include/deal.II/hp/dof_faces.h b/include/deal.II/hp/dof_faces.h index 70d5c263b3..445f49dc9d 100644 --- a/include/deal.II/hp/dof_faces.h +++ b/include/deal.II/hp/dof_faces.h @@ -24,13 +24,6 @@ DEAL_II_NAMESPACE_OPEN -namespace hp -{ - template - class FECollection; -} - - namespace internal { namespace hp diff --git a/include/deal.II/hp/dof_handler.h b/include/deal.II/hp/dof_handler.h index fd9547d766..263ce64be3 100644 --- a/include/deal.II/hp/dof_handler.h +++ b/include/deal.II/hp/dof_handler.h @@ -313,15 +313,8 @@ namespace hp * want a particular ordering, use the functions in namespace * DoFRenumbering. * - * @note In contrast to the dealii::DoFHandler::distribute_dofs() - * function, this function does not make a copy of the object - * given as argument. Rather, it stores a reference to the given - * object, and it is the responsibility of user code to ensure - * that the hp::FECollection given as argument lives at least as - * long as the hp::DoFHandler object. If you want to break this - * dependence by asking the hp::DoFHandler to release the - * reference to the hp::FECollection object, call the - * hp::DoFHandler::clear() function. + * @note In accordance with dealii::DoFHandler::distribute_dofs(), + * this function also makes a copy of the object given as argument. */ virtual void distribute_dofs (const hp::FECollection &fe); @@ -832,14 +825,9 @@ namespace hp SmartPointer,DoFHandler > tria; /** - * Store a pointer to the finite element set given latest for the - * distribution of dofs. In order to avoid destruction of the object - * before the lifetime of the DoF handler, we subscribe to the finite - * element object. To unlock the FE before the end of the lifetime of this - * DoF handler, use the clear() function (this clears all data of - * this object as well, though). + * Store a copy of the finite element set given latest to distribute_dofs(). */ - SmartPointer,hp::DoFHandler > finite_elements; + std::unique_ptr > fe_collection; /** * An object that describes how degrees of freedom should be distributed and @@ -1213,10 +1201,10 @@ namespace hp const hp::FECollection & DoFHandler::get_fe () const { - Assert (finite_elements != nullptr, + Assert (fe_collection, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return *finite_elements; + return *fe_collection; } @@ -1227,10 +1215,10 @@ namespace hp DoFHandler::get_fe (const unsigned int number) const { - Assert (finite_elements != nullptr, + Assert (fe_collection, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return (*finite_elements)[number]; + return (*fe_collection)[number]; } @@ -1240,10 +1228,10 @@ namespace hp const hp::FECollection & DoFHandler::get_fe_collection () const { - Assert (finite_elements != nullptr, + Assert (fe_collection, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return *finite_elements; + return *fe_collection; } diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index e43eac2ea9..a6812f0de7 100644 --- a/source/dofs/dof_handler_policy.cc +++ b/source/dofs/dof_handler_policy.cc @@ -772,8 +772,8 @@ namespace internal // and then only deal with those that are not identical of which we can // handle at most 2 dealii::Table<2,std::unique_ptr > - line_dof_identities (dof_handler.finite_elements->size(), - dof_handler.finite_elements->size()); + line_dof_identities (dof_handler.fe_collection->size(), + dof_handler.fe_collection->size()); for (typename hp::DoFHandler::active_cell_iterator cell=dof_handler.begin_active(); @@ -1004,8 +1004,8 @@ namespace internal // higher, and for quads only in 4d and higher, so this // isn't a particularly frequent case dealii::Table<2,std::unique_ptr > - quad_dof_identities (dof_handler.finite_elements->size(), - dof_handler.finite_elements->size()); + quad_dof_identities (dof_handler.fe_collection->size(), + dof_handler.fe_collection->size()); for (typename hp::DoFHandler::active_cell_iterator cell=dof_handler.begin_active(); diff --git a/source/hp/dof_handler.cc b/source/hp/dof_handler.cc index bb8f2ec2f2..07c15b4c7a 100644 --- a/source/hp/dof_handler.cc +++ b/source/hp/dof_handler.cc @@ -141,7 +141,7 @@ namespace internal locally_used_vertices[cell->vertex_index(v)] = true; std::vector > - vertex_fe_association (dof_handler.finite_elements->size(), + vertex_fe_association (dof_handler.fe_collection->size(), std::vector (dof_handler.tria->n_vertices(), false)); for (typename HpDoFHandler::active_cell_iterator @@ -162,10 +162,10 @@ namespace internal if (dof_handler.tria->vertex_used(v) == true) { unsigned int fe=0; - for (; fesize(); ++fe) + for (; fesize(); ++fe) if (vertex_fe_association[fe][v] == true) break; - Assert (fe != dof_handler.finite_elements->size(), ExcInternalError()); + Assert (fe != dof_handler.fe_collection->size(), ExcInternalError()); } #endif @@ -184,7 +184,7 @@ namespace internal { dof_handler.vertex_dof_offsets[v] = vertex_slots_needed; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fesize(); ++fe) if (vertex_fe_association[fe][v] == true) vertex_slots_needed += dof_handler.get_fe(fe).dofs_per_vertex + 1; @@ -201,7 +201,7 @@ namespace internal if (locally_used_vertices[v] == true) { unsigned int current_index = dof_handler.vertex_dof_offsets[v]; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fesize(); ++fe) if (vertex_fe_association[fe][v] == true) { // if this vertex uses this fe, then set the @@ -598,9 +598,9 @@ namespace internal void reserve_space (DoFHandler<1,spacedim> &dof_handler) { - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection, (typename DoFHandler<1,spacedim>::ExcNoFESelected())); - Assert (dof_handler.finite_elements->size() > 0, + Assert (dof_handler.fe_collection->size() > 0, (typename DoFHandler<1,spacedim>::ExcNoFESelected())); Assert (dof_handler.tria->n_levels() > 0, ExcMessage("The current Triangulation must not be empty.")); @@ -620,9 +620,9 @@ namespace internal void reserve_space (DoFHandler<2,spacedim> &dof_handler) { - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection, (typename DoFHandler<2,spacedim>::ExcNoFESelected())); - Assert (dof_handler.finite_elements->size() > 0, + Assert (dof_handler.fe_collection->size() > 0, (typename DoFHandler<2,spacedim>::ExcNoFESelected())); Assert (dof_handler.tria->n_levels() > 0, ExcMessage("The current Triangulation must not be empty.")); @@ -648,9 +648,9 @@ namespace internal { const unsigned int dim = 3; - Assert (dof_handler.finite_elements != nullptr, + Assert (dof_handler.fe_collection, (typename DoFHandler::ExcNoFESelected())); - Assert (dof_handler.finite_elements->size() > 0, + Assert (dof_handler.fe_collection->size() > 0, (typename DoFHandler::ExcNoFESelected())); Assert (dof_handler.tria->n_levels() > 0, ExcMessage("The current Triangulation must not be empty.")); @@ -679,7 +679,7 @@ namespace internal // given fe's, by setting a bit. in a later step, we // then actually allocate memory for the required dofs std::vector > - line_fe_association (dof_handler.finite_elements->size(), + line_fe_association (dof_handler.fe_collection->size(), std::vector (dof_handler.tria->n_raw_lines(), false)); @@ -696,7 +696,7 @@ namespace internal // case we do not have to allocate any memory at all std::vector line_is_used (dof_handler.tria->n_raw_lines(), false); for (unsigned int line=0; linen_raw_lines(); ++line) - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fesize(); ++fe) if (line_fe_association[fe][line] == true) { line_is_used[line] = true; @@ -718,7 +718,7 @@ namespace internal { dof_handler.faces->lines.dof_offsets[line] = line_slots_needed; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fesize(); ++fe) if (line_fe_association[fe][line] == true) line_slots_needed += dof_handler.get_fe(fe).dofs_per_line + 1; ++line_slots_needed; @@ -732,7 +732,7 @@ namespace internal if (line_is_used[line] == true) { unsigned int pointer = dof_handler.faces->lines.dof_offsets[line]; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fesize(); ++fe) if (line_fe_association[fe][line] == true) { // if this line uses this fe, then set the @@ -760,8 +760,8 @@ namespace internal max_couplings_between_dofs (const DoFHandler<1,spacedim> &dof_handler) { return std::min(static_cast (3* - dof_handler.finite_elements->max_dofs_per_vertex() + - 2*dof_handler.finite_elements->max_dofs_per_line()), + dof_handler.fe_collection->max_dofs_per_vertex() + + 2*dof_handler.fe_collection->max_dofs_per_line()), dof_handler.n_dofs()); } @@ -795,29 +795,29 @@ namespace internal switch (dof_handler.tria->max_adjacent_cells()) { case 4: - max_couplings=19*dof_handler.finite_elements->max_dofs_per_vertex() + - 28*dof_handler.finite_elements->max_dofs_per_line() + - 8*dof_handler.finite_elements->max_dofs_per_quad(); + max_couplings=19*dof_handler.fe_collection->max_dofs_per_vertex() + + 28*dof_handler.fe_collection->max_dofs_per_line() + + 8*dof_handler.fe_collection->max_dofs_per_quad(); break; case 5: - max_couplings=21*dof_handler.finite_elements->max_dofs_per_vertex() + - 31*dof_handler.finite_elements->max_dofs_per_line() + - 9*dof_handler.finite_elements->max_dofs_per_quad(); + max_couplings=21*dof_handler.fe_collection->max_dofs_per_vertex() + + 31*dof_handler.fe_collection->max_dofs_per_line() + + 9*dof_handler.fe_collection->max_dofs_per_quad(); break; case 6: - max_couplings=28*dof_handler.finite_elements->max_dofs_per_vertex() + - 42*dof_handler.finite_elements->max_dofs_per_line() + - 12*dof_handler.finite_elements->max_dofs_per_quad(); + max_couplings=28*dof_handler.fe_collection->max_dofs_per_vertex() + + 42*dof_handler.fe_collection->max_dofs_per_line() + + 12*dof_handler.fe_collection->max_dofs_per_quad(); break; case 7: - max_couplings=30*dof_handler.finite_elements->max_dofs_per_vertex() + - 45*dof_handler.finite_elements->max_dofs_per_line() + - 13*dof_handler.finite_elements->max_dofs_per_quad(); + max_couplings=30*dof_handler.fe_collection->max_dofs_per_vertex() + + 45*dof_handler.fe_collection->max_dofs_per_line() + + 13*dof_handler.fe_collection->max_dofs_per_quad(); break; case 8: - max_couplings=37*dof_handler.finite_elements->max_dofs_per_vertex() + - 56*dof_handler.finite_elements->max_dofs_per_line() + - 16*dof_handler.finite_elements->max_dofs_per_quad(); + max_couplings=37*dof_handler.fe_collection->max_dofs_per_vertex() + + 56*dof_handler.fe_collection->max_dofs_per_line() + + 16*dof_handler.fe_collection->max_dofs_per_quad(); break; default: Assert (false, ExcNotImplemented()); @@ -844,10 +844,10 @@ namespace internal types::global_dof_index max_couplings; if (max_adjacent_cells <= 8) - max_couplings=7*7*7*dof_handler.finite_elements->max_dofs_per_vertex() + - 7*6*7*3*dof_handler.finite_elements->max_dofs_per_line() + - 9*4*7*3*dof_handler.finite_elements->max_dofs_per_quad() + - 27*dof_handler.finite_elements->max_dofs_per_hex(); + max_couplings=7*7*7*dof_handler.fe_collection->max_dofs_per_vertex() + + 7*6*7*3*dof_handler.fe_collection->max_dofs_per_line() + + 9*4*7*3*dof_handler.fe_collection->max_dofs_per_quad() + + 27*dof_handler.fe_collection->max_dofs_per_hex(); else { Assert (false, ExcNotImplemented()); @@ -1106,7 +1106,7 @@ namespace hp template types::global_dof_index DoFHandler::n_boundary_dofs () const { - Assert (finite_elements != nullptr, ExcNoFESelected()); + Assert (fe_collection, ExcNoFESelected()); std::set boundary_dofs; std::vector dofs_on_face; @@ -1140,7 +1140,7 @@ namespace hp types::global_dof_index DoFHandler::n_boundary_dofs (const std::set &boundary_ids) const { - Assert (finite_elements != nullptr, ExcNoFESelected()); + Assert (fe_collection, ExcNoFESelected()); Assert (boundary_ids.find (numbers::internal_face_boundary_id) == boundary_ids.end(), ExcInvalidBoundaryIndicator()); @@ -1204,7 +1204,7 @@ namespace hp DoFHandler::memory_consumption () const { std::size_t mem = (MemoryConsumption::memory_consumption (tria) + - MemoryConsumption::memory_consumption (finite_elements) + + MemoryConsumption::memory_consumption (fe_collection) + MemoryConsumption::memory_consumption (tria) + MemoryConsumption::memory_consumption (levels) + MemoryConsumption::memory_consumption (*faces) + @@ -1288,7 +1288,8 @@ namespace hp Assert (tria->n_levels() > 0, ExcMessage("The Triangulation you are using is empty!")); - finite_elements = &ff; + if (fe_collection == nullptr || &(*fe_collection) != &ff) + fe_collection = std_cxx14::make_unique>(ff); // at the beginning, make sure every processor knows the // active_fe_indices on both its own cells and all ghost cells @@ -1327,9 +1328,9 @@ namespace hp // cover all fe indices presently in use on the mesh for (active_cell_iterator cell = begin_active(); cell != end(); ++cell) if (cell->is_locally_owned()) - Assert (cell->active_fe_index() < finite_elements->size(), + Assert (cell->active_fe_index() < fe_collection->size(), ExcInvalidFEIndex (cell->active_fe_index(), - finite_elements->size())); + fe_collection->size())); // then allocate space for all the other tables @@ -1371,7 +1372,7 @@ namespace hp Threads::TaskGroup<> tg; for (int level=levels.size()-1; level>=0; --level) tg += Threads::new_task (&dealii::internal::hp::DoFLevel::compress_data, - *levels[level], *finite_elements); + *levels[level], *fe_collection); tg.join_all (); } @@ -1410,7 +1411,7 @@ namespace hp void DoFHandler::clear () { // release lock to old fe - finite_elements = nullptr; + fe_collection = nullptr; // release memory clear_space (); @@ -1457,7 +1458,7 @@ namespace hp Threads::TaskGroup<> tg; for (int level=levels.size()-1; level>=0; --level) tg += Threads::new_task (&dealii::internal::hp::DoFLevel::uncompress_data, - *levels[level], *finite_elements); + *levels[level], *fe_collection); tg.join_all (); } @@ -1469,7 +1470,7 @@ namespace hp Threads::TaskGroup<> tg; for (int level=levels.size()-1; level>=0; --level) tg += Threads::new_task (&dealii::internal::hp::DoFLevel::compress_data, - *levels[level], *finite_elements); + *levels[level], *fe_collection); tg.join_all (); } } @@ -1480,7 +1481,7 @@ namespace hp unsigned int DoFHandler::max_couplings_between_dofs () const { - Assert (finite_elements != nullptr, ExcNoFESelected()); + Assert (fe_collection, ExcNoFESelected()); return dealii::internal::hp::DoFHandlerImplementation::Implementation::max_couplings_between_dofs (*this); } @@ -1490,16 +1491,16 @@ namespace hp unsigned int DoFHandler::max_couplings_between_boundary_dofs () const { - Assert (finite_elements != nullptr, ExcNoFESelected()); + Assert (fe_collection, ExcNoFESelected()); switch (dim) { case 1: - return finite_elements->max_dofs_per_vertex(); + return fe_collection->max_dofs_per_vertex(); case 2: - return (3*finite_elements->max_dofs_per_vertex() + return (3*fe_collection->max_dofs_per_vertex() + - 2*finite_elements->max_dofs_per_line()); + 2*fe_collection->max_dofs_per_line()); case 3: // we need to take refinement of one boundary face into // consideration here; in fact, this function returns what @@ -1510,9 +1511,9 @@ namespace hp // time. fortunately, omitting it for now does no harm since // the matrix will cry foul if its requirements are not // satisfied - return (19*finite_elements->max_dofs_per_vertex() + - 28*finite_elements->max_dofs_per_line() + - 8*finite_elements->max_dofs_per_quad()); + return (19*fe_collection->max_dofs_per_vertex() + + 28*fe_collection->max_dofs_per_line() + + 8*fe_collection->max_dofs_per_quad()); default: Assert (false, ExcNotImplemented()); return 0; @@ -1624,7 +1625,7 @@ namespace hp // indicators are zero anyway, and there is no point trying to set // anything (besides, we would trip over an assertion in // set_active_fe_index) - if (finite_elements != nullptr) + if (fe_collection) { cell_iterator cell = begin(), endc = end ();