Reference documentation for deal.II version Git ac854ef673 2021-03-06 10:09:45 +0100
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_dgp_nonparametric.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2002 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
18 
20 
21 #include <deal.II/fe/fe.h>
23 #include <deal.II/fe/fe_nothing.h>
24 #include <deal.II/fe/fe_values.h>
25 #include <deal.II/fe/mapping.h>
26 
27 #include <deal.II/grid/tria.h>
29 
30 #include <memory>
31 #include <sstream>
32 
33 
35 
36 template <int dim, int spacedim>
38  const unsigned int degree)
39  : FiniteElement<dim, spacedim>(
40  FiniteElementData<dim>(get_dpo_vector(degree),
41  1,
42  degree,
43  FiniteElementData<dim>::L2),
44  std::vector<bool>(
45  FiniteElementData<dim>(get_dpo_vector(degree), 1, degree)
46  .n_dofs_per_cell(),
47  true),
48  std::vector<ComponentMask>(
49  FiniteElementData<dim>(get_dpo_vector(degree), 1, degree)
50  .n_dofs_per_cell(),
51  std::vector<bool>(1, true)))
52  , polynomial_space(Polynomials::Legendre::generate_complete_basis(degree))
53 {
54  const unsigned int n_dofs = this->n_dofs_per_cell();
55  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
56  ref_case < RefinementCase<dim>::isotropic_refinement + 1;
57  ++ref_case)
58  {
59  if (dim != 2 && ref_case != RefinementCase<dim>::isotropic_refinement)
60  // do nothing, as anisotropic
61  // refinement is not
62  // implemented so far
63  continue;
64 
65  const unsigned int nc =
67  for (unsigned int i = 0; i < nc; ++i)
68  {
69  this->prolongation[ref_case - 1][i].reinit(n_dofs, n_dofs);
70  // Fill prolongation matrices with
71  // embedding operators
72  for (unsigned int j = 0; j < n_dofs; ++j)
73  this->prolongation[ref_case - 1][i](j, j) = 1.;
74  }
75  }
76 
77  // restriction can be defined
78  // through projection for
79  // discontinuous elements, but is
80  // presently not implemented for DGPNonparametric
81  // elements.
82  //
83  // if it were, then the following
84  // snippet would be the right code
85  // if ((degree < Matrices::n_projection_matrices) &&
86  // (Matrices::projection_matrices[degree] != 0))
87  // {
88  // restriction[0].fill (Matrices::projection_matrices[degree]);
89  // }
90  // else
91  // // matrix undefined, set size to zero
92  // for (unsigned int i=0;i<GeometryInfo<dim>::max_children_per_cell;++i)
93  // restriction[i].reinit(0, 0);
94  // since not implemented, set to
95  // "empty". however, that is done in the
96  // default constructor already, so do nothing
97  // for (unsigned int i=0;i<GeometryInfo<dim>::max_children_per_cell;++i)
98  // this->restriction[i].reinit(0, 0);
99 
100  // note further, that these
101  // elements have neither support
102  // nor face-support points, so
103  // leave these fields empty
104 }
105 
106 
107 
108 template <int dim, int spacedim>
109 std::string
111 {
112  // note that the
113  // FETools::get_fe_by_name
114  // function depends on the
115  // particular format of the string
116  // this function returns, so they
117  // have to be kept in synch
118 
119  std::ostringstream namebuf;
120  namebuf << "FE_DGPNonparametric<" << Utilities::dim_string(dim, spacedim)
121  << ">(" << this->degree << ")";
122 
123  return namebuf.str();
124 }
125 
126 
127 
128 template <int dim, int spacedim>
129 std::unique_ptr<FiniteElement<dim, spacedim>>
131 {
132  return std::make_unique<FE_DGPNonparametric<dim, spacedim>>(*this);
133 }
134 
135 
136 
137 template <int dim, int spacedim>
138 double
140  const Point<dim> & p) const
141 {
142  (void)i;
143  (void)p;
144  AssertIndexRange(i, this->n_dofs_per_cell());
145  AssertThrow(false,
147  return 0;
148 }
149 
150 
151 
152 template <int dim, int spacedim>
153 double
155  const unsigned int i,
156  const Point<dim> & p,
157  const unsigned int component) const
158 {
159  (void)i;
160  (void)p;
161  (void)component;
162  AssertIndexRange(i, this->n_dofs_per_cell());
163  AssertIndexRange(component, 1);
164  AssertThrow(false,
166  return 0;
167 }
168 
169 
170 
171 template <int dim, int spacedim>
174  const Point<dim> & p) const
175 {
176  (void)i;
177  (void)p;
178  AssertIndexRange(i, this->n_dofs_per_cell());
179  AssertThrow(false,
181  return Tensor<1, dim>();
182 }
183 
184 
185 template <int dim, int spacedim>
188  const unsigned int i,
189  const Point<dim> & p,
190  const unsigned int component) const
191 {
192  (void)i;
193  (void)p;
194  (void)component;
195  AssertIndexRange(i, this->n_dofs_per_cell());
196  AssertIndexRange(component, 1);
197  AssertThrow(false,
199  return Tensor<1, dim>();
200 }
201 
202 
203 
204 template <int dim, int spacedim>
207  const Point<dim> & p) const
208 {
209  (void)i;
210  (void)p;
211  AssertIndexRange(i, this->n_dofs_per_cell());
212  AssertThrow(false,
214  return Tensor<2, dim>();
215 }
216 
217 
218 
219 template <int dim, int spacedim>
222  const unsigned int i,
223  const Point<dim> & p,
224  const unsigned int component) const
225 {
226  (void)i;
227  (void)p;
228  (void)component;
229  AssertIndexRange(i, this->n_dofs_per_cell());
230  AssertIndexRange(component, 1);
231  AssertThrow(false,
233  return Tensor<2, dim>();
234 }
235 
236 
237 //---------------------------------------------------------------------------
238 // Auxiliary functions
239 //---------------------------------------------------------------------------
240 
241 
242 template <int dim, int spacedim>
243 std::vector<unsigned int>
245 {
246  std::vector<unsigned int> dpo(dim + 1, static_cast<unsigned int>(0));
247  dpo[dim] = deg + 1;
248  for (unsigned int i = 1; i < dim; ++i)
249  {
250  dpo[dim] *= deg + 1 + i;
251  dpo[dim] /= i + 1;
252  }
253  return dpo;
254 }
255 
256 
257 
258 template <int dim, int spacedim>
261  const UpdateFlags flags) const
262 {
263  UpdateFlags out = flags;
264 
267 
268  return out;
269 }
270 
271 
272 
273 //---------------------------------------------------------------------------
274 // Data field initialization
275 //---------------------------------------------------------------------------
276 
277 template <int dim, int spacedim>
278 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
280  const UpdateFlags update_flags,
281  const Mapping<dim, spacedim> &,
282  const Quadrature<dim> &,
284  spacedim>
285  & /*output_data*/) const
286 {
287  // generate a new data object
288  auto data_ptr =
289  std::make_unique<typename FiniteElement<dim, spacedim>::InternalDataBase>();
290  data_ptr->update_each = requires_update_flags(update_flags);
291 
292  // other than that, there is nothing we can add here as discussed
293  // in the general documentation of this class
294 
295  return data_ptr;
296 }
297 
298 
299 
300 //---------------------------------------------------------------------------
301 // Fill data of FEValues
302 //---------------------------------------------------------------------------
303 
304 template <int dim, int spacedim>
305 void
309  const Quadrature<dim> &,
310  const Mapping<dim, spacedim> &,
312  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
313  spacedim>
314  & mapping_data,
315  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
317  spacedim>
318  &output_data) const
319 {
321  ExcInternalError());
322 
323  const unsigned int n_q_points = mapping_data.quadrature_points.size();
324 
325  std::vector<double> values(
326  (fe_internal.update_each & update_values) ? this->n_dofs_per_cell() : 0);
327  std::vector<Tensor<1, dim>> grads(
328  (fe_internal.update_each & update_gradients) ? this->n_dofs_per_cell() : 0);
329  std::vector<Tensor<2, dim>> grad_grads(
330  (fe_internal.update_each & update_hessians) ? this->n_dofs_per_cell() : 0);
331  std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
332  std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
333 
334  if (fe_internal.update_each & (update_values | update_gradients))
335  for (unsigned int i = 0; i < n_q_points; ++i)
336  {
337  polynomial_space.evaluate(mapping_data.quadrature_points[i],
338  values,
339  grads,
340  grad_grads,
341  empty_vector_of_3rd_order_tensors,
342  empty_vector_of_4th_order_tensors);
343 
344  if (fe_internal.update_each & update_values)
345  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
346  output_data.shape_values[k][i] = values[k];
347 
348  if (fe_internal.update_each & update_gradients)
349  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
350  output_data.shape_gradients[k][i] = grads[k];
351 
352  if (fe_internal.update_each & update_hessians)
353  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
354  output_data.shape_hessians[k][i] = grad_grads[k];
355  }
356 }
357 
358 
359 
360 template <int dim, int spacedim>
361 void
364  const unsigned int,
365  const hp::QCollection<dim - 1> &,
366  const Mapping<dim, spacedim> &,
368  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
369  spacedim>
370  & mapping_data,
371  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
373  spacedim>
374  &output_data) const
375 {
377  ExcInternalError());
378 
379  const unsigned int n_q_points = mapping_data.quadrature_points.size();
380 
381  std::vector<double> values(
382  (fe_internal.update_each & update_values) ? this->n_dofs_per_cell() : 0);
383  std::vector<Tensor<1, dim>> grads(
384  (fe_internal.update_each & update_gradients) ? this->n_dofs_per_cell() : 0);
385  std::vector<Tensor<2, dim>> grad_grads(
386  (fe_internal.update_each & update_hessians) ? this->n_dofs_per_cell() : 0);
387  std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
388  std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
389 
390  if (fe_internal.update_each & (update_values | update_gradients))
391  for (unsigned int i = 0; i < n_q_points; ++i)
392  {
393  polynomial_space.evaluate(mapping_data.quadrature_points[i],
394  values,
395  grads,
396  grad_grads,
397  empty_vector_of_3rd_order_tensors,
398  empty_vector_of_4th_order_tensors);
399 
400  if (fe_internal.update_each & update_values)
401  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
402  output_data.shape_values[k][i] = values[k];
403 
404  if (fe_internal.update_each & update_gradients)
405  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
406  output_data.shape_gradients[k][i] = grads[k];
407 
408  if (fe_internal.update_each & update_hessians)
409  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
410  output_data.shape_hessians[k][i] = grad_grads[k];
411  }
412 }
413 
414 
415 
416 template <int dim, int spacedim>
417 void
420  const unsigned int,
421  const unsigned int,
422  const Quadrature<dim - 1> &,
423  const Mapping<dim, spacedim> &,
425  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
426  spacedim>
427  & mapping_data,
428  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
430  spacedim>
431  &output_data) const
432 {
434  ExcInternalError());
435 
436  const unsigned int n_q_points = mapping_data.quadrature_points.size();
437 
438  std::vector<double> values(
439  (fe_internal.update_each & update_values) ? this->n_dofs_per_cell() : 0);
440  std::vector<Tensor<1, dim>> grads(
441  (fe_internal.update_each & update_gradients) ? this->n_dofs_per_cell() : 0);
442  std::vector<Tensor<2, dim>> grad_grads(
443  (fe_internal.update_each & update_hessians) ? this->n_dofs_per_cell() : 0);
444  std::vector<Tensor<3, dim>> empty_vector_of_3rd_order_tensors;
445  std::vector<Tensor<4, dim>> empty_vector_of_4th_order_tensors;
446 
447  if (fe_internal.update_each & (update_values | update_gradients))
448  for (unsigned int i = 0; i < n_q_points; ++i)
449  {
450  polynomial_space.evaluate(mapping_data.quadrature_points[i],
451  values,
452  grads,
453  grad_grads,
454  empty_vector_of_3rd_order_tensors,
455  empty_vector_of_4th_order_tensors);
456 
457  if (fe_internal.update_each & update_values)
458  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
459  output_data.shape_values[k][i] = values[k];
460 
461  if (fe_internal.update_each & update_gradients)
462  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
463  output_data.shape_gradients[k][i] = grads[k];
464 
465  if (fe_internal.update_each & update_hessians)
466  for (unsigned int k = 0; k < this->n_dofs_per_cell(); ++k)
467  output_data.shape_hessians[k][i] = grad_grads[k];
468  }
469 }
470 
471 
472 
473 template <int dim, int spacedim>
474 void
476  const FiniteElement<dim, spacedim> &x_source_fe,
477  FullMatrix<double> & interpolation_matrix,
478  const unsigned int) const
479 {
480  // this is only implemented, if the source
481  // FE is also a DGPNonparametric element. in that case,
482  // both elements have no dofs on their
483  // faces and the face interpolation matrix
484  // is necessarily empty -- i.e. there isn't
485  // much we need to do here.
486  (void)interpolation_matrix;
487  AssertThrow(
488  (x_source_fe.get_name().find("FE_DGPNonparametric<") == 0) ||
489  (dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&x_source_fe) !=
490  nullptr),
492 
493  Assert(interpolation_matrix.m() == 0,
494  ExcDimensionMismatch(interpolation_matrix.m(), 0));
495  Assert(interpolation_matrix.n() == 0,
496  ExcDimensionMismatch(interpolation_matrix.n(), 0));
497 }
498 
499 
500 
501 template <int dim, int spacedim>
502 void
504  const FiniteElement<dim, spacedim> &x_source_fe,
505  const unsigned int,
506  FullMatrix<double> &interpolation_matrix,
507  const unsigned int) const
508 {
509  // this is only implemented, if the source
510  // FE is also a DGPNonparametric element. in that case,
511  // both elements have no dofs on their
512  // faces and the face interpolation matrix
513  // is necessarily empty -- i.e. there isn't
514  // much we need to do here.
515  (void)interpolation_matrix;
516  AssertThrow(
517  (x_source_fe.get_name().find("FE_DGPNonparametric<") == 0) ||
518  (dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&x_source_fe) !=
519  nullptr),
521 
522  Assert(interpolation_matrix.m() == 0,
523  ExcDimensionMismatch(interpolation_matrix.m(), 0));
524  Assert(interpolation_matrix.n() == 0,
525  ExcDimensionMismatch(interpolation_matrix.n(), 0));
526 }
527 
528 
529 
530 template <int dim, int spacedim>
531 bool
533 {
534  return true;
535 }
536 
537 
538 
539 template <int dim, int spacedim>
540 std::vector<std::pair<unsigned int, unsigned int>>
542  const FiniteElement<dim, spacedim> &fe_other) const
543 {
544  // there are no such constraints for DGPNonparametric
545  // elements at all
546  if (dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&fe_other) !=
547  nullptr)
548  return std::vector<std::pair<unsigned int, unsigned int>>();
549  else
550  {
551  Assert(false, ExcNotImplemented());
552  return std::vector<std::pair<unsigned int, unsigned int>>();
553  }
554 }
555 
556 
557 
558 template <int dim, int spacedim>
559 std::vector<std::pair<unsigned int, unsigned int>>
561  const FiniteElement<dim, spacedim> &fe_other) const
562 {
563  // there are no such constraints for DGPNonparametric
564  // elements at all
565  if (dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&fe_other) !=
566  nullptr)
567  return std::vector<std::pair<unsigned int, unsigned int>>();
568  else
569  {
570  Assert(false, ExcNotImplemented());
571  return std::vector<std::pair<unsigned int, unsigned int>>();
572  }
573 }
574 
575 
576 
577 template <int dim, int spacedim>
578 std::vector<std::pair<unsigned int, unsigned int>>
580  const FiniteElement<dim, spacedim> &fe_other,
581  const unsigned int) const
582 {
583  // there are no such constraints for DGPNonparametric
584  // elements at all
585  if (dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&fe_other) !=
586  nullptr)
587  return std::vector<std::pair<unsigned int, unsigned int>>();
588  else
589  {
590  Assert(false, ExcNotImplemented());
591  return std::vector<std::pair<unsigned int, unsigned int>>();
592  }
593 }
594 
595 
596 
597 template <int dim, int spacedim>
600  const FiniteElement<dim, spacedim> &fe_other,
601  const unsigned int codim) const
602 {
603  Assert(codim <= dim, ExcImpossibleInDim(dim));
604 
605  // vertex/line/face domination
606  // ---------------------------
607  if (codim > 0)
608  // this is a discontinuous element, so by definition there will
609  // be no constraints wherever this element comes together with
610  // any other kind of element
612 
613  // cell domination
614  // ---------------
615  if (const FE_DGPNonparametric<dim, spacedim> *fe_nonparametric_other =
616  dynamic_cast<const FE_DGPNonparametric<dim, spacedim> *>(&fe_other))
617  {
618  if (this->degree < fe_nonparametric_other->degree)
620  else if (this->degree == fe_nonparametric_other->degree)
622  else
624  }
625  else if (const FE_Nothing<dim> *fe_nothing =
626  dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
627  {
628  if (fe_nothing->is_dominating())
630  else
631  // the FE_Nothing has no degrees of freedom and it is typically used
632  // in a context where we don't require any continuity along the
633  // interface
635  }
636 
637  Assert(false, ExcNotImplemented());
639 }
640 
641 
642 
643 template <int dim, int spacedim>
644 bool
646  const unsigned int,
647  const unsigned int) const
648 {
649  return true;
650 }
651 
652 
653 
654 template <int dim, int spacedim>
655 std::size_t
657 {
658  Assert(false, ExcNotImplemented());
659  return 0;
660 }
661 
662 
663 
664 template <int dim, int spacedim>
665 unsigned int
667 {
668  return this->degree;
669 }
670 
671 
672 
673 // explicit instantiations
674 #include "fe_dgp_nonparametric.inst"
675 
676 
Shape function values.
size_type m() const
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)
Definition: exceptions.h:1691
const unsigned int degree
Definition: fe_base.h:435
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
const PolynomialSpace< dim > polynomial_space
STL namespace.
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)
Definition: exceptions.h:1576
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)
size_type n() const
unsigned int get_degree() const
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2418
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
#define Assert(cond, exc)
Definition: exceptions.h:1466
UpdateFlags
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
Definition: mapping.h:303
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
Definition: config.h:394
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)
Definition: utilities.cc:558
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
Definition: config.h:393
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.)
Definition: l2.h:160
UpdateFlags update_each
Definition: fe.h:710
static ::ExceptionBase & ExcInternalError()