From 173d0648d9c842c3f8fcb2c3bf014ae97709575b Mon Sep 17 00:00:00 2001
From: Timo Koch <timo.koch@iws.uni-stuttgart.de>
Date: Tue, 24 Aug 2021 19:21:24 +0200
Subject: [PATCH] [geometry] Add simple bounding sphere algorithm
---
dumux/geometry/boundingsphere.hh | 97 ++++++++++++++++++++++++++++++++
1 file changed, 97 insertions(+)
create mode 100644 dumux/geometry/boundingsphere.hh
diff --git a/dumux/geometry/boundingsphere.hh b/dumux/geometry/boundingsphere.hh
new file mode 100644
index 0000000000..993f02c87b
--- /dev/null
+++ b/dumux/geometry/boundingsphere.hh
@@ -0,0 +1,97 @@
+/*****************************************************************************
+ * See the file COPYING for full copying permissions. *
+ * *
+ * This program is free software: you can redistribute it and/or modify *
+ * it under the terms of the GNU General Public License as published by *
+ * the Free Software Foundation, either version 3 of the License, or *
+ * (at your option) any later version. *
+ * *
+ * This program is distributed in the hope that it will be useful, *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
+ * GNU General Public License for more details. *
+ * *
+ * You should have received a copy of the GNU General Public License *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>. *
+ *****************************************************************************/
+/*!
+ * \file
+ * \ingroup Geometry
+ * \brief A function to compute bounding spheres of points clouds or convex polytopes
+ */
+#ifndef DUMUX_GEOMETRY_BOUNDINGSPHERE_HH
+#define DUMUX_GEOMETRY_BOUNDINGSPHERE_HH
+
+#include <algorithm>
+#include <vector>
+
+#include <dumux/geometry/sphere.hh>
+
+namespace Dumux {
+
+/*!
+ * \ingroup Geometry
+ * \brief Computes a bounding sphere of a convex polytope geometry (Dune::Geometry interface)
+ * \note The bounding sphere is not necessarily the minimum enclosing sphere
+ * but computation is fast (AABB method)
+ */
+template<class ConvexGeometry>
+static inline Sphere<typename ConvexGeometry::ctype, ConvexGeometry::coorddimension>
+boundingSphere(const ConvexGeometry& geometry)
+{
+ constexpr int dimWorld = ConvexGeometry::coorddimension;
+ assert(geometry.corners() >= 1);
+
+ auto corner = geometry.corner(0);
+ auto xMin = corner, xMax = corner;
+ using std::max; using std::min;
+
+ // Compute the min and max over the remaining vertices
+ for (std::size_t i = 1; i < geometry.corners(); ++i)
+ {
+ corner = geometry.corner(i);
+ for (std::size_t dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
+ {
+ xMin[dimIdx] = min(xMin[dimIdx], corner[dimIdx]);
+ xMax[dimIdx] = max(xMax[dimIdx], corner[dimIdx]);
+ }
+ }
+
+ auto center = 0.5*xMax + xMin;
+ auto radius = (center-xMax).two_norm();
+ return { std::move(center), std::move(radius) };
+}
+
+namespace Detail {
+template<class Scalar, int dim>
+struct PointsToGeometryWrapper
+{
+ const std::vector<Dune::FieldVector<Scalar, dim>>& points_;
+ PointsToGeometryWrapper(const std::vector<Dune::FieldVector<Scalar, dim>>& points)
+ : points_(points) {}
+
+ static constexpr int coorddimension = dim;
+ using ctype = Scalar;
+
+ auto corners() const { return points_.size(); }
+ const auto& corner(std::size_t i) const { return points_[i]; }
+};
+} // end namespace Detail
+
+/*!
+ * \ingroup Geometry
+ * \brief Computes a bounding sphere of points
+ * \note The bounding sphere is not necessarily the minimum enclosing sphere
+ * but computation is fast (AABB method)
+ */
+template<class Scalar, int dim>
+static inline Sphere<Scalar, dim>
+boundingSphere(const std::vector<Dune::FieldVector<Scalar, dim>>& points)
+{
+ Detail::PointsToGeometryWrapper<Scalar, dim> geometry(points);
+ return boundingSphere(geometry);
+}
+
+} // end namespace Dumux
+
+#endif
--
GitLab