'power' method with integer exponent usable in compile-time expressions

Previously, the file 'power.hh' contained two template metaprograms
that computed integer powers and were usable in compile-time
expressions.  With the advent of constexpr these metaprograms can
be replaced by normal methods, leading to much code simplification.

CHANGELOG.md

@@ -11,6 +11,10 @@
   implementation in `dune-common` is automatically disabled, and the official
   implementation from the standard library is used instead.
 
+- There is now (finally!) a method `power` in the file `math.hh` that computes
+  powers with an integer exponent, and is usable in compile-time expressions.
+  The use of the old power methods in `power.hh` is henceforth discouraged.
+
 -   When run with an absolute build directory, `dunecontrol` now exposes the root build
     directory to CMake in the variable `DUNE_BUILD_DIRECTORY_ROOT_PATH`.
 

dune/common/math.hh

@@ -65,6 +65,26 @@ namespace Dune
   {};
 
 
+  /** \brief Power method for integer exponents
+   *
+   * \note Make sure that Mantissa is a non-integer type when using negative exponents!
+   */
+  template <class Mantissa, class Exponent>
+  constexpr Mantissa power(Mantissa m, Exponent p)
+  {
+    static_assert(std::numeric_limits<Exponent>::is_integer, "Exponent must be an integer type!");
+
+    auto result = Mantissa(1);
+    auto absp = (p<0) ? -p : p;   // This is simply abs, but std::abs is not constexpr
+    for (Exponent i = Exponent(0); i<absp; i++)
+      result *= m;
+
+    if (p<0)
+      result = Mantissa(1)/result;
+
+    return result;
+  }
+
   //! Calculates the factorial of m at compile time
   template <int m>
   struct Factorial

dune/common/test/powertest.cc

 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
 // vi: set et ts=4 sw=2 sts=2:
 
+#include <config.h>
+
 #include <iostream>
 
 #include <dune/common/exceptions.hh>
+#include <dune/common/math.hh>
 #include <dune/common/power.hh>
+#include <dune/common/unused.hh>
 
 using namespace Dune;
 
 int main (int argc, char** argv) try
 {
+  // Zero and positive powers
+  if (power(4,0) != 1)
+    DUNE_THROW(MathError, "power(4,0) does not compute the correct result");
+
+  if (power(4,1) != 4)
+    DUNE_THROW(MathError, "power(4,1) implementation does not compute the correct result");
+
+  if (power(4,2) != 16)
+    DUNE_THROW(MathError, "power(4,2) implementation does not compute the correct result");
+
+  if (power(4,3) != 64)
+    DUNE_THROW(MathError, "power(4,3) implementation does not compute the correct result");
+
+  // Negative powers
+  if (power(4.0,-1) != 0.25)
+    DUNE_THROW(MathError, "power(4,-1) implementation does not compute the correct result");
+
+  if (power(4.0,-2) != 0.0625)
+    DUNE_THROW(MathError, "power(4,-2) implementation does not compute the correct result");
+
+  if (power(4.0,-3) != 0.015625)
+    DUNE_THROW(MathError, "power(4,-3) implementation does not compute the correct result");
+
+  // Test whether the result can be used in a compile-time expression
+  enum { dummy = power(2,2) };
+
+  // Test legacy power implementation
   if (Power<0>::eval(4) != 1)
     DUNE_THROW(MathError, "Power implementation does not compute the correct result");