diff options
| -rw-r--r-- | ChangeLog | 10 | ||||
| -rw-r--r-- | NEWS | 2 | ||||
| -rw-r--r-- | math/libm-test.inc | 15 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/libm-test-ulps | 86 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/s_cexp.S | 253 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/s_cexpf.S | 257 | ||||
| -rw-r--r-- | sysdeps/i386/fpu/s_cexpl.S | 256 | ||||
| -rw-r--r-- | sysdeps/x86_64/fpu/libm-test-ulps | 22 |
8 files changed, 114 insertions, 787 deletions
diff --git a/ChangeLog b/ChangeLog index 75dbb36091..ebffe4f5d2 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,13 @@ +2012-03-21 Joseph Myers <joseph@codesourcery.com> + + [BZ #13883] + * sysdeps/i386/fpu/s_cexp.S: Remove. + * sysdeps/i386/fpu/s_cexpf.S: Likewise. + * sysdeps/i386/fpu/s_cexpl.S: Likewise. + * math/libm-test.inc (cexp_test): Add more tests. + * sysdeps/i386/fpu/libm-test-ulps: Update. + * sysdeps/x86_64/fpu/libm-test-ulps: Likewise. + 2012-03-21 Allan McRae <allan@archlinux.org> * timezone/Makefile: Do not install iso3166.tab and zone.tab diff --git a/NEWS b/NEWS index a8a3e5716c..5fc4444641 100644 --- a/NEWS +++ b/NEWS @@ -16,7 +16,7 @@ Version 2.16 13526, 13527, 13528, 13529, 13530, 13531, 13532, 13533, 13547, 13551, 13552, 13553, 13555, 13559, 13566, 13583, 13618, 13637, 13656, 13658, 13673, 13695, 13704, 13706, 13726, 13738, 13786, 13792, 13806, 13840, - 13841, 13844, 13846, 13851, 13852, 13854, 13871 + 13841, 13844, 13846, 13851, 13852, 13854, 13871, 13883 * ISO C11 support: diff --git a/math/libm-test.inc b/math/libm-test.inc index afc6f55d7b..05a000e0c1 100644 --- a/math/libm-test.inc +++ b/math/libm-test.inc @@ -1894,6 +1894,21 @@ cexp_test (void) TEST_c_c (cexp, 0.75L, 1.25L, 0.667537446429131586942201977015932112L, 2.00900045494094876258347228145863909L); TEST_c_c (cexp, -2.0, -3.0, -0.13398091492954261346140525546115575L, -0.019098516261135196432576240858800925L); + TEST_c_c (cexp, 0, 0x1p65, 0.99888622066058013610642172179340364209972L, -0.047183876212354673805106149805700013943218L); + TEST_c_c (cexp, 0, -0x1p65, 0.99888622066058013610642172179340364209972L, 0.047183876212354673805106149805700013943218L); + TEST_c_c (cexp, 50, 0x1p127, 4.053997150228616856622417636046265337193e21L, 3.232070315463388524466674772633810238819e21L); + +#ifndef TEST_FLOAT + TEST_c_c (cexp, 0, 1e22, 0.5232147853951389454975944733847094921409L, -0.8522008497671888017727058937530293682618L); + TEST_c_c (cexp, 0, 0x1p1023, -0.826369834614147994500785680811743734805L, 0.5631277798508840134529434079444683477104L); + TEST_c_c (cexp, 500, 0x1p1023, -1.159886268932754433233243794561351783426e217L, 7.904017694554466595359379965081774849708e216L); +#endif + +#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384 + TEST_c_c (cexp, 0, 0x1p16383L, 0.9210843909921906206874509522505756251609L, 0.3893629985894208126948115852610595405563L); + TEST_c_c (cexp, -10000, 0x1p16383L, 1.045876464564882298442774542991176546722e-4343L, 4.421154026488516836023811173959413420548e-4344L); +#endif + END (cexp, complex); } diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps index 1d87514e9b..4d61635f25 100644 --- a/sysdeps/i386/fpu/libm-test-ulps +++ b/sysdeps/i386/fpu/libm-test-ulps @@ -431,15 +431,44 @@ idouble: 1 ifloat: 1 # cexp +Test "Real part of: cexp (-10000 + 0x1p16383 i) == 1.045876464564882298442774542991176546722e-4343 + 4.421154026488516836023811173959413420548e-4344 i": +ildouble: 1 +ldouble: 1 Test "Real part of: cexp (-2.0 - 3.0 i) == -0.13398091492954261346140525546115575 - 0.019098516261135196432576240858800925 i": ildouble: 1 ldouble: 1 Test "Imaginary part of: cexp (-2.0 - 3.0 i) == -0.13398091492954261346140525546115575 - 0.019098516261135196432576240858800925 i": +float: 1 +ifloat: 1 ildouble: 1 ldouble: 1 +Test "Imaginary part of: cexp (0 + 0x1p65 i) == 0.99888622066058013610642172179340364209972 - 0.047183876212354673805106149805700013943218 i": +float: 1 +ifloat: 1 +Test "Imaginary part of: cexp (0 - 0x1p65 i) == 0.99888622066058013610642172179340364209972 + 0.047183876212354673805106149805700013943218 i": +float: 1 +ifloat: 1 Test "Real part of: cexp (0.75 + 1.25 i) == 0.667537446429131586942201977015932112 + 2.00900045494094876258347228145863909 i": +float: 1 +ifloat: 1 ildouble: 1 ldouble: 1 +Test "Imaginary part of: cexp (0.75 + 1.25 i) == 0.667537446429131586942201977015932112 + 2.00900045494094876258347228145863909 i": +ildouble: 1 +ldouble: 1 +Test "Real part of: cexp (50 + 0x1p127 i) == 4.053997150228616856622417636046265337193e21 + 3.232070315463388524466674772633810238819e21 i": +double: 2 +idouble: 2 +Test "Imaginary part of: cexp (50 + 0x1p127 i) == 4.053997150228616856622417636046265337193e21 + 3.232070315463388524466674772633810238819e21 i": +double: 1 +float: 1 +idouble: 1 +ifloat: 1 +ildouble: 1 +ldouble: 1 +Test "Real part of: cexp (500 + 0x1p1023 i) == -1.159886268932754433233243794561351783426e217 + 7.904017694554466595359379965081774849708e216 i": +double: 1 +idouble: 1 # clog Test "Real part of: clog (0.75 + 1.25 i) == 0.376885901188190075998919126749298416 + 1.03037682652431246378774332703115153 i": @@ -815,18 +844,22 @@ ifloat: 1 ildouble: 1 ldouble: 1 Test "Imaginary part of: cpow (0.75 + 1.25 i, 0.0 + 1.0 i) == 0.331825439177608832276067945276730566 + 0.131338600281188544930936345230903032 i": +float: 1 +ifloat: 1 ildouble: 1 ldouble: 1 Test "Real part of: cpow (0.75 + 1.25 i, 0.75 + 1.25 i) == 0.117506293914473555420279832210420483 + 0.346552747708338676483025352060418001 i": -float: 3 -ifloat: 3 +double: 1 +float: 4 +idouble: 1 +ifloat: 4 ildouble: 6 ldouble: 6 Test "Imaginary part of: cpow (0.75 + 1.25 i, 0.75 + 1.25 i) == 0.117506293914473555420279832210420483 + 0.346552747708338676483025352060418001 i": float: 1 ifloat: 1 -ildouble: 1 -ldouble: 1 +ildouble: 2 +ldouble: 2 Test "Real part of: cpow (0.75 + 1.25 i, 1.0 + 0.0 i) == 0.75 + 1.25 i": ildouble: 1 ldouble: 1 @@ -834,22 +867,27 @@ Test "Imaginary part of: cpow (0.75 + 1.25 i, 1.0 + 0.0 i) == 0.75 + 1.25 i": float: 1 ifloat: 1 Test "Real part of: cpow (0.75 + 1.25 i, 1.0 + 1.0 i) == 0.0846958290317209430433805274189191353 + 0.513285749182902449043287190519090481 i": -double: 1 +double: 2 float: 3 -idouble: 1 +idouble: 2 ifloat: 3 ildouble: 3 ldouble: 3 +Test "Real part of: cpow (2 + 0 i, 10 + 0 i) == 1024.0 + 0.0 i": +ildouble: 1 +ldouble: 1 Test "Real part of: cpow (2 + 3 i, 4 + 0 i) == -119.0 - 120.0 i": double: 1 -float: 4 +float: 5 idouble: 1 -ifloat: 4 +ifloat: 5 +ildouble: 1 +ldouble: 1 Test "Imaginary part of: cpow (2 + 3 i, 4 + 0 i) == -119.0 - 120.0 i": -float: 1 -ifloat: 1 -ildouble: 2 -ldouble: 2 +float: 2 +ifloat: 2 +ildouble: 4 +ldouble: 4 Test "Imaginary part of: cpow (e + 0 i, 0 + 2 * M_PIl i) == 1.0 + 0.0 i": double: 2 float: 3 @@ -2124,10 +2162,18 @@ ildouble: 1 ldouble: 1 Function: Real part of "cexp": +double: 2 +float: 1 +idouble: 2 +ifloat: 1 ildouble: 1 ldouble: 1 Function: Imaginary part of "cexp": +double: 1 +float: 1 +idouble: 1 +ifloat: 1 ildouble: 1 ldouble: 1 @@ -2212,20 +2258,20 @@ double: 1 ldouble: 1 Function: Real part of "cpow": -double: 1 -float: 4 -idouble: 1 -ifloat: 4 -ildouble: 6 -ldouble: 6 +double: 2 +float: 5 +idouble: 2 +ifloat: 5 +ildouble: 5 +ldouble: 5 Function: Imaginary part of "cpow": double: 2 float: 3 idouble: 2 ifloat: 3 -ildouble: 2 -ldouble: 2 +ildouble: 4 +ldouble: 4 Function: Real part of "csin": float: 1 diff --git a/sysdeps/i386/fpu/s_cexp.S b/sysdeps/i386/fpu/s_cexp.S deleted file mode 100644 index e5fdb7d735..0000000000 --- a/sysdeps/i386/fpu/s_cexp.S +++ /dev/null @@ -1,253 +0,0 @@ -/* ix87 specific implementation of complex exponential function for double. - Copyright (C) 1997, 2005, 2012 Free Software Foundation, Inc. - This file is part of the GNU C Library. - Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. - - The GNU C Library is free software; you can redistribute it and/or - modify it under the terms of the GNU Lesser General Public - License as published by the Free Software Foundation; either - version 2.1 of the License, or (at your option) any later version. - - The GNU C Library is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, see - <http://www.gnu.org/licenses/>. */ - -#include <sysdep.h> - - .section .rodata - - .align ALIGNARG(4) - ASM_TYPE_DIRECTIVE(huge_nan_null_null,@object) -huge_nan_null_null: - .byte 0, 0, 0, 0, 0, 0, 0xf0, 0x7f - .byte 0, 0, 0, 0, 0, 0, 0xff, 0x7f - .double 0.0 -zero: .double 0.0 -infinity: - .byte 0, 0, 0, 0, 0, 0, 0xf0, 0x7f - .byte 0, 0, 0, 0, 0, 0, 0xff, 0x7f - .double 0.0 - .byte 0, 0, 0, 0, 0, 0, 0, 0x80 - ASM_SIZE_DIRECTIVE(huge_nan_null_null) - - ASM_TYPE_DIRECTIVE(twopi,@object) -twopi: - .byte 0x35, 0xc2, 0x68, 0x21, 0xa2, 0xda, 0xf, 0xc9, 0x1, 0x40 - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(twopi) - - ASM_TYPE_DIRECTIVE(l2e,@object) -l2e: - .byte 0xbc, 0xf0, 0x17, 0x5c, 0x29, 0x3b, 0xaa, 0xb8, 0xff, 0x3f - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(l2e) - - ASM_TYPE_DIRECTIVE(one,@object) -one: .double 1.0 - ASM_SIZE_DIRECTIVE(one) - - -#ifdef PIC -#define MO(op) op##@GOTOFF(%ecx) -#define MOX(op,x,f) op##@GOTOFF(%ecx,x,f) -#else -#define MO(op) op -#define MOX(op,x,f) op(,x,f) -#endif - - .text -ENTRY(__cexp) - fldl 8(%esp) /* x */ - fxam - fnstsw - fldl 16(%esp) /* y : x */ -#ifdef PIC - LOAD_PIC_REG (cx) -#endif - movb %ah, %dh - andb $0x45, %ah - cmpb $0x05, %ah - je 1f /* Jump if real part is +-Inf */ - cmpb $0x01, %ah - je 2f /* Jump if real part is NaN */ - - fxam /* y : x */ - fnstsw - /* If the imaginary part is not finite we return NaN+i NaN, as - for the case when the real part is NaN. A test for +-Inf and - NaN would be necessary. But since we know the stack register - we applied `fxam' to is not empty we can simply use one test. - Check your FPU manual for more information. */ - andb $0x01, %ah - cmpb $0x01, %ah - je 20f - - /* We have finite numbers in the real and imaginary part. Do - the real work now. */ - fxch /* x : y */ - fldt MO(l2e) /* log2(e) : x : y */ - fmulp /* x * log2(e) : y */ - fld %st /* x * log2(e) : x * log2(e) : y */ - frndint /* int(x * log2(e)) : x * log2(e) : y */ - fsubr %st, %st(1) /* int(x * log2(e)) : frac(x * log2(e)) : y */ - fxch /* frac(x * log2(e)) : int(x * log2(e)) : y */ - f2xm1 /* 2^frac(x * log2(e))-1 : int(x * log2(e)) : y */ - faddl MO(one) /* 2^frac(x * log2(e)) : int(x * log2(e)) : y */ - fscale /* e^x : int(x * log2(e)) : y */ - fst %st(1) /* e^x : e^x : y */ - fxch %st(2) /* y : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 7f - fmulp %st, %st(3) /* sin(y) : e^x : e^x * cos(y) */ - fmulp %st, %st(1) /* e^x * sin(y) : e^x * cos(y) */ - movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpl 8(%eax) - fstpl (%eax) - ret $4 - - /* We have to reduce the argument to fsincos. */ - .align ALIGNARG(4) -7: fldt MO(twopi) /* 2*pi : y : e^x : e^x */ - fxch /* y : 2*pi : e^x : e^x */ -8: fprem1 /* y%(2*pi) : 2*pi : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 8b - fstp %st(1) /* y%(2*pi) : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fmulp %st, %st(3) - fmulp %st, %st(1) - movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpl 8(%eax) - fstpl (%eax) - ret $4 - - /* The real part is +-inf. We must make further differences. */ - .align ALIGNARG(4) -1: fxam /* y : x */ - fnstsw - movb %ah, %dl - testb $0x01, %ah /* See above why 0x01 is usable here. */ - jne 3f - - - /* The real part is +-Inf and the imaginary part is finite. */ - andl $0x245, %edx - cmpb $0x40, %dl /* Imaginary part == 0? */ - je 4f /* Yes -> */ - - fxch /* x : y */ - shrl $5, %edx - fstp %st(0) /* y */ /* Drop the real part. */ - andl $16, %edx /* This puts the sign bit of the real part - in bit 4. So we can use it to index a - small array to select 0 or Inf. */ - fsincos /* cos(y) : sin(y) */ - fnstsw - testl $0x0400, %eax - jnz 5f - fldl MOX(huge_nan_null_null,%edx,1) - movl 4(%esp), %edx /* Pointer to memory for result. */ - fstl 8(%edx) - fstpl (%edx) - ftst - fnstsw - shll $23, %eax - andl $0x80000000, %eax - orl %eax, 4(%edx) - fstp %st(0) - ftst - fnstsw - shll $23, %eax - andl $0x80000000, %eax - orl %eax, 12(%edx) - fstp %st(0) - ret $4 - /* We must reduce the argument to fsincos. */ - .align ALIGNARG(4) -5: fldt MO(twopi) - fxch -6: fprem1 - fnstsw - testl $0x400, %eax - jnz 6b - fstp %st(1) - fsincos - fldl MOX(huge_nan_null_null,%edx,1) - movl 4(%esp), %edx /* Pointer to memory for result. */ - fstl 8(%edx) - fstpl (%edx) - ftst - fnstsw - shll $23, %eax - andl $0x80000000, %eax - orl %eax, 4(%edx) - fstp %st(0) - ftst - fnstsw - shll $23, %eax - andl $0x80000000, %eax - orl %eax, 12(%edx) - fstp %st(0) - ret $4 - - /* The real part is +-Inf and the imaginary part is +-0. So return - +-Inf+-0i. */ - .align ALIGNARG(4) -4: movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpl 8(%eax) - shrl $5, %edx - fstp %st(0) - andl $16, %edx - fldl MOX(huge_nan_null_null,%edx,1) - fstpl (%eax) - ret $4 - - /* The real part is +-Inf, the imaginary is also is not finite. */ - .align ALIGNARG(4) -3: fstp %st(0) - fstp %st(0) /* <empty> */ - andb $0x45, %ah - andb $0x47, %dh - xorb %dh, %ah - jnz 30f - fldl MO(infinity) /* Raise invalid exception. */ - fmull MO(zero) - fstp %st(0) -30: movl %edx, %eax - shrl $5, %edx - shll $4, %eax - andl $16, %edx - andl $32, %eax - orl %eax, %edx - movl 4(%esp), %eax /* Pointer to memory for result. */ - - fldl MOX(huge_nan_null_null,%edx,1) - fldl MOX(huge_nan_null_null+8,%edx,1) - fxch - fstpl (%eax) - fstpl 8(%eax) - ret $4 - - /* The real part is NaN. */ - .align ALIGNARG(4) -20: fldl MO(infinity) /* Raise invalid exception. */ - fmull MO(zero) - fstp %st(0) -2: fstp %st(0) - fstp %st(0) - movl 4(%esp), %eax /* Pointer to memory for result. */ - fldl MO(huge_nan_null_null+8) - fstl (%eax) - fstpl 8(%eax) - ret $4 - -END(__cexp) -weak_alias (__cexp, cexp) diff --git a/sysdeps/i386/fpu/s_cexpf.S b/sysdeps/i386/fpu/s_cexpf.S deleted file mode 100644 index 6ed66e6d04..0000000000 --- a/sysdeps/i386/fpu/s_cexpf.S +++ /dev/null @@ -1,257 +0,0 @@ -/* ix87 specific implementation of complex exponential function for double. - Copyright (C) 1997, 2005, 2012 Free Software Foundation, Inc. - This file is part of the GNU C Library. - Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. - - The GNU C Library is free software; you can redistribute it and/or - modify it under the terms of the GNU Lesser General Public - License as published by the Free Software Foundation; either - version 2.1 of the License, or (at your option) any later version. - - The GNU C Library is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, see - <http://www.gnu.org/licenses/>. */ - -#include <sysdep.h> - - .section .rodata - - .align ALIGNARG(4) - ASM_TYPE_DIRECTIVE(huge_nan_null_null,@object) -huge_nan_null_null: - .byte 0, 0, 0x80, 0x7f - .byte 0, 0, 0xc0, 0x7f - .float 0.0 -zero: .float 0.0 -infinity: - .byte 0, 0, 0x80, 0x7f - .byte 0, 0, 0xc0, 0x7f - .float 0.0 - .byte 0, 0, 0, 0x80 - ASM_SIZE_DIRECTIVE(huge_nan_null_null) - - ASM_TYPE_DIRECTIVE(twopi,@object) -twopi: - .byte 0x35, 0xc2, 0x68, 0x21, 0xa2, 0xda, 0xf, 0xc9, 0x1, 0x40 - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(twopi) - - ASM_TYPE_DIRECTIVE(l2e,@object) -l2e: - .byte 0xbc, 0xf0, 0x17, 0x5c, 0x29, 0x3b, 0xaa, 0xb8, 0xff, 0x3f - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(l2e) - - ASM_TYPE_DIRECTIVE(one,@object) -one: .double 1.0 - ASM_SIZE_DIRECTIVE(one) - - -#ifdef PIC -#define MO(op) op##@GOTOFF(%ecx) -#define MOX(op,x,f) op##@GOTOFF(%ecx,x,f) -#else -#define MO(op) op -#define MOX(op,x,f) op(,x,f) -#endif - - .text -ENTRY(__cexpf) - flds 4(%esp) /* x */ - fxam - fnstsw - flds 8(%esp) /* y : x */ -#ifdef PIC - LOAD_PIC_REG (cx) -#endif - movb %ah, %dh - andb $0x45, %ah - cmpb $0x05, %ah - je 1f /* Jump if real part is +-Inf */ - cmpb $0x01, %ah - je 2f /* Jump if real part is NaN */ - - fxam /* y : x */ - fnstsw - /* If the imaginary part is not finite we return NaN+i NaN, as - for the case when the real part is NaN. A test for +-Inf and - NaN would be necessary. But since we know the stack register - we applied `fxam' to is not empty we can simply use one test. - Check your FPU manual for more information. */ - andb $0x01, %ah - cmpb $0x01, %ah - je 20f - - /* We have finite numbers in the real and imaginary part. Do - the real work now. */ - fxch /* x : y */ - fldt MO(l2e) /* log2(e) : x : y */ - fmulp /* x * log2(e) : y */ - fld %st /* x * log2(e) : x * log2(e) : y */ - frndint /* int(x * log2(e)) : x * log2(e) : y */ - fsubr %st, %st(1) /* int(x * log2(e)) : frac(x * log2(e)) : y */ - fxch /* frac(x * log2(e)) : int(x * log2(e)) : y */ - f2xm1 /* 2^frac(x * log2(e))-1 : int(x * log2(e)) : y */ - faddl MO(one) /* 2^frac(x * log2(e)) : int(x * log2(e)) : y */ - fscale /* e^x : int(x * log2(e)) : y */ - fst %st(1) /* e^x : e^x : y */ - fxch %st(2) /* y : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 7f - fmulp %st, %st(3) /* sin(y) : e^x : e^x * cos(y) */ - fmulp %st, %st(1) /* e^x * sin(y) : e^x * cos(y) */ - subl $8, %esp - cfi_adjust_cfa_offset (8) - fstps 4(%esp) - fstps (%esp) - popl %eax - cfi_adjust_cfa_offset (-4) - popl %edx - cfi_adjust_cfa_offset (-4) - ret - - /* We have to reduce the argument to fsincos. */ - .align ALIGNARG(4) -7: fldt MO(twopi) /* 2*pi : y : e^x : e^x */ - fxch /* y : 2*pi : e^x : e^x */ -8: fprem1 /* y%(2*pi) : 2*pi : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 8b - fstp %st(1) /* y%(2*pi) : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fmulp %st, %st(3) - fmulp %st, %st(1) - subl $8, %esp - cfi_adjust_cfa_offset (8) - fstps 4(%esp) - fstps (%esp) - popl %eax - cfi_adjust_cfa_offset (-4) - popl %edx - cfi_adjust_cfa_offset (-4) - ret - - /* The real part is +-inf. We must make further differences. */ - .align ALIGNARG(4) -1: fxam /* y : x */ - fnstsw - movb %ah, %dl - testb $0x01, %ah /* See above why 0x01 is usable here. */ - jne 3f - - - /* The real part is +-Inf and the imaginary part is finite. */ - andl $0x245, %edx - cmpb $0x40, %dl /* Imaginary part == 0? */ - je 4f /* Yes -> */ - - fxch /* x : y */ - shrl $6, %edx - fstp %st(0) /* y */ /* Drop the real part. */ - andl $8, %edx /* This puts the sign bit of the real part - in bit 3. So we can use it to index a - small array to select 0 or Inf. */ - fsincos /* cos(y) : sin(y) */ - fnstsw - testl $0x0400, %eax - jnz 5f - fxch - ftst - fnstsw - fstp %st(0) - shll $23, %eax - andl $0x80000000, %eax - orl MOX(huge_nan_null_null,%edx,1), %eax - movl MOX(huge_nan_null_null,%edx,1), %ecx - movl %eax, %edx - ftst - fnstsw - fstp %st(0) - shll $23, %eax - andl $0x80000000, %eax - orl %ecx, %eax - ret - /* We must reduce the argument to fsincos. */ - .align ALIGNARG(4) -5: fldt MO(twopi) - fxch -6: fprem1 - fnstsw - testl $0x400, %eax - jnz 6b - fstp %st(1) - fsincos - fxch - ftst - fnstsw - fstp %st(0) - shll $23, %eax - andl $0x80000000, %eax - orl MOX(huge_nan_null_null,%edx,1), %eax - movl MOX(huge_nan_null_null,%edx,1), %ecx - movl %eax, %edx - ftst - fnstsw - fstp %st(0) - shll $23, %eax - andl $0x80000000, %eax - orl %ecx, %eax - ret - - /* The real part is +-Inf and the imaginary part is +-0. So return - +-Inf+-0i. */ - .align ALIGNARG(4) -4: subl $4, %esp - cfi_adjust_cfa_offset (4) - fstps (%esp) - shrl $6, %edx - fstp %st(0) - andl $8, %edx - movl MOX(huge_nan_null_null,%edx,1), %eax - popl %edx - cfi_adjust_cfa_offset (-4) - ret - - /* The real part is +-Inf, the imaginary is also is not finite. */ - .align ALIGNARG(4) -3: fstp %st(0) - fstp %st(0) /* <empty> */ - andb $0x45, %ah - andb $0x47, %dh - xorb %dh, %ah - jnz 30f - flds MO(infinity) /* Raise invalid exception. */ - fmuls MO(zero) - fstp %st(0) -30: movl %edx, %eax - shrl $6, %edx - shll $3, %eax - andl $8, %edx - andl $16, %eax - orl %eax, %edx - - movl MOX(huge_nan_null_null,%edx,1), %eax - movl MOX(huge_nan_null_null+4,%edx,1), %edx - ret - - /* The real part is NaN. */ - .align ALIGNARG(4) -20: flds MO(infinity) /* Raise invalid exception. */ - fmuls MO(zero) - fstp %st(0) -2: fstp %st(0) - fstp %st(0) - movl MO(huge_nan_null_null+4), %eax - movl %eax, %edx - ret - -END(__cexpf) -weak_alias (__cexpf, cexpf) diff --git a/sysdeps/i386/fpu/s_cexpl.S b/sysdeps/i386/fpu/s_cexpl.S deleted file mode 100644 index ab02a172ad..0000000000 --- a/sysdeps/i386/fpu/s_cexpl.S +++ /dev/null @@ -1,256 +0,0 @@ -/* ix87 specific implementation of complex exponential function for double. - Copyright (C) 1997, 2005, 2012 Free Software Foundation, Inc. - This file is part of the GNU C Library. - Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997. - - The GNU C Library is free software; you can redistribute it and/or - modify it under the terms of the GNU Lesser General Public - License as published by the Free Software Foundation; either - version 2.1 of the License, or (at your option) any later version. - - The GNU C Library is distributed in the hope that it will be useful, - but WITHOUT ANY WARRANTY; without even the implied warranty of - MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU - Lesser General Public License for more details. - - You should have received a copy of the GNU Lesser General Public - License along with the GNU C Library; if not, see - <http://www.gnu.org/licenses/>. */ - -#include <sysdep.h> - - .section .rodata - - .align ALIGNARG(4) - ASM_TYPE_DIRECTIVE(huge_nan_null_null,@object) -huge_nan_null_null: - .byte 0, 0, 0, 0, 0, 0, 0xf0, 0x7f - .byte 0, 0, 0, 0, 0, 0, 0xff, 0x7f - .double 0.0 -zero: .double 0.0 -infinity: - .byte 0, 0, 0, 0, 0, 0, 0xf0, 0x7f - .byte 0, 0, 0, 0, 0, 0, 0xff, 0x7f - .double 0.0 - .byte 0, 0, 0, 0, 0, 0, 0, 0x80 - ASM_SIZE_DIRECTIVE(huge_nan_null_null) - - ASM_TYPE_DIRECTIVE(twopi,@object) -twopi: - .byte 0x35, 0xc2, 0x68, 0x21, 0xa2, 0xda, 0xf, 0xc9, 0x1, 0x40 - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(twopi) - - ASM_TYPE_DIRECTIVE(l2e,@object) -l2e: - .byte 0xbc, 0xf0, 0x17, 0x5c, 0x29, 0x3b, 0xaa, 0xb8, 0xff, 0x3f - .byte 0, 0, 0, 0, 0, 0 - ASM_SIZE_DIRECTIVE(l2e) - - ASM_TYPE_DIRECTIVE(one,@object) -one: .double 1.0 - ASM_SIZE_DIRECTIVE(one) - - -#ifdef PIC -#define MO(op) op##@GOTOFF(%ecx) -#define MOX(op,x,f) op##@GOTOFF(%ecx,x,f) -#else -#define MO(op) op -#define MOX(op,x,f) op(,x,f) -#endif - - .text -ENTRY(__cexpl) - fldt 8(%esp) /* x */ - fxam - fnstsw - fldt 20(%esp) /* y : x */ -#ifdef PIC - LOAD_PIC_REG (cx) -#endif - movb %ah, %dh - andb $0x45, %ah - cmpb $0x05, %ah - je 1f /* Jump if real part is +-Inf */ - cmpb $0x01, %ah - je 2f /* Jump if real part is NaN */ - - fxam /* y : x */ - fnstsw - /* If the imaginary part is not finite we return NaN+i NaN, as - for the case when the real part is NaN. A test for +-Inf and - NaN would be necessary. But since we know the stack register - we applied `fxam' to is not empty we can simply use one test. - Check your FPU manual for more information. */ - andb $0x01, %ah - cmpb $0x01, %ah - je 20f - - /* We have finite numbers in the real and imaginary part. Do - the real work now. */ - fxch /* x : y */ - fldt MO(l2e) /* log2(e) : x : y */ - fmulp /* x * log2(e) : y */ - fld %st /* x * log2(e) : x * log2(e) : y */ - frndint /* int(x * log2(e)) : x * log2(e) : y */ - fsubr %st, %st(1) /* int(x * log2(e)) : frac(x * log2(e)) : y */ - fxch /* frac(x * log2(e)) : int(x * log2(e)) : y */ - f2xm1 /* 2^frac(x * log2(e))-1 : int(x * log2(e)) : y */ - faddl MO(one) /* 2^frac(x * log2(e)) : int(x * log2(e)) : y */ - fscale /* e^x : int(x * log2(e)) : y */ - fst %st(1) /* e^x : e^x : y */ - fxch %st(2) /* y : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 7f - fmulp %st, %st(3) /* sin(y) : e^x : e^x * cos(y) */ - fmulp %st, %st(1) /* e^x * sin(y) : e^x * cos(y) */ - movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpt 12(%eax) - fstpt (%eax) - ret $4 - - /* We have to reduce the argument to fsincos. */ - .align ALIGNARG(4) -7: fldt MO(twopi) /* 2*pi : y : e^x : e^x */ - fxch /* y : 2*pi : e^x : e^x */ -8: fprem1 /* y%(2*pi) : 2*pi : e^x : e^x */ - fnstsw - testl $0x400, %eax - jnz 8b - fstp %st(1) /* y%(2*pi) : e^x : e^x */ - fsincos /* cos(y) : sin(y) : e^x : e^x */ - fmulp %st, %st(3) - fmulp %st, %st(1) - movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpt 12(%eax) - fstpt (%eax) - ret $4 - - /* The real part is +-inf. We must make further differences. */ - .align ALIGNARG(4) -1: fxam /* y : x */ - fnstsw - movb %ah, %dl - testb $0x01, %ah /* See above why 0x01 is usable here. */ - jne 3f - - - /* The real part is +-Inf and the imaginary part is finite. */ - andl $0x245, %edx - cmpb $0x40, %dl /* Imaginary part == 0? */ - je 4f /* Yes -> */ - - fxch /* x : y */ - shrl $5, %edx - fstp %st(0) /* y */ /* Drop the real part. */ - andl $16, %edx /* This puts the sign bit of the real part - in bit 4. So we can use it to index a - small array to select 0 or Inf. */ - fsincos /* cos(y) : sin(y) */ - fnstsw - testl $0x0400, %eax - jnz 5f - fldl MOX(huge_nan_null_null,%edx,1) - movl 4(%esp), %edx /* Pointer to memory for result. */ - fld %st - fstpt 12(%edx) - fstpt (%edx) - ftst - fnstsw - shll $7, %eax - andl $0x8000, %eax - orl %eax, 8(%edx) - fstp %st(0) - ftst - fnstsw - shll $7, %eax - andl $0x8000, %eax - orl %eax, 20(%edx) - fstp %st(0) - ret $4 - /* We must reduce the argument to fsincos. */ - .align ALIGNARG(4) -5: fldt MO(twopi) - fxch -6: fprem1 - fnstsw - testl $0x400, %eax - jnz 6b - fstp %st(1) - fsincos - fldl MOX(huge_nan_null_null,%edx,1) - movl 4(%esp), %edx /* Pointer to memory for result. */ - fld %st - fstpt 12(%edx) - fstpt (%edx) - ftst - fnstsw - shll $7, %eax - andl $0x8000, %eax - orl %eax, 8(%edx) - fstp %st(0) - ftst - fnstsw - shll $7, %eax - andl $0x8000, %eax - orl %eax, 20(%edx) - fstp %st(0) - ret $4 - - /* The real part is +-Inf and the imaginary part is +-0. So return - +-Inf+-0i. */ - .align ALIGNARG(4) -4: movl 4(%esp), %eax /* Pointer to memory for result. */ - fstpt 12(%eax) - shrl $5, %edx - fstp %st(0) - andl $16, %edx - fldl MOX(huge_nan_null_null,%edx,1) - fstpt (%eax) - ret $4 - - /* The real part is +-Inf, the imaginary is also is not finite. */ - .align ALIGNARG(4) -3: fstp %st(0) - fstp %st(0) /* <empty> */ - andb $0x45, %ah - andb $0x47, %dh - xorb %dh, %ah - jnz 30f - fldl MO(infinity) /* Raise invalid exception. */ - fmull MO(zero) - fstp %st(0) -30: movl %edx, %eax - shrl $5, %edx - shll $4, %eax - andl $16, %edx - andl $32, %eax - orl %eax, %edx - movl 4(%esp), %eax /* Pointer to memory for result. */ - - fldl MOX(huge_nan_null_null,%edx,1) - fldl MOX(huge_nan_null_null+8,%edx,1) - fxch - fstpt (%eax) - fstpt 12(%eax) - ret $4 - - /* The real part is NaN. */ - .align ALIGNARG(4) -20: fldl MO(infinity) /* Raise invalid exception. */ - fmull MO(zero) - fstp %st(0) -2: fstp %st(0) - fstp %st(0) - movl 4(%esp), %eax /* Pointer to memory for result. */ - fldl MO(huge_nan_null_null+8) - fld %st(0) - fstpt (%eax) - fstpt 12(%eax) - ret $4 - -END(__cexpl) -weak_alias (__cexpl, cexpl) diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps index fef6ea1a8a..e5eb8f9379 100644 --- a/sysdeps/x86_64/fpu/libm-test-ulps +++ b/sysdeps/x86_64/fpu/libm-test-ulps @@ -482,6 +482,9 @@ float: 1 ifloat: 1 # cexp +Test "Real part of: cexp (-10000 + 0x1p16383 i) == 1.045876464564882298442774542991176546722e-4343 + 4.421154026488516836023811173959413420548e-4344 i": +ildouble: 1 +ldouble: 1 Test "Imaginary part of: cexp (-2.0 - 3.0 i) == -0.13398091492954261346140525546115575 - 0.019098516261135196432576240858800925 i": float: 1 ifloat: 1 @@ -491,6 +494,19 @@ ifloat: 1 Test "Imaginary part of: cexp (0.75 + 1.25 i) == 0.667537446429131586942201977015932112 + 2.00900045494094876258347228145863909 i": ildouble: 1 ldouble: 1 +Test "Real part of: cexp (50 + 0x1p127 i) == 4.053997150228616856622417636046265337193e21 + 3.232070315463388524466674772633810238819e21 i": +double: 2 +float: 1 +idouble: 2 +ifloat: 1 +Test "Imaginary part of: cexp (50 + 0x1p127 i) == 4.053997150228616856622417636046265337193e21 + 3.232070315463388524466674772633810238819e21 i": +double: 1 +idouble: 1 +ildouble: 1 +ldouble: 1 +Test "Real part of: cexp (500 + 0x1p1023 i) == -1.159886268932754433233243794561351783426e217 + 7.904017694554466595359379965081774849708e216 i": +double: 1 +idouble: 1 # clog Test "Imaginary part of: clog (-2 - 3 i) == 1.2824746787307683680267437207826593 - 2.1587989303424641704769327722648368 i": @@ -2090,11 +2106,17 @@ ildouble: 1 ldouble: 1 Function: Real part of "cexp": +double: 2 float: 1 +idouble: 2 ifloat: 1 +ildouble: 1 +ldouble: 1 Function: Imaginary part of "cexp": +double: 1 float: 1 +idouble: 1 ifloat: 1 ildouble: 1 ldouble: 1 |