← Back to branch summary

~ubuntu-branches/ubuntu/utopic/eglibc/utopic

« back to all changes in this revision

Viewing changes to sysdeps/ia64/fpu/e_pow.S

  • Committer: Package Import Robot
  • Author(s): Adam Conrad
  • Date: 2012-10-26 05:14:58 UTC
  • mfrom: (1.5.1) (4.4.22 experimental)
  • Revision ID: package-import@ubuntu.com-20121026051458-oryotr4i03ob5pab
Tags: 2.16-0ubuntu1
* Merge with unreleased 2.16 in Debian experimental, remaining changes:
  - Drop the Breaks line from libc6, which refers to a Debian transition
  - Remove the libc6 recommends on libc6-i686, which we don't build
  - Enable libc6{,-dev}-armel on armhf and libc6{-dev}-armhf on armel
  - Ship update-locale and validlocale in /usr/sbin in libc-bin
  - Don't build locales or locales-all in Ubuntu, we rely on langpacks
  - Heavily mangle the way we do service restarting on major upgrades
  - Use different MIN_KERNEL_SUPPORTED versions than Debian, due to
    buildd needs.  This should be universally bumped to 3.2.0 once all
    our buildds (including the PPA guests) are running precise kernels
  - Build i386 variants as -march=i686, build amd64 with -O3, and build
    ppc64 variants (both 64-bit and 32-bit) with -O3 -fno-tree-vectorize
  - Re-enable unsubmitted-ldconfig-cache-abi.diff and rebuild the cache
    on upgrades from previous versions that used a different constant
  - debian/patches/any/local-CVE-2012-3406.diff: switch to malloc when
    array grows too large to handle via alloca extension (CVE-2012-3406)
  - Build generic i386/i686 flavour with -mno-tls-direct-seg-refs
* Changes added/dropped with this merge while reducing our delta:
  - Stop building glibc docs from the eglibc source, and instead make
    the glibc-docs stub have a hard dependency on glibc-doc-reference
  - Remove outdated conflicts against ancient versions of ia32-libs
  - Drop the tzdata dependency from libc6, it's in required and minimal
  - Use gcc-4.7/g++-4.7 by default on all our supported architectures
  - Save our historical changelog as changelog.ubuntu in the source
  - Drop nscd's libaudit build-dep for now, as libaudit is in universe
  - Drop the unnecessary Breaks from libc6 to locales and locales-all
  - Ship xen's ld.so.conf.d snippet as /etc/ld.so.conf.d/libc6-xen.conf
* Disable hard failures on the test suite for the first upload to raring

Show diffs side-by-side

added added

removed removed

Lines of Context:
sysdeps/ia64/fpu/e_pow.S
1
 
.file "pow.s"
2
 
 
3
 
 
4
 
// Copyright (c) 2000 - 2005, Intel Corporation
5
 
// All rights reserved.
6
 
//
7
 
// Contributed 2000 by the Intel Numerics Group, Intel Corporation
8
 
//
9
 
// Redistribution and use in source and binary forms, with or without
10
 
// modification, are permitted provided that the following conditions are
11
 
// met:
12
 
//
13
 
// * Redistributions of source code must retain the above copyright
14
 
// notice, this list of conditions and the following disclaimer.
15
 
//
16
 
// * Redistributions in binary form must reproduce the above copyright
17
 
// notice, this list of conditions and the following disclaimer in the
18
 
// documentation and/or other materials provided with the distribution.
19
 
//
20
 
// * The name of Intel Corporation may not be used to endorse or promote
21
 
// products derived from this software without specific prior written
22
 
// permission.
23
 
 
24
 
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
25
 
// "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
26
 
// LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
27
 
// A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL INTEL OR ITS
28
 
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
29
 
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
30
 
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
31
 
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
32
 
// OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY OR TORT (INCLUDING
33
 
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
34
 
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
35
 
//
36
 
// Intel Corporation is the author of this code, and requests that all
37
 
// problem reports or change requests be submitted to it directly at
38
 
// http://www.intel.com/software/products/opensource/libraries/num.htm.
39
 
//
40
 
// History
41
 
//==============================================================
42
 
// 02/02/00 Initial version
43
 
// 02/03/00 Added p12 to definite over/under path. With odd power we did not
44
 
//          maintain the sign of x in this path.
45
 
// 04/04/00 Unwind support added
46
 
// 04/19/00 pow(+-1,inf) now returns NaN
47
 
//          pow(+-val, +-inf) returns 0 or inf, but now does not call error
48
 
//          support
49
 
//          Added s1 to fcvt.fx because invalid flag was incorrectly set.
50
 
// 08/15/00 Bundle added after call to __libm_error_support to properly
51
 
//          set [the previously overwritten] GR_Parameter_RESULT.
52
 
// 09/07/00 Improved performance by eliminating bank conflicts and other stalls,
53
 
//          and tweaking the critical path
54
 
// 09/08/00 Per c99, pow(+-1,inf) now returns 1, and pow(+1,nan) returns 1
55
 
// 09/28/00 Updated NaN**0 path
56
 
// 01/20/01 Fixed denormal flag settings.
57
 
// 02/13/01 Improved speed.
58
 
// 03/19/01 Reordered exp polynomial to improve speed and eliminate monotonicity
59
 
//          problem in round up, down, and to zero modes.  Also corrected
60
 
//          overflow result when x negative, y odd in round up, down, zero.
61
 
// 06/14/01 Added brace missing from bundle
62
 
// 12/10/01 Corrected case where x negative, 2^52 <= |y| < 2^53, y odd integer.
63
 
// 12/20/01 Fixed monotonity problem in round to nearest.
64
 
// 02/08/02 Fixed overflow/underflow cases that were not calling error support.
65
 
// 05/20/02 Cleaned up namespace and sf0 syntax
66
 
// 08/29/02 Improved Itanium 2 performance
67
 
// 09/21/02 Added branch for |y*log(x)|<2^-11 to fix monotonicity problems.
68
 
// 02/10/03 Reordered header: .section, .global, .proc, .align
69
 
// 03/31/05 Reformatted delimiters between data tables
70
 
//
71
 
// API
72
 
//==============================================================
73
 
// double pow(double x, double y)
74
 
//
75
 
// Overview of operation
76
 
//==============================================================
77
 
//
78
 
// Three steps...
79
 
// 1. Log(x)
80
 
// 2. y Log(x)
81
 
// 3. exp(y log(x))
82
 
//
83
 
// This means we work with the absolute value of x and merge in the sign later.
84
 
//      Log(x) = G + delta + r -rsq/2 + p
85
 
// G,delta depend on the exponent of x and table entries. The table entries are
86
 
// indexed by the exponent of x, called K.
87
 
//
88
 
// The G and delta come out of the reduction; r is the reduced x.
89
 
//
90
 
// B = frcpa(x)
91
 
// xB-1 is small means that B is the approximate inverse of x.
92
 
//
93
 
//      Log(x) = Log( (1/B)(Bx) )
94
 
//             = Log(1/B) + Log(Bx)
95
 
//             = Log(1/B) + Log( 1 + (Bx-1))
96
 
//
97
 
//      x  = 2^K 1.x_1x_2.....x_52
98
 
//      B= frcpa(x) = 2^-k Cm
99
 
//      Log(1/B) = Log(1/(2^-K Cm))
100
 
//      Log(1/B) = Log((2^K/ Cm))
101
 
//      Log(1/B) = K Log(2) + Log(1/Cm)
102
 
//
103
 
//      Log(x)   = K Log(2) + Log(1/Cm) + Log( 1 + (Bx-1))
104
 
//
105
 
// If you take the significand of x, set the exponent to true 0, then Cm is
106
 
// the frcpa. We tabulate the Log(1/Cm) values. There are 256 of them.
107
 
// The frcpa table is indexed by 8 bits, the x_1 thru x_8.
108
 
// m = x_1x_2...x_8 is an 8-bit index.
109
 
//
110
 
//      Log(1/Cm) = log(1/frcpa(1+m/256)) where m goes from 0 to 255.
111
 
//
112
 
// We tabluate as two doubles, T and t, where T +t is the value itself.
113
 
//
114
 
//      Log(x)   = (K Log(2)_hi + T) + (Log(2)_hi + t) + Log( 1 + (Bx-1))
115
 
//      Log(x)   =  G + delta           + Log( 1 + (Bx-1))
116
 
//
117
 
// The Log( 1 + (Bx-1)) can be calculated as a series in r = Bx-1.
118
 
//
119
 
//      Log( 1 + (Bx-1)) = r - rsq/2 + p
120
 
//
121
 
// Then,
122
 
//
123
 
//      yLog(x) = yG + y delta + y(r-rsq/2) + yp
124
 
//      yLog(x) = Z1 + e3      + Z2         + Z3 + (e2 + e3)
125
 
//
126
 
//
127
 
//     exp(yLog(x)) = exp(Z1 + Z2 + Z3) exp(e1 + e2 + e3)
128
 
//
129
 
//
130
 
//       exp(Z3) is another series.
131
 
//       exp(e1 + e2 + e3) is approximated as f3 = 1 + (e1 + e2 + e3)
132
 
//
133
 
//       Z1 (128/log2) = number of log2/128 in Z1 is N1
134
 
//       Z2 (128/log2) = number of log2/128 in Z2 is N2
135
 
//
136
 
//       s1 = Z1 - N1 log2/128
137
 
//       s2 = Z2 - N2 log2/128
138
 
//
139
 
//       s = s1 + s2
140
 
//       N = N1 + N2
141
 
//
142
 
//       exp(Z1 + Z2) = exp(Z)
143
 
//       exp(Z)       = exp(s) exp(N log2/128)
144
 
//
145
 
//       exp(r)       = exp(Z - N log2/128)
146
 
//
147
 
//      r = s + d = (Z - N (log2/128)_hi) -N (log2/128)_lo
148
 
//                =  Z - N (log2/128)
149
 
//
150
 
//      Z         = s+d +N (log2/128)
151
 
//
152
 
//      exp(Z)    = exp(s) (1+d) exp(N log2/128)
153
 
//
154
 
//      N = M 128 + n
155
 
//
156
 
//      N log2/128 = M log2 + n log2/128
157
 
//
158
 
//      n is 8 binary digits = n_7n_6...n_1
159
 
//
160
 
//      n log2/128 = n_7n_6n_5 16 log2/128 + n_4n_3n_2n_1 log2/128
161
 
//      n log2/128 = n_7n_6n_5 log2/8 + n_4n_3n_2n_1 log2/128
162
 
//      n log2/128 = I2 log2/8 + I1 log2/128
163
 
//
164
 
//      N log2/128 = M log2 + I2 log2/8 + I1 log2/128
165
 
//
166
 
//      exp(Z)    = exp(s) (1+d) exp(log(2^M) + log(2^I2/8) + log(2^I1/128))
167
 
//      exp(Z)    = exp(s) (1+d1) (1+d2)(2^M) 2^I2/8 2^I1/128
168
 
//      exp(Z)    = exp(s) f1 f2 (2^M) 2^I2/8 2^I1/128
169
 
//
170
 
// I1, I2 are table indices. Use a series for exp(s).
171
 
// Then get exp(Z)
172
 
//
173
 
//     exp(yLog(x)) = exp(Z1 + Z2 + Z3) exp(e1 + e2 + e3)
174
 
//     exp(yLog(x)) = exp(Z) exp(Z3) f3
175
 
//     exp(yLog(x)) = exp(Z)f3 exp(Z3)
176
 
//     exp(yLog(x)) = A exp(Z3)
177
 
//
178
 
// We actually calculate exp(Z3) -1.
179
 
// Then,
180
 
//     exp(yLog(x)) = A + A( exp(Z3)   -1)
181
 
//
182
 
 
183
 
// Table Generation
184
 
//==============================================================
185
 
 
186
 
// The log values
187
 
// ==============
188
 
// The operation (K*log2_hi) must be exact. K is the true exponent of x.
189
 
// If we allow gradual underflow (denormals), K can be represented in 12 bits
190
 
// (as a two's complement number). We assume 13 bits as an engineering
191
 
// precaution.
192
 
//
193
 
//           +------------+----------------+-+
194
 
//           |  13 bits   | 50 bits        | |
195
 
//           +------------+----------------+-+
196
 
//           0            1                66
197
 
//                        2                34
198
 
//
199
 
// So we want the lsb(log2_hi) to be 2^-50
200
 
// We get log2 as a quad-extended (15-bit exponent, 128-bit significand)
201
 
//
202
 
//      0 fffe b17217f7d1cf79ab c9e3b39803f2f6af (4...)
203
 
//
204
 
// Consider numbering the bits left to right, starting at 0 thru 127.
205
 
// Bit 0 is the 2^-1 bit; bit 49 is the 2^-50 bit.
206
 
//
207
 
//  ...79ab
208
 
//     0111 1001 1010 1011
209
 
//     44
210
 
//     89
211
 
//
212
 
// So if we shift off the rightmost 14 bits, then (shift back only
213
 
// the top half) we get
214
 
//
215
 
//      0 fffe b17217f7d1cf4000 e6af278ece600fcb dabc000000000000
216
 
//
217
 
// Put the right 64-bit signficand in an FR register, convert to double;
218
 
// it is exact. Put the next 128 bits into a quad register and round to double.
219
 
// The true exponent of the low part is -51.
220
 
//
221
 
// hi is 0 fffe b17217f7d1cf4000
222
 
// lo is 0 ffcc e6af278ece601000
223
 
//
224
 
// Convert to double memory format and get
225
 
//
226
 
// hi is 0x3fe62e42fefa39e8
227
 
// lo is 0x3cccd5e4f1d9cc02
228
 
//
229
 
// log2_hi + log2_lo is an accurate value for log2.
230
 
//
231
 
//
232
 
// The T and t values
233
 
// ==================
234
 
// A similar method is used to generate the T and t values.
235
 
//
236
 
// K * log2_hi + T  must be exact.
237
 
//
238
 
// Smallest T,t
239
 
// ----------
240
 
// The smallest T,t is
241
 
//       T                   t
242
 
// 0x3f60040155d58800, 0x3c93bce0ce3ddd81  log(1/frcpa(1+0/256))=  +1.95503e-003
243
 
//
244
 
// The exponent is 0x3f6 (biased)  or -9 (true).
245
 
// For the smallest T value, what we want is to clip the significand such that
246
 
// when it is shifted right by 9, its lsb is in the bit for 2^-51. The 9 is the
247
 
// specific for the first entry. In general, it is 0xffff - (biased 15-bit
248
 
// exponent).
249
 
 
250
 
// Independently, what we have calculated is the table value as a quad
251
 
// precision number.
252
 
// Table entry 1 is
253
 
// 0 fff6 80200aaeac44ef38 338f77605fdf8000
254
 
//
255
 
// We store this quad precision number in a data structure that is
256
 
//    sign:           1
257
 
//    exponent:      15
258
 
//    signficand_hi: 64 (includes explicit bit)
259
 
//    signficand_lo: 49
260
 
// Because the explicit bit is included, the significand is 113 bits.
261
 
//
262
 
// Consider significand_hi for table entry 1.
263
 
//
264
 
//
265
 
// +-+--- ... -------+--------------------+
266
 
// | |
267
 
// +-+--- ... -------+--------------------+
268
 
// 0 1               4444444455555555556666
269
 
//                   2345678901234567890123
270
 
//
271
 
// Labeled as above, bit 0 is 2^0, bit 1 is 2^-1, etc.
272
 
// Bit 42 is 2^-42. If we shift to the right by 9, the bit in
273
 
// bit 42 goes in 51.
274
 
//
275
 
// So what we want to do is shift bits 43 thru 63 into significand_lo.
276
 
// This is shifting bit 42 into bit 63, taking care to retain shifted-off bits.
277
 
// Then shifting (just with signficaand_hi) back into bit 42.
278
 
//
279
 
// The shift_value is 63-42 = 21. In general, this is
280
 
//      63 - (51 -(0xffff - 0xfff6))
281
 
// For this example, it is
282
 
//      63 - (51 - 9) = 63 - 42  = 21
283
 
//
284
 
// This means we are shifting 21 bits into significand_lo. We must maintain more
285
 
// that a 128-bit signficand not to lose bits. So before the shift we put the
286
 
// 128-bit significand into a 256-bit signficand and then shift.
287
 
// The 256-bit significand has four parts: hh, hl, lh, and ll.
288
 
//
289
 
// Start off with
290
 
//      hh         hl         lh         ll
291
 
//      <64>       <49><15_0> <64_0>     <64_0>
292
 
//
293
 
// After shift by 21 (then return for significand_hi),
294
 
//      <43><21_0> <21><43>   <6><58_0>  <64_0>
295
 
//
296
 
// Take the hh part and convert to a double. There is no rounding here.
297
 
// The conversion is exact. The true exponent of the high part is the same as
298
 
// the true exponent of the input quad.
299
 
//
300
 
// We have some 64 plus significand bits for the low part. In this example, we
301
 
// have 70 bits. We want to round this to a double. Put them in a quad and then
302
 
// do a quad fnorm.
303
 
// For this example the true exponent of the low part is
304
 
//      true_exponent_of_high - 43 = true_exponent_of_high - (64-21)
305
 
// In general, this is
306
 
//      true_exponent_of_high - (64 - shift_value)
307
 
//
308
 
//
309
 
// Largest T,t
310
 
// ----------
311
 
// The largest T,t is
312
 
// 0x3fe62643fecf9742, 0x3c9e3147684bd37d  log(1/frcpa(1+255/256))=+6.92171e-001
313
 
//
314
 
// Table entry 256 is
315
 
// 0 fffe b1321ff67cba178c 51da12f4df5a0000
316
 
//
317
 
// The shift value is
318
 
//      63 - (51 -(0xffff - 0xfffe)) = 13
319
 
//
320
 
// The true exponent of the low part is
321
 
//      true_exponent_of_high - (64 - shift_value)
322
 
//      -1 - (64-13) = -52
323
 
// Biased as a double, this is 0x3cb
324
 
//
325
 
//
326
 
//
327
 
// So then lsb(T) must be >= 2^-51
328
 
// msb(Klog2_hi) <= 2^12
329
 
//
330
 
//              +--------+---------+
331
 
//              |       51 bits    | <== largest T
332
 
//              +--------+---------+
333
 
//              | 9 bits | 42 bits | <== smallest T
334
 
// +------------+----------------+-+
335
 
// |  13 bits   | 50 bits        | |
336
 
// +------------+----------------+-+
337
 
 
338
 
 
339
 
// Special Cases
340
 
//==============================================================
341
 
 
342
 
//                                   double     float
343
 
// overflow                          error 24   30
344
 
 
345
 
// underflow                         error 25   31
346
 
 
347
 
// X zero  Y zero
348
 
//  +0     +0                 +1     error 26   32
349
 
//  -0     +0                 +1     error 26   32
350
 
//  +0     -0                 +1     error 26   32
351
 
//  -0     -0                 +1     error 26   32
352
 
 
353
 
// X zero  Y negative
354
 
//  +0     -odd integer       +inf   error 27   33  divide-by-zero
355
 
//  -0     -odd integer       -inf   error 27   33  divide-by-zero
356
 
//  +0     !-odd integer      +inf   error 27   33  divide-by-zero
357
 
//  -0     !-odd integer      +inf   error 27   33  divide-by-zero
358
 
//  +0     -inf               +inf   error 27   33  divide-by-zero
359
 
//  -0     -inf               +inf   error 27   33  divide-by-zero
360
 
 
361
 
// X zero  Y positve
362
 
//  +0     +odd integer       +0
363
 
//  -0     +odd integer       -0
364
 
//  +0     !+odd integer      +0
365
 
//  -0     !+odd integer      +0
366
 
//  +0     +inf               +0
367
 
//  -0     +inf               +0
368
 
//  +0     Y NaN              quiet Y               invalid if Y SNaN
369
 
//  -0     Y NaN              quiet Y               invalid if Y SNaN
370
 
 
371
 
// X one
372
 
//  -1     Y inf              +1
373
 
//  -1     Y NaN              quiet Y               invalid if Y SNaN
374
 
//  +1     Y NaN              +1                    invalid if Y SNaN
375
 
//  +1     Y any else         +1
376
 
 
377
 
// X -     Y not integer      QNAN   error 28   34  invalid
378
 
 
379
 
// X NaN   Y 0                +1     error 29   35
380
 
// X NaN   Y NaN              quiet X               invalid if X or Y SNaN
381
 
// X NaN   Y any else         quiet X               invalid if X SNaN
382
 
// X !+1   Y NaN              quiet Y               invalid if Y SNaN
383
 
 
384
 
 
385
 
// X +inf  Y >0               +inf
386
 
// X -inf  Y >0, !odd integer +inf
387
 
// X -inf  Y >0, odd integer  -inf
388
 
 
389
 
// X +inf  Y <0               +0
390
 
// X -inf  Y <0, !odd integer +0
391
 
// X -inf  Y <0, odd integer  -0
392
 
 
393
 
// X +inf  Y =0               +1
394
 
// X -inf  Y =0               +1
395
 
 
396
 
// |X|<1   Y +inf             +0
397
 
// |X|<1   Y -inf             +inf
398
 
// |X|>1   Y +inf             +inf
399
 
// |X|>1   Y -inf             +0
400
 
 
401
 
// X any   Y =0               +1
402
 
 
403
 
// Assembly macros
404
 
//==============================================================
405
 
 
406
 
// integer registers used
407
 
 
408
 
pow_GR_signexp_X          = r14
409
 
pow_GR_17ones             = r15
410
 
pow_AD_P                  = r16
411
 
pow_GR_exp_2tom8          = r17
412
 
pow_GR_sig_X              = r18
413
 
pow_GR_10033              = r19
414
 
pow_GR_16ones             = r20
415
 
 
416
 
pow_AD_Tt                 = r21
417
 
pow_GR_exp_X              = r22
418
 
pow_AD_Q                  = r23
419
 
pow_GR_true_exp_X         = r24
420
 
pow_GR_y_zero             = r25
421
 
 
422
 
pow_GR_exp_Y              = r26
423
 
pow_AD_tbl1               = r27
424
 
pow_AD_tbl2               = r28
425
 
pow_GR_offset             = r29
426
 
pow_GR_exp_Xm1            = r30
427
 
pow_GR_xneg_yodd          = r31
428
 
 
429
 
pow_GR_signexp_Xm1        = r35
430
 
pow_GR_int_W1             = r36
431
 
pow_GR_int_W2             = r37
432
 
pow_GR_int_N              = r38
433
 
pow_GR_index1             = r39
434
 
pow_GR_index2             = r40
435
 
 
436
 
pow_AD_T1                 = r41
437
 
pow_AD_T2                 = r42
438
 
pow_int_GR_M              = r43
439
 
pow_GR_sig_int_Y          = r44
440
 
pow_GR_sign_Y_Gpr         = r45
441
 
 
442
 
pow_GR_17ones_m1          = r46
443
 
pow_GR_one                = r47
444
 
pow_GR_sign_Y             = r48
445
 
pow_GR_signexp_Y_Gpr      = r49
446
 
pow_GR_exp_Y_Gpr          = r50
447
 
 
448
 
pow_GR_true_exp_Y_Gpr     = r51
449
 
pow_GR_signexp_Y          = r52
450
 
pow_GR_x_one              = r53
451
 
pow_GR_exp_2toM63         = r54
452
 
pow_GR_big_pos            = r55
453
 
 
454
 
pow_GR_big_neg            = r56
455
 
 
456
 
GR_SAVE_B0                = r50
457
 
GR_SAVE_GP                = r51
458
 
GR_SAVE_PFS               = r52
459
 
 
460
 
GR_Parameter_X            = r53
461
 
GR_Parameter_Y            = r54
462
 
GR_Parameter_RESULT       = r55
463
 
pow_GR_tag                = r56
464
 
 
465
 
 
466
 
// floating point registers used
467
 
 
468
 
POW_B                     = f32
469
 
POW_NORM_X                = f33
470
 
POW_Xm1                   = f34
471
 
POW_r1                    = f34
472
 
POW_P4                    = f35
473
 
 
474
 
POW_P5                    = f36
475
 
POW_NORM_Y                = f37
476
 
POW_Q2                    = f38
477
 
POW_Q3                    = f39
478
 
POW_P2                    = f40
479
 
 
480
 
POW_P3                    = f41
481
 
POW_P0                    = f42
482
 
POW_log2_lo               = f43
483
 
POW_r                     = f44
484
 
POW_Q0_half               = f45
485
 
 
486
 
POW_Q1                    = f46
487
 
POW_tmp                   = f47
488
 
POW_log2_hi               = f48
489
 
POW_Q4                    = f49
490
 
POW_P1                    = f50
491
 
 
492
 
POW_log2_by_128_hi        = f51
493
 
POW_inv_log2_by_128       = f52
494
 
POW_rsq                   = f53
495
 
POW_Yrcub                 = f54
496
 
POW_log2_by_128_lo        = f55
497
 
 
498
 
POW_v6                    = f56
499
 
POW_xsq                   = f57
500
 
POW_v4                    = f58
501
 
POW_v2                    = f59
502
 
POW_T                     = f60
503
 
 
504
 
POW_Tt                    = f61
505
 
POW_RSHF                  = f62
506
 
POW_v21ps                 = f63
507
 
POW_s4                    = f64
508
 
POW_twoV                  = f65
509
 
 
510
 
POW_U                     = f66
511
 
POW_G                     = f67
512
 
POW_delta                 = f68
513
 
POW_v3                    = f69
514
 
POW_V                     = f70
515
 
 
516
 
POW_p                     = f71
517
 
POW_Z1                    = f72
518
 
POW_e3                    = f73
519
 
POW_e2                    = f74
520
 
POW_Z2                    = f75
521
 
 
522
 
POW_e1                    = f76
523
 
POW_W1                    = f77
524
 
POW_UmZ2                  = f78
525
 
POW_W2                    = f79
526
 
POW_Z3                    = f80
527
 
 
528
 
POW_int_W1                = f81
529
 
POW_e12                   = f82
530
 
POW_int_W2                = f83
531
 
POW_UmZ2pV                = f84
532
 
POW_Z3sq                  = f85
533
 
 
534
 
POW_e123                  = f86
535
 
POW_N1float               = f87
536
 
POW_N2float               = f88
537
 
POW_f3                    = f89
538
 
POW_q                     = f90
539
 
 
540
 
POW_s1                    = f91
541
 
POW_Nfloat                = f92
542
 
POW_s2                    = f93
543
 
POW_f2                    = f94
544
 
POW_f1                    = f95
545
 
 
546
 
POW_T1                    = f96
547
 
POW_T2                    = f97
548
 
POW_2M                    = f98
549
 
POW_s                     = f99
550
 
POW_f12                   = f100
551
 
 
552
 
POW_ssq                   = f101
553
 
POW_T1T2                  = f102
554
 
POW_1ps                   = f103
555
 
POW_A                     = f104
556
 
POW_es                    = f105
557
 
 
558
 
POW_Xp1                   = f106
559
 
POW_int_K                 = f107
560
 
POW_K                     = f108
561
 
POW_f123                  = f109
562
 
POW_Gpr                   = f110
563
 
 
564
 
POW_Y_Gpr                 = f111
565
 
POW_int_Y                 = f112
566
 
POW_abs_q                 = f114
567
 
POW_2toM63                = f115
568
 
 
569
 
POW_float_int_Y           = f116
570
 
POW_ftz_urm_f8            = f117
571
 
POW_wre_urm_f8            = f118
572
 
POW_big_neg               = f119
573
 
POW_big_pos               = f120
574
 
 
575
 
POW_GY_Z2                 = f121
576
 
POW_pYrcub_e3             = f122
577
 
POW_d                     = f123
578
 
POW_d2                    = f124
579
 
POW_poly_d_hi             = f121
580
 
POW_poly_d_lo             = f122
581
 
POW_poly_d                = f121
582
 
 
583
 
// Data tables
584
 
//==============================================================
585
 
 
586
 
RODATA
587
 
 
588
 
.align 16
589
 
 
590
 
LOCAL_OBJECT_START(pow_table_P)
591
 
data8 0x8000F7B249FF332D, 0x0000BFFC  // P_5
592
 
data8 0xAAAAAAA9E7902C7F, 0x0000BFFC  // P_3
593
 
data8 0x80000000000018E5, 0x0000BFFD  // P_1
594
 
data8 0xb8aa3b295c17f0bc, 0x00004006  // inv_ln2_by_128
595
 
//
596
 
//
597
 
data8 0x3FA5555555554A9E // Q_2
598
 
data8 0x3F8111124F4DD9F9 // Q_3
599
 
data8 0x3FE0000000000000 // Q_0
600
 
data8 0x3FC5555555554733 // Q_1
601
 
data8 0x3F56C16D9360FFA0 // Q_4
602
 
data8 0x43e8000000000000 // Right shift constant for exp
603
 
data8 0xc9e3b39803f2f6af, 0x00003fb7  // ln2_by_128_lo
604
 
data8 0x0000000000000000 // pad to eliminate bank conflicts with pow_table_Q
605
 
data8 0x0000000000000000 // pad to eliminate bank conflicts with pow_table_Q
606
 
LOCAL_OBJECT_END(pow_table_P)
607
 
 
608
 
LOCAL_OBJECT_START(pow_table_Q)
609
 
data8 0x9249FE7F0DC423CF, 0x00003FFC  // P_4
610
 
data8 0xCCCCCCCC4ED2BA7F, 0x00003FFC  // P_2
611
 
data8 0xAAAAAAAAAAAAB505, 0x00003FFD  // P_0
612
 
data8 0x3fe62e42fefa39e8, 0x3cccd5e4f1d9cc02 // log2 hi lo =  +6.93147e-001
613
 
data8 0xb17217f7d1cf79ab, 0x00003ff7  // ln2_by_128_hi
614
 
LOCAL_OBJECT_END(pow_table_Q)
615
 
 
616
 
 
617
 
LOCAL_OBJECT_START(pow_Tt)
618
 
data8 0x3f60040155d58800, 0x3c93bce0ce3ddd81 // log(1/frcpa(1+0/256))=  +1.95503e-003
619
 
data8 0x3f78121214586a00, 0x3cb540e0a5cfc9bc // log(1/frcpa(1+1/256))=  +5.87661e-003
620
 
data8 0x3f841929f9683200, 0x3cbdf1d57404da1f // log(1/frcpa(1+2/256))=  +9.81362e-003
621
 
data8 0x3f8c317384c75f00, 0x3c69806208c04c22 // log(1/frcpa(1+3/256))=  +1.37662e-002
622
 
data8 0x3f91a6b91ac73380, 0x3c7874daa716eb32 // log(1/frcpa(1+4/256))=  +1.72376e-002
623
 
data8 0x3f95ba9a5d9ac000, 0x3cacbb84e08d78ac // log(1/frcpa(1+5/256))=  +2.12196e-002
624
 
data8 0x3f99d2a807432580, 0x3cbcf80538b441e1 // log(1/frcpa(1+6/256))=  +2.52177e-002
625
 
data8 0x3f9d6b2725979800, 0x3c6095e5c8f8f359 // log(1/frcpa(1+7/256))=  +2.87291e-002
626
 
data8 0x3fa0c58fa19dfa80, 0x3cb4c5d4e9d0dda2 // log(1/frcpa(1+8/256))=  +3.27573e-002
627
 
data8 0x3fa2954c78cbce00, 0x3caa932b860ab8d6 // log(1/frcpa(1+9/256))=  +3.62953e-002
628
 
data8 0x3fa4a94d2da96c40, 0x3ca670452b76bbd5 // log(1/frcpa(1+10/256))=  +4.03542e-002
629
 
data8 0x3fa67c94f2d4bb40, 0x3ca84104f9941798 // log(1/frcpa(1+11/256))=  +4.39192e-002
630
 
data8 0x3fa85188b630f040, 0x3cb40a882cbf0153 // log(1/frcpa(1+12/256))=  +4.74971e-002
631
 
data8 0x3faa6b8abe73af40, 0x3c988d46e25c9059 // log(1/frcpa(1+13/256))=  +5.16017e-002
632
 
data8 0x3fac441e06f72a80, 0x3cae3e930a1a2a96 // log(1/frcpa(1+14/256))=  +5.52072e-002
633
 
data8 0x3fae1e6713606d00, 0x3c8a796f6283b580 // log(1/frcpa(1+15/256))=  +5.88257e-002
634
 
data8 0x3faffa6911ab9300, 0x3c5193070351e88a // log(1/frcpa(1+16/256))=  +6.24574e-002
635
 
data8 0x3fb0ec139c5da600, 0x3c623f2a75eb992d // log(1/frcpa(1+17/256))=  +6.61022e-002
636
 
data8 0x3fb1dbd2643d1900, 0x3ca649b2ef8927f0 // log(1/frcpa(1+18/256))=  +6.97605e-002
637
 
data8 0x3fb2cc7284fe5f00, 0x3cbc5e86599513e2 // log(1/frcpa(1+19/256))=  +7.34321e-002
638
 
data8 0x3fb3bdf5a7d1ee60, 0x3c90bd4bb69dada3 // log(1/frcpa(1+20/256))=  +7.71173e-002
639
 
data8 0x3fb4b05d7aa012e0, 0x3c54e377c9b8a54f // log(1/frcpa(1+21/256))=  +8.08161e-002
640
 
data8 0x3fb580db7ceb5700, 0x3c7fdb2f98354cde // log(1/frcpa(1+22/256))=  +8.39975e-002
641
 
data8 0x3fb674f089365a60, 0x3cb9994c9d3301c1 // log(1/frcpa(1+23/256))=  +8.77219e-002
642
 
data8 0x3fb769ef2c6b5680, 0x3caaec639db52a79 // log(1/frcpa(1+24/256))=  +9.14602e-002
643
 
data8 0x3fb85fd927506a40, 0x3c9f9f99a3cf8e25 // log(1/frcpa(1+25/256))=  +9.52125e-002
644
 
data8 0x3fb9335e5d594980, 0x3ca15c3abd47d99a // log(1/frcpa(1+26/256))=  +9.84401e-002
645
 
data8 0x3fba2b0220c8e5e0, 0x3cb4ca639adf6fc3 // log(1/frcpa(1+27/256))=  +1.02219e-001
646
 
data8 0x3fbb0004ac1a86a0, 0x3ca7cb81bf959a59 // log(1/frcpa(1+28/256))=  +1.05469e-001
647
 
data8 0x3fbbf968769fca00, 0x3cb0c646c121418e // log(1/frcpa(1+29/256))=  +1.09274e-001
648
 
data8 0x3fbccfedbfee13a0, 0x3ca0465fce24ab4b // log(1/frcpa(1+30/256))=  +1.12548e-001
649
 
data8 0x3fbda727638446a0, 0x3c82803f4e2e6603 // log(1/frcpa(1+31/256))=  +1.15832e-001
650
 
data8 0x3fbea3257fe10f60, 0x3cb986a3f2313d1a // log(1/frcpa(1+32/256))=  +1.19677e-001
651
 
data8 0x3fbf7be9fedbfde0, 0x3c97d16a6a621cf4 // log(1/frcpa(1+33/256))=  +1.22985e-001
652
 
data8 0x3fc02ab352ff25f0, 0x3c9cc6baad365600 // log(1/frcpa(1+34/256))=  +1.26303e-001
653
 
data8 0x3fc097ce579d2040, 0x3cb9ba16d329440b // log(1/frcpa(1+35/256))=  +1.29633e-001
654
 
data8 0x3fc1178e8227e470, 0x3cb7bc671683f8e6 // log(1/frcpa(1+36/256))=  +1.33531e-001
655
 
data8 0x3fc185747dbecf30, 0x3c9d1116f66d2345 // log(1/frcpa(1+37/256))=  +1.36885e-001
656
 
data8 0x3fc1f3b925f25d40, 0x3c8162c9ef939ac6 // log(1/frcpa(1+38/256))=  +1.40250e-001
657
 
data8 0x3fc2625d1e6ddf50, 0x3caad3a1ec384fc3 // log(1/frcpa(1+39/256))=  +1.43627e-001
658
 
data8 0x3fc2d1610c868130, 0x3cb3ad997036941b // log(1/frcpa(1+40/256))=  +1.47015e-001
659
 
data8 0x3fc340c597411420, 0x3cbc2308262c7998 // log(1/frcpa(1+41/256))=  +1.50414e-001
660
 
data8 0x3fc3b08b6757f2a0, 0x3cb2170d6cdf0526 // log(1/frcpa(1+42/256))=  +1.53825e-001
661
 
data8 0x3fc40dfb08378000, 0x3c9bb453c4f7b685 // log(1/frcpa(1+43/256))=  +1.56677e-001
662
 
data8 0x3fc47e74e8ca5f70, 0x3cb836a48fdfce9d // log(1/frcpa(1+44/256))=  +1.60109e-001
663
 
data8 0x3fc4ef51f6466de0, 0x3ca07a43919aa64b // log(1/frcpa(1+45/256))=  +1.63553e-001
664
 
data8 0x3fc56092e02ba510, 0x3ca85006899d97b0 // log(1/frcpa(1+46/256))=  +1.67010e-001
665
 
data8 0x3fc5d23857cd74d0, 0x3ca30a5ba6e7abbe // log(1/frcpa(1+47/256))=  +1.70478e-001
666
 
data8 0x3fc6313a37335d70, 0x3ca905586f0ac97e // log(1/frcpa(1+48/256))=  +1.73377e-001
667
 
data8 0x3fc6a399dabbd380, 0x3c9b2c6657a96684 // log(1/frcpa(1+49/256))=  +1.76868e-001
668
 
data8 0x3fc70337dd3ce410, 0x3cb50bc52f55cdd8 // log(1/frcpa(1+50/256))=  +1.79786e-001
669
 
data8 0x3fc77654128f6120, 0x3cad2eb7c9a39efe // log(1/frcpa(1+51/256))=  +1.83299e-001
670
 
data8 0x3fc7e9d82a0b0220, 0x3cba127e90393c01 // log(1/frcpa(1+52/256))=  +1.86824e-001
671
 
data8 0x3fc84a6b759f5120, 0x3cbd7fd52079f706 // log(1/frcpa(1+53/256))=  +1.89771e-001
672
 
data8 0x3fc8ab47d5f5a300, 0x3cbfae141751a3de // log(1/frcpa(1+54/256))=  +1.92727e-001
673
 
data8 0x3fc91fe490965810, 0x3cb69cf30a1c319e // log(1/frcpa(1+55/256))=  +1.96286e-001
674
 
data8 0x3fc981634011aa70, 0x3ca5bb3d208bc42a // log(1/frcpa(1+56/256))=  +1.99261e-001
675
 
data8 0x3fc9f6c407089660, 0x3ca04d68658179a0 // log(1/frcpa(1+57/256))=  +2.02843e-001
676
 
data8 0x3fca58e729348f40, 0x3c99f5411546c286 // log(1/frcpa(1+58/256))=  +2.05838e-001
677
 
data8 0x3fcabb55c31693a0, 0x3cb9a5350eb327d5 // log(1/frcpa(1+59/256))=  +2.08842e-001
678
 
data8 0x3fcb1e104919efd0, 0x3c18965fcce7c406 // log(1/frcpa(1+60/256))=  +2.11855e-001
679
 
data8 0x3fcb94ee93e367c0, 0x3cb503716da45184 // log(1/frcpa(1+61/256))=  +2.15483e-001
680
 
data8 0x3fcbf851c0675550, 0x3cbdf1b3f7ab5378 // log(1/frcpa(1+62/256))=  +2.18516e-001
681
 
data8 0x3fcc5c0254bf23a0, 0x3ca7aab9ed0b1d7b // log(1/frcpa(1+63/256))=  +2.21558e-001
682
 
data8 0x3fccc000c9db3c50, 0x3c92a7a2a850072a // log(1/frcpa(1+64/256))=  +2.24609e-001
683
 
data8 0x3fcd244d99c85670, 0x3c9f6019120edf4c // log(1/frcpa(1+65/256))=  +2.27670e-001
684
 
data8 0x3fcd88e93fb2f450, 0x3c6affb96815e081 // log(1/frcpa(1+66/256))=  +2.30741e-001
685
 
data8 0x3fcdedd437eaef00, 0x3c72553595897976 // log(1/frcpa(1+67/256))=  +2.33820e-001
686
 
data8 0x3fce530effe71010, 0x3c90913b020fa182 // log(1/frcpa(1+68/256))=  +2.36910e-001
687
 
data8 0x3fceb89a1648b970, 0x3c837ba4045bfd25 // log(1/frcpa(1+69/256))=  +2.40009e-001
688
 
data8 0x3fcf1e75fadf9bd0, 0x3cbcea6d13e0498d // log(1/frcpa(1+70/256))=  +2.43117e-001
689
 
data8 0x3fcf84a32ead7c30, 0x3ca5e3a67b3c6d77 // log(1/frcpa(1+71/256))=  +2.46235e-001
690
 
data8 0x3fcfeb2233ea07c0, 0x3cba0c6f0049c5a6 // log(1/frcpa(1+72/256))=  +2.49363e-001
691
 
data8 0x3fd028f9c7035c18, 0x3cb0a30b06677ff6 // log(1/frcpa(1+73/256))=  +2.52501e-001
692
 
data8 0x3fd05c8be0d96358, 0x3ca0f1c77ccb5865 // log(1/frcpa(1+74/256))=  +2.55649e-001
693
 
data8 0x3fd085eb8f8ae790, 0x3cbd513f45fe7a97 // log(1/frcpa(1+75/256))=  +2.58174e-001
694
 
data8 0x3fd0b9c8e32d1910, 0x3c927449047ca006 // log(1/frcpa(1+76/256))=  +2.61339e-001
695
 
data8 0x3fd0edd060b78080, 0x3c89b52d8435f53e // log(1/frcpa(1+77/256))=  +2.64515e-001
696
 
data8 0x3fd122024cf00638, 0x3cbdd976fabda4bd // log(1/frcpa(1+78/256))=  +2.67701e-001
697
 
data8 0x3fd14be2927aecd0, 0x3cb02f90ad0bc471 // log(1/frcpa(1+79/256))=  +2.70257e-001
698
 
data8 0x3fd180618ef18ad8, 0x3cbd003792c71a98 // log(1/frcpa(1+80/256))=  +2.73461e-001
699
 
data8 0x3fd1b50bbe2fc638, 0x3ca9ae64c6403ead // log(1/frcpa(1+81/256))=  +2.76675e-001
700
 
data8 0x3fd1df4cc7cf2428, 0x3cb43f0455f7e395 // log(1/frcpa(1+82/256))=  +2.79254e-001
701
 
data8 0x3fd214456d0eb8d0, 0x3cb0fbd748d75d30 // log(1/frcpa(1+83/256))=  +2.82487e-001
702
 
data8 0x3fd23ec5991eba48, 0x3c906edd746b77e2 // log(1/frcpa(1+84/256))=  +2.85081e-001
703
 
data8 0x3fd2740d9f870af8, 0x3ca9802e6a00a670 // log(1/frcpa(1+85/256))=  +2.88333e-001
704
 
data8 0x3fd29ecdabcdfa00, 0x3cacecef70890cfa // log(1/frcpa(1+86/256))=  +2.90943e-001
705
 
data8 0x3fd2d46602adcce8, 0x3cb97911955f3521 // log(1/frcpa(1+87/256))=  +2.94214e-001
706
 
data8 0x3fd2ff66b04ea9d0, 0x3cb12dabe191d1c9 // log(1/frcpa(1+88/256))=  +2.96838e-001
707
 
data8 0x3fd335504b355a30, 0x3cbdf9139df924ec // log(1/frcpa(1+89/256))=  +3.00129e-001
708
 
data8 0x3fd360925ec44f58, 0x3cb253e68977a1e3 // log(1/frcpa(1+90/256))=  +3.02769e-001
709
 
data8 0x3fd38bf1c3337e70, 0x3cb3d283d2a2da21 // log(1/frcpa(1+91/256))=  +3.05417e-001
710
 
data8 0x3fd3c25277333180, 0x3cadaa5b035eae27 // log(1/frcpa(1+92/256))=  +3.08735e-001
711
 
data8 0x3fd3edf463c16838, 0x3cb983d680d3c108 // log(1/frcpa(1+93/256))=  +3.11399e-001
712
 
data8 0x3fd419b423d5e8c0, 0x3cbc86dd921c139d // log(1/frcpa(1+94/256))=  +3.14069e-001
713
 
data8 0x3fd44591e0539f48, 0x3c86a76d6dc2782e // log(1/frcpa(1+95/256))=  +3.16746e-001
714
 
data8 0x3fd47c9175b6f0a8, 0x3cb59a2e013c6b5f // log(1/frcpa(1+96/256))=  +3.20103e-001
715
 
data8 0x3fd4a8b341552b08, 0x3c93f1e86e468694 // log(1/frcpa(1+97/256))=  +3.22797e-001
716
 
data8 0x3fd4d4f390890198, 0x3cbf5e4ea7c5105a // log(1/frcpa(1+98/256))=  +3.25498e-001
717
 
data8 0x3fd501528da1f960, 0x3cbf58da53e9ad10 // log(1/frcpa(1+99/256))=  +3.28206e-001
718
 
data8 0x3fd52dd06347d4f0, 0x3cb98a28cebf6eef // log(1/frcpa(1+100/256))=  +3.30921e-001
719
 
data8 0x3fd55a6d3c7b8a88, 0x3c9c76b67c2d1fd4 // log(1/frcpa(1+101/256))=  +3.33644e-001
720
 
data8 0x3fd5925d2b112a58, 0x3c9029616a4331b8 // log(1/frcpa(1+102/256))=  +3.37058e-001
721
 
data8 0x3fd5bf406b543db0, 0x3c9fb8292ecfc820 // log(1/frcpa(1+103/256))=  +3.39798e-001
722
 
data8 0x3fd5ec433d5c35a8, 0x3cb71a1229d17eec // log(1/frcpa(1+104/256))=  +3.42545e-001
723
 
data8 0x3fd61965cdb02c18, 0x3cbba94fe1dbb8d2 // log(1/frcpa(1+105/256))=  +3.45300e-001
724
 
data8 0x3fd646a84935b2a0, 0x3c9ee496d2c9ae57 // log(1/frcpa(1+106/256))=  +3.48063e-001
725
 
data8 0x3fd6740add31de90, 0x3cb1da3a6c7a9dfd // log(1/frcpa(1+107/256))=  +3.50833e-001
726
 
data8 0x3fd6a18db74a58c0, 0x3cb494c257add8dc // log(1/frcpa(1+108/256))=  +3.53610e-001
727
 
data8 0x3fd6cf31058670e8, 0x3cb0b244a70a8da9 // log(1/frcpa(1+109/256))=  +3.56396e-001
728
 
data8 0x3fd6f180e852f0b8, 0x3c9db7aefa866720 // log(1/frcpa(1+110/256))=  +3.58490e-001
729
 
data8 0x3fd71f5d71b894e8, 0x3cbe91c4bf324957 // log(1/frcpa(1+111/256))=  +3.61289e-001
730
 
data8 0x3fd74d5aefd66d58, 0x3cb06b3d9bfac023 // log(1/frcpa(1+112/256))=  +3.64096e-001
731
 
data8 0x3fd77b79922bd378, 0x3cb727d8804491f4 // log(1/frcpa(1+113/256))=  +3.66911e-001
732
 
data8 0x3fd7a9b9889f19e0, 0x3ca2ef22df5bc543 // log(1/frcpa(1+114/256))=  +3.69734e-001
733
 
data8 0x3fd7d81b037eb6a0, 0x3cb8fd3ba07a7ece // log(1/frcpa(1+115/256))=  +3.72565e-001
734
 
data8 0x3fd8069e33827230, 0x3c8bd1e25866e61a // log(1/frcpa(1+116/256))=  +3.75404e-001
735
 
data8 0x3fd82996d3ef8bc8, 0x3ca5aab9f5928928 // log(1/frcpa(1+117/256))=  +3.77538e-001
736
 
data8 0x3fd85855776dcbf8, 0x3ca56f33337789d6 // log(1/frcpa(1+118/256))=  +3.80391e-001
737
 
data8 0x3fd8873658327cc8, 0x3cbb8ef0401db49d // log(1/frcpa(1+119/256))=  +3.83253e-001
738
 
data8 0x3fd8aa75973ab8c8, 0x3cbb9961f509a680 // log(1/frcpa(1+120/256))=  +3.85404e-001
739
 
data8 0x3fd8d992dc8824e0, 0x3cb220512a53732d // log(1/frcpa(1+121/256))=  +3.88280e-001
740
 
data8 0x3fd908d2ea7d9510, 0x3c985f0e513bfb5c // log(1/frcpa(1+122/256))=  +3.91164e-001
741
 
data8 0x3fd92c59e79c0e50, 0x3cb82e073fd30d63 // log(1/frcpa(1+123/256))=  +3.93332e-001
742
 
data8 0x3fd95bd750ee3ed0, 0x3ca4aa7cdb6dd8a8 // log(1/frcpa(1+124/256))=  +3.96231e-001
743
 
data8 0x3fd98b7811a3ee58, 0x3caa93a5b660893e // log(1/frcpa(1+125/256))=  +3.99138e-001
744
 
data8 0x3fd9af47f33d4068, 0x3cac294b3b3190ba // log(1/frcpa(1+126/256))=  +4.01323e-001
745
 
data8 0x3fd9df270c1914a0, 0x3cbe1a58fd0cd67e // log(1/frcpa(1+127/256))=  +4.04245e-001
746
 
data8 0x3fda0325ed14fda0, 0x3cb1efa7950fb57e // log(1/frcpa(1+128/256))=  +4.06442e-001
747
 
data8 0x3fda33440224fa78, 0x3c8915fe75e7d477 // log(1/frcpa(1+129/256))=  +4.09379e-001
748
 
data8 0x3fda57725e80c380, 0x3ca72bd1062b1b7f // log(1/frcpa(1+130/256))=  +4.11587e-001
749
 
data8 0x3fda87d0165dd198, 0x3c91f7845f58dbad // log(1/frcpa(1+131/256))=  +4.14539e-001
750
 
data8 0x3fdaac2e6c03f890, 0x3cb6f237a911c509 // log(1/frcpa(1+132/256))=  +4.16759e-001
751
 
data8 0x3fdadccc6fdf6a80, 0x3c90ddc4b7687169 // log(1/frcpa(1+133/256))=  +4.19726e-001
752
 
data8 0x3fdb015b3eb1e790, 0x3c692dd7d90e1e8e // log(1/frcpa(1+134/256))=  +4.21958e-001
753
 
data8 0x3fdb323a3a635948, 0x3c6f85655cbe14de // log(1/frcpa(1+135/256))=  +4.24941e-001
754
 
data8 0x3fdb56fa04462908, 0x3c95252d841994de // log(1/frcpa(1+136/256))=  +4.27184e-001
755
 
data8 0x3fdb881aa659bc90, 0x3caa53a745a3642f // log(1/frcpa(1+137/256))=  +4.30182e-001
756
 
data8 0x3fdbad0bef3db160, 0x3cb32f2540dcc16a // log(1/frcpa(1+138/256))=  +4.32437e-001
757
 
data8 0x3fdbd21297781c28, 0x3cbd8e891e106f1d // log(1/frcpa(1+139/256))=  +4.34697e-001
758
 
data8 0x3fdc039236f08818, 0x3c809435af522ba7 // log(1/frcpa(1+140/256))=  +4.37718e-001
759
 
data8 0x3fdc28cb1e4d32f8, 0x3cb3944752fbd81e // log(1/frcpa(1+141/256))=  +4.39990e-001
760
 
data8 0x3fdc4e19b84723c0, 0x3c9a465260cd3fe5 // log(1/frcpa(1+142/256))=  +4.42267e-001
761
 
data8 0x3fdc7ff9c74554c8, 0x3c92447d5b6ca369 // log(1/frcpa(1+143/256))=  +4.45311e-001
762
 
data8 0x3fdca57b64e9db00, 0x3cb44344a8a00c82 // log(1/frcpa(1+144/256))=  +4.47600e-001
763
 
data8 0x3fdccb130a5ceba8, 0x3cbefaddfb97b73f // log(1/frcpa(1+145/256))=  +4.49895e-001
764
 
data8 0x3fdcf0c0d18f3268, 0x3cbd3e7bfee57898 // log(1/frcpa(1+146/256))=  +4.52194e-001
765
 
data8 0x3fdd232075b5a200, 0x3c9222599987447c // log(1/frcpa(1+147/256))=  +4.55269e-001
766
 
data8 0x3fdd490246defa68, 0x3cabafe9a767a80d // log(1/frcpa(1+148/256))=  +4.57581e-001
767
 
data8 0x3fdd6efa918d25c8, 0x3cb58a2624e1c6fd // log(1/frcpa(1+149/256))=  +4.59899e-001
768
 
data8 0x3fdd9509707ae528, 0x3cbdc3babce578e7 // log(1/frcpa(1+150/256))=  +4.62221e-001
769
 
data8 0x3fddbb2efe92c550, 0x3cb0ac0943c434a4 // log(1/frcpa(1+151/256))=  +4.64550e-001
770
 
data8 0x3fddee2f3445e4a8, 0x3cbba9d07ce820e8 // log(1/frcpa(1+152/256))=  +4.67663e-001
771
 
data8 0x3fde148a1a2726c8, 0x3cb6537e3375b205 // log(1/frcpa(1+153/256))=  +4.70004e-001
772
 
data8 0x3fde3afc0a49ff38, 0x3cbfed5518dbc20e // log(1/frcpa(1+154/256))=  +4.72350e-001
773
 
data8 0x3fde6185206d5168, 0x3cb6572601f73d5c // log(1/frcpa(1+155/256))=  +4.74702e-001
774
 
data8 0x3fde882578823d50, 0x3c9b24abd4584d1a // log(1/frcpa(1+156/256))=  +4.77060e-001
775
 
data8 0x3fdeaedd2eac9908, 0x3cb0ceb5e4d2c8f7 // log(1/frcpa(1+157/256))=  +4.79423e-001
776
 
data8 0x3fded5ac5f436be0, 0x3ca72f21f1f5238e // log(1/frcpa(1+158/256))=  +4.81792e-001
777
 
data8 0x3fdefc9326d16ab8, 0x3c85081a1639a45c // log(1/frcpa(1+159/256))=  +4.84166e-001
778
 
data8 0x3fdf2391a21575f8, 0x3cbf11015bdd297a // log(1/frcpa(1+160/256))=  +4.86546e-001
779
 
data8 0x3fdf4aa7ee031928, 0x3cb3795bc052a2d1 // log(1/frcpa(1+161/256))=  +4.88932e-001
780
 
data8 0x3fdf71d627c30bb0, 0x3c35c61f0f5a88f3 // log(1/frcpa(1+162/256))=  +4.91323e-001
781
 
data8 0x3fdf991c6cb3b378, 0x3c97d99419be6028 // log(1/frcpa(1+163/256))=  +4.93720e-001
782
 
data8 0x3fdfc07ada69a908, 0x3cbfe9341ded70b1 // log(1/frcpa(1+164/256))=  +4.96123e-001
783
 
data8 0x3fdfe7f18eb03d38, 0x3cb85718a640c33f // log(1/frcpa(1+165/256))=  +4.98532e-001
784
 
data8 0x3fe007c053c5002c, 0x3cb3addc9c065f09 // log(1/frcpa(1+166/256))=  +5.00946e-001
785
 
data8 0x3fe01b942198a5a0, 0x3c9d5aa4c77da6ac // log(1/frcpa(1+167/256))=  +5.03367e-001
786
 
data8 0x3fe02f74400c64e8, 0x3cb5a0ee4450ef52 // log(1/frcpa(1+168/256))=  +5.05793e-001
787
 
data8 0x3fe04360be7603ac, 0x3c9dd00c35630fe0 // log(1/frcpa(1+169/256))=  +5.08225e-001
788
 
data8 0x3fe05759ac47fe30, 0x3cbd063e1f0bd82c // log(1/frcpa(1+170/256))=  +5.10663e-001
789
 
data8 0x3fe06b5f1911cf50, 0x3cae8da674af5289 // log(1/frcpa(1+171/256))=  +5.13107e-001
790
 
data8 0x3fe078bf0533c568, 0x3c62241edf5fd1f7 // log(1/frcpa(1+172/256))=  +5.14740e-001
791
 
data8 0x3fe08cd9687e7b0c, 0x3cb3007febcca227 // log(1/frcpa(1+173/256))=  +5.17194e-001
792
 
data8 0x3fe0a10074cf9018, 0x3ca496e84603816b // log(1/frcpa(1+174/256))=  +5.19654e-001
793
 
data8 0x3fe0b5343a234474, 0x3cb46098d14fc90a // log(1/frcpa(1+175/256))=  +5.22120e-001
794
 
data8 0x3fe0c974c89431cc, 0x3cac0a7cdcbb86c6 // log(1/frcpa(1+176/256))=  +5.24592e-001
795
 
data8 0x3fe0ddc2305b9884, 0x3cb2f753210410ff // log(1/frcpa(1+177/256))=  +5.27070e-001
796
 
data8 0x3fe0eb524bafc918, 0x3c88affd6682229e // log(1/frcpa(1+178/256))=  +5.28726e-001
797
 
data8 0x3fe0ffb54213a474, 0x3cadeefbab9af993 // log(1/frcpa(1+179/256))=  +5.31214e-001
798
 
data8 0x3fe114253da97d9c, 0x3cbaf1c2b8bc160a // log(1/frcpa(1+180/256))=  +5.33709e-001
799
 
data8 0x3fe128a24f1d9afc, 0x3cb9cf4df375e650 // log(1/frcpa(1+181/256))=  +5.36210e-001
800
 
data8 0x3fe1365252bf0864, 0x3c985a621d4be111 // log(1/frcpa(1+182/256))=  +5.37881e-001
801
 
data8 0x3fe14ae558b4a92c, 0x3ca104c4aa8977d1 // log(1/frcpa(1+183/256))=  +5.40393e-001
802
 
data8 0x3fe15f85a19c7658, 0x3cbadf26e540f375 // log(1/frcpa(1+184/256))=  +5.42910e-001
803
 
data8 0x3fe16d4d38c119f8, 0x3cb3aea11caec416 // log(1/frcpa(1+185/256))=  +5.44592e-001
804
 
data8 0x3fe18203c20dd130, 0x3cba82d1211d1d6d // log(1/frcpa(1+186/256))=  +5.47121e-001
805
 
data8 0x3fe196c7bc4b1f38, 0x3cb6267acc4f4f4a // log(1/frcpa(1+187/256))=  +5.49656e-001
806
 
data8 0x3fe1a4a738b7a33c, 0x3c858930213c987d // log(1/frcpa(1+188/256))=  +5.51349e-001
807
 
data8 0x3fe1b981c0c9653c, 0x3c9bc2a4a30f697b // log(1/frcpa(1+189/256))=  +5.53895e-001
808
 
data8 0x3fe1ce69e8bb1068, 0x3cb7ae6199cf2a00 // log(1/frcpa(1+190/256))=  +5.56447e-001
809
 
data8 0x3fe1dc619de06944, 0x3c6b50bb38388177 // log(1/frcpa(1+191/256))=  +5.58152e-001
810
 
data8 0x3fe1f160a2ad0da0, 0x3cbd05b2778a5e1d // log(1/frcpa(1+192/256))=  +5.60715e-001
811
 
data8 0x3fe2066d7740737c, 0x3cb32e828f9c6bd6 // log(1/frcpa(1+193/256))=  +5.63285e-001
812
 
data8 0x3fe2147dba47a390, 0x3cbd579851b8b672 // log(1/frcpa(1+194/256))=  +5.65001e-001
813
 
data8 0x3fe229a1bc5ebac0, 0x3cbb321be5237ce8 // log(1/frcpa(1+195/256))=  +5.67582e-001
814
 
data8 0x3fe237c1841a502c, 0x3cb3b56e0915ea64 // log(1/frcpa(1+196/256))=  +5.69306e-001
815
 
data8 0x3fe24cfce6f80d98, 0x3cb34a4d1a422919 // log(1/frcpa(1+197/256))=  +5.71898e-001
816
 
data8 0x3fe25b2c55cd5760, 0x3cb237401ea5015e // log(1/frcpa(1+198/256))=  +5.73630e-001
817
 
data8 0x3fe2707f4d5f7c40, 0x3c9d30f20acc8341 // log(1/frcpa(1+199/256))=  +5.76233e-001
818
 
data8 0x3fe285e0842ca380, 0x3cbc4d866d5f21c0 // log(1/frcpa(1+200/256))=  +5.78842e-001
819
 
data8 0x3fe294294708b770, 0x3cb85e14d5dc54fa // log(1/frcpa(1+201/256))=  +5.80586e-001
820
 
data8 0x3fe2a9a2670aff0c, 0x3c7e6f8f468bbf91 // log(1/frcpa(1+202/256))=  +5.83207e-001
821
 
data8 0x3fe2b7fb2c8d1cc0, 0x3c930ffcf63c8b65 // log(1/frcpa(1+203/256))=  +5.84959e-001
822
 
data8 0x3fe2c65a6395f5f4, 0x3ca0afe20b53d2d2 // log(1/frcpa(1+204/256))=  +5.86713e-001
823
 
data8 0x3fe2dbf557b0df40, 0x3cb646be1188fbc9 // log(1/frcpa(1+205/256))=  +5.89350e-001
824
 
data8 0x3fe2ea64c3f97654, 0x3c96516fa8df33b2 // log(1/frcpa(1+206/256))=  +5.91113e-001
825
 
data8 0x3fe3001823684d70, 0x3cb96d64e16d1360 // log(1/frcpa(1+207/256))=  +5.93762e-001
826
 
data8 0x3fe30e97e9a8b5cc, 0x3c98ef96bc97cca0 // log(1/frcpa(1+208/256))=  +5.95531e-001
827
 
data8 0x3fe32463ebdd34e8, 0x3caef1dc9a56c1bf // log(1/frcpa(1+209/256))=  +5.98192e-001
828
 
data8 0x3fe332f4314ad794, 0x3caa4f0ac5d5fa11 // log(1/frcpa(1+210/256))=  +5.99970e-001
829
 
data8 0x3fe348d90e7464cc, 0x3cbe7889f0516acd // log(1/frcpa(1+211/256))=  +6.02643e-001
830
 
data8 0x3fe35779f8c43d6c, 0x3ca96bbab7245411 // log(1/frcpa(1+212/256))=  +6.04428e-001
831
 
data8 0x3fe36621961a6a98, 0x3ca31f32262db9fb // log(1/frcpa(1+213/256))=  +6.06217e-001
832
 
data8 0x3fe37c299f3c3668, 0x3cb15c72c107ee29 // log(1/frcpa(1+214/256))=  +6.08907e-001
833
 
data8 0x3fe38ae2171976e4, 0x3cba42a2554b2dd4 // log(1/frcpa(1+215/256))=  +6.10704e-001
834
 
data8 0x3fe399a157a603e4, 0x3cb99c62286d8919 // log(1/frcpa(1+216/256))=  +6.12504e-001
835
 
data8 0x3fe3afccfe77b9d0, 0x3ca11048f96a43bd // log(1/frcpa(1+217/256))=  +6.15210e-001
836
 
data8 0x3fe3be9d503533b4, 0x3ca4022f47588c3e // log(1/frcpa(1+218/256))=  +6.17018e-001
837
 
data8 0x3fe3cd7480b4a8a0, 0x3cb4ba7afc2dc56a // log(1/frcpa(1+219/256))=  +6.18830e-001
838
 
data8 0x3fe3e3c43918f76c, 0x3c859673d064b8ba // log(1/frcpa(1+220/256))=  +6.21554e-001
839
 
data8 0x3fe3f2acb27ed6c4, 0x3cb55c6b452a16a8 // log(1/frcpa(1+221/256))=  +6.23373e-001
840
 
data8 0x3fe4019c2125ca90, 0x3cb8c367879c5a31 // log(1/frcpa(1+222/256))=  +6.25197e-001
841
 
data8 0x3fe4181061389720, 0x3cb2c17a79c5cc6c // log(1/frcpa(1+223/256))=  +6.27937e-001
842
 
data8 0x3fe42711518df544, 0x3ca5f38d47012fc5 // log(1/frcpa(1+224/256))=  +6.29769e-001
843
 
data8 0x3fe436194e12b6bc, 0x3cb9854d65a9b426 // log(1/frcpa(1+225/256))=  +6.31604e-001
844
 
data8 0x3fe445285d68ea68, 0x3ca3ff9b3a81cd81 // log(1/frcpa(1+226/256))=  +6.33442e-001
845
 
data8 0x3fe45bcc464c8938, 0x3cb0a2d8011a6c05 // log(1/frcpa(1+227/256))=  +6.36206e-001
846
 
data8 0x3fe46aed21f117fc, 0x3c8a2be41f8e9f3d // log(1/frcpa(1+228/256))=  +6.38053e-001
847
 
data8 0x3fe47a1527e8a2d0, 0x3cba4a83594fab09 // log(1/frcpa(1+229/256))=  +6.39903e-001
848
 
data8 0x3fe489445efffcc8, 0x3cbf306a23dcbcde // log(1/frcpa(1+230/256))=  +6.41756e-001
849
 
data8 0x3fe4a018bcb69834, 0x3ca46c9285029fd1 // log(1/frcpa(1+231/256))=  +6.44543e-001
850
 
data8 0x3fe4af5a0c9d65d4, 0x3cbbc1db897580e3 // log(1/frcpa(1+232/256))=  +6.46405e-001
851
 
data8 0x3fe4bea2a5bdbe84, 0x3cb84d880d7ef775 // log(1/frcpa(1+233/256))=  +6.48271e-001
852
 
data8 0x3fe4cdf28f10ac44, 0x3cb3ec4b7893ce1f // log(1/frcpa(1+234/256))=  +6.50140e-001
853
 
data8 0x3fe4dd49cf994058, 0x3c897224d59d3408 // log(1/frcpa(1+235/256))=  +6.52013e-001
854
 
data8 0x3fe4eca86e64a680, 0x3cbccf620f24f0cd // log(1/frcpa(1+236/256))=  +6.53889e-001
855
 
data8 0x3fe503c43cd8eb68, 0x3c3f872c65971084 // log(1/frcpa(1+237/256))=  +6.56710e-001
856
 
data8 0x3fe513356667fc54, 0x3cb9ca64cc3d52c8 // log(1/frcpa(1+238/256))=  +6.58595e-001
857
 
data8 0x3fe522ae0738a3d4, 0x3cbe708164c75968 // log(1/frcpa(1+239/256))=  +6.60483e-001
858
 
data8 0x3fe5322e26867854, 0x3cb9988ba4aea615 // log(1/frcpa(1+240/256))=  +6.62376e-001
859
 
data8 0x3fe541b5cb979808, 0x3ca1662e3a6b95f5 // log(1/frcpa(1+241/256))=  +6.64271e-001
860
 
data8 0x3fe55144fdbcbd60, 0x3cb3acd4ca45c1e0 // log(1/frcpa(1+242/256))=  +6.66171e-001
861
 
data8 0x3fe560dbc45153c4, 0x3cb4988947959fed // log(1/frcpa(1+243/256))=  +6.68074e-001
862
 
data8 0x3fe5707a26bb8c64, 0x3cb3017fe6607ba9 // log(1/frcpa(1+244/256))=  +6.69980e-001
863
 
data8 0x3fe587f60ed5b8fc, 0x3cbe7a3266366ed4 // log(1/frcpa(1+245/256))=  +6.72847e-001
864
 
data8 0x3fe597a7977c8f30, 0x3ca1e12b9959a90e // log(1/frcpa(1+246/256))=  +6.74763e-001
865
 
data8 0x3fe5a760d634bb88, 0x3cb7c365e53d9602 // log(1/frcpa(1+247/256))=  +6.76682e-001
866
 
data8 0x3fe5b721d295f10c, 0x3cb716c2551ccbf0 // log(1/frcpa(1+248/256))=  +6.78605e-001
867
 
data8 0x3fe5c6ea94431ef8, 0x3ca02b2ed0e28261 // log(1/frcpa(1+249/256))=  +6.80532e-001
868
 
data8 0x3fe5d6bb22ea86f4, 0x3caf43a8bbb2f974 // log(1/frcpa(1+250/256))=  +6.82462e-001
869
 
data8 0x3fe5e6938645d38c, 0x3cbcedc98821b333 // log(1/frcpa(1+251/256))=  +6.84397e-001
870
 
data8 0x3fe5f673c61a2ed0, 0x3caa385eef5f2789 // log(1/frcpa(1+252/256))=  +6.86335e-001
871
 
data8 0x3fe6065bea385924, 0x3cb11624f165c5b4 // log(1/frcpa(1+253/256))=  +6.88276e-001
872
 
data8 0x3fe6164bfa7cc068, 0x3cbad884f87073fa // log(1/frcpa(1+254/256))=  +6.90222e-001
873
 
data8 0x3fe62643fecf9740, 0x3cb78c51da12f4df // log(1/frcpa(1+255/256))=  +6.92171e-001
874
 
LOCAL_OBJECT_END(pow_Tt)
875
 
 
876
 
 
877
 
// Table 1 is 2^(index_1/128) where
878
 
// index_1 goes from 0 to 15
879
 
LOCAL_OBJECT_START(pow_tbl1)
880
 
data8 0x8000000000000000 , 0x00003FFF
881
 
data8 0x80B1ED4FD999AB6C , 0x00003FFF
882
 
data8 0x8164D1F3BC030773 , 0x00003FFF
883
 
data8 0x8218AF4373FC25EC , 0x00003FFF
884
 
data8 0x82CD8698AC2BA1D7 , 0x00003FFF
885
 
data8 0x8383594EEFB6EE37 , 0x00003FFF
886
 
data8 0x843A28C3ACDE4046 , 0x00003FFF
887
 
data8 0x84F1F656379C1A29 , 0x00003FFF
888
 
data8 0x85AAC367CC487B15 , 0x00003FFF
889
 
data8 0x8664915B923FBA04 , 0x00003FFF
890
 
data8 0x871F61969E8D1010 , 0x00003FFF
891
 
data8 0x87DB357FF698D792 , 0x00003FFF
892
 
data8 0x88980E8092DA8527 , 0x00003FFF
893
 
data8 0x8955EE03618E5FDD , 0x00003FFF
894
 
data8 0x8A14D575496EFD9A , 0x00003FFF
895
 
data8 0x8AD4C6452C728924 , 0x00003FFF
896
 
LOCAL_OBJECT_END(pow_tbl1)
897
 
 
898
 
 
899
 
// Table 2 is 2^(index_1/8) where
900
 
// index_2 goes from 0 to 7
901
 
LOCAL_OBJECT_START(pow_tbl2)
902
 
data8 0x8000000000000000 , 0x00003FFF
903
 
data8 0x8B95C1E3EA8BD6E7 , 0x00003FFF
904
 
data8 0x9837F0518DB8A96F , 0x00003FFF
905
 
data8 0xA5FED6A9B15138EA , 0x00003FFF
906
 
data8 0xB504F333F9DE6484 , 0x00003FFF
907
 
data8 0xC5672A115506DADD , 0x00003FFF
908
 
data8 0xD744FCCAD69D6AF4 , 0x00003FFF
909
 
data8 0xEAC0C6E7DD24392F , 0x00003FFF
910
 
LOCAL_OBJECT_END(pow_tbl2)
911
 
 
912
 
.section .text
913
 
GLOBAL_LIBM_ENTRY(pow)
914
 
 
915
 
// Get exponent of x.  Will be used to calculate K.
916
 
{ .mfi
917
 
          getf.exp     pow_GR_signexp_X = f8
918
 
          fms.s1 POW_Xm1 = f8,f1,f1     // Will be used for r1 if x>0
919
 
          mov           pow_GR_17ones   = 0x1FFFF
920
 
}
921
 
{ .mfi
922
 
          addl          pow_AD_P        = @ltoff(pow_table_P), gp
923
 
          fma.s1 POW_Xp1 = f8,f1,f1     // Will be used for r1 if x<0
924
 
          nop.i 999
925
 
;;
926
 
}
927
 
 
928
 
// Get significand of x.  Will be used to get index to fetch T, Tt.
929
 
{ .mfi
930
 
          getf.sig      pow_GR_sig_X    = f8
931
 
          frcpa.s1      POW_B, p6       = f1,f8
932
 
          nop.i 999
933
 
}
934
 
{ .mfi
935
 
          ld8 pow_AD_P = [pow_AD_P]
936
 
          fma.s1        POW_NORM_X      = f8,f1,f0
937
 
          mov          pow_GR_exp_2tom8 = 0xFFF7
938
 
}
939
 
;;
940
 
 
941
 
// p13 = TRUE ==> X is unorm
942
 
// DOUBLE 0x10033  exponent limit at which y is an integer
943
 
{ .mfi
944
 
          nop.m 999
945
 
          fclass.m  p13,p0              = f8, 0x0b  // Test for x unorm
946
 
          addl pow_GR_10033             = 0x10033, r0
947
 
}
948
 
{ .mfi
949
 
          mov           pow_GR_16ones   = 0xFFFF
950
 
          fma.s1        POW_NORM_Y      = f9,f1,f0
951
 
          nop.i 999
952
 
}
953
 
;;
954
 
 
955
 
// p14 = TRUE ==> X is ZERO
956
 
{ .mfi
957
 
          adds          pow_AD_Tt       = pow_Tt - pow_table_P,  pow_AD_P
958
 
          fclass.m  p14,p0              = f8, 0x07
959
 
          and           pow_GR_exp_X    = pow_GR_signexp_X, pow_GR_17ones
960
 
}
961
 
{ .mfi
962
 
          adds          pow_AD_Q        = pow_table_Q - pow_table_P,  pow_AD_P
963
 
          nop.f 999
964
 
          nop.i 999
965
 
}
966
 
;;
967
 
 
968
 
{ .mfi
969
 
          ldfe          POW_P5          = [pow_AD_P], 16
970
 
          fcmp.lt.s1 p8,p9 = f8, f0     // Test for x<0
971
 
          nop.i 999
972
 
}
973
 
{ .mib
974
 
          ldfe          POW_P4          = [pow_AD_Q], 16
975
 
          sub       pow_GR_true_exp_X   = pow_GR_exp_X, pow_GR_16ones
976
 
(p13)     br.cond.spnt POW_X_DENORM
977
 
}
978
 
;;
979
 
 
980
 
// Continue normal and denormal paths here
981
 
POW_COMMON:
982
 
// p11 = TRUE ==> Y is a NAN
983
 
{ .mfi
984
 
          ldfe          POW_P3          = [pow_AD_P], 16
985
 
          fclass.m  p11,p0              = f9, 0xc3
986
 
          nop.i 999
987
 
}
988
 
{ .mfi
989
 
          ldfe          POW_P2          = [pow_AD_Q], 16
990
 
          nop.f 999
991
 
          mov pow_GR_y_zero = 0
992
 
}
993
 
;;
994
 
 
995
 
// Note POW_Xm1 and POW_r1 are used interchangably
996
 
{ .mfi
997
 
          alloc         r32=ar.pfs,2,19,4,0
998
 
          fms.s1        POW_r           = POW_B, POW_NORM_X,f1
999
 
          nop.i 999
1000
 
}
1001
 
{ .mfi
1002
 
          setf.sig POW_int_K            = pow_GR_true_exp_X
1003
 
(p8)      fnma.s1        POW_Xm1        = POW_Xp1,f1,f0
1004
 
          nop.i 999
1005
 
}
1006
 
;;
1007
 
 
1008
 
// p12 = TRUE if Y is ZERO
1009
 
// Compute xsq to decide later if |x|=1
1010
 
{ .mfi
1011
 
          ldfe          POW_P1          = [pow_AD_P], 16
1012
 
          fclass.m      p12,p0          = f9, 0x07
1013
 
          shl           pow_GR_offset   = pow_GR_sig_X, 1
1014
 
}
1015
 
{ .mfb
1016
 
          ldfe          POW_P0          = [pow_AD_Q], 16
1017
 
          fma.s1        POW_xsq = POW_NORM_X, POW_NORM_X, f0
1018
 
(p11)     br.cond.spnt  POW_Y_NAN       // Branch if y=nan
1019
 
}
1020
 
;;
1021
 
 
1022
 
// Get exponent of |x|-1 to use in comparison to 2^-8
1023
 
{ .mfi
1024
 
          getf.exp  pow_GR_signexp_Xm1  = POW_Xm1
1025
 
          fcvt.fx.s1   POW_int_Y        = POW_NORM_Y
1026
 
          shr.u     pow_GR_offset       = pow_GR_offset,56
1027
 
}
1028
 
;;
1029
 
 
1030
 
// p11 = TRUE ==> X is a NAN
1031
 
{ .mfi
1032
 
          ldfpd         POW_log2_hi, POW_log2_lo  = [pow_AD_Q], 16
1033
 
          fclass.m      p11,p0          = f8, 0xc3
1034
 
          shladd pow_AD_Tt = pow_GR_offset, 4, pow_AD_Tt
1035
 
}
1036
 
{ .mfi
1037
 
          ldfe          POW_inv_log2_by_128 = [pow_AD_P], 16
1038
 
          fma.s1 POW_delta              = f0,f0,f0 // delta=0 in case |x| near 1
1039
 
(p12)     mov pow_GR_y_zero = 1
1040
 
}
1041
 
;;
1042
 
 
1043
 
{ .mfi
1044
 
          ldfpd  POW_Q2, POW_Q3         = [pow_AD_P], 16
1045
 
          fma.s1 POW_G                  = f0,f0,f0  // G=0 in case |x| near 1
1046
 
          and       pow_GR_exp_Xm1      = pow_GR_signexp_Xm1, pow_GR_17ones
1047
 
}
1048
 
;;
1049
 
 
1050
 
// Determine if we will use the |x| near 1 path (p6) or normal path (p7)
1051
 
{ .mfi
1052
 
          getf.exp  pow_GR_signexp_Y    = POW_NORM_Y
1053
 
          nop.f 999
1054
 
          cmp.lt p6,p7                  = pow_GR_exp_Xm1, pow_GR_exp_2tom8
1055
 
}
1056
 
{ .mfb
1057
 
          ldfpd  POW_T, POW_Tt          = [pow_AD_Tt], 16
1058
 
          fma.s1        POW_rsq         = POW_r, POW_r,f0
1059
 
(p11)     br.cond.spnt  POW_X_NAN       // Branch if x=nan and y not nan
1060
 
}
1061
 
;;
1062
 
 
1063
 
// If on the x near 1 path, assign r1 to r and r1*r1 to rsq
1064
 
{ .mfi
1065
 
          ldfpd  POW_Q0_half, POW_Q1    = [pow_AD_P], 16
1066
 
(p6)      fma.s1    POW_r               = POW_r1, f1, f0
1067
 
          nop.i 999
1068
 
}
1069
 
{ .mfb
1070
 
          nop.m 999
1071
 
(p6)      fma.s1    POW_rsq             = POW_r1, POW_r1, f0
1072
 
(p14)     br.cond.spnt POW_X_0          // Branch if x zero and y not nan
1073
 
}
1074
 
;;
1075
 
 
1076
 
{ .mfi
1077
 
          ldfpd   POW_Q4, POW_RSHF      = [pow_AD_P], 16
1078
 
(p7)      fma.s1 POW_v6                 = POW_r,  POW_P5, POW_P4
1079
 
          nop.i 999
1080
 
}
1081
 
{ .mfi
1082
 
          mov pow_GR_exp_2toM63         = 0xffc0  // Exponent of 2^-63
1083
 
(p6)      fma.s1 POW_v6                 = POW_r1, POW_P5, POW_P4
1084
 
          nop.i 999
1085
 
}
1086
 
;;
1087
 
 
1088
 
{ .mfi
1089
 
          setf.exp POW_2toM63 = pow_GR_exp_2toM63  // Form 2^-63 for test of q
1090
 
(p7)      fma.s1 POW_v4                 = POW_P3, POW_r,  POW_P2
1091
 
          nop.i 999
1092
 
}
1093
 
{ .mfi
1094
 
          nop.m 999
1095
 
(p6)      fma.s1 POW_v4                 = POW_P3, POW_r1, POW_P2
1096
 
          nop.i 999
1097
 
}
1098
 
;;
1099
 
 
1100
 
{ .mfi
1101
 
          nop.m 999
1102
 
          fcvt.xf POW_K                 = POW_int_K
1103
 
          nop.i 999
1104
 
}
1105
 
;;
1106
 
 
1107
 
{ .mfi
1108
 
          getf.sig pow_GR_sig_int_Y     = POW_int_Y
1109
 
          fnma.s1 POW_twoV              = POW_NORM_Y, POW_rsq,f0
1110
 
          and pow_GR_exp_Y              = pow_GR_signexp_Y, pow_GR_17ones
1111
 
}
1112
 
{ .mfb
1113
 
          andcm pow_GR_sign_Y           = pow_GR_signexp_Y, pow_GR_17ones
1114
 
          fma.s1 POW_U                  = POW_NORM_Y,POW_r,f0
1115
 
(p12)     br.cond.spnt POW_Y_0   // Branch if y=zero, x not zero or nan
1116
 
}
1117
 
;;
1118
 
 
1119
 
// p11 = TRUE ==> X is NEGATIVE but not inf
1120
 
{ .mfi
1121
 
          ldfe      POW_log2_by_128_lo  = [pow_AD_P], 16
1122
 
          fclass.m  p11,p0              = f8, 0x1a
1123
 
          nop.i 999
1124
 
}
1125
 
{ .mfi
1126
 
          ldfe      POW_log2_by_128_hi  = [pow_AD_Q], 16
1127
 
          fma.s1 POW_v2                 = POW_P1, POW_r,  POW_P0
1128
 
          nop.i 999
1129
 
}
1130
 
;;
1131
 
 
1132
 
{ .mfi
1133
 
          nop.m 999
1134
 
          fcvt.xf   POW_float_int_Y     = POW_int_Y
1135
 
          nop.i 999
1136
 
}
1137
 
{ .mfi
1138
 
          nop.m 999
1139
 
          fma.s1 POW_v3                 = POW_v6, POW_rsq,  POW_v4
1140
 
          adds          pow_AD_tbl1     = pow_tbl1 - pow_Tt,  pow_AD_Q
1141
 
}
1142
 
;;
1143
 
 
1144
 
{ .mfi
1145
 
          nop.m 999
1146
 
(p7)      fma.s1 POW_delta              = POW_K, POW_log2_lo, POW_Tt
1147
 
          nop.i 999
1148
 
}
1149
 
{ .mfi
1150
 
          nop.m 999
1151
 
(p7)      fma.s1 POW_G                  = POW_K, POW_log2_hi, POW_T
1152
 
          adds pow_AD_tbl2              = pow_tbl2 - pow_tbl1,  pow_AD_tbl1
1153
 
}
1154
 
;;
1155
 
 
1156
 
{ .mfi
1157
 
          nop.m 999
1158
 
          fms.s1 POW_e2                 = POW_NORM_Y, POW_r, POW_U
1159
 
          nop.i 999
1160
 
}
1161
 
{ .mfi
1162
 
          nop.m 999
1163
 
          fma.s1 POW_Z2                 = POW_twoV, POW_Q0_half, POW_U
1164
 
          nop.i 999
1165
 
}
1166
 
;;
1167
 
 
1168
 
{ .mfi
1169
 
          nop.m 999
1170
 
          fma.s1 POW_Yrcub              = POW_rsq, POW_U, f0
1171
 
          nop.i 999
1172
 
}
1173
 
{ .mfi
1174
 
          nop.m 999
1175
 
          fma.s1 POW_p                  = POW_rsq, POW_v3, POW_v2
1176
 
          nop.i 999
1177
 
}
1178
 
;;
1179
 
 
1180
 
// p11 = TRUE ==> X is NEGATIVE but not inf
1181
 
//    p12 = TRUE ==> X is NEGATIVE  AND  Y  already even int
1182
 
//    p13 = TRUE ==> X is NEGATIVE  AND  Y possible int
1183
 
{ .mfi
1184
 
          nop.m 999
1185
 
          fma.s1 POW_Z1                 = POW_NORM_Y, POW_G, f0
1186
 
(p11)     cmp.gt.unc  p12,p13           = pow_GR_exp_Y, pow_GR_10033
1187
 
}
1188
 
{ .mfi
1189
 
          nop.m 999
1190
 
          fma.s1 POW_Gpr                = POW_G, f1, POW_r
1191
 
          nop.i 999
1192
 
}
1193
 
;;
1194
 
 
1195
 
// By adding RSHF (1.1000...*2^63) we put integer part in rightmost significand
1196
 
{ .mfi
1197
 
          nop.m 999
1198
 
          fma.s1 POW_W2  = POW_Z2, POW_inv_log2_by_128, POW_RSHF
1199
 
          nop.i 999
1200
 
}
1201
 
{ .mfi
1202
 
          nop.m 999
1203
 
          fms.s1 POW_UmZ2               = POW_U, f1, POW_Z2
1204
 
          nop.i 999
1205
 
}
1206
 
;;
1207
 
 
1208
 
{ .mfi
1209
 
          nop.m 999
1210
 
          fma.s1 POW_e3                 = POW_NORM_Y, POW_delta, f0
1211
 
          nop.i 999
1212
 
}
1213
 
;;
1214
 
 
1215
 
{ .mfi
1216
 
          nop.m 999
1217
 
          fma.s1 POW_Z3                 = POW_p, POW_Yrcub, f0
1218
 
          nop.i 999
1219
 
}
1220
 
{ .mfi
1221
 
          nop.m 999
1222
 
          fma.s1 POW_GY_Z2              = POW_G, POW_NORM_Y, POW_Z2
1223
 
          nop.i 999
1224
 
}
1225
 
;;
1226
 
 
1227
 
// By adding RSHF (1.1000...*2^63) we put integer part in rightmost significand
1228
 
{ .mfi
1229
 
          nop.m 999
1230
 
          fms.s1 POW_e1                 = POW_NORM_Y, POW_G, POW_Z1
1231
 
          nop.i 999
1232
 
}
1233
 
{ .mfi
1234
 
          nop.m 999
1235
 
          fma.s1 POW_W1  = POW_Z1, POW_inv_log2_by_128, POW_RSHF
1236
 
          nop.i 999
1237
 
}
1238
 
;;
1239
 
 
1240
 
// p13 = TRUE ==> X is NEGATIVE  AND  Y possible int
1241
 
//     p10 = TRUE ==> X is NEG and Y is an int
1242
 
//     p12 = TRUE ==> X is NEG and Y is not an int
1243
 
{ .mfi
1244
 
          nop.m 999
1245
 
(p13)     fcmp.eq.unc.s1 p10,p12        = POW_float_int_Y,  POW_NORM_Y
1246
 
          mov pow_GR_xneg_yodd = 0
1247
 
}
1248
 
{ .mfi
1249
 
          nop.m 999
1250
 
          fma.s1 POW_Y_Gpr              = POW_NORM_Y, POW_Gpr, f0
1251
 
          nop.i 999
1252
 
}
1253
 
;;
1254
 
 
1255
 
// By subtracting RSHF we get rounded integer POW_N2float
1256
 
{ .mfi
1257
 
          nop.m 999
1258
 
          fms.s1 POW_N2float  = POW_W2, f1, POW_RSHF
1259
 
          nop.i 999
1260
 
}
1261
 
{ .mfi
1262
 
          nop.m 999
1263
 
          fma.s1 POW_UmZ2pV             = POW_twoV,POW_Q0_half,POW_UmZ2
1264
 
          nop.i 999
1265
 
}
1266
 
;;
1267
 
 
1268
 
{ .mfi
1269
 
          nop.m 999
1270
 
          fma.s1 POW_Z3sq               = POW_Z3, POW_Z3, f0
1271
 
          nop.i 999
1272
 
}
1273
 
{ .mfi
1274
 
          nop.m 999
1275
 
          fma.s1 POW_v4                 = POW_Z3, POW_Q3, POW_Q2
1276
 
          nop.i 999
1277
 
}
1278
 
;;
1279
 
 
1280
 
// Extract rounded integer from rightmost significand of POW_W2
1281
 
// By subtracting RSHF we get rounded integer POW_N1float
1282
 
{ .mfi
1283
 
          getf.sig pow_GR_int_W2        = POW_W2
1284
 
          fms.s1 POW_N1float  = POW_W1, f1, POW_RSHF
1285
 
          nop.i 999
1286
 
}
1287
 
{ .mfi
1288
 
          nop.m 999
1289
 
          fma.s1 POW_v2                 = POW_Z3, POW_Q1, POW_Q0_half
1290
 
          nop.i 999
1291
 
}
1292
 
;;
1293
 
 
1294
 
{ .mfi
1295
 
          nop.m 999
1296
 
          fnma.s1 POW_s2 = POW_N2float, POW_log2_by_128_hi, POW_Z2
1297
 
          nop.i 999
1298
 
}
1299
 
{ .mfi
1300
 
          nop.m 999
1301
 
          fma.s1 POW_e2                 = POW_e2,f1,POW_UmZ2pV
1302
 
          nop.i 999
1303
 
}
1304
 
;;
1305
 
 
1306
 
// Extract rounded integer from rightmost significand of POW_W1
1307
 
// Test if x inf
1308
 
{ .mfi
1309
 
          getf.sig pow_GR_int_W1        = POW_W1
1310
 
          fclass.m p15,p0 = POW_NORM_X,  0x23
1311
 
          nop.i 999
1312
 
}
1313
 
{ .mfb
1314
 
          nop.m 999
1315
 
          fnma.s1 POW_f2  = POW_N2float, POW_log2_by_128_lo, f1
1316
 
(p12)     br.cond.spnt POW_X_NEG_Y_NONINT  // Branch if x neg, y not integer
1317
 
}
1318
 
;;
1319
 
 
1320
 
// p11 = TRUE ==> X is +1.0
1321
 
// p12 = TRUE ==> X is NEGATIVE  AND Y is an odd integer
1322
 
{ .mfi
1323
 
          getf.exp pow_GR_signexp_Y_Gpr = POW_Y_Gpr
1324
 
          fcmp.eq.s1 p11,p0 = POW_NORM_X, f1
1325
 
(p10)     tbit.nz.unc  p12,p0           = pow_GR_sig_int_Y,0
1326
 
}
1327
 
{ .mfi
1328
 
          nop.m 999
1329
 
          fma.s1 POW_v3                 = POW_Z3sq, POW_Q4, POW_v4
1330
 
          nop.i 999
1331
 
}
1332
 
;;
1333
 
 
1334
 
{ .mfi
1335
 
          nop.m 999
1336
 
          fnma.s1 POW_f1  = POW_N1float, POW_log2_by_128_lo, f1
1337
 
          nop.i 999
1338
 
}
1339
 
{ .mfb
1340
 
          nop.m 999
1341
 
          fnma.s1 POW_s1  = POW_N1float, POW_log2_by_128_hi, POW_Z1
1342
 
(p15)     br.cond.spnt POW_X_INF
1343
 
}
1344
 
;;
1345
 
 
1346
 
// Test x and y and flag denormal
1347
 
{ .mfi
1348
 
          nop.m 999
1349
 
          fcmp.eq.s0 p15,p0 = f8,f9
1350
 
          nop.i 999
1351
 
}
1352
 
{ .mfi
1353
 
          nop.m 999
1354
 
          fma.s1 POW_pYrcub_e3          = POW_p, POW_Yrcub, POW_e3
1355
 
          nop.i 999
1356
 
}
1357
 
;;
1358
 
 
1359
 
{ .mfi
1360
 
          nop.m 999
1361
 
          fcmp.eq.s1 p7,p0 = POW_NORM_Y, f1  // Test for y=1.0
1362
 
          nop.i 999
1363
 
}
1364
 
{ .mfi
1365
 
          nop.m 999
1366
 
          fma.s1  POW_e12               = POW_e1,f1,POW_e2
1367
 
          nop.i 999
1368
 
}
1369
 
;;
1370
 
 
1371
 
{ .mfi
1372
 
          add pow_GR_int_N              = pow_GR_int_W1, pow_GR_int_W2
1373
 
(p11)     fma.d.s0 f8 = f1,f1,f0    // If x=1, result is +1
1374
 
          nop.i 999
1375
 
}
1376
 
{ .mib
1377
 
(p12)     mov pow_GR_xneg_yodd = 1
1378
 
          nop.i 999
1379
 
(p11)     br.ret.spnt b0            // Early exit if x=1.0, result is +1
1380
 
}
1381
 
;;
1382
 
 
1383
 
{ .mfi
1384
 
          and pow_GR_index1             = 0x0f, pow_GR_int_N
1385
 
          fma.s1 POW_q                  = POW_Z3sq, POW_v3, POW_v2
1386
 
          shr pow_int_GR_M              = pow_GR_int_N, 7    // M = N/128
1387
 
}
1388
 
{ .mib
1389
 
          and pow_GR_index2             = 0x70, pow_GR_int_N
1390
 
          cmp.eq p6, p0                 = pow_GR_xneg_yodd, r0
1391
 
(p7)      br.ret.spnt b0        // Early exit if y=1.0, result is x
1392
 
}
1393
 
;;
1394
 
 
1395
 
{ .mfi
1396
 
          shladd pow_AD_T1              = pow_GR_index1, 4, pow_AD_tbl1
1397
 
          fma.s1 POW_s                  = POW_s1, f1, POW_s2
1398
 
          add pow_int_GR_M              = pow_GR_16ones, pow_int_GR_M
1399
 
}
1400
 
{ .mfi
1401
 
          add pow_AD_T2                 = pow_AD_tbl2, pow_GR_index2
1402
 
          fma.s1 POW_f12                = POW_f1, POW_f2,f0
1403
 
          and pow_GR_exp_Y_Gpr          = pow_GR_signexp_Y_Gpr, pow_GR_17ones
1404
 
}
1405
 
;;
1406
 
 
1407
 
{ .mmi
1408
 
          ldfe POW_T1                   = [pow_AD_T1]
1409
 
          ldfe POW_T2                   = [pow_AD_T2]
1410
 
          sub pow_GR_true_exp_Y_Gpr     = pow_GR_exp_Y_Gpr, pow_GR_16ones
1411
 
}
1412
 
;;
1413
 
 
1414
 
{ .mfi
1415
 
          setf.exp POW_2M               = pow_int_GR_M
1416
 
          fma.s1 POW_e123               = POW_e12, f1, POW_e3
1417
 
          nop.i 999
1418
 
}
1419
 
{ .mfb
1420
 
(p6)      cmp.gt p6, p0                 = -11, pow_GR_true_exp_Y_Gpr
1421
 
          fma.s1 POW_d                  = POW_GY_Z2, f1, POW_pYrcub_e3
1422
 
(p6)      br.cond.spnt POW_NEAR_ONE // branch if |y*log(x)| < 2^(-11)
1423
 
}
1424
 
;;
1425
 
 
1426
 
{ .mfi
1427
 
          nop.m 999
1428
 
          fma.s1 POW_q                  = POW_Z3sq, POW_q, POW_Z3
1429
 
          nop.i 999
1430
 
}
1431
 
;;
1432
 
 
1433
 
// p8 TRUE ==> |Y(G + r)| >= 10
1434
 
 
1435
 
// double
1436
 
//     -2^10  -2^9             2^9   2^10
1437
 
// -----+-----+----+ ... +-----+-----+-----
1438
 
//  p8  |             p9             |  p8
1439
 
//      |     |       p10      |     |
1440
 
 
1441
 
// Form signexp of constants to indicate overflow
1442
 
{ .mfi
1443
 
          mov         pow_GR_big_pos    = 0x103ff
1444
 
          fma.s1 POW_ssq                = POW_s, POW_s, f0
1445
 
          cmp.le p8,p9                  = 10, pow_GR_true_exp_Y_Gpr
1446
 
}
1447
 
{ .mfi
1448
 
          mov         pow_GR_big_neg    = 0x303ff
1449
 
          fma.s1 POW_v4                 = POW_s, POW_Q3, POW_Q2
1450
 
          andcm pow_GR_sign_Y_Gpr       = pow_GR_signexp_Y_Gpr, pow_GR_17ones
1451
 
}
1452
 
;;
1453
 
 
1454
 
// Form big positive and negative constants to test for possible overflow
1455
 
{ .mfi
1456
 
          setf.exp POW_big_pos          = pow_GR_big_pos
1457
 
          fma.s1 POW_v2                 = POW_s, POW_Q1, POW_Q0_half
1458
 
(p9)      cmp.le.unc p0,p10             = 9, pow_GR_true_exp_Y_Gpr
1459
 
}
1460
 
{ .mfb
1461
 
          setf.exp POW_big_neg          = pow_GR_big_neg
1462
 
          fma.s1 POW_1ps                = f1,f1,POW_s
1463
 
(p8)      br.cond.spnt POW_OVER_UNDER_X_NOT_INF
1464
 
}
1465
 
;;
1466
 
 
1467
 
// f123 = f12*(e123+1) = f12*e123+f12
1468
 
{ .mfi
1469
 
          nop.m 999
1470
 
          fma.s1 POW_f123               = POW_e123,POW_f12,POW_f12
1471
 
          nop.i 999
1472
 
}
1473
 
;;
1474
 
 
1475
 
{ .mfi
1476
 
          nop.m 999
1477
 
          fma.s1 POW_T1T2               = POW_T1, POW_T2, f0
1478
 
          nop.i 999
1479
 
}
1480
 
{ .mfi
1481
 
          nop.m 999
1482
 
          fma.s1 POW_v3                 = POW_ssq, POW_Q4, POW_v4
1483
 
          cmp.ne p12,p13 = pow_GR_xneg_yodd, r0
1484
 
}
1485
 
;;
1486
 
 
1487
 
{ .mfi
1488
 
          nop.m 999
1489
 
          fma.s1 POW_v21ps              = POW_ssq, POW_v2, POW_1ps
1490
 
          nop.i 999
1491
 
}
1492
 
{ .mfi
1493
 
          nop.m 999
1494
 
          fma.s1 POW_s4                 = POW_ssq, POW_ssq, f0
1495
 
          nop.i 999
1496
 
}
1497
 
;;
1498
 
 
1499
 
{ .mfi
1500
 
          nop.m 999
1501
 
(p12)     fnma.s1 POW_A                 =  POW_2M, POW_f123, f0
1502
 
          nop.i 999
1503
 
}
1504
 
{ .mfi
1505
 
          nop.m 999
1506
 
(p13)     fma.s1 POW_A                  =  POW_2M, POW_f123, f0
1507
 
          cmp.eq p14,p11 = r0,r0   // Initialize p14 on, p11 off
1508
 
}
1509
 
;;
1510
 
 
1511
 
{ .mfi
1512
 
          nop.m 999
1513
 
          fmerge.s POW_abs_q = f0, POW_q // Form |q| so can test its size
1514
 
          nop.i 999
1515
 
}
1516
 
;;
1517
 
 
1518
 
{ .mfi
1519
 
(p10)     cmp.eq p0,p14 = r0,r0    // Turn off p14 if no overflow
1520
 
          fma.s1 POW_es                 = POW_s4,  POW_v3, POW_v21ps
1521
 
          nop.i 999
1522
 
}
1523
 
{ .mfi
1524
 
          nop.m 999
1525
 
          fma.s1 POW_A                  = POW_A, POW_T1T2, f0
1526
 
          nop.i 999
1527
 
}
1528
 
;;
1529
 
 
1530
 
{ .mfi
1531
 
// Test for |q| < 2^-63.  If so then reverse last two steps of the result
1532
 
// to avoid monotonicity problems for results near 1.0 in round up/down/zero.
1533
 
// p11 will be set if need to reverse the order, p14 if not.
1534
 
          nop.m 999
1535
 
(p10)     fcmp.lt.s0 p11,p14 = POW_abs_q, POW_2toM63 // Test |q| <2^-63
1536
 
          nop.i 999
1537
 
}
1538
 
;;
1539
 
 
1540
 
.pred.rel "mutex",p11,p14
1541
 
{ .mfi
1542
 
          nop.m 999
1543
 
(p14)     fma.s1 POW_A                  = POW_A, POW_es, f0
1544
 
          nop.i 999
1545
 
}
1546
 
{ .mfi
1547
 
          nop.m 999
1548
 
(p11)     fma.s1 POW_A                  = POW_A, POW_q, POW_A
1549
 
          nop.i 999
1550
 
}
1551
 
;;
1552
 
 
1553
 
// Dummy op to set inexact if |q| < 2^-63
1554
 
{ .mfi
1555
 
          nop.m 999
1556
 
(p11)     fma.d.s0 POW_tmp              = POW_A, POW_q, POW_A
1557
 
          nop.i 999
1558
 
}
1559
 
;;
1560
 
 
1561
 
{ .mfi
1562
 
          nop.m 999
1563
 
(p14)     fma.d.s0 f8                   = POW_A, POW_q, POW_A
1564
 
          nop.i 999
1565
 
}
1566
 
{ .mfb
1567
 
          nop.m 999
1568
 
(p11)     fma.d.s0 f8                   = POW_A, POW_es, f0
1569
 
(p10)     br.ret.sptk     b0            // Exit main branch if no over/underflow
1570
 
}
1571
 
;;
1572
 
 
1573
 
// POSSIBLE_OVER_UNDER
1574
 
// p6 = TRUE ==> Y_Gpr negative
1575
 
// Result is already computed.  We just need to know if over/underflow occurred.
1576
 
 
1577
 
{ .mfb
1578
 
        cmp.eq p0,p6                    = pow_GR_sign_Y_Gpr, r0
1579
 
        nop.f 999
1580
 
(p6)    br.cond.spnt POW_POSSIBLE_UNDER
1581
 
}
1582
 
;;
1583
 
 
1584
 
// POSSIBLE_OVER
1585
 
// We got an answer.
1586
 
// overflow is a possibility, not a certainty
1587
 
 
1588
 
 
1589
 
// We define an overflow when the answer with
1590
 
//    WRE set
1591
 
//    user-defined rounding mode
1592
 
 
1593
 
// double
1594
 
// Largest double is 7FE (biased double)
1595
 
//                   7FE - 3FF + FFFF = 103FE
1596
 
// Create + largest_double_plus_ulp
1597
 
// Create - largest_double_plus_ulp
1598
 
// Calculate answer with WRE set.
1599
 
 
1600
 
// single
1601
 
// Largest single is FE (biased double)
1602
 
//                   FE - 7F + FFFF = 1007E
1603
 
// Create + largest_single_plus_ulp
1604
 
// Create - largest_single_plus_ulp
1605
 
// Calculate answer with WRE set.
1606
 
 
1607
 
// Cases when answer is ldn+1  are as follows:
1608
 
//  ldn                   ldn+1
1609
 
// --+----------|----------+------------
1610
 
//              |
1611
 
//    +inf          +inf      -inf
1612
 
//                  RN         RN
1613
 
//                             RZ
1614
 
 
1615
 
// Put in s2 (td set, wre set)
1616
 
{ .mfi
1617
 
        nop.m 999
1618
 
        fsetc.s2 0x7F,0x42
1619
 
        nop.i 999
1620
 
}
1621
 
;;
1622
 
 
1623
 
{ .mfi
1624
 
        nop.m 999
1625
 
        fma.d.s2 POW_wre_urm_f8         = POW_A, POW_q, POW_A
1626
 
        nop.i 999
1627
 
}
1628
 
;;
1629
 
 
1630
 
// Return s2 to default
1631
 
{ .mfi
1632
 
        nop.m 999
1633
 
        fsetc.s2 0x7F,0x40
1634
 
        nop.i 999
1635
 
}
1636
 
;;
1637
 
 
1638
 
// p7 = TRUE ==> yes, we have an overflow
1639
 
{ .mfi
1640
 
        nop.m 999
1641
 
        fcmp.ge.s1 p7, p8               =  POW_wre_urm_f8, POW_big_pos
1642
 
        nop.i 999
1643
 
}
1644
 
;;
1645
 
 
1646
 
{ .mfi
1647
 
        nop.m 999
1648
 
(p8)    fcmp.le.s1 p7, p0               =  POW_wre_urm_f8, POW_big_neg
1649
 
        nop.i 999
1650
 
}
1651
 
;;
1652
 
 
1653
 
{ .mbb
1654
 
(p7)   mov pow_GR_tag                   = 24
1655
 
(p7)   br.cond.spnt __libm_error_region // Branch if overflow
1656
 
       br.ret.sptk     b0               // Exit if did not overflow
1657
 
}
1658
 
;;
1659
 
 
1660
 
// Here if |y*log(x)| < 2^(-11)
1661
 
// pow(x,y) ~ exp(d) ~ 1 + d + 0.5*d^2 + Q1*d^3 + Q2*d^4, where d = y*log(x)
1662
 
.align 32
1663
 
POW_NEAR_ONE:
1664
 
 
1665
 
{ .mfi
1666
 
          nop.m 999
1667
 
          fma.s1 POW_d2                 = POW_d, POW_d, f0
1668
 
          nop.i 999
1669
 
}
1670
 
;;
1671
 
 
1672
 
{ .mfi
1673
 
          nop.m 999
1674
 
          fma.s1 POW_poly_d_hi          = POW_d, POW_Q0_half, f1
1675
 
          nop.i 999
1676
 
}
1677
 
{ .mfi
1678
 
          nop.m 999
1679
 
          fma.s1 POW_poly_d_lo          = POW_d, POW_Q2, POW_Q1
1680
 
          nop.i 999
1681
 
}
1682
 
;;
1683
 
 
1684
 
{ .mfi
1685
 
          nop.m 999
1686
 
          fma.s1 POW_poly_d             = POW_d2, POW_poly_d_lo, POW_poly_d_hi
1687
 
          nop.i 999
1688
 
}
1689
 
;;
1690
 
 
1691
 
{ .mfb
1692
 
          nop.m 999
1693
 
          fma.d.s0 f8                   = POW_d, POW_poly_d, f1
1694
 
          br.ret.sptk b0 // exit function for arguments |y*log(x)| < 2^(-11)
1695
 
}
1696
 
;;
1697
 
 
1698
 
POW_POSSIBLE_UNDER:
1699
 
// We got an answer. input was < -2^9 but > -2^10 (double)
1700
 
// We got an answer. input was < -2^6 but > -2^7  (float)
1701
 
// underflow is a possibility, not a certainty
1702
 
 
1703
 
// We define an underflow when the answer with
1704
 
//    ftz set
1705
 
// is zero (tiny numbers become zero)
1706
 
// Notice (from below) that if we have an unlimited exponent range,
1707
 
// then there is an extra machine number E between the largest denormal and
1708
 
// the smallest normal.
1709
 
// So if with unbounded exponent we round to E or below, then we are
1710
 
// tiny and underflow has occurred.
1711
 
// But notice that you can be in a situation where we are tiny, namely
1712
 
// rounded to E, but when the exponent is bounded we round to smallest
1713
 
// normal. So the answer can be the smallest normal with underflow.
1714
 
//                           E
1715
 
// -----+--------------------+--------------------+-----
1716
 
//      |                    |                    |
1717
 
//   1.1...10 2^-3fff    1.1...11 2^-3fff    1.0...00 2^-3ffe
1718
 
//   0.1...11 2^-3ffe                                   (biased, 1)
1719
 
//    largest dn                               smallest normal
1720
 
 
1721
 
// Put in s2 (td set, ftz set)
1722
 
{ .mfi
1723
 
        nop.m 999
1724
 
        fsetc.s2 0x7F,0x41
1725
 
        nop.i 999
1726
 
}
1727
 
;;
1728
 
 
1729
 
{ .mfi
1730
 
        nop.m 999
1731
 
        fma.d.s2 POW_ftz_urm_f8         = POW_A, POW_q, POW_A
1732
 
        nop.i 999
1733
 
}
1734
 
;;
1735
 
 
1736
 
// Return s2 to default
1737
 
{ .mfi
1738
 
        nop.m 999
1739
 
        fsetc.s2 0x7F,0x40
1740
 
        nop.i 999
1741
 
}
1742
 
;;
1743
 
 
1744
 
// p7 = TRUE ==> yes, we have an underflow
1745
 
{ .mfi
1746
 
        nop.m 999
1747
 
        fcmp.eq.s1 p7, p0               =  POW_ftz_urm_f8, f0
1748
 
        nop.i 999
1749
 
}
1750
 
;;
1751
 
 
1752
 
{ .mbb
1753
 
(p7)    mov pow_GR_tag                  = 25
1754
 
(p7)    br.cond.spnt __libm_error_region // Branch if underflow
1755
 
        br.ret.sptk     b0               // Exit if did not underflow
1756
 
}
1757
 
;;
1758
 
 
1759
 
POW_X_DENORM:
1760
 
// Here if x unorm. Use the NORM_X for getf instructions, and then back
1761
 
// to normal path
1762
 
{ .mfi
1763
 
        getf.exp      pow_GR_signexp_X  = POW_NORM_X
1764
 
        nop.f 999
1765
 
        nop.i 999
1766
 
}
1767
 
;;
1768
 
 
1769
 
{ .mmi
1770
 
        getf.sig      pow_GR_sig_X      = POW_NORM_X
1771
 
;;
1772
 
        and           pow_GR_exp_X      = pow_GR_signexp_X, pow_GR_17ones
1773
 
        nop.i 999
1774
 
}
1775
 
;;
1776
 
 
1777
 
{ .mib
1778
 
        sub       pow_GR_true_exp_X     = pow_GR_exp_X, pow_GR_16ones
1779
 
        nop.i 999
1780
 
        br.cond.sptk    POW_COMMON
1781
 
}
1782
 
;;
1783
 
 
1784
 
POW_X_0:
1785
 
// Here if x=0 and y not nan
1786
 
//
1787
 
// We have the following cases:
1788
 
//  p6  x=0  and  y>0 and is an integer (may be even or odd)
1789
 
//  p7  x=0  and  y>0 and is NOT an integer, return +0
1790
 
//  p8  x=0  and  y>0 and so big as to always be an even integer, return +0
1791
 
//  p9  x=0  and  y>0 and may not be integer
1792
 
//  p10 x=0  and  y>0 and is an odd  integer, return x
1793
 
//  p11 x=0  and  y>0 and is an even integer, return +0
1794
 
//  p12 used in dummy fcmp to set denormal flag if y=unorm
1795
 
//  p13 x=0  and  y>0
1796
 
//  p14 x=0  and  y=0, branch to code for calling error handling
1797
 
//  p15 x=0  and  y<0, branch to code for calling error handling
1798
 
//
1799
 
{ .mfi
1800
 
        getf.sig pow_GR_sig_int_Y = POW_int_Y // Get signif of int_Y
1801
 
        fcmp.lt.s1 p15,p13 = f9, f0           // Test for y<0
1802
 
        and pow_GR_exp_Y = pow_GR_signexp_Y, pow_GR_17ones
1803
 
}
1804
 
{ .mfb
1805
 
        cmp.ne p14,p0 = pow_GR_y_zero,r0      // Test for y=0
1806
 
        fcvt.xf   POW_float_int_Y = POW_int_Y
1807
 
(p14)   br.cond.spnt POW_X_0_Y_0              // Branch if x=0 and y=0
1808
 
}
1809
 
;;
1810
 
 
1811
 
// If x=0 and y>0, test y and flag denormal
1812
 
{ .mfb
1813
 
(p13)   cmp.gt.unc p8,p9 = pow_GR_exp_Y, pow_GR_10033 // Test y +big = even int
1814
 
(p13)   fcmp.eq.s0 p12,p0 = f9,f0    // If x=0, y>0 dummy op to flag denormal
1815
 
(p15)   br.cond.spnt POW_X_0_Y_NEG // Branch if x=0 and y<0
1816
 
}
1817
 
;;
1818
 
 
1819
 
// Here if x=0 and y>0
1820
 
{ .mfi
1821
 
        nop.m 999
1822
 
(p9)    fcmp.eq.unc.s1 p6,p7 = POW_float_int_Y,  POW_NORM_Y // Test y=int
1823
 
        nop.i 999
1824
 
}
1825
 
{ .mfi
1826
 
        nop.m 999
1827
 
(p8)    fma.d.s0 f8 = f0,f0,f0 // If x=0, y>0 and large even int, return +0
1828
 
        nop.i 999
1829
 
}
1830
 
;;
1831
 
 
1832
 
{ .mfi
1833
 
        nop.m 999
1834
 
(p7)    fma.d.s0 f8  = f0,f0,f0   // Result +0 if x=0 and y>0 and not integer
1835
 
(p6)    tbit.nz.unc p10,p11 = pow_GR_sig_int_Y,0 // If y>0 int, test y even/odd
1836
 
}
1837
 
;;
1838
 
 
1839
 
// Note if x=0, y>0 and odd integer, just return x
1840
 
{ .mfb
1841
 
        nop.m 999
1842
 
(p11)   fma.d.s0 f8  = f0,f0,f0   // Result +0 if x=0 and y even integer
1843
 
        br.ret.sptk b0            // Exit if x=0 and y>0
1844
 
}
1845
 
;;
1846
 
 
1847
 
POW_X_0_Y_0:
1848
 
// When X is +-0 and Y is +-0, IEEE returns 1.0
1849
 
// We call error support with this value
1850
 
 
1851
 
{ .mfb
1852
 
        mov pow_GR_tag                  = 26
1853
 
        fma.d.s0 f8                     = f1,f1,f0
1854
 
        br.cond.sptk __libm_error_region
1855
 
}
1856
 
;;
1857
 
 
1858
 
POW_X_0_Y_NEG:
1859
 
// When X is +-0 and Y is negative, IEEE returns
1860
 
// X     Y           answer
1861
 
// +0    -odd int    +inf
1862
 
// -0    -odd int    -inf
1863
 
 
1864
 
// +0    !-odd int   +inf
1865
 
// -0    !-odd int   +inf
1866
 
 
1867
 
// p6 == Y is a floating point number outside the integer.
1868
 
//       Hence it is an integer and is even.
1869
 
//       return +inf
1870
 
 
1871
 
// p7 == Y is a floating point number within the integer range.
1872
 
//      p9  == (int_Y = NORM_Y), Y is an integer, which may be odd or even.
1873
 
//           p11 odd
1874
 
//              return (sign_of_x)inf
1875
 
//           p12 even
1876
 
//              return +inf
1877
 
//      p10 == Y is not an integer
1878
 
//         return +inf
1879
 
//
1880
 
 
1881
 
{ .mfi
1882
 
          nop.m 999
1883
 
          nop.f 999
1884
 
          cmp.gt  p6,p7                 = pow_GR_exp_Y, pow_GR_10033
1885
 
}
1886
 
;;
1887
 
 
1888
 
{ .mfi
1889
 
          mov pow_GR_tag                = 27
1890
 
(p7)      fcmp.eq.unc.s1 p9,p10         = POW_float_int_Y,  POW_NORM_Y
1891
 
          nop.i 999
1892
 
}
1893
 
;;
1894
 
 
1895
 
{ .mfb
1896
 
          nop.m 999
1897
 
(p6)      frcpa.s0 f8,p13               = f1, f0
1898
 
(p6)      br.cond.sptk __libm_error_region   // x=0, y<0, y large neg int
1899
 
}
1900
 
;;
1901
 
 
1902
 
{ .mfb
1903
 
          nop.m 999
1904
 
(p10)     frcpa.s0 f8,p13               = f1, f0
1905
 
(p10)     br.cond.sptk __libm_error_region   // x=0, y<0, y not int
1906
 
}
1907
 
;;
1908
 
 
1909
 
// x=0, y<0, y an int
1910
 
{ .mib
1911
 
          nop.m 999
1912
 
(p9)      tbit.nz.unc p11,p12           = pow_GR_sig_int_Y,0
1913
 
          nop.b 999
1914
 
}
1915
 
;;
1916
 
 
1917
 
{ .mfi
1918
 
          nop.m 999
1919
 
(p12)     frcpa.s0 f8,p13               = f1,f0
1920
 
          nop.i 999
1921
 
}
1922
 
;;
1923
 
 
1924
 
{ .mfb
1925
 
          nop.m 999
1926
 
(p11)     frcpa.s0 f8,p13               = f1,f8
1927
 
          br.cond.sptk __libm_error_region
1928
 
}
1929
 
;;
1930
 
 
1931
 
 
1932
 
POW_Y_0:
1933
 
// Here for y zero, x anything but zero and nan
1934
 
// Set flag if x denormal
1935
 
// Result is +1.0
1936
 
{ .mfi
1937
 
        nop.m 999
1938
 
        fcmp.eq.s0 p6,p0 = f8,f0    // Sets flag if x denormal
1939
 
        nop.i 999
1940
 
}
1941
 
{ .mfb
1942
 
        nop.m 999
1943
 
        fma.d.s0 f8 = f1,f1,f0
1944
 
        br.ret.sptk b0
1945
 
}
1946
 
;;
1947
 
 
1948
 
 
1949
 
POW_X_INF:
1950
 
// Here when X is +-inf
1951
 
 
1952
 
// X +inf  Y +inf             +inf
1953
 
// X -inf  Y +inf             +inf
1954
 
 
1955
 
// X +inf  Y >0               +inf
1956
 
// X -inf  Y >0, !odd integer +inf     <== (-inf)^0.5 = +inf !!
1957
 
// X -inf  Y >0,  odd integer -inf
1958
 
 
1959
 
// X +inf  Y -inf             +0
1960
 
// X -inf  Y -inf             +0
1961
 
 
1962
 
// X +inf  Y <0               +0
1963
 
// X -inf  Y <0, !odd integer +0
1964
 
// X -inf  Y <0, odd integer  -0
1965
 
 
1966
 
// X + inf Y=+0                +1
1967
 
// X + inf Y=-0                +1
1968
 
// X - inf Y=+0                +1
1969
 
// X - inf Y=-0                +1
1970
 
 
1971
 
// p13 == Y negative
1972
 
// p14 == Y positive
1973
 
 
1974
 
// p6 == Y is a floating point number outside the integer.
1975
 
//       Hence it is an integer and is even.
1976
 
//       p13 == (Y negative)
1977
 
//          return +inf
1978
 
//       p14 == (Y positive)
1979
 
//          return +0
1980
 
 
1981
 
// p7 == Y is a floating point number within the integer range.
1982
 
//      p9  == (int_Y = NORM_Y), Y is an integer, which may be odd or even.
1983
 
//           p11 odd
1984
 
//              p13 == (Y negative)
1985
 
//                 return (sign_of_x)inf
1986
 
//              p14 == (Y positive)
1987
 
//                 return (sign_of_x)0
1988
 
//           pxx even
1989
 
//              p13 == (Y negative)
1990
 
//                 return +inf
1991
 
//              p14 == (Y positive)
1992
 
//                 return +0
1993
 
 
1994
 
//      pxx == Y is not an integer
1995
 
//           p13 == (Y negative)
1996
 
//                 return +inf
1997
 
//           p14 == (Y positive)
1998
 
//                 return +0
1999
 
//
2000
 
 
2001
 
// If x=inf, test y and flag denormal
2002
 
{ .mfi
2003
 
          nop.m 999
2004
 
          fcmp.eq.s0 p10,p11 = f9,f0
2005
 
          nop.i 999
2006
 
}
2007
 
;;
2008
 
 
2009
 
{ .mfi
2010
 
          nop.m 999
2011
 
          fcmp.lt.s0 p13,p14            = POW_NORM_Y,f0
2012
 
          cmp.gt  p6,p7                 = pow_GR_exp_Y, pow_GR_10033
2013
 
}
2014
 
{ .mfi
2015
 
          nop.m 999
2016
 
          fclass.m p12,p0               = f9, 0x23 //@inf
2017
 
          nop.i 999
2018
 
}
2019
 
;;
2020
 
 
2021
 
{ .mfi
2022
 
          nop.m 999
2023
 
          fclass.m p15,p0               = f9, 0x07 //@zero
2024
 
          nop.i 999
2025
 
}
2026
 
;;
2027
 
 
2028
 
{ .mfb
2029
 
          nop.m 999
2030
 
(p15)     fmerge.s f8 = f1,f1      // Return +1.0 if x=inf, y=0
2031
 
(p15)     br.ret.spnt b0           // Exit if x=inf, y=0
2032
 
}
2033
 
;;
2034
 
 
2035
 
{ .mfi
2036
 
          nop.m 999
2037
 
(p14)     frcpa.s1 f8,p10 = f1,f0  // If x=inf, y>0, assume result +inf
2038
 
          nop.i 999
2039
 
}
2040
 
{ .mfb
2041
 
          nop.m 999
2042
 
(p13)     fma.d.s0 f8 = f0,f0,f0   // If x=inf, y<0, assume result +0.0
2043
 
(p12)     br.ret.spnt b0           // Exit if x=inf, y=inf
2044
 
}
2045
 
;;
2046
 
 
2047
 
// Here if x=inf, and 0 < |y| < inf.  Need to correct results if y odd integer.
2048
 
{ .mfi
2049
 
          nop.m 999
2050
 
(p7)      fcmp.eq.unc.s1 p9,p0 = POW_float_int_Y,  POW_NORM_Y // Is y integer?
2051
 
          nop.i 999
2052
 
}
2053
 
;;
2054
 
 
2055
 
{ .mfi
2056
 
          nop.m 999
2057
 
          nop.f 999
2058
 
(p9)      tbit.nz.unc p11,p0 = pow_GR_sig_int_Y,0  // Test for y odd integer
2059
 
}
2060
 
;;
2061
 
 
2062
 
{ .mfb
2063
 
          nop.m 999
2064
 
(p11)     fmerge.s f8 = POW_NORM_X,f8    // If y odd integer use sign of x
2065
 
          br.ret.sptk b0                 // Exit for x=inf, 0 < |y| < inf
2066
 
}
2067
 
;;
2068
 
 
2069
 
 
2070
 
POW_X_NEG_Y_NONINT:
2071
 
// When X is negative and Y is a non-integer, IEEE
2072
 
// returns a qnan indefinite.
2073
 
// We call error support with this value
2074
 
 
2075
 
{ .mfb
2076
 
         mov pow_GR_tag                 = 28
2077
 
         frcpa.s0 f8,p6                 = f0,f0
2078
 
         br.cond.sptk __libm_error_region
2079
 
}
2080
 
;;
2081
 
 
2082
 
POW_X_NAN:
2083
 
// Here if x=nan, y not nan
2084
 
{ .mfi
2085
 
         nop.m 999
2086
 
         fclass.m  p9,p13 = f9, 0x07 // Test y=zero
2087
 
         nop.i 999
2088
 
}
2089
 
;;
2090
 
 
2091
 
{ .mfb
2092
 
         nop.m 999
2093
 
(p13)    fma.d.s0 f8 = f8,f1,f0
2094
 
(p13)    br.ret.sptk  b0            // Exit if x nan, y anything but zero or nan
2095
 
}
2096
 
;;
2097
 
 
2098
 
POW_X_NAN_Y_0:
2099
 
// When X is a NAN and Y is zero, IEEE returns 1.
2100
 
// We call error support with this value.
2101
 
{ .mfi
2102
 
         nop.m 999
2103
 
         fcmp.eq.s0 p6,p0 = f8,f0       // Dummy op to set invalid on snan
2104
 
         nop.i 999
2105
 
}
2106
 
{ .mfb
2107
 
         mov pow_GR_tag                 = 29
2108
 
         fma.d.s0 f8 = f0,f0,f1
2109
 
         br.cond.sptk __libm_error_region
2110
 
}
2111
 
;;
2112
 
 
2113
 
 
2114
 
POW_OVER_UNDER_X_NOT_INF:
2115
 
 
2116
 
// p8 is TRUE for overflow
2117
 
// p9 is TRUE for underflow
2118
 
 
2119
 
// if y is infinity, we should not over/underflow
2120
 
 
2121
 
{ .mfi
2122
 
          nop.m 999
2123
 
          fcmp.eq.s1     p14, p13       = POW_xsq,f1  // Test |x|=1
2124
 
          cmp.eq p8,p9                  = pow_GR_sign_Y_Gpr, r0
2125
 
}
2126
 
;;
2127
 
 
2128
 
{ .mfi
2129
 
          nop.m 999
2130
 
(p14)     fclass.m.unc       p15, p0    = f9, 0x23 // If |x|=1, test y=inf
2131
 
          nop.i 999
2132
 
}
2133
 
{ .mfi
2134
 
          nop.m 999
2135
 
(p13)     fclass.m.unc       p11,p0     = f9, 0x23 // If |x| not 1, test y=inf
2136
 
          nop.i 999
2137
 
}
2138
 
;;
2139
 
 
2140
 
// p15 = TRUE if |x|=1, y=inf, return +1
2141
 
{ .mfb
2142
 
          nop.m 999
2143
 
(p15)     fma.d.s0          f8          = f1,f1,f0 // If |x|=1, y=inf, result +1
2144
 
(p15)     br.ret.spnt b0                // Exit if |x|=1, y=inf
2145
 
}
2146
 
;;
2147
 
 
2148
 
.pred.rel "mutex",p8,p9
2149
 
{  .mfb
2150
 
(p8)      setf.exp           f8 = pow_GR_17ones // If exp(+big), result inf
2151
 
(p9)      fmerge.s           f8 = f0,f0         // If exp(-big), result 0
2152
 
(p11)     br.ret.sptk b0                // Exit if |x| not 1, y=inf
2153
 
}
2154
 
;;
2155
 
 
2156
 
{ .mfb
2157
 
          nop.m 999
2158
 
          nop.f 999
2159
 
          br.cond.sptk POW_OVER_UNDER_ERROR // Branch if y not inf
2160
 
}
2161
 
;;
2162
 
 
2163
 
 
2164
 
POW_Y_NAN:
2165
 
// Here if y=nan, x anything
2166
 
// If x = +1 then result is +1, else result is quiet Y
2167
 
{ .mfi
2168
 
       nop.m 999
2169
 
       fcmp.eq.s1         p10,p9        = POW_NORM_X, f1
2170
 
       nop.i 999
2171
 
}
2172
 
;;
2173
 
 
2174
 
{ .mfi
2175
 
       nop.m 999
2176
 
(p10)  fcmp.eq.s0 p6,p0 = f9,f1   // Set invalid, even if x=+1
2177
 
       nop.i 999
2178
 
}
2179
 
;;
2180
 
 
2181
 
{ .mfi
2182
 
       nop.m 999
2183
 
(p10)  fma.d.s0 f8 = f1,f1,f0
2184
 
       nop.i 999
2185
 
}
2186
 
{ .mfb
2187
 
       nop.m 999
2188
 
(p9)   fma.d.s0 f8 = f9,f8,f0
2189
 
       br.ret.sptk b0             // Exit y=nan
2190
 
}
2191
 
;;
2192
 
 
2193
 
 
2194
 
POW_OVER_UNDER_ERROR:
2195
 
// Here if we have overflow or underflow.
2196
 
// Enter with p12 true if x negative and y odd int to force -0 or -inf
2197
 
 
2198
 
{ .mfi
2199
 
         sub   pow_GR_17ones_m1         = pow_GR_17ones, r0, 1
2200
 
         nop.f 999
2201
 
         mov pow_GR_one                 = 0x1
2202
 
}
2203
 
;;
2204
 
 
2205
 
// overflow, force inf with O flag
2206
 
{ .mmb
2207
 
(p8)     mov pow_GR_tag                 = 24
2208
 
(p8)     setf.exp POW_tmp               = pow_GR_17ones_m1
2209
 
         nop.b 999
2210
 
}
2211
 
;;
2212
 
 
2213
 
// underflow, force zero with I, U flags
2214
 
{ .mmi
2215
 
(p9)    mov pow_GR_tag                  = 25
2216
 
(p9)    setf.exp POW_tmp                = pow_GR_one
2217
 
        nop.i 999
2218
 
}
2219
 
;;
2220
 
 
2221
 
{ .mfi
2222
 
        nop.m 999
2223
 
        fma.d.s0 f8                     = POW_tmp, POW_tmp, f0
2224
 
        nop.i 999
2225
 
}
2226
 
;;
2227
 
 
2228
 
// p12 x is negative and y is an odd integer, change sign of result
2229
 
{ .mfi
2230
 
        nop.m 999
2231
 
(p12)   fnma.d.s0 f8                    = POW_tmp, POW_tmp, f0
2232
 
        nop.i 999
2233
 
}
2234
 
;;
2235
 
 
2236
 
GLOBAL_LIBM_END(pow)
2237
 
 
2238
 
 
2239
 
LOCAL_LIBM_ENTRY(__libm_error_region)
2240
 
 
2241
 
.prologue
2242
 
{ .mfi
2243
 
        add   GR_Parameter_Y=-32,sp     // Parameter 2 value
2244
 
        nop.f 0
2245
 
.save   ar.pfs,GR_SAVE_PFS
2246
 
        mov  GR_SAVE_PFS=ar.pfs         // Save ar.pfs
2247
 
}
2248
 
{ .mfi
2249
 
.fframe 64
2250
 
        add sp=-64,sp                   // Create new stack
2251
 
        nop.f 0
2252
 
        mov GR_SAVE_GP=gp               // Save gp
2253
 
};;
2254
 
 
2255
 
{ .mmi
2256
 
        stfd [GR_Parameter_Y] = POW_NORM_Y,16 // STORE Parameter 2 on stack
2257
 
        add GR_Parameter_X = 16,sp      // Parameter 1 address
2258
 
.save   b0, GR_SAVE_B0
2259
 
        mov GR_SAVE_B0=b0               // Save b0
2260
 
};;
2261
 
 
2262
 
.body
2263
 
{ .mib
2264
 
        stfd [GR_Parameter_X] = POW_NORM_X // STORE Parameter 1 on stack
2265
 
        add   GR_Parameter_RESULT = 0,GR_Parameter_Y    // Parameter 3 address
2266
 
        nop.b 0
2267
 
}
2268
 
{ .mib
2269
 
        stfd [GR_Parameter_Y] = f8      // STORE Parameter 3 on stack
2270
 
        add   GR_Parameter_Y = -16,GR_Parameter_Y
2271
 
        br.call.sptk b0=__libm_error_support# // Call error handling function
2272
 
};;
2273
 
 
2274
 
{ .mmi
2275
 
        add   GR_Parameter_RESULT = 48,sp
2276
 
        nop.m 0
2277
 
        nop.i 0
2278
 
};;
2279
 
 
2280
 
{ .mmi
2281
 
        ldfd  f8 = [GR_Parameter_RESULT] // Get return result off stack
2282
 
.restore sp
2283
 
        add   sp = 64,sp                 // Restore stack pointer
2284
 
        mov   b0 = GR_SAVE_B0            // Restore return address
2285
 
};;
2286
 
 
2287
 
{ .mib
2288
 
        mov   gp = GR_SAVE_GP            // Restore gp
2289
 
        mov   ar.pfs = GR_SAVE_PFS       // Restore ar.pfs
2290
 
        br.ret.sptk     b0               // Return
2291
 
};;
2292
 
 
2293
 
LOCAL_LIBM_END(__libm_error_region)
2294
 
 
2295
 
.type   __libm_error_support#,@function
2296
 
.global __libm_error_support#
2297