36 template <
int dim,
int spacedim>
38 const unsigned int degree)
52 , polynomial_space(
Polynomials::Legendre::generate_complete_basis(degree))
56 ref_case < RefinementCase<dim>::isotropic_refinement + 1;
65 const unsigned int nc =
67 for (
unsigned int i = 0; i < nc; ++i)
69 this->
prolongation[ref_case - 1][i].reinit(n_dofs, n_dofs);
72 for (
unsigned int j = 0; j < n_dofs; ++j)
108 template <
int dim,
int spacedim>
119 std::ostringstream namebuf;
121 <<
">(" << this->
degree <<
")";
123 return namebuf.str();
128 template <
int dim,
int spacedim>
129 std::unique_ptr<FiniteElement<dim, spacedim>>
132 return std::make_unique<FE_DGPNonparametric<dim, spacedim>>(*this);
137 template <
int dim,
int spacedim>
152 template <
int dim,
int spacedim>
155 const unsigned int i,
157 const unsigned int component)
const 171 template <
int dim,
int spacedim>
185 template <
int dim,
int spacedim>
188 const unsigned int i,
190 const unsigned int component)
const 204 template <
int dim,
int spacedim>
219 template <
int dim,
int spacedim>
222 const unsigned int i,
224 const unsigned int component)
const 242 template <
int dim,
int spacedim>
243 std::vector<unsigned int>
246 std::vector<unsigned int> dpo(dim + 1, static_cast<unsigned int>(0));
248 for (
unsigned int i = 1; i < dim; ++i)
250 dpo[dim] *= deg + 1 + i;
258 template <
int dim,
int spacedim>
277 template <
int dim,
int spacedim>
278 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
289 std::make_unique<typename FiniteElement<dim, spacedim>::InternalDataBase>();
304 template <
int dim,
int spacedim>
312 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
323 const unsigned int n_q_points = mapping_data.quadrature_points.size();
325 std::vector<double>
values(
327 std::vector<Tensor<1, dim>> grads(
329 std::vector<Tensor<2, dim>> grad_grads(
331 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
332 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
335 for (
unsigned int i = 0; i < n_q_points; ++i)
341 empty_vector_of_3rd_order_tensors,
342 empty_vector_of_4th_order_tensors);
346 output_data.shape_values[k][i] =
values[k];
350 output_data.shape_gradients[k][i] = grads[k];
354 output_data.shape_hessians[k][i] = grad_grads[k];
360 template <
int dim,
int spacedim>
368 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
379 const unsigned int n_q_points = mapping_data.quadrature_points.size();
381 std::vector<double>
values(
383 std::vector<Tensor<1, dim>> grads(
385 std::vector<Tensor<2, dim>> grad_grads(
387 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
388 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
391 for (
unsigned int i = 0; i < n_q_points; ++i)
397 empty_vector_of_3rd_order_tensors,
398 empty_vector_of_4th_order_tensors);
402 output_data.shape_values[k][i] =
values[k];
406 output_data.shape_gradients[k][i] = grads[k];
410 output_data.shape_hessians[k][i] = grad_grads[k];
416 template <
int dim,
int spacedim>
425 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
436 const unsigned int n_q_points = mapping_data.quadrature_points.size();
438 std::vector<double>
values(
440 std::vector<Tensor<1, dim>> grads(
442 std::vector<Tensor<2, dim>> grad_grads(
444 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
445 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
448 for (
unsigned int i = 0; i < n_q_points; ++i)
454 empty_vector_of_3rd_order_tensors,
455 empty_vector_of_4th_order_tensors);
459 output_data.shape_values[k][i] =
values[k];
463 output_data.shape_gradients[k][i] = grads[k];
467 output_data.shape_hessians[k][i] = grad_grads[k];
473 template <
int dim,
int spacedim>
478 const unsigned int)
const 486 (void)interpolation_matrix;
488 (x_source_fe.
get_name().find(
"FE_DGPNonparametric<") == 0) ||
493 Assert(interpolation_matrix.
m() == 0,
495 Assert(interpolation_matrix.
n() == 0,
501 template <
int dim,
int spacedim>
507 const unsigned int)
const 515 (void)interpolation_matrix;
517 (x_source_fe.
get_name().find(
"FE_DGPNonparametric<") == 0) ||
522 Assert(interpolation_matrix.
m() == 0,
524 Assert(interpolation_matrix.
n() == 0,
530 template <
int dim,
int spacedim>
539 template <
int dim,
int spacedim>
540 std::vector<std::pair<unsigned int, unsigned int>>
548 return std::vector<std::pair<unsigned int, unsigned int>>();
552 return std::vector<std::pair<unsigned int, unsigned int>>();
558 template <
int dim,
int spacedim>
559 std::vector<std::pair<unsigned int, unsigned int>>
567 return std::vector<std::pair<unsigned int, unsigned int>>();
571 return std::vector<std::pair<unsigned int, unsigned int>>();
577 template <
int dim,
int spacedim>
578 std::vector<std::pair<unsigned int, unsigned int>>
581 const unsigned int)
const 587 return std::vector<std::pair<unsigned int, unsigned int>>();
591 return std::vector<std::pair<unsigned int, unsigned int>>();
597 template <
int dim,
int spacedim>
601 const unsigned int codim)
const 618 if (this->degree < fe_nonparametric_other->
degree)
620 else if (this->degree == fe_nonparametric_other->degree)
628 if (fe_nothing->is_dominating())
643 template <
int dim,
int spacedim>
647 const unsigned int)
const 654 template <
int dim,
int spacedim>
664 template <
int dim,
int spacedim>
674 #include "fe_dgp_nonparametric.inst"
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
friend class FE_DGPNonparametric
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
virtual std::size_t memory_consumption() 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 double shape_value(const unsigned int i, const Point< dim > &p) const override
#define AssertIndexRange(index, range)
const unsigned int degree
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
const PolynomialSpace< dim > polynomial_space
Transformed quadrature points.
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
#define AssertThrow(cond, exc)
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 UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual bool hp_constraints_are_implemented() const override
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const override
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
unsigned int get_degree() const
std::vector< std::vector< FullMatrix< double > > > prolongation
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
#define Assert(cond, exc)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
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
#define DEAL_II_NAMESPACE_CLOSE
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
virtual std::string get_name() const =0
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Second derivatives of shape functions.
std::string dim_string(const int dim, const int spacedim)
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const override
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
unsigned int n_dofs_per_cell() const
#define DEAL_II_NAMESPACE_OPEN
Shape function gradients.
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 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
static ::ExceptionBase & ExcNotImplemented()
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
static ::ExceptionBase & ExcInternalError()