45 template <
typename VectorType>
51 const VectorType &level_set,
54 , fe_collection(&fe_collection)
55 , q_collection_1D(quadrature)
56 , region_update_flags(region_update_flags)
57 , mesh_classifier(&mesh_classifier)
58 , quadrature_generator(q_collection_1D,
75 template <
typename VectorType>
83 const VectorType &level_set,
85 : mapping_collection(&mapping_collection)
86 , fe_collection(&fe_collection)
87 , q_collection_1D(q_collection_1D)
88 , region_update_flags(region_update_flags)
89 , mesh_classifier(&mesh_classifier)
90 , quadrature_generator(q_collection_1D,
107 Assert(fe_collection->size() > 0,
108 ExcMessage(
"Incoming hp::FECollection can not be empty."));
109 Assert(mapping_collection->size() == fe_collection->size() ||
110 mapping_collection->size() == 1,
111 ExcMessage(
"Size of hp::MappingCollection must be "
112 "the same as hp::FECollection or 1."));
113 Assert(q_collection.
size() == fe_collection->size() ||
114 q_collection.
size() == 1,
115 ExcMessage(
"Size of hp::QCollection<dim> must be the "
116 "same as hp::FECollection or 1."));
117 Assert(q_collection_1D.size() == fe_collection->size() ||
118 q_collection_1D.size() == 1,
119 ExcMessage(
"Size of hp::QCollection<1> must be the "
120 "same as hp::FECollection or 1."));
124 fe_values_inside_full_quadrature.resize(fe_collection->size());
125 fe_values_outside_full_quadrature.resize(fe_collection->size());
129 const unsigned int mapping_index =
130 mapping_collection->size() > 1 ?
fe_index : 0;
131 const unsigned int q_index = q_collection.
size() > 1 ?
fe_index : 0;
133 fe_values_inside_full_quadrature[
fe_index].emplace(
134 (*mapping_collection)[mapping_index],
136 q_collection[q_index],
137 region_update_flags.inside);
138 fe_values_outside_full_quadrature[
fe_index].emplace(
139 (*mapping_collection)[mapping_index],
141 q_collection[q_index],
142 region_update_flags.outside);
149 template <
bool level_dof_access>
154 current_cell_location = mesh_classifier->location_to_level_set(cell);
155 active_fe_index = cell->active_fe_index();
159 fe_values_inside.reset();
160 fe_values_surface.reset();
161 fe_values_outside.reset();
163 switch (current_cell_location)
167 fe_values_inside_full_quadrature.at(active_fe_index)->reinit(cell);
172 fe_values_outside_full_quadrature.at(active_fe_index)->reinit(cell);
177 const unsigned int mapping_index =
178 mapping_collection->size() > 1 ? active_fe_index : 0;
180 const unsigned int q1D_index =
181 q_collection_1D.size() > 1 ? active_fe_index : 0;
182 quadrature_generator.set_1D_quadrature(q1D_index);
183 quadrature_generator.generate(cell);
186 quadrature_generator.get_inside_quadrature();
188 quadrature_generator.get_outside_quadrature();
190 quadrature_generator.get_surface_quadrature();
195 if (inside_quadrature.
size() > 0)
197 fe_values_inside.emplace((*mapping_collection)[mapping_index],
198 (*fe_collection)[active_fe_index],
200 region_update_flags.inside);
202 fe_values_inside->reinit(cell);
205 if (outside_quadrature.
size() > 0)
207 fe_values_outside.emplace((*mapping_collection)[mapping_index],
208 (*fe_collection)[active_fe_index],
210 region_update_flags.outside);
212 fe_values_outside->reinit(cell);
215 if (surface_quadrature.
size() > 0)
217 fe_values_surface.emplace((*mapping_collection)[mapping_index],
218 (*fe_collection)[active_fe_index],
220 region_update_flags.surface);
221 fe_values_surface->reinit(cell);
237 const std::optional<::FEValues<dim>> &
241 return fe_values_inside_full_quadrature.at(active_fe_index);
243 return fe_values_inside;
249 const std::optional<::FEValues<dim>> &
253 return fe_values_outside_full_quadrature.at(active_fe_index);
255 return fe_values_outside;
261 const std::optional<FEImmersedSurfaceValues<dim>> &
264 return fe_values_surface;
270 template <
typename VectorType>
277 const VectorType &level_set,
280 , fe_collection(&fe_collection)
281 , q_collection_1D(quadrature)
282 , region_update_flags(region_update_flags)
283 , mesh_classifier(&mesh_classifier)
284 , face_quadrature_generator(q_collection_1D,
301 template <
typename VectorType>
310 const VectorType &level_set,
312 : mapping_collection(&mapping_collection)
313 , fe_collection(&fe_collection)
314 , q_collection_1D(q_collection_1D)
315 , region_update_flags(region_update_flags)
316 , mesh_classifier(&mesh_classifier)
317 , face_quadrature_generator(q_collection_1D,
335 Assert(fe_collection->size() > 0,
336 ExcMessage(
"Incoming hp::FECollection can not be empty."));
338 mapping_collection->size() == fe_collection->size() ||
339 mapping_collection->size() == 1,
341 "Size of hp::MappingCollection must be the same as hp::FECollection or 1."));
343 q_collection.
size() == fe_collection->size() || q_collection.
size() == 1,
345 "Size of hp::QCollection<dim> must be the same as hp::FECollection or 1."));
347 q_collection_1D.size() == fe_collection->size() ||
348 q_collection_1D.size() == 1,
350 "Size of hp::QCollection<1> must be the same as hp::FECollection or 1."));
354 fe_values_inside_full_quadrature.resize(fe_collection->size());
355 fe_values_outside_full_quadrature.resize(fe_collection->size());
359 const unsigned int mapping_index =
360 mapping_collection->size() > 1 ?
fe_index : 0;
361 const unsigned int q_index = q_collection.
size() > 1 ?
fe_index : 0;
363 fe_values_inside_full_quadrature[
fe_index].emplace(
364 (*mapping_collection)[mapping_index],
366 q_collection[q_index],
367 region_update_flags.inside);
368 fe_values_outside_full_quadrature[
fe_index].emplace(
369 (*mapping_collection)[mapping_index],
371 q_collection[q_index],
372 region_update_flags.outside);
379 template <
bool level_dof_access>
383 const unsigned int face_no,
386 current_face_location =
387 mesh_classifier->location_to_level_set(cell, face_no);
388 active_fe_index = cell->active_fe_index();
392 fe_values_inside.reset();
393 fe_values_outside.reset();
395 switch (current_face_location)
399 call_reinit(*fe_values_inside_full_quadrature.at(active_fe_index));
404 call_reinit(*fe_values_outside_full_quadrature.at(active_fe_index));
409 const unsigned int mapping_index =
410 mapping_collection->size() > 1 ? active_fe_index : 0;
411 const unsigned int q1D_index =
412 q_collection_1D.size() > 1 ? active_fe_index : 0;
414 face_quadrature_generator.set_1D_quadrature(q1D_index);
415 face_quadrature_generator.generate(cell, face_no);
417 const Quadrature<dim - 1> &inside_quadrature =
418 face_quadrature_generator.get_inside_quadrature();
419 const Quadrature<dim - 1> &outside_quadrature =
420 face_quadrature_generator.get_outside_quadrature();
425 if (inside_quadrature.size() > 0)
427 fe_values_inside.emplace((*mapping_collection)[mapping_index],
428 (*fe_collection)[active_fe_index],
430 region_update_flags.inside);
432 call_reinit(*fe_values_inside);
435 if (outside_quadrature.size() > 0)
437 fe_values_outside.emplace((*mapping_collection)[mapping_index],
438 (*fe_collection)[active_fe_index],
440 region_update_flags.outside);
442 call_reinit(*fe_values_outside);
457 const std::optional<::FEInterfaceValues<dim>> &
461 return fe_values_inside_full_quadrature.at(active_fe_index);
463 return fe_values_inside;
469 const std::optional<::FEInterfaceValues<dim>> &
473 return fe_values_outside_full_quadrature.at(active_fe_index);
475 return fe_values_outside;
479 #include "fe_values.inst"
const std::optional<::FEInterfaceValues< dim > > & get_outside_fe_values() const
void initialize(const hp::QCollection< dim - 1 > &q_collection)
const hp::QCollection< 1 > q_collection_1D
typename FaceQuadratureGenerator< dim >::AdditionalData AdditionalData
void do_reinit(const TriaIterator< DoFCellAccessor< dim, dim, level_dof_access >> &cell, const unsigned int face_no, const std::function< void(::FEInterfaceValues< dim > &)> &call_reinit)
const std::optional<::FEInterfaceValues< dim > > & get_inside_fe_values() const
FEInterfaceValues(const hp::FECollection< dim > &fe_collection, const Quadrature< 1 > &quadrature, const RegionUpdateFlags region_update_flags, const MeshClassifier< dim > &mesh_classifier, const DoFHandler< dim > &dof_handler, const VectorType &level_set, const AdditionalData &additional_data=AdditionalData())
void reinit(const TriaIterator< DoFCellAccessor< dim, dim, level_dof_access >> &cell)
const std::optional<::FEValues< dim > > & get_outside_fe_values() const
void initialize(const hp::QCollection< dim > &q_collection)
const std::optional< FEImmersedSurfaceValues< dim > > & get_surface_fe_values() const
const hp::QCollection< 1 > q_collection_1D
typename QuadratureGenerator< dim >::AdditionalData AdditionalData
FEValues(const hp::FECollection< dim > &fe_collection, const Quadrature< 1 > &quadrature, const RegionUpdateFlags region_update_flags, const MeshClassifier< dim > &mesh_classifier, const DoFHandler< dim > &dof_handler, const VectorType &level_set, const AdditionalData &additional_data=AdditionalData())
const std::optional<::FEValues< dim > > & get_inside_fe_values() const
unsigned int size() const
unsigned int size() const
void push_back(const Quadrature< dim_in > &new_quadrature)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
@ update_default
No update.
static const unsigned int invalid_unsigned_int
unsigned short int fe_index