preCICE v3.1.2
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
BoundingBox.cpp
Go to the documentation of this file.
1#include "BoundingBox.hpp"
2#include <algorithm>
3#include <limits>
4#include <ostream>
5#include <string>
6#include <utility>
8#include "logging/Logger.hpp"
9#include "mesh/Vertex.hpp"
10#include "utils/assertion.hpp"
11
12namespace precice::mesh {
13
14logging::Logger BoundingBox::_log{"mesh::BoundingBox"};
15
16BoundingBox::BoundingBox(Eigen::VectorXd boundMin, Eigen::VectorXd boundMax)
17{
18 PRECICE_ASSERT((boundMin.rows() == 2 && boundMax.rows() == 2) || (boundMin.rows() == 3 && boundMax.rows() == 3),
19 "Dimension of min {} and max {} vertices should be the same and both 2 or 3.", boundMin.rows(), boundMax.rows());
20 PRECICE_ASSERT((boundMin - boundMax).maxCoeff() <= 0, "Each component of min vertex {} must be <= max vertex {} in the same axis direction.", boundMin, boundMax);
21
22 _boundMin = std::move(boundMin);
23 _boundMax = std::move(boundMax);
24 _dimensions = _boundMin.rows();
25}
26
28{
29 PRECICE_ASSERT((int) bounds.size() == 4 || (int) bounds.size() == 6, "Dimension of a bounding box can only be 2 or 3. The dimension is {}.", bounds.size() / 2);
30
31 _dimensions = bounds.size() / 2;
32
33 // Eigen::Map maps an existing array to a Eigen's matrix or vector
34 // Here, its 1st input is the existing array's pointer (bounds.data()) and 2nd input is output VectorXd's size
35 // "bounds.data() + 1" means that output elements starts from the 2nd element of "bounds"
36 _boundMin = Eigen::Map<Eigen::VectorXd, 0, Eigen::InnerStride<2>>(bounds.data(), _dimensions);
37 _boundMax = Eigen::Map<Eigen::VectorXd, 0, Eigen::InnerStride<2>>(bounds.data() + 1, _dimensions);
38}
39
41 : _dimensions(dimension)
42{
43 PRECICE_ASSERT(dimension == 2 || dimension == 3, "Dimension of a bounding box can only be 2 or 3 (in this case, it is {}).", dimension);
44
45 _boundMin = Eigen::VectorXd(dimension);
46 _boundMax = Eigen::VectorXd(dimension);
47
48 // Define 'illegal' BoundingBox: _boundMin > _boundMax
51}
52
53bool BoundingBox::operator==(const BoundingBox &otherBB) const
54{
55 PRECICE_ASSERT(_dimensions == otherBB._dimensions, "Bounding boxes with different dimensions ({} and {}) cannot be compared.", _dimensions, otherBB._dimensions);
56
57 return _boundMin.isApprox(otherBB._boundMin) && _boundMax.isApprox(otherBB._boundMax);
58}
59
61{
62 return isDefault() || (_boundMin - _boundMax).isZero();
63}
64
66{
67 return _boundMin.isApproxToConstant(std::numeric_limits<double>::max()) && _boundMax.isApproxToConstant(std::numeric_limits<double>::lowest());
68}
69
70bool BoundingBox::contains(const mesh::Vertex &vertex) const
71{
72 PRECICE_ASSERT(_dimensions == vertex.getDimensions(), "Vertex with different dimensions than this bounding box cannot be checked.");
73
74 const auto &coords = vertex.rawCoords();
75 for (int d = 0; d < _dimensions; d++) {
76 if (coords[d] < _boundMin[d] || coords[d] > _boundMax[d]) {
77 return false;
78 }
79 }
80 return true;
81}
82
83Eigen::VectorXd BoundingBox::center() const
84{
85 PRECICE_ASSERT(!isDefault(), "Data of the bounding box is at default state.");
86
87 return (_boundMax + _boundMin) * 0.5;
88}
89
90Eigen::VectorXd BoundingBox::minCorner() const
91{
92 return _boundMin;
93}
94
95Eigen::VectorXd BoundingBox::maxCorner() const
96{
97 return _boundMax;
98}
99
100double BoundingBox::getEdgeLength(int direction) const
101{
102 PRECICE_ASSERT(direction < _dimensions);
103 PRECICE_ASSERT(_boundMax.size() > direction && _boundMin.size() > direction);
104 return std::abs(_boundMax(direction) - _boundMin(direction));
105}
106
108{
109 return (_boundMax - _boundMin).maxCoeff();
110}
111
113{
114 PRECICE_ASSERT(!isDefault(), "Data of the bounding box is at default state.");
115
116 double meshArea = 1.0;
117 for (int d = 0; d < _dimensions; d++)
118 if (not deadAxis[d])
119 meshArea *= _boundMax[d] - _boundMin[d];
120 return meshArea;
121}
122
124{
125 return _dimensions;
126}
127
129{
131 if (_dimensions == 3) {
132 bounds.insert(bounds.end(), {_boundMin[2], _boundMax[2]});
133 }
134 return bounds;
135}
136
138{
139 PRECICE_ASSERT(_dimensions == otherBB.getDimension(), "Other BoundingBox with different dimensions than this bounding box cannot be used to expand it.");
140
141 _boundMin = _boundMin.cwiseMin(otherBB._boundMin);
142 _boundMax = _boundMax.cwiseMax(otherBB._boundMax);
143}
144
145void BoundingBox::expandBy(const Vertex &vertices)
146{
147 PRECICE_ASSERT(_dimensions == vertices.getDimensions(), "Vertex with different dimensions than this bounding box cannot be used to expand it.");
148
149 const auto coords = vertices.getCoords();
150 _boundMin = _boundMin.cwiseMin(coords);
151 _boundMax = _boundMax.cwiseMax(coords);
152}
153
154void BoundingBox::expandBy(double value)
155{
156 if (!isDefault()) {
157 _boundMin.array() -= value;
158 _boundMax.array() += value;
159 }
160}
161
162void BoundingBox::scaleBy(double safetyFactor)
163{
164 if (!isDefault()) {
165 double maxSideLength = longestEdgeLength();
166 _boundMax.array() += safetyFactor * maxSideLength;
167 _boundMin.array() -= safetyFactor * maxSideLength;
168 PRECICE_DEBUG("Merged BoundingBox {}", *this);
169 }
170}
171
172bool BoundingBox::overlapping(const BoundingBox &otherBB) const
173{
174 for (int d = 0; d < _dimensions; d++) {
175 if ((_boundMin[d] < otherBB._boundMin[d] && _boundMax[d] < otherBB._boundMin[d]) ||
176 (otherBB._boundMin[d] < _boundMin[d] && otherBB._boundMax[d] < _boundMin[d])) {
177 return false;
178 }
179 }
180 return true;
181}
182
184{
185 out << "( ";
186 for (int d = 0; d < _dimensions; ++d) {
187 out << "[" << _boundMin[d] << " " << _boundMax[d] << "], ";
188 }
189 out << ")";
190}
191
193{
194 bb.print(os);
195 return os;
196}
197
198} // namespace precice::mesh
std::ostream & out
#define PRECICE_DEBUG(...)
Definition LogMacros.hpp:64
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
An axis-aligned bounding box around a (partition of a) mesh.
bool contains(const Vertex &vertex) const
Checks if vertex in contained in _bb.
Eigen::VectorXd maxCorner() const
the max corner of the bounding box
void print(std::ostream &out) const
Print bounds of bounding box, output operator overload.
Eigen::VectorXd minCorner() const
the min corner of the bounding box
bool empty() const
Check if every dimension's length is equal to zero.
double longestEdgeLength() const
returns the maximum length of the bounding box in any dimension
static logging::Logger _log
double getArea(std::vector< bool > deadAxis)
Calculate the area of bounding box.
void expandBy(const BoundingBox &otherBB)
Expand bounding box using another bounding box.
const std::vector< double > dataVector() const
Return data as std::vector.
bool overlapping(const BoundingBox &otherBB) const
Checks whether two bounding boxes are overlapping.
bool operator==(const BoundingBox &otherBB) const
Comparison Operator.
int getDimension() const
Getter dimension of the bounding box.
int _dimensions
Number of dimensions (2 or 3)
void scaleBy(double safetyFactor)
Increase the size of bounding box by safety margin.
Eigen::VectorXd center() const
Returns the Center Of Gravity of the mesh.
double getEdgeLength(int axis) const
returns the edge length of a specific axis
Vertex of a mesh.
Definition Vertex.hpp:16
const RawCoords & rawCoords() const
Direct access to the coordinates.
Definition Vertex.hpp:123
int getDimensions() const
Returns spatial dimensionality of vertex.
Definition Vertex.cpp:7
Eigen::VectorXd getCoords() const
Returns the coordinates of the vertex.
Definition Vertex.hpp:116
T data(T... args)
T insert(T... args)
provides Mesh, Data and primitives.
std::ostream & operator<<(std::ostream &os, const BoundingBox &bb)
T size(T... args)