1 files changed, 144 insertions, 0 deletions
diff --git a/src/utilities/benchmarkUtil.h b/src/utilities/benchmarkUtil.h
new file mode 100644
index 0000000..3c93b05
--- /dev/null
+++ b/src/utilities/benchmarkUtil.h
@@ -0,0 +1,144 @@
+#ifndef BENCHMARK_UTIL_H
+#define BENCHMARK_UTIL_H
+
+#include <deque>
+#include "io/ostreamManager.h"
+#include "functors/analytical/analyticalF.h"
+
+namespace olb {
+
+namespace util {
+
+/// Check time-convergence of a scalar.
+/** This class is useful, for example to check convergence of
+ * the velocity field for the simulation of a stationary flow.
+ * Convergence is claimed when the standard deviation of the
+ * monitored value is smaller than epsilon times the average.
+ * The statistics are taken over a macroscopic time scale of the
+ * system.
+ */
+template<typename T>
+class ValueTracer {
+public:
+  /// Ctor.
+  /** \param u The characteristic velocity of the system, for
+   *          computation of the characteristic time scale.
+   * \param L The characteristic length of the system, for
+   *          computation of the characteristic time scale.
+   * \param epsilon Precision of the convergence.
+   */
+  ValueTracer(T u, T L, T epsilon);
+  /// Ctor.
+  /** \param deltaT   corresponds to iteration steps (averaging)
+  \*  \param epsilon  allowed derivation of quantity
+  */
+  ValueTracer(int deltaT, T epsilon);
+  /// Change values of u and L to update characteristic scales of the system.
+  void resetScale(T u, T L);
+  /// reinitializes the values
+  void resetValues();
+  /// Get characteristic time scale.
+  int getDeltaT() const;
+  /// Feed the object with a new measured scalar.
+  void takeValue(T val, bool doPrint=false);
+  /// Test for convergence, with respect to stdDev.
+  bool hasConverged() const;
+  /// Test for convergence, with respect to difference between min and max value;
+  bool convergenceCheck() const;
+  /// Test for convergence, with respect to difference between min and max value;
+  bool hasConvergedMinMax() const;
+  T computeAverage() const;
+  T computeStdDev(T average) const;
+  void setEpsilon(T epsilon);
+private:
+  int    _deltaT;
+  T      _epsilon;
+  int    _t;
+  bool   _converged;
+  std::deque<T> _values;
+  mutable OstreamManager clout;
+};
+
+/// Propose successive test values of a scalar (e.g. Re) to check stability of a system.
+/** At first, the stability limit is explored by constant
+ * increments/decrements of the scalar, and then, by successive
+ * bisection.
+ */
+template<typename T>
+class BisectStepper {
+public:
+  /// The only constructor.
+  /** \param _iniVal Initial guess for the stability limit.
+   * \param _step   Step size at which the value is initially
+   *                incremented/decremented.
+   */
+  BisectStepper(T _iniVal, T _step=0.);
+  /// Get new value, and indicate if the previous value yielded a stable system or not.
+  T getVal(bool stable, bool doPrint=false);
+  /// Test for convergence.
+  bool hasConverged(T epsilon) const;
+private:
+  T iniVal, currentVal, lowerVal, upperVal;
+  T step;
+  enum {first, up, down, bisect} state;
+  mutable OstreamManager clout;
+};
+
+/// 1D Newton simple scheme
+template<typename T>
+class Newton1D {
+
+protected:
+  AnalyticalF1D<T,T>& _f;
+  AnalyticalDiffFD1D<T> _df;
+  T _yValue;
+  T _eps;
+  int _maxIterations;
+
+public:
+  Newton1D(AnalyticalF1D<T,T>& f, T yValue = T(), T eps = 1.e-8, int maxIterations = 100);
+
+  T solve(T startValue, bool print=false);
+};
+
+/// Trapezoidal rule
+template<typename T>
+class TrapezRuleInt1D {
+
+protected:
+  AnalyticalF1D<T,T>& _f;
+
+public:
+  TrapezRuleInt1D(AnalyticalF1D<T,T>& f);
+
+  T integrate(T min, T max, int nSteps);
+};
+
+/** Simple circular buffer to compute average and other quantities
+ *  over pre-defined temporal windows.
+ *  Works with every T supporting += (T or scalar), and /= (scalar) operations,
+ *  including double and Vector<double, size>.
+ */
+template<typename T>
+class CircularBuffer {
+public:
+  CircularBuffer(int size);
+  /// insert a new entry ed eventually erases the oldest one
+  void insert(T entry);
+  /// average over all the entries
+  T average();
+  /// get reference to the last entry for pos=0, the second-to-last for pos=1, and so on
+  T& get(int pos);
+  /// return size of the buffer
+  int getSize();
+  
+private:
+  int _size;
+  std::vector<T> _data;
+};
+
+} // namespace util
+
+} // namespace olb
+
+#endif