[BL] Added a class computing subdivisions for fixed number of sub divs.

BaseLib/Subdivision.cpp

+/**
+ * \copyright
+ * Copyright (c) 2012-2017, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "Subdivision.h"
+
+#include <BaseLib/Error.h>
+
+namespace BaseLib
+{
+GradualSubdivision::GradualSubdivision(const double L,
+                                       const double dL0,
+                                       const double max_dL,
+                                       const double multiplier)
+    : _length(L), _dL0(dL0), _max_dL(max_dL), _multiplier(multiplier)
+{
+    // Check if accumulated subdivisions can ever sum up to length.
+    // Cf. geometric series formula.
+    if (multiplier < 1.0 && dL0 / (1.0 - multiplier) < L) {
+        OGS_FATAL(
+            "Using dL0=%g and multiplier=%g the generated subdivisions can not "
+            "sum up to a total length of %g.",
+            dL0,
+            multiplier,
+            L);
+    }
+}
+
+GradualSubdivisionFixedNum::GradualSubdivisionFixedNum(
+    const double L, const std::size_t num_subdivisions, const double multiplier)
+    : _length{L}, _num_subdivisions{num_subdivisions}, _multiplier{multiplier}
+{
+}
+
+std::vector<double> GradualSubdivisionFixedNum::operator()() const
+{
+    std::vector<double> subdivisions;
+    subdivisions.reserve(_num_subdivisions + 1);
+    subdivisions.push_back(0.0);
+    auto const q = _multiplier;
+
+    if (q == 1.0) {
+        double const dx = _length / _num_subdivisions;
+
+        for (std::size_t i = 1; i < _num_subdivisions; ++i) {
+            subdivisions.push_back(dx * i);
+        }
+    } else {
+        // compute initial subdivision size
+        auto const a =
+            _length * (q - 1.0) / (std::pow(q, _num_subdivisions) - 1.0);
+
+        double qi = q;  // q^i
+        for (std::size_t i = 1; i < _num_subdivisions; ++i) {
+            subdivisions.push_back(a * (qi - 1.0) / (q - 1.0));
+            qi *= q;
+        }
+    }
+
+    subdivisions.push_back(_length);
+
+    return subdivisions;
+}
+
+}  // namespace BaseLib

BaseLib/Subdivision.h

@@ -10,6 +10,7 @@
 #pragma once
 
 #include <cmath>
+#include <vector>
 
 namespace BaseLib
 {
@@ -73,8 +74,7 @@ public:
             const double L,
             const double dL0,
             const double max_dL,
-            const double multiplier)
-    : _length(L), _dL0(dL0), _max_dL(max_dL), _multiplier(multiplier) {}
+            const double multiplier);
 
     /// Returns a vector of subdivided points
     std::vector<double> operator()() const override
@@ -106,4 +106,31 @@ private:
     const double _multiplier;
 };
 
+/**
+ * Gradual subdivision operator with a constant _multiplier.
+ *
+ * In this class the number of subdivisions is known a priori.
+ */
+class GradualSubdivisionFixedNum : public ISubdivision
+{
+public:
+    /**
+     * Constructor
+     * @param L                 total length to be subdivided
+     * @param num_subdivisions  number of subdivisions to generate
+     * @param multiplier        multiplier to cell length
+     */
+    GradualSubdivisionFixedNum(const double L,
+                               const std::size_t num_subdivisions,
+                               const double multiplier);
+
+    /// Returns a vector of subdivided points
+    std::vector<double> operator()() const override;
+
+private:
+    const double _length;
+    const std::size_t _num_subdivisions;
+    const double _multiplier;
+};
+
 } // BaseLib