cbrtl.c

Branch

Show log
Commit

Hash : a923e066
Author :
Date : 2025-09-10T22:55:08

math-h: Remove support for IRIX.

* lib/math.in.h (FP_ILOGB0, FP_ILOGBNAN): Remove code for IRIX.
(cbrtf): Test HAVE_CBRTF instead of HAVE_DECL_CBRTF.
(cbrtl): Test HAVE_CBRTL instead of HAVE_DECL_CBRTL.
(copysignf): Test HAVE_COPYSIGNF instead of HAVE_DECL_COPYSIGNF.
(exp2f): Test HAVE_EXP2F instead of HAVE_DECL_EXP2F.
(exp2): Test HAVE_EXP2 instead of HAVE_DECL_EXP2.
(exp2l): Test HAVE_EXP2L instead of HAVE_DECL_EXP2L.
(expm1l): Test HAVE_EXPM1L instead of HAVE_DECL_EXPM1L.
(log2f): Test HAVE_LOG2F instead of HAVE_DECL_LOG2F.
(log2): Test HAVE_LOG2 instead of HAVE_DECL_LOG2.
(log2l): Test HAVE_LOG2L instead of HAVE_DECL_LOG2L.
(remainder): Don't test HAVE_DECL_REMAINDER.
(remainderl): Test HAVE_REMAINDERL instead of HAVE_DECL_REMAINDERL.
(rintf): Test HAVE_RINTF instead of HAVE_DECL_RINTF.
(roundf): Test HAVE_ROUNDF instead of HAVE_DECL_ROUNDF.
(round): Test HAVE_ROUND instead of HAVE_DECL_ROUND.
(roundl): Test HAVE_ROUNDL instead of HAVE_DECL_ROUNDL.
(isnanf, isnand): Remove code for IRIX.
* m4/math_h.m4 (gl_MATH_H_DEFAULTS): Initialize HAVE_COPYSIGNF,
HAVE_EXP2, HAVE_EXP2F, HAVE_EXP2L, HAVE_EXPM1L, HAVE_LOG2, HAVE_LOG2F,
HAVE_LOG2L, HAVE_REMAINDERL, HAVE_RINTF, HAVE_ROUND, HAVE_ROUNDF,
HAVE_ROUNDL. Don't initialize HAVE_DECL_CBRTF, HAVE_DECL_CBRTL,
HAVE_DECL_COPYSIGNF, HAVE_DECL_EXP2, HAVE_DECL_EXP2F, HAVE_DECL_EXP2L,
HAVE_DECL_EXPM1L, HAVE_DECL_LOG2, HAVE_DECL_LOG2F, HAVE_DECL_LOG2L,
HAVE_DECL_REMAINDER, HAVE_DECL_REMAINDERL, HAVE_DECL_RINTF,
HAVE_DECL_ROUND, HAVE_DECL_ROUNDF, HAVE_DECL_ROUNDL.
* modules/math-h (Makefile.am): Substitute HAVE_COPYSIGNF, HAVE_EXP2,
HAVE_EXP2F, HAVE_EXP2L, HAVE_EXPM1L, HAVE_LOG2, HAVE_LOG2F, HAVE_LOG2L,
HAVE_REMAINDERL, HAVE_RINTF, HAVE_ROUND, HAVE_ROUNDF, HAVE_ROUNDL.
Don't substitute HAVE_DECL_CBRTF, HAVE_DECL_CBRTL, HAVE_DECL_COPYSIGNF,
HAVE_DECL_EXP2, HAVE_DECL_EXP2F, HAVE_DECL_EXP2L, HAVE_DECL_EXPM1L,
HAVE_DECL_LOG2, HAVE_DECL_LOG2F, HAVE_DECL_LOG2L, HAVE_DECL_REMAINDER,
HAVE_DECL_REMAINDERL, HAVE_DECL_RINTF, HAVE_DECL_ROUND,
HAVE_DECL_ROUNDF, HAVE_DECL_ROUNDL.
* m4/cbrtf.m4 (gl_FUNC_CBRTF_WORKS): Remove macro.
(gl_FUNC_CBRTF): Don't invoke it. Don't set HAVE_DECL_CBRTF.
* m4/cbrtl.m4 (gl_FUNC_CBRTL): Don't set HAVE_DECL_CBRTL. Remove test
regarding cbrtl (-0.0L).
* lib/cbrtl.c (cbrtl): Remove code for IRIX.
* m4/copysignf.m4 (gl_FUNC_COPYSIGNF): Don't set HAVE_DECL_COPYSIGNF.
* m4/exp2f.m4 (gl_FUNC_EXP2F): Don't set HAVE_DECL_EXP2F.
* m4/exp2.m4 (gl_FUNC_EXP2): Don't set HAVE_DECL_EXP2.
* m4/exp2l.m4 (gl_FUNC_EXP2L): Don't set HAVE_DECL_EXP2L. Remove code
for IRIX. Remove test regarding exp2l (-Inf).
* m4/expm1l.m4 (gl_FUNC_EXPM1L): Don't set HAVE_DECL_EXPM1L. Remove code
for IRIX.
* m4/log2f.m4 (gl_FUNC_LOG2F): Don't set HAVE_DECL_LOG2F.
* m4/log2.m4 (gl_FUNC_LOG2): Don't set HAVE_DECL_LOG2.
* lib/log2.c (log2): Remove code for IRIX.
* m4/log2l.m4 (gl_FUNC_LOG2L): Don't set HAVE_DECL_LOG2L. Remove code
for IRIX.
* lib/log2l.c (log2l): Remove code for IRIX.
* m4/remainder.m4 (gl_FUNC_REMAINDER): Don't set HAVE_DECL_REMAINDER.
* m4/remainderl.m4 (gl_FUNC_REMAINDERL): Don't set HAVE_DECL_REMAINDERL.
Remove code for IRIX.
* m4/rintf.m4 (gl_FUNC_RINTF): Don't set HAVE_DECL_RINTF.
* m4/roundf.m4 (gl_FUNC_ROUNDF): Don't set HAVE_DECL_ROUNDF.
* m4/round.m4 (gl_FUNC_ROUND): Don't set HAVE_DECL_ROUND.
* m4/roundl.m4 (gl_FUNC_ROUNDL): Don't set HAVE_DECL_ROUNDL.

Download

lib/cbrtl.c

/* Compute cubic root of long double value.
   Copyright (C) 2012-2025 Free Software Foundation, Inc.
   Cephes Math Library Release 2.2: January, 1991
   Copyright 1984, 1991 by Stephen L. Moshier
   Adapted for glibc October, 2001.

   This file 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 3 of the
   License, or (at your option) any later version.

   This file 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 this program.  If not, see <https://www.gnu.org/licenses/>.  */

#include <config.h>

/* Specification.  */
#include <math.h>

#if HAVE_SAME_LONG_DOUBLE_AS_DOUBLE

long double
cbrtl (long double x)
{
  return cbrt (x);
}

#else

/* Code based on glibc/sysdeps/ieee754/ldbl-128/s_cbrtl.c.  */

/*                                                      cbrtl.c
 *
 *      Cube root, long double precision
 *
 *
 *
 * SYNOPSIS:
 *
 * long double x, y, cbrtl();
 *
 * y = cbrtl( x );
 *
 *
 *
 * DESCRIPTION:
 *
 * Returns the cube root of the argument, which may be negative.
 *
 * Range reduction involves determining the power of 2 of
 * the argument.  A polynomial of degree 2 applied to the
 * mantissa, and multiplication by the cube root of 1, 2, or 4
 * approximates the root to within about 0.1%.  Then Newton's
 * iteration is used three times to converge to an accurate
 * result.
 *
 *
 *
 * ACCURACY:
 *
 *                      Relative error:
 * arithmetic   domain     # trials      peak         rms
 *    IEEE       -8,8       100000      1.3e-34     3.9e-35
 *    IEEE    exp(+-707)    100000      1.3e-34     4.3e-35
 *
 */

static const long double CBRT2 = 1.259921049894873164767210607278228350570251L;
static const long double CBRT4 = 1.587401051968199474751705639272308260391493L;
static const long double CBRT2I = 0.7937005259840997373758528196361541301957467L;
static const long double CBRT4I = 0.6299605249474365823836053036391141752851257L;

long double
cbrtl (long double x)
{
  if (isfinite (x) && x != 0.0L)
    {
      int e, rem, sign;
      long double z;

      if (x > 0)
        sign = 1;
      else
        {
          sign = -1;
          x = -x;
        }

      z = x;
      /* extract power of 2, leaving mantissa between 0.5 and 1  */
      x = frexpl (x, &e);

      /* Approximate cube root of number between .5 and 1,
         peak relative error = 1.2e-6  */
      x = ((((1.3584464340920900529734e-1L * x
              - 6.3986917220457538402318e-1L) * x
             + 1.2875551670318751538055e0L) * x
            - 1.4897083391357284957891e0L) * x
           + 1.3304961236013647092521e0L) * x + 3.7568280825958912391243e-1L;

      /* exponent divided by 3 */
      if (e >= 0)
        {
          rem = e;
          e /= 3;
          rem -= 3 * e;
          if (rem == 1)
            x *= CBRT2;
          else if (rem == 2)
            x *= CBRT4;
        }
      else
        {                           /* argument less than 1 */
          e = -e;
          rem = e;
          e /= 3;
          rem -= 3 * e;
          if (rem == 1)
            x *= CBRT2I;
          else if (rem == 2)
            x *= CBRT4I;
          e = -e;
        }

      /* multiply by power of 2 */
      x = ldexpl (x, e);

      /* Newton iteration */
      x -= (x - (z / (x * x))) * 0.3333333333333333333333333333333333333333L;
      x -= (x - (z / (x * x))) * 0.3333333333333333333333333333333333333333L;
      x -= (x - (z / (x * x))) * 0.3333333333333333333333333333333333333333L;

      if (sign < 0)
        x = -x;
      return x;
    }
  else
    return x + x;
}

#endif

Download

kc3-lang/gnulib/lib/cbrtl.c

Commit

lib/cbrtl.c