Reference documentation for deal.II version GIT 9c182271f7 2023-03-28 14:30: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 - 2021 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 
262  double
263  volume() const;
264 
269  center() const;
270 
274  Number
275  side_length(const unsigned int direction) const;
276 
280  Number
281  lower_bound(const unsigned int direction) const;
282 
286  Number
287  upper_bound(const unsigned int direction) const;
288 
293  bounds(const unsigned int direction) const;
294 
300  vertex(const unsigned int index) const;
301 
307  child(const unsigned int index) const;
308 
316  BoundingBox<spacedim - 1, Number>
317  cross_section(const unsigned int direction) const;
318 
329 
340 
346  template <class Archive>
347  void
348  serialize(Archive &ar, const unsigned int version);
349 
350 private:
351  std::pair<Point<spacedim, Number>, Point<spacedim, Number>> boundary_points;
352 };
353 
359 template <typename Number>
360 class BoundingBox<0, Number>
361 {
362 public:
367 
372 
376  template <class Container>
377  BoundingBox(const Container &);
378 };
379 
380 
386 template <int dim, typename Number = double>
389 
390 
391 namespace internal
392 {
425  template <int dim>
426  inline int
427  coordinate_to_one_dim_higher(const int locked_coordinate,
428  const int coordinate_in_dim)
429  {
430  AssertIndexRange(locked_coordinate, dim + 1);
431  AssertIndexRange(coordinate_in_dim, dim);
432  return (locked_coordinate + coordinate_in_dim + 1) % (dim + 1);
433  }
434 
435 } // namespace internal
436 
437 /*------------------------ Inline functions: BoundingBox --------------------*/
438 
439 #ifndef DOXYGEN
440 
441 
442 template <int spacedim, typename Number>
444  const Point<spacedim, Number> &p)
445  : BoundingBox({p, p})
446 {}
447 
448 
449 
450 template <int spacedim, typename Number>
453  &boundary_points)
454 {
455  // We check the Bounding Box is not degenerate
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!"));
460 
461  this->boundary_points = boundary_points;
462 }
463 
464 
465 
466 template <int spacedim, typename Number>
467 template <class Container>
468 inline BoundingBox<spacedim, Number>::BoundingBox(const Container &points)
469 {
470  // Use the default constructor in case points is empty instead of setting
471  // things to +oo and -oo
472  if (points.size() > 0)
473  {
474  auto &min = boundary_points.first;
475  auto &max = boundary_points.second;
476  for (unsigned int d = 0; d < spacedim; ++d)
477  {
478  min[d] = std::numeric_limits<Number>::infinity();
479  max[d] = -std::numeric_limits<Number>::infinity();
480  }
481 
482  for (const Point<spacedim, Number> &point : points)
483  for (unsigned int d = 0; d < spacedim; ++d)
484  {
485  min[d] = std::min(min[d], point[d]);
486  max[d] = std::max(max[d], point[d]);
487  }
488  }
489 }
490 
491 
492 
493 template <int spacedim, typename Number>
494 inline std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
496 {
497  return this->boundary_points;
498 }
499 
500 
501 
502 template <int spacedim, typename Number>
503 inline const std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
505 {
506  return this->boundary_points;
507 }
508 
509 
510 
511 template <int spacedim, typename Number>
512 inline bool
514  const BoundingBox<spacedim, Number> &box) const
515 {
516  return boundary_points == box.boundary_points;
517 }
518 
519 
520 
521 template <int spacedim, typename Number>
522 inline bool
524  const BoundingBox<spacedim, Number> &box) const
525 {
526  return boundary_points != box.boundary_points;
527 }
528 
529 
530 
531 template <int spacedim, typename Number>
532 inline void
533 BoundingBox<spacedim, Number>::extend(const Number amount)
534 {
535  for (unsigned int d = 0; d < spacedim; ++d)
536  {
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."));
542  }
543 }
544 
545 
546 
547 template <int spacedim, typename Number>
549 BoundingBox<spacedim, Number>::create_extended(const Number amount) const
550 {
551  // create and modify copy
552  auto bb = *this;
553  bb.extend(amount);
554 
555  return bb;
556 }
557 
558 
559 template <int spacedim, typename Number>
560 template <class Archive>
561 void
563  const unsigned int /*version*/)
564 {
565  ar &boundary_points;
566 }
567 
568 
569 
570 template <typename Number>
572 {
573  AssertThrow(false, ExcImpossibleInDim(0));
574 }
575 
576 
577 
578 template <typename Number>
580  const std::pair<Point<0, Number>, Point<0, Number>> &)
581 {
582  AssertThrow(false, ExcImpossibleInDim(0));
583 }
584 
585 
586 
587 template <typename Number>
588 template <class Container>
589 inline BoundingBox<0, Number>::BoundingBox(const Container &)
590 {
591  AssertThrow(false, ExcImpossibleInDim(0));
592 }
593 
594 
595 
596 #endif // DOXYGEN
598 
599 #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:427
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > boundary_points
Definition: bounding_box.h:351
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:61
BoundingBox< dim, Number > create_unit_bounding_box()
Point< spacedim, Number > center() const
BoundingBox()=default
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)
Definition: bounding_box.cc:45
bool point_inside(const Point< spacedim, Number > &p, const double tolerance=std::numeric_limits< Number >::epsilon()) const
Definition: bounding_box.cc:25
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
NeighborType get_neighbor_type(const BoundingBox< spacedim, Number > &other_bbox, const double tolerance=std::numeric_limits< Number >::epsilon()) const
Definition: bounding_box.cc:79
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:474
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
#define Assert(cond, exc)
Definition: exceptions.h:1586
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1827
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1675
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
Definition: utilities.cc:189
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)