Extract ComputeEnclosingBall and Add Helper Function Guideline (#99) …

…and Fixes (#98)

* Extract computeEnclosingBall Function to Helper File

* Improve Comment and Variable Names

* Introduce Coding Guideline about Helper Functions

Bembel/ClusterTree

@@ -1,6 +1,6 @@
 // This file is part of Bembel, the higher order C++ boundary element library.
 //
-// Copyright (C) 2022 see <http://www.bembel.eu>
+// Copyright (C) 2024 see <http://www.bembel.eu>
 //
 // It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz,
 // M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt,
@@ -37,6 +37,7 @@
 #include "Geometry"
 
 #include "src/util/Constants.hpp"
+#include "src/util/GeometryHelper.hpp"
 
 #include "src/ClusterTree/ElementTreeNode.hpp"
 #include "src/ClusterTree/ElementTree.hpp"

Bembel/src/ClusterTree/ElementTree.hpp

@@ -234,7 +234,7 @@ class ElementTree {
    *
    * This function stores all indices of an element in the column of the
    * returned matrix.
-   * 
+   *
    * \return A 4xN integer Matrix where N is the number elements.
    */
   Eigen::MatrixXi generateElementList() const {
@@ -372,23 +372,25 @@ class ElementTree {
     // compute point list
     if (!el.sons_.size()) {
       // assign enclosing balls to leafs
-      computeEnclosingBall(&mp1, &r1, P.col(el.vertices_[0]), 0,
-                           P.col(el.vertices_[2]), 0);
-      computeEnclosingBall(&mp2, &r2, P.col(el.vertices_[1]), 0,
-                           P.col(el.vertices_[3]), 0);
-      computeEnclosingBall(&(el.midpoint_), &(el.radius_), mp1, r1, mp2, r2);
+      util::computeEnclosingBall(&mp1, &r1, P.col(el.vertices_[0]), 0,
+                                 P.col(el.vertices_[2]), 0);
+      util::computeEnclosingBall(&mp2, &r2, P.col(el.vertices_[1]), 0,
+                                 P.col(el.vertices_[3]), 0);
+      util::computeEnclosingBall(&(el.midpoint_), &(el.radius_), mp1, r1, mp2,
+                                 r2);
     } else {
       // handle the four(!!!) children
       for (auto i = 0; i < 4; ++i)
         computeElementEnclosings_recursion(el.sons_[i], P);
       // assign enclosing balls to fathers bottom up
-      computeEnclosingBall(&mp1, &r1, el.sons_[0].midpoint_,
-                           el.sons_[0].radius_, el.sons_[2].midpoint_,
-                           el.sons_[2].radius_);
-      computeEnclosingBall(&mp2, &r2, el.sons_[1].midpoint_,
-                           el.sons_[1].radius_, el.sons_[3].midpoint_,
-                           el.sons_[3].radius_);
-      computeEnclosingBall(&(el.midpoint_), &(el.radius_), mp1, r1, mp2, r2);
+      util::computeEnclosingBall(&mp1, &r1, el.sons_[0].midpoint_,
+                                 el.sons_[0].radius_, el.sons_[2].midpoint_,
+                                 el.sons_[2].radius_);
+      util::computeEnclosingBall(&mp2, &r2, el.sons_[1].midpoint_,
+                                 el.sons_[1].radius_, el.sons_[3].midpoint_,
+                                 el.sons_[3].radius_);
+      util::computeEnclosingBall(&(el.midpoint_), &(el.radius_), mp1, r1, mp2,
+                                 r2);
     }
     return;
   }
@@ -498,35 +500,6 @@ class ElementTree {
     }
     return retval;
   }
-  //////////////////////////////////////////////////////////////////////////////
-  /// static members
-  //////////////////////////////////////////////////////////////////////////////
-  /**
-   *  \brief computes a ball enclosing the union of \f$B_r1(mp1)\f$ and \f$B_r2(mp2)\f$,
-   * i.e \f$B(mp,r)\supset B_r1(mp1) \cup B_r2(mp2)\f$.
-   */
-  static void computeEnclosingBall(Eigen::Vector3d *mp, double *r,
-                                   const Eigen::Vector3d &mp1, double r1,
-                                   const Eigen::Vector3d &mp2, double r2) {
-    // compute distance vector of the two spheres
-    auto z = mp1 - mp2;
-    auto norm = (mp1 - mp2).norm();
-    // B(d2,r2) subset B(d1,r1)
-    if (norm + r2 <= r1) {
-      *mp = mp1;
-      *r = r1;
-      // B(d1,r1) subset B(d2,r2)
-    } else if (norm + r1 <= r2) {
-      *mp = mp2;
-      *r = r2;
-      // the union is not a ball
-    } else {
-      *mp = 0.5 * (mp1 + mp2 + (r1 - r2) / norm * z);
-      *r = 0.5 * (r1 + r2 + norm);
-      *r = 0.5 * (r1 + r2 + norm);
-    }
-    return;
-  }
   /**
    * \brief Resolves neighborhood relations of the patches.
    *

Bembel/src/util/GeometryHelper.hpp

@@ -0,0 +1,55 @@
+// This file is part of Bembel, the higher order C++ boundary element library.
+//
+// Copyright (C) 2024 see <http://www.bembel.eu>
+//
+// It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz,
+// M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt,
+// Universitaet Basel, and Universita della Svizzera italiana, Lugano. This
+// source code is subject to the GNU General Public License version 3 and
+// provided WITHOUT ANY WARRANTY, see <http://www.bembel.eu> for further
+// information.
+
+#ifndef BEMBEL_SRC_UTIL_GEOMETRYHELPER_HPP_
+#define BEMBEL_SRC_UTIL_GEOMETRYHELPER_HPP_
+
+namespace Bembel {
+namespace util {
+
+/**
+ *  \brief computes a ball enclosing of two given enclosing balls.
+ *
+ * The computed enclosing ball contains for sure the two given balls with
+ * mit_point1 and radius1 and 2, respectively.
+ *
+ * \param *midpoint: Pointer to Eigen::Vector3d which is the new mid point of
+ * the enclosing ball
+ * \param *radius: Pointer where the radius of new enclosing ball is saved.
+ * \param midpoint1: Eigen::Vector3d midpoint of the first ball
+ * \param radius1: double radius of the first ball
+ * \param midpoint2: Eigen::Vector3d midpoint of the second ball
+ * \param radius2: double radius of the second ball
+ */
+void computeEnclosingBall(Eigen::Vector3d *midpoint, double *radius,
+                          const Eigen::Vector3d &midpoint1, double radius1,
+                          const Eigen::Vector3d &midpoint2, double radius2) {
+  // compute distance vector of the two spheres
+  auto z = midpoint1 - midpoint2;
+  auto norm = (midpoint1 - midpoint2).norm();
+  // B(d2,radius2) subset B(d1,radius1)
+  if (norm + radius2 <= radius1) {
+    *midpoint = midpoint1;
+    *radius = radius1;
+    // B(d1,radius1) subset B(d2,radius2)
+  } else if (norm + radius1 <= radius2) {
+    *midpoint = midpoint2;
+    *radius = radius2;
+    // the union is not a ball
+  } else {
+    *midpoint = 0.5 * (midpoint1 + midpoint2 + (radius1 - radius2) / norm * z);
+    *radius = 0.5 * (radius1 + radius2 + norm);
+  }
+  return;
+}
+}  // namespace util
+}  // namespace Bembel
+#endif  // BEMBEL_SRC_UTIL_GEOMETRYHELPER_HPP_

README.md

@@ -67,7 +67,7 @@ Please send bug reports and feature requests via [issue tracker](https://github.
 
 Any contribution to this project in fixing a bug or implementing a feature is welcome.
 Create a fork of this repository and create a [pull request](https://github.com/temf/bembel/pulls).
-To successfully merge your pull request you should follow our [Coding Guidelines](https://temf.github.io/bembel/Doxy_out/html/codingGuidelines.html)
+To successfully merge your pull request you should follow our [Coding Guidelines](https://temf.github.io/bembel/Doxy_out/html/coding_guidelines.html)
 
 
 ## 6. How to cite <a name="cite"></a>

doc/CodingGuidelines.dox

@@ -45,6 +45,15 @@ Citing Kent Beck: "You need to convey the purpose, type, and lifetime of the var
 
 \ref FAQNaming
 
+---
+\subsection ReuseCode Reuse Code
+
+### Helper Functions:
+
+Put code which might be useful outside of your feature into the `src/util` folder.
+Name the file topic + "Helper".
+Everything needs to be placed in the `util` namespace.
+
 ---
 \section git Git Workflow