41 "You are using MappingCartesian, but the incoming cell is not Cartesian.");
50template <
typename CellType>
54 if (!cell->reference_cell().is_hyper_cube())
60 const double abs_tol = 1e-30;
61 const double rel_tol = 1e-28;
62 const auto bounding_box = cell->bounding_box();
63 const auto &bounding_vertices = bounding_box.get_boundary_points();
64 const auto bb_diagonal_length_squared =
65 bounding_vertices.first.distance_square(bounding_vertices.second);
67 for (
const unsigned int v : cell->vertex_indices())
78 const double tolerance =
std::max(abs_tol * cell->vertex(v).norm_square(),
79 rel_tol * bb_diagonal_length_squared);
81 if (cell->vertex(v).distance_square(bounding_box.vertex(v)) > tolerance)
90template <
int dim,
int spacedim>
94 , inverse_cell_extents(
numbers::signaling_nan<
Tensor<1, dim>>())
95 , volume_element(
numbers::signaling_nan<double>())
96 , quadrature_points(q.get_points())
101template <
int dim,
int spacedim>
110 this->update_each = update_flags;
115template <
int dim,
int spacedim>
127template <
int dim,
int spacedim>
136template <
int dim,
int spacedim>
141 Assert(dim == reference_cell.get_dimension(),
142 ExcMessage(
"The dimension of your mapping (" +
144 ") and the reference cell cell_type (" +
146 " ) do not agree."));
148 return reference_cell.is_hyper_cube();
153template <
int dim,
int spacedim>
172template <
int dim,
int spacedim>
173std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
177 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
178 std::make_unique<InternalData>();
186template <
int dim,
int spacedim>
187std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
194 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
196 ReferenceCells::get_hypercube<dim>(), quadrature[0]));
206 data.update_each = update_flags;
213template <
int dim,
int spacedim>
214std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
219 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase> data_ptr =
221 ReferenceCells::get_hypercube<dim>(), quadrature));
231 data.update_each = update_flags;
238template <
int dim,
int spacedim>
250 for (
unsigned int d = 0; d < dim; ++d)
252 const double cell_extent_d = cell->vertex(1 << d)[d] - start[d];
253 data.cell_extents[d] = cell_extent_d;
254 Assert(cell_extent_d != 0.,
255 ExcMessage(
"Cell does not appear to be Cartesian!"));
256 data.inverse_cell_extents[d] = 1. / cell_extent_d;
267 transform_quadrature_points(
274 for (
unsigned int i = 0; i < quadrature_points.size(); ++i)
276 quadrature_points[i] = first_vertex;
277 for (
unsigned int d = 0; d < dim; ++d)
278 quadrature_points[i][d] +=
279 cell_extents[d] * unit_quadrature_points[i + offset][d];
286template <
int dim,
int spacedim>
292 std::vector<
Point<dim>> &quadrature_points)
const
298 transform_quadrature_points(cell->vertex(0),
300 unit_quadrature_points,
308template <
int dim,
int spacedim>
312 const unsigned int face_no,
314 std::vector<
Point<dim>> &quadrature_points)
const
321 ReferenceCells::get_hypercube<dim>(),
323 cell->combined_face_orientation(face_no),
324 quadrature_points.size());
327 transform_quadrature_points(cell->vertex(0),
337template <
int dim,
int spacedim>
341 const unsigned int face_no,
342 const unsigned int sub_no,
344 std::vector<
Point<dim>> &quadrature_points)
const
348 if (cell->face(face_no)->has_children())
356 ReferenceCells::get_hypercube<dim>(),
359 cell->combined_face_orientation(face_no),
360 quadrature_points.size(),
361 cell->subface_case(face_no));
363 transform_quadrature_points(cell->vertex(0),
373template <
int dim,
int spacedim>
376 const unsigned int face_no,
384 std::fill(normal_vectors.begin(),
385 normal_vectors.end(),
386 ReferenceCells::get_hypercube<dim>()
387 .template face_normal_vector<dim>(face_no));
393template <
int dim,
int spacedim>
404 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
408 for (
unsigned int i = 0;
414 for (
unsigned int i = 0;
421 for (
unsigned int i = 0;
428 for (
unsigned int i = 0;
435 for (
unsigned int i = 0;
445template <
int dim,
int spacedim>
452 double volume =
data.cell_extents[0];
453 for (
unsigned int d = 1; d < dim; ++d)
454 volume *=
data.cell_extents[d];
455 data.volume_element = volume;
461template <
int dim,
int spacedim>
473 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
476 for (
unsigned int j = 0; j < dim; ++j)
483template <
int dim,
int spacedim>
498 for (
unsigned int j = 0; j < dim; ++j)
500 data.inverse_cell_extents[j];
506template <
int dim,
int spacedim>
537 double J =
data.cell_extents[0];
538 for (
unsigned int d = 1; d < dim; ++d)
539 J *=
data.cell_extents[d];
540 data.volume_element = J;
542 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
551 return cell_similarity;
556template <
int dim,
int spacedim>
575 output_data.
initialize(unit_points.size(), update_flags);
578 data.update_each = update_flags;
593template <
int dim,
int spacedim>
597 const unsigned int face_no,
626 for (
unsigned int d = 0; d < dim; ++d)
628 J *=
data.cell_extents[d];
631 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
632 output_data.
JxW_values[i] = J * quadrature[0].weight(i);
635 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
646template <
int dim,
int spacedim>
650 const unsigned int face_no,
651 const unsigned int subface_no,
675 for (
unsigned int d = 0; d < dim; ++d)
677 J *=
data.cell_extents[d];
685 const unsigned int n_subfaces =
686 cell->face(face_no)->has_children() ?
687 cell->face(face_no)->n_children() :
689 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
694 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
705template <
int dim,
int spacedim>
728 quadrature.get_points(),
732 for (
unsigned int i = 0; i < output_data.
normal_vectors.size(); ++i)
737 for (
unsigned int d = 0; d < dim; ++d)
739 normal[d] = ref_space_normal[d] *
data.inverse_cell_extents[d];
741 normal /= normal.
norm();
746 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
752 double det_jacobian = 1.;
753 for (
unsigned int d = 0; d < dim; ++d)
755 det_jacobian *=
data.cell_extents[d];
757 ref_space_normal[d] *
data.inverse_cell_extents[d];
760 det_jacobian * invJTxNormal.
norm() * quadrature.weight(i);
771template <
int dim,
int spacedim>
784 switch (mapping_kind)
790 "update_covariant_transformation"));
792 for (
unsigned int i = 0; i < output.size(); ++i)
793 for (
unsigned int d = 0; d < dim; ++d)
794 output[i][d] = input[i][d] *
data.inverse_cell_extents[d];
802 "update_contravariant_transformation"));
804 for (
unsigned int i = 0; i < output.size(); ++i)
805 for (
unsigned int d = 0; d < dim; ++d)
806 output[i][d] = input[i][d] *
data.cell_extents[d];
813 "update_contravariant_transformation"));
816 "update_volume_elements"));
818 for (
unsigned int i = 0; i < output.size(); ++i)
819 for (
unsigned int d = 0; d < dim; ++d)
821 input[i][d] *
data.cell_extents[d] /
data.volume_element;
831template <
int dim,
int spacedim>
844 switch (mapping_kind)
850 "update_covariant_transformation"));
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)
856 input[i][d1][d2] *
data.inverse_cell_extents[d2];
864 "update_contravariant_transformation"));
866 for (
unsigned int i = 0; i < output.size(); ++i)
867 for (
unsigned int d1 = 0; d1 < dim; ++d1)
868 for (
unsigned int d2 = 0; d2 < dim; ++d2)
869 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2];
877 "update_covariant_transformation"));
879 for (
unsigned int i = 0; i < output.size(); ++i)
880 for (
unsigned int d1 = 0; d1 < dim; ++d1)
881 for (
unsigned int d2 = 0; d2 < dim; ++d2)
882 output[i][d1][d2] = input[i][d1][d2] *
883 data.inverse_cell_extents[d2] *
884 data.inverse_cell_extents[d1];
892 "update_contravariant_transformation"));
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] *
898 data.inverse_cell_extents[d1];
906 "update_contravariant_transformation"));
909 "update_volume_elements"));
911 for (
unsigned int i = 0; i < output.size(); ++i)
912 for (
unsigned int d1 = 0; d1 < dim; ++d1)
913 for (
unsigned int d2 = 0; d2 < dim; ++d2)
914 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2] /
923 "update_contravariant_transformation"));
926 "update_volume_elements"));
928 for (
unsigned int i = 0; i < output.size(); ++i)
929 for (
unsigned int d1 = 0; d1 < dim; ++d1)
930 for (
unsigned int d2 = 0; d2 < dim; ++d2)
931 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2] *
932 data.inverse_cell_extents[d1] /
944template <
int dim,
int spacedim>
957 switch (mapping_kind)
963 "update_covariant_transformation"));
965 for (
unsigned int i = 0; i < output.size(); ++i)
966 for (
unsigned int d1 = 0; d1 < dim; ++d1)
967 for (
unsigned int d2 = 0; d2 < dim; ++d2)
969 input[i][d1][d2] *
data.inverse_cell_extents[d2];
977 "update_contravariant_transformation"));
979 for (
unsigned int i = 0; i < output.size(); ++i)
980 for (
unsigned int d1 = 0; d1 < dim; ++d1)
981 for (
unsigned int d2 = 0; d2 < dim; ++d2)
982 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2];
990 "update_covariant_transformation"));
992 for (
unsigned int i = 0; i < output.size(); ++i)
993 for (
unsigned int d1 = 0; d1 < dim; ++d1)
994 for (
unsigned int d2 = 0; d2 < dim; ++d2)
995 output[i][d1][d2] = input[i][d1][d2] *
996 data.inverse_cell_extents[d2] *
997 data.inverse_cell_extents[d1];
1005 "update_contravariant_transformation"));
1007 for (
unsigned int i = 0; i < output.size(); ++i)
1008 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1009 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1010 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2] *
1011 data.inverse_cell_extents[d1];
1019 "update_contravariant_transformation"));
1022 "update_volume_elements"));
1024 for (
unsigned int i = 0; i < output.size(); ++i)
1025 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1026 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1027 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2] /
1028 data.volume_element;
1036 "update_contravariant_transformation"));
1039 "update_volume_elements"));
1041 for (
unsigned int i = 0; i < output.size(); ++i)
1042 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1043 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1044 output[i][d1][d2] = input[i][d1][d2] *
data.cell_extents[d2] *
1045 data.inverse_cell_extents[d1] /
1046 data.volume_element;
1057template <
int dim,
int spacedim>
1070 switch (mapping_kind)
1076 "update_covariant_transformation"));
1078 for (
unsigned int q = 0; q < output.size(); ++q)
1079 for (
unsigned int i = 0; i < spacedim; ++i)
1080 for (
unsigned int j = 0; j < spacedim; ++j)
1081 for (
unsigned int k = 0; k < spacedim; ++k)
1083 output[q][i][j][k] = input[q][i][j][k] *
1084 data.inverse_cell_extents[j] *
1085 data.inverse_cell_extents[k];
1096template <
int dim,
int spacedim>
1109 switch (mapping_kind)
1115 "update_covariant_transformation"));
1118 "update_contravariant_transformation"));
1120 for (
unsigned int q = 0; q < output.size(); ++q)
1121 for (
unsigned int i = 0; i < spacedim; ++i)
1122 for (
unsigned int j = 0; j < spacedim; ++j)
1123 for (
unsigned int k = 0; k < spacedim; ++k)
1125 output[q][i][j][k] = input[q][i][j][k] *
1126 data.cell_extents[i] *
1127 data.inverse_cell_extents[j] *
1128 data.inverse_cell_extents[k];
1137 "update_covariant_transformation"));
1139 for (
unsigned int q = 0; q < output.size(); ++q)
1140 for (
unsigned int i = 0; i < spacedim; ++i)
1141 for (
unsigned int j = 0; j < spacedim; ++j)
1142 for (
unsigned int k = 0; k < spacedim; ++k)
1144 output[q][i][j][k] = input[q][i][j][k] *
1145 (
data.inverse_cell_extents[i] *
1146 data.inverse_cell_extents[j]) *
1147 data.inverse_cell_extents[k];
1157 "update_covariant_transformation"));
1160 "update_contravariant_transformation"));
1163 "update_volume_elements"));
1165 for (
unsigned int q = 0; q < output.size(); ++q)
1166 for (
unsigned int i = 0; i < spacedim; ++i)
1167 for (
unsigned int j = 0; j < spacedim; ++j)
1168 for (
unsigned int k = 0; k < spacedim; ++k)
1170 output[q][i][j][k] =
1172 (
data.cell_extents[i] /
data.volume_element *
1173 data.inverse_cell_extents[j]) *
1174 data.inverse_cell_extents[k];
1187template <
int dim,
int spacedim>
1199 for (
unsigned int d = 0; d < dim; ++d)
1200 unit[d] += (cell->vertex(1 << d)[d] - unit[d]) * p[d];
1207template <
int dim,
int spacedim>
1220 for (
unsigned int d = 0; d < dim; ++d)
1221 real[d] = (real[d] - start[d]) / (cell->vertex(1 << d)[d] - start[d]);
1228template <
int dim,
int spacedim>
1238 if (dim != spacedim)
1244 std::array<double, dim> inverse_lengths;
1245 for (
unsigned int d = 0; d < dim; ++d)
1246 inverse_lengths[d] = 1. / (cell->vertex(1 << d)[d] - start[d]);
1248 for (
unsigned int i = 0; i < real_points.size(); ++i)
1249 for (
unsigned int d = 0; d < dim; ++d)
1250 unit_points[i][d] = (real_points[i][d] - start[d]) * inverse_lengths[d];
1255template <
int dim,
int spacedim>
1256std::unique_ptr<Mapping<dim, spacedim>>
1259 return std::make_unique<MappingCartesian<dim, spacedim>>(*this);
1265#include "fe/mapping_cartesian.inst"
ArrayView< std::remove_reference_t< typename std::iterator_traits< Iterator >::reference >, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
virtual void reinit(const UpdateFlags update_flags, const Quadrature< dim > &quadrature) override
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
void maybe_update_cell_quadrature_points(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const InternalData &data, const ArrayView< const Point< dim > > &unit_quadrature_points, std::vector< Point< dim > > &quadrature_points) const
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
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_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
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
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
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 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
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
Class storing the offset index into a Quadrature rule created by project_to_all_faces() or project_to...
static DataSetDescriptor subface(const ReferenceCell &reference_cell, 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 cell()
static DataSetDescriptor face(const ReferenceCell &reference_cell, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
Class which transforms dim - 1-dimensional quadrature rules to dim-dimensional face quadratures.
double weight(const unsigned int i) const
const std::vector< Point< dim > > & get_points() const
numbers::NumberTraits< Number >::real_type norm() const
unsigned int size() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
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)
@ 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
@ mapping_covariant_gradient
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
bool is_cartesian(const CellType &cell)
std::vector< index_type > data
std::enable_if_t< std::is_fundamental_v< T >, std::size_t > memory_consumption(const T &t)
std::string to_string(const number value, const unsigned int digits=numbers::invalid_unsigned_int)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)