42template <
int dim,
int spacedim>
52template <
int dim,
int spacedim>
64template <
int dim,
int spacedim>
73template <
int dim,
int spacedim>
81 ") and the reference cell cell_type (" +
90template <
int dim,
int spacedim>
109template <
int dim,
int spacedim>
110std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
114 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
115 std::make_unique<InternalData>(q);
128template <
int dim,
int spacedim>
129std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
136 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
138 ReferenceCells::get_hypercube<dim>(), quadrature[0]));
148 data.update_each = update_flags;
155template <
int dim,
int spacedim>
156std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
161 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
163 ReferenceCells::get_hypercube<dim>(), quadrature));
173 data.update_each = update_flags;
180template <
int dim,
int spacedim>
216template <
int dim,
int spacedim>
233template <
int dim,
int spacedim>
237 const unsigned int face_no,
246 ReferenceCells::get_hypercube<dim>(),
248 cell->face_orientation(face_no),
249 cell->face_flip(face_no),
250 cell->face_rotation(face_no),
260template <
int dim,
int spacedim>
264 const unsigned int face_no,
265 const unsigned int sub_no,
271 if (cell->face(face_no)->has_children())
279 ReferenceCells::get_hypercube<dim>(),
282 cell->face_orientation(face_no),
283 cell->face_flip(face_no),
284 cell->face_rotation(face_no),
286 cell->subface_case(face_no));
294template <
int dim,
int spacedim>
309 for (
unsigned int d = 0;
d < dim; ++
d)
317template <
int dim,
int spacedim>
320 const unsigned int face_no,
328 std::fill(normal_vectors.begin(),
329 normal_vectors.end(),
336template <
int dim,
int spacedim>
347 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
351 for (
unsigned int i = 0;
357 for (
unsigned int i = 0;
364 for (
unsigned int i = 0;
371 for (
unsigned int i = 0;
378 for (
unsigned int i = 0;
388template <
int dim,
int spacedim>
417 for (
unsigned int d = 1;
d < dim; ++
d)
421 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
428 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
431 for (
unsigned int j = 0; j < dim; ++j)
444 for (
unsigned int j = 0; j < dim; ++j)
448 return cell_similarity;
453template <
int dim,
int spacedim>
457 const unsigned int face_no,
485 for (
unsigned int d = 0;
d < dim; ++
d)
490 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
491 output_data.
JxW_values[i] = J * quadrature[0].weight(i);
494 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
500 for (
unsigned int d = 1;
d < dim; ++
d)
506 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
509 for (
unsigned int d = 0;
d < dim; ++
d)
519 for (
unsigned int d = 0;
d < dim; ++
d)
526template <
int dim,
int spacedim>
530 const unsigned int face_no,
531 const unsigned int subface_no,
553 for (
unsigned int d = 0;
d < dim; ++
d)
563 const unsigned int n_subfaces =
564 cell->face(face_no)->has_children() ?
565 cell->face(face_no)->n_children() :
567 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
572 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
578 for (
unsigned int d = 1;
d < dim; ++
d)
584 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
587 for (
unsigned int d = 0;
d < dim; ++
d)
597 for (
unsigned int d = 0;
d < dim; ++
d)
604template <
int dim,
int spacedim>
617 switch (mapping_kind)
623 "update_covariant_transformation"));
625 for (
unsigned int i = 0; i < output.size(); ++i)
626 for (
unsigned int d = 0;
d < dim; ++
d)
635 "update_contravariant_transformation"));
637 for (
unsigned int i = 0; i < output.size(); ++i)
638 for (
unsigned int d = 0;
d < dim; ++
d)
646 "update_contravariant_transformation"));
649 "update_volume_elements"));
651 for (
unsigned int i = 0; i < output.size(); ++i)
652 for (
unsigned int d = 0;
d < dim; ++
d)
664template <
int dim,
int spacedim>
677 switch (mapping_kind)
683 "update_covariant_transformation"));
685 for (
unsigned int i = 0; i < output.size(); ++i)
686 for (
unsigned int d1 = 0; d1 < dim; ++d1)
687 for (
unsigned int d2 = 0; d2 < dim; ++d2)
688 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
696 "update_contravariant_transformation"));
698 for (
unsigned int i = 0; i < output.size(); ++i)
699 for (
unsigned int d1 = 0; d1 < dim; ++d1)
700 for (
unsigned int d2 = 0; d2 < dim; ++d2)
701 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
709 "update_covariant_transformation"));
711 for (
unsigned int i = 0; i < output.size(); ++i)
712 for (
unsigned int d1 = 0; d1 < dim; ++d1)
713 for (
unsigned int d2 = 0; d2 < dim; ++d2)
714 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
723 "update_contravariant_transformation"));
725 for (
unsigned int i = 0; i < output.size(); ++i)
726 for (
unsigned int d1 = 0; d1 < dim; ++d1)
727 for (
unsigned int d2 = 0; d2 < dim; ++d2)
728 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
737 "update_contravariant_transformation"));
740 "update_volume_elements"));
742 for (
unsigned int i = 0; i < output.size(); ++i)
743 for (
unsigned int d1 = 0; d1 < dim; ++d1)
744 for (
unsigned int d2 = 0; d2 < dim; ++d2)
745 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
754 "update_contravariant_transformation"));
757 "update_volume_elements"));
759 for (
unsigned int i = 0; i < output.size(); ++i)
760 for (
unsigned int d1 = 0; d1 < dim; ++d1)
761 for (
unsigned int d2 = 0; d2 < dim; ++d2)
762 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
774template <
int dim,
int spacedim>
787 switch (mapping_kind)
793 "update_covariant_transformation"));
795 for (
unsigned int i = 0; i < output.size(); ++i)
796 for (
unsigned int d1 = 0; d1 < dim; ++d1)
797 for (
unsigned int d2 = 0; d2 < dim; ++d2)
798 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
806 "update_contravariant_transformation"));
808 for (
unsigned int i = 0; i < output.size(); ++i)
809 for (
unsigned int d1 = 0; d1 < dim; ++d1)
810 for (
unsigned int d2 = 0; d2 < dim; ++d2)
811 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
819 "update_covariant_transformation"));
821 for (
unsigned int i = 0; i < output.size(); ++i)
822 for (
unsigned int d1 = 0; d1 < dim; ++d1)
823 for (
unsigned int d2 = 0; d2 < dim; ++d2)
824 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
833 "update_contravariant_transformation"));
835 for (
unsigned int i = 0; i < output.size(); ++i)
836 for (
unsigned int d1 = 0; d1 < dim; ++d1)
837 for (
unsigned int d2 = 0; d2 < dim; ++d2)
838 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
847 "update_contravariant_transformation"));
850 "update_volume_elements"));
852 for (
unsigned int i = 0; i < output.size(); ++i)
853 for (
unsigned int d1 = 0; d1 < dim; ++d1)
854 for (
unsigned int d2 = 0; d2 < dim; ++d2)
855 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
864 "update_contravariant_transformation"));
867 "update_volume_elements"));
869 for (
unsigned int i = 0; i < output.size(); ++i)
870 for (
unsigned int d1 = 0; d1 < dim; ++d1)
871 for (
unsigned int d2 = 0; d2 < dim; ++d2)
872 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
884template <
int dim,
int spacedim>
897 switch (mapping_kind)
903 "update_covariant_transformation"));
905 for (
unsigned int q = 0; q < output.size(); ++q)
906 for (
unsigned int i = 0; i < spacedim; ++i)
907 for (
unsigned int j = 0; j < spacedim; ++j)
908 for (
unsigned int k = 0; k < spacedim; ++k)
910 output[q][i][j][k] = input[q][i][j][k] /
923template <
int dim,
int spacedim>
936 switch (mapping_kind)
942 "update_covariant_transformation"));
945 "update_contravariant_transformation"));
947 for (
unsigned int q = 0; q < output.size(); ++q)
948 for (
unsigned int i = 0; i < spacedim; ++i)
949 for (
unsigned int j = 0; j < spacedim; ++j)
950 for (
unsigned int k = 0; k < spacedim; ++k)
963 "update_covariant_transformation"));
965 for (
unsigned int q = 0; q < output.size(); ++q)
966 for (
unsigned int i = 0; i < spacedim; ++i)
967 for (
unsigned int j = 0; j < spacedim; ++j)
968 for (
unsigned int k = 0; k < spacedim; ++k)
982 "update_covariant_transformation"));
985 "update_contravariant_transformation"));
988 "update_volume_elements"));
990 for (
unsigned int q = 0; q < output.size(); ++q)
991 for (
unsigned int i = 0; i < spacedim; ++i)
992 for (
unsigned int j = 0; j < spacedim; ++j)
993 for (
unsigned int k = 0; k < spacedim; ++k)
1011template <
int dim,
int spacedim>
1022 length[0] = cell->vertex(1)(0) - start(0);
1025 length[0] = cell->vertex(1)(0) - start(0);
1026 length[1] = cell->vertex(2)(1) - start(1);
1029 length[0] = cell->vertex(1)(0) - start(0);
1030 length[1] = cell->vertex(2)(1) - start(1);
1031 length[2] = cell->vertex(4)(2) - start(2);
1038 for (
unsigned int d = 0;
d < dim; ++
d)
1039 p_real(
d) += length[
d] * p(
d);
1046template <
int dim,
int spacedim>
1052 if (dim != spacedim)
1061 real(0) /= cell->vertex(1)(0) - start(0);
1064 real(0) /= cell->vertex(1)(0) - start(0);
1065 real(1) /= cell->vertex(2)(1) - start(1);
1068 real(0) /= cell->vertex(1)(0) - start(0);
1069 real(1) /= cell->vertex(2)(1) - start(1);
1070 real(2) /= cell->vertex(4)(2) - start(2);
1080template <
int dim,
int spacedim>
1081std::unique_ptr<Mapping<dim, spacedim>>
1084 return std::make_unique<MappingCartesian<dim, spacedim>>(*this);
1090#include "mapping_cartesian.inst"
std::vector< Point< dim > > quadrature_points
Tensor< 1, dim > cell_extents
InternalData(const Quadrature< dim > &quadrature)
virtual std::size_t memory_consumption() 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 typename Mapping< dim, spacedim >::InternalDataBase &internal_data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const hp::QCollection< dim - 1 > &quadrature) const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim - 1 > &quadrature) const override
void update_cell_extents(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const InternalData &data) const
void maybe_update_subface_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
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
void maybe_update_face_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
virtual bool is_compatible_with(const ReferenceCell &reference_cell) const override
void maybe_update_cell_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, std::vector< Point< dim > > &quadrature_points) const
virtual bool preserves_vertex_locations() const override
virtual Point< spacedim > transform_unit_to_real_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< dim > &p) const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
void transform_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, const typename QProjector< dim >::DataSetDescriptor &offset, std::vector< Point< dim > > &quadrature_points) const
void maybe_update_jacobian_derivatives(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
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
void maybe_update_normal_vectors(const unsigned int face_no, const InternalData &data, std::vector< Tensor< 1, dim > > &normal_vectors) const
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
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
Abstract base class for mapping classes.
static DataSetDescriptor cell()
static DataSetDescriptor subface(const unsigned int face_no, const unsigned int subface_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points, const internal::SubfaceCase< dim > ref_case=internal::SubfaceCase< dim >::case_isotropic)
static DataSetDescriptor face(const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
double weight(const unsigned int i) const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_volume_elements
Determinant of the Jacobian.
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_jacobian_3rd_derivatives
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
@ mapping_covariant_gradient
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim > > > &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double > > &properties={})
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)