16#ifndef dealii_fe_poly_face_h
17#define dealii_fe_poly_face_h
58template <
class PolynomialType,
59 int dim = PolynomialType::dimension + 1,
90 virtual std::unique_ptr<
98 &output_data)
const override
104 return std::make_unique<InternalData>();
109 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
116 &output_data)
const override
124 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
125 data_ptr = std::make_unique<InternalData>();
129 const unsigned int n_q_points = quadrature[0].
size();
132 std::vector<double> values(0);
133 std::vector<
Tensor<1, dim - 1>> grads(0);
134 std::vector<
Tensor<2, dim - 1>> grad_grads(0);
135 std::vector<
Tensor<3, dim - 1>> empty_vector_of_3rd_order_tensors;
136 std::vector<
Tensor<4, dim - 1>> empty_vector_of_4th_order_tensors;
145 std::vector<double>(n_q_points));
146 for (
unsigned int i = 0; i < n_q_points; ++i)
152 empty_vector_of_3rd_order_tensors,
153 empty_vector_of_4th_order_tensors);
155 for (
unsigned int k = 0; k <
poly_space.n(); ++k)
156 data.shape_values[k][i] = values[k];
170 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
177 &output_data)
const override
183 ReferenceCells::get_hypercube<dim - 1>(), quadrature)),
199 &output_data)
const override;
206 const unsigned int face_no,
215 &output_data)
const override;
220 const unsigned int face_no,
221 const unsigned int sub_no,
230 &output_data)
const override;
std::vector< std::vector< double > > shape_values
PolynomialType poly_space
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
std::unique_ptr< typename FiniteElement< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual void fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
FE_PolyFace(const PolynomialType &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FiniteElement< dim, spacedim >::InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
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
unsigned int get_degree() const
const std::vector< bool > restriction_is_additive_flags
Abstract base class for mapping classes.
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_gradients
Shape function gradients.