16 #include <deal.II/base/array_view.h> 17 #include <deal.II/base/memory_consumption.h> 18 #include <deal.II/base/qprojector.h> 19 #include <deal.II/base/quadrature.h> 20 #include <deal.II/base/signaling_nan.h> 21 #include <deal.II/base/std_cxx14/memory.h> 22 #include <deal.II/base/tensor.h> 24 #include <deal.II/dofs/dof_accessor.h> 26 #include <deal.II/fe/fe_values.h> 27 #include <deal.II/fe/mapping_cartesian.h> 29 #include <deal.II/grid/tria.h> 30 #include <deal.II/grid/tria_iterator.h> 32 #include <deal.II/lac/full_matrix.h> 38 DEAL_II_NAMESPACE_OPEN
41 template <
int dim,
int spacedim>
46 template <
int dim,
int spacedim>
50 , volume_element(
numbers::signaling_nan<double>())
51 , quadrature_points(q.get_points())
56 template <
int dim,
int spacedim>
68 template <
int dim,
int spacedim>
77 template <
int dim,
int spacedim>
96 template <
int dim,
int spacedim>
97 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
101 auto data = std_cxx14::make_unique<InternalData>(q);
108 return std::move(data);
113 template <
int dim,
int spacedim>
114 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
119 auto data = std_cxx14::make_unique<InternalData>(
129 data->update_each = update_flags;
131 return std::move(data);
136 template <
int dim,
int spacedim>
137 std::unique_ptr<typename Mapping<dim, spacedim>::InternalDataBase>
142 auto data = std_cxx14::make_unique<InternalData>(
152 data->update_each = update_flags;
154 return std::move(data);
159 template <
int dim,
int spacedim>
163 const unsigned int face_no,
164 const unsigned int sub_no,
191 cell->face(face_no)->has_children() ||
197 !cell->face(face_no)->has_children() ||
198 (sub_no < cell->face(face_no)->n_children()),
199 ExcIndexRange(sub_no, 0, cell->face(face_no)->n_children()));
252 cell->face_orientation(face_no),
253 cell->face_flip(face_no),
254 cell->face_rotation(face_no),
260 cell->face_orientation(face_no),
261 cell->face_flip(face_no),
262 cell->face_rotation(face_no),
264 cell->subface_case(face_no))));
269 for (
unsigned int d = 0; d < dim; ++d)
293 std::fill(normal_vectors.begin(),
294 normal_vectors.end(),
306 std::fill(normal_vectors.begin(),
307 normal_vectors.end(),
321 std::fill(normal_vectors.begin(),
322 normal_vectors.end(),
335 template <
int dim,
int spacedim>
347 Assert(dynamic_cast<const InternalData *>(&internal_data) !=
nullptr,
351 std::vector<Tensor<1, dim>> dummy;
367 for (
unsigned int d = 1; d < dim; ++d)
371 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
378 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
381 for (
unsigned int j = 0; j < dim; ++j)
388 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
393 for (
unsigned int i = 0;
409 for (
unsigned int i = 0;
424 for (
unsigned int i = 0;
437 for (
unsigned int j = 0; j < dim; ++j)
441 return cell_similarity;
446 template <
int dim,
int spacedim>
450 const unsigned int face_no,
460 Assert(dynamic_cast<const InternalData *>(&internal_data) !=
nullptr,
475 for (
unsigned int d = 0; d < dim; ++d)
480 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
484 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
490 for (
unsigned int d = 1; d < dim; ++d)
496 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
499 for (
unsigned int d = 0; d < dim; ++d)
504 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
508 for (
unsigned int i = 0;
520 for (
unsigned int i = 0;
533 for (
unsigned int i = 0;
543 for (
unsigned int d = 0; d < dim; ++d)
550 template <
int dim,
int spacedim>
554 const unsigned int face_no,
555 const unsigned int subface_no,
563 Assert(dynamic_cast<const InternalData *>(&internal_data) !=
nullptr,
578 for (
unsigned int d = 0; d < dim; ++d)
588 const unsigned int n_subfaces =
589 cell->face(face_no)->has_children() ?
590 cell->face(face_no)->n_children() :
592 for (
unsigned int i = 0; i < output_data.
JxW_values.size(); ++i)
597 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
603 for (
unsigned int d = 1; d < dim; ++d)
609 for (
unsigned int i = 0; i < output_data.
jacobians.size(); ++i)
612 for (
unsigned int d = 0; d < dim; ++d)
617 for (
unsigned int i = 0; i < output_data.
jacobian_grads.size(); ++i)
621 for (
unsigned int i = 0;
633 for (
unsigned int i = 0;
646 for (
unsigned int i = 0;
656 for (
unsigned int d = 0; d < dim; ++d)
663 template <
int dim,
int spacedim>
672 Assert(dynamic_cast<const InternalData *>(&mapping_data) !=
nullptr,
676 switch (mapping_type)
682 "update_covariant_transformation"));
684 for (
unsigned int i = 0; i < output.size(); ++i)
685 for (
unsigned int d = 0; d < dim; ++d)
694 "update_contravariant_transformation"));
696 for (
unsigned int i = 0; i < output.size(); ++i)
697 for (
unsigned int d = 0; d < dim; ++d)
705 "update_contravariant_transformation"));
708 "update_volume_elements"));
710 for (
unsigned int i = 0; i < output.size(); ++i)
711 for (
unsigned int d = 0; d < dim; ++d)
723 template <
int dim,
int spacedim>
732 Assert(dynamic_cast<const InternalData *>(&mapping_data) !=
nullptr,
736 switch (mapping_type)
742 "update_covariant_transformation"));
744 for (
unsigned int i = 0; i < output.size(); ++i)
745 for (
unsigned int d1 = 0; d1 < dim; ++d1)
746 for (
unsigned int d2 = 0; d2 < dim; ++d2)
747 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
755 "update_contravariant_transformation"));
757 for (
unsigned int i = 0; i < output.size(); ++i)
758 for (
unsigned int d1 = 0; d1 < dim; ++d1)
759 for (
unsigned int d2 = 0; d2 < dim; ++d2)
760 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
768 "update_covariant_transformation"));
770 for (
unsigned int i = 0; i < output.size(); ++i)
771 for (
unsigned int d1 = 0; d1 < dim; ++d1)
772 for (
unsigned int d2 = 0; d2 < dim; ++d2)
773 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
782 "update_contravariant_transformation"));
784 for (
unsigned int i = 0; i < output.size(); ++i)
785 for (
unsigned int d1 = 0; d1 < dim; ++d1)
786 for (
unsigned int d2 = 0; d2 < dim; ++d2)
787 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
796 "update_contravariant_transformation"));
799 "update_volume_elements"));
801 for (
unsigned int i = 0; i < output.size(); ++i)
802 for (
unsigned int d1 = 0; d1 < dim; ++d1)
803 for (
unsigned int d2 = 0; d2 < dim; ++d2)
804 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
813 "update_contravariant_transformation"));
816 "update_volume_elements"));
818 for (
unsigned int i = 0; i < output.size(); ++i)
819 for (
unsigned int d1 = 0; d1 < dim; ++d1)
820 for (
unsigned int d2 = 0; d2 < dim; ++d2)
821 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2] /
833 template <
int dim,
int spacedim>
842 Assert(dynamic_cast<const InternalData *>(&mapping_data) !=
nullptr,
846 switch (mapping_type)
852 "update_covariant_transformation"));
854 for (
unsigned int i = 0; i < output.size(); ++i)
855 for (
unsigned int d1 = 0; d1 < dim; ++d1)
856 for (
unsigned int d2 = 0; d2 < dim; ++d2)
857 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2];
865 "update_contravariant_transformation"));
867 for (
unsigned int i = 0; i < output.size(); ++i)
868 for (
unsigned int d1 = 0; d1 < dim; ++d1)
869 for (
unsigned int d2 = 0; d2 < dim; ++d2)
870 output[i][d1][d2] = input[i][d1][d2] * data.
cell_extents[d2];
878 "update_covariant_transformation"));
880 for (
unsigned int i = 0; i < output.size(); ++i)
881 for (
unsigned int d1 = 0; d1 < dim; ++d1)
882 for (
unsigned int d2 = 0; d2 < dim; ++d2)
883 output[i][d1][d2] = input[i][d1][d2] / data.
cell_extents[d2] /
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] /
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] /
942 template <
int dim,
int spacedim>
951 Assert(dynamic_cast<const InternalData *>(&mapping_data) !=
nullptr,
955 switch (mapping_type)
961 "update_covariant_transformation"));
963 for (
unsigned int q = 0; q < output.size(); ++q)
964 for (
unsigned int i = 0; i < spacedim; ++i)
965 for (
unsigned int j = 0; j < spacedim; ++j)
966 for (
unsigned int k = 0; k < spacedim; ++k)
968 output[q][i][j][k] = input[q][i][j][k] /
979 template <
int dim,
int spacedim>
988 Assert(dynamic_cast<const InternalData *>(&mapping_data) !=
nullptr,
992 switch (mapping_type)
998 "update_covariant_transformation"));
1001 "update_contravariant_transformation"));
1003 for (
unsigned int q = 0; q < output.size(); ++q)
1004 for (
unsigned int i = 0; i < spacedim; ++i)
1005 for (
unsigned int j = 0; j < spacedim; ++j)
1006 for (
unsigned int k = 0; k < spacedim; ++k)
1008 output[q][i][j][k] =
1019 "update_covariant_transformation"));
1021 for (
unsigned int q = 0; q < output.size(); ++q)
1022 for (
unsigned int i = 0; i < spacedim; ++i)
1023 for (
unsigned int j = 0; j < spacedim; ++j)
1024 for (
unsigned int k = 0; k < spacedim; ++k)
1026 output[q][i][j][k] =
1038 "update_covariant_transformation"));
1041 "update_contravariant_transformation"));
1044 "update_volume_elements"));
1046 for (
unsigned int q = 0; q < output.size(); ++q)
1047 for (
unsigned int i = 0; i < spacedim; ++i)
1048 for (
unsigned int j = 0; j < spacedim; ++j)
1049 for (
unsigned int k = 0; k < spacedim; ++k)
1051 output[q][i][j][k] =
1066 template <
int dim,
int spacedim>
1077 length[0] = cell->vertex(1)(0) - start(0);
1080 length[0] = cell->vertex(1)(0) - start(0);
1081 length[1] = cell->vertex(2)(1) - start(1);
1084 length[0] = cell->vertex(1)(0) - start(0);
1085 length[1] = cell->vertex(2)(1) - start(1);
1086 length[2] = cell->vertex(4)(2) - start(2);
1093 for (
unsigned int d = 0; d < dim; ++d)
1094 p_real(d) += length[d] * p(d);
1101 template <
int dim,
int spacedim>
1107 if (dim != spacedim)
1116 real(0) /= cell->vertex(1)(0) - start(0);
1119 real(0) /= cell->vertex(1)(0) - start(0);
1120 real(1) /= cell->vertex(2)(1) - start(1);
1123 real(0) /= cell->vertex(1)(0) - start(0);
1124 real(1) /= cell->vertex(2)(1) - start(1);
1125 real(2) /= cell->vertex(4)(2) - start(2);
1134 template <
int dim,
int spacedim>
1135 std::unique_ptr<Mapping<dim, spacedim>>
1138 return std_cxx14::make_unique<MappingCartesian<dim, spacedim>>(*this);
1144 #include "mapping_cartesian.inst" 1147 DEAL_II_NAMESPACE_CLOSE
Transformed quadrature weights.
virtual Point< dim > transform_real_to_unit_cell(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const Point< spacedim > &p) const override
#define AssertDimension(dim1, dim2)
Contravariant transformation.
Tensor< 1, dim > cell_extents
Outer normal vector, not normalized.
std::vector< Point< dim > > quadrature_points
Determinant of the Jacobian.
Transformed quadrature points.
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_subface_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const override
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static DataSetDescriptor cell()
static Quadrature< dim > project_to_all_subfaces(const SubQuadrature &quadrature)
virtual std::size_t memory_consumption() const override
virtual std::unique_ptr< Mapping< dim, spacedim > > clone() const override
#define Assert(cond, exc)
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingType type, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const override
double weight(const unsigned int i) const
Abstract base class for mapping classes.
Gradient of volume element.
static const unsigned int invalid_face_number
static Quadrature< dim > project_to_all_faces(const SubQuadrature &quadrature)
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_face_data(const UpdateFlags flags, const Quadrature< dim-1 > &quadrature) const override
InternalData(const Quadrature< dim > &quadrature)
virtual std::unique_ptr< typename Mapping< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags, const Quadrature< dim > &quadrature) const override
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) 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
static::ExceptionBase & ExcNotImplemented()
void compute_fill(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const CellSimilarity::Similarity cell_similarity, const InternalData &data, std::vector< Point< dim >> &quadrature_points, std::vector< Tensor< 1, dim >> &normal_vectors) const
Covariant transformation.
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static::ExceptionBase & ExcInternalError()
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)