Split the implementation file per architecture, and use a BSD/Sun implementation of rint/rintf for arm, as the previous arm implementation used signed 32 bit integer intermediates.
Signed-off-by: Martin Storsjö <[email protected]> --- mingw-w64-crt/Makefile.am | 10 +++- mingw-w64-crt/math/arm/s_rint.c | 86 +++++++++++++++++++++++++++ mingw-w64-crt/math/arm/s_rintf.c | 51 ++++++++++++++++ mingw-w64-crt/math/arm64/rint.c | 12 ++++ mingw-w64-crt/math/arm64/rintf.c | 12 ++++ mingw-w64-crt/math/bsd_private_base.h | 22 +++++++ mingw-w64-crt/math/rint.c | 24 -------- mingw-w64-crt/math/rintf.c | 23 ------- mingw-w64-crt/math/x86/rint.c | 12 ++++ mingw-w64-crt/math/x86/rintf.c | 12 ++++ 10 files changed, 214 insertions(+), 50 deletions(-) create mode 100644 mingw-w64-crt/math/arm/s_rint.c create mode 100644 mingw-w64-crt/math/arm/s_rintf.c create mode 100644 mingw-w64-crt/math/arm64/rint.c create mode 100644 mingw-w64-crt/math/arm64/rintf.c delete mode 100644 mingw-w64-crt/math/rint.c delete mode 100644 mingw-w64-crt/math/rintf.c create mode 100644 mingw-w64-crt/math/x86/rint.c create mode 100644 mingw-w64-crt/math/x86/rintf.c diff --git a/mingw-w64-crt/Makefile.am b/mingw-w64-crt/Makefile.am index dc52cd6f0..e07feba74 100644 --- a/mingw-w64-crt/Makefile.am +++ b/mingw-w64-crt/Makefile.am @@ -412,7 +412,7 @@ src_libmingwex=\ math/lrint.c math/lrintf.c math/lrintl.c math/lround.c math/lroundf.c math/lroundl.c \ math/modf.c math/modff.c math/modfl.c math/nextafterf.c math/nextafterl.c math/nexttoward.c \ math/nexttowardf.c math/powf.c math/powi.c math/powif.c math/powil.c \ - math/rint.c math/rintf.c math/rintl.c math/round.c math/roundf.c \ + math/rintl.c math/round.c math/roundf.c \ math/roundl.c math/s_erf.c math/sf_erf.c math/signbit.c math/signbitf.c math/signbitl.c \ math/sinhf.c math/sinhl.c math/sqrt.c \ math/sqrtf.c math/sqrtl.c math/tanhf.c math/tanhl.c math/tgamma.c \ @@ -480,6 +480,7 @@ src_libmingwex_x86=\ math/x86/logbf.c math/x86/logbl.c math/x86/logl.c math/x86/nearbyint.S math/x86/nearbyintf.S \ math/x86/nearbyintl.S math/x86/pow.c math/x86/pow.def.h math/x86/powl.c math/x86/remainder.S \ math/x86/remainderf.S math/x86/remainderl.S math/x86/remquo.S math/x86/remquof.S math/x86/remquol.S \ + math/x86/rint.c math/x86/rintf.c \ math/x86/scalbn.S math/x86/scalbnf.S math/x86/scalbnl.S math/x86/sin.c math/x86/sin.def.h \ math/x86/sinf.c math/x86/sinl.c math/x86/sinl_internal.S math/x86/tanf.c math/x86/tanl.S \ math/x86/trunc.S math/x86/truncf.S @@ -490,7 +491,8 @@ src_libmingwex64=$(src_libmingwex_x86) # these only go into the ARM32 version: src_libmingwexarm32=\ - math/arm/_chgsignl.S + math/arm/_chgsignl.S \ + math/arm/s_rint.c math/arm/s_rintf.c if ENABLE_SOFTMATH src_libmingwexarm32+=\ @@ -517,7 +519,9 @@ endif # these only go into the ARM64 version: src_libmingwexarm64=\ - math/arm64/_chgsignl.S math/arm-common/ldexpl.c math/arm-common/sincos.c + math/arm64/_chgsignl.S \ + math/arm64/rint.c math/arm64/rintf.c \ + math/arm-common/ldexpl.c math/arm-common/sincos.c # These intrinsics are target independent: diff --git a/mingw-w64-crt/math/arm/s_rint.c b/mingw-w64-crt/math/arm/s_rint.c new file mode 100644 index 000000000..823a2e0c9 --- /dev/null +++ b/mingw-w64-crt/math/arm/s_rint.c @@ -0,0 +1,86 @@ +/* @(#)s_rint.c 5.1 93/09/24 */ +/* + * ==================================================== + * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. + * + * Developed at SunPro, a Sun Microsystems, Inc. business. + * Permission to use, copy, modify, and distribute this + * software is freely granted, provided that this notice + * is preserved. + * ==================================================== + */ + +#include <sys/cdefs.h> + +/* + * rint(x) + * Return x rounded to integral value according to the prevailing + * rounding mode. + * Method: + * Using floating addition. + * Exception: + * Inexact flag raised if x not equal to rint(x). + */ + +#include <float.h> + +#include "../bsd_private_base.h" + +static const double +TWO52[2]={ + 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */ + -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */ +}; + +double +rint(double x) +{ + int32_t i0,j0,sx; + u_int32_t i,i1; + double w,t; + EXTRACT_WORDS(i0,i1,x); + sx = (i0>>31)&1; + j0 = ((i0>>20)&0x7ff)-0x3ff; + if(j0<20) { + if(j0<0) { + if(((i0&0x7fffffff)|i1)==0) return x; + i1 |= (i0&0x0fffff); + i0 &= 0xfffe0000; + i0 |= ((i1|-i1)>>12)&0x80000; + SET_HIGH_WORD(x,i0); + STRICT_ASSIGN(double,w,TWO52[sx]+x); + t = w-TWO52[sx]; + GET_HIGH_WORD(i0,t); + SET_HIGH_WORD(t,(i0&0x7fffffff)|(sx<<31)); + return t; + } else { + i = (0x000fffff)>>j0; + if(((i0&i)|i1)==0) return x; /* x is integral */ + i>>=1; + if(((i0&i)|i1)!=0) { + /* + * Some bit is set after the 0.5 bit. To avoid the + * possibility of errors from double rounding in + * w = TWO52[sx]+x, adjust the 0.25 bit to a lower + * guard bit. We do this for all j0<=51. The + * adjustment is trickiest for j0==18 and j0==19 + * since then it spans the word boundary. + */ + if(j0==19) i1 = 0x40000000; else + if(j0==18) i1 = 0x80000000; else + i0 = (i0&(~i))|((0x20000)>>j0); + } + } + } else if (j0>51) { + if(j0==0x400) return x+x; /* inf or NaN */ + else return x; /* x is integral */ + } else { + i = ((u_int32_t)(0xffffffff))>>(j0-20); + if((i1&i)==0) return x; /* x is integral */ + i>>=1; + if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-20)); + } + INSERT_WORDS(x,i0,i1); + STRICT_ASSIGN(double,w,TWO52[sx]+x); + return w-TWO52[sx]; +} diff --git a/mingw-w64-crt/math/arm/s_rintf.c b/mingw-w64-crt/math/arm/s_rintf.c new file mode 100644 index 000000000..0ec01704a --- /dev/null +++ b/mingw-w64-crt/math/arm/s_rintf.c @@ -0,0 +1,51 @@ +/* s_rintf.c -- float version of s_rint.c. + * Conversion to float by Ian Lance Taylor, Cygnus Support, [email protected]. + */ + +/* + * ==================================================== + * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. + * + * Developed at SunPro, a Sun Microsystems, Inc. business. + * Permission to use, copy, modify, and distribute this + * software is freely granted, provided that this notice + * is preserved. + * ==================================================== + */ + +#include <sys/cdefs.h> + +#include <float.h> +#include <stdint.h> + +#include "../bsd_private_base.h" + +static const float +TWO23[2]={ + 8.3886080000e+06, /* 0x4b000000 */ + -8.3886080000e+06, /* 0xcb000000 */ +}; + +float +rintf(float x) +{ + int32_t i0,j0,sx; + float w,t; + GET_FLOAT_WORD(i0,x); + sx = (i0>>31)&1; + j0 = ((i0>>23)&0xff)-0x7f; + if(j0<23) { + if(j0<0) { + if((i0&0x7fffffff)==0) return x; + STRICT_ASSIGN(float,w,TWO23[sx]+x); + t = w-TWO23[sx]; + GET_FLOAT_WORD(i0,t); + SET_FLOAT_WORD(t,(i0&0x7fffffff)|(sx<<31)); + return t; + } + STRICT_ASSIGN(float,w,TWO23[sx]+x); + return w-TWO23[sx]; + } + if(j0==0x80) return x+x; /* inf or NaN */ + else return x; /* x is integral */ +} diff --git a/mingw-w64-crt/math/arm64/rint.c b/mingw-w64-crt/math/arm64/rint.c new file mode 100644 index 000000000..03a7ffa0a --- /dev/null +++ b/mingw-w64-crt/math/arm64/rint.c @@ -0,0 +1,12 @@ +/** + * This file has no copyright assigned and is placed in the Public Domain. + * This file is part of the mingw-w64 runtime package. + * No warranty is given; refer to the file DISCLAIMER.PD within this package. + */ +#include <math.h> + +double rint (double x) { + double retval = 0.0; + __asm__ __volatile__ ("frintx %d0, %d1\n\t" : "=w" (retval) : "w" (x)); + return retval; +} diff --git a/mingw-w64-crt/math/arm64/rintf.c b/mingw-w64-crt/math/arm64/rintf.c new file mode 100644 index 000000000..914769a1d --- /dev/null +++ b/mingw-w64-crt/math/arm64/rintf.c @@ -0,0 +1,12 @@ +/** + * This file has no copyright assigned and is placed in the Public Domain. + * This file is part of the mingw-w64 runtime package. + * No warranty is given; refer to the file DISCLAIMER.PD within this package. + */ +#include <math.h> + +float rintf (float x) { + float retval = 0.0F; + __asm__ __volatile__ ("frintx %s0, %s1\n\t" : "=w" (retval) : "w" (x)); + return retval; +} diff --git a/mingw-w64-crt/math/bsd_private_base.h b/mingw-w64-crt/math/bsd_private_base.h index f31c75d00..de5b93512 100644 --- a/mingw-w64-crt/math/bsd_private_base.h +++ b/mingw-w64-crt/math/bsd_private_base.h @@ -10,6 +10,7 @@ */ #include <inttypes.h> +#include <float.h> typedef unsigned int u_int32_t; @@ -99,3 +100,24 @@ do { \ gf_u.word = (i); \ (d) = gf_u.value; \ } while(0) + + +#ifdef FLT_EVAL_METHOD +/* + * Attempt to get strict C99 semantics for assignment with non-C99 compilers. + */ +#if FLT_EVAL_METHOD == 0 || __GNUC__ == 0 +#define STRICT_ASSIGN(type, lval, rval) ((lval) = (rval)) +#else +#define STRICT_ASSIGN(type, lval, rval) do { \ + volatile type __lval; \ + \ + if (sizeof(type) >= sizeof(long double)) \ + (lval) = (rval); \ + else { \ + __lval = (rval); \ + (lval) = __lval; \ + } \ +} while (0) +#endif +#endif /* FLT_EVAL_METHOD */ diff --git a/mingw-w64-crt/math/rint.c b/mingw-w64-crt/math/rint.c deleted file mode 100644 index d88331492..000000000 --- a/mingw-w64-crt/math/rint.c +++ /dev/null @@ -1,24 +0,0 @@ -/** - * This file has no copyright assigned and is placed in the Public Domain. - * This file is part of the mingw-w64 runtime package. - * No warranty is given; refer to the file DISCLAIMER.PD within this package. - */ -#include <math.h> - -double rint (double x) { - double retval = 0.0; -#if defined(_AMD64_) || defined(__x86_64__) || defined(_X86_) || defined(__i386__) - __asm__ __volatile__ ("frndint;" : "=t" (retval) : "0" (x)); -#elif defined(__arm__) || defined(_ARM_) - if (isnan(x) || isinf(x)) - return x; - float temp; - __asm__ __volatile__ ( - "vcvtr.s32.f64 %[tmp], %[src]\n\t" - "vcvt.f64.s32 %[dst], %[tmp]\n\t" - : [dst] "=w" (retval), [tmp] "=t" (temp) : [src] "w" (x)); -#elif defined(__aarch64__) || defined(_ARM64_) - __asm__ __volatile__ ("frintx %d0, %d1\n\t" : "=w" (retval) : "w" (x)); -#endif - return retval; -} diff --git a/mingw-w64-crt/math/rintf.c b/mingw-w64-crt/math/rintf.c deleted file mode 100644 index bfae51516..000000000 --- a/mingw-w64-crt/math/rintf.c +++ /dev/null @@ -1,23 +0,0 @@ -/** - * This file has no copyright assigned and is placed in the Public Domain. - * This file is part of the mingw-w64 runtime package. - * No warranty is given; refer to the file DISCLAIMER.PD within this package. - */ -#include <math.h> - -float rintf (float x) { - float retval = 0.0F; -#if defined(_AMD64_) || defined(__x86_64__) || defined(_X86_) || defined(__i386__) - __asm__ __volatile__ ("frndint;": "=t" (retval) : "0" (x)); -#elif defined(__arm__) || defined(_ARM_) - if (isnan(x) || isinf(x)) - return x; - __asm__ __volatile__ ( - "vcvtr.s32.f32 %[dst], %[src]\n\t" - "vcvt.f32.s32 %[dst], %[dst]\n\t" - : [dst] "=t" (retval) : [src] "w" (x)); -#elif defined(__aarch64__) || defined(_ARM64_) - __asm__ __volatile__ ("frintx %s0, %s1\n\t" : "=w" (retval) : "w" (x)); -#endif - return retval; -} diff --git a/mingw-w64-crt/math/x86/rint.c b/mingw-w64-crt/math/x86/rint.c new file mode 100644 index 000000000..55a73d123 --- /dev/null +++ b/mingw-w64-crt/math/x86/rint.c @@ -0,0 +1,12 @@ +/** + * This file has no copyright assigned and is placed in the Public Domain. + * This file is part of the mingw-w64 runtime package. + * No warranty is given; refer to the file DISCLAIMER.PD within this package. + */ +#include <math.h> + +double rint (double x) { + double retval = 0.0; + __asm__ __volatile__ ("frndint;" : "=t" (retval) : "0" (x)); + return retval; +} diff --git a/mingw-w64-crt/math/x86/rintf.c b/mingw-w64-crt/math/x86/rintf.c new file mode 100644 index 000000000..f5a6bf2b1 --- /dev/null +++ b/mingw-w64-crt/math/x86/rintf.c @@ -0,0 +1,12 @@ +/** + * This file has no copyright assigned and is placed in the Public Domain. + * This file is part of the mingw-w64 runtime package. + * No warranty is given; refer to the file DISCLAIMER.PD within this package. + */ +#include <math.h> + +float rintf (float x) { + float retval = 0.0F; + __asm__ __volatile__ ("frndint;": "=t" (retval) : "0" (x)); + return retval; +} -- 2.17.1 _______________________________________________ Mingw-w64-public mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/mingw-w64-public
