From 5f78adbe7e3761095cc3863205b9e0b647a56b1f Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Tue, 24 Apr 2018 01:32:38 +0200 Subject: [PATCH] Store the FECollection object and not just the pointer in DoFHandler and hp::DoFHandler --- include/deal.II/dofs/dof_accessor.templates.h | 22 ++-- include/deal.II/dofs/dof_handler.h | 8 +- include/deal.II/hp/dof_handler.h | 14 +-- source/dofs/dof_handler.cc | 11 +- source/dofs/dof_handler_policy.cc | 8 +- source/hp/dof_handler.cc | 111 ++++++++---------- 6 files changed, 80 insertions(+), 94 deletions(-) diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index 376a04bd18..65485e9aca 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.fe_collection != nullptr, + Assert (dof_handler.fe_collection.size() > 0, 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.fe_collection->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.fe_collection->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.fe_collection, + Assert (dof_handler.fe_collection.size() > 0, 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.fe_collection->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.fe_collection != nullptr, + Assert (dof_handler.fe_collection.size() > 0, 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.fe_collection != nullptr, + Assert (dof_handler.fe_collection.size() > 0, 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.fe_collection->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.fe_collection != nullptr, + Assert (dof_handler.fe_collection.size() > 0, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - Assert (fe_index < dof_handler.fe_collection->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.fe_collection->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/dofs/dof_handler.h b/include/deal.II/dofs/dof_handler.h index c4f8982b38..f988909cb5 100644 --- a/include/deal.II/dofs/dof_handler.h +++ b/include/deal.II/dofs/dof_handler.h @@ -1034,10 +1034,10 @@ private: /** - * Store a pointer to a hp::FECollection object containing the (one) + * Store a hp::FECollection object containing the (one) * FiniteElement this object is initialized with. */ - std::unique_ptr > fe_collection; + hp::FECollection fe_collection; /** * An object that describes how degrees of freedom should be distributed and @@ -1328,9 +1328,9 @@ inline const hp::FECollection & DoFHandler::get_fe_collection () const { - Assert(fe_collection != nullptr, + Assert(fe_collection.size() > 0, ExcMessage("You are trying to access the DoFHandler's FECollection object before it has been initialized.")); - return *fe_collection; + return fe_collection; } diff --git a/include/deal.II/hp/dof_handler.h b/include/deal.II/hp/dof_handler.h index 263ce64be3..aeb6024d01 100644 --- a/include/deal.II/hp/dof_handler.h +++ b/include/deal.II/hp/dof_handler.h @@ -827,7 +827,7 @@ namespace hp /** * Store a copy of the finite element set given latest to distribute_dofs(). */ - std::unique_ptr > fe_collection; + hp::FECollection fe_collection; /** * An object that describes how degrees of freedom should be distributed and @@ -1201,10 +1201,10 @@ namespace hp const hp::FECollection & DoFHandler::get_fe () const { - Assert (fe_collection, + Assert (fe_collection.size() > 0, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return *fe_collection; + return fe_collection; } @@ -1215,10 +1215,10 @@ namespace hp DoFHandler::get_fe (const unsigned int number) const { - Assert (fe_collection, + Assert (fe_collection.size() > 0, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return (*fe_collection)[number]; + return fe_collection[number]; } @@ -1228,10 +1228,10 @@ namespace hp const hp::FECollection & DoFHandler::get_fe_collection () const { - Assert (fe_collection, + Assert (fe_collection.size() > 0, ExcMessage ("No finite element collection is associated with " "this DoFHandler")); - return *fe_collection; + return fe_collection; } diff --git a/source/dofs/dof_handler.cc b/source/dofs/dof_handler.cc index 8bdd155f74..247bf6009d 100644 --- a/source/dofs/dof_handler.cc +++ b/source/dofs/dof_handler.cc @@ -653,7 +653,6 @@ template DoFHandler::DoFHandler (const Triangulation &tria) : tria(&tria, typeid(*this).name()), - fe_collection(nullptr), faces(nullptr), mg_faces (nullptr) { @@ -674,8 +673,7 @@ DoFHandler::DoFHandler (const Triangulation &tria) template DoFHandler::DoFHandler () : - tria(nullptr, typeid(*this).name()), - fe_collection(nullptr) + tria(nullptr, typeid(*this).name()) {} @@ -1002,8 +1000,8 @@ void DoFHandler::distribute_dofs (const FiniteElement>(ff); + if (fe_collection.size() == 0 || &(fe_collection[0]) != &ff) + fe_collection = hp::FECollection(ff); // delete all levels and set them // up newly. note that we still @@ -1094,9 +1092,6 @@ void DoFHandler::initialize_local_block_info () template void DoFHandler::clear () { - // release lock to old fe - fe_collection.reset(); - // release memory clear_space (); clear_mg_space (); diff --git a/source/dofs/dof_handler_policy.cc b/source/dofs/dof_handler_policy.cc index a6812f0de7..f43cf65e09 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.fe_collection->size(), - dof_handler.fe_collection->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.fe_collection->size(), - dof_handler.fe_collection->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 07c15b4c7a..1e741dd7bd 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.fe_collection->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(), 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) + for (unsigned int fe=0; fe &dof_handler) { - Assert (dof_handler.fe_collection, - (typename DoFHandler<1,spacedim>::ExcNoFESelected())); - Assert (dof_handler.fe_collection->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 +618,7 @@ namespace internal void reserve_space (DoFHandler<2,spacedim> &dof_handler) { - Assert (dof_handler.fe_collection, - (typename DoFHandler<2,spacedim>::ExcNoFESelected())); - Assert (dof_handler.fe_collection->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 +644,7 @@ namespace internal { const unsigned int dim = 3; - Assert (dof_handler.fe_collection, - (typename DoFHandler::ExcNoFESelected())); - Assert (dof_handler.fe_collection->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 +673,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.fe_collection->size(), + line_fe_association (dof_handler.fe_collection.size(), std::vector (dof_handler.tria->n_raw_lines(), false)); @@ -696,7 +690,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; felines.dof_offsets[line] = line_slots_needed; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; felines.dof_offsets[line]; - for (unsigned int fe=0; fesize(); ++fe) + for (unsigned int fe=0; fe &dof_handler) { return std::min(static_cast (3* - dof_handler.fe_collection->max_dofs_per_vertex() + - 2*dof_handler.fe_collection->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 +789,29 @@ namespace internal switch (dof_handler.tria->max_adjacent_cells()) { case 4: - 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(); + 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.fe_collection->max_dofs_per_vertex() + - 31*dof_handler.fe_collection->max_dofs_per_line() + - 9*dof_handler.fe_collection->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.fe_collection->max_dofs_per_vertex() + - 42*dof_handler.fe_collection->max_dofs_per_line() + - 12*dof_handler.fe_collection->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.fe_collection->max_dofs_per_vertex() + - 45*dof_handler.fe_collection->max_dofs_per_line() + - 13*dof_handler.fe_collection->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.fe_collection->max_dofs_per_vertex() + - 56*dof_handler.fe_collection->max_dofs_per_line() + - 16*dof_handler.fe_collection->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 +838,10 @@ namespace internal types::global_dof_index max_couplings; if (max_adjacent_cells <= 8) - 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(); + 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 +1100,7 @@ namespace hp template types::global_dof_index DoFHandler::n_boundary_dofs () const { - Assert (fe_collection, ExcNoFESelected()); + Assert(fe_collection.size() > 0, ExcNoFESelected()); std::set boundary_dofs; std::vector dofs_on_face; @@ -1140,7 +1134,7 @@ namespace hp types::global_dof_index DoFHandler::n_boundary_dofs (const std::set &boundary_ids) const { - Assert (fe_collection, ExcNoFESelected()); + Assert(fe_collection.size() > 0, ExcNoFESelected()); Assert (boundary_ids.find (numbers::internal_face_boundary_id) == boundary_ids.end(), ExcInvalidBoundaryIndicator()); @@ -1288,8 +1282,8 @@ namespace hp Assert (tria->n_levels() > 0, ExcMessage("The Triangulation you are using is empty!")); - if (fe_collection == nullptr || &(*fe_collection) != &ff) - fe_collection = std_cxx14::make_unique>(ff); + if (&fe_collection != &ff) + fe_collection = hp::FECollection(ff); // at the beginning, make sure every processor knows the // active_fe_indices on both its own cells and all ghost cells @@ -1328,9 +1322,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() < fe_collection->size(), + Assert (cell->active_fe_index() < fe_collection.size(), ExcInvalidFEIndex (cell->active_fe_index(), - fe_collection->size())); + fe_collection.size())); // then allocate space for all the other tables @@ -1372,7 +1366,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], *fe_collection); + *levels[level], fe_collection); tg.join_all (); } @@ -1410,9 +1404,6 @@ namespace hp template void DoFHandler::clear () { - // release lock to old fe - fe_collection = nullptr; - // release memory clear_space (); } @@ -1458,7 +1449,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], *fe_collection); + *levels[level], fe_collection); tg.join_all (); } @@ -1470,7 +1461,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], *fe_collection); + *levels[level], fe_collection); tg.join_all (); } } @@ -1481,7 +1472,7 @@ namespace hp unsigned int DoFHandler::max_couplings_between_dofs () const { - Assert (fe_collection, ExcNoFESelected()); + Assert(fe_collection.size() > 0, ExcNoFESelected()); return dealii::internal::hp::DoFHandlerImplementation::Implementation::max_couplings_between_dofs (*this); } @@ -1491,16 +1482,16 @@ namespace hp unsigned int DoFHandler::max_couplings_between_boundary_dofs () const { - Assert (fe_collection, ExcNoFESelected()); + Assert (fe_collection.size() > 0, ExcNoFESelected()); switch (dim) { case 1: - return fe_collection->max_dofs_per_vertex(); + return fe_collection.max_dofs_per_vertex(); case 2: - return (3*fe_collection->max_dofs_per_vertex() + return (3*fe_collection.max_dofs_per_vertex() + - 2*fe_collection->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 @@ -1511,9 +1502,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*fe_collection->max_dofs_per_vertex() + - 28*fe_collection->max_dofs_per_line() + - 8*fe_collection->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; @@ -1625,7 +1616,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 (fe_collection) + if (fe_collection.size() > 0) { cell_iterator cell = begin(), endc = end (); -- 2.39.5