1 /*################################################################################
2 ##
3 ## Copyright (C) 2016-2020 Keith O'Hara
4 ##
5 ## This file is part of the GCE-Math C++ library.
6 ##
7 ## Licensed under the Apache License, Version 2.0 (the "License");
8 ## you may not use this file except in compliance with the License.
9 ## You may obtain a copy of the License at
10 ##
11 ## http://www.apache.org/licenses/LICENSE-2.0
12 ##
13 ## Unless required by applicable law or agreed to in writing, software
14 ## distributed under the License is distributed on an "AS IS" BASIS,
15 ## WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
16 ## See the License for the specific language governing permissions and
17 ## limitations under the License.
18 ##
19 ################################################################################*/
20
21 /*
22 * compile-time power function
23 */
24
25 #ifndef _gcem_pow_integral_HPP
26 #define _gcem_pow_integral_HPP
27
28 namespace internal
29 {
30
31 template<typename T1, typename T2>
32 constexpr T1 pow_integral_compute(const T1 base, const T2 exp_term) noexcept;
33
34 // integral-valued powers using method described in
35 // https://en.wikipedia.org/wiki/Exponentiation_by_squaring
36
37 template<typename T1, typename T2>
38 constexpr
39 T1
pow_integral_compute_recur(const T1 base,const T1 val,const T2 exp_term)40 pow_integral_compute_recur(const T1 base, const T1 val, const T2 exp_term)
41 noexcept
42 {
43 return( exp_term > T2(1) ? \
44 (is_odd(exp_term) ? \
45 pow_integral_compute_recur(base*base,val*base,exp_term/2) :
46 pow_integral_compute_recur(base*base,val,exp_term/2)) :
47 (exp_term == T2(1) ? val*base : val) );
48 }
49
50 template<typename T1, typename T2, typename std::enable_if<std::is_signed<T2>::value>::type* = nullptr>
51 constexpr
52 T1
pow_integral_sgn_check(const T1 base,const T2 exp_term)53 pow_integral_sgn_check(const T1 base, const T2 exp_term)
54 noexcept
55 {
56 return( exp_term < T2(0) ? \
57 //
58 T1(1) / pow_integral_compute(base, - exp_term) :
59 //
60 pow_integral_compute_recur(base,T1(1),exp_term) );
61 }
62
63 template<typename T1, typename T2, typename std::enable_if<!std::is_signed<T2>::value>::type* = nullptr>
64 constexpr
65 T1
pow_integral_sgn_check(const T1 base,const T2 exp_term)66 pow_integral_sgn_check(const T1 base, const T2 exp_term)
67 noexcept
68 {
69 return( pow_integral_compute_recur(base,T1(1),exp_term) );
70 }
71
72 template<typename T1, typename T2>
73 constexpr
74 T1
pow_integral_compute(const T1 base,const T2 exp_term)75 pow_integral_compute(const T1 base, const T2 exp_term)
76 noexcept
77 {
78 return( exp_term == T2(3) ? \
79 base*base*base :
80 exp_term == T2(2) ? \
81 base*base :
82 exp_term == T2(1) ? \
83 base :
84 exp_term == T2(0) ? \
85 T1(1) :
86 // check for overflow
87 exp_term == GCLIM<T2>::min() ? \
88 T1(0) :
89 exp_term == GCLIM<T2>::max() ? \
90 GCLIM<T1>::infinity() :
91 // else
92 pow_integral_sgn_check(base,exp_term) );
93 }
94
95 template<typename T1, typename T2, typename std::enable_if<std::is_integral<T2>::value>::type* = nullptr>
96 constexpr
97 T1
pow_integral_type_check(const T1 base,const T2 exp_term)98 pow_integral_type_check(const T1 base, const T2 exp_term)
99 noexcept
100 {
101 return pow_integral_compute(base,exp_term);
102 }
103
104 template<typename T1, typename T2, typename std::enable_if<!std::is_integral<T2>::value>::type* = nullptr>
105 constexpr
106 T1
pow_integral_type_check(const T1 base,const T2 exp_term)107 pow_integral_type_check(const T1 base, const T2 exp_term)
108 noexcept
109 {
110 // return GCLIM<return_t<T1>>::quiet_NaN();
111 return pow_integral_compute(base,static_cast<llint_t>(exp_term));
112 }
113
114 //
115 // main function
116
117 template<typename T1, typename T2>
118 constexpr
119 T1
pow_integral(const T1 base,const T2 exp_term)120 pow_integral(const T1 base, const T2 exp_term)
121 noexcept
122 {
123 return internal::pow_integral_type_check(base,exp_term);
124 }
125
126 }
127
128 #endif
129