42 "You are using MappingCartesian, but the incoming cell is not Cartesian.");
51template <
class CellType>
55 if (!cell->reference_cell().is_hyper_cube())
58 const double tolerance = 1
e-14 * cell->diameter();
59 const auto bounding_box = cell->bounding_box();
61 for (
const unsigned int v : cell->vertex_indices())
62 if (cell->vertex(v).distance(bounding_box.vertex(v)) > tolerance)
70template <
int dim,
int spacedim>
80template <
int dim,
int spacedim>
92template <
int dim,
int spacedim>
101template <
int dim,
int spacedim>
107 ExcMessage(
"The dimension of your mapping (" +
109 ") and the reference cell cell_type (" +
111 " ) do not agree."));
118template <
int dim,
int spacedim>
137template <
int dim,
int spacedim>
138std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
142 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
143 std::make_unique<InternalData>(q);
156template <
int dim,
int spacedim>
157std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
164 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
166 ReferenceCells::get_hypercube<dim>(), quadrature[0]));
176 data.update_each = update_flags;
183template <
int dim,
int spacedim>
184std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
189 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
191 ReferenceCells::get_hypercube<dim>(), quadrature));
201 data.update_each = update_flags;
208template <
int dim,
int spacedim>
244template <
int dim,
int spacedim>
261template <
int dim,
int spacedim>
265 const unsigned int face_no,
274 ReferenceCells::get_hypercube<dim>(),
276 cell->face_orientation(face_no),
277 cell->face_flip(face_no),
278 cell->face_rotation(face_no),
288template <
int dim,
int spacedim>
292 const unsigned int face_no,
293 const unsigned int sub_no,
299 if (cell->face(face_no)->has_children())
307 ReferenceCells::get_hypercube<dim>(),
310 cell->face_orientation(face_no),
311 cell->face_flip(face_no),
312 cell->face_rotation(face_no),
314 cell->subface_case(face_no));
322template <
int dim,
int spacedim>
337 for (
unsigned int d = 0;
d < dim; ++
d)
345template <
int dim,
int spacedim>
348 const unsigned int face_no,
356 std::fill(normal_vectors.begin(),
357 normal_vectors.end(),
364template <
int dim,
int spacedim>
375 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
379 for (
unsigned int i = 0;
385 for (
unsigned int i = 0;
392 for (
unsigned int i = 0;
399 for (
unsigned int i = 0;
406 for (
unsigned int i = 0;
416template <
int dim,
int spacedim>
424 for (
unsigned int d = 1;
d < dim; ++
d)
432template <
int dim,
int spacedim>
444 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
447 for (
unsigned int j = 0; j < dim; ++j)
454template <
int dim,
int spacedim>
469 for (
unsigned int j = 0; j < dim; ++j)
476template <
int dim,
int spacedim>
507 for (
unsigned int d = 1;
d < dim; ++
d)
511 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
520 return cell_similarity;
525template <
int dim,
int spacedim>
544 output_data.
initialize(unit_points.size(), update_flags);
547 data.update_each = update_flags;
549 std::vector<Point<dim>>(unit_points.begin(), unit_points.end());
563template <
int dim,
int spacedim>
567 const unsigned int face_no,
596 for (
unsigned int d = 0;
d < dim; ++
d)
601 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
602 output_data.
JxW_values[i] = J * quadrature[0].weight(i);
605 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
616template <
int dim,
int spacedim>
620 const unsigned int face_no,
621 const unsigned int subface_no,
645 for (
unsigned int d = 0;
d < dim; ++
d)
655 const unsigned int n_subfaces =
656 cell->face(face_no)->has_children() ?
657 cell->face(face_no)->n_children() :
659 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
664 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
675template <
int dim,
int spacedim>
701 for (
unsigned int i = 0; i < output_data.
normal_vectors.size(); ++i)
706 for (
unsigned int d = 0;
d < dim; ++
d)
710 normal /= normal.
norm();
715 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
721 double det_jacobian = 1.;
722 for (
unsigned int d = 0;
d < dim; ++
d)
728 det_jacobian * invJTxNormal.
norm() * quadrature.
weight(i);
739template <
int dim,
int spacedim>
752 switch (mapping_kind)
758 "update_covariant_transformation"));
760 for (
unsigned int i = 0; i < output.size(); ++i)
761 for (
unsigned int d = 0;
d < dim; ++
d)
770 "update_contravariant_transformation"));
772 for (
unsigned int i = 0; i < output.size(); ++i)
773 for (
unsigned int d = 0;
d < dim; ++
d)
781 "update_contravariant_transformation"));
784 "update_volume_elements"));
786 for (
unsigned int i = 0; i < output.size(); ++i)
787 for (
unsigned int d = 0;
d < dim; ++
d)
799template <
int dim,
int spacedim>
812 switch (mapping_kind)
818 "update_covariant_transformation"));
820 for (
unsigned int i = 0; i < output.size(); ++i)
821 for (
unsigned int d1 = 0; d1 < dim; ++d1)
822 for (
unsigned int d2 = 0; d2 < dim; ++d2)
823 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
831 "update_contravariant_transformation"));
833 for (
unsigned int i = 0; i < output.size(); ++i)
834 for (
unsigned int d1 = 0; d1 < dim; ++d1)
835 for (
unsigned int d2 = 0; d2 < dim; ++d2)
836 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
844 "update_covariant_transformation"));
846 for (
unsigned int i = 0; i < output.size(); ++i)
847 for (
unsigned int d1 = 0; d1 < dim; ++d1)
848 for (
unsigned int d2 = 0; d2 < dim; ++d2)
849 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
858 "update_contravariant_transformation"));
860 for (
unsigned int i = 0; i < output.size(); ++i)
861 for (
unsigned int d1 = 0; d1 < dim; ++d1)
862 for (
unsigned int d2 = 0; d2 < dim; ++d2)
863 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
872 "update_contravariant_transformation"));
875 "update_volume_elements"));
877 for (
unsigned int i = 0; i < output.size(); ++i)
878 for (
unsigned int d1 = 0; d1 < dim; ++d1)
879 for (
unsigned int d2 = 0; d2 < dim; ++d2)
880 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
889 "update_contravariant_transformation"));
892 "update_volume_elements"));
894 for (
unsigned int i = 0; i < output.size(); ++i)
895 for (
unsigned int d1 = 0; d1 < dim; ++d1)
896 for (
unsigned int d2 = 0; d2 < dim; ++d2)
897 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
909template <
int dim,
int spacedim>
922 switch (mapping_kind)
928 "update_covariant_transformation"));
930 for (
unsigned int i = 0; i < output.size(); ++i)
931 for (
unsigned int d1 = 0; d1 < dim; ++d1)
932 for (
unsigned int d2 = 0; d2 < dim; ++d2)
933 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
941 "update_contravariant_transformation"));
943 for (
unsigned int i = 0; i < output.size(); ++i)
944 for (
unsigned int d1 = 0; d1 < dim; ++d1)
945 for (
unsigned int d2 = 0; d2 < dim; ++d2)
946 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
954 "update_covariant_transformation"));
956 for (
unsigned int i = 0; i < output.size(); ++i)
957 for (
unsigned int d1 = 0; d1 < dim; ++d1)
958 for (
unsigned int d2 = 0; d2 < dim; ++d2)
959 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
968 "update_contravariant_transformation"));
970 for (
unsigned int i = 0; i < output.size(); ++i)
971 for (
unsigned int d1 = 0; d1 < dim; ++d1)
972 for (
unsigned int d2 = 0; d2 < dim; ++d2)
973 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
982 "update_contravariant_transformation"));
985 "update_volume_elements"));
987 for (
unsigned int i = 0; i < output.size(); ++i)
988 for (
unsigned int d1 = 0; d1 < dim; ++d1)
989 for (
unsigned int d2 = 0; d2 < dim; ++d2)
990 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
999 "update_contravariant_transformation"));
1002 "update_volume_elements"));
1004 for (
unsigned int i = 0; i < output.size(); ++i)
1005 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1006 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1007 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
1019template <
int dim,
int spacedim>
1032 switch (mapping_kind)
1038 "update_covariant_transformation"));
1040 for (
unsigned int q = 0; q < output.size(); ++q)
1041 for (
unsigned int i = 0; i < spacedim; ++i)
1042 for (
unsigned int j = 0; j < spacedim; ++j)
1043 for (
unsigned int k = 0; k < spacedim; ++k)
1045 output[q][i][j][k] = input[q][i][j][k] /
1058template <
int dim,
int spacedim>
1071 switch (mapping_kind)
1077 "update_covariant_transformation"));
1080 "update_contravariant_transformation"));
1082 for (
unsigned int q = 0; q < output.size(); ++q)
1083 for (
unsigned int i = 0; i < spacedim; ++i)
1084 for (
unsigned int j = 0; j < spacedim; ++j)
1085 for (
unsigned int k = 0; k < spacedim; ++k)
1087 output[q][i][j][k] =
1098 "update_covariant_transformation"));
1100 for (
unsigned int q = 0; q < output.size(); ++q)
1101 for (
unsigned int i = 0; i < spacedim; ++i)
1102 for (
unsigned int j = 0; j < spacedim; ++j)
1103 for (
unsigned int k = 0; k < spacedim; ++k)
1105 output[q][i][j][k] =
1117 "update_covariant_transformation"));
1120 "update_contravariant_transformation"));
1123 "update_volume_elements"));
1125 for (
unsigned int q = 0; q < output.size(); ++q)
1126 for (
unsigned int i = 0; i < spacedim; ++i)
1127 for (
unsigned int j = 0; j < spacedim; ++j)
1128 for (
unsigned int k = 0; k < spacedim; ++k)
1130 output[q][i][j][k] =
1146template <
int dim,
int spacedim>
1159 length[0] = cell->vertex(1)(0) - start(0);
1162 length[0] = cell->vertex(1)(0) - start(0);
1163 length[1] = cell->vertex(2)(1) - start(1);
1166 length[0] = cell->vertex(1)(0) - start(0);
1167 length[1] = cell->vertex(2)(1) - start(1);
1168 length[2] = cell->vertex(4)(2) - start(2);
1175 for (
unsigned int d = 0;
d < dim; ++
d)
1176 p_real(
d) += length[
d] * p(
d);
1183template <
int dim,
int spacedim>
1191 if (dim != spacedim)
1200 real(0) /= cell->vertex(1)(0) - start(0);
1203 real(0) /= cell->vertex(1)(0) - start(0);
1204 real(1) /= cell->vertex(2)(1) - start(1);
1207 real(0) /= cell->vertex(1)(0) - start(0);
1208 real(1) /= cell->vertex(2)(1) - start(1);
1209 real(2) /= cell->vertex(4)(2) - start(2);
1219template <
int dim,
int spacedim>
1220std::unique_ptr<Mapping<dim, spacedim>>
1223 return std::make_unique<MappingCartesian<dim, spacedim>>(*this);
1229#include "mapping_cartesian.inst"
std::vector< Point< dim > > quadrature_points
Tensor< 1, dim > cell_extents
virtual std::size_t memory_consumption() const override
void fill_mapping_data_for_generic_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const ArrayView< const Point< dim > > &unit_points, const UpdateFlags update_flags, ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) const
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
void maybe_update_volume_elements(const InternalData &data) const
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_jacobians(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_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
virtual void fill_fe_immersed_surface_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const NonMatching::ImmersedSurfaceQuadrature< 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
void maybe_update_inverse_jacobians(const InternalData &data, const CellSimilarity::Similarity cell_similarity, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data) 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
Abstract base class for mapping classes.
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
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
numbers::NumberTraits< Number >::real_type norm() 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_default
No update.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
static ::ExceptionBase & ExcCellNotCartesian()
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertIndexRange(index, range)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
@ mapping_covariant_gradient
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
bool is_cartesian(const CellType &cell)
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 > e(const Tensor< 2, dim, Number > &F)
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)