Reference documentation for deal.II version 9.6.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
fe_values.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2021 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_non_matching_fe_values
16#define dealii_non_matching_fe_values
17
20
22
26
28
32
36
37#include <deque>
38#include <optional>
39
41
42namespace NonMatching
43{
79
80
142 template <int dim>
144 {
145 public:
147
169 template <typename VectorType>
171 const Quadrature<1> &quadrature,
174 const DoFHandler<dim> &dof_handler,
175 const VectorType &level_set,
176 const AdditionalData &additional_data = AdditionalData());
177
204 template <typename VectorType>
207 const hp::QCollection<dim> &q_collection,
208 const hp::QCollection<1> &q_collection_1d,
211 const DoFHandler<dim> &dof_handler,
212 const VectorType &level_set,
213 const AdditionalData &additional_data = AdditionalData());
214
246 template <bool level_dof_access>
247 void
248 reinit(
250 const unsigned int q_index = numbers::invalid_unsigned_int,
251 const unsigned int mapping_index = numbers::invalid_unsigned_int);
252
265 void
267 const unsigned int q_index = numbers::invalid_unsigned_int,
268 const unsigned int mapping_index = numbers::invalid_unsigned_int,
269 const unsigned int fe_index = numbers::invalid_unsigned_int);
270
279 const std::optional<::FEValues<dim>> &
280 get_inside_fe_values() const;
281
290 const std::optional<::FEValues<dim>> &
291 get_outside_fe_values() const;
292
300 const std::optional<FEImmersedSurfaceValues<dim>> &
301 get_surface_fe_values() const;
302
303 private:
307 template <typename CellIteratorType>
308 void
309 reinit_internal(const CellIteratorType &cell,
310 const unsigned int q_index,
311 const unsigned int mapping_index,
312 const unsigned int fe_index);
313
319 void
320 initialize(const hp::QCollection<dim> &q_collection);
321
326
331
337
342
346 unsigned int active_fe_index;
347
352
357
371 std::deque<std::optional<::FEValues<dim>>>
373
387 std::deque<std::optional<::FEValues<dim>>>
389
397 std::optional<::FEValues<dim>> fe_values_inside;
398
406 std::optional<::FEValues<dim>> fe_values_outside;
407
415 std::optional<NonMatching::FEImmersedSurfaceValues<dim>> fe_values_surface;
416
421 };
422
423
424
487 template <int dim>
489 {
490 public:
493
515 template <typename VectorType>
517 const Quadrature<1> &quadrature,
520 const DoFHandler<dim> &dof_handler,
521 const VectorType &level_set,
522 const AdditionalData &additional_data = AdditionalData());
523
550 template <typename VectorType>
553 const hp::QCollection<dim - 1> &q_collection,
554 const hp::QCollection<1> &q_collection_1d,
557 const DoFHandler<dim> &dof_handler,
558 const VectorType &level_set,
559 const AdditionalData &additional_data = AdditionalData());
560
574 template <typename CellIteratorType, typename CellNeighborIteratorType>
575 void
576 reinit(const CellIteratorType &cell,
577 const unsigned int face_no,
578 const unsigned int sub_face_no,
579 const CellNeighborIteratorType &cell_neighbor,
580 const unsigned int face_no_neighbor,
581 const unsigned int sub_face_no_neighbor,
582 const unsigned int q_index = numbers::invalid_unsigned_int,
583 const unsigned int mapping_index = numbers::invalid_unsigned_int,
584 const unsigned int fe_index = numbers::invalid_unsigned_int);
585
586
594 template <typename CellIteratorType>
595 void
596 reinit(const CellIteratorType &cell,
597 const unsigned int face_no,
598 const unsigned int q_index = numbers::invalid_unsigned_int,
599 const unsigned int mapping_index = numbers::invalid_unsigned_int,
600 const unsigned int fe_index = numbers::invalid_unsigned_int);
601
611 const std::optional<::FEInterfaceValues<dim>> &
612 get_inside_fe_values() const;
613
623 const std::optional<::FEInterfaceValues<dim>> &
624 get_outside_fe_values() const;
625
626 private:
632 void
633 initialize(const hp::QCollection<dim - 1> &q_collection);
634
641 template <typename CellAccessorType>
642 void
644 const unsigned int face_no,
645 const unsigned int q_index,
646 const unsigned int active_fe_index,
647 const std::function<void(::FEInterfaceValues<dim> &,
648 const unsigned int)> &call_reinit);
649
654
659
665
670
675
680
685 std::optional<::FEInterfaceValues<dim>>
687
692 std::optional<::FEInterfaceValues<dim>>
694
702 std::optional<::FEInterfaceValues<dim>> fe_values_inside;
703
711 std::optional<::FEInterfaceValues<dim>> fe_values_outside;
712
717 };
718
719
720#ifndef DOXYGEN
721
722 /*---------------------- Inline functions ---------------------*/
723
724 template <int dim>
725 template <typename CellIteratorType>
726 inline void
727 FEInterfaceValues<dim>::reinit(const CellIteratorType &cell,
728 const unsigned int face_no,
729 const unsigned int q_index,
730 const unsigned int mapping_index,
731 const unsigned int active_fe_index)
732 {
733 // Lambda describing how we should call reinit on a single
734 // ::FEInterfaceValues object.
735 const auto reinit_operation =
736 [&](::FEInterfaceValues<dim> &fe_interface_values,
737 const unsigned int q_index) {
738 fe_interface_values.reinit(
739 cell, face_no, q_index, mapping_index, active_fe_index);
740 };
741
742 do_reinit(cell, face_no, q_index, active_fe_index, reinit_operation);
743 }
744
745
746
747 template <int dim>
748 template <typename CellIteratorType, typename CellNeighborIteratorType>
749 inline void
750 FEInterfaceValues<dim>::reinit(const CellIteratorType &cell,
751 const unsigned int face_no,
752 const unsigned int sub_face_no,
753 const CellNeighborIteratorType &cell_neighbor,
754 const unsigned int face_no_neighbor,
755 const unsigned int sub_face_no_neighbor,
756 const unsigned int q_index,
757 const unsigned int mapping_index,
758 const unsigned int active_fe_index)
759 {
761 Assert(sub_face_no_neighbor == numbers::invalid_unsigned_int,
763
764 // Lambda describing how we should call reinit on a single
765 // ::FEInterfaceValues object.
766 const auto reinit_operation =
767 [&](::FEInterfaceValues<dim> &fe_interface_values,
768 const unsigned int q_index) {
769 fe_interface_values.reinit(cell,
770 face_no,
771 sub_face_no,
772 cell_neighbor,
773 face_no_neighbor,
774 sub_face_no_neighbor,
775 q_index,
776 mapping_index,
777 active_fe_index);
778 };
779
780 do_reinit(cell, face_no, q_index, active_fe_index, reinit_operation);
781 }
782
783#endif
784
785} // namespace NonMatching
787
788#endif
void reinit(const CellIteratorType &cell, const unsigned int face_no, const unsigned int sub_face_no, const CellNeighborIteratorType &cell_neighbor, const unsigned int face_no_neighbor, const unsigned int sub_face_no_neighbor, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int, const unsigned int fe_index_neighbor=numbers::invalid_unsigned_int)
const std::optional<::FEInterfaceValues< dim > > & get_outside_fe_values() const
Definition fe_values.cc:549
const SmartPointer< const hp::FECollection< dim > > fe_collection
Definition fe_values.h:658
void initialize(const hp::QCollection< dim - 1 > &q_collection)
Definition fe_values.cc:397
const SmartPointer< const hp::MappingCollection< dim > > mapping_collection
Definition fe_values.h:653
const hp::QCollection< 1 > q_collection_1D
Definition fe_values.h:664
void reinit(const CellIteratorType &cell, const unsigned int face_no, const unsigned int sub_face_no, const CellNeighborIteratorType &cell_neighbor, const unsigned int face_no_neighbor, const unsigned int sub_face_no_neighbor, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
const RegionUpdateFlags region_update_flags
Definition fe_values.h:674
void reinit(const CellIteratorType &cell, const unsigned int face_no, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
void do_reinit(const TriaIterator< CellAccessorType > &cell, const unsigned int face_no, const unsigned int q_index, const unsigned int active_fe_index, const std::function< void(::FEInterfaceValues< dim > &, const unsigned int)> &call_reinit)
Definition fe_values.cc:434
LocationToLevelSet current_face_location
Definition fe_values.h:669
std::optional<::FEInterfaceValues< dim > > fe_values_inside_full_quadrature
Definition fe_values.h:686
const std::optional<::FEInterfaceValues< dim > > & get_inside_fe_values() const
Definition fe_values.cc:537
const SmartPointer< const MeshClassifier< dim > > mesh_classifier
Definition fe_values.h:679
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())
Definition fe_values.cc:339
std::optional<::FEInterfaceValues< dim > > fe_values_outside_full_quadrature
Definition fe_values.h:693
std::optional<::FEInterfaceValues< dim > > fe_values_inside
Definition fe_values.h:702
DiscreteFaceQuadratureGenerator< dim > face_quadrature_generator
Definition fe_values.h:716
typename FaceQuadratureGenerator< dim >::AdditionalData AdditionalData
Definition fe_values.h:491
std::optional<::FEInterfaceValues< dim > > fe_values_outside
Definition fe_values.h:711
std::optional< NonMatching::FEImmersedSurfaceValues< dim > > fe_values_surface
Definition fe_values.h:415
std::deque< std::optional<::FEValues< dim > > > fe_values_inside_full_quadrature
Definition fe_values.h:372
std::optional<::FEValues< dim > > fe_values_outside
Definition fe_values.h:406
LocationToLevelSet current_cell_location
Definition fe_values.h:341
const SmartPointer< const hp::FECollection< dim > > fe_collection
Definition fe_values.h:330
const std::optional<::FEValues< dim > > & get_outside_fe_values() const
Definition fe_values.cc:318
void initialize(const hp::QCollection< dim > &q_collection)
Definition fe_values.cc:101
const std::optional< FEImmersedSurfaceValues< dim > > & get_surface_fe_values() const
Definition fe_values.cc:330
void reinit_internal(const CellIteratorType &cell, const unsigned int q_index, const unsigned int mapping_index, const unsigned int fe_index)
Definition fe_values.cc:178
const SmartPointer< const MeshClassifier< dim > > mesh_classifier
Definition fe_values.h:356
const hp::QCollection< 1 > q_collection_1D
Definition fe_values.h:336
const RegionUpdateFlags region_update_flags
Definition fe_values.h:351
void reinit(const TriaIterator< DoFCellAccessor< dim, dim, level_dof_access > > &cell, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int)
Definition fe_values.cc:150
std::optional<::FEValues< dim > > fe_values_inside
Definition fe_values.h:397
const SmartPointer< const hp::MappingCollection< dim > > mapping_collection
Definition fe_values.h:325
typename QuadratureGenerator< dim >::AdditionalData AdditionalData
Definition fe_values.h:146
DiscreteQuadratureGenerator< dim > quadrature_generator
Definition fe_values.h:420
std::deque< std::optional<::FEValues< dim > > > fe_values_outside_full_quadrature
Definition fe_values.h:388
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())
Definition fe_values.cc:45
unsigned int active_fe_index
Definition fe_values.h:346
const std::optional<::FEValues< dim > > & get_inside_fe_values() const
Definition fe_values.cc:306
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
UpdateFlags
static const unsigned int invalid_unsigned_int
Definition types.h:220