18 #include <deal.II/base/qprojector.h>
19 #include <deal.II/base/polynomial_space.h>
20 #include <deal.II/base/tensor_product_polynomials.h>
21 #include <deal.II/fe/fe_values.h>
22 #include <deal.II/fe/fe_poly.h>
27 template <
class POLY,
int dim,
int spacedim>
30 const std::vector<bool> &restriction_is_additive_flags,
31 const std::vector<ComponentMask> &nonzero_components):
33 restriction_is_additive_flags,
35 poly_space(poly_space)
41 template <
class POLY,
int dim,
int spacedim>
49 template <
class POLY,
int dim,
int spacedim>
55 return poly_space.compute_value(i, p);
59 template <
class POLY,
int dim,
int spacedim>
63 const unsigned int component)
const
67 return poly_space.compute_value(i, p);
72 template <
class POLY,
int dim,
int spacedim>
78 return poly_space.compute_grad(i, p);
83 template <
class POLY,
int dim,
int spacedim>
87 const unsigned int component)
const
91 return poly_space.compute_grad(i, p);
96 template <
class POLY,
int dim,
int spacedim>
102 return poly_space.compute_grad_grad(i, p);
107 template <
class POLY,
int dim,
int spacedim>
111 const unsigned int component)
const
115 return poly_space.compute_grad_grad(i, p);
127 template <
class POLY,
int dim,
int spacedim>
141 template <
class POLY,
int dim,
int spacedim>
163 template <
class POLY,
int dim,
int spacedim>
182 const unsigned int n_q_points = quadrature.
size();
185 std::vector<double> values(0);
186 std::vector<Tensor<1,dim> > grads(0);
187 std::vector<Tensor<2,dim> > grad_grads(0);
194 values.resize (this->dofs_per_cell);
196 std::vector<double> (n_q_points));
201 grads.resize (this->dofs_per_cell);
220 if (flags & (update_values | update_gradients))
221 for (
unsigned int i=0; i<n_q_points; ++i)
223 poly_space.compute(quadrature.
point(i),
224 values, grads, grad_grads);
226 if (flags & update_values)
227 for (
unsigned int k=0; k<this->dofs_per_cell; ++k)
230 if (flags & update_gradients)
231 for (
unsigned int k=0; k<this->dofs_per_cell; ++k)
245 template <
class POLY,
int dim,
int spacedim>
254 CellSimilarity::Similarity &cell_similarity)
const
265 for (
unsigned int k=0; k<this->dofs_per_cell; ++k)
268 for (
unsigned int i=0; i<quadrature.
size(); ++i)
271 if (flags &
update_gradients && cell_similarity != CellSimilarity::translation)
276 if (flags &
update_hessians && cell_similarity != CellSimilarity::translation)
278 mapping_data, fe_data, data);
283 template <
class POLY,
int dim,
int spacedim>
288 const unsigned int face,
307 cell->face_orientation(face),
308 cell->face_flip(face),
309 cell->face_rotation(face),
314 for (
unsigned int k=0; k<this->dofs_per_cell; ++k)
317 for (
unsigned int i=0; i<quadrature.
size(); ++i)
327 this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
398 template <
class POLY,
int dim,
int spacedim>
402 const unsigned int face,
403 const unsigned int subface,
422 cell->face_orientation(face),
423 cell->face_flip(face),
424 cell->face_rotation(face),
426 cell->subface_case(face));
430 for (
unsigned int k=0; k<this->dofs_per_cell; ++k)
433 for (
unsigned int i=0; i<quadrature.
size(); ++i)
443 this->compute_2nd (mapping, cell, offset, mapping_data, fe_data, data);
450 template <
class POLY>
452 std::vector<unsigned int>
453 get_poly_space_numbering (
const POLY &)
456 return std::vector<unsigned int>();
459 template <
class POLY>
461 std::vector<unsigned int>
462 get_poly_space_numbering_inverse (
const POLY &)
465 return std::vector<unsigned int>();
468 template <
int dim,
typename POLY>
470 std::vector<unsigned int>
476 template <
int dim,
typename POLY>
478 std::vector<unsigned int>
487 template <
class POLY,
int dim,
int spacedim>
488 std::vector<unsigned int>
491 return internal::get_poly_space_numbering (poly_space);
497 template <
class POLY,
int dim,
int spacedim>
498 std::vector<unsigned int>
501 return internal::get_poly_space_numbering_inverse (poly_space);
506 DEAL_II_NAMESPACE_CLOSE
Transformed quadrature weights.
#define AssertDimension(dim1, dim2)
void initialize_2nd(const FiniteElement< dim, spacedim > *element, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature)
virtual void fill_fe_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Quadrature< dim > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data, CellSimilarity::Similarity &cell_similarity) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual Mapping< dim, spacedim >::InternalDataBase * get_data(const UpdateFlags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature) const
virtual void fill_fe_subface_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
std::vector< unsigned int > get_poly_space_numbering() const
std::vector< std::vector< Tensor< 1, dim > > > shape_gradients
unsigned int get_degree() const
virtual void transform(const VectorSlice< const std::vector< Tensor< 1, dim > > > input, VectorSlice< std::vector< Tensor< 1, spacedim > > > output, const InternalDataBase &internal, const MappingType type) const =0
#define Assert(cond, exc)
UpdateFlags current_update_flags() const
virtual UpdateFlags update_once(const UpdateFlags flags) const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Second derivatives of shape functions.
std::vector< std::vector< double > > shape_values
GradientVector shape_gradients
unsigned int size() const
const std::vector< unsigned int > & get_numbering_inverse() const
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const std::vector< unsigned int > & get_numbering() const
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual void fill_fe_face_values(const Mapping< dim, spacedim > &mapping, const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, typename Mapping< dim, spacedim >::InternalDataBase &fe_internal, FEValuesData< dim, spacedim > &data) const
const VectorSlice< const VECTOR > make_slice(VECTOR &v)
std::vector< unsigned int > get_poly_space_numbering_inverse() const
Shape function gradients.
static DataSetDescriptor face(const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
Covariant mapping (see Mapping::transform() for details)
FE_Poly(const POLY &poly_space, const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
::ExceptionBase & ExcNotImplemented()
static DataSetDescriptor subface(const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points, const internal::SubfaceCase< dim > ref_case=internal::SubfaceCase< dim >::case_isotropic)
virtual UpdateFlags update_each(const UpdateFlags flags) const
::ExceptionBase & ExcInternalError()
const Point< dim > & point(const unsigned int i) const
Covariant transformation.