| // Copyright 2010 The Go Authors. All rights reserved. |
| // Use of this source code is governed by a BSD-style |
| // license that can be found in the LICENSE file. |
| |
| // The method is based on a paper by Naoki Shibata: "Efficient evaluation |
| // methods of elementary functions suitable for SIMD computation", Proc. |
| // of International Supercomputing Conference 2010 (ISC'10), pp. 25 -- 32 |
| // (May 2010). The paper is available at |
| // http://www.springerlink.com/content/340228x165742104/ |
| // |
| // The original code and the constants below are from the author's |
| // implementation available at http://freshmeat.net/projects/sleef. |
| // The README file says, "The software is in public domain. |
| // You can use the software without any obligation." |
| // |
| // This code is a simplified version of the original. |
| |
| #define LN2 0.6931471805599453094172321214581766 // log_e(2) |
| #define LOG2E 1.4426950408889634073599246810018920 // 1/LN2 |
| #define LN2U 0.69314718055966295651160180568695068359375 // upper half LN2 |
| #define LN2L 0.28235290563031577122588448175013436025525412068e-12 // lower half LN2 |
| |
| // func Exp(x float64) float64 |
| TEXT ·Exp(SB),7,$0 |
| // test bits for not-finite |
| MOVQ x+0(FP), AX |
| MOVQ $0x7ff0000000000000, BX |
| ANDQ BX, AX |
| CMPQ BX, AX |
| JEQ not_finite |
| MOVSD x+0(FP), X0 |
| MOVSD $LOG2E, X1 |
| MULSD X0, X1 |
| CVTTSD2SQ X1, BX // BX = exponent |
| CVTSQ2SD BX, X1 |
| MOVSD $LN2U, X2 |
| MULSD X1, X2 |
| SUBSD X2, X0 |
| MOVSD $LN2L, X2 |
| MULSD X1, X2 |
| SUBSD X2, X0 |
| // reduce argument |
| MOVSD $0.0625, X1 |
| MULSD X1, X0 |
| // Taylor series evaluation |
| MOVSD $2.4801587301587301587e-5, X1 |
| MULSD X0, X1 |
| MOVSD $1.9841269841269841270e-4, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $1.3888888888888888889e-3, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $8.3333333333333333333e-3, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $4.1666666666666666667e-2, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $1.6666666666666666667e-1, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $0.5, X2 |
| ADDSD X2, X1 |
| MULSD X0, X1 |
| MOVSD $1.0, X2 |
| ADDSD X2, X1 |
| MULSD X1, X0 |
| MOVSD $2.0, X1 |
| ADDSD X0, X1 |
| MULSD X1, X0 |
| MOVSD $2.0, X1 |
| ADDSD X0, X1 |
| MULSD X1, X0 |
| MOVSD $2.0, X1 |
| ADDSD X0, X1 |
| MULSD X1, X0 |
| MOVSD $2.0, X1 |
| ADDSD X0, X1 |
| MULSD X1, X0 |
| MOVSD $1.0, X1 |
| ADDSD X1, X0 |
| // return ldexp(fr, exp) |
| MOVQ $0x3ff, AX // bias + 1 |
| ADDQ AX, BX |
| MOVQ BX, X1 |
| MOVQ $52, AX // shift |
| MOVQ AX, X2 |
| PSLLQ X2, X1 |
| MULSD X1, X0 |
| MOVSD X0, r+8(FP) |
| RET |
| not_finite: |
| // test bits for -Inf |
| MOVQ x+0(FP), AX |
| MOVQ $0xfff0000000000000, BX |
| CMPQ BX, AX |
| JNE not_neginf |
| XORQ AX, AX |
| MOVQ AX, r+8(FP) |
| RET |
| not_neginf: |
| MOVQ AX, r+8(FP) |
| RET |
