16 #ifndef dealii_base_bounding_box_h
17 #define dealii_base_bounding_box_h
135 template <
int spacedim,
typename Number =
double>
176 template <
class Container>
309 bounds(
const unsigned int direction)
const;
316 vertex(
const unsigned int index)
const;
323 child(
const unsigned int index)
const;
365 const unsigned int direction)
const;
380 template <
class Archive>
393 template <
typename Number>
410 template <
class Container>
420 template <
int dim,
typename Number =
double>
462 const int coordinate_in_dim)
466 return (locked_coordinate + coordinate_in_dim + 1) % (dim + 1);
476 template <
int spacedim,
typename Number>
484 template <
int spacedim,
typename Number>
490 for (
unsigned int i = 0; i < spacedim; ++i)
491 Assert(boundary_points.first[i] <= boundary_points.second[i],
492 ExcMessage(
"Bounding Box can't be created: the points' "
493 "order should be bottom left, top right!"));
495 this->boundary_points = boundary_points;
500 template <
int spacedim,
typename Number>
501 template <
class Container>
506 if (points.size() > 0)
508 auto &
min = boundary_points.first;
509 auto &
max = boundary_points.second;
510 for (
unsigned int d = 0;
d < spacedim; ++
d)
512 min[
d] = std::numeric_limits<Number>::infinity();
513 max[
d] = -std::numeric_limits<Number>::infinity();
517 for (
unsigned int d = 0;
d < spacedim; ++
d)
527 template <
int spacedim,
typename Number>
531 return this->boundary_points;
536 template <
int spacedim,
typename Number>
540 return this->boundary_points;
545 template <
int spacedim,
typename Number>
555 template <
int spacedim,
typename Number>
565 template <
int spacedim,
typename Number>
569 for (
unsigned int d = 0;
d < spacedim; ++
d)
571 boundary_points.first[
d] -= amount;
572 boundary_points.second[
d] += amount;
573 Assert(boundary_points.first[
d] <= boundary_points.second[
d],
574 ExcMessage(
"Bounding Box can't be shrunk this much: the points' "
575 "order should remain bottom left, top right."));
581 template <
int spacedim,
typename Number>
594 template <
int spacedim,
typename Number>
597 const Number relative_amount)
const
602 for (
unsigned int d = 0;
d < spacedim; ++
d)
604 bb.boundary_points.first[
d] -= relative_amount * side_length(
d);
605 bb.boundary_points.second[
d] += relative_amount * side_length(
d);
606 Assert(bb.boundary_points.first[
d] <= bb.boundary_points.second[
d],
607 ExcMessage(
"Bounding Box can't be shrunk this much: the points' "
608 "order should remain bottom left, top right."));
616 template <
int spacedim,
typename Number>
617 template <
class Archive>
627 template <
typename Number>
635 template <
typename Number>
644 template <
typename Number>
645 template <
class Container>
BoundingBox(const std::pair< Point< 0, Number >, Point< 0, Number >> &)
BoundingBox(const Container &)
int coordinate_to_one_dim_higher(const int locked_coordinate, const int coordinate_in_dim)
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > boundary_points
BoundingBox< 1, Number > bounds(const unsigned int direction) const
const std::pair< Point< spacedim, Number >, Point< spacedim, Number > > & get_boundary_points() const
bool has_overlap_with(const BoundingBox< spacedim, Number > &other_bbox, const double tolerance=std::numeric_limits< Number >::epsilon()) const
BoundingBox< dim, Number > create_unit_bounding_box()
Point< spacedim, Number > center() const
Number lower_bound(const unsigned int direction) const
Number signed_distance(const Point< spacedim, Number > &point, const unsigned int direction) const
bool operator==(const BoundingBox< spacedim, Number > &box) const
void serialize(Archive &ar, const unsigned int version)
void merge_with(const BoundingBox< spacedim, Number > &other_bbox)
bool point_inside(const Point< spacedim, Number > &p, const double tolerance=std::numeric_limits< Number >::epsilon()) const
void extend(const Number amount)
BoundingBox(const Container &points)
BoundingBox< spacedim, Number > create_extended(const Number amount) const
Point< spacedim, Number > real_to_unit(const Point< spacedim, Number > &point) const
Number side_length(const unsigned int direction) const
BoundingBox(const std::pair< Point< spacedim, Number >, Point< spacedim, Number >> &boundary_points)
BoundingBox< spacedim, Number > child(const unsigned int index) const
bool operator!=(const BoundingBox< spacedim, Number > &box) const
BoundingBox< spacedim, Number > & operator=(const BoundingBox< spacedim, Number > &t)=default
BoundingBox< spacedim, Number > create_extended_relative(const Number relative_amount) const
NeighborType get_neighbor_type(const BoundingBox< spacedim, Number > &other_bbox, const double tolerance=std::numeric_limits< Number >::epsilon()) const
Point< spacedim, Number > vertex(const unsigned int index) const
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > & get_boundary_points()
BoundingBox(const Point< spacedim, Number > &point)
BoundingBox< spacedim - 1, Number > cross_section(const unsigned int direction) const
BoundingBox(const BoundingBox< spacedim, Number > &box)=default
Number upper_bound(const unsigned int direction) const
Point< spacedim, Number > unit_to_real(const Point< spacedim, Number > &point) const
VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &x, const ::VectorizedArray< Number, width > &y)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)