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) 2004 - 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_hp_fe_values_h
16#define dealii_hp_fe_values_h
17
18#include <deal.II/base/config.h>
19
21
22#include <deal.II/fe/fe.h>
24
25#include <deal.II/grid/tria.h>
27
31
32#include <memory>
33
35
36// Forward declaration
37#ifndef DOXYGEN
38template <int dim, int spacedim>
39class FiniteElement;
40#endif
41
42
43namespace hp
44{
66 template <int dim, int q_dim, typename FEValuesType>
68 {
69 public:
80
90 const std::vector<QCollection<q_dim>> &q_collection,
92
102
110 const std::vector<QCollection<q_dim>> &q_collection,
112
117
123 operator=(const FEValuesBase &) = delete;
124
134 void
135 precalculate_fe_values(const std::vector<unsigned int> &fe_indices,
136 const std::vector<unsigned int> &mapping_indices,
137 const std::vector<unsigned int> &q_indices);
138
151 void
153
160
166
170 const QCollection<q_dim> &
172
178
185 const FEValuesType &
187
188 protected:
196 FEValuesType &
197 select_fe_values(const unsigned int fe_index,
198 const unsigned int mapping_index,
199 const unsigned int q_index);
200
201 protected:
208
212 const SmartPointer<
216
221
230 const std::vector<QCollection<q_dim>> q_collections;
231
232 private:
245
251
256 };
257
258} // namespace hp
259
260
261namespace hp
262{
312 template <int dim, int spacedim = dim>
314 : public hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>
315 {
316 public:
317 static constexpr unsigned int dimension = dim;
318
319 static constexpr unsigned int space_dimension = spacedim;
320
328
329
338
339
385 template <bool lda>
386 void
388 const unsigned int q_index = numbers::invalid_unsigned_int,
389 const unsigned int mapping_index = numbers::invalid_unsigned_int,
390 const unsigned int fe_index = numbers::invalid_unsigned_int);
391
404 void
406 const unsigned int q_index = numbers::invalid_unsigned_int,
407 const unsigned int mapping_index = numbers::invalid_unsigned_int,
408 const unsigned int fe_index = numbers::invalid_unsigned_int);
409 };
410
411
412
436 template <int dim, int spacedim = dim>
438 : public hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>
439 {
440 public:
448
460 const std::vector<hp::QCollection<dim - 1>> &q_collections,
462
463
472
483 const std::vector<hp::QCollection<dim - 1>> &q_collections,
485
534 template <bool lda>
535 void
537 const unsigned int face_no,
538 const unsigned int q_index = numbers::invalid_unsigned_int,
539 const unsigned int mapping_index = numbers::invalid_unsigned_int,
540 const unsigned int fe_index = numbers::invalid_unsigned_int);
541
547 template <bool lda>
548 void
551 const unsigned int q_index = numbers::invalid_unsigned_int,
552 const unsigned int mapping_index = numbers::invalid_unsigned_int,
553 const unsigned int fe_index = numbers::invalid_unsigned_int);
554
567 void
569 const unsigned int face_no,
570 const unsigned int q_index = numbers::invalid_unsigned_int,
571 const unsigned int mapping_index = numbers::invalid_unsigned_int,
572 const unsigned int fe_index = numbers::invalid_unsigned_int);
573
579 void
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
587
588
595 template <int dim, int spacedim = dim>
597 : public hp::
598 FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>
599 {
600 public:
605 const hp::MappingCollection<dim, spacedim> &mapping_collection,
606 const hp::FECollection<dim, spacedim> &fe_collection,
607 const hp::QCollection<dim - 1> &q_collection,
608 const UpdateFlags update_flags);
609
610
617 const hp::QCollection<dim - 1> &q_collection,
618 const UpdateFlags update_flags);
619
655 template <bool lda>
656 void
658 const unsigned int face_no,
659 const unsigned int subface_no,
660 const unsigned int q_index = numbers::invalid_unsigned_int,
661 const unsigned int mapping_index = numbers::invalid_unsigned_int,
662 const unsigned int fe_index = numbers::invalid_unsigned_int);
663
676 void
678 const unsigned int face_no,
679 const unsigned int subface_no,
680 const unsigned int q_index = numbers::invalid_unsigned_int,
681 const unsigned int mapping_index = numbers::invalid_unsigned_int,
682 const unsigned int fe_index = numbers::invalid_unsigned_int);
683 };
684
685} // namespace hp
686
687
688// -------------- inline and template functions --------------
689
690namespace hp
691{
692 template <int dim, int q_dim, typename FEValuesType>
693 inline const FEValuesType &
695 {
696 Assert(
697 present_fe_values_index != TableIndices<3>(numbers::invalid_unsigned_int,
701 "The indices of the present FEValues object are all invalid. The reason"
702 " is likely that you have forgotten to call the reinit() function."));
703
704 return *fe_values_table(present_fe_values_index);
705 }
706
707
708
709 template <int dim, int q_dim, typename FEValuesType>
712 {
713 return *fe_collection;
714 }
715
716
717
718 template <int dim, int q_dim, typename FEValuesType>
721 {
722 return *mapping_collection;
723 }
724
725
726
727 template <int dim, int q_dim, typename FEValuesType>
728 inline const QCollection<q_dim> &
730 {
731 return q_collection;
732 }
733
734
735
736 template <int dim, int q_dim, typename FEValuesType>
737 inline UpdateFlags
739 {
740 return update_flags;
741 }
742} // namespace hp
743
745
746#endif
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda > > &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)
Definition fe_values.cc:423
FEFaceValues(const hp::MappingCollection< dim, spacedim > &mapping_collection, const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &q_collection, const UpdateFlags update_flags)
Definition fe_values.cc:370
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda > > &cell, const unsigned int face_no, const unsigned int subface_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)
Definition fe_values.cc:563
FESubfaceValues(const hp::MappingCollection< dim, spacedim > &mapping_collection, const hp::FECollection< dim, spacedim > &fe_collection, const hp::QCollection< dim - 1 > &q_collection, const UpdateFlags update_flags)
Definition fe_values.cc:535
const QCollection< q_dim > & get_quadrature_collection() const
Definition fe_values.h:729
FEValuesBase(const MappingCollection< dim, FEValuesType::space_dimension > &mapping_collection, const FECollection< dim, FEValuesType::space_dimension > &fe_collection, const QCollection< q_dim > &q_collection, const UpdateFlags update_flags)
Definition fe_values.cc:65
FEValuesBase & operator=(const FEValuesBase &)=delete
const FEValuesType & get_present_fe_values() const
Definition fe_values.h:694
const std::vector< QCollection< q_dim > > q_collections
Definition fe_values.h:230
FEValuesType & select_fe_values(const unsigned int fe_index, const unsigned int mapping_index, const unsigned int q_index)
Definition fe_values.cc:175
const QCollection< q_dim > q_collection
Definition fe_values.h:220
const SmartPointer< const FECollection< dim, FEValuesType::space_dimension >, FEValuesBase< dim, q_dim, FEValuesType > > fe_collection
Definition fe_values.h:207
UpdateFlags get_update_flags() const
Definition fe_values.h:738
const MappingCollection< dim, FEValuesType::space_dimension > & get_mapping_collection() const
Definition fe_values.h:720
Table< 3, std::unique_ptr< FEValuesType > > fe_values_table
Definition fe_values.h:244
const FECollection< dim, FEValuesType::space_dimension > & get_fe_collection() const
Definition fe_values.h:711
TableIndices< 3 > present_fe_values_index
Definition fe_values.h:250
void precalculate_fe_values()
Definition fe_values.cc:242
const SmartPointer< const MappingCollection< dim, FEValuesType::space_dimension >, FEValuesBase< dim, q_dim, FEValuesType > > mapping_collection
Definition fe_values.h:215
const UpdateFlags update_flags
Definition fe_values.h:255
static constexpr unsigned int dimension
Definition fe_values.h:317
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda > > &cell, 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)
Definition fe_values.cc:293
static constexpr unsigned int space_dimension
Definition fe_values.h:319
FEValues(const MappingCollection< dim, spacedim > &mapping_collection, const FECollection< dim, spacedim > &fe_collection, const QCollection< dim > &q_collection, const UpdateFlags update_flags)
Definition fe_values.cc:267
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
#define Assert(cond, exc)
static ::ExceptionBase & ExcMessage(std::string arg1)
UpdateFlags
Definition hp.h:117
static const unsigned int invalid_unsigned_int
Definition types.h:220