Eigen/src/Core/arch/AltiVec/MathFunctions.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2007 Julien Pommier
5 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
6 // Copyright (C) 2016 Konstantinos Margaritis <markos@freevec.org>
7 //
8 // This Source Code Form is subject to the terms of the Mozilla
9 // Public License v. 2.0. If a copy of the MPL was not distributed
10 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
11 
12 /* The sin, cos, exp, and log functions of this file come from
13  * Julien Pommier's sse math library: http://gruntthepeon.free.fr/ssemath/
14  */
15 
16 #ifndef EIGEN_MATH_FUNCTIONS_ALTIVEC_H
17 #define EIGEN_MATH_FUNCTIONS_ALTIVEC_H
18 
19 namespace Eigen {
20 
21 namespace internal {
22 
23 static _EIGEN_DECLARE_CONST_Packet4f(1 , 1.0f);
24 static _EIGEN_DECLARE_CONST_Packet4f(half, 0.5f);
25 static _EIGEN_DECLARE_CONST_Packet4i(0x7f, 0x7f);
26 static _EIGEN_DECLARE_CONST_Packet4i(23, 23);
27 
28 static _EIGEN_DECLARE_CONST_Packet4f_FROM_INT(inv_mant_mask, ~0x7f800000);
29 
30 /* the smallest non denormalized float number */
31 static _EIGEN_DECLARE_CONST_Packet4f_FROM_INT(min_norm_pos, 0x00800000);
32 static _EIGEN_DECLARE_CONST_Packet4f_FROM_INT(minus_inf, 0xff800000); // -1.f/0.f
33 static _EIGEN_DECLARE_CONST_Packet4f_FROM_INT(minus_nan, 0xffffffff);
34 
35 /* natural logarithm computed for 4 simultaneous float
36  return NaN for x <= 0
37 */
38 static _EIGEN_DECLARE_CONST_Packet4f(cephes_SQRTHF, 0.707106781186547524f);
39 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p0, 7.0376836292E-2f);
40 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p1, - 1.1514610310E-1f);
41 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p2, 1.1676998740E-1f);
42 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p3, - 1.2420140846E-1f);
43 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p4, + 1.4249322787E-1f);
44 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p5, - 1.6668057665E-1f);
45 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p6, + 2.0000714765E-1f);
46 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p7, - 2.4999993993E-1f);
47 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_p8, + 3.3333331174E-1f);
48 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_q1, -2.12194440e-4f);
49 static _EIGEN_DECLARE_CONST_Packet4f(cephes_log_q2, 0.693359375f);
50 
51 static _EIGEN_DECLARE_CONST_Packet4f(exp_hi, 88.3762626647950f);
52 static _EIGEN_DECLARE_CONST_Packet4f(exp_lo, -88.3762626647949f);
53 
54 static _EIGEN_DECLARE_CONST_Packet4f(cephes_LOG2EF, 1.44269504088896341f);
55 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_C1, 0.693359375f);
56 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_C2, -2.12194440e-4f);
57 
58 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p0, 1.9875691500E-4f);
59 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p1, 1.3981999507E-3f);
60 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p2, 8.3334519073E-3f);
61 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p3, 4.1665795894E-2f);
62 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p4, 1.6666665459E-1f);
63 static _EIGEN_DECLARE_CONST_Packet4f(cephes_exp_p5, 5.0000001201E-1f);
64 
65 #ifdef __VSX__
66 static _EIGEN_DECLARE_CONST_Packet2d(1 , 1.0);
67 static _EIGEN_DECLARE_CONST_Packet2d(2 , 2.0);
68 static _EIGEN_DECLARE_CONST_Packet2d(half, 0.5);
69 
70 static _EIGEN_DECLARE_CONST_Packet2d(exp_hi, 709.437);
71 static _EIGEN_DECLARE_CONST_Packet2d(exp_lo, -709.436139303);
72 
73 static _EIGEN_DECLARE_CONST_Packet2d(cephes_LOG2EF, 1.4426950408889634073599);
74 
75 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_p0, 1.26177193074810590878e-4);
76 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_p1, 3.02994407707441961300e-2);
77 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_p2, 9.99999999999999999910e-1);
78 
79 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_q0, 3.00198505138664455042e-6);
80 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_q1, 2.52448340349684104192e-3);
81 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_q2, 2.27265548208155028766e-1);
82 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_q3, 2.00000000000000000009e0);
83 
84 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_C1, 0.693145751953125);
85 static _EIGEN_DECLARE_CONST_Packet2d(cephes_exp_C2, 1.42860682030941723212e-6);
86 
87 #ifdef __POWER8_VECTOR__
88 static Packet2l p2l_1023 = { 1023, 1023 };
89 static Packet2ul p2ul_52 = { 52, 52 };
90 #endif
91 
92 #endif
93 
96 {
97  Packet4f x = _x;
98 
99  Packet4i emm0;
100 
101  /* isvalid_mask is 0 if x < 0 or x is NaN. */
102  Packet4ui isvalid_mask = reinterpret_cast<Packet4ui>(vec_cmpge(x, p4f_ZERO));
103  Packet4ui iszero_mask = reinterpret_cast<Packet4ui>(vec_cmpeq(x, p4f_ZERO));
104 
105  x = pmax(x, p4f_min_norm_pos); /* cut off denormalized stuff */
106  emm0 = vec_sr(reinterpret_cast<Packet4i>(x),
107  reinterpret_cast<Packet4ui>(p4i_23));
108 
109  /* keep only the fractional part */
110  x = pand(x, p4f_inv_mant_mask);
111  x = por(x, p4f_half);
112 
113  emm0 = psub(emm0, p4i_0x7f);
114  Packet4f e = padd(vec_ctf(emm0, 0), p4f_1);
115 
116  /* part2:
117  if( x < SQRTHF ) {
118  e -= 1;
119  x = x + x - 1.0;
120  } else { x = x - 1.0; }
121  */
122  Packet4f mask = reinterpret_cast<Packet4f>(vec_cmplt(x, p4f_cephes_SQRTHF));
123  Packet4f tmp = pand(x, mask);
124  x = psub(x, p4f_1);
125  e = psub(e, pand(p4f_1, mask));
126  x = padd(x, tmp);
127 
128  Packet4f x2 = pmul(x,x);
129  Packet4f x3 = pmul(x2,x);
130 
131  Packet4f y, y1, y2;
132  y = pmadd(p4f_cephes_log_p0, x, p4f_cephes_log_p1);
133  y1 = pmadd(p4f_cephes_log_p3, x, p4f_cephes_log_p4);
134  y2 = pmadd(p4f_cephes_log_p6, x, p4f_cephes_log_p7);
135  y = pmadd(y , x, p4f_cephes_log_p2);
136  y1 = pmadd(y1, x, p4f_cephes_log_p5);
137  y2 = pmadd(y2, x, p4f_cephes_log_p8);
138  y = pmadd(y, x3, y1);
139  y = pmadd(y, x3, y2);
140  y = pmul(y, x3);
141 
142  y1 = pmul(e, p4f_cephes_log_q1);
143  tmp = pmul(x2, p4f_half);
144  y = padd(y, y1);
145  x = psub(x, tmp);
146  y2 = pmul(e, p4f_cephes_log_q2);
147  x = padd(x, y);
148  x = padd(x, y2);
149  // negative arg will be NAN, 0 will be -INF
150  x = vec_sel(x, p4f_minus_inf, iszero_mask);
151  x = vec_sel(p4f_minus_nan, x, isvalid_mask);
152  return x;
153 }
154 
157 {
158  Packet4f x = _x;
159 
160  Packet4f tmp, fx;
161  Packet4i emm0;
162 
163  // clamp x
164  x = pmax(pmin(x, p4f_exp_hi), p4f_exp_lo);
165 
166  // express exp(x) as exp(g + n*log(2))
167  fx = pmadd(x, p4f_cephes_LOG2EF, p4f_half);
168 
169  fx = pfloor(fx);
170 
171  tmp = pmul(fx, p4f_cephes_exp_C1);
172  Packet4f z = pmul(fx, p4f_cephes_exp_C2);
173  x = psub(x, tmp);
174  x = psub(x, z);
175 
176  z = pmul(x,x);
177 
178  Packet4f y = p4f_cephes_exp_p0;
179  y = pmadd(y, x, p4f_cephes_exp_p1);
180  y = pmadd(y, x, p4f_cephes_exp_p2);
181  y = pmadd(y, x, p4f_cephes_exp_p3);
182  y = pmadd(y, x, p4f_cephes_exp_p4);
183  y = pmadd(y, x, p4f_cephes_exp_p5);
184  y = pmadd(y, z, x);
185  y = padd(y, p4f_1);
186 
187  // build 2^n
188  emm0 = vec_cts(fx, 0);
189  emm0 = vec_add(emm0, p4i_0x7f);
190  emm0 = vec_sl(emm0, reinterpret_cast<Packet4ui>(p4i_23));
191 
192  // Altivec's max & min operators just drop silent NaNs. Check NaNs in
193  // inputs and return them unmodified.
194  Packet4ui isnumber_mask = reinterpret_cast<Packet4ui>(vec_cmpeq(_x, _x));
195  return vec_sel(_x, pmax(pmul(y, reinterpret_cast<Packet4f>(emm0)), _x),
196  isnumber_mask);
197 }
198 
199 #ifndef EIGEN_COMP_CLANG
202 {
203  return vec_rsqrt(x);
204 }
205 #endif
206 
207 #ifdef __VSX__
208 #ifndef EIGEN_COMP_CLANG
211 {
212  return vec_rsqrt(x);
213 }
214 #endif
215 
218 {
219  return vec_sqrt(x);
220 }
221 
224 {
225  return vec_sqrt(x);
226 }
227 
228 // VSX support varies between different compilers and even different
229 // versions of the same compiler. For gcc version >= 4.9.3, we can use
230 // vec_cts to efficiently convert Packet2d to Packet2l. Otherwise, use
231 // a slow version that works with older compilers.
232 // Update: apparently vec_cts/vec_ctf intrinsics for 64-bit doubles
233 // are buggy, https://gcc.gnu.org/bugzilla/show_bug.cgi?id=70963
234 static inline Packet2l ConvertToPacket2l(const Packet2d& x) {
235 #if EIGEN_GNUC_AT_LEAST(5, 4) || \
236  (EIGEN_GNUC_AT(6, 1) && __GNUC_PATCHLEVEL__ >= 1)
237  return vec_cts(x, 0); // TODO: check clang version.
238 #else
239  double tmp[2];
240  memcpy(tmp, &x, sizeof(tmp));
241  Packet2l l = { static_cast<long long>(tmp[0]),
242  static_cast<long long>(tmp[1]) };
243  return l;
244 #endif
245 }
246 
249 {
250  Packet2d x = _x;
251 
252  Packet2d tmp, fx;
253  Packet2l emm0;
254 
255  // clamp x
256  x = pmax(pmin(x, p2d_exp_hi), p2d_exp_lo);
257 
258  /* express exp(x) as exp(g + n*log(2)) */
259  fx = pmadd(x, p2d_cephes_LOG2EF, p2d_half);
260 
261  fx = pfloor(fx);
262 
263  tmp = pmul(fx, p2d_cephes_exp_C1);
264  Packet2d z = pmul(fx, p2d_cephes_exp_C2);
265  x = psub(x, tmp);
266  x = psub(x, z);
267 
268  Packet2d x2 = pmul(x,x);
269 
270  Packet2d px = p2d_cephes_exp_p0;
271  px = pmadd(px, x2, p2d_cephes_exp_p1);
272  px = pmadd(px, x2, p2d_cephes_exp_p2);
273  px = pmul (px, x);
274 
275  Packet2d qx = p2d_cephes_exp_q0;
276  qx = pmadd(qx, x2, p2d_cephes_exp_q1);
277  qx = pmadd(qx, x2, p2d_cephes_exp_q2);
278  qx = pmadd(qx, x2, p2d_cephes_exp_q3);
279 
280  x = pdiv(px,psub(qx,px));
281  x = pmadd(p2d_2,x,p2d_1);
282 
283  // build 2^n
284  emm0 = ConvertToPacket2l(fx);
285 
286 #ifdef __POWER8_VECTOR__
287  emm0 = vec_add(emm0, p2l_1023);
288  emm0 = vec_sl(emm0, p2ul_52);
289 #else
290  // Code is a bit complex for POWER7. There is actually a
291  // vec_xxsldi intrinsic but it is not supported by some gcc versions.
292  // So we shift (52-32) bits and do a word swap with zeros.
293  _EIGEN_DECLARE_CONST_Packet4i(1023, 1023);
294  _EIGEN_DECLARE_CONST_Packet4i(20, 20); // 52 - 32
295 
296  Packet4i emm04i = reinterpret_cast<Packet4i>(emm0);
297  emm04i = vec_add(emm04i, p4i_1023);
298  emm04i = vec_sl(emm04i, reinterpret_cast<Packet4ui>(p4i_20));
299  static const Packet16uc perm = {
300  0x14, 0x15, 0x16, 0x17, 0x00, 0x01, 0x02, 0x03,
301  0x1c, 0x1d, 0x1e, 0x1f, 0x08, 0x09, 0x0a, 0x0b };
302 #ifdef _BIG_ENDIAN
303  emm0 = reinterpret_cast<Packet2l>(vec_perm(p4i_ZERO, emm04i, perm));
304 #else
305  emm0 = reinterpret_cast<Packet2l>(vec_perm(emm04i, p4i_ZERO, perm));
306 #endif
307 
308 #endif
309 
310  // Altivec's max & min operators just drop silent NaNs. Check NaNs in
311  // inputs and return them unmodified.
312  Packet2ul isnumber_mask = reinterpret_cast<Packet2ul>(vec_cmpeq(_x, _x));
313  return vec_sel(_x, pmax(pmul(x, reinterpret_cast<Packet2d>(emm0)), _x),
314  isnumber_mask);
315 }
316 #endif
317 
318 } // end namespace internal
319 
320 } // end namespace Eigen
321 
322 #endif // EIGEN_MATH_FUNCTIONS_ALTIVEC_H
static _EIGEN_DECLARE_CONST_Packet4i(0x7f, 0x7f)
__vector unsigned char Packet16uc
EIGEN_DEVICE_FUNC Packet por(const Packet &a, const Packet &b)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet2d psqrt< Packet2d >(const Packet2d &x)
static int f(const TensorMap< Tensor< int, 3 > > &tensor)
__vector unsigned int Packet4ui
Definition: LDLT.h:16
__vector int Packet4i
#define EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS
Definition: Macros.h:529
__vector long long Packet2l
EIGEN_DEVICE_FUNC Packet padd(const Packet &a, const Packet &b)
EIGEN_DEVICE_FUNC Packet pmin(const Packet &a, const Packet &b)
#define EIGEN_UNUSED
Definition: Macros.h:607
__vector unsigned long long Packet2ul
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet4f plog< Packet4f >(const Packet4f &_x)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet4f prsqrt< Packet4f >(const Packet4f &x)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet2d prsqrt< Packet2d >(const Packet2d &x)
TFSIMD_FORCE_INLINE const tfScalar & z() const
EIGEN_DEVICE_FUNC Packet pdiv(const Packet &a, const Packet &b)
EIGEN_DECLARE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS Packet pfloor(const Packet &a)
EIGEN_STRONG_INLINE Packet4f pmadd(const Packet4f &a, const Packet4f &b, const Packet4f &c)
EIGEN_DEVICE_FUNC Packet psub(const Packet &a, const Packet &b)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet4f psqrt< Packet4f >(const Packet4f &x)
static _EIGEN_DECLARE_CONST_Packet4f(1, 1.0f)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet4f pexp< Packet4f >(const Packet4f &_x)
EIGEN_DEFINE_FUNCTION_ALLOWING_MULTIPLE_DEFINITIONS EIGEN_UNUSED Packet2d pexp< Packet2d >(const Packet2d &_x)
EIGEN_DEVICE_FUNC Packet pmul(const Packet &a, const Packet &b)
static _EIGEN_DECLARE_CONST_Packet2d(1, 1.0)
EIGEN_DEVICE_FUNC Packet pmax(const Packet &a, const Packet &b)
EIGEN_DEVICE_FUNC Packet pand(const Packet &a, const Packet &b)
static _EIGEN_DECLARE_CONST_Packet4f_FROM_INT(inv_mant_mask,~0x7f800000)
const T & y


hebiros
Author(s): Xavier Artache , Matthew Tesch
autogenerated on Thu Sep 3 2020 04:08:24