Reference documentation for deal.II version Git ca0b05a7a7 2020-12-01 11:49:21 -0500
\(\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\}}\)
fe_values.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_hp_fe_values_h
17 #define dealii_hp_fe_values_h
18 
19 #include <deal.II/base/config.h>
20 
22 
23 #include <deal.II/fe/fe.h>
24 #include <deal.II/fe/fe_values.h>
25 
26 #include <deal.II/grid/tria.h>
28 
32 
33 #include <map>
34 #include <memory>
35 
37 
38 // Forward declaration
39 #ifndef DOXYGEN
40 template <int dim, int spacedim>
41 class FiniteElement;
42 #endif
43 
44 
45 namespace hp
46 {
68  template <int dim, int q_dim, class FEValuesType>
70  {
71  public:
82 
90  const QCollection<q_dim> & q_collection,
91  const UpdateFlags update_flags);
92 
97 
102  FEValuesBase &
103  operator=(const FEValuesBase &) = delete;
104 
114  void
115  precalculate_fe_values(const std::vector<unsigned int> &fe_indices,
116  const std::vector<unsigned int> &mapping_indices,
117  const std::vector<unsigned int> &q_indices);
118 
131  void
133 
139  get_fe_collection() const;
140 
145  get_mapping_collection() const;
146 
150  const QCollection<q_dim> &
152 
157  get_update_flags() const;
158 
165  const FEValuesType &
166  get_present_fe_values() const;
167 
168  protected:
176  FEValuesType &
177  select_fe_values(const unsigned int fe_index,
178  const unsigned int mapping_index,
179  const unsigned int q_index);
180 
181  protected:
188 
192  const SmartPointer<
196 
201 
202  private:
215 
221 
226  };
227 
228 } // namespace hp
229 
230 
231 namespace hp
232 {
282  template <int dim, int spacedim = dim>
283  class FEValues
284  : public hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>
285  {
286  public:
287  static const unsigned int dimension = dim;
288 
289  static const unsigned int space_dimension = spacedim;
290 
297  const UpdateFlags update_flags);
298 
299 
307  const UpdateFlags update_flags);
308 
309 
355  template <bool lda>
356  void
358  const unsigned int q_index = numbers::invalid_unsigned_int,
359  const unsigned int mapping_index = numbers::invalid_unsigned_int,
360  const unsigned int fe_index = numbers::invalid_unsigned_int);
361 
374  void
376  const unsigned int q_index = numbers::invalid_unsigned_int,
377  const unsigned int mapping_index = numbers::invalid_unsigned_int,
378  const unsigned int fe_index = numbers::invalid_unsigned_int);
379  };
380 
381 
382 
406  template <int dim, int spacedim = dim>
408  : public hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>
409  {
410  public:
417  const UpdateFlags update_flags);
418 
419 
425  FEFaceValues(const hp::FECollection<dim, spacedim> &fe_collection,
426  const hp::QCollection<dim - 1> & q_collection,
427  const UpdateFlags update_flags);
428 
474  template <bool lda>
475  void
477  const unsigned int face_no,
478  const unsigned int q_index = numbers::invalid_unsigned_int,
479  const unsigned int mapping_index = numbers::invalid_unsigned_int,
480  const unsigned int fe_index = numbers::invalid_unsigned_int);
481 
487  template <bool lda>
488  void
490  const typename Triangulation<dim, spacedim>::face_iterator &face,
491  const unsigned int q_index = numbers::invalid_unsigned_int,
492  const unsigned int mapping_index = numbers::invalid_unsigned_int,
493  const unsigned int fe_index = numbers::invalid_unsigned_int);
494 
507  void
509  const unsigned int face_no,
510  const unsigned int q_index = numbers::invalid_unsigned_int,
511  const unsigned int mapping_index = numbers::invalid_unsigned_int,
512  const unsigned int fe_index = numbers::invalid_unsigned_int);
513 
519  void
521  const typename Triangulation<dim, spacedim>::face_iterator &face,
522  const unsigned int q_index = numbers::invalid_unsigned_int,
523  const unsigned int mapping_index = numbers::invalid_unsigned_int,
524  const unsigned int fe_index = numbers::invalid_unsigned_int);
525  };
526 
527 
528 
535  template <int dim, int spacedim = dim>
537  : public hp::
538  FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>
539  {
540  public:
548  const UpdateFlags update_flags);
549 
550 
557  const hp::QCollection<dim - 1> & q_collection,
558  const UpdateFlags update_flags);
559 
595  template <bool lda>
596  void
598  const unsigned int face_no,
599  const unsigned int subface_no,
600  const unsigned int q_index = numbers::invalid_unsigned_int,
601  const unsigned int mapping_index = numbers::invalid_unsigned_int,
602  const unsigned int fe_index = numbers::invalid_unsigned_int);
603 
616  void
618  const unsigned int face_no,
619  const unsigned int subface_no,
620  const unsigned int q_index = numbers::invalid_unsigned_int,
621  const unsigned int mapping_index = numbers::invalid_unsigned_int,
622  const unsigned int fe_index = numbers::invalid_unsigned_int);
623  };
624 
625 } // namespace hp
626 
627 
628 // -------------- inline and template functions --------------
629 
630 namespace hp
631 {
632  template <int dim, int q_dim, class FEValuesType>
633  inline const FEValuesType &
635  {
637  }
638 
639 
640 
641  template <int dim, int q_dim, class FEValuesType>
644  {
645  return *fe_collection;
646  }
647 
648 
649 
650  template <int dim, int q_dim, class FEValuesType>
653  {
654  return *mapping_collection;
655  }
656 
657 
658 
659  template <int dim, int q_dim, class FEValuesType>
660  inline const QCollection<q_dim> &
662  {
663  return q_collection;
664  }
665 
666 
667 
668  template <int dim, int q_dim, class FEValuesType>
669  inline UpdateFlags
671  {
672  return update_flags;
673  }
674 } // namespace hp
675 
677 
678 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
const QCollection< q_dim > q_collection
Definition: fe_values.h:200
const SmartPointer< const MappingCollection< dim, FEValuesType::space_dimension >, FEValuesBase< dim, q_dim, FEValuesType > > mapping_collection
Definition: fe_values.h:195
TableIndices< 3 > present_fe_values_index
Definition: fe_values.h:220
UpdateFlags get_update_flags() const
Definition: fe_values.h:670
const UpdateFlags update_flags
Definition: fe_values.h:225
const MappingCollection< dim, FEValuesType::space_dimension > & get_mapping_collection() const
Definition: fe_values.h:652
UpdateFlags
void precalculate_fe_values()
Definition: fe_values.cc:177
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:372
FEValuesBase & operator=(const FEValuesBase &)=delete
Definition: hp.h:117
const SmartPointer< const FECollection< dim, FEValuesType::space_dimension >, FEValuesBase< dim, q_dim, FEValuesType > > fe_collection
Definition: fe_values.h:187
FEValuesType & select_fe_values(const unsigned int fe_index, const unsigned int mapping_index, const unsigned int q_index)
Definition: fe_values.cc:107
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:371
const FECollection< dim, FEValuesType::space_dimension > & get_fe_collection() const
Definition: fe_values.h:643
const QCollection< q_dim > & get_quadrature_collection() const
Definition: fe_values.h:661
Table< 3, std::unique_ptr< FEValuesType > > fe_values_table
Definition: fe_values.h:214
Definition: table.h:687
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:31
const FEValuesType & get_present_fe_values() const
Definition: fe_values.h:634