diff options
| -rw-r--r-- | lib/libm/Makefile | 4 | ||||
| -rw-r--r-- | lib/libm/arch/amd64/e_exp.S | 101 |
2 files changed, 103 insertions, 2 deletions
diff --git a/lib/libm/Makefile b/lib/libm/Makefile index fa0196ba871..cd0a0799563 100644 --- a/lib/libm/Makefile +++ b/lib/libm/Makefile @@ -1,5 +1,5 @@ # $NetBSD: Makefile,v 1.28 1995/11/20 22:06:19 jtc Exp $ -# $OpenBSD: Makefile,v 1.38 2006/03/18 20:13:24 brad Exp $ +# $OpenBSD: Makefile,v 1.39 2006/03/19 00:01:04 kettenis Exp $ # # @(#)Makefile 5.1beta 93/09/24 # @@ -66,7 +66,7 @@ ARCH_SRCS = e_acos.S e_asin.S e_atan2.S e_exp.S e_fmod.S e_log.S e_log10.S \ .elif (${MACHINE_ARCH} == "amd64") .PATH: ${.CURDIR}/arch/amd64 CPPFLAGS+=-I${.CURDIR}/arch/amd64 -ARCH_SRCS = e_acos.S e_asin.S e_atan2.S e_fmod.S e_log.S e_log10.S \ +ARCH_SRCS = e_acos.S e_asin.S e_atan2.S e_exp.S e_fmod.S e_log.S e_log10.S \ e_remainder.S e_remainderf.S e_scalb.S e_sqrt.S e_sqrtf.S \ s_atan.S s_atanf.S s_ceil.S s_ceilf.S s_copysign.S s_copysignf.S \ s_cos.S s_cosf.S s_finite.S s_finitef.S s_floor.S s_floorf.S \ diff --git a/lib/libm/arch/amd64/e_exp.S b/lib/libm/arch/amd64/e_exp.S new file mode 100644 index 00000000000..5ae6516e991 --- /dev/null +++ b/lib/libm/arch/amd64/e_exp.S @@ -0,0 +1,101 @@ +/* $OpenBSD: e_exp.S,v 1.4 2006/03/19 00:01:04 kettenis Exp $ */ +/* $NetBSD: e_exp.S,v 1.12 2002/02/27 16:32:46 christos Exp $ */ + +/* + * Copyright (c) 1993,94 Winning Strategies, Inc. + * All rights reserved. + * + * Redistribution and use in source and binary forms, with or without + * modification, are permitted provided that the following conditions + * are met: + * 1. Redistributions of source code must retain the above copyright + * notice, this list of conditions and the following disclaimer. + * 2. Redistributions in binary form must reproduce the above copyright + * notice, this list of conditions and the following disclaimer in the + * documentation and/or other materials provided with the distribution. + * 3. All advertising materials mentioning features or use of this software + * must display the following acknowledgement: + * This product includes software developed by Winning Strategies, Inc. + * 4. The name of the author may not be used to endorse or promote products + * derived from this software without specific prior written permission. + * + * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR + * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES + * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. + * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT, + * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT + * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, + * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY + * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT + * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF + * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + */ + +/* + * Written by: + * J.T. Conklin (jtc@wimsey.com), Winning Strategies, Inc. + */ + +#include <machine/asm.h> + +#include "abi.h" + +/* e^x = 2^(x * log2(e)) */ +ENTRY(__ieee754_exp) + XMM_ONE_ARG_DOUBLE_PROLOGUE + /* + * If x is +-Inf, then the subtraction would give Inf-Inf = NaN. + * Avoid this. Also avoid it if x is NaN for convenience. + */ + movl -4(%rsp),%eax + andl $0x7fffffff,%eax + cmpl $0x7ff00000,%eax + jae x_Inf_or_NaN + + fldl ARG_DOUBLE_ONE + + /* + * Ensure that the rounding mode is to nearest (to give the smallest + * possible fraction) and that the precision is as high as possible. + * We may as well mask interrupts if we switch the mode. + */ + fstcw -8(%rsp) + movl -8(%rsp),%eax + andl $0x0300,%eax + cmpl $0x0300,%eax /* RC == 0 && PC == 3? */ + je 1f /* jump if mode is good */ + movl $0x137f,-4(%rsp) + fldcw -4(%rsp) +1: + fldl2e + fmulp /* x * log2(e) */ + fst %st(1) + frndint /* int(x * log2(e)) */ + fst %st(2) + fsubrp /* fract(x * log2(e)) */ + f2xm1 /* 2^(fract(x * log2(e))) - 1 */ + fld1 + faddp /* 2^(fract(x * log2(e))) */ + fscale /* e^x */ + fstp %st(1) + je 1f + fldcw -8(%rsp) +1: + XMM_DOUBLE_EPILOGUE + ret + +x_Inf_or_NaN: + /* + * Return 0 if x is -Inf. Otherwise just return x, although the + * C version would return (x + x) (Real Indefinite) if x is a NaN. + */ + cmpl $0xfff00000,-4(%rsp) + jne x_not_minus_Inf + cmpl $0,-8(%rsp) + jne x_not_minus_Inf + xorpd %xmm0,%xmm0 + ret + +x_not_minus_Inf: + movsd ARG_DOUBLE_ONE,%xmm0 + ret |