Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
bounding_box.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_base_bounding_box_h
17 #define dealii_base_bounding_box_h
18 
19 
20 #include <deal.II/base/config.h>
21 
23 #include <deal.II/base/point.h>
24 
25 #include <limits>
26 
28 
33 enum class NeighborType
34 {
38  not_neighbors = 0,
39 
46  simple_neighbors = 1,
47 
55 
69 };
70 
135 template <int spacedim, typename Number = double>
137 {
138 public:
143  BoundingBox() = default;
144 
149 
155 
160 
167  &boundary_points);
168 
176  template <class Container>
177  BoundingBox(const Container &points);
178 
182  std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
184 
188  const std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
190 
194  bool
196 
200  bool
202 
207  bool
209  const BoundingBox<spacedim, Number> &other_bbox,
210  const double tolerance = std::numeric_limits<Number>::epsilon()) const;
211 
217  const BoundingBox<spacedim, Number> &other_bbox,
218  const double tolerance = std::numeric_limits<Number>::epsilon()) const;
219 
225  void
226  merge_with(const BoundingBox<spacedim, Number> &other_bbox);
227 
234  bool
235  point_inside(
236  const Point<spacedim, Number> &p,
237  const double tolerance = std::numeric_limits<Number>::epsilon()) const;
238 
249  void
250  extend(const Number amount);
251 
257  create_extended(const Number amount) const;
258 
273  create_extended_relative(const Number relative_amount) const;
274 
278  double
279  volume() const;
280 
285  center() const;
286 
290  Number
291  side_length(const unsigned int direction) const;
292 
296  Number
297  lower_bound(const unsigned int direction) const;
298 
302  Number
303  upper_bound(const unsigned int direction) const;
304 
309  bounds(const unsigned int direction) const;
310 
316  vertex(const unsigned int index) const;
317 
323  child(const unsigned int index) const;
324 
332  BoundingBox<spacedim - 1, Number>
333  cross_section(const unsigned int direction) const;
334 
345 
363  Number
365  const unsigned int direction) const;
366 
372  Number
374 
380  template <class Archive>
381  void
382  serialize(Archive &ar, const unsigned int version);
383 
384 private:
385  std::pair<Point<spacedim, Number>, Point<spacedim, Number>> boundary_points;
386 };
387 
393 template <typename Number>
394 class BoundingBox<0, Number>
395 {
396 public:
401 
406 
410  template <class Container>
411  BoundingBox(const Container &);
412 };
413 
414 
420 template <int dim, typename Number = double>
423 
424 
425 namespace internal
426 {
459  template <int dim>
460  inline int
461  coordinate_to_one_dim_higher(const int locked_coordinate,
462  const int coordinate_in_dim)
463  {
464  AssertIndexRange(locked_coordinate, dim + 1);
465  AssertIndexRange(coordinate_in_dim, dim);
466  return (locked_coordinate + coordinate_in_dim + 1) % (dim + 1);
467  }
468 
469 } // namespace internal
470 
471 /*------------------------ Inline functions: BoundingBox --------------------*/
472 
473 #ifndef DOXYGEN
474 
475 
476 template <int spacedim, typename Number>
478  const Point<spacedim, Number> &p)
479  : BoundingBox({p, p})
480 {}
481 
482 
483 
484 template <int spacedim, typename Number>
487  &boundary_points)
488 {
489  // We check the Bounding Box is not degenerate
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!"));
494 
495  this->boundary_points = boundary_points;
496 }
497 
498 
499 
500 template <int spacedim, typename Number>
501 template <class Container>
502 inline BoundingBox<spacedim, Number>::BoundingBox(const Container &points)
503 {
504  // Use the default constructor in case points is empty instead of setting
505  // things to +oo and -oo
506  if (points.size() > 0)
507  {
508  auto &min = boundary_points.first;
509  auto &max = boundary_points.second;
510  for (unsigned int d = 0; d < spacedim; ++d)
511  {
512  min[d] = std::numeric_limits<Number>::infinity();
513  max[d] = -std::numeric_limits<Number>::infinity();
514  }
515 
516  for (const Point<spacedim, Number> &point : points)
517  for (unsigned int d = 0; d < spacedim; ++d)
518  {
519  min[d] = std::min(min[d], point[d]);
520  max[d] = std::max(max[d], point[d]);
521  }
522  }
523 }
524 
525 
526 
527 template <int spacedim, typename Number>
528 inline std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
530 {
531  return this->boundary_points;
532 }
533 
534 
535 
536 template <int spacedim, typename Number>
537 inline const std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
539 {
540  return this->boundary_points;
541 }
542 
543 
544 
545 template <int spacedim, typename Number>
546 inline bool
548  const BoundingBox<spacedim, Number> &box) const
549 {
550  return boundary_points == box.boundary_points;
551 }
552 
553 
554 
555 template <int spacedim, typename Number>
556 inline bool
558  const BoundingBox<spacedim, Number> &box) const
559 {
560  return boundary_points != box.boundary_points;
561 }
562 
563 
564 
565 template <int spacedim, typename Number>
566 inline void
567 BoundingBox<spacedim, Number>::extend(const Number amount)
568 {
569  for (unsigned int d = 0; d < spacedim; ++d)
570  {
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."));
576  }
577 }
578 
579 
580 
581 template <int spacedim, typename Number>
583 BoundingBox<spacedim, Number>::create_extended(const Number amount) const
584 {
585  // create and modify copy
586  auto bb = *this;
587  bb.extend(amount);
588 
589  return bb;
590 }
591 
592 
593 
594 template <int spacedim, typename Number>
597  const Number relative_amount) const
598 {
599  // create and modify copy
600  auto bb = *this;
601 
602  for (unsigned int d = 0; d < spacedim; ++d)
603  {
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."));
609  }
610 
611  return bb;
612 }
613 
614 
615 
616 template <int spacedim, typename Number>
617 template <class Archive>
618 void
620  const unsigned int /*version*/)
621 {
622  ar &boundary_points;
623 }
624 
625 
626 
627 template <typename Number>
629 {
630  AssertThrow(false, ExcImpossibleInDim(0));
631 }
632 
633 
634 
635 template <typename Number>
637  const std::pair<Point<0, Number>, Point<0, Number>> &)
638 {
639  AssertThrow(false, ExcImpossibleInDim(0));
640 }
641 
642 
643 
644 template <typename Number>
645 template <class Container>
646 inline BoundingBox<0, Number>::BoundingBox(const Container &)
647 {
648  AssertThrow(false, ExcImpossibleInDim(0));
649 }
650 
651 
652 
653 #endif // DOXYGEN
655 
656 #endif
NeighborType
Definition: bounding_box.h:34
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)
Definition: bounding_box.h:461
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > boundary_points
Definition: bounding_box.h:385
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
Definition: bounding_box.cc:62
BoundingBox< dim, Number > create_unit_bounding_box()
Point< spacedim, Number > center() const
BoundingBox()=default
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)
Definition: bounding_box.cc:46
bool point_inside(const Point< spacedim, Number > &p, const double tolerance=std::numeric_limits< Number >::epsilon()) const
Definition: bounding_box.cc:26
double volume() 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
Definition: bounding_box.cc:80
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
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1732
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:192
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)