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>
293 bounds(
const unsigned int direction)
const;
300 vertex(
const unsigned int index)
const;
307 child(
const unsigned int index)
const;
346 template <
class Archive>
359 template <
typename Number>
376 template <
class Container>
386 template <
int dim,
typename Number =
double>
428 const int coordinate_in_dim)
432 return (locked_coordinate + coordinate_in_dim + 1) % (dim + 1);
442 template <
int spacedim,
typename Number>
450 template <
int spacedim,
typename Number>
456 for (
unsigned int i = 0; i < spacedim; ++i)
457 Assert(boundary_points.first[i] <= boundary_points.second[i],
458 ExcMessage(
"Bounding Box can't be created: the points' "
459 "order should be bottom left, top right!"));
461 this->boundary_points = boundary_points;
466 template <
int spacedim,
typename Number>
467 template <
class Container>
472 if (points.size() > 0)
474 auto &
min = boundary_points.first;
475 auto &
max = boundary_points.second;
476 for (
unsigned int d = 0;
d < spacedim; ++
d)
478 min[
d] = std::numeric_limits<Number>::infinity();
479 max[
d] = -std::numeric_limits<Number>::infinity();
483 for (
unsigned int d = 0;
d < spacedim; ++
d)
493 template <
int spacedim,
typename Number>
497 return this->boundary_points;
502 template <
int spacedim,
typename Number>
506 return this->boundary_points;
511 template <
int spacedim,
typename Number>
521 template <
int spacedim,
typename Number>
531 template <
int spacedim,
typename Number>
535 for (
unsigned int d = 0;
d < spacedim; ++
d)
537 boundary_points.first[
d] -= amount;
538 boundary_points.second[
d] += amount;
539 Assert(boundary_points.first[
d] <= boundary_points.second[
d],
540 ExcMessage(
"Bounding Box can't be shrunk this much: the points' "
541 "order should remain bottom left, top right."));
547 template <
int spacedim,
typename Number>
559 template <
int spacedim,
typename Number>
560 template <
class Archive>
570 template <
typename Number>
578 template <
typename Number>
587 template <
typename Number>
588 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
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
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)