Utilities::fixed_power - use exponentiation by squaring

diff --git a/include/deal.II/base/utilities.h b/include/deal.II/base/utilities.h
index 00442d9d9d63014804e6ae327f5af1445e557365..4c9ae69d84053fad60faf2535cd0e094cb746bc0 100644 (file)
--- a/include/deal.II/base/utilities.h
+++ b/include/deal.II/base/utilities.h
@@ -1030,27 +1030,15 @@ namespace Utilities
 {
   template <int N, typename T>
   inline T
-  fixed_power(const T n)
+  fixed_power(const T x)
   {
-    Assert(N >= 0, ExcNotImplemented());
-    switch (N)
-      {
-        case 0:
-          return dealii::internal::NumberType<T>::value(1);
-        case 1:
-          return n;
-        case 2:
-          return n * n;
-        case 3:
-          return n * n * n;
-        case 4:
-          return n * n * n * n;
-        default:
-          T result = n;
-          for (int d = 1; d < N; ++d)
-            result *= n;
-          return result;
-      }
+    if (N == 0)
+      return T(1.);
+    else if (N < 0)
+      return T(1.) / fixed_power<-N>(x);
+    else
+      // Use exponentiation by squaring:
+      return ((N % 2 == 1) ? x * fixed_power<N / 2>(x * x) : fixed_power<N / 2>(x * x));
   }