Reference documentation for deal.II version Git a4a8ba9cef 2021-01-17 11:42:49 -0500
\(\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\}}\)
mapping_q_generic.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 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 #ifndef dealii_mapping_q_generic_h
17 #define dealii_mapping_q_generic_h
18 
19 
20 #include <deal.II/base/config.h>
21 
25 #include <deal.II/base/table.h>
27 
28 #include <deal.II/fe/mapping.h>
29 
31 
34 
35 #include <array>
36 #include <cmath>
37 
39 
40 template <int, int>
41 class MappingQ;
42 
43 template <int, int>
44 class MappingQCache;
45 
46 
49 
50 
134 template <int dim, int spacedim = dim>
135 class MappingQGeneric : public Mapping<dim, spacedim>
136 {
137 public:
143  MappingQGeneric(const unsigned int polynomial_degree);
144 
149 
150  // for documentation, see the Mapping base class
151  virtual std::unique_ptr<Mapping<dim, spacedim>>
152  clone() const override;
153 
158  unsigned int
159  get_degree() const;
160 
165  virtual bool
166  preserves_vertex_locations() const override;
167 
173  // for documentation, see the Mapping base class
174  virtual Point<spacedim>
176  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
177  const Point<dim> &p) const override;
178 
179  // for documentation, see the Mapping base class
180  virtual Point<dim>
182  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
183  const Point<spacedim> &p) const override;
184 
185  // for documentation, see the Mapping base class
186  virtual void
188  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
189  const ArrayView<const Point<spacedim>> & real_points,
190  const ArrayView<Point<dim>> &unit_points) const override;
191 
201  // for documentation, see the Mapping base class
202  virtual void
203  transform(const ArrayView<const Tensor<1, dim>> & input,
204  const MappingKind kind,
205  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
206  const ArrayView<Tensor<1, spacedim>> &output) const override;
207 
208  // for documentation, see the Mapping base class
209  virtual void
211  const MappingKind kind,
212  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
213  const ArrayView<Tensor<2, spacedim>> &output) const override;
214 
215  // for documentation, see the Mapping base class
216  virtual void
217  transform(const ArrayView<const Tensor<2, dim>> & input,
218  const MappingKind kind,
219  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
220  const ArrayView<Tensor<2, spacedim>> &output) const override;
221 
222  // for documentation, see the Mapping base class
223  virtual void
225  const MappingKind kind,
226  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
227  const ArrayView<Tensor<3, spacedim>> &output) const override;
228 
229  // for documentation, see the Mapping base class
230  virtual void
231  transform(const ArrayView<const Tensor<3, dim>> & input,
232  const MappingKind kind,
233  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
234  const ArrayView<Tensor<3, spacedim>> &output) const override;
235 
244  template <typename Number, typename Number2>
245  void
247  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
248  const MappingKind kind,
249  const ArrayView<const Point<dim>> & unit_points,
250  const ArrayView<const Number> & input,
251  const ArrayView<Number2> & output) const;
252 
273  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
274  {
275  public:
280  InternalData(const unsigned int polynomial_degree);
281 
290  void
291  initialize(const UpdateFlags update_flags,
292  const Quadrature<dim> &quadrature,
293  const unsigned int n_original_q_points);
294 
300  void
301  initialize_face(const UpdateFlags update_flags,
302  const Quadrature<dim> &quadrature,
303  const unsigned int n_original_q_points);
304 
320  void
321  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
322 
327  const double &
328  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
329 
333  double &
334  shape(const unsigned int qpoint, const unsigned int shape_nr);
335 
339  const Tensor<1, dim> &
340  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
341 
346  derivative(const unsigned int qpoint, const unsigned int shape_nr);
347 
351  const Tensor<2, dim> &
352  second_derivative(const unsigned int qpoint,
353  const unsigned int shape_nr) const;
354 
359  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
360 
364  const Tensor<3, dim> &
365  third_derivative(const unsigned int qpoint,
366  const unsigned int shape_nr) const;
367 
372  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
373 
377  const Tensor<4, dim> &
378  fourth_derivative(const unsigned int qpoint,
379  const unsigned int shape_nr) const;
380 
385  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
386 
390  virtual std::size_t
391  memory_consumption() const override;
392 
398  std::vector<double> shape_values;
399 
405  std::vector<Tensor<1, dim>> shape_derivatives;
406 
413  std::vector<Tensor<2, dim>> shape_second_derivatives;
414 
421  std::vector<Tensor<3, dim>> shape_third_derivatives;
422 
429  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
430 
444  std::array<std::vector<Tensor<1, dim>>,
447 
452  const unsigned int polynomial_degree;
453 
463  const unsigned int n_shape_functions;
464 
465  /*
466  * The default line support points. Is used in when the shape function
467  * values are computed.
468  *
469  * The number of quadrature points depends on the degree of this
470  * class, and it matches the number of degrees of freedom of an
471  * FE_Q<1>(this->degree).
472  */
474 
482 
488 
494 
500 
506 
512 
517 
527  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
528 
536  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
537 
541  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
542 
547  mutable std::vector<Point<spacedim>> mapping_support_points;
548 
554 
559  mutable std::vector<double> volume_elements;
560  };
561 
562 
563  // documentation can be found in Mapping::requires_update_flags()
564  virtual UpdateFlags
565  requires_update_flags(const UpdateFlags update_flags) const override;
566 
567  // documentation can be found in Mapping::get_data()
568  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
569  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
570 
572 
573  // documentation can be found in Mapping::get_face_data()
574  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
575  get_face_data(const UpdateFlags flags,
576  const hp::QCollection<dim - 1> &quadrature) const override;
577 
578  // documentation can be found in Mapping::get_subface_data()
579  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
580  get_subface_data(const UpdateFlags flags,
581  const Quadrature<dim - 1> &quadrature) const override;
582 
583  // documentation can be found in Mapping::fill_fe_values()
586  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
587  const CellSimilarity::Similarity cell_similarity,
588  const Quadrature<dim> & quadrature,
589  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
591  &output_data) const override;
592 
594 
595  // documentation can be found in Mapping::fill_fe_face_values()
596  virtual void
598  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
599  const unsigned int face_no,
600  const hp::QCollection<dim - 1> & quadrature,
601  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
603  &output_data) const override;
604 
605  // documentation can be found in Mapping::fill_fe_subface_values()
606  virtual void
608  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
609  const unsigned int face_no,
610  const unsigned int subface_no,
611  const Quadrature<dim - 1> & quadrature,
612  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
614  &output_data) const override;
615 
620 protected:
625  const unsigned int polynomial_degree;
626 
627  /*
628  * The default line support points. These are used when computing the
629  * location in real space of the support points on lines and quads, which
630  * are needed by the Manifold<dim,spacedim> class.
631  *
632  * The number of points depends on the degree of this class, and it matches
633  * the number of degrees of freedom of an FE_Q<1>(this->degree).
634  */
635  const std::vector<Point<1>> line_support_points;
636 
637  /*
638  * The one-dimensional polynomials defined as Lagrange polynomials from the
639  * line support points. These are used for point evaluations and match the
640  * polynomial space of an FE_Q<1>(this->degree).
641  */
642  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
643 
644  /*
645  * The numbering from the lexicographic to the hierarchical ordering used
646  * when expanding the tensor product with the mapping support points (which
647  * come in hierarchical numbers).
648  */
649  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
650 
651  /*
652  * The support points in reference coordinates. These are used for
653  * constructing approximations of the output of
654  * compute_mapping_support_points() when evaluating the mapping on the fly,
655  * rather than going through the FEValues interface provided by
656  * InternalData.
657  *
658  * The number of points depends on the degree of this class, and it matches
659  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
660  */
661  const std::vector<Point<dim>> unit_cell_support_points;
662 
682  const std::vector<Table<2, double>>
684 
698 
728  virtual std::vector<Point<spacedim>>
730  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
731 
736  Point<dim>
738  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
739  const Point<spacedim> & p,
740  const Point<dim> &initial_p_unit) const;
741 
755  virtual void
757  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
758  std::vector<Point<spacedim>> & a) const;
759 
774  virtual void
776  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
777  std::vector<Point<spacedim>> & a) const;
778 
779  // Make MappingQ a friend since it needs to call the fill_fe_values()
780  // functions on its MappingQGeneric(1) sub-object.
781  template <int, int>
782  friend class MappingQ;
783 
784  // Make MappingQCache a friend since it needs to call the
785  // compute_mapping_support_points() function.
786  template <int, int>
787  friend class MappingQCache;
788 };
789 
790 
791 
794 /*----------------------------------------------------------------------*/
795 
796 #ifndef DOXYGEN
797 
798 template <int dim, int spacedim>
799 inline const double &
801  const unsigned int qpoint,
802  const unsigned int shape_nr) const
803 {
804  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
805  return shape_values[qpoint * n_shape_functions + shape_nr];
806 }
807 
808 
809 
810 template <int dim, int spacedim>
811 inline double &
813  const unsigned int shape_nr)
814 {
815  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
816  return shape_values[qpoint * n_shape_functions + shape_nr];
817 }
818 
819 
820 template <int dim, int spacedim>
821 inline const Tensor<1, dim> &
823  const unsigned int qpoint,
824  const unsigned int shape_nr) const
825 {
826  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
827  shape_derivatives.size());
828  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
829 }
830 
831 
832 
833 template <int dim, int spacedim>
834 inline Tensor<1, dim> &
836  const unsigned int qpoint,
837  const unsigned int shape_nr)
838 {
839  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
840  shape_derivatives.size());
841  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
842 }
843 
844 
845 template <int dim, int spacedim>
846 inline const Tensor<2, dim> &
848  const unsigned int qpoint,
849  const unsigned int shape_nr) const
850 {
851  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
852  shape_second_derivatives.size());
853  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
854 }
855 
856 
857 template <int dim, int spacedim>
858 inline Tensor<2, dim> &
860  const unsigned int qpoint,
861  const unsigned int shape_nr)
862 {
863  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
864  shape_second_derivatives.size());
865  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
866 }
867 
868 template <int dim, int spacedim>
869 inline const Tensor<3, dim> &
871  const unsigned int qpoint,
872  const unsigned int shape_nr) const
873 {
874  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
875  shape_third_derivatives.size());
876  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
877 }
878 
879 
880 template <int dim, int spacedim>
881 inline Tensor<3, dim> &
883  const unsigned int qpoint,
884  const unsigned int shape_nr)
885 {
886  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
887  shape_third_derivatives.size());
888  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
889 }
890 
891 
892 template <int dim, int spacedim>
893 inline const Tensor<4, dim> &
895  const unsigned int qpoint,
896  const unsigned int shape_nr) const
897 {
898  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
899  shape_fourth_derivatives.size());
900  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
901 }
902 
903 
904 template <int dim, int spacedim>
905 inline Tensor<4, dim> &
907  const unsigned int qpoint,
908  const unsigned int shape_nr)
909 {
910  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
911  shape_fourth_derivatives.size());
912  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
913 }
914 
915 
916 
917 template <int dim, int spacedim>
918 inline bool
920 {
921  return true;
922 }
923 
924 
925 
926 template <int dim, int spacedim>
927 template <typename Number, typename Number2>
928 inline void
930  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
931  const MappingKind kind,
932  const ArrayView<const Point<dim>> & unit_points,
933  const ArrayView<const Number> & input,
934  const ArrayView<Number2> & output) const
935 {
936  AssertDimension(unit_points.size(), output.size());
937  AssertDimension(unit_points.size(), input.size());
938 
939  const std::vector<Point<spacedim>> support_points =
940  this->compute_mapping_support_points(cell);
941 
942  const unsigned int n_points = unit_points.size();
943  const unsigned int n_lanes = VectorizedArray<double>::size();
944 
945  if (kind != mapping_covariant)
946  Assert(false, ExcNotImplemented());
947 
948  // Use the more heavy VectorizedArray code path if there is more than
949  // one point left to compute
950  for (unsigned int i = 0; i < n_points; i += n_lanes)
951  if (n_points - i > 1)
952  {
954  for (unsigned int j = 0; j < n_lanes; ++j)
955  if (i + j < n_points)
956  for (unsigned int d = 0; d < spacedim; ++d)
957  p_vec[d][j] = unit_points[i + j][d];
958  else
959  for (unsigned int d = 0; d < spacedim; ++d)
960  p_vec[d][j] = unit_points[i][d];
961 
965  support_points,
966  p_vec,
967  polynomial_degree == 1,
969  .second;
970 
972  grad.transpose().covariant_form();
973  for (unsigned int j = 0; j < n_lanes && i + j < n_points; ++j)
974  {
976  for (unsigned int d = 0; d < spacedim; ++d)
977  for (unsigned int e = 0; e < dim; ++e)
978  jac_j[d][e] = jac[d][e][j];
979  output[i + j] = apply_transformation(jac_j, input[i + j]);
980  }
981  }
982  else
983  {
987  support_points,
988  unit_points[i],
989  polynomial_degree == 1,
991  .second;
992  output[i] =
993  apply_transformation(grad.transpose().covariant_form(), input[i]);
994  }
995 }
996 
997 #endif // DOXYGEN
998 
999 /* -------------- declaration of explicit specializations ------------- */
1000 
1001 
1003 
1004 #endif
std::vector< Tensor< 2, dim > > shape_second_derivatives
Tensor< 1, spacedim, Number > apply_transformation(const DerivativeForm< 1, dim, spacedim, Number > &grad_F, const Tensor< 1, dim, Number > &d_x)
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
const unsigned int polynomial_degree
virtual void add_quad_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
AlignedVector< VectorizedArray< double > > values_quad
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
const Tensor< 1, dim > & derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual void add_line_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, std::vector< Point< spacedim >> &a) const
const Tensor< 3, dim > & third_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
const Table< 2, double > support_point_weights_cell
AlignedVector< VectorizedArray< double > > gradients_quad
virtual CellSimilarity::Similarity fill_fe_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
MappingQGeneric(const unsigned int polynomial_degree)
AlignedVector< VectorizedArray< double > > scratch
MappingKind
Definition: mapping.h:64
InternalData(const unsigned int polynomial_degree)
const unsigned int polynomial_degree
std::size_t size() const
Definition: array_view.h:542
const std::vector< Polynomials::Polynomial< double > > polynomials_1d
QGaussLobatto< 1 > line_support_points
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
virtual void transform_points_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< spacedim >> &real_points, const ArrayView< Point< dim >> &unit_points) const override
void compute_shape_function_values(const std::vector< Point< dim >> &unit_points)
std::vector< std::vector< Tensor< 1, spacedim > > > aux
#define Assert(cond, exc)
Definition: exceptions.h:1466
const std::vector< Point< dim > > unit_cell_support_points
UpdateFlags
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingKind kind, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const override
Abstract base class for mapping classes.
Definition: mapping.h:303
void transform_variable(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const MappingKind kind, const ArrayView< const Point< dim >> &unit_points, const ArrayView< const Number > &input, const ArrayView< Number2 > &output) const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:380
std::vector< Point< spacedim > > mapping_support_points
std::vector< Tensor< 3, dim > > shape_third_derivatives
DerivativeForm< 1, spacedim, dim, Number > transpose() const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
void initialize(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
std::vector< double > volume_elements
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
const unsigned int n_shape_functions
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
std::vector< Tensor< 1, dim > > shape_derivatives
const std::vector< Table< 2, double > > support_point_weights_perimeter_to_interior
Point< dim > transform_real_to_unit_cell_internal(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p, const Point< dim > &initial_p_unit) const
AlignedVector< VectorizedArray< double > > hessians_quad
AlignedVector< VectorizedArray< double > > values_dofs
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 typename Mapping< dim, spacedim >::InternalDataBase &internal_data, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
unsigned int get_degree() const
const double & shape(const unsigned int qpoint, const unsigned int shape_nr) const
std::array< std::vector< Tensor< 1, dim > >, GeometryInfo< dim >::faces_per_cell *(dim - 1)> unit_tangentials
const std::vector< Point< 1 > > line_support_points
std::pair< typename ProductTypeNoPoint< Number, Number2 >::type, Tensor< 1, dim, typename ProductTypeNoPoint< Number, Number2 >::type > > evaluate_tensor_product_value_and_gradient(const std::vector< Polynomials::Polynomial< double >> &poly, const std::vector< Number > &values, const Point< dim, Number2 > &p, const bool d_linear=false, const std::vector< unsigned int > &renumber={})
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:379
virtual bool preserves_vertex_locations() const override
const Tensor< 2, dim > & second_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
Triangulation< dim, spacedim >::cell_iterator cell_of_current_support_points
virtual std::size_t memory_consumption() const override
virtual std::vector< Point< spacedim > > compute_mapping_support_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const
const std::vector< unsigned int > renumber_lexicographic_to_hierarchic
static ::ExceptionBase & ExcNotImplemented()
std::vector< DerivativeForm< 1, dim, spacedim > > contravariant
std::vector< Tensor< 4, dim > > shape_fourth_derivatives
void initialize_face(const UpdateFlags update_flags, const Quadrature< dim > &quadrature, const unsigned int n_original_q_points)
std::vector< double > shape_values
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
const Tensor< 4, dim > & fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr) const
virtual void fill_fe_subface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const Quadrature< dim - 1 > &quadrature, const typename Mapping< dim, spacedim >::InternalDataBase &internal_data, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
std::vector< DerivativeForm< 1, dim, spacedim > > covariant
internal::MatrixFreeFunctions::ShapeInfo< VectorizedArray< double > > shape_info