16 #ifndef dealii__geometry_info_h 17 #define dealii__geometry_info_h 20 #include <deal.II/base/config.h> 21 #include <deal.II/base/exceptions.h> 22 #include <deal.II/base/point.h> 24 DEAL_II_NAMESPACE_OPEN
90 operator unsigned int ()
const;
162 isotropic_refinement =
static_cast<unsigned char>(-1)
228 isotropic_refinement = cut_x
299 cut_xy = cut_x | cut_y,
304 isotropic_refinement = cut_xy
375 cut_xy = cut_x | cut_y,
383 cut_xz = cut_x | cut_z,
387 cut_yz = cut_y | cut_z,
391 cut_xyz = cut_x | cut_y | cut_z,
396 isotropic_refinement = cut_xyz
446 operator unsigned char ()
const;
481 static std::size_t memory_consumption ();
487 template <
class Archive>
488 void serialize(Archive &ar,
489 const unsigned int version);
496 <<
"The refinement flags given (" << arg1 <<
") contain set bits that do not " 497 <<
"make sense for the space dimension of the object to which they are applied.");
504 unsigned char value :
549 case_isotropic =
static_cast<unsigned char>(-1)
581 case_isotropic = case_none
615 case_isotropic = case_none
654 case_isotropic = case_x
774 case_isotropic = case_xy
810 operator unsigned char ()
const;
815 static std::size_t memory_consumption ();
822 <<
"The subface case given (" << arg1 <<
") does not make sense " 823 <<
"for the space dimension of the object to which they are applied.");
830 unsigned char value :
871 static const unsigned int max_children_per_cell = 1;
876 static const unsigned int faces_per_cell = 0;
885 static const unsigned int max_children_per_face = 0;
897 static const unsigned int vertices_per_cell = 1;
905 static const unsigned int vertices_per_face = 0;
910 static const unsigned int lines_per_face = 0;
915 static const unsigned int quads_per_face = 0;
920 static const unsigned int lines_per_cell = 0;
925 static const unsigned int quads_per_cell = 0;
930 static const unsigned int hexes_per_cell = 0;
1472 static const unsigned int max_children_per_cell = 1 << dim;
1477 static const unsigned int faces_per_cell = 2 * dim;
1486 static const unsigned int max_children_per_face =
GeometryInfo<dim-1>::max_children_per_cell;
1491 static const unsigned int vertices_per_cell = 1 << dim;
1496 static const unsigned int vertices_per_face =
GeometryInfo<dim-1>::vertices_per_cell;
1501 static const unsigned int lines_per_face
1507 static const unsigned int quads_per_face
1519 static const unsigned int lines_per_cell
1530 static const unsigned int quads_per_cell
1537 static const unsigned int hexes_per_cell
1558 static const unsigned int ucd_to_deal[vertices_per_cell];
1573 static const unsigned int dx_to_deal[vertices_per_cell];
1585 static const unsigned int vertex_to_face[vertices_per_cell][dim];
1615 const unsigned int subface_no);
1625 const unsigned int face_no,
1626 const bool face_orientation =
true,
1627 const bool face_flip =
false,
1628 const bool face_rotation =
false);
1637 min_cell_refinement_case_for_face_refinement
1639 const unsigned int face_no,
1640 const bool face_orientation =
true,
1641 const bool face_flip =
false,
1642 const bool face_rotation =
false);
1651 const unsigned int line_no);
1659 min_cell_refinement_case_for_line_refinement(
const unsigned int line_no);
1710 const unsigned int face,
1711 const unsigned int subface,
1712 const bool face_orientation =
true,
1713 const bool face_flip =
false,
1714 const bool face_rotation =
false,
1733 line_to_cell_vertices (
const unsigned int line,
1734 const unsigned int vertex);
1758 face_to_cell_vertices (
const unsigned int face,
1759 const unsigned int vertex,
1760 const bool face_orientation =
true,
1761 const bool face_flip =
false,
1762 const bool face_rotation =
false);
1777 face_to_cell_lines (
const unsigned int face,
1778 const unsigned int line,
1779 const bool face_orientation =
true,
1780 const bool face_flip =
false,
1781 const bool face_rotation =
false);
1794 standard_to_real_face_vertex (
const unsigned int vertex,
1795 const bool face_orientation =
true,
1796 const bool face_flip =
false,
1797 const bool face_rotation =
false);
1810 real_to_standard_face_vertex (
const unsigned int vertex,
1811 const bool face_orientation =
true,
1812 const bool face_flip =
false,
1813 const bool face_rotation =
false);
1826 standard_to_real_face_line (
const unsigned int line,
1827 const bool face_orientation =
true,
1828 const bool face_flip =
false,
1829 const bool face_rotation =
false);
1842 real_to_standard_face_line (
const unsigned int line,
1843 const bool face_orientation =
true,
1844 const bool face_flip =
false,
1845 const bool face_rotation =
false);
1854 unit_cell_vertex (
const unsigned int vertex);
1878 cell_to_child_coordinates (
const Point<dim> &p,
1879 const unsigned int child_index,
1890 child_to_cell_coordinates (
const Point<dim> &p,
1891 const unsigned int child_index,
1942 d_linear_shape_function (
const Point<dim> &xi,
1943 const unsigned int i);
1951 d_linear_shape_function_gradient (
const Point<dim> &xi,
1952 const unsigned int i);
2005 template <
int spacedim>
2007 alternating_form_at_vertices
2008 #ifndef DEAL_II_CONSTEXPR_BUG 2013 Tensor<spacedim-dim,spacedim> *forms);
2025 static const unsigned int unit_normal_direction[faces_per_cell];
2043 static const int unit_normal_orientation[faces_per_cell];
2050 static const unsigned int opposite_face[faces_per_cell];
2058 <<
"The coordinates must satisfy 0 <= x_i <= 1, " 2059 <<
"but here we have x_i=" << arg1);
2066 <<
"RefinementCase<dim> " << arg1 <<
": face " << arg2
2067 <<
" has no subface " << arg3);
2078 #ifndef DEAL_II_MEMBER_ARRAY_SPECIALIZATION_BUG 2112 const unsigned int i);
2117 const unsigned int i);
2122 const unsigned int i);
2141 object (static_cast<Object>(object_dimension))
2146 GeometryPrimitive::operator
unsigned int ()
const 2148 return static_cast<unsigned int>(object);
2158 SubfaceCase<dim>::SubfaceCase (
const typename SubfacePossibilities<dim>::Possibilities subface_possibility)
2160 value (subface_possibility)
2166 SubfaceCase<dim>::operator
unsigned char ()
const 2181 return static_cast<unsigned char>(-1);
2190 const unsigned int dim = 1;
2204 const unsigned int dim = 2;
2220 const unsigned int dim = 3;
2246 value (refinement_case)
2254 ExcInvalidRefinementCase (refinement_case));
2263 value (refinement_case)
2271 ExcInvalidRefinementCase (refinement_case));
2328 template <
class Archive>
2334 unsigned char uchar_value = value;
2336 value = uchar_value;
2350 return Point<1>(
static_cast<double>(vertex));
2363 return Point<2>(vertex%2, vertex/2);
2376 return Point<3>(vertex%2, vertex/2%2, vertex/4);
2400 return (p[0] <= 0.5 ? 0 : 1);
2413 return (p[0] <= 0.5 ?
2414 (p[1] <= 0.5 ? 0 : 2) :
2415 (p[1] <= 0.5 ? 1 : 3));
2429 return (p[0] <= 0.5 ?
2431 (p[2] <= 0.5 ? 0 : 4) :
2432 (p[2] <= 0.5 ? 2 : 6)) :
2434 (p[2] <= 0.5 ? 1 : 5) :
2435 (p[2] <= 0.5 ? 3 : 7)));
2455 const unsigned int child_index,
2474 const unsigned int child_index,
2482 switch (refine_case)
2511 const unsigned int child_index,
2524 switch (refine_case)
2544 if (child_index%2==1)
2546 if (child_index/2==1)
2556 if (child_index/2==1)
2558 if (child_index%2==1)
2564 if (child_index%2==1)
2566 if (child_index/2==1)
2586 const unsigned int ,
2600 const unsigned int child_index,
2619 const unsigned int child_index,
2632 switch (refine_case)
2650 if (child_index%2==1)
2652 if (child_index/2==1)
2662 if (child_index/2==1)
2664 if (child_index%2==1)
2670 if (child_index%2==1)
2672 if (child_index/2==1)
2694 const unsigned int child_index,
2701 switch (refine_case)
2730 const unsigned int ,
2744 return (p[0] >= 0.) && (p[0] <= 1.);
2754 return (p[0] >= 0.) && (p[0] <= 1.) &&
2755 (p[1] >= 0.) && (p[1] <= 1.);
2765 return (p[0] >= 0.) && (p[0] <= 1.) &&
2766 (p[1] >= 0.) && (p[1] <= 1.) &&
2767 (p[2] >= 0.) && (p[2] <= 1.);
2776 return (p[0] >= -eps) && (p[0] <= 1.+eps);
2787 const double l = -eps, u = 1+eps;
2788 return (p[0] >= l) && (p[0] <= u) &&
2789 (p[1] >= l) && (p[1] <= u);
2800 const double l = -eps, u = 1.0+eps;
2801 return (p[0] >= l) && (p[0] <= u) &&
2802 (p[1] >= l) && (p[1] <= u) &&
2803 (p[2] >= l) && (p[2] <= u);
2808 DEAL_II_NAMESPACE_CLOSE
static Point< dim > child_to_cell_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
UpdateFlags operator&(UpdateFlags f1, UpdateFlags f2)
static unsigned int child_cell_from_point(const Point< dim > &p)
RefinementCase operator|(const RefinementCase &r) const
static bool is_inside_unit_cell(const Point< dim > &p)
GeometryPrimitive(const Object object)
void serialize(Archive &ar, const unsigned int version)
static Point< dim > cell_to_child_coordinates(const Point< dim > &p, const unsigned int child_index, const RefinementCase< dim > refine_case=RefinementCase< dim >::isotropic_refinement)
static Point< dim > unit_cell_vertex(const unsigned int vertex)
#define AssertThrow(cond, exc)
static::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static const unsigned int vertices_per_cell
RefinementCase operator~() const
#define DeclException1(Exception1, type1, outsequence)
#define Assert(cond, exc)
static::ExceptionBase & ExcInvalidCoordinate(double arg1)
UpdateFlags operator|(UpdateFlags f1, UpdateFlags f2)
static std::size_t memory_consumption()
RefinementCase operator&(const RefinementCase &r) const
static::ExceptionBase & ExcNotImplemented()
#define DeclException3(Exception3, type1, type2, type3, outsequence)
static Tensor< 1, dim > d_linear_shape_function_gradient(const Point< dim > &xi, const unsigned int i)
static RefinementCase cut_axis(const unsigned int i)
static::ExceptionBase & ExcInternalError()