From 74729ce14f481be69d246c3b965c224ac332cecb Mon Sep 17 00:00:00 2001
From: bolero-MURAKAMI <bolero.murakami@gmail.com>
Date: Sat, 8 Dec 2012 13:45:09 +0900
Subject: [PATCH] fix recursion depth math functions.

---
 sprout/detail/pow.hpp | 15 ++++++++++++---
 sprout/math/asin.hpp  | 35 ++++++++++-------------------------
 sprout/math/atan.hpp  | 30 +++++++++---------------------
 sprout/math/cos.hpp   | 10 +++++-----
 sprout/math/cosh.hpp  | 26 ++++++++++++--------------
 sprout/math/exp.hpp   | 26 ++++++++++++--------------
 sprout/math/log.hpp   | 27 ++++++++-------------------
 sprout/math/sinh.hpp  | 26 ++++++++++++--------------
 8 files changed, 80 insertions(+), 115 deletions(-)

diff --git a/sprout/detail/pow.hpp b/sprout/detail/pow.hpp
index 6afd2be8..020209cd 100644
--- a/sprout/detail/pow.hpp
+++ b/sprout/detail/pow.hpp
@@ -9,17 +9,26 @@ namespace sprout {
 		pow2(T const& x) {
 			return x * x;
 		}
+
 		template<typename T>
 		inline SPROUT_CONSTEXPR T
 		pow3(T const& x) {
 			return x * x * x;
 		}
+
+		template<typename T>
+		inline SPROUT_CONSTEXPR T
+		pow_n_impl(T const& x, int n) {
+			return n == 1 ? x
+				: n % 2 ? x * sprout::detail::pow2(sprout::detail::pow_n_impl(x, n / 2))
+				: sprout::detail::pow2(sprout::detail::pow_n_impl(x, n / 2))
+				;
+		}
 		template<typename T>
 		inline SPROUT_CONSTEXPR T
 		pow_n(T const& x, int n) {
-			return n == 1 ? x
-				: n % 2 ? x * sprout::detail::pow2(sprout::detail::pow_n(x, n / 2))
-				: sprout::detail::pow2(sprout::detail::pow_n(x, n / 2))
+			return n == 0 ? T(1)
+				: sprout::detail::pow_n_impl(x, n)
 				;
 		}
 	}	// namespace detail
diff --git a/sprout/math/asin.hpp b/sprout/math/asin.hpp
index 7c426ba6..2e88708f 100644
--- a/sprout/math/asin.hpp
+++ b/sprout/math/asin.hpp
@@ -5,6 +5,7 @@
 #include <limits>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/constants.hpp>
 #include <sprout/math/factorial.hpp>
@@ -16,37 +17,21 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			asin_impl_2(T x, T tmp, std::size_t n, T x2n1, T _4n) {
-				return 2 * n > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::asin_impl_2(
-						x,
-						tmp + sprout::math::factorial<T>(2 * n)
-							/ _4n / sprout::math::factorial<T>(n) / sprout::math::factorial<T>(n) / (2 * n + 1)
-							* x2n1
-							,
-						n + 1,
-						x2n1 * x * x,
-						_4n * 4
-						)
+			asin_impl_1(T x, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? sprout::math::factorial<T>(2 * n)
+						/ sprout::detail::pow_n(T(4), n) / sprout::detail::pow2(sprout::math::factorial<T>(n)) / (2 * n + 1)
+						* sprout::detail::pow_n(x, 2 * n + 1)
+					: sprout::math::detail::asin_impl_1(x, n, n + (last - n) / 2)
+						+ sprout::math::detail::asin_impl_1(x, n + (last - n) / 2, last)
 					;
 			}
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			asin_impl_1(T x) {
-				return sprout::math::detail::asin_impl_2(
-					x,
-					x,
-					1,
-					x * x * x,
-					T(4)
-					);
-			}
-			template<typename T>
-			inline SPROUT_CONSTEXPR T
 			asin_impl(T x) {
 				return x > sprout::math::half_root_two<T>()
-						? sprout::math::half_pi<T>() - sprout::math::detail::asin_impl_1(sprout::math::sqrt(1 - x * x))
-					: sprout::math::detail::asin_impl_1(x)
+					? sprout::math::half_pi<T>() - sprout::math::detail::asin_impl_1(sprout::math::sqrt(1 - x * x), 0, sprout::math::factorial_limit<T>() / 2 + 1)
+					: x + sprout::math::detail::asin_impl_1(x, 1, sprout::math::factorial_limit<T>() / 2 + 1)
 					;
 			}
 
diff --git a/sprout/math/atan.hpp b/sprout/math/atan.hpp
index 9a5893ad..ad4fa039 100644
--- a/sprout/math/atan.hpp
+++ b/sprout/math/atan.hpp
@@ -4,6 +4,7 @@
 #include <cstddef>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/constants.hpp>
 #include <sprout/math/factorial.hpp>
@@ -14,34 +15,21 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			atan_impl_2(T x, T tmp, std::size_t n, T x2n1) {
-				return n > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::atan_impl_2(
-						x,
-						tmp + (n % 2 ? -1 : 1) * x2n1 / (2 * n + 1),
-						n + 1,
-						x2n1 * x * x
-						)
+			atan_impl_1(T x, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? (n % 2 ? -1 : 1) * sprout::detail::pow_n(x, 2 * n + 1) / (2 * n + 1)
+					: sprout::math::detail::atan_impl_1(x, n, n + (last - n) / 2)
+						+ sprout::math::detail::atan_impl_1(x, n + (last - n) / 2, last)
 					;
 			}
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			atan_impl_1(T x) {
-				return sprout::math::detail::atan_impl_2(
-					x,
-					x,
-					1,
-					x * x * x
-					);
-			}
-			template<typename T>
-			inline SPROUT_CONSTEXPR T
 			atan_impl(T x) {
 				return x > sprout::math::root_two<T>() + 1
-						? sprout::math::half_pi<T>() - sprout::math::detail::atan_impl_1(1 / x)
+						? sprout::math::half_pi<T>() - sprout::math::detail::atan_impl_1(1 / x, 0, sprout::math::factorial_limit<T>() + 1)
 					: x > sprout::math::root_two<T>() - 1
-						? sprout::math::quarter_pi<T>() + sprout::math::detail::atan_impl_1((x - 1) / (x + 1))
-					: sprout::math::detail::atan_impl_1(x)
+						? sprout::math::quarter_pi<T>() + sprout::math::detail::atan_impl_1((x - 1) / (x + 1), 0, sprout::math::factorial_limit<T>() + 1)
+					: x + sprout::math::detail::atan_impl_1(x, 1, sprout::math::factorial_limit<T>() + 1)
 					;
 			}
 
diff --git a/sprout/math/cos.hpp b/sprout/math/cos.hpp
index c3707a8e..001e4822 100644
--- a/sprout/math/cos.hpp
+++ b/sprout/math/cos.hpp
@@ -17,11 +17,11 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			cos_impl_1(T x2, std::size_t first, std::size_t last) {
-				return last - first == 1
-					? (first % 2 ? -1 : 1) * sprout::detail::pow_n(x2, first) / sprout::math::factorial<T>(2 * first)
-					: sprout::math::detail::cos_impl_1(x2, first, first + (last - first) / 2)
-						+ sprout::math::detail::cos_impl_1(x2, first + (last - first) / 2, last)
+			cos_impl_1(T x2, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? (n % 2 ? -1 : 1) * sprout::detail::pow_n(x2, n) / sprout::math::factorial<T>(2 * n)
+					: sprout::math::detail::cos_impl_1(x2, n, n + (last - n) / 2)
+						+ sprout::math::detail::cos_impl_1(x2, n + (last - n) / 2, last)
 					;
 			}
 			template<typename FloatType>
diff --git a/sprout/math/cosh.hpp b/sprout/math/cosh.hpp
index 31b31608..c3c9d185 100644
--- a/sprout/math/cosh.hpp
+++ b/sprout/math/cosh.hpp
@@ -4,6 +4,7 @@
 #include <cstddef>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/factorial.hpp>
 #include <sprout/type_traits/enabler_if.hpp>
@@ -13,14 +14,11 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			cosh_impl(T x, T tmp, std::size_t n, T x2n) {
-				return 2 * n > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::cosh_impl(
-						x,
-						tmp + x2n / sprout::math::factorial<T>(2 * n),
-						n + 1,
-						x2n * x * x
-						)
+			cosh_impl(T x2, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? sprout::detail::pow_n(x2, n) / sprout::math::factorial<T>(2 * n)
+					: sprout::math::detail::cosh_impl(x2, n, n + (last - n) / 2)
+						+ sprout::math::detail::cosh_impl(x2, n + (last - n) / 2, last)
 					;
 			}
 
@@ -31,12 +29,12 @@ namespace sprout {
 			inline SPROUT_CONSTEXPR FloatType
 			cosh(FloatType x) {
 				typedef double type;
-				return static_cast<FloatType>(sprout::math::detail::cosh_impl(
-					static_cast<type>(x),
-					type(1),
-					1,
-					static_cast<type>(x) * static_cast<type>(x)
-					));
+				return static_cast<FloatType>(
+					type(1) + sprout::math::detail::cosh_impl(
+						static_cast<type>(x) * static_cast<type>(x),
+						1, sprout::math::factorial_limit<type>() / 2 + 1
+						)
+					);
 			}
 
 			template<
diff --git a/sprout/math/exp.hpp b/sprout/math/exp.hpp
index a529b012..148f41e3 100644
--- a/sprout/math/exp.hpp
+++ b/sprout/math/exp.hpp
@@ -4,6 +4,7 @@
 #include <cstddef>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/factorial.hpp>
 #include <sprout/type_traits/enabler_if.hpp>
@@ -13,14 +14,11 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			exp_impl(T x, T tmp, std::size_t n, T xn) {
-				return n > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::exp_impl(
-						x,
-						tmp + xn / sprout::math::factorial<T>(n),
-						n + 1,
-						xn * x
-						)
+			exp_impl(T x, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? sprout::detail::pow_n(x, n) / sprout::math::factorial<T>(n)
+					: sprout::math::detail::exp_impl(x, n, n + (last - n) / 2)
+						+ sprout::math::detail::exp_impl(x, n + (last - n) / 2, last)
 					;
 			}
 
@@ -31,12 +29,12 @@ namespace sprout {
 			inline SPROUT_CONSTEXPR FloatType
 			exp(FloatType x) {
 				typedef double type;
-				return static_cast<FloatType>(sprout::math::detail::exp_impl(
-					static_cast<type>(x),
-					type(1),
-					1,
-					static_cast<type>(x)
-					));
+				return static_cast<FloatType>(
+					type(1) + sprout::math::detail::exp_impl(
+						static_cast<type>(x),
+						1, sprout::math::factorial_limit<type>() + 1
+						)
+					);
 			}
 
 			template<
diff --git a/sprout/math/log.hpp b/sprout/math/log.hpp
index 2c4ead9c..09731056 100644
--- a/sprout/math/log.hpp
+++ b/sprout/math/log.hpp
@@ -5,6 +5,7 @@
 #include <limits>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/constants.hpp>
 #include <sprout/math/factorial.hpp>
@@ -16,30 +17,18 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			log_impl_2(T x, T tmp, std::size_t n, T xn) {
-				return n > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::log_impl_2(
-						x,
-						tmp + (n % 2 ? 1 : -1) * xn / n,
-						n + 1,
-						xn * x
-						)
+			log_impl_1(T x, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? (n % 2 ? 1 : -1) * sprout::detail::pow_n(x, n) / n
+					: sprout::math::detail::log_impl_1(x, n, n + (last - n) / 2)
+						+ sprout::math::detail::log_impl_1(x, n + (last - n) / 2, last)
 					;
 			}
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			log_impl_1(T x) {
-				return sprout::math::detail::log_impl_2(
-					x,
-					x,
-					2,
-					x * x
-					);
-			}
-			template<typename T>
-			inline SPROUT_CONSTEXPR T
 			log_impl(T x) {
-				return !(x > sprout::math::root_two<T>()) ? sprout::math::detail::log_impl_1(x - 1)
+				return !(x > sprout::math::root_two<T>())
+					? sprout::math::detail::log_impl_1(x - 1, 1, sprout::math::factorial_limit<T>() + 1)
 					: 2 * sprout::math::detail::log_impl(sprout::math::sqrt(x))
 					;
 			}
diff --git a/sprout/math/sinh.hpp b/sprout/math/sinh.hpp
index 4f276584..29211b66 100644
--- a/sprout/math/sinh.hpp
+++ b/sprout/math/sinh.hpp
@@ -4,6 +4,7 @@
 #include <cstddef>
 #include <type_traits>
 #include <sprout/config.hpp>
+#include <sprout/detail/pow.hpp>
 #include <sprout/math/detail/config.hpp>
 #include <sprout/math/factorial.hpp>
 #include <sprout/type_traits/enabler_if.hpp>
@@ -13,14 +14,11 @@ namespace sprout {
 		namespace detail {
 			template<typename T>
 			inline SPROUT_CONSTEXPR T
-			sinh_impl(T x, T tmp, std::size_t n, T x2n1) {
-				return 2 * n + 1 > sprout::math::factorial_limit<T>() ? tmp
-					: sprout::math::detail::sinh_impl(
-						x,
-						tmp + x2n1 / sprout::math::factorial<T>(2 * n + 1),
-						n + 1,
-						x2n1 * x * x
-						)
+			sinh_impl(T x, std::size_t n, std::size_t last) {
+				return last - n == 1
+					? sprout::detail::pow_n(x, 2 * n + 1) / sprout::math::factorial<T>(2 * n + 1)
+					: sprout::math::detail::sinh_impl(x, n, n + (last - n) / 2)
+						+ sprout::math::detail::sinh_impl(x, n + (last - n) / 2, last)
 					;
 			}
 
@@ -31,12 +29,12 @@ namespace sprout {
 			inline SPROUT_CONSTEXPR FloatType
 			sinh(FloatType x) {
 				typedef double type;
-				return static_cast<FloatType>(sprout::math::detail::sinh_impl(
-					static_cast<type>(x),
-					static_cast<type>(x),
-					1,
-					static_cast<type>(x) * static_cast<type>(x) * static_cast<type>(x)
-					));
+				return static_cast<FloatType>(
+					static_cast<type>(x) + sprout::math::detail::sinh_impl(
+						static_cast<type>(x),
+						1, (sprout::math::factorial_limit<type>() - 1) / 2 + 1
+						)
+					);
 			}
 
 			template<