Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
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
33enum class NeighborType
34{
38 not_neighbors = 0,
39
47
55
69};
70
135template <int spacedim, typename Number = double>
137{
138public:
143 BoundingBox() = default;
144
149
155
160
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
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
344 real_to_unit(const Point<spacedim, Number> &point) const;
345
355 unit_to_real(const Point<spacedim, Number> &point) const;
363 Number
365 const unsigned int direction) const;
366
372 Number
373 signed_distance(const Point<spacedim, Number> &point) const;
374
380 template <class Archive>
381 void
382 serialize(Archive &ar, const unsigned int version);
383
384private:
385 std::pair<Point<spacedim, Number>, Point<spacedim, Number>> boundary_points;
386};
387
393template <typename Number>
394class BoundingBox<0, Number>
395{
396public:
401
406
410 template <class Container>
411 BoundingBox(const Container &);
412};
413
414
420template <int dim, typename Number = double>
423
424
425namespace 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
476template <int spacedim, typename Number>
479 : BoundingBox({p, p})
480{}
481
482
483
484template <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
500template <int spacedim, typename Number>
501template <class Container>
502inline 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
527template <int spacedim, typename Number>
528inline std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
530{
531 return this->boundary_points;
532}
533
534
535
536template <int spacedim, typename Number>
537inline const std::pair<Point<spacedim, Number>, Point<spacedim, Number>> &
539{
540 return this->boundary_points;
541}
542
543
544
545template <int spacedim, typename Number>
546inline bool
548 const BoundingBox<spacedim, Number> &box) const
549{
550 return boundary_points == box.boundary_points;
551}
552
553
554
555template <int spacedim, typename Number>
556inline bool
558 const BoundingBox<spacedim, Number> &box) const
559{
560 return boundary_points != box.boundary_points;
561}
562
563
564
565template <int spacedim, typename Number>
566inline void
567BoundingBox<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
581template <int spacedim, typename Number>
583BoundingBox<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
594template <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
616template <int spacedim, typename Number>
617template <class Archive>
618void
620 const unsigned int /*version*/)
621{
622 ar &boundary_points;
623}
624
625
626
627template <typename Number>
629{
631}
632
633
634
635template <typename Number>
637 const std::pair<Point<0, Number>, Point<0, Number>> &)
638{
640}
641
642
643
644template <typename Number>
645template <class Container>
646inline BoundingBox<0, Number>::BoundingBox(const Container &)
647{
649}
650
651
652
653#endif // DOXYGEN
655
656#endif
NeighborType
BoundingBox< dim, Number > create_unit_bounding_box()
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
bool has_overlap_with(const BoundingBox< spacedim, Number > &other_bbox, const double tolerance=std::numeric_limits< Number >::epsilon()) const
BoundingBox(const std::pair< Point< spacedim, Number >, Point< spacedim, Number > > &boundary_points)
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)
std::pair< Point< spacedim, Number >, Point< spacedim, Number > > & get_boundary_points()
bool point_inside(const Point< spacedim, Number > &p, const double tolerance=std::numeric_limits< Number >::epsilon()) const
double volume() const
void extend(const Number amount)
BoundingBox< spacedim, Number > create_extended_relative(const Number relative_amount) const
BoundingBox(const Container &points)
Point< spacedim, Number > real_to_unit(const Point< spacedim, Number > &point) const
const std::pair< Point< spacedim, Number >, Point< spacedim, Number > > & get_boundary_points() const
Number side_length(const unsigned int direction) const
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
BoundingBox< spacedim, Number > create_extended(const Number amount) const
Point< spacedim, Number > vertex(const unsigned int index) const
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
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)