From 84bd78408bd9b73288c6aac389f76426d4eb7f3f Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Fri, 12 Nov 2021 23:54:21 +0100 Subject: [PATCH] FEValuesBase: move content around --- include/deal.II/fe/fe_values.h | 87 ++++++++++-- source/fe/fe_values.cc | 237 ++++++++++++++------------------- 2 files changed, 181 insertions(+), 143 deletions(-) diff --git a/include/deal.II/fe/fe_values.h b/include/deal.II/fe/fe_values.h index 2af6135b0d..3b6977c482 100644 --- a/include/deal.II/fe/fe_values.h +++ b/include/deal.II/fe/fe_values.h @@ -3795,6 +3795,17 @@ public: "only works for those. See FiniteElement::is_primitive() for more information."); protected: + /** + * Reset a unique_ptr. If we can, do not de-allocate the previously + * held memory but re-use it for the next item to avoid the repeated + * memory allocation. We do this because FEValues objects are heavily + * used in multithreaded contexts where memory allocations are evil. + */ + template + static void + reset_pointer_in_place_if_possible(std::unique_ptr &present_cell, + const Iterator & new_cell); + /** * Objects of the FEValues class need to store an iterator * to the present cell in order to be able to extract the values of the @@ -3825,15 +3836,75 @@ protected: * type * erasure design pattern. */ - class CellIteratorBase; + class CellIteratorBase + { + public: + DeclExceptionMsg( + ExcNeedsDoFHandler, + "You have previously called the FEValues::reinit function with a " + "cell iterator of type Triangulation::cell_iterator. However, " + "when you do this, you cannot call some functions in the FEValues " + "class, such as the get_function_values/gradients/hessians/third_derivatives " + "functions. If you need these functions, then you need to call " + "FEValues::reinit with an iterator type that allows to extract " + "degrees of freedom, such as DoFHandler::cell_iterator."); - /** - * Forward declaration of classes derived from CellIteratorBase. Their - * definition and implementation is given in the .cc file. - */ - template - class CellIterator; - class TriaCellIterator; + /** + * Constructor. + */ + template + CellIteratorBase( + const TriaIterator> &cell); + + /** + * Constructor. + */ + CellIteratorBase( + const typename Triangulation::cell_iterator &cell); + + /** + * Destructor. + */ + virtual ~CellIteratorBase() = default; + + /** + * Conversion operator to an iterator for triangulations. This + * conversion is implicit for the original iterators, since they are derived + * classes. However, since here we have kind of a parallel class hierarchy, + * we have to have a conversion operator. + */ + operator typename Triangulation::cell_iterator() const; + + /** + * Return the number of degrees of freedom the DoF + * handler object has to which the iterator belongs to. + */ + types::global_dof_index + n_dofs_for_dof_handler() const; + + /** + * Call @p get_interpolated_dof_values of the iterator with the + * given arguments. + */ + template + void + get_interpolated_dof_values( + const VectorType & in, + Vector &out) const; + + /** + * Call @p get_interpolated_dof_values of the iterator with the + * given arguments. + */ + void + get_interpolated_dof_values(const IndexSet & in, + Vector &out) const; + + private: + const typename Triangulation::cell_iterator cell; + const DoFHandler * dof_handler; + const bool level_dof_access; + }; /** * Store the cell selected last time the reinit() function was called. This diff --git a/source/fe/fe_values.cc b/source/fe/fe_values.cc index 0cec8de637..4a6b14bc7c 100644 --- a/source/fe/fe_values.cc +++ b/source/fe/fe_values.cc @@ -2592,119 +2592,91 @@ namespace internal /* ---------------- FEValuesBase::CellIteratorBase --------- */ template -class FEValuesBase::CellIteratorBase +template +FEValuesBase::CellIteratorBase::CellIteratorBase( + const TriaIterator> &cell) + : cell(cell) + , dof_handler(&cell->get_dof_handler()) + , level_dof_access(lda) +{} + + + +template +FEValuesBase::CellIteratorBase::CellIteratorBase( + const typename Triangulation::cell_iterator &cell) + : cell(cell) + , dof_handler(nullptr) + , level_dof_access(false) +{} + + + +template +FEValuesBase::CellIteratorBase:: +operator typename Triangulation::cell_iterator() const { -public: - DeclExceptionMsg( - ExcNeedsDoFHandler, - "You have previously called the FEValues::reinit function with a " - "cell iterator of type Triangulation::cell_iterator. However, " - "when you do this, you cannot call some functions in the FEValues " - "class, such as the get_function_values/gradients/hessians/third_derivatives " - "functions. If you need these functions, then you need to call " - "FEValues::reinit with an iterator type that allows to extract " - "degrees of freedom, such as DoFHandler::cell_iterator."); - - /** - * Constructor. - */ - template - CellIteratorBase( - const TriaIterator> &cell) - : cell(cell) - , dof_handler(&cell->get_dof_handler()) - , level_dof_access(lda) - {} + return cell; +} - /** - * Constructor. - */ - CellIteratorBase( - const typename Triangulation::cell_iterator &cell) - : cell(cell) - , dof_handler(nullptr) - , level_dof_access(false) - {} - /** - * Destructor. - */ - virtual ~CellIteratorBase() = default; - - /** - * Conversion operator to an iterator for triangulations. This - * conversion is implicit for the original iterators, since they are derived - * classes. However, since here we have kind of a parallel class hierarchy, - * we have to have a conversion operator. - */ - operator typename Triangulation::cell_iterator() const - { - return cell; - } - /** - * Return the number of degrees of freedom the DoF - * handler object has to which the iterator belongs to. - */ - types::global_dof_index - n_dofs_for_dof_handler() const - { - Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); - return dof_handler->n_dofs(); - } +template +types::global_dof_index +FEValuesBase::CellIteratorBase::n_dofs_for_dof_handler() const +{ + Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); + return dof_handler->n_dofs(); +} - /** - * Call @p get_interpolated_dof_values of the iterator with the - * given arguments. - */ - template - void - get_interpolated_dof_values( - const VectorType & in, - Vector &out) const - { - Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); - - if (level_dof_access) - DoFCellAccessor(&cell->get_triangulation(), - cell->level(), - cell->index(), - dof_handler) - .get_interpolated_dof_values(in, out); - else - DoFCellAccessor(&cell->get_triangulation(), - cell->level(), - cell->index(), - dof_handler) - .get_interpolated_dof_values(in, out); - } - /** - * Call @p get_interpolated_dof_values of the iterator with the - * given arguments. - */ - void - get_interpolated_dof_values(const IndexSet & in, - Vector &out) const - { - Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); - Assert(level_dof_access == false, ExcNotImplemented()); - const DoFCellAccessor cell_dofs( - &cell->get_triangulation(), cell->level(), cell->index(), dof_handler); +template +template +void +FEValuesBase::CellIteratorBase::get_interpolated_dof_values( + const VectorType & in, + Vector &out) const +{ + Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); + + if (level_dof_access) + DoFCellAccessor(&cell->get_triangulation(), + cell->level(), + cell->index(), + dof_handler) + .get_interpolated_dof_values(in, out); + else + DoFCellAccessor(&cell->get_triangulation(), + cell->level(), + cell->index(), + dof_handler) + .get_interpolated_dof_values(in, out); +} - std::vector dof_indices( - cell_dofs.get_fe().n_dofs_per_cell()); - cell_dofs.get_dof_indices(dof_indices); - for (unsigned int i = 0; i < cell_dofs.get_fe().n_dofs_per_cell(); ++i) - out[i] = (in.is_element(dof_indices[i]) ? 1 : 0); - } - const typename Triangulation::cell_iterator cell; - const DoFHandler * dof_handler; - const bool level_dof_access; -}; +template +void +FEValuesBase::CellIteratorBase::get_interpolated_dof_values( + const IndexSet & in, + Vector &out) const +{ + Assert(dof_handler != nullptr, ExcNeedsDoFHandler()); + Assert(level_dof_access == false, ExcNotImplemented()); + + const DoFCellAccessor cell_dofs( + &cell->get_triangulation(), cell->level(), cell->index(), dof_handler); + + std::vector dof_indices( + cell_dofs.get_fe().n_dofs_per_cell()); + cell_dofs.get_dof_indices(dof_indices); + + for (unsigned int i = 0; i < cell_dofs.get_fe().n_dofs_per_cell(); ++i) + out[i] = (in.is_element(dof_indices[i]) ? 1 : 0); +} + + namespace internal { @@ -4293,34 +4265,29 @@ FEValues::initialize(const UpdateFlags update_flags) -namespace +template +template +void +FEValuesBase::reset_pointer_in_place_if_possible( + std::unique_ptr &present_cell, + const Iterator & new_cell) { - // Reset a unique_ptr. If we can, do not de-allocate the previously - // held memory but re-use it for the next item to avoid the repeated - // memory allocation. We do this because FEValues objects are heavily - // used in multithreaded contexts where memory allocations are evil. - template - void - reset_pointer_in_place_if_possible(std::unique_ptr &present_cell, - const Iterator & new_cell) - { - // see if the existing pointer is non-null and if the type of - // the old object pointed to matches that of the one we'd - // like to create - if (present_cell.get() && (typeid(*present_cell.get()) == typeid(Type))) - { - // call destructor of the old object - static_cast(present_cell.get())->~Type(); + // see if the existing pointer is non-null and if the type of + // the old object pointed to matches that of the one we'd + // like to create + if (present_cell.get() && (typeid(*present_cell.get()) == typeid(Type))) + { + // call destructor of the old object + static_cast(present_cell.get())->~Type(); - // then construct a new object in-place - new (const_cast(static_cast(present_cell.get()))) - Type(new_cell); - } - else - // if the types don't match, there is nothing we can do here - present_cell = std::make_unique(new_cell); - } -} // namespace + // then construct a new object in-place + new (const_cast(static_cast(present_cell.get()))) + Type(new_cell); + } + else + // if the types don't match, there is nothing we can do here + present_cell = std::make_unique(new_cell); +} @@ -4341,7 +4308,7 @@ FEValues::reinit( this->maybe_invalidate_previous_present_cell(cell); this->check_cell_similarity(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); @@ -4375,7 +4342,7 @@ FEValues::reinit( this->maybe_invalidate_previous_present_cell(cell); this->check_cell_similarity(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); @@ -4642,7 +4609,7 @@ FEFaceValues::reinit( AssertIndexRange(face_no, GeometryInfo::faces_per_cell); this->maybe_invalidate_previous_present_cell(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); @@ -4676,7 +4643,7 @@ FEFaceValues::reinit( AssertIndexRange(face_no, GeometryInfo::faces_per_cell); this->maybe_invalidate_previous_present_cell(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); @@ -4889,7 +4856,7 @@ FESubfaceValues::reinit( "instead in these cases.")); this->maybe_invalidate_previous_present_cell(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); @@ -4936,7 +4903,7 @@ FESubfaceValues::reinit( cell->face(face_no)->n_children())); this->maybe_invalidate_previous_present_cell(cell); - reset_pointer_in_place_if_possible< + this->template reset_pointer_in_place_if_possible< typename FEValuesBase::CellIteratorBase>(this->present_cell, cell); -- 2.39.5