diff --git a/Bembel/ClusterTree b/Bembel/ClusterTree index 975201178..c4727d1ea 100644 --- a/Bembel/ClusterTree +++ b/Bembel/ClusterTree @@ -1,6 +1,6 @@ // This file is part of Bembel, the higher order C++ boundary element library. // -// Copyright (C) 2022 see +// Copyright (C) 2024 see // // 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" diff --git a/Bembel/src/ClusterTree/ElementTree.hpp b/Bembel/src/ClusterTree/ElementTree.hpp index 79628747c..6c7bb1b79 100644 --- a/Bembel/src/ClusterTree/ElementTree.hpp +++ b/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. * diff --git a/Bembel/src/util/GeometryHelper.hpp b/Bembel/src/util/GeometryHelper.hpp new file mode 100644 index 000000000..92b71a18e --- /dev/null +++ b/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 +// +// 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 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_ diff --git a/README.md b/README.md index 318eabfbb..9a85b6f71 100644 --- a/README.md +++ b/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 diff --git a/doc/CodingGuidelines.dox b/doc/CodingGuidelines.dox index c340ae93f..94211c5ee 100644 --- a/doc/CodingGuidelines.dox +++ b/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