Reference documentation for deal.II version Git 54662f1c23 2021-02-26 09:51:10 -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 
168  // for documentation, see the Mapping base class
169  virtual BoundingBox<spacedim>
171  &cell) const override;
172 
173  virtual bool
174  is_compatible_with(const ReferenceCell &reference_cell) const override;
175 
181  // for documentation, see the Mapping base class
182  virtual Point<spacedim>
184  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
185  const Point<dim> &p) const override;
186 
187  // for documentation, see the Mapping base class
188  virtual Point<dim>
190  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
191  const Point<spacedim> &p) const override;
192 
193  // for documentation, see the Mapping base class
194  virtual void
196  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
197  const ArrayView<const Point<spacedim>> & real_points,
198  const ArrayView<Point<dim>> &unit_points) const override;
199 
209  // for documentation, see the Mapping base class
210  virtual void
211  transform(const ArrayView<const Tensor<1, dim>> & input,
212  const MappingKind kind,
213  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
214  const ArrayView<Tensor<1, spacedim>> &output) const override;
215 
216  // for documentation, see the Mapping base class
217  virtual void
219  const MappingKind kind,
220  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
221  const ArrayView<Tensor<2, spacedim>> &output) const override;
222 
223  // for documentation, see the Mapping base class
224  virtual void
225  transform(const ArrayView<const Tensor<2, dim>> & input,
226  const MappingKind kind,
227  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
228  const ArrayView<Tensor<2, spacedim>> &output) const override;
229 
230  // for documentation, see the Mapping base class
231  virtual void
233  const MappingKind kind,
234  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
235  const ArrayView<Tensor<3, spacedim>> &output) const override;
236 
237  // for documentation, see the Mapping base class
238  virtual void
239  transform(const ArrayView<const Tensor<3, dim>> & input,
240  const MappingKind kind,
241  const typename Mapping<dim, spacedim>::InternalDataBase &internal,
242  const ArrayView<Tensor<3, spacedim>> &output) const override;
243 
252  template <typename Number, typename Number2>
253  void
255  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
256  const MappingKind kind,
257  const ArrayView<const Point<dim>> & unit_points,
258  const ArrayView<const Number> & input,
259  const ArrayView<Number2> & output) const;
260 
281  class InternalData : public Mapping<dim, spacedim>::InternalDataBase
282  {
283  public:
288  InternalData(const unsigned int polynomial_degree);
289 
298  void
299  initialize(const UpdateFlags update_flags,
300  const Quadrature<dim> &quadrature,
301  const unsigned int n_original_q_points);
302 
308  void
309  initialize_face(const UpdateFlags update_flags,
310  const Quadrature<dim> &quadrature,
311  const unsigned int n_original_q_points);
312 
328  void
329  compute_shape_function_values(const std::vector<Point<dim>> &unit_points);
330 
335  const double &
336  shape(const unsigned int qpoint, const unsigned int shape_nr) const;
337 
341  double &
342  shape(const unsigned int qpoint, const unsigned int shape_nr);
343 
347  const Tensor<1, dim> &
348  derivative(const unsigned int qpoint, const unsigned int shape_nr) const;
349 
354  derivative(const unsigned int qpoint, const unsigned int shape_nr);
355 
359  const Tensor<2, dim> &
360  second_derivative(const unsigned int qpoint,
361  const unsigned int shape_nr) const;
362 
367  second_derivative(const unsigned int qpoint, const unsigned int shape_nr);
368 
372  const Tensor<3, dim> &
373  third_derivative(const unsigned int qpoint,
374  const unsigned int shape_nr) const;
375 
380  third_derivative(const unsigned int qpoint, const unsigned int shape_nr);
381 
385  const Tensor<4, dim> &
386  fourth_derivative(const unsigned int qpoint,
387  const unsigned int shape_nr) const;
388 
393  fourth_derivative(const unsigned int qpoint, const unsigned int shape_nr);
394 
398  virtual std::size_t
399  memory_consumption() const override;
400 
406  std::vector<double> shape_values;
407 
413  std::vector<Tensor<1, dim>> shape_derivatives;
414 
421  std::vector<Tensor<2, dim>> shape_second_derivatives;
422 
429  std::vector<Tensor<3, dim>> shape_third_derivatives;
430 
437  std::vector<Tensor<4, dim>> shape_fourth_derivatives;
438 
452  std::array<std::vector<Tensor<1, dim>>,
455 
460  const unsigned int polynomial_degree;
461 
471  const unsigned int n_shape_functions;
472 
473  /*
474  * The default line support points. Is used in when the shape function
475  * values are computed.
476  *
477  * The number of quadrature points depends on the degree of this
478  * class, and it matches the number of degrees of freedom of an
479  * FE_Q<1>(this->degree).
480  */
482 
490 
496 
502 
508 
514 
520 
525 
535  mutable std::vector<DerivativeForm<1, dim, spacedim>> covariant;
536 
544  mutable std::vector<DerivativeForm<1, dim, spacedim>> contravariant;
545 
549  mutable std::vector<std::vector<Tensor<1, spacedim>>> aux;
550 
555  mutable std::vector<Point<spacedim>> mapping_support_points;
556 
562 
567  mutable std::vector<double> volume_elements;
568  };
569 
570 
571  // documentation can be found in Mapping::requires_update_flags()
572  virtual UpdateFlags
573  requires_update_flags(const UpdateFlags update_flags) const override;
574 
575  // documentation can be found in Mapping::get_data()
576  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
577  get_data(const UpdateFlags, const Quadrature<dim> &quadrature) const override;
578 
580 
581  // documentation can be found in Mapping::get_face_data()
582  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
583  get_face_data(const UpdateFlags flags,
584  const hp::QCollection<dim - 1> &quadrature) const override;
585 
586  // documentation can be found in Mapping::get_subface_data()
587  virtual std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
588  get_subface_data(const UpdateFlags flags,
589  const Quadrature<dim - 1> &quadrature) const override;
590 
591  // documentation can be found in Mapping::fill_fe_values()
594  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
595  const CellSimilarity::Similarity cell_similarity,
596  const Quadrature<dim> & quadrature,
597  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
599  &output_data) const override;
600 
602 
603  // documentation can be found in Mapping::fill_fe_face_values()
604  virtual void
606  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
607  const unsigned int face_no,
608  const hp::QCollection<dim - 1> & quadrature,
609  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
611  &output_data) const override;
612 
613  // documentation can be found in Mapping::fill_fe_subface_values()
614  virtual void
616  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
617  const unsigned int face_no,
618  const unsigned int subface_no,
619  const Quadrature<dim - 1> & quadrature,
620  const typename Mapping<dim, spacedim>::InternalDataBase & internal_data,
622  &output_data) const override;
623 
628 protected:
633  const unsigned int polynomial_degree;
634 
635  /*
636  * The default line support points. These are used when computing the
637  * location in real space of the support points on lines and quads, which
638  * are needed by the Manifold<dim,spacedim> class.
639  *
640  * The number of points depends on the degree of this class, and it matches
641  * the number of degrees of freedom of an FE_Q<1>(this->degree).
642  */
643  const std::vector<Point<1>> line_support_points;
644 
645  /*
646  * The one-dimensional polynomials defined as Lagrange polynomials from the
647  * line support points. These are used for point evaluations and match the
648  * polynomial space of an FE_Q<1>(this->degree).
649  */
650  const std::vector<Polynomials::Polynomial<double>> polynomials_1d;
651 
652  /*
653  * The numbering from the lexicographic to the hierarchical ordering used
654  * when expanding the tensor product with the mapping support points (which
655  * come in hierarchical numbers).
656  */
657  const std::vector<unsigned int> renumber_lexicographic_to_hierarchic;
658 
659  /*
660  * The support points in reference coordinates. These are used for
661  * constructing approximations of the output of
662  * compute_mapping_support_points() when evaluating the mapping on the fly,
663  * rather than going through the FEValues interface provided by
664  * InternalData.
665  *
666  * The number of points depends on the degree of this class, and it matches
667  * the number of degrees of freedom of an FE_Q<dim>(this->degree).
668  */
669  const std::vector<Point<dim>> unit_cell_support_points;
670 
690  const std::vector<Table<2, double>>
692 
706 
736  virtual std::vector<Point<spacedim>>
738  const typename Triangulation<dim, spacedim>::cell_iterator &cell) const;
739 
744  Point<dim>
746  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
747  const Point<spacedim> & p,
748  const Point<dim> &initial_p_unit) const;
749 
763  virtual void
765  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
766  std::vector<Point<spacedim>> & a) const;
767 
782  virtual void
784  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
785  std::vector<Point<spacedim>> & a) const;
786 
787  // Make MappingQ a friend since it needs to call the fill_fe_values()
788  // functions on its MappingQGeneric(1) sub-object.
789  template <int, int>
790  friend class MappingQ;
791 
792  // Make MappingQCache a friend since it needs to call the
793  // compute_mapping_support_points() function.
794  template <int, int>
795  friend class MappingQCache;
796 };
797 
798 
799 
802 /*----------------------------------------------------------------------*/
803 
804 #ifndef DOXYGEN
805 
806 template <int dim, int spacedim>
807 inline const double &
809  const unsigned int qpoint,
810  const unsigned int shape_nr) const
811 {
812  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
813  return shape_values[qpoint * n_shape_functions + shape_nr];
814 }
815 
816 
817 
818 template <int dim, int spacedim>
819 inline double &
821  const unsigned int shape_nr)
822 {
823  AssertIndexRange(qpoint * n_shape_functions + shape_nr, shape_values.size());
824  return shape_values[qpoint * n_shape_functions + shape_nr];
825 }
826 
827 
828 template <int dim, int spacedim>
829 inline const Tensor<1, dim> &
831  const unsigned int qpoint,
832  const unsigned int shape_nr) const
833 {
834  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
835  shape_derivatives.size());
836  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
837 }
838 
839 
840 
841 template <int dim, int spacedim>
842 inline Tensor<1, dim> &
844  const unsigned int qpoint,
845  const unsigned int shape_nr)
846 {
847  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
848  shape_derivatives.size());
849  return shape_derivatives[qpoint * n_shape_functions + shape_nr];
850 }
851 
852 
853 template <int dim, int spacedim>
854 inline const Tensor<2, dim> &
856  const unsigned int qpoint,
857  const unsigned int shape_nr) const
858 {
859  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
860  shape_second_derivatives.size());
861  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
862 }
863 
864 
865 template <int dim, int spacedim>
866 inline Tensor<2, dim> &
868  const unsigned int qpoint,
869  const unsigned int shape_nr)
870 {
871  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
872  shape_second_derivatives.size());
873  return shape_second_derivatives[qpoint * n_shape_functions + shape_nr];
874 }
875 
876 template <int dim, int spacedim>
877 inline const Tensor<3, dim> &
879  const unsigned int qpoint,
880  const unsigned int shape_nr) const
881 {
882  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
883  shape_third_derivatives.size());
884  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
885 }
886 
887 
888 template <int dim, int spacedim>
889 inline Tensor<3, dim> &
891  const unsigned int qpoint,
892  const unsigned int shape_nr)
893 {
894  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
895  shape_third_derivatives.size());
896  return shape_third_derivatives[qpoint * n_shape_functions + shape_nr];
897 }
898 
899 
900 template <int dim, int spacedim>
901 inline const Tensor<4, dim> &
903  const unsigned int qpoint,
904  const unsigned int shape_nr) const
905 {
906  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
907  shape_fourth_derivatives.size());
908  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
909 }
910 
911 
912 template <int dim, int spacedim>
913 inline Tensor<4, dim> &
915  const unsigned int qpoint,
916  const unsigned int shape_nr)
917 {
918  AssertIndexRange(qpoint * n_shape_functions + shape_nr,
919  shape_fourth_derivatives.size());
920  return shape_fourth_derivatives[qpoint * n_shape_functions + shape_nr];
921 }
922 
923 
924 
925 template <int dim, int spacedim>
926 inline bool
928 {
929  return true;
930 }
931 
932 
933 
934 template <int dim, int spacedim>
935 template <typename Number, typename Number2>
936 inline void
938  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
939  const MappingKind kind,
940  const ArrayView<const Point<dim>> & unit_points,
941  const ArrayView<const Number> & input,
942  const ArrayView<Number2> & output) const
943 {
944  AssertDimension(unit_points.size(), output.size());
945  AssertDimension(unit_points.size(), input.size());
946 
947  const std::vector<Point<spacedim>> support_points =
948  this->compute_mapping_support_points(cell);
949 
950  const unsigned int n_points = unit_points.size();
951  const unsigned int n_lanes = VectorizedArray<double>::size();
952 
953  if (kind != mapping_covariant)
954  Assert(false, ExcNotImplemented());
955 
956  // Use the more heavy VectorizedArray code path if there is more than
957  // one point left to compute
958  for (unsigned int i = 0; i < n_points; i += n_lanes)
959  if (n_points - i > 1)
960  {
962  for (unsigned int j = 0; j < n_lanes; ++j)
963  if (i + j < n_points)
964  for (unsigned int d = 0; d < spacedim; ++d)
965  p_vec[d][j] = unit_points[i + j][d];
966  else
967  for (unsigned int d = 0; d < spacedim; ++d)
968  p_vec[d][j] = unit_points[i][d];
969 
973  support_points,
974  p_vec,
975  polynomial_degree == 1,
977  .second;
978 
980  grad.transpose().covariant_form();
981  for (unsigned int j = 0; j < n_lanes && i + j < n_points; ++j)
982  {
984  for (unsigned int d = 0; d < spacedim; ++d)
985  for (unsigned int e = 0; e < dim; ++e)
986  jac_j[d][e] = jac[d][e][j];
987  output[i + j] = apply_transformation(jac_j, input[i + j]);
988  }
989  }
990  else
991  {
995  support_points,
996  unit_points[i],
997  polynomial_degree == 1,
999  .second;
1000  output[i] =
1001  apply_transformation(grad.transpose().covariant_form(), input[i]);
1002  }
1003 }
1004 
1005 #endif // DOXYGEN
1006 
1007 /* -------------- declaration of explicit specializations ------------- */
1008 
1009 
1011 
1012 #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:570
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
virtual bool is_compatible_with(const ReferenceCell &reference_cell) 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
void reference_cell(const ReferenceCell &reference_cell, Triangulation< dim, spacedim > &tria)
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:394
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:393
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
virtual BoundingBox< spacedim > get_bounding_box(const typename Triangulation< dim, spacedim >::cell_iterator &cell) const override
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