15#ifndef dealii_fe_poly_tensor_h
16#define dealii_fe_poly_tensor_h
140template <
int dim,
int spacedim = dim>
178 const unsigned int component)
const override;
193 const unsigned int component)
const override;
208 const unsigned int component)
const override;
223 const unsigned int n_q_points,
260 const unsigned int index,
261 const unsigned int face_no,
293 virtual std::unique_ptr<
301 &output_data)
const override
307 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
308 data_ptr = std::make_unique<InternalData>();
312 const unsigned int n_q_points = quadrature.
size();
315 std::vector<Tensor<1, dim>> values(0);
316 std::vector<Tensor<2, dim>> grads(0);
317 std::vector<Tensor<3, dim>> grad_grads(0);
318 std::vector<Tensor<4, dim>> third_derivatives(0);
319 std::vector<Tensor<5, dim>> fourth_derivatives(0);
328 const bool update_transformed_shape_values =
329 std::any_of(this->mapping_kind.begin(),
330 this->mapping_kind.end(),
331 [](
const MappingKind t) { return t != mapping_none; });
333 const bool update_transformed_shape_grads =
334 std::any_of(this->mapping_kind.begin(),
335 this->mapping_kind.end(),
337 return (t == mapping_raviart_thomas || t == mapping_piola ||
338 t == mapping_nedelec || t == mapping_contravariant);
341 const bool update_transformed_shape_hessian_tensors =
342 update_transformed_shape_values;
348 if (update_transformed_shape_values)
349 data.transformed_shape_values.resize(n_q_points);
356 data.transformed_shape_grads.resize(n_q_points);
358 if (update_transformed_shape_grads)
359 data.untransformed_shape_grads.resize(n_q_points);
366 data.transformed_shape_hessians.resize(n_q_points);
367 if (update_transformed_shape_hessian_tensors)
368 data.untransformed_shape_hessian_tensors.resize(n_q_points);
379 for (
unsigned int k = 0; k < n_q_points; ++k)
392 data.shape_values[i][k] = values[i];
399 data.shape_values[i][k] = add_values;
407 data.shape_grads[i][k] = grads[i];
414 data.shape_grads[i][k] = add_grads;
422 data.shape_grad_grads[i][k] = grad_grads[i];
430 data.shape_grad_grads[i][k] = add_grad_grads;
449 &output_data)
const override;
456 const unsigned int face_no,
465 &output_data)
const override;
470 const unsigned int face_no,
471 const unsigned int sub_no,
480 &output_data)
const override;
533 const std::unique_ptr<const TensorPolynomialsBase<dim>>
poly_space;
std::vector< Tensor< 3, spacedim > > transformed_shape_hessians
Table< 2, DerivativeForm< 2, dim, spacedim > > shape_grad_grads
Table< 2, DerivativeForm< 1, dim, spacedim > > shape_grads
std::vector< Tensor< 2, spacedim > > transformed_shape_grads
std::vector< Tensor< 3, dim > > untransformed_shape_hessian_tensors
std::vector< Tensor< 2, dim > > untransformed_shape_grads
Table< 2, Tensor< 1, dim > > shape_values
std::vector< Tensor< 1, spacedim > > transformed_shape_values
std::vector< double > dof_sign_change
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
std::vector< Table< 2, bool > > adjust_quad_dof_sign_for_face_orientation_table
std::vector< Tensor< 2, dim > > cached_grads
bool adjust_quad_dof_sign_for_face_orientation(const unsigned int index, const unsigned int face_no, const types::geometric_orientation combined_orientation) const
FullMatrix< double > inverse_node_matrix
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
FE_PolyTensor(const FE_PolyTensor &fe)
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Point< dim > cached_point
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
void compute_fill(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const typename QProjector< dim >::DataSetDescriptor &offset, const unsigned int n_q_points, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const typename FE_PolyTensor< dim, spacedim >::InternalData &fe_internal, internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
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
Threads::Mutex cache_mutex
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 double shape_value(const unsigned int i, const Point< dim > &p) 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
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
std::vector< MappingKind > mapping_kind
MappingKind get_mapping_kind(const unsigned int i) const
FE_PolyTensor(const TensorPolynomialsBase< dim > &polynomials, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
std::vector< Tensor< 1, dim > > cached_values
bool single_mapping_kind() const
const std::unique_ptr< const TensorPolynomialsBase< dim > > poly_space
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
std::vector< Tensor< 3, dim > > cached_grad_grads
unsigned int n_dofs_per_cell() const
const unsigned int dofs_per_cell
const std::vector< bool > restriction_is_additive_flags
const std::vector< ComponentMask > nonzero_components
Abstract base class for mapping classes.
Class storing the offset index into a Quadrature rule created by project_to_all_faces() or project_to...
const Point< dim > & point(const unsigned int i) const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_gradients
Shape function gradients.
std::vector< index_type > data
unsigned char geometric_orientation