On Thu, Oct 14, 2021 at 01:15:56AM +0200, Mark Kettenis wrote:
> Currently the lib/libm/msun/run-lrint_test regress fails on powerpc64
> and other platforms. Our implementation came from NetBSD, but NetBSD
> switched to the implementation from FreeBSD some time ago. That is
> the same implementation that we already use for lrintl(3) and
> llrintl(3).
>
> Diff below makes us use that implementation for lrint(3), lrintf(3),
> llrint(3) and llrintf(3) as well. This makes the regress test pass on
> powerpc64.
>
> ok?
libm Tests on powerpc64 pass now.
amd64 and i386 still passing.
arm64 and sparc64 still failing.
So it is an improvement.
OK bluhm@
> Index: lib/libm/src/s_llrint.c
> ===================================================================
> RCS file: /cvs/src/lib/libm/src/s_llrint.c,v
> retrieving revision 1.6
> diff -u -p -r1.6 s_llrint.c
> --- lib/libm/src/s_llrint.c 12 Sep 2016 19:47:02 -0000 1.6
> +++ lib/libm/src/s_llrint.c 13 Oct 2021 23:12:11 -0000
> @@ -1,14 +1,12 @@
> -/* $OpenBSD: s_llrint.c,v 1.6 2016/09/12 19:47:02 guenther Exp $ */
> -/* $NetBSD: llrint.c,v 1.2 2004/10/13 15:18:32 drochner Exp $ */
> +/* $OpenBSD$ */
>
> /*
> - * Written by Matthias Drochner <droch...@netbsd.org>.
> - * Public domain.
> + * Written by Martynas Venckus. Public domain
> */
>
> -#define LRINTNAME llrint
> -#define RESTYPE long long int
> -#define RESTYPE_MIN LLONG_MIN
> -#define RESTYPE_MAX LLONG_MAX
> +#define type double
> +#define roundit rint
> +#define dtype long long
> +#define fn llrint
>
> #include "s_lrint.c"
> Index: lib/libm/src/s_llrintf.c
> ===================================================================
> RCS file: /cvs/src/lib/libm/src/s_llrintf.c,v
> retrieving revision 1.2
> diff -u -p -r1.2 s_llrintf.c
> --- lib/libm/src/s_llrintf.c 25 Sep 2006 22:16:48 -0000 1.2
> +++ lib/libm/src/s_llrintf.c 13 Oct 2021 23:12:11 -0000
> @@ -1,14 +1,12 @@
> -/* $OpenBSD: s_llrintf.c,v 1.2 2006/09/25 22:16:48 kettenis Exp $ */
> -/* $NetBSD: llrintf.c,v 1.2 2004/10/13 15:18:32 drochner Exp $ */
> +/* $OpenBSD$ */
>
> /*
> - * Written by Matthias Drochner <droch...@netbsd.org>.
> - * Public domain.
> + * Written by Martynas Venckus. Public domain
> */
>
> -#define LRINTNAME llrintf
> -#define RESTYPE long long int
> -#define RESTYPE_MIN LLONG_MIN
> -#define RESTYPE_MAX LLONG_MAX
> +#define type float
> +#define roundit rintf
> +#define dtype long long
> +#define fn llrintf
>
> #include "s_lrintf.c"
> Index: lib/libm/src/s_lrint.c
> ===================================================================
> RCS file: /cvs/src/lib/libm/src/s_lrint.c,v
> retrieving revision 1.11
> diff -u -p -r1.11 s_lrint.c
> --- lib/libm/src/s_lrint.c 12 Sep 2016 19:47:02 -0000 1.11
> +++ lib/libm/src/s_lrint.c 13 Oct 2021 23:12:11 -0000
> @@ -1,9 +1,8 @@
> -/* $OpenBSD: s_lrint.c,v 1.11 2016/09/12 19:47:02 guenther Exp $ */
> -/* $NetBSD: lrint.c,v 1.3 2004/10/13 15:18:32 drochner Exp $ */
> +/* $OpenBSD$ */
>
> /*-
> - * Copyright (c) 2004
> - * Matthias Drochner. All rights reserved.
> + * Copyright (c) 2005 David Schultz <d...@freebsd.org>
> + * All rights reserved.
> *
> * Redistribution and use in source and binary forms, with or without
> * modification, are permitted provided that the following conditions
> @@ -27,75 +26,35 @@
> * SUCH DAMAGE.
> */
>
> -#include <sys/types.h>
> -#include <sys/limits.h>
> -#include <float.h>
> +#include <fenv.h>
> #include <math.h>
> -#include <ieeefp.h>
> -#include <machine/ieee.h>
>
> -#include "math_private.h"
> -
> -#ifndef LRINTNAME
> -#define LRINTNAME lrint
> -#define RESTYPE long int
> -#define RESTYPE_MIN LONG_MIN
> -#define RESTYPE_MAX LONG_MAX
> +#ifndef type
> +#define type double
> +#define roundit rint
> +#define dtype long
> +#define fn lrint
> #endif
>
> -#define RESTYPE_BITS (sizeof(RESTYPE) * 8)
> -
> -static const double
> -TWO52[2]={
> - 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
> - -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
> -};
> -
> -RESTYPE
> -LRINTNAME(double x)
> +/*
> + * C99 says we should not raise a spurious inexact exception when an
> + * invalid exception is raised. Unfortunately, the set of inputs
> + * that overflows depends on the rounding mode when 'dtype' has more
> + * significant bits than 'type'. Hence, we bend over backwards for the
> + * sake of correctness; an MD implementation could be more efficient.
> + */
> +dtype
> +fn(type x)
> {
> - u_int32_t i0, i1;
> - int e, s, shift;
> - RESTYPE res;
> -
> - GET_HIGH_WORD(i0, x);
> - e = i0 >> DBL_FRACHBITS;
> - s = e >> DBL_EXPBITS;
> - e = (e & 0x7ff) - DBL_EXP_BIAS;
> -
> - /* 1.0 x 2^31 (or 2^63) is already too large */
> - if (e >= (int)RESTYPE_BITS - 1)
> - return (s ? RESTYPE_MIN : RESTYPE_MAX); /* ??? unspecified */
> -
> - /* >= 2^52 is already an exact integer */
> - if (e < DBL_FRACBITS) {
> - volatile double t = x; /* clip extra precision */
> - /* round, using current direction */
> - t += TWO52[s];
> - t -= TWO52[s];
> - x = t;
> - }
> -
> - EXTRACT_WORDS(i0, i1, x);
> - e = ((i0 >> DBL_FRACHBITS) & 0x7ff) - DBL_EXP_BIAS;
> - i0 &= 0xfffff;
> - i0 |= (1 << DBL_FRACHBITS);
> -
> - if (e < 0)
> - return (0);
> -
> - shift = e - DBL_FRACBITS;
> - if (shift >=0)
> - res = (shift < RESTYPE_BITS ? (RESTYPE)i1 << shift : 0);
> - else
> - res = (shift > -RESTYPE_BITS ? (RESTYPE)i1 >> -shift : 0);
> - shift += 32;
> - if (shift >=0)
> - res |= (shift < RESTYPE_BITS ? (RESTYPE)i0 << shift : 0);
> - else
> - res |= (shift > -RESTYPE_BITS ? (RESTYPE)i0 >> -shift : 0);
> + fenv_t env;
> + dtype d;
>
> - return (s ? -res : res);
> + feholdexcept(&env);
> + d = (dtype)roundit(x);
> + if (fetestexcept(FE_INVALID))
> + feclearexcept(FE_INEXACT);
> + feupdateenv(&env);
> + return (d);
> }
> -DEF_STD(LRINTNAME);
> -LDBL_MAYBE_CLONE(LRINTNAME);
> +DEF_STD(fn);
> +LDBL_MAYBE_CLONE(fn);
> Index: lib/libm/src/s_lrintf.c
> ===================================================================
> RCS file: /cvs/src/lib/libm/src/s_lrintf.c,v
> retrieving revision 1.6
> diff -u -p -r1.6 s_lrintf.c
> --- lib/libm/src/s_lrintf.c 12 Sep 2016 19:47:02 -0000 1.6
> +++ lib/libm/src/s_lrintf.c 13 Oct 2021 23:12:11 -0000
> @@ -1,9 +1,8 @@
> -/* $OpenBSD: s_lrintf.c,v 1.6 2016/09/12 19:47:02 guenther Exp $ */
> -/* $NetBSD: lrintf.c,v 1.3 2004/10/13 15:18:32 drochner Exp $ */
> +/* $OpenBSD$ */
>
> /*-
> - * Copyright (c) 2004
> - * Matthias Drochner. All rights reserved.
> + * Copyright (c) 2005 David Schultz <d...@freebsd.org>
> + * All rights reserved.
> *
> * Redistribution and use in source and binary forms, with or without
> * modification, are permitted provided that the following conditions
> @@ -27,67 +26,34 @@
> * SUCH DAMAGE.
> */
>
> -#include <sys/types.h>
> -#include <sys/limits.h>
> +#include <fenv.h>
> #include <math.h>
> -#include <ieeefp.h>
> -#include <machine/ieee.h>
> -#include "math_private.h"
>
> -#ifndef LRINTNAME
> -#define LRINTNAME lrintf
> -#define RESTYPE long int
> -#define RESTYPE_MIN LONG_MIN
> -#define RESTYPE_MAX LONG_MAX
> +#ifndef type
> +#define type float
> +#define roundit rintf
> +#define dtype long
> +#define fn lrintf
> #endif
>
> -#define RESTYPE_BITS (sizeof(RESTYPE) * 8)
> -
> -static const float
> -TWO23[2]={
> - 8.3886080000e+06, /* 0x4b000000 */
> - -8.3886080000e+06, /* 0xcb000000 */
> -};
> -
> -RESTYPE
> -LRINTNAME(float x)
> +/*
> + * C99 says we should not raise a spurious inexact exception when an
> + * invalid exception is raised. Unfortunately, the set of inputs
> + * that overflows depends on the rounding mode when 'dtype' has more
> + * significant bits than 'type'. Hence, we bend over backwards for the
> + * sake of correctness; an MD implementation could be more efficient.
> + */
> +dtype
> +fn(type x)
> {
> - u_int32_t i0;
> - int e, s, shift;
> - RESTYPE res;
> -
> - GET_FLOAT_WORD(i0, x);
> - e = i0 >> SNG_FRACBITS;
> - s = e >> SNG_EXPBITS;
> - e = (e & 0xff) - SNG_EXP_BIAS;
> -
> - /* 1.0 x 2^31 (or 2^63) is already too large */
> - if (e >= (int)RESTYPE_BITS - 1)
> - return (s ? RESTYPE_MIN : RESTYPE_MAX); /* ??? unspecified */
> -
> - /* >= 2^23 is already an exact integer */
> - if (e < SNG_FRACBITS) {
> - volatile float t = x; /* clip extra precision */
> - /* round, using current direction */
> - t += TWO23[s];
> - t -= TWO23[s];
> - x = t;
> - }
> -
> - GET_FLOAT_WORD(i0, x);
> - e = ((i0 >> SNG_FRACBITS) & 0xff) - SNG_EXP_BIAS;
> - i0 &= 0x7fffff;
> - i0 |= (1 << SNG_FRACBITS);
> -
> - if (e < 0)
> - return (0);
> -
> - shift = e - SNG_FRACBITS;
> - if (shift >=0)
> - res = (shift < RESTYPE_BITS ? (RESTYPE)i0 << shift : 0);
> - else
> - res = (shift > -RESTYPE_BITS ? (RESTYPE)i0 >> -shift : 0);
> + fenv_t env;
> + dtype d;
>
> - return (s ? -res : res);
> + feholdexcept(&env);
> + d = (dtype)roundit(x);
> + if (fetestexcept(FE_INVALID))
> + feclearexcept(FE_INEXACT);
> + feupdateenv(&env);
> + return (d);
> }
> -DEF_STD(LRINTNAME);
> +DEF_STD(fn);
> Index: lib/libm/src/s_lrintl.c
> ===================================================================
> RCS file: /cvs/src/lib/libm/src/s_lrintl.c,v
> retrieving revision 1.3
> diff -u -p -r1.3 s_lrintl.c
> --- lib/libm/src/s_lrintl.c 15 Mar 2019 05:42:38 -0000 1.3
> +++ lib/libm/src/s_lrintl.c 13 Oct 2021 23:12:11 -0000
> @@ -31,9 +31,9 @@
>
> #ifndef type
> #define type long double
> -#define roundit rintl
> +#define roundit rintl
> #define dtype long
> -#define fn lrintl
> +#define fn lrintl
> #endif
>
> /*
