[PATCH v2 25/25] math: Use tanhf from CORE-MATH

DJ Delorie dj@redhat.com
Tue Dec 10 23:18:51 GMT 2024


Adhemerval Zanella <adhemerval.zanella@linaro.org> writes:
> The CORE-MATH implementation is correctly rounded (for any rounding mode)
> and shows slight better performance to the generic tanhf.

A comment about putting in the CORE-MATH url somewhere ;-)

A like the small comments on each line that clarify the values on/at
that line, but more thought needs to be put into using them
consistently, as they lack sufficient detail to always clarify their
meaning:

   return 0x1.921fb54442d18p+1;  /* +-inf */

Does the comment clarify what the hex number represents, or the
conditions under which we return it?

LGTM
Reviewed-by: DJ Delorie <dj@redhat.com>

> diff --git a/SHARED-FILES b/SHARED-FILES
> index 3bd4e7fb4a..032c407881 100644
> --- a/SHARED-FILES
> +++ b/SHARED-FILES
> @@ -330,3 +330,7 @@ sysdeps/ieee754/flt-32/e_sinhf.c:
>    (src/binary32/sinh/sinhf.c in CORE-MATH)
>    - the code was adapted to use glibc code style and internal
>      functions to handle errno, overflow, and underflow.
> +sysdeps/ieee754/flt-32/s_tanhf.c:
> +  (src/binary32/tanh/tanhf.c in CORE-MATH)
> +  - the code was adapted to use glibc code style and internal
> +    functions to handle errno, overflow, and underflow.

Ok.  As this is the last patch of the series, this is my last
opportunity to remind you to add an URL to CORE-MATH somewhere ;-)

> diff --git a/sysdeps/aarch64/libm-test-ulps b/sysdeps/aarch64/libm-test-ulps
> diff --git a/sysdeps/alpha/fpu/libm-test-ulps b/sysdeps/alpha/fpu/libm-test-ulps
> diff --git a/sysdeps/arc/fpu/libm-test-ulps b/sysdeps/arc/fpu/libm-test-ulps
> diff --git a/sysdeps/arc/nofpu/libm-test-ulps b/sysdeps/arc/nofpu/libm-test-ulps
> diff --git a/sysdeps/arm/libm-test-ulps b/sysdeps/arm/libm-test-ulps
> diff --git a/sysdeps/csky/fpu/libm-test-ulps b/sysdeps/csky/fpu/libm-test-ulps
> diff --git a/sysdeps/csky/nofpu/libm-test-ulps b/sysdeps/csky/nofpu/libm-test-ulps
> diff --git a/sysdeps/hppa/fpu/libm-test-ulps b/sysdeps/hppa/fpu/libm-test-ulps
> diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
> diff --git a/sysdeps/i386/i686/fpu/multiarch/libm-test-ulps b/sysdeps/i386/i686/fpu/multiarch/libm-test-ulps

Ok.

> diff --git a/sysdeps/ieee754/flt-32/s_tanhf.c b/sysdeps/ieee754/flt-32/s_tanhf.c
> index 2c12f04569..da03415e8a 100644
> --- a/sysdeps/ieee754/flt-32/s_tanhf.c
> +++ b/sysdeps/ieee754/flt-32/s_tanhf.c
> @@ -1,63 +1,90 @@
> -/* s_tanhf.c -- float version of s_tanh.c.
> - */
> +/* Correctly-rounded hyperbolic tangent function for binary32 value.
>  
> -/*
> - * ====================================================
> - * 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.
> - * ====================================================
> - */
> +Copyright (c) 2022-2024 Alexei Sibidanov.
>  
> -#if defined(LIBM_SCCS) && !defined(lint)
> -static char rcsid[] = "$NetBSD: s_tanhf.c,v 1.4 1995/05/10 20:48:24 jtc Exp $";
> -#endif
> +The original version of this file was copied from the CORE-MATH
> +project (file src/binary32/tanh/tanhf.c, revision bc385c2).
>  
> -#include <float.h>
> -#include <math.h>
> -#include <math_private.h>
> -#include <math-underflow.h>
> -#include <libm-alias-float.h>
> +Permission is hereby granted, free of charge, to any person obtaining a copy
> +of this software and associated documentation files (the "Software"), to deal
> +in the Software without restriction, including without limitation the rights
> +to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
> +copies of the Software, and to permit persons to whom the Software is
> +furnished to do so, subject to the following conditions:
>  
> -static const float one=1.0, two=2.0, tiny = 1.0e-30;
> -
> -float __tanhf(float x)
> -{
> -	float t,z;
> -	int32_t jx,ix;
> +The above copyright notice and this permission notice shall be included in all
> +copies or substantial portions of the Software.
>  
> -	GET_FLOAT_WORD(jx,x);
> -	ix = jx&0x7fffffff;
> +THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
> +IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
> +FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
> +AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
> +LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
> +OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
> +SOFTWARE.
> +*/
>  
> -    /* x is INF or NaN */
> -	if(ix>=0x7f800000) {
> -	    if (jx>=0) return one/x+one;    /* tanh(+-inf)=+-1 */
> -	    else       return one/x-one;    /* tanh(NaN) = NaN */
> -	}
> +#include <math.h>
> +#include <stdint.h>
> +#include <libm-alias-float.h>
> +#include "math_config.h"

Ok.

> -    /* |x| < 22 */
> -	if (ix < 0x41b00000) {		/* |x|<22 */
> -	    if (ix == 0)
> -		return x;		/* x == +-0 */
> -	    if (ix<0x24000000) 		/* |x|<2**-55 */
> -	      {
> -		math_check_force_underflow (x);
> -		return x*(one+x);    	/* tanh(small) = small */
> -	      }
> -	    if (ix>=0x3f800000) {	/* |x|>=1  */
> -		t = __expm1f(two*fabsf(x));
> -		z = one - two/(t+two);
> -	    } else {
> -	        t = __expm1f(-two*fabsf(x));
> -	        z= -t/(t+two);
> -	    }
> -    /* |x| > 22, return +-1 */
> -	} else {
> -	    z = one - tiny;		/* raised inexact flag */

Ok.

> +float
> +__tanhf (float x)
> +{
> +  double z = x;
> +  uint32_t ux = asuint (x);
> +  int e = (ux >> 23) & 0xff;
> +  if (__glibc_unlikely (e == 0xff))
> +    {
> +      if (ux << 9)
> +	return x + x; /* nan */
> +      static const float ir[] = { 1.0f, -1.0f };
> +      return ir[ux >> 31]; /* +-inf */
> +    }

Ok.  I'm a bit conflicted about the +-inf comment as the value being
returned is +-1.0f; I realize it means "input was +-inf" but these
comments are typically used on braces or if() lines,  Perhaps "/* x ==
+-inf */" would be more obvious?

(the previous one is not ambiguous, because both input and output are
NaNs)


> +  if (__glibc_unlikely (e < 115))
> +    {
> +      if (__glibc_unlikely (e < 102))
> +	{
> +	  if (__glibc_unlikely ((ux << 1) == 0))
> +	    return x;
> +	  return fmaf (-x, fabsf (x), x);
>  	}
> +      float x2 = x * x;
> +      return fmaf (x, -0x1.555556p-2f * x2, x);
> +    }

Ok.

> +  if ((ux << 1) > (0x41102cb3u << 1))
> +    return copysignf (1.0f, x) - copysignf (0x1p-25f, x);

Ok.  Obviously needed this way for rounding.

> +  double z2 = z * z;
> +  double z4 = z2 * z2;
> +  double z8 = z4 * z4;
> +  static const double cn[] =
> +    {
> +      0x1p+0,                0x1.30877b8b72d33p-3,  0x1.694aa09ae9e5ep-8,
> +      0x1.4101377abb729p-14, 0x1.e0392b1db0018p-22, 0x1.2533756e546f7p-30,
> +      0x1.d62e5abe6ae8ap-41, 0x1.b06be534182dep-54
> +    };
> +  static const double cd[] =
> +    {
> +      0x1p+0,                0x1.ed99131b0ebeap-2,  0x1.0d27ed6c95a69p-5,
> +      0x1.7cbdaca0e9fccp-11, 0x1.b4e60b892578ep-18, 0x1.a6f707c5c71abp-26,
> +      0x1.35a8b6e2cd94cp-35, 0x1.ca8230677aa01p-47
> +    };
> +  double n0 = cn[0] + z2 * cn[1];
> +  double n2 = cn[2] + z2 * cn[3];
> +  double n4 = cn[4] + z2 * cn[5];
> +  double n6 = cn[6] + z2 * cn[7];
> +  n0 += z4 * n2;
> +  n4 += z4 * n6;
> +  n0 += z8 * n4;
> +  double d0 = cd[0] + z2 * cd[1];
> +  double d2 = cd[2] + z2 * cd[3];
> +  double d4 = cd[4] + z2 * cd[5];
> +  double d6 = cd[6] + z2 * cd[7];
> +  d0 += z4 * d2;
> +  d4 += z4 * d6;
> +  d0 += z8 * d4;
> +  double r = z * n0 / d0;
> +  return r;
>  }
>  libm_alias_float (__tanh, tanh)

Ok.

> diff --git a/sysdeps/loongarch/lp64/libm-test-ulps b/sysdeps/loongarch/lp64/libm-test-ulps
> diff --git a/sysdeps/m68k/m680x0/fpu/libm-test-ulps b/sysdeps/m68k/m680x0/fpu/libm-test-ulps
> diff --git a/sysdeps/microblaze/libm-test-ulps b/sysdeps/microblaze/libm-test-ulps
> diff --git a/sysdeps/mips/mips32/libm-test-ulps b/sysdeps/mips/mips32/libm-test-ulps
> diff --git a/sysdeps/mips/mips64/libm-test-ulps b/sysdeps/mips/mips64/libm-test-ulps
> diff --git a/sysdeps/or1k/fpu/libm-test-ulps b/sysdeps/or1k/fpu/libm-test-ulps
> diff --git a/sysdeps/or1k/nofpu/libm-test-ulps b/sysdeps/or1k/nofpu/libm-test-ulps
> diff --git a/sysdeps/powerpc/fpu/libm-test-ulps b/sysdeps/powerpc/fpu/libm-test-ulps
> diff --git a/sysdeps/powerpc/nofpu/libm-test-ulps b/sysdeps/powerpc/nofpu/libm-test-ulps
> diff --git a/sysdeps/riscv/nofpu/libm-test-ulps b/sysdeps/riscv/nofpu/libm-test-ulps
> diff --git a/sysdeps/riscv/rvd/libm-test-ulps b/sysdeps/riscv/rvd/libm-test-ulps
> diff --git a/sysdeps/s390/fpu/libm-test-ulps b/sysdeps/s390/fpu/libm-test-ulps
> diff --git a/sysdeps/sh/libm-test-ulps b/sysdeps/sh/libm-test-ulps
> diff --git a/sysdeps/sparc/fpu/libm-test-ulps b/sysdeps/sparc/fpu/libm-test-ulps
> diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps

Ok.



More information about the Libc-alpha mailing list