16 #ifndef dealii_base_bounding_box_h 17 #define dealii_base_bounding_box_h 120 template <
int spacedim,
typename Number =
double>
145 template <
class Container>
152 get_boundary_points();
158 get_boundary_points()
const;
211 extend(
const Number &amount);
229 side_length(
const unsigned int direction)
const;
241 upper_bound(
const unsigned int direction)
const;
247 bounds(
const unsigned int direction)
const;
254 vertex(
const unsigned int index)
const;
261 child(
const unsigned int index)
const;
271 cross_section(
const unsigned int direction)
const;
300 template <
class Archive>
302 serialize(Archive &ar,
const unsigned int version);
313 template <
typename Number>
330 template <
class Container>
340 template <
int dim,
typename Number =
double>
382 const int coordinate_in_dim)
386 return (locked_coordinate + coordinate_in_dim + 1) % (dim + 1);
396 template <
int spacedim,
typename Number>
402 for (
unsigned int i = 0; i < spacedim; ++i)
403 Assert(boundary_points.first[i] <= boundary_points.second[i],
404 ExcMessage(
"Bounding Box can't be created: the points' " 405 "order should be bottom left, top right!"));
407 this->boundary_points = boundary_points;
412 template <
int spacedim,
typename Number>
413 template <
class Container>
418 if (points.size() > 0)
420 auto &
min = boundary_points.first;
421 auto &
max = boundary_points.second;
422 std::fill(
min.begin_raw(),
424 std::numeric_limits<Number>::infinity());
425 std::fill(
max.begin_raw(),
427 -std::numeric_limits<Number>::infinity());
430 for (
unsigned int d = 0;
d < spacedim; ++
d)
440 template <
int spacedim,
typename Number>
444 return this->boundary_points;
449 template <
int spacedim,
typename Number>
453 return this->boundary_points;
458 template <
int spacedim,
typename Number>
468 template <
int spacedim,
typename Number>
478 template <
int spacedim,
typename Number>
482 for (
unsigned int d = 0;
d < spacedim; ++
d)
484 boundary_points.first[
d] -= amount;
485 boundary_points.second[
d] += amount;
486 Assert(boundary_points.first[
d] <= boundary_points.second[
d],
487 ExcMessage(
"Bounding Box can't be shrunk this much: the points' " 488 "order should remain bottom left, top right."));
493 template <
int spacedim,
typename Number>
494 template <
class Archive>
504 template <
typename Number>
512 template <
typename Number>
521 template <
typename Number>
522 template <
class Container>
Iterator lower_bound(Iterator first, Iterator last, const T &val)
bool operator!=(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
void extend(const Number &amount)
#define AssertIndexRange(index, range)
#define AssertThrow(cond, exc)
int coordinate_to_one_dim_higher(const int locked_coordinate, const int coordinate_in_dim)
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
static ::ExceptionBase & ExcMessage(std::string arg1)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define Assert(cond, exc)
#define DEAL_II_NAMESPACE_CLOSE
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > boundary_points
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)
BoundingBox< dim, Number > create_unit_bounding_box()
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > & get_boundary_points()
bool operator!=(const BoundingBox< spacedim, Number > &box) const
#define DEAL_II_NAMESPACE_OPEN
T min(const T &t, const MPI_Comm &mpi_communicator)
bool operator==(const BoundingBox< spacedim, Number > &box) const
void serialize(Archive &ar, const unsigned int version)
T max(const T &t, const MPI_Comm &mpi_communicator)