Reference documentation for deal.II version Git 5b8b897cb2 2021-04-22 22:24:19 -0400
\(\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  &mapping_collection,
92  const std::vector<QCollection<q_dim>> & q_collection,
93  const UpdateFlags update_flags);
94 
100  FEValuesBase(
102  const QCollection<q_dim> & q_collection,
103  const UpdateFlags update_flags);
104 
110  FEValuesBase(
112  const std::vector<QCollection<q_dim>> & q_collection,
113  const UpdateFlags update_flags);
114 
119 
124  FEValuesBase &
125  operator=(const FEValuesBase &) = delete;
126 
136  void
137  precalculate_fe_values(const std::vector<unsigned int> &fe_indices,
138  const std::vector<unsigned int> &mapping_indices,
139  const std::vector<unsigned int> &q_indices);
140 
153  void
155 
161  get_fe_collection() const;
162 
167  get_mapping_collection() const;
168 
172  const QCollection<q_dim> &
174 
179  get_update_flags() const;
180 
187  const FEValuesType &
188  get_present_fe_values() const;
189 
190  protected:
198  FEValuesType &
199  select_fe_values(const unsigned int fe_index,
200  const unsigned int mapping_index,
201  const unsigned int q_index);
202 
203  protected:
210 
214  const SmartPointer<
218 
223 
232  const std::vector<QCollection<q_dim>> q_collections;
233 
234  private:
247 
253 
258  };
259 
260 } // namespace hp
261 
262 
263 namespace hp
264 {
314  template <int dim, int spacedim = dim>
315  class FEValues
316  : public hp::FEValuesBase<dim, dim, ::FEValues<dim, spacedim>>
317  {
318  public:
319  static const unsigned int dimension = dim;
320 
321  static const unsigned int space_dimension = spacedim;
322 
329  const UpdateFlags update_flags);
330 
331 
339  const UpdateFlags update_flags);
340 
341 
387  template <bool lda>
388  void
390  const unsigned int q_index = numbers::invalid_unsigned_int,
391  const unsigned int mapping_index = numbers::invalid_unsigned_int,
392  const unsigned int fe_index = numbers::invalid_unsigned_int);
393 
406  void
408  const unsigned int q_index = numbers::invalid_unsigned_int,
409  const unsigned int mapping_index = numbers::invalid_unsigned_int,
410  const unsigned int fe_index = numbers::invalid_unsigned_int);
411  };
412 
413 
414 
438  template <int dim, int spacedim = dim>
440  : public hp::FEValuesBase<dim, dim - 1, ::FEFaceValues<dim, spacedim>>
441  {
442  public:
449  const UpdateFlags update_flags);
450 
460  FEFaceValues(const hp::MappingCollection<dim, spacedim> &mapping_collection,
461  const hp::FECollection<dim, spacedim> & fe_collection,
462  const std::vector<hp::QCollection<dim - 1>> &q_collections,
463  const UpdateFlags update_flags);
464 
465 
471  FEFaceValues(const hp::FECollection<dim, spacedim> &fe_collection,
472  const hp::QCollection<dim - 1> & q_collection,
473  const UpdateFlags update_flags);
474 
484  FEFaceValues(const hp::FECollection<dim, spacedim> & fe_collection,
485  const std::vector<hp::QCollection<dim - 1>> &q_collections,
486  const UpdateFlags update_flags);
487 
533  template <bool lda>
534  void
536  const unsigned int face_no,
537  const unsigned int q_index = numbers::invalid_unsigned_int,
538  const unsigned int mapping_index = numbers::invalid_unsigned_int,
539  const unsigned int fe_index = numbers::invalid_unsigned_int);
540 
546  template <bool lda>
547  void
549  const typename Triangulation<dim, spacedim>::face_iterator &face,
550  const unsigned int q_index = numbers::invalid_unsigned_int,
551  const unsigned int mapping_index = numbers::invalid_unsigned_int,
552  const unsigned int fe_index = numbers::invalid_unsigned_int);
553 
566  void
568  const unsigned int face_no,
569  const unsigned int q_index = numbers::invalid_unsigned_int,
570  const unsigned int mapping_index = numbers::invalid_unsigned_int,
571  const unsigned int fe_index = numbers::invalid_unsigned_int);
572 
578  void
580  const typename Triangulation<dim, spacedim>::face_iterator &face,
581  const unsigned int q_index = numbers::invalid_unsigned_int,
582  const unsigned int mapping_index = numbers::invalid_unsigned_int,
583  const unsigned int fe_index = numbers::invalid_unsigned_int);
584  };
585 
586 
587 
594  template <int dim, int spacedim = dim>
596  : public hp::
597  FEValuesBase<dim, dim - 1, ::FESubfaceValues<dim, spacedim>>
598  {
599  public:
607  const UpdateFlags update_flags);
608 
609 
616  const hp::QCollection<dim - 1> & q_collection,
617  const UpdateFlags update_flags);
618 
654  template <bool lda>
655  void
657  const unsigned int face_no,
658  const unsigned int subface_no,
659  const unsigned int q_index = numbers::invalid_unsigned_int,
660  const unsigned int mapping_index = numbers::invalid_unsigned_int,
661  const unsigned int fe_index = numbers::invalid_unsigned_int);
662 
675  void
677  const unsigned int face_no,
678  const unsigned int subface_no,
679  const unsigned int q_index = numbers::invalid_unsigned_int,
680  const unsigned int mapping_index = numbers::invalid_unsigned_int,
681  const unsigned int fe_index = numbers::invalid_unsigned_int);
682  };
683 
684 } // namespace hp
685 
686 
687 // -------------- inline and template functions --------------
688 
689 namespace hp
690 {
691  template <int dim, int q_dim, class FEValuesType>
692  inline const FEValuesType &
694  {
696  }
697 
698 
699 
700  template <int dim, int q_dim, class FEValuesType>
703  {
704  return *fe_collection;
705  }
706 
707 
708 
709  template <int dim, int q_dim, class FEValuesType>
712  {
713  return *mapping_collection;
714  }
715 
716 
717 
718  template <int dim, int q_dim, class FEValuesType>
719  inline const QCollection<q_dim> &
721  {
722  return q_collection;
723  }
724 
725 
726 
727  template <int dim, int q_dim, class FEValuesType>
728  inline UpdateFlags
730  {
731  return update_flags;
732  }
733 } // namespace hp
734 
736 
737 #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:222
const std::vector< QCollection< q_dim > > q_collections
Definition: fe_values.h:232
const SmartPointer< const MappingCollection< dim, FEValuesType::space_dimension >, FEValuesBase< dim, q_dim, FEValuesType > > mapping_collection
Definition: fe_values.h:217
TableIndices< 3 > present_fe_values_index
Definition: fe_values.h:252
UpdateFlags get_update_flags() const
Definition: fe_values.h:729
const UpdateFlags update_flags
Definition: fe_values.h:257
const MappingCollection< dim, FEValuesType::space_dimension > & get_mapping_collection() const
Definition: fe_values.h:711
UpdateFlags
void precalculate_fe_values()
Definition: fe_values.cc:242
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
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:209
FEValuesType & select_fe_values(const unsigned int fe_index, const unsigned int mapping_index, const unsigned int q_index)
Definition: fe_values.cc:175
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
const FECollection< dim, FEValuesType::space_dimension > & get_fe_collection() const
Definition: fe_values.h:702
const QCollection< q_dim > & get_quadrature_collection() const
Definition: fe_values.h:720
Table< 3, std::unique_ptr< FEValuesType > > fe_values_table
Definition: fe_values.h:246
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:66
const FEValuesType & get_present_fe_values() const
Definition: fe_values.h:693