17 #include <deal.II/base/quadrature.h> 18 #include <deal.II/base/std_cxx14/memory.h> 20 #include <deal.II/dofs/dof_accessor.h> 22 #include <deal.II/fe/fe.h> 23 #include <deal.II/fe/fe_dgp_nonparametric.h> 24 #include <deal.II/fe/fe_values.h> 25 #include <deal.II/fe/mapping.h> 27 #include <deal.II/grid/tria.h> 28 #include <deal.II/grid/tria_iterator.h> 33 DEAL_II_NAMESPACE_OPEN
35 template <
int dim,
int spacedim>
37 const unsigned int degree)
48 std::vector<bool>(1, true)))
49 , polynomial_space(
Polynomials::Legendre::generate_complete_basis(degree))
53 ref_case < RefinementCase<dim>::isotropic_refinement + 1;
62 const unsigned int nc =
64 for (
unsigned int i = 0; i < nc; ++i)
66 this->
prolongation[ref_case - 1][i].reinit(n_dofs, n_dofs);
69 for (
unsigned int j = 0; j < n_dofs; ++j)
105 template <
int dim,
int spacedim>
116 std::ostringstream namebuf;
118 <<
">(" << this->
degree <<
")";
120 return namebuf.str();
125 template <
int dim,
int spacedim>
126 std::unique_ptr<FiniteElement<dim, spacedim>>
129 return std_cxx14::make_unique<FE_DGPNonparametric<dim, spacedim>>(*this);
134 template <
int dim,
int spacedim>
149 template <
int dim,
int spacedim>
152 const unsigned int i,
154 const unsigned int component)
const 168 template <
int dim,
int spacedim>
182 template <
int dim,
int spacedim>
185 const unsigned int i,
187 const unsigned int component)
const 201 template <
int dim,
int spacedim>
216 template <
int dim,
int spacedim>
219 const unsigned int i,
221 const unsigned int component)
const 239 template <
int dim,
int spacedim>
240 std::vector<unsigned int>
243 std::vector<unsigned int> dpo(dim + 1, static_cast<unsigned int>(0));
245 for (
unsigned int i = 1; i < dim; ++i)
247 dpo[dim] *= deg + 1 + i;
255 template <
int dim,
int spacedim>
274 template <
int dim,
int spacedim>
275 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
285 auto data = std_cxx14::make_unique<
292 return std::move(data);
301 template <
int dim,
int spacedim>
309 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
320 const unsigned int n_q_points = mapping_data.quadrature_points.size();
322 std::vector<double> values(
324 std::vector<Tensor<1, dim>> grads(
326 std::vector<Tensor<2, dim>> grad_grads(
328 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
329 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
332 for (
unsigned int i = 0; i < n_q_points; ++i)
338 empty_vector_of_3rd_order_tensors,
339 empty_vector_of_4th_order_tensors);
343 output_data.shape_values[k][i] = values[k];
347 output_data.shape_gradients[k][i] = grads[k];
351 output_data.shape_hessians[k][i] = grad_grads[k];
357 template <
int dim,
int spacedim>
365 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
376 const unsigned int n_q_points = mapping_data.quadrature_points.size();
378 std::vector<double> values(
380 std::vector<Tensor<1, dim>> grads(
382 std::vector<Tensor<2, dim>> grad_grads(
384 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
385 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
388 for (
unsigned int i = 0; i < n_q_points; ++i)
394 empty_vector_of_3rd_order_tensors,
395 empty_vector_of_4th_order_tensors);
399 output_data.shape_values[k][i] = values[k];
403 output_data.shape_gradients[k][i] = grads[k];
407 output_data.shape_hessians[k][i] = grad_grads[k];
413 template <
int dim,
int spacedim>
422 const ::internal::FEValuesImplementation::MappingRelatedData<dim,
433 const unsigned int n_q_points = mapping_data.quadrature_points.size();
435 std::vector<double> values(
437 std::vector<Tensor<1, dim>> grads(
439 std::vector<Tensor<2, dim>> grad_grads(
441 std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
442 std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
445 for (
unsigned int i = 0; i < n_q_points; ++i)
451 empty_vector_of_3rd_order_tensors,
452 empty_vector_of_4th_order_tensors);
456 output_data.shape_values[k][i] = values[k];
460 output_data.shape_gradients[k][i] = grads[k];
464 output_data.shape_hessians[k][i] = grad_grads[k];
470 template <
int dim,
int spacedim>
482 (void)interpolation_matrix;
484 (x_source_fe.
get_name().find(
"FE_DGPNonparametric<") == 0) ||
489 Assert(interpolation_matrix.
m() == 0,
491 Assert(interpolation_matrix.
n() == 0,
497 template <
int dim,
int spacedim>
510 (void)interpolation_matrix;
512 (x_source_fe.
get_name().find(
"FE_DGPNonparametric<") == 0) ||
517 Assert(interpolation_matrix.
m() == 0,
519 Assert(interpolation_matrix.
n() == 0,
525 template <
int dim,
int spacedim>
534 template <
int dim,
int spacedim>
535 std::vector<std::pair<unsigned int, unsigned int>>
543 return std::vector<std::pair<unsigned int, unsigned int>>();
547 return std::vector<std::pair<unsigned int, unsigned int>>();
553 template <
int dim,
int spacedim>
554 std::vector<std::pair<unsigned int, unsigned int>>
562 return std::vector<std::pair<unsigned int, unsigned int>>();
566 return std::vector<std::pair<unsigned int, unsigned int>>();
572 template <
int dim,
int spacedim>
573 std::vector<std::pair<unsigned int, unsigned int>>
581 return std::vector<std::pair<unsigned int, unsigned int>>();
585 return std::vector<std::pair<unsigned int, unsigned int>>();
591 template <
int dim,
int spacedim>
609 template <
int dim,
int spacedim>
613 const unsigned int)
const 620 template <
int dim,
int spacedim>
630 template <
int dim,
int spacedim>
640 #include "fe_dgp_nonparametric.inst" 643 DEAL_II_NAMESPACE_CLOSE
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 get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix) const override
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
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.
#define AssertThrow(cond, exc)
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual bool hp_constraints_are_implemented() const override
unsigned int get_degree() const
std::vector< std::vector< FullMatrix< double > > > prolongation
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) 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
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
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const override
Second derivatives of shape functions.
std::string dim_string(const int dim, const int spacedim)
const unsigned int dofs_per_cell
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
Shape function gradients.
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) 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 FiniteElementDomination::Domination compare_for_face_domination(const FiniteElement< dim, spacedim > &fe_other) 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
static::ExceptionBase & ExcInternalError()