https://gcc.gnu.org/bugzilla/show_bug.cgi?id=120993
Bug ID: 120993
Summary: powerpc64le LDBL_NORM_MAX does not conform to C23
Product: gcc
Version: 15.1.1
Status: UNCONFIRMED
Severity: normal
Priority: P3
Component: c
Assignee: unassigned at gcc dot gnu.org
Reporter: eggert at cs dot ucla.edu
Target Milestone: ---
ISO C23 §5.2.5.3.3¶32 says that LDBL_NORM_MAX must equal:
(1 - FLT_RADIX**-LDBL_MANT_DIG) * FLT_RADIX**LDBL_MAX_EXP
and on powerpc64le where FLT_RADIX=2, LDBL_MANT_DIG=106, LDBL_MAX_EXP=1024 this
must equal:
(1 - 2**-106) * 2**1024
= 2**1024 - 2**918
However, on this platform GCC defines LDBL_NORM_MAX to be
8.98846567431157953864652595394501288e+307L, which is 2**1023 - 2**917, i.e.,
half the value required by C23.
Presumably this discrepancy is due to the funky nature of the "double double"
arithmetic on this platform, where things get dicey above 2**1023 as overflow
can occur even before you get to 2**1024. However, the fact remains that
there's a disagreement between C23 and how GCC behaves. We ran across this
problem recently in some Gnulib test cases on GCC 15 powerpc64le.
One simple fix would be to change LDBL_MANT_DIG, LDBL_MAX_EXP, and
LDBL_NORM_MAX to equal DBL_MANT_DIG, DBL_MAX_EXP and DBL_NORM_MAX respectively.
There would be no need to change LDBL_MAX or any other part of the long double
implementation. This change would conform to C23, and would fix the problems we
ran into.
For an example of the problem, compile the following program test1.c with "gcc
-std=gnu23 test1.c":
#include <float.h>
#include <stdio.h>
int
main ()
{
long double epsilon = 1;
for (int i = 0; i < LDBL_MANT_DIG; i++)
epsilon /= FLT_RADIX;
long double one_minus_epsilon = 1 - epsilon;
long double ldbl_norm_max = one_minus_epsilon;
for (int i = 0; i < LDBL_MAX_EXP; i++)
ldbl_norm_max *= FLT_RADIX;
printf ("epsilon = %.27La\n", epsilon);
printf ("one_minus_epsilon = %.27La\n", one_minus_epsilon);
printf ("ldbl_norm_max = %.27La\n", ldbl_norm_max);
printf ("LDBL_NORM_MAX = %.27La\n", LDBL_NORM_MAX);
return ldbl_norm_max != LDBL_NORM_MAX;
}
This should exit with status 0 because all calculations are exact with no
rounding. However, on powerpc64le it exits with status 1 after outputting:
epsilon = 0x1.000000000000000000000000000p-106
one_minus_epsilon = 0x1.ffffffffffffffffffffffffff8p-1
ldbl_norm_max = inf
LDBL_NORM_MAX = 0x1.ffffffffffffffffffffffffff8p+1022
where the third output line is incorrect due to the funky nature of
calculations on this platform.
