Fix for powerpc64 long double complex divide failure
commitd91056851c5c60f226e3192fb955d018b53eb66f
authorPatrick McGehearty <patrick.mcgehearty@oracle.com>
Sun, 3 Oct 2021 22:07:06 +0000 (3 18:07 -0400)
committerJeff Law <jeffreyalaw@gmail.com>
Sun, 3 Oct 2021 22:07:06 +0000 (3 18:07 -0400)
tree85fed2a9aa6ccb6611b56a96e31c915cc2fc5bbc
parentfdf0b6ce6c1cfa1c328c0c40473c71ca11fd8303
Fix for powerpc64 long double complex divide failure

- - - -

New in version 6: Due to an oversight (i.e. coding error), version 5
changed the use of __LIBGCC_TF_EPSILON__ to __LIBGCC_DF_EPSILON__ but
not the other LIBGCC_TF values. For correct execution of the long
double test case it is necessary to also switch to using
__LIBGCC_DF_MIN__. For consistency we also switch to using
__LIBGCC_DF_MAX__. LDBL_MIN is 2**53 times as larger than DBL_MIN.
The larger value causes the code to switch the order of computation
when it is not optimal, resulting in failure for one of the values
in the cdivchk_ld.c test. Using DBL_MIN does not cause that failure..

There may be opportunity for further refinement of IBM128 format
Long Double complex divide, but that's beyond the scope of this
patch.

- - - -

This revision adds a test in libgcc/libgcc2.c for when
"__LIBGCC_TF_MANT_DIG__ == 106" to use __LIBGCC_DF_EPSILON__ instead
of __LIBGCC_TF_EPSILON__. That is specific to IBM 128-bit format long
doubles where EPSILON is very, very small and 1/EPSILON oveflows to
infinity. This change avoids the overflow without affecting any other
platform. Discussion in the patch is adjusted to reflect this
limitation.

It does not make any changes to .../rs6000/_divkc3.c, leaving it to
use __LIBGCC_KF__*. That means the upstream gcc will not build in
older IBM environments that do not recognize the KF floating point
mode properly. Environments that do not need IBM longdouble support
do build cleanly.

- - - -
This patch addresses the failure of powerpc64 long double complex divide
in native ibm long double format after the patch "Practical improvement
to libgcc complex divide".

The new code uses the following macros which are intended to be mapped
to appropriate values according to the underlying hardware representation.
See https://gcc.gnu.org/bugzilla/show_bug.cgi?id=101104

RBIG     a value near the maximum representation
RMIN     a value near the minimum representation
         (but not in the subnormal range)
RMIN2    a value moderately less than 1
RMINSCAL the inverse of RMIN2
RMAX2    RBIG * RMIN2  - a value to limit scaling to not overflow

When "long double" values were not using the IEEE 128-bit format but
the traditional IBM 128-bit, the previous code used the LDBL values
which caused overflow for RMINSCAL. The new code uses the DBL values.

RBIG  LDBL_MAX = 0x1.fffffffffffff800p+1022
      DBL_MAX  = 0x1.fffffffffffff000p+1022

RMIN  LDBL_MIN = 0x1.0000000000000000p-969
RMIN  DBL_MIN  = 0x1.0000000000000000p-1022

RMIN2 LDBL_EPSILON = 0x0.0000000000001000p-1022 = 0x1.0p-1074
RMIN2 DBL_EPSILON  = 0x1.0000000000000000p-52

\e[ORMINSCAL 1/LDBL_EPSILON = inf (1.0p+1074 does not fit in IBM 128-bit).
         1/DBL_EPSILON  = 0x1.0000000000000000p+52

RMAX2 = RBIG * RMIN2 = 0x1.fffffffffffff800p-52
        RBIG * RMIN2 = 0x1.fffffffffffff000p+970

The MAX and MIN values have only modest changes since the maximum and
minimum values are about the same as for double precision.  The
EPSILON field is considerably different. Due to how very small values
can be represented in the lower 64 bits of the IBM 128-bit floating
point, EPSILON is extremely small, so far beyond the desired value
that inversion of the value overflows and even without the overflow,
the RMAX2 is so small as to eliminate most usage of the test.

The change has been tested on gcc135.fsffrance.org and gains the
expected improvements in accuracy for long double complex divide.

libgcc/
PR target/101104
* libgcc2.c (RMIN2, RMINSCAL, RMAX2):
Use more correct values for native IBM 128-bit.
libgcc/libgcc2.c