From 530c6f4fcd90a6efa4e86808cdb38536d281424d Mon Sep 17 00:00:00 2001
From: Guilhem Lavaux <guilhem.lavaux@gmail.com>
Date: Mon, 21 May 2012 19:11:14 -0400
Subject: [PATCH] Few clever algorithms for fast power x^N

---
 sample/testAlgo.cpp | 20 +++++++++++++++
 src/algo.hpp        | 59 +++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 79 insertions(+)
 create mode 100644 sample/testAlgo.cpp
 create mode 100644 src/algo.hpp

diff --git a/sample/testAlgo.cpp b/sample/testAlgo.cpp
new file mode 100644
index 0000000..c7dba56
--- /dev/null
+++ b/sample/testAlgo.cpp
@@ -0,0 +1,20 @@
+#include <iostream>
+#include <stdio.h>
+#include "algo.hpp"
+
+using namespace CosmoTool;
+using namespace std;
+
+int main(int argc, char **argv)
+{
+  cout << square(2) << endl;
+  cout << cube(2) << endl;
+  cout << square(2.1f) << endl;
+  cout << cube(2.1f) << endl;
+  
+  cout << spower<2>(2.1f) << endl;
+  cout << spower<3>(2.1f) << endl;
+  cout << spower<4>(2.1f) << endl;
+  cout << spower<5>(2.1f) << endl;
+  return 0;
+}
diff --git a/src/algo.hpp b/src/algo.hpp
new file mode 100644
index 0000000..9b6e892
--- /dev/null
+++ b/src/algo.hpp
@@ -0,0 +1,59 @@
+#ifndef __COSMOTOOL_ALGO_HPP
+#define __COSMOTOOL_ALGO_HPP
+
+#include "config.hpp"
+
+namespace  CosmoTool
+{
+
+  template<typename T>
+  T cube(T a)
+  {
+    return a*a*a;
+  }
+
+  template<typename T>
+  T square(T a)
+  {
+    return a*a;
+  }
+
+  template<int N>
+  class SPowerBase
+  {
+  public:
+    template<typename T>
+    static T spower(T a)
+    {
+      if ((N%2)==0)
+	{
+	  T b = SPowerBase<N/2>::spower(a);
+	  return b*b;
+	}
+      T b = SPowerBase<(N-1)/2>::spower(a);
+      return a*b*b;
+    }
+  };
+
+  template<>
+  class SPowerBase<0>
+  {
+  public:
+    template<typename T>
+    static T spower(T a)
+    {
+      return T(1);
+    }
+  };
+
+  template<int N, typename T>
+  T spower(T a)
+  {
+    return SPowerBase<N>::spower(a);
+  }
+
+
+
+};
+
+#endif