From fab172f2c12619d56dce3596fc31d97f573454b6 Mon Sep 17 00:00:00 2001 From: Jean-Paul Pelteret Date: Thu, 6 Jan 2022 20:52:17 +0100 Subject: [PATCH] Add hp-capabilities to MeshWorker::ScratchData --- include/deal.II/meshworker/scratch_data.h | 266 +++++++++++-- source/meshworker/scratch_data.cc | 436 ++++++++++++++++++---- 2 files changed, 616 insertions(+), 86 deletions(-) diff --git a/include/deal.II/meshworker/scratch_data.h b/include/deal.II/meshworker/scratch_data.h index 261c4df13b..0c97cadd73 100644 --- a/include/deal.II/meshworker/scratch_data.h +++ b/include/deal.II/meshworker/scratch_data.h @@ -24,7 +24,10 @@ #include #include -#include + +#include +#include +#include #include @@ -274,7 +277,7 @@ namespace MeshWorker const UpdateFlags & neighbor_face_update_flags = update_default); /** - * Same as the other constructor, using the default MappingQ1. + * Same as the other constructor, using the default linear mapping. * * @param fe The finite element * @param quadrature The cell quadrature @@ -293,7 +296,7 @@ namespace MeshWorker const UpdateFlags & face_update_flags = update_default); /** - * Same as the other constructor, using the default MappingQ1. + * Same as the other constructor, using the default linear mapping. * * @param fe The finite element * @param quadrature The cell quadrature @@ -315,6 +318,101 @@ namespace MeshWorker const UpdateFlags & face_update_flags = update_default, const UpdateFlags & neighbor_face_update_flags = update_default); + /** + * Create an empty ScratchData object. A SmartPointer pointing to + * @p mapping_collection and @p fe_collection is stored internally. Make sure they live longer + * than this class instance. + * + * The constructor does not initialize any of the internal hp::FEValues + * objects. These are initialized the first time one of the reinit() + * functions is called, using the arguments passed here. + * + * @param mapping_collection The mapping collection to use in the internal hp::FEValues objects + * @param fe_collection The finite element collection + * @param cell_quadrature_collection The cell quadrature collection + * @param cell_update_flags UpdateFlags for the current cell hp::FEValues and + * neighbor cell hp::FEValues + * @param face_quadrature_collection The face quadrature collection, used for hp::FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + * @param face_update_flags UpdateFlags used for FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + */ + ScratchData(const hp::MappingCollection &mapping_collection, + const hp::FECollection & fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const hp::QCollection &face_quadrature_collection = + hp::QCollection(), + const UpdateFlags &face_update_flags = update_default); + + /** + * Similar to the other constructor, but this one allows to specify + * different flags for neighbor cells and faces. + * + * @param mapping_collection The mapping collection to use in the internal hp::FEValues objects + * @param fe_collection The finite element collection + * @param cell_quadrature_collection The cell quadrature collection + * @param cell_update_flags UpdateFlags for the current cell hp::FEValues + * @param neighbor_cell_update_flags UpdateFlags for the neighbor cell hp::FEValues + * @param face_quadrature_collection The face quadrature collection, used for hp::FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + * @param face_update_flags UpdateFlags used for FEFaceValues and + * hp::FESubfaceValues for the current cell + * @param neighbor_face_update_flags UpdateFlags used for hp::FEFaceValues and + * hp::FESubfaceValues for the neighbor cell + */ + ScratchData(const hp::MappingCollection &mapping_collection, + const hp::FECollection & fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const UpdateFlags & neighbor_cell_update_flags, + const hp::QCollection &face_quadrature_collection = + hp::QCollection(), + const UpdateFlags &face_update_flags = update_default, + const UpdateFlags &neighbor_face_update_flags = update_default); + + /** + * Same as the other constructor, using the default linear mapping. + * + * @param fe_collection The finite element collection + * @param cell_quadrature_collection The cell quadrature collection + * @param cell_update_flags UpdateFlags for the current cell hp::FEValues and + * neighbor cell hp::FEValues + * @param face_quadrature_collection The face quadrature collection, used for hp::FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + * @param face_update_flags UpdateFlags used for FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + */ + ScratchData(const hp::FECollection &fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const hp::QCollection &face_quadrature_collection = + hp::QCollection(), + const UpdateFlags &face_update_flags = update_default); + + /** + * Same as the other constructor, using the default linear mapping. + * + * @param fe_collection The finite element collection + * @param cell_quadrature_collection The cell quadrature collection + * @param cell_update_flags UpdateFlags for the current cell hp::FEValues + * @param neighbor_cell_update_flags UpdateFlags for the neighbor cell hp::FEValues + * @param face_quadrature_collection The face quadrature collection, used for hp::FEFaceValues and + * hp::FESubfaceValues for both the current cell and the neighbor cell + * @param face_update_flags UpdateFlags used for FEFaceValues and + * hp::FESubfaceValues for the current cell + * @param neighbor_face_update_flags UpdateFlags used for hp::FEFaceValues and + * hp::FESubfaceValues for the neighbor cell + */ + ScratchData(const hp::FECollection &fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const UpdateFlags & neighbor_cell_update_flags, + const hp::QCollection &face_quadrature_collection = + hp::QCollection(), + const UpdateFlags &face_update_flags = update_default, + const UpdateFlags &neighbor_face_update_flags = update_default); + /** * Deep copy constructor. FEValues objects are not copied. */ @@ -550,6 +648,37 @@ namespace MeshWorker const Quadrature & get_face_quadrature() const; + /** + * Return a reference to the used mapping. + */ + const hp::MappingCollection & + get_mapping_collection() const; + + /** + * Return a reference to the selected finite element object. + */ + const hp::FECollection & + get_fe_collection() const; + + /** + * Return a reference to the cell quadrature object in use. + */ + const hp::QCollection & + get_cell_quadrature_collection() const; + + /** + * Return a reference to the face quadrature object in use. + */ + const hp::QCollection & + get_face_quadrature_collection() const; + + /** + * Returns a boolean indicating whether or not this ScratchData object has + * hp-capabilities enabled. + */ + bool + has_hp_capabilities() const; + /** * Return the cell update flags set. */ @@ -1100,6 +1229,11 @@ namespace MeshWorker const unsigned int size, const Number & exemplar_number) const; + /** + * @name Data that supports the standard FE implementation + */ + /** @{ */ // non-hp data + /** * The mapping used by the internal FEValues. Make sure it lives * longer than this class. @@ -1125,61 +1259,151 @@ namespace MeshWorker Quadrature face_quadrature; /** - * UpdateFlags to use when initializing the cell FEValues object. + * Finite element values on the current cell. */ - UpdateFlags cell_update_flags; + std::unique_ptr> fe_values; /** - * UpdateFlags to use when initializing the neighbor cell FEValues objects. + * Finite element values on the current face. */ - UpdateFlags neighbor_cell_update_flags; + std::unique_ptr> fe_face_values; /** - * UpdateFlags to use when initializing FEFaceValues and FESubfaceValues - * objects. + * Finite element values on the current subface. */ - UpdateFlags face_update_flags; + std::unique_ptr> fe_subface_values; /** - * UpdateFlags to use when initializing neighbor FEFaceValues and - * FESubfaceValues objects. + * Finite element values on the neighbor cell. */ - UpdateFlags neighbor_face_update_flags; + std::unique_ptr> neighbor_fe_values; + + /** + * Finite element values on the neighbor face. + */ + std::unique_ptr> neighbor_fe_face_values; + + /** + * Finite element values on the neighbor subface. + */ + std::unique_ptr> neighbor_fe_subface_values; + + /** + * Interface values on facets. + */ + std::unique_ptr> interface_fe_values; + + /** @} */ // non-hp data + + /** + * @name Data that supports the hp-FE implementation + */ + /** @{ */ // hp data + + /** + * The mapping collection used by the internal hp::FEValues. Make sure it + * lives longer than this class. + */ + SmartPointer> mapping_collection; + + /** + * The finite element used by the internal FEValues. Make sure it lives + * longer than this class. + */ + SmartPointer> fe_collection; + + /** + * Quadrature formula used to integrate on the current cell, and on its + * neighbor. + */ + hp::QCollection cell_quadrature_collection; + + /** + * Quadrature formula used to integrate on faces, subfaces, and neighbor + * faces and subfaces. + */ + hp::QCollection face_quadrature_collection; + + /** + * Boolean indicating whether or not the current ScratchData has hp- + * capabilities. + */ + bool hp_capability_enabled; /** * Finite element values on the current cell. */ - std::unique_ptr> fe_values; + std::unique_ptr> hp_fe_values; /** * Finite element values on the current face. */ - std::unique_ptr> fe_face_values; + std::unique_ptr> hp_fe_face_values; /** * Finite element values on the current subface. */ - std::unique_ptr> fe_subface_values; + std::unique_ptr> hp_fe_subface_values; /** * Finite element values on the neighbor cell. */ - std::unique_ptr> neighbor_fe_values; + std::unique_ptr> neighbor_hp_fe_values; /** * Finite element values on the neighbor face. */ - std::unique_ptr> neighbor_fe_face_values; + std::unique_ptr> neighbor_hp_fe_face_values; /** * Finite element values on the neighbor subface. */ - std::unique_ptr> neighbor_fe_subface_values; + std::unique_ptr> + neighbor_hp_fe_subface_values; /** - * Interface values on facets. + * Exception used when a certain feature doesn't make sense when + * ScratchData does has hp-capabilities enabled. + * + * @ingroup Exceptions */ - std::unique_ptr> interface_fe_values; + DeclExceptionMsg(ExcOnlyAvailableWithoutHP, + "The current function doesn't make sense when used with a " + "ScratchData object with hp-capabilities."); + + /** + * Exception used when a certain feature doesn't make sense when + * ScratchData does not have hp-capabilities enabled. + * + * @ingroup Exceptions + */ + DeclExceptionMsg(ExcOnlyAvailableWithHP, + "The current function doesn't make sense when used with a " + "ScratchData object without hp-capabilities."); + + /** @} */ // hp data + + /** + * UpdateFlags to use when initializing the cell FEValues object. + */ + UpdateFlags cell_update_flags; + + /** + * UpdateFlags to use when initializing the neighbor cell FEValues objects. + */ + UpdateFlags neighbor_cell_update_flags; + + /** + * UpdateFlags to use when initializing FEFaceValues and FESubfaceValues + * objects. + */ + UpdateFlags face_update_flags; + + /** + * UpdateFlags to use when initializing neighbor FEFaceValues and + * FESubfaceValues objects. + */ + UpdateFlags neighbor_face_update_flags; /** * Dof indices on the current cell. diff --git a/source/meshworker/scratch_data.cc b/source/meshworker/scratch_data.cc index d8844ad2b8..208acece0b 100644 --- a/source/meshworker/scratch_data.cc +++ b/source/meshworker/scratch_data.cc @@ -33,6 +33,7 @@ namespace MeshWorker , fe(&fe) , cell_quadrature(quadrature) , face_quadrature(face_quadrature) + , hp_capability_enabled(false) , cell_update_flags(update_flags) , neighbor_cell_update_flags(update_flags) , face_update_flags(face_update_flags) @@ -57,6 +58,7 @@ namespace MeshWorker , fe(&fe) , cell_quadrature(quadrature) , face_quadrature(face_quadrature) + , hp_capability_enabled(false) , cell_update_flags(update_flags) , neighbor_cell_update_flags(neighbor_update_flags) , face_update_flags(face_update_flags) @@ -107,6 +109,94 @@ namespace MeshWorker + template + ScratchData::ScratchData( + const hp::MappingCollection &mapping_collection, + const hp::FECollection & fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const hp::QCollection & face_quadrature_collection, + const UpdateFlags & face_update_flags) + : mapping_collection(&mapping_collection) + , fe_collection(&fe_collection) + , cell_quadrature_collection(cell_quadrature_collection) + , face_quadrature_collection(face_quadrature_collection) + , hp_capability_enabled(true) + , cell_update_flags(cell_update_flags) + , neighbor_cell_update_flags(cell_update_flags) + , face_update_flags(face_update_flags) + , neighbor_face_update_flags(face_update_flags) + { + local_dof_indices.reserve(fe_collection.max_dofs_per_cell()); + neighbor_dof_indices.reserve(fe_collection.max_dofs_per_cell()); + } + + + + template + ScratchData::ScratchData( + const hp::MappingCollection &mapping_collection, + const hp::FECollection & fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const UpdateFlags & neighbor_cell_update_flags, + const hp::QCollection & face_quadrature_collection, + const UpdateFlags & face_update_flags, + const UpdateFlags & neighbor_face_update_flags) + : mapping_collection(&mapping_collection) + , fe_collection(&fe_collection) + , cell_quadrature_collection(cell_quadrature_collection) + , face_quadrature_collection(face_quadrature_collection) + , hp_capability_enabled(true) + , cell_update_flags(cell_update_flags) + , neighbor_cell_update_flags(neighbor_cell_update_flags) + , face_update_flags(face_update_flags) + , neighbor_face_update_flags(neighbor_face_update_flags) + { + local_dof_indices.reserve(fe_collection.max_dofs_per_cell()); + neighbor_dof_indices.reserve(fe_collection.max_dofs_per_cell()); + } + + + + template + ScratchData::ScratchData( + const hp::FECollection &fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const hp::QCollection & face_quadrature_collection, + const UpdateFlags & face_update_flags) + : ScratchData(fe_collection.get_reference_cell_default_linear_mapping(), + fe_collection, + cell_quadrature_collection, + cell_update_flags, + face_quadrature_collection, + face_update_flags) + {} + + + + template + ScratchData::ScratchData( + const hp::FECollection &fe_collection, + const hp::QCollection & cell_quadrature_collection, + const UpdateFlags & cell_update_flags, + const UpdateFlags & neighbor_cell_update_flags, + const hp::QCollection & face_quadrature_collection, + const UpdateFlags & face_update_flags, + const UpdateFlags & neighbor_face_update_flags) + : ScratchData(fe_collection.get_reference_cell_default_linear_mapping(), + fe_collection, + cell_quadrature_collection, + cell_update_flags, + neighbor_cell_update_flags, + face_quadrature_collection, + face_update_flags, + neighbor_face_update_flags) + {} + + + template ScratchData::ScratchData( const ScratchData &scratch) @@ -114,6 +204,11 @@ namespace MeshWorker , fe(scratch.fe) , cell_quadrature(scratch.cell_quadrature) , face_quadrature(scratch.face_quadrature) + , mapping_collection(scratch.mapping_collection) + , fe_collection(scratch.fe_collection) + , cell_quadrature_collection(scratch.cell_quadrature_collection) + , face_quadrature_collection(scratch.face_quadrature_collection) + , hp_capability_enabled(scratch.hp_capability_enabled) , cell_update_flags(scratch.cell_update_flags) , neighbor_cell_update_flags(scratch.neighbor_cell_update_flags) , face_update_flags(scratch.face_update_flags) @@ -131,17 +226,39 @@ namespace MeshWorker ScratchData::reinit( const typename DoFHandler::active_cell_iterator &cell) { - if (!fe_values) - fe_values = std::make_unique>(*mapping, - *fe, - cell_quadrature, - cell_update_flags); + if (hp_capability_enabled == false) + { + if (!fe_values) + fe_values = std::make_unique>( + *mapping, *fe, cell_quadrature, cell_update_flags); + + fe_values->reinit(cell); + local_dof_indices.resize(fe_values->dofs_per_cell); + cell->get_dof_indices(local_dof_indices); + current_fe_values = fe_values.get(); + return *fe_values; + } + else + { + if (!hp_fe_values) + hp_fe_values = std::make_unique>( + *mapping_collection, + *fe_collection, + cell_quadrature_collection, + cell_update_flags); + + hp_fe_values->reinit(cell); + const auto &fe_values = hp_fe_values->get_present_fe_values(); + + AssertDimension( + (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(), + fe_values.dofs_per_cell); + local_dof_indices.resize(fe_values.dofs_per_cell); + cell->get_dof_indices(local_dof_indices); - fe_values->reinit(cell); - local_dof_indices.resize(fe_values->dofs_per_cell); - cell->get_dof_indices(local_dof_indices); - current_fe_values = fe_values.get(); - return *fe_values; + current_fe_values = &fe_values; + return fe_values; + } } @@ -152,15 +269,37 @@ namespace MeshWorker const typename DoFHandler::active_cell_iterator &cell, const unsigned int face_no) { - if (!fe_face_values) - fe_face_values = std::make_unique>( - *mapping, *fe, face_quadrature, face_update_flags); + if (hp_capability_enabled == false) + { + if (!fe_face_values) + fe_face_values = std::make_unique>( + *mapping, *fe, face_quadrature, face_update_flags); - fe_face_values->reinit(cell, face_no); - local_dof_indices.resize(fe->n_dofs_per_cell()); - cell->get_dof_indices(local_dof_indices); - current_fe_values = fe_face_values.get(); - return *fe_face_values; + fe_face_values->reinit(cell, face_no); + local_dof_indices.resize(fe->n_dofs_per_cell()); + cell->get_dof_indices(local_dof_indices); + current_fe_values = fe_face_values.get(); + return *fe_face_values; + } + else + { + if (!hp_fe_face_values) + hp_fe_face_values = std::make_unique>( + *mapping_collection, + *fe_collection, + face_quadrature_collection, + face_update_flags); + + hp_fe_face_values->reinit(cell, face_no); + const auto &fe_face_values = hp_fe_face_values->get_present_fe_values(); + const auto &fe = (*fe_collection)[cell->active_fe_index()]; + + local_dof_indices.resize(fe.n_dofs_per_cell()); + cell->get_dof_indices(local_dof_indices); + + current_fe_values = &fe_face_values; + return fe_face_values; + } } @@ -174,15 +313,40 @@ namespace MeshWorker { if (subface_no != numbers::invalid_unsigned_int) { - if (!fe_subface_values) - fe_subface_values = std::make_unique>( - *mapping, *fe, face_quadrature, face_update_flags); - fe_subface_values->reinit(cell, face_no, subface_no); - local_dof_indices.resize(fe->n_dofs_per_cell()); - cell->get_dof_indices(local_dof_indices); - - current_fe_values = fe_subface_values.get(); - return *fe_subface_values; + if (hp_capability_enabled == false) + { + if (!fe_subface_values) + fe_subface_values = + std::make_unique>( + *mapping, *fe, face_quadrature, face_update_flags); + fe_subface_values->reinit(cell, face_no, subface_no); + local_dof_indices.resize(fe->n_dofs_per_cell()); + cell->get_dof_indices(local_dof_indices); + + current_fe_values = fe_subface_values.get(); + return *fe_subface_values; + } + else + { + if (!hp_fe_subface_values) + hp_fe_subface_values = + std::make_unique>( + *mapping_collection, + *fe_collection, + face_quadrature_collection, + face_update_flags); + + hp_fe_subface_values->reinit(cell, face_no, subface_no); + const auto &fe_subface_values = + hp_fe_subface_values->get_present_fe_values(); + const auto &fe = (*fe_collection)[cell->active_fe_index()]; + + local_dof_indices.resize(fe.n_dofs_per_cell()); + cell->get_dof_indices(local_dof_indices); + + current_fe_values = &fe_subface_values; + return fe_subface_values; + } } else return reinit(cell, face_no); @@ -201,22 +365,32 @@ namespace MeshWorker const unsigned int face_no_neighbor, const unsigned int sub_face_no_neighbor) { - if (!interface_fe_values) - interface_fe_values = std::make_unique>( - *mapping, *fe, face_quadrature, face_update_flags); - interface_fe_values->reinit(cell, - face_no, - sub_face_no, - cell_neighbor, - face_no_neighbor, - sub_face_no_neighbor); - - current_fe_values = &interface_fe_values->get_fe_face_values(0); - current_neighbor_fe_values = &interface_fe_values->get_fe_face_values(1); - - cell_neighbor->get_dof_indices(neighbor_dof_indices); - local_dof_indices = interface_fe_values->get_interface_dof_indices(); - return *interface_fe_values; + if (hp_capability_enabled == false) + { + if (!interface_fe_values) + interface_fe_values = + std::make_unique>( + *mapping, *fe, face_quadrature, face_update_flags); + interface_fe_values->reinit(cell, + face_no, + sub_face_no, + cell_neighbor, + face_no_neighbor, + sub_face_no_neighbor); + + current_fe_values = &interface_fe_values->get_fe_face_values(0); + current_neighbor_fe_values = + &interface_fe_values->get_fe_face_values(1); + + cell_neighbor->get_dof_indices(neighbor_dof_indices); + local_dof_indices = interface_fe_values->get_interface_dof_indices(); + return *interface_fe_values; + } + else + { + AssertThrow(false, ExcOnlyAvailableWithoutHP()); + return *interface_fe_values; + } } @@ -226,14 +400,39 @@ namespace MeshWorker ScratchData::reinit_neighbor( const typename DoFHandler::active_cell_iterator &cell) { - if (!neighbor_fe_values) - neighbor_fe_values = std::make_unique>( - *mapping, *fe, cell_quadrature, neighbor_cell_update_flags); + if (hp_capability_enabled == false) + { + if (!neighbor_fe_values) + neighbor_fe_values = std::make_unique>( + *mapping, *fe, cell_quadrature, neighbor_cell_update_flags); - neighbor_fe_values->reinit(cell); - cell->get_dof_indices(neighbor_dof_indices); - current_neighbor_fe_values = neighbor_fe_values.get(); - return *neighbor_fe_values; + neighbor_fe_values->reinit(cell); + cell->get_dof_indices(neighbor_dof_indices); + current_neighbor_fe_values = neighbor_fe_values.get(); + return *neighbor_fe_values; + } + else + { + if (!neighbor_hp_fe_values) + neighbor_hp_fe_values = std::make_unique>( + *mapping_collection, + *fe_collection, + cell_quadrature_collection, + neighbor_cell_update_flags); + + neighbor_hp_fe_values->reinit(cell); + const auto &neighbor_fe_values = + neighbor_hp_fe_values->get_present_fe_values(); + + AssertDimension( + (*fe_collection)[cell->active_fe_index()].n_dofs_per_cell(), + neighbor_fe_values.dofs_per_cell); + neighbor_dof_indices.resize(neighbor_fe_values.dofs_per_cell); + cell->get_dof_indices(neighbor_dof_indices); + + current_neighbor_fe_values = &neighbor_fe_values; + return neighbor_fe_values; + } } @@ -244,13 +443,38 @@ namespace MeshWorker const typename DoFHandler::active_cell_iterator &cell, const unsigned int face_no) { - if (!neighbor_fe_face_values) - neighbor_fe_face_values = std::make_unique>( - *mapping, *fe, face_quadrature, neighbor_face_update_flags); - neighbor_fe_face_values->reinit(cell, face_no); - cell->get_dof_indices(neighbor_dof_indices); - current_neighbor_fe_values = neighbor_fe_face_values.get(); - return *neighbor_fe_face_values; + if (hp_capability_enabled == false) + { + if (!neighbor_fe_face_values) + neighbor_fe_face_values = + std::make_unique>( + *mapping, *fe, face_quadrature, neighbor_face_update_flags); + neighbor_fe_face_values->reinit(cell, face_no); + cell->get_dof_indices(neighbor_dof_indices); + current_neighbor_fe_values = neighbor_fe_face_values.get(); + return *neighbor_fe_face_values; + } + else + { + if (!neighbor_hp_fe_face_values) + neighbor_hp_fe_face_values = + std::make_unique>( + *mapping_collection, + *fe_collection, + face_quadrature_collection, + neighbor_face_update_flags); + + neighbor_hp_fe_face_values->reinit(cell, face_no); + const auto &neighbor_fe_face_values = + neighbor_hp_fe_face_values->get_present_fe_values(); + const auto &neighbor_fe = (*fe_collection)[cell->active_fe_index()]; + + neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell()); + cell->get_dof_indices(neighbor_dof_indices); + + current_neighbor_fe_values = &neighbor_fe_face_values; + return neighbor_fe_face_values; + } } @@ -264,14 +488,39 @@ namespace MeshWorker { if (subface_no != numbers::invalid_unsigned_int) { - if (!neighbor_fe_subface_values) - neighbor_fe_subface_values = - std::make_unique>( - *mapping, *fe, face_quadrature, neighbor_face_update_flags); - neighbor_fe_subface_values->reinit(cell, face_no, subface_no); - cell->get_dof_indices(neighbor_dof_indices); - current_neighbor_fe_values = neighbor_fe_subface_values.get(); - return *neighbor_fe_subface_values; + if (hp_capability_enabled == false) + { + if (!neighbor_fe_subface_values) + neighbor_fe_subface_values = + std::make_unique>( + *mapping, *fe, face_quadrature, neighbor_face_update_flags); + neighbor_fe_subface_values->reinit(cell, face_no, subface_no); + cell->get_dof_indices(neighbor_dof_indices); + current_neighbor_fe_values = neighbor_fe_subface_values.get(); + return *neighbor_fe_subface_values; + } + else + { + if (!neighbor_hp_fe_subface_values) + neighbor_hp_fe_subface_values = + std::make_unique>( + *mapping_collection, + *fe_collection, + face_quadrature_collection, + neighbor_face_update_flags); + + neighbor_hp_fe_subface_values->reinit(cell, face_no, subface_no); + const auto &neighbor_fe_subface_values = + neighbor_hp_fe_subface_values->get_present_fe_values(); + const auto &neighbor_fe = (*fe_collection)[cell->active_fe_index()]; + + neighbor_dof_indices.resize(neighbor_fe.n_dofs_per_cell()); + cell->get_dof_indices(neighbor_dof_indices); + cell->get_dof_indices(neighbor_dof_indices); + + current_neighbor_fe_values = &neighbor_fe_subface_values; + return neighbor_fe_subface_values; + } } else return reinit_neighbor(cell, face_no); @@ -400,6 +649,8 @@ namespace MeshWorker const Mapping & ScratchData::get_mapping() const { + Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP()); + Assert(mapping, ExcNotInitialized()); return *mapping; } @@ -409,6 +660,8 @@ namespace MeshWorker const FiniteElement & ScratchData::get_fe() const { + Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP()); + Assert(fe, ExcNotInitialized()); return *fe; } @@ -418,6 +671,7 @@ namespace MeshWorker const Quadrature & ScratchData::get_cell_quadrature() const { + Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP()); return cell_quadrature; } @@ -427,11 +681,63 @@ namespace MeshWorker const Quadrature & ScratchData::get_face_quadrature() const { + Assert(hp_capability_enabled == false, ExcOnlyAvailableWithoutHP()); return face_quadrature; } + template + const hp::MappingCollection & + ScratchData::get_mapping_collection() const + { + Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP()); + Assert(mapping_collection, ExcNotInitialized()); + return *mapping_collection; + } + + + + template + const hp::FECollection & + ScratchData::get_fe_collection() const + { + Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP()); + Assert(fe_collection, ExcNotInitialized()); + return *fe_collection; + } + + + + template + const hp::QCollection & + ScratchData::get_cell_quadrature_collection() const + { + Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP()); + return cell_quadrature_collection; + } + + + + template + const hp::QCollection & + ScratchData::get_face_quadrature_collection() const + { + Assert(hp_capability_enabled == true, ExcOnlyAvailableWithHP()); + return face_quadrature_collection; + } + + + + template + bool + ScratchData::has_hp_capabilities() const + { + return hp_capability_enabled; + } + + + template UpdateFlags ScratchData::get_cell_update_flags() const -- 2.39.5