Reference documentation for deal.II version GIT relicensing-233-g802318d791 2024-03-28 20:20:02+00:00
\(\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_nedelec_sz.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2018 - 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_fe_nedelec_sz_h
16#define dealii_fe_nedelec_sz_h
17
18#include <deal.II/base/config.h>
19
23
24#include <deal.II/fe/fe.h>
26#include <deal.II/fe/mapping.h>
27
29
82template <int dim, int spacedim = dim>
83class FE_NedelecSZ : public FiniteElement<dim, dim>
84{
85public:
86 static_assert(dim == spacedim,
87 "FE_NedelecSZ is only implemented for dim==spacedim!");
88
103 FE_NedelecSZ(const unsigned int order);
104
105 virtual UpdateFlags
106 requires_update_flags(const UpdateFlags update_flags) const override;
107
108 virtual std::string
109 get_name() const override;
110
111 virtual std::unique_ptr<FiniteElement<dim, dim>>
112 clone() const override;
113
118 virtual double
119 shape_value(const unsigned int i, const Point<dim> &p) const override;
120
124 virtual double
125 shape_value_component(const unsigned int i,
126 const Point<dim> &p,
127 const unsigned int component) const override;
128
133 virtual Tensor<1, dim>
134 shape_grad(const unsigned int i, const Point<dim> &p) const override;
135
139 virtual Tensor<1, dim>
140 shape_grad_component(const unsigned int i,
141 const Point<dim> &p,
142 const unsigned int component) const override;
143
148 virtual Tensor<2, dim>
149 shape_grad_grad(const unsigned int i, const Point<dim> &p) const override;
150
154 virtual Tensor<2, dim>
155 shape_grad_grad_component(const unsigned int i,
156 const Point<dim> &p,
157 const unsigned int component) const override;
158
159protected:
165
166 virtual std::unique_ptr<
167 typename ::FiniteElement<dim, spacedim>::InternalDataBase>
168 get_data(
169 const UpdateFlags update_flags,
170 const Mapping<dim, spacedim> &mapping,
171 const Quadrature<dim> &quadrature,
173 spacedim>
174 &output_data) const override;
175
181 virtual void
183 const typename Triangulation<dim, dim>::cell_iterator &cell,
184 const CellSimilarity::Similarity cell_similarity,
185 const Quadrature<dim> &quadrature,
186 const Mapping<dim, dim> &mapping,
187 const typename Mapping<dim, dim>::InternalDataBase &mapping_internal,
189 &mapping_data,
190 const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
192 &data) const override;
193
194 using FiniteElement<dim, spacedim>::fill_fe_face_values;
195
201 virtual void
203 const typename Triangulation<dim, dim>::cell_iterator &cell,
204 const unsigned int face_no,
205 const hp::QCollection<dim - 1> &quadrature,
206 const Mapping<dim, dim> &mapping,
207 const typename Mapping<dim, dim>::InternalDataBase &mapping_internal,
209 &mapping_data,
210 const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
212 &data) const override;
213
217 virtual void
219 const typename Triangulation<dim, dim>::cell_iterator &cell,
220 const unsigned int face_no,
221 const unsigned int sub_no,
222 const Quadrature<dim - 1> &quadrature,
223 const Mapping<dim, dim> &mapping,
224 const typename Mapping<dim, dim>::InternalDataBase &mapping_internal,
226 &mapping_data,
227 const typename FiniteElement<dim, dim>::InternalDataBase &fedata,
229 &data) const override;
230
251 {
252 public:
259 mutable std::vector<std::vector<Tensor<1, dim>>> shape_values;
260
268 mutable std::vector<std::vector<DerivativeForm<1, dim, dim>>> shape_grads;
269
277 mutable std::vector<std::vector<DerivativeForm<2, dim, dim>>>
279
295 std::vector<std::vector<std::vector<double>>> sigma_imj_values;
296
312 std::vector<std::vector<std::vector<double>>> sigma_imj_grads;
313
325 std::vector<std::vector<double>> edge_sigma_values;
326
338 std::vector<std::vector<double>> edge_sigma_grads;
339
355 std::vector<std::vector<double>> edge_lambda_values;
356
366 std::vector<std::vector<double>> edge_lambda_grads_2d;
367
377 std::vector<std::vector<std::vector<double>>> edge_lambda_grads_3d;
378
389 std::vector<std::vector<std::vector<double>>> edge_lambda_gradgrads_3d;
390
407 std::vector<std::vector<double>> face_lambda_values;
408
417 std::vector<std::vector<double>> face_lambda_grads;
418 };
419
420private:
429 static std::vector<unsigned int>
430 get_dpo_vector(const unsigned int degree);
431
435 std::vector<Polynomials::Polynomial<double>> IntegratedLegendrePolynomials;
436
441 void
442 create_polynomials(const unsigned int degree);
443
447 unsigned int
448 compute_num_dofs(const unsigned int degree) const;
449
454 void
456 const Quadrature<dim> &quadrature,
457 const InternalData &fedata) const;
458
463 void
465 const Quadrature<dim> &quadrature,
466 const InternalData &fedata) const;
467};
468
469
470
474
475#endif
std::vector< std::vector< double > > edge_lambda_grads_2d
std::vector< std::vector< Tensor< 1, dim > > > shape_values
std::vector< std::vector< std::vector< double > > > sigma_imj_values
std::vector< std::vector< double > > edge_sigma_grads
std::vector< std::vector< double > > edge_sigma_values
std::vector< std::vector< DerivativeForm< 2, dim, dim > > > shape_hessians
std::vector< std::vector< std::vector< double > > > sigma_imj_grads
std::vector< std::vector< double > > face_lambda_values
std::vector< std::vector< DerivativeForm< 1, dim, dim > > > shape_grads
std::vector< std::vector< std::vector< double > > > edge_lambda_grads_3d
std::vector< std::vector< double > > edge_lambda_values
std::vector< std::vector< std::vector< double > > > edge_lambda_gradgrads_3d
std::vector< std::vector< double > > face_lambda_grads
unsigned int compute_num_dofs(const unsigned int degree) const
virtual void fill_fe_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
virtual void fill_fe_subface_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
void create_polynomials(const unsigned int degree)
virtual std::unique_ptr< typename ::FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual void fill_fe_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
virtual std::string get_name() const override
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
std::vector< Polynomials::Polynomial< double > > IntegratedLegendrePolynomials
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
void fill_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
void fill_edge_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
MappingKind mapping_kind
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
const unsigned int degree
Definition fe_data.h:452
friend class InternalDataBase
Definition fe.h:3061
Abstract base class for mapping classes.
Definition mapping.h:316
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
UpdateFlags
MappingKind
Definition mapping.h:77